Reputation: 21
I am working with some atmospheric model outputs available in NetCDF format. I want to calculate Temperature Advection by 850 hPa winds. I tried using the metpy function https://unidata.github.io/MetPy/latest/api/generated/metpy.calc.advection.html#metpy.calc.advection for the same. I am attaching a part of my code. Please help to resolve this error. All required libraries are updated with their latest versions.
import xarray as xr
import numpy as np
import matplotlib.pyplot as plt
from metpy.units import units
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import metpy.calc as mpcalc
temp = xr.open_dataset('ts.nc')
u10 = xr.open_dataset('u10.nc')
v10 = xr.open_dataset('v10.nc')
lat = temp['lat']
lon = temp['lon']
lonn, latt = np.meshgrid(lon,lat)
dx, dy = mpcalc.lat_lon_grid_deltas(lon, lat, )
time = np.array(temp.time[0:9])
T1 = temp['ts'][0:9,:,:]
u10_1 = u10['u10'][0:9,:,:]
v10_1 = v10['v10'][0:9,:,:]
T11 = xr.DataArray(data = T1.values, dims = ['time', 'lat', 'lon'], coords=dict(lat=lat, lon = lon, time = time))
u10_11 = xr.DataArray(data = u10_1.values, dims = ['time', 'lat', 'lon'], coords=dict(lat=lat, lon = lon, time = time))
v10_11 = xr.DataArray(data = v10_1.values, dims = ['time', 'lat', 'lon'], coords=dict(lat=lat, lon = lon, time = time))
advec = mpcalc.advection(T1, [u10_1, v10_1], (dx, dy))
[code][1]
Error: AttributeError: crs attribute is not available.
Please let me know if any further references are required.
Upvotes: 1
Views: 289
Reputation: 5843
If there are CF projection metadata present in your dataset, MetPy may be able to automatically identify the proper CRS, by calling the parse_cf()
method:
v10 = xr.open_dataset('v10.nc')
v10 = v10.metpy.parse_cf()
Upvotes: 0
Reputation: 11
Mentioned error (AttributeError: crs attribute is not available) speaks itself the problem in the code.
This error is because of the absence of cartopy projection attribute in the data coordinates. You need to manually add the projection in the data.
T1=T1.metpy.assign_crs(grid_mapping_name='latitude_longitude', earth_radius=6371229.0)
You can pick any one of the projections from this list.
Upvotes: 1