[example] cmor xwrf example

[larsbuntemeyer]

This is my example for running cmorization with WRF datasets (ignore the REMO specific global attributes please):

from cordex.tables import cordex_cmor_table
from cordex import cmor
import pandas as pd
import xarray as xr
import xwrf

table_prefix = "CORDEX-CMIP6"

# get example dataset from xwrf 
ds_old = xwrf.tutorial.open_dataset("wrfout")
ds_new = ds_old.xwrf.postprocess()
# we use the lowest layer as tas dummy
ds = (
    ds_new[
        ["x", "y", "air_potential_temperature", "XLAT", "XLONG", "wrf_projection"]
    ]
    .isel(z=0)
    .rename(Time="time", air_potential_temperature="tas")
)
# extend time axis and make dummy monthly dataset
time = pd.date_range("2025-01-01", periods=12, freq="MS")
ds = (
    xr.decode_cf(ds, decode_coords="all")
    .squeeze(drop=True)
    .expand_dims(time=time, axis=0)
)
table_id = "mon"
filename = cmor.cmorize_variable(
    ds,
    "tas",
    cmor_table=cordex_cmor_table(f"{table_prefix}_{table_id}"),
    dataset_table=cordex_cmor_table(f"{table_prefix}_remo_example"),
    grids_table=cordex_cmor_table(f"{table_prefix}_grids"),
    domain_id="West-Coast",
    crop=False,
)
ds_out = xr.open_dataset(filename)

def plot(ds):
    from cartopy import crs as ccrs
    import cartopy.feature as cf
    import matplotlib.pyplot as plt
    grid_mapping = ds.cf['grid_mapping']
    transform = ccrs.LambertConformal(central_longitude=grid_mapping.longitude_of_central_meridian,
                                      central_latitude=grid_mapping.latitude_of_projection_origin,
                                      false_easting=grid_mapping.false_easting,
                                      false_northing=grid_mapping.false_northing,
                                      standard_parallels=grid_mapping.standard_parallel,)
                                    
    plt.figure(figsize=(10,8))
    ax = plt.axes(projection=transform)
    ax.gridlines(draw_labels=True, linewidth=0.5, color='gray',
                xlocs=range(-180,180,10), ylocs=range(-90,90,5))
    ds.tas.plot(ax=ax, transform=transform)
    ax.coastlines(resolution='50m', color='black', linewidth=1)
    ax.add_feature(cf.BORDERS)

plot(ds_out.isel(time=0))

[larsbuntemeyer]

pinging @alessandroanav