C file PCPROF2.FOR
C ----------------------------------------------------------------------!-------
CONTAINS PLOTAL, PLTINI, PSCALE
C ----------------------------------------------------------------------!-------
C 10 Nov 2009 - Added:     double precision d1mach
C 26 Aug 2008 - Replaced using rmach() with direct calls to 
C                I1mach(), R1mach(), or D1mach()
C 13 Mar 2008 - Corrected length in common block /HDR/
C 31 Mar 2006 - Added phase integral to SUBROUTINE PLOTAL
C 24 Mar 2006 - Moved HEIGHT to S6
C 24 Mar 2006 - Moved GET to S5
C 24 Mar 2006 - Moved GRACE, GRAFB to S8
C 24 Mar 2006 - Changed COMMON block /RAYDEV/
C 17 Mar 2006 - Added phase speed to SUBROUTINE PLOTAL
C 17 Mar 2006 - Added wavelength to SUBROUTINE PLOTAL
C 14 Mar 2006 - Added vertical wavenumber to SUBROUTINE PLOTAL
C  1 Mar 2006 - Removed ITEST, ITOC
C  1 Mar 2006 - Added square of Brunt Vaisala frequency to SUBROUTINE PLOTAL
C 28 Feb 2006 - Added ACF, IACF, BV, IBV to EXTERNAL
C 24 Feb 2006 - Added acoustic cutoff frequency to SUBROUTINE PLOTAL
C 24 Feb 2006 - Added Brunt-Vaisala frequency to SUBROUTINE PLOTAL
C 24 June 2004 - added wind/current direction to SUBROUTINE PLOTAL
C 4 Jun 2004 - added an argument to SETGET
C Note: changes marked as "10 Jun 2003" were not made to this file on that date.
C They are the changes in the 10 Jun 2003 version of the raytracing program.
C 3 May 2004 - /CRMACH/RMACH(5) --> (9)
C 10 Jun 2003 - lines 243, 249, 263, (FUNCTION GET) added ,PK(10)
C 3 May 2004 - (FUNCTION GET) NZ, etc. from replacements, not IMOVE
C ----------------------------------------------------------------------!-------
      SUBROUTINE PLOTAL(PLTMODE,PLTCNT,TLTS,NP,MP,WK,FG,XN,YN,MARK)
C ----------------------------------------------------------------------!-------
C PLOT ONE PROFILE 
C ----------------------------------------------------------------------!-------
      CHARACTER*(*) TLTS(4)
      DOUBLE PRECISION WK(NP,2)
      INTEGER PLTCNT
C
C*****"PLTMODE" IS A LOGICAL ARRAY.  SEE SUBROUTINE PRFPLT.  JAH.
      LOGICAL PLTMODE(12)
C
      EXTERNAL REARTH,TOPV,COLAT,SPEEDV,ISPEDV,ITEMPV,TEMPV
      EXTERNAL WNMAGV,WNVRTV,WNSTHV,WNESTV,IWINDV
      EXTERNAL PTPLTX,PTPLTY,IPTPLT,ITOPV
     1 ,WNANGV                                                ! 24 June 2004
     2, ACF, IACF, BV, IBV                                    ! 28 Feb 2006
     3, BV2, IBV2                                             !  1 Mar 2006
     4, VWN, IVWN                                             ! 14 Mar 2006
     5, WL, IWL, PHSPEED, IPHSPEED                            ! 17 Mar 2006
     6, VRT, IVRT                                             ! 31 Mar 2006

C.....Added the following entry points in PRFSUP for east, south, vertical
C     components of current velocity. (and direction 24 June 2004)

      external ieastv, isoutv, ivertv, iangv ! 24 June 2004
C
C.....Just to re-iterate (see the blank input data forms for profiles
C     for more information),
C
C     PLTCNT (aka W(90)) Value        Type of Profile
C
C           1                         Sound Speed
C           2                         Current Speed
C           3,4,5                     Vertical,Southward,Eastward
C                                        Current Components
C           6                         Bottom Topography
C           7                         Temperature
C           8                         Wind/Current Direction         24 Jun 2004
C           9                         Acoustic cutoff frequency      24 Feb 2006
C          10                         Brunt-Vaisala frequency        24 Feb 2006
C          11                         Static stability (square of N)  1 Mar 2006
C          12                         Vertical wavenumber            14 Mar 2006
C          13                         wavelength                     17 Mar 2006
C          14                         phase speed                    17 Mar 2006
C          15                         phase integral                 31 Mar 2006
C

      GO TO (1,2,3,4,5,6,7,8,9,10,11,12,13,14,15) PLTCNT           ! 31 Mar 2006
C
      CALL PRFPLT(PLTMODE,IPTPLT,PTPLTX,PTPLTY,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
C
1     CALL PRFPLT(PLTMODE,ISPEDV,REARTH,SPEEDV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
2     CALL PRFPLT(PLTMODE,IWINDV,REARTH,WNMAGV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
3     CALL PRFPLT(PLTMODE,IVERTV,REARTH,WNVRTV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
4     CALL PRFPLT(PLTMODE,ISOUTV,REARTH,WNSTHV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
5     CALL PRFPLT(PLTMODE,IEASTV,REARTH,WNESTV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
6     CALL PRFPLT(PLTMODE,ITOPV,COLAT,TOPV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
7     CALL PRFPLT(PLTMODE,ITEMPV,REARTH,TEMPV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK)
      RETURN
8     CALL PRFPLT(PLTMODE,IANGV,REARTH,WNANGV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK) ! 24 June 2004
      RETURN
9     CALL PRFPLT(PLTMODE,IACF,REARTH,ACF,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK) ! 24 Feb 2006
      RETURN
10    CALL PRFPLT(PLTMODE,IBV,REARTH,BV,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK) ! 24 Feb 2006
      RETURN
11    CALL PRFPLT(PLTMODE,IBV2,REARTH,BV2,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK) !  1 Mar 2006
      RETURN
12    CALL PRFPLT(PLTMODE,IVWN,REARTH,VWN,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK) ! 14 Mar 2006
      RETURN
13    CALL PRFPLT(PLTMODE,IWL,REARTH,WL,TLTS,NP,MP,WK,FG,XN,YN,
     &            MARK) ! 17 Mar 2006
      RETURN
14    CALL PRFPLT(PLTMODE,IPHSPEED,REARTH,PHSPEED,TLTS,
     & NP,MP,WK,FG,XN,YN, MARK) ! 17 Mar 2006
      RETURN
15    CALL PRFPLT(PLTMODE,IVRT,REARTH,VRT,TLTS,
     & NP,MP,WK,FG,XN,YN, MARK)                            ! 31 Mar 2006
      RETURN
C ----------------------------------------------------------------------!-------
      END  SUBROUTINE PLOTAL
C ----------------------------------------------------------------------!-------
C ----------------------------------------------------------------------!-------
      SUBROUTINE PLOTINI_test(DHITE,SYMRED)
C ----------------------------------------------------------------------!-------
C+P   (PLOTTING) INITIALIZE PLOTTING FOR EIGEN PROGRAM
C
      implicit double precision (a-h,o-z)

C     COMMON DECK "HDR" INSERTED HERE
      CHARACTER*10 INITID*80,DAT,TOD
         COMMON/HDR/SEC,stsec                         !Rev 7/13/90 ! 13 Mar 2008
      COMMON/HDRC/INITID,DAT,TOD
C
      CALL COMPRS
      CALL GRACE(0.)
      CALL PHYSOR(2.,1.3)
      CALL PAGE(12.,10.)
      CALL HEIGHT(DHITE)
      CALL SETANN
      S=SYMRED
C  DEFAULT TO 15% REDUCTION ON POINT SIZE OVER DISSPLAY'S DEFAULT
      IF(S.EQ.0.0) S=.15
      CALL SCLPIC(S)
C ----------------------------------------------------------------------!-------
      END SUBROUTINE PLOTINI_test
C ----------------------------------------------------------------------!-------
C ----------------------------------------------------------------------!-------
      SUBROUTINE PLTINI
C ----------------------------------------------------------------------!-------
C PLOT INITIALIZATION

      implicit double precision (a-h,o-z)

C     COMMON DECK "WWE" INSERTED HERE
      PARAMETER (NWARSZ=1000)
      COMMON/WW/ dum(10), MAXW,W(NWARSZ)
      DOUBLE PRECISION NRSLMI
      EQUIVALENCE (W(80),APOG) ,(W(401),PTYP) ,(W(403),ELMIN)
     1   ,(W(404),ELMAX) ,(W(405),ELVINC) ,(W(406),RGINC)
     2   ,(W(407),HOPLME) ,(W(408),NRSLME)
     3   ,(W(409),GAPLO), (W(410),GAPHI), (W(411),SYMRED)
     4   ,(W(428),GRPMIN) ,(W(429),GRPMAX) ,(W(430),GRPVINC)
     5   ,(W(432),HOPLMG) ,(W(433),NRSLMG)
     6   ,(W(434),RFSPEED) ,(W(435),TMFT) ,(W(436),RMAX)
     7   ,(W(453),RCRNG) ,(W(454),SLPCUT)
     8   ,(W(456),HOPLMI) ,(W(457),NRSLMI) ,(W(458),ESET)
C
C Added "d0" to signify double precision constant.

      CALL PLOTINI(.3d0,SYMRED)
C
C ----------------------------------------------------------------------!-------
      END SUBROUTINE PLTINI
C ----------------------------------------------------------------------!-------
C ----------------------------------------------------------------------!-------
      SUBROUTINE PSCALE(AM,AX,AI,NTICKS,AXIS)
C ----------------------------------------------------------------------!-------
C THIS SUBROUTINE CHOOSES THE SCALE FOR AN AXIS BASED ON THE
C MINIMUM (AM) AND MAXIMUM (AX) VALUES FURNISHED.

      implicit double precision (a-h,o-z)

      DOUBLE PRECISION S(5)
      INTEGER NTM(5)
      DATA NPTS/5/
C
      DATA S/1.0,2.0,5.0,10.0,20.0/
      DATA NTM/5,4,5,5,4/
      DATA SLOP/0.01/
c      print *, 'AX = ', AX, 'AM = ', AM
C
C If the maximum is less than the minimum, then things are hosed.

      if (ax .lt. am) then
	   PRINT*, 'ERROR: The minimum value exceeds the maximum value'
	   STOP    '       in WORK array in subroutine PRFPLT.'
      endif

C If the minimum value = the maximum value, then all values are the same
C for this coordinate of the plot (i.e. a straight line)
C Also, we set the increment to 1, if it was zero, in order to avoid
C division by zero in PRFPLT.

      IF(AX.eq.AM) then
	  ax = ax + 1
	  am = am - 1
	  if (ai .eq. 0) ai = 1
	  NTICKS = 3
	  RETURN
      endif
C
      NTICKS=1
      IF(NCLOSE(AM,AX,AI,NTICKS,AXIS).EQ.1) RETURN
C
c      print *, 'AX = ', AX, 'AM = ', AM
      D=10.0**(dint(dlog10(AX-AM))-1.0)
c      print *, 'AX = ', AX, 'AM = ', AM
      AI=1.E30
      U=AI
      DO 100 I=1,NPTS
      T=S(I)*D
      AM1=dINT(AM/T+SLOP)*T
      IF(AM.LT.0.0) AM1=AM1-T
      AX1=dINT(AX/T+1.0-SLOP)*T
      IF(AX.LT.0.0) AX1=AX1-T
      TMP=(AX1-AM1)/T
      IF( (I.LT.NPTS.OR.AI.NE.U).AND.(T.GT.AI .OR. TMP.GT.6.5))
     1     GO TO 100
      NTICKS=NTM(I)
      AI=T
      AX2=AX1
      AM2=AM1
100   CONTINUE
      IF(AX.GT.AX2+AI/10.) AX2=AX2+AI
      AX=AX2
      IF(AM.LT.AM2-AI/10.) AM2=AM2-AI
      AM=AM2
      IF(NTICKS.GT.0) RETURN
C
      N=NCLOSE(AM,AX,AI,NTICKS,AXIS)
C ----------------------------------------------------------------------!-------
      END SUBROUTINE PSCALE
C ----------------------------------------------------------------------!-------
C ----------------------------------------------------------------------!-------
      FUNCTION GET_test(Z)
C ----------------------------------------------------------------------!-------

C     Note:  The only reason this function is included for the PROFIL
C	     program is for the SETGET entry point (which is called by
C	     CONBLK).
C     'GET' AND ENTRY 'DOT' PROVIDE A CONTROL METHOD FOR AVOIDING
C     REDUNDANT CALLS TO THE TERRAIN AND RECEIVER MODELS. THE VALUES
C     RETURNED ARE THE FUNCTION VALUES FOR 'F' OR 'G' OR THEIR DERIVATI
C     (VIA 'DOT' ENTRY). UNNECESSARY CALLS ARE ELIMINATED THROUGH USE
C     OF 'TIME OF CALL' VARIABLES WHICH ARE COMPARED WITH THE CURRENT 
C     LAST CALL TIME MAINTAINED BY THE 'RKAM' PROGRAM. WHEN VALUES ARE
C     NOT CURRENT THEY ARE UPDATED BY CALLS TO THE APPROPRIATE ROUTINES
C

      implicit double precision (a-h,o-z)

      double precision d1mach                                     ! 10 Nov 2009
      DOUBLE PRECISION Z(*),PF(10),PG(10),pu(10) ! 3 May 2004
C     COMMON DECK "SS" INSERTED HERE
      DOUBLE PRECISION MODSURF
      COMMON/SS/ MODSURF(4)
      COMMON/SS/U,PUR,PURR,PURTH,PURPH
      COMMON/SS/PUTH,PUPH,PUTHTH,PUPHPH,PUTHPH,USELECT,UTIME
C     COMMON DECK "RKAM" INSERTED HERE
      DOUBLE PRECISION KR,KTH,KPH
      COMMON//R,TH,PH,KR,KTH,KPH,RKVARS(14),TPULSE,CSTEP,DRDT(20)
C     COMMON DECK "FNDER" INSERTED HERE 
      COMMON/FNDER/GY(20),NZ,NPZR,NPZRR,NPZRTH,NPZRPH,NPZTH 
      COMMON/FNDER/NPZPH,NPZTHTH,NPZPHPH,NPZTHPH,NSELECT,NTIME
C     COMMON DECK "RR" INSERTED HERE
      DOUBLE PRECISION MODREC
      COMMON/RR/ MODREC(4)
      COMMON/RR/F,PFR,PFRR,PFRTH,PFRPH
      COMMON/RR/PFTH,PFPH,PFTHTH,PFPHPH,PFTHPH,FSELECT,FTIME
C     COMMON DECK "GG" INSERTED HERE
      DOUBLE PRECISION MODG
      COMMON/GG/MODG(4)
      COMMON/GG/G,PGR,PGRR,PGRTH,PGRPH
      COMMON/GG/PGTH,PGPH,PGTHTH,PGPHPH,PGTHPH,GSELECT,GTIME
C     COMMON DECK "RKTIME" INSERTED HERE
      COMMON/CRKTIME/RKTIME
C
C     COMMON DECK "RMACH" INSERTED HERE 
C      COMMON/CRMACH/RMACH(5) ! 3 May 2004
      COMMON/CRMACH/RMACH(9) ! 3 May 2004
      EQUIVALENCE(RKTIME,IRKTIME)
C
C     COMMON DECK "WWR" INSERTED HERE
      PARAMETER (NWARSZ=1000) 
      COMMON/WW/ dum(10), MAXW,W(NWARSZ)
      DOUBLE PRECISION MAXSTP,MAXERR,INTYP,LLAT,LLON
      EQUIVALENCE (EARTHR,W(1)),(RAY,W(2)),(XMTRH,W(3)),(TLAT,W(4)),
     1 (TLON,W(5)),(OW,W(6)),(FBEG,W(7)),(FEND,W(8)),(FSTEP,W(9)),
     2 (AZ1,W(10)),(AZBEG,W(11)),(AZEND,W(12)),(AZSTEP,W(13)),
     3 (BETA,W(14)),(ELBEG,W(15)),(ELEND,W(16)),(ELSTEP,W(17)),
     8   (RUNSUP,W(18)),(RCVRH,W(20)),
     4 (ONLY,W(21)),(HOP,W(22)),(MAXSTP,W(23)),(PLAT,W(24)),(PLON,W(25))
     5,(HMAX,W(26)),(RAYFNC,W(29)),(EXTINC,W(33)),
     6 (HMIN,W(27)),(RGMAX,W(28)),
     8 (INTYP,W(41)),(MAXERR,W(42)),(ERATIO,W(43)),
     6 (STEP1,W(44)),(STPMAX,W(45)),(STPMIN,W(46)),(FACTR,W(47)),
     7 (SKIP,W(71)),(RAYSET,W(72)),(PRTSRP,W(74)),(HITLET,W(75))
     9 ,(BINRAY,W(76)),(PAGLN,W(77)),(PLT,W(81)),(PFACTR,W(82)),
     1 (LLAT,W(83)),(LLON,W(84)),(RLAT,W(85)),(RLON,W(86))
     2,(TIC,W(87)),(HB,W(88)),(HT,W(89)),(TICV,W(96))
C
      COMMON/CGET/ZERO
C
      DOUBLE PRECISION PI(10),PJ(10), PK(10) ! 3 May 2004
      COMMON/TMP1/NC,NSZR    ,NSZRR    ,NSZRTH
     1 ,NSZRPH    ,NSZTH    ,NSZPH    ,NSZTHTH    ,NSZPHPH    ,NSZTHPH
     2 ,NVELECT    ,NWIME
C
C      DATA PI,PJ/1.0,9*0.0,1.0,9*0.0/ ! 3 May 2004
      DATA PI,PJ,pk/1.d0,9*0.d0,1.d0,9*0.d0,1.d0,9*0.d0/ ! 3 May 2004
      DATA NC,NSZR,NSZRR,NSZRTH,NSZRPH,NSZTH,NSZPH,NSZTHTH,NSZPHPH
     1   ,NSZTHPH,NVELECT,NWIME/1,2,3,4,5,6,7,8,9,10,11,12/ 
      DATA FWIME,GWIME,FVELECT,GVELECT/2*-1.0,8HRECEIVER ,7HTERRAIN / 
      DATA UWIME,UVELECT/-1.0,7HSURFACE /
      DATA ZHRO/0.0/
C
      IENTRY=1
      GO TO 5
C
      ENTRY SETGET_test (z) ! 4 Jun 2004
C
      CALL RMOVE(PF,PI,10)
      CALL RMOVE(PG,PJ,10)
      call rmove(pu,pk,10) ! 3 May 2004
C      CALL IMOVE(NZ,NC,12) ! 3 May 2004
      FTIME=FWIME
      GTIME=GWIME
      FSELECT=FVELECT
      GSELECT=GVELECT
      UTIME=UWIME
      USELECT=UVELECT
      ZERO=ZHRO
      nz=1 ! 3 May 2004
      npzr=2 ! 3 May 2004
      npzrr=3 ! 3 May 2004
      npzrth=4 ! 3 May 2004
      npzrph=5 ! 3 May 2004
      npzth=6 ! 3 May 2004
      npzph=7 ! 3 May 2004
      npzthth=8 ! 3 May 2004
      npzphph=9 ! 3 May 2004
      npzthph=10 ! 3 May 2004
      nselect=11 ! 3 May 2004
      ntime=12 ! 3 May 2004
C
      RETURN
C
      ENTRY DOT_test(Z)
      IENTRY=2
5     CONTINUE
C
C     THE FOLLOWING THREE ASSIGNMENT STATEMENTS ARE NEEDED TO REPAIR
C     ARRAY ELEMENTS OVERWRITTEN ON THE CRAY VERSION OF THIS PROGRAM. 
C     THE ADDITIONAL ARRAY 'GY' HAS ALSO BEEN ADDED TO THE /FNDER/
C     LABELED COMMON TO ACT AS A BUFFER AREA.
C     NO CAUSE FOR THIS PROBLEM HAS YET BEEN FOUND. THESE STATEMENTS
C     ARE NOT NEEDED ON OTHER MACHINES AS FAR AS IS KNOWN.
      NZ=1
      NPZR=2
      NPZRR=3
C 26 Aug 2008      IF(ZERO.EQ.0.0) ZERO=EARTHR*RMACH(3)*2.0
      IF(ZERO.EQ.0.0) ZERO=EARTHR*D1MACH(3)*2.0                    ! 26 Aug 2008
C
      IF(ITEST(Z(NTIME)).EQ.IRKTIME) GO TO 10
C
      IF(Z(NSELECT).EQ.USELECT) CALL SURFACE
      IF(Z(NSELECT).EQ.FSELECT) CALL RECEVER
      IF(Z(NSELECT).EQ.GSELECT) CALL TOPOG
      Z(NTIME)=RKTIME
C     REMOVE MACHINE ROUND OFF NOISE FROM EXACT RECEIVER LOCATIONS
      IF(ABS(Z(NZ)).LE.ZERO) Z(NZ)=0.0
C
10    IF(IENTRY.NE.1) GO TO 20
      GET_test=Z(NZ)
      RETURN
C
20    GET_test=Z(NPZR)*DRDT(1)+Z(NPZTH)*DRDT(2)+Z(NPZPH)*DRDT(3) 
      RETURN
C ----------------------------------------------------------------------!-------
      END FUNCTION GET_test
C ----------------------------------------------------------------------!-------

C ----------------------------------------------------------------------!-------
      SUBROUTINE HEIGHT_test(H)
C ----------------------------------------------------------------------!-------
C     SETS REFERENCE CHARACTER HEIGHT IN INCHES(DISSPLA)
C     WRITE PARAMETERS TO GRAPHICS CALLS FILE
C
      implicit double precision (a-h,o-z)

c     COMMON/RAYDEV/NRYIND,NDEVTMP,NFRMAT,NDEVGRP,NDEVBIN,NDEVTXT ! 24 Mar 2006
      COMMON/RAYDEV/NRYIND,NDEVTMP,NFRMAT,NDEVTXT,NDEVBIN         ! 24 Mar 2006
C
      WRITE(NDEVTXT,100) 13,H
100   FORMAT(I2,F4.2)
C ----------------------------------------------------------------------!-------
      END SUBROUTINE HEIGHT_test
C ----------------------------------------------------------------------!-------
C ----------------------------------------------------------------------!-------
      FUNCTION ITEST_test(I)
C ----------------------------------------------------------------------!-------
C     USED TO PASS INTEGER VALUES THROUGH FOR VARIABLES TYPED REAL
C     (AS IN MIXED REAL/INTEGER ARRAYS)
      ITEST=I
C ----------------------------------------------------------------------!-------
      END FUNCTION ITEST_test
C ----------------------------------------------------------------------!-------
C ----------------------------------------------------------------------!-------
      FUNCTION ITOC_test(N)
C ----------------------------------------------------------------------!-------
C     RETURN 7 CHARACTER STRING REPRESENTATION OF INTEGER N
C     IF NUMBER IS TOO LARGE OR SMALL USE FLOATING POINT FORMAT
      CHARACTER ITOC*7
      IF(N.LT.-9999.OR.N.GT.99999) GO TO 100
      ITOC=' '
      WRITE(ITOC,'(I7)',ERR=100) N
      RETURN
100   WRITE(ITOC,'(2PG7.0)') FLOAT(N)
C ----------------------------------------------------------------------!-------
      END FUNCTION ITOC_test
C ----------------------------------------------------------------------!-------
