Apply filtering to CM2.6 data

[paigem]

This issue is to track any progress made on applying GCM-filters to CM2.6 data.

Pertinent information:

• CM2.6 uses MOM5 ocean model, which is on a B-grid
• The vector Laplacian filtering is still in dev (see issue here)
• For scalar quantities, we can use the implemented MOM5T Laplacian for T-grid variables and MOM5U for U-grid variables

How to access CM2.6 data (more info at #6):

# Load atmosphere data: daily fields
import xarray as xr
import gcsfs
gcs = gcsfs.GCSFileSystem()

path = f'gs://cmip6/GFDL_CM2_6/control/atmos_daily.zarr'
ds_atm = xr.open_zarr(gcs.get_mapper(path))

# Load ocean data: monthly air-sea flux output
from intake import open_catalog
cat = open_catalog("https://raw.githubusercontent.com/pangeo-data/pangeo-datastore/master/intake-catalogs/ocean/GFDL_CM2.6.yaml")
ds_oc  = cat["GFDL_CM2_6_control_ocean_boundary_flux"].to_dask()

# Load ocean model grid data
ds_oc_grid  = cat["GFDL_CM2_6_grid"].to_dask()

Tasks:

• Proof of concept:
  • Can we apply a filter from GCM-filters to (at least a subset of) CM2.6 data?
  • Does the computation work with Dask cluster?
• What variables do we want to filter?
• What scales do we want to filter out?

@dhruvbalwada and I will meet to kick this off tomorrow.

[paigem]

Great chat with @dhruvbalwada today! We discussed a few things:

• We successfully applied the Regular_with_land_area_weighted filter to some test ocean CM2.6 data (did not try using Dask cluster yet)
• We need to better understand the different use cases for the Laplacians in GCM-filters
  • e.g. do we want to use fixed factor filtering (e.g. filtering out scales less than 10x the grid spacing) or fixed length scale filtering (e.g. filter out scales smaller than 100km)?
  • the 2 MOM-specific filters (MOM5T and MOM5U) appear to only work for fixed length scale filtering (see list of all grid types in GCM-filters)
  • we may be able to use the Regular_with_land_area_weighted Laplacian for fixed factor filtering of MOM fields
• The numerical stability and error in the Gaussian filter computation can depend rather sensitively on the filter scale (compared to the grid scale), the number of steps (n_steps), and the precision of the data (float32 vs. float64)
  • @dhruvbalwada has done some diagnostics looking into this

Link to notebook

Question for @jbusecke: For the CM2.6 ocean data, there is a separate zarr store that has the ocean grid variables. Is there a similar dataset that has atmospheric grid variables?

[jbusecke]

Question for @jbusecke: For the CM2.6 ocean data, there is a separate zarr store that has the ocean grid variables. Is there a similar dataset that has atmospheric grid variables?

Good question. Short answer: I did not find anything like that in the original repository. What information would you need exactly? Cell area?

[jbusecke]

Also great progress! This looks awesome. I personally do not understand that much about the intricacies of filters. It might be worth having some sort of synchronous meeting about it?

[paigem]

Good question. Short answer: I did not find anything like that in the original repository. What information would you need exactly? Cell area?

Also a good question. I guess we will need to look at exactly which variables we will need to filter, but if we do want to filter atmospheric fields then yes, we will likely need cell area and potentially other grid variables (e.g. cell width) as well.

But, thinking about this more, I wonder how we plan to handle the different ocean and atmosphere grids (as we discussed in the Miro brainstorming session). The atmospheric variables are already on the coarser atmosphere grid, so would we interpolate first to the high-res ocean grid and then filter as for the ocean fields? Or would we filter the atmospheric fields directly. In the first case, we wouldn't need the atmosphere grid info (at least for the filtering stage - we would need grid info for the interpolation to the ocean grid).

I personally do not understand that much about the intricacies of filters. It might be worth having some sort of synchronous meeting about it?

I agree - it would probably be worthwhile to have a meeting about filtering at some point! I personally find that most of my questions will arise as I code, so maybe we can make it a hacking/discussion session. 🙂

[paigem]

Also @dhruvbalwada - I just noticed a new PR in GCM-filters that addresses the units/dimensionality issue that we were confused about: https://github.com/ocean-eddy-cpt/gcm-filters/pull/146. The PR will clarify whether dx_min should be 1 (non dimensional unit) or a physical length scale.

[jbusecke]

But, thinking about this more, I wonder how we plan to handle the different ocean and atmosphere grids (as we discussed in the Miro brainstorming session). The atmospheric variables are already on the coarser atmosphere grid, so would we interpolate first to the high-res ocean grid and then filter as for the ocean fields? Or would we filter the atmospheric fields directly. In the first case, we wouldn't need the atmosphere grid info (at least for the filtering stage - we would need grid info for the interpolation to the ocean grid).

Uh this is another great question. I would tend to first interpolate, then filter, just to avoid any ambiguity when working with filter details on different grids. If we interpolate first, we know that the filter is doing exactly the same thing on both the ocean and atmospheric fields. Curious what @rabernat thinks?

[paigem]

@jbusecke and I made some progress on this today! 🎉 You can check out the notebook I just added via PR #22.

A few take-aways from today's session where we used the MOM5T Laplacian implemented in GCM-filters to filter a single time slice of MOM5 sensible heat data:

• It took us a bit to understand why this Laplacian is necessary in addition to the 4 listed in this example for the tripolar POP grid (since both POP and MOM5 are B-grid, tripolar grids). We realized the difference is due to the fact that MOM5 uses a northeast grid placement, compared to POP's southwest.
• We ended up masking out the Arctic to avoid having issues dealing with very small grid lengths compared to the 100km filter we applied. It still works if we don't mask out the Arctic, but seems to have a lot of numerical noise in that (and the Antarctic) region. We expect that the wet mask somehow does not line up perfectly, as we can see just a few grid points with length scales smaller than ~1 km that we would expect the wet mask to take care of. Or perhaps this could be dealt with by accounting for the tripolar boundaries? I plan to open up an issue in GCM-filters to discuss ways to get around this.

Also, in parallel, Julius was working on the interpolation code to interpolate atmospheric data onto the ocean grid, and appeared to be successful with that too! Productive day! 😄 💪

[paigem]

Closing since we are not using gcm-filters and have a working workflow already (see notebook in #32 )