Trouble with format of fieldset dataset & out of bounds error

[kaylapng]

Question

Hi, I am trying to track particles using ocean currents and drop them in by delaying their start times. What should be the format of my .nc file for my fieldset? Right now, my .nc file downloaded from OSCAR has variables (time, lat, lon, u, v), where u and v have dimensions (time, latitude, and longitude). However, the currents aren't applying, and when I print pset for the final positions, I'm getting the same coordinates as the start.

I'm also getting a field sampled out of bounds error which I'm not understanding because I did not set bounds. If anyone has any ideas it would be much appreciated--i've been stuck for days. Thank you very much for your time and assistance and apologies for the messiness.

Supporting code & Error Message

class MicroplasticParticle(JITParticle):
    concentration = Variable('concentration', dtype=np.float32, initial=0.0)

def convert_to_days_since_0001(unix_seconds):
    # Define the offset in days between 0001-01-01 and 1970-01-01
    offset_days = 719162
    # Convert Unix time (seconds since 1970-01-01) to days
    days_since_epoch = unix_seconds / 86400
    # Add offset to get days since 0001-01-01
    days_since_0001 = days_since_epoch + offset_days
    return days_since_0001
    
# Get latitude and longitude data sorted by time.
def get_sorted_ocean_data():
    u_data = xr.open_dataset('currents_u.nc')
    v_data = xr.open_dataset('currents_v.nc')

    u_data_transposed = u_data['u'].transpose('time', 'latitude', 'longitude')
    v_data_transposed = v_data['v'].transpose('time', 'latitude', 'longitude')
    time_values = u_data['time'].values.astype('datetime64[s]').astype('float')
    time_values2 = convert_to_days_since_0001(time_values)
    print (time_values2)   
   
    lats = u_data['lat'].values
    lons = u_data['lon'].values
    
    return u_data_transposed,v_data_transposed,time_values2,lats,lons

    
# Get mpaa data
def get_mpaa_data():
    mpaa = pd.read_csv('mpaa.csv')
    mpaa['Day'] = 1
    mpaa = mpaa[mpaa['Year'] >= 1993]
    times = pd.to_datetime(mpaa[['Year', 'Month', 'Day']])
    time_as_datetime = [cftime.DatetimeJulian(t.year, t.month, t.day) for t in times]
    # Convert the cftime object to a floating-point timestamp.
    time_as_float = [cftime.date2num(dtJulian, units='days since 0001-01-01', calendar='julian') for dtJulian in time_as_datetime]
    initial_lons = mpaa['Longitude'].values
    initial_lats = mpaa['Latitude'].values
    
    return time_as_float,initial_lons,initial_lats

def main():
    u_data_transposed,v_data_transposed,time_values2,lats,lons = get_sorted_ocean_data();
      
    data = {
        'U': u_data_transposed.values,
        'V': v_data_transposed.values
    }
    dimensions = {
        'time': time_values2,
        'lat': lats,
        'lon': lons
    }

    time_as_float,initial_lons,initial_lats = get_mpaa_data();
    
    # create fieldset
    fieldset = FieldSet.from_data(data, dimensions, allow_time_extrapolation=True)
    # Create particle set
    pset = ParticleSet.from_list(fieldset=fieldset, pclass=MicroplasticParticle,time=time_as_float, lon=initial_lons, lat=initial_lats)
    # Define kernel
    kernels = pset.Kernel(AdvectionRK4)
    
    # Simulation particle tracing
    output_file = pset.ParticleFile(name="microplastics_output.zarr", outputdt=timedelta(days=30))
    
    pset.execute(
        kernels,
        runtime=timedelta(days=365 * 10),  # Simulate 10 years
        dt=timedelta(days=30),              # Timestep of 1 month
        output_file=output_file,
        verbose_progress=True
    )
    
    print("Final particle positions:")
    print(pset.lon, pset.lat)
    
if __name__ == "__main__":
    main()```

## The Error Message

FieldOutOfBoundError                      Traceback (most recent call last)
Cell In[2], line 113
    110     print(pset.lon, pset.lat)
    112 if __name__ == "__main__":
--> 113     main()

Cell In[2], line 100, in main()
     97 # Simulation particle tracing
     98 output_file = pset.ParticleFile(name="microplastics_output.zarr", outputdt=timedelta(days=30))
--> 100 pset.execute(
    101     kernels,
    102     runtime=timedelta(days=365 * 10),  # Simulate 10 years
    103     dt=timedelta(days=30),              # Timestep of 1 month
    104     output_file=output_file,
    105     verbose_progress=True
    106 )
    108 # Output the final position of the particle
    109 print("Final particle positions:")

File ~\AppData\Local\Programs\Python\Python310\lib\site-packages\parcels\particleset.py:1171, in ParticleSet.execute(self, pyfunc, pyfunc_inter, endtime, runtime, dt, output_file, verbose_progress, postIterationCallbacks, callbackdt, delete_cfiles)
   1169 # If we don't perform interaction, only execute the normal kernel efficiently.
   1170 if self._interaction_kernel is None:
-> 1171     res = self._kernel.execute(self, endtime=next_time, dt=dt)
   1172     if res == StatusCode.StopAllExecution:
   1173         return StatusCode.StopAllExecution

File ~\AppData\Local\Programs\Python\Python310\lib\site-packages\parcels\kernel.py:651, in Kernel.execute(self, pset, endtime, dt)
    649     raise TimeExtrapolationError(p.time)
    650 elif p.state == StatusCode.ErrorOutOfBounds:
--> 651     raise FieldOutOfBoundError(p.lon, p.lat, p.depth)
    652 elif p.state == StatusCode.ErrorThroughSurface:
    653     raise FieldOutOfBoundSurfaceError(p.lon, p.lat, p.depth)

FieldOutOfBoundError: Field sampled out-of-bound, at (-71.12999725341797, 41.290000915527344, 0.0)```

[erikvansebille]

Hi @kaylapng, thanks for using Parcels! Your error seems to be because your particles are started in a region of the ocean (71W, 41N) where you don't have velocity data for. What domain is your hydrodynamic data? Could it be that your hydrodynamic data is between 0E and 360E? What happens if you add 360 to your initial longitudes for the particles?

initial_lons = [lon+360 for lon in mpaa['Longitude'].values]