Scheduled Downtime
On Friday 21 April 2023 @ 5pm MT, this website will be down for maintenance and expected to return online the morning of 24 April 2023 at the latest

difference of z and z_agl in wrf-python

P

prick

New member
I tried to use the built-in scripts in wrf-python to get the elevation at Denver. My understanding is z represents height relative to MSL, z_agl is relative to ground level. However, The minimum z that I got at Denver is ~26m, and the minumum z_agl I got at Denver is ~63m. The difference is only <40m, while the ground level at Denver is supposed to be around 5000ft (~1700m) above MSL.
Is there any mistake with my interpretation with the difference of z and z_agl?
Thank you.
 
I think a moderator transferred your question to a separate thread! I was wondering where the reply went.

The "z" coordinate is above MSL. Not sure which function you used, did you use vinterp()? That one has an option to interpolate to ght_agl (above ground level), maybe you used that?

I cannot answer your question without specifics.
 
Here is the complete code. Thank you!


from __future__ import print_function, division

from netCDF4 import Dataset
from wrf import getvar, interpline, CoordPair, xy_to_ll, ll_to_xy, latlon_coords
import scipy
from scipy.interpolate import interp1d
import matplotlib.pyplot as plt

ncfile = Dataset("./wrfout/wrfout_d01_2019-11-26_21_00_00")

z = getvar(ncfile, "z", units="m", timeidx=0)
z_agl = getvar(ncfile, "height_agl", units="m", timeidx=0)
wspd_wdir_earth = getvar(ncfile, 'uvmet_wspd_wdir', timeidx=0)
geopt=getvar(ncfile,'geopt', timeidx=0)


# denver, 39.87 N, 104.67 W
xy_wrfpy = ll_to_xy(ncfile, 39.87, -104.67)
wspd_1d=wspd_wdir_earth[0,:,xy_wrfpy[0], xy_wrfpy[1]]
zagl_1d=z_agl[:,xy_wrfpy[0],xy_wrfpy[1]]
z_1d=z_agl[:,xy_wrfpy[0],xy_wrfpy[1]]

print(zagl_1d)
print(z_1d)
 
Here is what I got after running:


<xarray.DataArray 'height' (bottom_top: 49)> Size: 196B
array([ 63.682693, 123.44999 , 199.56142 , 295.9557 ,
417.1872 , 568.3548 , 755.0117 , 982.8348 ,
1256.7452 , 1580.3383 , 1955.5077 , 2381.7334 ,
2856.307 , 3375.887 , 3939.921 , 4523.442 ,
5100.6646 , 5671.849 , 6236.5776 , 6793.988 ,
7343.4565 , 7884.6074 , 8417.316 , 8941.73 ,
9457.925 , 9965.789 , 10465.27 , 10956.482 ,
11439.354 , 11913.822 , 12380.541 , 12840.4795 ,
13294.357 , 13742.679 , 14186.403 , 14626.831 ,
15064.389 , 15498.923 , 15930.86 , 16361.131 ,
16790.205 , 17218.424 , 17647.033 , 18077.584 ,
18511.283 , 18948.877 , 19390.133 , 19834.543 ,
20282.416 ], dtype=float32)
Coordinates:
XLONG float32 4B -93.02
XLAT float32 4B 31.57
XTIME float32 4B 180.0
Time datetime64[ns] 8B 2019-11-26T21:00:00
latlon_coord object 8B CoordPair(lat=39.87, lon=-104.67)
Dimensions without coordinates: bottom_top
Attributes:
FieldType: 104
MemoryOrder: XYZ
description: model height - [MSL] (mass grid)
units: m
stagger:
coordinates: XLONG XLAT XTIME
projection: LambertConformal(stand_lon=-92.5, moad_cen_lat=37.489994049...



<xarray.DataArray 'height_agl' (bottom_top: 49)> Size: 196B
array([ 26.25071 , 86.018005, 162.12943 , 258.5237 ,
379.75522 , 530.9228 , 717.5797 , 945.4028 ,
1219.3132 , 1542.9062 , 1918.0757 , 2344.3015 ,
2818.875 , 3338.455 , 3902.489 , 4486.01 ,
5063.2324 , 5634.417 , 6199.1455 , 6756.5557 ,
7306.0244 , 7847.1753 , 8379.885 , 8904.299 ,
9420.493 , 9928.357 , 10427.838 , 10919.051 ,
11401.922 , 11876.391 , 12343.109 , 12803.048 ,
13256.926 , 13705.247 , 14148.972 , 14589.399 ,
15026.957 , 15461.491 , 15893.429 , 16323.699 ,
16752.773 , 17180.992 , 17609.602 , 18040.152 ,
18473.852 , 18911.445 , 19352.701 , 19797.111 ,
20244.984 ], dtype=float32)
Coordinates:
XLONG float32 4B -93.02
XLAT float32 4B 31.57
XTIME float32 4B 180.0
Time datetime64[ns] 8B 2019-11-26T21:00:00
latlon_coord object 8B CoordPair(lat=39.87, lon=-104.67)
Dimensions without coordinates: bottom_top
Attributes:
FieldType: 104
MemoryOrder: XYZ
description: model height - [AGL] (mass grid)
units: m
stagger:
coordinates: XLONG XLAT XTIME
projection: LambertConformal(stand_lon=-92.5, moad_cen_lat=37.489994049...
 
Huh, no idea why it does that. Are there some staggering shenanigans going on? Geopotential is vertically staggered.

You'll have to ask someone who is an actual programmer for this, but that might be a problem since (from what I've heard) wrf-python is not maintained anymore.

EDIT: The first height grid should be above ground at all times, even for "z" since the vertical coordinate is quasi-terrain following at the bottom, so it does not go below the ground. No idea, sorry I can't help anymore.
 
Do you think there is any problem the way how I use ll_to_xy?

# denver, 39.87 N, 104.67 W
xy_wrfpy = ll_to_xy(ncfile, 39.87, -104.67)

I noticed that I have to put "-104.67" instead of "104.67" to make the result look reasonable. Is it because it's W?
 
That makes sense. Thanks.
I just want to make sure that I used the right functions to plot the vertical wind speed profile.
I used the following functions to plot wind speed vs. elevation above MSL: (instead of z_agl, I used z)
z = getvar(ncfile, "z", units="m", timeidx=0)
wspd_wdir_earth = getvar(ncfile, 'uvmet_wspd_wdir', timeidx=0)

Do you see any problems with these?
 
Well, the fields seem fine and I don't see any problem with the procedure. However, one problem that you may wish to fix is that with the way you are doing this now, your heights will not be constant with time (since these are the heights of the eta surfaces).

If you wanted to look at profiles at constant height levels with time, you'd need to interpolate them to some fixes levels.
 
You must log in or register to reply here.
Top