1

I want to make a cross section of frontogenesis my code is

  domain_1 = os.path.abspath(filenames_in)
info = os.path.join(domain_1, filename)
data = xr.open_dataset(info)
data = data.metpy.parse_cf().squeeze()
print(data)
p1 = data['level'].values[:]* units.hPa
start = (56, -9)
end = (56.5, -6)
cross = cross_section(data, start, end).set_coords(('latitude', 'longitude'))
r,le,t,w =xr.broadcast(cross['r'],cross['level'],cross['t'],cross['w'])
# print('r',r.shape)
# print(le)
# print('w',w.shape)
# print(w)
theta=mpcalc.potential_temperature(level,t)
cross['Potential_temperature'] = xr.DataArray(theta, coords=t.coords, dims=t.dims)
ptemp=cross['Potential_temperature']
cross['u'] = cross['u']
u = cross['u']
cross['v'] = cross['v']
v = cross['v']
# Set subset slice for the geographic extent of data to limit download
lon_slice = slice(-40, 20)
lat_slice = slice(70, 40)
# Grab lat/lon values (GFS will be 1D)
lats = data.latitude.sel(latitude=lat_slice).values
lons = data.longitude.sel(longitude=lon_slice).values
# Compute dx and dy spacing for use in vorticity calculation
dx, dy = mpcalc.lat_lon_grid_deltas(lons,lats)
front = mpcalc.frontogenesis(ptemp, u, v, dx[None,:, :], dy[None,:, :], x_dim=-1, y_dim=-2)
tmpk = data.t.metpy.sel(
    latitude=lat_slice, longitude=lon_slice).metpy.unit_array.squeeze()
uwnd = data['u'].metpy.sel(
    latitude=lat_slice, longitude=lon_slice).metpy.unit_array.squeeze()
vwnd = data['v'].metpy.sel(
    latitude=lat_slice, longitude=lon_slice).metpy.unit_array.squeeze()
potent_temp=mpcalc.potential_temperature(p1[:,None,None], tmpk)
fronto = mpcalc.frontogenesis(potent_temp, uwnd, vwnd, dx[None, :, :], dy[None, :, :], x_dim=-1, y_dim=-2)
# f=cross_section(fronto, start, end).set_coords(('latitude', 'longitude'))

there are always some bugs. how to fix it I do not know if i need to calculate the frontogenesis first and then use cross-function or use cross function first and then calculate frontogenesis. I tries both of them but I cannot get the cross section.

ValueError: operands could not be broadcast together with shapes (0,120,241) (19,100)

if i calculate the frontogenesis first and then use cross-function

lon_slice = slice(-40, 20) 
lat_slice = slice(70, 40) 
subset = data.metpy.sel(latitude=lat_slice, longitude=lon_slice)
subset['potential_temperature'] = mpcalc.potential_temperature(
subset['level'],
subset['t'])
subset['frontogenesis'] = mpcalc.frontogensis(
subset['potential_temperature'],
subset['u'],
subset['v'])
start = (56, -9) end = (56.5, -6)
cross = cross_section(subset, start, end).set_coords(('latitude', 'longitude'))

i always meet this error

ValueError: Data missing required coordinate information. Verify that your data have been parsed by MetPy with proper x and y dimension coordinates and added crs coordinate of the correct projection for each variable.
tianhang zhang
  • 13
  • 3

1 Answers1

1

Given that MetPy's frontogensis calculation operates on a 2D horizontal grid, you will indeed need to calculate frontogensis prior to taking the cross section. Similarly, subsetting should be done prior to calculations, as MetPy does not yet support delayed computations. Additionally, instead of manually handling grid deltas and unit arrays, you will typically run into less errors when you let MetPy handle those for you, like so:

domain_1 = os.path.abspath(filenames_in)
info = os.path.join(domain_1, filename)
data = xr.open_dataset(info)
data = data.metpy.parse_cf().squeeze()
print(data)

lon_slice = slice(-40, 20)
lat_slice = slice(70, 40)
subset = data.metpy.sel(latitude=lat_slice, longitude=lon_slice)

subset['potential_temperature'] = mpcalc.potential_temperature(
    subset['level'],
    subset['t']
)
subset['frontogenesis'] = mpcalc.frontogensis(
    subset['potential_temperature'],
    subset['u'],
    subset['v']
)

start = (56, -9)
end = (56.5, -6)
cross = cross_section(subset, start, end).set_coords(('latitude', 'longitude'))

(This assumes that your data have proper 'units' attributes and aligned grids in horizontal and vertical. If this is not the case, some extra data cleaning will be needed after opening.)

Jon Thielen
  • 479
  • 2
  • 7
  • Hi Jon, i always meet this error `code`AttributeError: crs attribute is not available. You may need to use the `parse_cf` or `assign_crs` methods. Consult the "xarray with MetPy Tutorial" for more details. During handling of the above exception, another exception occurred: ValueError Traceback (most recent call last) ~\AppData\Local\Temp/ipykernel_5372/1239301948.py in 24 start = (56, -9) 25 end = (56.5, -6) ---> 26 cross = cross_section(subset, start, end).set_coords(('latitude', 'longitude')) 27 – tianhang zhang May 24 '22 at 20:56
  • ValueError: Data missing required coordinate information. Verify that your data have been parsed by MetPy with proper x and y dimension coordinates and added crs coordinate of the correct projection for each variable. @Jon Thielen – tianhang zhang May 24 '22 at 21:00
  • the data what i am using is ERA5 data. thank yo for your answer! – tianhang zhang May 24 '22 at 21:10
  • That message indicates that the data you have do not have the proper metadata for MetPy to understand the data's coordinate system, which is needed to properly calculate the cross-section coordinates. The message indicates to consult MetPy's docs about how to add this information: https://unidata.github.io/MetPy/latest/tutorials/xarray_tutorial.html – DopplerShift May 25 '22 at 00:07