diff --git a/parcels/particleset.py b/parcels/particleset.py index 8ff352322..e13ba81ee 100644 --- a/parcels/particleset.py +++ b/parcels/particleset.py @@ -347,221 +347,7 @@ def populate_indices(self): self._data["ei"][:, i] = idx # assumes that we are in the surface layer (zi=0) @classmethod - def from_list( - cls, fieldset, pclass, lon, lat, depth=None, time=None, repeatdt=None, lonlatdepth_dtype=None, **kwargs - ): - """Initialise the ParticleSet from lists of lon and lat. - - Parameters - ---------- - fieldset : - mod:`parcels.fieldset.FieldSet` object from which to sample velocity - pclass : - Particle class. May be a parcels.particle.Particle class as defined in parcels, or a subclass defining a custom particle. - lon : - List of initial longitude values for particles - lat : - List of initial latitude values for particles - depth : - Optional list of initial depth values for particles. Default is 0m - time : - Optional list of start time values for particles. Default is fieldset.U.time[0] - repeatdt : - Optional interval (in seconds) on which to repeat the release of the ParticleSet (Default value = None) - lonlatdepth_dtype : - Floating precision for lon, lat, depth particle coordinates. - It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear' - and np.float64 if the interpolation method is 'cgrid_velocity' - Other Variables can be initialised using further arguments (e.g. v=... for a Variable named 'v') - **kwargs : - Keyword arguments passed to the particleset constructor. - """ - return cls( - fieldset=fieldset, - pclass=pclass, - lon=lon, - lat=lat, - depth=depth, - time=time, - repeatdt=repeatdt, - lonlatdepth_dtype=lonlatdepth_dtype, - **kwargs, - ) - - @classmethod - def from_line( - cls, - fieldset, - pclass, - start, - finish, - size, - depth=None, - time=None, - repeatdt=None, - lonlatdepth_dtype=None, - **kwargs, - ): - """Create a particleset in the shape of a line (according to a cartesian grid). - - Initialise the ParticleSet from start/finish coordinates with equidistant spacing - Note that this method uses simple numpy.linspace calls and does not take into account - great circles, so may not be a exact on a globe - - Parameters - ---------- - fieldset : - mod:`parcels.fieldset.FieldSet` object from which to sample velocity - pclass : - Particle class. May be a parcels.particle.Particle as defined in parcels, or a subclass defining a custom particle. - start : - Start point (longitude, latitude) for initialisation of particles on a straight line. - finish : - End point (longitude, latitude) for initialisation of particles on a straight line. - size : - Initial size of particle set - depth : - Optional list of initial depth values for particles. Default is 0m - time : - Optional start time value for particles. Default is fieldset.U.time[0] - repeatdt : - Optional interval (in seconds) on which to repeat the release of the ParticleSet (Default value = None) - lonlatdepth_dtype : - Floating precision for lon, lat, depth particle coordinates. - It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear' - and np.float64 if the interpolation method is 'cgrid_velocity' - """ - lon = np.linspace(start[0], finish[0], size) - lat = np.linspace(start[1], finish[1], size) - if type(depth) in [int, float]: - depth = [depth] * size - return cls( - fieldset=fieldset, - pclass=pclass, - lon=lon, - lat=lat, - depth=depth, - time=time, - repeatdt=repeatdt, - lonlatdepth_dtype=lonlatdepth_dtype, - **kwargs, - ) - - @classmethod - def _monte_carlo_sample(cls, start_field, size, mode="monte_carlo"): - """Converts a starting field into a monte-carlo sample of lons and lats. - - Parameters - ---------- - start_field : parcels.field.Field - mod:`parcels.fieldset.Field` object for initialising particles stochastically (horizontally) according to the presented density field. - size : - - mode : - (Default value = 'monte_carlo') - - Returns - ------- - list of float - A list of longitude values. - list of float - A list of latitude values. - """ - if mode == "monte_carlo": - data = start_field.data if isinstance(start_field.data, np.ndarray) else np.array(start_field.data) - if start_field.interp_method == "cgrid_tracer": - p_interior = np.squeeze(data[0, 1:, 1:]) - else: # if A-grid - d = data - p_interior = (d[0, :-1, :-1] + d[0, 1:, :-1] + d[0, :-1, 1:] + d[0, 1:, 1:]) / 4.0 - p_interior = np.where(d[0, :-1, :-1] == 0, 0, p_interior) - p_interior = np.where(d[0, 1:, :-1] == 0, 0, p_interior) - p_interior = np.where(d[0, 1:, 1:] == 0, 0, p_interior) - p_interior = np.where(d[0, :-1, 1:] == 0, 0, p_interior) - p = np.reshape(p_interior, (1, p_interior.size)) - inds = np.random.choice(p_interior.size, size, replace=True, p=p[0] / np.sum(p)) - xsi = np.random.uniform(size=len(inds)) - eta = np.random.uniform(size=len(inds)) - j, i = np.unravel_index(inds, p_interior.shape) - grid = start_field.grid - lon, lat = ([], []) - if grid._gtype in [GridType.RectilinearZGrid, GridType.RectilinearSGrid]: - lon = grid.lon[i] + xsi * (grid.lon[i + 1] - grid.lon[i]) - lat = grid.lat[j] + eta * (grid.lat[j + 1] - grid.lat[j]) - else: - lons = np.array([grid.lon[j, i], grid.lon[j, i + 1], grid.lon[j + 1, i + 1], grid.lon[j + 1, i]]) - if grid.mesh == "spherical": - lons[1:] = np.where(lons[1:] - lons[0] > 180, lons[1:] - 360, lons[1:]) - lons[1:] = np.where(-lons[1:] + lons[0] > 180, lons[1:] + 360, lons[1:]) - lon = ( - (1 - xsi) * (1 - eta) * lons[0] - + xsi * (1 - eta) * lons[1] - + xsi * eta * lons[2] - + (1 - xsi) * eta * lons[3] - ) - lat = ( - (1 - xsi) * (1 - eta) * grid.lat[j, i] - + xsi * (1 - eta) * grid.lat[j, i + 1] - + xsi * eta * grid.lat[j + 1, i + 1] - + (1 - xsi) * eta * grid.lat[j + 1, i] - ) - return list(lat), list(lon) - else: - raise NotImplementedError(f'Mode {mode} not implemented. Please use "monte carlo" algorithm instead.') - - @classmethod - def from_field( - cls, - fieldset, - pclass, - start_field, - size, - mode="monte_carlo", - depth=None, - time=None, - repeatdt=None, - lonlatdepth_dtype=None, - ): - """Initialise the ParticleSet randomly drawn according to distribution from a field. - - Parameters - ---------- - fieldset : parcels.fieldset.FieldSet - mod:`parcels.fieldset.FieldSet` object from which to sample velocity - pclass : - Particle class. May be a parcels.particle.Particle class as defined in parcels, or a subclass defining a custom particle. - start_field : parcels.field.Field - Field for initialising particles stochastically (horizontally) according to the presented density field. - size : - Initial size of particle set - mode : - Type of random sampling. Currently only 'monte_carlo' is implemented (Default value = 'monte_carlo') - depth : - Optional list of initial depth values for particles. Default is 0m - time : - Optional start time value for particles. Default is fieldset.U.time[0] - repeatdt : - Optional interval (in seconds) on which to repeat the release of the ParticleSet (Default value = None) - lonlatdepth_dtype : - Floating precision for lon, lat, depth particle coordinates. - It is either np.float32 or np.float64. Default is np.float32 if fieldset.U.interp_method is 'linear' - and np.float64 if the interpolation method is 'cgrid_velocity' - """ - lat, lon = cls._monte_carlo_sample(start_field, size, mode) - - return cls( - fieldset=fieldset, - pclass=pclass, - lon=lon, - lat=lat, - depth=depth, - time=time, - lonlatdepth_dtype=lonlatdepth_dtype, - repeatdt=repeatdt, - ) - - @classmethod - def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime=None, repeatdt=None, **kwargs): + def from_particlefile(cls, fieldset, pclass, filename, restart=True, restarttime=None, **kwargs): """Initialise the ParticleSet from a zarr ParticleFile. This creates a new ParticleSet based on locations of all particles written in a zarr ParticleFile at a certain time. Particle IDs are preserved if restart=True diff --git a/tests/test_fieldset.py b/tests/test_fieldset.py index 97c632d9d..fd83a28be 100644 --- a/tests/test_fieldset.py +++ b/tests/test_fieldset.py @@ -343,7 +343,7 @@ def sampleTemp(particle, fieldset, time): # pragma: no cover ] ) - pset = ParticleSet.from_list(fieldset, pclass=MyParticle, lon=[0.5], lat=[0.5], depth=[0.5]) + pset = ParticleSet(fieldset, pclass=MyParticle, lon=[0.5], lat=[0.5], depth=[0.5]) pset.execute( AdvectionRK4_3D + pset.Kernel(sampleTemp), runtime=timedelta(hours=51), dt=timedelta(hours=dt_sign * 1) ) diff --git a/tests/test_particlesets.py b/tests/test_particlesets.py index e1579123c..fc186519d 100644 --- a/tests/test_particlesets.py +++ b/tests/test_particlesets.py @@ -32,7 +32,7 @@ def test_pset_create_list_with_customvariable(fieldset): MyParticle = Particle.add_variable("v") v_vals = np.arange(npart) - pset = ParticleSet.from_list(fieldset, lon=lon, lat=lat, v=v_vals, pclass=MyParticle) + pset = ParticleSet(fieldset, lon=lon, lat=lat, v=v_vals, pclass=MyParticle) assert np.allclose([p.lon for p in pset], lon, rtol=1e-12) assert np.allclose([p.lat for p in pset], lat, rtol=1e-12) assert np.allclose([p.v for p in pset], v_vals, rtol=1e-12) @@ -127,7 +127,7 @@ def test_pset_create_with_time(fieldset): time = 5.0 pset = ParticleSet(fieldset, lon=lon, lat=lat, pclass=Particle, time=time) assert np.allclose([p.time for p in pset], time, rtol=1e-12) - pset = ParticleSet.from_list(fieldset, lon=lon, lat=lat, pclass=Particle, time=[time] * npart) + pset = ParticleSet(fieldset, lon=lon, lat=lat, pclass=Particle, time=[time] * npart) assert np.allclose([p.time for p in pset], time, rtol=1e-12) pset = ParticleSet.from_line(fieldset, size=npart, start=(0, 1), finish=(1, 0), pclass=Particle, time=time) assert np.allclose([p.time for p in pset], time, rtol=1e-12) diff --git a/tests/v4/test_particleset.py b/tests/v4/test_particleset.py index 714574d91..d21322790 100644 --- a/tests/v4/test_particleset.py +++ b/tests/v4/test_particleset.py @@ -37,19 +37,6 @@ def test_pset_create_lon_lat(fieldset): assert np.allclose([p.lat for p in pset], lat, rtol=1e-12) -@pytest.mark.parametrize("lonlatdepth_dtype", [np.float64, np.float32]) -def test_pset_create_line(fieldset, lonlatdepth_dtype): - npart = 100 - lon = np.linspace(0, 1, npart, dtype=lonlatdepth_dtype) - lat = np.linspace(1, 0, npart, dtype=lonlatdepth_dtype) - pset = ParticleSet.from_line( - fieldset, size=npart, start=(0, 1), finish=(1, 0), pclass=Particle, lonlatdepth_dtype=lonlatdepth_dtype - ) - assert np.allclose([p.lon for p in pset], lon, rtol=1e-12) - assert np.allclose([p.lat for p in pset], lat, rtol=1e-12) - assert isinstance(pset[0].lat, lonlatdepth_dtype) - - def test_create_empty_pset(fieldset): pset = ParticleSet(fieldset, pclass=Particle) assert pset.size == 0 diff --git a/tests/v4/test_structured_gcm.py b/tests/v4/test_structured_gcm.py index fa8eee40a..f724cff9b 100644 --- a/tests/v4/test_structured_gcm.py +++ b/tests/v4/test_structured_gcm.py @@ -12,6 +12,6 @@ def test_fieldset_frompop(): # dimensions = {"lon": "lon", "lat": "lat", "time": "time"} # fieldset = FieldSet.from_pop(filenames, variables, dimensions, mesh="flat") - # pset = ParticleSet.from_list(fieldset, Particle, lon=[3, 5, 1], lat=[3, 5, 1]) + # pset = ParticleSet(fieldset, Particle, lon=[3, 5, 1], lat=[3, 5, 1]) # pset.execute(AdvectionRK4, runtime=3, dt=1) pass diff --git a/v3to4-breaking-changes.md b/v3to4-breaking-changes.md index 1dc1094ae..c6e6979f4 100644 --- a/v3to4-breaking-changes.md +++ b/v3to4-breaking-changes.md @@ -12,3 +12,4 @@ FieldSet ParticleSet - ParticleSet.execute() expects `numpy.datetime64`/`numpy.timedelta.64` for `runtime`, `endtime` and `dt` +- `ParticleSet.from_field()`, `ParticleSet.from_line()`, `ParticleSet.from_list()` has been removed