мои спектры мощности правдоподобны? сравнение между lomb-scargle и fft (scipy.signal и numpy.fft)
Может ли кто-нибудь любезно указать, почему я получаю совсем другие результаты?
Есть много пиков, которые не должны появляться. Фактически, должен быть только один пик.
Я новичок в Python, и все комментарии о моем коде ниже приветствуются.
Тестовые данные здесь.введите описание ссылки здесь
Вы можете напрямуюwget https://clbin.com/YJkwr
.
Его первый столбец - время прибытия серии фотонов. Источник света испускает фотоны случайным образом. Общее время 55736 секунд и 67398 фотонов. Я пытаюсь обнаружить некоторую периодичность интенсивности света.
Мы можем сопоставить время и интенсивность света пропорционально количеству фотонов в каждом отсчете времени.
Я попробовал numpy.fft и lomb-scargle of scipy.signal, чтобы получить спектр мощности, но получил очень разные результаты.
быстрое преобразование Фурье import pylab as pl
import numpy as np
timepoints=np.loadtxt('timesequence',usecols=(0,),unpack=True,delimiter=",")
binshu=50000
interd=54425./binshu
t=np.histogram(timepoints,bins=binshu)[0]
sp = np.fft.fft(t)
freq = np.fft.fftfreq(len(t),d=interd)
freqnum = np.fft.fftfreq(len(t),d=interd).argsort()
pl.xlabel("frequency(Hz)")
pl.plot(freq[freqnum],np.abs(sp)[freqnum])
Lomb-scargle timepoints=np.loadtxt('timesequence',usecols=(0,),unpack=True,delimiter=",")
binshu=50000
intensity=np.histogram(timepoints,bins=binshu)[0].astype('float64')
middletime=np.histogram(timepoints,bins=binshu)[1][1:binshu+1]-np.histogram(timepoints,bins=binshu)[1][3]*0.5
freq1=1./(timepoints.max()-timepoints.min())
freq2=freq1*len(timepoints)
freqs=np.linspace(freq1,freq2,1000)
pm1=spectral.lombscargle(middletime,intensity,freqs)
pl.xlabel("frequency(Hz)")
pl.plot(freqs,pm1)
pl.xlim(0.5*freq1,2*freq2)
pl.ylim(0,250)
pl.show()
*********************************
Спасибо, Уоррен, я ценю это. Спасибо за подробный ответ.
Вы правы, время интеграции есть, здесь оно около 1.7с.
Есть много одиночных воздействий. Каждая экспозиция стоит 1,7 с.
В одной экспозиции мы не можем точно сказать время ее прибытия.
Если временные ряды похожи на:
0 1,7 3,4 8,5 8,5
Время интегрирования последних двух фотонов1.7s
,не(8.5-3.4)s
Так что я буду пересматривать часть вашего кода.
ОДНАКО мой вопрос все еще остается. Вы настраиваете несколько параметров, чтобы получить0.024Hz
Пик в пике ломаного пика до некоторой степени. И вы используете это, чтобы вести ваши параметры в FFT.
Если вы не знаете номер0.024
, может быть, вы можете использовать разные параметры, чтобы получить другой самый высокий пик?
Как гарантировать каждый раз, когда мы можем получить правоnum_ls_freqs
? Вы можете увидеть, если мы выбираем разныеnum_ls_freqs
, самый высокий пик сдвигов.
Если у меня много временных рядов, каждый раз мне приходится указывать разные параметры? И как их получить?