diff --git a/docs/Demo.ipynb b/docs/Demo-AK.ipynb similarity index 99% rename from docs/Demo.ipynb rename to docs/Demo-AK.ipynb index 4b4dc51..c1f37d0 100644 --- a/docs/Demo.ipynb +++ b/docs/Demo-AK.ipynb @@ -5,7 +5,7 @@ "id": "a9da2fa8-9a25-497e-8ec6-d08609ad38ca", "metadata": {}, "source": [ - "# Demonstration of `ocean-model-skill-assessor`\n", + "# Demonstration of `ocean-model-skill-assessor`: Alaska example\n", "\n", "Here we demonstrate how to use `ocean-model-skill-assessor` as a Python package." ] diff --git a/notebooks/Demo.ipynb b/docs/Demo-CA.ipynb similarity index 99% rename from notebooks/Demo.ipynb rename to docs/Demo-CA.ipynb index 4b4dc51..343094e 100644 --- a/notebooks/Demo.ipynb +++ b/docs/Demo-CA.ipynb @@ -5,7 +5,7 @@ "id": "a9da2fa8-9a25-497e-8ec6-d08609ad38ca", "metadata": {}, "source": [ - "# Demonstration of `ocean-model-skill-assessor`\n", + "# Demonstration of `ocean-model-skill-assessor`: California example\n", "\n", "Here we demonstrate how to use `ocean-model-skill-assessor` as a Python package." ] diff --git a/docs/index.rst b/docs/index.rst index bd03781..e0d1dc5 100644 --- a/docs/index.rst +++ b/docs/index.rst @@ -9,7 +9,8 @@ Welcome to ocean-model-skill-assessor's documentation! .. toctree:: :maxdepth: 3 - Demo + Demo-AK + Demo-CA Demo_workflows api diff --git a/ocean_model_skill_assessor/main.py b/ocean_model_skill_assessor/main.py index 76f60d5..4f0921c 100644 --- a/ocean_model_skill_assessor/main.py +++ b/ocean_model_skill_assessor/main.py @@ -52,11 +52,13 @@ def find_bbox(ds): lon = ds.cf["longitude"].values lat = ds.cf["latitude"].values - except KeyError as e: + except KeyError: # In case there are multiple grids, just take first one; # they are close enough - lon = list(ds.cf[["longitude"]].coords.keys())[0].values - lat = list(ds.cf[["latitude"]].coords.keys())[0].values + lon = ds[ds.cf.coordinates["longitude"][0]] + lat = ds[ds.cf.coordinates["latitude"][0]] + # lon = list(ds.cf[["longitude"]].coords.keys())[0].values + # lat = list(ds.cf[["latitude"]].coords.keys())[0].values min_lon = lon.min() max_lon = lon.max() @@ -97,7 +99,10 @@ def read_model(loc_model, xarray_kwargs, time_range=None): xarray Dataset containing model output. """ - dsm = xr.open_dataset(loc_model, **xarray_kwargs) + if isinstance(loc_model, list): + dsm = xr.open_mfdataset(loc_model, **xarray_kwargs) + else: + dsm = xr.open_dataset(loc_model, **xarray_kwargs) # add more cf-xarray info dsm = dsm.cf.guess_coord_axis() @@ -112,10 +117,12 @@ def read_model(loc_model, xarray_kwargs, time_range=None): dsm = dsm.cf.isel(T=index).cf.sel(T=slice(time_range[0], time_range[1])) # force longitude to be from -180 to 180 - lkey = dsm.cf["longitude"].name - dsm[lkey] = dsm.cf["longitude"].where( - dsm.cf["longitude"] < 180, dsm.cf["longitude"] - 360 - ) + for lkey in dsm.cf.coordinates["longitude"]: + dsm[lkey] = dsm[lkey].where(dsm[lkey] < 180, dsm[lkey] - 360) + # lkey = dsm.cf["longitude"].name + # dsm[lkey] = dsm.cf["longitude"].where( + # dsm.cf["longitude"] < 180, dsm.cf["longitude"] - 360 + # ) return dsm @@ -208,6 +215,7 @@ def run( only_search=False, only_searchplot=False, parallel=True, + proj=None, readers=None, run_qc=False, skip_units=False, @@ -251,6 +259,8 @@ def run( parallel : boolean, optional Whether to run in parallel with `multiprocessing` library where possible. Default is True. + proj: proj instance + Projection from cartopy. Example: `cartopy.crs.Mercator()`. readers : odg reader or list of readers, optional Can specify which of the available readers to use in your search. Options are odg.erddap, odg.axds, and odg.local. Default is to use all. @@ -336,6 +346,7 @@ def run( # Plot discovered datasets lls_stations, names_stations, lls_box, names_boxes = prep_plot(search) + # import pdb; pdb.set_trace() omsa.map.plot( lls_stations=lls_stations, names_stations=names_stations, @@ -344,6 +355,7 @@ def run( boundary=boundary, res="10m", figname=figname_map, + proj=proj, ) if only_searchplot: @@ -392,10 +404,18 @@ def run( kwargs["T"] = T # xoak doesn't work for 1D lon/lat coords - if dsm.cf["longitude"].ndim == dsm.cf["latitude"].ndim == 1: + if ( + dsm.cf[variable].cf["longitude"].ndim + == dsm.cf[variable].cf["latitude"].ndim + == 1 + ): model_var = dsm.cf[variable].cf.sel(**kwargs).to_dataset() - elif dsm.cf["longitude"].ndim == dsm.cf["latitude"].ndim == 2: + elif ( + dsm.cf[variable].cf["longitude"].ndim + == dsm.cf[variable].cf["latitude"].ndim + == 2 + ): model_var = dsm.cf[variable].em.sel2dcf(**kwargs).to_dataset() # Combine and align the two time series of variable diff --git a/ocean_model_skill_assessor/plot/map.py b/ocean_model_skill_assessor/plot/map.py index bd7c88e..9cbc939 100644 --- a/ocean_model_skill_assessor/plot/map.py +++ b/ocean_model_skill_assessor/plot/map.py @@ -40,7 +40,7 @@ def plot( boundary: array-like, Nx2, optional Model boundary locations. boundary[:,0] longitudes and boundary[:,1] latitudes. - proj: proj instance + proj: proj instance, optional Projection from cartopy. Example: `cartopy.crs.Mercator()`. res: str Resolution of Natural Earth features. Options: '110m', '50m', '10m'. @@ -68,7 +68,7 @@ def plot( names_stations = [names_stations] if not proj: - proj = cartopy.crs.Mercator(central_longitude=central_longitude) + proj = cartopy.crs.Mercator(central_longitude=float(central_longitude)) fig = plt.figure(figsize=(8, 7), dpi=100) ax = fig.add_axes([0.06, 0.01, 0.93, 0.95], projection=proj) # ax.set_frame_on(False) # kind of like it without the box