周波数があり、Pythonでプロットシヌス波が必要です

Question

私はちょうど洞の波の変調で戦っています。 私は周波数を取得しました（時間が経つにつれて変化するデータ）から、対応する周波数の副波をプロットする必要があります。

青い線は、本物のデータ点をプロットし、緑色のは、私が今までやったことですが、それがすべてで、実際のデータとcorespondしません。

正弦波をプロットするためのコードは、底部である。

def plotmodulsin(): 
    n = 530 
    f1, f2 = 16, 50 # frequency 

    t = linspace(6.94,8.2,530) 
    dt = t[1] - t[0] # needed for integration 
    print t[1] 
    print t[0] 
    f_inst = logspace(log10(f1), log10(f2), n) 
    phi = 2 * pi * cumsum(f_inst) * dt # integrate to get phase 
    pylab.plot(t, 5*sin(phi)) 

振幅ベクトル：

[2.64、-2.64、6.14、-6.14、9.56、-9.56、12.57、-12.57、 15.55、-15.55、18.04、-18.04、21.17、-21.17、23.34、-23.34、25.86、-25.86、28.03、-28.03、30.49、-30.49、33.28、-33.28、35.36、-35.36、36.47、-36.47、 38.86、-38.86、41.49、-41.49、42.91、-42.91、44.41、-44.41、45.98、-45.98、47.63、-47.63、47.63、-47.63、51.23、-51.23、51.23、-51.23、53.18、-53.18、 55.24、-55.24、55.24、-55.24、55.24、-55.24、57.43、-57.43、57.43、-57 0.43、59.75、-59.75、59.75、-59.75、59.75、-59.75、59.75、-59.75、62.22、-62.22、59.75、-59.75、62.22、-62.22、59.75、-59.75、62.22、-62.22、62.22、 -62.22、59.75、-59.75、62.22、-62.22、62.22、-62.22、59.75、-59.75、62.22、-62.22、62.22、-62.22、62.22、-62.22、59.75、-59.75、62.22、-62.22、59.75、実際のデータのため-59.75、62.22、-62.22、59.75、-59.75、59.75]

時間ベクトル：

[6.954、6.985、7.016、7.041、7.066、7.088、7.11、7.13、 7.149、7.167、7.186、7.202、7.219、7.235、7.251、7.266、7.282、7.296、7.311、7.325、7.339、7.352、7.366、7.379、7.392、7.404、7.417、7.43、7.442、7.454、7.466、7.478、7.49、 7.501、7.513、7.524、7.536、7.547、7.558、 7.569、7.58、7.591、7.602、7.613、7.624、7.634、7.645、7.655、7.666、7.676、7.686、7.697、7.707、7.717、7.728、7.738、7.748、7.758、7.768、7.778、7.788、7.798、7.808、7.818、 7.828、7.838、7.848、7.858、7.868、7.877、7.887、7.897、7.907、7.917、7.927、7.937、7.946、7.956、7.966、7.976、7.986、7.996、8.006、8.016、8.026、8.035、8.045、8.055、8.065、 8.075、8.084、8.094、8.104、8.114、8.124、8.134、8.144、8.154、8.164、8.174、8.184、8.194、8.20]

だから私は一定の振幅と、以下の周波数で洞を生成する必要があります。

[10.5、16.03、20.0、22.94、25.51、27.47、29.76、31.25、32.89、34.25、35.71、37.31、38.46、39.06、40.32、41.67、42.37、4 3.1、43.86、44.64、44.64、46.3、46.3、47.17、48.08、48.08、48.08、49.02、49.02、50.0、50.0、50.0、50.0]

Answer 1

あなたが何かで機能と一致するように試みることができる正弦かデータから周波数と振幅の推定値を抽出することによって、実際に余弦のようなものになります。私があなたを正しく理解していれば、あなたのデータは最大値と最小値なので、それと似た三角関数が必要です。データが2つの配列timeとvalueに保存されている場合、振幅の見積もりは単にnp.abs(value)で与えられます。周波数は、最大値と最小値との時間差の2倍の逆数として与えられる。 freq = 0.5/(time[1:]-time[:-1])は、各時間間隔の中点の頻度推定値を提供します。したがって、対応する時間は、freqTimes = (time[1:]+time[:-1])/2.として与えられる。

より滑らかなカーブを得るために、振幅値と周波数値を補間して、その間の値の推定値を得ることができます。これを行う簡単な方法は、np.interpの使用です。単純な線形補間を行います。補間するタイミングを指定する必要があります。私たちは、そのための配列を作成し、その後で補間します：

phi = (2*np.pi*np.cumsum(freqInterpolated) 
     *(timesToInterpolate[1]-timesToInterpolate[0])) 

そして今、あなたはプロットできます：今、あなたはあなたがすでに実行して、あなたの例であったこと、統合を行うことができます

n = 10000 
timesToInterpolate = np.linspace(time[0], time[-1], n, endpoint=True) 
freqInterpolated = np.interp(timesToInterpolate, freqTimes, freq) 
amplInterpolated = np.interp(timesToInterpolate, time, np.abs(value)) 

。だから、一緒にすべてを置くことは、あなたに与えられる：（あなたが完全な配列が含まれている場合）

import numpy as np 
import matplotlib.pyplot as plt 

time = np.array([6.954, 6.985, 7.016, 7.041, 7.066, 7.088, 7.11, 7.13]) #... 
value = np.array([2.64, -2.64, 6.14, -6.14, 9.56, -9.56, 12.57, -12.57]) #... 

freq = 0.5/(time[1:]-time[:-1]) 
freqTimes = (time[1:]+time[:-1])/2. 

n = 10000 
timesToInterpolate = np.linspace(time[0], time[-1], n, endpoint=True) 
freqInterpolated = np.interp(timesToInterpolate, freqTimes, freq) 
amplInterpolated = np.interp(timesToInterpolate, time, np.abs(value)) 

phi = (2*np.pi*np.cumsum(freqInterpolated) 
     *(timesToInterpolate[1]-timesToInterpolate[0])) 

plt.plot(time, value) 
plt.plot(timesToInterpolate, amplInterpolated*np.cos(phi)) #or np.sin(phi+np.pi/2) 
plt.show() 

結果は次のようになります。