PG2PLplot
Aid in the conversion from PGPlot to PLplot in Fortran programs
pg2plplot.f90
Go to the documentation of this file.
1 !> \file pg2plplot.f90 Contains pgplot-to-plplot bindings (i.e., call PLplot from PGplot commands)
2 
3 !
4 ! LICENCE:
5 !
6 ! Copyright 2010-2015 AstroFloyd, joequant
7 !
8 ! This file is part of the PG2PLplot package.
9 !
10 ! This is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published
11 ! by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
12 !
13 ! This software is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
14 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
15 !
16 ! You should have received a copy of the GNU General Public License along with this code (LICENCE). If not, see
17 ! <http://www.gnu.org/licenses/>.
18 !
19 !
20 ! PG2PLplot can be found on http://pg2plplot.sourceforge.net
21 ! Some routines are taken from libSUFR, http://libsufr.sourceforge.net
22 
23 
24 
25 !***********************************************************************************************************************************
26 !> \brief PG2PLplot module
27 
28 module pg2plplot
29  use plplot, only: plflt
30  implicit none
31  private plflt
32  save
33 
34  !> Conversion factor for the character height
35  real(plflt), parameter :: ch_fac = 0.35_plflt
36  !logical, parameter :: compatibility_warnings = .false. ! Don't warn
37  logical, parameter :: compatibility_warnings = .true. ! Do warn
38 
39  real(plflt) :: xcur=0., ycur=0.
40  integer :: save_level = 0
41  integer, parameter :: max_level = 20
42  integer, parameter :: max_open = 100
43  real(plflt), parameter :: mm_per_inch = 25.4
44  integer :: save_ffamily(max_level)
45  integer :: save_fstyle(max_level)
46  integer :: save_fweight(max_level)
47  integer :: save_lwidth(max_level)
48  integer :: save_lstyle(max_level)
49  integer :: save_color(max_level)
50  integer :: cur_lwidth=1, cur_color=1, cur_lstyle=1
51  integer :: devid=0
52  logical :: is_init
53  real(plflt) :: paper_width, paper_ratio
54 
55 contains
56 
57  !*********************************************************************************************************************************
58  subroutine do_init()
59  implicit none
60  if(.not. is_init) then
61  call plinit()
62  call plbop()
63  is_init = .true.
64  end if
65  end subroutine do_init
66  !*********************************************************************************************************************************
67 
68 end module pg2plplot
69 !***********************************************************************************************************************************
70 
71 
72 
73 !***********************************************************************************************************************************
74 !> \brief Set line style
75 !!
76 !! \param ls Line style
77 
78 subroutine pgsls(ls)
79  use pg2plplot, only : cur_lstyle
80  implicit none
81  integer, intent(in) :: ls
82  integer :: ls1,ls2, styles(5)
83 
84  cur_lstyle = ls
85  styles = (/1,4,5,2,6/)
86 
87  ls1 = ls ! 1: solid, 2: dashes, 3: dash-dotted, 4: dotted, 5: dash-dot-dot-dot
88  if(ls1.lt.1.or.ls1.gt.5) ls1 = 1
89 
90 
91  ls2 = styles(ls1) ! 1: solid, 2: short dashes, 3: long dashes, 4: long dashes, short gaps, 5: long-short dashes,
92  ! 6: long-short dashes, long-short gaps, 7: ?, 8: ?
93 
94  call pllsty(ls2)
95 
96  !print*,'pgsls: ',ls,ls1,ls2
97 
98 end subroutine pgsls
99 !***********************************************************************************************************************************
100 
101 !***********************************************************************************************************************************
102 !> \brief Query line style
103 !!
104 !! \param ls Line style
105 
106 subroutine pgqls(ls)
107  use pg2plplot, only : cur_lstyle
108  implicit none
109  integer, intent(out) :: ls
110  ls = cur_lstyle
111 end subroutine pgqls
112 !***********************************************************************************************************************************
113 
114 !***********************************************************************************************************************************
115 !> \brief Set line width
116 !!
117 !! \param lw Line width
118 
119 subroutine pgslw(lw)
120  use plplot, only: plflt
121  use pg2plplot, only : cur_lwidth
122  implicit none
123  integer, intent(in) :: lw
124  integer :: lw1
125  real(kind=plflt) :: lw2
126 
127  cur_lwidth = lw
128  lw1 = max(min(lw,201),1)
129  lw2 = real(lw1 - 1)
130 
131  call plwidth(lw2)
132 
133  !print*,'pgslw: ',lw,lw1,lw2
134 
135 end subroutine pgslw
136 !***********************************************************************************************************************************
137 
138 !***********************************************************************************************************************************
139 !> \brief Query line width
140 !!
141 !! \param lw Line width
142 
143 subroutine pgqlw(lw)
145  implicit none
146  integer, intent(out) :: lw
147 
148  lw = max(1,cur_lwidth)
149 
150 end subroutine pgqlw
151 !***********************************************************************************************************************************
152 
153 !***********************************************************************************************************************************
154 !> \brief Set character font
155 !!
156 !! \param cf Character font
157 
158 subroutine pgscf(cf)
159  implicit none
160  integer, intent(in) :: cf
161 
162  call plfont(cf)
163 
164 end subroutine pgscf
165 !***********************************************************************************************************************************
166 
167 !***********************************************************************************************************************************
168 !> \brief Query character font
169 !!
170 !! \retval cf Character font
171 
172 subroutine pgqcf(cf)
173  implicit none
174  integer, intent(out) :: cf
175  integer :: family, style, weight
176  call plgfont(family, style, weight)
177  if(style .eq. 1) then
178  cf = 3
179  return
180  else if(family .eq. 1) then
181  cf = 2
182  return
183  else if(family .eq. 3) then
184  cf = 4
185  return
186  end if
187  cf = 1
188 end subroutine pgqcf
189 !***********************************************************************************************************************************
190 
191 !***********************************************************************************************************************************
192 !> \brief Inquire PGPLOT general information - dummy routine
193 !!
194 !! \param item
195 !! \retval value
196 !! \retval length
197 
198 subroutine pgqinf(item, value, length)
199  implicit none
200  character, intent(in) :: item*(*)
201  character, intent(out) :: value*(*)
202  integer, intent(out) :: length
203  value = item(1:1)
204  value = 'Dummy'
205  length = 5
206 end subroutine pgqinf
207 !***********************************************************************************************************************************
208 
209 !***********************************************************************************************************************************
210 !> \brief Set colour index
211 !!
212 !! \param ci1 Colour index
213 
214 subroutine pgsci(ci1)
215  use pg2plplot, only: cur_color
216  implicit none
217  integer, intent(in) :: ci1
218  integer :: ci2,colours(0:15)
219  cur_color = ci1
220  ci2 = 15 ! White
221  ci2 = ci1
222  colours = (/0,15,1,3,9,7,13,2,8,12,4,11,10,5,7,7/)
223  if(ci1.ge.0.and.ci1.le.15) ci2 = colours(ci1)
224 
225  call plcol0(ci2)
226 
227  !write(6,'(A,2I6)')' pgsci: ',ci1,ci2
228 
229 end subroutine pgsci
230 !***********************************************************************************************************************************
231 
232 !***********************************************************************************************************************************
233 !> \brief Set colour representation
234 !!
235 !! \param ci1 Colour index
236 !! \param r Red colour (0-1)
237 !! \param g Green colour (0-1)
238 !! \param b Blue colour (0-1)
239 
240 subroutine pgscr(ci1, r,g,b)
241  implicit none
242  integer, intent(in) :: ci1
243  real, intent(in) :: r,g,b
244  integer :: ci2,ri,gi,bi,colours(0:15)
245 
246  ri = nint(r*255)
247  gi = nint(g*255)
248  bi = nint(b*255)
249 
250  colours = (/0,15,1,3,9,7,13,2,8,12,4,11,10,5,7,7/)
251  ci2 = ci1
252  if(ci1.ge.0.and.ci1.le.15) ci2 = colours(ci1)
253  if(ci2.gt.15) call plscmap0n(256) ! Allows indices 0-255
254 
255  call plscol0(ci2, ri,gi,bi)
256 
257  !write(6,'(A,2I6, 5x,3F10.3, 5x,3I6)')' pgscr: ',ci1,ci2, r,g,b, ri,gi,bi
258 
259 end subroutine pgscr
260 !***********************************************************************************************************************************
261 
262 !***********************************************************************************************************************************
263 !> \brief Query colour representation
264 !!
265 !! \param ci1 Colour index
266 !! \retval r Red colour (0-1)
267 !! \retval g Green colour (0-1)
268 !! \retval b Blue colour (0-1)
269 
270 subroutine pgqcr(ci1, r,g,b)
271  implicit none
272  integer, intent(in) :: ci1
273  real, intent(out) :: r,g,b
274  integer :: ci2,ri,gi,bi,colours(0:15)
275 
276  ci2 = ci1
277  colours = (/0,15,1,3,9,7,13,2,8,12,4,11,10,5,7,7/)
278  if(ci1.ge.0.and.ci1.le.15) ci2 = colours(ci1)
279 
280  call plgcol0(ci2,ri,gi,bi)
281 
282  r = real(ri)/255.
283  g = real(gi)/255.
284  b = real(bi)/255.
285 
286  !write(6,'(A,2I6,5x,3I6,5x,3F10.3)')' pgqcr: ',ci1,ci2, ri,gi,bi, r,g,b
287 
288 end subroutine pgqcr
289 !***********************************************************************************************************************************
290 
291 !***********************************************************************************************************************************
292 !> \brief Set colour-index range for pggray and pgimag - dummy routine!
293 !!
294 !! \param ci1 Lower colour index
295 !! \param ci2 Upper colour index
296 
297 subroutine pgscir(ci1,ci2)
299  implicit none
300  integer, intent(in) :: ci1,ci2
301  integer :: tmp
302  integer, save :: warn
303 
304  tmp = ci1
305  tmp = ci2
306  tmp = tmp ! Avoid 'variable is set but not used' warnings from compiler for dummy variable
307 
308  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
309  if(warn.ne.123) call warn_dummy_routine('pgscir', '')
310  warn = 123
311 
312 end subroutine pgscir
313 !***********************************************************************************************************************************
314 
315 
316 !***********************************************************************************************************************************
317 !> \brief Query colour-index range for pggray and pgimag - dummy routine!
318 !!
319 !! \param ci1 Lower colour index
320 !! \param ci2 Upper colour index
321 
322 subroutine pgqcir(ci1,ci2)
324  implicit none
325  integer, intent(out) :: ci1,ci2
326  integer, save :: warn
327 
328  ci1 = 0
329  ci2 = 255
330 
331  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
332  if(warn.ne.123) call warn_dummy_routine('pgqcir', '')
333  warn = 123
334 
335 end subroutine pgqcir
336 !***********************************************************************************************************************************
337 
338 
339 
340 !***********************************************************************************************************************************
341 !> \brief Set fill style
342 !!
343 !! \param fs Fill style (1-4, no match for 2: outline)
344 
345 subroutine pgsfs(fs)
346  implicit none
347  integer, intent(in) :: fs
348  integer :: fs1,fs2, styles(4)
349 
350  fs1 = fs
351  if(fs.lt.1.or.fs.gt.4) fs1 = 1
352 
353  ! fs1: 1: solid, 2: outline, 3-hatched, 4-cross-hatched
354  ! fs2: 0: solid, 1: H lines, 2: V lines, 3: hatches 45d up, 4: hatches 45d down, 5: hatches 30d up, 6: hatches 30d down
355  ! 7: upright cross-hatces, 8: 45d cross-hatces
356 
357  styles = (/0,0,3,8/)
358 
359  fs2 = styles(fs1)
360 
361  call plpsty(fs2)
362 
363  !print*,'pgsfs: ',fs1,fs2
364 
365 end subroutine pgsfs
366 !***********************************************************************************************************************************
367 
368 !***********************************************************************************************************************************
369 !> \brief Set hash style
370 !!
371 !! \param ang Angle of the lines (deg)
372 !! \param sep Spacing (in % of view-surface size): >0!
373 !! \param ph Phase of hatching. Use e.g. 0.0 and 0.5 for double hatching - dummy variable: not used
374 
375 subroutine pgshs(ang, sep, ph)
376  implicit none
377  real, intent(in) :: ang, sep, ph
378  integer :: inc, del, tmp
379 
380  inc = nint(ang*10.) ! Tenths of a degree
381  del = nint(sep*1000) ! Spacing in micrometers(!)
382  tmp = nint(ph)
383  tmp = tmp ! Avoid 'variable is set but not used' warnings from compiler for dummy variable
384 
385  call plpat(1, inc, del)
386 
387 end subroutine pgshs
388 !***********************************************************************************************************************************
389 
390 
391 !***********************************************************************************************************************************
392 !> \brief Set character height
393 !!
394 !! \param ch Character height
395 
396 subroutine pgsch(ch)
397  use plplot, only: plflt
398  use pg2plplot, only: ch_fac
399 
400  implicit none
401  real, intent(in) :: ch
402  real(kind=plflt) :: ch1,ch2
403 
404  ch1 = 0.0_plflt ! 0: don't change
405  ch2 = ch * ch_fac
406 
407  call plschr(ch1,ch2)
408 
409  !print*
410  !print*,'pgsch1: ',ch,ch1,ch2
411  !call plgchr(ch1,ch2)
412  !print*,'pgsch2: ',ch,ch1,ch2
413 
414 end subroutine pgsch
415 !***********************************************************************************************************************************
416 
417 !***********************************************************************************************************************************
418 !> \brief Inquire character height in a variety of units
419 !!
420 !! \param unit 0: normalized device coordinates, 1: in, 2: mm, 3: pixels, 4: world coordinates
421 !! \retval xch The character height for text written with a vertical baseline
422 !! \retval ych The character height for text written with a horizontal baseline (the usual case)
423 
424 subroutine pgqcs(unit, xch, ych)
425  use plplot, only: plflt
426  use pg2plplot, only: mm_per_inch
427 
428  implicit none
429  integer, intent(in) :: unit
430  real, intent(out) :: xch, ych
431  real(kind=plflt) :: ch1,ch2
432  real(kind=plflt) xp, yp, xleng, yleng, xoff, yoff
433  call plgchr(ch1,ch2)
434  if(unit.eq.2) then
435  xch = real(ch2)
436  ych = real(ch2)
437  else if(unit .eq. 1) then
438  xch = real(ch2 / mm_per_inch)
439  ych = real(ch2 / mm_per_inch)
440  else if(unit .eq. 0) then
441  call plgpage(xp, yp, xleng, yleng, xoff, yoff)
442  xch = real(ch2 / yleng)
443  ych = real(ch2 / yleng)
444  else
445  print *, 'unknown unit in pgqcs', unit
446  end if
447 
448 end subroutine pgqcs
449 !***********************************************************************************************************************************
450 
451 !***********************************************************************************************************************************
452 !> \brief Query character height
453 !!
454 !! \param ch Character height
455 
456 subroutine pgqch(ch)
457  use plplot, only: plflt
458  use pg2plplot, only: ch_fac
459 
460  implicit none
461  real, intent(out) :: ch
462  real(kind=plflt) :: ch1,ch2
463 
464  call plgchr(ch1,ch2)
465  ch = real(ch2/(ch1*ch_fac))
466 
467  !print*,'pgqch: ',ch,ch1,ch2
468 
469 end subroutine pgqch
470 !***********************************************************************************************************************************
471 
472 
473 !***********************************************************************************************************************************
474 !> \brief Set arrow head - dummy routine!
475 !!
476 !! \param fs Fill style
477 !! \param angle Arrow-point angle
478 !! \param barb Arrow back cut
479 
480 subroutine pgsah(fs, angle, barb)
482  implicit none
483  integer, intent(in) :: fs
484  real, intent(in) :: angle, barb
485  integer :: tmp
486  integer, save :: warn
487 
488  tmp = fs
489  tmp = nint(angle)
490  tmp = nint(barb)
491  tmp = tmp ! Avoid 'variable is set but not used' warnings from compiler for dummy variable
492 
493  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
494  if(warn.ne.123) call warn_dummy_routine('pgsah', '')
495  warn = 123
496 
497 end subroutine pgsah
498 !***********************************************************************************************************************************
499 
500 
501 !***********************************************************************************************************************************
502 !> \brief Draw a line
503 !!
504 !! \param n Number of points
505 !! \param x1 X-values of points
506 !! \param y1 Y-values of points
507 
508 subroutine pgline(n,x1,y1)
509  use plplot, only: plflt, plline
510 
511  implicit none
512  integer, intent(in) :: n
513  real, intent(in) :: x1(n),y1(n)
514  real(kind=plflt) :: x2(n),y2(n)
515 
516  x2 = x1
517  y2 = y1
518 
519  call plline(x2,y2)
520 
521  !write(6,'(A,99ES16.9)')' pgline1: ',x1(1:min(n,9)),y1(1:min(n,9))
522  !write(6,'(A,99ES16.9)')' pgline2: ',x2(1:min(n,9)),y2(1:min(n,9))
523 
524 end subroutine pgline
525 !***********************************************************************************************************************************
526 
527 !***********************************************************************************************************************************
528 !> \brief Draw an arrow - only a line is drawn, arrows not supported in PLplot!
529 !!
530 !! \param x1 X-value of start point
531 !! \param y1 Y-value of start point
532 !! \param x2 X-value of end point
533 !! \param y2 Y-value of end point
534 
535 subroutine pgarro(x1,y1, x2,y2)
536  use plplot, only: plflt, plline, plpoin
537 
538  implicit none
539  real, intent(in) :: x1,x2, y1,y2
540  real(kind=plflt) :: x(2),y(2)
541 
542  x = (/x1,x2/)
543  y = (/y1,y2/)
544 
545  call plline(x,y)
546  call plpoin((/x(2)/),(/y(2)/),2)
547 
548  !print*,'pgarro: ',x1,x2,y1,y2,' ',x,y
549 
550 end subroutine pgarro
551 !***********************************************************************************************************************************
552 
553 
554 
555 !***********************************************************************************************************************************
556 !> \brief Draw points
557 !!
558 !! \param n Number of points
559 !! \param x1 X-values of points
560 !! \param y1 Y-values of points
561 !! \param s Plot symbol to use
562 
563 subroutine pgpoint(n,x1,y1,s)
564  use plplot, only: plflt, plpoin
565 
566  implicit none
567  integer, intent(in) :: n,s
568  real, intent(in) :: x1(n),y1(n)
569  real(kind=plflt) :: x2(n),y2(n)
570 
571  x2 = x1
572  y2 = y1
573 
574  call plpoin(x2,y2,s)
575  !call plsym(1,x2,y2,143) ! Produces Hershey -> many letters, etc
576 
577 end subroutine pgpoint
578 !***********************************************************************************************************************************
579 
580 !***********************************************************************************************************************************
581 !> \brief Draw a polygone
582 !!
583 !! \param n Number of points
584 !! \param x1 X-values of points
585 !! \param y1 Y-values of points
586 
587 subroutine pgpoly(n,x1,y1)
588  use plplot, only: plflt, plfill
589 
590  implicit none
591  integer, intent(in) :: n
592  real, intent(in) :: x1(n),y1(n)
593  real(kind=plflt) :: x2(n),y2(n)
594 
595  x2 = x1
596  y2 = y1
597 
598  call plfill(x2,y2)
599 
600  !print*,'pgpoly: ',n,x1(1),x1(n),y1(1),y1(n)
601 
602 end subroutine pgpoly
603 !***********************************************************************************************************************************
604 
605 !***********************************************************************************************************************************
606 !> \brief Draw a rectangle
607 !!
608 !! \param x1 Lower x value
609 !! \param x2 Upper x value
610 !! \param y1 Lower y value
611 !! \param y2 Upper y value
612 
613 subroutine pgrect(x1,x2,y1,y2)
614  use plplot, only: plflt, plfill
615 
616  implicit none
617  real, intent(in) :: x1,x2,y1,y2
618  real(kind=plflt) :: x(4),y(4)
619 
620  x = (/x1,x1,x2,x2/)
621  y = (/y1,y2,y2,y1/)
622 
623  call plfill(x,y)
624 
625 end subroutine pgrect
626 !***********************************************************************************************************************************
627 
628 
629 !***********************************************************************************************************************************
630 !> \brief Draw a circle
631 !!
632 !! \param xc X-value of centre
633 !! \param yc Y-value of centre
634 !! \param r Radius
635 
636 subroutine pgcirc(xc, yc, r)
637  use plplot, only: plflt, plfill
638 
639  implicit none
640  real, intent(in) :: xc, yc, r
641 
642  integer, parameter :: n=100
643  integer :: i
644  real(kind=plflt) :: x(n),y(n),twopi
645 
646  twopi = 8.*atan(1.)
647 
648  do i=1,n
649  x(i) = xc + r * cos(twopi/real(n-1))
650  y(i) = yc + r * sin(twopi/real(n-1))
651  end do
652 
653  call plfill(x,y)
654 
655 end subroutine pgcirc
656 !***********************************************************************************************************************************
657 
658 
659 
660 
661 !***********************************************************************************************************************************
662 !> \brief Make a contour plot
663 !!
664 !! \param arr Data array
665 !! \param nx Dimension 1 of data array
666 !! \param ny Dimension 2 of data array
667 !! \param ix1 Start index range in dimension 1 to use from data array
668 !! \param ix2 End index range in dimension 1 to use from data array
669 !! \param iy1 Start index range in dimension 2 to use from data array
670 !! \param iy2 End index range in dimension 2 to use from data array
671 !! \param c Array with heights to draw contours for
672 !! \param nc Dimension of c
673 !! \param tr Affine transformation elements
674 
675 subroutine pgcont(arr, nx,ny, ix1,ix2, iy1,iy2, c, nc, tr)
676  use plplot, only: plflt, plcont
677 
678  implicit none
679  integer, intent(in) :: nx,ny, ix1,ix2, iy1,iy2, nc
680  real, intent(in) :: arr(nx,ny), c(*), tr(6)
681  real(kind=plflt) :: arr1(nx,ny), clevel(nc), tr1(6)
682 
683  arr1 = arr
684  clevel = c(1:nc)
685  tr1 = tr
686 
687  call plcont(arr1, ix1,ix2, iy1,iy2, clevel, tr1)
688 
689 end subroutine pgcont
690 !***********************************************************************************************************************************
691 
692 
693 !***********************************************************************************************************************************
694 !> \brief Shade a region (between contours/heights) - dummy routine!
695 !!
696 !! \param arr Data array
697 !! \param nx Dimension 1 of data array
698 !! \param ny Dimension 2 of data array
699 !! \param ix1 Start index range in dimension 1 to use from data array
700 !! \param ix2 End index range in dimension 1 to use from data array
701 !! \param iy1 Start index range in dimension 2 to use from data array
702 !! \param iy2 End index range in dimension 2 to use from data array
703 !! \param c1 Lower limit in height/contour to fill
704 !! \param c2 Upper limit in height/contour to fill
705 !! \param tr Affine transformation elements
706 
707 subroutine pgconf(arr, nx,ny, ix1,ix2, iy1,iy2, c1, c2, tr)
709  implicit none
710  integer, intent(in) :: nx,ny, ix1,ix2, iy1,iy2
711  real, intent(in) :: arr(nx,ny), c1,c2, tr(6)
712  !real(kind=plflt) :: arr1(nx,ny), clevel(nc), tr1(6)
713  integer :: tmp
714  integer, save :: warn
715 
716  tmp = nint(arr(1,1))
717  tmp = nx
718  tmp = ny
719  tmp = ix1
720  tmp = ix2
721  tmp = iy1
722  tmp = iy2
723  tmp = nint(c1)
724  tmp = nint(c2)
725  tmp = nint(tr(1))
726  tmp = tmp ! Avoid 'variable is set but not used' warnings from compiler for dummy variable
727 
728 
729  !arr1 = arr
730  !clevel = c(1:nc)
731  !tr1 = tr
732  !
733  !call plshade1(arr1, ix1,ix2, iy1,iy2, clevel, tr1)
734 
735  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
736  if(warn.ne.123) call warn_dummy_routine('pgconf', '')
737  warn = 123
738 
739 end subroutine pgconf
740 !***********************************************************************************************************************************
741 
742 
743 
744 
745 
746 
747 
748 
749 
750 
751 
752 
753 
754 !***********************************************************************************************************************************
755 !> \brief Print text with arbitrary angle and justification
756 !!
757 !! \param x1 X-coordinate of text
758 !! \param y1 Y-coordinate of text
759 !! \param ang Angle
760 !! \param just1 Justification
761 !! \param text Text to print
762 !!
763 !! \note Angle only correct for 0,90,180,270deg or square viewport
764 
765 subroutine pgptxt(x1,y1,ang,just1,text)
766  use plplot, only: plflt, plptex
767 
768  implicit none
769  real, intent(in) :: x1,y1,ang,just1
770  character, intent(in) :: text*(*)
771  real :: d2r
772  real(kind=plflt) :: x2,y2,just2,dx,dy,xmin,xmax,ymin,ymax
773  character :: text1*(len(text))
774 
775  d2r = atan(1.)/45.
776  call plgvpw(xmin,xmax,ymin,ymax)
777 
778  ! Convert angle -> dy/dx
779  dx = (xmax-xmin)*0.1
780 
781  if(abs(ang).lt.1.e-5) then ! ang ~ 0 deg
782  dy = 0.0
783  else if(abs(mod(ang-90.,180.)).lt.1.e-5) then ! ang = +/-90deg
784  dx = 0.0
785  dy = -1.0 ! ang = -90deg
786  if(abs(mod(ang-90.,360.)).lt.1.e-5) dy = 1. ! ang = +90deg
787  else
788  dx = 1.0
789  dy = dx*tan(ang*d2r) * (ymax-ymin)/(xmax-xmin)
790  if(ang.gt.90. .and. ang.lt.270. .or. ang.lt.-90. .and. ang.gt.-270.) then
791  dx = -dx
792  dy = -dy
793  end if
794  end if
795 
796  x2 = x1
797  y2 = y1
798  just2 = just1
799 
800  text1 = text
801  call pg2pltext(text1)
802 
803  call plptex(x2,y2,dx,dy,just2,trim(text1))
804 
805  !write(6,'(A,4F10.3,A)')' pgptxt: ',x1,y1,ang,just1,trim(text)
806  !write(6,'(A,5F10.3,A)')' pgptxt: ',x2,y2,dx,dy,just2,trim(text1)
807 
808 end subroutine pgptxt
809 !***********************************************************************************************************************************
810 
811 !***********************************************************************************************************************************
812 !> \brief Non-standard alias for pgptxt()
813 !!
814 !! \param x X-coordinate of text
815 !! \param y Y-coordinate of text
816 !! \param ang Angle
817 !! \param just Justification
818 !! \param text Text to print
819 !!
820 !! \note Angle only right for 0,90,180,270deg or square viewport
821 
822 subroutine pgptext(x,y,ang,just,text)
823  implicit none
824  real, intent(in) :: x,y,ang,just
825  character, intent(in) :: text*(*)
826 
827  call pgptxt(x,y,ang,just,text)
828 
829 end subroutine pgptext
830 !***********************************************************************************************************************************
831 
832 
833 !***********************************************************************************************************************************
834 !> \brief Print text with default angle and justification
835 !!
836 !! \param x1 X-coordinate of text
837 !! \param y1 Y-coordinate of text
838 !! \param text Text to print
839 
840 subroutine pgtext(x1,y1,text)
841  use plplot, only: plflt, plptex
842 
843  implicit none
844  real, intent(in) :: x1,y1
845  character, intent(in) :: text*(*)
846  real(kind=plflt) :: x2,y2,just,dx,dy
847  character :: text1*(len(text))
848 
849  ! Convert angle=0deg -> dy/dx
850  dx = 1.
851  dy = 0.
852  just = 0. !Left-adjusted
853 
854  x2 = x1
855  y2 = y1
856 
857  text1 = text
858  call pg2pltext(text1)
859  call plptex(x2,y2,dx,dy,just,text1)
860 
861 end subroutine pgtext
862 !***********************************************************************************************************************************
863 
864 
865 !***********************************************************************************************************************************
866 !> \brief Print text in the margin
867 !!
868 !! \param side Side to print text on ('L','R','T','B')
869 !! \param disp1 Distance from axis
870 !! \param pos1 Position along axis
871 !! \param just1 Justification
872 !! \param text Text to print
873 
874 subroutine pgmtxt(side, disp1, pos1, just1, text)
875  use plplot, only: plflt, plmtex
876 
877  implicit none
878  real, intent(in) :: disp1,pos1,just1
879  real(kind=plflt) :: disp2,pos2,just2
880  character, intent(in) :: side*(*)
881  character, intent(in) :: text*(*)
882  character :: text1*(len(text))
883 
884  disp2 = disp1
885  pos2 = pos1
886  just2 = just1
887 
888  !write(6,'(2A,2(3F10.3,5x),A)')' pgmtxt: ',trim(side),disp1,pos1,just1,disp2,pos2,just2,trim(text)
889 
890  text1 = text
891  call pg2pltext(text1)
892  call plmtex(side, disp2, pos2, just2, text1)
893 
894 end subroutine pgmtxt
895 !***********************************************************************************************************************************
896 
897 !***********************************************************************************************************************************
898 !> \brief Alias for pgmtxt()
899 !!
900 !! \param side Side to print text on ('L','R','T','B')
901 !! \param disp
902 !! \param pos Position
903 !! \param just Justification
904 !! \param text Text to print
905 
906 subroutine pgmtext(side, disp, pos, just, text)
907  implicit none
908  real, intent(in) :: disp,pos,just
909  character, intent(in) :: side*(*)
910  character, intent(in) :: text*(*)
911  character :: text1*(len(text))
912 
913  text1 = text
914  call pg2pltext(text1)
915  call pgmtxt(side, disp, pos, just, text1)
916 
917 end subroutine pgmtext
918 !***********************************************************************************************************************************
919 
920 
921 !***********************************************************************************************************************************
922 !> \brief Interface for pglab
923 !!
924 !! \param xlbl Label for the horizontal axis
925 !! \param ylbl Label for the vertical axis
926 !! \param toplbl Plot title
927 
928 subroutine pglab(xlbl, ylbl, toplbl)
929  implicit none
930  character, intent(in) :: xlbl*(*), ylbl*(*), toplbl*(*)
931  call pgbbuf
932  call pgmtxt('T', 2.0, 0.5, 0.5, toplbl)
933  call pgmtxt('B', 3.2, 0.5, 0.5, xlbl)
934  call pgmtxt('L', 2.2, 0.5, 0.5, ylbl)
935  call pgebuf
936 end subroutine pglab
937 !***********************************************************************************************************************************
938 
939 !***********************************************************************************************************************************
940 !> \brief Alias for pglab
941 !!
942 !! \param xlbl Label for the horizontal axis
943 !! \param ylbl Label for the vertical axis
944 !! \param toplbl Plot title
945 
946 subroutine pglabel(xlbl, ylbl, toplbl)
947  implicit none
948  character, intent(in) :: xlbl*(*), ylbl*(*), toplbl*(*)
949  call pglab(xlbl, ylbl, toplbl)
950 end subroutine pglabel
951 !***********************************************************************************************************************************
952 
953 !***********************************************************************************************************************************
954 !> \brief Start a new plot
955 !!
956 !! \param pgdev PGplot device
957 !!
958 !! This function creates a new stream for each xwin request, but uses the same stream (0) for each file stream request.
959 !! Xwin streams are treated differently, since we want bufferring for smooth animations whereas the other streams are non-buffered.
960 
961 function pgopen(pgdev)
962  use plplot, only: plspause, plsfnam, plsdev
963  use plplot, only: plmkstrm, plsetopt
964  use pg2plplot, only : is_init, devid
965  implicit none
966  integer :: pgopen
967  character, intent(in) :: pgdev*(*)
968 
969  integer :: cur_stream=0, check_error
970  character :: pldev*(99),filename*(99)
971 
972  filename = 'plot_temp.png'
973 
974  call pg2pldev(pgdev, pldev,filename) ! Extract pldev and filename from pgdev
975  call plmkstrm(cur_stream)
976  call plssub(1, 1)
977  call plsdev(trim(pldev))
978  if(trim(pldev).ne.'xwin') then
979  if(check_error(trim(filename)).ne. 0) then
980  pgopen = -1
981  return
982  end if
983  call plsfnam(trim(filename)) ! Set output file name
984  else
985  call plsetopt("db", "")
986  end if
987 
988  !call plscolbg(255,255,255) ! Set background colour to white
989  call plfontld(1) ! Load extended character set(?)
990  call plspause(.false.) ! Pause at plend()
991 
992  pgopen = cur_stream + 1
993  devid = pgopen
994  is_init = .false.
995 end function pgopen
996 !***********************************************************************************************************************************
997 
998 !***********************************************************************************************************************************
999 !> \brief Inquire current device identifier
1000 !!
1001 !! \retval Device identifier
1002 
1003 subroutine pgqid(id)
1004  use pg2plplot, only : devid
1005  implicit none
1006  integer, intent(out) :: id
1007  id = devid
1008 end subroutine pgqid
1009 !***********************************************************************************************************************************
1010 
1011 
1012 !***********************************************************************************************************************************
1013 !> \brief Begin a new plot
1014 !!
1015 !! \param i Display ID
1016 !! \param pgdev PGplot device
1017 !! \param nx Number of frames in the x-direction
1018 !! \param ny Number of frames in the y-direction
1019 !!
1020 !! This is a bit simpler than pgopen(), since I don't need to create a new stream for each request.
1021 
1022 subroutine pgbegin(i,pgdev,nx,ny)
1023  use plplot, only: plspause, plsfnam, plsetopt, plssub, plsdev
1024  use pg2plplot, only : is_init
1025  implicit none
1026  integer, intent(in) :: i,nx,ny
1027  character, intent(in) :: pgdev*(*)
1028  integer :: i1, check_error
1029  character :: pldev*(99),filename*(99)
1030 
1031  ! These are ignored by pgbegin. Can't be self, since calling arguments are likely constants
1032  i1 = i
1033  i1 = nx
1034  i1 = ny
1035  i1 = i1
1036 
1037  call pg2pldev(pgdev, pldev,filename) ! Extract pldev and filename from pgdev
1038 
1039  if(trim(pldev).ne.'xwin') then
1040  if(check_error(trim(filename)).ne.0) return
1041  call plsfnam(trim(filename)) ! Set output file name
1042  else
1043  call plsetopt("db", "")
1044  end if
1045 
1046  call plfontld(1) ! Load extended character set(?)
1047  call plssub(1, 1)
1048  call plsdev(trim(pldev))
1049  is_init = .false.
1050  call plspause(.false.) ! Pause at plend()
1051 end subroutine pgbegin
1052 !***********************************************************************************************************************************
1053 
1054 !***********************************************************************************************************************************
1055 !> \brief Alias for pgbegin
1056 !!
1057 !! \param i Display ID
1058 !! \param pgdev PGplot device
1059 !! \param nx Number of frames in the x-direction
1060 !! \param ny Number of frames in the y-direction
1061 
1062 subroutine pgbeg(i, pgdev, nx, ny)
1063  implicit none
1064  integer, intent(in) :: i,nx,ny
1065  character, intent(in) :: pgdev*(*)
1066  call pgbegin(i, pgdev, nx, ny)
1067 end subroutine pgbeg
1068 !***********************************************************************************************************************************
1069 
1070 
1071 !***********************************************************************************************************************************
1072 !> \brief End a plot
1073 !!
1074 
1075 subroutine pgend()
1076  implicit none
1077 
1078  call plflush()
1079  call plend1()
1080 end subroutine pgend
1081 !***********************************************************************************************************************************
1082 
1083 
1084 !***********************************************************************************************************************************
1085 !> \brief Set paper size
1086 !!
1087 !! \param width Paper width
1088 !! \param ratio Paper aspect ratio
1089 
1090 subroutine pgpap(width,ratio)
1092  use plplot, only : plflt
1093  implicit none
1094  real, intent(in) :: width, ratio
1095  integer :: xlen,ylen,xoff,yoff
1096  real(kind=plflt) :: xp,yp
1097 
1098  xp = 300. ! DPI
1099  yp = 300.
1100  paper_width = width
1101  paper_ratio = ratio
1102 
1103  xlen = nint(paper_width*xp)
1104  ylen = nint(paper_width*xp*paper_ratio)
1105 
1106  xoff = 0 ! Offset
1107  yoff = 0
1108 
1109  call plspage(xp,yp,xlen,ylen,xoff,yoff) ! Must be called before plinit()!
1110  call do_init()
1111 
1112 end subroutine pgpap
1113 !***********************************************************************************************************************************
1114 
1115 
1116 !***********************************************************************************************************************************
1117 !> \brief Set view port
1118 !!
1119 !! \param xl1 Left side of the x-axis
1120 !! \param xr1 Right side of the x-axis
1121 !! \param yb1 Bottom side of the y-axis
1122 !! \param yt1 Top side of the y-axis
1123 
1124 subroutine pgsvp(xl1,xr1,yb1,yt1)
1125  use plplot, only: plflt
1126 
1127  implicit none
1128  real, intent(in) :: xl1,xr1,yb1,yt1
1129  real(kind=plflt) :: xl2,xr2,yb2,yt2
1130 
1131  xl2 = xl1
1132  xr2 = xr1
1133  yb2 = yb1
1134  yt2 = yt1
1135  !write(6,'(A,2(4F10.3,5x))')' pgsvp: ',xl1,xr1,yb1,yt1,xl2,xr2,yb2,yt2
1136 
1137  call plvpor(xl2,xr2,yb2,yt2)
1138 
1139 end subroutine pgsvp
1140 !***********************************************************************************************************************************
1141 
1142 
1143 !***********************************************************************************************************************************
1144 !> \brief Set window
1145 !!
1146 !! \param xmin1 Left
1147 !! \param xmax1 Right
1148 !! \param ymin1 Top
1149 !! \param ymax1 Bottom
1150 
1151 subroutine pgswin(xmin1,xmax1,ymin1,ymax1)
1152  use plplot, only: plflt
1153 
1154  implicit none
1155  real, intent(in) :: xmin1,xmax1,ymin1,ymax1
1156  real(kind=plflt) :: xmin2,xmax2,ymin2,ymax2
1157 
1158  xmin2 = xmin1
1159  xmax2 = xmax1
1160  ymin2 = ymin1
1161  ymax2 = ymax1
1162  !write(6,'(A,2(4F10.3,5x))')' pgswin: ',xmin1,xmax1,ymin1,ymax1,xmin2,xmax2,ymin2,ymax2
1163 
1164  call plwind(xmin2,xmax2,ymin2,ymax2)
1165 
1166 end subroutine pgswin
1167 !***********************************************************************************************************************************
1168 
1169 !***********************************************************************************************************************************
1170 !> \brief alias for pgswin
1171 !!
1172 !! \param xmin1 Left
1173 !! \param xmax1 Right
1174 !! \param ymin1 Top
1175 !! \param ymax1 Bottom
1176 
1177 subroutine pgwindow(xmin1,xmax1,ymin1,ymax1)
1178  use plplot, only: plflt
1179 
1180  implicit none
1181  real, intent(in) :: xmin1,xmax1,ymin1,ymax1
1182  real(kind=plflt) :: xmin2,xmax2,ymin2,ymax2
1183 
1184  xmin2 = dble(xmin1)
1185  xmax2 = dble(xmax1)
1186  ymin2 = dble(ymin1)
1187  ymax2 = dble(ymax1)
1188  call plwind(xmin2,xmax2,ymin2,ymax2)
1189 
1190 end subroutine pgwindow
1191 !***********************************************************************************************************************************
1192 
1193 
1194 !***********************************************************************************************************************************
1195 !> \brief Subdivide view surface into panels
1196 !!
1197 !! \param nxsub Number of subticks on the x-axis
1198 !! \param nysub Number of subticks on the y-axis
1199 
1200 subroutine pgsubp(nxsub, nysub)
1201  implicit none
1202  integer, intent(in) :: nxsub,nysub
1203 
1204  call plssub(nxsub, nysub)
1205 
1206 end subroutine pgsubp
1207 !***********************************************************************************************************************************
1208 
1209 
1210 !***********************************************************************************************************************************
1211 !> \brief Advance to the next (sub-)page
1212 !!
1213 
1214 subroutine pgpage()
1215  use pg2plplot, only: do_init
1216  implicit none
1217 
1218  call pladv(0)
1219  call do_init()
1220 
1221 end subroutine pgpage
1222 !***********************************************************************************************************************************
1223 
1224 !***********************************************************************************************************************************
1225 !> Alias for pgpage
1226 subroutine pgadvance()
1227  call pgpage()
1228 end subroutine pgadvance
1229 !***********************************************************************************************************************************
1230 
1231 
1232 !***********************************************************************************************************************************
1233 !> \brief Start buffering output - dummy routine!
1234 !!
1235 
1236 subroutine pgbbuf()
1238  implicit none
1239  integer, save :: warn
1240 
1241  call do_init()
1242 
1243  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
1244  if(warn.ne.123) call warn_dummy_routine('pgbbuf', '')
1245  warn = 123
1246 
1247 end subroutine pgbbuf
1248 !***********************************************************************************************************************************
1249 
1250 
1251 !***********************************************************************************************************************************
1252 !> \brief End buffering output - dummy routine!
1253 !!
1254 
1255 subroutine pgebuf()
1257  implicit none
1258  integer, save :: warn
1259 
1260  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
1261  if(warn.ne.123) call warn_dummy_routine('pgebuf', '')
1262  warn = 123
1263 
1264 end subroutine pgebuf
1265 !***********************************************************************************************************************************
1266 
1267 
1268 
1269 
1270 !***********************************************************************************************************************************
1271 !> \brief Draw a box (+axes) around a plot
1272 !!
1273 !! \param xopt Options for the x-axis
1274 !! \param xtick1 Distance between ticks on the x-axis
1275 !! \param nxsub Number of subticks on the x-axis
1276 !! \param yopt Options for the y-axis
1277 !! \param ytick1 Distance between ticks on the y-axis
1278 !! \param nysub Number of subticks on the y-axis
1279 
1280 subroutine pgbox(xopt, xtick1, nxsub, yopt, ytick1, nysub)
1281  use plplot, only: plflt, plbox
1282 
1283  implicit none
1284  integer, intent(in) :: nxsub,nysub
1285  real, intent(in) :: xtick1,ytick1
1286  character, intent(in) :: xopt*(*),yopt*(*)
1287  real(kind=plflt) :: xtick2,ytick2
1288 
1289  xtick2 = xtick1
1290  ytick2 = ytick1
1291  !write(6,'(A,2(2F10.3,5x))')' pgbox: ',xtick1,ytick1,xtick2,ytick2
1292 
1293  call plbox(xopt, xtick2, nxsub, yopt, ytick2, nysub)
1294 
1295 end subroutine pgbox
1296 !***********************************************************************************************************************************
1297 
1298 
1299 
1300 !***********************************************************************************************************************************
1301 !> \brief Draw a single tick mark, optionally with label - no PLplot routine found, using an ad-hoc routine
1302 !!
1303 !! \param x1 world coordinates of one endpoint of the axis
1304 !! \param y1 world coordinates of one endpoint of the axis
1305 !! \param x2 world coordinates of the other endpoint of the axis
1306 !! \param y2 world coordinates of the other endpoint of the axis
1307 !!
1308 !! \param pos draw the tick mark at fraction pos (0<=pos<=1) along the line from (X1,Y1) to (X2,Y2)
1309 !! \param tikl length of tick mark drawn to the left or bottom of the axis - drawing ticks outside the box may not be supported
1310 !! \param tikr length of major tick marks drawn to the right or top of the axis - drawing outside the box may not be supported
1311 !!
1312 !! \param disp displacement of label text from the axis
1313 !! \param orient orientation of label text, in degrees
1314 !! \param lbl text of label (may be blank)
1315 
1316 subroutine pgtick(x1, y1, x2, y2, pos, tikl, tikr, disp, orient, lbl)
1317  use plplot, only: plflt, plmtex
1318 
1319  implicit none
1320  real, intent(in) :: x1, y1, x2, y2, pos, tikl, tikr, disp, orient
1321  character, intent(in) :: lbl*(*)
1322 
1323  real :: x,y,dx,dy, dpx,dpy, reldiff
1324  real(plflt) :: p_xmin,p_xmax,p_ymin,p_ymax, plx1,plx2,ply1,ply2, tlen, disp1,pos1,just
1325  character :: lbl1*(len(lbl))
1326  logical :: seq,deq
1327 
1328  disp1 = disp
1329  x = orient ! Unused
1330  x = x ! Avoid 'variable is set but not used' warnings from compiler for dummy variable
1331  lbl1 = lbl
1332 
1333  dx = x2 - x1
1334  dy = y2 - y1
1335  x = x1 + pos*dx
1336  y = y1 + pos*dy
1337 
1338  call plgvpw(p_xmin, p_xmax, p_ymin, p_ymax)
1339 
1340  dpx = real(abs(p_xmax-p_xmin))
1341  dpy = real(abs(p_ymax-p_ymin))
1342 
1343 
1344  tlen = 0.03 ! default size of a tickmark
1345  pos1 = pos
1346  just = 0.5
1347 
1348  if(reldiff(y1,y2) .le. epsilon(y1)*10 .or. (seq(y1,0.) .and. seq(y2,0.))) then ! Want horizontal axis
1349 
1350  plx1 = x
1351  plx2 = x
1352  ply1 = y
1353 
1354  ! Plot ticks below the axis:
1355  if(tikl.gt.tiny(tikl)) then
1356  ply2 = y - dpy * tikl * tlen
1357  call pljoin(plx1,ply1,plx2,ply2)
1358  end if
1359 
1360  ! Plot ticks above the axis:
1361  if(tikr.gt.tiny(tikr)) then
1362  ply2 = y + dpy * tikr * tlen
1363  call pljoin(plx1,ply1,plx2,ply2)
1364  end if
1365 
1366  ! Print labels:
1367  if(abs(reldiff(y,real(p_ymin))) .le. epsilon(y)*10 .or. (seq(y,0.) .and. deq(p_ymin,0.0_plflt))) then
1368  call plmtex('B', disp1, pos1, just, trim(lbl1)) ! Print label below the bottom axis
1369  else
1370  call plmtex('T', disp1, pos1, just, trim(lbl1)) ! Print label above the top axis
1371  end if
1372 
1373  else if(reldiff(x1,x2) .le. epsilon(x1)*10 .or. (seq(x1,0.) .and. seq(x2,0.))) then ! Want vertical axis
1374 
1375  ply1 = y
1376  ply2 = y
1377  plx1 = x
1378 
1379  ! Plot ticks to the left of the axis:
1380  if(tikl.gt.tiny(tikl)) then
1381  plx2 = x - dpx * tikl * tlen
1382  call pljoin(plx1,ply1, plx2,ply2)
1383  end if
1384 
1385  ! Plot ticks to the right of the axis:
1386  if(tikr.gt.tiny(tikr)) then
1387  plx2 = x + dpx * tikr * tlen
1388  call pljoin(plx1,ply1, plx2,ply2)
1389  end if
1390 
1391  ! Print labels:
1392  if(abs(reldiff(x,real(p_xmin))) .le. epsilon(x)*10 .or. (seq(x,0.) .and. deq(p_xmin,0.0_plflt))) then
1393  call plmtex('L', disp1, pos1, just, trim(lbl1)) ! Print label to the left of the left-hand axis
1394  else
1395  call plmtex('R', disp1, pos1, just, trim(lbl1)) ! Print label to the right of the right-hand axis
1396  end if
1397 
1398  else
1399  write(0, '(A)') "PG2PLplot, pgtick(): x1!=x2 and y1!=y2 - I don't know which axis you want to use..."
1400  print*,'x:',x1,x2,reldiff(x1,x2)
1401  print*,'y:',y1,y2,reldiff(y1,y2)
1402  end if
1403 
1404 end subroutine pgtick
1405 !***********************************************************************************************************************************
1406 
1407 
1408 
1409 !***********************************************************************************************************************************
1410 !> \brief Read data from screen - no PLplot equivalent (yet) - dummy routine!
1411 !!
1412 !! \param maxpt Maximum number of points that may be accepted
1413 !! \param npt Number of points entered (zero on first call)
1414 !! \param x Array of x-coordinates
1415 !! \param y Array of y-coordinates
1416 !! \param symbol Code number of symbol to use for markingentered points
1417 !!
1418 !! \todo No plplot routine found yet - using dummy
1419 
1420 subroutine pgolin(maxpt, npt, x, y, symbol)
1422 
1423  implicit none
1424  integer, intent(in) :: maxpt,symbol
1425  integer, intent(out) :: npt
1426  real, intent(out) :: x(maxpt),y(maxpt)
1427 
1428  integer :: symbol1
1429  integer, save :: warn
1430 
1431  npt = maxpt
1432  x = 0.d0
1433  y = 0.d0
1434  symbol1 = symbol
1435  symbol1 = symbol1 ! Avoid 'variable is set but not used' warnings from compiler for dummy variable
1436 
1437  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
1438  if(warn.ne.123) call warn_dummy_routine('pgolin', '')
1439  warn = 123
1440 
1441 end subroutine pgolin
1442 !***********************************************************************************************************************************
1443 
1444 
1445 !***********************************************************************************************************************************
1446 !> \brief Erase screen
1447 !!
1448 
1449 subroutine pgeras()
1450  use plplot, only: plclear
1451  implicit none
1452  call plclear()
1453 end subroutine pgeras
1454 !***********************************************************************************************************************************
1455 
1456 
1457 
1458 !***********************************************************************************************************************************
1459 !> \brief Replace the PGPlot escape character '\' with the PLplot escape character '#'
1460 !!
1461 !! \param string Text string to convert
1462 
1463 subroutine pg2pltext(string)
1464  implicit none
1465  character, intent(inout) :: string*(*)
1466 
1467  !print*,trim(string)
1468 
1469  call replace_substring(string, '\', '#') ! Replace the PGPlot escape character \ with the PLplot escape character # '
1470  call replace_substring(string, '#(0248)', '#(2246)') ! \approx -> \sim
1471  call replace_substring(string, '#(062', '#(212') ! alpha - gamma
1472  call replace_substring(string, '#(063', '#(213') ! delta - nu
1473  call replace_substring(string, '#(064', '#(214') ! xi - psi
1474  call replace_substring(string, '#(0650', '#(2150') ! omega
1475  call replace_substring(string, '#(0685)', '#(2185)') ! (var)theta
1476 
1477  !print*,trim(string)
1478 
1479 end subroutine pg2pltext
1480 !***********************************************************************************************************************************
1481 
1482 !***********************************************************************************************************************************
1483 !> \brief Converts PGplot device ID to PLplot device ID + file name
1484 !!
1485 !! \param pgdev PGplot device ID (includes filename: 'file.name/dev')
1486 !! \retval pldev PLplot device ID ('dev')
1487 !! \retval filename PLplot file name ('file.name')
1488 !!
1489 !!
1490 !! \note Possible values for pldev:
1491 !! - xwin:
1492 !! - wxwidgets:
1493 !! - xcairo:
1494 !! - qtwidget:
1495 !!
1496 !! - png: no anti-aliasing in lines
1497 !! - wxpng: no extended characters
1498 !! - pngcairo:
1499 !! - pngqt:
1500 
1501 subroutine pg2pldev(pgdev, pldev,filename)
1502  implicit none
1503  character, intent(in) :: pgdev*(*)
1504  character, intent(out) :: pldev*(*), filename*(*)
1505  integer :: i
1506 
1507 
1508  i = index(pgdev, '/', back=.true.)
1509 
1510  filename = ' '
1511  if(i.ne.1) filename = pgdev(1:i-1)
1512 
1513  pldev = ' '
1514  pldev = pgdev(i+1:)
1515 
1516  ! Change PGPlot ps to PLplot ps:
1517  call replace_substring(pldev,'vcps','psc')
1518  call replace_substring(pldev,'cvps','psc')
1519  call replace_substring(pldev,'cps','psc')
1520  call replace_substring(pldev,'vps','ps')
1521 
1522  ! Change PGplot ppm output tot png:
1523  call replace_substring(pldev,'ppm','png')
1524  call replace_substring(filename,'ppm','png')
1525 
1526  ! Use pngqt rather than png:
1527  !call replace_substring(pldev,'png','pngqt')
1528 
1529  ! Use pngcairo rather than png, since I get segfaults if outputing to both X11 and png which pulls in pngqt (joequant):
1530  if(pldev .eq. 'png') pldev = 'pngcairo'
1531 
1532 end subroutine pg2pldev
1533 !***********************************************************************************************************************************
1534 
1535 
1536 !***********************************************************************************************************************************
1537 !> \brief Select output stream
1538 !!
1539 !! \param pgdev stream number
1540 
1541 subroutine pgslct(pgdev)
1542  use pg2plplot, only: do_init
1543  use plplot, only : plspause, plgdev
1544 
1545  implicit none
1546  integer, intent(in) :: pgdev
1547  character :: pldev*(99)
1548 
1549  call do_init()
1550  call plgdev(pldev)
1551  if(trim(pldev).eq.'xwin') call pleop()
1552  call plsstrm(pgdev-1)
1553 
1554 end subroutine pgslct
1555 !***********************************************************************************************************************************
1556 
1557 
1558 !***********************************************************************************************************************************
1559 !> \brief Save parameters
1560 !!
1561 
1562 subroutine pgsave()
1566  use pg2plplot, only: cur_lstyle, cur_lwidth
1567  use plplot, only : plgfont
1568  implicit none
1569 
1570  save_level = save_level + 1
1571  if(save_level.gt.max_level) then
1572  write(0,'(/,A,/)') '*** PG2PLplot WARNING: too many save calls in pgsave()'
1573  return
1574  end if
1575 
1576  call plgfont(save_ffamily(save_level), save_fstyle(save_level), save_fweight(save_level))
1577  save_lwidth(save_level) = cur_lwidth
1578  save_color(save_level) = cur_color
1579  save_lstyle(save_level) = cur_lstyle
1580 end subroutine pgsave
1581 !***********************************************************************************************************************************
1582 
1583 
1584 !***********************************************************************************************************************************
1585 !> \brief Unsave parameters
1586 !!
1587 
1588 subroutine pgunsa()
1592  use plplot, only : plsfont
1593  implicit none
1594 
1595  if(save_level.eq.0) then
1596  write(0,'(/,A,/)') '*** PG2PLplot WARNING: no save call in stack in pgunsa()'
1597  return
1598  end if
1599 
1600  if(save_level.gt.max_level) then
1601  write(0,'(/,A,/)') '*** PG2PLplot WARNING: unsave greater than stack in pgunsa()'
1602  save_level = save_level - 1
1603  return
1604  end if
1605 
1606  call plsfont(save_ffamily(save_level), save_fstyle(save_level), save_fweight(save_level))
1607  call pgslw(save_lwidth(save_level))
1608  call pgsci(save_color(save_level))
1609  call pgsls(save_lstyle(save_level))
1610  save_level = save_level - 1
1611 
1612 end subroutine pgunsa
1613 !***********************************************************************************************************************************
1614 
1615 
1616 !***********************************************************************************************************************************
1617 !> \brief Move pen to location - for use with pgdraw() only!
1618 !!
1619 !! \param x Horizontal location
1620 !! \param y Vertical location
1621 
1622 subroutine pgmove(x, y)
1623  use pg2plplot, only: xcur, ycur
1624  implicit none
1625  real, intent(in) :: x, y
1626 
1627  xcur = x
1628  ycur = y
1629 
1630 end subroutine pgmove
1631 !***********************************************************************************************************************************
1632 
1633 
1634 !***********************************************************************************************************************************
1635 !> \brief Draw line to location
1636 !!
1637 !! \param x Horizontal location
1638 !! \param y Vertical location
1639 
1640 subroutine pgdraw(x, y)
1641  use pg2plplot, only: xcur, ycur
1642  use plplot, only: pljoin, plflt
1643  implicit none
1644  real, intent(in) :: x, y
1645  real(plflt) :: x1, y1
1646 
1647  x1 = x
1648  y1 = y
1649 
1650  call pljoin(xcur, ycur, x1, y1)
1651  call pgmove(x, y)
1652 
1653 end subroutine pgdraw
1654 !***********************************************************************************************************************************
1655 
1656 
1657 !***********************************************************************************************************************************
1658 !> \brief Draw a point
1659 !!
1660 !! \param n Number of data points to draw
1661 !! \param x Horizontal locations
1662 !! \param y Vertical locations
1663 !! \param ncode Plot symbol
1664 
1665 subroutine pgpt(n, x, y, ncode)
1666  use plplot, only: plpoin, plflt, plsym
1667 
1668  implicit none
1669  integer, intent(in) :: n, ncode
1670  real(plflt), intent(in), dimension(n) :: x, y
1671 
1672  integer :: code
1673 
1674 
1675  if(ncode == -3) then
1676  code = 7
1677  else if(ncode == -4) then
1678  code = 6
1679  else if(ncode < 0) then
1680  code = 22
1681  else
1682  code = ncode
1683  end if
1684 
1685  call plsym(x, y, code)
1686 
1687 end subroutine pgpt
1688 !***********************************************************************************************************************************
1689 
1690 
1691 !***********************************************************************************************************************************
1692 !> \brief Draw a single point
1693 !!
1694 !! \param x Horizontal location
1695 !! \param y Vertical location
1696 !! \param ncode Plot symbol
1697 
1698 subroutine pgpt1(x, y, ncode)
1699  use plplot, only: plpoin, plflt
1700  implicit none
1701  integer, intent(in) :: ncode
1702  real(plflt), intent(in) :: x, y
1703  real(plflt), dimension(1) :: xin, yin
1704 
1705  xin(1) = x
1706  yin(1) = y
1707 
1708  call pgpt(1, xin, yin, ncode)
1709 
1710 end subroutine pgpt1
1711 !***********************************************************************************************************************************
1712 
1713 
1714 !***********************************************************************************************************************************
1715 !> \brief Close stream
1716 !!
1717 
1718 subroutine pgclos()
1719  implicit none
1720  call plend1()
1721 end subroutine pgclos
1722 !***********************************************************************************************************************************
1723 
1724 
1725 !***********************************************************************************************************************************
1726 !> \brief Get cursor location - dummy routine!
1727 !!
1728 !! \param mode Display mode (0-7)
1729 !! \param posn If POSN=1, PGBAND tries to place the cursor at (X,Y); if POSN=0, it leaves the cursor at its current position
1730 !! \param xref World x-coordinate of the anchor point
1731 !! \param yref World y-coordinate of the anchor point
1732 !! \param x World x-coordinate of the cursor
1733 !! \retval y World y-coordinate of the cursor
1734 !! \retval ch Character typed by the user
1735 
1736 subroutine pgband(mode, posn, xref, yref, x, y, ch)
1738 
1739  implicit none
1740  integer, intent(in) :: mode, posn
1741  real, intent(in) :: xref, yref
1742  real, intent(out) :: x, y
1743  character, intent(out) :: ch*(*)
1744 
1745  integer, save :: warn
1746  integer :: i1
1747  real :: r1
1748 
1749 
1750  i1 = mode
1751  i1 = posn
1752  i1 = i1 ! Remove 'unused variable error from g95'
1753 
1754  r1 = xref
1755  r1 = yref
1756  r1 = r1 ! Remove 'unused variable error from g95'
1757 
1758  x = 0.0
1759  y = 0.0
1760  ch = char(0)
1761 
1762  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
1763  if(warn.ne.123) call warn_dummy_routine('pgband', '')
1764  warn = 123
1765 
1766 end subroutine pgband
1767 !***********************************************************************************************************************************
1768 
1769 
1770 !***********************************************************************************************************************************
1771 !> \brief Get colour range - dummy routine!
1772 !!
1773 !! \param c1 Lower colour index
1774 !! \param c2 Upper colour index
1775 
1776 subroutine pgqcol(c1, c2)
1778 
1779  implicit none
1780  integer, intent(out) :: c1, c2
1781  integer, save :: warn
1782 
1783 
1784  c1 = 0
1785  c2 = 255
1786 
1787  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
1788  if(warn.ne.123) call warn_dummy_routine('pgqcir','using a default colour range')
1789  warn = 123
1790 
1791 end subroutine pgqcol
1792 !***********************************************************************************************************************************
1793 
1794 
1795 
1796 !***********************************************************************************************************************************
1797 !> \brief Get current colour index - dummy routine!
1798 !!
1799 !! \retval ci Colour index
1800 
1801 subroutine pgqci(ci)
1803 
1804  implicit none
1805  integer, intent(out) :: ci
1806  integer, save :: warn
1807 
1808  ci = 0
1809 
1810  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
1811  if(warn.ne.123) call warn_dummy_routine('pgqcir','using a default colour index')
1812  warn = 123
1813 
1814 end subroutine pgqci
1815 !***********************************************************************************************************************************
1816 
1817 
1818 !***********************************************************************************************************************************
1819 !> \brief Get viewport size
1820 !!
1821 !! \param units specify the units of the output parameters:
1822 !! - 0 : normalized device coordinates
1823 !! - 1 : inches
1824 !! - 2 : millimeters
1825 !! - 3 : pixels
1826 !!
1827 !! \retval x1 Left
1828 !! \retval x2 Right
1829 !! \retval y1 Bottom
1830 !! \retval y2 Top
1831 
1832 subroutine pgqvp(units, x1, x2, y1, y2)
1833  use plplot, only: plflt
1834  use pg2plplot, only : mm_per_inch
1835  implicit none
1836  integer, intent(in) :: units
1837  real, intent(out) :: x1, x2, y1, y2
1838  real(kind=plflt) x1d, x2d, y1d, y2d
1839  real(kind=plflt) xp, yp, xleng, yleng, xoff, yoff
1840 
1841  call plgvpd(x1d, x2d, y1d, y2d)
1842  if(units .ne. 0) then
1843  call plgpage(xp, yp, xleng, yleng, xoff, yoff)
1844  x1 = real(x1d * xleng - xoff)
1845  x2 = real(x2d * xleng - xoff)
1846  y1 = real(y1d * yleng - yoff)
1847  y2 = real(y2d * yleng - yoff)
1848  if(units .eq. 3) then
1849  return
1850  else if(units .eq. 2) then
1851  x1 = x1 / real(xp * mm_per_inch)
1852  x2 = x2 / real(xp* mm_per_inch)
1853  y1 = y1 / real(yp * mm_per_inch)
1854  y2 = y2 / real(yp* mm_per_inch)
1855  return
1856  else if(units .eq. 1) then
1857  x1 = x1 / real(xp)
1858  x2 = x2 / real(xp)
1859  y1 = y1 / real(yp)
1860  y2 = y2 / real(yp)
1861  return
1862  end if
1863  else
1864  print *, "unknown units in pgqvp", units
1865  end if
1866  x1 = real(x1d)
1867  x2 = real(x2d)
1868  y1 = real(y1d)
1869  y2 = real(y2d)
1870 end subroutine pgqvp
1871 
1872 !> \brief get view surface
1873 subroutine pgqvsz(units, x1, x2, y1, y2)
1874  use plplot, only: plflt
1875  use pg2plplot, only : mm_per_inch
1876  implicit none
1877  integer, intent(in) :: units
1878  real, intent(out) :: x1, x2, y1, y2
1879  real(kind=plflt) xp, yp, xleng, yleng, xoff, yoff
1880 
1881  x1 = 0.0
1882  y1 = 0.0
1883 
1884  if(units .eq. 0) then
1885  x2 = 1.0
1886  y2 = 1.0
1887  return
1888  end if
1889 
1890  call plgpage(xp, yp, xleng, yleng, xoff, yoff)
1891  if(units .eq. 3) then
1892  x2 = real(xleng)
1893  y2 = real(yleng)
1894  else if(units .eq. 2) then
1895  x2 = real(xleng / xp* mm_per_inch)
1896  y2 = real(yleng / yp* mm_per_inch)
1897  else if(units .eq. 1) then
1898  x2 = x2 / real(xp)
1899  y2 = y2 / real(yp)
1900  return
1901  else
1902  print *, 'undefined units in pgqvsz'
1903  return
1904  end if
1905 end subroutine pgqvsz
1906 !***********************************************************************************************************************************
1907 
1908 
1909 !***********************************************************************************************************************************
1910 !> \brief Print user name and date in plot - dummy routine!
1911 !!
1912 
1913 subroutine pgiden()
1915 
1916  implicit none
1917  integer, save :: warn
1918 
1919  if(.not.compatibility_warnings) warn = 123 ! Don't warn about compatibility between PGPlot and PLplot
1920  if(warn.ne.123) call warn_dummy_routine('pgiden','ignoring')
1921  warn = 123
1922 
1923 end subroutine pgiden
1924 !***********************************************************************************************************************************
1925 
1926 
1927 
1928 !***********************************************************************************************************************************
1929 !> \brief get window size
1930 !!
1931 !! \retval x1 Left
1932 !! \retval x2 Right
1933 !! \retval y1 Bottom
1934 !! \retval y2 Top
1935 
1936 subroutine pgqwin(x1, x2, y1, y2)
1937  use plplot, only: plflt
1938  implicit none
1939  real, intent(out) :: x1, x2, y1, y2
1940  real(kind=plflt) :: xmin, xmax, ymin, ymax
1941 
1942  call plgspa(xmin, xmax, ymin, ymax)
1943  x1 = real(xmin)
1944  x2 = real(xmax)
1945  y1 = real(ymin)
1946  y2 = real(ymax)
1947 
1948 end subroutine pgqwin
1949 !***********************************************************************************************************************************
1950 
1951 
1952 !***********************************************************************************************************************************
1953 !> \brief Flush - dummy routine
1954 !!
1955 
1956 subroutine pgupdt()
1957  return
1958 end subroutine pgupdt
1959 !***********************************************************************************************************************************
1960 
1961 
1962 !***********************************************************************************************************************************
1963 !> \brief Search and replace occurences of a substring in a string, taken from libSUFR
1964 !!
1965 !! \param string Original string to replace in
1966 !! \param str_in Search string
1967 !! \param str_out Replacement string
1968 
1969 subroutine replace_substring(string, str_in, str_out)
1970  implicit none
1971  character, intent(inout) :: string*(*)
1972  character, intent(in) :: str_in*(*),str_out*(*)
1973  integer :: is,l
1974 
1975  l = len_trim(str_in)
1976  is = huge(is)
1977  do
1978  is = index(string, str_in, back=.false.)
1979  if(is.le.0) exit
1980  string = string(1:is-1)//trim(str_out)//trim(string(is+l:))
1981  end do
1982 
1983 end subroutine replace_substring
1984 !***********************************************************************************************************************************
1985 
1986 !***********************************************************************************************************************************
1987 !> \brief Inquire fill-area style - dummy routine
1988 !!
1989 !! \retval fs The current fill-area style
1990 
1991 subroutine pgqfs(fs)
1992  implicit none
1993  integer, intent(out) :: fs
1994  fs = 1
1995 end subroutine pgqfs
1996 !***********************************************************************************************************************************
1997 
1998 !***********************************************************************************************************************************
1999 !> \brief Set text background color index - dummy routine
2000 !!
2001 !! \param b Background color index
2002 
2003 subroutine pgstbg(b)
2004  implicit none
2005  integer, intent(in) :: b
2006  integer :: bb
2007  bb = b ! Use the dummy variable to prevent compiler complaints
2008  bb = bb ! Use the local variable after setting it to prevent compiler complaints
2009 end subroutine pgstbg
2010 !***********************************************************************************************************************************
2011 
2012 !***********************************************************************************************************************************
2013 !> \brief Query text background color index - dummy routine
2014 !!
2015 !! \param b Background color index
2016 
2017 subroutine pgqtbg(b)
2018  implicit none
2019  integer, intent(out) :: b
2020  b = 0
2021 end subroutine pgqtbg
2022 !***********************************************************************************************************************************
2023 
2024 
2025 !***********************************************************************************************************************************
2026 !> \brief Control new page prompting - dummy routine
2027 !!
2028 !! \param prompt Set prompt state to ON for interactive devices
2029 
2030 subroutine pgask(prompt)
2031  implicit none
2032  logical, intent(in) :: prompt
2033  logical :: prompt1
2034  prompt1 = prompt ! Use the dummy variable to prevent compiler complaints
2035  prompt1 = prompt1 ! Use the local variable after setting it to prevent compiler complaints
2036 end subroutine pgask
2037 !***********************************************************************************************************************************
2038 
2039 !***********************************************************************************************************************************
2040 !> \brief Set window and viewport and draw labeled frame
2041 !!
2042 !! \param xmin Left
2043 !! \param xmax Right
2044 !! \param ymin Top (?)
2045 !! \param ymax Bottom (?)
2046 !! \param just JUST=1: scales of x and y axes are equal, otherwise independent
2047 !! \param axis Controls the plotting of axes, tick marks, etc
2048 
2049 subroutine pgenv(xmin, xmax, ymin, ymax, just, axis)
2050  use plplot, only: plflt
2051  implicit none
2052  real, intent(in) :: xmin, xmax, ymin, ymax
2053  integer, intent(in) :: just, axis
2054  real(kind=plflt) :: xminp, xmaxp, yminp, ymaxp
2055  xminp = xmin
2056  xmaxp = xmax
2057  yminp = ymin
2058  ymaxp = ymax
2059  call plenv(xminp, xmaxp, yminp, ymaxp, just, axis)
2060 end subroutine pgenv
2061 !***********************************************************************************************************************************
2062 
2063 !***********************************************************************************************************************************
2064 !> \brief Plot a histogram of binned data
2065 !!
2066 !! \param nbin Number of bins
2067 !! \param x Abscissae of the bins
2068 !! \param data Data values of bins
2069 !! \param center If .TRUE., the X values denote the center of the bin, else the lower edge
2070 
2071 subroutine pgbin(nbin, x, data, center)
2072  use plplot, only: plflt, plbin
2073  implicit none
2074  integer, intent(in) :: nbin
2075  real, intent(in) :: x(*), data(*)
2076  logical, intent(in) :: center
2077  real(kind=plflt) :: xp(nbin), datap(nbin)
2078  integer :: opt
2079  xp = x(1:nbin)
2080  datap = data(1:nbin)
2081  if(center) opt = 1
2082  call plbin(xp, datap, opt)
2083 end subroutine pgbin
2084 !***********************************************************************************************************************************
2085 
2086 !***********************************************************************************************************************************
2087 !> \brief Vertical error bar
2088 !!
2089 !! \param n Number of error bars to plot
2090 !! \param x World horizontal coordinates of the data
2091 !! \param ymax World vertical coordinates of top end of the error bars
2092 !! \param ymin World vertical coordinates of bottom end of the error bars
2093 !! \param t Length of terminals to be drawn at the ends of the error bar
2094 
2095 subroutine pgerry(n, x, ymax, ymin, t)
2096  use plplot, only: plflt, plerry
2097  implicit none
2098  integer, intent(in) :: n
2099  real, intent(in) :: x(*), ymin(*), ymax(*), t
2100  real(kind=plflt) :: xp(n), yminp(n), ymaxp(n)
2101 
2102  xp = x(1:n) + 0*t ! Use dummy variable t somewhere to prevent compiler complaints
2103  yminp = ymin(1:n)
2104  ymaxp = ymax(1:n)
2105  call plerry(xp, yminp, ymaxp)
2106 
2107 end subroutine pgerry
2108 !***********************************************************************************************************************************
2109 
2110 !***********************************************************************************************************************************
2111 !> \brief Set viewport (inches)
2112 !!
2113 !! \param xleft Horizontal coordinate of left-hand edge of viewport, in inches from left edge of view surface
2114 !! \param xright Horizontal coordinate of right-hand edge of viewport, in inches from left edge of view surface
2115 !! \param ybot Vertical coordinate of bottom edge of viewport, in inches from bottom of view surface
2116 !! \param ytop Vertical coordinate of top edge of viewport, in inches from bottom of view surface
2117 
2118 subroutine pgvsiz(xleft, xright, ybot, ytop)
2119  use plplot, only: plflt
2120  implicit none
2121  real, intent(in) :: xleft, xright, ybot, ytop
2122  real(kind=plflt) :: x1p, x2p, y1p, y2p
2123  x1p=xleft
2124  x2p=xright
2125  y1p=ybot
2126  y2p=ytop
2127  call plsvpa(x1p, x2p, y1p, y2p)
2128 end subroutine pgvsiz
2129 !***********************************************************************************************************************************
2130 
2131 !***********************************************************************************************************************************
2132 !> \brief Non-standard alias for PGVSIZ
2133 !!
2134 !! \param xleft Horizontal coordinate of left-hand edge of viewport, in inches from left edge of view surface
2135 !! \param xright Horizontal coordinate of right-hand edge of viewport, in inches from left edge of view surface
2136 !! \param ybot Vertical coordinate of bottom edge of viewport, in inches from bottom of view surface
2137 !! \param ytop Vertical coordinate of top edge of viewport, in inches from bottom of view surface
2138 
2139 subroutine pgvsize(xleft, xright, ybot, ytop)
2140  implicit none
2141  real, intent(in) :: xleft, xright, ybot, ytop
2142  call pgvsiz(xleft, xright, ybot, ytop)
2143 end subroutine pgvsize
2144 !***********************************************************************************************************************************
2145 
2146 !***********************************************************************************************************************************
2147 !> \brief Find length of a string in a variety of units - dummy variable
2148 !!
2149 !! \param units 0: normalized device coordinates, 1: in, 2: mm, 3: absolute device coordinates (dots), 4: world coordinates,
2150 !! 5: fraction of the current viewport size
2151 !! \param string String of interest
2152 !! \param xl Length of string in horizontal direction
2153 !! \param yl Length of string in vertical direction
2154 
2155 subroutine pglen(units, string, xl, yl)
2156  implicit none
2157  character, intent(in) :: string*(*)
2158  integer, intent(in) :: units
2159  real, intent(out) :: xl, yl
2160  character :: dummy
2161 
2162  print *, "pglen not implemented"
2163  dummy = string(1:1) ! Use the dummy variable string somewhere to prevent compiler complaints
2164  dummy = dummy ! Use the local variable after setting it to prevent compiler complaints
2165  xl = 1.0 + real(0*units) ! Use the dummy variable units somewhere to prevent compiler complaints
2166  yl = 1.0
2167 
2168 end subroutine pglen
2169 !***********************************************************************************************************************************
2170 
2171 !***********************************************************************************************************************************
2172 !> \brief Find bounding box of text string - dummy variable
2173 !!
2174 !! \param x World x-coordinate
2175 !! \param y World y-coordinate
2176 !! \param angle Angle between baseline and horizontal (degrees)
2177 !! \param fjust Horizontal justification
2178 !! \param text Text string that would be written
2179 !!
2180 !! \retval xbox Horizontal limits of the bounding box
2181 !! \retval ybox Vertical limits of the bounding box
2182 
2183 
2184 subroutine pgqtxt(x, y, angle, fjust, text, xbox, ybox)
2185  implicit none
2186  character, intent(in) :: text*(*)
2187  integer, intent(in) :: x, y, angle, fjust
2188  real, intent(out) :: xbox(4), ybox(4)
2189  integer :: dumInt
2190  character :: dumStr
2191 
2192  ! Use the dummy input variables somewhere to prevent compiler complaints:
2193  dumstr = text(1:1)
2194  dumint = x + y + angle + fjust
2195  ! Use the local variable after setting it to prevent compiler complaints:
2196  dumstr = dumstr
2197  dumint = dumint
2198 
2199  print *, "pgqtxt not implemented"
2200 
2201  xbox(1) = 0.0
2202  xbox(2) = 1.0
2203  xbox(3) = 0.0
2204  xbox(4) = 1.0
2205 
2206  ybox(1) = 0.0
2207  ybox(2) = 1.0
2208  ybox(3) = 0.0
2209  ybox(4) = 1.0
2210 
2211 end subroutine pgqtxt
2212 !***********************************************************************************************************************************
2213 
2214 
2215 
2216 
2217 !***********************************************************************************************************************************
2218 !> \brief Return the relative difference between two numbers: dx/<x> - taken from libSUFR, turned into single precision
2219 !!
2220 !! \param x1 First number
2221 !! \param x2 Second number
2222 
2223 function reldiff(x1,x2)
2224  implicit none
2225  real, intent(in) :: x1,x2
2226  real :: reldiff, xsum,xdiff
2227 
2228  xsum = x1+x2
2229  xdiff = x2-x1
2230 
2231  if(abs(xsum).gt.tiny(xsum)) then
2232  reldiff = xdiff / (xsum*0.5)
2233  else ! Can't divide by zero
2234  if(abs(xdiff).gt.tiny(xdiff)) then
2235  reldiff = xdiff
2236  else
2237  reldiff = 1.0 ! 0/0
2238  end if
2239  end if
2240 
2241 end function reldiff
2242 !***********************************************************************************************************************************
2243 
2244 
2245 !***********************************************************************************************************************************
2246 !> \brief Check whether opening a file gives an error
2247 !!
2248 !! \param fname File name
2249 
2250 function check_error(fname)
2251  implicit none
2252  character, intent(in) :: fname*(*)
2253  integer :: err, check_error
2254 
2255  open(9999, file=trim(fname), err=450, iostat=err)
2256  close(9999)
2257  check_error = 0
2258  return
2259 
2260 450 continue
2261  check_error = err
2262 
2263 end function check_error
2264 !***********************************************************************************************************************************
2265 
2266 
2267 !***********************************************************************************************************************************
2268 !> \brief Print a warning when no (proper) PLplot equivalent is defined and a dummy routine is used instead
2269 !!
2270 !! \param routine Name of the undefined routine
2271 !! \param message Message to be postponed
2272 
2273 subroutine warn_dummy_routine(routine, message)
2274  implicit none
2275  character, intent(in) :: routine*(*),message*(*)
2276 
2277  write(0,'(/,A)') '*** PG2PLplot WARNING: no PLplot equivalent was found for the PGplot routine '//trim(routine)//'() '// &
2278  trim(message)//' ***'
2279 
2280 end subroutine warn_dummy_routine
2281 !***********************************************************************************************************************************
2282 
2283 
2284 !***********************************************************************************************************************************
2285 !> \brief Test whether two double-precision variables are equal to better than twice the machine precision, taken from libSUFR
2286 !!
2287 !! \param x1 First number
2288 !! \param x2 Second number
2289 
2290 function deq(x1,x2)
2291  use plplot, only: plflt
2292 
2293  implicit none
2294  real(plflt), intent(in) :: x1,x2
2295  real(plflt) :: eps
2296  logical :: deq
2297 
2298  eps = 2*tiny(x1)
2299  if(abs(x1-x2).le.eps) then
2300  deq = .true.
2301  else
2302  deq = .false.
2303  end if
2304 
2305 end function deq
2306 !***********************************************************************************************************************************
2307 
2308 
2309 !***********************************************************************************************************************************
2310 !> \brief Test whether two single-precision variables are equal to better than twice the machine precision, taken from libSUFR
2311 !!
2312 !! \param x1 First number
2313 !! \param x2 Second number
2314 
2315 function seq(x1,x2)
2316  implicit none
2317  real, intent(in) :: x1,x2
2318  real :: eps
2319  logical :: seq
2320 
2321  eps = 2*tiny(x1)
2322  if(abs(x1-x2).le.eps) then
2323  seq = .true.
2324  else
2325  seq = .false.
2326  end if
2327 
2328 end function seq
2329 !***********************************************************************************************************************************
2330 
2331 
subroutine pgscf(cf)
Set character font.
Definition: pg2plplot.f90:159
logical, parameter compatibility_warnings
Definition: pg2plplot.f90:37
subroutine pgebuf()
End buffering output - dummy routine!
Definition: pg2plplot.f90:1256
subroutine pgqlw(lw)
Query line width.
Definition: pg2plplot.f90:144
subroutine pgqcir(ci1, ci2)
Query colour-index range for pggray and pgimag - dummy routine!
Definition: pg2plplot.f90:323
subroutine pgupdt()
Flush - dummy routine.
Definition: pg2plplot.f90:1957
subroutine pgsah(fs, angle, barb)
Set arrow head - dummy routine!
Definition: pg2plplot.f90:481
subroutine pgptxt(x1, y1, ang, just1, text)
Print text with arbitrary angle and justification.
Definition: pg2plplot.f90:766
subroutine pglen(units, string, xl, yl)
Find length of a string in a variety of units - dummy variable.
Definition: pg2plplot.f90:2156
integer function check_error(fname)
Check whether opening a file gives an error.
Definition: pg2plplot.f90:2251
integer cur_lstyle
Definition: pg2plplot.f90:50
subroutine pgsubp(nxsub, nysub)
Subdivide view surface into panels.
Definition: pg2plplot.f90:1201
subroutine pgband(mode, posn, xref, yref, x, y, ch)
Get cursor location - dummy routine!
Definition: pg2plplot.f90:1737
subroutine pg2pldev(pgdev, pldev, filename)
Converts PGplot device ID to PLplot device ID + file name.
Definition: pg2plplot.f90:1502
subroutine pgqtxt(x, y, angle, fjust, text, xbox, ybox)
Find bounding box of text string - dummy variable.
Definition: pg2plplot.f90:2185
subroutine pgolin(maxpt, npt, x, y, symbol)
Read data from screen - no PLplot equivalent (yet) - dummy routine!
Definition: pg2plplot.f90:1421
subroutine pgscir(ci1, ci2)
Set colour-index range for pggray and pgimag - dummy routine!
Definition: pg2plplot.f90:298
subroutine pgdraw(x, y)
Draw line to location.
Definition: pg2plplot.f90:1641
subroutine pgbbuf()
Start buffering output - dummy routine!
Definition: pg2plplot.f90:1237
integer, dimension(max_level) save_color
Definition: pg2plplot.f90:49
subroutine pgiden()
Print user name and date in plot - dummy routine!
Definition: pg2plplot.f90:1914
integer, dimension(max_level) save_ffamily
Definition: pg2plplot.f90:44
integer, parameter max_open
Definition: pg2plplot.f90:42
subroutine pgconf(arr, nx, ny, ix1, ix2, iy1, iy2, c1, c2, tr)
Shade a region (between contours/heights) - dummy routine!
Definition: pg2plplot.f90:708
subroutine pgsave()
Save parameters.
Definition: pg2plplot.f90:1563
subroutine pgptext(x, y, ang, just, text)
Non-standard alias for pgptxt()
Definition: pg2plplot.f90:823
subroutine pgqch(ch)
Query character height.
Definition: pg2plplot.f90:457
subroutine pgqtbg(b)
Query text background color index - dummy routine.
Definition: pg2plplot.f90:2018
integer, dimension(max_level) save_fweight
Definition: pg2plplot.f90:46
subroutine warn_dummy_routine(routine, message)
Print a warning when no (proper) PLplot equivalent is defined and a dummy routine is used instead...
Definition: pg2plplot.f90:2274
real(plflt) xcur
Definition: pg2plplot.f90:39
subroutine pgstbg(b)
Set text background color index - dummy routine.
Definition: pg2plplot.f90:2004
real(plflt) paper_ratio
Definition: pg2plplot.f90:53
logical function deq(x1, x2)
Test whether two double-precision variables are equal to better than twice the machine precision...
Definition: pg2plplot.f90:2291
subroutine pgbox(xopt, xtick1, nxsub, yopt, ytick1, nysub)
Draw a box (+axes) around a plot.
Definition: pg2plplot.f90:1281
subroutine pgqcr(ci1, r, g, b)
Query colour representation.
Definition: pg2plplot.f90:271
subroutine pgqwin(x1, x2, y1, y2)
get window size
Definition: pg2plplot.f90:1937
subroutine pgtick(x1, y1, x2, y2, pos, tikl, tikr, disp, orient, lbl)
Draw a single tick mark, optionally with label - no PLplot routine found, using an ad-hoc routine...
Definition: pg2plplot.f90:1317
subroutine pgline(n, x1, y1)
Draw a line.
Definition: pg2plplot.f90:509
subroutine pgvsize(xleft, xright, ybot, ytop)
Non-standard alias for PGVSIZ.
Definition: pg2plplot.f90:2140
subroutine pgclos()
Close stream.
Definition: pg2plplot.f90:1719
subroutine pglabel(xlbl, ylbl, toplbl)
Alias for pglab.
Definition: pg2plplot.f90:947
subroutine pgarro(x1, y1, x2, y2)
Draw an arrow - only a line is drawn, arrows not supported in PLplot!
Definition: pg2plplot.f90:536
real(plflt) ycur
Definition: pg2plplot.f90:39
subroutine pgadvance()
Alias for pgpage.
Definition: pg2plplot.f90:1227
subroutine pgenv(xmin, xmax, ymin, ymax, just, axis)
Set window and viewport and draw labeled frame.
Definition: pg2plplot.f90:2050
real(plflt), parameter mm_per_inch
Definition: pg2plplot.f90:43
subroutine pgunsa()
Unsave parameters.
Definition: pg2plplot.f90:1589
subroutine pgwindow(xmin1, xmax1, ymin1, ymax1)
alias for pgswin
Definition: pg2plplot.f90:1178
subroutine pgend()
End a plot.
Definition: pg2plplot.f90:1076
subroutine pgqvp(units, x1, x2, y1, y2)
Get viewport size.
Definition: pg2plplot.f90:1833
subroutine pgask(prompt)
Control new page prompting - dummy routine.
Definition: pg2plplot.f90:2031
logical is_init
Definition: pg2plplot.f90:52
subroutine pgqvsz(units, x1, x2, y1, y2)
get view surface
Definition: pg2plplot.f90:1874
subroutine pglab(xlbl, ylbl, toplbl)
Interface for pglab.
Definition: pg2plplot.f90:929
subroutine pgqls(ls)
Query line style.
Definition: pg2plplot.f90:107
subroutine pgmtxt(side, disp1, pos1, just1, text)
Print text in the margin.
Definition: pg2plplot.f90:875
subroutine pgsch(ch)
Set character height.
Definition: pg2plplot.f90:397
integer, parameter max_level
Definition: pg2plplot.f90:41
integer, dimension(max_level) save_lwidth
Definition: pg2plplot.f90:47
subroutine pgqci(ci)
Get current colour index - dummy routine!
Definition: pg2plplot.f90:1802
subroutine pgerry(n, x, ymax, ymin, t)
Vertical error bar.
Definition: pg2plplot.f90:2096
subroutine pgtext(x1, y1, text)
Print text with default angle and justification.
Definition: pg2plplot.f90:841
subroutine pgvsiz(xleft, xright, ybot, ytop)
Set viewport (inches)
Definition: pg2plplot.f90:2119
real function reldiff(x1, x2)
Return the relative difference between two numbers: dx/ - taken from libSUFR, turned into single p...
Definition: pg2plplot.f90:2224
subroutine pgqcol(c1, c2)
Get colour range - dummy routine!
Definition: pg2plplot.f90:1777
integer cur_color
Definition: pg2plplot.f90:50
subroutine pgswin(xmin1, xmax1, ymin1, ymax1)
Set window.
Definition: pg2plplot.f90:1152
real(plflt) paper_width
Definition: pg2plplot.f90:53
integer, dimension(max_level) save_fstyle
Definition: pg2plplot.f90:45
subroutine pgscr(ci1, r, g, b)
Set colour representation.
Definition: pg2plplot.f90:241
integer devid
Definition: pg2plplot.f90:51
subroutine pgcirc(xc, yc, r)
Draw a circle.
Definition: pg2plplot.f90:637
subroutine pgpage()
Advance to the next (sub-)page.
Definition: pg2plplot.f90:1215
subroutine pgrect(x1, x2, y1, y2)
Draw a rectangle.
Definition: pg2plplot.f90:614
integer cur_lwidth
Definition: pg2plplot.f90:50
subroutine pgsci(ci1)
Set colour index.
Definition: pg2plplot.f90:215
subroutine pgqcs(unit, xch, ych)
Inquire character height in a variety of units.
Definition: pg2plplot.f90:425
subroutine pgbeg(i, pgdev, nx, ny)
Alias for pgbegin.
Definition: pg2plplot.f90:1063
subroutine pgslct(pgdev)
Select output stream.
Definition: pg2plplot.f90:1542
subroutine pgcont(arr, nx, ny, ix1, ix2, iy1, iy2, c, nc, tr)
Make a contour plot.
Definition: pg2plplot.f90:676
subroutine pgqinf(item, value, length)
Inquire PGPLOT general information - dummy routine.
Definition: pg2plplot.f90:199
subroutine pgsvp(xl1, xr1, yb1, yt1)
Set view port.
Definition: pg2plplot.f90:1125
subroutine pgqid(id)
Inquire current device identifier.
Definition: pg2plplot.f90:1004
subroutine pgmove(x, y)
Move pen to location - for use with pgdraw() only!
Definition: pg2plplot.f90:1623
subroutine pgbegin(i, pgdev, nx, ny)
Begin a new plot.
Definition: pg2plplot.f90:1023
integer function pgopen(pgdev)
Start a new plot.
Definition: pg2plplot.f90:962
subroutine pgpoint(n, x1, y1, s)
Draw points.
Definition: pg2plplot.f90:564
subroutine pgpt1(x, y, ncode)
Draw a single point.
Definition: pg2plplot.f90:1699
subroutine pgbin(nbin, x, data, center)
Plot a histogram of binned data.
Definition: pg2plplot.f90:2072
PG2PLplot module.
Definition: pg2plplot.f90:28
integer, dimension(max_level) save_lstyle
Definition: pg2plplot.f90:48
subroutine pgqfs(fs)
Inquire fill-area style - dummy routine.
Definition: pg2plplot.f90:1992
subroutine pgpap(width, ratio)
Set paper size.
Definition: pg2plplot.f90:1091
subroutine pgmtext(side, disp, pos, just, text)
Alias for pgmtxt()
Definition: pg2plplot.f90:907
subroutine pgpoly(n, x1, y1)
Draw a polygone.
Definition: pg2plplot.f90:588
subroutine pg2pltext(string)
Replace the PGPlot escape character '\' with the PLplot escape character '#'.
Definition: pg2plplot.f90:1464
real(plflt), parameter ch_fac
Conversion factor for the character height.
Definition: pg2plplot.f90:35
subroutine pgsls(ls)
Set line style.
Definition: pg2plplot.f90:79
subroutine replace_substring(string, str_in, str_out)
Search and replace occurences of a substring in a string, taken from libSUFR.
Definition: pg2plplot.f90:1970
subroutine pgshs(ang, sep, ph)
Set hash style.
Definition: pg2plplot.f90:376
integer save_level
Definition: pg2plplot.f90:40
subroutine pgpt(n, x, y, ncode)
Draw a point.
Definition: pg2plplot.f90:1666
subroutine pgslw(lw)
Set line width.
Definition: pg2plplot.f90:120
logical function seq(x1, x2)
Test whether two single-precision variables are equal to better than twice the machine precision...
Definition: pg2plplot.f90:2316
subroutine pgeras()
Erase screen.
Definition: pg2plplot.f90:1450
subroutine pgsfs(fs)
Set fill style.
Definition: pg2plplot.f90:346
subroutine do_init()
Definition: pg2plplot.f90:59
subroutine pgqcf(cf)
Query character font.
Definition: pg2plplot.f90:173