Я пытаюсь найти градиент, под которым строится мой график, в то время как линия соответствует двойной логарифмической шкале. Поэтому я написал следующую функцию:
def calc_coefficients_signal_shift(n: int, N: int, num: int, shift: int, operations: int):
wnss = white_noise_signal_shift(n, N, num, shift, operations)
whitenoise_filtered = np.abs(wnss[1][0:n // 2, :]).copy()
shifted_whitenoise = np.abs(wnss[4][0:n // 2, :]).copy()
x = np.linspace(1, n, n // 2)
original_coefficients = []
shifted_coefficients = []
for i in range(num):
original_coefficients.append(np.polyfit(np.log10(x), np.log10(whitenoise_filtered[:, i]), 1))
shifted_coefficients.append(np.polyfit(np.log10(x), np.log10(shifted_whitenoise[:, i]), 1))
original_coefficients, shifted_coefficients = \
np.asarray((original_coefficients, shifted_coefficients))
fig, (ax1, ax2) = plt.subplots(2, 1, sharey=True, figsize=(10, 7.5))
ax1.plot(whitenoise_filtered)
ax1.loglog(10 ** (original_coefficients.mean() * np.log10(x) + 1), 'r')
ax1.set_title('Original noise')
ax2.loglog(shifted_whitenoise)
ax2.loglog(10 ** (shifted_coefficients.mean() * np.log10(x) + 1), 'r')
ax2.set_title('Shifted noise')
print(original_coefficients.mean(), shifted_coefficients.mean())
Цель функции calc_coefficients_signal_shift - определить, изменяется ли градиент графика после смещения сигнала. Я ожидаю, что он будет где-то около -5/3., поскольку это значение я применил после импорта в функцию white_noise с переменной slope_loglog (фильтрация шума по коэффициенту наклона). Если для числа operations введено значение 0, смещение выполняется, что должно привести к одинаковым массивам для обоих коэффициентов. Однако в результате получается nan для исходного шума и реальное значение для смещенного шума (который в данном случае вообще не смещается, следовательно, первоначальный шум).
Может кто-нибудь сказать мне, кто я что-то не так?
ПРИМЕЧАНИЕ: вы можете предположить, что операция сдвига работает правильно, так как я уже некоторое время отлаживаю эту операцию, и она ведет себя так, как должна. Этот вопрос касается только моей функции подгонки линии.
import matplotlib.pyplot as plt
import numpy as np
import numpy.fft as fft
import numpy.random as rnd
# slope to use for every function
grad = -5 / 3.
rnd.seed(10)
def white_noise(n: int, N: int, slope: int = grad):
x = np.linspace(1, n, n // 2)
slope_loglog = 10 ** (slope * np.log10(x) + 1)
whitenoise = rnd.randn(n // 2, N) + 1j * rnd.randn(n // 2, N)
whitenoise[0, :] = 0 # zero-mean noise
whitenoise_filtered = whitenoise * slope_loglog[:, np.newaxis]
whitenoise = np.concatenate((whitenoise, whitenoise[0:1, :], np.conj(whitenoise[-1:0:-1, :])), axis=0)
whitenoise_filtered = np.concatenate(
(whitenoise_filtered, whitenoise_filtered[0:1, :], np.conj(whitenoise_filtered[-1:0:-1, :])), axis=0)
whitenoise_signal = fft.ifft(whitenoise_filtered, axis=0)
whitenoise_signal = np.real_if_close(whitenoise_signal)
if np.iscomplex(whitenoise_signal).any():
print('Warning! whitenoise_signal is complex-valued!')
whitenoise_retransformed = fft.fft(whitenoise_signal, axis=0)
return whitenoise, whitenoise_filtered, whitenoise_signal, whitenoise_retransformed, slope_loglog
def random_arrays(N: int, num: int):
res = np.asarray(range(N))
rnd.shuffle(res)
return res[:num]
def modified_roll(inp, shift: int, operations: int):
count = 0
array = inp[:]
array_rolled = array.copy()
for k in range(operations):
count += shift
array = np.roll(array, shift, axis=0)
array[:count] = 0
array_rolled += array
out = array_rolled / operations
return out
def white_noise_signal_shift(n: int, N: int, num: int, shift: int, operations: int):
whitenoise, whitenoise_filtered, whitenoise_signal = white_noise(n, N)[:3]
# only showing the selected arrays
arrays_to_select = random_arrays(N, num)
selected_whitenoise = whitenoise[:, arrays_to_select].copy()
selected_whitenoise_filtered = whitenoise_filtered[:, arrays_to_select].copy()
selected_whitenoise_signal = whitenoise_signal[:, arrays_to_select].copy()
# shifting the signal as a field of different refractive index would do
if operations == 0:
shifted_signal = selected_whitenoise_signal
else:
shifted_signal = modified_roll(selected_whitenoise_signal.copy(), shift, operations)
# fourier transform back to the power frequency domain
shifted_whitenoise = fft.fft(shifted_signal, axis=0)
return selected_whitenoise, selected_whitenoise_filtered, selected_whitenoise_signal, shifted_signal, \
shifted_whitenoise
def calc_coefficients_signal_shift(n: int, N: int, num: int, shift: int, operations: int):
wnss = white_noise_signal_shift(n, N, num, shift, operations)
whitenoise_filtered = np.abs(wnss[1][0:n // 2, :]).copy()
shifted_whitenoise = np.abs(wnss[4][0:n // 2, :]).copy()
x = np.linspace(1, n, n // 2)
original_coefficients = []
shifted_coefficients = []
for i in range(num):
original_coefficients.append(np.polyfit(np.log10(x), np.log10(whitenoise_filtered[:, i]), 1))
shifted_coefficients.append(np.polyfit(np.log10(x), np.log10(shifted_whitenoise[:, i]), 1))
original_coefficients, shifted_coefficients = \
np.asarray((original_coefficients, shifted_coefficients))
fig, (ax1, ax2) = plt.subplots(2, 1, sharey=True, figsize=(10, 7.5))
ax1.loglog(whitenoise_filtered)
ax1.loglog(10 ** (original_coefficients.mean() * np.log10(x) + 1), 'r')
ax1.set_title('Original noise')
ax2.loglog(shifted_whitenoise)
ax2.loglog(10 ** (shifted_coefficients.mean() * np.log10(x) + 1), 'r')
ax2.set_title('Shifted noise')
print(original_coefficients.mean(), shifted_coefficients.mean())
calc_coefficients_signal_shift(200, 1, 1, 0, 0)