; ; This program creates a skewt diagram from the sounding data ; pro plot_arm_skewt,t_range,temp_wrpn,dp_wrpn,pres_wrpn, graphic_path, graphic_prefix ;;site,date COMMON RANGES, trange, prange common mycolors, my_pmain,my_d_red,my_d_dgreen,my_d_blue,my_color ; ; Make sure everything passed in comes in as a string ; ;site=strcompress(site,/remove_all) ;date=strcompress(date,/remove_all) ; ; Establish the location of the data based on the site ; ;if site eq 'nsa' then begin ; armdatadir='/data/mace2/arm_data/nsa/' ; datastream='nsaisssonde10sC1.a1' ; t_range = [-30, 35] ;endif else if site eq 'sgp' then begin ; armdatadir='/data/mace3/arm_data/sgp/' ; if long(date) lt 20001124 then $ ; datastream='sgpsondewrpnC1.a1' $ ; else datastream='sgpsondewnpnC1.a1' ; t_range = [-30, 35] ;endif else if site eq 'twpc1' then begin ; armdatadir='/data/mace/arm_data/twp/Manus_C1/' ; datastream='twpsondewnpnC1.a1' ; t_range = [-30, 35] ;endif else if site eq 'twpc2' then begin ; armdatadir='/data/mace/arm_data/twp/Nauru_C2/' ; datastream='twpsondewnpnC2.a1' ; t_range = [-30, 35] ;endif else if site eq 'twpc3' then begin ; armdatadir='/data/mace/arm_data/twp/Darwin_C3/' ; datastream='twpsondewnpnC3.a1' ; t_range = [-30, 35] ; p_range = [1050,50] ;endif ;filename=armdatadir+datastream+'/'+datastream+'.'+date+'*.cdf*' ; ; Use z-buffer to generate gifs ; ; ; set z buffer image to white ; ;tempimage=make_array(700,800,/byte,value=255) ;tv,tempimage ; ; Set up the colors ; common colors,r_orig,g_orig,b_orig,r_curr,g_curr,b_curr loadct,12, ncolors=32 my_pmain=0 my_d_red=24 ;red temperature lines my_d_dgreen=3 ;green saturation adiabats my_d_blue=15 ;purple temp lines my_color=12 ;sounding r_curr[250:255]=255 & r_orig=r_curr g_curr[250:255]=255 & g_orig=g_curr b_curr[250:255]=255 & b_orig=b_curr !p.background=!d.n_colors set_plot,'Z' device,set_resolution=[700,800] ; ; Stores the names of the files that match our date ; ;unix_command='gunzip '+filename ;spawn,unix_command ;files = findfile(filename, count=count) ; ; Found matching files ; ;if count gt 0 then begin ; ; Do for all the matching files in the directory ; ; for j=0, count-1 do begin ; ; Open the netcdf file ; ; cdfid = ncdf_open(files[j]) ;************************************************ ; Read in data and cut array to only good values, which are ; defined to be from surface pressure to pressure minimum ;************************************************ ; ; Read in the data ; ; ncdf_varget,cdfid,'base_time',base_time ;base time ; ncdf_varget,cdfid,'time_offset',time_offset ;time offset ;; Pressure can be named press or pres ; pres_id=ncdf_varid(cdfid,'press') ; if pres_id eq -1 then pres_id=ncdf_varid(cdfid,'pres') ; ncdf_varget,cdfid,pres_id,pres ;pressure ; ncdf_varget,cdfid,'tdry',temp ;temperature ; ncdf_varget,cdfid,'dp',dp ;dewpoint ; ncdf_varget,cdfid,'wspd',wspd ;wind speed ; ncdf_varget,cdfid,'deg',wdir ;wind direction ; ncdf_varget,cdfid,'rh', rh ;relative humidity ;; ; Find the minimum pressure value ; ; presmin=min(pres,min_index) ;minimum in the pressure data ; ; Cut the arrays to be from the surface to the pressure minimum ; ; time_offset=time_offset[0:min_index] ; pres=pres[0:min_index] ; temp=temp[0:min_index] ; dp=dp[0:min_index] ; wspd=wspd[0:min_index] ; wdir=wdir[0:min_index] ; rh=rh[0:min_index] ; ; Close the netcdf file ; ; ncdf_close,cdfid ; ; Restrict the arrays to only good values ; ; result=where(temp ne 999 and dp ne 999,count) ; time_offset=time_offset[result] ; pres=pres[result] ; dp=dp[result] ; temp=temp[result] ; wspd=wspd[result] ; wdir=wdir[result] ; ; Read hhmmss ; ; sepstring=strsplit(files[j],'.',/extract) ; hms = sepstring[4] ; ; If it is the first file ; ; if (j eq 0) then begin ; btto = base_time+time_offset[0] ; num_elem_wrpn = n_elements(time_offset) ; hms_wrpn = hms ; pres_wrpn = pres ; dp_wrpn = dp ; temp_wrpn = temp ; wspd_wrpn = wspd ; wdir_wrpn = wdir ; endif $ ;end of the first file ; ; If it isn't the first file. Put everything from the different ; files into one array ; ; else begin ;if it is not the first file ; btto = [btto, base_time+time_offset[0]] ; num_elem_wrpn = [num_elem_wrpn, n_elements(time_offset)] ; hms_wrpn = [hms_wrpn, hms] ; pres_wrpn = [pres_wrpn, pres] ; dp_wrpn = [dp_wrpn, dp] ; temp_wrpn = [temp_wrpn, temp] ; wspd_wrpn = [wspd_wrpn, wspd] ; wdir_wrpn = [wdir_wrpn, wdir] ; endelse ;end of not the first file ; endfor ;end of loop through all matching files ;***************************************** ; Plot skew-T diagrams for the wrpn file ;***************************************** mnth = [' Jan ',' Feb ',' Mar ',' Apr ',' May ',' Jun ', $ ' Jul ',' Aug ',' Sep ',' Oct ',' Nov ',' Dec '] ; for i=0, n_elements(hms_wrpn)-1 do begin ; ; Prepare some strings for the title ; ; yyyy=strmid(date,0,4) ; mm=strmid(date,4,2) ; dd=strmid(date,6,2) ; hh=strmid(hms_wrpn[i],0,2) ; nn=strmid(hms_wrpn[i],2,2) ; ss=strmid(hms_wrpn[i],4,2) ; ; Call skewt subroutine ; skewt, t_range ; ; Make a date and time string ; ; datestring = string(format='(I0)',dd(0)) + mnth(mm(0)-1) + $ ; string(format='(I0)',yyyy(0)) ; timestring = hh+':'+nn+':'+ss ; ; Title the plot ; ;title=datastream+' sonde data for '+datestring+', '+timestring title='Capricorn Voyage Sounding '+graphic_prefix xyouts,0.5,0.95,title,color=my_pmain,align=0.5,charsize=1.5,/normal ; ; Plot WRPN data ; ; if (i eq 0) then first = 0 $ ; else first = fix(total(num_elem_wrpn[0:i-1])) ; if (i eq 0) then last=fix(total(num_elem_wrpn[0:i]-1)) $ ; else last = fix(total(num_elem_wrpn[0:i])-1) ; ; Call plot_skewt subroutine ; ; plot_skewt,temp_wrpn[first:last],dp_wrpn[first:last],pres_wrpn[first:last] plot_skewt,temp_wrpn,dp_wrpn,pres_wrpn ; ; Plot sonde winds ; ; windbarb, wdir_wrpn[first:last], wspd_wrpn[first:last], pres_wrpn[first:last] ; ; Create the gif ; ;outname = datastream+'.'+date+'.'+hms_wrpn[i] ;gname='/data/mace4/arm_data/temp_image/'+outname+'.gif' gname='skewt.gif' ; write_gif,gname,tvrd() write_gif, graphic_path+graphic_prefix+'.gif', tvrd() ;unix_command='rm '+'/data/mace4/arm_data/temp_image/'+outname+'.png' ;spawn,unix_command ; endfor ;end of loop through all the wrpn data ; ; Close the zbuffer ; device,/close ; ; Reset these variables to zero ; !p.multi=0 !p.region=0 !p.position=0 ; ; Zip up the files ; ; unix_command='gzip '+filename ; spawn,unix_command ; ; Fix permissions ; ; unix_command='chmod 664 '+filename ; spawn,unix_command ; ; Fix group ; ; unix_command='chgrp met-mace '+filename ; spawn,unix_command ;endif ;found files end ;------------------------------------------------------------------------------------------------------------------------------------------------------------ ;------------------------------------------------------------------------------------------------------------------------------------------------------------ ;------------------------------------------------------------------------------------------------------------------------------------------------------------ ; HERE IS WHERE THE SKEWT STUFF STARTS ;------------------------------------------------------------------------------------------------------------------------------------------------------------ ;------------------------------------------------------------------------------------------------------------------------------------------------------------ ;------------------------------------------------------------------------------------------------------------------------------------------------------------ ; $Id: skewt.pro,v 1.1 1999/11/23 21:24:27 halter Release_1_0 $ ;======================================================================== ; SKEWT.PRO (IDL CODE) ; ; Draw a Skew-T, Log(P) diagram given a temperature range for your data. ; ; Originator: Andrew F. Loughe (afl@cdc.noaa.gov) ; CIRES/NOAA ; Boulder, CO USA ; This code carries no warranty or claim ; as to its usefulness or accuracy! ; ; A Number of the functions found in this file were converted from ; FORTRAN code that was received from NCAR in Boulder, CO USA. ; The original source of the equations is thought to be: ; "Algorithms for Generating a Skew-T, Log P Diagram ; and Computing Selected Meteorological Quantities" ; by G.S. Stipanuk, White Sands Missle Range, Report ECOM-5515. ; ;======================================================================== ; FUNCTION TO COMPUTE SATURATION VAPOR PRESSURE GIVEN TEMP IN KELVIN. ; ESAT_sk(MILLIBARS), T(KELVIN) FUNCTION ESAT_sk, T TC = T - 273.16 ESAT_sk = 6.1078 * EXP( (17.2693882 * TC) / (TC + 237.3) ) RETURN, ESAT_sk END ;======================================================================== ; FUNCTION TO COMPUTE MIXING RATIO GIVEN TEMP. AND PRESS. ; W(GRAMS WATER VAPOR/KILOGRAM DRY AIR), P(MILLIBAR) FUNCTION W_sk, T, P IF (T GE 999.) THEN GOTO, JUMP10 X = ESAT_sk(T) W_sk = 621.97 * X / (P - X) RETURN, W_sk JUMP10: W_sk = 0.0 RETURN, W_sk END ;======================================================================== ; FUNCTION TO COMPUTE SATURATION ADIABATIC TEMP AT 1000 MB GIVEN T & P. ; OS AND T (KELVIN), P (MILLIBARS ) FUNCTION OS_sk, T, P IF (T LT 100.) THEN T1 = T + 273.16 IF (T GE 100.) THEN T1 = T OS_sk = T1 * ((1000./P)^.286) / (EXP( -2.6518986*W_sk(T1,P) / T1) ) RETURN, OS_sk END ;======================================================================== ; FUNCTION TO COMPUTE THE TEMPERATURE (KELVIN) OF AIR AT A GIVEN ; PRESSURE AND WITH A GIVEN MIXING RATIO. ; TMR(KELVIN), W(GRAMS WATER VAPOR/KILOGRAM DRY AIR), P(MILLIBAR) FUNCTION TMR_sk, W, P X = ALOG10 ( W * P / (622.+ W) ) TMR_sk = 10. ^ ( .0498646455 * X + 2.4082965 ) - 7.07475 + $ 38.9114 * ( (10.^( .0915 * X ) - 1.2035 )^2 ) RETURN, TMR_sk END ;======================================================================== ; FUNCTION TO COMPUTE TEMPERATUE (KELVIN) OF A MOIST ADIABAT GIVEN ; OS(KELVIN), P(MILLIBARS) ; SIGN(A,B) REPLACES THE ALGEBRAIC SIGN OF A WITH THE SIGN OF B FUNCTION TSA_sk, OS, P A = OS TQ = 253.16 D = 120 FOR I = 1, 12 DO BEGIN D = D/2. ; IF THE TEMPERATURE DIFFERENCE, X, IS SMALL, EXIT THIS LOOP. X = A * EXP (-2.6518986*W_sk(TQ,P)/TQ)-TQ*((1000./P)^.286) IF ( ABS(X) LT 0.01 ) THEN GOTO, JUMP2 ; TQ = TQ + SIGN(D,X) IF (X LT 0) THEN D = -ABS(D) IF (X GT 0) THEN D = ABS(D) TQ = TQ + D ENDFOR JUMP2: TSA_sk=TQ RETURN, TSA_sk END ;======================================================================== ; Function to determine osition (temp, press) when the isotherms ; in the diagram are rotated (skewed) 45 degrees to the right. ; Used for finding the points needed to connect the dots when ; drawing ALL of the lines (except pressure). ; Originator: Andrew F. Loughe Function Tnew, T, P COMMON RANGES, trange, prange P0 = prange(0) xy1 = convert_coord( [T, P0], /data, /to_device) xy2 = convert_coord( [T, P], /data, /to_device) dy = xy2(1) - xy1(1) dx = dy ; dx = dy for this 45-45-90 triangle xy = convert_coord( [xy2(0)+dx, xy2(1)], /device, /to_data) Tnew = xy(0) return, Tnew end ;======================================================================== ; Function to determine osition (temp, press) in the unskewed ; coordinate system (Opposite of Tnew). ; Used only when placing the labels on various lines. ; Originator: Andrew F. Loughe Function Told, T, METHOD COMMON RANGES, trange, prange P0 = prange(0) P1 = prange(1) T0 = trange(0) T1 = trange(1) if (method eq 1) then begin xy1 = convert_coord( [T, P0], /data, /to_device ) xy2 = convert_coord( [T0, P0], /data, /to_device ) dx = xy2(0) - xy1(0) xy = convert_coord( [xy2(0), xy2(1)+dx], /device, /to_data ) xy1 = convert_coord( [xy(0), xy(1)], /data, /to_device ) xy2 = convert_coord( [xy(0), P1], /data, /to_device ) dy = xy2(1) - xy1(1) xy = convert_coord([xy1(0)+(dy/2.), xy1(1)+(dy/2.)],$ /device, /to_data) endif if (method eq 2) then begin xy1 = convert_coord( [T, P0], /data, /to_device ) xy2 = convert_coord( [T1, P0], /data, /to_device ) dx = xy2(0) - xy1(0) xy = convert_coord( [xy1(0)+dx/2. , xy1(1)+dx/2.], $ /device, /to_data) endif return, xy end ;======================================================================== ; Function to determine the necessary trange for plotting the sounding. Function T_RANGE, p_range, t, td, p COMMON RANGES, trange, prange trange = [-40, 40] ; Default which can be changed if (p_range eq 0) then $ prange = [1050, 100] ; Default which can be changed ; Find number of data levels szd = size(t) nlevels = szd(1) ; Ensure that temperatures are in Celsius. if ( (total(t)/nlevels) gt 100. ) then t = t - 273.16 if ( (total(td)/nlevels) gt 100. ) then td = td - 273.16 ; Set up dummy plot space. if (!d.name eq 'PS') then device, /inch, xsize=7, ysize=7 daspect = FLOAT(!D.Y_SIZE)/FLOAT(!D.X_SIZE) * $ (trange(1)-trange(0))/80. margin = 0.1 aspect = 1.0 ; A square x0 = 0.50 - (0.5 - margin)*(daspect/aspect) y0 = margin x1 = 0.50 + (0.5 - margin)*(daspect/aspect) y1 = 1.0 - margin !P.position=[x0,y0,x1,y1] plot_io, trange, prange, yrange=prange, /nodata, /xs, /ys, $ color=!p.background, /noerase !p.multi(0) = !p.multi(0) + 1 ; Determine necessary temperature range for the diagram. xx0 = fltarr(nlevels) & yy0=xx0 & xx1=xx0 & yy1=yy0 for i = 0, nlevels-1 do begin xx0(i) = tnew( t(i), p(i) ) yy0(i) = p(i) xx1(i) = tnew( td(i), p(i) ) yy1(i) = p(i) endfor xbegin = fix( (min(xx1)-10.) / 10. ) * 10. xend = fix( (max(xx0)+10.) / 10. ) * 10. return, [xbegin, xend] end ;======================================================================== ; Routine to plot the temperature and dew point temperature sounding ; on top of a skew-T, Log(P) diagram. ; ; Keyword(s): ; COLOR: If set, the color of the line to be plotted over the data line ; for purposes of data visualization pro plot_skewt, t, td, p, color=color common mycolors, my_pmain,my_d_red,my_d_dgreen,my_d_blue,my_color pmain=my_pmain d_red=my_d_red d_dgreen=my_d_dgreen d_blue=my_d_blue color=my_color COMMON RANGES, trange, prange ; Find number of data levels szd = size(t) nlevels = szd(1) ; Define clipping space. clip=[trange(0),prange(0),trange(1),prange(1)] ; Ensure that temperatures are in Celsius. if ( (total(t)/nlevels) gt 100. ) then t = t - 273.16 if ( (total(td)/nlevels) gt 100. ) then td = td - 273.16 ; Overplot the data onto the digram. for i = 0, nlevels-2 do begin ; Plot temperature sounding data. x0 = tnew( t(i), p(i) ) y0 = p(i) x1 = tnew( t(i+1), p(i+1) ) y1 = p(i+1) if (keyword_set(color) eq 1) then plots, [x0, x1], [y0, y1], clip=clip, noclip=0, $ thick=2, color=color $ else plots, [x0, x1], [y0, y1], clip=clip, noclip=0, thick=2, color=pmain ; Plot dew point temperature sounding data. x0 = tnew( td(i), p(i) ) y0 = p(i) x1 = tnew( td(i+1), p(i+1) ) y1 = p(i+1) if (keyword_set(color) eq 1) then plots, [x0, x1], [y0, y1], clip=clip, noclip=0, $ thick=2, color=color $ else plots, [x0, x1], [y0, y1], clip=clip, noclip=0, thick=2, color=pmain endfor ; Close Postscript device, rename output file. ;tdh if ( !d.name eq 'PS') then begin ;tdh device, /close ;tdh spawn, 'mv idl.ps skewt.ps' ;tdh set_plot, 'X' ;tdh endif end ;======================================================================== ; ; PROCEDURE TO DRAW A SKEW-T, Log(P) DIAGRAM GIVEN A DESIRED ; TEMPERATURE RANGE FOR THE DATA. ; ; Originator: Andrew F. Loughe ; PRO SKEWT, t_range, everyT=everyT, everyDA=everyDA, $ everySA=everySA, everyW=everyW, title=title, notitle=notitle on_error, 2 COMMON RANGES, trange, prange common mycolors, my_pmain,my_d_red,my_d_dgreen,my_d_blue,my_color pmain=my_pmain d_red=my_d_red d_dgreen=my_d_dgreen d_blue=my_d_blue if (n_elements(everyT) le 0) then everyT = 10 ; T = Temperature if (n_elements(everyDA) le 0) then everyDA = 10 ; DA = Dry adiabat if (n_elements(everySA) le 0) then everySA = 1 ; SA = Saturated adiabat if (n_elements(everyW) le 0) then everyW = 1 ; W = Mixing ratio if (not keyword_set(title)) then title='Skew-T, Log(P) Diagram' if (keyword_set(notitle)) then title=' ' if (n_elements(prange)) eq 0 then prange = [1050., 100] if (N_params() eq 0) then $ message,$ 'EXAMPLE: skewt, [-20, 20], everyT=10, everyDA=10, everySA=2, everyW=2' if (n_elements(t_range)) eq 1 then t_range=[-40., 40.] ; Set some defaults trange = t_range charsize=1 ;charsize = 1.25 ; Set default character size ;window, 0, xsize=800, ysize=800, /pixmap !p.multi = [0,1,1] ; Make plot square for arbitrarily chosen trange of 80 degrees. ; Code from Ken Bowman if (!d.name eq 'PS') then device, /inch, xsize=7, ysize=7 daspect = FLOAT(!D.Y_SIZE)/FLOAT(!D.X_SIZE) * (trange(1)-trange(0))/80. margin = 0.1 aspect = 1.0 ; A square x0 = 0.50 - (0.5 - margin)*(daspect/aspect) y0 = margin x1 = 0.50 + (0.5 - margin)*(daspect/aspect) y1 = 1.0 - margin !P.POSITION = [x0, y0, x1, y1] ; Set value of sytem variable. ; Determine character height and width. Apply charsize. char_ht = convert_coord([0, !d.y_ch_size], /device, /to_norm) char_ht = char_ht(1) * 1.0 if (!d.name ne 'X' and charsize gt 1.) then $ char_ht = char_ht * charsize char_wd = convert_coord([0, !d.x_ch_size], /device, /to_norm) char_wd = char_wd(1) ; Create the plot space. plot_io, trange, prange, yrange=prange, /nodata, xs=5, /ys, $ xticklen=.01, ytickname=replicate(' ',30), charsize=charsize, $ title=title, color=pmain, background=!d.n_colors-1 ; Print PRESSURE title along the y-axis. lnt=alog(prange(1)) & lnb=alog(prange(0)) & avg=exp(.5*(lnt+lnb)) xy = convert_coord([trange(0), avg],/data,/to_norm) xyouts, xy(0)-(5.*char_wd), xy(1), 'PRESSURE (hPa)', orient=90, $ /norm, align=.5, color=pmain ; Print TEMPERATURE title along the x-axis. xy = convert_coord([.5*(trange(0)+trange(1)), prange(0)], /data, /to_norm) ;xyouts, xy(0), xy(1)-(3.*char_ht), 'TEMPERATURE (!uo!nC)', align=.5, /norm, color=pmain ; Draw Pressure labels next to tick marks along the y-axis. pressures = [1050,1000,900,800,700,600,500,400,300,200,100] for i = 0, 10 do begin ytick = pressures(i) if (ytick ge prange(1)) then begin xy = convert_coord( [trange(0), ytick], /data, /to_norm ) xyouts, xy(0)-(.2*char_wd), xy(1)-(.25*char_ht), $ strcompress(string(ytick),/remove_all), align=1, $ charsize=charsize, /norm, color=pmain plots, [trange(0), trange(1)], [ytick, ytick], color=pmain ; Horizontal line endif endfor clip=[trange(0),prange(0),trange(1),prange(1)] ; Define clipping space. ;======================================================================== ; Draw skewed isotherms every "everyT (10C)" (Lines are straight). for temp = trange(0)-100, trange(1)+5, everyT do begin x0 = temp y0 = prange(0) x1 = temp y1 = prange(1) ; Draw the line. newx0 = tnew(x0, y0) ; Find rotated temperature position newx1 = tnew(x1, y1) ; Find rotated temperature position plots, [newx0, newx1], [y0, y1], color=pmain, clip=clip, noclip=0 ; Draw line labels ; Use method #1 in xy function to determine a place for the label. drew_label = 'no' xy = Told(temp, 1) if ( xy(0) gt trange(0) and xy(0) lt trange(1) and $ xy(1) gt prange(1) and xy(1) lt prange(0) ) then begin drew_label = 'yes' xyouts, xy(0), xy(1), strcompress(string(fix(temp)), /rem),$ color=pmain, orient=45, align=.5, charsize=charsize endif ; Use method #2 in xy function to determine a place for the label. if (drew_label eq 'no') then xy = Told(temp, 2) if ( xy(0) gt trange(0) and xy(0) lt trange(1) and $ xy(1) gt prange(1) and xy(1) lt prange(0) and $ drew_label eq 'no') then begin xyouts, xy(0), xy(1), strcompress(string(fix(temp)), /rem),$ color=pmain, orient=45, align=.5, charsize=charsize endif endfor ;======================================================================== ; Draw dry adiabats every "everyDA (10C)" (Lines are curved). for temp = trange(0), trange(0)+220, everyDA do begin x1 = float(temp) y1 = 1050. inc = -2. ; Lines will be curved, so use a small press. increment. drew_label='no' icount = 0 ; Dry adiabats from 1050mb up to prange(1). ; For a given temperature and pressure, compute theta and plot a line. for press = y1, prange(1), inc do begin icount = icount + 1 x0 = float(x1) ; Orig Temp y0 = float(press + inc) ; Orig Press y1 = float(y0 + inc) ; New Press x1 = (temp+273.16) * ( y1 / 1000. ) ^ (287./1004.) ; New Temp x1 = x1 - 273.16 newx0 = tnew(x0, y0) ; Find rotated temperature position newx1 = tnew(x1, y1) ; Find rotated temperature position ; Draw the labels. if (fix(x1) eq fix(trange(0)) and drew_label eq 'no') then begin drew_label='yes' if ( newx1 gt trange(0) and newx1 lt trange(1) and $ y1 gt prange(1) and y1 lt prange(0) ) then $ xyouts,newx1,y1,strcompress(string(fix(temp)),/remove),$ align=.5, color=d_red, charsize=charsize, orientation=-45 endif ; Draw the line. if (icount gt 1) then $ plots, [newx0, newx1], [y0, y1], color=d_red, clip=clip, noclip=0 if (newx1 lt trange(0)) then goto, jump2 endfor jump2: dummy=0 endfor ;======================================================================== ; Draw saturated adiabats. Begin at 40C and step backwards by 4C. ; These lines are curved. TS = 40. FOR TS = 40, -64, -everySA*4 DO BEGIN P = 1060. TK = TS + 273.16 AOS = OS_sk(TK, 1000.) ATSA = TSA_sk(AOS, P) - 273.16 FOR J = 0, 85 DO BEGIN P0 = P T0 = ATSA P = P - 10. ATSA = TSA_sk(AOS, P) - 273.16 if (j gt 0) then begin newx0=tnew(T0,P0) ; Find rotated temperature position newx1=tnew(ATSA,P) ; Find rotated temperature position ; Leave a space for the labels and draw them. if (P gt 730 or P lt 700) then $ plots, [newx0, newx1], [P0, P], $ color=d_dgreen, clip=clip, noclip=0 if ( P eq 730 ) then begin if (newx1 gt trange(0) and newx1 lt trange(1)) then $ xyouts,newx1,P,strcompress(string(fix(TS)),/remove),align=.5,$ color=d_dgreen, charsize=charsize endif endif ENDFOR ENDFOR ;======================================================================== ; Draw mixing ratio lines (Lines are straight). ; Find temperature for a given Ws (g/kg) and Press (mb). Ws=[ .1,.2,.4,.6,.8,1.,1.5,2,3,4,5,6,7,8,9,10,12, $ 14,16,18,20,24,28,32,36,40,44,48,52,56,60,68,76,84 ] for i = 0, N_elements(Ws)-1, everyW do begin press1 = prange(0) - 20 tmr1 = tmr_sk(Ws(i), press1) - 273.16 press2 = 200. tmr2 = tmr_sk(Ws(i), press2) - 273.16 newx0=tnew(tmr1,press1) ; Find rotated temperature position newx1=tnew(tmr2,press2) ; Find rotated temperature position ; Draw the line. plots, [newx0, newx1], [press1, press2], linestyle=0, $ clip=clip, noclip=0 , color=d_blue ; Draw the line label. drew_label='no' if (newx0 gt trange(0) and newx0 lt trange(1)) then begin drew_label='yes' if (Ws(i) ge 3.0) then $ xyouts, newx0, press1+70, strcompress(string(fix(Ws(i))),/remove),$ align=.5,color=d_blue, charsize=charsize-.3 if (Ws(i) lt 3.0) then $ xyouts, newx0, press1+70, string(Ws(i),format='(f3.1)'), align=.5,$ color=d_blue, charsize=charsize-.3 endif if (newx1 gt trange(0) and newx1 lt trange(1)) then begin if (Ws(i) ge 2.0) then $ xyouts, newx1, press2-2, strcompress(string(fix(Ws(i))),/remove),$ align=.5, color=d_blue, charsize=charsize-.2 if (Ws(i) lt 2.0) then $ xyouts, newx1, press2-2, string(Ws(i),format='(f3.1)'), align=.5,$ color=d_blue, charsize=charsize-.2 endif endfor ;======================================================================== ; Redraw the plot boundary. plots, [trange(0),trange(1),trange(1),trange(0),trange(0)], $ [prange(0),prange(0),prange(1),prange(1),prange(0)], thick=2, color=pmain !p.position = [.05, .05, .95, .95] ; Reset position parameter END