WRF-H domain is different from WRF domain

[Erin Dougherty]

Hi there, 

I am plotting my WRF-H output and it seems shifted too far south compared to the geo_em_d01.nc file I provided (attached here). The original gridded streamflow output was plotted using the x and y dimensions output from the file. I also tried converting the x,y dimensions to lat, lon using pyproj: 

source_crs = 'esri:102009' # Coordinate system of the file
target_crs = 'epsg:4326' # Global lat-lon coordinate system

pyproj.Transformer.from_crs(source_crs, target_crs)

However, that resulted in a domain that is too far east compared to the original WRF domain. Do you have any idea as to what is going on here?

Thank you, 

Erin

[Kevin]

Erin,

It is unlikely that the 'esri:102009' standard coordinate system (https://epsg.io/102009 ) actually matches the coordinate system of your geogrid file. Typically, geogrid files are created in WPS and the coordinate system will be centered on the outermost domain center. Also, the datum will be different because WRF is based on a sphere and WGS84 is based on a spheroidal datum. You can use the Esri WKT or the Proj4 string provided in the Fulldom_hires netCDF file to define the correct coordinate system for your domain. No need to project the coordinates to lat/lon if your mapping application supports projected coordinates.

Kevin

[Erin Dougherty]

Thank you for the insight, Kevin. 

We figured out the issue, and it came down to the projection of the WRF-H data used in conjection with the cartopy projections in python. To resolve this issue, the WRF-H native x,y grid data was projected to lat/lon using the pyproj package: 

from pyproj import Proj
X, Y = np.meshgrid(ctrl_wrfh_chrtout_gr.x.values, ctrl_wrfh_chrtout_gr.y.values)

source = pyproj.Proj(ctrl_wrfh_chrtout_gr.proj4)

lon, lat = source(X.flatten(), Y.flatten(), inverse=True)

The lat/lons were then added as coordinates to the gridded streamflow xarray dataset:

mean_stmflow_grd.coords['lat'] = (mean_stmflow_grd.dims, lat_2d)
mean_stmflow_grd.coords['lon'] = (mean_stmflow_grd.dims, lon_2d)

Then the streamflow could be plotted over shapefiles of U.S. states using the transform function in cartopy:

fig = plt.figure(figsize = (8,8))
v = np.linspace(0,1200,11,endpoint = True)
cmap = 'BuPu'
ax = plt.subplot(1, 1, 1, projection=ccrs.LambertConformal(central_longitude=-93))
ax.add_geometries(state, crs=ccrs.PlateCarree(), facecolor='none', edgecolor='0.1', linewidth = 0.6, zorder=6)
c1 = plt.pcolormesh(mean_stmflow_grd.lon.values, mean_stmflow_grd.lat.values, mean_stmflow_grd.values, vmax=100, cmap = cmap, transform=ccrs.PlateCarree(),)

plt.show()

I hope this helps if others have the same issue! 

- Erin