From a254db0b1122bce9f4238a29b3b429a348f90016 Mon Sep 17 00:00:00 2001 From: Tyler Sutterley Date: Tue, 26 Nov 2024 15:57:56 -0800 Subject: [PATCH] feat: expose buffer distance to crop tide model data for #367 --- pyTMD/compute.py | 13 ++++++++++--- pyTMD/io/ATLAS.py | 17 ++++++++++++++--- pyTMD/io/FES.py | 14 ++++++++++++-- pyTMD/io/GOT.py | 14 ++++++++++++-- pyTMD/io/OTIS.py | 22 ++++++++++++++++------ scripts/compute_tidal_currents.py | 13 ++++++++++--- scripts/compute_tidal_elevations.py | 15 ++++++++++++--- 7 files changed, 86 insertions(+), 22 deletions(-) diff --git a/pyTMD/compute.py b/pyTMD/compute.py index 97246b79..8fcf6b39 100644 --- a/pyTMD/compute.py +++ b/pyTMD/compute.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" compute.py -Written by Tyler Sutterley (10/2024) +Written by Tyler Sutterley (11/2024) Calculates tidal elevations for correcting elevation or imagery data Calculates tidal currents at locations and times @@ -60,6 +60,7 @@ interpolate.py: interpolation routines for spatial data UPDATE HISTORY: + Updated 11/2024: expose buffer distance for cropping tide model data Updated 10/2024: compute delta times based on corrections type simplify by using wrapper functions to read and interpolate constants added option to append equilibrium amplitudes for node tides @@ -216,6 +217,7 @@ def tide_elevations( DEFINITION_FILE: str | pathlib.Path | IOBase | None = None, CROP: bool = False, BOUNDS: list | np.ndarray | None = None, + BUFFER: int | float | None = None, EPSG: str | int = 3031, EPOCH: list | tuple = (2000, 1, 1, 0, 0, 0), TYPE: str | None = 'drift', @@ -255,6 +257,8 @@ def tide_elevations( Crop tide model data to (buffered) bounds BOUNDS: list, np.ndarray or NoneType, default None Boundaries for cropping tide model data + BUFFER: int, float or NoneType, default None + Buffer distance for cropping tide model data EPSG: int, default: 3031 (Polar Stereographic South, WGS84) Input coordinate system EPOCH: tuple, default (2000,1,1,0,0,0) @@ -359,7 +363,7 @@ def tide_elevations( # read tidal constants and interpolate to grid points amp, ph, c = model.extract_constants(lon, lat, type=model.type, - crop=CROP, bounds=BOUNDS, method=METHOD, + crop=CROP, bounds=BOUNDS, buffer=BUFFER, method=METHOD, extrapolate=EXTRAPOLATE, cutoff=CUTOFF, append_node=APPEND_NODE, apply_flexure=APPLY_FLEXURE) # calculate complex phase in radians for Euler's @@ -441,6 +445,7 @@ def tide_currents( DEFINITION_FILE: str | pathlib.Path | IOBase | None = None, CROP: bool = False, BOUNDS: list | np.ndarray | None = None, + BUFFER: int | float | None = None, EPSG: str | int = 3031, EPOCH: list | tuple = (2000, 1, 1, 0, 0, 0), TYPE: str | None = 'drift', @@ -478,6 +483,8 @@ def tide_currents( Crop tide model data to (buffered) bounds BOUNDS: list, np.ndarray or NoneType, default None Boundaries for cropping tide model data + BUFFER: int, float or NoneType, default None + Buffer distance for cropping tide model data EPSG: int, default: 3031 (Polar Stereographic South, WGS84) Input coordinate system EPOCH: tuple, default (2000,1,1,0,0,0) @@ -580,7 +587,7 @@ def tide_currents( for t in model.type: # read tidal constants and interpolate to grid points amp, ph, c = model.extract_constants(lon, lat, type=t, - crop=CROP, bounds=BOUNDS, method=METHOD, + crop=CROP, bounds=BOUNDS, buffer=BUFFER, method=METHOD, extrapolate=EXTRAPOLATE, cutoff=CUTOFF) # calculate complex phase in radians for Euler's cph = -1j*ph*np.pi/180.0 diff --git a/pyTMD/io/ATLAS.py b/pyTMD/io/ATLAS.py index 96d385f3..ab6aba7b 100644 --- a/pyTMD/io/ATLAS.py +++ b/pyTMD/io/ATLAS.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" ATLAS.py -Written by Tyler Sutterley (10/2024) +Written by Tyler Sutterley (11/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 @@ -55,6 +55,7 @@ interpolate.py: interpolation routines for spatial data UPDATE HISTORY: + Updated 11/2024: expose buffer distance for cropping tide model data Updated 10/2024: fix error when using default bounds in extract_constants Updated 07/2024: added crop and bounds keywords for trimming model data Updated 02/2024: changed variable for setting global grid flag to is_global @@ -164,6 +165,8 @@ def extract_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: int, float or NoneType, default None + Buffer angle for cropping tide model data method: str, default 'spline' Interpolation method @@ -196,6 +199,7 @@ def extract_constants( kwargs.setdefault('type', 'z') kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', None) kwargs.setdefault('method', 'spline') kwargs.setdefault('extrapolate', False) kwargs.setdefault('cutoff', 10.0) @@ -235,6 +239,8 @@ def extract_constants( bounds = kwargs['bounds'] or [xmin, xmax, ymin, ymax] # grid step size of tide model dlon = lon[1] - lon[0] + # default buffer if cropping data + buffer = kwargs['buffer'] or 4*dlon # if global: extend limits is_global = False @@ -244,7 +250,7 @@ def extract_constants( mlon, mlat = np.copy(lon), np.copy(lat) bathymetry, lon, lat = _crop(bathymetry, mlon, mlat, bounds=bounds, - buffer=4*dlon + buffer=buffer ) elif (np.min(ilon) < 0.0) & (np.max(lon) > 180.0): # input points convention (-180:180) @@ -320,7 +326,7 @@ def extract_constants( if kwargs['crop'] and np.any(bounds): hc, _, _ = _crop(hc, mlon, mlat, bounds=bounds, - buffer=4*dlon + buffer=buffer ) # replace original values with extend matrices if is_global: @@ -412,6 +418,8 @@ def read_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: int or float, default 0 + Buffer angle for cropping tide model data Returns ------- @@ -423,6 +431,7 @@ def read_constants( kwargs.setdefault('compressed', True) kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', 0) # raise warning if model files are entered as a string or path if isinstance(model_files, (str, pathlib.Path)): @@ -444,6 +453,7 @@ def read_constants( mlon, mlat = np.copy(lon), np.copy(lat) bathymetry, lon, lat = _crop(bathymetry, mlon, mlat, bounds=kwargs['bounds'], + buffer=kwargs['buffer'], ) # grid step size of tide model dlon = lon[1] - lon[0] @@ -473,6 +483,7 @@ def read_constants( if kwargs['crop'] and np.any(kwargs['bounds']): hc, lon, lat = _crop(hc, mlon, mlat, bounds=kwargs['bounds'], + buffer=kwargs['buffer'], ) # replace original values with extend matrices if is_global: diff --git a/pyTMD/io/FES.py b/pyTMD/io/FES.py index 91370464..8a5b66c3 100644 --- a/pyTMD/io/FES.py +++ b/pyTMD/io/FES.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" FES.py -Written by Tyler Sutterley (10/2024) +Written by Tyler Sutterley (11/2024) Reads files for a tidal model and makes initial calculations to run tide program Includes functions to extract tidal harmonic constants from the @@ -57,6 +57,7 @@ interpolate.py: interpolation routines for spatial data UPDATE HISTORY: + Updated 11/2024: expose buffer distance for cropping tide model data Updated 10/2024: fix error when using default bounds in extract_constants Updated 07/2024: added new FES2022 to available known model versions FES2022 have masked longitudes, only extract longitude data @@ -175,6 +176,8 @@ def extract_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: float or NoneType, default None + Buffer angle for cropping tide model data method: str, default 'spline' Interpolation method @@ -203,6 +206,7 @@ def extract_constants( kwargs.setdefault('compressed', False) kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', None) kwargs.setdefault('method', 'spline') kwargs.setdefault('extrapolate', False) kwargs.setdefault('cutoff', 10.0) @@ -255,12 +259,14 @@ def extract_constants( hc, lon, lat = read_netcdf_file(model_file, **kwargs) # grid step size of tide model dlon = lon[1] - lon[0] + # default buffer if cropping data + buffer = kwargs['buffer'] or 4*dlon # crop tide model data to (buffered) bounds # or adjust longitudinal convention to fit tide model if kwargs['crop'] and np.any(bounds): hc, lon, lat = _crop(hc, lon, lat, bounds=bounds, - buffer=4*dlon + buffer=buffer, ) elif (np.min(ilon) < 0.0) & (np.max(lon) > 180.0): # input points convention (-180:180) @@ -372,6 +378,8 @@ def read_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: int or float, default 0 + Buffer angle for cropping tide model data Returns ------- @@ -384,6 +392,7 @@ def read_constants( kwargs.setdefault('compressed', False) kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', 0) # raise warning if model files are entered as a string or path if isinstance(model_files, (str, pathlib.Path)): @@ -414,6 +423,7 @@ def read_constants( if kwargs['crop'] and np.any(kwargs['bounds']): hc, lon, lat = _crop(hc, lon, lat, bounds=kwargs['bounds'], + buffer=kwargs['buffer'], ) # grid step size of tide model dlon = lon[1] - lon[0] diff --git a/pyTMD/io/GOT.py b/pyTMD/io/GOT.py index 72b2a720..7ce8eabc 100644 --- a/pyTMD/io/GOT.py +++ b/pyTMD/io/GOT.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" GOT.py -Written by Tyler Sutterley (10/2024) +Written by Tyler Sutterley (11/2024) Reads files for Richard Ray's Global Ocean Tide (GOT) models and makes initial calculations to run the tide program @@ -42,6 +42,7 @@ interpolate.py: interpolation routines for spatial data UPDATE HISTORY: + Updated 11/2024: expose buffer distance for cropping tide model data Updated 10/2024: fix error when using default bounds in extract_constants Updated 07/2024: added crop and bounds keywords for trimming model data use parse function from constituents class to extract names @@ -147,6 +148,8 @@ def extract_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: int or float, default None + Buffer angle for cropping tide model data method: str, default 'spline' Interpolation method @@ -176,6 +179,7 @@ def extract_constants( kwargs.setdefault('compressed', False) kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', None) kwargs.setdefault('method', 'spline') kwargs.setdefault('extrapolate', False) kwargs.setdefault('cutoff', 10.0) @@ -232,12 +236,14 @@ def extract_constants( constituents.append(cons) # grid step size of tide model dlon = np.abs(lon[1] - lon[0]) + # default buffer if cropping data + buffer = kwargs['buffer'] or 4*dlon # crop tide model data to (buffered) bounds # or adjust longitudinal convention to fit tide model if kwargs['crop'] and np.any(bounds): hc, lon, lat = _crop(hc, lon, lat, bounds=bounds, - buffer=4*dlon + buffer=buffer, ) elif (np.min(ilon) < 0.0) & (np.max(lon) > 180.0): # input points convention (-180:180) @@ -326,6 +332,8 @@ def read_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: int or float, default 0 + Buffer angle for cropping tide model data Returns ------- @@ -337,6 +345,7 @@ def read_constants( kwargs.setdefault('compressed', False) kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', 0) # raise warning if model files are entered as a string if isinstance(model_files, (str, pathlib.Path)): @@ -362,6 +371,7 @@ def read_constants( if kwargs['crop'] and np.any(kwargs['bounds']): hc, lon, lat = _crop(hc, lon, lat, bounds=kwargs['bounds'], + buffer=kwargs['buffer'], ) # grid step size of tide model dlon = np.abs(lon[1] - lon[0]) diff --git a/pyTMD/io/OTIS.py b/pyTMD/io/OTIS.py index e86c1466..3c29aa83 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 (10/2024) +Written by Tyler Sutterley (11/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 @@ -58,6 +58,7 @@ interpolate.py: interpolation routines for spatial data UPDATE HISTORY: + Updated 11/2024: expose buffer distance for cropping tide model data Updated 10/2024: save latitude and longitude to output constituent object fix error when using default bounds in extract_constants Updated 09/2024: using new JSON dictionary format for model projections @@ -212,6 +213,8 @@ def extract_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: int, float or NoneType, default None + Buffer angle or distance for cropping tide model data method: str, default 'spline' Interpolation method @@ -243,6 +246,7 @@ def extract_constants( kwargs.setdefault('grid', 'OTIS') kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', None) kwargs.setdefault('method', 'spline') kwargs.setdefault('extrapolate', False) kwargs.setdefault('cutoff', 10.0) @@ -291,15 +295,17 @@ def extract_constants( xmin, xmax = np.min(x), np.max(x) ymin, ymax = np.min(y), np.max(y) bounds = kwargs['bounds'] or [xmin, xmax, ymin, ymax] + # default buffer if cropping data + buffer = kwargs['buffer'] or 4*dx # crop mask and bathymetry data to (buffered) bounds # or adjust longitudinal convention to fit tide model if kwargs['crop'] and np.any(bounds): mx, my = np.copy(xi), np.copy(yi) mz, xi, yi = _crop(mz, mx, my, bounds=bounds, - buffer=4*dx, is_geographic=is_geographic) + buffer=buffer, is_geographic=is_geographic) hz, xi, yi = _crop(hz, mx, my, bounds=bounds, - buffer=4*dx, is_geographic=is_geographic) + buffer=buffer, is_geographic=is_geographic) elif (np.min(x) < np.min(xi)) & is_geographic: # input points convention (-180:180) # tide model convention (0:360) @@ -445,7 +451,7 @@ def extract_constants( # crop tide model data to (buffered) bounds if kwargs['crop'] and np.any(bounds): hc, _, _ = _crop(hc, mx, my, - bounds=bounds, buffer=4*dx, + bounds=bounds, buffer=buffer, is_geographic=is_geographic) # replace original values with extend matrices if is_global: @@ -549,6 +555,8 @@ def read_constants( Crop tide model data to (buffered) bounds bounds: list or NoneType, default None Boundaries for cropping tide model data + buffer: int or float, default 0 + Buffer angle or distance for cropping tide model data apply_flexure: bool, default False Apply ice flexure scaling factor to height values @@ -562,6 +570,7 @@ def read_constants( kwargs.setdefault('grid', 'OTIS') kwargs.setdefault('crop', False) kwargs.setdefault('bounds', None) + kwargs.setdefault('buffer', 0) kwargs.setdefault('apply_flexure', False) # check that grid file is accessible @@ -599,9 +608,9 @@ def read_constants( # crop tide model data mx, my = np.copy(xi), np.copy(yi) mz, xi, yi = _crop(mz, mx, my, bounds=kwargs['bounds'], - is_geographic=is_geographic) + buffer=kwargs['buffer'], is_geographic=is_geographic) hz, xi, yi = _crop(hz, mx, my, bounds=kwargs['bounds'], - is_geographic=is_geographic) + buffer=kwargs['buffer'], is_geographic=is_geographic) # replace original values with extend arrays/matrices if ((xi[-1] - xi[0]) == (360.0 - dx)) & is_geographic: @@ -710,6 +719,7 @@ def read_constants( if kwargs['crop'] and np.any(kwargs['bounds']): hc, _, _ = _crop(hc, mx, my, bounds=kwargs['bounds'], + buffer=kwargs['buffer'], is_geographic=is_geographic) # replace original values with extend matrices if is_global: diff --git a/scripts/compute_tidal_currents.py b/scripts/compute_tidal_currents.py index 452f03af..8b003737 100755 --- a/scripts/compute_tidal_currents.py +++ b/scripts/compute_tidal_currents.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" compute_tidal_currents.py -Written by Tyler Sutterley (10/2024) +Written by Tyler Sutterley (11/2024) Calculates zonal and meridional tidal currents for an input file Uses OTIS format tidal solutions provided by Oregon State University and ESR @@ -23,6 +23,7 @@ --gzip, -G: Tide model files are gzip compressed --definition-file X: Model definition file for use in calculating currents -C, --crop: Crop tide model to (buffered) bounds of data + -B X, --buffer X: Buffer for cropping tide model --format X: input and output data format csv (default) netCDF4 @@ -103,6 +104,7 @@ predict.py: predict tidal values using harmonic constants UPDATE HISTORY: + Updated 11/2024: add option for buffer distance to crop tide model data Updated 10/2024: compute delta times based on corrections type simplify by using wrapper functions to read and interpolate constants Updated 09/2024: use JSON database for known model parameters @@ -221,6 +223,7 @@ def compute_tidal_currents(tide_dir, input_file, output_file, GZIP=True, DEFINITION_FILE=None, CROP=False, + BUFFER=None, FORMAT='csv', VARIABLES=[], HEADER=0, @@ -311,8 +314,8 @@ def compute_tidal_currents(tide_dir, input_file, output_file, for t in model.type: # read tidal constants and interpolate to grid points amp, ph, c = model.extract_constants(np.ravel(lon), np.ravel(lat), - type=t, crop=CROP, method=METHOD, extrapolate=EXTRAPOLATE, - cutoff=CUTOFF) + type=t, crop=CROP, buffer=BUFFER, method=METHOD, + extrapolate=EXTRAPOLATE, cutoff=CUTOFF) # calculate complex phase in radians for Euler's cph = -1j*ph*np.pi/180.0 # calculate constituent oscillation @@ -490,6 +493,9 @@ def arguments(): parser.add_argument('--crop', '-C', default=False, action='store_true', help='Crop tide model to bounds of data') + parser.add_argument('--buffer', '-B', + type=float, default=None, + help='Buffer for cropping tide model') # input and output data format parser.add_argument('--format','-F', type=str, default='csv', @@ -604,6 +610,7 @@ def main(): GZIP=args.gzip, DEFINITION_FILE=args.definition_file, CROP=args.crop, + BUFFER=args.buffer, FORMAT=args.format, VARIABLES=args.variables, HEADER=args.header, diff --git a/scripts/compute_tidal_elevations.py b/scripts/compute_tidal_elevations.py index 2c9ef0e9..7c7bd6d7 100755 --- a/scripts/compute_tidal_elevations.py +++ b/scripts/compute_tidal_elevations.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" compute_tidal_elevations.py -Written by Tyler Sutterley (10/2024) +Written by Tyler Sutterley (11/2024) Calculates tidal elevations for an input file Uses OTIS format tidal solutions provided by Oregon State University and ESR @@ -24,6 +24,7 @@ --gzip, -G: Tide model files are gzip compressed --definition-file X: Model definition file for use as correction -C, --crop: Crop tide model to (buffered) bounds of data + -B X, --buffer X: Buffer for cropping tide model --format X: input and output data format csv (default) netCDF4 @@ -107,6 +108,7 @@ predict.py: predict tidal values using harmonic constants UPDATE HISTORY: + Updated 11/2024: add option for buffer distance to crop tide model data Updated 10/2024: compute delta times based on corrections type simplify by using wrapper functions to read and interpolate constants added option to append equilibrium amplitudes for node tides @@ -227,6 +229,7 @@ def compute_tidal_elevations(tide_dir, input_file, output_file, GZIP=True, DEFINITION_FILE=None, CROP=False, + BUFFER=None, FORMAT='csv', VARIABLES=[], HEADER=0, @@ -315,8 +318,10 @@ def compute_tidal_elevations(tide_dir, input_file, output_file, # read tidal constants and interpolate to grid points amp, ph, c = model.extract_constants(np.ravel(lon), np.ravel(lat), - type=model.type, crop=CROP, method=METHOD, extrapolate=EXTRAPOLATE, - cutoff=CUTOFF, append_node=APPEND_NODE, apply_flexure=APPLY_FLEXURE) + type=model.type, crop=CROP, buffer=BUFFER, + method=METHOD, extrapolate=EXTRAPOLATE, + cutoff=CUTOFF, append_node=APPEND_NODE, + apply_flexure=APPLY_FLEXURE) # calculate complex phase in radians for Euler's cph = -1j*ph*np.pi/180.0 # calculate constituent oscillation @@ -483,6 +488,9 @@ def arguments(): parser.add_argument('--crop', '-C', default=False, action='store_true', help='Crop tide model to bounds of data') + parser.add_argument('--buffer', '-B', + type=float, default=None, + help='Buffer for cropping tide model') # input and output data format parser.add_argument('--format','-F', type=str, default='csv', @@ -606,6 +614,7 @@ def main(): GZIP=args.gzip, DEFINITION_FILE=args.definition_file, CROP=args.crop, + BUFFER=args.buffer, FORMAT=args.format, VARIABLES=args.variables, HEADER=args.header,