Re: pt&ept

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Thu May 27 2010 - 10:52:12 MDT

Hello,

Please ... do not address emails directly to me. When you do address me
directly, no one else will respond. I am very busy and can not keep
up with personal emails.

[1]
As I have noted many times on ncl-talk, the golden rule of data
processing is to *look at your data before processing*. The air and shum
variables have different level sizes:

netcdf air.2001 {
dimensions:
         lon = 144 ;
         lat = 73 ;
         level = 17 ; <*********************************
         time = UNLIMITED ; // (365 currently)

netcdf shum.2001 {
dimensions:
         lon = 144 ;
         lat = 73 ;
         level = 8 ; <*********************************
         time = UNLIMITED ; // (365 currently)

This means, they do not "conform".

http://www.ncl.ucar.edu/Document/Functions/Built-in/conform.shtml
http://www.ncl.ucar.edu/Document/Functions/Built-in/conform_dims.shtml

[2]
It looks to me like you took pieces from several scripts and stuck
them into a script. You read in the data twice using different variable
names. The second time you did not unpack the values. Really, you must
think about what you want to do and how the programming steps flow
into each other. Use printVarSummary" to look at your variables.

[3]
Instead of doing a computation you set a string value:

    T_e = "/T+(L_v/c_pd)*X/" ; common approximation
    theta_e = "/T_e*(p0/P)^Rcpd/"

Why did you use quotes[ "..."]? This just made T_e and theta_e string
variables. No computation were performed.

[4]
I have no idea why the following was being used:
    THETAE =
zonal_mpsi(aavg(::-1,:,{30:70}),aavg&lat,lev_p(::-1),shavg(:,{30:70}))

[5]
Attached is a script that seems to do what you want.

Good luck

On 05/26/2010 11:39 AM, wei huang wrote:
> Dear Respected Dr.Dennis Shea,
>
> Thanks for your help. I am interesting to plot the figure Press/height
> vs. latitude (like the ncl example (h_lat_2.ncl). My code is attached
> herewith, I'm confused about the function used in the code. Thanks
> advance for your help.
>
> regrads,
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
>
> function theta_eqv(T:numeric, X:numeric, p0[1], p[*]:numeric)
>
> local T_e, L_v, c_pd, R, Rcpd, dimT, dimX, P
>
> begin
>
> L_v = 2400. ; Latent heat of evaporation [approx 2400 kJ/kg at 25C]
>
> c_pd = 1004. ; Specific heat at constant pressure for air [approx 1004
> J/(kg-K)]
>
> R = 287. ; Specific gas constant for air [J/(kg-K)]
>
> Rcpd = R/c_pd
>
> dimT = dimsizes(T)
>
> dimX = dimsizes(X)
>
> if (dimsizes(dimT) .ne. dimsizes(dimX)) then
>
> print("theta_eqv: rank of T .ne. rank of X")
>
> exit
>
> end if
>
> P = conform(T,p,1) ; make p same shape/rank/size as T
>
> T_e = "/T+(L_v/c_pd)*X/" ; common approximation
>
> theta_e = "/T_e*(p0/P)^Rcpd/"
>
> theta_e@long_name = "Equivalent potential temperature"
>
> theta_e@units = "K"
>
> ;copy_VarCoords(T, theta_e) ; assign coordinates
>
> return(theta_e)
>
> end
>
> begin
>
> diri = "/data/ncep_data/pressure/"
>
> f = addfile(diri+"air.2001.nc <http://air.2001.nc>", "r")
>
> a = short2flt(f->air(0:29,:,:,:))
>
> f1 = addfile(diri+"shum.2001.nc <http://shum.2001.nc>","r")
>
> sh = short2flt(f1->shum(0:29,:,:,:))
>
> aavg = dim_avg_Wrap(a(level|:,lat|:,lon|:,time|:))
>
> shavg = dim_avg_Wrap(sh(level|:,lat|:,lon|:,time|:))
>
> printVarSummary(shavg)
>
> p = f->level
>
> T = f->air
>
> x = f1->shum ; spc humidity
>
> X = x/(1-x) ; convert spc humidity to mixing ratio
>
> X@long_name = "mixing ratio"
>
> p0 = 1000
>
> p0@units = "hPa" ; must match units of p
>
> THETAE = theta_eqv(T, X, p0, p)
>
> printVarSummary( THETAE )
>
> ; printMinMax(THETAE, True )
>
> lat = f->lat
>
> lev_p = 100*(f->level)
>
> lev_p@units = "Pa"
>
> THETAE =
> zonal_mpsi(aavg(::-1,:,{30:70}),aavg&lat,lev_p(::-1),shavg(:,{30:70}))
>
> THETAE!0 = "lev"
>
> THETAE!1 = "lat"
>
> THETAE&lev = lev_p
>
> THETAE&lat = lat
>
> THETAE@long_name = "ept"
>
> wks = gsn_open_wks ("x11", "pot" )
>
> plot = gsn_csm_pres_hgt(wks,THETAE(:,{10:20}), False )
>
> end
>
> Variable: THETAE
>
> Type: string
>
> Total Size: 4 bytes
>
> 1 values
>
> Number of Dimensions: 1
>
> Dimensions and sizes: [1]
>
> Coordinates:
>
> Number Of Attributes: 2
>
> long_name : Equivalent potential temperature
>
> units : K
>
> fatal:Number of subscripts do not match number of dimensions of
> variable,(2) Subscripts used, (3) Subscripts expected
>
> fatal:Execute: Error occurred at or near line 56 in file pot1.ncl
>
>
>
> On Wed, May 26, 2010 at 11:15 PM, Dennis Shea <shea@ucar.edu
> <mailto:shea@ucar.edu>> wrote:
>
> This question is ambiguous. Please be specific.
>
> [1] Do you mean a plot over a map at a specific level?
>
> See Example 2 at: http://www.ncl.ucar.edu/Applications/iso.shtml
>
> This also includes an example of calculating potential temperature
> and equivalent potential temperature.
>
> [2] Do you mean pres/height vs lat or lon or time or ?
>
>
> http://www.ncl.ucar.edu/Applications/
>
> Click:
> Press/height vs. latitude
> Press/height vs. longitude
> Press/height vs. time
>
> =============
>
>
>
>
> On 5/25/10 7:57 PM, wei huang wrote:
>
> Dear All,
> Does anyone could share me a script or guided me for plotting the
> potential temperature at a height level (from surface to 200mb)
> and also
> equivalent potential temperature? Your help will be appreciated.
>
> Thanks in advance
>
> Wei
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>

-- 
======================================================
Dennis J. Shea                  tel: 303-497-1361    |
P.O. Box 3000                   fax: 303-497-1333    |
Climate Analysis Section                             |
Climate & Global Dynamics Div.                       |
National Center for Atmospheric Research             |
Boulder, CO  80307                                   |
USA                        email: shea 'at' ucar.edu |
======================================================


_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk

Received on Thu May 27 10:52:18 2010

This archive was generated by hypermail 2.1.8 : Tue Jun 01 2010 - 09:12:20 MDT