diff --git a/providers/ESR.json b/providers/ESR.json index 30cbb32..6510fa5 100644 --- a/providers/ESR.json +++ b/providers/ESR.json @@ -7,7 +7,16 @@ "u": "aodtm5_tmd/UV0_Arc5km" }, "name": "AODTM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aodtm-5/", "type": [ "u", @@ -21,7 +30,16 @@ "u": "aotim5_tmd/UV_Arc5km" }, "name": "AOTIM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": [ "u", @@ -35,7 +53,16 @@ "u": "Arc5km2018/UV_Arc5km2018" }, "name": "AOTIM-5-2018", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": [ "u", @@ -50,7 +77,17 @@ "u": "Arc2kmTM/UV_Arc2kmTM_v1" }, "name": "Arc2kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2D21RK6K", "type": [ "u", @@ -65,7 +102,7 @@ "u": "cats0201_tmd/UV0_CATS02_01" }, "name": "CATS0201", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://mail.esr.org/polar_tide_models/Model_CATS0201.html", "type": [ "u", @@ -79,7 +116,17 @@ "u": "CATS2008/uv.CATS2008.out" }, "name": "CATS2008", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601235", "type": [ "u", @@ -93,7 +140,17 @@ "u": "greenlandTMD_v2/u_Greenland8_rot.v2" }, "name": "Gr1km-v2", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.1002/2016RG000546", "type": [ "u", @@ -108,7 +165,17 @@ "u": "Gr1kmTM/UV_Gr1kmTM_v1" }, "name": "Gr1kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2B853K18", "type": [ "u", @@ -123,7 +190,16 @@ "grid_file": "aodtm5_tmd/grid_Arc5km", "model_file": "aodtm5_tmd/h0_Arc5km.oce", "name": "AODTM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aodtm-5/", "type": "z", "variable": "tide_ocean" @@ -133,7 +209,16 @@ "grid_file": "aotim5_tmd/grid_Arc5km", "model_file": "aotim5_tmd/h_Arc5km.oce", "name": "AOTIM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": "z", "variable": "tide_ocean" @@ -143,7 +228,16 @@ "grid_file": "Arc5km2018/grid_Arc5km2018", "model_file": "Arc5km2018/h_Arc5km2018", "name": "AOTIM-5-2018", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": "z", "variable": "tide_ocean", @@ -154,7 +248,17 @@ "grid_file": "Arc2kmTM/grid_Arc2kmTM_v1", "model_file": "Arc2kmTM/h_Arc2kmTM_v1", "name": "Arc2kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2D21RK6K", "type": "z", "variable": "tide_ocean", @@ -165,7 +269,7 @@ "grid_file": "cats0201_tmd/grid_CATS", "model_file": "cats0201_tmd/h0_CATS02_01", "name": "CATS0201", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://mail.esr.org/polar_tide_models/Model_CATS0201.html", "type": "z", "variable": "tide_ocean" @@ -175,7 +279,17 @@ "grid_file": "CATS2008/grid_CATS2008", "model_file": "CATS2008/hf.CATS2008.out", "name": "CATS2008", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601235", "type": "z", "variable": "tide_ocean" @@ -185,7 +299,17 @@ "grid_file": "CATS2008a_SPOTL_Load/grid_CATS2008a_opt", "model_file": "CATS2008a_SPOTL_Load/h_CATS2008a_SPOTL_load", "name": "CATS2008_load", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601235", "type": "z", "variable": "tide_load" @@ -195,7 +319,17 @@ "grid_file": "greenlandTMD_v2/grid_Greenland8.v2", "model_file": "greenlandTMD_v2/h_Greenland8.v2", "name": "Gr1km-v2", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.1002/2016RG000546", "type": "z", "variable": "tide_ocean", @@ -206,11 +340,21 @@ "grid_file": "Gr1kmTM/grid_Gr1kmTM_v1", "model_file": "Gr1kmTM/h_Gr1kmTM_v1", "name": "Gr1kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2B853K18", "type": "z", "variable": "tide_ocean", "version": "v1" } } -} \ No newline at end of file +} diff --git a/providers/TPXO.json b/providers/TPXO.json index 1c8b4a9..1360e9f 100644 --- a/providers/TPXO.json +++ b/providers/TPXO.json @@ -23,7 +23,7 @@ ] }, "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.0001, "type": [ @@ -72,7 +72,7 @@ ] }, "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.0001, "type": [ @@ -88,7 +88,7 @@ "u": "TPXO7.2_tmd/u_tpxo7.2" }, "name": "TPXO7.2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": [ "u", @@ -103,7 +103,7 @@ "u": "tpxo8_atlas/uv.tpxo8_atlas_30_v1" }, "name": "TPXO8-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo8_atlas.html", "type": [ "u", @@ -131,7 +131,7 @@ ] }, "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.0001, "type": [ @@ -174,7 +174,7 @@ ] }, "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.0001, "type": [ @@ -203,7 +203,7 @@ ] }, "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -246,7 +246,7 @@ ] }, "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -277,7 +277,7 @@ ] }, "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -324,7 +324,7 @@ ] }, "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -355,7 +355,7 @@ ] }, "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -402,7 +402,7 @@ ] }, "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -434,7 +434,7 @@ ] }, "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -483,7 +483,7 @@ ] }, "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -499,7 +499,7 @@ "u": "TPXO9.1/u_tpxo9.v1" }, "name": "TPXO9.1", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": [ "u", @@ -530,7 +530,7 @@ "TPXO10_atlas_v2/h_s2_tpxo10_atlas_30_v2" ], "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.001, "type": "z", @@ -558,7 +558,7 @@ "TPXO10_atlas_v2/h_s2_tpxo10_atlas_30_v2.nc" ], "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.001, "type": "z", @@ -570,7 +570,7 @@ "grid_file": "TPXO7.2_tmd/grid_tpxo7.2", "model_file": "TPXO7.2_tmd/h_tpxo7.2", "name": "TPXO7.2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": "z", "variable": "tide_ocean", @@ -581,7 +581,7 @@ "grid_file": "TPXO7.2_load/grid_tpxo6.2", "model_file": "TPXO7.2_load/h_tpxo7.2_load", "name": "TPXO7.2_load", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": "z", "variable": "tide_load", @@ -592,7 +592,7 @@ "grid_file": "tpxo8_atlas/grid_tpxo8atlas_30_v1", "model_file": "tpxo8_atlas/hf.tpxo8_atlas_30_v1", "name": "TPXO8-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo8_atlas.html", "type": "z", "variable": "tide_ocean", @@ -600,7 +600,7 @@ }, "TPXO8-atlas-nc": { "format": "ATLAS-netcdf", - "name": "TPXO8-atlas-nc", + "grid_file": "TPXO8_atlas_v1/grid_tpxo8atlas_30_v1.nc", "model_file": [ "TPXO8_atlas_v1/hf.k1_tpxo8_atlas_30c_v1.nc", "TPXO8_atlas_v1/hf.n2_tpxo8_atlas_30c_v1.nc", @@ -612,12 +612,12 @@ "TPXO8_atlas_v1/hf.m4_tpxo8_atlas_30c_v1.nc", "TPXO8_atlas_v1/hf.m2_tpxo8_atlas_30c_v1.nc" ], - "grid_file": "TPXO8_atlas_v1/grid_tpxo8atlas_30_v1.nc", + "name": "TPXO8-atlas-nc", + "reference": "http://volkov.oce.orst.edu/tides/tpxo8_atlas.html", + "scale": 0.001, "type": "z", - "version": "v1", "variable": "tide_ocean", - "scale": 0.001, - "reference": "http://volkov.oce.orst.edu/tides/tpxo8_atlas.html" + "version": "v1" }, "TPXO9-atlas": { "format": "OTIS", @@ -637,7 +637,7 @@ "TPXO9_atlas/h_2n2_tpxo9_atlas_30" ], "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.001, "type": "z", @@ -662,7 +662,7 @@ "TPXO9_atlas/h_2n2_tpxo9_atlas_30.nc" ], "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.001, "type": "z", @@ -687,7 +687,7 @@ "TPXO9_atlas_v2/h_2n2_tpxo9_atlas_30_v2" ], "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -712,7 +712,7 @@ "TPXO9_atlas_v2/h_2n2_tpxo9_atlas_30_v2.nc" ], "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -739,7 +739,7 @@ "TPXO9_atlas_v3/h_mm_tpxo9_atlas_30_v3" ], "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -766,7 +766,7 @@ "TPXO9_atlas_v3/h_mm_tpxo9_atlas_30_v3.nc" ], "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -793,7 +793,7 @@ "TPXO9_atlas_v4/h_mm_tpxo9_atlas_30_v4" ], "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -820,7 +820,7 @@ "TPXO9_atlas_v4/h_mm_tpxo9_atlas_30_v4.nc" ], "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -848,7 +848,7 @@ "TPXO9_atlas_v5/h_mm_tpxo9_atlas_30_v5" ], "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -876,7 +876,7 @@ "TPXO9_atlas_v5/h_mm_tpxo9_atlas_30_v5.nc" ], "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -888,7 +888,7 @@ "grid_file": "TPXO9.1/DATA/grid_tpxo9", "model_file": "TPXO9.1/DATA/h_tpxo9.v1", "name": "TPXO9.1", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": "z", "variable": "tide_ocean", diff --git a/providers/_update_providers.py b/providers/_update_providers.py new file mode 100644 index 0000000..b6fb3dc --- /dev/null +++ b/providers/_update_providers.py @@ -0,0 +1,159 @@ +""" +_update_providers.py (08/2024) +Update the projection variables in the providers +""" +import json +import pyproj +import inspect +import pathlib +import argparse +import numpy as np +import pyTMD.utilities + +# current file path +filename = inspect.getframeinfo(inspect.currentframe()).filename +filepath = pathlib.Path(filename).absolute().parent + +# PURPOSE: create argument parser +def arguments(): + parser = argparse.ArgumentParser( + description="""Update the projection variables in the providers" + """, + fromfile_prefix_chars="@" + ) + # command line parameters + parser.add_argument('--pretty', '-p', + action='store_true', + help='Pretty print the json file') + parser.add_argument('--verbose', '-v', + action='store_true', + help='Verbose output') + return parser + +def get_projection(PROJ): + # previously defined coordinate reference systems + # python dictionary with named conversion functions + CRS = {} + CRS['3031'] = _EPSG3031 + CRS['3413'] = _EPSG3413 + CRS['CATS2008'] = _CATS2008 + CRS['3976'] = _EPSG3976 + CRS['AEDNorth'] = _AEDNorth + CRS['4326'] = _EPSG4326 + # check that PROJ for conversion was entered correctly + # run named conversion program and return values + try: + crs = CRS[PROJ].__call__() + except KeyError as exc: + pass + else: + # return the output variables + return crs + +def _EPSG3031(): + """ + CRS for models in EPSG:3031 (Antarctic Polar Stereographic) + """ + # coordinate reference system information + crs = pyproj.CRS.from_user_input({'proj':'stere', + 'lat_0':-90, 'lat_ts':-71, 'lon_0':0, 'x_0':0., 'y_0':0., + 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) + return crs + +# wrapper function for models in EPSG 3413 (Sea Ice Polar Stereographic North) +def _EPSG3413(): + """ + CRS for models in EPSG:3413 (Sea Ice Polar Stereographic North) + """ + # coordinate reference system information + crs = pyproj.CRS.from_user_input({'proj':'stere', + 'lat_0':90, 'lat_ts':70, 'lon_0':-45, 'x_0':0., 'y_0':0., + 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) + return crs + +# wrapper function for CATS2008 tide models +def _CATS2008(): + """ + CRS for Circum-Antarctic Tidal Simulation models + """ + # coordinate reference system information + crs = pyproj.CRS.from_user_input({'proj':'stere', + 'lat_0':-90, 'lat_ts':-71, 'lon_0':-70, 'x_0':0., 'y_0':0., + 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) + return crs + +# wrapper function for models in EPSG 3976 (NSIDC Sea Ice Stereographic South) +def _EPSG3976(): + """ + CRS for models in EPSG:3976 (Sea Ice Polar Stereographic South) + """ + # coordinate reference system information + crs = pyproj.CRS.from_user_input({'proj':'stere', + 'lat_0':-90, 'lat_ts':-70, 'lon_0':0, 'x_0':0., 'y_0':0., + 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) + return crs + +# function for models in (idealized) Azimuthal Equidistant projection +def _AEDNorth(): + """ + CRS for models in idealized Azimuthal Equidistant projections + """ + # projections for converting to and from input EPSG + R = 111700.0*180.0/np.pi + crs = pyproj.CRS.from_user_input({'proj':'aeqd','lat_0':90, + 'lon_0':270,'x_0':0.,'y_0':0.,'ellps':'sphere', + 'R':R,'units':'km'}) + return crs + +# wrapper function to pass lat/lon values or convert if EPSG +def _EPSG4326(): + """ + CRS for models in EPSG:4326 (WGS84 Latitude/Longitude) + """ + crs = pyproj.CRS.from_epsg(4326) + return crs + +# wrapper function for using custom projections +def _custom(self, PROJ: int | str): + """ + CRS for models in a custom projection + + Parameters + ---------- + PROJ: int or str + Spatial reference system code for coordinate transformations + """ + # coordinate reference system information + crs = self.from_input(PROJ) + return crs + +def main(): + # Read the system arguments listed after the program + parser = arguments() + args,_ = parser.parse_known_args() + + # find providers + providers = [f for f in filepath.iterdir() if (f.suffix == '.json')] + # for each provider + for provider in providers: + with provider.open('r', encoding='utf-8') as fid: + p = json.load(fid) + for key, val in p.items(): + # update projections + for model, parameters in val.items(): + if 'projection' in parameters: + crs = get_projection(parameters['projection']) + if crs.crs.to_epsg() is not None: + p[key][model]['projection'] = crs.crs.to_string() + else: + d = crs.crs.to_dict() + d.pop('no_defs') + p[key][model]['projection'] = d.copy() + + # writing model parameters to JSON database file + with provider.open('w', encoding='utf-8') as fid: + indent = 4 if args.pretty else None + json.dump(p, fid, indent=indent, sort_keys=True) + +if __name__ == '__main__': + main() \ No newline at end of file diff --git a/providers/providers.json b/providers/providers.json index 326fa3d..84a9fdb 100644 --- a/providers/providers.json +++ b/providers/providers.json @@ -7,7 +7,17 @@ "u": "CATS2008_v2023/CATS2008_v2023.nc" }, "name": "CATS2008-v2023", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601772", "type": [ "u", @@ -56,7 +66,17 @@ "grid_file": "CATS2008_v2023/CATS2008_v2023.nc", "model_file": "CATS2008_v2023/CATS2008_v2023.nc", "name": "CATS2008-v2023", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601772", "type": "z", "variable": "tide_ocean" @@ -138,4 +158,4 @@ "version": "HAMTIDE11" } } -} \ No newline at end of file +} diff --git a/pyTMD/crs.py b/pyTMD/crs.py index 9ab549e..3e3ba55 100644 --- a/pyTMD/crs.py +++ b/pyTMD/crs.py @@ -11,7 +11,7 @@ INPUTS: i1: longitude ('F') or projection easting x ('B') i2: latitude ('F') or projection northing y ('B') - PROJ: spatial reference system code for coordinate transformations + PROJ: spatial reference system for coordinate transformations BF: backwards ('B') or forward ('F') translations OPTIONS: @@ -31,6 +31,7 @@ UPDATE HISTORY: Updated 09/2024: added function for idealized Arctic Azimuthal projection + complete refactor to use JSON dictionary format for model projections Updated 07/2024: added function to get the CRS transform Updated 05/2024: make subscriptable and allow item assignment Updated 04/2024: use wrapper to importlib for optional dependencies @@ -83,7 +84,7 @@ def __init__(self): def convert(self, i1: np.ndarray, i2: np.ndarray, - PROJ: str, + PROJ: str | dict, BF: str, EPSG: int | str = 4326 ): @@ -96,8 +97,8 @@ def convert(self, Input x-coordinates i2: np.ndarray Input y-coordinates - PROJ: str - Spatial reference system code for coordinate transformations + PROJ: str or dict + Spatial reference system for coordinate transformations BF: str Direction of transformation @@ -115,65 +116,37 @@ def convert(self, """ # name of the projection self.name = PROJ - # set the transform and transform direction - self.get(PROJ, EPSG=EPSG) + # get the CRS and transform direction + self.get(PROJ) self._direction = BF[0].upper() # run conversion program and return values - return self.transform(i1, i2) + return self.transform(i1, i2, EPSG=EPSG) # PURPOSE: try to get the projection information - def get(self, PROJ: str, EPSG: int | str = 4326): + def get(self, PROJ: str | dict): """ - Tries to get the CRS transformer for given values + Tries to get the coordinate reference system Parameters ---------- - PROJ: str - Spatial reference system code for coordinate transformations - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - - Returns - ------- - o1: np.ndarray - Output transformed x-coordinates - o2: np.ndarray - Output transformed y-coordinates + PROJ: str or dict + Spatial reference system for coordinate transformations """ - # name of the projection - self.name = PROJ - # python dictionary with named conversion functions - transforms = {} - transforms['3031'] = self._EPSG3031 - transforms['3413'] = self._EPSG3413 - transforms['CATS2008'] = self._CATS2008 - transforms['3976'] = self._EPSG3976 - transforms['PSNorth'] = self._PSNorth - transforms['AEDNorth'] = self._AEDNorth - transforms['4326'] = self._EPSG4326 - # check that PROJ for conversion was entered correctly - # run named conversion program and return values + # get the coordinate reference system try: - transforms[PROJ](EPSG) - except KeyError as exc: - pass - else: - # return the output variables - return self - # try changing the projection using a custom projection - # run custom conversion program and return values - try: - self._custom(PROJ, EPSG=EPSG) + self.crs = self.from_input(PROJ) + self.name = self.crs.name except Exception as exc: pass else: return self # projection not found or available - raise Exception(f'PROJ: {PROJ} conversion function not found') + raise pyproj.exceptions.CRSError def transform(self, i1: np.ndarray, i2: np.ndarray, + EPSG: int | str = 4326, **kwargs): """ Performs Coordinates Reference System (CRS) transformations @@ -184,6 +157,8 @@ def transform(self, Input x-coordinates i2: np.ndarray Input y-coordinates + EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) + input (``'F'``) or output (``'B'``) coordinate system kwargs: dict Keyword arguments for the transformation @@ -196,38 +171,32 @@ def transform(self, """ # set the direction of the transformation kwargs.setdefault('direction', self.direction) - if (self.name == 'PSNorth') and (self.direction.name == 'FORWARD'): - # convert input coordinate reference system to lat/lon - lon, lat = self.transformer.transform(i1, i2, **kwargs) - # convert lat/lon to (idealized) Polar-Stereographic x/y - o1 = (90.0 - lat)*111.7*np.cos(lon/180.0*np.pi) - o2 = (90.0 - lat)*111.7*np.sin(lon/180.0*np.pi) - elif (self.name == 'PSNorth') and (self.direction.name == 'INVERSE'): - # convert (idealized) Polar-Stereographic x/y to lat/lon - lon = np.arctan2(i2, i1)*180.0/np.pi - lat = 90.0 - np.sqrt(i1**2 + i2**2)/111.7 - # adjust longitudes to be -180:180 - ii, = np.nonzero(lon < 0) - lon[ii] += 360.0 - # convert to output coordinate reference system - o1, o2 = self.transformer.transform(lon, lat, **kwargs) - else: - # convert coordinate reference system - o1, o2 = self.transformer.transform(i1, i2, **kwargs) + # get the coordinate reference system and transform + source_crs = self.from_input(EPSG) + self.transformer = pyproj.Transformer.from_crs( + source_crs, self.crs, always_xy=True) + # convert coordinate reference system + o1, o2 = self.transformer.transform(i1, i2, **kwargs) # return the transformed coordinates return (o1, o2) # PURPOSE: try to get the projection information - def from_input(self, PROJECTION: int | str): + def from_input(self, PROJECTION: int | str | dict): """ Attempt to retrieve the Coordinate Reference System - for an input code Parameters ---------- - PROJECTION: int or str - Coordinate Reference System code + PROJECTION: int, str or dict + Coordinate Reference System """ + # coordinate reference system dictoinary + try: + CRS = pyproj.CRS.from_user_input(PROJECTION) + except (ValueError, pyproj.exceptions.CRSError): + pass + else: + return CRS # EPSG projection code try: CRS = pyproj.CRS.from_epsg(int(PROJECTION)) @@ -245,153 +214,6 @@ def from_input(self, PROJECTION: int | str): # no projection can be made raise pyproj.exceptions.CRSError - def _EPSG3031(self, EPSG: int | str = 4326): - """ - Transform for models in EPSG:3031 (Antarctic Polar Stereographic) - - Parameters - ---------- - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - # coordinate reference system information - crs1 = self.from_input(EPSG) - crs2 = pyproj.CRS.from_user_input({'proj':'stere', - 'lat_0':-90, 'lat_ts':-71, 'lon_0':0, 'x_0':0., 'y_0':0., - 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - - # wrapper function for models in EPSG 3413 (Sea Ice Polar Stereographic North) - def _EPSG3413(self, EPSG: int | str = 4326): - """ - Transform for models in EPSG:3413 (Sea Ice Polar Stereographic North) - - Parameters - ---------- - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - # coordinate reference system information - crs1 = self.from_input(EPSG) - crs2 = pyproj.CRS.from_user_input({'proj':'stere', - 'lat_0':90, 'lat_ts':70, 'lon_0':-45, 'x_0':0., 'y_0':0., - 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - - # wrapper function for CATS2008 tide models - def _CATS2008(self, EPSG: int | str = 4326): - """ - Transform for Circum-Antarctic Tidal Simulation models - - Parameters - ---------- - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - # coordinate reference system information - crs1 = self.from_input(EPSG) - crs2 = pyproj.CRS.from_user_input({'proj':'stere', - 'lat_0':-90, 'lat_ts':-71, 'lon_0':-70, 'x_0':0., 'y_0':0., - 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - - # wrapper function for models in EPSG 3976 (NSIDC Sea Ice Stereographic South) - def _EPSG3976(self, EPSG: int | str = 4326): - """ - Transform for models in EPSG:3976 (Sea Ice Polar Stereographic South) - - Parameters - ---------- - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - # coordinate reference system information - crs1 = self.from_input(EPSG) - crs2 = pyproj.CRS.from_user_input({'proj':'stere', - 'lat_0':-90, 'lat_ts':-70, 'lon_0':0, 'x_0':0., 'y_0':0., - 'ellps':'WGS84', 'datum':'WGS84', 'units':'km'}) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - - # function for models in (idealized) PSNorth projection - def _PSNorth(self, EPSG: int | str = 4326): - """ - Transform for idealized Arctic Polar Stereographic models - - Parameters - ---------- - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - # projections for converting to and from input EPSG - crs1 = self.from_input(EPSG) - crs2 = pyproj.CRS.from_epsg(4326) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - - # function for models in (idealized) Azimuthal Equidistant projection - def _AEDNorth(self, EPSG: int | str = 4326): - """ - Transform for models in idealized Azimuthal Equidistant projections - - Parameters - ---------- - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - # projections for converting to and from input EPSG - crs1 = self.from_input(EPSG) - R = 111700.0*180.0/np.pi - crs2 = pyproj.CRS.from_user_input({'proj':'aeqd','lat_0':90, - 'lon_0':270,'x_0':0.,'y_0':0.,'ellps':'sphere', - 'R':R,'units':'km'}) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - - # wrapper function to pass lat/lon values or convert if EPSG - def _EPSG4326(self, EPSG: int | str = 4326): - """ - Transform for models in EPSG:4326 (WGS84 Latitude/Longitude) - - Parameters - ---------- - EPSG: int, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - crs1 = self.from_input(EPSG) - crs2 = pyproj.CRS.from_epsg(4326) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - - # wrapper function for using custom projections - def _custom(self, PROJ: int | str, EPSG: int | str = 4326): - """ - Transform for models in a custom projection - - Parameters - ---------- - PROJ: int or str - Spatial reference system code for coordinate transformations - EPSG: int or str, default 4326 (WGS84 Latitude/Longitude) - input (``'F'``) or output (``'B'``) coordinate system - """ - # coordinate reference system information - crs1 = self.from_input(EPSG) - crs2 = self.from_input(PROJ) - self.transformer = pyproj.Transformer.from_crs(crs1, crs2, - always_xy=True) - return self - @property def direction(self): """ @@ -409,10 +231,7 @@ def is_geographic(self): """ Check if the coordinate reference system is geographic """ - if (self.name == 'PSNorth'): - return False - else: - return self.transformer.target_crs.is_geographic + return self.crs.is_geographic def __str__(self): """String representation of the ``crs`` object diff --git a/pyTMD/data/database.json b/pyTMD/data/database.json index 8c68a23..d195374 100644 --- a/pyTMD/data/database.json +++ b/pyTMD/data/database.json @@ -7,7 +7,16 @@ "u": "aodtm5_tmd/UV0_Arc5km" }, "name": "AODTM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aodtm-5/", "type": [ "u", @@ -21,7 +30,16 @@ "u": "aotim5_tmd/UV_Arc5km" }, "name": "AOTIM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": [ "u", @@ -35,7 +53,16 @@ "u": "Arc5km2018/UV_Arc5km2018" }, "name": "AOTIM-5-2018", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": [ "u", @@ -50,7 +77,17 @@ "u": "Arc2kmTM/UV_Arc2kmTM_v1" }, "name": "Arc2kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2D21RK6K", "type": [ "u", @@ -65,7 +102,7 @@ "u": "cats0201_tmd/UV0_CATS02_01" }, "name": "CATS0201", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://mail.esr.org/polar_tide_models/Model_CATS0201.html", "type": [ "u", @@ -79,7 +116,17 @@ "u": "CATS2008/uv.CATS2008.out" }, "name": "CATS2008", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601235", "type": [ "u", @@ -93,7 +140,17 @@ "u": "CATS2008_v2023/CATS2008_v2023.nc" }, "name": "CATS2008-v2023", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601772", "type": [ "u", @@ -192,7 +249,17 @@ "u": "greenlandTMD_v2/u_Greenland8_rot.v2" }, "name": "Gr1km-v2", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.1002/2016RG000546", "type": [ "u", @@ -207,7 +274,17 @@ "u": "Gr1kmTM/UV_Gr1kmTM_v1" }, "name": "Gr1kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2B853K18", "type": [ "u", @@ -273,7 +350,7 @@ ] }, "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.0001, "type": [ @@ -322,7 +399,7 @@ ] }, "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.0001, "type": [ @@ -338,7 +415,7 @@ "u": "TPXO7.2_tmd/u_tpxo7.2" }, "name": "TPXO7.2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": [ "u", @@ -353,7 +430,7 @@ "u": "tpxo8_atlas/uv.tpxo8_atlas_30_v1" }, "name": "TPXO8-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo8_atlas.html", "type": [ "u", @@ -381,7 +458,7 @@ ] }, "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.0001, "type": [ @@ -424,7 +501,7 @@ ] }, "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.0001, "type": [ @@ -453,7 +530,7 @@ ] }, "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -496,7 +573,7 @@ ] }, "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -527,7 +604,7 @@ ] }, "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -574,7 +651,7 @@ ] }, "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -605,7 +682,7 @@ ] }, "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -652,7 +729,7 @@ ] }, "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -684,7 +761,7 @@ ] }, "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -733,7 +810,7 @@ ] }, "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.0001, "type": [ @@ -749,7 +826,7 @@ "u": "TPXO9.1/u_tpxo9.v1" }, "name": "TPXO9.1", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": [ "u", @@ -764,7 +841,16 @@ "grid_file": "aodtm5_tmd/grid_Arc5km", "model_file": "aodtm5_tmd/h0_Arc5km.oce", "name": "AODTM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aodtm-5/", "type": "z", "variable": "tide_ocean" @@ -774,7 +860,16 @@ "grid_file": "aotim5_tmd/grid_Arc5km", "model_file": "aotim5_tmd/h_Arc5km.oce", "name": "AOTIM-5", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": "z", "variable": "tide_ocean" @@ -784,7 +879,16 @@ "grid_file": "Arc5km2018/grid_Arc5km2018", "model_file": "Arc5km2018/h_Arc5km2018", "name": "AOTIM-5-2018", - "projection": "AEDNorth", + "projection": { + "R": 6399938.5716113, + "lat_0": 90, + "lon_0": 270, + "proj": "aeqd", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://www.esr.org/research/polar-tide-models/list-of-polar-tide-models/aotim-5/", "type": "z", "variable": "tide_ocean", @@ -795,7 +899,17 @@ "grid_file": "Arc2kmTM/grid_Arc2kmTM_v1", "model_file": "Arc2kmTM/h_Arc2kmTM_v1", "name": "Arc2kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2D21RK6K", "type": "z", "variable": "tide_ocean", @@ -806,7 +920,7 @@ "grid_file": "cats0201_tmd/grid_CATS", "model_file": "cats0201_tmd/h0_CATS02_01", "name": "CATS0201", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://mail.esr.org/polar_tide_models/Model_CATS0201.html", "type": "z", "variable": "tide_ocean" @@ -816,7 +930,17 @@ "grid_file": "CATS2008/grid_CATS2008", "model_file": "CATS2008/hf.CATS2008.out", "name": "CATS2008", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601235", "type": "z", "variable": "tide_ocean" @@ -826,7 +950,17 @@ "grid_file": "CATS2008_v2023/CATS2008_v2023.nc", "model_file": "CATS2008_v2023/CATS2008_v2023.nc", "name": "CATS2008-v2023", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601772", "type": "z", "variable": "tide_ocean" @@ -836,7 +970,17 @@ "grid_file": "CATS2008a_SPOTL_Load/grid_CATS2008a_opt", "model_file": "CATS2008a_SPOTL_Load/h_CATS2008a_SPOTL_load", "name": "CATS2008_load", - "projection": "CATS2008", + "projection": { + "datum": "WGS84", + "lat_0": -90, + "lat_ts": -71, + "lon_0": -70, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.15784/601235", "type": "z", "variable": "tide_load" @@ -1536,7 +1680,17 @@ "grid_file": "greenlandTMD_v2/grid_Greenland8.v2", "model_file": "greenlandTMD_v2/h_Greenland8.v2", "name": "Gr1km-v2", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.1002/2016RG000546", "type": "z", "variable": "tide_ocean", @@ -1547,7 +1701,17 @@ "grid_file": "Gr1kmTM/grid_Gr1kmTM_v1", "model_file": "Gr1kmTM/h_Gr1kmTM_v1", "name": "Gr1kmTM", - "projection": "3413", + "projection": { + "datum": "WGS84", + "lat_0": 90, + "lat_ts": 70, + "lon_0": -45, + "proj": "stere", + "type": "crs", + "units": "km", + "x_0": 0, + "y_0": 0 + }, "reference": "https://doi.org/10.18739/A2B853K18", "type": "z", "variable": "tide_ocean", @@ -1611,7 +1775,7 @@ "TPXO10_atlas_v2/h_s2_tpxo10_atlas_30_v2" ], "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.001, "type": "z", @@ -1639,7 +1803,7 @@ "TPXO10_atlas_v2/h_s2_tpxo10_atlas_30_v2.nc" ], "name": "TPXO10-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo10-atlas", "scale": 0.001, "type": "z", @@ -1651,7 +1815,7 @@ "grid_file": "TPXO7.2_tmd/grid_tpxo7.2", "model_file": "TPXO7.2_tmd/h_tpxo7.2", "name": "TPXO7.2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": "z", "variable": "tide_ocean", @@ -1662,7 +1826,7 @@ "grid_file": "TPXO7.2_load/grid_tpxo6.2", "model_file": "TPXO7.2_load/h_tpxo7.2_load", "name": "TPXO7.2_load", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": "z", "variable": "tide_load", @@ -1673,7 +1837,7 @@ "grid_file": "tpxo8_atlas/grid_tpxo8atlas_30_v1", "model_file": "tpxo8_atlas/hf.tpxo8_atlas_30_v1", "name": "TPXO8-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo8_atlas.html", "type": "z", "variable": "tide_ocean", @@ -1718,7 +1882,7 @@ "TPXO9_atlas/h_2n2_tpxo9_atlas_30" ], "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.001, "type": "z", @@ -1743,7 +1907,7 @@ "TPXO9_atlas/h_2n2_tpxo9_atlas_30.nc" ], "name": "TPXO9-atlas", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/tpxo9_atlas.html", "scale": 0.001, "type": "z", @@ -1768,7 +1932,7 @@ "TPXO9_atlas_v2/h_2n2_tpxo9_atlas_30_v2" ], "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1793,7 +1957,7 @@ "TPXO9_atlas_v2/h_2n2_tpxo9_atlas_30_v2.nc" ], "name": "TPXO9-atlas-v2", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1820,7 +1984,7 @@ "TPXO9_atlas_v3/h_mm_tpxo9_atlas_30_v3" ], "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1847,7 +2011,7 @@ "TPXO9_atlas_v3/h_mm_tpxo9_atlas_30_v3.nc" ], "name": "TPXO9-atlas-v3", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1874,7 +2038,7 @@ "TPXO9_atlas_v4/h_mm_tpxo9_atlas_30_v4" ], "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1901,7 +2065,7 @@ "TPXO9_atlas_v4/h_mm_tpxo9_atlas_30_v4.nc" ], "name": "TPXO9-atlas-v4", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1929,7 +2093,7 @@ "TPXO9_atlas_v5/h_mm_tpxo9_atlas_30_v5" ], "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1957,7 +2121,7 @@ "TPXO9_atlas_v5/h_mm_tpxo9_atlas_30_v5.nc" ], "name": "TPXO9-atlas-v5", - "projection": "4326", + "projection": "EPSG:4326", "reference": "https://www.tpxo.net/global/tpxo9-atlas", "scale": 0.001, "type": "z", @@ -1969,7 +2133,7 @@ "grid_file": "TPXO9.1/DATA/grid_tpxo9", "model_file": "TPXO9.1/DATA/h_tpxo9.v1", "name": "TPXO9.1", - "projection": "4326", + "projection": "EPSG:4326", "reference": "http://volkov.oce.orst.edu/tides/global.html", "type": "z", "variable": "tide_ocean", diff --git a/pyTMD/io/OTIS.py b/pyTMD/io/OTIS.py index 8b52033..cbe5280 100644 --- a/pyTMD/io/OTIS.py +++ b/pyTMD/io/OTIS.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" OTIS.py -Written by Tyler Sutterley (08/2024) +Written by Tyler Sutterley (09/2024) Reads files for a tidal model and makes initial calculations to run tide program Includes functions to extract tidal harmonic constants from OTIS tide models for @@ -17,7 +17,6 @@ ilat: latitude to interpolate grid_file: grid file for model model_file: model file containing each constituent - EPSG: projection of tide model data OPTIONS: type: tidal variable to run @@ -59,6 +58,7 @@ interpolate.py: interpolation routines for spatial data UPDATE HISTORY: + Updated 09/2024: using new JSON dictionary format for model projections Updated 08/2024: revert change and assume crop bounds are projected Updated 07/2024: added crop and bounds keywords for trimming model data convert the crs of bounds when cropping model data @@ -126,7 +126,6 @@ import struct import logging import pathlib -import warnings import numpy as np import scipy.interpolate import pyTMD.crs @@ -171,7 +170,7 @@ def extract_constants( ilat: np.ndarray, grid_file: str | pathlib.Path | None = None, model_file: str | pathlib.Path | list | None = None, - EPSG: str | int | None = None, + projection: dict | str | int | None = None, **kwargs ): """ @@ -191,7 +190,7 @@ def extract_constants( grid file for model model_file: str, pathlib.Path, list or NoneType, default None model file containing each constituent - EPSG: str or NoneType, default None, + projection: str or NoneType, default None, projection of tide model data type: str, default 'z' Tidal variable to read @@ -250,8 +249,8 @@ def extract_constants( EXTRAPOLATE='extrapolate',CUTOFF='cutoff',GRID='grid') for old,new in deprecated_keywords.items(): if old in kwargs.keys(): - warnings.warn(f"""Deprecated keyword argument {old}. - Changed to '{new}'""", DeprecationWarning) + logging.warning(f"""Deprecated keyword argument {old}. + Changed to '{new}'""") # set renamed argument to not break workflows kwargs[new] = copy.copy(kwargs[old]) @@ -279,9 +278,9 @@ def extract_constants( ilon = np.atleast_1d(np.copy(ilon)) ilat = np.atleast_1d(np.copy(ilat)) # run wrapper function to convert coordinate systems of input lat/lon - transformer = pyTMD.crs().get(EPSG) - x,y = transformer.transform(ilon, ilat, direction='FORWARD') - is_geographic = transformer.is_geographic + crs = pyTMD.crs().get(projection) + x,y = crs.transform(ilon, ilat, direction='FORWARD') + is_geographic = crs.is_geographic # grid step size of tide model dx = xi[1] - xi[0] dy = yi[1] - yi[0] @@ -515,7 +514,7 @@ def extract_constants( def read_constants( grid_file: str | pathlib.Path | None = None, model_file: str | pathlib.Path | list | None = None, - EPSG: str | int | None = None, + projection: dict | str | int | None = None, **kwargs ): """ @@ -527,7 +526,7 @@ def read_constants( grid file for model model_file: str, pathlib.Path, list or NoneType, default None model file containing each constituent - EPSG: str or NoneType, default None, + projection: str, dict or NoneType, default None, projection of tide model data type: str, default 'z' Tidal variable to read @@ -586,9 +585,9 @@ def read_constants( dy = yi[1] - yi[0] # run wrapper function to convert coordinate systems - transformer = pyTMD.crs().get(EPSG) + crs = pyTMD.crs().get(projection) # if global: extend limits - is_geographic = transformer.is_geographic + is_geographic = crs.is_geographic is_global = False # crop mask and bathymetry data to (buffered) bounds @@ -653,9 +652,9 @@ def read_constants( cons = [read_constituents(m)[0].pop() for m in model_file] else: cons,_ = read_constituents(model_file, grid=kwargs['grid']) - # save output constituents + # save output constituents and coordinate reference system constituents = pyTMD.io.constituents(x=xi, y=yi, - bathymetry=bathymetry.data, mask=mask) + bathymetry=bathymetry.data, mask=mask, crs=crs) # read each model constituent for i,c in enumerate(cons): @@ -720,7 +719,6 @@ def interpolate_constants( ilon: np.ndarray, ilat: np.ndarray, constituents, - EPSG: str | int | None = None, **kwargs ): """ @@ -737,8 +735,6 @@ def interpolate_constants( latitude to interpolate constituents: obj Tide model constituents (complex form) - EPSG: str or NoneType, default None, - projection of tide model data type: str, default 'z' Tidal variable to read @@ -782,10 +778,9 @@ def interpolate_constants( # adjust dimensions of input coordinates to be iterable ilon = np.atleast_1d(np.copy(ilon)) ilat = np.atleast_1d(np.copy(ilat)) - # run wrapper function to convert coordinate systems of input lat/lon - transformer = pyTMD.crs().get(EPSG) - x,y = transformer.transform(ilon, ilat) - is_geographic = transformer.is_geographic + # convert coordinate systems of input lat/lon + x,y = constituents.crs.transform(ilon, ilat) + is_geographic = constituents.crs.is_geographic # adjust longitudinal convention of input latitude and longitude # to fit tide model convention if (np.min(x) < np.min(xi)) & is_geographic: diff --git a/scripts/compute_tidal_currents.py b/scripts/compute_tidal_currents.py index bf0d40c..750dcf2 100755 --- a/scripts/compute_tidal_currents.py +++ b/scripts/compute_tidal_currents.py @@ -103,7 +103,6 @@ drop support for the ascii definition file format use model class attributes for file format and corrections add command line option to select nodal corrections type - use attribute for inferring minor long period constituents Updated 08/2024: allow inferring only specific minor constituents added option to try automatic detection of definition file format changed from 'geotiff' to 'GTiff' and 'cog' formats diff --git a/scripts/compute_tidal_elevations.py b/scripts/compute_tidal_elevations.py index 89b38eb..e6e010b 100755 --- a/scripts/compute_tidal_elevations.py +++ b/scripts/compute_tidal_elevations.py @@ -106,7 +106,6 @@ drop support for the ascii definition file format use model class attributes for file format and corrections add command line option to select nodal corrections type - use attribute for inferring minor long period constituents Updated 08/2024: allow inferring only specific minor constituents added option to try automatic detection of definition file format changed from 'geotiff' to 'GTiff' and 'cog' formats diff --git a/test/test_coordinates.py b/test/test_coordinates.py index c3dd569..5c5e36a 100644 --- a/test/test_coordinates.py +++ b/test/test_coordinates.py @@ -5,6 +5,7 @@ UPDATE HISTORY: Updated 09/2024: add test for Arctic regions with new projection + using new JSON dictionary format for model projections Updated 07/2024: add check for if projections are geographic Updated 12/2023: use new crs class for coordinate reprojection Written 08/2020 @@ -14,19 +15,21 @@ import pyTMD.crs # parameterize projections -@pytest.mark.parametrize("PROJ", ['3031','CATS2008','3976','AEDNorth','4326']) +models = ['AOTIM-5-2018','Gr1km-v2','CATS2008','TPXO9-atlas-v5'] +@pytest.mark.parametrize("MODEL", models) # PURPOSE: verify forward and backwards coordinate conversions -def test_coordinates(PROJ): - startlat = {'3031':-60,'CATS2008':-60,'3976':-60,'AEDNorth':60,'4326':90} - endlat = {'3031':-70,'CATS2008':-70,'3976':-70,'AEDNorth':70,'4326':-90} - is_geographic = {'3031':False,'CATS2008':False,'3976':False, - 'AEDNorth':False,'4326':True} +def test_coordinates(MODEL): + startlat = {'Gr1km-v2':60,'CATS2008':-60,'AOTIM-5-2018':60,'TPXO9-atlas-v5':90} + endlat = {'Gr1km-v2':70,'CATS2008':-70,'AOTIM-5-2018':70,'TPXO9-atlas-v5':-90} + is_geographic = {'Gr1km-v2':False,'CATS2008':False, + 'AOTIM-5-2018':False,'TPXO9-atlas-v5':True} i1 = np.arange(-180,180+1,1) - i2 = np.linspace(startlat[PROJ],endlat[PROJ],len(i1)) + i2 = np.linspace(startlat[MODEL],endlat[MODEL],len(i1)) # convert latitude and longitude to and from projection - transform = pyTMD.crs().get(PROJ) - o1, o2 = pyTMD.crs().convert(i1,i2,PROJ,'F') - lon, lat = pyTMD.crs().convert(o1,o2,PROJ,'B') + model = pyTMD.models.elevation[MODEL] + crs = pyTMD.crs().get(model['projection']) + o1, o2 = crs.transform(i1, i2, direction='FORWARD') + lon, lat = crs.transform(o1, o2, direction='INVERSE') # calculate great circle distance between inputs and outputs cdist = np.arccos(np.sin(i2*np.pi/180.0)*np.sin(lat*np.pi/180.0) + np.cos(i2*np.pi/180.0)*np.cos(lat*np.pi/180.0)* @@ -34,7 +37,18 @@ def test_coordinates(PROJ): # test that forward and backwards conversions are within tolerance eps = np.finfo(np.float32).eps assert np.all(cdist < eps) - assert transform.is_geographic == is_geographic[PROJ] + # convert latitude and longitude to and from projection + # using the convert function + o1, o2 = pyTMD.crs().convert(i1, i2, model['projection'], 'F') + lon, lat = pyTMD.crs().convert(o1, o2, model['projection'], 'B') + # calculate great circle distance between inputs and outputs + cdist = np.arccos(np.sin(i2*np.pi/180.0)*np.sin(lat*np.pi/180.0) + + np.cos(i2*np.pi/180.0)*np.cos(lat*np.pi/180.0)* + np.cos((lon-i1)*np.pi/180.0),dtype=np.float32) + # test that forward and backwards conversions are within tolerance + eps = np.finfo(np.float32).eps + assert np.all(cdist < eps) + assert crs.is_geographic == is_geographic[MODEL] # PURPOSE: verify coordinate conversions are close for Arctic regions def test_arctic_projection(): @@ -43,8 +57,10 @@ def test_arctic_projection(): i1 = -180.0 + 360.0*np.random.rand(N) i2 = 60.0 + 30.0*np.random.rand(N) # convert latitude and longitude to and from projection - o1, o2 = pyTMD.crs().convert(i1,i2,'AEDNorth','F') - lon, lat = pyTMD.crs().convert(o1,o2,'AEDNorth','B') + model = pyTMD.models.elevation['AOTIM-5-2018'] + crs = pyTMD.crs().get(model['projection']) + o1, o2 = crs.transform(i1, i2, direction='FORWARD') + lon, lat = crs.transform(o1, o2, direction='INVERSE') # calculate great circle distance between inputs and outputs cdist = np.arccos(np.sin(i2*np.pi/180.0)*np.sin(lat*np.pi/180.0) + np.cos(i2*np.pi/180.0)*np.cos(lat*np.pi/180.0)* diff --git a/test/test_download_and_read.py b/test/test_download_and_read.py index b5edab7..c7a54b6 100644 --- a/test/test_download_and_read.py +++ b/test/test_download_and_read.py @@ -21,6 +21,7 @@ UPDATE HISTORY: Updated 09/2024: drop support for the ascii definition file format use model class attributes for file format and corrections + using new JSON dictionary format for model projections Updated 07/2024: add parametrize over cropping the model fields Updated 04/2024: use timescale for temporal operations Updated 01/2024: refactored compute functions into compute.py @@ -244,12 +245,9 @@ def test_check_CATS2008(self): # PURPOSE: Tests that interpolated results are comparable to AntTG database def test_compare_CATS2008(self): # model parameters for CATS2008 - modelpath = filepath.joinpath('CATS2008') - grid_file = modelpath.joinpath('grid_CATS2008') - model_file = modelpath.joinpath('hf.CATS2008.out') - GRID = 'OTIS' - EPSG = 'CATS2008' - TYPE = 'z' + model = pyTMD.io.model(filepath).elevation('CATS2008') + grid_file = model.grid_file + model_file = model.model_file # open Antarctic Tide Gauge (AntTG) database AntTG = filepath.joinpath('AntTG_ocean_height_v1.txt') @@ -295,8 +293,8 @@ def test_compare_CATS2008(self): # extract amplitude and phase from tide model amp,ph,D,cons = pyTMD.io.OTIS.extract_constants(station_lon, - station_lat, grid_file, model_file, EPSG, type=TYPE, - method='spline', grid=GRID) + station_lat, grid_file, model_file, model['projection'], + type=model['type'], method='spline', grid=model['format']) # reorder constituents of model and convert amplitudes to cm model_amp = np.ma.zeros((antarctic_stations,len(constituents))) model_ph = np.ma.zeros((antarctic_stations,len(constituents))) @@ -333,19 +331,18 @@ def test_compare_CATS2008(self): # parameterize type: heights versus currents parameters = [] - parameters.append(dict(type='z',model='hf.CATS2008.out',grid='grid_CATS2008')) - parameters.append(dict(type='U',model='uv.CATS2008.out',grid='grid_CATS2008')) - parameters.append(dict(type='V',model='uv.CATS2008.out',grid='grid_CATS2008')) + parameters.append(dict(type='z',model='hf.CATS2008.out')) + parameters.append(dict(type='U',model='uv.CATS2008.out')) + parameters.append(dict(type='V',model='uv.CATS2008.out')) @pytest.mark.parametrize("parameters", parameters) # PURPOSE: Tests that interpolated results are comparable to Matlab program def test_verify_CATS2008(self, parameters): # model parameters for CATS2008 - modelpath = filepath.joinpath('CATS2008') - grid_file = modelpath.joinpath(parameters['grid']) + model = pyTMD.io.model(filepath).current('CATS2008') + modelpath = model.grid_file.parent + grid_file = model.grid_file model_file = modelpath.joinpath(parameters['model']) TYPE = parameters['type'] - GRID = 'OTIS' - EPSG = 'CATS2008' # open Antarctic Tide Gauge (AntTG) database AntTG = filepath.joinpath('AntTG_ocean_height_v1.txt') @@ -392,8 +389,8 @@ def test_verify_CATS2008(self, parameters): # extract amplitude and phase from tide model amp,ph,D,c = pyTMD.io.OTIS.extract_constants(station_lon, - station_lat, grid_file, model_file, EPSG, type=TYPE, - method='spline', grid=GRID) + station_lat, grid_file, model_file, model['projection'], + type=TYPE, method='spline', grid=model['format']) # calculate complex phase in radians for Euler's cph = -1j*ph*np.pi/180.0 @@ -436,9 +433,9 @@ def test_verify_CATS2008(self, parameters): # predict tidal elevations at time and infer minor corrections tide.mask[:] = np.any(hc.mask) tide.data[:] = pyTMD.predict.time_series(tide_time, hc, c, - deltat=deltat, corrections=GRID) + deltat=deltat, corrections=model['corrections']) minor = pyTMD.predict.infer_minor(tide_time, hc, c, - deltat=deltat, corrections=GRID) + deltat=deltat, corrections=model['corrections']) tide.data[:] += minor.data[:] # calculate differences between matlab and python version @@ -452,12 +449,11 @@ def test_verify_CATS2008(self, parameters): # PURPOSE: Tests that tidal ellipse results are comparable to Matlab program def test_tidal_ellipse(self): # model parameters for CATS2008 - modelpath = filepath.joinpath('CATS2008') - grid_file = modelpath.joinpath('grid_CATS2008') - model_file = modelpath.joinpath('uv.CATS2008.out') + model = pyTMD.io.model(filepath).current('CATS2008') + modelpath = model.grid_file.parent + grid_file = model.grid_file + model_file = model.model_file['u'] TYPES = ['U','V'] - GRID = 'OTIS' - EPSG = 'CATS2008' # open Antarctic Tide Gauge (AntTG) database AntTG = filepath.joinpath('AntTG_ocean_height_v1.txt') @@ -521,8 +517,9 @@ def test_tidal_ellipse(self): for TYPE in TYPES: # extract amplitude and phase from tide model amp,ph,D,c = pyTMD.io.OTIS.extract_constants(station_lon[i], - station_lat[i], grid_file, model_file, EPSG, type=TYPE, - method='spline', grid=GRID) + station_lat[i], grid_file, model_file, + model['projection'], type=TYPE, + method='spline', grid=model.format) # calculate complex phase in radians for Euler's cph = -1j*ph*np.pi/180.0 # calculate constituent oscillation for station @@ -615,20 +612,19 @@ def test_solve(self, SOLVER): # parameterize type: heights versus currents # parameterize interpolation method parameters = [] - parameters.append(dict(type='z',model='hf.CATS2008.out',grid='grid_CATS2008')) - parameters.append(dict(type='U',model='uv.CATS2008.out',grid='grid_CATS2008')) - parameters.append(dict(type='V',model='uv.CATS2008.out',grid='grid_CATS2008')) + parameters.append(dict(type='z',model='hf.CATS2008.out')) + parameters.append(dict(type='U',model='uv.CATS2008.out')) + parameters.append(dict(type='V',model='uv.CATS2008.out')) @pytest.mark.parametrize("parameters", parameters) @pytest.mark.parametrize("METHOD", ['spline']) # PURPOSE: Tests that interpolated results are comparable def test_compare_constituents(self, parameters, METHOD): # model parameters for CATS2008 - modelpath = filepath.joinpath('CATS2008') - grid_file = modelpath.joinpath(parameters['grid']) + model = pyTMD.io.model(filepath).current('CATS2008') + modelpath = model.grid_file.parent + grid_file = model.grid_file model_file = modelpath.joinpath(parameters['model']) TYPE = parameters['type'] - GRID = 'OTIS' - EPSG = 'CATS2008' # open Antarctic Tide Gauge (AntTG) database AntTG = filepath.joinpath('AntTG_ocean_height_v1.txt') @@ -676,18 +672,17 @@ def test_compare_constituents(self, parameters, METHOD): # extract amplitude and phase from tide model amp1, ph1, D1, c = pyTMD.io.OTIS.extract_constants(station_lon[i], - station_lat[i], grid_file, model_file, EPSG, type=TYPE, - method=METHOD, grid=GRID) + station_lat[i], grid_file, model_file, model['projection'], + type=TYPE, method=METHOD, grid=model.format) # calculate complex form of constituent oscillation hc1 = amp1*np.exp(-1j*ph1*np.pi/180.0) # read complex constituents from tide model constituents = pyTMD.io.OTIS.read_constants(grid_file, model_file, - EPSG, type=TYPE, grid=GRID) + model['projection'], type=TYPE, grid=model.format) # interpolate constituents to station coordinates amp2, ph2, D2 = pyTMD.io.OTIS.interpolate_constants(station_lon[i], - station_lat[i], constituents, EPSG, type=TYPE, - method=METHOD) + station_lat[i], constituents, type=TYPE, method=METHOD) # calculate complex form of constituent oscillation hc2 = amp2*np.exp(-1j*ph2*np.pi/180.0) @@ -822,19 +817,18 @@ def download_Arctic_Tide_Atlas(self): # parameterize type: heights versus currents parameters = [] - parameters.append(dict(type='z',model='h_Arc5km2018',grid='grid_Arc5km2018')) - parameters.append(dict(type='U',model='UV_Arc5km2018',grid='grid_Arc5km2018')) - parameters.append(dict(type='V',model='UV_Arc5km2018',grid='grid_Arc5km2018')) + parameters.append(dict(type='z',model='h_Arc5km2018')) + parameters.append(dict(type='U',model='UV_Arc5km2018')) + parameters.append(dict(type='V',model='UV_Arc5km2018')) @pytest.mark.parametrize("parameters", parameters) # PURPOSE: Tests that interpolated results are comparable to Matlab program def test_verify_AOTIM5_2018(self, parameters): # model parameters for AOTIM-5-2018 - modelpath = filepath.joinpath('Arc5km2018') - grid_file = modelpath.joinpath(parameters['grid']) + model = pyTMD.io.model(filepath).current('AOTIM-5-2018') + modelpath = model.grid_file.parent + grid_file = model.grid_file model_file = modelpath.joinpath(parameters['model']) TYPE = parameters['type'] - GRID = 'OTIS' - EPSG = 'AEDNorth' # open Arctic Tidal Current Atlas list of records ATLAS = filepath.joinpath('List_of_records.txt') @@ -869,8 +863,8 @@ def test_verify_AOTIM5_2018(self, parameters): # extract amplitude and phase from tide model amp,ph,D,c = pyTMD.io.OTIS.extract_constants(station_lon, - station_lat, grid_file, model_file, EPSG, type=TYPE, - method='spline', grid=GRID) + station_lat, grid_file, model_file, model['projection'], + type=TYPE, method='spline', grid=model['format']) # calculate complex phase in radians for Euler's cph = -1j*ph*np.pi/180.0 # will verify differences between model outputs are within tolerance @@ -908,9 +902,9 @@ def test_verify_AOTIM5_2018(self, parameters): # predict tidal elevations at time and infer minor corrections tide.mask[:] = np.any(hc.mask) tide.data[:] = pyTMD.predict.time_series(tide_time, hc, c, - deltat=deltat, corrections=GRID) + deltat=deltat, corrections=model['corrections']) minor = pyTMD.predict.infer_minor(tide_time, hc, c, - deltat=deltat, corrections=GRID) + deltat=deltat, corrections=model['corrections']) tide.data[:] += minor.data[:] # calculate differences between matlab and python version @@ -924,20 +918,19 @@ def test_verify_AOTIM5_2018(self, parameters): # parameterize type: heights versus currents # parameterize interpolation method parameters = [] - parameters.append(dict(type='z',model='h_Arc5km2018',grid='grid_Arc5km2018')) - parameters.append(dict(type='U',model='UV_Arc5km2018',grid='grid_Arc5km2018')) - parameters.append(dict(type='V',model='UV_Arc5km2018',grid='grid_Arc5km2018')) + parameters.append(dict(type='z',model='h_Arc5km2018')) + parameters.append(dict(type='U',model='UV_Arc5km2018')) + parameters.append(dict(type='V',model='UV_Arc5km2018')) @pytest.mark.parametrize("parameters", parameters) @pytest.mark.parametrize("METHOD", ['spline']) # PURPOSE: Tests that interpolated results are comparable def test_compare_constituents(self, parameters, METHOD): # model parameters for AOTIM-5-2018 - modelpath = filepath.joinpath('Arc5km2018') - grid_file = modelpath.joinpath(parameters['grid']) + model = pyTMD.io.model(filepath).current('AOTIM-5-2018') + modelpath = model.grid_file.parent + grid_file = model.grid_file model_file = modelpath.joinpath(parameters['model']) TYPE = parameters['type'] - GRID = 'OTIS' - EPSG = 'AEDNorth' # open Arctic Tidal Current Atlas list of records ATLAS = filepath.joinpath('List_of_records.txt') @@ -968,18 +961,17 @@ def test_compare_constituents(self, parameters, METHOD): # extract amplitude and phase from tide model amp1,ph1,D1,c = pyTMD.io.OTIS.extract_constants(station_lon[i], - station_lat[i], grid_file, model_file, EPSG, type=TYPE, - method=METHOD, grid=GRID) + station_lat[i], grid_file, model_file, model['projection'], + type=TYPE, method=METHOD, grid=model['format']) # calculate complex form of constituent oscillation hc1 = amp1*np.exp(-1j*ph1*np.pi/180.0) # read complex constituents from tide model constituents = pyTMD.io.OTIS.read_constants(grid_file, model_file, - EPSG, type=TYPE, grid=GRID) + model['projection'], type=TYPE, grid=model['format']) # interpolate constituents to station coordinates amp2,ph2,D2 = pyTMD.io.OTIS.interpolate_constants(station_lon[i], - station_lat[i], constituents, EPSG, type=TYPE, - method=METHOD) + station_lat[i], constituents, type=TYPE, method=METHOD) # calculate complex form of constituent oscillation hc2 = amp2*np.exp(-1j*ph2*np.pi/180.0) diff --git a/test/test_model.py b/test/test_model.py index 753b1c9..e6b5cac 100644 --- a/test/test_model.py +++ b/test/test_model.py @@ -5,6 +5,7 @@ UPDATE HISTORY: Updated 09/2024: drop support for the ascii definition file format fix parsing of TPXO8-atlas-nc constituents + using new JSON dictionary format for model projections Updated 08/2024: add automatic detection of definition file format Updated 07/2024: add new JSON format definition file format Written 04/2024 @@ -34,7 +35,9 @@ def test_definition_CATS2008(): assert m.name == 'CATS2008' assert m.model_file == pathlib.Path('CATS2008/hf.CATS2008.out') assert m.grid_file == pathlib.Path('CATS2008/grid_CATS2008') - assert m.projection == 'CATS2008' + assert m.projection == {'datum': 'WGS84', 'lat_0': -90, 'lat_ts': -71, + 'lon_0': -70, 'proj': 'stere', 'type': 'crs', 'units': 'km', + 'x_0': 0, 'y_0': 0} assert m.type == 'z' assert m.variable == 'tide_ocean' # test properties