!################################################################## !################################################################## !###### ###### !###### SUBROUTINE CTR3D ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE ctr3d(b,x,y,z, x1,x2,dx,y1,y2,dy,z1,z2,dz, & 39,29 nx,ibgn,iend, ny,jbgn,jend, nz,kbgn,kend, & label,time,slicopt, kslice, jslice, islice, & n,xp,yp,axy2d,av2d,zp, runname, factor,tem1,tem2,tem3, & tem4,bb,tem5,hterain,pltopt) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set-up 2-d slices of a 3-d data array to contour with ! subroutine ctr2d. ! ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! 6/08/92 Added full documentation (K. Brewster) ! ! 12/25/1992 M. Xue and H. Jin ! Added capability to plot arbitary cross sections. ! ! 8/28/1994 M. Zou ! Added color shader to contour plot,add full documentation ! ! 3/25/96 (K. Brewster) ! Added variables isize,jsize,ksize ! !----------------------------------------------------------------------- ! ! INPUT: ! ! b 3-dimension array of variable ! x x-coord of scalar point (km) ! y y-coord of scalar point (km) ! z z-coord of scalar point in computation space (km) ! label character string describing the contents of a plot ! time model runing time step ! slicopt slice orientation indicator ! = 1, x-y slice of at k=kslice is plotted. ! = 2, x-z slice of at j=jslice is plotted. ! = 3, y-z slice of at i=islice is plotted. ! = 4, horizontal slice at z index islice is plotted. ! = 5, xy-z cross section of wind islice is plotted. ! = 6, data field on constant p-level is plotted. ! = 0, all of the three slices above are plotted. ! axy2d 2d x-y array ! av2d 2D array for the vertical slice ! xp x-coordinate of grid points on arbitary vertical ! cross-section ! yp y-coordinate of grid points on arbitary vertical ! cross-section ! zp z-coordinate of grid points on arbitary vertical ! cross-section ! runname character string describing the model run ! factor scaling factor ! hterain 2-D terrain data for contour ! trnplt flag to plot terrain (0/1) ! WORK ARRAY: ! tem1 Temporary work array. ! tem2 Temporary work array. ! tem3 Temporary work array. ! tem4 Temporary work array. ! tem5 Temporary work array. ! ! (These arrays are defined and used locally (i.e. inside this ! subroutine), they may also be passed into routines called by ! this one. Exiting the call to this subroutine, these temporary ! work arrays may be used for other purposes therefore their ! contents overwritten. Please examine the usage of work arrays ! before you alter the code.) ! ! pp01 The pressure (mb) value at the specific p-level ! ercpl reciprocal of exponent ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! INTEGER :: nx,ny,nz INTEGER :: n REAL :: b(nx,ny,nz) REAL :: x(nx,ny,nz) REAL :: y(nx,ny,nz) REAL :: z(nx,ny,nz) ! REAL :: axy2d(nx,ny) REAL :: av2d(n,nz),zp(n,nz) REAL :: xp(n),yp(n) REAL :: x1,x2,dx,y1,y2,dy,z1,z2,dz INTEGER :: ibgn,iend,jbgn,jend,kbgn,kend,length ! CHARACTER (LEN=6) :: timhms CHARACTER (LEN=*) :: label ! REAL :: time ! INTEGER :: slicopt INTEGER :: kslice,jslice,islice ! CHARACTER (LEN=*) :: runname ! REAL :: factor ! INTEGER :: trnplt ! plot terrain option (0/1/2/3) INTEGER :: pltopt ! plot variable option (0/1/2/3) REAL :: hterain(nx,ny) ! The height of the terrain. ! REAL :: tem1(*) REAL :: tem2(*) REAL :: tem3(*) REAL :: tem4(*) REAL :: bb(nx,ny,nz) REAL :: tem5(*) ! !----------------------------------------------------------------------- ! ! Some constants ! !----------------------------------------------------------------------- ! REAL :: pp01, ercpl PARAMETER (ercpl=0.3678794) ! exp(-1.0) ! !----------------------------------------------------------------------- ! ! Common blocks for plotting control parameters ! !----------------------------------------------------------------------- ! REAL :: x01,y01 ! the first point of interpolation REAL :: x02,y02 ! the second point of interpolation REAL :: zlevel ! the given height of the slice REAL :: sinaf,cosaf,dist,sqrtdxy COMMON /slicev/x01,y01,x02,y02,sinaf,cosaf,dist,sqrtdxy COMMON /sliceh/zlevel INTEGER :: ovrtrn ! overlay terrain option (0/1) REAL :: trninc,trnmin, trnmax ! terrain interval minimum, maximum REAL :: ztmin,ztmax COMMON /trnpar/ trnplt,ovrtrn,trninc,trnmin, trnmax,ztmin,ztmax INTEGER :: smooth COMMON /smoothopt/smooth ! INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat CHARACTER (LEN=4) :: stem2 CHARACTER (LEN=1) :: stem1 REAL :: x_tmp COMMON /tmphc2/ x_tmp REAL :: tmpx, tmpy CHARACTER (LEN=20) :: distc REAL :: x101, y101, x102,y102 COMMON /slicev1/x101, y101, x102,y102 ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! INTEGER :: i,j,k,ij,ik,jk,isize,jsize,ksize, llabel CHARACTER (LEN=120) :: label_copy CHARACTER (LEN=120) :: title INTEGER :: wrtflag CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab CHARACTER (LEN=25) :: timestring COMMON /timelev/wrtflag, timelab, levlab, timestring ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! isize = (iend-ibgn)+1 jsize = (jend-jbgn)+1 ksize = (kend-kbgn)+1 label_copy = label llabel = 120 CALL xstrlnth(label_copy, llabel) CALL xpscmnt('Start plotting '//label_copy(1:llabel)) ! !----------------------------------------------------------------------- ! ! Set up terrain, if needed. ! !----------------------------------------------------------------------- ! IF(trnplt == 1 .OR.trnplt == 2 .OR. ovrtrn == 1) THEN ! IF(ovrtrn.eq.1) THEN DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem4(ij)=hterain(i,j) END DO END DO END IF CALL cvttim( time, timhms ) IF( timhms(1:1) == '0' ) timhms(1:1)=' ' WRITE(timelab,'(''T='',F8.1,A)') time, & ' s ('//timhms(1:2)//':'//timhms(3:4)//':'//timhms(5:6)//')' CALL get_time_string ( time, timestring) IF ( slicopt == 2 .OR. slicopt == 3 .OR. slicopt == 5 ) THEN CALL cal_dist(haxisu,dx,dy,x01,y01,x02,y02, & slicopt,tmpx,tmpy,distc) END IF IF(slicopt == 1 .OR. slicopt == 0 ) THEN k = kslice DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij) = -9999.0 IF( tem2(ij)=x(i,j,k) tem3(ij)=y(i,j,k) END DO END DO IF (k /= 2) THEN WRITE(levlab,'(''GRID LEVEL='',I3)')k WRITE(title,'(a)') label ELSE WRITE(levlab,'(''FIRST LEVEL ABOVE GROUND (SURFACE)'')') WRITE(title,'(a)') label END IF length = 120 CALL strlnth( title, length) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem5) END DO CALL ctr2d(tem1,tem2,tem3, x1,x2,dx, y1,y2,dy, & isize,jsize,title(1:length),runname, & tem4,slicopt,pltopt) ! !----------------------------------------------------------------------- ! ! slicopt=2 Plot x-z cross-section ! !----------------------------------------------------------------------- ! ELSE IF (slicopt == 2 .OR. slicopt == 0 ) THEN x_tmp = y(1,jslice,1) j = jslice DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 IF( tem2(ik)=x(i,j,k) tem3(ik)=z(i,j,k) END DO END DO dist = (j-1.5)*tmpy length=20 CALL strmin ( distc, length) WRITE(levlab,'(''X-Z PLANE AT Y='',F8.1,A)')dist,distc(1:length) WRITE(title,'( a )') label length = 120 CALL strlnth( title, length) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,ksize,1,isize,1,ksize,tem5) END DO CALL ctr2d(tem1,tem2,tem3, x1,x2,dx, z1,z2,dz, & isize,ksize,title(1:length),runname, & tem4,slicopt,pltopt) ! !----------------------------------------------------------------------- ! ! slicopt=3 Plot y-z cross-section ! !----------------------------------------------------------------------- ! ELSE IF ( slicopt == 3 .OR. slicopt == 0) THEN x_tmp = x(islice,1,1) i = islice DO k=kbgn,kend DO j=jbgn,jend jk = j-jbgn+1 + (k-kbgn)*jsize tem1(jk) = -9999.0 IF( tem2(jk)=y(i,j,k) tem3(jk)=z(i,j,k) END DO END DO dist = (i-1.5)*tmpx length=20 CALL strmin ( distc, length) WRITE(levlab,'(''Y-Z PLANE AT X='',F8.1,A)')dist,distc(1:length) WRITE(title,'( a )' ) label length = 120 CALL strlnth( title, length) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,jsize,ksize,1,jsize,1,ksize,tem5) END DO CALL ctr2d(tem1,tem2,tem3, y1,y2,dy, z1,z2,dz, & jsize,ksize,title(1:length),runname, & tem4,slicopt,pltopt) ! !----------------------------------------------------------------------- ! ! slicopt=4 Plot horizontal slice at given height ! slicopt=6 Plot constant pressure slice at given pressure(mb) ! slicopt=7 Plot isentropic surfaces ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 4.OR.slicopt == 6.OR.slicopt == 7 ) THEN DO k=kbgn,kend DO j=jbgn,jend DO i=ibgn,iend bb(i,j,k) = -9999.0 IF( END DO END DO END DO CALL hintrp1(nx,ny,nz,kbgn,kend,bb,z,zlevel,axy2d) DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij)=axy2d(i,j) tem2(ij)=x(i,j,2) tem3(ij)=y(i,j,2) END DO END DO IF(slicopt == 4) THEN WRITE(levlab,'(''Z='',F7.3,'' KM MSL'')') & zlevel ELSE IF(slicopt == 6) THEN pp01 = 0.01*ercpl**zlevel WRITE(levlab,'(''P='',F7.2,A)') pp01, ' MB' ELSE WRITE(levlab,'(''THETA='',F5.1,A)') zlevel, ' (K)' END IF WRITE(title,'(a)') label length = 120 CALL strlnth( title, length) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem5) END DO CALL ctr2d(tem1,tem2,tem3, x1,x2,dx, y1,y2,dy, & isize,jsize,title(1:length),runname, & tem4,slicopt,pltopt) ! !----------------------------------------------------------------------- ! ! slicopt=5 Plot vectical slice through two given points ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 5 ) THEN CALL sectvrt(nx,ny,nz,b,x,y,z,av2d,zp,n,xp,yp) DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 IF(av2d(i,k) /= -9999.0) tem1(ik)=av2d(i,k)*factor tem2(ik)=x1+(i-ibgn)* sqrtdxy tem3(ik)=zp(i,k) END DO END DO IF(xfmat == -1) THEN length=20 CALL strmin ( distc, length) WRITE(levlab, & '(''VERTICAL PLANE FROM '',4(A,F8.1),A,A)') & '(',x101,',',y101,') through (',x102,',',y102,') ', & distc(1:length) WRITE(title,'(a)') label WRITE(title,'(a)') label ELSE IF(xfmat == 0) THEN length=20 CALL strmin ( distc, length) WRITE(levlab, & '(''VERTICAL PLANE FROM '',4(A,I5),A,A)') & '(',INT(x101),',',INT(y101),') through (',INT(x102),',' & ,INT(y102),') ', distc(1:length) WRITE(title,'(a)') label ELSE WRITE(title,'(a)') label WRITE(stem1,'(i1)')xfmat WRITE(stem2,43)stem1 43 FORMAT('f8.',a1) WRITE(title,'(a)') label END IF length = 120 CALL strlnth( title, length) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,ksize,1,isize,1,ksize,tem5) END DO CALL ctr2d(tem1,tem2,tem3, x1,x2,sqrtdxy, z1,z2,dz, & isize,ksize,title(1:length),runname, & tem4,slicopt,pltopt) ! END IF CALL xpscmnt('End plotting '//label_copy(1:llabel)) RETURN END SUBROUTINE ctr3d ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE CTR2D ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE ctr2d(a,x,y,xl,xr,dx,yb,yt,dy,m,n,title,runname, & 6,26 hterain,slicopt,pltopt) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Generate contour plots of 2-d field A given its coordinates ! using ZXPLOT package.. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! ! 6/08/92 (K. Brewster) ! Added full documentation. ! ! 8/28/94 (M. Zou) ! Added color routing , overlay terrain. ! ! 1/24/96 (J. Zong and M. Xue) ! Fixed a problem related to finding the minimum and maximum of the ! 2D array, a, when there exist missing data. Initial min. and max. ! should be set to values other than the missing value, -9999.0. ! !----------------------------------------------------------------------- ! ! INPUT: ! ! a 2-dimensional slice of data to contour ! ! x x coordinate of grid points in plot space (over on page) ! y y coordinate of grid points in plot space (up on page) ! ! ! xl Left bound of the physical domain ! xr Right bound of the physical domain ! dx Spacing between x-axis tick marks ! yb Bottom bound of the physical domain. ! yt Top bound of the physical domain. ! dy Spacing between y-axis tick marks ! ! m first dimension of a ! n second dimension of a ! ! title character string describing the contents of a ! runname character string describing the model run ! ! hterain 2-D terrain data to contour ! slicopt slice orientation indicator ! = 1, x-y slice of at k=kslice is plotted. ! = 2, x-z slice of at j=jslice is plotted. ! = 3, y-z slice of at i=islice is plotted. ! = 4, horizontal slice at z index islice is plotted. ! = 5, xy-z cross section of wind islice is plotted. ! = 6, data field on constant p-level is plotted. ! = 0, all of the three slices above are plotted. ! plot variable plot option (0/1/2/3) ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INCLUDE 'arpsplt.inc' ! ! INTEGER :: m,n ! REAL :: a(m,n) REAL :: x(m,n) REAL :: y(m,n) REAL :: xl,xr,dx,yb,yt,dy REAL :: hterain(m,n) ! CHARACTER (LEN=*) :: runname CHARACTER (LEN=*) :: title ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: layover COMMON /laypar/ layover INTEGER :: ovrobs,obsset,obscol,obs_marktyp REAL :: obs_marksz COMMON /obspar/ ovrobs,obsset,obscol,obs_marktyp, obs_marksz ! INTEGER :: ovrstaopt INTEGER :: ovrstam,staset,ovrstan,ovrstav,stacol,markprio,wrtstax INTEGER :: nsta_typ,sta_typ(10),sta_marktyp(10),sta_markcol(10) REAL :: sta_marksz(10),wrtstad COMMON /sta_par/ ovrstaopt,ovrstam,staset,ovrstan,ovrstav,stacol, & markprio,nsta_typ,sta_typ,sta_marktyp, & sta_markcol,sta_marksz,wrtstax,wrtstad ! REAL :: ctinc,ctmin,ctmax,vtunt ! contour interval and vector unit COMMON /incunt/ ctinc,ctmin,ctmax,vtunt ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! INTEGER :: flag INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat ! REAL :: yxratio COMMON /yratio/ yxratio ! the scaling factor the y/x ratio. ! INTEGER :: ntitle,titcol, nxpic, nypic, wpltime REAL :: titsiz CHARACTER (LEN=132) :: ptitle(3), footer_l, footer_c, footer_r COMMON /titpar1/ptitle, footer_l, footer_c, footer_r COMMON /titpar2/ntitle,titcol,wpltime, nxpic, nypic COMMON /titpar3/titsiz ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! CHARACTER (LEN=150) :: ch1 CHARACTER (LEN=150) :: ch INTEGER, ALLOCATABLE :: iwrk(:) REAL , ALLOCATABLE :: xwk(:),ywk(:) INTEGER :: istatus INTEGER :: i,j REAL :: cl(500) ! contour levels REAL :: pl,pr,pb,pt ! plot space left, right, bottom, top coordinate REAL :: px,py ! plot space left-right length and up-down height REAL :: pxc,pyc ! plot space left-right center and ! up-down center REAL :: xs,ys ! real space left-right length and up-down height REAL :: zinc ! contour interval REAL :: zmin,zmax ! max and min of data array INTEGER :: ncl,slicopt,mode1 REAL :: zlevel COMMON/sliceh/zlevel INTEGER :: pltopt ! variavle plot option (0/1/2/3) INTEGER :: timeovr COMMON /timover/ timeovr REAL :: lblmag, ctrlbsiz, axlbsiz COMMON /labmag/ lblmag, ctrlbsiz, axlbsiz REAL :: xfinc ! INTEGER :: col_table,pcolbar COMMON /coltable/col_table,pcolbar ! INTEGER :: LEN,len1 ! CHARACTER (LEN=12) :: varname COMMON /varplt1/ varname ! CHARACTER (LEN=150) :: f_ch ! INTEGER :: setcontopt, setcontnum REAL :: setconts(maxunevm,maxuneva) COMMON /setcon_par/setcontopt,setcontnum,setconts INTEGER :: ncont REAL :: tcont(maxunevm) INTEGER :: wrtflag CHARACTER (LEN=25) :: timestring CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab COMMON /timelev/wrtflag,timelab, levlab, timestring CHARACTER (LEN=80) :: prestr INTEGER :: preflag COMMON /preinfo/ prestr,preflag REAL :: x101, y101, x102,y102 COMMON /slicev1/x101, y101, x102,y102 REAL :: ytmp !! local temporary variable REAL :: f_cputime,f_walltime,cpu1, cpu2,second1, second2 REAL :: hatch_angle ! INTEGER :: missval_colind, missfill_opt ! miss value color index COMMON /multi_value/ missfill_opt,missval_colind INTEGER :: missfill DATA missfill/0/ INTEGER :: mxset INTEGER :: xnwpic_called COMMON /callnwpic/xnwpic_called INTEGER :: iclfrq INTEGER :: ctrlbfrq COMMON /clb_frq/ ctrlbfrq ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! !----------------------------------------------------------------------- ! ! Check for adequate room in work array ! !----------------------------------------------------------------------- ! second1= f_walltime() cpu1 = f_cputime() ncont = 0 ncl = 1 ALLOCATE(iwrk(m*n),STAT=istatus) ALLOCATE(xwk(m*n),ywk(m*n),STAT=istatus) WRITE(6,'(//a,a)') 'Plotting ',title IF( layover == 0 .OR. xnwpic_called == 0) THEN CALL xnwpic xnwpic_called = 1 timeovr=0 ! set overlayer terrain agian wrtflag = 0 ! preflag = 0 prestr = levlab len1=80 CALL strmin(prestr,len1) ELSE timeovr=1 wrtflag = wrtflag + 1 END IF ! !----------------------------------------------------------------------- ! ! Get plotting space variables ! !----------------------------------------------------------------------- ! CALL xqpspc( pl, pr, pb, pt) px = pr - pl py = pt - pb pxc = (pr+pl)/2 pyc = (pb+pt)/2 xs = xr-xl ys = yt-yb ! !----------------------------------------------------------------------- ! ! Let the longest lenth determine size scaling of plot ! !----------------------------------------------------------------------- ! IF( py/px >= (ys*yxratio)/xs ) THEN py = (ys*yxratio)/xs*px CALL xpspac(pl, pr, pyc-py/2, pyc+py/2 ) ELSE px = xs/(ys*yxratio)*py CALL xpspac(pxc-px/2, pxc+px/2, pb, pt) END IF ! !----------------------------------------------------------------------- ! ! Set the real distance to plot distance scaling ! !----------------------------------------------------------------------- ! CALL xmap( xl, xr, yb, yt ) ! CALL xlbsiz( ctrlbsiz*(yt-yb)*lblmag ) ! !----------------------------------------------------------------------- ! ! Find max and min of data array ! !----------------------------------------------------------------------- ! mxset = 0 missfill = 0 DO j=1,n DO i=1,m IF(ABS( missfill = 1 CYCLE END IF IF( mxset == 0) THEN zmax= a(i,j) zmin= a(i,j) mxset = 1 ELSE zmax= MAX (zmax,a(i,j)) zmin= MIN (zmin,a(i,j)) END IF END DO END DO IF (missfill == 1 .AND. missfill_opt == 1) & CALL fillmissval ( m,n,xl, xr, yb,yt ) CALL xcolor(lbcolor) ! !----------------------------------------------------------------------- ! ! Find proper contour interval and then contour field ! using ZXPLOT routine xconta ! !----------------------------------------------------------------------- cl(1)=0.0 IF( zmax-zmin > 1.0E-20 ) THEN ! !----------------------------------------------------------------------- ! ! Check to see if user defined contour levels is available for the ! current variable. ! !----------------------------------------------------------------------- CALL get_contour (ncont, tcont) iclfrq = ctrlbfrq IF(setcontopt > 0 .AND.ncont > 0) THEN ch1(1:11)=' contours: ' LEN=11 DO i =1,ncont CALL xrch1(tcont(i),f_ch,len1) WRITE(ch1,'(a,a,'' '')')ch1(1:LEN), f_ch(1:len1) LEN=LEN+len1+1 END DO DO i=1,ncont cl(i)=tcont(i) END DO ncl = ncont mode1 = 4 iclfrq = 1 GO TO 150 END IF IF( ctinc == 0.0) THEN cl(2)=cl(1)+ xfinc(zmax-zmin)/2 IF(cl(2)-cl(1) == 0.0) cl(2)=cl(1)+1.0 mode1=1 CALL xnctrs( 8,20) ELSE IF ( ctinc == -9999.) THEN CALL xnctrs( 8,20) CALL set_interval(a, m,n,zmin,zmax,ctmin,ctmax,cl) ctinc = cl(2)-cl(1) zinc = ctinc mode1=1 ELSE cl(2)=cl(1)+ctinc CALL xnctrs(1,300) mode1=1 END IF 150 CONTINUE zinc = cl(2)-cl(1) !----------------------------------------------------------------------- ! ! Plot contour or color filled contour fields ! !----------------------------------------------------------------------- CALL xwindw(xl, xr, yb, yt) CALL xctrlim(ctmin,ctmax) CALL xclfrq(iclfrq) IF(pltopt == 1) THEN CALL xctrclr(icolor, icolor) CALL xconta(a,x,y,iwrk,m,m,n,cl,ncl,mode1) ELSE IF( pltopt == 2) THEN CALL xctrclr(icolor, icolor1) CALL xcolfil(a,x,y,iwrk,xwk,ywk,m,m,n,cl,ncl,mode1) CALL xchsiz(0.025*(yt-yb)) CALL xcpalet(pcolbar) ELSE IF(pltopt == 4) THEN CALL xctrclr(icolor, icolor1) CALL xconta(a,x,y,iwrk,m,m,n,cl,ncl,mode1) ELSE IF(pltopt == 5) THEN CALL xctrclr(icolor, icolor1) CALL xcolfil(a,x,y,iwrk,xwk,ywk,m,m,n,cl,ncl,mode1) CALL xchsiz(0.025*(yt-yb)) CALL xcpalet(pcolbar) CALL xctrclr(lbcolor, lbcolor) CALL xconta(a,x,y,iwrk,m,m,n,cl,ncl,mode1) ELSE IF(pltopt == 6) THEN CALL xctrclr(icolor, icolor) CALL xdhtch(0.003) CALL xctrclr(icolor, icolor) ncl = 2 mode1 = 4 cl(1) = ctmin cl(2) = ctmax CALL xclfrq(1) CALL xhilit(0) CALL xconta(a,x,y,iwrk,m,m,n,cl,ncl,mode1) CALL xhilit(1) hatch_angle = 45.0 CALL xdhtch(0.004) CALL xhatcha(a,x,y,xwk,ywk,m,m,n,ctmin,1.0E10,hatch_angle) CALL xdhtch(0.002) CALL xhatcha(a,x,y,xwk,ywk,m,m,n,ctmax,1.0E10,hatch_angle) END IF CALL xclfrq(2) ELSE cl(2)=1.0 ncl=2 END IF IF(ctinc == 0.0) THEN zinc = cl(2) - cl(1) ELSE zinc = ctinc END IF ! !----------------------------------------------------------------------- ! ! Plot map, boxes and polygons. ! !----------------------------------------------------------------------- ! CALL pltextra(slicopt, pltopt) !----------------------------------------------------------------------- ! ! Terrain outline in vertical slices. ! !----------------------------------------------------------------------- IF(slicopt == 2 .OR. slicopt == 3 .OR.slicopt == 5) THEN CALL xcolor(trcolor) CALL xthick(3) CALL xpenup( x(1,1), y(1,1)-0.5*(y(1,2)-y(1,1)) ) DO i=2,m CALL xpendn( x(i,1), y(i,1)-0.5*(y(i,2)-y(i,1)) ) END DO CALL xthick(1) END IF ! !----------------------------------------------------------------------- ! ! Overlay terrain contour if required in x-y level ! or Plot terrain outline in slice zlevel ! !----------------------------------------------------------------------- ! IF ( timeovr == 0 ) CALL plttrn(hterain,x,y,m,n,slicopt) ! !----------------------------------------------------------------------- ! ! Plot observations ! !----------------------------------------------------------------------- ! IF(obsset == 1 .AND. ovrobs == 1 .AND. (slicopt == 1.OR.slicopt == 4.OR. & slicopt == 6.OR.slicopt == 7.OR.slicopt == 8) ) THEN CALL xchsiz(0.025*ys * lblmag) CALL pltobs(1) obsset=0 END IF ! !----------------------------------------------------------------------- ! ! Plot station labels ! !----------------------------------------------------------------------- ! IF( ovrstaopt == 1 .AND. wrtstax == 1 & .AND. (timeovr == 0 .OR. (timeovr == 1 .AND. pltopt == 2)) & .AND.(slicopt == 2.OR.slicopt == 3.OR. slicopt == 5) ) THEN CALL xchsiz(0.025*ys * lblmag) flag=1 CALL pltsta(a,a,x,y,m,n,flag,slicopt) END IF IF(ovrstaopt == 1 .AND. staset == 1 .AND. & (ovrstam == 1.OR.ovrstan == 1.OR.ovrstav == 1) & .AND.(slicopt == 1.OR.slicopt == 4 & .OR.slicopt == 6.OR.slicopt == 7.OR.slicopt == 8) & .AND.(timeovr == 0.OR.(timeovr == 1.AND.pltopt == 2))) THEN CALL xchsiz(0.025*ys * lblmag) flag=0 CALL pltsta(a,a,x,y,m,n,flag,slicopt) ! staset=0 END IF CALL xwdwof ! !----------------------------------------------------------------------- ! ! Plot axes with tick marks ! !----------------------------------------------------------------------- ! CALL pltaxes(slicopt,dx,dy) ! IF(ntitle > 0 .AND. nxpic == 1 .AND. nypic == 1 & .AND. timeovr == 0 ) THEN CALL xcolor(titcol) CALL xchsiz(0.025*ys * titsiz) DO i=1,ntitle LEN=132 CALL strlnth(ptitle(i),LEN) CALL xchori(0.) CALL xcharc( xl+xs/2,yt+(0.25-(i-1)*0.06)*ys,ptitle(i)(1:LEN)) END DO CALL xcolor(lbcolor) END IF CALL xchsiz( 0.030*ys * lblmag ) ! plot time and level label IF ( layover < 1 ) THEN IF(levlab /= ' ') THEN len1=80 CALL strmin(levlab,len1) CALL xcharc((xl+xr)*0.5,yt+0.015*ys, levlab(1:len1)) preflag = 1 END IF len1=50 CALL strmin(timelab,len1) CALL xcharc((xl+xr)*0.5,yt+0.06*ys, & timestring(1:25)//' '//timelab(1:len1)) END IF IF(preflag == 0 .AND. levlab /= ' ') THEN len1=80 CALL strmin(levlab,len1) CALL xcharc((xl+xr)*0.5,yt+0.015*ys, levlab(1:len1)) preflag = 1 END IF LEN = 120 CALL strmin(title, LEN) IF( title(LEN-1:LEN-1) == ')' ) THEN title(LEN-1:LEN-1)=',' ELSE title(LEN:LEN+1)=' (' LEN = LEN+1 END IF IF(pltopt == 2) THEN WRITE(f_ch, '(a,''SHADED)'')')title(1:LEN) ELSE IF( pltopt == 5 ) THEN WRITE(f_ch, '(a,''SHADED/CONTOUR)'')')title(1:LEN) ELSE WRITE(f_ch, '(a,''CONTOUR)'')')title(1:LEN) END IF ! if first levlab is not equal second levlab then attatch levlab on f_ch LEN=120 CALL strmin(f_ch, LEN) len1=80 CALL strmin(levlab,len1) ! !mx ! IF (preflag.eq.1 .and. prestr(1:len1).ne.levlab(1:len1) ! : .and. prestr(1:1).ne.' ' ! : .and.layover.ne.0 .and. levlab(1:1).ne.' ') THEN ! write(f_ch,'(a,a)') f_ch(1:len),levlab(1:len1) ! ENDIF ! IF(pltopt.eq.1 .and. layover.ge.1 ) CALL xcolor(icolor) IF(pltopt == 1) CALL xcolor(icolor) ! plot variable name IF (preflag == 1 .AND. prestr /= levlab .AND. prestr /= ' ' & .AND.layover /= 0 .AND. levlab /= ' ') THEN CALL xchsiz( 0.018*ys * lblmag ) ELSE CALL xchsiz( 0.028*ys * lblmag ) END IF ! save for next time use IF(prestr(1:1) == ' ' .AND.layover /= 0 ) prestr=levlab IF(lbaxis == 1 ) THEN IF( wrtstax == 0) THEN ytmp = 0.08 ELSE ytmp = 0.14 END IF ELSE ytmp = 0.12 END IF LEN=150 CALL strmin(f_ch,LEN) CALL xchsiz( 0.025*ys * lblmag ) !mx CALL xcolor(lbcolor) CALL xcharl(xl-0.20*(xr-xl), yb-(ytmp+wrtflag*0.030)*ys , & f_ch(1:LEN)) IF (pltopt == 1.OR.pltopt == 3.OR.pltopt == 4.OR.pltopt == 5)THEN IF(ABS(zmin-zmax) <= 1.e-15 .OR. ncont > 0) THEN WRITE(ch,'(''MIN='',G9.3E2,'' MAX='',G9.3E2)')zmin,zmax ELSE WRITE(ch,'(''MIN='',G10.4E2,'' MAX='',G10.4E2, & & '' inc='',g10.4E2)')zmin,zmax,zinc END IF ELSE IF( pltopt == 2 ) THEN WRITE(ch,'(''MIN='',G9.3E2,'' MAX='',G9.3E2)')zmin,zmax END IF LEN=150 CALL strmin(ch,LEN) CALL xcharr(xr+0.20*(xr-xl), yb-(ytmp+wrtflag*0.030)*ys , & ch(1:LEN)) IF (ncont > 1 .AND. (pltopt == 1 .OR. pltopt == 4) ) THEN wrtflag = wrtflag+1 len1=150 CALL strmin(ch1,len1) CALL xcharr(xr+0.20*(xr-xl),yb-(ytmp+wrtflag*0.030)*ys, & ch1(1:len1)) END IF !----------------------------------------------------------------------- ! ! Plot additional text below the figure ! !----------------------------------------------------------------------- CALL label2d(runname) CALL xpspac(pl, pr, pb, pt) ! set frame back cpu2 = f_cputime() second2 = f_walltime() ! write(6,*)'!!!! total cpu time for one CTR2D :', ! : cpu2-cpu1,' PLOT:',varname DEALLOCATE(iwrk) DEALLOCATE(xwk,ywk) RETURN END SUBROUTINE ctr2d ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE CTRINC ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE ctrinc( ctinc0, ctmin0, ctmax0 ) 3 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set the contour interval for field to plotted by CTR2D. ! ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! 6/08/92 Added full documentation (K. Brewster) ! !----------------------------------------------------------------------- ! ! INPUT: ! ! ctinc0 Contour interval ! If CTINC0 = 0.0, the interval is internally determined. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE REAL :: ctinc0,ctmin0,ctmax0 ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! REAL :: ctinc,ctmin,ctmax,vtunt ! contour interval and vector unit COMMON /incunt/ ctinc,ctmin,ctmax,vtunt ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ctinc = ctinc0 ctmin = ctmin0 ctmax = ctmax0 RETURN END SUBROUTINE ctrinc ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE LABEL2D ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE label2d(runname) 2,3 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! Plot certain text labels for VTR2D and CTR2d. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! Taked from former CTR2D. ! !----------------------------------------------------------------------- ! ! INPUT: ! ! xl Left bound of the physical domain ! xr Right bound of the physical domain ! yb Bottom bound of the physical domain. ! yt Top bound of the physical domain. ! ! runname character string describing the model run ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE CHARACTER (LEN=*) :: runname INTEGER :: layover COMMON /laypar/ layover INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat ! INTEGER :: ntitle,titcol, nxpic, nypic, wpltime REAL :: titsiz CHARACTER (LEN=132) :: ptitle(3), footer_l, footer_c, footer_r COMMON /titpar1/ptitle, footer_l, footer_c, footer_r COMMON /titpar2/ntitle,titcol,wpltime, nxpic, nypic COMMON /titpar3/titsiz REAL :: xl,xr,yb,yt ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! INTEGER :: nopic REAL :: xlimit, ylimit, rotang INTEGER :: nhpic, nvpic,ifont INTEGER :: ovrtrn ,trnplt ! overlay terrain option (0/1) REAL :: trninc,trnmin, trnmax ! terrain interval minimum, maximum REAL :: ztmin,ztmax COMMON /trnpar/ trnplt,ovrtrn,trninc,trnmin, trnmax,ztmin,ztmax INTEGER :: timeovr COMMON /timover/ timeovr REAL :: lblmag, ctrlbsiz, axlbsiz COMMON /labmag/ lblmag, ctrlbsiz, axlbsiz ! INTEGER :: col_table,pcolbar COMMON /coltable/col_table,pcolbar ! CHARACTER (LEN=24) :: tzstring CHARACTER (LEN=24) :: tz CHARACTER (LEN=132) :: datetimestr INTEGER :: lnblnk, len1, len2, len3 CHARACTER (LEN=120) :: string_l, string_c, string_r CHARACTER (LEN=8) :: tzone CHARACTER (LEN=10) :: cur_time CHARACTER (LEN=8) :: cur_date INTEGER :: t_values(8) REAL :: ytmp, hch ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! CALL xqmap(xl,xr,yb,yt) CALL xcolor(lbcolor) CALL xqnpic(nopic) CALL xqspac(nhpic, nvpic, rotang, xlimit, ylimit) CALL xchsiz( 0.021*(yt-yb) * lblmag ) IF(timeovr == 0) THEN IF(nopic == nhpic*(nvpic-1)+1 ) THEN IF ( wpltime == 1) THEN CALL date_and_time(cur_date,cur_time,tzone,t_values) IF(t_values(4) == 0) THEN tzstring = ' UTC' ELSE tzstring = ' Local Time' END IF WRITE (datetimestr,999) 'Plotted ', & t_values(1),t_values(2),t_values(3), & t_values(5),t_values(6),tzstring 999 FORMAT (a, i4.4,'/',i2.2,'/',i2.2,' ',i2.2,':',i2.2,a) END IF IF ( footer_l == ' ') THEN string_l = 'ARPSPLT/ZXPLOT ' ELSE string_l = footer_l END IF IF( footer_c == ' ') THEN string_c = runname ELSE string_c = footer_c END IF IF(wpltime == 1 ) THEN string_r = datetimestr(:lnblnk(datetimestr)) ELSE string_r = footer_r END IF CALL xqcfnt(ifont) CALL xcfont(xfont) ytmp = 0.29 CALL xqchsz(hch) IF ( layover < 1) THEN len1=120 CALL strmin(string_l, len1) len2=120 CALL strmin(string_c, len2) len3=120 CALL strmin(string_r, len3) CALL xcharc(xl+0.5*(xr-xl), & yb-(ytmp+layover*0.03)*(yt-yb), & string_l(1:len1)//' '//string_c(1:len2)//' '// & string_r(1:len3)) END IF CALL xcfont(ifont) END IF timeovr=1 END IF RETURN END SUBROUTINE label2d ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE VTR3D ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE vtr3d(u,v,w, x,y,z, xw,xe,dx, ys,yn,dy, zb,zt,dz, & 3,38 nx,ibgn,iend,ist, ny,jbgn,jend,jst, nz,kbgn,kend,kst, & kslice, jslice, islice, label,time, runname, factor, & slicopt,n,xp,yp,zp,u1,v1,u2,v2,w2, & tem1,tem2,tem3,tem4, & tem5,tem6,hterain) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Plot vector fields in 2-d slices ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! 6/08/92 Added full documentation (K. Brewster) ! ! 3/25/96 (K. Brewster) ! Added variables isize,jsize,ksize ! !----------------------------------------------------------------------- ! ! INPUT: ! ! u 3-dimensional array of u wind components (m/s) ! v 3-dimensional array of v wind components (m/s) ! w 3-dimensional array of w wind components (m/s) ! ! x x coordinate of grid points in physical/comp. space (m) ! y y coordinate of grid points in physical/comp. space (m) ! z z coordinate of grid points in physical space (m) ! ! xw value of x for first i grid point to plot ! xe value of x for last i grid point to plot ! ys value of y for first j grid point to plot ! yn value of y for last j grid point to plot ! zb value of z for first k grid point to plot ! zt value of z for last k grid point to plot ! ! nx first dimension of b ! ibgn index of first i grid point to plot ! iend index of last i grid point to plot ! ! ny second dimension of b ! jbgn index of first j grid point to plot ! jend index of last j grid point to plot ! ! nz third dimension of b ! kbgn index of first k grid point to plot ! kend index of last k grid point to plot ! ! ist step size in x direction ! jst step size in y direction ! kst step size in z direction ! ! time time of data in seconds ! ! kslice k index of plane for slicopt=1 x-y slice ! jslice j index of plane for slicopt=2 x-z slice ! islice i index of plane for slicopt=1 y-z slice ! ! runname character string decribing run ! ! factor scaling factor for winds ! V*factor wind vectors are plotted ! ! slicopt slice orientation indicator ! = 1, x-y slice of at k=kslice is plotted. ! = 2, x-z slice of at j=jslice is plotted. ! = 3, y-z slice of at i=islice is plotted. ! = 4, horizontal slice at z index islice is plotted. ! = 5, xy-z cross section of wind islice is plotted. ! = 6, data field on constant p-level is plotted. ! = 0, all of the three slices above are plotted. ! ! WORK ARRAYS: ! ! tem1 Temporary work array. ! tem2 Temporary work array. ! tem3 Temporary work array. ! tem4 Temporary work array. ! tem5 Temporary work array. ! ! (These arrays are defined and used locally (i.e. inside this ! subroutine), they may also be passed into routines called by ! this one. Exiting the call to this subroutine, these temporary ! work arrays may be used for other purposes therefore their ! contents overwritten. Please examine the usage of work arrays ! before you alter the code.) ! ! pp01 The pressure (mb) value at the specific p-level ! ercpl reciprocal of exponent ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! INTEGER :: nx,ny,nz INTEGER :: n ! REAL :: u(nx,ny,nz) REAL :: v(nx,ny,nz) REAL :: w(nx,ny,nz) REAL :: x(nx,ny,nz) REAL :: y(nx,ny,nz) REAL :: z(nx,ny,nz) ! REAL :: u1(nx,ny),v1(nx,ny) REAL :: u2(n,nz),v2(n,nz),w2(n,nz),zp(n,nz) REAL :: xp(n),yp(n) INTEGER :: kslice,jslice,islice CHARACTER (LEN=*) :: runname CHARACTER (LEN=*) :: label REAL :: xw,xe,dx,ys,yn,dy,zb,zt,dz INTEGER :: ibgn,iend,ist, jbgn,jend,jst, kbgn,kend,kst REAL :: time,factor INTEGER :: slicopt INTEGER :: iunits, itype COMMON /windvtr/iunits, itype CHARACTER (LEN=12) :: varname COMMON /varplt1/ varname REAL :: xw1,xe1,ys1,yn1 COMMON /xuvpar/xw1,xe1,ys1,yn1 ! !----------------------------------------------------------------------- ! ! Some constants ! !----------------------------------------------------------------------- ! REAL :: pp01, ercpl PARAMETER (ercpl=0.3678794) ! exp(-1.0) ! !----------------------------------------------------------------------- ! ! Work arrays: tem1,tem2,tem3,tem4,tem5 of size at least ! max( nx*ny, nx*nz, ny*nz). ! !----------------------------------------------------------------------- ! REAL :: tem1(*),tem2(*),tem3(*),tem4(*),tem5(*) REAL :: tem6(*) ! !----------------------------------------------------------------------- ! ! Common blocks for plotting control parameters ! !----------------------------------------------------------------------- ! REAL :: x01,y01 ! the first point of interpolation REAL :: x02,y02 ! the second point of interpolation REAL :: zlevel ! the given height of the slice REAL :: sinaf,cosaf,dist,sqrtdxy COMMON /slicev/x01,y01,x02,y02,sinaf,cosaf,dist,sqrtdxy COMMON /sliceh/zlevel INTEGER :: ovrobs,obsset,obscol,obs_marktyp REAL :: obs_marksz COMMON /obspar/ ovrobs,obsset,obscol,obs_marktyp, obs_marksz ! !----------------------------------------------------------------------- ! ! Misc. local Variables ! !----------------------------------------------------------------------- ! INTEGER :: i,j,k,ij,ik,jk,istep,jstep,length,isize,jsize,ksize REAL :: uunit CHARACTER (LEN=6) :: timhms CHARACTER (LEN=120) :: title ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! INTEGER :: trnplt ! flag to plot terain (1 or 0) REAL :: hterain(nx,ny) ! The height of the terrain. INTEGER :: ovrtrn ! overlay terrain option (0/1) REAL :: trninc,trnmin, trnmax ! terrain interval minimum, maximum REAL :: ztmin,ztmax COMMON /trnpar/ trnplt,ovrtrn,trninc,trnmin, trnmax,ztmin,ztmax ! INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat CHARACTER (LEN=6) :: stem2 CHARACTER (LEN=1) :: stem1 INTEGER :: smooth COMMON /smoothopt/smooth INTEGER :: id REAL :: x_tmp COMMON /tmphc2/ x_tmp INTEGER :: wrtflag CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab CHARACTER (LEN=25) :: timestring COMMON /timelev/wrtflag, timelab, levlab, timestring REAL :: tmpx, tmpy CHARACTER (LEN=20) :: distc REAL :: x101, y101, x102,y102 COMMON /slicev1/x101, y101, x102,y102 INTEGER :: llabel CHARACTER (LEN=120) :: label_copy ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! isize=(iend-ibgn)+1 jsize=(jend-jbgn)+1 ksize=(kend-kbgn)+1 label_copy = label llabel = 120 CALL xstrlnth(label_copy, llabel) CALL xpscmnt('Start plotting '//label_copy(1:llabel)) ! !----------------------------------------------------------------------- ! ! slicopt=1 Plot u-v field ! !----------------------------------------------------------------------- ! CALL cvttim ( time, timhms) IF( timhms(1:1) == '0' ) timhms(1:1)=' ' WRITE(timelab,'(''T='',F8.1,A)') time, & ' s ('//timhms(1:2)//':'//timhms(3:4)//':'//timhms(5:6)//')' CALL get_time_string ( time, timestring) IF ( slicopt == 2 .OR. slicopt == 3 .OR. slicopt == 5) THEN CALL cal_dist(haxisu,dx,dy,x01,y01,x02,y02,slicopt, & tmpx,tmpy,distc) END IF ! !----------------------------------------------------------------------- ! ! Set up terrain, if needed. ! !----------------------------------------------------------------------- ! IF(trnplt == 1 .OR.trnplt == 2 .OR. ovrtrn == 1) THEN DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem5(ij)=hterain(i,j) END DO END DO END IF IF( slicopt == 1 .OR. slicopt == 0 ) THEN k = kslice DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij) = -9999.0 tem2(ij) = -9999.0 IF( IF( tem3(ij)=x(i,j,k) tem4(ij)=y(i,j,k) END DO END DO IF (k /= 2) THEN WRITE(levlab,'(''GRID LEVEL='',I3)')k WRITE(title,'(''U-V '',A)')label ELSE WRITE(levlab,'(''FIRST LEVEL ABOVE GROUND (SURFACE)'')') WRITE(title,'(''U-V '',A)') label END IF length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) uunit = 10.0 CALL xvmode(1) istep = ist jstep = jst DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem6) CALL smooth9pmv(tem2,isize,jsize,1,isize,1,jsize,tem6) END DO CALL vtr2d(tem1,tem2,tem3,tem4, uunit, xw,xe,dx,ys,yn,dy, & isize,istep,jsize,jstep, & title(1:length),runname, 1, & tem5,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=2 Plot u-w field ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 2 .OR. slicopt == 0 ) THEN x_tmp = y(1,jslice,1) j = jslice dist = (j-1.5)*tmpy length=20 CALL strmin ( distc, length) WRITE(levlab,'(''X-Z PLANE AT Y='',F8.1,A)')dist,distc(1:length) IF(varname(1:6) == 'xuvplt') THEN xw1=xw xe1=xe ys1=ys yn1=yn id=4 DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 tem2(ik) = -9999.0 IF( IF( tem3(ik)=x(i,j,k) tem4(ik)=z(i,j,k) END DO END DO WRITE(title,'(''U-V '',A)')label ELSE id=2 DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 tem2(ik) = -9999.0 IF( IF( tem3(ik)=x(i,j,k) tem4(ik)=z(i,j,k) END DO END DO WRITE(title,'(''U-W '',A)')label END IF length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,ksize,1,isize,1,ksize,tem6) CALL smooth9pmv(tem2,isize,ksize,1,isize,1,ksize,tem6) END DO uunit = 10.0 CALL xvmode(1) istep = ist jstep = kst CALL vtr2d(tem1,tem2,tem3,tem4,uunit, xw,xe,dx,zb,zt,dz, & isize,istep,ksize,jstep, & title(1:length),runname, id, & tem5,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=3 Plot v-w field ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 3 .OR. slicopt == 0 ) THEN x_tmp = y(1,jslice,1) i = islice dist = (i-1.5)*tmpx length=20 CALL strmin ( distc, length) WRITE(levlab,'(''Y-Z PLANE AT X='',F8.1,A)')dist,distc(1:length) IF(varname(1:6) == 'xuvplt') THEN xw1=xw xe1=xe ys1=ys yn1=yn id=4 DO k=kbgn,kend DO j=jbgn,jend jk = j-jbgn+1 + (k-kbgn)*jsize tem1(jk) = -9999.0 tem2(jk) = -9999.0 IF( IF( tem3(jk)=y(i,j,k) tem4(jk)=z(i,j,k) END DO END DO WRITE(title,'(''U-V '',A)')label ELSE id=3 DO k=kbgn,kend DO j=jbgn,jend jk = j-jbgn+1 + (k-kbgn)*jsize tem1(jk) = -9999.0 tem2(jk) = -9999.0 IF( IF( tem3(jk)=y(i,j,k) tem4(jk)=z(i,j,k) END DO END DO WRITE(title,'(''V-W '',A)')label END IF length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,jsize,ksize,1,jsize,1,ksize,tem6) CALL smooth9pmv(tem2,jsize,ksize,1,jsize,1,ksize,tem6) END DO uunit = 10.0 CALL xvmode(1) istep = jst jstep = kst CALL vtr2d(tem1,tem2,tem3,tem4,uunit, ys,yn,dy,zb,zt,dz, & jsize,istep,ksize,jstep, & title(1:length),runname, id, & tem5,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=4 Plot u-v field on constant z levels ! slicopt=6 Plot u-v field on constant pressure levels ! slicopt=7 Plot u-v field on constant PT levels ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 4.OR.slicopt == 6.OR.slicopt == 7 ) THEN ! CALL hintrp(nx,ny,nz,u,z,zlevel,u1) ! CALL hintrp(nx,ny,nz,v,z,zlevel,v1) CALL hintrp1(nx,ny,nz,kbgn,kend,u,z,zlevel,u1) CALL hintrp1(nx,ny,nz,kbgn,kend,v,z,zlevel,v1) DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij) = -9999.0 tem2(ij) = -9999.0 IF(u1(i,j) /= -9999.0) tem1(ij)=u1(i,j)*factor IF(v1(i,j) /= -9999.0) tem2(ij)=v1(i,j)*factor tem3(ij)=x(i,j,2) tem4(ij)=y(i,j,2) END DO END DO IF( slicopt == 4) THEN WRITE(levlab,'(''Z='',F7.3,'' KM MSL'')') & zlevel ELSE IF( slicopt == 6) THEN pp01 = 0.01*ercpl**zlevel WRITE(levlab,'(''P='',F7.2,A)') pp01, ' MB' ELSE WRITE(levlab,'(''THETA='',F5.1,A)') zlevel, ' (K)' END IF WRITE(title,'(''U-V '',A)') label length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) uunit = 10.0 CALL xvmode(1) istep = ist jstep = jst DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem6) CALL smooth9pmv(tem2,isize,jsize,1,isize,1,jsize,tem6) END DO CALL vtr2d(tem1,tem2,tem3,tem4, uunit, xw,xe,dx,ys,yn,dy, & isize,istep,jsize,jstep, & title(1:length),runname, 1, & tem5,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=5 Plot u-v field ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 5 ) THEN CALL sectvrt(nx,ny,nz,u,x,y,z,u2,zp,n,xp,yp) CALL sectvrt(nx,ny,nz,v,x,y,z,v2,zp,n,xp,yp) CALL sectvrt(nx,ny,nz,w,x,y,z,w2,zp,n,xp,yp) IF(varname(1:6) == 'xuvplt') THEN xw1=xw xe1=xe ys1=ys yn1=yn id=4 DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 tem2(ik) = -9999.0 IF(u2(i,k) /= -9999.0) tem1(ik)= u2(i,k)*factor IF(v2(i,k) /= -9999.0) tem2(ik)= v2(i,k)*factor tem3(ik)=xw+(i-ibgn)* sqrtdxy tem4(ik)=zp(i,k) END DO END DO ELSE id=2 DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 tem2(ik) = -9999.0 IF(u2(i,k) /= -9999.0 .AND. v2(i,k) /= -9999.0) & tem1(ik)=(u2(i,k)*cosaf+v2(i,k)*sinaf)*factor IF(w2(i,k) /= -9999.0) tem2(ik)=w2(i,k)*factor tem3(ik)=xw+(i-ibgn)* sqrtdxy tem4(ik)=zp(i,k) END DO END DO END IF IF(xfmat == -1) THEN length=20 CALL strmin(distc,length) IF(varname(1:6) == 'xuvplt') THEN length=20 CALL strmin(distc,length) WRITE(title,'(''U-V '',A)') label WRITE(levlab,'(''XY-Z PLOT FROM '',4(A,F5.1),A,A)') & '(',x101,',',y101,') through (',x102,',',y102,') ', & distc(1:length) ELSE WRITE(title,'(''V-W '',A)') label WRITE(levlab, & '(''VERTICAL PLANE FROM '',4(A,F8.1),A,A)') & '(',x101,',',y101,') through (',x102,',',y102,') ', & distc(1:length) END IF ELSE IF(xfmat == 0) THEN length=20 CALL strmin(distc,length) IF(varname(1:6) == 'xuvplt') THEN WRITE(title,'(''U-V '',A)') label WRITE(levlab,'(''XY-Z PLOT FROM '',4(A,F5.1),A,A)') & '(',x101,',',y101,') through (',x102,',',y102,') ', & distc(1:length) ELSE WRITE(title,'(''V-W '',A)') label WRITE(levlab, & '(''VERTICAL PLANE FROM '',4(A,I5),A,A)') & '(',INT(x101),',',INT(y101),') through (',INT(x102),',', & INT(y102),') ', distc(1:length) END IF ELSE WRITE(stem1,'(i1)')xfmat WRITE(stem2,43)stem1 43 FORMAT('f8.',a1) END IF length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,ksize,1,isize,1,ksize,tem6) CALL smooth9pmv(tem2,isize,ksize,1,isize,1,ksize,tem6) END DO uunit = 10.0 CALL xvmode(1) istep = ist jstep = kst CALL vtr2d(tem1,tem2,tem3,tem4,uunit, xw,xe,sqrtdxy,zb,zt,dz, & isize,istep,ksize,jstep, & title(1:length),runname, id, & tem5,slicopt) END IF CALL xpscmnt('End plotting '//label_copy(1:llabel)) RETURN END SUBROUTINE vtr3d ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE VTR2D ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE vtr2d(u,v,x,y,uunit1, xl,xr,dx,yb,yt,dy, & 6,20 m,istep,n,jstep,char1,char2, vpltmod, & hterain,slicopt) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Plot 2-d wind (u1,u2) vector field defined on grid points (x,y) ! using ZXPLOT package.. ! !----------------------------------------------------------------------- ! ! AUTHOR: ! ! MODIFICATION HISTORY: ! ! 1/24/96 (J. Zong and M. Xue) ! Fixed a problem related to finding the minima and maxima of u & v ! when there exist missing data. The initial min. and max. should be ! set to values other than the missing value, -9999.0. ! !----------------------------------------------------------------------- ! ! ! INPUT: ! u m by n 2-dimensional array of u (left-to-right) ! wind components (m/s) ! v m by n 2-dimensional array of v (down-to-up) ! wind components (m/s) ! ! x x coordinate of grid points in physical/comp. space (m) ! y y coordinate of grid points in physical/comp. space (m) ! ! uunit1 ! ! xl,xr The left and right bound of the physical domain. ! dx Spacing between the x-axis tick marks ! yb,yt Bottom and top bound of the physical domain. ! dy Spacing between the y-axis tick marks ! ! m First dimension of vector component array ! istep Step increment for plotting in x direction ! ! n Second dimension of vector component array ! jstep Step increment for plotting in y direction ! ! char1 First character string to plot (title) ! char2 Second character string to plot (runname) ! ! vpltmod vpltmod = 1 for u-v vector (u=u, v=v in model space) ! vpltmod = 2 for u-w vector (u=u, v=w in model space) ! vpltmod = 3 for v-w vector (u=v, v=w in model space) ! hterain the height of terrain ! slicopt slice orientation indicator ! = 1, x-y slice of at k=kslice is plotted. ! = 2, x-z slice of at j=jslice is plotted. ! = 3, y-z slice of at i=islice is plotted. ! = 4, horizontal slice at z index islice is plotted. ! = 5, xy-z cross section of wind islice is plotted. ! = 6, data field on constant p-level is plotted. ! = 0, all of the three slices above are plotted. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INCLUDE 'arpsplt.inc' ! INTEGER :: m,n ! REAL :: u(m,n) REAL :: v(m,n) REAL :: x(m,n) REAL :: y(m,n) REAL :: uunit1 REAL :: xl,xr,dx,yb,yt,dy INTEGER :: istep,jstep ! CHARACTER (LEN=*) :: char2 CHARACTER (LEN=*) :: char1 ! INTEGER :: slicopt,vpltmod ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: layover COMMON /laypar/ layover REAL :: ctinc,ctmin,ctmax,vtunt !contour interval and vector unit COMMON /incunt/ ctinc,ctmin,ctmax,vtunt REAL :: xleng,vunit COMMON /vecscl/ xleng,vunit ! REAL :: yxratio !the scaling factor the y/x ratio. COMMON /yratio/ yxratio ! INTEGER :: iunits, itype COMMON /windvtr/iunits, itype ! INTEGER :: ovrstaopt INTEGER :: ovrstam,staset,ovrstan,ovrstav,stacol,markprio,wrtstax INTEGER :: nsta_typ,sta_typ(10),sta_marktyp(10),sta_markcol(10) REAL :: sta_marksz(10),wrtstad COMMON /sta_par/ ovrstaopt,ovrstam,staset,ovrstan,ovrstav,stacol, & markprio,nsta_typ,sta_typ,sta_marktyp, & sta_markcol,sta_marksz,wrtstax,wrtstad ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor INTEGER :: ovrobs,obsset,obscol,obs_marktyp REAL :: obs_marksz COMMON /obspar/ ovrobs,obsset,obscol,obs_marktyp, obs_marksz REAL :: lblmag, ctrlbsiz, axlbsiz COMMON /labmag/ lblmag, ctrlbsiz, axlbsiz ! INTEGER :: flag INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat ! REAL :: ubarb(200,200), vbarb(200,200) COMMON /windtmp/ubarb, vbarb ! ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! INTEGER :: i,j,key REAL :: pl,pr,pb,pt ! plot space left, right, bottom, top coordinate REAL :: px,py ! plot space left-right length and up-down height REAL :: xs,ys ! real space left-right length and up-down height REAL :: pxc,pyc ! plot space left-right center and ! up-down center REAL :: x0,y0 REAL :: umax,umin ! max and min of u component REAL :: vmax,vmin ! max and min of v component REAL :: uunit REAL :: am ! REAL :: zlevel COMMON/sliceh/zlevel INTEGER, ALLOCATABLE :: iwrk(:) REAL , ALLOCATABLE :: xwk(:),ywk(:) INTEGER :: istatus ! REAL :: hterain(m,n) ! The height of the terrain. INTEGER :: timeovr COMMON /timover/ timeovr ! INTEGER :: LEN, len1 REAL :: xleng0,istand INTEGER :: iunits0 CHARACTER (LEN=15) :: ichar2 CHARACTER (LEN=150) :: f_char1 CHARACTER (LEN=150) :: ch INTEGER :: ntitle,titcol, nxpic, nypic, wpltime REAL :: titsiz CHARACTER (LEN=132) :: ptitle(3), footer_l, footer_c, footer_r COMMON /titpar1/ptitle, footer_l, footer_c, footer_r COMMON /titpar2/ntitle,titcol,wpltime, nxpic, nypic COMMON /titpar3/titsiz INTEGER :: col_table,pcolbar COMMON /coltable/col_table,pcolbar ! CHARACTER (LEN=12) :: varname COMMON /varplt1/ varname REAL :: xw1,xe1,ys1,yn1 COMMON /xuvpar/xw1,xe1,ys1,yn1 INTEGER :: wrtflag CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab CHARACTER (LEN=25) :: timestring COMMON /timelev/wrtflag, timelab, levlab, timestring CHARACTER (LEN=80) :: prestr INTEGER :: preflag COMMON /preinfo/ prestr,preflag REAL :: x101, y101, x102,y102 COMMON /slicev1/x101, y101, x102,y102 REAL :: ytmp !!local temporary variable REAL :: f_cputime,f_walltime,cpu1, cpu2,second1, second2 INTEGER :: xnwpic_called COMMON /callnwpic/xnwpic_called ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ALLOCATE(iwrk(m*n),STAT=istatus) ALLOCATE(xwk(m*n),ywk(m*n),STAT=istatus) second1= f_walltime() cpu1 = f_cputime() WRITE(6,'(/1x,a,a)') ' Plotting ',char1 IF( layover == 0 .OR. xnwpic_called == 0) THEN CALL xnwpic xnwpic_called=1 timeovr=0 wrtflag = 0 preflag = 0 prestr = levlab len1=80 CALL strmin(prestr,len1) ELSE timeovr=1 wrtflag = wrtflag + 1 END IF ! !----------------------------------------------------------------------- ! ! Get plotting space variables ! !----------------------------------------------------------------------- ! CALL xqpspc( pl, pr, pb, pt) px = pr - pl py = pt - pb xs = xr-xl ys = yt-yb pxc = (pr+pl)/2 pyc = (pb+pt)/2 ! !----------------------------------------------------------------------- ! ! Let the longest lenth determine size scaling of plot ! !----------------------------------------------------------------------- ! IF( py/px >= ys*yxratio/xs ) THEN py = ys*yxratio/xs*px CALL xpspac(pl, pr, pyc-py/2, pyc+py/2 ) ELSE px = xs/(ys*yxratio)*py CALL xpspac(pxc-px/2, pxc+px/2, pb, pt) END IF ! !----------------------------------------------------------------------- ! ! Set the real distance to plot distance scaling ! !----------------------------------------------------------------------- ! CALL xmap( xl, xr, yb,yt) ! !----------------------------------------------------------------------- ! ! Plot maps, boxes, and polygons ! !----------------------------------------------------------------------- ! CALL xcolor(lbcolor) CALL pltextra(slicopt, 1 ) ! !----------------------------------------------------------------------- ! ! Find max and min of data array ! !----------------------------------------------------------------------- ! DO j=1,n DO i=1,m IF( umin = u(i,j) vmin = v(i,j) GO TO 110 END DO END DO 110 CONTINUE umax=umin vmax=vmin DO j=1,n DO i=1,m IF( IF( IF( IF( END DO END DO ! !----------------------------------------------------------------------- ! ! Fill various labels ! !----------------------------------------------------------------------- ! CALL xchsiz( 0.030*(yt-yb) * lblmag ) CALL xcolor(lbcolor) IF ( layover < 1 ) THEN len1=50 CALL strmin(timelab,len1) CALL xcharl(xl,yt+0.07*ys, timestring(1:25)) CALL xcharr(xr+0.05*(xr-xl),yt+0.07*ys, timelab(1:len1)) IF(levlab /= ' ') THEN len1=80 CALL strmin(levlab,len1) CALL xcharc(xl+xs*0.5,yt+0.015*ys, levlab(1:len1)) preflag = 1 END IF ! len1=80 ! CALL strmin(levlab,len1) ! CALL xcharc(xl+xs*0.5,yt+0.015*ys, levlab(1:len1)) END IF IF(preflag == 0 .AND. levlab /= ' ') THEN len1=80 CALL strmin(levlab,len1) CALL xcharc(xl+xs*0.5,yt+0.015*ys, levlab(1:len1)) preflag = 1 END IF IF( vpltmod == 1 .OR. vpltmod == 4 ) THEN WRITE(ch,'(''Umin='',F7.2,'' Umax='',F7.2, & & '' Vmin='',f7.2,'' Vmax='',f7.2)') umin,umax,vmin,vmax ELSE IF( vpltmod == 2 ) THEN WRITE(ch,'(''Umin='',F7.2,'' Umax='',F7.2, & & '' Wmin='',f7.2,'' Wmax='',f7.2)') umin,umax,vmin,vmax ELSE WRITE(ch,'(''Vmin='',F7.2,'' Vmax='',F7.2, & & '' Wmin='',f7.2,'' Wmax='',f7.2)') umin,umax,vmin,vmax END IF LEN=150 CALL strmin(char1,LEN) IF( char1(LEN-1:LEN-1) == ')' ) char1(LEN-1:LEN-1)=',' IF(itype == 1) THEN WRITE(f_char1, '(a,''VECTOR)'')')char1(1:LEN) ELSE IF(itype == 2) THEN WRITE(f_char1, '(a, ''BARB)'')')char1(1:LEN) END IF ! if first levlab is not equal second levlab then attatch levlab on f_ch ! !mx LEN=150 CALL strmin(f_char1,LEN) len1=80 CALL strmin(levlab,len1) ! IF (preflag.eq.1 .and. prestr(1:len1).ne.levlab(1:len1) ! : .and. prestr(1:1).ne.' ' ! : .and.layover.ne.0 .and. levlab(1:1).ne.' ') THEN ! write(f_char1,'(a,a)') f_char1(1:len),levlab(1:len1) ! ENDIF WRITE(6,'(1x,a51)') ch(1:51) CALL xcolor(icolor) IF(lbaxis == 1) THEN IF(wrtstax == 0) THEN ytmp = 0.08 ELSE ytmp =0.14 END IF ELSE ytmp = 0.12 END IF LEN=150 CALL strmin(f_char1,LEN) CALL xchsiz(0.025*ys * lblmag ) CALL xcharl(xl-0.20*(xr-xl), yb-(yt-yb)*(ytmp+wrtflag*0.030), & f_char1(1:LEN)) len1=150 CALL strmin(ch,len1) CALL xcharr(xr+0.20*(xr-xl), yb-(yt-yb)*(ytmp+wrtflag*0.030), & ch(1:len1)) ! !----------------------------------------------------------------------- ! ! Set vector unit and plot vectors. ! !----------------------------------------------------------------------- ! uunit=uunit1 IF( vtunt /= 0.0 ) THEN uunit=vtunt CALL xvmode(2) END IF ! Set parameter for barb xleng0 = (pr-pl)/(m-1) * istep * 0.65 IF(iunits == 1 .AND. itype == 2) THEN iunits0=1 istand = 5. WRITE(ichar2,'(a15)')'5 m/s' ELSE IF(iunits == 2 .AND. itype == 2) THEN iunits0=2 istand = 10. WRITE(ichar2,'(a15)')'10 knots' ELSE IF (iunits == 3 .AND. itype == 2) THEN iunits0=2 istand = 10. WRITE(ichar2,'(a15)')'10 MPH' END IF IF(layover >= 1) CALL xcolor(icolor) CALL xmap(xl,xr, yb,yt) CALL xvectu(u,v,m,m,istep,n,jstep,xleng,uunit) CALL xcolor(icolor) CALL xwindw(xl, xr, yb, yt) IF(itype == 1) THEN CALL xvectr(u,v,x,y,m,m,istep,n,jstep,xleng,uunit) ELSE IF(itype == 2) THEN CALL xbarbs(u,v,x,y,m,m,istep,n,jstep,iunits0,xleng*0.65,2) END IF CALL xwdwof ! !----------------------------------------------------------------------- ! ! Plot axes with tick marks ! !----------------------------------------------------------------------- ! CALL pltaxes(slicopt,dx,dy) vunit=uunit x0=xl-(xr-xl)*0.08 y0=yb-(yt-yb)*0.07 key=0 am=0.5 IF( ((m-1)/istep) > 30 ) am=1.0 IF(itype == 1) THEN IF(varname(1:6) == 'xuvplt') CALL xmap(xl, xr,yb,yt) CALL xvectk(x0,y0,xleng*am,uunit*am, key) CALL xmap(xl, xr, yb, yt) END IF ! !----------------------------------------------------------------------- ! ! Plot terrain profile in vertical slices ! !----------------------------------------------------------------------- ! IF(slicopt == 2 .OR. slicopt == 3 .OR.slicopt == 5) THEN CALL xcolor(trcolor) CALL xthick(2) CALL xpenup( x(1,1), y(1,1)-0.5*(y(1,2)-y(1,1)) ) DO i=2,m CALL xpendn( x(i,1), y(i,1)-0.5*(y(i,2)-y(i,1)) ) END DO CALL xthick(1) END IF ! !----------------------------------------------------------------------- ! ! Overlay terrain contour if required in x-y level ! or Plot terrain outline in this slice zlevel . ! !----------------------------------------------------------------------- ! IF(timeovr == 0) CALL plttrn(hterain,x,y,m,n,slicopt) CALL xcolor(lbcolor) CALL xwindw(xl, xr, yb, yt) ! !----------------------------------------------------------------------- ! ! Plot observations ! !----------------------------------------------------------------------- ! IF(ovrobs == 1 .AND. obsset == 1.AND. (slicopt == 1.OR.slicopt == 4.OR. & slicopt == 6.OR.slicopt == 7.OR.slicopt == 8)) THEN CALL pltobs(3) obsset=0 END IF ! !----------------------------------------------------------------------- ! ! Plot station labels ! !----------------------------------------------------------------------- ! CALL xcolor(lbcolor) IF(ovrstaopt == 1 .AND. staset == 1 .AND. & (ovrstam == 1.OR.ovrstan == 1.OR.ovrstav == 1).AND. & (slicopt == 1.OR.slicopt == 4.OR.slicopt == 6 & .OR.slicopt == 7.OR.slicopt == 8) & .AND.timeovr == 0 ) THEN CALL xchsiz(0.025*ys * lblmag) CALL pltsta(u,v,x,y,m,n,0,slicopt) ! staset=0 END IF IF( ovrstaopt == 1 .AND. wrtstax == 1 .AND. timeovr == 0 & .AND.(slicopt == 2.OR.slicopt == 3.OR. slicopt == 5) ) THEN CALL xchsiz(0.025*ys * lblmag) flag=1 CALL pltsta(u,v,x,y,m,n,flag,slicopt) END IF CALL xwdwof !----------------------------------------------------------------------- ! ! Plot additional text below the figure ! !----------------------------------------------------------------------- CALL label2d(char2) cpu2 = f_cputime() second2 = f_walltime() ! write(6,*)'!!!! total cpu time for one VTR2D :', ! : cpu2-cpu1,' PLOT:',varname DEALLOCATE(iwrk) DEALLOCATE(xwk,ywk) RETURN END SUBROUTINE vtr2d ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE VTRUNT ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE vtrunt( vtunt0 ) 6 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set the wind vector unit for wind field to be plotted by VTR2D. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! 6/08/92 Added full documentation (K. Brewster) ! !----------------------------------------------------------------------- ! ! INPUT: ! ! vtunt0 Unit vector ! If VTUNT0 = 0.0, the unit is internally determined. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE REAL :: vtunt0 ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! REAL :: ctinc,ctmin,ctmax,vtunt ! contour interval and vector unit COMMON /incunt/ ctinc,ctmin,ctmax,vtunt ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! vtunt = vtunt0 RETURN END SUBROUTINE vtrunt ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE STRM3D ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE strm3d(u,v,w, x,y,z, xw,xe,dx, ys,yn,dy, zb,zt,dz, & 2,43 nx,ibgn,iend,ist, ny,jbgn,jend,jst, nz,kbgn,kend,kst, & kslice, jslice, islice, time, runname,factor,slicopt, & n,xp,yp,zp,u1,v1,u2,v2,w2, & tem1,tem2,tem3,tem4,tem5, & tem6,hterain) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Plot a streamline field in 2-d slices ! ! AUTHOR: Ming Xue ! 1/16/1992 ! ! MODIFICATION HISTORY: ! ! 3/25/96 (K. Brewster) ! Added variables isize,jsize,ksize ! !----------------------------------------------------------------------- ! ! INPUT: ! ! u 3-dimensional array of u wind components (m/s) ! v 3-dimensional array of v wind components (m/s) ! w 3-dimensional array of w wind components (m/s) ! ! x x coordinate of grid points in physical/comp. space (m) ! y y coordinate of grid points in physical/comp. space (m) ! z z coordinate of grid points in physcal space (m) ! ! xw value of x for first i grid point to plot ! xe value of x for last i grid point to plot ! ys value of y for first j grid point to plot ! yn value of y for last j grid point to plot ! zb value of z for first k grid point to plot ! zt value of z for last k grid point to plot ! ! nx first dimension of b ! ibgn index of first i grid point to plot ! iend index of last i grid point to plot ! ! ny second dimension of b ! jbgn index of first j grid point to plot ! jend index of last j grid point to plot ! ! nz third dimension of b ! kbgn index of first k grid point to plot ! kend index of last k grid point to plot ! ! ist step size in x direction ! jst step size in y direction ! kst step size in z direction ! ! time time of data in seconds ! ! kslice k index of plane for slicopt=1 x-y slice ! jslice j index of plane for slicopt=2 x-z slice ! islice i index of plane for slicopt=1 y-z slice ! ! runname character string decribing run ! ! factor scaling factor for winds ! V*factor wind vectors are plotted ! ! slicopt slice orientation indicator ! = 1, x-y slice of at k=kslice is plotted. ! = 2, x-z slice of at j=jslice is plotted. ! = 3, y-z slice of at i=islice is plotted. ! = 4, horizontal slice at z index islice is plotted. ! = 5, xy-z cross section of wind islice is plotted. ! = 6, data field on constant p-level is plotted. ! = 0, all of the three slices above are plotted. ! ! WORK ARRAYS: ! ! tem1 Temporary work array. ! tem2 Temporary work array. ! tem3 Temporary work array. ! tem4 Temporary work array. ! tem5 Temporary work array. ! tem6 Temporary work array. ! ! (These arrays are defined and used locally (i.e. inside this ! subroutine), they may also be passed into routines called by ! this one. Exiting the call to this subroutine, these temporary ! work arrays may be used for other purposes therefore their ! contents overwritten. Please examine the usage of work arrays ! before you alter the code.) ! ! pp01 The pressure (mb) value at the specific p-level ! ercpl reciprocal of exponent ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! INTEGER :: nx,ny,nz INTEGER :: n ! REAL :: u(nx,ny,nz) REAL :: v(nx,ny,nz) REAL :: w(nx,ny,nz) REAL :: x(nx,ny,nz) REAL :: y(nx,ny,nz) REAL :: z(nx,ny,nz) ! REAL :: u1(nx,ny),v1(nx,ny) REAL :: u2(n,nz),v2(n,nz),w2(n,nz),zp(n,nz) REAL :: xp(n),yp(n) INTEGER :: kslice,jslice,islice CHARACTER (LEN=*) :: runname REAL :: xw,xe,dx,ys,yn,dy,zb,zt,dz INTEGER :: ibgn,iend,ist, jbgn,jend,jst, kbgn,kend,kst REAL :: time,factor INTEGER :: slicopt REAL :: x_tmp COMMON /tmphc2/ x_tmp ! !----------------------------------------------------------------------- ! ! Some constants ! !----------------------------------------------------------------------- ! REAL :: pp01, ercpl PARAMETER (ercpl=0.3678794) ! exp(-1.0) ! !----------------------------------------------------------------------- ! ! Work arrays: tem1,tem2,tem3 of size at least ! max( nx*ny, nx*nz, ny*nz). ! !----------------------------------------------------------------------- ! REAL :: tem1(*),tem2(*),tem3(*),tem4(*),tem5(*),tem6(*) INTEGER :: nzmax PARAMETER (nzmax = 300) REAL :: fdata(nzmax),zdata(nzmax),fprof(nzmax),zprof(nzmax) ! !----------------------------------------------------------------------- ! ! Common blocks for plotting control parameters ! !----------------------------------------------------------------------- ! REAL :: x01,y01 ! the first point of interpolation REAL :: x02,y02 ! the second point of interpolation REAL :: zlevel ! the given height of the slice REAL :: sinaf,cosaf,dist,sqrtdxy COMMON /slicev/x01,y01,x02,y02,sinaf,cosaf,dist,sqrtdxy COMMON /sliceh/zlevel ! !----------------------------------------------------------------------- ! ! Misc. local Variables ! !----------------------------------------------------------------------- ! INTEGER :: i,j,k,ij,ik,jk,length,isize,jsize,ksize CHARACTER (LEN=6) :: timhms CHARACTER (LEN=120) :: title ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! INTEGER :: trnplt ! flag to plot terain (1 or 0) REAL :: hterain(nx,ny) ! The height of the terrain. INTEGER :: ovrtrn ! overlay terrain option (0/1) REAL :: trninc,trnmin, trnmax ! terrain interval minimum, maximum REAL :: ztmin,ztmax COMMON /trnpar/ trnplt,ovrtrn,trninc,trnmin, trnmax,ztmin,ztmax ! INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat CHARACTER (LEN=4) :: stem2 CHARACTER (LEN=1) :: stem1 INTEGER :: smooth COMMON /smoothopt/smooth INTEGER :: wrtflag CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab CHARACTER (LEN=25) :: timestring COMMON /timelev/wrtflag, timelab, levlab, timestring REAL :: tmpx, tmpy CHARACTER (LEN=20) :: distc REAL :: x101, y101, x102,y102 COMMON /slicev1/x101, y101, x102,y102 ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! isize=(iend-ibgn)+1 jsize=(jend-jbgn)+1 ksize=(kend-kbgn)+1 ! !----------------------------------------------------------------------- ! ! setup time label ! !----------------------------------------------------------------------- ! CALL cvttim ( time, timhms) IF( timhms(1:1) == '0' ) timhms(1:1)=' ' WRITE(timelab,'(''T='',F8.1,A)') time, & ' s ('//timhms(1:2)//':'//timhms(3:4)//':'//timhms(5:6)//')' CALL get_time_string ( time, timestring) ! !----------------------------------------------------------------------- ! ! Set up terrain, if needed. ! !----------------------------------------------------------------------- ! IF(trnplt == 1 .OR.trnplt == 2 .OR. ovrtrn == 1) THEN DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem6(ij)=hterain(i,j) END DO END DO END IF IF ( slicopt == 2 .OR. slicopt == 3 .OR. slicopt == 5) THEN CALL cal_dist(haxisu,dx,dy,x01,y01,x02,y02,slicopt, & tmpx,tmpy,distc) END IF ! !----------------------------------------------------------------------- ! ! slicopt=1 Plot u-v field ! !----------------------------------------------------------------------- ! IF( slicopt == 1 .OR. slicopt == 0 ) THEN k = kslice DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij) = -9999.0 tem2(ij) = -9999.0 IF( IF( tem3(ij)=x(i,j,k) tem5(ij)=y(i,j,k) END DO END DO IF (k /= 2) THEN WRITE(title,'(''U-V STREAMLINE'')') WRITE(levlab,'(''X-Y CROSS SECTION THROUGH K='',I3)')k ELSE WRITE(title,'(''U-V STREAMLINE'')') WRITE(levlab,'(''X-Y CROSS SECTION THROUGH K=2 (SURFACE)'')') END IF length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem4) CALL smooth9pmv(tem2,isize,jsize,1,isize,1,jsize,tem4) END DO CALL strm2d(tem1,tem2, xw,xe,ys,yn, dx, dy, & isize,jsize, & title(1:length),runname, tem3,tem5, & tem6,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=2 Plot u-w streamline ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 2 .OR. slicopt == 0 ) THEN x_tmp = y(1,jslice,1) j = jslice DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 tem2(ik) = -9999.0 IF( IF( tem4(ik)=z(i,j,k) tem3(ik)=x(i,j,k) END DO END DO IF( nzmax < ksize) THEN WRITE(6,'(1x,a)') & 'nzmax given in STRM3D too small. Job stopped.' STOP END IF DO k=1,ksize zprof(k)= zb+(k-1)*(zt-zb)/(kend-kbgn) END DO CALL unigrid(isize,ksize,tem1,tem4, & fdata,zdata,fprof,zprof) CALL unigrid(isize,ksize,tem2,tem4, & fdata,zdata,fprof,zprof) WRITE(title,'(''U-W STREAMLINE '')') dist = (j-1)*tmpy length=20 CALL strmin ( distc, length) WRITE(levlab,'(''X-Z CROSS SECTION THROUGH J='',I3, & & '' (y = '',f8.1,a)')j,dist,distc(1:length) length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,ksize,1,isize,1,ksize,tem4) CALL smooth9pmv(tem2,isize,ksize,1,isize,1,ksize,tem4) END DO CALL strm2d(tem1,tem2, xw,xe,zb,zt, dx, dz, & isize,ksize, & title(1:length),runname, tem3 ,tem4, & tem6,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=3 Plot v-w field ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 3 .OR. slicopt == 0 ) THEN x_tmp = y(1,jslice,1) i = islice DO k=kbgn,kend DO j=jbgn,jend jk = j-jbgn+1 + (k-kbgn)*jsize tem1(jk) = -9999.0 tem2(jk) = -9999.0 IF( IF( tem4(jk)=z(i,j,k) tem5(jk)=y(i,j,k) END DO END DO IF( nzmax < ksize) THEN WRITE(6,'(1x,a)') & 'nzmax given in STRM3D too small. Job stopped.' STOP END IF DO k=1,ksize zprof(k)= zb+(k-1)*(zt-zb)/(kend-kbgn) END DO CALL unigrid(jsize,ksize,tem1,tem4, & fdata,zdata,fprof,zprof) CALL unigrid(jsize,ksize,tem2,tem4, & fdata,zdata,fprof,zprof) dist = (i-1)*tmpx length=20 CALL strmin ( distc, length) WRITE(levlab,'(''Y-Z CROSS SECTION THROUGH I='',I3, & & '' ( x='',f8.1,a )')i,dist, distc(1:length) WRITE(title,'(''V-W STREAMLINE'')') length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,jsize,ksize,1,jsize,1,ksize,tem4) CALL smooth9pmv(tem2,jsize,ksize,1,jsize,1,ksize,tem4) END DO CALL strm2d(tem1,tem2, ys,yn,zb,zt, dy, dz, & jsize,ksize, & title(1:length),runname, tem5 ,tem4, & tem6,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=4 Plot u-v streamlines on constant z levels ! slicopt=6 Plot u-v streamlines on constant pressure levels ! slicopt=7 Plot u-v streamlines on constant PT levels ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 4.OR.slicopt == 6.OR.slicopt == 7 ) THEN ! CALL hintrp(nx,ny,nz,u,z,zlevel,u1) ! CALL hintrp(nx,ny,nz,v,z,zlevel,v1) CALL hintrp1(nx,ny,nz,kbgn,kend,u,z,zlevel,u1) CALL hintrp1(nx,ny,nz,kbgn,kend,v,z,zlevel,v1) DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij) = -9999.0 tem2(ij) = -9999.0 IF(u1(i,j) /= -9999.0) tem1(ij)=u1(i,j)*factor IF(v1(i,j) /= -9999.0) tem2(ij)=v1(i,j)*factor tem3(ij)=x(i,j,2) tem5(ij)=y(i,j,2) END DO END DO IF( slicopt == 4) THEN WRITE(levlab,'(''Z='',F7.3,A,'' MSL'')')zlevel,' KM' ELSE IF( slicopt == 6) THEN pp01=0.01*ercpl**zlevel WRITE(levlab,'(''P='',F7.2,A)') pp01, ' MB' ELSE WRITE(levlab,'(''THETA='',F5.1,A)') zlevel, ' (K)' END IF WRITE(title,'(''U-V STREAMLINE'')') length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem4) CALL smooth9pmv(tem2,isize,jsize,1,isize,1,jsize,tem4) END DO CALL strm2d(tem1,tem2, xw,xe,ys,yn, dx, dy, & isize,jsize, & title(1:length),runname, tem3 ,tem5, & tem6,slicopt) ! !----------------------------------------------------------------------- ! ! slicopt=5 Plot V-w field ! !----------------------------------------------------------------------- ! ELSE IF( slicopt == 5 ) THEN CALL sectvrt(nx,ny,nz,u,x,y,z,u2,zp,n,xp,yp) CALL sectvrt(nx,ny,nz,v,x,y,z,v2,zp,n,xp,yp) CALL sectvrt(nx,ny,nz,w,x,y,z,w2,zp,n,xp,yp) DO k=kbgn,kend DO i=ibgn,iend ik = i-ibgn+1 + (k-kbgn)*isize tem1(ik) = -9999.0 tem2(ik) = -9999.0 IF(u2(i,k) /= -9999.0 .AND. v2(i,k) /= -9999.0) & tem1(ik)=(u2(i,k)*cosaf+v2(i,k)*sinaf)*factor IF(w2(i,k) /= -9999.0) tem2(ik)=w2(i,k)*factor tem3(ik)=xw+(i-ibgn)*sqrtdxy tem4(ik)=zp(i,k) END DO END DO IF( nzmax < ksize) THEN WRITE(6,'(1x,a)') & 'nzmax given in STRM3D too small. Job stopped.' STOP END IF DO k=1,ksize zprof(k)= zb+(k-1)*(zt-zb)/(kend-kbgn) END DO CALL unigrid(isize,ksize,tem1,tem4, & fdata,zdata,fprof,zprof) CALL unigrid(isize,ksize,tem2,tem4, & fdata,zdata,fprof,zprof) IF(xfmat == -1) THEN length=20 CALL strmin(distc,length) WRITE(title,'(''V-W STREAMLINE'')') WRITE(levlab, & '(''VERT CROSS SECTION THROUGH '',4(A,F8.1),A,A)') & '(',x101,',',y101,') (',x102,',',y102,')',distc(1:length) ELSE IF(xfmat == 0 ) THEN length=20 CALL strmin(distc,length) WRITE(title,'(''V-W STREAMLINE'')') WRITE(levlab, & '(''VERT CROSS SECTION THROUGH '',4(A,I5),A,A)') & '(',INT(x101),',',INT(y101),') (',INT(x102),',',INT(y102), & ')',distc(1:length) ELSE WRITE(stem1,'(i1)')xfmat WRITE(stem2,43)stem1 43 FORMAT('f8.',a1) ! write(title,'(''V-w streamline at t='',f8.1,a, ! : '' through '',4(a,stem2),a)') ! : time,' s ('//timhms(1:2)//':'//timhms(3:4)//':'// ! : timhms(5:6)//')','(',x01,',',y01,') (',x02,',',y02,')' END IF length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,ksize,1,isize,1,ksize,tem4) CALL smooth9pmv(tem2,isize,ksize,1,isize,1,ksize,tem4) END DO CALL strm2d(tem1,tem2, xw,xe,zb,zt, sqrtdxy, dz, & isize,ksize, & title(1:length),runname, tem3,tem4, & tem6,slicopt) END IF RETURN END SUBROUTINE strm3d ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE STRM2D ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE strm2d(u,v,xl,xr,yb,yt,dx,dy,m,n,char1,char2, x,y, & 5,6 hterain,slicopt) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Plot streamlines of a 2-d wind (u1,u2) field using ncargraphic ! subroutine strmln ! ! INPUT: ! u m by n 2-dimensional array of u (left-to-right) ! wind components (m/s) ! v m by n 2-dimensional array of v (down-to-up) ! wind components (m/s) ! ! xl,xr The left and right bound of the physical domain. ! yb,yt Bottom and top bound of the physical domain. ! dx,dy Grid interval in x and y direction (km) ! ! m First dimension of vector component array ! n Second dimension of vector component array ! ! char1 First character string to plot (title) ! char2 Second character string to plot (runname) ! ! x x coordinate of grid points in plot space (over on page) ! y y coordinate of grid points in plot space (up on page) ! ! hterain the height of terrain ! slicopt slice orientation indicator ! slicopt = 1, x-y slice of u,v at z index kslice is plotted. ! slicopt = 2, x-z slice of u,w at y index jslice is plotted. ! slicopt = 3, y-z slice of v,w at x index islice is plotted. ! slicopt = 4, x-y slice of u,v at z index islice is plotted. ! slicopt = 5, xy-z cross section of wind islice is plotted. ! slicopt = 6, data field on constant p-level is plotted. ! slicopt = 0, all of the three slices above are plotted. ! ! WORK ARRAY ! tem1 A work array of size at least m*n*2 ! tem2 A work array of size at least m*n*2 ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INCLUDE 'arpsplt.inc' ! INTEGER :: m,n INTEGER :: i ! REAL :: u(m,n) REAL :: v(m,n) REAL :: x(m,n) REAL :: y(m,n) REAL :: xl,xr,yb,yt,dx,dy ! CHARACTER (LEN=*) :: char2 CHARACTER (LEN=*) :: char1 INTEGER :: ierror ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: layover COMMON /laypar/ layover REAL :: ctinc,ctmin,ctmax,vtunt !contour interval and vector unit COMMON /incunt/ ctinc,ctmin,ctmax,vtunt ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! INTEGER :: flag INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat ! INTEGER :: ovrstaopt INTEGER :: ovrstam,staset,ovrstan,ovrstav,stacol,markprio,wrtstax INTEGER :: nsta_typ,sta_typ(10),sta_marktyp(10),sta_markcol(10) REAL :: sta_marksz(10),wrtstad COMMON /sta_par/ ovrstaopt,ovrstam,staset,ovrstan,ovrstav,stacol, & markprio,nsta_typ,sta_typ,sta_marktyp, & sta_markcol,sta_marksz,wrtstax,wrtstad ! REAL :: yxratio !the scaling factor the y/x ratio. COMMON /yratio/ yxratio INTEGER :: col_table,pcolbar COMMON /coltable/col_table,pcolbar ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! INTEGER :: nopic,nhpic,nvpic,ifont REAL :: pl,pr,pb,pt ! plot space left, right, bottom, top coordinate REAL :: px,py ! plot space left-right length and up-down height REAL :: xs,ys ! real space left-right length and up-down height REAL :: pxc,pyc ! plot space left-right center and ! up-down center REAL :: xlimit,ylimit REAL :: rotang REAL :: xp1,xp2,yp1,yp2 REAL :: xd1,xd2,yd1,yd2,xpos1,xpos2,ypos1,ypos2 ! REAL :: zlevel COMMON/sliceh/zlevel ! REAL :: hterain(m,n) ! The height of the terrain. ! INTEGER :: slicopt INTEGER :: timeovr COMMON /timover/ timeovr REAL :: lblmag, ctrlbsiz, axlbsiz COMMON /labmag/ lblmag, ctrlbsiz, axlbsiz ! INTEGER :: len1 INTEGER :: wrtflag CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab CHARACTER (LEN=25) :: timestring COMMON /timelev/wrtflag, timelab, levlab, timestring REAL :: x101, y101, x102,y102 COMMON /slicev1/x101, y101, x102,y102 INTEGER :: xnwpic_called COMMON /callnwpic/xnwpic_called REAL :: ytmp !! local temporary variable ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! WRITE(6,'(/1x,a,a)') ' Plotting ',char1 IF( layover == 0 .OR. xnwpic_called == 0) THEN CALL xnwpic xnwpic_called=1 timeovr=0 wrtflag = 0 ELSE timeovr=1 wrtflag = wrtflag + 1 END IF ! !----------------------------------------------------------------------- ! ! Get plotting space variables ! !----------------------------------------------------------------------- ! CALL xqpspc( pl, pr, pb, pt) px = pr - pl py = pt - pb xs = xr-xl ys = yt-yb pxc = (pr+pl)/2 pyc = (pb+pt)/2 ! !----------------------------------------------------------------------- ! ! Let the longest lenth determine size scaling of plot ! !----------------------------------------------------------------------- ! IF( py/px >= ys*yxratio/xs ) THEN py = ys*yxratio/xs*px CALL xpspac(pl, pr, pyc-py/2, pyc+py/2 ) ELSE px = xs/(ys*yxratio)*py CALL xpspac(pxc-px/2, pxc+px/2, pb, pt) END IF ! !----------------------------------------------------------------------- ! ! Set the real distance to plot distance scaling ! !----------------------------------------------------------------------- ! CALL xmap( xl, xr, yb,yt) ! !----------------------------------------------------------------------- ! ! Plot map, boxes and polygons. ! !----------------------------------------------------------------------- ! CALL xcolor(lbcolor) CALL pltextra(slicopt, 1 ) xpos1 = xl xpos2 = xr ypos1 = yb ypos2 = yt CALL xtrans(xpos1,ypos1) CALL xtrans(xpos2,ypos2) CALL xzx2ncar(xpos1,ypos1) CALL xzx2ncar(xpos2,ypos2) ! IF(slicopt == 2 .OR. slicopt == 3 .OR.slicopt == 5) THEN CALL xcolor(trcolor) CALL xthick(3) CALL xpenup( x(1,1), y(1,1)-0.5*(y(1,2)-y(1,1)) ) DO i=2,m CALL xpendn( x(i,1), y(i,1)-0.5*(y(i,2)-y(i,1)) ) END DO CALL xthick(1) END IF ! !----------------------------------------------------------------------- ! ! Overlay terrain contour if required in x-y level ! or Plot terrain outline in this slice zlevel . ! !----------------------------------------------------------------------- ! IF( timeovr == 0)CALL plttrn(hterain,x,y,m,n,slicopt) CALL xcolor(lbcolor) ! !----------------------------------------------------------------------- ! ! Plot station labels ! !----------------------------------------------------------------------- ! IF(ovrstaopt == 1 .AND. staset == 1 .AND. & (ovrstam == 1.OR.ovrstan == 1.OR.ovrstav == 1).AND. & (slicopt == 1.OR.slicopt == 4.OR.slicopt == 6 & .OR.slicopt == 7.OR.slicopt == 8) & .AND.timeovr == 0 ) THEN CALL xchsiz(0.025*ys * lblmag) CALL pltsta(u,v,x,y,m,n,0,slicopt) ! staset=0 END IF ! !----------------------------------------------------------------------- ! ! Plot observations ! !----------------------------------------------------------------------- ! IF( ovrstaopt == 1 .AND. wrtstax == 1 .AND. timeovr == 0 & .AND.(slicopt == 2.OR.slicopt == 3.OR. slicopt == 5) ) THEN CALL xchsiz(0.025*ys * lblmag) flag=1 CALL pltsta(u,v,x,y,m,n,flag,slicopt) END IF ! !----------------------------------------------------------------------- ! ! Plot streamlines ! !----------------------------------------------------------------------- ! CALL xcolor(lbcolor) ! CALL xqset(xp1,xp2,yp1,yp2, xd1,xd2,yd1,yd2) CALL set(xpos1,xpos2,ypos1,ypos2, 1.0, FLOAT(m), 1.0, FLOAT(n),1) CALL xcolor(icolor) CALL strmln(u,v,y ,m,m,n, 1, ierror) CALL set(xp1,xp2,yp1,yp2, xd1,xd2,yd1,yd2, 1) ! !----------------------------------------------------------------------- ! ! Plot axes with tick marks ! !----------------------------------------------------------------------- ! CALL pltaxes(slicopt,dx,dy) !----------------------------------------------------------------------- ! ! Plot labels ! !----------------------------------------------------------------------- CALL xcolor(lbcolor) CALL xqnpic(nopic) CALL xqspac(nhpic, nvpic, rotang, xlimit, ylimit) ! write time and level CALL xchsiz( 0.030*ys * lblmag ) IF ( layover < 1 ) THEN len1=50 CALL strmin(timelab,len1) CALL xcharl(xl,yt+0.07*ys, timestring(1:25)) CALL xcharr(xr+0.05*(xr-xl),yt+0.07*ys, timelab(1:len1)) len1=80 CALL strmin(levlab,len1) CALL xcharc(xl+xs*0.5,yt+0.015*ys, levlab(1:len1)) END IF ! write variable label CALL xcolor(icolor) IF(lbaxis == 1) THEN IF(wrtstax == 0) THEN ytmp = 0.08 ELSE ytmp = 0.14 END IF ELSE ytmp =0.12 END IF CALL xchsiz( 0.025*ys * lblmag ) CALL xcharl(xl-0.20*(xr-xl), yb-(yt-yb)*(ytmp+layover*0.030), & char1) CALL xcolor(lbcolor) IF (timeovr == 0) THEN IF(nopic == nhpic*(nvpic-1)+1 ) THEN ytmp =0.25 IF(layover < 1) CALL xcharl(xl,yb-(ytmp+layover*0.03)*(yt-yb), char2 ) CALL xqcfnt(ifont) CALL xcfont(xfont) ytmp = 0.20 IF(layover < 1) CALL xcharl(xl,yb-(0.20+layover*0.03)*(yt-yb), & 'CAPS/ARPS ' ) ! : 'Project Hub-CAPS Experimental ' ) !Hub-CAPS+ CALL xcfont(ifont) END IF timeovr=1 END IF RETURN END SUBROUTINE strm2d ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE CTRSFC ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE ctrsfc(a,x,y,x1,x2,dx,y1,y2,dy, & 45,6 nx,ibgn,iend, ny,jbgn,jend, & label,time, runname, factor,tem1,tem2,tem3, & tem4,tem5,hterain,slicopt,pltopt) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! To plot a contour map for a 2-d surface array. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! 4/20/1994 ! ! MODIFICATION HISTORY: ! ! 9/27/95 (Yuhe Liu) ! Fixed a bug in call of smth. Added the temporary array tem5 to ! the argument list. ! ! 3/25/96 (Keith Brewster) ! Added variables isize,jsize and replaced smth with smooth9p ! !----------------------------------------------------------------------- ! ! INPUT : ! ! a 2-d surface array. ! ! x x coordinate of grid points in physical/comp. space (m) ! y y coordinate of grid points in physical/comp. space (m) ! ! x1 value of x for first i grid point to plot ! x2 value of x for last i grid point to plot ! dx ! y1 value of y for first j grid point to plot ! y2 value of y for last j grid point to plot ! dy ! ! nx first dimension of a ! ibgn index of first i grid point to plot ! iend index of last i grid point to plot ! ! ny second dimension of a ! jbgn index of first j grid point to plot ! jend index of last j grid point to plot ! ! label character string describing the contents of a ! ! time time of data in seconds ! ! runname character string decribing run ! ! factor scaling factor for data ! contours are labelled a*factor ! slicopt slice orientation indicator ! slicopt = 1, x-y slice of u,v at z index kslice is plotted. ! slicopt = 2, x-z slice of u,w at y index jslice is plotted. ! slicopt = 3, y-z slice of v,w at x index islice is plotted. ! slicopt = 4, x-y slice of u,v at z index islice is plotted. ! slicopt = 5, xy-z cross section of wind islice is plotted. ! slicopt = 6, data field on constant p-level is plotted. ! slicopt = 0, all of the three slices above are plotted. ! plot variable plot option (0/1/2/3) ! ! WORK ARRAYS: ! ! tem1 Temporary work array. ! tem2 Temporary work array. ! tem3 Temporary work array. ! tem4 Temporary work array. ! tem5 Temporary work array. ! ! ! hterain The height of the terrain. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! INTEGER :: nx,ny REAL :: a(nx,ny) REAL :: x(nx,ny) REAL :: y(nx,ny) REAL :: x1,x2,dx,y1,y2,dy INTEGER :: ibgn,iend,jbgn,jend,length CHARACTER (LEN=6) :: timhms CHARACTER (LEN=*) :: label CHARACTER (LEN=*) :: runname REAL :: time REAL :: factor REAL :: tem1(*) REAL :: tem2(*) REAL :: tem3(*) REAL :: tem4(*) REAL :: tem5(*) REAL :: hterain(nx,ny) INTEGER :: slicopt INTEGER :: pltopt ! variable plot option (0/1/2/3) INTEGER :: ovrtrn,trnplt ! overlay terrain option (0/1) REAL :: trninc,trnmin, trnmax ! terrain interval minimum, maximum REAL :: ztmin,ztmax COMMON /trnpar/ trnplt,ovrtrn,trninc,trnmin, trnmax,ztmin,ztmax INTEGER :: smooth COMMON /smoothopt/smooth INTEGER :: wrtflag CHARACTER (LEN=120) :: label_copy CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab CHARACTER (LEN=25) :: timestring COMMON /timelev/wrtflag, timelab, levlab, timestring ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! INTEGER :: i,j,ij,isize,jsize,llabel CHARACTER (LEN=120) :: title ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! isize=(iend-ibgn)+1 jsize=(jend-jbgn)+1 label_copy = label llabel = 120 CALL xstrlnth(label_copy, llabel) CALL xpscmnt('Start plotting '//label_copy(1:llabel)) ! !----------------------------------------------------------------------- ! ! Set up terrain, if needed. ! !----------------------------------------------------------------------- ! IF(ovrtrn == 1) THEN DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem4(ij)=hterain(i,j) END DO END DO END IF ! CALL cvttim( time, timhms ) IF( timhms(1:1) == '0' ) timhms(1:1)=' ' WRITE(timelab,'(''T='',F8.1,A)') time, & ' s ('//timhms(1:2)//':'//timhms(3:4)//':'//timhms(5:6)//')' CALL get_time_string ( time, timestring) DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij) = -9999.0 IF( tem2(ij)=x(i,j) tem3(ij)=y(i,j) END DO END DO levlab=' ' WRITE(title,'(a)') label length = 120 CALL strlnth( title, length) CALL strmin ( title, length) DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem5) END DO CALL ctr2d(tem1,tem2,tem3, x1,x2,dx, y1,y2,dy, & isize,jsize,title(1:length),runname, & tem4,slicopt,pltopt) CALL xpscmnt('End plotting '//label_copy(1:llabel)) RETURN END SUBROUTINE ctrsfc ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE OVERLAY ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE overlay (layovr) 14 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set the layover counter parameter in the laypar common block ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! 6/08/92 Added full documentation (K. Brewster) ! ! 8/08/93 (MX) ! Automatically set the overlay parameter when input is not zero. ! !----------------------------------------------------------------------- ! ! INPUT: ! ! layovr The 'overlay' parameter. ! If layover .ne. 0, the following 2-d contour plot will be ! superimposed on the previous plot. ! layover =1, 2, ... indicating this is the ! layover'th (1st or 2nd ...) plot to be overlayed ! on the previous one. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! INTEGER :: layovr ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: layover, first_frame COMMON /laypar/ layover DATA first_frame /1/ ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! IF( first_frame == 1 .OR. layovr == 0 ) THEN layover = 0 ELSE layover = layover + 1 END IF first_frame = 0 RETURN END SUBROUTINE overlay ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE STYXRT ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE styxrt( yxrt ) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set the scaling factor of the y/x ratio of the plot. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! 6/08/92 Added full documentation (K. Brewster) ! !----------------------------------------------------------------------- ! ! INPUT: ! ! yxrt Ratio of height to length of plot space ! Default is set in the main program to 1.0 ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE REAL :: yxrt ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! REAL :: yxratio COMMON /yratio/ yxratio ! the scaling factor the y/x ratio. ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! yxratio = yxrt RETURN END SUBROUTINE styxrt ! !################################################################## !################################################################## !###### ###### !###### FUNCTION XFINC ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! REAL FUNCTION xfinc(x) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Automatically divide domain (0,x) to a number of subdomain ! with interval xfinc which is >=4 and =<16 for fold=1.0 ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! sometime ! ! MODIFICATIONS: ! 6/09/92 Added full documentation (K. Brewster) ! !----------------------------------------------------------------------- ! ! INPUT: ! x not sure ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! REAL :: x ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! INTEGER :: ipower REAL :: d,fold ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ipower= INT( ALOG10(x) ) d= INT(x/(10.0**ipower)) fold=1.0 xfinc=0.1*x IF( d >= 0.0 .AND. d < 3.0 ) THEN xfinc=2.0*10.0**(ipower-1) ELSE IF( d >= 3.0 .AND. d < 7.0 ) THEN xfinc=5.0*10.0**(ipower-1)*fold ELSE IF( d >= 7.0 .AND. d < 10. ) THEN xfinc=1.0*10.0** ipower*fold END IF IF(xfinc == 0.0) xfinc=x*0.1 RETURN END FUNCTION xfinc ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE CLIPWD ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE clipwd(x1,y1,x2,y2,idispl) 1,2 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Return the portion of a line that is within a given window ! (xw1,xw2,yw1,yw2) ! ! If the given line is completely outside the window, ! idispl=0, otherwise, idispl=1. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! 3/6/93 ! !----------------------------------------------------------------------- ! ! INPUT: ! ! x1 value of x for first i grid point to plot ! x2 value of x for last i grid point to plot ! y1 value of y for first j grid point to plot ! y2 value of y for last j grid point to plot ! ! idispl line orientation indicator ! idispl = 0, the given line is completely outside the window ! idispl = 1, the given line is partly inside the window ! !----------------------------------------------------------------------- ! ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! REAL :: x1,x2,y1,y2 INTEGER :: idispl ! !----------------------------------------------------------------------- ! ! Common blocks for plotting control parameters ! !----------------------------------------------------------------------- ! REAL :: xw1,xw2,yw1,yw2 COMMON /pltwdw/ xw1,xw2,yw1,yw2 INTEGER :: ic1(4),ic2(4) ! !----------------------------------------------------------------------- ! ! Misc. local Variables ! !----------------------------------------------------------------------- ! INTEGER :: i,knt,isw REAL :: x0,y0 REAL :: isum1,isum2,ic01,ic02,ic03,ic04 ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! knt = 0 5 knt = knt+1 CALL encodwd(x1,y1,ic1) CALL encodwd(x2,y2,ic2) isum1=ic1(1)+ic1(2)+ic1(3)+ic1(4) isum2=ic2(1)+ic2(2)+ic2(3)+ic2(4) idispl=1 IF(isum1+isum2 == 0) GO TO 999 idispl=0 DO i=1,4 20 IF(ic1(i)+ic2(i) == 2) GO TO 999 END DO ! !----------------------------------------------------------------------- ! ! Make sure (x1,y1) is outside the window ! !----------------------------------------------------------------------- ! isw=0 15 IF(isum1 == 0) THEN ic01=ic1(1) ic02=ic1(2) ic03=ic1(3) ic04=ic1(4) DO i=1,4 ic1(i)=ic2(i) END DO ic2(1)=ic01 ic2(2)=ic02 ic2(3)=ic03 ic2(4)=ic04 x0=x1 y0=y1 x1=x2 y1=y2 x2=x0 y2=y0 isw=1 END IF IF(ic1(1) == 1) THEN y1=y1+(xw1-x1)*(y2-y1)/(x2-x1) x1=xw1 ELSE IF(ic1(2) == 1) THEN y1=y1+(xw2-x1)*(y2-y1)/(x2-x1) x1=xw2 ELSE IF(ic1(3) == 1) THEN x1=x1+(yw1-y1)*(x2-x1)/(y2-y1) y1=yw1 ELSE IF(ic1(4) == 1) THEN x1=x1+(yw2-y1)*(x2-x1)/(y2-y1) y1=yw2 END IF IF(isw == 1) THEN x0=x1 y0=y1 x1=x2 y1=y2 x2=x0 y2=y0 END IF idispl=1 IF(knt > 10) THEN WRITE(6,*)'Dead loop encountered in CLIPWD, job stopped.' STOP 991 END IF GO TO 5 999 RETURN END SUBROUTINE clipwd ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE ENCODWD ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE encodwd(x,y,ic) 2 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Encode a line section for window clipping purpose. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! 3/6/93 ! !----------------------------------------------------------------------- ! ! INPUT: ! ! x value of x ! y value of y ! ic ! !----------------------------------------------------------------------- ! ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! REAL :: x,y INTEGER :: ic(4) ! !----------------------------------------------------------------------- ! ! Common blocks for plotting control parameters ! !----------------------------------------------------------------------- ! REAL :: xw1,xw2,yw1,yw2 COMMON /pltwdw/ xw1,xw2,yw1,yw2 ! !----------------------------------------------------------------------- ! ! Misc. local Variables ! !----------------------------------------------------------------------- ! INTEGER :: i ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! DO i=1,4 ic(i)=0 END DO IF(x < xw1) ic(1)=1 IF(x > xw2) ic(2)=1 IF(y < yw1) ic(3)=1 IF(y > yw2) ic(4)=1 RETURN END SUBROUTINE encodwd ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE CTRCOL ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE ctrcol (icol,icol0) 9 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set the color for field to plotted by CTR2D. ! ! !----------------------------------------------------------------------- ! ! AUTHOR: Min Zou ! ! !----------------------------------------------------------------------- ! ! INPUT: ! ! icol begin color index ! icol0 end color index ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: icol,icol0 ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: icolor,icolor1,lbcolor,trcolor COMMON /recolor/icolor,icolor1,lbcolor,trcolor !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! icolor=icol icolor1=icol0 RETURN END SUBROUTINE ctrcol ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE CTRVTR ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE ctrvtr (units0,type0) 6 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set the units and type for plot wind ! ! !----------------------------------------------------------------------- ! ! AUTHOR: Min Zou ! ! !----------------------------------------------------------------------- ! ! INPUT: ! ! units0 the units of wind ! type0 the type of wind ! wcolor0 the color index ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: units0,type0 ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: iunits, itype COMMON /windvtr/iunits, itype !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! iunits = units0 itype = type0 RETURN END SUBROUTINE ctrvtr !################################################################## !################################################################## !###### ###### !###### SUBROUTINE VARPLT ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE varplt( var ) 13 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Set the variable plot name for xconta. ! ! !----------------------------------------------------------------------- ! ! AUTHOR: Min Zou ! ! MODIFICATION HISTORY: ! 3/28/96 ! !----------------------------------------------------------------------- ! ! INPUT: ! ! var variable plot name ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! ! !----------------------------------------------------------------------- ! IMPLICIT NONE CHARACTER (LEN=*) :: var CHARACTER (LEN=12) :: varname COMMON /varplt1/ varname varname=var RETURN END SUBROUTINE varplt ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE VTRSFC ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE vtrsfc(u,v, x,y, xw,xe,dx, ys,yn,dy, & 3,7 nx,ibgn,iend,ist, ny,jbgn,jend,jst, & label,time, runname, factor, slicopt, & tem1,tem2,tem3,tem4, & tem5,tem6,hterain) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Plot vector fields in 2-d array ! ! AUTHOR: Min Zou ! 4/28/97 ! !----------------------------------------------------------------------- ! ! INPUT: ! ! u 2-dimensional array of u wind components (m/s) ! v 2-dimensional array of v wind components (m/s) ! ! x x coordinate of grid points in physical/comp. space (m) ! y y coordinate of grid points in physical/comp. space (m) ! z z coordinate of grid points in physical space (m) ! ! xw value of x for first i grid point to plot ! xe value of x for last i grid point to plot ! ys value of y for first j grid point to plot ! yn value of y for last j grid point to plot ! ! nx first dimension of b ! ibgn index of first i grid point to plot ! iend index of last i grid point to plot ! ! ny second dimension of b ! jbgn index of first j grid point to plot ! jend index of last j grid point to plot ! ! ! time time of data in seconds ! ! runname character string decribing run ! ! factor scaling factor for winds ! V*factor wind vectors are plotted ! slicopt slice orientation indicator ! slicopt = 1, x-y slice of u,v at z index kslice is plotted. ! slicopt = 2, x-z slice of u,w at y index jslice is plotted. ! slicopt = 3, y-z slice of v,w at x index islice is plotted. ! slicopt = 4, x-y slice of u,v at z index islice is plotted. ! slicopt = 5, xy-z cross section of wind islice is plotted. ! slicopt = 6, data field on constant p-level is plotted. ! slicopt = 0, all of the three slices above are plotted. ! ! WORK ARRAYS: ! ! tem1 Temporary work array. ! tem2 Temporary work array. ! tem3 Temporary work array. ! tem4 Temporary work array. ! tem5 Temporary work array. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! INTEGER :: nx,ny ! REAL :: u(nx,ny) REAL :: v(nx,ny) REAL :: x(nx,ny) REAL :: y(nx,ny) ! CHARACTER (LEN=*) :: runname CHARACTER (LEN=*) :: label REAL :: xw,xe,dx,ys,yn,dy INTEGER :: ibgn,iend,ist, jbgn,jend,jst REAL :: time,factor INTEGER :: slicopt INTEGER :: iunits, itype COMMON /windvtr/iunits, itype CHARACTER (LEN=12) :: varname COMMON /varplt1/ varname REAL :: xw1,xe1,ys1,yn1 COMMON /xuvpar/xw1,xe1,ys1,yn1 ! !----------------------------------------------------------------------- ! ! Work arrays: tem1,tem2,tem3,tem4,tem5 of size at least ! max( nx*ny, nx*nz, ny*nz). ! !----------------------------------------------------------------------- ! REAL :: tem1(*),tem2(*),tem3(*),tem4(*),tem5(*) REAL :: tem6(*) ! !----------------------------------------------------------------------- ! ! Common blocks for plotting control parameters ! !----------------------------------------------------------------------- ! REAL :: x01,y01 ! the first point of interpolation REAL :: x02,y02 ! the second point of interpolation REAL :: zlevel ! the given height of the slice REAL :: sinaf,cosaf,dist,sqrtdxy COMMON /slicev/x01,y01,x02,y02,sinaf,cosaf,dist,sqrtdxy COMMON /sliceh/zlevel INTEGER :: ovrobs,obsset,obscol,obs_marktyp REAL :: obs_marksz COMMON /obspar/ ovrobs,obsset,obscol,obs_marktyp, obs_marksz ! !----------------------------------------------------------------------- ! ! Misc. local Variables ! !----------------------------------------------------------------------- ! INTEGER :: i,j,ij,istep,jstep,length,isize,jsize REAL :: uunit CHARACTER (LEN=6) :: timhms CHARACTER (LEN=120) :: title ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! INTEGER :: trnplt ! flag to plot terain (1 or 0) REAL :: hterain(nx,ny) ! The height of the terrain. INTEGER :: ovrtrn ! overlay terrain option (0/1) REAL :: trninc,trnmin, trnmax ! terrain interval minimum, maximum REAL :: ztmin,ztmax COMMON /trnpar/ trnplt,ovrtrn,trninc,trnmin, trnmax,ztmin,ztmax ! INTEGER :: smooth COMMON /smoothopt/smooth INTEGER :: wrtflag, llabel CHARACTER (LEN=80) :: levlab CHARACTER (LEN=50) :: timelab CHARACTER (LEN=25) :: timestring COMMON /timelev/wrtflag,timelab, levlab, timestring CHARACTER (LEN=120) :: label_copy ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! isize=(iend-ibgn)+1 jsize=(jend-jbgn)+1 ! ! !----------------------------------------------------------------------- ! ! Set up terrain, if needed. ! !----------------------------------------------------------------------- ! label_copy = label llabel = 120 CALL xstrlnth(label_copy, llabel) CALL xpscmnt('Start plotting '//label_copy(1:llabel)) IF(trnplt == 1 .OR.trnplt == 2 .OR. ovrtrn == 1) THEN DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem5(ij)=hterain(i,j) END DO END DO END IF CALL cvttim ( time, timhms) IF( timhms(1:1) == '0' ) timhms(1:1)=' ' WRITE(timelab,'(''T='',F8.1,A)') time, & ' s ('//timhms(1:2)//':'//timhms(3:4)//':'//timhms(5:6)//')' CALL get_time_string ( time, timestring) ! length=50 ! CALL strmin(timelab,length) ! write(timelab,'(a,'' '',a)') timestring(1:21), timelab(1:length) ! print*,'in vtrsfc', timelab DO j=jbgn,jend DO i=ibgn,iend ij = i-ibgn+1 + (j-jbgn)*isize tem1(ij) = -9999.0 tem2(ij) = -9999.0 IF( IF( tem3(ij)=x(i,j) tem4(ij)=y(i,j) END DO END DO levlab = 'First level above ground (surface)' WRITE(title,'(''U-V '',A)') label length = 120 CALL strlnth( title, length ) CALL strmin ( title, length) uunit = 10.0 CALL xvmode(1) istep = ist jstep = jst DO i=1,smooth CALL smooth9pmv(tem1,isize,jsize,1,isize,1,jsize,tem6) CALL smooth9pmv(tem2,isize,jsize,1,isize,1,jsize,tem6) END DO CALL vtr2d(tem1,tem2,tem3,tem4, uunit, xw,xe,dx,ys,yn,dy, & isize,istep,jsize,jstep, & title(1:length),runname, 1, & tem5,slicopt) CALL xpscmnt('End plotting '//label_copy(1:llabel)) RETURN END SUBROUTINE vtrsfc ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE SET_INTERVAL ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE set_interval(z, m,n,zmin1,zmax1,ctmin,ctmax,cl) 1 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Limited contour interval when uinc = -9999. ! !----------------------------------------------------------------------- ! ! AUTHOR: Min Zou ! ! !----------------------------------------------------------------------- ! ! INPUT: ! ! z 2-D array ! m,n dimension of 2-d array ! zmin1 the minimum value of 2-D array ! zmax1 The maximum value of 2-D array ! ctmin the input minimum value ! ctmax the input maximum value ! cl the intervals ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! REAL :: z(m,n),cl(*) REAL :: zmin, zmax REAL :: zinc COMMON /xclm19/ nmin, nmax COMMON /xcrf17/clref,lcptn,labtyp,iclf,lhilit,ihlf,kct0 COMMON /zchole/ nhole,specia,nvtrbadv COMMON /xoutch/ nch IF(ctmin == 0.0 .AND. ctmax == 0.0) THEN cl(2)=cl(1)+ xfinc(zmax1-zmin1)/2 ELSE cl(2)=cl(1)+ xfinc(ctmax-ctmin)/2 END IF IF(cl(2)-cl(1) == 0.0) cl(2)=cl(1)+1.0 zinc = cl(2)-cl(1) ncmin=nmin ncmax=nmax diff=ctmax-ctmin IF(diff-ABS( zinc)*1.0E-6 > 0.0) THEN GO TO 4 END IF 999 WRITE(nch,'(a,a)') & ' Bad first guess of contour increment or field is constant', & ', number of contours is one.' ncnt=1 cl(1)= ctmin 899 RETURN 4 kcount=0 1 CONTINUE eps=0.1*zinc kcount=kcount+1 IF( kcount > 20) GO TO 998 kzinc=(ctmin-clref)/zinc zmin=kzinc*zinc+clref kzinc=(ctmax-clref)/zinc zmax=kzinc*zinc+clref IF(ctmin-clref > 0.0) zmin=zmin+zinc IF(ctmax-clref < 0.0) zmax=zmax-zinc ! clv=zmin-zinc ncnt=0 6 clv=clv+zinc IF(clv-zmax-eps > 0.0) THEN GO TO 8 END IF 3 ncnt=ncnt+1 IF(ncnt > ncmax) THEN zinc=zinc*2 WRITE(nch,1000) ncmax, zinc 1000 FORMAT(' Number of contours > ',i3,' ,Zinc is doubled. Zinc=' & ,e10.3) GO TO 1 END IF IF( ABS( clv-clref ) < eps ) clv=clref cl(ncnt)=clv GO TO 6 8 CONTINUE IF( ncnt < ncmin) THEN zinc=zinc/2 WRITE(nch,2000) ncmin,zinc 2000 FORMAT(' Number of contours < ',i3,' ,Zinc is halved. Zinc=' & ,e10.3) GO TO 1 END IF WRITE(nch,'('' * NUMBER OF CONTOURS= '',I5,'' MIN='',E12.4, & & '' MAX='', e12.4,'' inc='',e12.5 )') & ncnt,ctmin,ctmax,zinc IF( zmin1 >= ctmin .AND. zmax1 <= ctmax) THEN zinc = cl(2) - cl(1) WRITE(nch,'(''SET MINIMUM CONTOUR INTERVAL IS'',E12.4, & & '' ctmin='',e12.4,'' ctmax='',e12.4 )')zinc,ctmin,ctmax CALL xctref(zinc) CALL xnctrs( 1,300) ELSE WRITE(nch,'(''NO NEED SET MINIMUM CONTOUR INTERVAL'')' ) WRITE(nch,'(''CNTOUR INTERVAL IS SET AUTOMATICALLY'')' ) cl(2)=cl(1)+ xfinc(zmax1-zmin1)/2 IF(cl(2)-cl(1) == 0.0) cl(2)=cl(1)+1.0 END IF RETURN 998 WRITE(nch,*)' Contour levels can not be selected by XCNTLV.' WRITE(nch,*) & ' Plz alter input contour interval or limits of contour number' END SUBROUTINE set_interval ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE XRCH1 ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE xrch1( r,ch,lch) 1,3 ! Return real number R as a character string in automatically set format REAL :: r CHARACTER (LEN=20) :: str CHARACTER (LEN=20) :: ch CALL get_format(r,str) IF(ABS(r-0.0) < 1.e-20) THEN WRITE(ch,'(F3.1)') r ELSE WRITE(ch,str) r END IF lch=20 CALL strlnth( ch, lch) CALL strmin ( ch, lch) RETURN END SUBROUTINE xrch1 ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE GET_FORMAT ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE get_format(r,ch) 1 INTEGER :: npoz CHARACTER (LEN=1) :: FORM,ndrob CHARACTER (LEN=20) :: ch WRITE(ch,10)r 10 FORMAT(g11.4) DO i=20,1,-1 IF(ch(i:i) == '0'.OR.ch(i:i) == ' ') THEN ch(i:i)=' ' ELSE GO TO 1 END IF END DO 1 CONTINUE npoz=0 ndot=0 nmant=0 ndrob=' ' FORM='F' DO i = 1,20 IF(ch(i:i) /= ' ' ) npoz=npoz+1 IF(ch(i:i) == 'E') FORM='E' IF(ndrob == '.'.AND.ch(i:i) /= ' ') ndot=ndot+1 IF(ch(i:i) == '.') ndrob='.' IF(FORM /= 'E') nmant=npoz END DO npoz=npoz IF(FORM == 'F') THEN IF(ndot /= 0) THEN WRITE(ch,20) '(',FORM,npoz,'.',ndot,')' ELSE WRITE(ch,20) '(',FORM,npoz,'.',ndot,')' END IF ELSE IF(FORM == 'E') THEN ch = '(1PE20.2)' ELSE WRITE(ch,20) '(',FORM,npoz,'.',nmant,')' END IF 20 FORMAT(a1,a1,i1,a1,i1,a1) RETURN END SUBROUTINE get_format !################################################################## !################################################################## !###### ###### !###### SUBROUTINE DRAWMAP ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE drawmap(nunit) 1 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! This subroutine will plot the map ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! 6/2/97 Min Zou ! Read multiple mapfile only once. Using differnt line style to ! plot mapdata. ! !----------------------------------------------------------------------- ! ! INPUT: ! ! nunit the channel of the mapfile data ! mapfile character of map file name ! !----------------------------------------------------------------------- ! ! Variable Declarations ! !----------------------------------------------------------------------- ! IMPLICIT NONE INCLUDE 'arpsplt.inc' INTEGER :: nunit,i INTEGER :: lmapfile CHARACTER (LEN=132) :: mapfile(maxmap) INTEGER :: mapgrid,mapgridcol, kolor REAL :: latgrid,longrid INTEGER :: nmapfile,mapcol(maxmap),mapline_style(maxmap) COMMON /mappar1/nmapfile,mapcol,mapline_style,mapfile COMMON /mappar2/mapgrid,mapgridcol,latgrid,longrid REAL :: x1,x2,y1,y2 CALL xpscmnt('Start of map plotting ') CALL xqmap (x1,x2,y1,y2) CALL xwindw(x1,x2,y1,y2) CALL xqcolor(kolor) DO i=1,nmapfile CALL xcolor(mapcol(i)) IF(mapline_style(i) == 1) THEN CALL xthick(1) CALL xbrokn(6,3,6,3) ELSE IF(mapline_style(i) == 2) THEN CALL xthick(1) ELSE IF(mapline_style(i) == 3) THEN CALL xthick(3) CALL xfull END IF lmapfile=132 CALL xstrlnth(mapfile(i), lmapfile) ! print*,'plotting ',mapfile(i)(1:lmapfile) CALL xdrawmap(nunit,mapfile(i)(1:lmapfile),latgrid,longrid) END DO CALL xcolor(kolor) CALL xfull CALL xwdwof CALL xpscmnt('End of map plotting ') RETURN END SUBROUTINE drawmap ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE PLTOBS ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE pltobs(obopt) 2 ! ! !----------------------------------------------------------------------- ! ! PURPOSE: ! Plots observations on an arpsplt contour map. ! !----------------------------------------------------------------------- ! ! INPUT: ! obopt Plotting option ! 1 Plot data in obs1 as characters ! 2 Plot data in obs1 and obs2 as characters ! 3 Plot wind arrows with obs1 as u and obs2 as v. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INCLUDE 'arpsplt.inc' ! ! Arguments ! INTEGER :: obopt ! INTEGER :: nobs REAL :: latob(mxsfcob) REAL :: lonob(mxsfcob) REAL :: obs1(mxsfcob) REAL :: obs2(mxsfcob) ! COMMON /sfc_obs1/ nobs COMMON /sfc_obs2/ latob,lonob,obs1,obs2 ! ! Plotting parameters ! CHARACTER (LEN=1) :: cross PARAMETER(cross='+') ! INTEGER :: ovrobs,obsset,obscol,obs_marktyp REAL :: obs_marksz COMMON /obspar/ ovrobs,obsset,obscol,obs_marktyp, obs_marksz ! ! Plotting common blocks ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor REAL :: ctinc,ctmin,ctmax,vtunt ! contour interval and vector unit COMMON /incunt/ ctinc,ctmin,ctmax,vtunt REAL :: xleng,vunit COMMON /vecscl/ xleng,vunit INTEGER :: iunits, itype COMMON /windvtr/iunits, itype ! ! Misc local variables ! INTEGER :: iob REAL :: orgmag,obmag,yoff,yoff2 REAL :: x1,x2,y1,y2 REAL :: xob,yob CHARACTER (LEN=4) :: chplot INTEGER :: imkrfil ! ! Set-up plotting space and zxplot variables ! CALL xcolor(lbcolor) CALL xqmap (x1,x2,y1,y2) CALL xwindw(x1,x2,y1,y2) CALL xchori(0.0) CALL xqchsz(orgmag) obmag=0.8*orgmag yoff=0.5*orgmag yoff2=2.*yoff CALL xchsiz(obmag) IF(obopt == 1) THEN CALL xcolor(obscol) CALL xmrksz(obs_marksz) DO iob=1,nobs IF(obs1(iob) > -98. .AND. obs1(iob) < 500.) THEN CALL xlltoxy(1,1,latob(iob),lonob(iob),xob,yob) WRITE(chplot,810) nint(obs1(iob)) 810 FORMAT(i4) CALL xcharc((0.001*xob),(0.001*yob+yoff),chplot) ! call XCHARC((0.001*xob),(0.001*yob),cross) CALL xmarker((0.001*xob),(0.001*yob),obs_marktyp) IF(obs_marktyp > 5) THEN CALL xqmkrfil(imkrfil) CALL xmkrfil(1) CALL xmarker((0.001*xob),(0.001*yob),MOD(obs_marktyp,5)) CALL xmkrfil(imkrfil) END IF END IF END DO ELSE IF(obopt == 2) THEN CALL xcolor(obscol) CALL xmrksz(obs_marksz) DO iob=1,nobs CALL xlltoxy(1,1,latob(iob),lonob(iob),xob,yob) IF(obs1(iob) > -98. .AND. obs1(iob) < 500.) THEN WRITE(chplot,810) nint(obs1(iob)) CALL xcharc((0.001*xob),(0.001*yob+yoff),chplot) ! call XCHARC((0.001*xob),(0.001*yob),cross) IF(obs_marktyp > 5) THEN CALL xqmkrfil(imkrfil) CALL xmkrfil(1) CALL xmarker((0.001*xob),(0.001*yob),MOD(obs_marktyp,5)) CALL xmkrfil(imkrfil) END IF END IF IF(obs2(iob) > -98. .AND. obs2(iob) < 500.) THEN WRITE(chplot,810) nint(obs2(iob)) CALL xcharc((0.001*xob),(0.001*yob+yoff),chplot) ! call XCHARC((0.001*xob),(0.001*yob),cross) CALL xmarker((0.001*xob),(0.001*yob),obs_marktyp) IF(obs_marktyp > 5) THEN CALL xqmkrfil(imkrfil) CALL xmkrfil(1) CALL xmarker((0.001*xob),(0.001*yob),MOD(obs_marktyp,5)) CALL xmkrfil(imkrfil) END IF END IF END DO ELSE IF(obopt == 3) THEN CALL xcolor(obscol) CALL xmrksz(obs_marksz) DO iob=1,nobs IF(obs1(iob) > -98. .AND. obs1(iob) < 500. .AND. & obs2(iob) > -98. .AND. obs2(iob) < 500.) THEN CALL xlltoxy(1,1,latob(iob),lonob(iob),xob,yob) xob=0.001*xob yob=0.001*yob IF( xob > x1 .AND. xob < x2 .AND. yob > y1 .AND. yob < y2 ) THEN CALL xarrow(obs1(iob),obs2(iob),xob,yob,xleng,vunit) CALL xmarker(xob,yob,obs_marktyp) IF(obs_marktyp > 5) THEN CALL xqmkrfil(imkrfil) CALL xmkrfil(1) CALL xmarker(xob,yob,MOD(obs_marktyp,5)) CALL xmkrfil(imkrfil) END IF END IF END IF END DO END IF CALL xcolor(lbcolor) CALL xchsiz(orgmag) CALL xfull CALL xwdwof RETURN END SUBROUTINE pltobs ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE PLTSTA ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE pltsta(a,b,x,y,m,n,flag,slicopt) 6,5 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! This subroutine will plot some station information. ! !----------------------------------------------------------------------- ! ! !----------------------------------------------------------------------- ! ! AUTHOR: ! Min Zou (6/1/97) ! ! Modification history: ! !----------------------------------------------------------------------- ! ! INPUT: ! ! a, b 2-dimension array of Variable ! m, n Array dimensions ! x, y x-coord and y-coord of the staions ! flag a flag for different plot ! slicopt slice orientation indicator ! !----------------------------------------------------------------------- ! ! Variable Declarations ! !----------------------------------------------------------------------- ! IMPLICIT NONE INCLUDE 'arpsplt.inc' INTEGER :: m,n REAL :: a(m,n) REAL :: b(m,n) REAL :: x(m,n) REAL :: y(m,n) INTEGER :: nsta,nstapro(mxstalo),nstatyp(mxstalo) REAL :: latsta(mxstalo), lonsta(mxstalo) CHARACTER (LEN=5) :: s_name(mxstalo) INTEGER :: ovrstaopt INTEGER :: ovrstam,staset,ovrstan,ovrstav,stacol,markprio,wrtstax INTEGER :: nsta_typ,sta_typ(10),sta_marktyp(10), sta_markcol(10) REAL :: sta_marksz(10) REAL :: wrtstad COMMON /sta_par/ ovrstaopt,ovrstam,staset,ovrstan,ovrstav,stacol, & markprio, nsta_typ,sta_typ,sta_marktyp, & sta_markcol,sta_marksz,wrtstax,wrtstad COMMON /sta_loc/latsta,lonsta,nstatyp,nstapro,nsta COMMON /sta_loc1/s_name REAL :: xob(mxstalo), yob(mxstalo),aob(mxstalo),bob(mxstalo) COMMON /xob_yob/xob, yob INTEGER :: LEN,i,j ! REAL :: x01,x02,y01,y02 REAL :: sinaf,cosaf,dist,sqrtdxy COMMON /slicev/x01,y01,x02,y02,sinaf,cosaf,dist,sqrtdxy ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! INTEGER :: layover COMMON /laypar/ layover ! CHARACTER (LEN=12) :: varname COMMON /varplt1/ varname ! REAL :: xori1,xori2,yori1,yori2,zori1,zori2 COMMON /tmphc1/xori1,xori2,yori1,yori2,zori1,zori2 ! REAL :: xleng,vunit COMMON /vecscl/ xleng,vunit INTEGER :: iunits, itype COMMON /windvtr/iunits, itype REAL :: x_tmp COMMON /tmphc2/ x_tmp ! REAL :: x1,x2,y1,y2 REAL :: orgmag,obmag,yoff,yoff2,xoff, xoff2 CHARACTER (LEN=30) :: ctmp ! INTEGER :: flag,slicopt,fg REAL :: xdist, ydist,xd0,yd0,xa,xb SAVE fg REAL :: xleng0, spd, dir, istand INTEGER :: iunits0, imkrfil ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! calculate xob and yob on the xy plane CALL xwindw(xori1,xori2,yori1,yori2) IF(fg == 0) THEN CALL xlltoxy(nsta,1,latsta,lonsta,xob,yob) DO i= 1,nsta xob(i) = xob(i)*0.001 yob(i) = yob(i)*0.001 END DO fg=1 END IF ! CALL xqmap (x1,x2,y1,y2) CALL xwindw(x1,x2,y1,y2) CALL xchori(0.0) CALL xqchsz(orgmag) obmag=0.8*orgmag yoff=0.5*orgmag yoff2=3.*yoff xoff = 0.001*(x2-x1) xoff2 = 4.*xoff CALL xchsiz(obmag) IF(ovrstav == 1 ) THEN ! interpolation CALL intepo (nsta,xob,yob,aob,m,n,x,y,a) IF(varname(1:6) == 'vtrplt' .OR. varname(1:6) == 'vtpplt') & CALL intepo (nsta,xob,yob,bob,m,n,x,y,b) END IF IF (flag == 1) THEN CALL xcolor(stacol) IF(slicopt == 5) THEN xa = (y02-y01)/(x02-x01) xb = y01 - xa*x01 END IF DO i = 1,nsta IF(nstapro(i) <= markprio) THEN IF(wrtstax == 1 )THEN CALL xwindw(x1,x2-0.005*(x2-x1), & y1-0.3*(y2-y1),y1) LEN=5 CALL strlnth(s_name(i),LEN) CALL xchori(90.) IF(slicopt == 2) THEN CALL xwindw(xori1,xori2-0.005*(xori2-xori1), & y1-0.3*(y2-y1),y1) IF ( (xob(i) <= xori2.AND.xob(i) >= xori1) & .AND. (yob(i) <= yori2.AND.yob(i) >= yori1) ) THEN IF( ABS(yob(i)-x_tmp) <= wrtstad ) THEN ! CALL XCHARR((xob(i)),y1-1.75*yoff2, CALL xcharr((xob(i)),y1-1.50*yoff2, & s_name(i)(1:LEN)) END IF END IF ELSE IF(slicopt == 3) THEN CALL xwindw(yori1,yori2-0.005*(yori2-yori1), & y1-0.3*(y2-y1),y1) IF ( (xob(i) <= xori2.AND.xob(i) >= xori1) & .AND. (yob(i) <= yori2.AND.yob(i) >= yori1) ) THEN IF( ABS(xob(i)-x_tmp) <= wrtstad ) THEN CALL xcharr((yob(i)),y1-1.75*yoff2, & s_name(i)(1:LEN)) END IF END IF ELSE IF( slicopt == 5) THEN CALL xwindw(x1,x2-0.005*(x2-x1), & y1-0.3*(y2-y1),y1) IF ( (xob(i) <= xori2.AND.xob(i) >= xori1) & .AND. (yob(i) <= yori2.AND.yob(i) >= yori1) ) THEN xd0 = 1./(xa*xa+1.0)*((yob(i)-xb)*xa+xob(i)) yd0 = xa*xd0+xb xdist = SQRT((x01-xd0)*(x01-xd0) + (y01-yd0)*(y01-yd0)) ydist = SQRT((xob(i)-xd0)*(xob(i)-xd0)+ & (yob(i)-yd0)*(yob(i)-yd0)) xdist = xdist+x1 IF(ABS(ydist) <= wrtstad ) THEN CALL xcharr(xdist,y1-1.75*yoff2, & s_name(i)(1:LEN)) END IF END IF END IF CALL xchori(0.) END IF END IF END DO ELSE IF(flag == 0) THEN CALL xwindw(x1,x2,y1,y2) DO i = 1,nsta IF( (xob(i) >= x1.AND.xob(i) <= x2) .AND. & (yob(i) >= y1.AND.yob(i) <= y2) ) THEN IF(nstapro(i) <= markprio) THEN IF(ovrstan == 1) THEN LEN=5 CALL strlnth(s_name(i),LEN) CALL xcharc((xob(i)),(yob(i)-yoff2), & s_name(i)(1:LEN)) END IF IF(ovrstam == 1) THEN DO j=1,nsta_typ CALL xmrksz(sta_marksz(j)) CALL xcolor(sta_markcol(j)) IF(nstatyp(i) == sta_typ(j)) THEN CALL xmarker((xob(i)),(yob(i)), & sta_marktyp(j)) IF(sta_marktyp(j) > 5) THEN CALL xqmkrfil(imkrfil) CALL xmkrfil(1) CALL xmarker((xob(i)),(yob(i)), & MOD(sta_marktyp(j) ,5)) CALL xmkrfil(imkrfil) END IF IF(ovrstan == 1) THEN LEN=5 CALL strlnth(s_name(i),LEN) CALL xcharc((xob(i)),(yob(i)-yoff2), & s_name(i)(1:LEN)) END IF END IF END DO END IF IF(ovrstav == 1) THEN CALL xcolor(stacol) IF(varname(1:6) == 'vtrplt' .OR. varname(1:6) == 'vtpplt') THEN ! IF(i.eq.1) THEN xleng0=xleng*0.0004 IF(iunits == 2 ) THEN iunits0=2 istand = 10. WRITE(ctmp,'(a30)')'10 knots' ELSE IF (iunits == 3) THEN iunits0=2 istand = 10. WRITE(ctmp,'(a30)')'10 MPH' ELSE IF(iunits == 1) THEN iunits0=1 istand = 5. WRITE(ctmp,'(a30)')'5 m/s' END IF ! ENDIF IF(aob(i) /= -9999. .AND. bob(i) /= -9999.) THEN spd = SQRT(aob(i)*aob(i)+bob(i)*bob(i)) dir = ATAN2(-1.*aob(i),-1.*bob(i))*180./3.1415926 IF(dir <= 0.) dir = 360.+dir ! CALL barb((xob(i)), ! : (yob(i)),dir,spd,iunits0-1, xleng0) CALL xbarb(aob(i),bob(i),xob(i),yob(i), & iunits0,xleng*0.65,2) END IF ELSE IF(aob(i) /= -9999.) THEN CALL xrch(aob(i),ctmp,LEN) IF(layover == 0) CALL xcharr((xob(i)-xoff2), & (yob(i)+yoff),ctmp(1:LEN)) IF(layover == 1) CALL xcharl((xob(i)+xoff2), & (yob(i)+yoff) ,ctmp(1:LEN)) IF(layover == 2) CALL xcharc((xob(i)+xoff2), & (yob(i)-yoff),ctmp(1:LEN)) IF(layover == 3) CALL xcharl((xob(i))+xoff2, & (yob(i)-yoff) ,ctmp(1:LEN)) END IF END IF END IF END IF END IF END DO END IF CALL xcolor(lbcolor) CALL xchsiz(orgmag) CALL xfull CALL xwdwof RETURN END SUBROUTINE pltsta ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE RUNLAB ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE runlab(runname) 36 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Plot a run label at the lower left cornor of the picture frame. ! !----------------------------------------------------------------------- ! ! INPUT: ! ! runname character string of run label ! !----------------------------------------------------------------------- ! ! Variables Declarations ! !----------------------------------------------------------------------- ! IMPLICIT NONE CHARACTER (LEN=*) :: runname REAL :: xl, xr, yb, yt, rotang, xlimit, ylimit INTEGER :: nopic, nxpic, nypic ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! CALL xqmap(xl,xr,yb,yt) CALL xqnpic(nopic) CALL xqspac(nxpic, nypic, rotang, xlimit, ylimit) IF( rotang == 0.0 ) THEN IF(nopic == nxpic*nypic -(nxpic-1)) THEN CALL xcharl( xl, yb-0.15*(yt-yb), runname ) END IF ELSE IF(nopic == nypic*nxpic -(nypic-1)) THEN CALL xcharl( xl, yb-0.15*(yt-yb), runname ) END IF END IF RETURN END SUBROUTINE runlab ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE VPROFIL ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE vprofil(nx,ny,nz,var,xc,yc,zpc,plwr,pupr, & 25,2 xpnt,ypnt,npoints,zlwr,zupr,xcaptn,ycaptn,npicprof, & profil,height) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! This subroutine will plot the vertical profiles of a given ! variable through points (xpnt(i),ypnt(i),i=1,npoints). ! !----------------------------------------------------------------------- ! ! !----------------------------------------------------------------------- ! ! AUTHOR: ! Adwait Sathye (2/28/94) ! ! Modification history: ! ! 4/18/94, (Ming Xue) ! Major overhaul. Many temporary arrays removed. New frame option ! added. ! ! 9/18/1995 (Ming Xue) ! Fixed a problem in the code that determines kbgn and kend. ! ! 10/8/1996 (Y. Richardson) ! Corrected a bug in the interpolation weights. ! !----------------------------------------------------------------------- ! ! !----------------------------------------------------------------------- ! ! INPUT: ! ! nx, ny, nz Array dimensions ! var Variable data array ! xc,yc,zpc The coordinate of input data var. ! plwr,pupr Lower and upper bounds for the horiz. axis of profile ! xpnt, ypnt Arrays containing the X and Y locations of the ! mulitple profiles to be plotted ! npoints Number of profile points to be plotted ! zlwr, zupr Bounds in the vertical direction ! xcaptn Caption for the X axis ! ycaptn Caption for the Y axis ! ! Work arrays: ! ! profil,height Temporary arrays ! !----------------------------------------------------------------------- ! ! !----------------------------------------------------------------------- ! ! Variable Declarations ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! !----------------------------------------------------------------------- ! ! Variables passed in ! !----------------------------------------------------------------------- ! INTEGER :: nx, ny, nz REAL :: var(nx,ny,nz) REAL :: xc(nx,ny,nz), yc(nx,ny,nz), zpc(nx,ny,nz) REAL :: plwr,pupr REAL :: zlwr, zupr INTEGER :: npoints REAL :: xpnt(npoints), ypnt(npoints) CHARACTER (LEN=*) :: xcaptn CHARACTER (LEN=*) :: ycaptn INTEGER :: npicprof REAL :: profil(nz,npoints) REAL :: height(nz,npoints) LOGICAL :: multiprof ! !----------------------------------------------------------------------- ! ! Temporary local variables ! !----------------------------------------------------------------------- ! REAL :: lower, upper, zmin, zmax REAL :: x1, x2, y1, y2 REAL :: a1, a2, a3, a4 INTEGER :: i, j, k, ix, jy,kbgn,kend,ip,lchar REAL :: dx,dy,temp,hmaxk,hmink CHARACTER (LEN=80) :: ch REAL :: lblmag, ctrlbsiz, axlbsiz COMMON /labmag/ lblmag, ctrlbsiz, axlbsiz ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! multiprof = .false. IF( npicprof == 0 .AND. npoints > 1 ) multiprof = .true. ! !----------------------------------------------------------------------- ! ! If lower boundary is bigger, then swap boundaries. ! !----------------------------------------------------------------------- ! IF (plwr > pupr) THEN lower = pupr upper = plwr ELSE lower = plwr upper = pupr END IF ! !----------------------------------------------------------------------- ! ! Find corresponding coordinates for the boundaries in the Z ! dimension. IF both have been set to 0, use the boundary values. ! Else, loop through the zpc array and find the location of the ! point. ! !----------------------------------------------------------------------- ! dx = xc(2,1,1)-xc(1,1,1) dy = yc(1,2,1)-yc(1,1,1) zmin = 0.5*(zpc(1,1,1)+zpc(1,1,2)) zmax = 0.5*(zpc(1,1,nz-1)+zpc(1,1,nz-2)) DO j=1,ny-1 DO i=1,nx-1 zmin=MIN( zmin, 0.5*(zpc(i,j,1)+zpc(i,j,2)) ) zmax=MAX( zmax, 0.5*(zpc(i,j,nz-1)+zpc(i,j,nz-2)) ) END DO END DO IF( zlwr /= zupr ) THEN zmin = MAX(zmin, zlwr) zmax = MIN(zmax, zupr) END IF IF( zmax < zmin) WRITE(6,'(a, f10.3, a, f10.3)') & 'Warning: zmax is less then zmin. Check input data zprofbgn', & zlwr , 'and zprofbgn',zupr DO ip = 1, npoints ix = INT ( (xpnt(ip) - xc(1,1,1))/dx ) + 1 jy = INT ( (ypnt(ip) - yc(1,1,1))/dy ) + 1 ix = MIN(MAX(1,ix), nx-2) jy = MIN(MAX(1,jy), ny-2) WRITE(6,'(1x,2a,2(a,f10.3),2(a,i4) )') & 'Plotting ',xcaptn,' profile through (', & xpnt(ip),',',ypnt(ip),') km, at i=',ix,' j=',jy ! !----------------------------------------------------------------------- ! ! Interpolate the data value and its height to the specified point. ! !----------------------------------------------------------------------- ! a1 = ABS ((xc(ix ,jy ,1)-xpnt(ip))*(yc(ix ,jy ,1)-ypnt(ip))) a2 = ABS ((xc(ix+1,jy ,1)-xpnt(ip))*(yc(ix+1,jy ,1)-ypnt(ip))) a3 = ABS ((xc(ix+1,jy+1,1)-xpnt(ip))*(yc(ix+1,jy+1,1)-ypnt(ip))) a4 = ABS ((xc(ix ,jy+1,1)-xpnt(ip))*(yc(ix ,jy+1,1)-ypnt(ip))) DO k = 1,nz-1 profil(k,ip)= (a3*var(ix ,jy,k)+a2*var(ix ,jy+1,k)+ & a4*var(ix+1,jy,k)+a1*var(ix+1,jy+1,k)) & /(a1 + a2 + a3 + a4) height(k,ip)= (a3*zpc(ix ,jy,k)+a2*zpc(ix ,jy+1,k)+ & a4*zpc(ix+1,jy,k)+a1*zpc(ix+1,jy+1,k)) & /(a1 + a2 + a3 + a4) END DO END DO kbgn = nz-1 DO k=nz-1,1,-1 hmaxk = height(k,1) DO ip=1,npoints hmaxk = MAX(hmaxk,height(k,ip)) END DO IF( hmaxk >= zmin) kbgn = k END DO kend = 1 DO k=1,nz-1 hmink = height(k,1) DO ip=1,npoints hmink = MIN(hmink,height(k,ip)) END DO IF( hmink <= zmax) kend=k END DO ! !----------------------------------------------------------------------- ! ! If input bounds for the profile are zero, use the min. and max. ! in the profile as the lower and upper bounds for the horizontal ! axis. ! !----------------------------------------------------------------------- ! IF( plwr == 0.0 .AND. pupr == 0.0 ) THEN lower = profil(kbgn,1) upper = profil(kend,1) DO ip=1,npoints DO k = kbgn,kend lower = MIN(lower, profil(k,ip)) upper = MAX(upper, profil(k,ip)) END DO END DO ELSE lower = plwr upper = pupr END IF ! !----------------------------------------------------------------------- ! ! If the lower and upper bounds are equal, set the horizontal ! axis scale to 1.0. ! !----------------------------------------------------------------------- ! IF ((lower == 0.0 .AND. upper == 0.0).OR.upper == lower) upper = lower+1.0 ! !----------------------------------------------------------------------- ! ! Start to plot the profile... ! !----------------------------------------------------------------------- ! DO ip=1,npoints IF( (.NOT.multiprof) .OR. (multiprof.AND.ip == 1) ) THEN CALL xnwpic CALL xaxtik(1, 1) CALL xaxant(-1, -1) CALL xmap (lower, upper, zmin, zmax) CALL xaxnsz ( axlbsiz*(zmax-zmin)*lblmag ) CALL xqmap(x1,x2,y1,y2) CALL xchsiz(0.03*(y2-y1)*lblmag) CALL xchori(0.0) IF( .NOT.multiprof ) THEN lchar = LEN( xcaptn) ch = xcaptn WRITE(ch(lchar+1:lchar+33), '(a,f13.3,a,f13.3,a)') & ' at (',xpnt(ip),',',ypnt(ip),')' lchar = lchar+33 CALL strmin(ch(1:lchar), lchar) CALL xcharc((x1+x2)*0.5, y1-(y2-y1)*0.10, ch(1:lchar)) ELSE CALL xcharc((x1+x2)*0.5, y1-(y2-y1)*0.10, xcaptn ) END IF ! !----------------------------------------------------------------------- ! ! Check if the points lie on one side of the axis. If the points are ! all positive, draw the y-axis on the left border, if all points are ! negative, draw the y-axis on the right border. If points lie on both ! sides, draw the y-axis through x=0.0. ! !----------------------------------------------------------------------- ! temp = lower * upper IF (temp > 0.0) THEN IF (lower > 0.0) THEN CALL xaxes(lower,0.0,zmin,0.0) CALL xchori(90.0) CALL xcharc(x1-0.12*(x2-x1), (y1+y2)*0.5, ycaptn) ELSE CALL xaxant(-1, 1) CALL xaxtik(1, -1) CALL xaxes(upper,0.0,zmin,0.0) CALL xchori(90.0) CALL xcharc(x1-0.05*(x2-x1), (y1+y2)*0.5, ycaptn) END IF ELSE CALL xaxes(0.0,0.0,zmin,0.0) CALL xchori(90.0) CALL xcharc(x1-0.10*(x2-x1), (y1+y2)*0.5, ycaptn) END IF END IF CALL xchori(0.0) CALL xbordr CALL xfull ! !----------------------------------------------------------------------- ! ! The first plot is labeled `A'. The subsequent plots will be `B'... ! !----------------------------------------------------------------------- ! IF( multiprof ) THEN CALL xlbon CALL xlabel(CHAR(64+ip)) ch(1:1) = CHAR(64+ip) ch(2:2) = ' ' lchar = 2 WRITE(ch(lchar+1:lchar+33), '(a,f13.2,a,f13.2,a)') & ' at (',xpnt(ip),',',ypnt(ip),')' lchar = lchar+33 CALL strmin(ch(1:lchar), lchar) CALL xqmap(x1,x2,y1,y2) CALL xchsiz(0.025*(y2-y1)*lblmag) CALL xchori(0.0) CALL xcharl(x1+(x2-x1)*0.03, y2-(y2-y1)*(0.03+0.035*ip), & ch(1:lchar)) CALL xlbsiz( ctrlbsiz*(y2-y1)*lblmag ) ELSE CALL xlboff END IF CALL xwindw(lower, upper, zmin, zmax) CALL xqmap(x1,x2,y1,y2) CALL xcurve(profil(kbgn,ip),height(kbgn,ip),kend-kbgn+1,0) CALL xwdwof END DO RETURN END SUBROUTINE vprofil ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE SPLTPARA ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE spltpara(inc,MIN,MAX,ovr,hlf,zro,col1,col2,pltvar),4 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! Set some parameters for one plot. ! !----------------------------------------------------------------------- ! ! AUTHOR: ! Min Zou (3/2/98) ! ! Modification history: ! !----------------------------------------------------------------------- ! ! INPUT: ! inc interval of the contour ! min the minimum value for the contour ! max the maximum valur foj the contour ! ovr overlay option ! hlf the contour highlight frequency ! zro define the attributes of zero contours ! col1 the start color index for contour ! col2 the end color index for contour ! pltvar the plot name ! len the length of pltvar ! !----------------------------------------------------------------------- ! INTEGER :: ovr,hlf,zro,col1,col2 REAL :: inc, MIN, MAX CHARACTER (LEN=12) :: pltvar !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! CALL ctrinc( inc, MIN, MAX ) CALL overlay(ovr) CALL xhlfrq (hlf) CALL xczero (zro) CALL ctrcol( col1,col2) CALL varplt( pltvar ) RETURN END SUBROUTINE spltpara SUBROUTINE fillmissval (m, n, xl, xr, yb,yt ) 1 REAL :: x1,x2,y1,y2 REAL :: xra(4), yra(4) INTEGER :: missval_colind, missfill_opt ! miss value color index COMMON /multi_value/ missfill_opt,missval_colind x1 = xl + (xr-xl)/REAL(m)*0.5 x2 = xr - (xr-xl)/REAL(m)*0.5 y1 = yb + (yt-yb)/REAL(n)*0.5 y2 = yt - (yt-yb)/REAL(n)*0.5 xra(1) = x1 xra(2) = x2 xra(3) = x2 xra(4) = x1 yra(1) = y1 yra(2) = y1 yra(3) = y2 yra(4) = y2 CALL xcolor(missval_colind) CALL xfilarea(xra, yra, 4) RETURN END SUBROUTINE fillmissval ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE HINTRP ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE hintrp(nx,ny,nz,a3din,z3d,zlevel, a2dout) 1 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Interpolate a 3-D array to horizontal level z=zlevel. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! Based on original SECTHRZ. ! 12/10/98. ! ! MODIFICATION HISTORY: ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! ! INPUT : ! ! nx Number of grid points in the x-direction (east/west) ! ny Number of grid points in the y-direction (north/south) ! nz Number of grid points in the vertical ! ! a3din 3-d input array ! z3d z-coordinate of data in a3din ! zlevel Level to which data is interpolated. ! ! OUTPUT: ! a2dout 2-d output array interpolated to zlevel ! !----------------------------------------------------------------------- ! ! Parameters of output ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: nx,ny,nz REAL :: a3din(nx,ny,nz) ! 3-d input array REAL :: z3d (nx,ny,nz) ! z-coordinate of data in a3din REAL :: zlevel ! Level to which data is interpolated. REAL :: a2dout(nx,ny) ! 2-d output array interpolated to zlevel INTEGER :: i,j,k ! !----------------------------------------------------------------------- ! ! Find index for interpolation ! !----------------------------------------------------------------------- ! DO i=1,nx-1 DO j=1,ny-1 IF(zlevel <= z3d(i,j,1)) GO TO 11 IF(zlevel >= z3d(i,j,nz-1)) GO TO 12 DO k=2,nz-2 IF(zlevel >= z3d(i,j,k).AND.zlevel < z3d(i,j,k+1)) GO TO 15 END DO 11 k=1 GO TO 15 12 k=nz-1 GO TO 15 15 a2dout(i,j)=a3din(i,j,k)+(a3din(i,j,k+1)-a3din(i,j,k))* & (zlevel-z3d(i,j,k))/(z3d(i,j,k+1)-z3d(i,j,k)) !----------------------------------------------------------------------- ! ! If the data point is below the ground level, set the ! data value to the missing value. ! !----------------------------------------------------------------------- IF( zlevel < z3d(i,j,2) ) a2dout(i,j) = -9999.0 IF( zlevel > z3d(i,j,nz-1)) a2dout(i,j) = -9999.0 END DO END DO RETURN END SUBROUTINE hintrp ! ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE HINTRP1 ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE hintrp1(nx,ny,nz, kbgn,kend,a3din,z3d,zlevel, a2dout) 8 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Interpolate a 3-D array to horizontal level z=zlevel. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! Based on original SECTHRZ. ! 12/10/98. ! ! MODIFICATION HISTORY: ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! ! INPUT : ! ! nx Number of grid points in the x-direction (east/west) ! ny Number of grid points in the y-direction (north/south) ! nz Number of grid points in the vertical ! kbgn ! kend ! ! a3din 3-d input array ! z3d z-coordinate of data in a3din ! zlevel Level to which data is interpolated. ! ! OUTPUT: ! a2dout 2-d output array interpolated to zlevel ! !----------------------------------------------------------------------- ! ! Parameters of output ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: nx,ny,nz INTEGER :: kbgn, kend REAL :: a3din(nx,ny,nz) ! 3-d input array REAL :: z3d (nx,ny,nz) ! z-coordinate of data in a3din REAL :: zlevel ! Level to which data is interpolated. REAL :: a2dout(nx,ny) ! 2-d output array interpolated to zlevel INTEGER :: i,j,k ! !----------------------------------------------------------------------- ! ! Find index for interpolation ! !----------------------------------------------------------------------- ! DO i=1,nx-1 DO j=1,ny-1 IF(zlevel <= z3d(i,j,kbgn)) GO TO 11 IF(zlevel >= z3d(i,j,kend)) GO TO 12 DO k=kbgn,kend-1 IF(zlevel >= z3d(i,j,k).AND.zlevel < z3d(i,j,k+1)) GO TO 15 END DO 11 k=kbgn GO TO 15 12 k=kend-1 GO TO 15 15 a2dout(i,j)=a3din(i,j,k)+(a3din(i,j,k+1)-a3din(i,j,k))* & (zlevel-z3d(i,j,k))/(z3d(i,j,k+1)-z3d(i,j,k)) !----------------------------------------------------------------------- ! ! If the data point is below the ground level, set the ! data value to the missing value. ! !----------------------------------------------------------------------- IF( zlevel < z3d(i,j,kbgn) ) a2dout(i,j) = -9999.0 IF( zlevel > z3d(i,j,kend) ) a2dout(i,j) = -9999.0 END DO END DO RETURN END SUBROUTINE hintrp1 SUBROUTINE indxbnds(xc,yc,zpc,nx,ny,nz, & 1 xbgn,xend,ybgn,yend,zbgn,zend, & ibgn,iend,jbgn,jend,kbgn,kend) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! Return index bounds of the domain to be plotted ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: nx,ny,nz REAL :: xc (nx,ny,nz) ! x-coor of sacalar point (km) REAL :: yc (nx,ny,nz) ! y-coor of sacalar point (km) REAL :: zpc (nx,ny,nz) ! z-coor of sacalar point in physical ! space (km) REAL :: xbgn,xend,ybgn,yend,zbgn,zend INTEGER :: ibgn,iend,jbgn,jend,kbgn,kend INTEGER :: i,j,k IF(xbgn /= xend) THEN DO i = 1,nx-1 IF(xc(i,2,2) >= xbgn) EXIT END DO ibgn=MAX(i,1) DO i = nx-1,1,-1 IF(xc(i,2,2) <= xend ) EXIT END DO iend=MIN(i,nx-1) ELSE ibgn = 2 iend = nx- 2 xbgn = xc(ibgn,2,2) xend = xc(iend,2,2) END IF IF(ybgn /= yend) THEN DO j = 1,ny-1 IF(yc(2,j,2) >= ybgn) EXIT END DO jbgn=MAX(j,1) DO j = ny-1,1,-1 IF(yc(2,j,2) <= yend) EXIT END DO jend=MIN(j,ny-1) ELSE jbgn = 2 jend = ny-2 ybgn = yc(2,jbgn,2) yend = yc(2,jend,2) END IF IF(zbgn /= zend) THEN kend = 2 DO k = 2,nz-1 DO j = jbgn,jend DO i = ibgn,iend IF(zpc(i,j,k) < zend) THEN kend=k GO TO 225 END IF END DO END DO GO TO 235 225 CONTINUE END DO 235 kend = MIN(kend, nz-1) kbgn= nz-1 DO k = nz-1,2,-1 DO j = jbgn,jend DO i = ibgn,iend IF(zpc(i,j,k) > zbgn) THEN kbgn=k GO TO 250 END IF END DO END DO GO TO 245 250 CONTINUE END DO 245 kbgn = MAX(kbgn,2) ELSE kbgn = 2 kend = nz-2 END IF WRITE(6,'(/1x,a,i3,a,i3)') 'ibgn =',ibgn,', iend =',iend IF(iend < ibgn) THEN WRITE(6,'(1x,a,/1x,a)') & 'iend was found smaller than ibgn. Check the input', & 'domain bounds in x direction. Program stopped.' STOP END IF WRITE(6,'(1x,a,i3,a,i3)') 'jbgn =',jbgn,', jend =',jend IF(jend < jbgn) THEN WRITE(6,'(1x,a,/1x,a)') & 'jend was found smaller than jbgn. Check the input', & 'domain bounds in y direction. Program stopped.' STOP END IF WRITE(6,'(1x,a,i3,a,i3)') 'kbgn =',kbgn,', kend =',kend IF(kend < kbgn) THEN WRITE(6,'(1x,a,/1x,a)') & 'kend was found smaller than kbgn. Check the input', & 'domain bounds in z direction. Program stopped.' STOP END IF RETURN END SUBROUTINE indxbnds SUBROUTINE ctrsetup(zinc,zminc,zmaxc, & 84,4 zovr,zhlf,zzro,zcol1,zcol2,zlabel) IMPLICIT NONE REAL :: zinc,zminc,zmaxc INTEGER :: zovr,zhlf,zzro,zcol1,zcol2 CHARACTER (LEN=*) :: zlabel CALL ctrinc (zinc,zminc,zmaxc ) CALL overlay(zovr) CALL xhlfrq (zhlf) CALL xczero (zzro) CALL ctrcol (zcol1,zcol2 ) CALL varplt (zlabel) RETURN END SUBROUTINE ctrsetup ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE PLTTRN ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE plttrn(hterain,x,y,m,n,slicopt) 2,1 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Generate terrain contours ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! !----------------------------------------------------------------------- ! ! INPUT: ! ! hterain 2-D terrain data to contour ! x x coordinate of grid points in plot space (over on page) ! y y coordinate of grid points in plot space (up on page) ! m first dimension ! n second dimension ! ! slicopt slice orientation indicator ! slicopt = 1, x-y slice of u,v at z index kslice is plotted. ! slicopt = 2, x-z slice of u,w at y index jslice is plotted. ! slicopt = 3, y-z slice of v,w at x index islice is plotted. ! slicopt = 4, x-y slice of u,v at z index islice is plotted. ! slicopt = 5, xy-z cross section of wind islice is plotted. ! slicopt = 6, data field on constant p-level is plotted. ! slicopt = 0, all of the three slices above are plotted. ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: m,n REAL :: hterain(m,n) REAL :: x(m,n) REAL :: y(m,n) INTEGER :: slicopt ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! REAL :: ctinc,ctmin,ctmax,vtunt ! contour interval and vector unit COMMON /incunt/ ctinc,ctmin,ctmax,vtunt ! INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! REAL :: ztmin,ztmax INTEGER :: ovrtrn ,trnplt ! overlay terrain option (0/1) REAL :: trninc,trnmin, trnmax ! terrain interval minimum, maximum COMMON /trnpar/ trnplt,ovrtrn,trninc,trnmin, trnmax,ztmin,ztmax REAL :: zlevel COMMON /sliceh/zlevel INTEGER :: col_table,pcolbar COMMON /coltable/col_table,pcolbar ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- REAL :: cl(500) INTEGER :: ncl REAL :: z02,xl,xr,yt,yb,xfinc INTEGER :: mode1 INTEGER, ALLOCATABLE :: iwrk(:) REAL , ALLOCATABLE :: xwk(:),ywk(:) INTEGER :: istatus ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ALLOCATE(iwrk(m*n),STAT=istatus) ALLOCATE(xwk(m*n),ywk(m*n),STAT=istatus) IF( ovrtrn == 0 .OR. ztmax-ztmin < 1.0E-20 ) RETURN ! !----------------------------------------------------------------------- ! ! Overlay terrain contour if required in x-y level ! or Plot terrain outline in slice zlevel ! !----------------------------------------------------------------------- ! CALL xqmap(xl,xr,yb,yt) cl(1)=0.0 IF(slicopt == 1.OR.slicopt == 8) THEN CALL ctrinc( trninc,trnmin, trnmax ) IF( trninc == 0.0) THEN cl(2)=cl(1)+ xfinc(ztmax-ztmin)/2 IF(cl(2)-cl(1) == 0.0) cl(2)=cl(1)+1.0 mode1=1 CALL xnctrs(6,18) ELSE cl(2)=cl(1)+trninc CALL xnctrs(1,300) IF(ztmin == 0.0 .AND. ztmax == 0.0) THEN mode1=1 ELSE mode1=3 END IF END IF CALL xctrlim(ctmin,ctmax) IF (trnplt == 1) THEN CALL xthick(2) CALL xctrclr(trcolor, trcolor) IF(mode1 == 3) THEN ncl=(ztmax-ztmin)/trninc+1 cl(1)=ztmin cl(2)=cl(1)+trninc END IF CALL xconta(hterain,x,y,iwrk,m,m,n,cl,ncl,mode1) ELSE IF (trnplt == 2) THEN CALL xctrclr(icolor, icolor1) CALL xcolfil(hterain,x,y,iwrk,xwk,ywk,m,m,n,cl,ncl,mode1) CALL xchsiz(0.025*(yt-yb)) CALL xcpalet(pcolbar) ELSE IF (trnplt == 4) THEN CALL xctrclr(icolor, icolor1) CALL xconta(hterain,x,y,iwrk,m,m,n,cl,ncl,mode1) END IF ELSE IF(slicopt == 4.OR.slicopt == 6.OR.slicopt == 7) THEN CALL xcolor(trcolor) z02=zlevel*1000. CALL xthick(2) CALL xcontr(hterain,x,y,iwrk,m,m,n,z02) CALL xthick(1) END IF DEALLOCATE(iwrk) DEALLOCATE(xwk,ywk) RETURN END SUBROUTINE plttrn !################################################################## !################################################################## !###### ###### !###### SUBROUTINE PLTAXES ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE pltaxes(slicopt,dx,dy) 3,3 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! !----------------------------------------------------------------------- ! ! AUTHOR: M. Xue ! ! MODIFICATION HISTORY: ! !----------------------------------------------------------------------- ! ! ! INPUT: ! slicopt slice orientation indicator ! = 1, x-y slice of at k=kslice is plotted. ! = 2, x-z slice of at j=jslice is plotted. ! = 3, y-z slice of at i=islice is plotted. ! = 4, horizontal slice at z index islice is plotted. ! = 5, xy-z cross section of wind islice is plotted. ! = 6, data field on constant p-level is plotted. ! = 0, all of the three slices above are plotted. ! dx Spacing between the x-axis tick marks ! dy Spacing between the y-axis tick marks ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE REAL :: dx,dy INTEGER :: slicopt ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: layover COMMON /laypar/ layover INTEGER :: presaxis_no REAL :: pres_val(20), pres_z(20) COMMON /pressbar_par/presaxis_no,pres_val,pres_z REAL :: lblmag, ctrlbsiz, axlbsiz COMMON /labmag/ lblmag, ctrlbsiz, axlbsiz INTEGER :: timeovr COMMON /timover/ timeovr REAL :: x101, y101, x102,y102 COMMON /slicev1/x101, y101, x102,y102 INTEGER :: xfont ! the font of character INTEGER :: haxisu, vaxisu INTEGER :: lbaxis INTEGER :: tickopt INTEGER :: xfmat REAL :: hmintick,vmajtick,vmintick,hmajtick COMMON /var_par/ xfont,haxisu,vaxisu,lbaxis,tickopt,hmintick, & vmajtick, vmintick,hmajtick,xfmat INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! REAL :: axnmag REAL :: xl,xr,yb,yt,pl,pr,pb,pt REAL :: xtem1, xtem2 !local temporary variable REAL :: x1,x2, y1,y2, xstep, ystep, xmstep, ymstep CHARACTER (LEN=15) :: ylabel CHARACTER (LEN=15) :: xlabel INTEGER :: LEN ! !----------------------------------------------------------------------- ! ! Set-up variables for tick marks and draw axes ! !----------------------------------------------------------------------- ! CALL xqmap(xl,xr, yb,yt) CALL setcords(xl,xr,yb,yt,dx,dy, slicopt, & x1,x2,y1,y2,xlabel,ylabel,xstep,ystep,xmstep,ymstep) CALL xqpspc( pl, pr, pb, pt) axnmag = axlbsiz*MIN(pt-pb, pr-pl)*lblmag CALL xaxnmg( axnmag ) IF(slicopt == 5) THEN IF( ABS(y101-y102) <= 1.0E-3 ) THEN xtem1 = x101 xtem2 = x102 ELSE IF(ABS(x101-x102) <= 1.0E-3 ) THEN xtem1 = y101 xtem2 = y102 ELSE xtem1 = SQRT(x101*x101 + y101*y101) xtem2 = xtem1 + SQRT( (y102-y101)*(y102-y101) + & (x102-x101)*(x102-x101) ) END IF ELSE xtem1 = x1 xtem2 = x2 END IF CALL xmap(xtem1,xtem2,y1,y2) IF( layover == 0) THEN CALL xaxsor(0.0, 0.0) ! call xthick(2) CALL xaxsca1( xtem1,xtem2,xstep,xmstep, y1,y2,ystep,ymstep ) CALL xthick(1) END IF ! ! Plot pressure axis ! CALL xqpspc(pl, pr, pb, pt) IF(presaxis_no > 0 .AND.timeovr == 0 .AND. & (slicopt == 2 .OR. slicopt == 3 .OR. slicopt == 5) ) THEN x1 = pl - (pr-pl)*0.25 x2 = pl y1 = pb y2 = pt CALL xpspac(x1,x2,y1,y2) y1 = yb y2 = yt CALL xmap(x1,x2,y1,y2) CALL xaxfmt('(I4)') CALL xyaxis(x1+0.40*(x2-x1),pres_z,pres_val,presaxis_no) CALL xchori(90.) CALL xcharc(x1,yb+(yt-yb)*0.5 ,'Pressure(mb)') CALL xchori(0.) END IF ! ! Restore the original plotting scape ! CALL xpspac( pl, pr, pb, pt) CALL xmap(xl,xr, yb,yt) IF(layover > 1) THEN CALL xchsiz( 0.018*(yt-yb) * lblmag ) ELSE CALL xchsiz( 0.020*(yt-yb) * lblmag ) END IF IF(lbaxis == 1 .AND. timeovr == 0) THEN CALL xcolor(lbcolor) LEN=15 CALL strmin(xlabel,LEN) CALL xcharc( xl+(xr-xl)*0.5,yb-0.08*(yt-yb),xlabel(1:LEN)) LEN=15 CALL strmin(ylabel,LEN) CALL xchori(90.) CALL xcharc(xl-0.10*(xr-xl),yb+(yt-yb)*0.5,ylabel(1:LEN)) CALL xchori(0.) END IF RETURN END SUBROUTINE pltaxes ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE PLTEXTRA ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! SUBROUTINE pltextra(slicopt, pltopt) 3,2 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Plot extra things such as map, boxes, polygons and stations ! in a 2D plot ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! ! MODIFICATION HISTORY: ! ! 6/08/92 (K. Brewster) ! Added full documentation. ! ! 8/28/94 (M. Zou) ! Added color routing , overlay terrain. ! ! 1/24/96 (J. Zong and M. Xue) ! Fixed a problem related to finding the minimum and maximum of the ! 2D array, a, when there exist missing data. Initial min. and max. ! should be set to values other than the missing value, -9999.0. ! !----------------------------------------------------------------------- ! ! INPUT: ! slicopt slice orientation indicator ! = 1, x-y slice of at k=kslice is plotted. ! = 2, x-z slice of at j=jslice is plotted. ! = 3, y-z slice of at i=islice is plotted. ! = 4, horizontal slice at z index islice is plotted. ! = 5, xy-z cross section of wind islice is plotted. ! = 6, data field on constant p-level is plotted. ! = 0, all of the three slices above are plotted. ! plot ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: slicopt, pltopt INCLUDE 'arpsplt.inc' ! !----------------------------------------------------------------------- ! ! Plotting control common blocks ! !----------------------------------------------------------------------- ! INTEGER :: ovrstaopt INTEGER :: ovrstam,staset,ovrstan,ovrstav,stacol,markprio,wrtstax INTEGER :: nsta_typ,sta_typ(10),sta_marktyp(10),sta_markcol(10) REAL :: sta_marksz(10),wrtstad COMMON /sta_par/ ovrstaopt,ovrstam,staset,ovrstan,ovrstav,stacol, & markprio,nsta_typ,sta_typ,sta_marktyp, & sta_markcol,sta_marksz,wrtstax,wrtstad INTEGER :: icolor,icolor1,lbcolor,trcolor ! required color COMMON /recolor/icolor,icolor1,lbcolor,trcolor ! !----------------------------------------------------------------------- ! ! Misc. local variables ! !----------------------------------------------------------------------- ! INTEGER :: timeovr COMMON /timover/ timeovr INTEGER :: ovrmap COMMON /mappar / ovrmap INTEGER :: number_of_boxes,boxcol REAL :: bx1(10),bx2(10),by1(10),by2(10) COMMON /boxesopt/number_of_boxes,boxcol,bx1,bx2,by1,by2 INTEGER :: num_of_verts INTEGER :: number_of_polys,polycol REAL :: vertx(max_verts,max_polys),verty(max_verts,max_polys) COMMON /polysopt/number_of_polys,polycol,vertx,verty REAL :: lblmag, ctrlbsiz, axlbsiz COMMON /labmag/ lblmag, ctrlbsiz, axlbsiz INTEGER :: nunit INTEGER :: i,j ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! !----------------------------------------------------------------------- ! ! Plot map boundaries. ! !----------------------------------------------------------------------- ! CALL xcolor(lbcolor) IF((slicopt == 1.OR.slicopt == 4.OR.slicopt == 6 & .OR.slicopt == 7.OR.slicopt == 8) & .AND.ovrmap == 1 & .AND.(timeovr == 0 .OR. (timeovr == 1 .AND. pltopt == 2) ))THEN CALL getunit(nunit) CALL drawmap(nunit) CLOSE (UNIT=nunit) CALL retunit(nunit) CALL xthick(1) END IF ! !----------------------------------------------------------------------- ! ! Draw boxes ! !----------------------------------------------------------------------- ! IF(number_of_boxes /= 0 .AND. (slicopt == 1.OR.slicopt == 4.OR.slicopt == 6 & .OR.slicopt == 7.OR.slicopt == 8) & .AND. timeovr == 0 ) THEN CALL xthick(1) CALL xcolor(boxcol) CALL xbrokn(6,3,6,3) DO i=1,number_of_boxes CALL xbox(bx1(i),bx2(i),by1(i),by2(i)) END DO CALL xthick(1) CALL xfull END IF ! !----------------------------------------------------------------------- ! ! Draw polylines ! !----------------------------------------------------------------------- ! IF(number_of_polys /= 0 .AND. (slicopt == 1.OR.slicopt == 4.OR.slicopt == 6 & .OR.slicopt == 7.OR.slicopt == 8) & .AND. timeovr == 0 ) THEN CALL xthick(2) CALL xcolor(polycol) ! CALL xbrokn(6,3,6,3) DO j=1,number_of_polys num_of_verts=0 DO i=1,max_verts IF(vertx(i,j) /= -9999. .AND. verty(i,j) /= -9999.) & num_of_verts = num_of_verts +1 END DO IF(num_of_verts /= 0 ) CALL xcurve(vertx(1,j),verty(1,j),num_of_verts, 0) END DO CALL xthick(1) CALL xfull END IF RETURN END SUBROUTINE pltextra ! !################################################################## !################################################################## !###### ###### !###### SUBROUTINE SMOOTH9PMV ###### !###### ###### !###### Developed by ###### !###### Center for Analysis and Prediction of Storms ###### !###### University of Oklahoma ###### !###### ###### !################################################################## !################################################################## ! ! SUBROUTINE smooth9pmv( arr, nx,ny,ibgn,iend,jbgn,jend, tem1 ) 28 ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! 1 2 1 ! Smooth a 2-D array by the filter of { 2 4 2 } ! 1 2 1 ! !----------------------------------------------------------------------- ! ! AUTHOR: Yuhe Liu ! ! 5/3/94 ! ! Modification History ! 8/20/1995 (M. Xue) ! Fixed errors in the index bound of loops 100 and 200. ! !----------------------------------------------------------------------- ! ! INPUT: ! ! nx Number of grid points in the x-direction ! ny Number of grid points in the y-direction ! ibgn First index in x-direction in the soomthing region. ! iend Last index in x-direction in the soomthing region. ! jbgn First index in j-direction in the soomthing region. ! jend Last index in j-direction in the soomthing region. ! ! arr 2-D array ! ! OUTPUT: ! ! arr 2-D array ! ! TEMPORARY: ! ! tem1 Temporary 2-D array ! !----------------------------------------------------------------------- ! ! Variable Declarations. ! !----------------------------------------------------------------------- ! IMPLICIT NONE ! INTEGER :: nx ! Number of grid points in the x-direction INTEGER :: ny ! Number of grid points in the y-direction INTEGER :: ibgn INTEGER :: iend INTEGER :: jbgn INTEGER :: jend ! REAL :: arr (nx,ny) ! 2-D array ! REAL :: tem1(nx,ny) ! Temporary array ! !----------------------------------------------------------------------- ! ! Misc. local variables: ! !----------------------------------------------------------------------- ! INTEGER :: i,j,ip,im,jp,jm REAL :: wtf,mv ! !----------------------------------------------------------------------- ! ! Include files: ! !----------------------------------------------------------------------- ! ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! wtf = 1.0/16.0 mv = -9999.0 ! missing value flag DO i=1,nx DO j=1,ny IF( ABS(arr(i,j)-mv) <= 0.1) arr(i,j)=mv END DO END DO DO j = jbgn,jend DO i = ibgn,iend ip=MIN(nx,i+1) im=MAX( 1,i-1) jp=MIN(ny,j+1) jm=MAX( 1,j-1) tem1(i,j) = wtf & * ( arr(im,jm) + 2.*arr(i,jm) + arr(ip,jm) & + 2.*arr(im,j ) + 4.*arr(i,j ) + 2.*arr(ip,j ) & + arr(im,jp) + 2.*arr(i,jp) + arr(ip,jp)) IF(arr(im,jm) == mv.OR.arr(i,jm) == mv.OR.arr(ip,jm) == mv.OR. & arr(im,j ) == mv.OR.arr(i,j ) == mv.OR.arr(ip,j ) == mv.OR. & arr(im,jp) == mv.OR.arr(i,jp) == mv.OR.arr(ip,jp) == mv)THEN tem1(i,j)=mv END IF END DO END DO DO j = jbgn,jend DO i = ibgn,iend arr(i,j) = tem1(i,j) END DO END DO RETURN END SUBROUTINE smooth9pmv SUBROUTINE buoycy_plt(nx,ny,nz,ptprt,pprt,qv,qc,qr,qi,qs,qh, & 1 ptbar,pbar,rhobar,qvbar, wbuoy, tem1) ! !----------------------------------------------------------------------- ! ! PURPOSE: ! ! Calculate the total buoyancy including liquid and solid water ! loading. ! !----------------------------------------------------------------------- ! ! AUTHOR: Ming Xue ! 10/10/91. ! ! MODIFICATION HISTORY: ! ! 5/05/92 (M. Xue) ! Added full documentation. ! ! 3/10/93 (M. Xue) ! The buoyancy term is reformulated. The previous formula was ! in error. The water loading was calculated wrong, resulting in ! a value of the water loading that is typically an order of ! magnitude too small. ! ! 3/25/94 (G. Bassett & M. Xue) ! The buoyancy terms are reformulated for better numerical accuracy. ! Instead of storing numbers which had the form (1+eps)*(1+eps1) ! (eps << 1 and eps1 <<1), terms were expanded out, and most of the ! high order terms neglected, except for the second order terms ! in ptprt, pprt and qvbar. ! ! 9/10/94 (D. Weber & Y. Lu) ! Cleaned up documentation. ! ! 6/21/95 (Alan Shapiro) ! Fixed bug involving missing qvpert term in buoyancy formulation. ! ! 10/15/97 (Donghai Wang) ! Added a new option for including the second order terms. ! ! 11/05/97 (D. Weber) ! Changed lower loop bounds in DO LOOP 400 for computing the ! buoyancy term from k=3,nz-2 to k=2,nz-1. Level k=2 data will be ! used in the hydrostatic pprt lower boundary condition (removed ! DO LOOP 410 used to set wbuoy = 0.0 for k= 2 and nz-1). ! !----------------------------------------------------------------------- ! ! INPUT : ! ! nx Number of grid points in the x-direction (east/west) ! ny Number of grid points in the y-direction (north/south) ! nz Number of grid points in the vertical direction. ! ! ptprt Perturbation potential temperature at a time level (K) ! pprt Perturbation pressure at a given time level (Pascal) ! qv Water vapor specific humidity at a given time level ! (kg/kg) ! qc Cloud water mixing ratio at a given time level (kg/kg) ! qr Rainwater mixing ratio at a given time level (kg/kg) ! qi Cloud ice mixing ratio at a given time level (kg/kg) ! qs Snow mixing ratio at a given time level (kg/kg) ! qh Hail mixing ratio at a given time level (kg/kg) ! ! ptbar Base state potential temperature (K) ! pbar Base state pressure (Pascal) ! rhobar Base state density rhobar ! qvbar Base state water vapor specific humidity (kg/kg) ! ! OUTPUT: ! ! wbuoy The total buoyancy force (kg/(m*s)**2) ! ! WORK ARRAYS: ! ! tem1 Temporary work array. ! !----------------------------------------------------------------------- ! ! Variable Declarations ! !----------------------------------------------------------------------- ! IMPLICIT NONE INTEGER :: nx,ny,nz ! Number of grid points in 3 directions REAL :: ptprt (nx,ny,nz) ! Perturbation potential temperature ! at a given time level (K) REAL :: pprt (nx,ny,nz) ! Perturbation pressure at a given time ! level (Pascal) REAL :: qv (nx,ny,nz) ! Water vapor specific humidity (kg/kg) REAL :: qc (nx,ny,nz) ! Cloud water mixing ratio (kg/kg) REAL :: qr (nx,ny,nz) ! Rain water mixing ratio (kg/kg) REAL :: qi (nx,ny,nz) ! Cloud ice mixing ratio (kg/kg) REAL :: qs (nx,ny,nz) ! Snow mixing ratio (kg/kg) REAL :: qh (nx,ny,nz) ! Hail mixing ratio (kg/kg) REAL :: ptbar (nx,ny,nz) ! Base state potential temperature (K) REAL :: pbar (nx,ny,nz) ! Base state pressure (Pascal). REAL :: rhobar(nx,ny,nz) ! Base state density rhobar REAL :: qvbar (nx,ny,nz) ! Base state water vapor specific ! humidity(kg/kg) REAL :: wbuoy(nx,ny,nz) ! Total buoyancy in w-eq. (kg/(m*s)**2) REAL :: tem1 (nx,ny,nz) ! Temporary work array. ! !----------------------------------------------------------------------- ! ! Misc. local variables: ! !----------------------------------------------------------------------- ! INTEGER :: i,j,k REAL :: g5 REAL :: pttem,tema ! !----------------------------------------------------------------------- ! ! Include files: ! !----------------------------------------------------------------------- ! INCLUDE 'globcst.inc' ! Global model control parameters INCLUDE 'phycst.inc' ! Physical constants ! !----------------------------------------------------------------------- ! ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Beginning of executable code... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! !----------------------------------------------------------------------- ! ! The total buoyancy ! ! wbuoy = rhobar*g ( ptprt/ptbar-pprt/(rhobar*csndsq)+ ! qvprt/(rddrv+qvbar)-(qvprt+qc+qr+qs+qi+qh)/(1+qvbar) ! -(ptprt*ptprt)/(ptbar*ptbar) !2nd-order ! +0.5*(ptprt*pprt)/(cpdcv*ptbar*pbar)) !2nd-order ! ! and rddrv=rd/rv, cp, cv, rd and rv are defined in phycst.inc. ! ! Here, the contribution from pprt (i.e., term pprt/(rhobar*csndsq)) ! is evaluated inside the small time steps, therefore wbuoy ! does not include this part. ! ! The contribution from ptprt is calculated inside the small time ! steps if the potential temperature equation is solved inside ! small time steps, i.e., if ptsmlstp=1. ! !----------------------------------------------------------------------- ! tema = 1.0/cpdcv DO k=1,nz-1 DO j=1,ny-1 DO i=1,nx-1 pttem = ptprt(i,j,k)/ptbar(i,j,k) tem1(i,j,k) = pttem* & (1.0-pttem+0.5*pprt(i,j,k)*(tema/pbar(i,j,k))) END DO END DO END DO ! !----------------------------------------------------------------------- ! ! Add on the contributions to the buoyancy from the water vapor ! content and the liquid and ice water loading. ! !----------------------------------------------------------------------- ! DO k=1,nz-1 DO j=1,ny-1 DO i=1,nx-1 tem1(i,j,k) = tem1(i,j,k) & + (qv(i,j,k) - qvbar(i,j,k))/(rddrv + qvbar(i,j,k)) & - (qv(i,j,k) - qvbar(i,j,k) + qc(i,j,k) + qr(i,j,k) + & qs(i,j,k) + qi(i,j,k) + qh(i,j,k))/(1 + qvbar(i,j,k)) END DO END DO END DO ! !----------------------------------------------------------------------- ! ! Then the total buoyancy: ! ! wbuoy = tem1 * rhobar * g ! ! averged to the w-point on the staggered grid. ! !----------------------------------------------------------------------- ! DO k=1,nz-1 DO j=1,ny-1 DO i=1,nx-1 wbuoy(i,j,k)= tem1(i,j, k )*rhobar(i,j, k )*g END DO END DO END DO RETURN END SUBROUTINE buoycy_plt