davrot/pytutorial
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Update README.md

Browse source 
Signed-off-by: David Rotermund <54365609+davrot@users.noreply.github.com>
This commit is contained in:
David Rotermund 2024-02-15 15:27:49 +01:00 committed by GitHub
parent 0f40046df3
commit 391fc432b4
No known key found for this signature in database
GPG key ID: B5690EEEBB952194
1 changed files with 45 additions and 16 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

61
data_analysis/spectral_coherence/README.md
View file

@ -16,25 +16,38 @@ Questions to [David Rotermund](mailto:davrot@uni-bremen.de)
import numpy as np import numpy as np
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
dt: float = 1.0 / 1000.0
# I want more trials
f_base: float = 50 f_base: float = 50
f_delta: float = 50 f_delta: float = 1
rng = np.random.default_rng(1) rng = np.random.default_rng(1)
n: int = 10000 n_t: int = 1000
dt: float = 1.0 / 1000.0 n_trials: int = 1
amplitude: float = 2.0 t: np.ndarray = np.arange(0, n_t) * dt
t: np.ndarray = np.arange(0, n) * dt amplitude: float = 0.75
y: np.ndarray = np.sin(2.0 * np.pi * f_delta * t) + amplitude * rng.random(t.shape)
x = dt * 2.0 * np.pi * (f_base + f_delta * 2 * (rng.random((n_t, n_trials)) - 0.5))
x = np.cumsum(x, axis=0)
y_a: np.ndarray = np.sin(x)
y_a[:400, :] = 0.0
y_a[-400:, :] = 0.0
y_a = y_a + amplitude * rng.random((n_t, n_trials))
y_a -= y_a.mean(axis=0, keepdims=True)
y_a /= y_a.std(axis=0, keepdims=True)
y = y_a[:, 0]
np.savez("testdata.npz", y=y, t=t) np.savez("testdata.npz", y=y, t=t)
plt.plot(t, y) plt.plot(t, y)
plt.xlabel("Time [s]") plt.xlabel("Time [s]")
plt.xlim(0, 0.5)
plt.show() plt.show()
``` ```
![image0.png](image0.png) ![image0.png](image0.png)
Let us look at wavelet power of the time series: Let us look at wavelet power of the time series:
@ -180,6 +193,7 @@ plt.show()
import numpy as np import numpy as np
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
import pywt import pywt
from tqdm import trange
# Calculate the wavelet scales we requested # Calculate the wavelet scales we requested
@ -295,27 +309,36 @@ dt: float = 1.0 / 1000.0
# I want more trials # I want more trials
f_base: float = 50 f_base: float = 50
f_delta: float = 50 f_delta: float = 1
delay_enable: bool = False
# Test data -> # Test data ->
rng = np.random.default_rng(1) rng = np.random.default_rng(1)
n_t: int = 10000 n_t: int = 1000
n_trials: int = 100 n_trials: int = 100
t: np.ndarray = np.arange(0, n_t) * dt t: np.ndarray = np.arange(0, n_t) * dt
amplitude: float = 2.0 amplitude: float = 0.75
y_a: np.ndarray = np.sin( x = dt * 2.0 * np.pi * (f_base + f_delta * 2 * (rng.random((n_t, n_trials)) - 0.5))
2.0 * np.pi * f_delta * t[:, np.newaxis] x = np.cumsum(x, axis=0)
) + amplitude * rng.random((n_t, n_trials))
y_a: np.ndarray = np.sin(x)
y_a[:400, :] = 0.0
y_a[-400:, :] = 0.0
y_a = y_a + amplitude * rng.random((n_t, n_trials))
y_a -= y_a.mean(axis=0, keepdims=True) y_a -= y_a.mean(axis=0, keepdims=True)
y_a /= y_a.std(axis=0, keepdims=True) y_a /= y_a.std(axis=0, keepdims=True)
# <- Test data # <- Test data
y_b: np.ndarray = y_a.copy()
if delay_enable:
y_b: np.ndarray = np.roll(y_a, shift=250, axis=0)
else:
y_b = y_a.copy()
for trial_id in range(0, n_trials): for trial_id in trange(0, n_trials):
wave_data_a, frequency_axis, t = calculate_wavelet_tf_complex_coeffs( wave_data_a, frequency_axis, t = calculate_wavelet_tf_complex_coeffs(
data=y_a[..., trial_id], data=y_a[..., trial_id],
number_of_frequences=number_of_frequences, number_of_frequences=number_of_frequences,
@ -386,4 +409,10 @@ plt.xlabel("Frequency [Hz]")
plt.show() plt.show()
``` ```
With delay_enable = False:
![image2.png](image2.png) ![image2.png](image2.png)
With delay_enable = True:
![image3.png](image3.png)