how to make maximum wind plot along with level showing maximum wind for each grid
sushant puranik
sushantpuranik at GMAIL.COM
Fri Nov 27 05:17:44 EST 2009
Sorry to disturb you again.
I am not able to solve my problem till yet.
Since i want to study the structure of a jet stream.I am trying to find
level of maximum wind.
At present i am able to print maximum wind speed with the help of
d max(mag(u,v),lev=700,lev=1000). (Since my data is arranged in reverse
order)
But i need alongwith wind speed it should print level to which maximum value
corresponds.
Because d maxloc(mag(u,v),lev=700,lev=1000) gives location of grid to which
maximum values belongs. But instead of that i need pressure level.
how to solve this problem.
My single file contain data of one month (in netcdf format) and data is at 6
hour interval.
If you have any solution please send me, i am trying hard for solving this
particular problem.
Thanking you
On Wed, Nov 25, 2009 at 11:49 PM, Charles Seman <Charles.Seman at noaa.gov>wrote:
> Sushant,
>
> Sorry for not explaining further... The script is designed to plot time
> and zonally-averaged mass fluxes and estimated tropopause heights for
> two of GFDL's models, and uses GrADS ctl files to access datasets that
> you do not have... The datasets that are accessed are quite large (think
> too large to send via email). Maybe the script could be useful to you
> as an example of how to do some calculations in GrADS? What is the
> problem you are trying to solve? Or what diagnostic are you trying to
> calculate? If it is to find the level of maximum wind, the code (see
> earlier email) to find the level of minimum temperature could be
> modified to find the level of maximum wind (the code was offered as an
> example of how to use GrADS to do this class of problems and not to
> solve your specific problem...).
>
>
> Hope this helps,
> Chuck
>
> sushant puranik wrote:
>
>> Thank you for immediate reply i am testing those scripts. but while
>> running the script it shows error in the information of GrADS files
>>
>> level_dir = '/net2/cjs/ljd/lima/ch3i1/pp/atmos_level/av/monthly_16yr'
>> don_dir = '/net2/cjs/ljd/lima/ch3i1/pp/atmosdon/av/monthly_16yr'
>> ras_dir =
>> '/net2/cjs/lwh/fms/lima/m45_am2p14_ch3i/pp/atmos_level/av/monthly_16yr'
>>
>> sdir_months = 'Aug Sep'
>> month_days = '31 30'
>>
>> /home/Sushant/Desktop/level_ctl = 'mc.o_zl.ctl' ; level_var = 'mc'
>> don_ctl1 = 'uceml_deep.o_zl.ctl' ; don_var1 = 'uceml_deep'
>> don_ctl2 = 'umeml_deep.o_zl.ctl' ; don_var2 = 'umeml_deep'
>> don_ctl3 = 'dmeml_deep.o_zl.ctl' ; don_var3 = 'dmeml_deep'
>>
>> ras_ctl = 'mc.o_zl.ctl' ; ras_var = 'mc'
>>
>> temp_ctl = 'temp.o_zl.ctl' ; temp_var = 'temp'
>>
>> So can you please send me those files.
>>
>> Also i am using ECMWF data in the netcdf format. So there is no need
>> to use ctl files. It is possible to run the script without using such
>> ctl file info.
>>
>> Sushant
>>
>>
>> On Wed, Nov 25, 2009 at 1:01 AM, Charles Seman <Charles.Seman at noaa.gov
>> <mailto:Charles.Seman at noaa.gov>> wrote:
>>
>> Hi Phillip,
>>
>> Please find attached a script "fig3ab.gs <http://fig3ab.gs>" and
>>
>> called script functions
>> "pause.gsf" and "scale_zdef_zlevs.gsf". "scale_zdef_zlevs.gsf"
>> makes a
>> local (temporary) ctl file with scaled vertical levels for plotting
>> purposes: here specifically to convert the original zdef zlevs values
>> from meters to kilometers (there may be a better way to do this; maybe
>> using commands "set ylab off", followed by "set yaxis ..."). I
>> hope all
>> of the scripts needed for the fig3ab.gs <http://fig3ab.gs> script
>>
>> are attached; if not,
>> please let me know. You may have questions about these scripts;
>> please
>> contact me if you do...
>>
>>
>> Hope this helps,
>> Chuck
>>
>> Phillip Lin wrote:
>>
>> Hi Chuck,
>> Thanks for your dedicated input to this problem we have.
>> And could you be kind to send me a script you did before for
>> the sane or
>> similar question? Then, maybe it is much helpful for me to
>> solve this
>> problem!
>> Thanks much!
>> Phillip
>>
>>
>> Hi Phillip,
>>
>> Sorry I didn't explain further: the code was meant as an
>> example of how
>> to find and plot extrema in a field of data. To do the
>> max wind, the
>> code would have to be modified to search in a wind speed
>> field using
>> "maxloc"... and the exact code would also depend on the
>> domain...
>>
>> Hope this helps,
>> Chuck
>>
>> Phillip Lin wrote:
>>
>> Hi Chuck,
>> I also met this problem. And there are something
>> unclear to me (only to
>> me!)in your last email, because:
>> (1) the level that has minimum temp isn't the same as
>> the level where
>> there is maximum zonal wind speed;
>> (2) what does the meaning of the step (3) you said;
>> sorry how can plot
>> the
>> values (xloc,j,yloc,j) on a (y,z)plot (this is the
>> cross section in the
>> height-latitude plane?
>> (3)could you be kind to tell me the basic ideas to get
>> the maximum jet
>> stream, or... this means the contours of zonal wind
>> speed in a specific
>> level as 250mb or ...?
>> (4)could you be kind to send me the script to do this
>> plotting, and data
>> file too for a test.
>> I am quite confused about this! Thank you so much in
>> advance!
>> Phillip
>>
>>
>>
>>
>> Sushant,
>>
>> Here's some code that:
>> 1) finds the height level "zmin" of the minimum
>> temperature
>> (represented
>> by variable ztav = zonal average of weighted time
>> average temperature)
>> at each horizontal location "j" in a (y,z)
>> physical domain
>> 2) uses GrADS function "gr2xy"
>> (
>> http://grads.iges.org/grads/gadoc/script.html#commands)
>> to transform
>> the grid point location (j,zmin) to plot location
>> (xloc.j,yloc.j)
>> 3) plots the values (xloc.j,yloc.j) on a (y,z) plot
>>
>> say ' find and save location of minimum
>> temperature...'
>>
>> j=1
>> while ( j <= ny )
>> 'set y 'j
>> 'set gxout print'
>> 'd minloc(ztav,z=1,z='nz')' ; rec3 =
>> sublin(result,3) ; zmin =
>> subwrd(rec3,1)
>> * say result
>> say ' ...j, zmin = 'j', 'zmin
>> 'query gr2xy 'j' 'zmin
>> rec = sublin(result,1) ; xloc.j =
>> subwrd(rec,3) ; yloc.j =
>> subwrd(rec,6)
>> * say result
>> say ' xloc, yloc = 'xloc.j', 'yloc.j
>> j=j+1
>> endwhile
>> ...
>> if( plot_min_temp = 'yes' )
>> say ' plot level of minimum temperature... '
>> j=1
>> while ( j <= ny )
>> 'set y 'j
>> * say ' xloc, yloc = 'xloc.j', 'yloc.j
>> 'set line 1'
>> 'draw mark 3 'xloc.j' 'yloc.j' 0.03'
>> j=j+1
>> endwhile
>> endif
>>
>> Please contact me if you have any questions about
>> this code, or would
>> like a script that does a plot.
>>
>> Hope this helps,
>> Chuck
>>
>>
>> sushant puranik wrote:
>>
>>
>> Hello Sir
>> I want to study the structure of a jet stream
>> for that purpose i wish
>> to make a plot showing maximum wind speed for
>> each grid box alongwith
>> the level in which the speed is maximum. Means
>> it should print level
>> with maximum wind speed value. I am trying to
>> solve this problem since
>> last 15 days but i am not able to solve it.
>> Does anyone has any idea
>> how i can make such plot using GrADS. I am
>> using ECMWF interim
>> reanalysis data for this purpose.
>>
>> Any suggestion is appreciated
>>
>> Thanking you.
>>
>> --
>> Sushant Puranik
>> Junior Research Fellow
>> Dept. of Atmospheric & Space Sciences,
>> University of Pune,
>> Pune-07,
>> India.
>>
>>
>> --
>>
>> Please note that Charles.Seman at noaa.gov
>> <mailto:Charles.Seman at noaa.gov> should be
>>
>> considered my NOAA
>> email address, not cjs at gfdl.noaa.gov
>> <mailto:cjs at gfdl.noaa.gov>.
>>
>>
>>
>> ********************************************************************
>> Charles Seman
>> Charles.Seman at noaa.gov
>> <mailto:Charles.Seman at noaa.gov>
>>
>> U.S. Department of Commerce / NOAA / OAR
>> Geophysical Fluid Dynamics Laboratory
>> voice: (609) 452-6547
>> 201 Forrestal Road
>> fax: (609) 987-5063
>> Princeton, NJ 08540-6649
>> http://www.gfdl.noaa.gov/~cjs/<http://www.gfdl.noaa.gov/%7Ecjs/>
>> <http://www.gfdl.noaa.gov/%7Ecjs/>
>>
>>
>> ********************************************************************
>>
>> "The contents of this message are mine personally
>> and do not
>> necessarily
>> reflect any position of the Government or NOAA."
>>
>>
>>
>> --
>>
>> Please note that Charles.Seman at noaa.gov
>> <mailto:Charles.Seman at noaa.gov> should be considered my NOAA
>>
>> email address, not cjs at gfdl.noaa.gov
>> <mailto:cjs at gfdl.noaa.gov>.
>>
>>
>> ********************************************************************
>> Charles Seman
>> Charles.Seman at noaa.gov <mailto:Charles.Seman at noaa.gov>
>>
>> U.S. Department of Commerce / NOAA / OAR
>> Geophysical Fluid Dynamics Laboratory voice:
>> (609) 452-6547
>> 201 Forrestal Road fax:
>> (609) 987-5063
>> Princeton, NJ 08540-6649
>> http://www.gfdl.noaa.gov/~cjs/<http://www.gfdl.noaa.gov/%7Ecjs/>
>> <http://www.gfdl.noaa.gov/%7Ecjs/>
>>
>>
>> ********************************************************************
>>
>> "The contents of this message are mine personally and do
>> not necessarily
>> reflect any position of the Government or NOAA."
>>
>>
>>
>> --
>>
>> Please note that Charles.Seman at noaa.gov
>> <mailto:Charles.Seman at noaa.gov> should be considered my NOAA
>> email address, not cjs at gfdl.noaa.gov <mailto:cjs at gfdl.noaa.gov>.
>>
>> ********************************************************************
>> Charles Seman
>> Charles.Seman at noaa.gov <mailto:Charles.Seman at noaa.gov>
>>
>> U.S. Department of Commerce / NOAA / OAR
>> Geophysical Fluid Dynamics Laboratory voice: (609) 452-6547
>> 201 Forrestal Road fax: (609) 987-5063
>> Princeton, NJ 08540-6649
>> http://www.gfdl.noaa.gov/~cjs/ <http://www.gfdl.noaa.gov/%7Ecjs/> <
>> http://www.gfdl.noaa.gov/%7Ecjs/>
>>
>> ********************************************************************
>>
>> "The contents of this message are mine personally and do not
>> necessarily
>> reflect any position of the Government or NOAA."
>>
>>
>>
>> ************************************************************************
>> *
>> * To use: "grads -p"
>> *
>> * Script calculates and plots zonally averaged cross-sections of the
>> * time averaged height interpolated "mass flux" from the GrADS files
>> * specified below (with an option to overlay the level of minimum
>> * temperature):
>> *
>>
>> ************************************************************************
>> function mass_flux (args)
>>
>> 'reinit'
>>
>>
>> ************************************************************************
>> * Allow external GrADS functions
>>
>> ************************************************************************
>>
>> rc = gsfallow("on")
>>
>>
>> ************************************************************************
>> * Pause between plots?
>>
>> ************************************************************************
>>
>> do_pause = 'no'
>>
>>
>> ************************************************************************
>> * GrADS meta file output?
>>
>> ************************************************************************
>>
>> do_print = 'yes'
>>
>>
>> ************************************************************************
>> * Information for GrADS files
>>
>> ************************************************************************
>>
>> level_dir = '/net2/cjs/ljd/lima/ch3i1/pp/atmos_level/av/monthly_16yr'
>> don_dir = '/net2/cjs/ljd/lima/ch3i1/pp/atmosdon/av/monthly_16yr'
>> ras_dir =
>> '/net2/cjs/lwh/fms/lima/m45_am2p14_ch3i/pp/atmos_level/av/monthly_16yr'
>>
>> sdir_months = 'Aug Sep'
>> month_days = '31 30'
>>
>> level_ctl = 'mc.o_zl.ctl' ; level_var = 'mc'
>> don_ctl1 = 'uceml_deep.o_zl.ctl' ; don_var1 = 'uceml_deep'
>> don_ctl2 = 'umeml_deep.o_zl.ctl' ; don_var2 = 'umeml_deep'
>> don_ctl3 = 'dmeml_deep.o_zl.ctl' ; don_var3 = 'dmeml_deep'
>>
>> ras_ctl = 'mc.o_zl.ctl' ; ras_var = 'mc'
>>
>> temp_ctl = 'temp.o_zl.ctl' ; temp_var = 'temp'
>>
>>
>> ************************************************************************
>> * Information for the plots
>>
>> ************************************************************************
>>
>> plot_min_max = 'no'
>> *
>> * Mass flux variable scale and units
>> *
>> var_scale = 1000
>> var_units = '(g m`a-2`n s`a-1`n)'
>> *
>> * Height levels (scale GrADS ctl file levels from "m" to "km")
>> *
>> zscale = 1./1000.
>> zlevs = '0 18'
>> ylabel_int = 2
>> ylabel = 'Height (km)'
>> *
>> * Define the contour levels...
>> *
>> *contour_levels = '1 3 5 7 9 11 13'
>> *contour_levels = '0 2 4 6 8 10 12 14 16'
>> contour_levels = '0 0.01 0.1 1 3 5 7 10 15'
>> *
>> * Plot level of minimum temperature?
>> *
>> plot_min_temp = 'yes'
>>
>>
>> ************************************************************************
>> * Plot variables' page limits
>>
>> ************************************************************************
>>
>> pvar_vpage.1 = '0 8.5 4.8 9.3'
>> pvar_vpage.2 = '0 8.5 0.3 4.8'
>> pvar_parea = '1 7.0 0.5 4.0'
>>
>> xl_pvar_parea = subwrd(pvar_parea,1)
>> yt_pvar_parea = subwrd(pvar_parea,4)
>>
>>
>> ************************************************************************
>> * Title information for the plots
>>
>> ************************************************************************
>>
>> mod_title = 'FMS lima 1983-1998'
>> var_title = 'Mass Flux 'var_units
>> avg_title = 'Zonal average of weighted time average'
>>
>> figure_panel.1 = '(a)' ; model.1 = 'AM2-D'
>> figure_panel.2 = '(b)' ; model.2 = 'AM2'
>> *
>> * Define names for the "sdir_months" subdirectories...
>> *
>> define_months()
>>
>>
>> ************************************************************************
>> * GrADS metafile for output
>>
>> ************************************************************************
>>
>> if( do_print = 'yes' )
>> 'enable print fig3ab.gx'
>> endif
>>
>>
>> *-----------------------------------------------------------------------
>> * Determine number of months to average...
>>
>> *-----------------------------------------------------------------------
>>
>> n=1
>> month = subwrd(sdir_months,n)
>> if( month = '' )
>> say
>> say 'No months specified in variable "sdir_months"...'
>> say 'Exiting script'
>> say
>> exit
>> else
>> while ( month != '' )
>> n=n+1
>> month = subwrd(sdir_months,n)
>> endwhile
>> nmonths=n-1
>> endif
>>
>> first = subwrd(sdir_months,1) ; last = subwrd(sdir_months,nmonths)
>> months = ''_name.first'-'_name.last''
>>
>>
>> *-----------------------------------------------------------------------
>> * Verify height levels...
>>
>> *-----------------------------------------------------------------------
>>
>> l=1
>> zlv.l = subwrd(zlevs,l)
>> if( zlv.l = '' )
>> say
>> say 'No pressure levels specified in variable "zlevs"...'
>> say 'Exiting script'
>> say
>> exit
>> else
>> while ( zlv.l != '' )
>> l=l+1
>> zlv.l = subwrd(zlevs,l)
>> endwhile
>> nl = l-1
>> if( nl > 2 )
>> say
>> say 'More than two height levels ("nl" = 'nl') specified in
>> variable "zlevs"...'
>> say 'Exiting script'
>> say
>> exit
>> else
>> if( zlv.2 <= zlv.1 )
>> say
>> say 'Height level "2" (= 'zlv.2') not greater than height
>> level "1" (= 'plv.1')...'
>> say 'Exiting script'
>> say
>> exit
>> endif
>> endif
>> endif
>>
>> say
>> say
>> '------------------------------------------------------------------'
>> say 'Calculate zonal averages of weighted time average "mass
>> fluxes"...'
>> say
>> '------------------------------------------------------------------'
>> say
>> say 'Initialize time sum variables for time average...'
>> say
>> *
>> * Use GrADS ctl files for the first "month" to define time sum
>> variables...
>> *
>> month1 = subwrd(sdir_months,1)
>> say
>> say '--------------------'
>> say 'Model: 'model.1
>> say '--------------------'
>> say
>>
>> ctl_file = ''level_dir'/'month1'/'level_ctl
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file_scale = 'mass_flux.scale.ctl'
>> scale_zdef_zlevs(ctl_file,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>>
>> 'query file'
>> rec5 = sublin(result,5)
>> nx = subwrd(rec5,3) ; ny = subwrd(rec5,6) ; nz = subwrd(rec5,9)
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> say ' initialize: wtsum1m'
>> 'define wtsum1m = const('level_var',0,-a)'
>> 'close 1'
>> '!rm 'ctl_file_scale
>>
>> if( plot_min_temp = 'yes' )
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file = ''level_dir'/'month1'/'temp_ctl
>> ctl_file_scale = 'temp.scale.ctl'
>> scale_zdef_zlevs(ctl_file,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> 'query file'
>> rec5 = sublin(result,5)
>> nxt = subwrd(rec5,3) ; nyt = subwrd(rec5,6) ; nzt = subwrd(rec5,9)
>> if( nxt != nx )
>> say '...nxt = 'nxt
>> say '...nx = 'nx
>> exit
>> endif
>> if( nyt != ny )
>> say '...nyt = 'nyt
>> say '...ny = 'ny
>> exit
>> endif
>> if( nzt != nz )
>> say '...nzt = 'nzt
>> say '...nz = 'nz
>> exit
>> endif
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> say ' initialize: wtsum1t'
>> 'define wtsum1t = const('temp_var',0,-a)'
>> 'close 1'
>> '!rm 'ctl_file_scale
>> endif
>>
>> say
>> say '--------------------'
>> say 'Model: 'model.2
>> say '--------------------'
>> say
>>
>> ctl_file = ''ras_dir'/'month1'/'ras_ctl
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file_scale = 'mass_flux.scale.ctl'
>> scale_zdef_zlevs(ctl_file,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> 'query file'
>> rec5 = sublin(result,5)
>> nx = subwrd(rec5,3) ; ny = subwrd(rec5,6) ; nz = subwrd(rec5,9)
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> say ' initialize: wtsum2m'
>> 'define wtsum2m = const('ras_var',0,-a)'
>> 'close 1'
>> '!rm 'ctl_file_scale
>>
>> if( plot_min_temp = 'yes' )
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file = ''ras_dir'/'month1'/'temp_ctl
>> ctl_file_scale = 'temp.scale.ctl'
>> scale_zdef_zlevs(ctl_file,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> 'query file'
>> rec5 = sublin(result,5)
>> nxt = subwrd(rec5,3) ; nyt = subwrd(rec5,6) ; nzt = subwrd(rec5,9)
>> if( nxt != nx )
>> say '...nxt = 'nxt
>> say '...nx = 'nx
>> exit
>> endif
>> if( nyt != ny )
>> say '...nyt = 'nyt
>> say '...ny = 'ny
>> exit
>> endif
>> if( nzt != nz )
>> say '...nzt = 'nzt
>> say '...nz = 'nz
>> exit
>> endif
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> say ' initialize: wtsum2t'
>> 'define wtsum2t = const('temp_var',0,-a)'
>> 'close 1'
>> '!rm 'ctl_file_scale
>> endif
>>
>> if( do_pause = 'yes' )
>> pause()
>> endif
>> *
>> * Calculate running time sum variables for "mass flux" variable...
>> *
>> n=1
>> ndays=0
>>
>> *-----------------------------------------------------------------------
>> while ( n <= nmonths )
>>
>> *-----------------------------------------------------------------------
>> *
>> * Name of the month...
>> *
>> month = subwrd(sdir_months,n)
>> *
>> * Days in the month...
>> *
>> days = subwrd(month_days,n)
>>
>> say
>> say '------------------------'
>> say 'Month = 'month', days = 'days
>> say '------------------------'
>> say
>> say
>> say '--------------------'
>> say 'Model: 'model.1
>> say '--------------------'
>> say
>>
>> ctl_file.1 = ''level_dir'/'month'/'level_ctl
>> ctl_file.2 = ''don_dir'/'month'/'don_ctl1
>> ctl_file.3 = ''don_dir'/'month'/'don_ctl2
>> ctl_file.4 = ''don_dir'/'month'/'don_ctl3
>> *
>> * Make local ctl files with scaled vertical levels...
>> *
>> ctl_file_scale.1 = 'mass_flux.1.scale.ctl'
>> scale_zdef_zlevs(ctl_file.1,zscale,ctl_file_scale.1)
>> ctl_file_scale.2 = 'mass_flux.2.scale.ctl'
>> scale_zdef_zlevs(ctl_file.2,zscale,ctl_file_scale.2)
>> ctl_file_scale.3 = 'mass_flux.3.scale.ctl'
>> scale_zdef_zlevs(ctl_file.3,zscale,ctl_file_scale.3)
>> ctl_file_scale.4 = 'mass_flux.4.scale.ctl'
>> scale_zdef_zlevs(ctl_file.4,zscale,ctl_file_scale.4)
>>
>> say 'Open GrADS file 1:'
>> say '...'ctl_file_scale.1
>> 'open 'ctl_file_scale.1
>> say 'Open GrADS file 2:'
>> say '...'ctl_file_scale.2
>> 'open 'ctl_file_scale.2
>> say 'Open GrADS file 3:'
>> say '...'ctl_file_scale.3
>> 'open 'ctl_file_scale.3
>> say 'Open GrADS file 4:'
>> say '...'ctl_file_scale.4
>> 'open 'ctl_file_scale.4
>> say 'GrADS files ready'
>>
>> say
>> say ' calculate weighted (using days = 'days') "mass flux" time
>> sum'
>> say ' ('level_var'.1+'don_var1'.2+'don_var2'.3+'don_var3'.4)'
>> say ' for current month and add it to the time sum variable...'
>> say
>> 'set dfile 1'
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define wtsum1m = wtsum1m +
>> 'days'*('level_var'.1+'don_var1'.2+'don_var2'.3+'don_var3'.4)'
>>
>> 'close 4' ; 'close 3' ; 'close 2' ; 'close 1'
>> '!rm 'ctl_file_scale.1
>> '!rm 'ctl_file_scale.2
>> '!rm 'ctl_file_scale.3
>> '!rm 'ctl_file_scale.4
>>
>> if( plot_min_temp = 'yes' )
>>
>> ctl_file = ''level_dir'/'month'/'temp_ctl
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file_scale = 'temp.scale.ctl'
>> scale_zdef_zlevs(ctl_file,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> say 'GrADS file ready'
>> say
>> say ' calculate weighted (using days = 'days') temperature
>> time sum'
>> say ' for current month and add it to the time sum variable...'
>> say
>>
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define wtsum1t = wtsum1t + 'days'*'temp_var
>>
>> 'close 1'
>> '!rm 'ctl_file_scale
>>
>> endif
>>
>> say
>> say '--------------------'
>> say 'Model: 'model.2
>> say '--------------------'
>> say
>>
>> ctl_file = ''ras_dir'/'month'/'ras_ctl
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file_scale = 'mass_flux.scale.ctl'
>> scale_zdef_zlevs(ctl_file,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> say 'GrADS file ready'
>>
>> say
>> say ' calculate weighted (using days = 'days') "mass flux" time
>> sum'
>> say ' ('ras_var'.1)'
>> say ' for current month and add it to the time sum variable...'
>> say
>> 'query file'
>> rec5 = sublin(result,5)
>> nx = subwrd(rec5,3) ; ny = subwrd(rec5,6) ; nz = subwrd(rec5,9)
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define wtsum2m = wtsum2m + 'days'*'ras_var
>>
>> 'close 1'
>> '!rm 'ctl_file_scale
>>
>> if( plot_min_temp = 'yes' )
>>
>> ctl_file = ''ras_dir'/'month'/'temp_ctl
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file_scale = 'temp.scale.ctl'
>> scale_zdef_zlevs(ctl_file,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> say 'GrADS file ready'
>> say
>> say ' calculate weighted (using days = 'days') temperature
>> time sum'
>> say ' for current month and add it to the time sum variable...'
>> say
>>
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define wtsum2t = wtsum2t + 'days'*'temp_var
>>
>> 'close 1'
>> '!rm 'ctl_file_scale
>>
>> endif
>>
>> say
>> say '...finished month'
>> say
>>
>> if( do_pause = 'yes' )
>> pause()
>> endif
>>
>> n=n+1
>> ndays=ndays+days
>>
>> *-----------------------------------------------------------------------
>> endwhile
>>
>> *-----------------------------------------------------------------------
>> say
>> say '---------------------------'
>> say 'Finished weighted time sums'
>> say '---------------------------'
>> say
>> say
>> say
>> '------------------------------------------------------------------'
>> say 'Calculate zonal averages of the weighted time averages'
>> say 'of the weighted time sums, and plot the data...'
>> say
>> '------------------------------------------------------------------'
>> say
>> say 'Use GrADS ctl files for the first "month" to define time
>> averaged variables for plotting...'
>> *
>> * Model #1
>> *
>> ctl_file.1 = ''level_dir'/'month1'/'level_ctl
>> ctl_file_temp.1 = ''level_dir'/'month1'/'temp_ctl
>> *
>> * Model #2
>> *
>> ctl_file.2 = ''ras_dir'/'month1'/'ras_ctl
>> ctl_file_temp.2 = ''ras_dir'/'month1'/'temp_ctl
>> *
>> * Initialize min/max variables...
>> *
>> min = 1e6 ; max = -1e6
>> *
>> * Start a new frame, and plot a title at top of page...
>> *
>> ptitle_top(var_title,avg_title,months,pvar_parea)
>>
>> m=1
>> while ( m <= 2 )
>>
>> say
>> say '--------------------'
>> say 'Model: 'model.m
>> say '--------------------'
>> say
>>
>> if( plot_min_temp = 'yes' )
>> say
>> say 'Find and save level of minimum temperature for the'
>> say 'zonal and weighted time average temperature...'
>> say
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file_scale = 'temp.scale.ctl'
>> scale_zdef_zlevs(ctl_file_temp.m,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> say 'GrADS file ready'
>> say
>> 'query file'
>> rec5 = sublin(result,5)
>> nx = subwrd(rec5,3) ; ny = subwrd(rec5,6) ; nz = subwrd(rec5,9)
>>
>> say ' calculate weighted time average (using ndays = 'ndays')...'
>>
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define tav = wtsum'm't/'ndays
>>
>> say ' calculate zonal averages of the weighted time average...'
>>
>> 'set x 1'
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define ztav = ave(tav,x=1,x='nx')'
>>
>> say ' plot zonal average of the weighted time average...'
>>
>> 'set vpage 'pvar_vpage.m
>> 'set parea 'pvar_parea
>> 'set grads off'
>> 'set gxout contour'
>> 'set clab off'
>> 'set ccolor 0'
>> 'set xlab off'
>> 'set ylab off'
>> 'd ztav'
>>
>> say ' find and save location of minimum temperature...'
>>
>> j=1
>> while ( j <= ny )
>> 'set y 'j
>> 'set gxout print'
>> 'd minloc(ztav,z=1,z='nz')' ; rec3 = sublin(result,3) ; zmin
>> = subwrd(rec3,1)
>> * say result
>> say ' ...j, zmin = 'j', 'zmin
>> 'query gr2xy 'j' 'zmin
>> rec = sublin(result,1) ; xloc.j = subwrd(rec,3) ; yloc.j =
>> subwrd(rec,6)
>> * say result
>> say ' xloc, yloc = 'xloc.j', 'yloc.j
>> j=j+1
>> endwhile
>>
>> 'undefine wtsum'm't' ; 'undefine tav' ; 'undefine ztav'
>> 'close 1'
>> '!rm 'ctl_file_scale
>>
>> endif
>>
>> say
>> say 'Calculate and plot zonally averaged "mass flux"...'
>> say
>> *
>> * Make local ctl file with scaled vertical levels...
>> *
>> ctl_file_scale = 'mass_flux.scale.ctl'
>> scale_zdef_zlevs(ctl_file.m,zscale,ctl_file_scale)
>>
>> say 'Open GrADS file:'
>> say '...'ctl_file_scale
>> 'open 'ctl_file_scale
>> say 'GrADS file ready'
>> say
>> say
>> 'query file'
>> rec5 = sublin(result,5)
>> nx = subwrd(rec5,3) ; ny = subwrd(rec5,6) ; nz = subwrd(rec5,9)
>>
>> say ' calculate weighted "mass flux" time average (using ndays
>> = 'ndays')'
>>
>> 'set x 1 'nx
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define tav = wtsum'm'm/'ndays
>>
>> say ' calculate zonal averages of the weighted time average...'
>>
>> 'set x 1'
>> 'set y 1 'ny
>> 'set z 1 'nz
>> 'set t 1'
>> 'define ztav = 'var_scale'*ave(tav,x=1,x='nx')'
>>
>> say ' plot zonal averages of time averaged "mass flux"'
>> 'set vpage 'pvar_vpage.m
>> 'set parea 'pvar_parea
>> 'set grads off'
>> 'set clevs 'contour_levels
>> 'set gxout shaded'
>> 'set xlab on'
>> 'set ylab on'
>> 'set ylint 'ylabel_int
>> 'set x 1'
>> 'set y 1 'ny
>> 'set lev 'zlevs
>> 'set t 1'
>> 'd ztav'
>> 'run cbar.gs <http://cbar.gs> 0.7'
>>
>> 'draw title \\\ 'model.m
>> * 'draw xlab Latitude \\\'
>> 'draw ylab \\\ 'ylabel
>> xfp = xl_pvar_parea + 0.1
>> yfp = yt_pvar_parea + 0.2
>> 'set string 1 l 6'
>> 'set strsiz 0.13'
>> 'draw string 'xfp' 'yfp' 'figure_panel.m
>> 'set gxout stat'
>> 'd ztav'
>> min_max = sublin(result,8)
>> say ' 'min_max
>> min = subwrd(min_max,4)
>> max = subwrd(min_max,5)
>> if( plot_min_max = 'yes' )
>> xmm = xl_pvar_parea + 0.1
>> ymm = yt_pvar_parea - 0.1
>> 'set string 1 l 4'
>> 'set strsiz 0.08'
>> 'draw string 'xmm' 'ymm' 'min_max' '
>> endif
>>
>> if( plot_min_temp = 'yes' )
>> say ' plot level of minimum temperature... '
>> j=1
>> while ( j <= ny )
>> 'set y 'j
>> * say ' xloc, yloc = 'xloc.j', 'yloc.j
>> 'set line 1'
>> 'draw mark 3 'xloc.j' 'yloc.j' 0.03'
>> j=j+1
>> endwhile
>> endif
>>
>> 'undefine wtsum'm'm' ; 'undefine tav' ; 'undefine ztav'
>> 'close 1'
>> '!rm 'ctl_file_scale
>>
>> say
>> say '...finished plot for: 'model.m
>> say
>> m=m+1
>> endwhile
>>
>> say
>> say '--------------'
>> say 'Finished plots'
>> say '--------------'
>> say
>>
>> if( do_print = 'yes' )
>> 'print'
>> 'disable print'
>> endif
>>
>> say
>> say '--------------------------'
>> say '"Mass flux" min/max values'
>> say '--------------------------'
>> say
>> say 'min = 'min' max = 'max
>>
>> say
>> say '************************************************'
>> say ' Finished with this script.'
>> say '************************************************'
>> say
>>
>> if( batch_mode = 'yes' )
>> quit
>> endif
>>
>> return
>>
>>
>> ************************************************************************
>> *
>> function define_months()
>> *
>> * Define "names" for all months...
>> *
>> _name.Jan = 'January'
>> _name.Feb = 'February'
>> _name.Mar = 'March'
>> _name.Apr = 'April'
>> _name.May = 'May'
>> _name.Jun = 'June'
>> _name.Jul = 'July'
>> _name.Aug = 'August'
>> _name.Sep = 'September'
>> _name.Oct = 'October'
>> _name.Nov = 'November'
>> _name.Dec = 'December'
>> *
>> return
>> *
>>
>> ************************************************************************
>>
>>
>> ************************************************************************
>> *
>> function ptitle_top(title,subtitle,months,plim)
>> *
>> * Plots a title at the top of the page
>> *
>> x1 = subwrd(plim,1)
>> x2 = subwrd(plim,2)
>> xc = x1 + (x2-x1)/2
>> 'c'
>> 'set vpage 0 8.5 0 11'
>> 'set string 1 c 6'
>> 'set strsiz 0.14'
>> 'draw string 'xc' 10.2 'title
>> *'draw string 'xc' 10.2 'subtitle
>> 'draw string 'xc' 9.8 'months
>> *
>> return
>> *
>>
>> ************************************************************************
>>
>>
>>
>>
>> --
>> Sushant Puranik
>> Junior Research Fellow
>> Dept. of Atmospheric & Space Sciences,
>> University of Pune,
>> Pune-07,
>> India.
>>
>
> --
>
> Please note that Charles.Seman at noaa.gov should be considered my NOAA
> email address, not cjs at gfdl.noaa.gov.
>
> ********************************************************************
> Charles Seman Charles.Seman at noaa.gov
> U.S. Department of Commerce / NOAA / OAR
> Geophysical Fluid Dynamics Laboratory voice: (609) 452-6547
> 201 Forrestal Road fax: (609) 987-5063
> Princeton, NJ 08540-6649 http://www.gfdl.noaa.gov/~cjs/<http://www.gfdl.noaa.gov/%7Ecjs/>
> ********************************************************************
>
> "The contents of this message are mine personally and do not necessarily
> reflect any position of the Government or NOAA."
>
--
Sushant Puranik
Junior Research Fellow
Dept. of Atmospheric & Space Sciences,
University of Pune,
Pune-07,
India.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://gradsusr.org/pipermail/gradsusr/attachments/20091127/c828cc4c/attachment.html
More information about the gradsusr
mailing list