'patch' interpolation for salt_sfc_restore.nc

[aekiss]

'patch' interpolation is needed for salt_sfc_restore.nc to avoid blocky fluxes, as discussed here:
https://arccss.slack.com/archives/C6PP0GU9Y/p1515557852000122
also see
#52 (comment)
#71

[aekiss]

@aidanheerdegen there's no need to be conservative as this isn't a flux. So I'd say 'patch' is better because it is smoother.

[aekiss]

@aidanheerdegen Both the 0.1 deg and 0.25 deg need smoothing - e.g. see
/short/v45/aek156/access-om2/input/mom_01deg/salt_sfc_restore.nc
and
/short/v45/aek156/access-om2/input/mom_025deg/salt_sfc_restore.nc
but the 1deg is coarse enough that these artefacts seem (by eye) to be negligible, even near the tripoles.

[aekiss]

There are two separate issues:

1. the underlying WOA13v2 1/4 deg SSS monthly climatology has weird discontinuities and glitches (see Noise and discontinuities in WOA #103) that need to be dealt with first
2. we need to interpolate the corrected WOA13v2 1/4 deg data more smoothly onto the 1/10 and 1/4 degree tripolar grids to avoid blocky fluxes

I've just heard back from NOAA - they are currently investigating 1. I'm happy to delay fixing this until late Feb when I'm back from conferences and we will have hopefully heard more from NOAA by then.

But so I don't forget, this is what I currently think would be the best approach:
1. interpolate WOA13v2 1/4 deg SSS monthly climatology in time (piecewise linear?) so that discontinuities don't persist at full strength for the whole month (using the fact that the glitches change location from month to month). This also avoids jolting the ocean at each change of month (since restoring timescale is currently 10 days). Perhaps we could even match JRA55 temporal sampling (3 hourly?)
2. fill in land mask with interpolation from ocean so land fill values can't contaminate ocean
3. spatially filter, eg 9-point median or mean, to smooth remaining WOA artefacts
4. interpolate onto MOM tripolar grid with patch

[russfiedler]

@aekiss MOM interpolates the surface restoring fields in time on the fly

if (id_restore(index_salt) > 0 .and. salt_damp_factor > 0.0) then call time_interp_external(id_restore(index_salt), Time%model_time, data)
So step 1 is done automatically for you. There shouldn't be any jolts.

[aidanheerdegen]

Good point @russfiedler.

I was going to caution against time interpolation as the first step in any case, as it would propagate a much larger smooth/fill task further down the stack of operations.

[aekiss]

OK great, so we can skip step 1 then.

[aekiss]

This slipped off my radar, but since we are starting IAF runs this would be a good time to resolve it.

[aekiss]

Runoff at ocean points also seems to need smoother interpolation.

e.g. near Canadian tripole in /g/data3/hh5/tmp/cosima/access-om2-01/01deg_jra55v13_iaf/output194/ocean/ocean_month.nc (with stretched colour scale to exaggerate small values):

[russfiedler]

I had a look in Ferret and I don't see some of the features in this plot. For instance that stuff off Melville Island. It looks like it could be aliasing in ncview. The runoff is a bit "blocky" in the sense there is a cutoff to zero but that is really how it should be. Some of the holes in the way it is spread out look a bit odd.
What do the original fields from JRA look like?

[aekiss]

I'm not sure what you mean by aliasing - I am displaying the model output at 1:1 scale.

I've stretched the colourmap to highlight very low values, which might be why they aren't showing up in Ferret?

I think we are seeing two different features here.

Runoff that hits land is mapped to the nearest coastal point and then spread on the model grid via @nichannah's KDTree method until it falls below a prescribed cap. This produces the red semicircle at the Mackenzie River mouth at lower left.

Runoff that hits the ocean and is below the cap seems to be put on the ocean grid by a nearest-value method (i.e. ocean grid points get the runoff data from the nearest JRA55 grid point), so the coarser resolution of JRA55 shows up as blocky features on the finer model grid (distorted by the tripole of course). I think this would be more nicely done by a smooth conservative interpolation.

[aekiss]

This is /g/data/ua8/JRA55-do/latest/runoff_all.2017.nc. The white dots are just ncview's coastline overlay.

[russfiedler]

I first had a look at the fields using an increment of 5e-8 between 0 and 5e-6kg/m^2/s (which is equivalent to 0.4mm/day) and I don't see the blue/cyan jitter that you get. I reckon it's an artifact in ncview. Here's a plot of Melville island with increment 1e-6 (<0.1mm/day) below 2e-5 and 1e-5 up to 1e-4. The graticules indicate index space. It all looks fine to me.

[aekiss]

Hmmm, I think we have focused in on different features of these plots. I'm not sure what you mean by the blue/cyan jitter.
I am concerned by the rectangular steps (of latitude scale ~0.25deg in your plot), which presumably result from JRA55 having coarser resolution than ACCESS-OM2-01.
I was suggesting we use a smoother conservative interpolation method for runoff, which is consistent with what we've done for the other forcing fields.

[russfiedler]

Ok, I see your concern but the difference between runoff/rivers is that we are attempting to distribute what are point sources. We also need to make sure that the maximum occurs by the coast.

I had a look at some areas in the regular portion of the grid and the initial interpolation is definitely NOT doing nearest neighbour interpolation. It looks like it is first order conservative. I think this might be the remapping file /short/public/access-om2/input_236a3011/yatm_01deg/rmp_jrar_to_cict_CONSERV.nc.

Since we have to guarantee positivity and conservation I'm not sure we have other options unless a smoother gets applied in YATM.

[nichannah]

Yes, it is a two step process: 1) conservative 1st order interpolation from JRA55 river grid to mom grid, 2) now on the mom grid moving all 'landed' runoff to the nearest ocean neighbour.

I have recently moved to 2nd order conservative for all flux fields except runoff. I will switch over the runoff now too and we can see whether that makes a difference. I'm guessing it will look a lot smoother.

[nichannah]

There is about output that uses 2nd order conservative for the runoff here: /short/x77/nah599/access-om2/control/minimal_01deg_jra55_ryf/archive/output013

[russfiedler]

@nichannah Can you put the runoff output file somewhere under v45 or in a public area? I can't see x77.

[nichannah]

I've moved this to /short/public/nah599/ocean_month.nc

[russfiedler]

Ta, but as I expected it goes negative.

yes? stats runoff

         mass flux of liquid river runoff entering ocean
         LONGITUDE: 80E(-280) to 80E(80)
         LATITUDE: 81.1S to 90N
         Z:  N/A
         TIME: 16-AUG-1900 12:00 NOLEAP
         E:  N/A
         F:  N/A
         DATA SET: ./ocean_month.nc

Total # of data points: 9720000 (36002700111*1)
# flagged as bad data: 3644761
Minimum value: -0.013164
Maximum value: 0.03
Mean value: 0.0000014966 (unweighted average)
Standard deviation: 0.000090934

It goes negative at 170K points!

yes? let negrun= if runoff lt 0 then runoff
yes? list negrun[i=@ngd,j=@ngd]
!-> list 170693
VARIABLE : constant
170693.

[russfiedler]

We see that there is negative offshore runoff pretty well everywhere.

[aekiss]

Just saving Russ' comment from the Slack event horizon: https://arccss.slack.com/archives/C9Q7Y1400/p1547423551089000
"I'm just reading a paper that that may adress my concerns about the non positive definite character of second order interpolation. It seems to introduce a "repair" step that fixes any local problems like negative values.
https://cnls.lanl.gov/~shashkov/papers/jcp_report.ps
Section 2.3 is the bit of interest.
Ok, nic what they do is similar to how you remap excess runoff. They check if things are wrong and then "steal" from the surrounding 8 neighbours. If there isn't enough then they expand the domain to be borrowed from."

[aekiss]

For salt_sfc_restore.nc I think we should stick with WOA13v2 for the new RYF 0.1deg config https://github.com/OceansAus/access-om2/projects/1 for consistency with previous runs.

I propose we

1. start with /g/data1/v45/aph502/jra55/sss/mean_salinity_of_0_and_10m.nc
   (this is the basis of our current salt_sfc_restore.nc; it was created using Ferret (ostensibly to fill in the continents) from /g/data/v45/akm157/data/WOA13v2/averaged_decades/woa13_decav_sss_clim_04v2_0m_10m.nc which Adele created and contains an average of 0m and 10m data from WOA13v2... I think)
2. spatially filter this with (say) a 3x3-point running mean, to smooth remaining WOA artefacts (see Noise and discontinuities in WOA #103)
3. interpolate this onto MOM tripolar grid with patch

[aekiss]

Thanks @nichannah - is there a new salt_sfc_restore.nc we can try?

[aekiss]

@nichannah's new salt_sfc_restore.nc seems to do the trick.

Here's what the restoring flux used to look like in /g/data3/hh5/tmp/cosima/access-om2-01/01deg_jra55v13_ryf8485_spinup6_newtopog/output436/ocean/ocean_month.nc which used the old /short/public/access-om2/input_38570c62/mom_01deg/salt_sfc_restore.nc

And here's what it looks like for the same month in /g/data3/hh5/tmp/cosima/access-om2-01/01deg_jra55v13_ryf8485_spinup6_newtopog_scalewind/output435/ocean/ocean_month.nc which used the new /short/public/access-om2/input_rc/mom_01deg/salt_sfc_restore.nc

[StephenGriffies]

Great job!

Are there remaining issues to resolve for the MOM01, or are we now ready to rock and roll?

[nichannah]

Can we close this issue now?

[aekiss]

Not quite, there's one last thing: /short/public/access-om2/input_08022019/mom_025deg/salt_sfc_restore.nc should (ideally) be re-done with patch interpolation as it is slightly blocky near the tripoles.

(We'll have to live with the blocky runoff until we have a positive-definite conservative smooth interpolation method.)

[nichannah]

A new salt_sfc_restore.nc has been placed at /short/public/access-om2/input_rc/mom_025deg/salt_sfc_restore.nc

The attached image shows some differences between the old (on the left) and the new (right).

@aekiss if you can confirm that this a fix then we will close the issue.

[aekiss]

Thanks @nichannah that looks great. I'll leave this issue open as a reminder until I've set up the default 0.25 deg configurations to use this.

[russfiedler]

It's a curse I tell you!

time:units = "months since 0000-01-01 00:00:00" ;

The year 0000 doesn't exist and results in hundreds->thousands of redundant warning messages as MOM/FMS changes it to 0001

[aekiss]

Would 1900 be safer? #117

[russfiedler]

It's specified on a NOLEAP calendar so it shouldn't matter. FMS/MOM handles this sort of thing internally if there's a calendar mismatch in certain circumstances and it's possibly safer to avoid -ve dates. For instance if some wants to start a run over the historical period starting in the 19th century.

I like the 0001 starting date because it feels like a 'natural' start :)

[nichannah]

Sorry, the meta-data is copied over. I'll fix it and replace the file.

[aekiss]

I've fixed this by doing

ncatted -O -a units,time,m,char,"months since 0001-01-01 00:00:00" /short/public/access-om2/input_rc/mom_025deg/salt_sfc_restore.nc

(the original is in /short/public/access-om2/input_rc/mom_025deg/salt_sfc_restore-bkp.nc)

@nichannah, @russfiedler - is this all that was needed?

also see #49

[russfiedler]

That looks right.