Jeff,<br><br>I liked the idea, but that didn't help. I tried the following, both resulting in a black background for the map area:<br><br>printim precip_NE.png png x570 y435 -b /home/jeff/public_html/images/NE.png -t 00<br>
<br>printim precip_NE.png png x570 y435 -t 00 -b /home/jeff/public_html/images/NE.png<br><br>Example: <a href="http://jeffsweatherservice.com/grads/nam/precip_NE_2.png">http://jeffsweatherservice.com/grads/nam/precip_NE_2.png</a><br>
<br>I even tried using "00" in the color table, but that had the same impact.<br><br>Now, if I use the following command, the black background becomes white:<br><br>printim precip_NE.png png x570 y435 -b /home/jeff/public_html/images/NE.png white<br>
<br>Thanks for the idea though.<br><br><br>Jeff<br><br><br><br><br><br><br><br><br><div class="gmail_quote">On Mon, Nov 8, 2010 at 10:59 AM, <span dir="ltr"><<a href="mailto:gradsusr-request@gradsusr.org">gradsusr-request@gradsusr.org</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">Send gradsusr mailing list submissions to<br>
<a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
<br>
To subscribe or unsubscribe via the World Wide Web, visit<br>
<a href="http://gradsusr.org/mailman/listinfo/gradsusr" target="_blank">http://gradsusr.org/mailman/listinfo/gradsusr</a><br>
or, via email, send a message with subject or body 'help' to<br>
<a href="mailto:gradsusr-request@gradsusr.org">gradsusr-request@gradsusr.org</a><br>
<br>
You can reach the person managing the list at<br>
<a href="mailto:gradsusr-owner@gradsusr.org">gradsusr-owner@gradsusr.org</a><br>
<br>
When replying, please edit your Subject line so it is more specific<br>
than "Re: Contents of gradsusr digest..."<br>
<br>
<br>
Today's Topics:<br>
<br>
1. Re: GrADS help, plotting shaded contour precip/ptype products<br>
over a basemap/background (Jeffrey Duda)<br>
<br>
<br>
----------------------------------------------------------------------<br>
<br>
Message: 1<br>
Date: Mon, 8 Nov 2010 14:04:13 -0600<br>
From: Jeffrey Duda <<a href="mailto:jdduda@iastate.edu">jdduda@iastate.edu</a>><br>
Subject: Re: [gradsusr] GrADS help, plotting shaded contour<br>
precip/ptype products over a basemap/background<br>
To: GrADS Users Forum <<a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a>><br>
Message-ID:<br>
<AANLkTinO=D7vgkJzwo=+60vLO5VgOZjp=<a href="mailto:4Hhr7P93qo%2B@mail.gmail.com">4Hhr7P93qo+@mail.gmail.com</a>><br>
Content-Type: text/plain; charset="iso-8859-1"<br>
<br>
Jeff,<br>
Try using a two-digit color number (i.e., 00) in the printim command as per<br>
the recommendation on the index page for the printim command:<br>
<a href="http://www.iges.org/grads/gadoc/gradcomdprintim.html" target="_blank">http://www.iges.org/grads/gadoc/gradcomdprintim.html</a><br>
<br>
Jeff<br>
<br>
On Mon, Nov 8, 2010 at 1:41 PM, Jeff Chabot <<a href="mailto:jsc219@gmail.com">jsc219@gmail.com</a>> wrote:<br>
<br>
> I am using two separate color tables. For precip:<br>
><br>
> *Categorical Rain<br>
> set rgb 30 210 180 140<br>
> set rgb 16 107 159 91<br>
> set rgb 17 90 142 75<br>
> set rgb 18 74 125 60<br>
> set rgb 19 58 109 45<br>
> set rgb 20 41 92 19<br>
> set rgb 21 25 75 14<br>
> set rgb 22 253 248 2<br>
> set rgb 23 229 188 0<br>
> set rgb 24 253 149 0<br>
> set rgb 25 253 0 0<br>
> set rgb 26 212 0 0<br>
> set rgb 27 255 0 153<br>
><br>
> set clevs 0.01 0.05 0.10 0.25 0.50 1.00 2.00 3.00 4.00 5.00<br>
> set ccols 0 16 18 19 20 22 23 24 25 26 27<br>
><br>
> For precip type:<br>
><br>
> *Categorical Rain<br>
> set rgb 16 107 159 91<br>
> set rgb 17 90 142 75<br>
> set rgb 18 74 125 60<br>
> set rgb 19 58 109 45<br>
> set rgb 20 41 92 19<br>
> set rgb 21 25 75 14<br>
> set rgb 22 253 248 2<br>
> set rgb 23 229 188 0<br>
> set rgb 24 253 149 0<br>
> set rgb 25 253 0 0<br>
> set rgb 26 212 0 0<br>
> set rgb 27 255 0 153<br>
><br>
> *Categorical Freezing Rain<br>
> set rgb 28 255 188 188<br>
> set rgb 29 255 172 172<br>
> set rgb 30 255 156 156<br>
> set rgb 31 253 140 140<br>
> set rgb 32 255 124 124<br>
> set rgb 33 255 108 108<br>
> set rgb 34 255 96 96<br>
> set rgb 35 255 80 80<br>
> set rgb 36 255 64 56<br>
> set rgb 37 240 32 128<br>
> set rgb 38 240 16 255<br>
><br>
> *Categorical Ice Pellets / Mix<br>
> set rgb 39 255 200 0<br>
> set rgb 40 255 180 0<br>
> set rgb 41 255 160 0<br>
> set rgb 42 255 140 0<br>
> set rgb 43 255 124 0<br>
> set rgb 44 255 108 0<br>
> set rgb 45 255 96 0<br>
> set rgb 46 255 80 0<br>
> set rgb 47 255 64 0<br>
> set rgb 48 255 32 0<br>
> set rgb 49 255 16 0<br>
><br>
> *Categorical Snow<br>
> set rgb 50 4 233 231<br>
> set rgb 51 0 173 255<br>
> set rgb 52 0 148 255<br>
> set rgb 53 0 123 255<br>
> set rgb 54 0 104 255<br>
> set rgb 55 0 85 255<br>
> set rgb 56 4 67 245<br>
> set rgb 57 0 38 255<br>
> set rgb 58 0 14 255<br>
> set rgb 59 0 0 255<br>
> set rgb 60 0 0 223<br>
> set rgb 61 255 255 255<br>
><br>
> set clevs 0.01 0.05 0.1 0.25 0.5 1 2 3 4 5 10 10.01 10.05 10.1 10.25 10.5<br>
> 11 12<br>
> 13 14 15 20 20.01 20.05 20.1 20.25 20.5 21 22 23 24 25 30 30.01 30.05 30.1<br>
> 30.25<br>
> 30.5 31 32 33 34 35<br>
> set ccols 0 16 18 19 20 22 23 24 25 26 27 0 39 40 41 42 43 44 45 46 47 0 0<br>
> 28 29<br>
> 30 31 32 33 34 35 36 0 0 50 51 52 53 54 55 56 57 58 59<br>
><br>
> So, Jeff, I believe that I am using #0 correctly. But, when I do, it just<br>
> uses white or black as the background as in this example:<br>
><br>
> <a href="http://jeffsweatherservice.com/grads/nam/ptype_NE_2.gif" target="_blank">http://jeffsweatherservice.com/grads/nam/ptype_NE_2.gif</a><br>
><br>
> Thanks,<br>
><br>
><br>
> Jeff<br>
><br>
><br>
><br>
> On Mon, Nov 8, 2010 at 10:18 AM, <<a href="mailto:gradsusr-request@gradsusr.org">gradsusr-request@gradsusr.org</a>> wrote:<br>
><br>
>> Send gradsusr mailing list submissions to<br>
>> <a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
>><br>
>> To subscribe or unsubscribe via the World Wide Web, visit<br>
>><br>
>> <a href="http://gradsusr.org/mailman/listinfo/gradsusr" target="_blank">http://gradsusr.org/mailman/listinfo/gradsusr</a><br>
>> or, via email, send a message with subject or body 'help' to<br>
>> <a href="mailto:gradsusr-request@gradsusr.org">gradsusr-request@gradsusr.org</a><br>
>><br>
>> You can reach the person managing the list at<br>
>> <a href="mailto:gradsusr-owner@gradsusr.org">gradsusr-owner@gradsusr.org</a><br>
>><br>
>> When replying, please edit your Subject line so it is more specific<br>
>> than "Re: Contents of gradsusr digest..."<br>
>><br>
>><br>
>> Today's Topics:<br>
>><br>
>> 1. GrADS help, plotting shaded contour precip/ptype products<br>
>> over a basemap/background (Jeff Chabot)<br>
>> 2. Re: GrADS help, plotting shaded contour precip/ptype products<br>
>> over a basemap/background (Jeffrey Duda)<br>
>> 3. Problem with skewt plot using WRF output (Cristian ...)<br>
>><br>
>><br>
>> ----------------------------------------------------------------------<br>
>><br>
>> Message: 1<br>
>> Date: Mon, 8 Nov 2010 13:51:33 -0500<br>
>> From: Jeff Chabot <<a href="mailto:jsc219@gmail.com">jsc219@gmail.com</a>><br>
>> Subject: [gradsusr] GrADS help, plotting shaded contour precip/ptype<br>
>> products over a basemap/background<br>
>> To: <a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
>> Message-ID:<br>
>> <AANLkTinKn=<a href="mailto:eTr6myE00LdycRWRxdyy5e8hpBtf%2BFgLPz@mail.gmail.com">eTr6myE00LdycRWRxdyy5e8hpBtf+FgLPz@mail.gmail.com</a><<a href="mailto:eTr6myE00LdycRWRxdyy5e8hpBtf%252BFgLPz@mail.gmail.com">eTr6myE00LdycRWRxdyy5e8hpBtf%2BFgLPz@mail.gmail.com</a>><br>
>> ><br>
>> Content-Type: text/plain; charset="iso-8859-1"<br>
>><br>
>><br>
>> Dear GrADS Users,<br>
>> I have been trying to apply a basemap to my precip/ptype shaded contour<br>
>> products on my personal website for some time now. I have brainstormed<br>
>> and<br>
>> searched the GrADS online documents and have come up with ideas, but they<br>
>> have never worked. My latest idea was to use Google Maps terrain data as<br>
>> a<br>
>> background to the map. I thought it was a good idea, but I discovered<br>
>> that<br>
>> I still cannot resolve the issue in that precip = 0 shows either all white<br>
>> or all black depending on what I choose for the background, covering up my<br>
>> basemap. I really thought the following command would have resolved that<br>
>> issue, but it did not:<br>
>><br>
>> printim precip_NE_1.png png x570 y435 -b NE.png -t 0 Versions:GrADS: 2.0a9<br>
>> OS: Fedora 13<br>
>> Sample: <a href="http://jeffsweatherservice.com/grads/nam/ptype_NE_2.png" target="_blank">http://jeffsweatherservice.com/grads/nam/ptype_NE_2.png</a><br>
>><br>
>> I know this can be done because of the example:<br>
>><br>
>> <a href="http://wxmaps.org/pix/ez.east.html" target="_blank">http://wxmaps.org/pix/ez.east.html</a><br>
>><br>
>> Any assistance here would be much appreciated.Sincerely,<br>
>><br>
>><br>
>> Jeff Chabot<br>
>> Email: <a href="mailto:jsc219@gmail.com">jsc219@gmail.com</a><br>
>> Web: <a href="http://jeffsweatherservice.com" target="_blank">http://jeffsweatherservice.com</a><br>
>> -------------- next part --------------<br>
>> An HTML attachment was scrubbed...<br>
>> URL:<br>
>> <a href="http://gradsusr.org/pipermail/gradsusr/attachments/20101108/93f97642/attachment-0001.html" target="_blank">http://gradsusr.org/pipermail/gradsusr/attachments/20101108/93f97642/attachment-0001.html</a><br>
>><br>
>> ------------------------------<br>
>><br>
>> Message: 2<br>
>> Date: Mon, 8 Nov 2010 12:57:30 -0600<br>
>> From: Jeffrey Duda <<a href="mailto:jdduda@iastate.edu">jdduda@iastate.edu</a>><br>
>> Subject: Re: [gradsusr] GrADS help, plotting shaded contour<br>
>> precip/ptype products over a basemap/background<br>
>> To: GrADS Users Forum <<a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a>><br>
>> Message-ID:<br>
>> <<a href="mailto:AANLkTikPmH-Og0-6uNAi4tpUMCoZbAZ3YTqm_HXkndpV@mail.gmail.com">AANLkTikPmH-Og0-6uNAi4tpUMCoZbAZ3YTqm_HXkndpV@mail.gmail.com</a>><br>
>> Content-Type: text/plain; charset="iso-8859-1"<br>
>><br>
>><br>
>> Jeff,<br>
>> What color table are you using? Are you sure you have a #0 color? I know<br>
>> that's supposed to be the background, but perhaps you are overwriting it<br>
>> with some other color that has the same color as your background but a<br>
>> different color number.<br>
>><br>
>> Jeff Duda<br>
>><br>
><br>
><br>
><br>
>><br>
>> On Mon, Nov 8, 2010 at 12:51 PM, Jeff Chabot <<a href="mailto:jsc219@gmail.com">jsc219@gmail.com</a>> wrote:<br>
>><br>
>> > Dear GrADS Users,<br>
>> > I have been trying to apply a basemap to my precip/ptype shaded contour<br>
>> > products on my personal website for some time now. I have brainstormed<br>
>> and<br>
>> > searched the GrADS online documents and have come up with ideas, but<br>
>> they<br>
>> > have never worked. My latest idea was to use Google Maps terrain data<br>
>> as a<br>
>> > background to the map. I thought it was a good idea, but I discovered<br>
>> that<br>
>> > I still cannot resolve the issue in that precip = 0 shows either all<br>
>> white<br>
>> > or all black depending on what I choose for the background, covering up<br>
>> my<br>
>> > basemap. I really thought the following command would have resolved<br>
>> that<br>
>> > issue, but it did not:<br>
>> ><br>
>> > printim precip_NE_1.png png x570 y435 -b NE.png -t 0 Versions:GrADS:<br>
>> 2.0a9<br>
>> > OS: Fedora 13<br>
>> > Sample: <a href="http://jeffsweatherservice.com/grads/nam/ptype_NE_2.png" target="_blank">http://jeffsweatherservice.com/grads/nam/ptype_NE_2.png</a><br>
>> ><br>
>> > I know this can be done because of the example:<br>
>> ><br>
>> > <a href="http://wxmaps.org/pix/ez.east.html" target="_blank">http://wxmaps.org/pix/ez.east.html</a><br>
>> ><br>
>> > Any assistance here would be much appreciated. Sincerely,<br>
>> ><br>
>> ><br>
>> > Jeff Chabot<br>
>> > Email: <a href="mailto:jsc219@gmail.com">jsc219@gmail.com</a><br>
>> > Web: <a href="http://jeffsweatherservice.com" target="_blank">http://jeffsweatherservice.com</a><br>
>> ><br>
>> ><br>
>> > _______________________________________________<br>
>> > gradsusr mailing list<br>
>> > <a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
>> > <a href="http://gradsusr.org/mailman/listinfo/gradsusr" target="_blank">http://gradsusr.org/mailman/listinfo/gradsusr</a><br>
>> ><br>
>> ><br>
>><br>
>><br>
>> --<br>
>> Jeff Duda<br>
>> Iowa State University<br>
>> Meteorology Graduate Student<br>
>> 3134 Agronomy Hall<br>
>> <a href="http://www.meteor.iastate.edu/%7Ejdduda" target="_blank">www.meteor.iastate.edu/~jdduda</a> <<a href="http://www.meteor.iastate.edu/%7Ejdduda" target="_blank">http://www.meteor.iastate.edu/%7Ejdduda</a>><br>
>> -------------- next part --------------<br>
>> An HTML attachment was scrubbed...<br>
>> URL:<br>
>> <a href="http://gradsusr.org/pipermail/gradsusr/attachments/20101108/159827b2/attachment-0001.html" target="_blank">http://gradsusr.org/pipermail/gradsusr/attachments/20101108/159827b2/attachment-0001.html</a><br>
>><br>
>> ------------------------------<br>
>><br>
>> Message: 3<br>
>> Date: Mon, 8 Nov 2010 16:23:09 -0300<br>
>> From: "Cristian ..." <<a href="mailto:galenso85@hotmail.com">galenso85@hotmail.com</a>><br>
>> Subject: [gradsusr] Problem with skewt plot using WRF output<br>
>> To: <<a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a>><br>
>> Message-ID: <SNT112-W63DF5D09C452E1C6FB44FBF4F0@phx.gbl><br>
>> Content-Type: text/plain; charset="iso-8859-1"<br>
>><br>
>><br>
>> Hello,<br>
>><br>
>> My name is Cristian. I am having some troubles to perform some skewt<br>
>> graphics using WRF output.<br>
>> I have created a binary and .ctl archives using ARWpost.<br>
>> The file is OK, and no troubles are expected.<br>
>><br>
>> The problem starts when I want to create a skewt graphic.<br>
>> The temperature and dewpoint profiles are OK, but no information about<br>
>> wind is shown. Neither the Hodograpf, nor the wind speed and direction<br>
>> profiles are shown.<br>
>><br>
>> No message of an error apears, and i have no idea what is happening.<br>
>> Even the undef value is th same as the .ctl archive.<br>
>><br>
>> I am using Grads version 1.a8, in UBUNTU 9.04.<br>
>><br>
>> Here is the script. I only changed parameters at the top of the archive.<br>
>><br>
>> I really appreciate your help.<br>
>><br>
>> Thanks<br>
>><br>
>><br>
>> -----------------------------------------------------------------------------------------------------------------------------------------<br>
>> * Script para construir radiosondeos<br>
>><br>
>> 'reinit'<br>
>> 'run /home/cristian/Escritorio/Experimentos/lib/jaecol'<br>
>> 'set display color white'<br>
>> 'reset'<br>
>><br>
>> * Defino el archivo a utilizar<br>
>><br>
>> 'open 20081004_YSU.ctl'<br>
>> *path=''<br>
>><br>
>> * Latitud Longitud y tiempo<br>
>><br>
>> lat=-46<br>
>> lon=-68<br>
>> time=12z04Oct2008<br>
>><br>
>><br>
>> 'set time 'time<br>
>> 'set lon 'lon<br>
>> 'set lat 'lat<br>
>><br>
>> *Entre que niveles hace el perfil hacerlo entre 1000 y 600<br>
>><br>
>> levmax=925<br>
>> levmin=500<br>
>><br>
>> 'set lev 'levmax ' ' levmin<br>
>><br>
>> *calculo Td a partir de HR, P y T<br>
>> *'define et=(-2937.4/(t))-(4.9283*log10(t))+22.5518'<br>
>> *'define es=pow(10,et)*10'<br>
>> *'define e=es*0.01*hr'<br>
>> *'define q=0.622*e/(lev-e)'<br>
>> *'define aux1=q*1000*(1.+0.81*q)'<br>
>> *'define aux2=0.2854*(1.-0.28*q)'<br>
>> *'define aux3=(1000/lev)'<br>
>> *'define td=243.5/((17.67/(log(e/6.11)))-1)'<br>
>> *'define tc=t-273'<br>
>><br>
>> 'set ylopts 1 5 0.18'<br>
>><br>
>> plotskew(tc,td,u,v)<br>
>><br>
>> 'q dims'<br>
>> line1=sublin(result,4)<br>
>> line2=sublin(result,5)<br>
>> itime1=subwrd(line1,6)<br>
>> itit=substr(itime1,1,6)<br>
>> itime1=subwrd(line2,6)<br>
>> itime=substr(itime1,1,12)<br>
>> d=substr(itime1,4,2)<br>
>> m=substr(itime1,6,3)<br>
>> y=substr(itime1,9,4)<br>
>> h=substr(itime1,1,2)<br>
>> lati=-1*lat<br>
>> loni=-1*lon<br>
>><br>
>> 'draw title Perfil vertical de T y Td en 'lati'S y 'loni'W\'itime<br>
>> 'printim 'd'.'h'sondeo.png png x600 y800 white'<br>
>><br>
>> function plotskew(sndtemp,snddewp,u,v)<br>
>><br>
>> *El valor Undef que se coloca ac? debe coincidir con el valor undef del<br>
>> CTL, sin? hay problemas<br>
>> *en la construcci?n de la hod?grafa y del perfil vertical de viento.<br>
>> undef=1.e30<br>
>> *************************************************************************<br>
>> *<br>
>> * GrADS Script to Plot a SkewT/LogP Diagram<br>
>> *<br>
>> * Bob Hart<br>
>> * Penn State University / Dept of Meteorology<br>
>> * Last Update: January 23, 2001<br>
>> *<br>
>> * Recent Changes:<br>
>> *<br>
>> * 01/23/01 - Fixed a small bug in the theta-e calculation.<br>
>> * Errors averaged 0.5-3K. Thank you George Bryan.<br>
>> *<br>
>> * 11/10/99 - Change in calculation method for CAPE/CIN. Trapezoid<br>
>> * integration method is now used. Speeds up execution<br>
>> * by 25%, and increases accuracy by 5-10%.<br>
>> *<br>
>> * 10/18/99 - Minor glitch fixed that occasionally caused crash.<br>
>> *<br>
>> * 8/26/99 - Datasets with missing data can now be used.<br>
>> *<br>
>> * Features:<br>
>> * - All features of standard skewt/logp plot<br>
>> * - RH sounding<br>
>> * - LCL location<br>
>> * - Parcel trajectory for both sfc based convection and elevated from<br>
>> * most unstable level (highest theta-e level reported)<br>
>> * - Stability indices and precipitable water calculations<br>
>> * - CAPE & CIN Calculations<br>
>> * - Wind Profile<br>
>> * - Hodograph / Hodograph scaling<br>
>> * - Helicity and SR Helicity Calculations and Display<br>
>> * - Color aspects of output<br>
>> * - Line Thickness, style aspects of output<br>
>> * - Can be run in either PORTRAIT or LANDSCAPE mode.<br>
>> *<br>
>> * There are numerous tunable parameters below to change the structure<br>
>> * and output for the diagram.<br>
>> *<br>
>> * Function Arguments:<br>
>> * sndtemp - temperature data (Celsius) as a function of pressure<br>
>> * snddewp - dewpoint data (Celsius) as a function of pressure<br>
>> * sndspd - wind speed data (knots) as a function of pressure<br>
>> * snddir - wind direction data as a function of pressure<br>
>> *<br>
>> * Use '-1' for any of the above 4 arguments to indicate that you<br>
>> * are not passing that variable. The appropriate options will<br>
>> * be ignored based on your specifying '-1' for that variable.<br>
>> *<br>
>> * NOTE: Make sure to set the vertical range of the plot before running.<br>
>> * I.e., "SET LEV 1050 150", for example. This does not have to<br>
>> * be limited to the pressure range of your data.<br>
>> *<br>
>> * Labelling: Pressure/Height is labelled along left side. Temperature is<br>
>> * labelled along bottom. Mixing ratio is labelled along right<br>
>> * side/top.<br>
>> *<br>
>> *<br>
>> * PROBLEMS: First check out the web page for the script (which also<br>
>> * has a link to a FAQ with answers to many common questions<br>
>> * about using the script):<br>
>> * <a href="http://www.ems.psu.edu/%7Ehart/skew.html" target="_blank">http://www.ems.psu.edu/~hart/skew.html</a><<a href="http://www.ems.psu.edu/%7Ehart/skew.html" target="_blank">http://www.ems.psu.edu/%7Ehart/skew.html</a>><br>
>> *<br>
>> * Please send any further problems, comments, or suggestions to<br>
>> * <<a href="mailto:hart@ems.psu.edu">hart@ems.psu.edu</a>><br>
>> *<br>
>> * ACKNOWLEDGMENTS: Thanks go to the innumerable users who have helped<br>
>> * fine tune the script from the horrible mess from which it began.<br>
>> * In particular, thanks go out to Steve Lord (NCEP), Mike Fiorino (ECMWF),<br>
>> * George Bryan (PSU), Davide Sacchetti (CMIRL), and Enrico Minguzzi<br>
>> (CMIRL).<br>
>> *<br>
>> **************************************************************************<br>
>> * !!!!! BEGINNING OF USER-SPECIFIED OPTIONS !!!!!!<br>
>> **************************************************************************<br>
>> *<br>
>> * --------------------- Initialization options ----------------------<br>
>> *<br>
>> * ClrScrn = Whether to clear the screen before drawing diagram<br>
>> * [1 = yes, 0 = no]<br>
>><br>
>> ClrScrn = 1<br>
>><br>
>> *<br>
>> * ------------------- Define Skew-T Diagram Shape/Slope-----------------<br>
>> *<br>
>> * (P1,T1) = Pres, Temp of some point on left-most side<br>
>> * (P2,T2) = Pres, Temp of some point on right-most side<br>
>> * (P3,T3) = Pres, Temp of some point in diagram which is mid-point<br>
>> * in the horizontal between 1 and 2.<br>
>> *<br>
>> * P1, P2, P3 are in mb ; T1, T2, T3 are in Celsius<br>
>> *<br>
>> * These define the SLOPE and WIDTH of the diagram as you see it but DO NOT<br>
>> * DEFINE THE HEIGHT of the diagram as you see it. In other words,<br>
>> * 1 and 2 do NOT necessarily need to be at the bottom of the diagram and<br>
>> * 3 does NOT necessarily need to be at the top. THE VERTICAL PRESSURE<br>
>> * RANGE OF THE SKEWT AS YOU SEE IT IS DETERMINED BY YOUR 'SET Z ...'<br>
>> * COMMAND OR THE 'SET LEV ...' COMMAND BEFORE RUNNING THIS SCRIPT.<br>
>> *<br>
>> * _______________________<br>
>> * | |<br>
>> * | |<br>
>> * | 3 |<br>
>> * | |<br>
>> * | |<br>
>> * | |<br>
>> * | |<br>
>> * | |<br>
>> * | |<br>
>> * | |<br>
>> * | |<br>
>> * |1 2|<br>
>> * | |<br>
>> * |_______________________|<br>
>> *<br>
>> *<br>
>> * A good set of defining points are given below. Feel free<br>
>> * to experiment with variations.<br>
>><br>
>><br>
>> *P1 = 1000<br>
>> *T1 = -40<br>
>><br>
>> *P2 = 1000<br>
>> *T2 = 40<br>
>><br>
>> *P3 = 200<br>
>> *T3 = 0<br>
>><br>
>> * Another good set of defining points suggested by Juan Ruiz (Emagrama)<br>
>> * are:<br>
>> *<br>
>> P1 = 1000<br>
>> T1 = -80<br>
>> *<br>
>> P2 = 1000<br>
>> T2 = 40<br>
>> *<br>
>> P3 = 500<br>
>> T3 = -20<br>
>><br>
>> * ------------------- Contour Intervals / Levels<br>
>> --------------------------<br>
>> *<br>
>> * All variables below are contour intervals/levels for diagram<br>
>> *<br>
>> * Thetaint = interval for potential temperature lines<br>
>> * Thetwint = interval for moist pseudo adiabats<br>
>> * tempint = interval for temperature lines<br>
>> * wsclevs = contour LEVELS for mixing ratio lines<br>
>> *<br>
>> *<br>
>> thetaint= 20<br>
>> thetwint= 10<br>
>> tempint = 10<br>
>> wsclevs = ".1 .5 1 3 6 10 15 20 25 30"<br>
>> *<br>
>> *<br>
>> * ------------------------ Output Options --------------------------------<br>
>> *<br>
>> * All variables below are logical .. 1=yes, 0=no, unless otherwise<br>
>> * specified.<br>
>> *<br>
>> * DrawBarb = Draw wind barbs along right side of plot<br>
>> * DrawThet = Draw dry adiabats<br>
>> * DrawThtw = Draw moist pseudo-adiabats<br>
>> * DrawTemp = Draw temperature lines<br>
>> * DrawMix = Draw mixing ratio lines<br>
>> * DrawTSnd = Draw temperature sounding<br>
>> * DrawDSnd = Draw dewpoint sounding<br>
>> * DrawRH = Draw relative humidity sounding<br>
>> * DrawPrcl = Draw parcel path from surface upward<br>
>> * DrawPMax = Draw parcel path from most unstable level upward<br>
>> * DrawIndx = Display stability indices & CAPE<br>
>> * DrawHeli = Calculate and display absolute and storm-relative helicity<br>
>> * DrawHodo = Draw hodograph<br>
>> * DrawPLev = Draw Pressure Levels<br>
>> * DrawZLev = Draw height levels and lines<br>
>> * 0 = no lines<br>
>> * 1 = above ground level (AGL)<br>
>> * 2 = above sea level (ASL)<br>
>> * DrawZSTD = Draw Height levels using standard atm lapse rate<br>
>> * LblAxes = Label the x,y axes (temperature, pressure,mixing ratio)<br>
>> *<br>
>> * ThtwStop = Pressure level at which to stop drawing Theta-w lines<br>
>> * MixStop = Pressure level at which to stop drawing Mixratio lines<br>
>><br>
>> DrawBarb= 1<br>
>> DrawThet= 1<br>
>> DrawThtw= 1<br>
>> DrawTemp= 1<br>
>> DrawMix = 1<br>
>> DrawTSnd= 1<br>
>> DrawDSnd= 1<br>
>> DrawRH = 0<br>
>> DrawPrcl= 1<br>
>> DrawPMax= 1<br>
>> DrawIndx= 1<br>
>> DrawHeli= 0<br>
>> DrawHodo= 1<br>
>> DrawPLev= 1<br>
>> DrawZLev= 0<br>
>> DrawZSTD= 0<br>
>> LblAxes = 1<br>
>><br>
>> ThtwStop = 200<br>
>> MixStop = 400<br>
>><br>
>><br>
>> *<br>
>> * ----------------- Sounding Geography options ------------------------<br>
>> *<br>
>> * SfcElev = Elevation above sea-level (meters) of lowest level reported<br>
>> * in sounding. Used only if DrawZLev = 2<br>
>><br>
>> SfcElev = 0<br>
>><br>
>><br>
>> *<br>
>> * ------------------ Thermodynamic Index Options --------------------<br>
>> *<br>
>> * All variables here are in inches. Use -1 for the default values.<br>
>> *<br>
>> * Text1XC = X-location of midpoint of K,TT,PW output box<br>
>> * Text1YC = Y-location of midpoint of K,TT,PW output box<br>
>> * Text2XC = X-Location of midpoint of surface indices output box<br>
>> * Text2YC = Y-location of midpoint of surface indices output box<br>
>> * Text3XC = X-Location of midpoint of most unstable level-based indices<br>
>> * output box<br>
>> * Text3YC = Y-location of midpoint of most unstable level-based indices<br>
>> * output box<br>
>><br>
>> Text1XC = -1<br>
>> Text1YC = -1<br>
>> Text2XC = -1<br>
>> Text2YC = -1<br>
>> Text3XC = -1<br>
>> Text3YC = -1<br>
>><br>
>> *<br>
>> * ----------------- Wind Barb Profile Options ----------------------------<br>
>> *<br>
>> * All variables here are in units of inches, unless otherwise specified<br>
>> *<br>
>> * barbint = Interval for plotting barbs (in units of levels)<br>
>> * poleloc = X-Location of profile. Choose -1 for the default.<br>
>> * polelen = Length of wind-barb pole<br>
>> * Len05 = Length of each 5-knot barb<br>
>> * Len10 = Length of each 10-knot barb<br>
>> * Len50 = Length of each 50-knot flag<br>
>> * Wid50 = Width of base of 50-knot flag<br>
>> * Spac50 = Spacing between 50-knot flag and next ror occurred on<br>
>> libarb/flag<br>
>> * Spac10 = Spacing between 10-knot flag and next flag<br>
>> * Spac05 = Spacing between 5-knot flag and next flag<br>
>> * Flagbase= Draw flagbase (filled circle) for each windbarb [1=yes, 0<br>
>> =no]<br>
>> * Fill50 = Solid-fill 50-knot flag [1=yes, 0=no]<br>
>> * barbline= Draw a vertical line connecting all the wind barbs [1=yes,<br>
>> 0=no]<br>
>> *<br>
>> barbint = 1<br>
>> poleloc = -1<br>
>> polelen = 0.35<br>
>> len05 = 0.07<br>
>> len10 = 0.15<br>
>> len50 = 0.15<br>
>> wid50 = 0.06<br>
>> spac50 = 0.07<br>
>> spac10 = 0.05<br>
>> spac05 = 0.05<br>
>> Fill50 = 1<br>
>> flagbase= 1<br>
>> barbline= 1<br>
>><br>
>> *<br>
>> *<br>
>> *---------------- Hodograph Options -------------------------------------<br>
>> *<br>
>> * All variables here are in units of inches, unless otherwise specified<br>
>> *<br>
>> * HodXcent= x-location of hodograph center. Use -1 for default location.<br>
>> * HodYcent= y-location of hodograph center. Use -1 for default location.<br>
>> * HodSize = Size of hodograph in inches<br>
>> * NumRing = Number of rings to place in hodograph (must be at least 1)<br>
>> * HodRing = Wind speed increment of each hodograph ring<br>
>> * HodoDep = Depth (above lowest level in mb) of end of hodograph trace<br>
>> * TickInt = Interval (in kts) at which tick marks are drawn along the axes<br>
>> * Use 0 for no tick marks.<br>
>> * TickSize= Size of tick mark in inches<br>
>> * Text4XC = X-location of midpoint of hodograph text output. Use -1 for<br>
>> default.<br>
>> * Text4YC = Y-location of midpoint of hodograph text output. Use -1 for<br>
>> default.<br>
>><br>
>> HodXcent= 6<br>
>> HodYcent= 9<br>
>> HodSize = 2<br>
>> NumRing = 3<br>
>> HodRing = 10<br>
>> HodoDep = 500<br>
>> TickInt = 5<br>
>> TickSize= 0.1<br>
>> Text4XC = -1<br>
>> Text4YC = -1<br>
>><br>
>> *--------------- Helicity Options ---------------------------------------<br>
>> *<br>
>> * MeanVTop = Top pressure level (mb) of mean-wind calculation<br>
>> * MeanVBot = Bottom pressure level (mb) of mean-wind calculation<br>
>> * HelicDep = Depth in mb (above ground) of helicity integration<br>
>> * StormMot = Type of storm motion estimation scheme. Use following:<br>
>> * 0 = No departure from mean wind.<br>
>> * 1 = Davies-Jones (1990) approach<br>
>> * FillArrw = Whether to fill the arrowhead of the storm motion vector<br>
>> * [1 = yes, 0 = no]<br>
>><br>
>> MeanVTop= 300<br>
>> MeanVBot= 850<br>
>> HelicDep= 300<br>
>> StormMot= 0<br>
>> FillArrw= 1<br>
>><br>
>> *<br>
>> *---------------- Color Options ------------------------------------------<br>
>> *<br>
>> * ThetCol = Color of dry adiabats<br>
>> * TempCol = Color of temperature lines<br>
>> * MixCol = Color of mixing ratio lines<br>
>> * ThtwCol = Color of moist adiabats<br>
>> * TSndCol = Color of Temperature Sounding<br>
>> * DSndCol = Color of Dewpoint Sounding<br>
>> * RHCol = Color of RH Sounding<br>
>> * PrclCol = Color of parcel trace<br>
>> * BarbCol = Color of wind barbs (choose -1 for color according to speed)<br>
>> * HodoCol = Color of hodograph trace<br>
>><br>
>> ThetCol = 23<br>
>> TempCol = 79<br>
>> MixCol = 38<br>
>> ThtwCol = 39<br>
>> TSndCol = 29<br>
>> DSndCol = 49<br>
>> RHCol = 3<br>
>> PrclCol = 1<br>
>> BarbCol = 1<br>
>> HodoCol = 2<br>
>><br>
>> *<br>
>> *-------------------- Line Style Options<br>
>> ------------------------------------<br>
>> *<br>
>> * GrADS Styles: 1=solid;2=long dash;3=short dash;4=long,short dashed;<br>
>> * 5=dotted;6=dot dash;7=dot dot dash<br>
>> *<br>
>> * ThetLine = Line Style of dry adiabats<br>
>> * TempLine = Line Style of temperature lines<br>
>> * MixLine = Line Style of mixing ratio lines<br>
>> * ThtwLine = Line Style of moist adiabats<br>
>> * TSndLine = Line Style of Temperature Sounding<br>
>> * DSndLine = Line Style of Dewpoint Sounding<br>
>> * RHLine = Line Style of RH sounding<br>
>> * PrclLine = Line Style of parcel trace<br>
>> * HodoLine = Line Style of hodograph trace<br>
>> *<br>
>><br>
>> ThetLine = 1<br>
>> TempLine = 1<br>
>> MixLine = 5<br>
>> ThtwLine = 2<br>
>> TSndLine = 1<br>
>> DSndLine = 2<br>
>> RHLine = 1<br>
>> PrclLine = 3<br>
>> HodoLine = 1<br>
>><br>
>> *<br>
>> *------------------- Line Thickness<br>
>> Options---------------------------------<br>
>> * GrADS Line Thickness: increases with increasing number. Influences<br>
>> * hardcopy output more strongly than screen output.<br>
>> *<br>
>> *<br>
>> * ThetThk = Line Thickness of dry adiabats<br>
>> * TempThk = Line Thickness of temperature lines<br>
>> * MixThk = Line Thickness of mixing ratio lines<br>
>> * ThtwThk = Line Thickness of moist adiabats<br>
>> * TSndThk = Line Thickness of temperature sounding<br>
>> * DSndThk = Line thickness of dewpoint sounding<br>
>> * RHThk = Line thickness of RH sounding<br>
>> * PrclThk = Line thickness of parcel trace<br>
>> * HodoThk = Line thickness of hodograph trace<br>
>> * BarbThk = Line thickness of wind barbs<br>
>><br>
>> ThetThk = 3<br>
>> TempThk = 1<br>
>> MixThk = 7<br>
>> ThtwThk = 3<br>
>> TSndThk = 12<br>
>> DSndThk = 12<br>
>> RHThk = 8<br>
>> PrclThk = 6<br>
>> HodoThk = 6<br>
>> BarbThk = 2<br>
>><br>
>> *<br>
>> *------------------- Data Point Marker Options<br>
>> -----------------------------<br>
>> * GrADS Marker Types: 0 = none ; 1 = cross ; 2 = open circle ;<br>
>> * 3 = closed circle ; 4 = open square ; 5 = closed<br>
>> square<br>
>> * 6 = X ; 7 = diamond ; 8 = triangle ; 9 = none<br>
>> * 10 = open circle with vertical line ; 11 = open oval<br>
>> *<br>
>> * TSndMrk = Mark type of data point marker for temperature sounding<br>
>> * DSndMrk = Mark type of data point marker for dewpoint sounding<br>
>> * RHMrk = Mark type of data point marker for relative humidity sounding<br>
>> * MrkSize = Mark size (inches) of each data marker<br>
>><br>
>> TSndMrk = 3<br>
>> DSndMrk = 3<br>
>> RHMrk = 0<br>
>> MrkSize = 0.1<br>
>><br>
>><br>
>> * !!!!! YOU SHOULD NOT NEED TO CHANGE ANYTHING BELOW HERE !!!!!<br>
>><br>
>> ****************************************************************************<br>
>><br>
>> *-------------------------------------------<br>
>> * grab user-specified environment dimensions<br>
>> *-------------------------------------------<br>
>><br>
>> "q dims"<br>
>> rec=sublin(result,2)<br>
>> _xtype=subwrd(rec,3)<br>
>> _xval=subwrd(rec,9)<br>
>> rec=sublin(result,3)<br>
>> _yval=subwrd(rec,9)<br>
>> _ytype=subwrd(rec,3)<br>
>> rec=sublin(result,4)<br>
>> _ptype=subwrd(rec,3)<br>
>> _pmax=subwrd(rec,6)<br>
>> _pmin=subwrd(rec,8)<br>
>> _zmin=subwrd(rec,11)<br>
>> _zmax=subwrd(rec,13)<br>
>> rec=sublin(result,5)<br>
>> _ttype=subwrd(rec,3)<br>
>> _tval=subwrd(rec,9)<br>
>><br>
>> "q file"<br>
>> rec=sublin(result,5)<br>
>> _zmaxfile=subwrd(rec,9)<br>
>><br>
>> *-------------------------------------------------------------<br>
>> * Check to ensure that dimensions are valid. Warn & exit if not.<br>
>> *--------------------------------------------------------------<br>
>><br>
>> dimrc=0<br>
>> If (_xtype != "fixed")<br>
>> say "X-Dims Error: Not fixed. Use 'set lon' or 'set x' to specify a<br>
>> value."<br>
>> dimrc=-1<br>
>> Endif<br>
>><br>
>> If (_ytype != "fixed")<br>
>> say "Y-Dims Error: Not fixed. Use 'set lat' or 'set y' to specify a<br>
>> value"<br>
>> dimrc=-1<br>
>> Endif<br>
>><br>
>> If (_ptype != "varying")<br>
>> say "Z-Dims Error: Not varying. Use 'set lev' or 'set z' to specify a<br>
>> range."<br>
>> dimrc=-1<br>
>> Endif<br>
>><br>
>> If (_ttype != "fixed")<br>
>> say "Time Error: Not fixed. Use 'set time' or 'set t' to specify a<br>
>> value"<br>
>> dimrc=-1<br>
>> Endif<br>
>><br>
>><br>
>> If (dimrc < 0)<br>
>> Return(-1)<br>
>> Endif<br>
>><br>
>><br>
>> *<br>
>> * A few global variables used in units conversion<br>
>> *<br>
>><br>
>> _pi=3.14159265<br>
>> _dtr=_pi/180<br>
>> _rtd=1/_dtr<br>
>> _ktm=0.514444<br>
>> _mtk=1/_ktm<br>
>><br>
>> * A few global constants used in thermo calcs<br>
>><br>
>> _C0=0.99999683<br>
>> _C1=-0.90826951/100<br>
>> _C2= 0.78736169/10000<br>
>> _C3=-0.61117958/1000000<br>
>> _C4= 0.43884187/pow(10,8)<br>
>> _C5=-0.29883885/pow(10,10)<br>
>> _C6= 0.21874425/pow(10,12)<br>
>> _C7=-0.17892321/pow(10,14)<br>
>> _C8= 0.11112018/pow(10,16)<br>
>> _C9=-0.30994571/pow(10,19)<br>
>><br>
>> *Calculo la direcci?n e intensidad del viento usando las funciones GetWdir<br>
>> y GetWsp<br>
>><br>
>><br>
>><br>
>><br>
>><br>
>><br>
>> * A pressure array of power calculations which should be performed<br>
>> * only once to reduce execution time.<br>
>><br>
>> zz=1100<br>
>> while (zz > 10)<br>
>> subscr=zz/10<br>
>> _powpres.subscr=pow(zz,0.286)<br>
>> zz=zz-10<br>
>> endwhile<br>
>><br>
>> *<br>
>> * Turn off options not available due to user data limitations<br>
>> *<br>
>><br>
>> If (ClrScrn = 1)<br>
>> "clear"<br>
>> Endif<br>
>><br>
>> If (sndspd = -1 | snddir = -1)<br>
>> DrawBarb = 0<br>
>> DrawHodo = 0<br>
>> DrawHeli = 0<br>
>> Endif<br>
>><br>
>> If (snddewp = -1)<br>
>> DrawDSnd = 0<br>
>> DrawRH = 0<br>
>> DrawPrcl = 0<br>
>> DrawPMax = 0<br>
>> DrawIndx = 0<br>
>> Endif<br>
>><br>
>> If (sndtemp = -1)<br>
>> DrawTSnd = 0<br>
>> DrawRH = 0<br>
>> DrawPrcl = 0<br>
>> DrawPMax = 0<br>
>> DrawIndx = 0<br>
>> DrawZLev = 0<br>
>> Endif<br>
>><br>
>> If (NumRing < 1)<br>
>> DrawHodo = 0<br>
>> Endif<br>
>><br>
>> "q gxinfo"<br>
>> rec=sublin(result,2)<br>
>> xsize=subwrd(rec,4)<br>
>><br>
>> If (xsize = 11)<br>
>> PageType = "Landscape"<br>
>> Else<br>
>> PageType = "Portrait"<br>
>> Endif<br>
>><br>
>> *------------------------------------------------------<br>
>> * calculate constants determining slope/shape of diagram<br>
>> * based on temp/pressure values given by user<br>
>> *-------------------------------------------------------<br>
>><br>
>> "set x 1"<br>
>> "set y 1"<br>
>> "set z 1"<br>
>> "set t 1"<br>
>> _m1=(T1+T2-2*T3)/(2*log10(P2/P3))<br>
>> _m2=(T2-T3-_m1*log10(P2/P3))/50<br>
>> _m3=(T1-_m1*log10(P1))<br>
>><br>
>> "set z "_zmin" "_zmax<br>
>> "set zlog on"<br>
>> "set xlab off"<br>
>><br>
>> *-------------------------------------------------<br>
>> * perform coordinate transformation to Skew-T/LogP<br>
>> *-------------------------------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set t "_tval<br>
>> "define tempx=("sndtemp"-"_m1"*log10(lev)-"_m3")/"_m2<br>
>> "define dewpx=("snddewp"-"_m1"*log10(lev)-"_m3")/"_m2<br>
>><br>
>> If (PageType = "Portrait")<br>
>> "set parea 0.7 7 0.75 10"<br>
>> Else<br>
>> "set parea 0.7 6.5 0.5 8"<br>
>> Endif<br>
>><br>
>> "set axlim 0 100"<br>
>> "set lon 0 100"<br>
>> "set grid on 1 1"<br>
>><br>
>> "set z "_zmin " " _zmax<br>
>> "set lon 0 100"<br>
>> "set clevs -900"<br>
>> "set gxout contour"<br>
>><br>
>> *-------------------------------------<br>
>> * Draw pressure lines<br>
>> *-------------------------------------<br>
>><br>
>> If (DrawPLev = 0)<br>
>> "set ylab off"<br>
>> Else<br>
>> "set ylab on"<br>
>> "set ylopts 1 5 0.18"<br>
>> "set xlopts 1 3 0.18"<br>
>> Endif<br>
>><br>
>> "d lon"<br>
>><br>
>> *--------------------------------------<br>
>> * Determine corners of skewt/logp frame<br>
>> *--------------------------------------<br>
>><br>
>> "q w2xy 100 "_pmin<br>
>> rxloc=subwrd(result,3)<br>
>> tyloc=subwrd(result,6)<br>
>> "q w2xy 0 "_pmax<br>
>> lxloc=subwrd(result,3)<br>
>> byloc=subwrd(result,6)<br>
>><br>
>> If (DrawPLev = 1 & LblAxes = 1)<br>
>> "set strsiz 0.15"<br>
>> "set string 1 c 3 0"<br>
>> If (PageType = "Portrait")<br>
>> * "draw string 0.5 10.5 hPa."<br>
>> Else<br>
>> "draw string 0.5 8.35 hPa."<br>
>> Endif<br>
>> Endif<br>
>><br>
>> *---------------------------------------------------<br>
>> * Calculate & draw actual height lines using temp data<br>
>> *---------------------------------------------------<br>
>><br>
>> If (DrawZLev > 0)<br>
>> say "Calculating observed height levels from temp/pressure data."<br>
>> zz=1<br>
>> "set gxout stat"<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set t "_tval<br>
>> count=0<br>
>> while (zz < _zmax)<br>
>> "set z "zz<br>
>> pp.zz=subwrd(result,4)<br>
>> lpp.zz=log(pp.zz)<br>
>> "d "sndtemp<br>
>> rec=sublin(result,8)<br>
>> tt=subwrd(rec,4)<br>
>> if (tt > -900)<br>
>> tk=tt+273.15<br>
>> count=count+1<br>
>> zzm=zz-1<br>
>> If (count = 1)<br>
>> If (DrawZLev = 2)<br>
>> htlb="ASL"<br>
>> height.zz=SfcElev<br>
>> Else<br>
>> htlb="AGL"<br>
>> height.zz=0<br>
>> Endif<br>
>> sfcz=height.zz<br>
>> Else<br>
>><br>
>> DZ=29.2857*(lpp.zzm-lpp.zz)*(lpp.zz*tk+lpp.zzm*tkold)/(lpp.zz+lpp.zzm)<br>
>> height.zz=height.zzm+DZ<br>
>> highz=height.zz<br>
>> Endif<br>
>> else<br>
>> height.zz = -9999<br>
>> endif<br>
>> tkold=tk<br>
>> zz=zz+1<br>
>> endwhile<br>
>><br>
>> maxht=int(highz/1000)<br>
>> if (int(sfcz/1000) = sfcz/1000)<br>
>> minht=int(sfcz/1000)<br>
>> else<br>
>> minht=1+int(sfcz/1000)<br>
>> endif<br>
>><br>
>> ht=minht<br>
>> "set line 1 3 1"<br>
>> "set strsiz 0.10"<br>
>> "set string 1 l 3 0"<br>
>> while (ht <= maxht)<br>
>> zz=1<br>
>> while (height.zz/1000 <= ht)<br>
>> zz=zz+1<br>
>> endwhile<br>
>> zzm=zz-1<br>
>> PBelow=pp.zzm<br>
>> PAbove=pp.zz<br>
>> HBelow=height.zzm<br>
>> HAbove=height.zz<br>
>> DZ=HAbove-HBelow<br>
>> DP=PAbove-PBelow<br>
>> Del=ht*1000-HBelow<br>
>> Est=PBelow+Del*DP/DZ<br>
>> If (Est >= _pmin & Est <= _pmax)<br>
>> "q w2xy 1 " Est<br>
>> yloc=subwrd(result,6)<br>
>> "draw line " lxloc " " yloc " " rxloc " " yloc<br>
>> "draw string 0.22 "yloc-0.05" "ht<br>
>> Endif<br>
>> ht=ht+1<br>
>> endwhile<br>
>> "set strsiz 0.10"<br>
>> "set string 1"<br>
>> If (LblAxes = 1)<br>
>> If (PageType = "Portrait")<br>
>> "draw string 0.25 10.85 km"<br>
>> "draw string 0.25 10.75 "htlb<br>
>> "draw string 0.25 10.65 OBS"<br>
>> Else<br>
>> "draw string 0.25 8.35 km"<br>
>> "draw string 0.25 8.25 "htlb<br>
>> "draw string 0.25 8.15 OBS"<br>
>> Endif<br>
>> Endif<br>
>> Endif<br>
>><br>
>><br>
>> *---------------------------------------------------<br>
>> * Draw height levels (height above MSL using Std Atm)<br>
>> *---------------------------------------------------<br>
>><br>
>> If (DrawZSTD = 1)<br>
>> "set strsiz 0.10"<br>
>> minht=30.735*(1-pow(_pmax/1013.26,0.287))<br>
>> minht=int(minht+0.5)<br>
>> maxht=30.735*(1-pow(_pmin/1013.26,0.287))<br>
>> maxht=int(maxht)<br>
>> "set gxout stat"<br>
>> zcount=minht<br>
>> while (zcount <= maxht)<br>
>> plev=1013.26*pow((1-zcount/30.735),3.4843)<br>
>> "q w2xy 0 "plev<br>
>> yloc=subwrd(result,6)<br>
>> "draw string 0 "yloc-0.05" "zcount<br>
>> zcount=zcount+1<br>
>> endwhile<br>
>> "set strsiz 0.10"<br>
>> If (LblAxes = 1)<br>
>> If (PageType = "Portrait")<br>
>> "draw string 0 10.85 km"<br>
>> "draw string 0 10.75 ASL"<br>
>> "draw string 0 10.65 STD"<br>
>> Else<br>
>> "draw string 0 8.35 km"<br>
>> "draw string 0 8.25 ASL"<br>
>> "draw string 0 8.15 STD"<br>
>> Endif<br>
>> Endif<br>
>> Endif<br>
>><br>
>><br>
>> *-----------------------<br>
>> * Plot temperature lines<br>
>> *-----------------------<br>
>><br>
>> If (DrawTemp = 1)<br>
>> "set strsiz 0.1"<br>
>> "set z "_zmin " " _zmax<br>
>> "set line "TempCol " " TempLine " "TempThk<br>
>> "set string 1 c 3 0"<br>
>> "set gxout stat"<br>
>> maxtline=GetTemp(100,_pmax)<br>
>> mintline=GetTemp(0,_pmin)<br>
>><br>
>> maxtline=tempint*int(maxtline/tempint)<br>
>> mintline=tempint*int(mintline/tempint)<br>
>><br>
>> tloop=mintline<br>
>> While (tloop <= maxtline)<br>
>> Botxtemp=GetXLoc(tloop,_pmax)<br>
>> "q w2xy "Botxtemp " " _pmax<br>
>> Botxloc=subwrd(result,3)<br>
>> Botyloc=byloc<br>
>> Topxtemp=GetXLoc(tloop,_pmin)<br>
>> "q w2xy "Topxtemp " " _pmin<br>
>> Topxloc=subwrd(result,3)<br>
>> Topyloc=tyloc<br>
>> If (Botxtemp <= 100 | Topxtemp <= 100)<br>
>> If (Topxtemp > 100)<br>
>> Slope=(Topyloc-Botyloc)/(Topxtemp-Botxtemp)<br>
>> b=Topyloc-Slope*Topxtemp<br>
>> Topyloc=Slope*100+b<br>
>> Topxloc=rxloc<br>
>> Endif<br>
>> If (Botxtemp < 0)<br>
>> Slope=(Topyloc-Botyloc)/(Topxtemp-Botxtemp)<br>
>> b=Botyloc-Slope*Botxtemp<br>
>> Botyloc=b<br>
>> Botxloc=lxloc<br>
>> Else<br>
>> "set strsiz 0.15"<br>
>> "draw string " Botxloc-0.05 " " Botyloc-0.15 " " tloop<br>
>> Endif<br>
>> "draw line "Botxloc " " Botyloc " " Topxloc " " Topyloc<br>
>> Endif<br>
>> tloop=tloop+tempint<br>
>> EndWhile<br>
>> If (LblAxes = 1)<br>
>> "set strsiz 0.15"<br>
>> "set string 1 c"<br>
>> If (PageType = "Portrait")<br>
>> "draw string 4.0 0.35 Temperatura (`3.`0C)"<br>
>> "draw string 7.7 0.35 Viento m/s "<br>
>> Else<br>
>> "draw string 3.5 0.15 Temperatura (`3.`0C)"<br>
>> Endif<br>
>> Endif<br>
>> Endif<br>
>><br>
>><br>
>> *------------------<br>
>> * Plot dry adiabats<br>
>> *------------------<br>
>><br>
>> If (DrawThet = 1)<br>
>> temp=GetTemp(100,_pmin)<br>
>> maxtheta=GetThet2(temp,-100,_pmin)<br>
>> maxtheta=thetaint*int(maxtheta/thetaint)<br>
>> temp=GetTemp(0,_pmax)<br>
>> mintheta=GetThet2(temp,-100,_pmax)<br>
>> mintheta=thetaint*int(mintheta/thetaint)<br>
>><br>
>> "set lon 0 100"<br>
>> "set y 1"<br>
>> "set z 1"<br>
>> tloop=mintheta<br>
>> "set line "ThetCol" "ThetLine " "ThetThk<br>
>> While (tloop <= maxtheta)<br>
>> PTemp=LiftDry(tloop,1000,_pmin,1,_pmin,_pmax)<br>
>> tloop=tloop+thetaint<br>
>> Endwhile<br>
>> Endif<br>
>><br>
>> *------------------------<br>
>> * Plot mixing ratio lines<br>
>> *------------------------<br>
>><br>
>> If (DrawMix = 1)<br>
>> If (MixStop < _pmin)<br>
>> MixStop = _pmin<br>
>> Endif<br>
>> "set string 1 l"<br>
>> "set z "_zmin " " _zmax<br>
>> "set cint 1"<br>
>> "set line "MixCol" " MixLine " "MixThk<br>
>> cont = 1<br>
>> mloop=subwrd(wsclevs,1)<br>
>> count = 1<br>
>> While (cont = 1)<br>
>> BotCoef=log(mloop*_pmax/3801.66)<br>
>> BotTval=-245.5*BotCoef/(BotCoef-17.67)<br>
>> Botxtemp=GetXLoc(BotTval,_pmax)<br>
>> "q w2xy "Botxtemp " " _pmax<br>
>> Botxloc=subwrd(result,3)<br>
>> Botyloc=byloc<br>
>> TopCoef=log(mloop*MixStop/3801.66)<br>
>> TopTval=-245.5*TopCoef/(TopCoef-17.67)<br>
>> Topxtemp=GetXLoc(TopTval,MixStop)<br>
>> "q w2xy "Topxtemp " " MixStop<br>
>> Topxloc=subwrd(result,3)<br>
>> Topyloc=subwrd(result,6)<br>
>> "set string "MixCol" l 3"<br>
>> "set strsiz 0.09"<br>
>> If (Botxtemp <= 100 | Topxtemp <= 100)<br>
>> If (Topxtemp > 100)<br>
>> Slope=(Topyloc-Botyloc)/(Topxtemp-Botxtemp)<br>
>> b=Topyloc-Slope*Topxtemp<br>
>> Topyloc=Slope*100+b<br>
>> Topxloc=rxloc<br>
>> "draw string " Topxloc+0.05 " " Topyloc " " mloop<br>
>> Else<br>
>> "draw string " Topxloc " " Topyloc+0.1 " " mloop<br>
>> Endif<br>
>> If (Botxtemp < 0)<br>
>> Slope=(Topyloc-Botyloc)/(Topxtemp-Botxtemp)<br>
>> b=Botyloc-Slope*Botxtemp<br>
>> Botyloc=b<br>
>> Botxloc=lxloc<br>
>> Endif<br>
>> "draw line "Botxloc " " Botyloc " " Topxloc " " Topyloc<br>
>> Endif<br>
>> count=count+1<br>
>> mloop=subwrd(wsclevs,count)<br>
>> If (mloop = "" | count > 50)<br>
>> cont = 0<br>
>> Endif<br>
>> EndWhile<br>
>> If (LblAxes = 1)<br>
>> "set strsiz 0.15"<br>
>> "set string 1 c 3 90"<br>
>> If (PageType = "Portrait")<br>
>> * "draw string 7.40 4.75 Relacion de Mezcla (g/kg)"<br>
>> Else<br>
>> * "draw string 6.90 4.25 Mixing Ratio (g/kg)"<br>
>> Endif<br>
>> "set string 1 c 3 0"<br>
>> Endif<br>
>> Endif<br>
>><br>
>> *-----------------------------<br>
>> * Plot moist (pseudo) adiabats<br>
>> *-----------------------------<br>
>><br>
>> If (DrawThtw = 1)<br>
>> "set lon 0 100"<br>
>> "set y 1"<br>
>> "set z 1"<br>
>> "set gxout stat"<br>
>> tloop=80<br>
>> "set line "ThtwCol" "ThtwLine " "ThtwThk<br>
>> While (tloop > -80)<br>
>> PTemp=LiftWet(tloop,1000,ThtwStop,1,_pmin,_pmax)<br>
>> tloop=tloop-thetwint<br>
>> Endwhile<br>
>> Endif<br>
>><br>
>><br>
>> *-----------------------------------------------------<br>
>> * Plot transformed user-specified temperature sounding<br>
>> *-----------------------------------------------------<br>
>><br>
>> If (DrawTSnd = 1)<br>
>> say "Drawing temperature sounding."<br>
>> "set gxout line"<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set z "_zmin" "_zmax<br>
>> "set ccolor "TSndCol<br>
>> "set cstyle "TSndLine<br>
>> "set cmark "TSndMrk<br>
>> "set digsiz "MrkSize<br>
>> "set cthick "TSndThk<br>
>> "set missconn on"<br>
>> "d tempx"<br>
>> Endif<br>
>><br>
>> *---------------------------------------------------<br>
>> * Plot transformed user-specified dewpoint sounding<br>
>> *---------------------------------------------------<br>
>><br>
>> If (DrawDSnd = 1)<br>
>> say "Drawing dewpoint sounding."<br>
>> "set gxout line"<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set z "_zmin" "_zmax<br>
>> "set cmark "DSndMrk<br>
>> "set digsiz "MrkSize<br>
>> "set ccolor "DSndCol<br>
>> "set cstyle "DSndLine<br>
>> "set cthick "DSndThk<br>
>> "set missconn on"<br>
>> "d dewpx"<br>
>> Endif<br>
>><br>
>> *----------------------------------------<br>
>> * Determine lowest level of reported data<br>
>> *----------------------------------------<br>
>><br>
>> If (DrawTSnd = 1)<br>
>> field=sndtemp<br>
>> Else<br>
>> field=sndspd<br>
>> Endif<br>
>><br>
>> "set gxout stat"<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set t "_tval<br>
>> "set lev " _pmax " " _pmin<br>
>> "d maskout(lev,"field"+300)"<br>
>> rec=sublin(result,1)<br>
>> check=substr(rec,1,6)<br>
>> If (check = "Notice")<br>
>> rec=sublin(result,9)<br>
>> Else<br>
>> rec=sublin(result,8)<br>
>> Endif<br>
>> SfcPlev=subwrd(rec,5)<br>
>><br>
>> If (DrawTSnd = 1 & DrawDSnd = 1)<br>
>> "set lev "SfcPlev<br>
>> "d "sndtemp<br>
>> rec=sublin(result,8)<br>
>> Sfctemp=subwrd(rec,4)<br>
>> "d "snddewp<br>
>> rec=sublin(result,8)<br>
>> Sfcdewp=subwrd(rec,4)<br>
>> SfcThee=Thetae(Sfctemp,Sfcdewp,SfcPlev)<br>
>><br>
>> *------------------------------------------<br>
>> * Calculate temperature and pressure of LCL<br>
>> *------------------------------------------<br>
>><br>
>> TLcl=Templcl(Sfctemp,Sfcdewp)<br>
>> PLcl=Preslcl(Sfctemp,Sfcdewp,SfcPlev)<br>
>> Endif<br>
>><br>
>> *----------------------------------------------------------<br>
>> * Plot parcel path from surface to LCL and up moist adiabat<br>
>> *----------------------------------------------------------<br>
>><br>
>> If (DrawPrcl = 1)<br>
>> say "Drawing parcel path from surface upward."<br>
>> If (PageType = "Portrait")<br>
>> xloc=7.15<br>
>> Else<br>
>> xloc=6.65<br>
>> Endif<br>
>> "q w2xy 1 "PLcl<br>
>> rec=sublin(result,1)<br>
>> yloc=subwrd(rec,6)<br>
>> "set strsiz 0.1"<br>
>> If (PLcl < _pmax)<br>
>> "set string 1 l"<br>
>> "draw string "xloc" "yloc" LCL"<br>
>> "set line 1 1 1"<br>
>> "draw line "xloc-0.15" "yloc" "xloc-0.05" "yloc<br>
>> Endif<br>
>> "set lon 0 100"<br>
>> "set gxout stat"<br>
>> "set line "PrclCol" "PrclLine " " PrclThk<br>
>> PTemp=LiftDry(Sfctemp,SfcPlev,PLcl,1,_pmin,_pmax)<br>
>> Ptemp=LiftWet(TLcl,PLcl,_pmin,1,_pmin,_pmax)<br>
>> Endif<br>
>><br>
>> *-------------------------------------------------------<br>
>> * Determine level within lowest 250mb having highest<br>
>> * theta-e value<br>
>> *-------------------------------------------------------<br>
>><br>
>> If (DrawTSnd = 1 & DrawDSnd = 1)<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set t "_tval<br>
>> zz=1<br>
>> MaxThee=-999<br>
>> "set gxout stat"<br>
>> while (zz <= _zmax & pp > _pmax-250)<br>
>> "set z "zz<br>
>> pp=subwrd(result,4)<br>
>> "d "sndtemp<br>
>> rec=sublin(result,8)<br>
>> tt=subwrd(rec,4)<br>
>> "d "snddewp<br>
>> rec=sublin(result,8)<br>
>> dd=subwrd(rec,4)<br>
>> If (abs(tt) < 130 & abs(dd) < 130)<br>
>> Thee=Thetae(tt,dd,pp)<br>
>> If (Thee > MaxThee)<br>
>> MaxThee=Thee<br>
>> TMaxThee=tt<br>
>> DMaxThee=dd<br>
>> PMaxThee=pp<br>
>> Endif<br>
>> endif<br>
>> zz=zz+1<br>
>> Endwhile<br>
>> If (PMaxThee = SfcPlev-250)<br>
>> PMaxThee = SfcPlev<br>
>> Endif<br>
>> *------------------------------------------------------<br>
>> * Calculate temperature and pressure of LCL from highest<br>
>> * theta-e level<br>
>> *------------------------------------------------------<br>
>> If (SfcPlev != PMaxThee)<br>
>> TLclMax=Templcl(TMaxThee,DMaxThee)<br>
>> PLclMax=Preslcl(TMaxThee,DMaxThee,PMaxThee)<br>
>> Endif<br>
>> Endif<br>
>><br>
>> *----------------------------------------------------------<br>
>> * Plot parcel path from highest theta-e level to LCL and up<br>
>> * moist adiabat<br>
>> *----------------------------------------------------------<br>
>><br>
>> If (DrawPMax = 1 & SfcPlev != PMaxThee)<br>
>> say "Drawing parcel path from most unstable level upward."<br>
>> If (PageType = "Portrait")<br>
>> xloc=7.15<br>
>> Else<br>
>> xloc=6.65<br>
>> Endif<br>
>> "q w2xy 1 "PLclMax<br>
>> rec=sublin(result,1)<br>
>> yloc=subwrd(rec,6)<br>
>> "set strsiz 0.1"<br>
>> If (PLclMax < _pmax)<br>
>> "set string 1 l"<br>
>> "draw string "xloc" "yloc" LCL"<br>
>> "set line 1 1 1"<br>
>> "draw line "xloc-0.15" "yloc" "xloc-0.05" "yloc<br>
>> Endif<br>
>> "set lon 0 100"<br>
>> "set gxout stat"<br>
>> "set line "PrclCol" "PrclLine " " PrclThk<br>
>> PTemp=LiftDry(TMaxThee,PMaxThee,PLclMax,1,_pmin,_pmax)<br>
>> Ptemp=LiftWet(TLclMax,PLclMax,_pmin,1,_pmin,_pmax)<br>
>> Endif<br>
>><br>
>> *--------------------------------<br>
>> * Draw thermodynamic indices<br>
>> *--------------------------------<br>
>><br>
>> If (DrawIndx = 1)<br>
>> "set string 1 l"<br>
>> "set strsiz 0.10"<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set t "_tval<br>
>> say "Calculating precipitable water."<br>
>> pw=precipw(sndtemp,snddewp,_pmax,_pmin)<br>
>> say "Calculating thermodynamic indices."<br>
>> Temp850=interp(sndtemp,850)<br>
>> Temp700=interp(sndtemp,700)<br>
>> Temp500=interp(sndtemp,500)<br>
>> Dewp850=interp(snddewp,850)<br>
>> Dewp700=interp(snddewp,700)<br>
>> Dewp500=interp(snddewp,500)<br>
>> If (Temp850>-900 & Dewp850>-900 & Dewp700>-900 & Temp700>-900 &<br>
>> Temp500>-900)<br>
>> K=Temp850+Dewp850+Dewp700-Temp700-Temp500<br>
>> Else<br>
>> K=-999<br>
>> Endif<br>
>> If (Temp850 > -900 & Dewp850 > -900 & Temp500 > -900)<br>
>> tt=Temp850+Dewp850-2*Temp500<br>
>> Else<br>
>> tt=-999<br>
>> Endif<br>
>> Temp500V=virtual2(Temp500+273.15,Dewp500+273.15,500)-273.15<br>
>> PclTemp=LiftWet(TLcl,PLcl,500,0)<br>
>> PclTempV=virtual2(PclTemp+273.15,PclTemp+273.15,500)-273.15<br>
>> SLI=Temp500V-PclTempV<br>
>> rec=CAPE(TLcl,PLcl,200,sndtemp,snddewp)<br>
>> Pos=subwrd(rec,1)<br>
>> CIN=subwrd(rec,2)<br>
>><br>
>> If (SfcPlev != PMaxThee)<br>
>> PclTemp=LiftWet(TLclMax,PLclMax,500,0)<br>
>> PclTempV=virtual2(PclTemp+273.15,PclTemp+273.15,500)-273.15<br>
>> LIMax=Temp500V-PclTempV<br>
>> rec=CAPE(TLclMax,PLclMax,100,sndtemp,snddewp)<br>
>> PosMax=subwrd(rec,1)<br>
>> CINMax=subwrd(rec,2)<br>
>> Else<br>
>> LIMax=SLI<br>
>> PosMax=Pos<br>
>> CINMax=CIN<br>
>> MaxThee=SfcThee<br>
>> Endif<br>
>><br>
>> If (PageType = "Portrait")<br>
>> If (Text1XC = -1)<br>
>> Text1XC=rxloc-0.75<br>
>> Endif<br>
>> If (Text1YC = -1)<br>
>> Text1YC=tyloc-2.25<br>
>> Endif<br>
>> If (Text2XC = -1)<br>
>> Text2XC=rxloc-0.75<br>
>> Endif<br>
>> If (Text2YC = -1)<br>
>> Text2YC=tyloc-3.25<br>
>> Endif<br>
>> If (Text3XC = -1)<br>
>> Text3XC=rxloc-0.75<br>
>> Endif<br>
>> If (Text3YC = -1)<br>
>> Text3YC=tyloc-4.40<br>
>> Endif<br>
>> Else<br>
>> If (Text1XC = -1)<br>
>> Text1XC=rxloc+2.50<br>
>> Endif<br>
>> If (Text1YC = -1)<br>
>> Text1YC=tyloc-3.00<br>
>> Endif<br>
>> If (Text2XC = -1)<br>
>> Text2XC=rxloc+2.50<br>
>> Endif<br>
>> If (Text2YC = -1)<br>
>> Text2YC=tyloc-4.00<br>
>> Endif<br>
>> If (Text3XC = -1)<br>
>> Text3XC=rxloc+2.50<br>
>> Endif<br>
>> If (Text3YC = -1)<br>
>> Text3YC=tyloc-5.10<br>
>> Endif<br>
>> Endif<br>
>> "set string 1 l 3"<br>
>> "set line 0 1 3"<br>
>> "draw recf "Text1XC-0.75 " " Text1YC-0.40 " " Text1XC+0.75 " "<br>
>> Text1YC+0.25<br>
>> "set line 1 1 3"<br>
>> "draw rec "Text1XC-0.75 " " Text1YC-0.40 " " Text1XC+0.75 " "<br>
>> Text1YC+0.25<br>
>> "draw string "Text1XC-0.70 " " Text1YC+0.10" K"<br>
>> "draw string "Text1XC+0.25 " " Text1YC+0.10" " int(K)<br>
>> "draw string "Text1XC-0.70 " " Text1YC-0.10 " TT"<br>
>> "draw string "Text1XC+0.25 " " Text1YC-0.10 " " int(tt)<br>
>> "draw string "Text1XC-0.70 " " Text1YC-0.25 " AP(cm)"<br>
>> "draw string "Text1XC+0.25 " " Text1YC-0.25 " " int(pw*100)/100<br>
>> "set line 0 1 3"<br>
>> "draw recf "Text2XC-0.75 " " Text2YC-0.60 " " Text2XC+0.75 " "<br>
>> Text2YC+0.60<br>
>> "set line 1 1 3"<br>
>> "draw rec "Text2XC-0.75 " " Text2YC-0.60 " " Text2XC+0.75 " "<br>
>> Text2YC+0.60<br>
>> "draw string "Text2XC-0.35 " " Text2YC+0.50 " Superficie"<br>
>> "draw string "Text2XC-0.70 " " Text2YC+0.30 " Temp(`3.`0C)"<br>
>> "draw string "Text2XC+0.25 " " Text2YC+0.30 " " int(Sfctemp*10)/10<br>
>> "draw string "Text2XC-0.70 " " Text2YC+0.15 " Dewp(`3.`0C)"<br>
>> "draw string "Text2XC+0.25 " " Text2YC+0.15 " " int(Sfcdewp*10)/10<br>
>> "draw string "Text2XC-0.70 " " Text2YC " `3z`0`bE`n(K)"<br>
>> "draw string "Text2XC+0.25 " " Text2YC " " int(SfcThee)<br>
>> "draw string "Text2XC-0.70 " " Text2YC-0.15 " LI"<br>
>> "draw string "Text2XC+0.25 " " Text2YC-0.15 " " round(SLI)<br>
>> "draw string "Text2XC-0.70 " " Text2YC-0.30 " CAPE(J)"<br>
>> "draw string "Text2XC+0.25 " " Text2YC-0.30 " " int(Pos)<br>
>> "draw string "Text2XC-0.70 " " Text2YC-0.45 " CIN(J)"<br>
>> "draw string "Text2XC+0.25 " " Text2YC-0.45 " " int(CIN)<br>
>> "set line 0 1 3"<br>
>> "draw recf "Text3XC-0.75 " " Text3YC-0.55 " " Text3XC+0.75 " "<br>
>> Text3YC+0.55<br>
>> "set line 1 1 3"<br>
>> "draw rec "Text3XC-0.75 " " Text3YC-0.55 " " Text3XC+0.75 " "<br>
>> Text3YC+0.55<br>
>> "draw string "Text3XC-0.60 " " Text3YC+0.45 " Mas inestable"<br>
>> "draw string "Text3XC-0.70 " " Text3YC+0.20 " Pres(hPa)"<br>
>> "draw string "Text3XC+0.25 " " Text3YC+0.20 " " int(PMaxThee)<br>
>> "draw string "Text3XC-0.70 " " Text3YC+0.05 " `3z`0`bE`n(K)"<br>
>> "draw string "Text3XC+0.25 " " Text3YC+0.05 " " int(MaxThee)<br>
>> "draw string "Text3XC-0.70 " " Text3YC-0.10 " LI"<br>
>> "draw string "Text3XC+0.25 " " Text3YC-0.10 " "round(LIMax)<br>
>> "draw string "Text3XC-0.70 " " Text3YC-0.25 " CAPE(J)"<br>
>> "draw string "Text3XC+0.25 " " Text3YC-0.25 " "int(PosMax)<br>
>> "draw string "Text3XC-0.70 " " Text3YC-0.40 " CIN(J)"<br>
>> "draw string "Text3XC+0.25 " " Text3YC-0.40 " " int(CINMax)<br>
>> Endif<br>
>><br>
>> *-----------------------------<br>
>> * Draw wind profile along side<br>
>> *-----------------------------<br>
>><br>
>> If (DrawBarb = 1)<br>
>> say "Drawing Wind Profile."<br>
>> If (poleloc = -1)<br>
>> If (PageType = "Portrait")<br>
>> poleloc = 8.0<br>
>> Else<br>
>> poleloc = 7.5<br>
>> Endif<br>
>> Endif<br>
>> If (barbline = 1)<br>
>> "set line 1 1 3"<br>
>> "draw line "poleloc " " byloc " " poleloc " " tyloc<br>
>> Endif<br>
>> If (BarbCol = -1)<br>
>> 'set rgb 41 255 0 132'<br>
>> 'set rgb 42 255 0 168'<br>
>> 'set rgb 43 255 0 204'<br>
>> 'set rgb 44 255 0 240'<br>
>> 'set rgb 45 255 0 255'<br>
>> 'set rgb 46 204 0 255'<br>
>> 'set rgb 47 174 0 255'<br>
>> 'set rgb 48 138 0 255'<br>
>> 'set rgb 49 108 0 255'<br>
>> 'set rgb 50 84 0 255'<br>
>> 'set rgb 51 40 0 255'<br>
>> 'set rgb 52 0 0 255'<br>
>> 'set rgb 53 0 42 255'<br>
>> 'set rgb 54 0 84 255'<br>
>> 'set rgb 55 0 120 255'<br>
>> 'set rgb 56 0 150 255'<br>
>> 'set rgb 57 0 192 255'<br>
>> 'set rgb 58 0 240 255'<br>
>> 'set rgb 59 0 255 210'<br>
>> 'set rgb 60 0 255 160'<br>
>> 'set rgb 61 0 255 126'<br>
>> 'set rgb 62 0 255 78'<br>
>> 'set rgb 63 84 255 0'<br>
>> 'set rgb 64 138 255 0'<br>
>> 'set rgb 65 188 255 0'<br>
>> 'set rgb 66 236 255 0'<br>
>> 'set rgb 67 255 255 0'<br>
>> 'set rgb 68 255 222 0'<br>
>> 'set rgb 69 255 192 0'<br>
>> 'set rgb 70 255 162 0'<br>
>> 'set rgb 71 255 138 0'<br>
>> 'set rgb 72 255 108 0'<br>
>> 'set rgb 73 255 84 0'<br>
>> 'set rgb 74 255 54 0'<br>
>> 'set rgb 75 255 12 0'<br>
>> 'set rgb 76 255 0 34'<br>
>> 'set rgb 77 255 0 70'<br>
>> 'set rgb 78 255 0 105'<br>
>> 'set rgb 79 255 0 140'<br>
>> 'set rgb 80 255 0 175'<br>
>> 'set rgb 81 255 0 215'<br>
>> 'set rgb 82 255 0 255'<br>
>> 'set rgb 83 255 255 255'<br>
>><br>
>> col1='83 83 83 83 83 83 83 83 83 83 82 81 80 79 78'<br>
>> col2='77 76 75 74 73 72 71 70 69 68 67 66 65 64 63'<br>
>> col3='62 61 60 59 58 57 56 55 54 53 52 51 50 49 48'<br>
>> 'set rbcols 'col1' 'col2' 'col3<br>
>> Endif<br>
>> "set z "_zmin" "_zmax<br>
>> "set gxout stat"<br>
>> zz=1<br>
>> wspd=-999<br>
>> cont=1<br>
>> While (cont = 1 & zz < _zmax)<br>
>> "set z "zz<br>
>> pres=subwrd(result,4)<br>
>> "d "u<br>
>> rec=sublin(result,8)<br>
>> uwnd=subwrd(rec,4)<br>
>><br>
>> "d "v<br>
>> rec=sublin(result,8)<br>
>> vwnd=subwrd(rec,4)<br>
>><br>
>> wspd=GetWSpd(uwnd,vwnd)<br>
>><br>
>> if (wspd < 0 | pres > _pmax)<br>
>> zz=zz+1<br>
>> else<br>
>> cont=0<br>
>> Endif<br>
>> Endwhile<br>
>> While (zz <= _zmax)<br>
>><br>
>> "d "u"(z="zz")"<br>
>> rec=sublin(result,8)<br>
>> uwnd=subwrd(rec,4)<br>
>><br>
>> "d "v"(z="zz")"<br>
>> rec=sublin(result,8)<br>
>> vwnd=subwrd(rec,4)<br>
>><br>
>> wspd=GetWSpd(uwnd,vwnd)<br>
>><br>
>> If (BarbCol >= 0)<br>
>> "set line "BarbCol " 1 "BarbThk<br>
>> Else<br>
>> tempbcol=55+wspd/5<br>
>> If (tempbcol > 83)<br>
>> tempbcol = 83<br>
>> Endif<br>
>> "set line "tempbcol " 1 "BarbThk<br>
>> Endif<br>
>><br>
>> "d "u"(z="zz")"<br>
>> rec=sublin(result,8)<br>
>> uwnd=subwrd(rec,4)<br>
>><br>
>> "d "v"(z="zz")"<br>
>> rec=sublin(result,8)<br>
>> vwnd=subwrd(rec,4)<br>
>> if(vwnd != undef)<br>
>> wspd=GetWSpd(uwnd,vwnd)<br>
>> wdir=GetWDir(uwnd,vwnd)<br>
>> Else<br>
>> wspd=0<br>
>> wdir=0<br>
>> Endif<br>
>><br>
>> xwind=GetUWnd(wspd,wdir)<br>
>> ywind=GetVWnd(wspd,wdir)<br>
>> "query gr2xy 5 "zz<br>
>> y1=subwrd(result,6)<br>
>> if (wspd > 0)<br>
>> cc=polelen/wspd<br>
>> xendpole=poleloc-xwind*cc<br>
>> yendpole=y1-ywind*cc<br>
>> endif<br>
>> if (xendpole>0 & wspd >= 0.5)<br>
>> if (flagbase = 1)<br>
>> "draw mark 3 "poleloc " " y1 " 0.05"<br>
>> endif<br>
>> "draw line " poleloc " " y1 " " xendpole " " yendpole<br>
>> flagloop=wspd/10<br>
>> windcount=wspd<br>
>> flagcount=0<br>
>> xflagstart=xendpole<br>
>> yflagstart=yendpole<br>
>> dx=cos((180-wdir)*_dtr)<br>
>> dy=sin((180-wdir)*_dtr)<br>
>> while (windcount > 47.5)<br>
>> flagcount=flagcount+1<br>
>> dxflag=-len50*dx<br>
>> dyflag=-len50*dy<br>
>> xflagend=xflagstart+dxflag<br>
>> yflagend=yflagstart+dyflag<br>
>> windcount=windcount-50<br>
>> x1=xflagstart+0.5*wid50*xwind/wspd<br>
>> y1=yflagstart+0.5*wid50*ywind/wspd<br>
>> x2=xflagstart-0.5*wid50*xwind/wspd<br>
>> y2=yflagstart-0.5*wid50*ywind/wspd<br>
>> If (Fill50 = 1)<br>
>> "draw polyf "x1" "y1" "x2" "y2" "xflagend" "yflagend" "x1"<br>
>> "y1<br>
>> Else<br>
>> "draw line "x1 " "y1 " " xflagend " " yflagend " "<br>
>> "draw line "x2 " "y2 " " xflagend " " yflagend<br>
>> "draw line "x1 " "y1 " " x2 " " y2<br>
>> Endif<br>
>> xflagstart=xflagstart+spac50*xwind/wspd<br>
>> yflagstart=yflagstart+spac50*ywind/wspd<br>
>> endwhile<br>
>> while (windcount > 7.5 )<br>
>> flagcount=flagcount+1<br>
>> dxflag=-len10*dx<br>
>> dyflag=-len10*dy<br>
>> xflagend=xflagstart+dxflag<br>
>> yflagend=yflagstart+dyflag<br>
>> windcount=windcount-10<br>
>> "draw line " xflagstart " " yflagstart " " xflagend " " yflagend<br>
>> xflagstart=xflagstart+spac10*xwind/wspd<br>
>> yflagstart=yflagstart+spac10*ywind/wspd<br>
>> endwhile<br>
>> if (windcount > 2.5)<br>
>> flagcount=flagcount+1<br>
>> if (flagcount = 1)<br>
>> xflagstart=xflagstart+spac05*xwind/wspd<br>
>> yflagstart=yflagstart+spac05*ywind/wspd<br>
>> endif<br>
>> dxflag=-len05*dx<br>
>> dyflag=-len05*dy<br>
>> xflagend=xflagstart+dxflag<br>
>> yflagend=yflagstart+dyflag<br>
>> windcount=windcount-5<br>
>> "draw line " xflagstart " " yflagstart " " xflagend " " yflagend<br>
>> endif<br>
>> else<br>
>> if (wspd < 0.5 & wspd >= 0)<br>
>> "draw mark 2 " poleloc " " y1 " 0.08"<br>
>> endif<br>
>> endif<br>
>> zz=zz+barbint<br>
>> endwhile<br>
>> Endif<br>
>><br>
>> *----------------<br>
>> * Draw Hodograph<br>
>> *----------------<br>
>><br>
>> If (DrawHodo = 1)<br>
>> say "Drawing Hodograph."<br>
>><br>
>> If (HodXcent = -1 | HodYcent = -1)<br>
>> If (PageType = "Portrait")<br>
>> HodXcent=6<br>
>> HodYcent=9.5<br>
>> rec=sublin(result,8)<br>
>> uwnd=subwrd(rec,4)<br>
>> Else<br>
>> HodXcent=9<br>
>> HodYcent=7.0<br>
>> Endif<br>
>> Endif<br>
>> HodL=HodXcent-HodSize/2.0<br>
>> HodR=HodXcent+HodSize/2.0<br>
>> HodT=HodYcent+HodSize/2.0<br>
>> HodB=HodYcent-HodSize/2.0<br>
>> RingSpac=HodSize/(NumRing*2)<br>
>> "set line 0"<br>
>> "draw recf "HodL" "HodB" "HodR" "HodT<br>
>> "set line "HodoCol" 1 6"<br>
>> "draw rec "HodL" "HodB" "HodR" "HodT<br>
>> "set line 1 1 3"<br>
>> "set string 1 c"<br>
>> "draw mark 1 "HodXcent " "HodYcent " " HodSize<br>
>> i=1<br>
>> While (i <= NumRing)<br>
>> "set strsiz 0.10"<br>
>> "set string 1 c 3 45"<br>
>> uwnd=-i*HodRing*cos(45*_dtr)<br>
>> xloc=HodXcent+uwnd*RingSpac/HodRing<br>
>> yloc=HodYcent+uwnd*RingSpac/HodRing<br>
>><br>
>> "draw mark 2 "HodXcent " " HodYcent " " i*HodSize/NumRing<br>
>> "draw string "xloc " " yloc " " HodRing*i<br>
>> i=i+1<br>
>> Endwhile<br>
>> "set string 1 l 3 0"<br>
>> If (TickInt > 0)<br>
>> i=0<br>
>> while (i < HodRing*NumRing)<br>
>> dist=i*RingSpac/HodRing<br>
>> hrxloc=HodXcent+dist<br>
>> hlxloc=HodXcent-dist<br>
>> htyloc=HodYcent+dist<br>
>> hbyloc=HodYcent-dist<br>
>> "set line 1 1 3"<br>
>> "draw line "hrxloc " " HodYcent-TickSize/2 " " hrxloc " "<br>
>> HodYcent+TickSize/2<br>
>> "draw line "hlxloc " " HodYcent-TickSize/2 " " hlxloc " "<br>
>> HodYcent+TickSize/2<br>
>> "draw line "HodXcent+TickSize/2 " " htyloc " " HodXcent-TickSize/2<br>
>> " " htyloc<br>
>> "draw line "HodXcent+TickSize/2 " " hbyloc " " HodXcent-TickSize/2<br>
>> " " hbyloc<br>
>> i=i+TickInt<br>
>> endwhile<br>
>> Endif<br>
>> "set line "HodoCol " " HodoLine " "HodoThk<br>
>> "draw string "HodL+0.05 " " HodT-0.1 " m/s"<br>
>> zloop=_zmin<br>
>> xold=-999<br>
>> yold=-999<br>
>> count=0<br>
>> Depth=0<br>
>> While (zloop < _zmax & Depth < HodoDep)<br>
>> "set z "zloop<br>
>> pres=subwrd(result,4)<br>
>><br>
>> "d "u<br>
>> rec=sublin(result,8)<br>
>> uwnd=subwrd(rec,4)<br>
>> "d "v<br>
>> rec=sublin(result,8)<br>
>> vwnd=subwrd(rec,4)<br>
>> if (uwnd = undef)<br>
>> uwnd=0<br>
>> vwnd=0<br>
>> endif<br>
>><br>
>> If (wspd >= 0)<br>
>> xloc=HodXcent+uwnd*RingSpac/HodRing<br>
>> yloc=HodYcent+vwnd*RingSpac/HodRing<br>
>> If (xloc > 0 & yloc > 0 & xold > 0 & yold > 0)<br>
>> Depth=Depth+pold-pres<br>
>> count=count+1<br>
>> If (count = 1)<br>
>> "draw mark 3 "xold " " yold " 0.05"<br>
>> Endif<br>
>> "draw line "xold" "yold" "xloc" "yloc<br>
>> Endif<br>
>> xold=xloc<br>
>> yold=yloc<br>
>> Endif<br>
>> zloop=zloop+1<br>
>> pold=pres<br>
>> EndWhile<br>
>><br>
>> If (count > 0)<br>
>> "draw mark 3 "xold " " yold " 0.05"<br>
>> Endif<br>
>> Endif<br>
>><br>
>> *----------------------------------------------<br>
>> * Calculate and Display Absolute & S-R Helicity<br>
>> *----------------------------------------------<br>
>><br>
>> If (DrawHeli = 1)<br>
>> say "Calculating Helicity & SR Helicity."<br>
>> delp=10<br>
>> UTotal=0<br>
>> VTotal=0<br>
>><br>
>> * First, calculate mass-weighted mean wind<br>
>> * Since delp is a constant, and mass is proportional to<br>
>> * delp, this is a simple sum.<br>
>><br>
>> "set lev "_pmax " " _pmin<br>
>> "define uwndarr="sndspd"*cos((270-"snddir")*"_dtr")"<br>
>> "define vwndarr="sndspd"*sin((270-"snddir")*"_dtr")"<br>
>><br>
>> pres=MeanVBot<br>
>> While (pres >= MeanVTop)<br>
>> uwnd=interp(uwndarr,pres)*_ktm<br>
>> vwnd=interp(vwndarr,pres)*_ktm<br>
>> If (uwnd > -900 & vwnd > -900)<br>
>> UTotal=UTotal+uwnd<br>
>> VTotal=VTotal+vwnd<br>
>> Endif<br>
>> pres=pres-delp<br>
>> EndWhile<br>
>> vcount=1+(MeanVBot-MeanVTop)/delp<br>
>> Umean=UTotal/vcount<br>
>> Vmean=VTotal/vcount<br>
>> Spdmean=GetWSpd(Umean,Vmean)<br>
>> MeanDir=GetWDir(Umean,Vmean)<br>
>><br>
>> * Now, rotate and reduce mean wind to get storm motion<br>
>><br>
>> If (StormMot = 1)<br>
>> If (Spdmean < 15)<br>
>> Reduct=0.25<br>
>> Rotate=30<br>
>> Else<br>
>> Reduct=0.20<br>
>> Rotate=20<br>
>> Endif<br>
>> Else<br>
>> Reduct=0.0<br>
>> Rotate=0.0<br>
>> Endif<br>
>><br>
>> UReduce=(1-Reduct)*Umean<br>
>> VReduce=(1-Reduct)*Vmean<br>
>> StormSpd=GetWSpd(UReduce,VReduce)<br>
>><br>
>> StormDir=GetWDir(UReduce,VReduce)+Rotate<br>
>> If (StormDir >= 360)<br>
>> StormDir=StormDir-360<br>
>> Endif<br>
>><br>
>> StormU=GetUWnd(StormSpd,StormDir)<br>
>> StormV=GetVWnd(StormSpd,StormDir)<br>
>><br>
>> * Draw Storm Motion Vector<br>
>><br>
>> xloc=HodXcent+_mtk*StormU*RingSpac/HodRing<br>
>> yloc=HodYcent+_mtk*StormV*RingSpac/HodRing<br>
>><br>
>> "set line 1 1 4"<br>
>> "draw line "HodXcent " " HodYcent " " xloc " " yloc<br>
>> Arr1U=GetUWnd(HodRing/10,StormDir+30)<br>
>> Arr1V=GetVWnd(HodRing/10,StormDir+30)<br>
>> Arr2U=GetUWnd(HodRing/10,StormDir-30)<br>
>> Arr2V=GetVWnd(HodRing/10,StormDir-30)<br>
>><br>
>> xloc2=xloc-Arr1U/HodRing<br>
>> xloc3=xloc-Arr2U/HodRing<br>
>> yloc2=yloc-Arr1V/HodRing<br>
>> yloc3=yloc-Arr2V/HodRing<br>
>><br>
>> "set line 1 1 3"<br>
>><br>
>> If (FillArrw = 0)<br>
>> "draw line "xloc" "yloc" "xloc2" "yloc2<br>
>> "draw line "xloc" "yloc" "xloc3" "yloc3<br>
>> Else<br>
>> "draw polyf "xloc" "yloc" "xloc2" "yloc2" "xloc3" "yloc3" "xloc"<br>
>> "yloc<br>
>> Endif<br>
>><br>
>><br>
>> * Now, calculate SR and Environmental Helicity<br>
>><br>
>> helic=0<br>
>> SRhelic=0<br>
>> MinP=SfcPlev-HelicDep<br>
>> pres=SfcPlev<br>
>> uwndold=-999<br>
>> vwndold=-999<br>
>> While (pres >= MinP)<br>
>> uwnd=interp(uwndarr,pres)*_ktm<br>
>> vwnd=interp(vwndarr,pres)*_ktm<br>
>> If (uwnd > -900 & uwndold > -900)<br>
>> du=uwnd-uwndold<br>
>> dv=vwnd-vwndold<br>
>> ubar=0.5*(uwnd+uwndold)<br>
>> vbar=0.5*(vwnd+vwndold)<br>
>> uhelic=-dv*ubar<br>
>> vhelic=du*vbar<br>
>> SRuhelic=-dv*(ubar-StormU)<br>
>> SRvhelic=du*(vbar-StormV)<br>
>> SRhelic=SRhelic+SRuhelic+SRvhelic<br>
>> helic=helic+uhelic+vhelic<br>
>> Endif<br>
>> uwndold=uwnd<br>
>> vwndold=vwnd<br>
>> pres=pres-delp<br>
>> EndWhile<br>
>><br>
>> "set strsiz 0.1"<br>
>> "set string 1 l 3"<br>
>> If (PageType = "Portrait")<br>
>> If (Text4XC = -1)<br>
>> Text4XC=rxloc-0.75<br>
>> Endif<br>
>> If (Text4YC = -1)<br>
>> Text4YC=tyloc-5.45<br>
>> Endif<br>
>> Else<br>
>> If (Text4XC = -1)<br>
>> Text4XC=rxloc+2.50<br>
>> Endif<br>
>> If (Text4YC = -1)<br>
>> Text4YC=tyloc-6.10<br>
>> Endif<br>
>> Endif<br>
>> "set line 0 1 3"<br>
>> "draw recf "Text4XC-0.75 " "Text4YC-0.5 " " Text4XC+0.75 " "<br>
>> Text4YC+0.5<br>
>> "set line 1 1 3"<br>
>> "draw rec "Text4XC-0.75 " "Text4YC-0.5 " " Text4XC+0.75 " " Text4YC+0.5<br>
>> "draw string "Text4XC-0.45 " " Text4YC+0.40 " Hodograph"<br>
>> "draw string "Text4XC-0.70 " " Text4YC+0.20 " EH"<br>
>> "draw string "Text4XC+0.25 " " Text4YC+0.20 " "int(helic)<br>
>> "draw string "Text4XC-0.70 " " Text4YC+0.05 " SREH"<br>
>> "draw string "Text4XC+0.25 " " Text4YC+0.05 " " int(SRhelic)<br>
>> "draw string "Text4XC-0.70 " " Text4YC-0.20 " StmDir"<br>
>> "draw string "Text4XC+0.25 " " Text4YC-0.20 " " int(StormDir)"`3.`0"<br>
>> "draw string "Text4XC-0.70 " " Text4YC-0.35 " StmSpd(kt)"<br>
>> "draw string "Text4XC+0.25 " " Text4YC-0.35 " " int(_mtk*StormSpd)<br>
>> Endif<br>
>><br>
>> *---------------------------------------------------<br>
>> * Plot RH profile.<br>
>> *---------------------------------------------------<br>
>><br>
>> If (DrawRH = 1)<br>
>> "set z "_zmin" "_zmax<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set t "_tval<br>
>> "set zlog on"<br>
>> "set gxout line"<br>
>> "set ccolor "RHCol<br>
>> "set cstyle "RHLine<br>
>> "set cmark "RHMrk<br>
>> "set digsiz "MrkSize<br>
>> "set missconn on"<br>
>> "set xlab on"<br>
>> "set frame off"<br>
>> "set vrange 0 350"<br>
>> "set xlpos 0 t"<br>
>> "set xlevs 25 50 75 100"<br>
>> "set grid vertical 5"<br>
>> "define<br>
>> rh=100*exp((17.2694*"snddewp")/("snddewp"+237.3)-(17.2694*"sndtemp")/("sndtemp"+237.3))"<br>
>> "d rh"<br>
>> If (LblAxes = 1)<br>
>> "set string 1 c 3 0"<br>
>> "set strsiz 0.125"<br>
>> If (PageType = "Portrait")<br>
>> "draw string 1.5 10.85 RH (%)"<br>
>> Else<br>
>> "draw string 1.75 8.35 RH (%)"<br>
>> Endif<br>
>> Endif<br>
>> Endif<br>
>><br>
>> *------------------------------------------<br>
>> * Reset environment to original dimensions<br>
>> *------------------------------------------<br>
>><br>
>> "set t "_tval<br>
>> "set x "_xval<br>
>> "set y "_yval<br>
>> "set z "_zmin " "_zmax<br>
>><br>
>> say "Done."<br>
>><br>
>> *Return(0)<br>
>><br>
>> *************************************************************************<br>
>> function Templcl(temp,dewp)<br>
>><br>
>> *------------------------------------------------------<br>
>> * Calculate the temp at the LCL given temp & dewp in C<br>
>> *------------------------------------------------------<br>
>><br>
>> tempk=temp+273.15<br>
>> dewpk=dewp+273.15<br>
>> Parta=1/(dewpk-56)<br>
>> Partb=log(tempk/dewpk)/800<br>
>> Tlcl=1/(Parta+Partb)+56<br>
>> return(Tlcl-273.15)<br>
>><br>
>> **************************************************************************<br>
>><br>
>> function Preslcl(temp,dewp,pres)<br>
>><br>
>> *-------------------------------------------------------<br>
>> * Calculate press of LCL given temp & dewp in C and pressure<br>
>> *-------------------------------------------------------<br>
>><br>
>> Tlclk=Templcl(temp,dewp)+273.15<br>
>> tempk=temp+273.15<br>
>> theta=tempk*pow(1000/pres,0.286)<br>
>> plcl=1000*pow(Tlclk/theta,3.48)<br>
>> return(plcl)<br>
>><br>
>> **************************************************************************<br>
>> function LiftWet(startt,startp,endp,display,Pmin,Pmax)<br>
>><br>
>> *--------------------------------------------------------------------<br>
>> * Lift a parcel moist adiabatically from startp to endp.<br>
>> * Init temp is startt in C. If you wish to see the parcel's<br>
>> * path plotted, display should be 1. Returns temp of parcel at endp.<br>
>> *--------------------------------------------------------------------<br>
>><br>
>> temp=startt<br>
>> pres=startp<br>
>> cont = 1<br>
>> delp=10<br>
>> While (pres >= endp & cont = 1)<br>
>> If (display = 1)<br>
>> xtemp=GetXLoc(temp,pres)<br>
>> "q w2xy "xtemp" "pres<br>
>> xloc=subwrd(result,3)<br>
>> yloc=subwrd(result,6)<br>
>> If (xtemp < 0 | xtemp > 100)<br>
>> cont=0<br>
>> Else<br>
>> If (pres >= Pmin & pres < Pmax & pres < startp)<br>
>> "draw line "xold" "yold" "xloc" "yloc<br>
>> Endif<br>
>> Endif<br>
>> Endif<br>
>> xold=xloc<br>
>> yold=yloc<br>
>> temp=temp-100*delp*gammaw(temp,pres-delp/2,100)<br>
>> pres=pres-delp<br>
>> EndWhile<br>
>> return(temp)<br>
>><br>
>><br>
>> **************************************************************************<br>
>> function LiftDry(startt,startp,endp,display,Pmin,Pmax)<br>
>><br>
>> *--------------------------------------------------------------------<br>
>> * Lift a parcel dry adiabatically from startp to endp.<br>
>> * Init temp is startt in C. If you wish to see the parcel's<br>
>> * path plotted, display should be 1. Returns temp of parcel at endp.<br>
>> *--------------------------------------------------------------------<br>
>><br>
>> starttk=startt+273.15<br>
>> cont = 1<br>
>> delp=10<br>
>> round=int(startp/10)*10<br>
>> subscr=0.1*round<br>
>> powstart=pow(startp,-0.286)<br>
>> temp=starttk*_powpres.subscr*powstart-273.15<br>
>> pres=round-10<br>
>> While (pres >= endp & cont = 1)<br>
>> subscr=0.1*pres<br>
>> temp=starttk*_powpres.subscr*powstart-273.15<br>
>> If (display = 1)<br>
>> xtemp=GetXLoc(temp,pres)<br>
>> "q w2xy "xtemp" "pres<br>
>> xloc=subwrd(result,3)<br>
>> yloc=subwrd(result,6)<br>
>> If (xtempold > 0 & xtempold < 100 & xtemp > 0 & xtemp < 100)<br>
>> If (pres >= Pmin & pres < Pmax & pres < startp)<br>
>> "draw line "xold" "yold" "xloc" "yloc<br>
>> Endif<br>
>> Endif<br>
>> Endif<br>
>> xold=xloc<br>
>> xtempold=xtemp<br>
>> yold=yloc<br>
>> pres=pres-delp<br>
>> EndWhile<br>
>> return(temp)<br>
>><br>
>> **************************************************************************<br>
>> function CAPE(startt,startp,endp,sndtemp,snddewp)<br>
>><br>
>> *---------------------------------------------------------------------<br>
>> * Returns all postive area and convective inhibition above LCL.<br>
>> * Parcel is lifted from LCL at startt,startp and is halted<br>
>> * at endp. Integration method used is trapezoid method.<br>
>> *---------------------------------------------------------------------<br>
>><br>
>> pres=startp<br>
>> PclTemp=startt<br>
>> PclTempV=virtual2(PclTemp+273.15,PclTemp+273.15,pres)-273.15<br>
>> delp=10<br>
>> Pos=0<br>
>> Neg=0<br>
>> Neg2=0<br>
>><br>
>> count=0<br>
>> While (pres >= endp)<br>
>> EnvTemp=interp(sndtemp,pres)<br>
>> EnvDewp=interp(snddewp,pres)<br>
>> EnvTempV=virtual2(EnvTemp+273.15,EnvDewp+273.15,pres)-273.15<br>
>> DelT=PclTempV-EnvTempV<br>
>> If (abs(EnvTempV) < 130 & abs(PclTempV) < 130)<br>
>> count=count+1<br>
>> If (count > 1)<br>
>> Val=DelT/pres+Prev<br>
>> If (Val > 0)<br>
>> Pos=Pos+Val<br>
>> Neg2=0<br>
>> Else<br>
>> Neg=Neg+abs(Val)<br>
>> Neg2=Neg2+abs(Val)<br>
>> Endif<br>
>> Endif<br>
>> Prev=DelT/pres<br>
>> Endif<br>
>> pres=pres-delp<br>
>> PclTemp=PclTemp-100*delp*gammaw(PclTemp,pres,100)<br>
>> PclTempV=virtual2(PclTemp+273.15,PclTemp+273.15,pres)-273.15<br>
>> Endwhile<br>
>><br>
>> Pos=0.5*Pos*287*delp<br>
>> CIN=0.5*(Neg-Neg2)*287*delp<br>
>><br>
>> return(Pos" "CIN)<br>
>><br>
>><br>
>> ***************************************************************************<br>
>> function gammaw(tempc,pres,rh)<br>
>><br>
>> *-----------------------------------------------------------------------<br>
>> * Function to calculate the moist adiabatic lapse rate (deg C/Pa) based<br>
>> * on the temperature, pressure, and rh of the environment.<br>
>> *----------------------------------------------------------------------<br>
>><br>
>> tempk=tempc+273.15<br>
>> es=satvap2(tempc)<br>
>> ws=mixratio(es,pres)<br>
>> w=rh*ws/100<br>
>> tempv=virtual(tempk,w)<br>
>> latent=latentc(tempc)<br>
>><br>
>> A=1.0+latent*ws/(287*tempk)<br>
>> B=1.0+0.622*latent*latent*ws/(1005*287*tempk*tempk)<br>
>> Density=100*pres/(287*tempv)<br>
>> lapse=(A/B)/(1005*Density)<br>
>> return(lapse)<br>
>><br>
>> *************************************************************************<br>
>> function latentc(tempc)<br>
>><br>
>> *-----------------------------------------------------------------------<br>
>> * Function to return the latent heat of condensation in J/kg given<br>
>> * temperature in degrees Celsius.<br>
>> *-----------------------------------------------------------------------<br>
>><br>
>> val=2502.2-2.43089*tempc<br>
>><br>
>> return(val*1000)<br>
>><br>
>> *************************************************************************<br>
>> function precipw(sndtemp,snddewp,startp,endp)<br>
>><br>
>> *-----------------------------------------------------------------------<br>
>> * Function to calculate the precipitable water (cm) in a sounding<br>
>> * starting at pressure level startp and ending at pressure level endp.<br>
>> *-----------------------------------------------------------------------<br>
>><br>
>> ppold=-999<br>
>> ttold=-999<br>
>> ddold=-999<br>
>> delp=10<br>
>> Int=0<br>
>> mix=0<br>
>> pres=startp<br>
>> logpp=log(pres)<br>
>> logppm=log(pres-delp)<br>
>> while (pres >= endp)<br>
>> tt=interp(sndtemp,pres)<br>
>> dd=interp(snddewp,pres)<br>
>> if (tt>-900 & ttold>-900 & dd>-900 & ddold>-900)<br>
>> e=satvap2(dd)<br>
>> mix=mixratio(e,pres)<br>
>> mixavg=(logpp*mix+logppm*mixold)/(logpp+logppm)<br>
>> Int=Int+1.020408*mixavg*delp<br>
>> endif<br>
>> ttold=tt<br>
>> ddold=dd<br>
>> ppold=pp<br>
>> mixold=mix<br>
>> pres=pres-delp<br>
>> logpp=logppm<br>
>> logppm=log(pres-delp)<br>
>> endwhile<br>
>><br>
>> return(Int)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function virtual(temp,mix)<br>
>><br>
>> *------------------------------------------------------------<br>
>> * Function to return virtual temperature given temperature in<br>
>> * kelvin and mixing ratio in g/g.<br>
>> *-------------------------------------------------------------<br>
>><br>
>> tempv=temp*(1.0+0.6*mix)<br>
>><br>
>> return (tempv)<br>
>><br>
>> ************************************************************************<br>
>><br>
>> function virtual2(temp,dewp,pres)<br>
>><br>
>> *------------------------------------------------------------<br>
>> * Function to return virtual temperature in kelvin given temperature in<br>
>> * kelvin and dewpoint in kelvin and pressure in mb<br>
>> *-------------------------------------------------------------<br>
>><br>
>> if (temp > 0 & dewp > 0)<br>
>> vap=satvap2(dewp-273.15)<br>
>> mix=mixratio(vap,pres)<br>
>> tempv=virtual(temp,mix)<br>
>> else<br>
>> tempv=-9999<br>
>> endif<br>
>><br>
>> return (tempv)<br>
>><br>
>> ************************************************************************<br>
>><br>
>> function satvapor(temp)<br>
>><br>
>> *---------------------------------------------------------------<br>
>> * Given temp in Celsius, returns saturation vapor pressure in mb<br>
>> *---------------------------------------------------------------<br>
>><br>
>><br>
>> pol=_C0+temp*(_C1+temp*(_C2+temp*(_C3+temp*(_C4+temp*(_C5+temp*(_C6+temp*(_C7+temp*(_C8+temp*(_C9)))))))))<br>
>><br>
>> return(6.1078/pow(pol,8))<br>
>><br>
>> ************************************************************************<br>
>><br>
>> function satvap2(temp)<br>
>><br>
>> *---------------------------------------------------------------<br>
>> * Given temp in Celsius, returns saturation vapor pressure in mb<br>
>> *---------------------------------------------------------------<br>
>><br>
>> es=6.112*exp(17.67*temp/(temp+243.5))<br>
>><br>
>> return(es)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function mixratio(e,p)<br>
>><br>
>> *------------------------------------------------------<br>
>> * Given vapor pressure and pressure, return mixing ratio<br>
>> *-------------------------------------------------------<br>
>><br>
>> mix=0.622*e/(p-e)<br>
>><br>
>> return(mix)<br>
>><br>
>> ************************************************************************<br>
>><br>
>> function getrh(temp,dewp,pres)<br>
>><br>
>> tempk=temp+273.15<br>
>> dewpk=dewp+273.15<br>
>><br>
>> es=satvap2(temp)<br>
>><br>
>> If (temp > 0)<br>
>> A=2.53*pow(10,9)<br>
>> B=5420<br>
>> Else<br>
>> A=3.41*pow(10,10)<br>
>> B=6130<br>
>> Endif<br>
>><br>
>> w=A*0.622*exp(-B/dewpk)/pres<br>
>> ws=mixratio(es,pres)<br>
>><br>
>> return(100*w/ws)<br>
>><br>
>> ************************************************************************<br>
>> function interp(array,pres)<br>
>><br>
>> *------------------------------------------------------------------------<br>
>> * Interpolate inside array for pressure level pres.<br>
>> * Returns estimated value of array at pressure pres.<br>
>> *------------------------------------------------------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "set lev "pres<br>
>> altpres=subwrd(result,4)<br>
>> "q dims"<br>
>> rec=sublin(result,4)<br>
>> zlev=subwrd(rec,9)<br>
>><br>
>> If (zlev < 2 | zlev > _zmaxfile)<br>
>> Vest = -9999.0<br>
>> Else<br>
>> If (altpres > pres)<br>
>> zlev=zlev+1<br>
>> Endif<br>
>> "set z "zlev<br>
>> PAbove=subwrd(result,4)<br>
>> "d "array"(lev="PAbove")"<br>
>> rec=sublin(result,8)<br>
>> VAbove=subwrd(rec,4)<br>
>> "set z "zlev-1<br>
>> PBelow=subwrd(result,4)<br>
>> "d "array"(lev="PBelow")"<br>
>> rec=sublin(result,8)<br>
>> VBelow=subwrd(rec,4)<br>
>><br>
>> * Now if we are in a region of missing data, find next good level.<br>
>><br>
>> If (abs(VAbove) > 130 & zlev > 1 & zlev < _zmaxfile)<br>
>> zz=zlev<br>
>> While (abs(VAbove) > 130 & zz < _zmaxfile)<br>
>> zz=zz+1<br>
>> "set z "zz<br>
>> PAbove=subwrd(result,4)<br>
>> "d "array"(lev="PAbove")"<br>
>> rec=sublin(result,8)<br>
>> VAbove=subwrd(rec,4)<br>
>> EndWhile<br>
>> Endif<br>
>><br>
>> If (abs(VBelow) > 130 & zlev > 1 & zlev < _zmaxfile)<br>
>> zz=zlev-1<br>
>> While (abs(VBelow) > 130 & zz > 1)<br>
>> zz=zz-1<br>
>> "set z "zz<br>
>> PBelow=subwrd(result,4)<br>
>> "d "array"(lev="PBelow")"<br>
>> rec=sublin(result,8)<br>
>> VBelow=subwrd(rec,4)<br>
>> EndWhile<br>
>> Endif<br>
>><br>
>> If (abs(VAbove) < 130 & abs(VBelow) < 130)<br>
>> Vest=VBelow+log(PBelow/pres)*(VAbove-VBelow)/log(PBelow/PAbove)<br>
>> Else<br>
>> Vest=-9999.0<br>
>> Endif<br>
>><br>
>> Endif<br>
>><br>
>> Return(Vest)<br>
>><br>
>><br>
>> ****************************************************************************<br>
>><br>
>> function GetUWnd(wspd,wdir)<br>
>><br>
>> *------------------------<br>
>> * Get x-component of wind.<br>
>> *------------------------<br>
>><br>
>><br>
>> If (wspd >= 0)<br>
>> xwind=wspd*cos((270-wdir)*_dtr)<br>
>> Else<br>
>> xwind = -9999.0<br>
>> Endif<br>
>> return(xwind)<br>
>><br>
>> **************************************************************************<br>
>><br>
>> function GetVWnd(wspd,wdir)<br>
>><br>
>> *-----------------------<br>
>> * Get y-component of wind<br>
>> *------------------------<br>
>><br>
>> If (wspd >= 0)<br>
>> ywind=wspd*sin((270-wdir)*_dtr)<br>
>> Else<br>
>> ywind = -9999.0<br>
>> Endif<br>
>> return(ywind)<br>
>><br>
>><br>
>> *************************************************************************<br>
>><br>
>> function GetWSpd(xwind,ywind)<br>
>><br>
>><br>
>> "set gxout stat"<br>
>> "d mag("xwind","ywind")"<br>
>> rec=sublin(result,8)<br>
>> val=subwrd(rec,4)<br>
>><br>
>> return (val)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function GetWDir(xwind,ywind)<br>
>><br>
>> * Return wind direction given x and y components<br>
>><br>
>> "set gxout stat"<br>
>> "define theta=270-"_rtd"*atan2("ywind","xwind")"<br>
>> "d theta"<br>
>> rec=sublin(result,8)<br>
>> Dir=subwrd(rec,4)<br>
>><br>
>> If (Dir > 360)<br>
>> Dir=Dir-360<br>
>> Endif<br>
>><br>
>> If (Dir < 0)<br>
>> Dir=360+Dir<br>
>> Endif<br>
>><br>
>> return(Dir)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function GetXLoc(temp,pres)<br>
>><br>
>> *-------------------------------------------------<br>
>> * Get x-location on skew-t based on temp, pressure<br>
>> *-------------------------------------------------<br>
>><br>
>> xloc=(temp-_m1*log10(pres)-_m3)/_m2<br>
>> return(xloc)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function GetTemp(xloc,pres)<br>
>><br>
>> *-------------------------------------------------<br>
>> * Return temperature at location given by xloc,pres<br>
>> *-------------------------------------------------<br>
>><br>
>> tempval=_m1*log10(pres)+_m2*xloc+_m3<br>
>> return(tempval)<br>
>><br>
>> **************************************************************************<br>
>><br>
>> function GetTheta(temp,pres)<br>
>><br>
>> *---------------------------------------------------<br>
>> * Calculate theta for a given temperature and pressure<br>
>> *---------------------------------------------------<br>
>><br>
>> theta=(temp+273.15)*pow(1000/pres,0.286)-273.15<br>
>> return(theta)<br>
>><br>
>><br>
>> *************************************************************************<br>
>><br>
>> function GetThet2(temp,dewp,pres)<br>
>><br>
>> *---------------------------------------------------<br>
>> * Calculate theta for a given temperature,dewp, and pressure<br>
>> *---------------------------------------------------<br>
>><br>
>> tempk=273.15+temp<br>
>> dewpk=273.15+dewp<br>
>><br>
>> es=satvap2(temp)<br>
>> ws=mixratio(es,pres)<br>
>><br>
>> mix=10*getrh(temp,dewp,pres)*ws<br>
>><br>
>> exponent=0.2854*(1.0-0.00028*mix)<br>
>> theta=(temp+273.15)*pow(1000/pres,exponent)-273.15<br>
>> return(theta)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function Thetae(temp,dewp,pres)<br>
>><br>
>> *--------------------------------------------------------------<br>
>> * Return equiv. pot. temp in Kelvin given temp, dewp in celsius<br>
>> * From Bolton (1980) Mon Wea Rev<br>
>> *--------------------------------------------------------------<br>
>><br>
>> es=satvap2(temp)<br>
>> ws=mixratio(es,pres)<br>
>> mix=10*getrh(temp,dewp,pres)*ws<br>
>> theta=GetThet2(temp,dewp,pres)+273.15<br>
>> TLcl=Templcl(temp,dewp)+273.15<br>
>> thetae=theta*exp((3.376/TLcl-0.00254)*mix*(1.0+0.00081*mix))<br>
>><br>
>> return(thetae)<br>
>><br>
>> **************************************************************************<br>
>><br>
>><br>
>> function int(i0)<br>
>><br>
>> *--------------------------<br>
>> * Return integer of i0<br>
>> *--------------------------<br>
>> i=0<br>
>> while(i<12)<br>
>> i=i+1<br>
>> if(substr(i0,i,1)='.')<br>
>> i0=substr(i0,1,i-1)<br>
>> break<br>
>> endif<br>
>> endwhile<br>
>> return(i0)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function abs(i)<br>
>><br>
>> *----------------------------<br>
>> * return absolute value of i<br>
>> *----------------------------<br>
>><br>
>> if (i < 0)<br>
>> absval=-i<br>
>> else<br>
>> absval=i<br>
>> endif<br>
>><br>
>> return(absval)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function log(i)<br>
>><br>
>> *---------------------------<br>
>> * return natural log of i<br>
>> *---------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "d log("i")"<br>
>> rec=sublin(result,8)<br>
>> val=subwrd(rec,4)<br>
>> return(val)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function log10(i)<br>
>><br>
>> *--------------------------<br>
>> * return log base 10 of i<br>
>> *--------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "d log10("i")"<br>
>> rec=sublin(result,8)<br>
>> val=subwrd(rec,4)<br>
>> return(val)<br>
>><br>
>> *************************************************************************<br>
>><br>
>> function pow(i,j)<br>
>><br>
>> *-------------------------------<br>
>> * return power of i raised to j<br>
>> *-------------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "d pow("i","j")"<br>
>> rec=sublin(result,8)<br>
>> val=subwrd(rec,4)<br>
>> return(val)<br>
>><br>
>> ************************************************************************<br>
>><br>
>> function cos(i)<br>
>><br>
>> *-----------------------------------------<br>
>> * return cosine of i, where i is in radians<br>
>> *------------------------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "d cos("i")"<br>
>> rec=sublin(result,8)<br>
>> val=subwrd(rec,4)<br>
>> return(val)<br>
>><br>
>> ************************************************************************<br>
>><br>
>> function sin(i)<br>
>><br>
>> *------------------------------------------<br>
>> * return sine of i, where i is in radians<br>
>> *------------------------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "d sin("i")"<br>
>> rec=sublin(result,8)<br>
>> val=subwrd(rec,4)<br>
>> return(val)<br>
>><br>
>> ************************************************************************<br>
>><br>
>> function exp(i)<br>
>><br>
>> *------------------------------------------<br>
>> * return exponential of i<br>
>> *------------------------------------------<br>
>><br>
>> "set gxout stat"<br>
>> "d exp("i")"<br>
>> rec=sublin(result,8)<br>
>> val=subwrd(rec,4)<br>
>> return(val)<br>
>><br>
>> ***********************************************************************<br>
>> function round(i)<br>
>><br>
>> rr=abs(1.0*i)<br>
>> rr=int(rr+0.5)<br>
>> if (i < 0)<br>
>> rr=-1*rr<br>
>> endif<br>
>> return(rr)<br>
>><br>
>><br>
>><br>
>> -------------- next part --------------<br>
>> An HTML attachment was scrubbed...<br>
>> URL:<br>
>> <a href="http://gradsusr.org/pipermail/gradsusr/attachments/20101108/1c1328b3/attachment.html" target="_blank">http://gradsusr.org/pipermail/gradsusr/attachments/20101108/1c1328b3/attachment.html</a><br>
>><br>
>> ------------------------------<br>
>><br>
>><br>
>> _______________________________________________<br>
>> gradsusr mailing list<br>
>> <a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
>> <a href="http://gradsusr.org/mailman/listinfo/gradsusr" target="_blank">http://gradsusr.org/mailman/listinfo/gradsusr</a><br>
>><br>
>><br>
>> End of gradsusr Digest, Vol 9, Issue 13<br>
>> ***************************************<br>
>><br>
><br>
><br>
> _______________________________________________<br>
> gradsusr mailing list<br>
> <a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
> <a href="http://gradsusr.org/mailman/listinfo/gradsusr" target="_blank">http://gradsusr.org/mailman/listinfo/gradsusr</a><br>
><br>
><br>
<br>
<br>
--<br>
Jeff Duda<br>
Iowa State University<br>
Meteorology Graduate Student<br>
3134 Agronomy Hall<br>
<a href="http://www.meteor.iastate.edu/%7Ejdduda" target="_blank">www.meteor.iastate.edu/~jdduda</a><br>
-------------- next part --------------<br>
An HTML attachment was scrubbed...<br>
URL: <a href="http://gradsusr.org/pipermail/gradsusr/attachments/20101108/f8b27adf/attachment.html" target="_blank">http://gradsusr.org/pipermail/gradsusr/attachments/20101108/f8b27adf/attachment.html</a><br>
<br>
------------------------------<br>
<br>
_______________________________________________<br>
gradsusr mailing list<br>
<a href="mailto:gradsusr@gradsusr.org">gradsusr@gradsusr.org</a><br>
<a href="http://gradsusr.org/mailman/listinfo/gradsusr" target="_blank">http://gradsusr.org/mailman/listinfo/gradsusr</a><br>
<br>
<br>
End of gradsusr Digest, Vol 9, Issue 16<br>
***************************************<br>
</blockquote></div><br>