! ! Richard Carpenter ! Univ. of Oklahoma ! August 1992 ! ! ! ######################################## ! ######################################## ! ######################################## ! ######## ######## ! ######## SKEWT ######## ! ######## ######## ! ######################################## ! ######################################## ! ######################################## PROGRAM skewt,102 !----------------------------------------------------------------------- ! ! HISTORY: ! 09/13/93 Modified somewhat for parallel ARPSTools version. ! 12/08/93 Added option for plotting Tv with or without water ldng. ! 1995/01/24 -tvnowl now also plots Tv of environment. ! 1995/12/07 In Derive added code not to bomb when Td = 0. ! 1997/03/06 Ported to throttle (Sun) from rossby (Sun). ! Change istnm (INTEGER) to stnm (CHAR). ! 1997/03/21 Hodograph now in upper left corner of main plot. ! Using color table input file (hardwired). ! 1997/04/22 Heights are now MSL, not AGL. ! 1997/04/23 Add a 2nd color scheme - all lines of each snd are 1 color. ! 1997/05/22 Improve handling of multiple soundings. ! 1999/10/11 Merged in COMET-Tinker version ("skewti"). Last COMET mod ! was 1999/01/31. ! Changed default tv_wl to False. [RLC] ! Converted to Fortran-90 using TO_F90 by Alan Miller ! 1999/10/14 Numerous changes, merged version checks out. ! Color file may be specified by -colorfile or ! setenv SKEWT_COLORFILE. [RLC ! ! Skewed temperatures are indicated by _sk. ! ! Theta-E lines calculated using Bolton (1980), eq. 43. The slightly less ! accurate method, eq. 35, dIFfers by no more that 0.2 K. Overall, they seem ! to by within 0.5 K of the DOD Skew-T chart. ! ! UNITS: All quantities are SI units unless labeled otherwise (e.g., ! presmb, tempc). Mixing ratios are kg/kg. ! ! m/s = mph / 2.237 = kts / 1.944 ! mph = kts / 1.151 ! !----------------------------------------------------------------------= IMPLICIT NONE INCLUDE 'thermo.consts' LOGICAL, PARAMETER :: thetalines=.true., thetaelines=.true., & mixinglines=.true., do_legend=.true. INTEGER, PARAMETER :: nmax=1000, npmx=200, ntmx=40, nticmx=80, nsmax=8 REAL, PARAMETER :: km2hft=32.8084 !----------------------------------------------------------------------= LOGICAL :: & plot_sfc_parcel, plot_cb_parcel, plot_frame_2, plot_tv, & has_wind , orig , gempak_fmt, & weightq , cape_use_t , cape_use_irrev, write_modsnd , & has_pres , has_zz , arpstools = .false., & do_hodo = .false., hodo_denmwind = .false., helcontrs = .false., & verbose = .false., modify_snd = .false., tv_wl = .false., & wrt_indcs = .true., logfiles = .false., do_indices = .false., & plot_wind = .false., print_info = .false. INTEGER :: i , j , k , n , & ifcb , nxtic , nptic , npres , & isnd , nztic , nc=0 , kc , & np , klcl , jgray ,lcolors(nsmax), jhunits , & kp , irev , ltypes(nsmax), kk , nargs , jdefcolor, & ksc , n_td_sk , kmod , nsnd=0 , kmix , nitems , & iyr , imon , iday , ihr , imin , ifont , & jtcolor , jtdcolor , jparcolor, jtvcolor , jhodoclr , & jcolor(nsmax), jscheme,ii,time(6,4) INTEGER :: preciptype_me, preciptype_laps(nmax),meprectypeskewt REAL :: tw REAL :: & tmin=0.0 , tmax=0.0 , tinc=10.0, pmin=0.0 , pref=1050.0 , & skewfac , tskew , dtskew , t1 , t2 , p1 , & p2 , thetae , pinc=100.0,zcb , & dz , stnelev , px , py , pres , es , & plogmin , plogmax , theta , dp , tc ,qqsgkg(14), & pmbcb , tccb , esmb , preslcl , templcl , tlcl , & tvlcl , zlcl , deg2rad , pi , xx , yy , & psfc , px1 , px2 , py1 , py2 , temp , & siz0 , qs , & zmin , zmax , qmin , qmax , siz08 , pyoffset , & yykm1 , dlnp , sfctf , sfctdf ,rad2deg , pmax , & pspacing , pydelta , tvavg , y1 , y2 , & sfctf0 , sfctdf0 , pmb_mod , pxoffset , zlnb , preslnb , & capetvwl , cintvwl ,capetvnowl, cintvnowl, capet , cint , & capetvwlc, cintvwlc,capetvnowlc,cintvnowlc, capetc , cintc , & pxhodo1 , pxhodo2 , pyhodo1 , pyhodo2 , x1 , x2 REAL :: presmb_s(nmax), tempc_s(nmax), tdewc_s(nmax) , & temp_c(nmax) , preslog(nmax) , pp(npmx) ,th_sk(npmx,ntmx), & xtic(nticmx) , ptic(nticmx) , pticlab(nticmx),pplog(npmx) , & the_sk(npmx,ntmx),tq(npmx,ntmx) , alwcpr(nmax) , alwcpa(nmax), & pres_c(nmax) , zzkm(nmax) , & plog_c(nmax) , ztic(nticmx) , zticlab(nticmx),zz_s(nmax) , & rh_s(nmax) , & theta_s(nmax) , alwc_cr(nmax) , tempp(nmax) , plogp(nmax) , & tvp(nmax) , tvwlp(nmax) , thete(nmax) , thetq(nmax) , & tv_c(nmax) , tvwl_c(nmax) , & tempc_sk(nmax), tdewc_sk(nmax), tvc_sk(nmax) , tvcp_sk(nmax) , & tcp_sk(nmax) ,tvcwlp_sk(nmax), tcc_sk(nmax) ,tvcwlc_sk(nmax), & tvcnowlc_sk(nmax), & pres_s(nmax) , temp_s(nmax) , tdew_s(nmax) , qv_s(nmax) , & qs_s(nmax) , tv_s(nmax) , qvgkg_s(nmax) , & uu_s(nmax) , vv_s(nmax) , dir_s(nmax) , spd_s(nmax) , & alwcgkgpa(nmax),alwcgkgpr(nmax), & tempc_sk_orig(nmax), tdewc_sk_orig(nmax), presmb_s_orig(nmax), & preslog_orig(nmax), wbt(nmax) ! Array wbt(nmax) added for precip type algorithms. 1/18/98 EMK REAL :: pres500,ht500,temp500,tdewp500,spd500,dir500,pres700,ht700, & temp700,tdewp700,spd700,dir700,pres850,ht850,temp850,tdewp850, & spd850,dir850,showalter,kindex,liftedindex,totaltotals,sweat, & brn,brnshear,wetbulbzero,convtemp,zlfc,preslfc,mixrat, & mixrattot,prestopbl,xmxrat,precipwat,pccl,slat,slon, & lidstrength,dirmean,spmean,helicity REAL :: indx(6,27),lat(6),lon(6) CHARACTER (LEN=255) :: sndfile(nsmax), strings(nsmax),file1,file2, indoutfile CHARACTER (LEN=3) :: stid CHARACTER (LEN=6) :: stnm CHARACTER (LEN=2) :: chr,cmin,cmoi,cdayi,chri,cmini CHARACTER (LEN=3) :: cstni CHARACTER (LEN=132) :: string, string2 CHARACTER (LEN=100) :: capestring, color_file CHARACTER (LEN=19) :: nm(33) DATA nm / & 'ht LCL (m AGL) ', 'ht LFC (m AGL) ', 'ht LNB (m AGL) ',& 'pres LCL (mb) ', 'pres LFC (mb) ', 'pres LNB (mb) ',& 'CAPE-Tv-w/loading ', 'CAPE-Tv-w/o loading', 'CAPE-using T only ',& 'CIN-Tv-w/loading ', 'CIN-Tv-w/o loading ', 'CIN-using T only ',& 'Showalter Index ', 'K-Index ', 'Lifted Index ',& 'Total-totals Index ', 'SWEAT Index ', 'Bulk Rich Shear ',& 'Precip. Water (cm) ', 'Wet-bulb zero (m) ', 'Convective Temp ',& 'BL mr (sfc- 50mb) ', 'pres CCL (mb) ', 'Lidstrength Index ',& 'Storm Direction ', 'Storm Speed (m/s) ', 'SR Helicity ',& 'Place-Time: ', 'Date: ', 'Latitude: ',& 'Longitude: ', 'MEta Sfc Prcp Type ', 'LAPS Sfc Precip '/ CHARACTER (LEN=14) :: precindex_me(6),precindex_laps(6) CHARACTER (LEN=3) :: place(6) INTEGER :: tstart, tstop, tstep INTEGER :: istatus ! jtdcolor was 12, jt was 20, jpar was 22 DATA jdefcolor/01/, jgray/04/, jtcolor/39/, jtdcolor/33/, & jparcolor/22/, jtvcolor/22/, jhodoclr/21/ DATA jcolor /30, 01, 20, 33, 22, 21, 03, 05/ !best for jscheme=2 !DATA jcolor /01, 30, 20, 33, 22, 21, 03, 05/ !best for jscheme=2 !DATA jcolor /01, 02, 30, 33, 22, 21, 03, 05/ !best for jscheme=1 ! jscheme: 1 = first snd is multi-color ! jscheme: 2 = all snd are single color. Do not plot sfc parcel for snd>2 ! jscheme: 3 = all snd are single color. Plot sfc parcel etc for snd>2 DATA jscheme /2/ !color scheme -- 2 for single color per snd ! Update the following for each version DATA color_file /'../data/arpsplt/skewt.pltcbar'/ DATA plot_sfc_parcel, plot_cb_parcel, plot_frame_2, plot_tv & /.false., .false., .false., .false./ DATA weightq, cape_use_t, cape_use_irrev /.false., .false., .TRUE./ DATA nargs /0/, siz0/0.02/, ifont /2/ DATA pspacing/50.0/ ! normalized vert spacing of lvls; ! larger number -> closer spacing DATA qqsgkg /0.05, 0.1, 0.2, 0.5, 1.0, 2.0, 5.0, & 10.0, 15.0, 20.0, 25.0, 30.0, 35.0, 40.0/ !Mixing ratio ref lines INCLUDE 'thermo.stfunc' !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ ! ! Begin of executable code ... ! !@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ !----------------------------------------------------------------------= ! Set up !----------------------------------------------------------------------= pi = 4. * ATAN (1.) deg2rad = pi/180. rad2deg = 1. / deg2rad ! Get input ! PRINT *, 'Call ParseInput' CALL parseinput ( & nsmax , nsnd , sndfile , has_wind , & print_info , verbose , logfiles , do_indices , & color_file , & do_hodo , hodo_denmwind , jhunits , helcontrs , & plot_cb_parcel,plot_sfc_parcel, plot_tv , & modify_snd , weightq , write_modsnd , & tv_wl , cape_use_t , cape_use_irrev, & tmin , tmax , pmin , pmax ) ! Confirm options PRINT '(/,a)', ' === Confirm user options ===' PRINT *, 'Done ParseInput, color file is ',TRIM(color_file) IF(plot_tv) PRINT *, 'Plotting Tv' IF(plot_frame_2) PRINT *, 'Plotting frame 2' IF(plot_cb_parcel) PRINT *, 'Plotting parcel lifted from cloud base' IF(plot_sfc_parcel) PRINT *, 'Plotting parcel lifted from surface' IF(weightq) PRINT *, 'Weighting mixing ratio when modifying sounding' IF(cape_use_t) PRINT *, 'Using T instead of Tv when computing CAPE' IF(cape_use_irrev) PRINT *, 'Using pseudoadiabat' ! string describing CAPE capestring = 'Parcel CAPE computed by integ. ' capestring = 'CAPE: ' string = capestring IF (cape_use_t) THEN ! string2 = ' T (no virt. correction)' string2 = ' T' ELSE IF (tv_wl) THEN ! string2 = ' Tv (incl. condensate loading)' string2 = ' Tv*' ELSE ! string2 = ' Tv (no condensate loading)' string2 = ' Tv' END IF capestring = trim(string) // trim(string2) string = capestring IF (cape_use_irrev) THEN ! string2 = ' psuedoadiabat.' string2 = ' psuedo' ELSE ! string2 = ' reversible adiabat.' string2 = ' rev' END IF capestring = trim(string) // trim(string2) PRINT *, TRIM(capestring) PRINT *, '=== End of user options ===' PRINT *, ' ' IF (nsnd == 1) jscheme = 1 DO i=1,6,1 place(i)='---' END DO !----------------------------------------------------------------------= ! LOOP THROUGH SOUNDINGS !----------------------------------------------------------------------= DO isnd=1,nsnd !...Read the sounding ! ReadSound returns P, Z, Theta, Qv, U, V. ! P or Z is computed hydrostatically, IF needed. PRINT *, 'Calling ReadSound, file = ', trim(sndfile(isnd)) CALL readsound (pres_s, zz_s, theta_s, qv_s, uu_s, vv_s, & has_pres, has_zz, pmbcb, tccb, dz, stnelev, psfc, & sndfile(isnd), nmax, n, ifcb, plot_wind, gempak_fmt, has_wind, & stid,stnm,slat,slon, iyr,imon,iday,ihr,imin,istatus) IF (istatus /= 0) CYCLE !...Compute other quantities from P,Z,Th,Qv,U,V. CALL derive (n, pres_s, zz_s, theta_s, qv_s, uu_s, vv_s, & temp_s, tempc_s, tdew_s, tdewc_s, presmb_s, spd_s, dir_s, qs_s, rh_s, tv_s) ! Graphics stuff -- first sounding only IF (isnd == 1) THEN ! Figure out coords based on pres at top of sounding IF (tmin == 0. .AND. tmax == 0.) THEN tmin = -20. tmax = 40. pmin = 300. pmax = pref IF (presmb_s(n) < pmin) THEN pmin = pmin - pinc END IF IF (presmb_s(n) < pmin) THEN pmin = pmin - pinc tmin = tmin - tinc END IF IF (pmin == 300.) tmin = -10. IF (presmb_s(1) > pmax) pmax = pref + 50. END IF plogmin = LOG(pmin) plogmax = LOG(pmax) dtskew = tmax - tmin IF (plot_wind) dtskew = dtskew + tinc IF (pmin >= 200.) dtskew = dtskew - tinc IF (pmin >= 300.) dtskew = dtskew - tinc skewfac = dtskew/(plogmax - plogmin) END IF !isnd !----------------------------------------------------------------------= ! ANALYSIS !----------------------------------------------------------------------= PRINT *, 'Writing sounding info to 1.skewt, 2.skewt' ! Save original values; truncate sounding above 100 mb DO k=1,n presmb_s_orig(k) = presmb_s(k) tempc_sk_orig(k) = tempc_s(k) tdewc_sk_orig(k) = tdewc_s(k) ! print *, k,pmin,presmb_S(k) IF (presmb_s(k) <= pmin) THEN n = k PRINT *, 'Truncating sounding to n,p(n): ', n, presmb_s(n) GO TO 8100 END IF END DO 8100 CONTINUE ! :::::::::::::::::::::::::: ! :::::::: SOUNDING :::::::: ! :::::::::::::::::::::::::: ! Note: both alwcPA and alwcPR are no longer obtained (01/28/93) irev = 1 IF (cape_use_irrev) irev = 0 PRINT *, 'CALLing AnalyzSnd2, irev = ', irev CALL analyzsnd2 (n,pres_s,temp_s,tdew_s, irev, & zz_s,qv_s,qs_s,rh_s,tv_s,theta_s, preslcl,templcl,zlcl,preslnb,zlnb, & alwcpr,tempp,tvp,tvwlp, & cintvwl,capetvwl, cintvnowl,capetvnowl, cint,capet, 0) DO k=1,n qvgkg_s(k)= qv_s(k) * 1000. END DO !----------------------------------------------------------------------= ! Modify sounding (first sounding only) !----------------------------------------------------------------------= IF (isnd == 1 .AND. modify_snd) THEN PRINT * CALL sfcupd (presmb_s,zz_s,tempc_s,tdewc_s,qvgkg_s,theta_s, & uu_s,vv_s, sfctf,sfctdf,orig,n,nmax, kmod,pmb_mod,sfctf0,sfctdf0,weightq) PRINT *, 'Done modifying sounding. n,kmod: ', n, kmod PRINT '(1x,a,4f8.2)', 'sfctf,sfctdf,sfctf0,sfctdf0: ', & sfctf,sfctdf,sfctf0,sfctdf0 PRINT * ! ...debrief the modified values (obtain p,Z,Th,Qv,U,V) DO k=1,n pres_s(k) = presmb_s(k) * 100. qv_s(k) = qvgkg_s(k) * 1.e-3 END DO CALL derive (n, pres_s, zz_s, theta_s, qv_s, uu_s, vv_s, & temp_s, tempc_s, tdew_s, tdewc_s, presmb_s, spd_s, dir_s, & qs_s, rh_s, tv_s) ! Compute heights hydrostatiCALLy ! zz_s(1) = 0. DO k=2,kmod tvavg = .5 * (tv_s(k-1) + tv_s(k)) dz = rd/grav*tvavg * LOG(pres_s(k-1)/pres_s(k)) zz_s(k) = zz_s(k-1) + dz END DO PRINT *, 'Calling AnalyzSnd2 after modifying sounding. irev=',irev CALL analyzsnd2(n,pres_s,temp_s,tdew_s, irev, & zz_s,qv_s,qs_s,rh_s,tv_s,theta_s, & preslcl,templcl,zlcl,preslnb,zlnb, & alwcpr,tempp,tvp,tvwlp, & cintvwl,capetvwl,cintvnowl,capetvnowl,cint,capet,0) END IF ! modify_snd !----------------------------------------------------------------------= tvlcl = ftvirtnowl(templcl, fmixrat(preslcl,fsvpres(templcl-tfrz))) DO k=1,n zzkm(k) = zz_s(k) * 1.e-3 END DO !----------------------------------------------------------------------= ! WRITE SOUNDING INFO TO FILE !----------------------------------------------------------------------= PRINT * PRINT *, 'SOUNDING PARAMETERS' PRINT * IF (logfiles) THEN WRITE (file1,'(A,I1,A)') 'a-', isnd, '.skewt' PRINT *, 'Writing to ', trim(file1) OPEN (UNIT=1, FILE=file1) WRITE (1,'(2a)') '# sounding file: ', sndfile(isnd) WRITE (1,9903) WRITE (1,9909) END IF PRINT * 9904 FORMAT (4(a,f8.2)) PRINT 9904, 'CAPE: ', capetvwl, ' CIN: ', cintvwl, & ' Max updraft:', SQRT(2.*capetvwl) IF (cape_use_t) PRINT *, 'Values using T, not Tv*:' PRINT 9904, 'CAPE: ', capet, ' CIN: ', cint, & ' Max updraft:', SQRT(2.*capet) IF (verbose) THEN PRINT * PRINT 9903 PRINT 9909 END IF 9903 FORMAT ('# Pres Z T[C] Td[C] Qv[g/kg] Qs ', & ' RH Tv[C] Theta U V') 9909 FORMAT ('#',77('-')) 9901 FORMAT (f7.2,f7.0,12F7.2) 9902 FORMAT (4F8.2,2X,8F8.2) DO k=1,n IF (logfiles) WRITE (1,9901) & pres_s(k)*1.e-2,zz_s(k),temp_s(k)-tfrz,tdew_s(k)-tfrz, & qv_s(k)*1000.,qs_s(k)*1000.,rh_s(k),tv_s(k)-tfrz,theta_s(k), & uu_s(k), vv_s(k) IF (verbose) THEN PRINT * PRINT 9901, pres_s(k)*1.e-2,zz_s(k),temp_s(k)-tfrz, tdew_s(k)-tfrz, & qv_s(k)*1000.,qs_s(k)*1000.,rh_s(k),tv_s(k)-tfrz,theta_s(k), & uu_s(k), vv_s(k) END IF END DO IF (logfiles) CLOSE (1) ! WRITE to ARPS input format file IF (logfiles) THEN WRITE (file1,'(A,I1,A)') 'a-', isnd, '.skewt' CALL writearpssnd (file1, stnelev, psfc, n, & pres_s, zz_s, theta_s, temp_s, qv_s, rh_s, tdew_s, uu_s, vv_s) END IF 9950 FORMAT (i3,2X,f8.2,f7.0,2X,2F8.2,2X,2F8.2) IF (plot_wind .AND. verbose) THEN PRINT * DO k=1,n PRINT 9950, k,pres_s(k)*1.e-2,zz_s(k), dir_s(k),spd_s(k),uu_s(k),vv_s(k) END DO PRINT * END IF ! :::::::::::::::::::::::::::::::::::::::::::: ! :::::::: PARCEL LIFTED FROM SURFACE :::::::: ! :::::::::::::::::::::::::::::::::::::::::::: PRINT * PRINT *, 'PARAMETERS FOR PARCEL LIFTED FROM SFC ' PRINT * IF (logfiles) THEN WRITE (file1,'(A,I1,A)') 'b-', isnd, '.skewt' PRINT *, 'Writing ', trim(file1) OPEN (UNIT=1, FILE=file1) WRITE (1,'(2a)') '# sounding file: ', sndfile(isnd) WRITE (1,9910) preslcl*1.e-2,templcl-tfrz, zlcl WRITE (1,9911) WRITE (1,9909) END IF PRINT * PRINT 9910, preslcl*1.e-2,templcl-tfrz, zlcl IF (verbose) THEN PRINT * PRINT 9911 PRINT 9909 END IF 9910 FORMAT ('# Parcel lifted from surface: LCL P,t,z: ',2F7.2,f7.0) 9911 FORMAT ('# Pres Z ALWC[g/kg] TP[C] dTvP[C] dTv*P[C]', & ' dT') DO k=1,n IF (pres_s(k) <= preslcl) THEN IF (verbose) THEN PRINT * PRINT 9902, presmb_s(k),zz_s(k),alwcpr(k)*1000., tempp(k)-tfrz, & tvp(k)-tv_s(k),tvwlp(k)-tv_s(k), tempp(k)-temp_s(k) END IF IF (logfiles) WRITE (1,9902) & presmb_s(k),zz_s(k),alwcpr(k)*1000.,tempp(k)-tfrz, & tvp(k)-tv_s(k),tvwlp(k)-tv_s(k), tempp(k)-temp_s(k) END IF END DO IF (logfiles) CLOSE (1) ! ::::::::::::::::::::::::::::::::::::::::::::::::::::: ! :::::::: PARCEL LIFTED FROM KNOWN CLOUD BASE :::::::: ! ::::::::::::::::::::::::::::::::::::::::::::::::::::: IF (ifcb == 1) THEN PRINT * PRINT *, 'PARAMETERS FOR PARCEL LIFTED FROM CLOUD BASE ' PRINT * IF (logfiles) THEN WRITE ( file1,'(A,I1,A)') 'c-', isnd, '.skewt' WRITE (file2,'(A,I1,A)') 'd-', isnd, '.skewt' PRINT *, 'Writing ', trim(file1), ', ', trim(file2) OPEN (UNIT=1, FILE='3.skewt') OPEN (UNIT=2, FILE='4.skewt') WRITE (1,'(2a)') '# sounding file: ', sndfile(isnd) WRITE (2,'(2a)') '# sounding file: ', sndfile(isnd) WRITE (1,9911) DO k=1,n WRITE (2,9902) pres_s(k)/100., zz_s(k),1.00 END DO END IF ! point kC = 1 (cloud base) kc = 1 pres_c(kc) = pmbcb * 100. temp_c(kc) = tccb + tfrz qs = fmixrat(pres_c(kc),fsvpres(tccb)) tv_c(kc) = ftvirt(temp_c(kc),qs,0.) tvwl_c(kc)= ftvirt(temp_c(kc),qs,0.) capetvwlc = 0. capetvnowlc = 0. capetc = 0. ! points kC = 2:nC ! point kSC is the first sounding point below the cloud base. ! note that k=kSC will not be used. DO k=1,n IF (presmb_s(k) < pmbcb) THEN kc = kc + 1 pres_c(kc)= pres_s(k) END IF END DO nc = kc ksc = n - nc + 1 PRINT *, 'nC,kSC: ', nc,ksc PRINT 9902, pres_c(1)/100., zcb,alwc_cr(1)*1000.,temp_c(1)-tfrz PRINT *, 'CALLing AnalyzSnd2 (cloud base), irev = ', irev CALL analyzsnd2 (nc,pres_s(ksc),temp_s(ksc),tdew_s(ksc), irev, & zz_s(ksc),qv_s(ksc),qs_s(ksc),rh_s(ksc),tv_s(ksc),theta_s(ksc), & pres_c,temp_c,zcb,xx,xx, alwc_cr,temp_c,tv_c,tvwl_c, & cintvwlc,capetvwlc, cintvnowlc,capetvnowlc, cint,capet, 1) IF (logfiles) THEN WRITE (1,9922) pmbcb, tccb, zcb WRITE (1,9923) WRITE (1,9909) END IF PRINT * PRINT 9922, pmbcb, tccb, zcb IF (verbose) THEN PRINT * PRINT 9923 PRINT 9909 END IF 9922 FORMAT ('# Cloud base p,T,z: ',2F7.2,f9.0) 9923 FORMAT ('# Pres Z ALWC[g/kg] TP[C]') PRINT 9904, 'CAPE: ', capetvwlc, ' CIN: ', cintvwlc kc = 1 IF (verbose) THEN PRINT * PRINT 9902, pres_c(kc)/100., -99.,alwc_cr(kc)*1000.,temp_c(kc)-tfrz, & tv_c(kc),tvwl_c(kc), temp_c(kc) END IF IF (logfiles) WRITE (1,9902) pres_c(kc)/100., & -99.,alwc_cr(kc)*1000.,temp_c(kc)-tfrz, tv_c(kc),tvwl_c(kc), temp_c(kc) DO kc=2,nc k = kc + ksc - 1 !k=kC will be bogus. IF (verbose) THEN PRINT 9902, pres_c(kc)/100., & zz_s(k),alwc_cr(kc)*1000.,temp_c(kc)-tfrz, & tv_c(kc)-tv_s(k),tvwl_c(kc)-tv_s(k), temp_c(kc)-temp_s(k) END IF IF (logfiles) WRITE (1,9902) pres_c(kc)/100., & zz_s(k),alwc_cr(kc)*1000.,temp_c(kc)-tfrz, & tv_c(kc)-tv_s(k),tvwl_c(kc)-tv_s(k), temp_c(kc)-temp_s(k) IF (logfiles) WRITE (2,9902) pres_c(kc)/100., & zz_s(k),alwc_cr(kc)*1000. END DO IF (logfiles) THEN CLOSE (1) CLOSE (2) END IF END IF ! Write the modified sounding in model format. ! Remove point k=kmod, which had been added to the sounding. IF (isnd == 1 .AND. write_modsnd .AND. modify_snd) THEN file1 = trim(sndfile(isnd)) // '.mod' PRINT *, 'Writing modified sounding to: ', file1 OPEN (UNIT=2, FILE=file1) WRITE (2,*) 'Modified sounding:', trim(file1) WRITE (2,8801) (n-1)/2, 2.*dz, stnelev, psfc, ifcb, pmbcb, tccb WRITE (2,8802) (theta_s(k),k=1,kmod-1), (theta_s(k),k=kmod+1,n) WRITE (2,8802) (qv_s(k),k=1,kmod-1), (qv_s(k),k=kmod+1,n) WRITE (2,8802) (uu_s(k),k=1,kmod-1), (uu_s(k),k=kmod+1,n) WRITE (2,8802) (vv_s(k),k=1,kmod-1), (vv_s(k),k=kmod+1,n) 8801 FORMAT (//,i4, 3F9.2, i8,2F7.2) 8802 FORMAT (/(1P,5G15.6)) CLOSE (UNIT=2) END IF !----------------------------------------------------------------------= IF (isnd == 1) THEN ! Height axis kk = 1 yykm1 = 10000. !dlnp = (LOG(pref)-LOG(pmin))/70. !spacing between labels !dlnp = (LOG(pref)-LOG(pmin))/30. !RLC 1995/12/07 dlnp = (LOG(pref)-LOG(pmin))/pspacing DO k=1,n yy = LOG (presmb_s(k)) IF (yykm1-yy > dlnp) THEN ztic(kk)= yy zticlab(kk)= zz_s(k) * 1.0E-3 IF (jhunits == 1) zticlab(kk)= zticlab(kk) * km2hft yykm1 = yy kk = kk + 1 END IF END DO nztic = kk - 1 END IF !isnd !----------------------------------------------------------------------= ! SKEW THE TEMPERATURES !----------------------------------------------------------------------= ! Sounding n_td_sk = n DO k=1,n preslog(k)= LOG (presmb_s(k)) preslog_orig(k) = LOG (presmb_s_orig(k)) tskew = skewfac * (plogmax - preslog(k)) tvc_sk(k) = ftvirtnowl(temp_s(k),qv_s(k)) + tskew - tfrz tempc_sk(k) = tempc_s(k) + tskew tdewc_sk(k) = tdewc_s(k) + tskew IF (n_td_sk == n .AND. tdewc_s(k) < -90.) n_td_sk = k - 1 tskew = skewfac * (plogmax - preslog_orig(k)) tempc_sk_orig(k) = tempc_sk_orig(k) + tskew tdewc_sk_orig(k) = tdewc_sk_orig(k) + tskew END DO ! Parcel ! Insert an extra point at the LCL np = n + 1 ! ...below LCL klcl = 0 DO kp=1,np k = kp IF (klcl == 0 .AND. preslcl > pres_s(k)) THEN klcl = k GO TO 8800 END IF IF (klcl > 0) k = kp - 1 plogp(kp) = preslog(k) tskew = skewfac * (plogmax - plogp(kp)) tvcwlp_sk(kp) = tvwlp(k) - tfrz + tskew tvcp_sk(kp) = tvp(k) - tfrz + tskew tcp_sk(kp)= tempp(k) - tfrz + tskew 8800 CONTINUE END DO IF (klcl == 0) klcl = np ! EMK (1/31/98) ! ...at LCL kp = klcl PRINT *, 'kLCL', klcl plogp(kp) = LOG (preslcl*1.e-2) tskew = skewfac * (plogmax - plogp(kp)) tvcwlp_sk(kp)= tvlcl - tfrz + tskew tvcp_sk(kp)= tvlcl - tfrz + tskew tcp_sk(kp)= templcl - tfrz + tskew ! Known Cloud Base ! the option tvnowl, CB, is not available RLC 1994/03/04 DO kc=1,nc plog_c(kc) = LOG (pres_c(kc)*1.e-2) tskew = skewfac * (plogmax - plog_c(kc)) tvcwlc_sk(kc) = tvwl_c(kc) - tfrz + tskew !tvcnowlC_sk(kC) = tvnowl_C(kC) - tfrz + tskew tvcnowlc_sk(kc) = tvwl_c(kc) - tfrz + tskew tcc_sk(kc) = temp_c(kc) - tfrz + tskew END DO !----------------------------------------------------------------------- ! EXTENDED CALCULATION OF STABILITY PARAMETERS !----------------------------------------------------------------------- IF (do_indices) THEN !----------------------------------------------------------------------= ! INTERPOLATE THE VARIABLES TO 850, 700 and 500 mbs. JJM !----------------------------------------------------------------------= CALL interpress (n,pres_s,tempc_s,tdewc_s,spd_s,zz_s, & dir_s,pres850,ht850,temp850,tdewp850,spd850,dir850, & pres700,ht700,temp700,tdewp700,spd700,dir700,pres500, & ht500,temp500,tdewp500,spd500,dir500,nmax) !----------------------------------------------------------------------= ! COMPUTE AN AVERAGE MIXING RATIO IN THE LOWEST 50 mb FOR USE AS ! A MIXED SFC PARCEL !----------------------------------------------------------------------= mixrattot=0. prestopbl=0. ii=0 DO i=1,n,1 IF (pres_s(i) > (pres_s(1)-5000.0)) THEN mixrat=xmxrat((pres_s(i)/100.),(tdewc_s(i)+273.16)) mixrattot=mixrattot+mixrat ii=ii+1 END IF END DO IF (ii > 1) THEN mixrat=mixrattot/(ii) prestopbl=pres_s(ii) ELSE mixrat=mixrattot prestopbl = pres_s(1) ! EMK 1/31/98 ! prestopBL=pres_S(1)+1000. ! Original. Why set this to be below ! ! the lowest sounding level? END IF !----------------------------------------------------------------------= ! CALCULATE THE PERTINENT INDICES (JJM) !----------------------------------------------------------------------= CALL indices(pres850,ht850,temp850,tdewp850,spd850,dir850, & pres700,ht700,temp700,tdewp700,spd700,dir700,pres500,ht500, & temp500,tdewp500,spd500,dir500,showalter,kindex,liftedindex, & totaltotals,sweat,preslcl,zlcl,templcl,mixrat,prestopbl, & pres_s,tempc_s,zz_s,tdewc_s,n,nmax,convtemp,precipwat, & wetbulbzero,pccl,zlfc,preslfc,zlnb,uu_s,vv_s,brnshear, & lidstrength,dirmean,spmean,helicity) IF ( pres_s(1) < 85000. ) THEN ! cms showalter = 1000000000000000. kindex = 1000000000000000. totaltotals = 1000000000000000. sweat = 1000000000000000. END IF ! Determine the wet bulb temperature profile (1/18/98) EMK k = 1 DO WHILE (k <= n) wbt(k) = tw((tempc_s(k)),(tdewc_s(k)),(pres_s(k)/100.)) wbt(k) = wbt(k) + 273.15 ! print*,'k = ',k ! print*,'Wet Bulb Temperature (K) = ',wbt(k) k = k + 1 END DO ! Call precip. type algorithms (10/24/97, 1/18/98) EMK PRINT*, 'Calling meprectypeskewt...' preciptype_me = meprectypeskewt(nmax,n,pres_s,zz_s,tempc_s, tdewc_s,wbt) PRINT*,'Finished with meprectypeskewt' ! print*,'MesoEta Precip. Type = ',preciptype_me PRINT*, 'Calling prectype_laps_skewt...' CALL prectype_laps_skewt(nmax,n,tempc_s,pres_s,wbt, preciptype_laps) PRINT*,'Finished with prectype_laps_skewt' ! print*,'LAPS surface Precip. Type = ',preciptype_laps(1) !----------------------------------------------------------------------= ! FORMATTED PRINTOUTS (USED BY POST-PROCESSING PROGRAMS) !----------------------------------------------------------------------= indx(isnd,1)=zlcl 2001 FORMAT(a19,4X,f6.0,4X,f6.0,4X,f6.0,4X,f6.0,4X,f6.0,4X,f6.0) indx(isnd,2)=zlfc 2002 FORMAT(a19,4X,f6.1,4X,f6.1,4X,f6.1,4X,f6.1,4X,f6.1,4X,f6.1) indx(isnd,3)=zlnb indx(isnd,4)=preslcl/100.0 indx(isnd,5)=preslfc/100.0 indx(isnd,6)=preslnb/100.0 indx(isnd,7)=capetvwl 2003 FORMAT(a19,5X,f5.0,5X,f5.0,5X,f5.0,5X,f5.0,5X,f5.0,5X,f5.0) indx(isnd,8)=capetvnowl indx(isnd,9)=capet indx(isnd,10)=cintvwl indx(isnd,11)=cintvnowl indx(isnd,12)=cint indx(isnd,13)=showalter 2004 FORMAT(a19,5X,f5.1,5X,f5.1,5X,f5.1,5X,f5.1,5X,f5.1,5X,f5.1) indx(isnd,14)=kindex indx(isnd,15)=liftedindex indx(isnd,16)=totaltotals indx(isnd,17)=sweat indx(isnd,18)=brnshear indx(isnd,19)=precipwat indx(isnd,20)=wetbulbzero indx(isnd,21)=convtemp indx(isnd,22)=mixrat*1000 indx(isnd,23)=pccl indx(isnd,24)=lidstrength indx(isnd,25)=dirmean indx(isnd,26)=spmean indx(isnd,27)=helicity ! Adding precindex_me and precindex_laps for precip. type ! algorithms. (EMK 10/23/1997, 1/18/1998) IF (preciptype_me == 1) THEN precindex_me(isnd) = 'SN' ELSE IF (preciptype_me == 2) THEN precindex_me(isnd) = 'IP' ELSE IF (preciptype_me == 3) THEN precindex_me(isnd) = 'FZ RN' ELSE precindex_me(isnd) = 'RN' END IF IF (preciptype_laps(1) == 1) THEN precindex_laps(isnd) = 'SN' ELSE IF (preciptype_laps(1) == 2) THEN precindex_laps(isnd) = 'IP' ELSE IF (preciptype_laps(1) == 3) THEN precindex_laps(isnd) = 'FZ RN' ELSE precindex_laps(isnd) = 'RN' END IF place(isnd)=stid time(isnd,1)=iday time(isnd,2)=ihr time(isnd,3)=imin time(isnd,4)=imon lat(isnd)=slat lon(isnd)=slon 2005 FORMAT(a19,a4,1X,i2.2,':',i2.2,1X,a3,1X,i2.2,':',i2.2,1X,a3,1X, & i2.2,':',i2.2,1X,a3,1X,i2.2,':',i2.2,1X,a3,1X,i2.2,':',i2.2, & 1X,a3,1X,i2.2,':',i2.2) 2006 FORMAT(a19) 2007 FORMAT(a19,2X,i2.2,'/',i2.2,5X,i2.2,'/',i2.2,5X,i2.2,'/',i2.2, & 5X,i2.2,'/',i2.2,5X,i2.2,'/',i2.2,5X,i2.2,'/',i2.2) 2010 FORMAT(a19,3X,f7.2,3X,f7.2,3X,f7.2,3X,f7.2,3X,f7.2,3X,f7.2) 2011 FORMAT(a,8X,a5,5X,a5,5X,a5,5X,a5,5X,a5,5X,a5) IF (isnd == nsnd) THEN PRINT * PRINT * PRINT 2007, nm(29), (time(i,4),time(i,1), i=1,6) PRINT 2005, nm(28), (place(i),time(i,2),time(i,3), i=1,6) PRINT 2010, nm(30), (lat(i), i=1,6) PRINT 2010, nm(31), (lon(i), i=1,6) PRINT * DO j=1,3 PRINT 2001, nm(j), (indx(i,j), i=1,6) END DO DO j=4,6 PRINT 2002, nm(j), (indx(i,j), i=1,6) END DO DO j=7,9 PRINT 2003, nm(j), (indx(i,j), i=1,6) END DO DO j=10,12 PRINT 2002, nm(j), (indx(i,j), i=1,6) END DO DO j=13,16 PRINT 2004, nm(j), (indx(i,j), i=1,6) END DO PRINT 2002, nm(17), (indx(i,17), i=1,6) PRINT 2010, nm(18), (indx(i,18), i=1,6) PRINT 2004, nm(19), (indx(i,19), i=1,6) PRINT 2001, nm(20), (indx(i,20), i=1,6) PRINT 2004, nm(21), (indx(i,21), i=1,6) PRINT 2010, nm(22), (indx(i,22), i=1,6) PRINT 2002, nm(23), (indx(i,23), i=1,6) PRINT 2010, nm(24), (indx(i,24), i=1,6) PRINT 2003, nm(25), (indx(i,25), i=1,6) PRINT 2010, nm(26), (indx(i,26), i=1,6) PRINT 2010, nm(27), (indx(i,27), i=1,6) PRINT 2011, nm(32), (precindex_me(i), i=1,6) PRINT 2011, nm(33), (precindex_laps(i), i=1,6) END IF IF (wrt_indcs) THEN indoutfile=stid//'_ind'//'.txt' OPEN (UNIT=7,FILE=TRIM(indoutfile),STATUS='unknown') IF (nsnd == 1) THEN WRITE (7,2057) nm(29), time(1,4), time(1,1), iyr WRITE (7,2055) nm(28), place(1), time(1,2), time(1,3) WRITE (7,2058) nm(30), lat(1) WRITE (7,2060) nm(31), lon(1) WRITE (7,*) DO j=1,27 WRITE (7,2051) nm(j), indx(1,j) END DO WRITE (7,2061) nm(32), precindex_me(1) WRITE (7,2061) nm(33), precindex_laps(1) 2051 FORMAT(a19,4X,f6.0) 2052 FORMAT(a19,4X,f6.1) 2053 FORMAT(a19,5X,f5.0) 2054 FORMAT(a19,5X,f5.1) 2055 FORMAT(a19,a4,1X,i2.2,':',i2.2,'Z') 2057 FORMAT(a19,2X,i2.2,'/',i2.2,'/',i2.2) 2058 FORMAT(a19,5X,f5.2) 2059 FORMAT(a19,4X,f6.2) 2060 FORMAT(a19,3X,f7.2) 2061 FORMAT(a19,8X,a5) ELSE IF (isnd == nsnd) THEN WRITE (7,2007) nm(29), (time(i,4),time(i,1), i=1,6) WRITE (7,2005) nm(28), (place(i),time(i,2),time(i,3), i=1,6) WRITE (7,2010) nm(30), (lat(i), i=1,6) WRITE (7,2010) nm(31), (lon(i), i=1,6) WRITE (7,*) DO j=1,3 WRITE (7,2001) nm(j), (indx(i,j), i=1,6) END DO DO j=4,6 WRITE (7,2002) nm(j), (indx(i,j), i=1,6) END DO DO j=7,9 WRITE (7,2003) nm(j), (indx(i,j), i=1,6) END DO DO j=10,12 WRITE (7,2002) nm(j), (indx(i,j), i=1,6) END DO DO j=13,16 WRITE (7,2004) nm(j), (indx(i,j), i=1,6) END DO WRITE (7,2002) nm(17), (indx(i,17), i=1,6) WRITE (7,2010) nm(18), (indx(i,18), i=1,6) WRITE (7,2004) nm(19), (indx(i,19), i=1,6) WRITE (7,2001) nm(20), (indx(i,20), i=1,6) WRITE (7,2004) nm(21), (indx(i,21), i=1,6) WRITE (7,2010) nm(22), (indx(i,22), i=1,6) WRITE (7,2002) nm(23), (indx(i,23), i=1,6) WRITE (7,2010) nm(24), (indx(i,24), i=1,6) WRITE (7,2003) nm(25), (indx(i,25), i=1,6) WRITE (7,2010) nm(26), (indx(i,26), i=1,6) WRITE (7,2010) nm(27), (indx(i,27), i=1,6) WRITE (7,2011) nm(32), (precindex_me(i), i=1,6) WRITE (7,2011) nm(33), (precindex_laps(i), i=1,6) !----------------------------------------------------------------------= ! END FORMATTED PRINTOUTS (JJM) !----------------------------------------------------------------------= CLOSE(7) END IF END IF END IF ! do_indices !----------------------------------------------------------------------= ! PLOTTING !----------------------------------------------------------------------= PRINT *, 'Begin plotting' IF (isnd == 1) THEN ! X-tics nxtic = nint ((tmax-tmin)/tinc) + 1 DO i=1,nxtic xtic(i) = (i-1)*tinc + tmin END DO ! P-tics nptic = INT((pmax-pmin)/pinc + 1.0) DO i=1,nptic pticlab(i)= pmin + (i-1)*pinc ptic(i) = LOG(pticlab(i)) END DO ! Z-tics px1 = 0.1 IF (do_indices) px1 = 0.065 ! EMK (1/18/98) px2 = 0.94 py1 = 0.10 py2 = 0.97 IF (plot_wind) THEN px2 = px2 - 0.12 IF (nsnd > 1) px2 = px2 - 0.05 END IF PRINT *, px1, px2, py1, py2 IF (jscheme > 1) THEN jtcolor = jcolor(isnd) jtvcolor = jcolor(isnd) jtdcolor = jcolor(isnd) jparcolor = jcolor(isnd) END IF CALL xdevic CALL xstctfn (color_file) CALL setcolors (-1) ! Use -1 for file, or >=1 for predefined ! reverse black/white ! CALL gscr (1, 0, 1.0, 1.0, 1.0) ! uncomment them for NCARG ! CALL gscr (1, 1, 0.0, 0.0, 0.0) IF (arpstools) ifont = 2 CALL xcfont (ifont) CALL color (jdefcolor) CALL xpspac (px1,px2,py1,py2) CALL xmap (tmin,tmax,plogmax,plogmin) ! Protected region for hodograph pxhodo1 = px1 pxhodo2 = pxhodo1 + 0.3 pyhodo2 = py2 pyhodo1 = pyhodo2 - 0.3 IF (do_hodo) THEN PRINT *, 'Protected region for hodograph' CALL plt2wrld (pxhodo1,pyhodo1,x1,y1) CALL plt2wrld (pxhodo2,pyhodo2,x2,y2) !print *, tmin,tmax,plogmax,plogmin !print *, x1,x2,y1,y2 CALL xmask (x1,x2,y2,y1) !must reverse ys; see CALL XMAP END IF CALL xchmag (0.9*siz0) CALL xaxnmg (0.9*siz0) CALL xaxtik (0,0) CALL xaxfmt ('(I3)') CALL xxaxis (xtic, xtic, nxtic, plogmax) ! CALL XAXISX (tmin, plogmax, tinc) CALL xaxfmt ('(I4)') CALL xyaxis (tmin, ptic, pticlab, nptic) !pres ! Height axis along right edge CALL xaxnmg (0.65*siz0) CALL xaxnmg (0.90*siz0) !RLC 1995/12/07 CALL xaxtik (1,-1) CALL xaxant (1,1) CALL xaxfmt ('(F4.1)') string = 'km' IF (jhunits == 1) THEN CALL xaxfmt ('(I3.3)') string = 'hft' END IF CALL xyaxis (tmax, ztic, zticlab, nztic) !hgt CALL wrld2plt (tmax,ztic(1),px,py) CALL plt2wrld (px+0.035,py-0.03,xx,yy) ! Original CALL xcharc (xx,yy,trim(string)) CALL xbordr pyoffset = 0.03 pxoffset = 0.01 pydelta = 0.02 ! Clip lines outside window CALL xwindw (tmin,tmax,plogmin,plogmax) ! P-lines CALL color (jgray) DO i=2,nptic CALL myline (tmin,ptic(i),tmax,ptic(i)) END DO ! T-lines #ifdef INTLOOP tstart = NINT(tmin-4*tinc) tstop = NINT(tmax-tinc) tstep = NINT(tinc) DO i=tstart,tstop,tstep #else DO i=tmin-4.*tinc,tmax-tinc,tinc #endif t1 = FLOAT(i) t2 = FLOAT(i) + dtskew p1 = plogmax p2 = plogmin ! IF (t1.LT.tmin) THEN ! p1 = p1 + (tmin-t1)/(t2-t1) * (p2 - p1) ! t1 = tmin ! END IF ! IF (t2.GT.tmax) THEN ! p2 = p1 + (tmax-t1)/(t2-t1) * (p2 - p1) ! t2 = tmax ! END IF CALL myline (t1,p1,t2,p2) END DO ! Theta-lines ! Note that th,the are the (skewed) values of T[C] corresp to theta,thetae. PRINT *, '... Plotting theta lines' ! CALL XBROKN (1,6,5,6) dp = 25. npres = nint ((pref-pmin)/dp) + 1 ! PRINT *, 'npres ', npres DO k=1,npres pp(k) = pmin + (k-1) * dp pplog(k) = LOG (pp(k)) END DO IF (thetalines) THEN DO i=1,nxtic+4 theta = tfrz + tmin + i*tinc DO k=npres,1,-1 tskew = skewfac * (plogmax - pplog(k)) tc = theta * (pp(k)/p00mb) ** rcp - tfrz th_sk(k,i) = tc + tskew END DO CALL xcurve (th_sk(1,i), pplog, npres, 0) END DO END IF ! ThetaE-lines ! Procedure: given theta(1000mb), set a reference Theta-E. THEN loop through ! other pres levels, changing T,Theta, until we find the T at which Th-E = the ! reference value. ! See Bolton (1980). IF (thetaelines) THEN CALL xbrokn (1,6,1,6) PRINT *, '... Plotting theta-e lines' ! DO i=1,2*nxtic-3 DO i=1,nxtic-1 ! lowest point (reference) k = npres pres = pp(k) * 100. ! theta = tfrz + tmin + i*tinc*.5 theta = tfrz + tmin + i*tinc temp = theta * (pres/p00) ** rcp tc = temp - tfrz es = fsvpres(tc) qs = fmixrat(pres,es) thetae = fthetae(pres,temp,temp,qs) ! PRINT *, 'Theta-E line:', i, theta-tfrz, thetae-tfrz tskew = skewfac * (plogmax - pplog(k)) the_sk(k,i)= tc + tskew ! other points (iterate) DO k=npres-1,1,-1 pres = pp(k) * 100. CALL the2t (temp, thetae, pres) tskew = skewfac * (plogmax - pplog(k)) the_sk(k,i) = temp - tfrz + tskew END DO CALL xcurve (the_sk(1,i), pplog, npres, 0) END DO !CALL Color (jdefcolor) END IF !----------------------------------------------------------------------- ! ! Mixing ratio lines ! !----------------------------------------------------------------------- IF (mixinglines) THEN CALL xchmag (0.6*siz0) !CALL Color (3) PRINT *, '... Plotting mixing ratio lines' !CALL XDASH CALL xbrokn (6,16,6,16) DO i=1,14 DO k=1,npres tskew = skewfac * (plogmax - pplog(k)) esmb = pp(k) * qqsgkg(i) / (1000.*cp622 + qqsgkg(i)) tc = 243.5 / (17.67/LOG(esmb/6.112) - 1.) tq(k,i) = tc + tskew IF ( pp(k) <= 490.0 ) kmix = k END DO ! kmix is the uppermost point to plot !kmix = MIN (kmix,npres/2) CALL xcurve (tq(kmix,i), pplog(kmix), npres-kmix+1, 0) IF (qqsgkg(i) >= 1.0) THEN WRITE (string,'(I2)') INT(qqsgkg(i)) ELSE IF (qqsgkg(i) >= 0.1) THEN WRITE (string,'(F3.1)') qqsgkg(i) ELSE WRITE (string,'(F4.2)') qqsgkg(i) END IF CALL xcharc( tq(kmix,i),pplog(kmix),trim(string) ) END DO END IF CALL xwdwof !Turn off window clipping IF (do_hodo) CALL xunmsk (1) ! Un-protect hodograph region END IF !isnd=1 ! Label the plot CALL color (jcolor(isnd)) CALL plot_strings (sndfile(isnd), stid, do_indices, & cape_use_t, gempak_fmt, modify_snd, plot_cb_parcel, & plot_sfc_parcel, print_info, tv_wl, & isnd , nsnd , iyr , imon , iday , ihr , & imin , px1 , px2 , py1 , py2 , & capet , capetc ,capetvnowl,capetvnowlc,capetvwl , capetvwlc, & cint , cintc , cintvnowl,cintvnowlc, cintvwl , cintvwlc , & pmb_mod , pmbcb , preslcl , preslnb , & sfctf , sfctf0 , sfctdf , sfctdf0 , zcb, zlcl, zlnb , & preslfc , zlfc , convtemp , precipwat , liftedindex , & totaltotals , kindex , sweat , wetbulbzero , brnshear , & lidstrength , dirmean , spmean , helicity , preciptype_me, & preciptype_laps,nmax,n,jcolor(isnd)) ! Linetypes: 0,1 full; 2 dash; 3 dot; 4 my dot; 5 my dash; 6 my dash-dot; ! negative bold. !----------------------------------------------------------------------= ! ! FINALLY, PLOT THE CURVES! ! !----------------------------------------------------------------------= CALL color (jcolor(isnd)) ! Hodograph IF (do_hodo) THEN IF (isnd <= 2) THEN CALL hodograph (pres_s, zz_s, theta_s, qv_s, uu_s, vv_s, & n,isnd, sndfile(isnd), helcontrs,hodo_denmwind, & jdefcolor, jgray, jcolor(isnd), & pxhodo1,pxhodo2,pyhodo1,pyhodo2) END IF CALL xpspac (px1,px2,py1,py2) CALL xmap (tmin,tmax,plogmax,plogmin) END IF !----------------------------------------------------------------------= CALL xwindw (tmin,tmax,plogmin,plogmax) !----------------------------------------------------------------------- ! ! Plot T, Td, Tparcel ! !----------------------------------------------------------------------- CALL color (jcolor(isnd)) CALL xthick(2) !thickness of lines, 1 or 2 IF (isnd > 2) CALL xthick(1) ! T CALL xfull IF (isnd == 1 .AND. jscheme == 1) CALL color (jtcolor) CALL xcurve (tempc_sk, preslog, n, 0) ! Td CALL xdash IF (isnd == 1 .AND. jscheme == 1) CALL color (jtdcolor) CALL xcurve (tdewc_sk, preslog, n_td_sk, 0) ! Tv IF (plot_tv .AND. (isnd == 1 .OR. jscheme == 3) ) THEN IF (isnd == 1 .AND. jscheme == 1) CALL color (jtvcolor) CALL xbrokn (10,5,1,5) CALL xcurve (tvc_sk, preslog, n, 0) END IF ! Original part of modified sounding IF (modify_snd .AND. (isnd == 1 .OR. jscheme == 3) ) THEN CALL xthick (1) CALL xfull IF (isnd == 1 .AND. jscheme == 1) CALL color (jtcolor) CALL xcurve (tempc_sk_orig, preslog_orig, kmod, 0) CALL xdash IF (isnd == 1 .AND. jscheme == 1) CALL color (jtdcolor) CALL xcurve (tdewc_sk_orig, preslog_orig, kmod, 0) IF (isnd <= 2) CALL xthick (2) END IF ! Parcel lifted from surface IF (plot_sfc_parcel .AND. (isnd == 1 .OR. jscheme == 3) ) THEN IF (isnd == 1 .AND. jscheme == 1) CALL color (jparcolor) CALL lintyp (03) IF (plot_tv) THEN IF (tv_wl) THEN CALL xcurve (tvcwlp_sk, plogp, np, 0) CALL xscatter (tvcwlp_sk(klcl),plogp(klcl),1,4, 0.005) ELSE CALL xcurve (tvcp_sk, plogp, np, 0) CALL xscatter (tvcp_sk(klcl),plogp(klcl),1,4, 0.005) END IF ELSE CALL xcurve (tcp_sk, plogp, np, 0) CALL xscatter (tcp_sk(klcl),plogp(klcl),1,4, 0.005) END IF END IF ! Parcel lifted from cloud base IF (ifcb == 1 .AND. plot_cb_parcel .AND. & (isnd == 1 .OR. jscheme == 3) ) THEN CALL xdot IF (isnd == 1 .AND. jscheme == 1) CALL color (jparcolor) IF (plot_tv) THEN IF (tv_wl) THEN CALL xcurve (tvcwlc_sk, plog_c, nc, 0) CALL xfull CALL xscatter (tvcwlc_sk(1),plog_c(1),1,4, 0.005) ELSE CALL xcurve (tvcnowlc_sk, plog_c, nc, 0) CALL xfull CALL xscatter (tvcnowlc_sk(1),plog_c(1),1,4, 0.005) END IF ELSE CALL xcurve (tcc_sk, plog_c, nc, 0) CALL xfull CALL xscatter (tcc_sk(1),plog_c(1),1,4, 0.005) END IF END IF CALL xwdwof !Turn off window clipping !...Plot wind CALL color (jcolor(isnd)) CALL xfull CALL xthick(1) IF (isnd <= 2 .AND. plot_wind) THEN PRINT *, '...Plotting wind barbs' string = ' ' IF (nsnd > 1) WRITE (string,'(I1)') isnd CALL plotwind (presmb_s, dir_s,spd_s, n, pspacing, & jscheme, string, do_indices) END IF !...Plot LAPS Precip (EMK, 1/18/98) IF (do_indices) THEN CALL color (jcolor(isnd)) CALL xfull CALL xthick(1) IF (isnd <= 2) THEN PRINT *, '...Plotting LAPS Precip. Type' string = ' ' IF (nsnd > 1) WRITE (string,'(I1)') isnd CALL plotprecip (presmb_s, n, pspacing, jscheme, string,preciptype_laps) END IF END IF !...Plot legend IF (isnd == 1 .AND. do_legend) THEN nitems = 2 strings(1)= 'T' strings(2)= 'Td' ltypes(1) = -1 ltypes(2) = -2 lcolors(1)= jtcolor lcolors(2)= jtdcolor IF (modify_snd) THEN nitems = 4 strings(1) = 'T (modified)' strings(2) = 'Td (modified)' strings(3) = 'T (original)' strings(4) = 'Td (original)' ltypes(1) = -1 ltypes(2) = -2 ltypes(3) = 1 ltypes(4) = 2 lcolors(1) = jtcolor lcolors(2) = jtdcolor lcolors(3) = jtcolor lcolors(4) = jtdcolor END IF IF (plot_tv) THEN nitems = nitems + 1 strings(nitems) = 'Tv' ltypes(nitems) = -6 lcolors(nitems) = jtvcolor END IF IF (plot_sfc_parcel .OR. plot_cb_parcel) THEN nitems = nitems + 1 strings(nitems) = 'T (parcel)' IF (plot_tv) strings(nitems) = 'Tv* (parcel)' ltypes(nitems) = 03 lcolors(nitems) = jparcolor END IF CALL legend (nitems, strings, ltypes,lcolors, & px1+0.01, 0.5*(py1+py2), 0, 0.9*siz0) END IF ! CALL sflush ! uncomment this for NCARG END DO !isnd !----------------------------------------------------------------------= ! FRAME 2 !----------------------------------------------------------------------= IF (.NOT. plot_frame_2) GO TO 9000 CALL xframe PRINT *, '... Plotting frame 2' siz08 = 0.8 * siz0 CALL xchmag (0.9*siz0) CALL xaxnmg (0.9*siz0) ! plot ALWC DO k=1,n alwcgkgpa(k) = alwcpa(k) * 1000. alwcgkgpr(k) = alwcpr(k) * 1000. qvgkg_s(k) = qv_s(k) * 1000. END DO px1 = 0.10 px2 = 0.95 py2 = 0.90 py1 = 0.50 zmin = 0. zmax = 4. DO i=1,10 IF (zzkm(n) > zmax) zmax = zmax + 1. END DO qmin = 0. qmax = 10. CALL pcurves (px1,px2,py1,py2, zmin,zmax,qmin,qmax, 0, & 2,'Height [km]','*',0., 1,'Mixing ratio [g/kg]','(I2)',0., & n,3,zzkm, qvgkg_s,alwcgkgpa,alwcgkgpr,xx, 2,1,4,0, 2,1,4,0) DO k=1,n WRITE (2,*) zz_s(k), alwcpr(k)*1000. END DO strings(1)= 'Water vapor' strings(2)= 'LWC (pseudoadiab)' strings(3)= 'LWC (reversible)' ltypes(1) = 2 ltypes(2) = 1 ltypes(3) = 4 CALL legend (3,strings,ltypes,ltypes, px1,py2-0.02,0,siz08) ! compute potential temperatures py2 = 0.40 py1 = 0.10 DO k=1,n tlcl = ftlcl(temp_s(k),tdew_s(k)) thete(k) = fthetae(pres_s(k),temp_s(k),tlcl,qv_s(k)) thetq(k) = fthetaq(pres_s(k),temp_s(k),qv_s(k), fthq_cwcptrm(qv_s(k))) END DO qmin = 320. qmax = 320. DO i=1,10 DO k=1,n IF (thete(k) > qmax) qmax = qmax + 10. IF (thetq(k) < qmin) qmin = qmin - 10. END DO END DO CALL pcurves (px1,px2,py1,py2, zmin,zmax,qmin,qmax, 0, & 1,'Height [km]','*',0., 1,'Potential Temperature [K]','(I3)',0., & n,2,zzkm, thete,thetq,xx,xx, 1,4,0,0, 1,4,0,0) strings(1)= 'Theta-E' strings(2)= 'Theta-Q' ltypes(1) = 1 ltypes(2) = 4 CALL legend (2,strings,ltypes,ltypes, px1,py1+.10,0,siz08) string = 'File: ' // sndfile(isnd) CALL xchmag (0.6*siz0) CALL plt2wrld (px1,0.03,xx,yy) CALL xcharl (xx,yy, trim(string)) ! CALL sflush ! uncommented this for NCARG 9000 CONTINUE ! CALL sflush ! uncommented this for NCARG CALL xgrend PRINT *, 'SKEWT: Normal completion' STOP END PROGRAM skewt ! ! ######################################## ! ######################################## ! ######################################## ! ######## ######## ! ######## PLOTWIND ######## ! ######## ######## ! ######################################## ! ######################################## ! ######################################## SUBROUTINE plotwind (presmb, dir, spd, n, pspacing, jscheme, & 1,7 label, do_indices) IMPLICIT NONE INTEGER, INTENT(IN) :: n REAL, INTENT(IN OUT) :: presmb(n) REAL, INTENT(IN OUT) :: dir(n) REAL, INTENT(IN OUT) :: spd(n) REAL, INTENT(IN) :: pspacing INTEGER, INTENT(IN OUT) :: jscheme CHARACTER (LEN=*), INTENT(IN) :: label LOGICAL, INTENT(IN OUT) :: do_indices INTEGER :: k REAL :: a,xx,yy,yykm1,xmin,xmax,ymin,ymax, dlnp, px1,px2,py1,py2 REAL, SAVE :: px=0.0 ! spacing between printed levels CALL xqmap (xmin,xmax,ymin,ymax) dlnp = ABS(ymax-ymin)/pspacing yykm1 = 10000. ! location of wind barbs IF (px == 0.0) THEN CALL xqpspc (px1,px2,py1,py2) a = 0.11 IF (do_indices) a = 0.12 px = px2 + a CALL plt2wrld (px,0.0,xx,yy) ELSE a = 0.05 IF (do_indices) a = 0.08 px = px + a CALL plt2wrld (px,0.0,xx,yy) END IF DO k=1,n IF (ABS(dir(k)) > 360.0 .OR. ABS(spd(k)) > 250.0) GO TO 99 yy = LOG (presmb(k)) IF ( yykm1-yy > dlnp ) THEN yykm1 = yy ! IF (jscheme.NE.2) THEN ! IF (spd(k).GT.117.5/2.) THEN ! CALL Color (04) ! Else IF (spd(k).GT.77.5/2.) THEN ! CALL Color (24) ! Else IF (spd(k).GT.37.5/2.) THEN ! CALL Color (23) ! Else ! CALL Color (01) ! END IF ! END IF CALL barb (xx,yy,dir(k),spd(k), 0,0.035) END IF 99 CONTINUE END DO CALL plt2wrld (px,py1-0.02,xx,yy) CALL xcharc (xx,yy,trim(label)) IF (jscheme /= 2) CALL color (01) END SUBROUTINE plotwind ! ######################################## ! ######################################## ! ######################################## ! ######## ######## ! ######## DERIVE ######## ! ######## ######## ! ######################################## ! ######################################## ! ######################################## SUBROUTINE derive (n, pres, zz, theta, qv, uu, vv, & 2 temp, tempc, tdew, tdewc, presmb, spd, dir, qs, rh, tv) IMPLICIT NONE INTEGER, INTENT(IN) :: n REAL, INTENT(IN), DIMENSION(n) :: zz, pres, theta, qv, uu, vv REAL, INTENT(OUT), DIMENSION(n) :: & temp, tempc, tdew, tdewc, presmb, spd, dir, qs, rh, tv INCLUDE 'thermo.consts' INTEGER :: k REAL, PARAMETER :: rad2deg=180/3.14159265 REAL :: vpres INCLUDE 'thermo.stfunc' ! Input: n, pres, zz, theta, qv, uu, vv ! Output: temp, tempc, tdew, tdewc, presmb, spd, dir, qs, rh, tv DO k=1,n temp(k) = theta(k) * (pres(k) * p00inv) ** rcp tempc(k)= temp(k) - tfrz IF (qv(k) <= 0.) THEN tdewc(k) = -243.5 ELSE vpres = fvpres(pres(k),qv(k)) IF (vpres == 611.2) vpres = 611.3 !rlc 1995/12/07 tdewc(k) = ftdewc(vpres) END IF tdew(k) = tdewc(k) + tfrz presmb(k) = pres(k) * 1.e-2 spd(k) = SQRT (uu(k)**2 + vv(k)**2 ) dir(k) = 0. IF (spd(k) /= 0.) dir(k) = 270.- ATAN2(vv(k),uu(k)) * rad2deg dir(k) = MOD (dir(k)+360.,360.) qs(k) = fmixrat(pres(k),fsvpres(temp(k)-tfrz)) rh(k) = frh_q(qv(k),qs(k)) tv(k) = ftvirtnowl(temp(k),qv(k)) END DO END SUBROUTINE derive ! ! ! ! ! ######################################## ! ######################################## ! ######################################## ! ######## ######## ! ######## PLOT_STRINGS ######## ! ######## ######## ! ######################################## ! ######################################## ! ######################################## ! ! SUBROUTINE plot_strings (sndfile, stid, do_indices, & 1,44 cape_use_t, gempak_fmt, modify_snd, plot_cb_parcel, & plot_sfc_parcel, print_info, tv_wl, & isnd , nsnd , iyr , imon , iday , ihr , & imin , px1 , px2 , py1 , py2 , & capet , capetc ,capetvnowl,capetvnowlc,capetvwl , capetvwlc, & cint , cintc , cintvnowl,cintvnowlc, cintvwl , cintvwlc , & pmb_mod , pmbcb , preslcl , preslnb , & sfctf , sfctf0 , sfctdf , sfctdf0 , zcb, zlcl, zlnb , & preslfc , zlfc , convtemp , precipwat , liftedindex , & totaltotals , kindex , sweat , wetbulbzero , brnshear , & lidstrength , dirmean , spmean , helicity , preciptype_me, & preciptype_laps,nmax,n,sndcolor) ! !======================================================================- ! &&&& G E N E R A L I N F O R M A T I O N &&&& !======================================================================- ! ! PURPOSE: This subroutine does x. ! ! AUTHOR: Richard Carpenter, Univ. of Oklahoma (rcarpenter@ou.edu) ! ! HISTORY: ! 1997/04/24 First written ! ! INPUT ARGUMENTS: ! x() : Array ! ! INPUT/OUTPUT ARGUMENTS: ! <none> ! ! OUTPUT ARGUMENTS: ! <none> ! ! I/O: ! <none> ! ! SPECIAL REQUIREMENTS: ! <none> ! ! OTHER INFORMATION: ! <none> ! !======================================================================- ! %%%% D E C L A R A T I O N S %%%% !======================================================================- ! IMPLICIT NONE CHARACTER (LEN=*), INTENT(IN) :: sndfile, stid LOGICAL, INTENT(IN) :: & do_indices, cape_use_t, gempak_fmt, modify_snd, plot_cb_parcel, & plot_sfc_parcel, print_info, tv_wl INTEGER, INTENT(IN) :: isnd, nsnd, iyr, imon, iday, ihr, imin REAL, INTENT(IN) :: & px1, px2, py1, py2, capet, capetc, capetvnowl, capetvnowlc, capetvwl, & capetvwlc, cint, cintc, cintvnowl, cintvnowlc, cintvwl, cintvwlc, & pmb_mod, pmbcb, preslcl, preslnb, sfctf, sfctf0, sfctdf, sfctdf0, & zcb, zlcl, zlnb, preslfc, zlfc, convtemp, precipwat, liftedindex, & totaltotals, kindex, sweat, wetbulbzero, brnshear, lidstrength, & dirmean, spmean, helicity INTEGER, INTENT(IN) :: preciptype_me, nmax, preciptype_laps(nmax), n, sndcolor ! - - - - - - - - - INCLUDE files - - - - - - - - - - - - - - - - - - - ! - - - - - - - - - Constant declarations - - - - - - - - - - - - - - - ! - - - - - - - - - Argument declarations - - - - - - - - - - - - - - - CHARACTER (LEN=85) :: key_laps REAL :: spmkts, a, b ! - - - - - - - - - Global/External declarations - - - - - - - - - - - ! - - - - - - - - - Local declarations - - - - - - - - - - - - - - - - INTEGER :: i REAL :: xx, yy, pxoffset, pyoffset, pydelta, siz0 CHARACTER :: runtype*20, file0*40 CHARACTER (LEN=132) :: string, string2 CHARACTER (LEN=2) :: prcp_string(4) DATA prcp_string / 'SN', 'IP', 'FZ RN', 'RN' / ! - - - - - - - - - DATA statements - - - - - - - - - - - - - - - - - - DATA pyoffset/0.03/, pxoffset/0.01/, pydelta/0.02/, siz0/0.02/ SAVE pyoffset, pxoffset ! - - - - - - - - - Statement functions - - - - - - - - - - - - - - - - !6 & 12345678 , 12345678 , 12345678 , 12345678 , 12345678 , 12345678 , !======================================================================- ! @@@@ E X E C U T A B L E C O D E @@@@ !======================================================================- ! PRINT *, 'PLOT_STRINGS' file0 = sndfile(INDEX(sndfile,'/',back=.true.)+1:) ! filename minus path ! Determine run type string (GEMPAK) IF (gempak_fmt) THEN IF (file0(:2) == 'ar') THEN i = INDEX(file0,'+') IF (i == 0) THEN string = '00' ELSE string = file0(INDEX(file0,'+')+1:INDEX(file0,'_',back=.true.)-1) END IF runtype = trim(string) // 'h FCST' ELSE IF (file0(:2) == 'ad') THEN runtype = 'ADAS' ELSE IF (file0(:2) == 'rc') THEN runtype = 'RUC' ELSE runtype = ' ' END IF END IF ! Label for first file IF (isnd == 1) THEN CALL color (01) string = 'File: ' // file0 IF (gempak_fmt) THEN WRITE (string,9988) stid, iyr,imon,iday, ihr,imin, trim(runtype) 9988 FORMAT (a,1X, i4.2,2('/',i2.2), 1X,2I2.2, 'Z', 1X,a) END IF CALL plt2wrld (0.5,0.02,xx,yy) IF (do_indices) THEN a = 1.0 b = 0.03 ELSE a = 1.2 b = 0.0 END IF CALL xchmag (a*siz0) CALL xcharc (xx,yy+b,trim(string)) PRINT *, 'Comment string = ',string END IF ! Key for LAPS Precip. Type Symbols (EMK, 1/18/98) IF (isnd==1 .AND. do_indices) THEN CALL xchmag (.75*siz0) ! (EMK) key_laps= 'LAPS Potnl Prcp Type: R=. ZR=~ IP=p S=*' CALL xcharc(xx,yy-0.075,trim(key_laps)) END IF ! other labels CALL color (sndcolor) CALL xchmag (0.8*siz0) IF (gempak_fmt) THEN WRITE (string,9987) isnd, stid, iyr,imon,iday, ihr,imin, trim(runtype) 9987 FORMAT ('(',i1,') ',a,1X,i4.2,2('/',i2.2),1X,2I2.2,'Z',1X,a,1X,a) !CALL Plt2Wrld (px1+pxoffset,py1+0.088-0.02*isnd,xx,yy) CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) PRINT *, string pyoffset = pyoffset + pydelta END IF string = 'File: ' // trim(file0) CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta IF (plot_sfc_parcel .OR. print_info) THEN WRITE (string,'(a,i4,a,i5,a)') & ' LCL:', nint(preslcl*1.e-2),' mb,', nint(zlcl),' m' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset= pyoffset + pydelta WRITE (string,'(a,i4,a,i5,a)') & 'LFC:', nint(preslfc*1.e-2),' mb,', nint(zlfc),' m' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset= pyoffset + pydelta WRITE (string,'(a,i4,a,i5,a)') & 'LNB:', nint(preslnb*1.e-2),' mb,', nint(zlnb),' m' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta IF (do_indices) THEN WRITE (string,'(a,f4.1,a,f4.1,a)') & 'Conv.T:',convtemp,' C, Prec.Wat.:',precipwat,' cm' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta WRITE (string,'(a,f5.1,a,i2,a,i2)') & 'LI: ',liftedindex,'C, TT: ',nint(totaltotals),', KI: ' ,nint(kindex) CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta WRITE (string,'(a,i4,a,i4,a)') & 'SWEAT: ',nint(sweat),', W/B zero: ',nint(wetbulbzero), ' m' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta WRITE (string,'(a,f4.1,a,f6.1,a)') & 'BRNshear: ',brnshear,' m/s, LSI:',lidstrength,'K' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta spmkts = spmean * 1.944 WRITE (string,'(a,i3,a,i2,a)') & 'Mean Storm Motion: ',nint(dirmean),' at ',nint(spmean), ' m/s' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta WRITE (string,'(a,i3,a)') 'Storm Rel. Hel.:',nint(helicity),' m2/s2' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta ! Precipitation Types Added to Skew-T. EMK 10/24/1997, 1/18/98 WRITE (string,'(A)') 'MEta/LAPS Sfc Cond Prcp: ' // & TRIM(prcp_string(preciptype_me)) // ' / ' // & TRIM(prcp_string(preciptype_laps(1))) CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset= pyoffset + pydelta END IF ! ORIGINAL PLOT LABELS CMS ! WRITE (string,'(a,i4,a,F5.1,a)') ! > 'LCL:', NINT(presLCL*1.e-2),' mb,', zLCL/1000.0,' km' ! CALL Plt2Wrld (px2-pxoffset,py2-pyoffset,xx,yy) ! CALL XCHARR (xx,yy,TRIM(string)) ! pyoffset = pyoffset + pydelta !c ! WRITE (string,'(a,i4,a,F5.1,a)') ! > 'LNB:', NINT(presLNB*1.e-2),' mb,', zLNB/1000.0,' km' ! CALL Plt2Wrld (px2-pxoffset,py2-pyoffset,xx,yy) ! CALL XCHARR (xx,yy,TRIM(string)) ! pyoffset = pyoffset + pydelta !c ! CAPE string string = ' ' string2 = ' ' IF (cape_use_t) THEN WRITE (string,9951) nint(capet), nint(cint) WRITE (string2,9951) nint(capetc), nint(cintc) ELSE IF (tv_wl) THEN WRITE (string,9951) nint(capetvwl), nint(cintvwl) WRITE (string2,9951) nint(capetvwlc), nint(cintvwlc) ELSE WRITE (string,9951) nint(capetvnowl), nint(cintvnowl) WRITE (string2,9951) nint(capetvnowlc), nint(cintvnowlc) END IF 9951 FORMAT ('CAPE:', i5, ' J/kg CIN:', i4, ' J/kg') CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta END IF ! CAPE of CB parcel IF (plot_cb_parcel) THEN WRITE (string,'(a,i4,a,i5,a)') & 'Cloud base:', nint(pmbcb), ' mb,', nint(zcb), ' m' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta IF (cape_use_t) THEN WRITE (string2,9951) nint(capetc), nint(cintc) ELSE IF (tv_wl) THEN WRITE (string2,9951) nint(capetvwlc), nint(cintvwlc) ELSE WRITE (string2,9951) nint(capetvnowlc), nint(cintvnowlc) END IF CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string2)) pyoffset = pyoffset + pydelta END IF ! Modified snd IF (modify_snd) THEN 9011 FORMAT (a,f5.1,'/',f4.1,'F',1X, f5.1,'/',f4.1,'C') WRITE (string,9011) 'Init sfc T/Td:', & sfctf0, sfctdf0, (sfctf0-32.)*5./9., (sfctdf0-32.)*5./9. CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta WRITE (string,9011) 'Mod sfc T/Td:', & sfctf, sfctdf, (sfctf-32.)*5./9., (sfctdf-32.)*5./9. CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta WRITE (string,'(a,i4,a)') 'Top of modified layer:', nint(pmb_mod), ' mb' CALL plt2wrld (px2-pxoffset,py2-pyoffset,xx,yy) yy = yy - 0.035 ! EMK CALL xcharr (xx,yy,trim(string)) pyoffset = pyoffset + pydelta END IF pyoffset = pyoffset + pydelta ! ready for next file END SUBROUTINE plot_strings