[ncl-talk] Vertical cross section plot of WRF_Reflectivity (10 Cm) (REFL_10CM) through various longitude corresponding to a single lat
Bill Ladwig
ladwig at ucar.edu
Mon Sep 24 15:18:12 MDT 2018
Hi Trisanu,
See the attached script. This uses only the WRF routines for the
computations. Since I don't have your data, here is a cross section that it
created for a hurricane. It uses a start point and an end point for the
cross section line, rather than a pivot point and angle. Hopefully there
aren't too many typos.
Regards,
Bill
On Mon, Sep 24, 2018 at 7:50 AM Trisanu Banik <baniktrisanu at gmail.com>
wrote:
> Sir
>
> Thanks for your feedback.
>
> The problem is still persisting in a different manner. I have attached
> three plots for your convenience. The plot REFL_10CM.000002.jpg is vertical
> cross section of Reflectivity and model level (36) is in Y axis. I just
> wish to plot same reflectivity keeping height in Y axis instead of model
> level. I have utilized two function 1st "int2p_n_Wrap" and plotted cross
> sectional plot attached in the email as a name of cross.000001.jpg using
> plot function "gsn_csm_contour". Secondly I have utilized
> "wrf_user_intrp3d" and plotted cross.jpg utilizing the same plot function.
>
> But interpolated plot does not resemble with the REFL_10CM. cross.jpg
> represents reflectivity well but it provide the grid info in X and Y axis
> which is undesirable for me. As the function "wrf_user_intrp3d" generates
> two dimensional variable with vertical and horizontal grid info in two
> axis. Please guide how generate plot like REFL_10CM.000002.jpg with height
> in Y axis.
>
> Regards
> Trisanu
>
>
> *Trisanu Banik,PhD*
>
> *Research Scientist*
>
> *North Eastern Space Applications Centre (NESAC)*
> *Government of India*
> *Department of Space*
>
> *Umiam-793103, Meghalaya*
> *Mobile-9774837581*
>
>
>
> On Thu, Sep 20, 2018 at 2:02 AM Bill Ladwig <ladwig at ucar.edu> wrote:
>
>> Is there a reason why you can't use the cross section capability of
>> wrf_user_intrp3d? See Example 1 here:
>>
>> https://www.ncl.ucar.edu/Document/Functions/WRF_arw/wrf_user_intrp3d.shtml
>>
>> If you need specific height levels, rather than the 1% increments the
>> routine does by default, I can show you how to accomplish that (see code
>> below, or do a search, since I've guided others through this in the past).
>>
>> Also, WRF uses PH + PHB to get the geopotential, not ZNW, so to get the
>> height at each grid point, you do (PH + PHB) / 9.81 (or just use
>> wrf_user_getvar with 'z' as the product type). (If you don't have both of
>> those variables in your WRF file, you're going to be in trouble). Also, the
>> right two dimensions in a WRF file aren't lat and lon, they're south_north,
>> and west_east, which is a bit of a misnomer since it's really used for
>> "bottom to top" and "left to right" for a grid. WRF uses curvilinear
>> coordinates, so the latitude and longitude values require two dimensions to
>> determine them (this is why XLAT and XLONG are two dimensions).
>>
>> For reference, the code that computes the cross section in
>> wrf_user_intrp3d is here (removed the horizontal interpolation stuff):
>>
>> undef("wrf_user_intrp3d")
>>> function wrf_user_intrp3d( var3d:numeric, z_in:numeric, \
>>> plot_type:string, \
>>> loc_param:numeric, angle:numeric,
>>> opts:logical )
>>> ; var3d - 3d field to interpolate (all input fields must be
>>> unstaggered)
>>> ; z_in - interpolate to this field (either p/z)
>>> ; plot_type - interpolate horizontally "h", or vertically "v"
>>> ; loc_param - level(s) for horizontal plots (eg. 500hPa ; 3000m -
>>> scalar),
>>> ; plane for vertical plots (2 values representing an xy
>>> point
>>> ; on the model domain through which the vertical plane
>>> will pass
>>> ; OR 4 values specifying start and end values
>>> ; angle - 0.0 for horizontal plots, and
>>> ; an angle for vertical plots - 90 represent a WE cross
>>> section
>>> ; opts Used IF opts is TRUE, else use loc_param and angle to
>>> determine crosssection
>>> begin
>>>
>>> if(plot_type .eq. "v" ) then ; vertical cross section needed
>>> dims = dimsizes(var3d)
>>> nd = dimsizes(dims)
>>> dimX = dims(nd-1)
>>> dimY = dims(nd-2)
>>> dimZ = dims(nd-3)
>>> if ( nd .eq. 4 ) then
>>> if ( z_in(0,dimZ-1,0,0) .lt. z_in(0,dimZ-2,0,0) ) then
>>> ; We must be interpolating to pressure
>>> ; This routine needs input field and level in hPa - lets
>>> make sure of this
>>> if ( max(z_in) .gt. 2000. ) then
>>> ; looks like we have Pa as input - make this hPa
>>> z_in = z_in * 0.01
>>> end if
>>> end if
>>> z = z_in(0,:,:,:)
>>> else
>>> if ( z_in(dimZ-1,0,0) .lt. z_in(dimZ-2,0,0) ) then
>>> ; We must be interpolating to pressure
>>> ; This routine needs input field and level in hPa - lets
>>> make sure of this
>>> if ( z_in(0,0,0) .gt. 2000. ) then
>>> ; looks like we have Pa as input - make this hPa
>>> z_in = z_in * 0.01
>>> end if
>>> end if
>>> z = z_in
>>> end if
>>>
>>> ; set vertical cross section
>>> if (opts) then
>>> xy = wrf_user_set_xy( z, loc_param(0), loc_param(1), \ ;
>>> assumes 0-based indexing in v6.5.0
>>> loc_param(2), loc_param(3), \
>>> angle, opts )
>>> else
>>> xy = wrf_user_set_xy( z, loc_param(0), loc_param(1), \ ;
>>> 0.0, 0.0, angle, opts )
>>> end if
>>> xp = dimsizes(xy)
>>>
>>>
>>> ; first we interp z
>>> var2dz = wrf_interp_2d_xy( z, xy)
>>> ; interp to constant z grid
>>> if(var2dz(0,0) .gt. var2dz(1,0) ) then ; monotonically decreasing
>>> coordinate
>>> z_max = floor(max(z)/10)*10 ; bottom value
>>> z_min = ceil(min(z)/10)*10 ; top value
>>> dz = 10
>>> nlevels = tointeger( (z_max-z_min)/dz)
>>> z_var2d = new( (/nlevels/), typeof(z))
>>> z_var2d(0) = z_max
>>> dz = -dz
>>> else
>>> z_max = max(z)
>>> z_min = 0.
>>> dz = 0.01 * z_max
>>> nlevels = tointeger( z_max/dz )
>>> z_var2d = new( (/nlevels/), typeof(z))
>>> z_var2d(0) = z_min
>>> end if
>>>
>>> do i=1, nlevels-1
>>> z_var2d(i) = z_var2d(0)+i*dz
>>> end do
>>>
>>> ; interp the variable
>>> if ( dimsizes(dims) .eq. 4 ) then
>>> var2d = new( (/dims(0), nlevels, xp(0)/), typeof(var2dz))
>>> do it = 0,dims(0)-1
>>> var2dtmp = wrf_interp_2d_xy( var3d(it,:,:,:), xy)
>>> do i=0,xp(0)-1
>>> var2d(it,:,i) = wrf_interp_1d( var2dtmp(:,i), var2dz(:,i),
>>> z_var2d)
>>> end do
>>> end do
>>> var2d!0 = var3d!0
>>> var2d!1 = "Vertical"
>>> var2d!2 = "Horizontal"
>>> else
>>> var2d = new( (/nlevels, xp(0)/), typeof(var2dz))
>>> var2dtmp = wrf_interp_2d_xy( var3d, xy)
>>> do i=0,xp(0)-1
>>> var2d(:,i) = wrf_interp_1d( var2dtmp(:,i), var2dz(:,i),
>>> z_var2d)
>>> end do
>>> var2d!0 = "Vertical"
>>> var2d!1 = "Horizontal"
>>> end if
>>>
>>>
>>> st_x = tointeger(xy(0,0)) ; + 1 (removed 1-based indexing in 6.5.0)
>>>
>>> st_y = tointeger(xy(0,1)) ; + 1
>>> ed_x = tointeger(xy(xp(0)-1,0)) ; + 1
>>> ed_y = tointeger(xy(xp(0)-1,1)) ; + 1
>>> if (opts) then
>>> var2d at Orientation = "Cross-Section: (" + \
>>> st_x + "," + st_y + ") to (" + \
>>> ed_x + "," + ed_y + ")"
>>> else
>>> var2d at Orientation = "Cross-Section: (" + \
>>> st_x + "," + st_y + ") to (" + \
>>> ed_x + "," + ed_y + ") ; center=(" + \
>>> loc_param(0) + "," + loc_param(1) + \
>>> ") ; angle=" + angle
>>> end if
>>> return(var2d)
>>> end if
>>>
>>>
>>> end
>>
>>
>> If you want to set your own levels, set z_var2d to the levels you want
>> and comment out the section of code that computes the 1% increments. You'll
>> want to rename this function and copy this code in to your script to use it.
>>
>> Hope this helps,
>>
>> Bill
>>
>> On Wed, Sep 19, 2018 at 7:44 AM Trisanu Banik <baniktrisanu at gmail.com>
>> wrote:
>>
>>> Hi
>>>
>>> I am still struggling with the plot of vertical cross section of
>>> reflectivity with respect to height and various lat keeping long constant
>>> and vice versa.
>>> During the search of cross sectional plots in web, I got few functions
>>> that might be helpful for generating vertical cross sectional plot with
>>> respect to height.
>>>
>>> But my data has four dimensions time, eta level, lat, long. If I
>>> consider the values for particular time it will become 3D i.e. level, lat,
>>> long. Now to convert eta level to pressure I have utilized "q_isobaric =
>>> int2p_n_Wrap(P_HYD,REFL_10CM,level,linlog,0) " this. Presently I am able to
>>> plot reflectivity w..r.t pressure for various lat long.
>>>
>>> Now the challenge is to convert vertical level to height. Various
>>> inbuilt functions are there in NCL that convert pressure to height like "wrf_user_vert_interp"
>>> and gsn_csm_pres_hgt etc. But my data dimension does not satisfy the
>>> functions requirement. In wrf_out file, the necessary variables related to
>>> this are, one dimentional vertical eta level (ZNW), Hydrostatic Pressure
>>> level (P_HYD) and reflectivity (REFL_10CM). No height related 3D parameter
>>> is available in the wrf_out file. Only terrain height is available which is
>>> 2D.
>>>
>>> For your convenience, I have attached my tentative code. Please guide me
>>> to plot the vertical cross section of reflectivity with Height in Y axis
>>> for various lat long.
>>>
>>> Thanks & Regards
>>> Trisanu
>>>
>>>
>>>
>>> *Trisanu Banik,PhD*
>>>
>>> *Research Scientist*
>>>
>>> *North Eastern Space Applications Centre (NESAC)*
>>> *Government of India*
>>> *Department of Space*
>>>
>>> *Umiam-793103, Meghalaya*
>>> *Mobile-9774837581*
>>>
>>>
>>>
>>> On Thu, Sep 13, 2018 at 7:58 PM Jim Means <jim at weatherextreme.com>
>>> wrote:
>>>
>>>> In addition to Soma's suggestion, you might also look at the cross
>>>> section plot examples here:
>>>>
>>>>
>>>> http://www2.mmm.ucar.edu/wrf/OnLineTutorial/Graphics/NCL/NCL_examples.htm
>>>>
>>>> On 9/13/2018 7:09 AM, Soma Roy wrote:
>>>>
>>>> See the different examples at the below link:
>>>> https://www.ncl.ucar.edu/Applications/height_long.shtml
>>>>
>>>> Best,
>>>> Soma
>>>>
>>>> On Thu, Sep 13, 2018, 19:10 Trisanu Banik <baniktrisanu at gmail.com>
>>>> wrote:
>>>>
>>>>> Hi,
>>>>>
>>>>> I need help to plot a vertical cross section plot of WRF REFL_10CM
>>>>> (Reflectivity). The parameter is four dimensional (Time, eta level, lat,
>>>>> lon). I tried to plot the vertical cross section of reflectivity with
>>>>> respect to height for a particular time. But the REFL_10CM parameter
>>>>> contain only eta level. Therefore whenever I generated plot, it always
>>>>> provide vertical cross section w.r.t eta level only. i need suggestion how
>>>>> I can convert eta level to pressure level/Height level. So that, I can plot
>>>>> a vertical cross section of WRF reflectivity for various longitude
>>>>> corresponding to a single latitude. Here I have added a sample plot. I need
>>>>> to plot vertical cross section w.r.t height (not w.r.t eta level)
>>>>>
>>>>> Please help me to resolve the issue.
>>>>>
>>>>>
>>>>> *Trisanu Banik,PhD*
>>>>>
>>>>> *Research Scientist *
>>>>>
>>>>> *North Eastern Space Applications Centre (NESAC) *
>>>>> *Government of India*
>>>>> *Department of Space*
>>>>>
>>>>> *Umiam-793103, Meghalaya *
>>>>> *Mobile-9774837581*
>>>>>
>>>>> _______________________________________________
>>>>> ncl-talk mailing list
>>>>> ncl-talk at ucar.edu
>>>>> List instructions, subscriber options, unsubscribe:
>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>
>>>>
>>>>
>>>> _______________________________________________
>>>> ncl-talk mailing listncl-talk at ucar.edu
>>>> List instructions, subscriber options, unsubscribe:http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>
>>>>
>>>> --
>>>>
>>>> James D. Means, Ph.D.
>>>> Senior Atmospheric & Climate Scientist
>>>> California Office
>>>> Tele: 619-495-1638 | Fax: 775-636-8430
>>>> 930 Tahoe Blvd., Suite 802-560
>>>> Incline Village, Nevada 89451
>>>> jim at weatherextreme.com | vcard
>>>> <http://www.weatherextreme.com/vcards/James%20Means.vcf>
>>>> www.weatherextreme.com
>>>> <http://www.weatherextreme.com/>
>>>> _______________________________________________
>>>> ncl-talk mailing list
>>>> ncl-talk at ucar.edu
>>>> List instructions, subscriber options, unsubscribe:
>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>
>>> _______________________________________________
>>> ncl-talk mailing list
>>> ncl-talk at ucar.edu
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20180924/c12d5b3f/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: refl10_cross.ncl
Type: application/octet-stream
Size: 2525 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20180924/c12d5b3f/attachment.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: cross.png
Type: image/png
Size: 72749 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20180924/c12d5b3f/attachment.png>
More information about the ncl-talk
mailing list