Skip to content

GitLab

Commit 14ffd7d3 authored by thopri's avatar thopri
Browse files

updated fes tide code

parent a6d44375
Hide whitespace changes
Inline Side-by-side
Showing with 113 additions and 56 deletions
inputs/namelist_remote.bdy
View file @ 14ffd7d3
......@@ -83,9 +83,12 @@
nn_month_end = 11 ! month end (default = 12 is years>1)
sn_dst_calendar = 'gregorian' ! output calendar format
nn_base_year = 1960 ! base year for time counter
sn_tide_grid = './src_data/tide/grid_tpxo7.2.nc'
sn_tide_h = './src_data/tide/h_tpxo7.2.nc'
sn_tide_u = './src_data/tide/u_tpxo7.2.nc'
! TPXO file locations
sn_tide_grid = '/Users/thopri/Projects/PyNEMO/DATA/TPXO/grid_tpxo7.2.nc'
sn_tide_h = '/Users/thopri/Projects/PyNEMO/DATA/TPXO/h_tpxo7.2.nc'
sn_tide_u = './Users/thopri/Projects/PyNEMO/DATA/TPXO/u_tpxo7.2.nc'
! location of FES data
sn_tide_fes = '/Users/thopri/Projects/PyNEMO/DATA/FES/'
!------------------------------------------------------------------------------
! Additional parameters
......
pynemo/profile.py
View file @ 14ffd7d3
......@@ -438,13 +438,15 @@ def process_bdy(setup_filepath=0, mask_gui=False):
if settings['tide']:
if settings['tide_model']=='tpxo':
cons = tide.nemo_bdy_tpx7p2_rot(
cons = tide.nemo_bdy_tide_rot(
Setup, DstCoord, bdy_ind['t'], bdy_ind['u'], bdy_ind['v'],
settings['clname'])
settings['clname'], settings['tide_model'])
elif settings['tide_model']=='fes':
logger.error('Tidal model: %s, not yet implimented',
cons = tide.nemo_bdy_tide_rot(
Setup, DstCoord, bdy_ind['t'], bdy_ind['u'], bdy_ind['v'],
settings['clname'],settings['tide_model'])
logger.warning('Tidal model: %s, not yet properly implimented',
settings['tide_model'])
return
else:
logger.error('Tidal model: %s, not recognised',
settings['tide_model'])
......
pynemo/tide/fes_extract_HC.py
View file @ 14ffd7d3
......@@ -12,8 +12,8 @@ from netCDF4 import Dataset
from scipy import interpolate
import numpy as np
class FESExtract(object):
""" This is TPXO model extract_hc.c implementation in python"""
class FesExtract(object):
""" This is FES model extract_hc.c implementation in python adapted from tpxo_extract_HC.py"""
def __init__(self, settings, lat, lon, grid_type):
"""initialises the Extract of tide information from the netcdf
Tidal files"""
......@@ -21,37 +21,69 @@ class FESExtract(object):
tide_model = 'fes'
if tide_model == 'fes': # Define stuff to generalise Tide model
hRe_name = 'hRe'
hIm_name = 'hIm'
lon_z_name = 'lon_z'
lat_z_name = 'lat_z'
URe_name = 'URe'
UIm_name = 'UIm'
lon_u_name = 'lon_u'
lat_u_name = 'lat_u'
VRe_name = 'VRe'
VIm_name = 'VIm'
lon_v_name = 'lon_v'
lat_v_name = 'lat_v'
mz_name = 'mz'
mu_name = 'mu'
mv_name = 'mv'
self.grid = Dataset(settings['tide_grid'])#../data/tide/grid_tpxo7.2.nc')
#read the height_dataset file
self.height_dataset = Dataset(settings['tide_h'])#../data/tide/h_tpxo7.2.nc')
#read the velocity_dataset file
self.velocity_dataset = Dataset(settings['tide_u'])#../data/tide/u_tpxo7.2.nc')
height_z = self.grid.variables['hz']
mask_z = self.grid.variables['mz']
lon_z = self.height_dataset.variables[lon_z_name][:, 0]
lat_z = self.height_dataset.variables[lat_z_name][0, :]
lon_resolution = lon_z[1] - lon_z[0]
data_in_km = 0 # added to maintain the reference to matlab tmd code
# Pull out the constituents that are avaibable
self.cons = []
for ncon in range(self.height_dataset.variables['con'].shape[0]):
self.cons.append(self.height_dataset.variables['con'][ncon, :].tostring().strip())
hRe_name = 'hRe'
hIm_name = 'hIm'
lon_z_name = 'lon'
lat_z_name = 'lat'
URe_name = 'URe'
UIm_name = 'UIm'
lon_u_name = 'lon'
lat_u_name = 'lat'
VRe_name = 'VRe'
VIm_name = 'VIm'
lon_v_name = 'lon'
lat_v_name = 'lat'
mz_name = 'mz'
mu_name = 'mu'
mv_name = 'mv'
#constituents = ['2N2', 'EPS2', 'J1', 'K1', 'K2', 'L2', 'LA2', 'M2', 'M3', 'M4', 'M6', 'M8', 'MF', 'MKS2',
#'MM', 'MN4', 'MS4', 'MSF', 'MSQM', 'MTM', 'MU2', 'N2', 'N4', 'NU2', 'O1', 'P1', 'Q1', 'R2',
#'S1', 'S2', 'S4', 'SA', 'SSA', 'T2']
constituents = ['M2','S2']
# extract lon and lat z data
lon_z = Dataset[settings['tide_fes']+constituents[0]+'_Z.nc'].variables['lon']
lat_z = Dataset[settings['tide_fes']+constituents[0]+'_Z.nc'].variabels['lat']
lon_resolution = lon_z[1] - lon_z[0]
data_in_km = 0 # added to maintain the reference to matlab tmd code
self.cons = constituents
#read and convert the height_dataset file to complex
for ncon in constituents:
amp = Dataset(settings['tide_fes']+constituents[ncon]+'_Z.nc').variables['amplitude'][:]
phase = Dataset(settings['tide_fes']+constituents[ncon]+'_Z.nc').variables['phase'][:]
hRe = amp*np.sin(phase)
hIm = amp*np.cos(phase)
self.height_dataset = {'hRe':hRe,'hIm':hIm}
#read and convert the velocity_dataset files to complex
for ncon in constituents:
amp = Dataset(settings['tide_fes']+constituents[ncon]+'_U.nc').variables['Ua'][:]
phase = Dataset(settings['tide_fes']+constituents[ncon]+'_U.nc').variables['Ug'][:]
URe = amp*np.sin(phase)
UIm = amp*np.cos(phase)
self.Uvelocity_dataset = {'URe':URe,'UIm':UIm}
for ncon in constituents:
amp = Dataset(settings['tide_fes']+constituents[ncon]+'_V.nc').variables['Va'][:]
phase = Dataset(settings['tide_fes']+constituents[ncon]+'_V.nc').variables['Vg'][:]
VRe = amp*np.sin(phase)
VIm = amp*np.cos(phase)
self.Vvelocity_dataset = {'VRe':VRe,'VIm':VIm}
# extract example amplitude grid and change NaNs to 0 (for land) and other values to 1 (for water)
mask_z = Dataset(settings['tide_fes']+constituents[ncon]+'_Z.nc').variables['amplitude'][:]
where_are_NaNs = np.isnan(mask_z)
where_no_NaNs = np.invert(np.isnan(mask_z))
mask_z[where_no_NaNs] = 1
mask_z[where_are_NaNs] = 0
# need to sort for FES
self.grid = Dataset(settings['tide_grid'])#../data/tide/grid_tpxo7.2.nc')
height_z = self.grid.variables['hz']
else:
print('Don''t know that tide model')
......
pynemo/tide/nemo_bdy_tide3.py
View file @ 14ffd7d3
......@@ -8,6 +8,7 @@ Module to extract constituents for the input grid mapped onto output grid
# pylint: disable=no-name-in-module
import copy
from . import tpxo_extract_HC
from . import fes_extract_HC
import numpy as np
from netCDF4 import Dataset
from pynemo import nemo_bdy_grid_angle
......@@ -16,8 +17,8 @@ from pynemo.reader.factory import GetFile
import logging
def nemo_bdy_tpx7p2_rot(setup, DstCoord, Grid_T, Grid_U, Grid_V, comp):
""" TPXO Global Tidal model interpolation including rotation grid"""
def nemo_bdy_tide_rot(setup, DstCoord, Grid_T, Grid_U, Grid_V, comp,tide_model):
""" Global Tidal model interpolation including rotation grid"""
key_transport = 0 # compute the velocities from transport
numharm = len(comp)
logger = logging.getLogger(__name__)
......@@ -52,16 +53,25 @@ def nemo_bdy_tpx7p2_rot(setup, DstCoord, Grid_T, Grid_U, Grid_V, comp):
#convert the U-longitudes into the TMD conventions (0/360E)
dst_lon[dst_lon < 0.0] = dst_lon[dst_lon < 0.0]+360.0
if tide_model == 'tpxo':
tpxo_ux = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_U.grid_type)
tpxo_vx = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_V.grid_type)
tpxo_ux = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_U.grid_type)
tpxo_vx = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_V.grid_type)
ampuX = tpxo_ux.amp
phauX = tpxo_ux.gph
ampvX = tpxo_vx.amp
phavX = tpxo_vx.gph
ampuX = tpxo_ux.amp
phauX = tpxo_ux.gph
ampvX = tpxo_vx.amp
phavX = tpxo_vx.gph
if tide_model == 'fes':
fes_uy = fes_extract_HC.FesExtract(setup.settings, dst_lat, dst_lon, Grid_U.grid_type)
fes_vy = fes_extract_HC.FesExtract(setup.settings, dst_lat, dst_lon, Grid_V.grid_type)
#check if ux data are missing
ampuY = fes_uy.amp
phauY = fes_uy.gph
ampvY = fes_vy.amp
phavY = fes_vy.gph
#check if ux data are missing
ind = np.where((np.isnan(ampuX)) | (np.isnan(phauX)))
if ind[0].size > 0:
logger.warning('Missing zonal velocity along the x open boundary')
......@@ -83,13 +93,23 @@ def nemo_bdy_tpx7p2_rot(setup, DstCoord, Grid_T, Grid_U, Grid_V, comp):
#convert the U-longitudes into the TMD conventions (0/360E)
dst_lon[dst_lon < 0.0] = dst_lon[dst_lon < 0.0]+360.0
tpxo_uy = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_U.grid_type)
tpxo_vy = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_V.grid_type)
ampuY = tpxo_uy.amp
phauY = tpxo_uy.gph
ampvY = tpxo_vy.amp
phavY = tpxo_vy.gph
if tide_model == 'tpxo':
tpxo_uy = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_U.grid_type)
tpxo_vy = tpxo_extract_HC.TpxoExtract(setup.settings, dst_lat, dst_lon, Grid_V.grid_type)
ampuY = tpxo_uy.amp
phauY = tpxo_uy.gph
ampvY = tpxo_vy.amp
phavY = tpxo_vy.gph
if tide_model == 'fes':
fes_uy = fes_extract_HC.FesExtract(setup.settings, dst_lat, dst_lon, Grid_U.grid_type)
fes_vy = fes_extract_HC.FesExtract(setup.settings, dst_lat, dst_lon, Grid_V.grid_type)
ampuY = fes_uy.amp
phauY = fes_uy.gph
ampvY = fes_vy.amp
phavY = fes_vy.gph
#check if ux data are missing
ind = np.where((np.isnan(ampuY)) | (np.isnan(phauY)))
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment