Update our analysis code to only smooth winds or tracers

[paigem]

To isolate the largest contributors to the small scale differences, we discussed smoothing only the winds or only the tracers. This may take some thought on how to implement efficiently.

This was based on discussions on Friday, Sep. 30th (notes here).

[jbusecke]

I have thought more about this. I currently think that the 'cleanest' way of doing this would be to 'degrade' the effective resolution (while leaving the actual resolution unchanged) would be a two step process:

1. Coarsen the data (either with a naive .coarsen or by conservative xesmf interpolation to e.g. a 2x2 grid)
2. Linearly interpolate these values back to the initial grid using xesmf.

We could also consider a nearest neighbor interpolation in step 2, but I would be worried about the strong edges we would create with this.

[jbusecke]

Ok I might need some help from the gcm-filters pro users here:
I set up a basic regular filter which should respect land like this:

import gcm_filters
wet_mask = ds_oc_grid.wet
simple_filter = gcm_filters.Filter(
    filter_scale=40,
    dx_min=1,
    filter_shape=gcm_filters.FilterShape.TAPER,
    grid_type=gcm_filters.GridType.REGULAR_WITH_LAND,
    grid_vars={'wet_mask': wet_mask}
)

And then I am trying to filter the following field:

with

da_filtered = simple_filter.apply(ds_merged.surface_temp, dims=['yt_ocean', 'xt_ocean'])

This gives me a smoothed field as desired, but there are tons of artifacts near the coasts.

This is totally undesirable for our purposes, since this seems to influence the signal into the WBC. Am I making some obvious mistake here?

[rabernat]

Probably a dumb question, but are we sure that the wet mask is correct?

[jbusecke]

This is the wet mask (which looks ok to me):

Could it be that this behavior is caused by the large values (K not C)? If I removed a constant value (300) from the data and then smooth the effect seems to be much smaller:

[jbusecke]

Hmm upon further inspection there are a few data points that do not perfectly overlap with the map. Let me investigate a bit more...

[jbusecke]

Oh wow, that was it!
Rebuilding the mask from np.insnan() gives this result:

Much better! And a lesson learned about gcm-filters haha.

[jbusecke]

I have worked a fair bit on the smoothing procedure over the past two days, and I think I have a satisfying solution for now:

As @paigem experienced in her other project, a filter scale of e.g. 10 grid boxes, does not at all produce smoothed fields that would resemble the amount of detail a 1 deg model would produce (eddies are still very much visible, see below). So in order to separate the mesoscale and smaller we need to use much larger filter_scales.
We do however want to avoid too large* filter scales too, since that could pick up larger 'regional' features which might hinder physical interpretability of the results.

To visually find the sweet-spot I have made some comparison plots with a range of filter scales:

The top left image is the raw model SST output. Below each row shows the smoothed output for a particular filter scale (see title) on the left, and the difference from the raw field on the right.

I think that filter_scale=50 is the sweet spot, where the eddies are not visible anymore in the smoothed field.

@ocean-transport/collab_team please let me know if you think different.

Also @TomNicholas mentioned that for a publication we might want to rerun our results with slightly larger/smaller filter scales to test the robustness of the results.

Here is the notebook to reproduce the results.

[dhruvbalwada]

The general rule of thumb is that the resolution of the model is about 5times the grid size. So if dx=10km, then the smallest features resolved are roughly 50km in size. This is due to the intrinsic dissipation in the model.

It seems like a filter scale of 50 would be reasonable.