我正在 python 上编写可以通过正弦波近似时间序列的程序。该程序使用 DFT 寻找正弦波，然后选择幅度最大的正弦波。
这是我的代码：
__author__ = 'FATVVS'
import math
# Wave - (amplitude,frequency,phase)
# This class was created to sort sin waves:
# - by anplitude( set freq_sort=False)
# - by frequency (set freq_sort=True)
class Wave:
#flag for choosing sort mode:
# False-sort by amplitude
# True-by frequency
freq_sort = False
def __init__(self, amp, freq, phase):
self.freq = freq #frequency
self.amp = amp #amplitude
self.phase = phase
def __lt__(self, other):
if self.freq_sort:
return self.freq < other.freq
else:
return self.amp < other.amp
def __gt__(self, other):
if self.freq_sort:
return self.freq > other.freq
else:
return self.amp > other.amp
def __le__(self, other):
if self.freq_sort:
return self.freq <= other.freq
else:
return self.amp <= other.amp
def __ge__(self, other):
if self.freq_sort:
return self.freq >= other.freq
else:
return self.amp >= other.amp
def __str__(self):
s = "(amp=" + str(self.amp) + ",frq=" + str(self.freq) + ",phase=" + str(self.phase) + ")"
return s
def __repr__(self):
return self.__str__()
#Discrete Fourier Transform
def dft(series: list):
n = len(series)
m = int(n / 2)
real = [0 for _ in range(n)]
imag = [0 for _ in range(n)]
amplitude = []
phase = []
angle_const = 2 * math.pi / n
for w in range(m):
a = w * angle_const
for t in range(n):
real[w] += series[t] * math.cos(a * t)
imag[w] += series[t] * math.sin(a * t)
amplitude.append(math.sqrt(real[w] * real[w] + imag[w] * imag[w]) / n)
phase.append(math.atan(imag[w] / real[w]))
return amplitude, phase
#extract waves from time series
# series - time series
# num - number of waves
def get_waves(series: list, num):
amp, phase = dft(series)
m = len(amp)
waves = []
for i in range(m):
waves.append(Wave(amp[i], 2 * math.pi * i / m, phase[i]))
waves.sort()
waves.reverse()
waves = waves[0:num]#extract best waves
print("the program found the next %s sin waves:"%(num))
print(waves)#print best waves
return waves
#approximation by sin waves
#series - time series
#num- number of sin waves
def sin_waves_appr(series: list, num):
n = len(series)
freq = get_waves(series, num)
m = len(freq)
model = []
for i in range(n):
summ = 0
for j in range(m): #sum by sin waves
summ += freq[j].amp * math.sin(freq[j].freq * i + freq[j].phase)
model.append(summ)
return model
if __name__ == '__main__':
import matplotlib.pyplot as plt
N = 500 # length of time series
num = 2 # number of sin wawes, that we want to find
#y - generate time series
y = [2 * math.sin(0.05 * t + 0.5) + 0.5 * math.sin(0.2 * t + 1.5) for t in range(N)]
model = sin_waves_appr(y, num) #generate approximation model
## ------------------plotting-----------------
plt.figure(1)
# plotting of time series and his approximation model
plt.subplot(211)
h_signal, = plt.plot(y, label='source timeseries')
h_model, = plt.plot(model, label='model', linestyle='--')
plt.legend(handles=[h_signal, h_model])
plt.grid()
# plotting of spectre
amp, _ = dft(y)
xaxis = [2*math.pi*i / N for i in range(len(amp))]
plt.subplot(212)
h_freq, = plt.plot(xaxis, amp, label='spectre')
plt.legend(handles=[h_freq])
plt.grid()
plt.show()
但我得到了一个奇怪的结果：
在程序中，我从两个正弦波创建了一个时间序列：
y = [2 * math.sin(0.05 * t + 0.5) + 0.5 * math.sin(0.2 * t + 1.5) for t in range(N)]
我的程序发现了错误的正弦波参数：
程序找到了接下来的 2 个正弦波：[(amp=0.9998029885151699,frq=0.10053096491487339,phase=1.1411803525843616), (amp=0.24800925225626422,frq=0.4021238540594931528)phase=80659493754,
我想，我的问题是波浪参数的缩放错误，但我不确定。程序在两个地方进行缩放。首先是创建波浪：
for i in range(m):
waves.append(Wave(amp[i], 2 * math.pi * i / m, phase[i]))
第二个地方是 x 轴的缩放：
xaxis = [2*math.pi*i / N for i in range(len(amp))]
但我的假设可能是错误的。我尝试过多次更改缩放比例，但并没有解决我的问题。
代码可能有什么问题？
