Я работаю с данными об аномалиях уровня моря, которые представляют собой трехмерную матрицу, состоящую из времени, широты и долготы. И я сделал 1D дискретное преобразование Фурье оси времени, чтобы создать спектр мощности. Из теоремы Парсеваля дисперсия должна быть равна спектру мощности, суммируемому по всем возможным частотам.
Дисперсия определяется как:
var = (1 / # образцов) * сумма (данные ** 2)
В пространственном отношении я получаю те же результаты, но по величине моя дисперсия не равна величине спектра мощности ... прикрепив картинку и мой код ниже. Почему это может быть так, что величины отличаются друг от друга?
import numpy as np
from matplotlib import pyplot as plt
import netCDF4 as s
from scipy.fftpack import fft, ifft
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from scipy import signal
from scipy import stats
Yearmin = 2018
Yearmax = 2019
year_len = Yearmax - Yearmin + 1.0 # number of years
direcInput = "filepath"
a = s.Dataset(direcInput+"test.nc", mode='r')
#creating arrays
lat = a.variables["latitude"][:]
lon = a.variables["longitude"][:]
time1 = a.variables["time"][:] #DAYS SINCE JAN 1ST 1950
sla = a.variables["sla"][:,:,:] #t, y, x
time = Yearmin + (year_len * (time1 - np.min(time1)) / ( np.max(time1) - np.min(time1))) #TIME IN 2018,2019 DECIMALS
#get fft over first axis
sla_fft = np.fft.rfft(tapered_sla, axis=0)
spectrum = np.abs(sla_fft)**2
t = time1
freq_array = np.fft.rfftfreq(len(t), d=(t[1]-t[0])) #d = 1 days, step size spacing
sla_taper_square = np.square(tapered_sla)
variance_tapered = (np.sum(sla_taper_square, axis=0))/sla.shape[0]#/sla.shape[0] #/tapered_sla.shape[0] #dividing by M
spec_sum = (1/2*np.pi)*np.sum(spectrum, axis=0)
log_spec_sum = np.log(spec_sum)
fig = plt.figure(figsize=(15,5))
ax = plt.axes(projection=ccrs.PlateCarree())
x = ax.contourf(lon,lat,variance_tapered, 20, extend = 'both', cmap = 'RdBu_r')
ax.coastlines()
gridlines = ax.gridlines(draw_labels=True)
cbar = plt.colorbar(x, fraction=.046, pad=0.04)
cbar.set_label('Variance (m$^2$/#samples)', labelpad=15, y=.5, rotation=270)
plt.title('Variance of tapered SLA Data', y=1.1)
plt.show()
fig = plt.figure(figsize=(15,5))
ax = plt.axes(projection=ccrs.PlateCarree())
x = ax.contourf(lon,lat,spec_sum, 20, extend = 'both', cmap = 'RdBu_r')
ax.coastlines()
gridlines = ax.gridlines(draw_labels=True)
cbar = plt.colorbar(x, fraction=.046, pad=0.04)
cbar.set_label('Power spectrum [m$^2$ per day]', labelpad=15, y=.5, rotation=270)
plt.title('Sum of power spectrum over frequencies', y=1.1)
plt.show()
