From c4a6313995363a6de3176e383c6f43d139c55e7b Mon Sep 17 00:00:00 2001
From: MaNyh <nyhuis@tfd.uni-hannover.de>
Date: Mon, 28 Mar 2022 17:59:11 +0200
Subject: [PATCH] better fit

---
 .../timeseries/stationary_signal_check.py     | 29 +++++++-----
 tests/test_ntrfc_postprocessing.py            | 45 +++++++++----------
 2 files changed, 40 insertions(+), 34 deletions(-)

diff --git a/ntrfc/postprocessing/timeseries/stationary_signal_check.py b/ntrfc/postprocessing/timeseries/stationary_signal_check.py
index e9db44f0..2aca1150 100644
--- a/ntrfc/postprocessing/timeseries/stationary_signal_check.py
+++ b/ntrfc/postprocessing/timeseries/stationary_signal_check.py
@@ -1,5 +1,7 @@
 import numpy as np
 from matplotlib import pyplot as plt
+import matplotlib
+#matplotlib.use('TkAgg')
 
 from ntrfc.postprocessing.timeseries.integral_scales import integralscales
 
@@ -22,7 +24,7 @@ def parsed_timeseries_analysis(timesteps, signal, resolvechunks=20, verbose=True
     if stationarity==True:
         csig=signal
         ctime=timesteps
-        plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps)
+        plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps,signal_type)
 
         return stationarity, timescale, stationarity_timestep
 
@@ -30,12 +32,13 @@ def parsed_timeseries_analysis(timesteps, signal, resolvechunks=20, verbose=True
         csig = np.concatenate([*checksigchunks[resolvechunks - i:]])
         ctime = np.concatenate([*checktimechunks[resolvechunks - i:]])
 
-        signal_type,newstationarity,newstationarity_timestep, newtimescale,newlengthscale = check_signal_stationarity(resolvechunks, csig, ctime)
+        new_signal_type,newstationarity,newstationarity_timestep, newtimescale,newlengthscale = check_signal_stationarity(resolvechunks, csig, ctime)
 
         if newstationarity:
             newscales = (timescale, lengthscale)
             stationarity=True
             scales = newscales
+            signal_type=new_signal_type
             stationarity_timestep = newstationarity_timestep
 
 
@@ -43,15 +46,15 @@ def parsed_timeseries_analysis(timesteps, signal, resolvechunks=20, verbose=True
             # when no further stationarity found, return status
             # when done, return last status
 
-            plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps)
+            plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps,signal_type)
 
             return stationarity, scales, stationarity_timestep
 
-    plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps)
+    plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps,signal_type)
     return stationarity, scales, stationarity_timestep
 
 
-def plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps):
+def plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, timesteps, signaltype):
     plt.figure()
     plt.plot(timesteps, signal)
     plt.plot(ctime, csig, color="black", linewidth=4)
@@ -61,6 +64,10 @@ def plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, times
     if ymax-ymin<0.01:
         ymax=0.5+np.mean(signal)
         ymin =-0.5+np.mean(signal)
+    else:
+        ymin -=1
+        ymax += 1
+
     if sts==0.0 or sts:
         plt.vlines(sts, ymin=ymin, ymax=ymax,
                    linewidth=4, color="k", linestyles="dashed")
@@ -69,9 +76,9 @@ def plot_stationarity_analisys(csig, ctime, signal, stationarity_timestep, times
     else:
         plt.axvspan(0, timesteps[-1],
                     facecolor='red', alpha=0.5)
-
+    plt.title(signaltype)
     plt.ylim(ymin, ymax)
-    plt.show()  #
+    plt.show()
 
 
 def check_signal_stationarity(resolvechunks, signal, timesteps, verbose = True):
@@ -93,9 +100,9 @@ def check_signal_stationarity(resolvechunks, signal, timesteps, verbose = True):
     # todo: now it is only mean, val and var that is being investigated.
     # it makes sense to also investigate the behaviour of the autocorrelation
     # but as the signal is divided into chunks, one has to
-    const_mean = np.allclose(mean, means,rtol=0.05)
-    const_val = np.allclose(mean, signal,rtol=0.05)
-    const_var = np.allclose(var,vars,rtol=0.4)
+    const_mean = np.allclose(mean, means,rtol=0.02)
+    const_val = np.allclose(mean, signal,rtol=0.02)
+    const_var = np.allclose(var,vars,rtol=0.02)
     #
     if const_mean and const_var:
         timescale, lengthscale = integralscales(signal, timesteps)
@@ -107,7 +114,7 @@ def check_signal_stationarity(resolvechunks, signal, timesteps, verbose = True):
         timescale, lengthscale = 0,0
         timescales, lengthscales = np.zeros(resolvechunks),np.zeros(resolvechunks)
 
-    const_tscales = np.allclose(timescales, timescale,rtol=0.1)
+    const_tscales = np.allclose(timescales, timescale,rtol=0.02)
 
     """
     from itertools import product
diff --git a/tests/test_ntrfc_postprocessing.py b/tests/test_ntrfc_postprocessing.py
index 9c4fe577..bd6f2abd 100644
--- a/tests/test_ntrfc_postprocessing.py
+++ b/tests/test_ntrfc_postprocessing.py
@@ -31,7 +31,7 @@ def sine_signal(numtimesteps, time):
 
 
 def noise_signal(numtimesteps, time):
-    mu, sigma = 0, 0.2  # mean and standard deviation
+    mu, sigma = 0, 0.02  # mean and standard deviation
     s = np.random.normal(mu, sigma, size=numtimesteps)
 
     Dt = time / numtimesteps
@@ -74,14 +74,14 @@ def tanh_sin_signal(numtimesteps, time):
 
 def tanh_sin_noise_signal(numtimesteps, time):
     timesteps = np.linspace(0, time, numtimesteps)
-    stationary = time/4
+    stationary = time/10
     tanhsignal = np.tanh(timesteps * stationary ** -1)
 
     timesteps = np.linspace(0, 1, numtimesteps) * time
-    omega = 1000
-    signal = np.sin(timesteps * omega) + 1
+    omega = 100
+    signal = np.sin(timesteps * omega)*0.3 + 1
 
-    mu, sigma = 0, 0.2  # mean and standard deviation
+    mu, sigma = 0, 0.1  # mean and standard deviation
     s = np.random.normal(mu, sigma, size=numtimesteps)
 
     Dt = time / numtimesteps
@@ -89,7 +89,7 @@ def tanh_sin_noise_signal(numtimesteps, time):
     dt = time / 1000
     timescale = int(dt / Dt)
     weights = np.repeat(1.0, timescale) / timescale
-    noise = np.convolve(s, weights, 'same') + 1
+    noise = np.convolve(s, weights, 'same')
 
     omega = 1000
     timesteps = np.linspace(0, time, numtimesteps)
@@ -111,11 +111,11 @@ def sine_abate_signal(numtimesteps, time):
 
 def complex_signal(numtimesteps, time):
     timesteps = np.linspace(0, time, numtimesteps)
-    stationary = time/4
+    stationary = time/10
     tanhsignal = np.tanh(timesteps * stationary ** -1)
 
     timesteps = np.linspace(0, 1, numtimesteps) * time
-    omega = 1000
+    omega = 10
     signal = np.sin(timesteps * omega) + 1
 
     mu, sigma = 0, 0.2  # mean and standard deviation
@@ -145,32 +145,31 @@ def complex_signal(numtimesteps, time):
 
 
 def test_analize_stationarity(verbose=True):
-    tight_tolerance = 1e-2
     loose_tolerance = 1
     signals = {"constant": constant_signal(1000, 1),
                "ramp": ramp_signal(1000, 1),
                "rampconst": ramptoconst_signal(10000, 1),
                "sine": sine_signal(10000, 1),
-               "noise": noise_signal(10000, 1),
-               "convergentsine_signal": convergentsine_signal(10000, 10),
-               "tanh":tanh_signal(10000,20),
-               "tanh_sin":tanh_sin_signal(10000,10),
-               "tanh_sin_noise": tanh_sin_noise_signal(10000, 2.8),
-               "sine_abate":sine_abate_signal(20000, 4),
-               "complex":complex_signal(40000,60)
+               "noise": noise_signal(100000, 10),
+               "convergentsine_signal": convergentsine_signal(10000, 40),
+               "tanh":tanh_signal(20000,20),
+               "tanh_sin":tanh_sin_signal(20000,32),
+               "tanh_sin_noise": tanh_sin_noise_signal(40000, 480),
+               "sine_abate":sine_abate_signal(30000, 4),
+               "complex":complex_signal(50000,480)
                }
 
     expections = {"constant": (True, (0, 0), [0,0]),
                   "ramp": (False, (0, 0), [-1,-1]),
                   "rampconst": (True, (0, 0), [0.45,0.55]),
-                  "sine": (True, (0.0005, 0.0005), [0,0]),
-                  "noise": (True, (0.0005, 0.0005), [0,0]),
-                  "convergentsine_signal": (True, (0.0005, 0.0005), [5,5.5]),
-                  "tanh":(True,(0,0),[3.5,4.5]),
-                  "tanh_sin":(True,(0.0005,0.0005),[4,4.75]),
-                  "tanh_sin_noise":(True,(0.0005,0.0005),[0.8,1.2]),
+                  "sine": (True, (0.005, 0.005), [0,0]),
+                  "noise": (True, (0.005, 0.005), [0,0]),
+                  "convergentsine_signal": (True, (0.0005, 0.0005), [15.5,16.5]),
+                  "tanh":(True,(0,0),[4.5,5.5]),
+                  "tanh_sin":(True,(0.0005,0.0005),[15.5,16.5]),
+                  "tanh_sin_noise":(True,(0.0005,0.0005),[110,130]),
                   "sine_abate":(True,(0.00025,0.00025),[0.35,1]),
-                  "complex":(True,(0.00025,0.0025),[20,26])
+                  "complex":(True,(0.00025,0.0025),[90,130])
                   }
 
     for name, sig in signals.items():
-- 
GitLab