[ncl-talk] Fwd: hyi2hyo: rightmost dimension of ps and xi must be the same

Kevin Hallock hallock at ucar.edu
Wed Nov 1 17:28:57 MDT 2017


Hi Jayant,

I found several lines of code that I think might be contributing to this issue.

First of all, the values that have been assigned in this snippet of code might not be what you expect them to be:
  dNames = getvardims(xo)
  time  = dNames(0)
  lev   = dNames(1)
  lat   = dNames(2)
  lon   = dNames(3)

If you print() any of the variables “time”, “lev”, “lat”, or “lon”, the result will just be the scalar string “time”, “lev”, “lat”, or “lon”, respectively. Based on the rest of the code, I suspect we really want those variables to contain the coordinate arrays from the variable “xo”.

When NCL reaches the line:
filevardef(fout, "time" ,typeof(time),getvardims(time))

“typeof(time)” will be “string” as mentioned above, and “getvardims(time)” will actually be a missing value. Really what we want is the coordinate array for “time”, which should show “typeof(time)” as “double” and “getvardims(time)” as “time”.

Since “time”, “lev”, “lat”, and “lon” are all coordinate arrays associated with the dimensions of the variable “xo”, we can instead access them using NCL’s coordinate array operator ‘&’:
  time  = xo&time
  lev   = xo&lev
  lat   = xo&lat
  lon   = xo&lon

After making this change, the script should run without showing any error messages. However, there is one other issue I found near the end of your script when assigning variables in the output data file:
  fout->time   = (/time/)
  fout->lev    = (/lev/)
  fout->lat    = (/lat/)
  fout->lon    = (/lon/)
  fout->xo     = (/xo/)

In this case, “fout->xo” will be assigned only the data array of the “xo” script variable. Previously in the script, the metadata from the script variable “xo” was assigned to the output file variable named “OX_VMR_avrg”:
  filevarattdef(fout,"OX_VMR_avrg",xo)

Essentially what this means is that the output file will have a variable named “xo” containing the data and a variable named “OX_VMR_avrg” containing only the metadata.

Changing the line:
  fout->xo     = (/xo/)
to instead be:
  fout->OX_VMR_avrg     = (/xo/)

should eliminate the unnecessary “xo” variable in the output file.

I hope this helps!
Kevin

> Begin forwarded message:
> From: Jayant <jayantkp2979 at gmail.com <mailto:jayantkp2979 at gmail.com>>
> Date: Tue, Oct 31, 2017 at 9:06 PM
> Subject: Re: [ncl-talk] hyi2hyo: rightmost dimension of ps and xi must be the same
> To: Dennis Shea <shea at ucar.edu <mailto:shea at ucar.edu>>
> Cc: "ncl-talk at ucar.edu <mailto:ncl-talk at ucar.edu>" <ncl-talk at ucar.edu <mailto:ncl-talk at ucar.edu>>
> 
> 
> Hi Denis,
> Thank you for the corrections. I am expanding the script to write the output in a netcdf file.
> I am getting error with 'filevardef' function. Can you help me fix teh error? Below I paste the complete code.
> Thanks in advance,
> Jayant
> ----------------------------
> the error code
> fatal:["Execute.c":8575]:Execute: Error occurred at or near line 103 in file scr_vregrid_ham.ncl
> ----------------------------
> 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"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/esmf/ESMF_regridding.ncl"
> ;
> begin
> ;
> 
>   f    = addfile("ham_oxi_aps_T63L31.nc","r")
>   hyai = f->hyam          ; input A(k)
>   hybi = f->hybm          ; input B(k)
>   xi   = f->OX_VMR_avrg   ; variable to be interpolated
>   printVarSummary(xi)     ; (time,lev,lon,lat);  (0,1,2,3)
>   print("---")
> 
>   xi  := xi(time|:,lev|:,lat|:,lon|:)  ; (time,lev,lat,lon); REORDER
>   printVarSummary(xi)     ; (time,lev,lat,lon);  (0,1,2,3)
>   print("---")
> 
>   xi   = dim_rmvmean_Wrap(xi)        ; remove zonal means
>  ;xi   = dim_rmvmean_n_Wrap(xi, 0)   ; remove time mean
>   printVarSummary(xi)          ; rmvmean_op_NCL:  dim_rmvmean over dimension(s): lon
>   print("---")
> 
>   fP   = addfile("surf_pres_new_exp_echam_exp.nc <http://surf_pres_new_exp_echam_exp.nc/>","r")
>   psfc = fP->aps          ; surface pressure (Pa)
>   printVarSummary(psfc)   ; (time,lat,lon)
>   print("---")
> 
>   PSFC = conform( xi(:,0,:,:), psfc(0,:,:), (/1,2/) )
>   copy_VarMeta(psfc(0,:,:), PSFC(0,:,:))
>   printVarSummary(PSFC)
> 
>   sigma = (/0.00657, 0.01386, 0.02309, 0.03469, 0.04920, 0.06723, 0.08945, 0.11654, 0.14916, \
>          0.18783, 0.23286, 0.28421, 0.34137, 0.40334, 0.46860, 0.53529, 0.60135, 0.66482,    \
>          0.72401, 0.77773, 0.82527, 0.86642, 0.90135, 0.93054, 0.95459, 0.97418, 0.99000,    \
>          1.00000/)
> 
>   hybo = sigma            ; *your* B(k) [ sigma levels ]
>   nsig = dimsizes(hybo)
>   hyao = new(nsig,typeof(hybo))
>   hyao = 0.0             ; set output A(k) to 0.0
> 
>   p0   = 101325           ; reference pressure (Pa)
>   p0 at units = "Pa"
> 
>   xo   = hyi2hyo_Wrap(p0,hyai,hybi,PSFC,xi,hyao,hybo,0) ; contributed.ncl
>   printVarSummary(xo)
>   ;
>   ; writing output in netdf
>   ;
>   dNames = getvardims(xo)
>   time  = dNames(0)
>   lev   = dNames(1)
>   lat   = dNames(2)
>   lon   = dNames(3)
> 
>   ntim  = dimsizes(time)                 ; get dimension sizes  
>   klev  = dimsizes(lev)                                               
>   nlat  = dimsizes(lat)  
>   nlon  = dimsizes(lon)      
> 
>   diro = "./"                     ; Output directory
>   filo = "example.nc <http://example.nc/>"             ; Output file
>   system("/bin/rm -f " + diro + filo)    ; remove if exists
>   fout  = addfile (diro + filo, "c")  ; open output file
> 
>   ; explicitly declare file definition mode. Improve efficiency.
>   setfileoption(fout,"DefineMode",True)
> 
>   ; create global attributes of the file
>   fAtt               = True            ; assign file attributes
>   fAtt at title         = "NCL Efficient Approach to netCDF Creation"  
>   fAtt at source_file   = "ham_oxi_aps_T63L31.nc"
>   fAtt at Conventions   = "None"   
>   fAtt at creation_date = systemfunc ("date")        
>   fileattdef( fout, fAtt )            ; copy file attributes
> 
>   ;===================================================================
>   ; predefine the coordinate variables and their dimensionality
>   ; Note: to get an UNLIMITED record dimension, we set the dimensionality
>   ; to -1 (or the actual size) and set the dimension name to True.
>   ;===================================================================
>   dimNames = (/"time", "lat", "lon", "lev"/)  
>   dimSizes = (/ ntim ,  nlat,  nlon, klev /) 
>   dimUnlim = (/ True , False, False, False/)   
>   filedimdef(fout,dimNames,dimSizes,dimUnlim)
> 
>   ;===================================================================
>   ; predefine the the dimensionality of the variables to be written out
>   ;===================================================================
>   ; Here we are using NCL functions to facilitate defining 
>   ; each variable's dimension name(s) and type. 
>   ; The following could be replaced with explicit, user defined dimension 
>   ; names different from those associated with the variable in memory. 
>   ; Say, PS(time,lat,lon) in the NCL script. They could be redefined for the file via: 
>   ; filevardef(fout, "PS"   ,typeof(PS) ,(/"TIME","latitude","longitude"/)) 
>   ;===================================================================
>   filevardef(fout, "time" ,typeof(time),getvardims(time)) 
>   filevardef(fout, "lev"  ,typeof(lev),getvardims(lev) )                           
>   filevardef(fout, "lat"  ,typeof(lat),getvardims(lat))                          
>   filevardef(fout, "lon"  ,typeof(lon),getvardims(lon))                          
>   filevardef(fout, "OX_VMR_avrg"    ,typeof(xo)  ,getvardims(xo))
> 
>   ;===================================================================
>   ; Copy attributes associated with each variable to the file
>   ; All attributes associated with each variable will be copied.
>   ;====================================================================
>   filevarattdef(fout,"OX_VMR_avrg",xo)                ; copy T attributes
>   filevarattdef(fout,"time" ,time)                    ; copy time attributes
>   filevarattdef(fout,"lev"  ,lev)                     ; copy lev attributes
>   filevarattdef(fout,"lat"  ,lat)                     ; copy lat attributes
>   filevarattdef(fout,"lon"  ,lon)                     ; copy lon attributes    
>   
>   ;===================================================================
>   ; explicitly exit file definition mode. **NOT REQUIRED**
>   ;===================================================================
>   setfileoption(fout,"DefineMode",False)
> 
>   ;===================================================================
>   ; output only the data values since the dimensionality and such have
>   ; been predefined. The "(/", "/)" syntax tells NCL to only output the
>   ; data values to the predefined locations on the file.
>   ;====================================================================
>   fout->time   = (/time/)     
>   fout->lev    = (/lev/)
>   fout->lat    = (/lat/)
>   fout->lon    = (/lon/) 
>   fout->xo      = (/xo/)
> 
> end
> 
> 
> On Tue, Oct 31, 2017 at 6:17 PM, Dennis Shea <shea at ucar.edu <mailto:shea at ucar.edu>> wrote:
> I think the refrain 'look at your data' applies here.
> 
> ===
> There are several dimension related issues with your data.
> 
> [A] Variable xi
> 
> xi:  [time | 12] x [lev | 31] x [lon | 192] x [lat | 96]
> 
> Note the two rightmost 'spatial' dimensions: names and sizes
> 
> [B] Variable psfc
> 
> psfc: [time | 1] x [lat | 96] x [lon | 192]
> 
> ----------------------
> 
> spatial dimensions:
> xi has order (...,lon,lat)     ==> (...,192,96)
> psfc gas order (...,lat,lon)   ==> (...,96,192)
> 
> temporal dimension size
> xi has time=12, psfc=1
> 
> -----------------------
> 
> You could expand psfc 12 size 12 but I am not sure that is wht you want to do.
> 
> [C]
> 
> P0 must have a units attribute
> 
> [D]
> 
> Did you want to remove the zonal means or the time mena???
> 
> ++++++++++++++++++++++++++++++++++++++++++++
> 
>   f    = addfile("ham_oxi_aps_T63L31.nc","r")
>   hyai = f->hyam          ; input A(k)
>   hybi = f->hybm          ; input B(k)
>   xi   = f->OX_VMR_avrg   ; variable to be interpolated
>   printVarSummary(xi)     ; (time,lev,lon,lat);  (0,1,2,3)
>   print("---")
> 
>   xi  := xi(time|:,lev|:,lat|:,lon|:)  ; (time,lev,lat,lon); REORDER
>   printVarSummary(xi)     ; (time,lev,lat,lon);  (0,1,2,3)
>   print("---")
> 
>   xi   = dim_rmvmean_Wrap(xi)        ; remove zonal means
>  ;xi   = dim_rmvmean_n_Wrap(xi, 0)   ; remove time mean
>   printVarSummary(xi)          ; rmvmean_op_NCL:  dim_rmvmean over dimension(s): lon
>   print("---")
> 
>   fP   = addfile("surf_pres_new_exp_echam_exp.nc <http://surf_pres_new_exp_echam_exp.nc/>","r")
>   psfc = fP->aps          ; surface pressure (Pa)
>   printVarSummary(psfc)   ; (time,lat,lon)
>   print("---")
> 
>   PSFC = conform( xi(:,0,:,:), psfc(0,:,:), (/1,2/) )
>   copy_VarMeta(psfc(0,:,:), PSFC(0,:,:))
>   printVarSummary(PSFC)
> 
>   sigma = (/0.00657, 0.01386, 0.02309, 0.03469, 0.04920, 0.06723, 0.08945, 0.11654, 0.14916, \
>          0.18783, 0.23286, 0.28421, 0.34137, 0.40334, 0.46860, 0.53529, 0.60135, 0.66482,    \
>          0.72401, 0.77773, 0.82527, 0.86642, 0.90135, 0.93054, 0.95459, 0.97418, 0.99000,    \
>          1.00000/)
> 
>   hybo = sigma            ; *your* B(k) [ sigma levels ]
>   nsig = dimsizes(hybo)
>   hyao = new(nsig,typeof(hybo))
>   hyao = 0.0             ; set output A(k) to 0.0
> 
>   p0   = 101325           ; reference pressure (Pa)
>   p0 at units = "Pa"
> 
>   xo   = hyi2hyo_Wrap(p0,hyai,hybi,PSFC,xi,hyao,hybo,0) ; contributed.ncl
>   printVarSummary(xo)
> 
> 
> 
> On Tue, Oct 31, 2017 at 11:01 AM, Jayant <jayantkp2979 at gmail.com <mailto:jayantkp2979 at gmail.com>> wrote:
> Dear,
> 
> I am following the steps given in thread (https://www.ncl.ucar.edu/Support/talk_archives/2006/0001.html <https://www.ncl.ucar.edu/Support/talk_archives/2006/0001.html>) to convert a file from hybrid level to sigma levels. I paste below my script and the outputs. I am not able to understand where I am going wrong. Please advice. Also, let me know if I need to send my netcdf files as well.
> Thanks in advance,
> Jayant
> 
>   f=addfile("ham_oxi_aps_T63L31.nc","r")
>   hyai = f->hyam         ; input A(k)
>   hybi = f->hybm         ; input B(k) 
>   xi   = f->OX_VMR_avrg        ; (0,:,:,:)     ; variable to be interpolated
>   dims = getvardims(xi)
>   xin  = dim_rmvmean_Wrap(xi($dims(0)$|:,$dims(1)$|:,$dims(3)$|:,$dims(2)$|:))
>   printVarSummary(xin)
>   p0   = 101325            ; reference pressure (Pa)
> 
>   fP   = addfile("surf_pres_new_exp_echam_exp.nc <http://surf_pres_new_exp_echam_exp.nc/>","r")
>   psfc = fP->aps            ; surface pressure (Pa)
>   printVarSummary(psfc)
>   ; dims = getvardims(psfc)
>   ; psfcr = dim_rmvmean_Wrap(psfc($dims(0)$|:,$dims(2)$|:,$dims(1)$|:))
>   ; psfc2d = psfcr(0,:,:)
> 
>   sigma = (/0.00657, 0.01386, 0.02309, 0.03469, 0.04920, 0.06723, 0.08945, 0.11654, 0.14916, \
>          0.18783, 0.23286, 0.28421, 0.34137, 0.40334, 0.46860, 0.53529, 0.60135, 0.66482,    \
>          0.72401, 0.77773, 0.82527, 0.86642, 0.90135, 0.93054, 0.95459, 0.97418, 0.99000,    \
>          1.00000/)
> 
>   hybo = sigma            ; *your* B(k) [ sigma levels ] 
>   nsig = dimsizes(hybo)
>   hyao = new(nsig,typeof(hybo))
>   hyao = 0.0             ; set output A(k) to 0.0 
> 
>   ; xo = hyi2hyo_Wrap(p0,hyai,hybi,psfcr,xi,hyao,hybo,0) ; contributed.ncl 
>   xo = hyi2hyo_Wrap(p0,hyai,hybi,psfc,xin,hyao,hybo,0) ; contributed.ncl 
>   printVarSummary(xo)
> 
> On running, I am getting the following error:
> 
> Variable: xin
> Type: float
> Total Size: 27426816 bytes
>             6856704 values
> Number of Dimensions: 4
> Dimensions and sizes:    [time | 12] x [lev | 31] x [lat | 96] x [lon | 192]
> Coordinates: 
>             time: [   0..  11]
>             lev: [ 1..31]
>             lat: [88.5721685140073..-88.5721685140073]
>             lon: [   0..358.125]
> Number Of Attributes: 5
>   rmvmean_op_NCL :    dim_rmvmean over dimension(s): lon
>   long_name :    Deviation from mean
>   grid_type :    gaussian
>   units :    VMR
>   _FillValue :    9.96921e+36
> 
> Variable: psfc
> Type: float
> Total Size: 73728 bytes
>             18432 values
> Number of Dimensions: 3
> Dimensions and sizes:    [time | 1] x [lat | 96] x [lon | 192]
> Coordinates: 
>             time: [30.99166666666667..30.99166666666667]
>             lat: [88.57216851400727..-88.57216851400727]
>             lon: [   0..358.125]
> Number Of Attributes: 5
>   grid_type :    gaussian
>   table :    128
>   code :    134
>   units :    Pa
>   long_name :    surface pressure
> fatal:hyi2hyo: The rightmost dimensions of 'ps' and 'xi' must be the same
> fatal:["Execute.c":8575]:Execute: Error occurred at or near line 9139 in file $NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl
> 
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu <mailto:ncl-talk at ucar.edu>
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk <http://mailman.ucar.edu/mailman/listinfo/ncl-talk>
> 
> 
> 

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171101/0979b9ce/attachment.html>


More information about the ncl-talk mailing list