COSIMA / cosima-recipes

A cookbook of recipes (i.e., examples) for analysing ocean and sea ice model output
https://cosima-recipes.readthedocs.io
Apache License 2.0
46 stars 66 forks source link

Plotting maps with correct tripole region -- method in `Making_Maps_with_Cartopy.ipynb` *not* working for 1 degree output #133

Open navidcy opened 3 years ago

navidcy commented 3 years ago

When one tries to follow the tips within the tutorial Making_Maps_with_Cartopy.ipynb for fixing the map polewards of 65 degrees latitude for 1 degree model output it fails miserably.

Here is a copy of the Making_Maps_with_Cartopy.ipynb tutorial a the with the only difference that loads SST from a 1 degree model + the `grid

https://gist.github.com/navidcy/fee1a84f7fa4a05b548191ae36f142d3 or view in nbviewer

Perhaps it's because the 1 degree grid in /g/data/ik11/grids/ includes masked arrays for geolon_t/geolat_t? (I don't know how to check that -- can anybody help me out here?)

(For the record, this was discovered by @fzhang0.)

As a general remark though, the grids in/g/data/ik11/grids/ are inconsistent across the three resolutions and this create lots of trouble in studies where quantities from model output across resolutions are being compared. I would happily regenerate those grid .nc files if I knew how to do that... :(

navidcy commented 3 years ago

btw, @aekiss suggested a fix/hack for when working with 1 degree model output which is:

geolon_t = xr.open_dataset(‘/g/data/ik11/grids/ocean_grid_10.nc’).geolon_t
geolat_t = xr.open_dataset(‘/g/data/ik11/grids/ocean_grid_10.nc’).geolat_t
SST[‘xt_ocean’] = geolon_t[‘xt_ocean’].data
SST[‘yt_ocean’] = geolon_t[‘yt_ocean’].data
SST = SST.assign_coords({‘latitude’: geolat_t, ‘longitude’: geolon_t})

This works but why? Ideally we'd like to have a method that works the same across resolutions...

navidcy commented 3 years ago

cc: @AndyHoggANU, @aidanheerdegen, @angus-g in case they have any ideas.

aidanheerdegen commented 3 years ago

To be clear, are you saying this still works without the two extra lines

SST[‘xt_ocean’] = geolon_t[‘xt_ocean’].data
SST[‘yt_ocean’] = geolon_t[‘yt_ocean’].data

with a 0.25 degree dataset?

navidcy commented 3 years ago

To be clear, are you saying this still works without the two extra lines

SST[‘xt_ocean’] = geolon_t[‘xt_ocean’].data
SST[‘yt_ocean’] = geolon_t[‘yt_ocean’].data

with a 0.25 degree dataset?

Exactly! The current tutorial uses 0.25 degree output and does not include these two lines and everything works OK. But for 1 degree output one gets the errors you see here: https://gist.github.com/navidcy/fee1a84f7fa4a05b548191ae36f142d3 or view in nbviewer

navidcy commented 1 year ago

this is related to #212