NaN-to-zero conversion during interpolation in coastal domain of MOM6

[BoBOcean]

Question

Question

Particles sample unrealistic salinity values during backward tracking with MOM6 dataset. The model data contains only NaN for land cells, but Parcels converts NaN to 0.00 during interpolation when particles approach coastlines. How can I prevent this NaN-to-zero conversion or handle it properly for coastal applications where zero salinity is scientifically invalid?

Supporting code/error messages

import parcels
import numpy as np
from datetime import timedelta
import math
fieldset = parcels.FieldSet.from_netcdf(
    "mom6_hourly_surface_1994.nc",
    {"U": "ssu", "V": "ssv", "S": "sos"},
    {"U": {"lon": "xq", "lat": "yh", "time": "time"},
     "V": {"lon": "xh", "lat": "yq", "time": "time"},
     "S": {"lon": "xh", "lat": "yh", "time": "time"}},
    allow_time_extrapolation=True
)
fieldset.U.interp_method = 'cgrid_velocity'
fieldset.V.interp_method = 'cgrid_velocity'
SampleParticle = parcels.JITParticle.add_variable("salinity")
def SampleS(particle, fieldset, time):
    salt_val = fieldset.S[time, particle.depth, particle.lat, particle.lon]
    if math.isnan(salt_val):
        particle.delete()
    else:
        particle.salinity = salt_val
pset = parcels.ParticleSet.from_list(
    fieldset, pclass=SampleParticle,
    lon=[particle_lon], lat=[particle_lat], depth=[0]
)
pset.execute([parcels.AdvectionRK4, SampleS], runtime=timedelta(days=60), dt=-timedelta(hours=1))

# Model data verification shows no zeros exist:
Original salinity zeros in model: 0
Original NaN values in model: 557845560
# Parcels output shows interpolation artifacts:
Zero salinities in output: 930 points
Salinity range: 0.00 to 32.84
# Tried approaches still produce artifacts:
- cgrid_velocity interpolation for C-grid
- linear_invdist_land_tracer interpolation
- nearest neighbor interpolation  
- landmask-based particle deletion
- NaN deletion during sampling

[michaeldenes]

Hi @BoBOcean, you could try the linear_invdist_land_tracer interpolation method.
See more details here: https://docs.oceanparcels.org/en/latest/examples/tutorial_interpolation.html

[erikvansebille]

Hi @BoBOcean, just for a bit more context: NaNs are indeed converted to zero in Parcels. This is because the interpolators otherwise don't work, so there is no easy way to turn it off. If you want to know, this is implemented at
https://github.com/OceanParcels/Parcels/blob/4073545a8ac4a0a49469277eb7603dd6a4c6cb0a/parcels/field.py#L312

Note that this will change in Parcels v4, which we will release as an alpha-version in the coming weeks.

What is the reason you want to maintain the NaNs for land cells?