Re: vertical Interpolation

Dear Dennis,

Thanks for your suggestions.
Now it is interpolating reasonably for the 1st output level (100 hPa), at least.
But for rest of the levels, all values are ZERO!!
Can you please point out the mistake?

"data" & "ps" are inputs from files.

plevs = (/100.,150.,200.,250.,300.,400.,500.,600.,700,750.,800.,850.,900.,925.,950.,975.,1000./)
sigma = (/.99733, .99165, .98521, .97794, .96973, .96049, .95011, .93849,\
                .92551, .91107, .89506, .87739, .85797, .83674, .81366, .78872,\
                .76194, .73340, .70319, .67148, .63846, .60437, .56950, .53415,\
                .49863, .46329, .42843, .39437, .36138, .32970, .29953, .27102,\
                .24429, .21939, .19635, .17517, .15579, .13815, .12218, .10778,\
                .09483, .08324, .07289, .06367, .05549, .04823, .04182, .03615,\
                .03115, .02675, .02288, .01948, .01649, .01387, .01157, .00956,\
                .00780, .00625, .00490, .00372, .00269, .00179, .00099, .00027 \
              /)

hybm = sigma(::-1)
hyam = hybm*0.0

x = vinth2p(data,hyam,hybm,plevs,ps,1,1000,1,True)

Thanks again.

With regards,
Abhik
___________________________________________

Abhik Santra
CSIR Research Fellow,
Indian Institute of Tropical Meteorology,
Pune - 411008.
India
___________________________________________

----- Original Message -----
From: "Dennis Shea" <shea@ucar.edu>
To: "R Phani" <rphani@tropmet.res.in>
Cc: "NCL" <ncl-talk@ucar.edu>, "Abhik Santra" <abhiksantra@tropmet.res.in>
Sent: Saturday, May 10, 2014 7:16:55 PM
Subject: Re: vertical Interpolation

I have no idea why you included
> datai = (/data(:,::-1,:,:)/)
> [SNIP]
> datai&lon = lon
>
=========================
The following should work.

    latS = 0.0 ; whatever you specify
    latN = 50.0

    plevs = (/100.,150.,200.,250.,300.,400.,500.,600.\
              ,700, 750.,800.,850.,900.,925.,950.,975.,1000./)

  ;;plevs!0 = "plevs" ; optional; add meta data; recommended
  ;;plevs@lon_name = "pressure level"
  ;;plevs@units = "hPa"

    sigma = (/.99733, .99165, .98521, .97794, .96973, .96049,\
               .95011, .93849, .92551, .91107, .89506, .87739,\
               .85797, .83674, .81366, .78872, .76194, .73340,\
               .70319, .67148, .63846, .60437, .56950, .53415,\
               .49863, .46329, .42843, .39437, .36138, .32970,\
               .29953, .27102, .24429, .21939, .19635, .17517,\
               .15579, .13815, .12218, .10778, .09483, .08324,\
               .07289, .06367, .05549, .04823, .04182, .03615,\
               .03115, .02675, .02288, .01948, .01649, .01387,\
               .01157, .00956, .00780, .00625, .00490, .00372,\
               .00269, .00179, .00099, .00027 \
             /)
        
    f = addfile ("CLW_T126_2015_orig.grb", "r")
    data = f->MIXR_GDS4_HYBL(:,:,{latS:latN},:)

    printVarSummary(data)
    print("data: min="+min(data)+" max="+max(data))

    hybm = sigma(::-1) ; want top-to-bottom order
    hyam = hybm*0.0 ; same size but filled with 0.0

; file Surf_pres_384_190_T126_2015.grb is not available
; create a bogus variable

  ;;f = addfile ("Surf_pres_384_190_T126_2015.grb", "r")
  ;;ps = f->PRES_GDS4_SFC(:,{latS:latN},:)

    ps = data(:,0,:,:) ; (:,:,:)
    ps = 100000.
    ps@long_name = "Surface Pressure"
    ps@units = "Pa"
    printVarSummary(ps)

    x = vinth2p(data,hyam,hybm,plevs,ps,1,1000,1,True)
    x@long_name = data@long_name
    x@units = data@units

    printVarSummary(x)
    print("x: min="+min(x)+" max="+max(x))

; following *NOT* necessary; rename dimensions

    x!0 = "time"
    x!2 = "lat"
    x!3 = "lon"
    printVarSummary(x)

On 5/9/14, 6:54 AM, R Phani wrote:
> Dear NCL forum,
>
> I am trying to do vertical interpolation for a variable with the the information given below. The levels are in hybrid as integer, so converted the levels to take the sigma levels. There was a communication by Abhik Santra and Dennis has suggested that to reverse the levels. I followed the same but still the maximum and minimum values of the input are 0.0002147481 and 0 respectively but when we use vinth2p function, the max, min values are 0.03422477, -0.0336985 respectively.
>
> The script is..
> -----------------------------
> plevs = (/100.,150.,200.,250.,300.,400.,500.,600.,700,750.,800.,850.,900.,925.,950.,975.,1000./)
> sigma = (/.99733, .99165, .98521, .97794, .96973, .96049, .95011, .93849,\
> .92551, .91107, .89506, .87739, .85797, .83674, .81366, .78872,\
> .76194, .73340, .70319, .67148, .63846, .60437, .56950, .53415,\
> .49863, .46329, .42843, .39437, .36138, .32970, .29953, .27102,\
> .24429, .21939, .19635, .17517, .15579, .13815, .12218, .10778,\
> .09483, .08324, .07289, .06367, .05549, .04823, .04182, .03615,\
> .03115, .02675, .02288, .01948, .01649, .01387, .01157, .00956,\
> .00780, .00625, .00490, .00372, .00269, .00179, .00099, .00027 \
> /)
> f = addfile ("CLW_T126_2015.grb", "r")
> data = f->MIXR_GDS4_HYBL(:,:,{latS:latN},:)
> datai = new(dimsizes(data),typeof(data),data@_FillValue)
> datai = (/data(:,::-1,:,:)/)
> lat = f->g4_lat_2({latS:latN})
> lon = f->g4_lon_3
> hybm = sigma(::-1)
> dimx = dimsizes(data)
> ntim =dimx(0)
> delete(data)
>
> datai!0 = "time"
> datai&time = ispan(0,ntim-1,1)
> datai!1 = "lev"
> datai&lev = hybm
> datai&lev@units = "hPa"
> datai!2 = "lat"
> datai&lat = lat
> datai!3 = "lon"
> datai&lon = lon
> printVarSummary(datai)
>
> hyam = hybm*0.0
>
> f = addfile ("Surf_pres_384_190_T126_2015.grb", "r")
> ps = f->PRES_GDS4_SFC(:,{latS:latN},:)
>
> x = vinth2p(datai,hyam,hybm,plevs,ps,1,1000,1,True)
> --------------------------------------------------------------
> The input "datai" has the dimensions
> Number of Dimensions: 4
> Dimensions and sizes: [time | 365] x [lev | 64] x [lat | 74] x [lon | 384]
> Coordinates:
> time: [0..364]
> lev: [0.00027..0.99733]
> lat: [-24.0944..44.88169]
> lon: [ 0..359.254]
> Number Of Attributes: 1
> _FillValue : 1e+20
>
>
> The file information is
> -------------------------------------------------------
> filename: CLW_T126_2015
> path: CLW_T126_2015.grb
> file global attributes:
> dimensions:
> initial_time0_hours = 365
> lv_HYBL1 = 64
> g4_lat_2 = 190
> g4_lon_3 = 384
> variables:
> float MIXR_GDS4_HYBL ( initial_time0_hours, lv_HYBL1, g4_lat_2, g4_lon_3 )
> center : US National Weather Service - NCEP (WMC)
> long_name : Humidity mixing ratio
> units : kg/kg
> _FillValue : 1e+20
> level_indicator : 109
> gds_grid_type : 4
> parameter_table_version : 2
> parameter_number : 53
> model : Climate Data Assimilation System (CDAS)
> forecast_time : 0
> forecast_time_units : days
>
> double initial_time0_hours ( initial_time0_hours )
> long_name : initial time
> units : hours since 1800-01-01 00:00
>
> double initial_time0_encoded ( initial_time0_hours )
> long_name : initial time encoded as double
> units : yyyymmddhh.hh_frac
>
> float g4_lat_2 ( g4_lat_2 )
> long_name : latitude
> GridType : Gaussian Latitude/Longitude Grid
> units : degrees_north
> N : 95
> Di : 0.938
> Lo2 : 359.254
> La2 : -89.277
> Lo1 : 0
> La1 : 89.277
>
> float g4_lon_3 ( g4_lon_3 )
> long_name : longitude
> GridType : Gaussian Latitude/Longitude Grid
> units : degrees_east
> N : 95
> Di : 0.938
> Lo2 : 359.254
> La2 : -89.277
> Lo1 : 0
> La1 : 89.277
>
> integer lv_HYBL1 ( lv_HYBL1 )
> long_name : Hybrid level
> units : number
>
> string initial_time0 ( initial_time0_hours )
> long_name : Initial time of first record
> units : mm/dd/yyyy (hh:mm)
> ------------------------------------------------------
>
> Phani
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Sat May 10 13:54:56 2014