Udpt README.md

.gitignore

@@ -0,0 +1,15 @@
+# FOLDERS 
+./.vscode
+tests/data
+env
+./src/_pycache__
+
+# FILES
+.env
+*.pyc
+*.log
+*.csv
+*.nc
+*.log
+!.vscode/tasks.json
+!.vscode/settings.json

.vscode/settings.json

@@ -0,0 +1,3 @@
+{
+    "workbench.colorTheme": "Dracula Soft"
+}

IHSetShoreFor.egg-info/PKG-INFO

@@ -0,0 +1,23 @@
+Metadata-Version: 2.1
+Name: IHSetShoreFor
+Version: 1.1.0
+Summary: IH-SET Miller and Dean (2004)
+Home-page: https://github.com/defreitasL/IHSetShoreFor
+Author: Lucas de Freitas Pereira
+Author-email: [email protected]
+Classifier: Development Status :: 3 - Alpha
+Classifier: Intended Audience :: Developers
+Classifier: License :: OSI Approved :: MIT License
+Classifier: Programming Language :: Python :: 3
+Classifier: Programming Language :: Python :: 3.6
+Classifier: Programming Language :: Python :: 3.7
+Classifier: Programming Language :: Python :: 3.8
+Classifier: Programming Language :: Python :: 3.9
+License-File: LICENSE
+Requires-Dist: numpy
+Requires-Dist: xarray
+Requires-Dist: numba
+Requires-Dist: datetime
+Requires-Dist: spotpy
+Requires-Dist: IHSetCalibration@ git+https://github.com/defreitasL/IHSetCalibration.git
+Requires-Dist: IHSetUtils@ git+https://github.com/IHCantabria/IHSetUtils.git

IHSetShoreFor.egg-info/SOURCES.txt

@@ -0,0 +1,11 @@
+LICENSE
+README.md
+setup.py
+IHSetShoreFor/__init__.py
+IHSetShoreFor/calibration.py
+IHSetShoreFor/shoreFor.py
+IHSetShoreFor.egg-info/PKG-INFO
+IHSetShoreFor.egg-info/SOURCES.txt
+IHSetShoreFor.egg-info/dependency_links.txt
+IHSetShoreFor.egg-info/requires.txt
+IHSetShoreFor.egg-info/top_level.txt

IHSetShoreFor.egg-info/dependency_links.txt

@@ -0,0 +1 @@
+

IHSetShoreFor.egg-info/requires.txt

@@ -0,0 +1,7 @@
+numpy
+xarray
+numba
+datetime
+spotpy
+IHSetCalibration@ git+https://github.com/defreitasL/IHSetCalibration.git
+IHSetUtils@ git+https://github.com/IHCantabria/IHSetUtils.git

IHSetShoreFor.egg-info/top_level.txt

@@ -0,0 +1 @@
+IHSetShoreFor

IHSetShoreFor/__init__.py

@@ -0,0 +1,5 @@
+# src/__init__.py
+
+# Import modules and functions from your package here
+from .shoreFor import shoreFor
+from .calibration import cal_ShoreFor

IHSetShoreFor/calibration.py

@@ -0,0 +1,193 @@
+import numpy as np
+import xarray as xr
+from datetime import datetime
+from spotpy.parameter import Uniform
+from .shoreFor import shoreFor
+from IHSetCalibration import objective_functions
+from IHSetUtils import wMOORE
+
+class cal_ShoreFor(object):
+    """
+    cal_ShoreFor
+    
+    Configuration to calibrate and run the Davidson et al. (2013) Shoreline Evolution Model.
+    
+    This class reads input datasets, performs its calibration.
+    """
+
+    def __init__(self, path):
+
+        self.path = path
+
+
+        mkTime = np.vectorize(lambda Y, M, D, h: datetime(int(Y), int(M), int(D), int(h), 0, 0))
+
+        cfg = xr.open_dataset(path+'config.nc')
+        wav = xr.open_dataset(path+'wav.nc')
+        ens = xr.open_dataset(path+'ens.nc')
+        slv = xr.open_dataset(path+'sl.nc')
+
+        self.cal_alg = cfg['cal_alg'].values
+        self.metrics = cfg['metrics'].values
+        self.dt = cfg['dt'].values
+        self.switch_Yini = cfg['switch_Yini'].values
+        self.switch_D = cfg['switch_D'].values
+        self.D50 = cfg['D50'].values
+
+        if self.cal_alg == 'NSGAII':
+            self.n_pop = cfg['n_pop'].values
+            self.generations = cfg['generations'].values
+            self.n_obj = cfg['n_obj'].values
+            self.cal_obj = objective_functions(self.cal_alg, self.metrics, n_pop=self.n_pop, generations=self.generations, n_obj=self.n_obj)
+        else:
+            self.repetitions = cfg['repetitions'].values
+            self.cal_obj = objective_functions(self.cal_alg, self.metrics, repetitions=self.repetitions)
+
+        self.Hs = wav['Hs'].values
+        self.Tp = wav['Tp'].values
+        self.Dir = wav['Dirb'].values
+
+        self.time = mkTime(wav['Y'].values, wav['M'].values, wav['D'].values, wav['h'].values)
+
+        self.Obs = ens['Obs'].values
+        self.time_obs = mkTime(ens['Y'].values, ens['M'].values, ens['D'].values, ens['h'].values)
+
+        self.start_date = datetime(int(cfg['Ysi'].values), int(cfg['Msi'].values), int(cfg['Dsi'].values))
+        self.end_date = datetime(int(cfg['Ysf'].values), int(cfg['Msf'].values), int(cfg['Dsf'].values))
+
+        self.P = self.Hs ** 2 * self.Tp
+        self.ws = wMOORE(self.D50)
+        self.Omega = self.Hb / (self.ws * self.Tp)
+
+        self.split_data()
+
+        if self.switch_Yini == 0:
+            self.Yini = self.Obs_splited[0]
+
+        cfg.close()
+        wav.close()
+        ens.close()
+        mkIdx = np.vectorize(lambda t: np.argmin(np.abs(self.time - t)))
+        self.idx_obs = mkIdx(self.time_obs)
+
+        if self.switch_Yini == 0 and self.switch_D == 0:
+            def model_simulation(par):
+                phi = par['phi']
+                c = par['c']
+                D = 2 * phi
+
+                Ymd, _ = shoreFor(self.P_splited,
+                                    self.Omega_splited,
+                                    self.dt,
+                                    phi,
+                                    c,
+                                    D,
+                                    self.Yini)
+                return Ymd[self.idx_obs_splited]
+
+            self.params = [
+                Uniform('phi', 3, 365),
+                Uniform('c', 1e-6, 1e-2)
+            ]
+            self.model_sim = model_simulation
+
+        elif self.switch_Yini == 1 and self.switch_D == 0:
+            def model_simulation(par):
+                phi = par['phi']
+                c = par['c']
+                D = 2 * phi
+                Yini = par['Yini']
+
+                Ymd, _ = shoreFor(self.P_splited,
+                                    self.Omega_splited,
+                                    self.dt,
+                                    phi,
+                                    c,
+                                    D,
+                                    Yini)
+                return Ymd[self.idx_obs_splited]
+
+            self.params = [
+                Uniform('phi', 3, 365),
+                Uniform('c', 1e-6, 1e-2),
+                Uniform('Yini', np.min(self.Obs), (self.Obs))
+            ]
+            self.model_sim = model_simulation
+
+        elif self.switch_Yini == 0 and self.switch_D == 1:
+            def model_simulation(par):
+                phi = par['phi']
+                c = par['c']
+                D = par['D']
+
+                Ymd, _ = shoreFor(self.P_splited,
+                                    self.Omega_splited,
+                                    self.dt,
+                                    phi,
+                                    c,
+                                    D,
+                                    self.Yini)
+                return Ymd[self.idx_obs_splited]
+
+            self.params = [
+                Uniform('phi', 3, 365),
+                Uniform('c', 1e-6, 1e-2),
+                Uniform('D', 6, 730)
+            ]
+            self.model_sim = model_simulation
+
+        elif self.switch_Yini == 1 and self.switch_D == 1:
+            def model_simulation(par):
+                phi = par['phi']
+                c = par['c']
+                D = par['D']
+                Yini = par['Yini']
+
+                Ymd, _ = shoreFor(self.P_splited,
+                                    self.Omega_splited,
+                                    self.dt,
+                                    phi,
+                                    c,
+                                    D,
+                                    Yini)
+                return Ymd[self.idx_obs_splited]
+
+            self.params = [
+                Uniform('phi', 3, 365),
+                Uniform('c', 1e-6, 1e-2),
+                Uniform('D', 6, 730),
+                Uniform('Yini', np.min(self.Obs), (self.Obs))
+            ]
+            self.model_sim = model_simulation
+
+
+    def split_data(self):
+        """
+        Split the data into calibration and validation datasets.
+        """ 
+
+        idx = np.where((self.time < self.start_date) | (self.time > self.end_date))
+        self.idx_validation = idx
+
+        idx = np.where((self.time >= self.start_date) & (self.time <= self.end_date))
+        self.idx_calibration = idx
+        self.P_splited = self.P[idx]
+        self.Omega_splited = self.Omega[idx]
+        self.wast_splited = self.wast[idx]
+        self.time_splited = self.time[idx]
+
+        idx = np.where((self.time_obs >= self.start_date) & (self.time_obs <= self.end_date))
+
+        self.Obs_splited = self.Obs[idx]
+        self.time_obs_splited = self.time_obs[idx]
+
+        mkIdx = np.vectorize(lambda t: np.argmin(np.abs(self.time_splited - t)))
+        self.idx_obs_splited = mkIdx(self.time_obs_splited)
+        self.observations = self.Obs_splited
+
+        # Validation    
+        idx = np.where((self.time_obs < self.start_date) | (self.time_obs > self.end_date))
+        self.idx_validation_obs = idx
+        mkIdx = np.vectorize(lambda t: np.argmin(np.abs(self.time[self.idx_validation] - t)))
+        self.idx_validation_for_obs = mkIdx(self.time_obs[idx])
+

IHSetShoreFor/shoreFor.py

@@ -0,0 +1,39 @@
+import numpy as np
+from numba import jit
+
+@jit
+def shoreFor(P, Omega, dt, phi=0, c=0, D=0, Sini=0):
+    '''
+    This function calculates the equilibrium profile and the profile evolution   
+    
+    '''
+
+    ii = np.arange(0, (D-1)*24, dt)
+    phivecP = 10 ** (-np.abs(ii) / (phi * 24))
+    IDX = len(phivecP)
+
+    phivecP = np.hstack((np.zeros(IDX), phivecP))
+
+    window = phivecP / np.sum(phivecP)
+
+    OmegaEQ = np.convolve(Omega - np.mean(Omega), window, mode='same') + np.mean(Omega)
+
+    F = (P ** 0.5) * (OmegaEQ - Omega) / np.std(OmegaEQ)
+
+    F[:IDX - 1] = 0
+
+    S = np.full(len(Omega), np.nan)
+
+    rero = F < 0
+    racr = F >= 0
+    S[0] = Sini
+
+    r = np.abs(np.sum(F[racr]) / np.sum(F[rero]))
+
+    r_rero_F = r * rero[1:] * F[1:]
+    racr_F = racr[1:] * F[1:]
+    r_rero_F_prev = r * rero[:-1] * F[:-1]
+    racr_F_prev = racr[:-1] * F[:-1]
+    S[1:] = 0.5 * dt * c * np.cumsum(r_rero_F + racr_F + r_rero_F_prev + racr_F_prev) + S[0]
+
+    return S, OmegaEQ