From 58fca42e0fb2968208c13de35eb82ee3e1d55a5e Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 28 May 2024 16:33:17 +0200 Subject: [PATCH 01/40] wip drifters --- virtual_ship/instruments/drifter.py | 85 +++++++++++++++++++++++++++++ 1 file changed, 85 insertions(+) create mode 100644 virtual_ship/instruments/drifter.py diff --git a/virtual_ship/instruments/drifter.py b/virtual_ship/instruments/drifter.py new file mode 100644 index 00000000..adb48e12 --- /dev/null +++ b/virtual_ship/instruments/drifter.py @@ -0,0 +1,85 @@ +"""Argo float instrument.""" + +import math +from dataclasses import dataclass +from datetime import timedelta + +import numpy as np +from parcels import ( + AdvectionRK4, + FieldSet, + JITParticle, + ParticleSet, + StatusCode, + Variable, +) + +from .location import Location + + +@dataclass +class Drifter: + """Configuration for a single Argo float.""" + + location: Location + deployment_time: float + + +class _ArgoParticle(JITParticle): + temperature = Variable( + "temperature", + initial=np.nan, + ) + + +def simulate_drifters( + drifters: list[Drifter], + fieldset: FieldSet, + out_file_name: str, + outputdt: timedelta, +) -> None: + """ + Use parcels to simulate a set of Argo floats in a fieldset. + + :param argo_floats: A list of Argo floats to simulate. + :param fieldset: The fieldset to simulate the Argo floats in. + :param out_file_name: The file to write the results to. + :param outputdt: Interval which dictates the update frequency of file output during simulation + """ + lon = [argo.location.lon for argo in argo_floats] + lat = [argo.location.lat for argo in argo_floats] + time = [argo.deployment_time for argo in argo_floats] + + # define the parcels simulation + argo_float_fieldset = ParticleSet( + fieldset=fieldset, + pclass=_ArgoParticle, + lon=lon, + lat=lat, + depth=[argo.min_depth for argo in argo_floats], + time=time, + min_depth=[argo.min_depth for argo in argo_floats], + max_depth=[argo.max_depth for argo in argo_floats], + drift_depth=[argo.drift_depth for argo in argo_floats], + vertical_speed=[argo.vertical_speed for argo in argo_floats], + cycle_days=[argo.cycle_days for argo in argo_floats], + drift_days=[argo.drift_days for argo in argo_floats], + ) + + # define output file for the simulation + out_file = argo_float_fieldset.ParticleFile( + name=out_file_name, + outputdt=timedelta(minutes=5), + chunks=(1, 500), + ) + + # get time when the fieldset ends + fieldset_endtime = fieldset.time_origin.fulltime(fieldset.U.grid.time_full[-1]) + + # execute simulation + drifter_fieldset.execute( + [AdvectionRK4, _sample_temperature, _check_error], + endtime=fieldset_endtime, + dt=outputdt, + output_file=out_file, + ) From 3524ab327d078fc1548e899e50bbee7405f1e76b Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 28 May 2024 16:44:45 +0200 Subject: [PATCH 02/40] drifters instrument. just missing test --- virtual_ship/instruments/__init__.py | 4 +-- virtual_ship/instruments/drifter.py | 54 +++++++++++++--------------- 2 files changed, 27 insertions(+), 31 deletions(-) diff --git a/virtual_ship/instruments/__init__.py b/virtual_ship/instruments/__init__.py index 283e7cba..261bd7a8 100644 --- a/virtual_ship/instruments/__init__.py +++ b/virtual_ship/instruments/__init__.py @@ -1,6 +1,6 @@ """Measurement instrument that can be used with Parcels.""" -from . import argo_float +from . import argo_float, drifter from .location import Location -__all__ = ["Location", "argo_float"] +__all__ = ["Location", "argo_float", "drifter"] diff --git a/virtual_ship/instruments/drifter.py b/virtual_ship/instruments/drifter.py index adb48e12..c12215fa 100644 --- a/virtual_ship/instruments/drifter.py +++ b/virtual_ship/instruments/drifter.py @@ -1,37 +1,39 @@ -"""Argo float instrument.""" +"""Drifter instrument.""" -import math from dataclasses import dataclass from datetime import timedelta import numpy as np -from parcels import ( - AdvectionRK4, - FieldSet, - JITParticle, - ParticleSet, - StatusCode, - Variable, -) +from parcels import AdvectionRK4, FieldSet, JITParticle, ParticleSet, Variable from .location import Location @dataclass class Drifter: - """Configuration for a single Argo float.""" + """Configuration for a single Drifter.""" location: Location deployment_time: float + min_depth: float -class _ArgoParticle(JITParticle): +class _DrifterParticle(JITParticle): temperature = Variable( "temperature", initial=np.nan, ) +def _sample_temperature(particle, fieldset, time): + particle.temperature = fieldset.T[time, particle.depth, particle.lat, particle.lon] + + +def _check_error(particle, fieldset, time): + if particle.state >= 50: # This captures all Errors + particle.delete() + + def simulate_drifters( drifters: list[Drifter], fieldset: FieldSet, @@ -39,35 +41,29 @@ def simulate_drifters( outputdt: timedelta, ) -> None: """ - Use parcels to simulate a set of Argo floats in a fieldset. + Use parcels to simulate a set of drifters in a fieldset. - :param argo_floats: A list of Argo floats to simulate. - :param fieldset: The fieldset to simulate the Argo floats in. + :param drifters: A list of drifters to simulate. + :param fieldset: The fieldset to simulate the drifters in. :param out_file_name: The file to write the results to. :param outputdt: Interval which dictates the update frequency of file output during simulation """ - lon = [argo.location.lon for argo in argo_floats] - lat = [argo.location.lat for argo in argo_floats] - time = [argo.deployment_time for argo in argo_floats] + lon = [drifter.location.lon for drifter in drifters] + lat = [drifter.location.lat for drifter in drifters] + time = [drifter.deployment_time for drifter in drifters] # define the parcels simulation - argo_float_fieldset = ParticleSet( + drifter_particleset = ParticleSet( fieldset=fieldset, - pclass=_ArgoParticle, + pclass=_DrifterParticle, lon=lon, lat=lat, - depth=[argo.min_depth for argo in argo_floats], + depth=[drifter.min_depth for drifter in drifters], time=time, - min_depth=[argo.min_depth for argo in argo_floats], - max_depth=[argo.max_depth for argo in argo_floats], - drift_depth=[argo.drift_depth for argo in argo_floats], - vertical_speed=[argo.vertical_speed for argo in argo_floats], - cycle_days=[argo.cycle_days for argo in argo_floats], - drift_days=[argo.drift_days for argo in argo_floats], ) # define output file for the simulation - out_file = argo_float_fieldset.ParticleFile( + out_file = drifter_particleset.ParticleFile( name=out_file_name, outputdt=timedelta(minutes=5), chunks=(1, 500), @@ -77,7 +73,7 @@ def simulate_drifters( fieldset_endtime = fieldset.time_origin.fulltime(fieldset.U.grid.time_full[-1]) # execute simulation - drifter_fieldset.execute( + drifter_particleset.execute( [AdvectionRK4, _sample_temperature, _check_error], endtime=fieldset_endtime, dt=outputdt, From 9c00fcf334e0a70ebc8c45695066e88a4f957ac3 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 28 May 2024 17:03:34 +0200 Subject: [PATCH 03/40] drifter test --- tests/instruments/test_drifters.py | 35 ++++++++++++++++++++++++++++++ 1 file changed, 35 insertions(+) create mode 100644 tests/instruments/test_drifters.py diff --git a/tests/instruments/test_drifters.py b/tests/instruments/test_drifters.py new file mode 100644 index 00000000..0af64143 --- /dev/null +++ b/tests/instruments/test_drifters.py @@ -0,0 +1,35 @@ +"""Test the simulation of drifters.""" + +from virtual_ship.instruments.drifter import simulate_drifters, Drifter +from virtual_ship.instruments import Location +from parcels import FieldSet +import numpy as np +from datetime import timedelta + + +def test_simulate_drifters() -> None: + fieldset = FieldSet.from_data( + {"U": 0, "V": 0, "T": 0}, + { + "lon": 0, + "lat": 0, + "time": [np.datetime64("1950-01-01") + np.timedelta64(632160, "h")], + }, + ) + + min_depth = -fieldset.U.depth[0] + + drifters = [ + Drifter( + location=Location(latitude=0, longitude=0), + deployment_time=0, + min_depth=min_depth, + ) + ] + + simulate_drifters( + drifters=drifters, + fieldset=fieldset, + out_file_name="test", + outputdt=timedelta(minutes=5), + ) From 695c73eb73db529a98d094587231ccd9293dc3a1 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 28 May 2024 17:24:06 +0200 Subject: [PATCH 04/40] Remove old drifter file and fix errors --- tests/sailship_config.json | 2 +- tests/test_sailship.py | 7 +- virtual_ship/drifter_deployments.py | 109 ---------------------------- virtual_ship/sailship.py | 33 ++++++--- 4 files changed, 29 insertions(+), 122 deletions(-) delete mode 100644 virtual_ship/drifter_deployments.py diff --git a/tests/sailship_config.json b/tests/sailship_config.json index 23cf4c1a..d7b90b05 100644 --- a/tests/sailship_config.json +++ b/tests/sailship_config.json @@ -25,7 +25,7 @@ "max_depth": "max" }, "drifter_deploylocations": [ - + [-23.081289, 63.743631] ], "argo_deploylocations": [ [-23.081289, 63.743631] diff --git a/tests/test_sailship.py b/tests/test_sailship.py index afd8a95f..33f1e290 100644 --- a/tests/test_sailship.py +++ b/tests/test_sailship.py @@ -17,11 +17,12 @@ def test_sailship() -> None: ctd_fieldset.add_constant("max_depth", -ctd_fieldset.U.depth[-1]) drifter_fieldset = FieldSet.from_data( + {"U": 0, "V": 0, "T": 0}, { - "U": 0, - "V": 0, + "lon": 0, + "lat": 0, + "time": [np.datetime64("1950-01-01") + np.timedelta64(632160, "h")], }, - {"lon": 0, "lat": 0}, ) argo_float_fieldset = FieldSet.from_data( diff --git a/virtual_ship/drifter_deployments.py b/virtual_ship/drifter_deployments.py deleted file mode 100644 index 69836bc3..00000000 --- a/virtual_ship/drifter_deployments.py +++ /dev/null @@ -1,109 +0,0 @@ -"""drifter_deployments function.""" - -import datetime -import os -from datetime import timedelta - -import numpy as np -from parcels import AdvectionRK4, FieldSet, JITParticle, ParticleSet, Variable - -from .virtual_ship_configuration import VirtualShipConfiguration - - -def drifter_deployments(config: VirtualShipConfiguration, drifter_time): - """ - Deploy drifters. - - :param config: The cruise configuration. - :param drifter_time: TODO - """ - if len(config.drifter_deploylocations) > 0: - - fieldset = config.drifter_fieldset - - # Create particle to sample water underway - class DrifterParticle(JITParticle): - """Define a new particle class that samples Temperature as a surface drifter.""" - - temperature = Variable("temperature", initial=np.nan) - - # define function sampling Temperature - def SampleT(particle, fieldset, time): - particle.temperature = fieldset.T[ - time, particle.depth, particle.lat, particle.lon - ] - - def CheckError(particle, fieldset, time): - if particle.state >= 50: # This captures all Errors - particle.delete() - - # initialize drifters - lon = [] - lat = [] - for i in range(len(config.drifter_deploylocations)): - lon.append(config.drifter_deploylocations[i][0]) - lat.append(config.drifter_deploylocations[i][1]) - time = drifter_time - - # Create and execute drifter particles - pset = ParticleSet( - fieldset=fieldset, - pclass=DrifterParticle, - lon=lon, - lat=lat, - depth=np.repeat(fieldset.mindepth, len(time)), - time=time, - ) - output_file = pset.ParticleFile( - name=os.path.join("results", "Drifters.zarr"), outputdt=timedelta(hours=1) - ) - - fieldset_endtime = fieldset.time_origin.fulltime(fieldset.U.grid.time_full[-1]) - drifter_endtime = np.array( - ( - datetime.datetime.strptime( - config.requested_ship_time["start"], "%Y-%m-%dT%H:%M:%S" - ) - + timedelta(weeks=6) - ) - ).astype("datetime64[ms]") - - pset.execute( - [AdvectionRK4, SampleT, CheckError], - endtime=min(fieldset_endtime, drifter_endtime), - dt=timedelta(minutes=5), - output_file=output_file, - ) - - -def create_drifter_fieldset(config, data_dir: str): - """ - Create a fieldset from netcdf files for drifters, returns fieldset with negative depth values. - - :param config: The cruise configuration. - :param data_dir: TODO - :returns: The fieldset. - """ - filenames = { - "U": os.path.join(data_dir, "drifterdata_UV.nc"), - "V": os.path.join(data_dir, "drifterdata_UV.nc"), - "T": os.path.join(data_dir, "drifterdata_T.nc"), - } - variables = {"U": "uo", "V": "vo", "T": "thetao"} - dimensions = { - "lon": "longitude", - "lat": "latitude", - "time": "time", - "depth": "depth", - } - - # create the fieldset and set interpolation methods - fieldset = FieldSet.from_netcdf( - filenames, variables, dimensions, allow_time_extrapolation=False - ) - fieldset.T.interp_method = "linear_invdist_land_tracer" - for g in fieldset.gridset.grids: - if max(g.depth) > 0: - g.depth = -g.depth # make depth negative - fieldset.mindepth = -fieldset.U.depth[0] # uppermost layer in the hydrodynamic data - return fieldset diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 51238069..ed0bc082 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -9,8 +9,8 @@ from shapely.geometry import Point, Polygon from .costs import costs -from .drifter_deployments import drifter_deployments from .instruments.argo_float import ArgoFloat, simulate_argo_floats +from .instruments.drifter import Drifter, simulate_drifters from .instruments.location import Location from .postprocess import postprocess from .virtual_ship_configuration import VirtualShipConfiguration @@ -137,11 +137,12 @@ def SampleT(particle, fieldset, time): # initialize drifters and argo floats drifter = 0 - drifter_time = [] + drifters: list[Drifter] = [] argo = 0 argo_floats: list[ArgoFloat] = [] - ARGO_MIN_DEPTH = -config.argo_float_fieldset.U.depth[0] + argo_min_depth = -config.argo_float_fieldset.U.depth[0] + drifter_min_depth = -config.drifter_fieldset.U.depth[0] # run the model for the length of the sample_lons list for i in range(len(sample_lons) - 1): @@ -203,14 +204,23 @@ def SampleT(particle, fieldset, time): pset_CTD.time[0] + timedelta(minutes=20).total_seconds() ) # add CTD time and 20 minutes for deployment - # check if we are at a drifter deployment location + # check if we are at a `drifter` deployment location if drifter < len(config.drifter_deploylocations): while ( abs(sample_lons[i] - config.drifter_deploylocations[drifter][0]) < 0.01 and abs(sample_lats[i] - config.drifter_deploylocations[drifter][1]) < 0.01 ): - drifter_time.append(total_time) + drifters.append( + Drifter( + location=Location( + latitude=config.drifter_deploylocations[drifter][0], + longitude=config.drifter_deploylocations[drifter][1], + ), + deployment_time=total_time, + min_depth=drifter_min_depth, + ) + ) drifter += 1 print( f"Drifter {drifter} deployed at {sample_lons[i]}, {sample_lats[i]}" @@ -231,7 +241,7 @@ def SampleT(particle, fieldset, time): longitude=config.argo_deploylocations[argo][1], ), deployment_time=total_time, - min_depth=ARGO_MIN_DEPTH, + min_depth=argo_min_depth, max_depth=config.argo_characteristics["maxdepth"], drift_depth=config.argo_characteristics["driftdepth"], vertical_speed=config.argo_characteristics["vertical_speed"], @@ -267,10 +277,15 @@ def SampleT(particle, fieldset, time): ) print("Cruise has ended. Please wait for drifters and/or Argo floats to finish.") - # simulate drifter deployments - drifter_deployments(config, drifter_time) + # simulate drifters + simulate_drifters( + drifters=drifters, + fieldset=config.drifter_fieldset, + out_file_name=os.path.join("results", "drifters.zarr"), + outputdt=timedelta(minutes=5), + ) - # simulate argo deployments + # simulate argo floats simulate_argo_floats( argo_floats=argo_floats, fieldset=config.argo_float_fieldset, From 42cd06e8345370af9117c81544a0a69389610c3d Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 28 May 2024 17:30:40 +0200 Subject: [PATCH 05/40] codetools --- tests/instruments/test_drifters.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/tests/instruments/test_drifters.py b/tests/instruments/test_drifters.py index 0af64143..dafbb934 100644 --- a/tests/instruments/test_drifters.py +++ b/tests/instruments/test_drifters.py @@ -1,11 +1,13 @@ """Test the simulation of drifters.""" -from virtual_ship.instruments.drifter import simulate_drifters, Drifter -from virtual_ship.instruments import Location -from parcels import FieldSet -import numpy as np from datetime import timedelta +import numpy as np +from parcels import FieldSet + +from virtual_ship.instruments import Location +from virtual_ship.instruments.drifter import Drifter, simulate_drifters + def test_simulate_drifters() -> None: fieldset = FieldSet.from_data( From ffbc702111e4497ae31c1c643a086cb7821af3a3 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 17:25:44 +0200 Subject: [PATCH 06/40] outputdf wherever i forgot to use it --- virtual_ship/instruments/argo_float.py | 2 +- virtual_ship/instruments/drifter.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/virtual_ship/instruments/argo_float.py b/virtual_ship/instruments/argo_float.py index d9e7d50e..4918ddec 100644 --- a/virtual_ship/instruments/argo_float.py +++ b/virtual_ship/instruments/argo_float.py @@ -152,7 +152,7 @@ def simulate_argo_floats( # define output file for the simulation out_file = argo_float_particleset.ParticleFile( name=out_file_name, - outputdt=timedelta(minutes=5), + outputdt=outputdt, chunks=(1, 500), ) diff --git a/virtual_ship/instruments/drifter.py b/virtual_ship/instruments/drifter.py index c12215fa..6a603f62 100644 --- a/virtual_ship/instruments/drifter.py +++ b/virtual_ship/instruments/drifter.py @@ -65,7 +65,7 @@ def simulate_drifters( # define output file for the simulation out_file = drifter_particleset.ParticleFile( name=out_file_name, - outputdt=timedelta(minutes=5), + outputdt=outputdt, chunks=(1, 500), ) From 02a8b3c8b4e9a148077355eb3b69835e7c285de5 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 17:25:50 +0200 Subject: [PATCH 07/40] added todo list --- todo | 2 ++ 1 file changed, 2 insertions(+) create mode 100644 todo diff --git a/todo b/todo new file mode 100644 index 00000000..56ed7f74 --- /dev/null +++ b/todo @@ -0,0 +1,2 @@ +rename outputdt +end time for instruments \ No newline at end of file From f7917b34aca97d50e6220f57a4111f584fbb341f Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 17:26:01 +0200 Subject: [PATCH 08/40] ctd instrument --- virtual_ship/instruments/ctd.py | 109 ++++++++++++++++++++++++++++++++ 1 file changed, 109 insertions(+) create mode 100644 virtual_ship/instruments/ctd.py diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py new file mode 100644 index 00000000..a195f070 --- /dev/null +++ b/virtual_ship/instruments/ctd.py @@ -0,0 +1,109 @@ +"""CTD instrument.""" + +from dataclasses import dataclass +from datetime import timedelta + +import numpy as np +from parcels import FieldSet, JITParticle, ParticleSet, Variable + +from .location import Location + + +@dataclass +class CTDInstrument: + """Configuration for a single CTD instrument.""" + + location: Location + deployment_time: float + min_depth: float + + +class _CTDParticle(JITParticle): + salinity = Variable("salinity", initial=np.nan) + temperature = Variable("temperature", initial=np.nan) + raising = Variable("raising", dtype=np.int32, initial=0.0) + + +def _sample_temperature(particle, fieldset, time): + particle.temperature = fieldset.T[time, particle.depth, particle.lat, particle.lon] + + +def _sample_salinity(particle, fieldset, time): + particle.salinity = fieldset.S[time, particle.depth, particle.lat, particle.lon] + + +def _ctd_cast(particle, fieldset, time): + # Lowering and raising CTD + # TODO question: if is executed every time... move outside function? Not if "drifting" now possible + if ( + -fieldset.bathymetry[time, particle.depth, particle.lat, particle.lon] + > fieldset.max_depth + ): + maxdepth = ( + -fieldset.bathymetry[time, particle.depth, particle.lat, particle.lon] + 20 + ) + else: + maxdepth = fieldset.max_depth + winch_speed = -1.0 # sink and rise speed in m/s + + if particle.raising == 0: + # Sinking with winch_speed until near seafloor + particle_ddepth = winch_speed * particle.dt + if particle.depth <= maxdepth: + particle.raising = 1 + + if particle.raising == 1: + # Rising with winch_speed until depth is -2 m + if particle.depth < -2: + particle_ddepth = -winch_speed * particle.dt + if particle.depth + particle_ddepth >= -2: + # to break the loop ... + particle.state = 41 + print("CTD cast finished.") + + +def simulate_drifters( + ctd_instruments: list[CTDInstrument], + fieldset: FieldSet, + out_file_name: str, + outputdt: timedelta, +) -> None: + """ + Use parcels to simulate a set of CTD instruments in a fieldset. + + :param ctd_instruments: A list of CTD instruments to simulate. + :param fieldset: The fieldset to simulate the CTD instruments in. + :param out_file_name: The file to write the results to. + :param outputdt: Interval which dictates the update frequency of file output during simulation + """ + lon = [ctd.location.lon for ctd in ctd_instruments] + lat = [ctd.location.lat for ctd in ctd_instruments] + time = [ctd.deployment_time for ctd in ctd_instruments] + + # release CTD particle + ctd_particleset = ParticleSet( + fieldset=fieldset, + pclass=_CTDParticle, + lon=lon, + lat=lat, + depth=[ctd.min_depth for ctd in ctd_instruments], + time=time, + ) + + # define output file for the simulation + out_file = ctd_particleset.ParticleFile( + name=out_file_name, + outputdt=outputdt, + chunks=(1, 500), + ) + + # get time when the fieldset ends + fieldset_endtime = fieldset.time_origin.fulltime(fieldset.U.grid.time_full[-1]) + + # execute simulation + ctd_particleset.execute( + [_sample_salinity, _sample_temperature, _ctd_cast], + endtime=fieldset_endtime, + dt=outputdt, + output_file=out_file, + ) From 82093968af30d15e5eaef3ae61965ddef3f4eccf Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 17:35:43 +0200 Subject: [PATCH 09/40] CTDtest --- tests/instruments/test_ctd.py | 39 +++++++++++++++++++++++++++++++++ todo | 3 ++- virtual_ship/instruments/ctd.py | 9 +++++--- 3 files changed, 47 insertions(+), 4 deletions(-) create mode 100644 tests/instruments/test_ctd.py diff --git a/tests/instruments/test_ctd.py b/tests/instruments/test_ctd.py new file mode 100644 index 00000000..e37da6d3 --- /dev/null +++ b/tests/instruments/test_ctd.py @@ -0,0 +1,39 @@ +"""Test the simulation of CTD instruments.""" + +from datetime import timedelta + +import numpy as np +from parcels import FieldSet + +from virtual_ship.instruments import Location +from virtual_ship.instruments.ctd import CTDInstrument, simulate_ctd + + +def test_simulate_drifters() -> None: + fieldset = FieldSet.from_data( + {"U": 0, "V": 0, "T": 0, "S": 0, "bathymetry": 100}, + { + "lon": 0, + "lat": 0, + "time": [np.datetime64("1950-01-01") + np.timedelta64(632160, "h")], + }, + ) + + min_depth = -fieldset.U.depth[0] + max_depth = -fieldset.U.depth[-1] + + ctd_instruments = [ + CTDInstrument( + location=Location(latitude=0, longitude=0), + deployment_time=0, + min_depth=min_depth, + max_depth=max_depth, + ) + ] + + simulate_ctd( + ctd_instruments=ctd_instruments, + fieldset=fieldset, + out_file_name="test", + outputdt=timedelta(minutes=5), + ) diff --git a/todo b/todo index 56ed7f74..511fc884 100644 --- a/todo +++ b/todo @@ -1,2 +1,3 @@ rename outputdt -end time for instruments \ No newline at end of file +end time for instruments +check dtypes for particles \ No newline at end of file diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index a195f070..155a4951 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -16,12 +16,14 @@ class CTDInstrument: location: Location deployment_time: float min_depth: float + max_depth: float class _CTDParticle(JITParticle): salinity = Variable("salinity", initial=np.nan) temperature = Variable("temperature", initial=np.nan) raising = Variable("raising", dtype=np.int32, initial=0.0) + max_depth = Variable("max_depth", dtype=np.float32) def _sample_temperature(particle, fieldset, time): @@ -37,13 +39,13 @@ def _ctd_cast(particle, fieldset, time): # TODO question: if is executed every time... move outside function? Not if "drifting" now possible if ( -fieldset.bathymetry[time, particle.depth, particle.lat, particle.lon] - > fieldset.max_depth + > particle.max_depth ): maxdepth = ( -fieldset.bathymetry[time, particle.depth, particle.lat, particle.lon] + 20 ) else: - maxdepth = fieldset.max_depth + maxdepth = particle.max_depth winch_speed = -1.0 # sink and rise speed in m/s if particle.raising == 0: @@ -62,7 +64,7 @@ def _ctd_cast(particle, fieldset, time): print("CTD cast finished.") -def simulate_drifters( +def simulate_ctd( ctd_instruments: list[CTDInstrument], fieldset: FieldSet, out_file_name: str, @@ -88,6 +90,7 @@ def simulate_drifters( lat=lat, depth=[ctd.min_depth for ctd in ctd_instruments], time=time, + max_depth=[ctd.max_depth for ctd in ctd_instruments], ) # define output file for the simulation From 1b30b0ded3e5f92d782935c5494ea24a7202e88a Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 17:59:16 +0200 Subject: [PATCH 10/40] use new ctd instrument in sailship --- tests/instruments/test_ctd.py | 2 +- tests/sailship_config.json | 2 +- virtual_ship/sailship.py | 103 ++++++++++------------------------ 3 files changed, 32 insertions(+), 75 deletions(-) diff --git a/tests/instruments/test_ctd.py b/tests/instruments/test_ctd.py index e37da6d3..0c3784cd 100644 --- a/tests/instruments/test_ctd.py +++ b/tests/instruments/test_ctd.py @@ -35,5 +35,5 @@ def test_simulate_drifters() -> None: ctd_instruments=ctd_instruments, fieldset=fieldset, out_file_name="test", - outputdt=timedelta(minutes=5), + outputdt=timedelta(seconds=10), ) diff --git a/tests/sailship_config.json b/tests/sailship_config.json index d7b90b05..1f3805dc 100644 --- a/tests/sailship_config.json +++ b/tests/sailship_config.json @@ -19,7 +19,7 @@ "bin_size_m": 24 }, "CTD_locations": [ - [-23.071289, 63.743631] + [-23.081289, 63.743631] ], "CTD_settings": { "max_depth": "max" diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index ed0bc082..e41b52f7 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -11,6 +11,7 @@ from .costs import costs from .instruments.argo_float import ArgoFloat, simulate_argo_floats from .instruments.drifter import Drifter, simulate_drifters +from .instruments.ctd import CTDInstrument, simulate_ctd from .instruments.location import Location from .postprocess import postprocess from .virtual_ship_configuration import VirtualShipConfiguration @@ -18,7 +19,7 @@ def sailship(config: VirtualShipConfiguration): """ - Use parcels to simulate the ship, take CTDs and measure ADCP and underwaydata. + Use parcels to simulate the ship, take ctd_instruments and measure ADCP and underwaydata. :param config: The cruise configuration. :raises NotImplementedError: TODO @@ -43,50 +44,12 @@ class UnderwayDataParticle(JITParticle): salinity = Variable("salinity", initial=np.nan) temperature = Variable("temperature", initial=np.nan) - # Create CTD like particles to sample the ocean - class CTDParticle(JITParticle): - """Define a new particle class that does CTD like measurements.""" - - salinity = Variable("salinity", initial=np.nan) - temperature = Variable("temperature", initial=np.nan) - raising = Variable("raising", dtype=np.int32, initial=0.0) - # define ADCP sampling function without conversion (because of A grid) def SampleVel(particle, fieldset, time): particle.U, particle.V = fieldset.UV.eval( time, particle.depth, particle.lat, particle.lon, applyConversion=False ) - # define function lowering and raising CTD - def CTDcast(particle, fieldset, time): - # TODO question: if is executed every time... move outside function? Not if "drifting" now possible - if ( - -fieldset.bathymetry[time, particle.depth, particle.lat, particle.lon] - > fieldset.max_depth - ): - maxdepth = ( - -fieldset.bathymetry[time, particle.depth, particle.lat, particle.lon] - + 20 - ) - else: - maxdepth = fieldset.max_depth - winch_speed = -1.0 # sink and rise speed in m/s - - if particle.raising == 0: - # Sinking with winch_speed until near seafloor - particle_ddepth = winch_speed * particle.dt - if particle.depth <= maxdepth: - particle.raising = 1 - - if particle.raising == 1: - # Rising with winch_speed until depth is -2 m - if particle.depth < -2: - particle_ddepth = -winch_speed * particle.dt - if particle.depth + particle_ddepth >= -2: - # to break the loop ... - particle.state = 41 - print("CTD cast finished.") - # define function sampling Salinity def SampleS(particle, fieldset, time): particle.salinity = fieldset.S[time, particle.depth, particle.lat, particle.lon] @@ -128,19 +91,18 @@ def SampleT(particle, fieldset, time): name=os.path.join("results", "UnderwayData.zarr") ) - # initialize CTD station number and time + # initialize time total_time = timedelta(hours=0).total_seconds() - ctd = 0 - ctd_dt = timedelta( - seconds=10 - ) # timestep of CTD output reflecting post-process binning into 10m bins # initialize drifters and argo floats drifter = 0 drifters: list[Drifter] = [] argo = 0 argo_floats: list[ArgoFloat] = [] + ctd = 0 + ctd_instruments: list[CTDInstrument] = [] + ctd_min_depth = -config.ctd_fieldset.U.depth[0] argo_min_depth = -config.argo_float_fieldset.U.depth[0] drifter_min_depth = -config.drifter_fieldset.U.depth[0] @@ -171,38 +133,25 @@ def SampleT(particle, fieldset, time): # check if virtual ship is at a CTD station if ctd < len(config.CTD_locations): if ( - abs(sample_lons[i] - config.CTD_locations[ctd][0]) < 0.001 - and abs(sample_lats[i] - config.CTD_locations[ctd][1]) < 0.001 + abs(sample_lons[i] - config.CTD_locations[ctd][0]) < 0.01 + and abs(sample_lats[i] - config.CTD_locations[ctd][1]) < 0.01 ): - ctd += 1 - - # release CTD particle - pset_CTD = ParticleSet( - fieldset=fieldset, - pclass=CTDParticle, - lon=sample_lons[i], - lat=sample_lats[i], - depth=fieldset.mindepth, - time=total_time, - ) - - # create a ParticleFile to store the CTD output - ctd_output_file = pset_CTD.ParticleFile( - name=f"{os.path.join('results', 'CTDs', 'CTD_')}{ctd:03d}.zarr", - outputdt=ctd_dt, + ctd_instruments.append( + CTDInstrument( + location=Location( + latitude=config.CTD_locations[ctd][0], + longitude=config.CTD_locations[ctd][1], + ), + deployment_time=total_time, + min_depth=ctd_min_depth, + max_depth=-config.ctd_fieldset.U.depth[-1], + ) ) + ctd += 1 - # record the temperature and salinity of the particle - pset_CTD.execute( - [SampleS, SampleT, CTDcast], - runtime=timedelta(hours=8), - dt=ctd_dt, - output_file=ctd_output_file, - verbose_progress=False, - ) - total_time = ( - pset_CTD.time[0] + timedelta(minutes=20).total_seconds() - ) # add CTD time and 20 minutes for deployment + total_time += timedelta( + minutes=20 + ).total_seconds() # Add 20 minutes for deployment # check if we are at a `drifter` deployment location if drifter < len(config.drifter_deploylocations): @@ -277,6 +226,14 @@ def SampleT(particle, fieldset, time): ) print("Cruise has ended. Please wait for drifters and/or Argo floats to finish.") + # simulate ctd + simulate_ctd( + ctd_instruments=ctd_instruments, + fieldset=config.ctd_fieldset, + out_file_name=os.path.join("results", "ctd_instruments.zarr"), + outputdt=timedelta(seconds=10), + ) + # simulate drifters simulate_drifters( drifters=drifters, From af9536e6dc44ab017a416b55eb7e964aad5090af Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 17:59:52 +0200 Subject: [PATCH 11/40] fix typo --- tests/instruments/test_argo.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/instruments/test_argo.py b/tests/instruments/test_argo.py index 349ed051..83bbb261 100644 --- a/tests/instruments/test_argo.py +++ b/tests/instruments/test_argo.py @@ -12,7 +12,7 @@ def test_simulate_argo_floats() -> None: DRIFT_DEPTH = -1000 MAX_DEPTH = -2000 - VERTICLE_SPEED = -0.10 + VERTICAL_SPEED = -0.10 CYCLE_DAYS = 10 DRIFT_DAYS = 9 @@ -34,7 +34,7 @@ def test_simulate_argo_floats() -> None: min_depth=min_depth, max_depth=MAX_DEPTH, drift_depth=DRIFT_DEPTH, - vertical_speed=VERTICLE_SPEED, + vertical_speed=VERTICAL_SPEED, cycle_days=CYCLE_DAYS, drift_days=DRIFT_DAYS, ) From fe47b0c70df407747eb0a11a6776693aa9ea77cd Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 18:00:15 +0200 Subject: [PATCH 12/40] codetools --- virtual_ship/sailship.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index e41b52f7..6e9cb22d 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -10,8 +10,8 @@ from .costs import costs from .instruments.argo_float import ArgoFloat, simulate_argo_floats -from .instruments.drifter import Drifter, simulate_drifters from .instruments.ctd import CTDInstrument, simulate_ctd +from .instruments.drifter import Drifter, simulate_drifters from .instruments.location import Location from .postprocess import postprocess from .virtual_ship_configuration import VirtualShipConfiguration From 49564949001a4cb7f94707aacd80cf9536e31c3c Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 18:02:19 +0200 Subject: [PATCH 13/40] use new parcels api for particle class --- virtual_ship/instruments/ctd.py | 13 ++++++++----- virtual_ship/instruments/drifter.py | 10 +++++----- 2 files changed, 13 insertions(+), 10 deletions(-) diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index 155a4951..cc034c5a 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -19,11 +19,14 @@ class CTDInstrument: max_depth: float -class _CTDParticle(JITParticle): - salinity = Variable("salinity", initial=np.nan) - temperature = Variable("temperature", initial=np.nan) - raising = Variable("raising", dtype=np.int32, initial=0.0) - max_depth = Variable("max_depth", dtype=np.float32) +_CTDParticle = JITParticle.add_variables( + [ + Variable("salinity", initial=np.nan), + Variable("temperature", initial=np.nan), + Variable("raising", dtype=np.int32, initial=0.0), + Variable("max_depth", dtype=np.float32), + ] +) def _sample_temperature(particle, fieldset, time): diff --git a/virtual_ship/instruments/drifter.py b/virtual_ship/instruments/drifter.py index 6a603f62..a8f8c86f 100644 --- a/virtual_ship/instruments/drifter.py +++ b/virtual_ship/instruments/drifter.py @@ -18,11 +18,11 @@ class Drifter: min_depth: float -class _DrifterParticle(JITParticle): - temperature = Variable( - "temperature", - initial=np.nan, - ) +_DrifterParticle = JITParticle.add_variables( + [ + Variable("temperature", initial=np.nan), + ] +) def _sample_temperature(particle, fieldset, time): From e33c04597aa32d66c9119d3a3b3b5d42d5df878e Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 29 May 2024 18:09:06 +0200 Subject: [PATCH 14/40] add ctd to init --- virtual_ship/instruments/__init__.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/virtual_ship/instruments/__init__.py b/virtual_ship/instruments/__init__.py index 261bd7a8..f6d89a5a 100644 --- a/virtual_ship/instruments/__init__.py +++ b/virtual_ship/instruments/__init__.py @@ -1,6 +1,6 @@ """Measurement instrument that can be used with Parcels.""" -from . import argo_float, drifter +from . import argo_float, ctd, drifter from .location import Location -__all__ = ["Location", "argo_float", "drifter"] +__all__ = ["Location", "argo_float", "ctd", "drifter"] From 2fd016926ddf914bb1af1a51337985de71967107 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Sun, 2 Jun 2024 12:47:28 +0200 Subject: [PATCH 15/40] comments --- virtual_ship/instruments/argo_float.py | 2 +- virtual_ship/instruments/ctd.py | 2 +- virtual_ship/instruments/drifter.py | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/virtual_ship/instruments/argo_float.py b/virtual_ship/instruments/argo_float.py index 4918ddec..2b02a611 100644 --- a/virtual_ship/instruments/argo_float.py +++ b/virtual_ship/instruments/argo_float.py @@ -133,7 +133,7 @@ def simulate_argo_floats( lat = [argo.location.lat for argo in argo_floats] time = [argo.deployment_time for argo in argo_floats] - # define the parcels simulation + # define parcel particles argo_float_particleset = ParticleSet( fieldset=fieldset, pclass=_ArgoParticle, diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index cc034c5a..d5e14f6a 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -85,7 +85,7 @@ def simulate_ctd( lat = [ctd.location.lat for ctd in ctd_instruments] time = [ctd.deployment_time for ctd in ctd_instruments] - # release CTD particle + # define parcel particles ctd_particleset = ParticleSet( fieldset=fieldset, pclass=_CTDParticle, diff --git a/virtual_ship/instruments/drifter.py b/virtual_ship/instruments/drifter.py index a8f8c86f..7c4212a3 100644 --- a/virtual_ship/instruments/drifter.py +++ b/virtual_ship/instruments/drifter.py @@ -52,7 +52,7 @@ def simulate_drifters( lat = [drifter.location.lat for drifter in drifters] time = [drifter.deployment_time for drifter in drifters] - # define the parcels simulation + # define parcel particles drifter_particleset = ParticleSet( fieldset=fieldset, pclass=_DrifterParticle, From 98b17c7291f53f417802b80dbaf5005a0160b58a Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Sun, 2 Jun 2024 12:47:37 +0200 Subject: [PATCH 16/40] first try adcp --- virtual_ship/instruments/adcp.py | 65 ++++++++++++++++++++++++++++++++ 1 file changed, 65 insertions(+) create mode 100644 virtual_ship/instruments/adcp.py diff --git a/virtual_ship/instruments/adcp.py b/virtual_ship/instruments/adcp.py new file mode 100644 index 00000000..8c1c687e --- /dev/null +++ b/virtual_ship/instruments/adcp.py @@ -0,0 +1,65 @@ +"""ADCP.""" + +from dataclasses import dataclass +from datetime import timedelta + +import numpy as np +from parcels import FieldSet, JITParticle, ParticleSet, Variable + +from .location import Location + + +@dataclass +class SamplePoint: + location: Location + time: float + + +class _ADCPParticle(JITParticle): + """ADCP particle sampling a single point at some depth.""" + + U = Variable("U", dtype=np.float32, initial=0.0) + V = Variable("V", dtype=np.float32, initial=0.0) + + +def _sample_velocity(particle, fieldset, time): + particle.U, particle.V = fieldset.UV.eval( + time, particle.depth, particle.lat, particle.lon, applyConversion=False + ) + + +def simulate_adcp( + fieldset: FieldSet, + out_file_name: str, + max_depth: float, + min_depth: float, + bin_size: float, + sample_points: list[SamplePoint], +) -> None: + sample_points.sort(key=lambda p: p.time) + + bins = np.arange(max_depth, min_depth, bin_size) + num_particles = len(bins) + particleset = ParticleSet.from_list( + fieldset=fieldset, + pclass=_ADCPParticle, + lon=np.full( + num_particles, 0.0 + ), # initial lat/lon are irrelevant and will be overruled later. + lat=np.full(num_particles, 0.0), + depth=bins, + time=0, # same for time + ) + + # define output file for the simulation + out_file = particleset.ParticleFile( + name=out_file_name, + ) + + for point in sample_points: + particleset.lon_nextloop[:] = point.location.lon + particleset.lat_nextloop[:] = point.location.lat + particleset.time_nextloop[:] = point.time + + particleset.execute([_sample_velocity], dt=1, runtime=1, verbose_progress=False) + out_file.write(particleset, time=particleset[0].time) From f9df8bcde3caf192ffb23f5464b8a3149b3844fa Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Sun, 2 Jun 2024 12:58:15 +0200 Subject: [PATCH 17/40] adcp test --- tests/instruments/test_adcp.py | 35 ++++++++++++++++++++++++++++++++ virtual_ship/instruments/adcp.py | 3 +-- 2 files changed, 36 insertions(+), 2 deletions(-) create mode 100644 tests/instruments/test_adcp.py diff --git a/tests/instruments/test_adcp.py b/tests/instruments/test_adcp.py new file mode 100644 index 00000000..9360164f --- /dev/null +++ b/tests/instruments/test_adcp.py @@ -0,0 +1,35 @@ +"""Test the simulation of Argo floats.""" + +from datetime import timedelta + +import numpy as np +from parcels import FieldSet + +from virtual_ship.instruments import Location +from virtual_ship.instruments.adcp import simulate_adcp, SamplePoint + + +def test_simulate_argo_floats() -> None: + MAX_DEPTH = -1000 + MIN_DEPTH = -5 + BIN_SIZE = 24 + + fieldset = FieldSet.from_data( + {"U": 0, "V": 0}, + { + "lon": 0, + "lat": 0, + "time": [np.datetime64("1950-01-01") + np.timedelta64(632160, "h")], + }, + ) + + sample_points = [SamplePoint(Location(0, 0), 0)] + + simulate_adcp( + fieldset=fieldset, + out_file_name="test", + max_depth=MAX_DEPTH, + min_depth=MIN_DEPTH, + bin_size=BIN_SIZE, + sample_points=sample_points, + ) diff --git a/virtual_ship/instruments/adcp.py b/virtual_ship/instruments/adcp.py index 8c1c687e..e7fb9bee 100644 --- a/virtual_ship/instruments/adcp.py +++ b/virtual_ship/instruments/adcp.py @@ -1,7 +1,6 @@ -"""ADCP.""" +"""ADCP instrument.""" from dataclasses import dataclass -from datetime import timedelta import numpy as np from parcels import FieldSet, JITParticle, ParticleSet, Variable From 9ded5e244aaf7ef59489e20ec44a2cebb1225b13 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Mon, 3 Jun 2024 15:37:14 +0200 Subject: [PATCH 18/40] use adcp instrument in sailship --- virtual_ship/sailship.py | 210 +++++++++++++++++---------------------- 1 file changed, 92 insertions(+), 118 deletions(-) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 6e9cb22d..2ef92060 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -5,13 +5,14 @@ import numpy as np import pyproj -from parcels import Field, FieldSet, JITParticle, ParticleSet, Variable +from parcels import JITParticle, ParticleSet, Variable from shapely.geometry import Point, Polygon from .costs import costs from .instruments.argo_float import ArgoFloat, simulate_argo_floats from .instruments.ctd import CTDInstrument, simulate_ctd from .instruments.drifter import Drifter, simulate_drifters +from .instruments.adcp import simulate_adcp, SamplePoint as ADCPSamplePoint from .instruments.location import Location from .postprocess import postprocess from .virtual_ship_configuration import VirtualShipConfiguration @@ -25,18 +26,12 @@ def sailship(config: VirtualShipConfiguration): :raises NotImplementedError: TODO """ # Create fieldset and retreive final schip route as sample_lons and sample_lats - fieldset = config.ctd_fieldset # create_fieldset(config, data_dir) + fieldset = config.ctd_fieldset + adcp_fieldset = config.ctd_fieldset sample_lons, sample_lats = shiproute(config) print("Arrived in region of interest, starting to gather data.") - # Create Vessel Mounted ADCP like particles to sample the ocean - class VM_ADCPParticle(JITParticle): - """Define a new particle class that does Vessel Mounted ADCP like measurements.""" - - U = Variable("U", dtype=np.float32, initial=0.0) - V = Variable("V", dtype=np.float32, initial=0.0) - # Create particle to sample water underway class UnderwayDataParticle(JITParticle): """Define a new particle class that samples water directly under the hull.""" @@ -44,12 +39,6 @@ class UnderwayDataParticle(JITParticle): salinity = Variable("salinity", initial=np.nan) temperature = Variable("temperature", initial=np.nan) - # define ADCP sampling function without conversion (because of A grid) - def SampleVel(particle, fieldset, time): - particle.U, particle.V = fieldset.UV.eval( - time, particle.depth, particle.lat, particle.lon, applyConversion=False - ) - # define function sampling Salinity def SampleS(particle, fieldset, time): particle.salinity = fieldset.S[time, particle.depth, particle.lat, particle.lon] @@ -60,24 +49,6 @@ def SampleT(particle, fieldset, time): time, particle.depth, particle.lat, particle.lon ] - # Create ADCP like particleset and output file - ADCP_bins = np.arange( - config.ADCP_settings["max_depth"], -5, config.ADCP_settings["bin_size_m"] - ) - vert_particles = len(ADCP_bins) - pset_ADCP = ParticleSet.from_list( - fieldset=fieldset, - pclass=VM_ADCPParticle, - lon=np.full(vert_particles, sample_lons[0]), - lat=np.full(vert_particles, sample_lats[0]), - depth=ADCP_bins, - time=0, - ) - adcp_output_file = pset_ADCP.ParticleFile(name=os.path.join("results", "ADCP.zarr")) - adcp_dt = timedelta( - minutes=5 - ).total_seconds() # timestep of ADCP output, every 5 min - # Create underway particle pset_UnderwayData = ParticleSet.from_list( fieldset=fieldset, @@ -102,6 +73,12 @@ def SampleT(particle, fieldset, time): ctd = 0 ctd_instruments: list[CTDInstrument] = [] + adcp_dt = timedelta(minutes=5).total_seconds() + adcp_sample_points = [ + ADCPSamplePoint(Location(latitude=lat, longitude=lon), n * adcp_dt) + for n, (lat, lon) in enumerate(zip(sample_lats, sample_lons)) + ] + ctd_min_depth = -config.ctd_fieldset.U.depth[0] argo_min_depth = -config.argo_float_fieldset.U.depth[0] drifter_min_depth = -config.drifter_fieldset.U.depth[0] @@ -112,23 +89,13 @@ def SampleT(particle, fieldset, time): if i % 96 == 0: print(f"Gathered data {timedelta(seconds=total_time)} hours since start.") - # execute the ADCP kernels to sample U and V and underway T and S - if config.ADCP_data: - pset_ADCP.execute( - [SampleVel], dt=adcp_dt, runtime=1, verbose_progress=False - ) - adcp_output_file.write(pset_ADCP, time=pset_ADCP[0].time) if config.underway_data: pset_UnderwayData.execute( [SampleS, SampleT], dt=adcp_dt, runtime=1, verbose_progress=False ) - UnderwayData_output_file.write(pset_UnderwayData, time=pset_ADCP[0].time) - if pset_ADCP[0].time > fieldset.maxtime: - print( - "Ship time is over, waiting for drifters and/or Argo floats to finish." + UnderwayData_output_file.write( + pset_UnderwayData, time=pset_UnderwayData[0].time ) - raise NotImplementedError() - # return drifter_time, argo_time # check if virtual ship is at a CTD station if ctd < len(config.CTD_locations): @@ -204,19 +171,16 @@ def SampleT(particle, fieldset, time): break # update the particle time and location - pset_ADCP.lon_nextloop[:] = sample_lons[i + 1] - pset_ADCP.lat_nextloop[:] = sample_lats[i + 1] pset_UnderwayData.lon_nextloop[:] = sample_lons[i + 1] pset_UnderwayData.lat_nextloop[:] = sample_lats[i + 1] total_time += adcp_dt - pset_ADCP.time_nextloop[:] = total_time pset_UnderwayData.time_nextloop[:] = total_time # write the final locations of the ADCP and Underway data particles - if config.ADCP_data: - pset_ADCP.execute(SampleVel, dt=adcp_dt, runtime=1, verbose_progress=False) - adcp_output_file.write_latest_locations(pset_ADCP, time=total_time) + # if config.ADCP_data: # TODO what is this + # pset_ADCP.execute(SampleVel, dt=adcp_dt, runtime=1, verbose_progress=False) + # adcp_output_file.write_latest_locations(pset_ADCP, time=total_time) if config.underway_data: pset_UnderwayData.execute( [SampleS, SampleT], dt=adcp_dt, runtime=1, verbose_progress=False @@ -226,6 +190,16 @@ def SampleT(particle, fieldset, time): ) print("Cruise has ended. Please wait for drifters and/or Argo floats to finish.") + # simulate adcp + simulate_adcp( + fieldset=adcp_fieldset, + out_file_name=os.path.join("results", "adcp.zarr"), + max_depth=config.ADCP_settings["max_depth"], + min_depth=-5, + bin_size=config.ADCP_settings["bin_size_m"], + sample_points=adcp_sample_points, + ) + # simulate ctd simulate_ctd( ctd_instruments=ctd_instruments, @@ -261,73 +235,73 @@ def SampleT(particle, fieldset, time): ) -def create_fieldset(config, data_dir: str): - """ - Create fieldset from netcdf files and adds bathymetry data for CTD cast, returns fieldset with negative depth values. - - :param config: The cruise configuration. - :param data_dir: TODO - :returns: The fieldset. - :raises ValueError: If downloaded data is not as expected. - """ - filenames = { - "U": os.path.join(data_dir, "studentdata_UV.nc"), - "V": os.path.join(data_dir, "studentdata_UV.nc"), - "S": os.path.join(data_dir, "studentdata_S.nc"), - "T": os.path.join(data_dir, "studentdata_T.nc"), - } - variables = {"U": "uo", "V": "vo", "S": "so", "T": "thetao"} - dimensions = { - "lon": "longitude", - "lat": "latitude", - "time": "time", - "depth": "depth", - } - - # create the fieldset and set interpolation methods - fieldset = FieldSet.from_netcdf( - filenames, variables, dimensions, allow_time_extrapolation=True - ) - fieldset.T.interp_method = "linear_invdist_land_tracer" - fieldset.S.interp_method = "linear_invdist_land_tracer" - for g in fieldset.gridset.grids: - if max(g.depth) > 0: - g.depth = -g.depth # make depth negative - fieldset.mindepth = -fieldset.U.depth[0] # uppermost layer in the hydrodynamic data - if config.CTD_settings["max_depth"] == "max": - fieldset.add_constant("max_depth", -fieldset.U.depth[-1]) - else: - fieldset.add_constant("max_depth", config.CTD_settings["max_depth"]) - fieldset.add_constant("maxtime", fieldset.U.grid.time_full[-1]) - - # add bathymetry data to the fieldset for CTD cast - bathymetry_file = os.path.join(data_dir, "GLO-MFC_001_024_mask_bathy.nc") - bathymetry_variables = ("bathymetry", "deptho") - bathymetry_dimensions = {"lon": "longitude", "lat": "latitude"} - bathymetry_field = Field.from_netcdf( - bathymetry_file, bathymetry_variables, bathymetry_dimensions - ) - fieldset.add_field(bathymetry_field) - # read in data already - fieldset.computeTimeChunk(0, 1) - - if fieldset.U.lon.min() > config.region_of_interest["West"]: - raise ValueError( - "FieldSet western boundary is outside region of interest. Please run download_data.py again." - ) - if fieldset.U.lon.max() < config.region_of_interest["East"]: - raise ValueError( - "FieldSet eastern boundary is outside region of interest. Please run download_data.py again." - ) - if fieldset.U.lat.min() > config.region_of_interest["South"]: - raise ValueError( - "FieldSet southern boundary is outside region of interest. Please run download_data.py again." - ) - if fieldset.U.lat.max() < config.region_of_interest["North"]: - raise ValueError( - "FieldSet northern boundary is outside region of interest. Please run download_data.py again." - ) - return fieldset +# def create_fieldset(config, data_dir: str): +# """ +# Create fieldset from netcdf files and adds bathymetry data for CTD cast, returns fieldset with negative depth values. + +# :param config: The cruise configuration. +# :param data_dir: TODO +# :returns: The fieldset. +# :raises ValueError: If downloaded data is not as expected. +# """ +# filenames = { +# "U": os.path.join(data_dir, "studentdata_UV.nc"), +# "V": os.path.join(data_dir, "studentdata_UV.nc"), +# "S": os.path.join(data_dir, "studentdata_S.nc"), +# "T": os.path.join(data_dir, "studentdata_T.nc"), +# } +# variables = {"U": "uo", "V": "vo", "S": "so", "T": "thetao"} +# dimensions = { +# "lon": "longitude", +# "lat": "latitude", +# "time": "time", +# "depth": "depth", +# } + +# # create the fieldset and set interpolation methods +# fieldset = FieldSet.from_netcdf( +# filenames, variables, dimensions, allow_time_extrapolation=True +# ) +# fieldset.T.interp_method = "linear_invdist_land_tracer" +# fieldset.S.interp_method = "linear_invdist_land_tracer" +# for g in fieldset.gridset.grids: +# if max(g.depth) > 0: +# g.depth = -g.depth # make depth negative +# fieldset.mindepth = -fieldset.U.depth[0] # uppermost layer in the hydrodynamic data +# if config.CTD_settings["max_depth"] == "max": +# fieldset.add_constant("max_depth", -fieldset.U.depth[-1]) +# else: +# fieldset.add_constant("max_depth", config.CTD_settings["max_depth"]) +# fieldset.add_constant("maxtime", fieldset.U.grid.time_full[-1]) + +# # add bathymetry data to the fieldset for CTD cast +# bathymetry_file = os.path.join(data_dir, "GLO-MFC_001_024_mask_bathy.nc") +# bathymetry_variables = ("bathymetry", "deptho") +# bathymetry_dimensions = {"lon": "longitude", "lat": "latitude"} +# bathymetry_field = Field.from_netcdf( +# bathymetry_file, bathymetry_variables, bathymetry_dimensions +# ) +# fieldset.add_field(bathymetry_field) +# # read in data already +# fieldset.computeTimeChunk(0, 1) + +# if fieldset.U.lon.min() > config.region_of_interest["West"]: +# raise ValueError( +# "FieldSet western boundary is outside region of interest. Please run download_data.py again." +# ) +# if fieldset.U.lon.max() < config.region_of_interest["East"]: +# raise ValueError( +# "FieldSet eastern boundary is outside region of interest. Please run download_data.py again." +# ) +# if fieldset.U.lat.min() > config.region_of_interest["South"]: +# raise ValueError( +# "FieldSet southern boundary is outside region of interest. Please run download_data.py again." +# ) +# if fieldset.U.lat.max() < config.region_of_interest["North"]: +# raise ValueError( +# "FieldSet northern boundary is outside region of interest. Please run download_data.py again." +# ) +# return fieldset def shiproute(config): From 4e571173832c3ae047c49a541f42a12735bef320 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Mon, 3 Jun 2024 15:51:23 +0200 Subject: [PATCH 19/40] cleanup --- tests/instruments/test_adcp.py | 2 -- virtual_ship/instruments/adcp.py | 11 ++++++----- virtual_ship/instruments/argo_float.py | 4 ++-- virtual_ship/instruments/ctd.py | 4 ++-- virtual_ship/instruments/drifter.py | 2 +- 5 files changed, 11 insertions(+), 12 deletions(-) diff --git a/tests/instruments/test_adcp.py b/tests/instruments/test_adcp.py index 9360164f..84540c2e 100644 --- a/tests/instruments/test_adcp.py +++ b/tests/instruments/test_adcp.py @@ -1,7 +1,5 @@ """Test the simulation of Argo floats.""" -from datetime import timedelta - import numpy as np from parcels import FieldSet diff --git a/virtual_ship/instruments/adcp.py b/virtual_ship/instruments/adcp.py index e7fb9bee..7868aa7b 100644 --- a/virtual_ship/instruments/adcp.py +++ b/virtual_ship/instruments/adcp.py @@ -14,11 +14,12 @@ class SamplePoint: time: float -class _ADCPParticle(JITParticle): - """ADCP particle sampling a single point at some depth.""" - - U = Variable("U", dtype=np.float32, initial=0.0) - V = Variable("V", dtype=np.float32, initial=0.0) +_ADCPParticle = JITParticle.add_variables( + [ + Variable("U", dtype=np.float32, initial=np.nan), + Variable("V", dtype=np.float32, initial=np.nan), + ] +) def _sample_velocity(particle, fieldset, time): diff --git a/virtual_ship/instruments/argo_float.py b/virtual_ship/instruments/argo_float.py index 2b02a611..1218d965 100644 --- a/virtual_ship/instruments/argo_float.py +++ b/virtual_ship/instruments/argo_float.py @@ -36,8 +36,8 @@ class ArgoFloat: Variable("cycle_phase", dtype=np.int32, initial=0.0), Variable("cycle_age", dtype=np.float32, initial=0.0), Variable("drift_age", dtype=np.float32, initial=0.0), - Variable("salinity", initial=np.nan), - Variable("temperature", initial=np.nan), + Variable("salinity", dtype=np.float32, initial=np.nan), + Variable("temperature", dtype=np.float32, initial=np.nan), Variable("min_depth", dtype=np.float32), Variable("max_depth", dtype=np.float32), Variable("drift_depth", dtype=np.float32), diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index d5e14f6a..53b91c8a 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -21,8 +21,8 @@ class CTDInstrument: _CTDParticle = JITParticle.add_variables( [ - Variable("salinity", initial=np.nan), - Variable("temperature", initial=np.nan), + Variable("salinity", dtype=np.float32, initial=np.nan), + Variable("temperature", dtype=np.float32, initial=np.nan), Variable("raising", dtype=np.int32, initial=0.0), Variable("max_depth", dtype=np.float32), ] diff --git a/virtual_ship/instruments/drifter.py b/virtual_ship/instruments/drifter.py index 7c4212a3..fd38782c 100644 --- a/virtual_ship/instruments/drifter.py +++ b/virtual_ship/instruments/drifter.py @@ -20,7 +20,7 @@ class Drifter: _DrifterParticle = JITParticle.add_variables( [ - Variable("temperature", initial=np.nan), + Variable("temperature", dtype=np.float32, initial=np.nan), ] ) From 0eaaff1aeb52649498b061fcebb717d0b0506a16 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Mon, 3 Jun 2024 15:51:41 +0200 Subject: [PATCH 20/40] instrument ship s t --- tests/instruments/test_ship_st.py | 29 ++++++++++++ virtual_ship/instruments/ship_st.py | 68 +++++++++++++++++++++++++++++ 2 files changed, 97 insertions(+) create mode 100644 tests/instruments/test_ship_st.py create mode 100644 virtual_ship/instruments/ship_st.py diff --git a/tests/instruments/test_ship_st.py b/tests/instruments/test_ship_st.py new file mode 100644 index 00000000..f4871859 --- /dev/null +++ b/tests/instruments/test_ship_st.py @@ -0,0 +1,29 @@ +"""Test the simulation of ship salinity temperature measurements.""" + +import numpy as np +from parcels import FieldSet + +from virtual_ship.instruments import Location +from virtual_ship.instruments.ship_st import simulate_ship_st, SamplePoint + + +def test_simulate_argo_floats() -> None: + DEPTH = -2 + + fieldset = FieldSet.from_data( + {"U": 0, "V": 0, "S": 0, "T": 0}, + { + "lon": 0, + "lat": 0, + "time": [np.datetime64("1950-01-01") + np.timedelta64(632160, "h")], + }, + ) + + sample_points = [SamplePoint(Location(0, 0), 0)] + + simulate_ship_st( + fieldset=fieldset, + out_file_name="test", + depth=DEPTH, + sample_points=sample_points, + ) diff --git a/virtual_ship/instruments/ship_st.py b/virtual_ship/instruments/ship_st.py new file mode 100644 index 00000000..80004799 --- /dev/null +++ b/virtual_ship/instruments/ship_st.py @@ -0,0 +1,68 @@ +"""Ship salinity and temperature.""" + +from dataclasses import dataclass + +import numpy as np +from parcels import FieldSet, JITParticle, ParticleSet, Variable + +from .location import Location + + +@dataclass +class SamplePoint: + location: Location + time: float + + +_ShipSTParticle = JITParticle.add_variables( + [ + Variable("salinity", dtype=np.float32, initial=np.nan), + Variable("temperature", dtype=np.float32, initial=np.nan), + ] +) + + +# define function sampling Salinity +def _sample_salinity(particle, fieldset, time): + particle.salinity = fieldset.S[time, particle.depth, particle.lat, particle.lon] + + +# define function sampling Temperature +def _sample_temperature(particle, fieldset, time): + particle.temperature = fieldset.T[time, particle.depth, particle.lat, particle.lon] + + +def simulate_ship_st( + fieldset: FieldSet, + out_file_name: str, + depth: float, + sample_points: list[SamplePoint], +) -> None: + sample_points.sort(key=lambda p: p.time) + + particleset = ParticleSet.from_list( + fieldset=fieldset, + pclass=_ShipSTParticle, + lon=0.0, # initial lat/lon are irrelevant and will be overruled later. + lat=0.0, + depth=depth, + time=0, # same for time + ) + + # define output file for the simulation + out_file = particleset.ParticleFile( + name=out_file_name, + ) + + for point in sample_points: + particleset.lon_nextloop[:] = point.location.lon + particleset.lat_nextloop[:] = point.location.lat + particleset.time_nextloop[:] = point.time + + particleset.execute( + [_sample_salinity, _sample_temperature], + dt=1, + runtime=1, + verbose_progress=False, + ) + out_file.write(particleset, time=particleset[0].time) From 0ae077601ad84ec256e89bb29c6f0db9adf6d3cf Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Mon, 3 Jun 2024 16:05:12 +0200 Subject: [PATCH 21/40] use ship st in sailship --- virtual_ship/sailship.py | 93 +++++++++++----------------------------- 1 file changed, 25 insertions(+), 68 deletions(-) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 2ef92060..79ede8b4 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -5,7 +5,6 @@ import numpy as np import pyproj -from parcels import JITParticle, ParticleSet, Variable from shapely.geometry import Point, Polygon from .costs import costs @@ -13,6 +12,7 @@ from .instruments.ctd import CTDInstrument, simulate_ctd from .instruments.drifter import Drifter, simulate_drifters from .instruments.adcp import simulate_adcp, SamplePoint as ADCPSamplePoint +from .instruments.ship_st import simulate_ship_st, SamplePoint as ShipSTSamplePoint from .instruments.location import Location from .postprocess import postprocess from .virtual_ship_configuration import VirtualShipConfiguration @@ -23,48 +23,14 @@ def sailship(config: VirtualShipConfiguration): Use parcels to simulate the ship, take ctd_instruments and measure ADCP and underwaydata. :param config: The cruise configuration. - :raises NotImplementedError: TODO """ # Create fieldset and retreive final schip route as sample_lons and sample_lats - fieldset = config.ctd_fieldset adcp_fieldset = config.ctd_fieldset + ship_st_fieldset = config.ctd_fieldset sample_lons, sample_lats = shiproute(config) print("Arrived in region of interest, starting to gather data.") - # Create particle to sample water underway - class UnderwayDataParticle(JITParticle): - """Define a new particle class that samples water directly under the hull.""" - - salinity = Variable("salinity", initial=np.nan) - temperature = Variable("temperature", initial=np.nan) - - # define function sampling Salinity - def SampleS(particle, fieldset, time): - particle.salinity = fieldset.S[time, particle.depth, particle.lat, particle.lon] - - # define function sampling Temperature - def SampleT(particle, fieldset, time): - particle.temperature = fieldset.T[ - time, particle.depth, particle.lat, particle.lon - ] - - # Create underway particle - pset_UnderwayData = ParticleSet.from_list( - fieldset=fieldset, - pclass=UnderwayDataParticle, - lon=sample_lons[0], - lat=sample_lats[0], - depth=-2, - time=0, - ) - UnderwayData_output_file = pset_UnderwayData.ParticleFile( - name=os.path.join("results", "UnderwayData.zarr") - ) - - # initialize time - total_time = timedelta(hours=0).total_seconds() - # initialize drifters and argo floats drifter = 0 drifters: list[Drifter] = [] @@ -73,9 +39,13 @@ def SampleT(particle, fieldset, time): ctd = 0 ctd_instruments: list[CTDInstrument] = [] - adcp_dt = timedelta(minutes=5).total_seconds() + route_points_dt = timedelta(minutes=5).total_seconds() adcp_sample_points = [ - ADCPSamplePoint(Location(latitude=lat, longitude=lon), n * adcp_dt) + ADCPSamplePoint(Location(latitude=lat, longitude=lon), n * route_points_dt) + for n, (lat, lon) in enumerate(zip(sample_lats, sample_lons)) + ] + ship_st_sample_points = [ + ShipSTSamplePoint(Location(latitude=lat, longitude=lon), n * route_points_dt) for n, (lat, lon) in enumerate(zip(sample_lats, sample_lons)) ] @@ -84,19 +54,12 @@ def SampleT(particle, fieldset, time): drifter_min_depth = -config.drifter_fieldset.U.depth[0] # run the model for the length of the sample_lons list + total_time = timedelta(hours=0).total_seconds() for i in range(len(sample_lons) - 1): if i % 96 == 0: print(f"Gathered data {timedelta(seconds=total_time)} hours since start.") - if config.underway_data: - pset_UnderwayData.execute( - [SampleS, SampleT], dt=adcp_dt, runtime=1, verbose_progress=False - ) - UnderwayData_output_file.write( - pset_UnderwayData, time=pset_UnderwayData[0].time - ) - # check if virtual ship is at a CTD station if ctd < len(config.CTD_locations): if ( @@ -170,27 +133,20 @@ def SampleT(particle, fieldset, time): if argo == len(config.argo_deploylocations): break - # update the particle time and location - pset_UnderwayData.lon_nextloop[:] = sample_lons[i + 1] - pset_UnderwayData.lat_nextloop[:] = sample_lats[i + 1] - - total_time += adcp_dt - pset_UnderwayData.time_nextloop[:] = total_time - - # write the final locations of the ADCP and Underway data particles - # if config.ADCP_data: # TODO what is this - # pset_ADCP.execute(SampleVel, dt=adcp_dt, runtime=1, verbose_progress=False) - # adcp_output_file.write_latest_locations(pset_ADCP, time=total_time) - if config.underway_data: - pset_UnderwayData.execute( - [SampleS, SampleT], dt=adcp_dt, runtime=1, verbose_progress=False - ) - UnderwayData_output_file.write_latest_locations( - pset_UnderwayData, time=total_time - ) + # update timey + total_time += route_points_dt + print("Cruise has ended. Please wait for drifters and/or Argo floats to finish.") - # simulate adcp + print("Simulating onboard salinity and temperature measurements.") + simulate_ship_st( + fieldset=ship_st_fieldset, + out_file_name=os.path.join("results", "ship_st.zarr"), + depth=-2, + sample_points=ship_st_sample_points, + ) + + print("Simulating onboard ADCP.") simulate_adcp( fieldset=adcp_fieldset, out_file_name=os.path.join("results", "adcp.zarr"), @@ -200,7 +156,7 @@ def SampleT(particle, fieldset, time): sample_points=adcp_sample_points, ) - # simulate ctd + print("Simulating CTD casts.") simulate_ctd( ctd_instruments=ctd_instruments, fieldset=config.ctd_fieldset, @@ -208,7 +164,7 @@ def SampleT(particle, fieldset, time): outputdt=timedelta(seconds=10), ) - # simulate drifters + print("Simulating drifters") simulate_drifters( drifters=drifters, fieldset=config.drifter_fieldset, @@ -216,7 +172,7 @@ def SampleT(particle, fieldset, time): outputdt=timedelta(minutes=5), ) - # simulate argo floats + print("Simulating argo floats") simulate_argo_floats( argo_floats=argo_floats, fieldset=config.argo_float_fieldset, @@ -225,6 +181,7 @@ def SampleT(particle, fieldset, time): ) # convert CTD data to CSV + print("Postprocessing..") postprocess() print("All data has been gathered and postprocessed, returning home.") From 9d6b0c59c686ccc0559f6b6dd2812a0c22f82a5a Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Mon, 3 Jun 2024 16:17:41 +0200 Subject: [PATCH 22/40] cleanup --- tests/instruments/test_ctd.py | 8 ++++---- virtual_ship/instruments/ctd.py | 22 +++++++++++----------- virtual_ship/sailship.py | 12 ++++++------ 3 files changed, 21 insertions(+), 21 deletions(-) diff --git a/tests/instruments/test_ctd.py b/tests/instruments/test_ctd.py index 0c3784cd..ba6aa81a 100644 --- a/tests/instruments/test_ctd.py +++ b/tests/instruments/test_ctd.py @@ -6,7 +6,7 @@ from parcels import FieldSet from virtual_ship.instruments import Location -from virtual_ship.instruments.ctd import CTDInstrument, simulate_ctd +from virtual_ship.instruments.ctd import CTD, simulate_ctd def test_simulate_drifters() -> None: @@ -22,8 +22,8 @@ def test_simulate_drifters() -> None: min_depth = -fieldset.U.depth[0] max_depth = -fieldset.U.depth[-1] - ctd_instruments = [ - CTDInstrument( + ctds = [ + CTD( location=Location(latitude=0, longitude=0), deployment_time=0, min_depth=min_depth, @@ -32,7 +32,7 @@ def test_simulate_drifters() -> None: ] simulate_ctd( - ctd_instruments=ctd_instruments, + ctds=ctds, fieldset=fieldset, out_file_name="test", outputdt=timedelta(seconds=10), diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index 53b91c8a..4271d63c 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -10,8 +10,8 @@ @dataclass -class CTDInstrument: - """Configuration for a single CTD instrument.""" +class CTD: + """Configuration for a single CTD.""" location: Location deployment_time: float @@ -68,22 +68,22 @@ def _ctd_cast(particle, fieldset, time): def simulate_ctd( - ctd_instruments: list[CTDInstrument], + ctds: list[CTD], fieldset: FieldSet, out_file_name: str, outputdt: timedelta, ) -> None: """ - Use parcels to simulate a set of CTD instruments in a fieldset. + Use parcels to simulate a set of CTDs in a fieldset. - :param ctd_instruments: A list of CTD instruments to simulate. - :param fieldset: The fieldset to simulate the CTD instruments in. + :param ctd: A list of CTDs to simulate. + :param fieldset: The fieldset to simulate the CTDs in. :param out_file_name: The file to write the results to. :param outputdt: Interval which dictates the update frequency of file output during simulation """ - lon = [ctd.location.lon for ctd in ctd_instruments] - lat = [ctd.location.lat for ctd in ctd_instruments] - time = [ctd.deployment_time for ctd in ctd_instruments] + lon = [ctd.location.lon for ctd in ctds] + lat = [ctd.location.lat for ctd in ctds] + time = [ctd.deployment_time for ctd in ctds] # define parcel particles ctd_particleset = ParticleSet( @@ -91,9 +91,9 @@ def simulate_ctd( pclass=_CTDParticle, lon=lon, lat=lat, - depth=[ctd.min_depth for ctd in ctd_instruments], + depth=[ctd.min_depth for ctd in ctds], time=time, - max_depth=[ctd.max_depth for ctd in ctd_instruments], + max_depth=[ctd.max_depth for ctd in ctds], ) # define output file for the simulation diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 79ede8b4..7178c9f7 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -9,7 +9,7 @@ from .costs import costs from .instruments.argo_float import ArgoFloat, simulate_argo_floats -from .instruments.ctd import CTDInstrument, simulate_ctd +from .instruments.ctd import CTD, simulate_ctd from .instruments.drifter import Drifter, simulate_drifters from .instruments.adcp import simulate_adcp, SamplePoint as ADCPSamplePoint from .instruments.ship_st import simulate_ship_st, SamplePoint as ShipSTSamplePoint @@ -37,7 +37,7 @@ def sailship(config: VirtualShipConfiguration): argo = 0 argo_floats: list[ArgoFloat] = [] ctd = 0 - ctd_instruments: list[CTDInstrument] = [] + ctds: list[CTD] = [] route_points_dt = timedelta(minutes=5).total_seconds() adcp_sample_points = [ @@ -66,8 +66,8 @@ def sailship(config: VirtualShipConfiguration): abs(sample_lons[i] - config.CTD_locations[ctd][0]) < 0.01 and abs(sample_lats[i] - config.CTD_locations[ctd][1]) < 0.01 ): - ctd_instruments.append( - CTDInstrument( + ctds.append( + CTD( location=Location( latitude=config.CTD_locations[ctd][0], longitude=config.CTD_locations[ctd][1], @@ -158,9 +158,9 @@ def sailship(config: VirtualShipConfiguration): print("Simulating CTD casts.") simulate_ctd( - ctd_instruments=ctd_instruments, + ctds=ctds, fieldset=config.ctd_fieldset, - out_file_name=os.path.join("results", "ctd_instruments.zarr"), + out_file_name=os.path.join("results", "ctd.zarr"), outputdt=timedelta(seconds=10), ) From b06ce7dd478596cbb3f79db4c909db752514f355 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Mon, 3 Jun 2024 17:24:55 +0200 Subject: [PATCH 23/40] refactor sailship ctd cast --- virtual_ship/sailship.py | 72 +++++++++++++++++++++++++--------------- 1 file changed, 45 insertions(+), 27 deletions(-) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 7178c9f7..4752b828 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -31,12 +31,19 @@ def sailship(config: VirtualShipConfiguration): sample_lons, sample_lats = shiproute(config) print("Arrived in region of interest, starting to gather data.") - # initialize drifters and argo floats drifter = 0 drifters: list[Drifter] = [] argo = 0 argo_floats: list[ArgoFloat] = [] - ctd = 0 + + # the preferred locations of all ctd casts yet to be done + all_ctd_locations = [ + Location(latitude=ctd[1], longitude=ctd[0]) for ctd in config.CTD_locations + ] + remaining_ctd_locations = set(all_ctd_locations) + if len(remaining_ctd_locations) != len(all_ctd_locations): + print("WARN: Some CTD locations are identical and will be combined.") + # ctd cast objects to be used in ctd simulation ctds: list[CTD] = [] route_points_dt = timedelta(minutes=5).total_seconds() @@ -55,33 +62,38 @@ def sailship(config: VirtualShipConfiguration): # run the model for the length of the sample_lons list total_time = timedelta(hours=0).total_seconds() - for i in range(len(sample_lons) - 1): + for i, (sample_lat, sample_lon) in enumerate( + zip(sample_lats, sample_lons, strict=True) + ): + location_here = Location(latitude=sample_lat, longitude=sample_lon) if i % 96 == 0: print(f"Gathered data {timedelta(seconds=total_time)} hours since start.") - # check if virtual ship is at a CTD station - if ctd < len(config.CTD_locations): - if ( - abs(sample_lons[i] - config.CTD_locations[ctd][0]) < 0.01 - and abs(sample_lats[i] - config.CTD_locations[ctd][1]) < 0.01 - ): - ctds.append( - CTD( - location=Location( - latitude=config.CTD_locations[ctd][0], - longitude=config.CTD_locations[ctd][1], - ), - deployment_time=total_time, - min_depth=ctd_min_depth, - max_depth=-config.ctd_fieldset.U.depth[-1], - ) - ) - ctd += 1 - - total_time += timedelta( - minutes=20 - ).total_seconds() # Add 20 minutes for deployment + # find CTD casts to be done at this location + ctds_here = set( + [ + ctd + for ctd in remaining_ctd_locations + if all(np.isclose([ctd.lat, ctd.lon], [sample_lat, sample_lon])) + ] + ) + if len(ctds_here) > 1: + print( + "WARN: Multiple CTD casts match the current location. Only a single cast will be performed." + ) + ctds.append( + CTD( + location=location_here, + deployment_time=total_time, + min_depth=ctd_min_depth, + max_depth=-config.ctd_fieldset.U.depth[-1], + ) + ) + remaining_ctd_locations -= ctds_here + # add 20 minutes to sailing time for deployment + if len(ctds_here) != 0: + total_time += timedelta(minutes=20).total_seconds() # check if we are at a `drifter` deployment location if drifter < len(config.drifter_deploylocations): @@ -133,8 +145,14 @@ def sailship(config: VirtualShipConfiguration): if argo == len(config.argo_deploylocations): break - # update timey + # add time it takes to move to the next route point total_time += route_points_dt + total_time -= route_points_dt + + if len(remaining_ctd_locations) != 0: + print( + "WARN: some CTD casts were not along the route and have not been performed." + ) print("Cruise has ended. Please wait for drifters and/or Argo floats to finish.") @@ -261,7 +279,7 @@ def sailship(config: VirtualShipConfiguration): # return fieldset -def shiproute(config): +def shiproute(config) -> tuple[list[float], list[float]]: """ Take in route coordinates and return lat and lon points within region of interest to sample. From 1716df20bc3d1ccd44c01cb7322fdb97686ac439 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Mon, 3 Jun 2024 17:41:22 +0200 Subject: [PATCH 24/40] cleanup --- virtual_ship/sailship.py | 25 +++++++++++++++---------- 1 file changed, 15 insertions(+), 10 deletions(-) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 4752b828..ee81d040 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -24,6 +24,16 @@ def sailship(config: VirtualShipConfiguration): :param config: The cruise configuration. """ + # TODO this will be in the config later, but for now we don't change the config structure + # from here ----- + ctd_locations_list = [ + Location(latitude=ctd[1], longitude=ctd[0]) for ctd in config.CTD_locations + ] + ctd_locations = set(ctd_locations_list) + if len(ctd_locations) != len(ctd_locations_list): + print("WARN: Some CTD locations are identical and have been combined.") + # until here ---- + # Create fieldset and retreive final schip route as sample_lons and sample_lats adcp_fieldset = config.ctd_fieldset ship_st_fieldset = config.ctd_fieldset @@ -36,13 +46,8 @@ def sailship(config: VirtualShipConfiguration): argo = 0 argo_floats: list[ArgoFloat] = [] - # the preferred locations of all ctd casts yet to be done - all_ctd_locations = [ - Location(latitude=ctd[1], longitude=ctd[0]) for ctd in config.CTD_locations - ] - remaining_ctd_locations = set(all_ctd_locations) - if len(remaining_ctd_locations) != len(all_ctd_locations): - print("WARN: Some CTD locations are identical and will be combined.") + # ctd cast locations that have been visited + ctd_locations_visited: set[Location] = set() # ctd cast objects to be used in ctd simulation ctds: list[CTD] = [] @@ -74,7 +79,7 @@ def sailship(config: VirtualShipConfiguration): ctds_here = set( [ ctd - for ctd in remaining_ctd_locations + for ctd in ctd_locations - ctd_locations_visited if all(np.isclose([ctd.lat, ctd.lon], [sample_lat, sample_lon])) ] ) @@ -90,7 +95,7 @@ def sailship(config: VirtualShipConfiguration): max_depth=-config.ctd_fieldset.U.depth[-1], ) ) - remaining_ctd_locations -= ctds_here + ctd_locations_visited = ctd_locations_visited.union(ctds_here) # add 20 minutes to sailing time for deployment if len(ctds_here) != 0: total_time += timedelta(minutes=20).total_seconds() @@ -149,7 +154,7 @@ def sailship(config: VirtualShipConfiguration): total_time += route_points_dt total_time -= route_points_dt - if len(remaining_ctd_locations) != 0: + if len(ctd_locations_visited) != len(ctd_locations): print( "WARN: some CTD casts were not along the route and have not been performed." ) From 78872ed322cd1b30ca87b9a2473ccd9282bd3467 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 10:35:20 +0200 Subject: [PATCH 25/40] large cleanup --- tests/test_sailship.py | 15 +- virtual_ship/costs.py | 9 +- virtual_ship/sailship.py | 213 +++++++++++---------- virtual_ship/virtual_ship_configuration.py | 12 +- 4 files changed, 144 insertions(+), 105 deletions(-) diff --git a/tests/test_sailship.py b/tests/test_sailship.py index 33f1e290..70b3d0eb 100644 --- a/tests/test_sailship.py +++ b/tests/test_sailship.py @@ -8,13 +8,20 @@ def test_sailship() -> None: + adcp_fieldset = FieldSet.from_data( + {"U": 0, "V": 0}, + {"lon": 0, "lat": 0}, + ) + + ship_st_fieldset = FieldSet.from_data( + {"U": 0, "V": 0, "S": 0, "T": 0}, + {"lon": 0, "lat": 0}, + ) + ctd_fieldset = FieldSet.from_data( {"U": 0, "V": 0, "S": 0, "T": 0, "bathymetry": 0}, {"lon": 0, "lat": 0}, ) - ctd_fieldset.add_constant("maxtime", ctd_fieldset.U.grid.time_full[-1]) - ctd_fieldset.add_constant("mindepth", -ctd_fieldset.U.depth[0]) - ctd_fieldset.add_constant("max_depth", -ctd_fieldset.U.depth[-1]) drifter_fieldset = FieldSet.from_data( {"U": 0, "V": 0, "T": 0}, @@ -36,6 +43,8 @@ def test_sailship() -> None: config = VirtualShipConfiguration( "sailship_config.json", + adcp_fieldset=adcp_fieldset, + ship_st_fieldset=ship_st_fieldset, ctd_fieldset=ctd_fieldset, drifter_fieldset=drifter_fieldset, argo_float_fieldset=argo_float_fieldset, diff --git a/virtual_ship/costs.py b/virtual_ship/costs.py index 6c27f81f..02ee4803 100644 --- a/virtual_ship/costs.py +++ b/virtual_ship/costs.py @@ -1,19 +1,22 @@ """costs function.""" +from datetime import timedelta +from .virtual_ship_configuration import VirtualShipConfiguration -def costs(config, total_time): + +def costs(config: VirtualShipConfiguration, total_time: timedelta): """ Calculate the cost of the virtual ship (in US$). :param config: The cruise configuration. - :param total_time: Time cruised in seconds. + :param total_time: Time cruised. :returns: The calculated cost of the cruise. """ ship_cost_per_day = 30000 drifter_deploy_cost = 2500 argo_deploy_cost = 15000 - ship_cost = ship_cost_per_day / 24 * total_time // 3600 + ship_cost = ship_cost_per_day / 24 * total_time.total_seconds() // 3600 argo_cost = len(config.argo_deploylocations) * argo_deploy_cost drifter_cost = len(config.drifter_deploylocations) * drifter_deploy_cost diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index ee81d040..c0d74739 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -26,6 +26,26 @@ def sailship(config: VirtualShipConfiguration): """ # TODO this will be in the config later, but for now we don't change the config structure # from here ----- + argo_locations_list = [ + Location(latitude=argo[1], longitude=argo[0]) + for argo in config.argo_deploylocations + ] + argo_locations = set(argo_locations_list) + if len(argo_locations) != len(argo_locations_list): + print( + "WARN: Some argo float deployment locations are identical and have been combined." + ) + + drifter_locations_list = [ + Location(latitude=drifter[1], longitude=drifter[0]) + for drifter in config.drifter_deploylocations + ] + drifter_locations = set(drifter_locations_list) + if len(drifter_locations) != len(drifter_locations_list): + print( + "WARN: Some drifter deployment locations are identical and have been combined." + ) + ctd_locations_list = [ Location(latitude=ctd[1], longitude=ctd[0]) for ctd in config.CTD_locations ] @@ -34,53 +54,98 @@ def sailship(config: VirtualShipConfiguration): print("WARN: Some CTD locations are identical and have been combined.") # until here ---- - # Create fieldset and retreive final schip route as sample_lons and sample_lats - adcp_fieldset = config.ctd_fieldset - ship_st_fieldset = config.ctd_fieldset + # get discrete points along the ships route were sampling and deployments will be performed + route_dt = timedelta(minutes=5) + route_points = shiproute(config=config, dt=route_dt) - sample_lons, sample_lats = shiproute(config) - print("Arrived in region of interest, starting to gather data.") + # adcp objects to be used in simulation + adcps: list[ADCPSamplePoint] = [] - drifter = 0 - drifters: list[Drifter] = [] - argo = 0 + # ship st objects to be used in simulation + ship_sts: list[ShipSTSamplePoint] = [] + + # argo float deployment locations that have been visited + argo_locations_visited: set[Location] = set() + # argo float objects to be used in simulation argo_floats: list[ArgoFloat] = [] + # drifter deployment locations that have been visited + drifter_locations_visited: set[Location] = set() + # drifter objects to be used in simulation + drifters: list[Drifter] = [] + # ctd cast locations that have been visited ctd_locations_visited: set[Location] = set() - # ctd cast objects to be used in ctd simulation + # ctd cast objects to be used in simulation ctds: list[CTD] = [] - route_points_dt = timedelta(minutes=5).total_seconds() - adcp_sample_points = [ - ADCPSamplePoint(Location(latitude=lat, longitude=lon), n * route_points_dt) - for n, (lat, lon) in enumerate(zip(sample_lats, sample_lons)) - ] - ship_st_sample_points = [ - ShipSTSamplePoint(Location(latitude=lat, longitude=lon), n * route_points_dt) - for n, (lat, lon) in enumerate(zip(sample_lats, sample_lons)) - ] - - ctd_min_depth = -config.ctd_fieldset.U.depth[0] - argo_min_depth = -config.argo_float_fieldset.U.depth[0] - drifter_min_depth = -config.drifter_fieldset.U.depth[0] + # iterate over each discrete route point, find deployment and measurement locations and times, and measure how much time this took + time_past = timedelta() + for i, route_point in enumerate(route_points): + if i % 96 == 0: + print(f"Gathered data {time_past} hours since start.") - # run the model for the length of the sample_lons list - total_time = timedelta(hours=0).total_seconds() - for i, (sample_lat, sample_lon) in enumerate( - zip(sample_lats, sample_lons, strict=True) - ): - location_here = Location(latitude=sample_lat, longitude=sample_lon) + # find drifter deployments to be done at this location + drifters_here = set( + [ + drifter + for drifter in drifter_locations - drifter_locations_visited + if all( + np.isclose( + [drifter.lat, drifter.lon], [route_point.lat, route_point.lon] + ) + ) + ] + ) + if len(drifters_here) > 1: + print( + "WARN: Multiple drifter deployments match the current location. Only a single deployment will be performed." + ) + drifters.append( + Drifter( + location=route_point, + deployment_time=time_past.total_seconds(), + min_depth=-config.drifter_fieldset.U.depth[0], + ) + ) + drifter_locations_visited = drifter_locations_visited.union(drifters_here) - if i % 96 == 0: - print(f"Gathered data {timedelta(seconds=total_time)} hours since start.") + # find argo float deployments to be done at this location + argos_here = set( + [ + argo + for argo in argo_locations - argo_locations_visited + if all( + np.isclose([argo.lat, argo.lon], [route_point.lat, route_point.lon]) + ) + ] + ) + if len(argos_here) > 1: + print( + "WARN: Multiple argo float deployments match the current location. Only a single deployment will be performed." + ) + argo_floats.append( + ArgoFloat( + location=route_point, + deployment_time=time_past.total_seconds(), + min_depth=-config.argo_float_fieldset.U.depth[0], + max_depth=config.argo_characteristics["maxdepth"], + drift_depth=config.argo_characteristics["driftdepth"], + vertical_speed=config.argo_characteristics["vertical_speed"], + cycle_days=config.argo_characteristics["cycle_days"], + drift_days=config.argo_characteristics["drift_days"], + ) + ) + argo_locations_visited = argo_locations_visited.union(argos_here) # find CTD casts to be done at this location ctds_here = set( [ ctd for ctd in ctd_locations - ctd_locations_visited - if all(np.isclose([ctd.lat, ctd.lon], [sample_lat, sample_lon])) + if all( + np.isclose([ctd.lat, ctd.lon], [route_point.lat, route_point.lon]) + ) ] ) if len(ctds_here) > 1: @@ -89,70 +154,20 @@ def sailship(config: VirtualShipConfiguration): ) ctds.append( CTD( - location=location_here, - deployment_time=total_time, - min_depth=ctd_min_depth, + location=route_point, + deployment_time=time_past.total_seconds(), + min_depth=-config.ctd_fieldset.U.depth[0], max_depth=-config.ctd_fieldset.U.depth[-1], ) ) ctd_locations_visited = ctd_locations_visited.union(ctds_here) # add 20 minutes to sailing time for deployment if len(ctds_here) != 0: - total_time += timedelta(minutes=20).total_seconds() - - # check if we are at a `drifter` deployment location - if drifter < len(config.drifter_deploylocations): - while ( - abs(sample_lons[i] - config.drifter_deploylocations[drifter][0]) < 0.01 - and abs(sample_lats[i] - config.drifter_deploylocations[drifter][1]) - < 0.01 - ): - drifters.append( - Drifter( - location=Location( - latitude=config.drifter_deploylocations[drifter][0], - longitude=config.drifter_deploylocations[drifter][1], - ), - deployment_time=total_time, - min_depth=drifter_min_depth, - ) - ) - drifter += 1 - print( - f"Drifter {drifter} deployed at {sample_lons[i]}, {sample_lats[i]}" - ) - if drifter == len(config.drifter_deploylocations): - break - - # check if we are at a argo deployment location - if argo < len(config.argo_deploylocations): - while ( - abs(sample_lons[i] - config.argo_deploylocations[argo][0]) < 0.01 - and abs(sample_lats[i] - config.argo_deploylocations[argo][1]) < 0.01 - ): - argo_floats.append( - ArgoFloat( - location=Location( - latitude=config.argo_deploylocations[argo][0], - longitude=config.argo_deploylocations[argo][1], - ), - deployment_time=total_time, - min_depth=argo_min_depth, - max_depth=config.argo_characteristics["maxdepth"], - drift_depth=config.argo_characteristics["driftdepth"], - vertical_speed=config.argo_characteristics["vertical_speed"], - cycle_days=config.argo_characteristics["cycle_days"], - drift_days=config.argo_characteristics["drift_days"], - ) - ) - argo += 1 - print(f"Argo {argo} deployed at {sample_lons[i]}, {sample_lats[i]}") - if argo == len(config.argo_deploylocations): - break + time_past += timedelta(minutes=20) # add time it takes to move to the next route point - total_time += route_points_dt - total_time -= route_points_dt + time_past += route_dt + time_past -= route_dt if len(ctd_locations_visited) != len(ctd_locations): print( @@ -163,20 +178,20 @@ def sailship(config: VirtualShipConfiguration): print("Simulating onboard salinity and temperature measurements.") simulate_ship_st( - fieldset=ship_st_fieldset, + fieldset=config.ship_st_fieldset, out_file_name=os.path.join("results", "ship_st.zarr"), depth=-2, - sample_points=ship_st_sample_points, + sample_points=ship_sts, ) print("Simulating onboard ADCP.") simulate_adcp( - fieldset=adcp_fieldset, + fieldset=config.adcp_fieldset, out_file_name=os.path.join("results", "adcp.zarr"), max_depth=config.ADCP_settings["max_depth"], min_depth=-5, bin_size=config.ADCP_settings["bin_size_m"], - sample_points=adcp_sample_points, + sample_points=adcps, ) print("Simulating CTD casts.") @@ -209,10 +224,8 @@ def sailship(config: VirtualShipConfiguration): print("All data has been gathered and postprocessed, returning home.") - cost = costs(config, total_time) - print( - f"This cruise took {timedelta(seconds=total_time)} and would have cost {cost:,.0f} euros." - ) + cost = costs(config, time_past) + print(f"This cruise took {time_past} and would have cost {cost:,.0f} euros.") # def create_fieldset(config, data_dir: str): @@ -284,7 +297,7 @@ def sailship(config: VirtualShipConfiguration): # return fieldset -def shiproute(config) -> tuple[list[float], list[float]]: +def shiproute(config: VirtualShipConfiguration, dt: timedelta) -> list[Location]: """ Take in route coordinates and return lat and lon points within region of interest to sample. @@ -308,7 +321,7 @@ def shiproute(config) -> tuple[list[float], list[float]]: cruise_speed = 5.14 geod = pyproj.Geod(ellps="WGS84") azimuth1, azimuth2, distance = geod.inv(startlong, startlat, endlong, endlat) - if distance > (cruise_speed * 60 * 5): + if distance > (cruise_speed * dt.total_seconds()): r = geod.inv_intermediate( startlong, startlat, @@ -340,4 +353,8 @@ def shiproute(config) -> tuple[list[float], list[float]]: if poly.contains(Point(lons[i], lats[i])): sample_lons.append(lons[i]) sample_lats.append(lats[i]) - return sample_lons, sample_lats + points = [ + Location(latitude=lat, longitude=lon) + for lat, lon in zip(sample_lats, sample_lons, strict=True) + ] + return points diff --git a/virtual_ship/virtual_ship_configuration.py b/virtual_ship/virtual_ship_configuration.py index f74251b7..b3300d53 100644 --- a/virtual_ship/virtual_ship_configuration.py +++ b/virtual_ship/virtual_ship_configuration.py @@ -11,6 +11,10 @@ class VirtualShipConfiguration: """Configuration of the virtual ship, initialized from a json file.""" + adcp_fieldset: FieldSet + ship_st_fieldset: FieldSet + argo_float_fieldset: FieldSet + drifter_fieldset: FieldSet ctd_fieldset: FieldSet drifter_fieldset: FieldSet argo_float_fieldset: FieldSet @@ -18,6 +22,8 @@ class VirtualShipConfiguration: def __init__( self, json_file, + adcp_fieldset: FieldSet, + ship_st_fieldset: FieldSet, ctd_fieldset: FieldSet, drifter_fieldset: FieldSet, argo_float_fieldset: FieldSet, @@ -26,11 +32,15 @@ def __init__( Initialize this object. :param json_file: Path to the JSON file to init from. + :param adcp_fieldset: Fieldset for ADCP measurements. + :param ship_st_fieldset: Fieldset for ship salinity temperature measurements. :param ctd_fieldset: Fieldset for CTD measurements. :param drifter_fieldset: Fieldset for CTD measurements. - :param argo_float_fieldset: FIeldset for argo float measurements. + :param argo_float_fieldset: Fieldset for argo float measurements. :raises ValueError: If JSON file not valid. """ + self.adcp_fieldset = adcp_fieldset + self.ship_st_fieldset = ship_st_fieldset self.ctd_fieldset = ctd_fieldset self.drifter_fieldset = drifter_fieldset self.argo_float_fieldset = argo_float_fieldset From 087371d138a5eac178e28cc9309092f641fef6d3 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 10:47:55 +0200 Subject: [PATCH 26/40] add checks that all drifter and argo floats are deployed --- virtual_ship/sailship.py | 17 ++++++++++++++--- 1 file changed, 14 insertions(+), 3 deletions(-) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index c0d74739..00d51f38 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -80,6 +80,7 @@ def sailship(config: VirtualShipConfiguration): ctds: list[CTD] = [] # iterate over each discrete route point, find deployment and measurement locations and times, and measure how much time this took + print("Traversing ship route.") time_past = timedelta() for i, route_point in enumerate(route_points): if i % 96 == 0: @@ -167,14 +168,24 @@ def sailship(config: VirtualShipConfiguration): # add time it takes to move to the next route point time_past += route_dt + # remove the last one, because no sailing to the next point was needed time_past -= route_dt - if len(ctd_locations_visited) != len(ctd_locations): + # check if all drifter, argo float, ctd locations were visited + if len(drifter_locations_visited) != len(drifter_locations): + print( + "WARN: some drifter deployments were not planned along the route and have not been performed." + ) + + if len(argo_locations_visited) != len(argo_locations): print( - "WARN: some CTD casts were not along the route and have not been performed." + "WARN: some argo float deployments were not planned along the route and have not been performed." ) - print("Cruise has ended. Please wait for drifters and/or Argo floats to finish.") + if len(ctd_locations_visited) != len(ctd_locations): + print( + "WARN: some CTD casts were not planned along the route and have not been performed." + ) print("Simulating onboard salinity and temperature measurements.") simulate_ship_st( From 3e10d46fb5602caece28eef2c29b24f73520858f Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 10:49:48 +0200 Subject: [PATCH 27/40] comments --- virtual_ship/sailship.py | 1 + 1 file changed, 1 insertion(+) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 00d51f38..c6d34d90 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -80,6 +80,7 @@ def sailship(config: VirtualShipConfiguration): ctds: list[CTD] = [] # iterate over each discrete route point, find deployment and measurement locations and times, and measure how much time this took + # TODO between drifters, argo floats, ctd there is quite some repetition print("Traversing ship route.") time_past = timedelta() for i, route_point in enumerate(route_points): From 47c7e05faf4c959ef99e6bf53b9426666578d938 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 11:05:00 +0200 Subject: [PATCH 28/40] docstrings and other codetools fixes --- tests/instruments/test_adcp.py | 6 +++--- tests/instruments/test_ship_st.py | 6 +++--- virtual_ship/costs.py | 1 + virtual_ship/instruments/__init__.py | 3 ++- virtual_ship/instruments/adcp.py | 23 ++++++++++++----------- virtual_ship/instruments/ctd.py | 2 +- virtual_ship/instruments/ship_st.py | 21 ++++++++++----------- virtual_ship/instruments/spacetime.py | 14 ++++++++++++++ virtual_ship/sailship.py | 11 ++++++----- 9 files changed, 52 insertions(+), 35 deletions(-) create mode 100644 virtual_ship/instruments/spacetime.py diff --git a/tests/instruments/test_adcp.py b/tests/instruments/test_adcp.py index 84540c2e..06808879 100644 --- a/tests/instruments/test_adcp.py +++ b/tests/instruments/test_adcp.py @@ -3,8 +3,8 @@ import numpy as np from parcels import FieldSet -from virtual_ship.instruments import Location -from virtual_ship.instruments.adcp import simulate_adcp, SamplePoint +from virtual_ship.instruments import Location, Spacetime +from virtual_ship.instruments.adcp import simulate_adcp def test_simulate_argo_floats() -> None: @@ -21,7 +21,7 @@ def test_simulate_argo_floats() -> None: }, ) - sample_points = [SamplePoint(Location(0, 0), 0)] + sample_points = [Spacetime(Location(0, 0), 0)] simulate_adcp( fieldset=fieldset, diff --git a/tests/instruments/test_ship_st.py b/tests/instruments/test_ship_st.py index f4871859..0556e586 100644 --- a/tests/instruments/test_ship_st.py +++ b/tests/instruments/test_ship_st.py @@ -3,8 +3,8 @@ import numpy as np from parcels import FieldSet -from virtual_ship.instruments import Location -from virtual_ship.instruments.ship_st import simulate_ship_st, SamplePoint +from virtual_ship.instruments import Location, Spacetime +from virtual_ship.instruments.ship_st import simulate_ship_st def test_simulate_argo_floats() -> None: @@ -19,7 +19,7 @@ def test_simulate_argo_floats() -> None: }, ) - sample_points = [SamplePoint(Location(0, 0), 0)] + sample_points = [Spacetime(Location(0, 0), 0)] simulate_ship_st( fieldset=fieldset, diff --git a/virtual_ship/costs.py b/virtual_ship/costs.py index 02ee4803..9ba1d970 100644 --- a/virtual_ship/costs.py +++ b/virtual_ship/costs.py @@ -1,6 +1,7 @@ """costs function.""" from datetime import timedelta + from .virtual_ship_configuration import VirtualShipConfiguration diff --git a/virtual_ship/instruments/__init__.py b/virtual_ship/instruments/__init__.py index f6d89a5a..b70f85ba 100644 --- a/virtual_ship/instruments/__init__.py +++ b/virtual_ship/instruments/__init__.py @@ -2,5 +2,6 @@ from . import argo_float, ctd, drifter from .location import Location +from .spacetime import Spacetime -__all__ = ["Location", "argo_float", "ctd", "drifter"] +__all__ = ["Location", "Spacetime", "argo_float", "ctd", "drifter"] diff --git a/virtual_ship/instruments/adcp.py b/virtual_ship/instruments/adcp.py index 7868aa7b..b6f5ef25 100644 --- a/virtual_ship/instruments/adcp.py +++ b/virtual_ship/instruments/adcp.py @@ -1,18 +1,9 @@ """ADCP instrument.""" -from dataclasses import dataclass - import numpy as np from parcels import FieldSet, JITParticle, ParticleSet, Variable -from .location import Location - - -@dataclass -class SamplePoint: - location: Location - time: float - +from .spacetime import Spacetime _ADCPParticle = JITParticle.add_variables( [ @@ -34,8 +25,18 @@ def simulate_adcp( max_depth: float, min_depth: float, bin_size: float, - sample_points: list[SamplePoint], + sample_points: list[Spacetime], ) -> None: + """ + Use parcels to simulate an ADCP in a fieldset. + + :param fieldset: The fieldset to simulate the ADCP in. + :param out_file_name: The file to write the results to. + :param max_depth: Maximum depth the ADCP can measure. + :param min_depth: Minimum depth the ADCP can measure. + :param bin_size: How many samples to take in the complete range between max_depth and min_depth. + :param sample_points: The places and times to sample at. + """ sample_points.sort(key=lambda p: p.time) bins = np.arange(max_depth, min_depth, bin_size) diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index 4271d63c..281ae6d3 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -76,7 +76,7 @@ def simulate_ctd( """ Use parcels to simulate a set of CTDs in a fieldset. - :param ctd: A list of CTDs to simulate. + :param ctds: A list of CTDs to simulate. :param fieldset: The fieldset to simulate the CTDs in. :param out_file_name: The file to write the results to. :param outputdt: Interval which dictates the update frequency of file output during simulation diff --git a/virtual_ship/instruments/ship_st.py b/virtual_ship/instruments/ship_st.py index 80004799..eb40fbe1 100644 --- a/virtual_ship/instruments/ship_st.py +++ b/virtual_ship/instruments/ship_st.py @@ -1,18 +1,9 @@ """Ship salinity and temperature.""" -from dataclasses import dataclass - import numpy as np from parcels import FieldSet, JITParticle, ParticleSet, Variable -from .location import Location - - -@dataclass -class SamplePoint: - location: Location - time: float - +from .spacetime import Spacetime _ShipSTParticle = JITParticle.add_variables( [ @@ -36,8 +27,16 @@ def simulate_ship_st( fieldset: FieldSet, out_file_name: str, depth: float, - sample_points: list[SamplePoint], + sample_points: list[Spacetime], ) -> None: + """ + Use parcels to simulate the measurement of salinity and temperature aboard a ship in a fieldset. + + :param fieldset: The fieldset to simulate the Argo floats in. + :param out_file_name: The file to write the results to. + :param depth: The depth at which to measure. 0 is water surface, negative is into the water. + :param sample_points: The places and times to sample at. + """ sample_points.sort(key=lambda p: p.time) particleset = ParticleSet.from_list( diff --git a/virtual_ship/instruments/spacetime.py b/virtual_ship/instruments/spacetime.py new file mode 100644 index 00000000..c277e642 --- /dev/null +++ b/virtual_ship/instruments/spacetime.py @@ -0,0 +1,14 @@ +"""Location class. See class description.""" + +from dataclasses import dataclass + +from .location import Location + + +@dataclass +# TODO I take suggestions for a better name +class Spacetime: + """A location and time.""" + + location: Location + time: float diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index c6d34d90..6d012b29 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -8,12 +8,12 @@ from shapely.geometry import Point, Polygon from .costs import costs +from .instruments import Location, Spacetime +from .instruments.adcp import simulate_adcp from .instruments.argo_float import ArgoFloat, simulate_argo_floats from .instruments.ctd import CTD, simulate_ctd from .instruments.drifter import Drifter, simulate_drifters -from .instruments.adcp import simulate_adcp, SamplePoint as ADCPSamplePoint -from .instruments.ship_st import simulate_ship_st, SamplePoint as ShipSTSamplePoint -from .instruments.location import Location +from .instruments.ship_st import simulate_ship_st from .postprocess import postprocess from .virtual_ship_configuration import VirtualShipConfiguration @@ -59,10 +59,10 @@ def sailship(config: VirtualShipConfiguration): route_points = shiproute(config=config, dt=route_dt) # adcp objects to be used in simulation - adcps: list[ADCPSamplePoint] = [] + adcps: list[Spacetime] = [] # ship st objects to be used in simulation - ship_sts: list[ShipSTSamplePoint] = [] + ship_sts: list[Spacetime] = [] # argo float deployment locations that have been visited argo_locations_visited: set[Location] = set() @@ -314,6 +314,7 @@ def shiproute(config: VirtualShipConfiguration, dt: timedelta) -> list[Location] Take in route coordinates and return lat and lon points within region of interest to sample. :param config: The cruise configuration. + :param dt: Sailing time between each discrete route point. :returns: lat and lon points within region of interest to sample. """ # Initialize lists to store intermediate points From 072437cf149f28386612808f092e9f3e294f2607 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 11:15:40 +0200 Subject: [PATCH 29/40] rename SamplePoint to Spacetime and move them one dir up --- tests/instruments/test_adcp.py | 2 +- tests/instruments/test_argo.py | 5 ++--- tests/instruments/test_ctd.py | 5 ++--- tests/instruments/test_drifters.py | 5 ++--- tests/instruments/test_ship_st.py | 2 +- virtual_ship/__init__.py | 4 +++- virtual_ship/instruments/__init__.py | 6 ++---- virtual_ship/instruments/adcp.py | 2 +- virtual_ship/instruments/argo_float.py | 11 +++++------ virtual_ship/instruments/ctd.py | 11 +++++------ virtual_ship/instruments/drifter.py | 11 +++++------ virtual_ship/instruments/ship_st.py | 2 +- virtual_ship/{instruments => }/location.py | 0 virtual_ship/sailship.py | 18 +++++++++++------- virtual_ship/{instruments => }/spacetime.py | 0 15 files changed, 41 insertions(+), 43 deletions(-) rename virtual_ship/{instruments => }/location.py (100%) rename virtual_ship/{instruments => }/spacetime.py (100%) diff --git a/tests/instruments/test_adcp.py b/tests/instruments/test_adcp.py index 06808879..bb03a90d 100644 --- a/tests/instruments/test_adcp.py +++ b/tests/instruments/test_adcp.py @@ -3,7 +3,7 @@ import numpy as np from parcels import FieldSet -from virtual_ship.instruments import Location, Spacetime +from virtual_ship import Location, Spacetime from virtual_ship.instruments.adcp import simulate_adcp diff --git a/tests/instruments/test_argo.py b/tests/instruments/test_argo.py index 83bbb261..1b6493ea 100644 --- a/tests/instruments/test_argo.py +++ b/tests/instruments/test_argo.py @@ -5,7 +5,7 @@ import numpy as np from parcels import FieldSet -from virtual_ship.instruments import Location +from virtual_ship import Location, Spacetime from virtual_ship.instruments.argo_float import ArgoFloat, simulate_argo_floats @@ -29,8 +29,7 @@ def test_simulate_argo_floats() -> None: argo_floats = [ ArgoFloat( - location=Location(latitude=0, longitude=0), - deployment_time=0, + spacetime=Spacetime(location=Location(latitude=0, longitude=0), time=0), min_depth=min_depth, max_depth=MAX_DEPTH, drift_depth=DRIFT_DEPTH, diff --git a/tests/instruments/test_ctd.py b/tests/instruments/test_ctd.py index ba6aa81a..793e2f6c 100644 --- a/tests/instruments/test_ctd.py +++ b/tests/instruments/test_ctd.py @@ -5,7 +5,7 @@ import numpy as np from parcels import FieldSet -from virtual_ship.instruments import Location +from virtual_ship import Location, Spacetime from virtual_ship.instruments.ctd import CTD, simulate_ctd @@ -24,8 +24,7 @@ def test_simulate_drifters() -> None: ctds = [ CTD( - location=Location(latitude=0, longitude=0), - deployment_time=0, + spacetime=Spacetime(location=Location(latitude=0, longitude=0), time=0), min_depth=min_depth, max_depth=max_depth, ) diff --git a/tests/instruments/test_drifters.py b/tests/instruments/test_drifters.py index dafbb934..421499d6 100644 --- a/tests/instruments/test_drifters.py +++ b/tests/instruments/test_drifters.py @@ -5,7 +5,7 @@ import numpy as np from parcels import FieldSet -from virtual_ship.instruments import Location +from virtual_ship import Location, Spacetime from virtual_ship.instruments.drifter import Drifter, simulate_drifters @@ -23,8 +23,7 @@ def test_simulate_drifters() -> None: drifters = [ Drifter( - location=Location(latitude=0, longitude=0), - deployment_time=0, + spacetime=Spacetime(location=Location(latitude=0, longitude=0), time=0), min_depth=min_depth, ) ] diff --git a/tests/instruments/test_ship_st.py b/tests/instruments/test_ship_st.py index 0556e586..7764364b 100644 --- a/tests/instruments/test_ship_st.py +++ b/tests/instruments/test_ship_st.py @@ -3,7 +3,7 @@ import numpy as np from parcels import FieldSet -from virtual_ship.instruments import Location, Spacetime +from virtual_ship import Location, Spacetime from virtual_ship.instruments.ship_st import simulate_ship_st diff --git a/virtual_ship/__init__.py b/virtual_ship/__init__.py index bf152389..dbf175e3 100644 --- a/virtual_ship/__init__.py +++ b/virtual_ship/__init__.py @@ -1,5 +1,7 @@ """Code for the Virtual Ship Classroom, where Marine Scientists can combine Copernicus Marine Data with an OceanParcels ship to go on a virtual expedition.""" from . import instruments, sailship +from .location import Location +from .spacetime import Spacetime -__all__ = ["instruments", "sailship"] +__all__ = ["Location", "Spacetime", "instruments", "sailship"] diff --git a/virtual_ship/instruments/__init__.py b/virtual_ship/instruments/__init__.py index b70f85ba..9ca81420 100644 --- a/virtual_ship/instruments/__init__.py +++ b/virtual_ship/instruments/__init__.py @@ -1,7 +1,5 @@ """Measurement instrument that can be used with Parcels.""" -from . import argo_float, ctd, drifter -from .location import Location -from .spacetime import Spacetime +from . import adcp, argo_float, ctd, drifter, ship_st -__all__ = ["Location", "Spacetime", "argo_float", "ctd", "drifter"] +__all__ = ["adcp", "argo_float", "ctd", "drifter", "ship_st"] diff --git a/virtual_ship/instruments/adcp.py b/virtual_ship/instruments/adcp.py index b6f5ef25..2d7d032c 100644 --- a/virtual_ship/instruments/adcp.py +++ b/virtual_ship/instruments/adcp.py @@ -3,7 +3,7 @@ import numpy as np from parcels import FieldSet, JITParticle, ParticleSet, Variable -from .spacetime import Spacetime +from ..spacetime import Spacetime _ADCPParticle = JITParticle.add_variables( [ diff --git a/virtual_ship/instruments/argo_float.py b/virtual_ship/instruments/argo_float.py index 1218d965..cd014671 100644 --- a/virtual_ship/instruments/argo_float.py +++ b/virtual_ship/instruments/argo_float.py @@ -14,15 +14,14 @@ Variable, ) -from .location import Location +from ..spacetime import Spacetime @dataclass class ArgoFloat: """Configuration for a single Argo float.""" - location: Location - deployment_time: float + spacetime: Spacetime min_depth: float max_depth: float drift_depth: float @@ -129,9 +128,9 @@ def simulate_argo_floats( :param out_file_name: The file to write the results to. :param outputdt: Interval which dictates the update frequency of file output during simulation """ - lon = [argo.location.lon for argo in argo_floats] - lat = [argo.location.lat for argo in argo_floats] - time = [argo.deployment_time for argo in argo_floats] + lon = [argo.spacetime.location.lon for argo in argo_floats] + lat = [argo.spacetime.location.lat for argo in argo_floats] + time = [argo.spacetime.time for argo in argo_floats] # define parcel particles argo_float_particleset = ParticleSet( diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index 281ae6d3..fe0c909d 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -6,15 +6,14 @@ import numpy as np from parcels import FieldSet, JITParticle, ParticleSet, Variable -from .location import Location +from ..spacetime import Spacetime @dataclass class CTD: """Configuration for a single CTD.""" - location: Location - deployment_time: float + spacetime: Spacetime min_depth: float max_depth: float @@ -81,9 +80,9 @@ def simulate_ctd( :param out_file_name: The file to write the results to. :param outputdt: Interval which dictates the update frequency of file output during simulation """ - lon = [ctd.location.lon for ctd in ctds] - lat = [ctd.location.lat for ctd in ctds] - time = [ctd.deployment_time for ctd in ctds] + lon = [ctd.spacetime.location.lon for ctd in ctds] + lat = [ctd.spacetime.location.lat for ctd in ctds] + time = [ctd.spacetime.time for ctd in ctds] # define parcel particles ctd_particleset = ParticleSet( diff --git a/virtual_ship/instruments/drifter.py b/virtual_ship/instruments/drifter.py index fd38782c..5bed859f 100644 --- a/virtual_ship/instruments/drifter.py +++ b/virtual_ship/instruments/drifter.py @@ -6,15 +6,14 @@ import numpy as np from parcels import AdvectionRK4, FieldSet, JITParticle, ParticleSet, Variable -from .location import Location +from ..spacetime import Spacetime @dataclass class Drifter: """Configuration for a single Drifter.""" - location: Location - deployment_time: float + spacetime: Spacetime min_depth: float @@ -48,9 +47,9 @@ def simulate_drifters( :param out_file_name: The file to write the results to. :param outputdt: Interval which dictates the update frequency of file output during simulation """ - lon = [drifter.location.lon for drifter in drifters] - lat = [drifter.location.lat for drifter in drifters] - time = [drifter.deployment_time for drifter in drifters] + lon = [drifter.spacetime.location.lon for drifter in drifters] + lat = [drifter.spacetime.location.lat for drifter in drifters] + time = [drifter.spacetime.time for drifter in drifters] # define parcel particles drifter_particleset = ParticleSet( diff --git a/virtual_ship/instruments/ship_st.py b/virtual_ship/instruments/ship_st.py index eb40fbe1..6aa733a9 100644 --- a/virtual_ship/instruments/ship_st.py +++ b/virtual_ship/instruments/ship_st.py @@ -3,7 +3,7 @@ import numpy as np from parcels import FieldSet, JITParticle, ParticleSet, Variable -from .spacetime import Spacetime +from ..spacetime import Spacetime _ShipSTParticle = JITParticle.add_variables( [ diff --git a/virtual_ship/instruments/location.py b/virtual_ship/location.py similarity index 100% rename from virtual_ship/instruments/location.py rename to virtual_ship/location.py diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 6d012b29..cbc36372 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -8,13 +8,14 @@ from shapely.geometry import Point, Polygon from .costs import costs -from .instruments import Location, Spacetime from .instruments.adcp import simulate_adcp from .instruments.argo_float import ArgoFloat, simulate_argo_floats from .instruments.ctd import CTD, simulate_ctd from .instruments.drifter import Drifter, simulate_drifters from .instruments.ship_st import simulate_ship_st +from .location import Location from .postprocess import postprocess +from .spacetime import Spacetime from .virtual_ship_configuration import VirtualShipConfiguration @@ -105,8 +106,9 @@ def sailship(config: VirtualShipConfiguration): ) drifters.append( Drifter( - location=route_point, - deployment_time=time_past.total_seconds(), + spacetime=Spacetime( + location=route_point, time=time_past.total_seconds() + ), min_depth=-config.drifter_fieldset.U.depth[0], ) ) @@ -128,8 +130,9 @@ def sailship(config: VirtualShipConfiguration): ) argo_floats.append( ArgoFloat( - location=route_point, - deployment_time=time_past.total_seconds(), + spacetime=Spacetime( + location=route_point, time=time_past.total_seconds() + ), min_depth=-config.argo_float_fieldset.U.depth[0], max_depth=config.argo_characteristics["maxdepth"], drift_depth=config.argo_characteristics["driftdepth"], @@ -156,8 +159,9 @@ def sailship(config: VirtualShipConfiguration): ) ctds.append( CTD( - location=route_point, - deployment_time=time_past.total_seconds(), + spacetime=Spacetime( + location=route_point, time=time_past.total_seconds() + ), min_depth=-config.ctd_fieldset.U.depth[0], max_depth=-config.ctd_fieldset.U.depth[-1], ) diff --git a/virtual_ship/instruments/spacetime.py b/virtual_ship/spacetime.py similarity index 100% rename from virtual_ship/instruments/spacetime.py rename to virtual_ship/spacetime.py From c717f5c252fdf985c38dae3a5487f0697839f953 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 11:16:10 +0200 Subject: [PATCH 30/40] rm todo --- todo | 3 --- 1 file changed, 3 deletions(-) delete mode 100644 todo diff --git a/todo b/todo deleted file mode 100644 index 511fc884..00000000 --- a/todo +++ /dev/null @@ -1,3 +0,0 @@ -rename outputdt -end time for instruments -check dtypes for particles \ No newline at end of file From 1241f7040cd529ef543965cfccc0ee927124d36d Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 11:25:02 +0200 Subject: [PATCH 31/40] run pydocstyle on tests --- codetools.sh | 2 +- tests/conftest.py | 9 ++++++++- tests/instruments/test_adcp.py | 2 +- 3 files changed, 10 insertions(+), 3 deletions(-) diff --git a/codetools.sh b/codetools.sh index a359c5d6..8e0a4acc 100755 --- a/codetools.sh +++ b/codetools.sh @@ -20,7 +20,7 @@ flake8 ./$PACKAGE ./$TESTS echo "--------------" echo "pydocstyle" echo "--------------" -pydocstyle ./$PACKAGE +pydocstyle ./$PACKAGE ./$TESTS echo "--------------" echo "sort-all" diff --git a/tests/conftest.py b/tests/conftest.py index 01e91653..f656c38f 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,7 +1,14 @@ +"""Test configuration that is ran for every test.""" + import pytest -# Set the working directory for each test to the directory of that test. @pytest.fixture(autouse=True) def change_test_dir(request, monkeypatch): + """ + Set the working directory for each test to the directory of that test. + + :param request: - + :param monkeypatch: - + """ monkeypatch.chdir(request.fspath.dirname) diff --git a/tests/instruments/test_adcp.py b/tests/instruments/test_adcp.py index bb03a90d..6b472f32 100644 --- a/tests/instruments/test_adcp.py +++ b/tests/instruments/test_adcp.py @@ -1,4 +1,4 @@ -"""Test the simulation of Argo floats.""" +"""Test the simulation of ADCP instruments.""" import numpy as np from parcels import FieldSet From 55c290db4951e11c3928b8e07bbf78dc31219978 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Tue, 4 Jun 2024 14:48:06 +0200 Subject: [PATCH 32/40] fix some names --- tests/instruments/test_adcp.py | 2 +- tests/instruments/test_ctd.py | 2 +- tests/instruments/test_ship_st.py | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/tests/instruments/test_adcp.py b/tests/instruments/test_adcp.py index 6b472f32..fdbd50e3 100644 --- a/tests/instruments/test_adcp.py +++ b/tests/instruments/test_adcp.py @@ -7,7 +7,7 @@ from virtual_ship.instruments.adcp import simulate_adcp -def test_simulate_argo_floats() -> None: +def test_simulate_adcp() -> None: MAX_DEPTH = -1000 MIN_DEPTH = -5 BIN_SIZE = 24 diff --git a/tests/instruments/test_ctd.py b/tests/instruments/test_ctd.py index 793e2f6c..92c456cd 100644 --- a/tests/instruments/test_ctd.py +++ b/tests/instruments/test_ctd.py @@ -9,7 +9,7 @@ from virtual_ship.instruments.ctd import CTD, simulate_ctd -def test_simulate_drifters() -> None: +def test_simulate_ctds() -> None: fieldset = FieldSet.from_data( {"U": 0, "V": 0, "T": 0, "S": 0, "bathymetry": 100}, { diff --git a/tests/instruments/test_ship_st.py b/tests/instruments/test_ship_st.py index 7764364b..e8c0d513 100644 --- a/tests/instruments/test_ship_st.py +++ b/tests/instruments/test_ship_st.py @@ -7,7 +7,7 @@ from virtual_ship.instruments.ship_st import simulate_ship_st -def test_simulate_argo_floats() -> None: +def test_simulate_ship_st() -> None: DEPTH = -2 fieldset = FieldSet.from_data( From ae077e631c6657fd6cb8035682da375249cff516 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 16:35:32 +0200 Subject: [PATCH 33/40] Update virtual_ship/instruments/ctd.py Co-authored-by: Erik van Sebille --- virtual_ship/instruments/ctd.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index fe0c909d..c076f3d9 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -22,7 +22,7 @@ class CTD: [ Variable("salinity", dtype=np.float32, initial=np.nan), Variable("temperature", dtype=np.float32, initial=np.nan), - Variable("raising", dtype=np.int32, initial=0.0), + Variable("raising", dtype=np.int32, initial=0), Variable("max_depth", dtype=np.float32), ] ) From b7ee4ff151598f4170e8fbbf8d5fb78ac5e5ed93 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 16:36:28 +0200 Subject: [PATCH 34/40] removed comment about ctd --- virtual_ship/instruments/ctd.py | 1 - 1 file changed, 1 deletion(-) diff --git a/virtual_ship/instruments/ctd.py b/virtual_ship/instruments/ctd.py index c076f3d9..89990970 100644 --- a/virtual_ship/instruments/ctd.py +++ b/virtual_ship/instruments/ctd.py @@ -38,7 +38,6 @@ def _sample_salinity(particle, fieldset, time): def _ctd_cast(particle, fieldset, time): # Lowering and raising CTD - # TODO question: if is executed every time... move outside function? Not if "drifting" now possible if ( -fieldset.bathymetry[time, particle.depth, particle.lat, particle.lon] > particle.max_depth From dfe461d6cbb97ed4c991b1f8f4cf0957fc1b5c23 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 16:38:22 +0200 Subject: [PATCH 35/40] fix incorrect docstring --- virtual_ship/instruments/ship_st.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/virtual_ship/instruments/ship_st.py b/virtual_ship/instruments/ship_st.py index 6aa733a9..1932a3a4 100644 --- a/virtual_ship/instruments/ship_st.py +++ b/virtual_ship/instruments/ship_st.py @@ -32,7 +32,7 @@ def simulate_ship_st( """ Use parcels to simulate the measurement of salinity and temperature aboard a ship in a fieldset. - :param fieldset: The fieldset to simulate the Argo floats in. + :param fieldset: The fieldset to simulate the sampling in. :param out_file_name: The file to write the results to. :param depth: The depth at which to measure. 0 is water surface, negative is into the water. :param sample_points: The places and times to sample at. From 4b065a582bb024d953786c47045e7242c209f59d Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 16:40:54 +0200 Subject: [PATCH 36/40] renamed ship_st to ship_underwater_st --- .../{test_ship_st.py => test_ship_underwater_st.py} | 6 +++--- tests/test_sailship.py | 4 ++-- virtual_ship/instruments/__init__.py | 4 ++-- .../{ship_st.py => ship_underwater_st.py} | 2 +- virtual_ship/sailship.py | 12 ++++++------ virtual_ship/virtual_ship_configuration.py | 8 ++++---- 6 files changed, 18 insertions(+), 18 deletions(-) rename tests/instruments/{test_ship_st.py => test_ship_underwater_st.py} (78%) rename virtual_ship/instruments/{ship_st.py => ship_underwater_st.py} (98%) diff --git a/tests/instruments/test_ship_st.py b/tests/instruments/test_ship_underwater_st.py similarity index 78% rename from tests/instruments/test_ship_st.py rename to tests/instruments/test_ship_underwater_st.py index e8c0d513..05de1d49 100644 --- a/tests/instruments/test_ship_st.py +++ b/tests/instruments/test_ship_underwater_st.py @@ -4,10 +4,10 @@ from parcels import FieldSet from virtual_ship import Location, Spacetime -from virtual_ship.instruments.ship_st import simulate_ship_st +from virtual_ship.instruments.ship_underwater_st import simulate_ship_underwater_st -def test_simulate_ship_st() -> None: +def test_simulate_ship_underwater_st() -> None: DEPTH = -2 fieldset = FieldSet.from_data( @@ -21,7 +21,7 @@ def test_simulate_ship_st() -> None: sample_points = [Spacetime(Location(0, 0), 0)] - simulate_ship_st( + simulate_ship_underwater_st( fieldset=fieldset, out_file_name="test", depth=DEPTH, diff --git a/tests/test_sailship.py b/tests/test_sailship.py index 70b3d0eb..6307af92 100644 --- a/tests/test_sailship.py +++ b/tests/test_sailship.py @@ -13,7 +13,7 @@ def test_sailship() -> None: {"lon": 0, "lat": 0}, ) - ship_st_fieldset = FieldSet.from_data( + ship_underwater_st_fieldset = FieldSet.from_data( {"U": 0, "V": 0, "S": 0, "T": 0}, {"lon": 0, "lat": 0}, ) @@ -44,7 +44,7 @@ def test_sailship() -> None: config = VirtualShipConfiguration( "sailship_config.json", adcp_fieldset=adcp_fieldset, - ship_st_fieldset=ship_st_fieldset, + ship_underwater_st_fieldset=ship_underwater_st_fieldset, ctd_fieldset=ctd_fieldset, drifter_fieldset=drifter_fieldset, argo_float_fieldset=argo_float_fieldset, diff --git a/virtual_ship/instruments/__init__.py b/virtual_ship/instruments/__init__.py index 9ca81420..edcaa8c1 100644 --- a/virtual_ship/instruments/__init__.py +++ b/virtual_ship/instruments/__init__.py @@ -1,5 +1,5 @@ """Measurement instrument that can be used with Parcels.""" -from . import adcp, argo_float, ctd, drifter, ship_st +from . import adcp, argo_float, ctd, drifter, ship_underwater_st -__all__ = ["adcp", "argo_float", "ctd", "drifter", "ship_st"] +__all__ = ["adcp", "argo_float", "ctd", "drifter", "ship_underwater_st"] diff --git a/virtual_ship/instruments/ship_st.py b/virtual_ship/instruments/ship_underwater_st.py similarity index 98% rename from virtual_ship/instruments/ship_st.py rename to virtual_ship/instruments/ship_underwater_st.py index 1932a3a4..4c8eade5 100644 --- a/virtual_ship/instruments/ship_st.py +++ b/virtual_ship/instruments/ship_underwater_st.py @@ -23,7 +23,7 @@ def _sample_temperature(particle, fieldset, time): particle.temperature = fieldset.T[time, particle.depth, particle.lat, particle.lon] -def simulate_ship_st( +def simulate_ship_underwater_st( fieldset: FieldSet, out_file_name: str, depth: float, diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index cbc36372..781af2ca 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -12,7 +12,7 @@ from .instruments.argo_float import ArgoFloat, simulate_argo_floats from .instruments.ctd import CTD, simulate_ctd from .instruments.drifter import Drifter, simulate_drifters -from .instruments.ship_st import simulate_ship_st +from .instruments.ship_underwater_st import simulate_ship_underwater_st from .location import Location from .postprocess import postprocess from .spacetime import Spacetime @@ -63,7 +63,7 @@ def sailship(config: VirtualShipConfiguration): adcps: list[Spacetime] = [] # ship st objects to be used in simulation - ship_sts: list[Spacetime] = [] + ship_underwater_sts: list[Spacetime] = [] # argo float deployment locations that have been visited argo_locations_visited: set[Location] = set() @@ -193,11 +193,11 @@ def sailship(config: VirtualShipConfiguration): ) print("Simulating onboard salinity and temperature measurements.") - simulate_ship_st( - fieldset=config.ship_st_fieldset, - out_file_name=os.path.join("results", "ship_st.zarr"), + simulate_ship_underwater_st( + fieldset=config.ship_underwater_st_fieldset, + out_file_name=os.path.join("results", "ship_underwater_st.zarr"), depth=-2, - sample_points=ship_sts, + sample_points=ship_underwater_sts, ) print("Simulating onboard ADCP.") diff --git a/virtual_ship/virtual_ship_configuration.py b/virtual_ship/virtual_ship_configuration.py index b3300d53..ca2b8547 100644 --- a/virtual_ship/virtual_ship_configuration.py +++ b/virtual_ship/virtual_ship_configuration.py @@ -12,7 +12,7 @@ class VirtualShipConfiguration: """Configuration of the virtual ship, initialized from a json file.""" adcp_fieldset: FieldSet - ship_st_fieldset: FieldSet + ship_underwater_st_fieldset: FieldSet argo_float_fieldset: FieldSet drifter_fieldset: FieldSet ctd_fieldset: FieldSet @@ -23,7 +23,7 @@ def __init__( self, json_file, adcp_fieldset: FieldSet, - ship_st_fieldset: FieldSet, + ship_underwater_st_fieldset: FieldSet, ctd_fieldset: FieldSet, drifter_fieldset: FieldSet, argo_float_fieldset: FieldSet, @@ -33,14 +33,14 @@ def __init__( :param json_file: Path to the JSON file to init from. :param adcp_fieldset: Fieldset for ADCP measurements. - :param ship_st_fieldset: Fieldset for ship salinity temperature measurements. + :param ship_underwater_st_fieldset: Fieldset for ship salinity temperature measurements. :param ctd_fieldset: Fieldset for CTD measurements. :param drifter_fieldset: Fieldset for CTD measurements. :param argo_float_fieldset: Fieldset for argo float measurements. :raises ValueError: If JSON file not valid. """ self.adcp_fieldset = adcp_fieldset - self.ship_st_fieldset = ship_st_fieldset + self.ship_underwater_st_fieldset = ship_underwater_st_fieldset self.ctd_fieldset = ctd_fieldset self.drifter_fieldset = drifter_fieldset self.argo_float_fieldset = argo_float_fieldset From 23a0b312af271f2fab60554d5e6a0d887839cbc7 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 16:41:56 +0200 Subject: [PATCH 37/40] remove sailship create fieldset commented out function --- virtual_ship/sailship.py | 69 ---------------------------------------- 1 file changed, 69 deletions(-) diff --git a/virtual_ship/sailship.py b/virtual_ship/sailship.py index 781af2ca..62bf5767 100644 --- a/virtual_ship/sailship.py +++ b/virtual_ship/sailship.py @@ -244,75 +244,6 @@ def sailship(config: VirtualShipConfiguration): print(f"This cruise took {time_past} and would have cost {cost:,.0f} euros.") -# def create_fieldset(config, data_dir: str): -# """ -# Create fieldset from netcdf files and adds bathymetry data for CTD cast, returns fieldset with negative depth values. - -# :param config: The cruise configuration. -# :param data_dir: TODO -# :returns: The fieldset. -# :raises ValueError: If downloaded data is not as expected. -# """ -# filenames = { -# "U": os.path.join(data_dir, "studentdata_UV.nc"), -# "V": os.path.join(data_dir, "studentdata_UV.nc"), -# "S": os.path.join(data_dir, "studentdata_S.nc"), -# "T": os.path.join(data_dir, "studentdata_T.nc"), -# } -# variables = {"U": "uo", "V": "vo", "S": "so", "T": "thetao"} -# dimensions = { -# "lon": "longitude", -# "lat": "latitude", -# "time": "time", -# "depth": "depth", -# } - -# # create the fieldset and set interpolation methods -# fieldset = FieldSet.from_netcdf( -# filenames, variables, dimensions, allow_time_extrapolation=True -# ) -# fieldset.T.interp_method = "linear_invdist_land_tracer" -# fieldset.S.interp_method = "linear_invdist_land_tracer" -# for g in fieldset.gridset.grids: -# if max(g.depth) > 0: -# g.depth = -g.depth # make depth negative -# fieldset.mindepth = -fieldset.U.depth[0] # uppermost layer in the hydrodynamic data -# if config.CTD_settings["max_depth"] == "max": -# fieldset.add_constant("max_depth", -fieldset.U.depth[-1]) -# else: -# fieldset.add_constant("max_depth", config.CTD_settings["max_depth"]) -# fieldset.add_constant("maxtime", fieldset.U.grid.time_full[-1]) - -# # add bathymetry data to the fieldset for CTD cast -# bathymetry_file = os.path.join(data_dir, "GLO-MFC_001_024_mask_bathy.nc") -# bathymetry_variables = ("bathymetry", "deptho") -# bathymetry_dimensions = {"lon": "longitude", "lat": "latitude"} -# bathymetry_field = Field.from_netcdf( -# bathymetry_file, bathymetry_variables, bathymetry_dimensions -# ) -# fieldset.add_field(bathymetry_field) -# # read in data already -# fieldset.computeTimeChunk(0, 1) - -# if fieldset.U.lon.min() > config.region_of_interest["West"]: -# raise ValueError( -# "FieldSet western boundary is outside region of interest. Please run download_data.py again." -# ) -# if fieldset.U.lon.max() < config.region_of_interest["East"]: -# raise ValueError( -# "FieldSet eastern boundary is outside region of interest. Please run download_data.py again." -# ) -# if fieldset.U.lat.min() > config.region_of_interest["South"]: -# raise ValueError( -# "FieldSet southern boundary is outside region of interest. Please run download_data.py again." -# ) -# if fieldset.U.lat.max() < config.region_of_interest["North"]: -# raise ValueError( -# "FieldSet northern boundary is outside region of interest. Please run download_data.py again." -# ) -# return fieldset - - def shiproute(config: VirtualShipConfiguration, dt: timedelta) -> list[Location]: """ Take in route coordinates and return lat and lon points within region of interest to sample. From 06092b00abbdc2fd0492220772b9496f1f14ca33 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 16:48:27 +0200 Subject: [PATCH 38/40] ship underway st docstring improvement --- virtual_ship/instruments/ship_underwater_st.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/virtual_ship/instruments/ship_underwater_st.py b/virtual_ship/instruments/ship_underwater_st.py index 4c8eade5..53545e45 100644 --- a/virtual_ship/instruments/ship_underwater_st.py +++ b/virtual_ship/instruments/ship_underwater_st.py @@ -30,7 +30,7 @@ def simulate_ship_underwater_st( sample_points: list[Spacetime], ) -> None: """ - Use parcels to simulate the measurement of salinity and temperature aboard a ship in a fieldset. + Use parcels to simulate underway data, measuring salinity and temperature at the given depth in a fieldset. :param fieldset: The fieldset to simulate the sampling in. :param out_file_name: The file to write the results to. From ced9221ed46b1c2b026d68f383811466b273c054 Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 17:09:53 +0200 Subject: [PATCH 39/40] minor docstring change --- virtual_ship/instruments/ship_underwater_st.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/virtual_ship/instruments/ship_underwater_st.py b/virtual_ship/instruments/ship_underwater_st.py index 53545e45..b093b62f 100644 --- a/virtual_ship/instruments/ship_underwater_st.py +++ b/virtual_ship/instruments/ship_underwater_st.py @@ -30,7 +30,7 @@ def simulate_ship_underwater_st( sample_points: list[Spacetime], ) -> None: """ - Use parcels to simulate underway data, measuring salinity and temperature at the given depth in a fieldset. + Use parcels to simulate underway data, measuring salinity and temperature at the given depth as the ship moves in a fieldset. :param fieldset: The fieldset to simulate the sampling in. :param out_file_name: The file to write the results to. From e50ca96060db15c5417a842eb8287082ca83cc4d Mon Sep 17 00:00:00 2001 From: Aart Stuurman Date: Wed, 5 Jun 2024 17:10:33 +0200 Subject: [PATCH 40/40] minor docstring change --- virtual_ship/instruments/ship_underwater_st.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/virtual_ship/instruments/ship_underwater_st.py b/virtual_ship/instruments/ship_underwater_st.py index b093b62f..389a1e8e 100644 --- a/virtual_ship/instruments/ship_underwater_st.py +++ b/virtual_ship/instruments/ship_underwater_st.py @@ -30,7 +30,7 @@ def simulate_ship_underwater_st( sample_points: list[Spacetime], ) -> None: """ - Use parcels to simulate underway data, measuring salinity and temperature at the given depth as the ship moves in a fieldset. + Use parcels to simulate underway data, measuring salinity and temperature at the given depth along the ship track in a fieldset. :param fieldset: The fieldset to simulate the sampling in. :param out_file_name: The file to write the results to.