2016-05-05 11 views
1

私は時間と周波数領域で自家製クワッドの振動をプロットしようとしました。 周波数領域プロットで最高ピークの値を印刷するにはどうすればよいですか?周波数領域プロットの最高ピーク値を出力します

コード:

import matplotlib.pyplot as plt 
import numpy as np 
from scipy import fft, arange 

csv = np.genfromtxt ('/Users/shaunbarney/Desktop/Results/quadOscillations.csv', delimiter=",",dtype=float) 
x = csv[:,0] 
y = csv[:,1] 
x = x - 6318  #Remove start offset 
av=0 
for i in xrange(1097):  #Calculate average sampling time in seconds oscillations 
    if i == 1076: 
     avSampleTime = av/1097000  # 
     break 
    av = av + (x[i+1]-x[i]) 

Fs = 1/avSampleTime #Average sampling freq. 
n = 1079    #no.Samples 
k = arange(n) 
Ts = n/Fs 
frq = k/Ts   #Frequency range two sided 
frq = frq[range(n/2)] #Frequency range one sided 
Y = fft(y)/n   #Fast fourier transfors 
Y = Y[range(n/2)]  #Normalise 

#  PLOTS 

plt.subplot(2,1,1) 
plt.plot(frq,abs(Y),'r') # plotting the spectrum 
plt.xlabel('Freq (Hz)') 
plt.ylabel('|Y(freq)|') 
plt.grid('on') 
plt.subplot(2,1,2) 
plt.plot(x,y) 
plt.xlabel('Time (ms)') 
plt.ylabel('Angle (degrees)') 
plt.grid('on') 
plt.show() 

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

enter image description here

おかげで、 ショーン

答えて

4

あなたがnumpyを使用しているので、ただ単に決定するためにnumpy.maxnumpy.argmaxを使用ピークとピークの位置が表示されます彼はあなたの画面に出ています。この位置を見つけたら、周波数配列にインデックスを付けて最終座標を取得します。あなたのコードを実行するときに、あなたのすべての変数が作成されたと仮定すると

、単純に次の操作を行います。

mY = np.abs(Y) # Find magnitude 
peakY = np.max(mY) # Find max peak 
locY = np.argmax(mY) # Find its location 
frqY = frq[locY] # Get the actual frequency value 

peakYがあなたのグラフとfrqYで最大である大きさの値が含まれていますが、この、その周波数が含まれています最大値(すなわちピーク)は、に位置する。ボーナスとして、あなたはグラフ上のそれを別の色でプロットし、大きいマーカを使ってそれを主要なマグニチュードグラフと区別することができます。 plotの複数の呼び出しを呼び出すと、現在のフォーカスの図の上に追加されます。したがって、スペクトルをプロットし、スペクトルの上にこの点をプロットします。ポイントのサイズをプロットの太さよりも大きくし、ポイントを別の色でマーキングします。この最大のピーク値と対応する場所を反映したタイトルを作成することもできます。使用アレイスライシング:

#  PLOTS 
# New - Find peak of spectrum - Code from above 
mY = np.abs(Y) # Find magnitude 
peakY = np.max(mY) # Find max peak 
locY = np.argmax(mY) # Find its location 
frqY = frq[locY] # Get the actual frequency value 

# Code from before 
plt.subplot(2,1,1) 
plt.plot(frq,abs(Y),'r') # plotting the spectrum 

# New - Plot the max point 
plt.plot(frqY, peakY, 'b.', markersize=18) 

# New - Make title reflecting peak information 
plt.title('Peak value: %f, Location: %f Hz' % (peakY, frqY)) 

# Rest of the code is the same 
plt.xlabel('Freq (Hz)') 
plt.ylabel('|Y(freq)|') 
plt.grid('on') 
plt.subplot(2,1,2) 
plt.plot(x,y) 
plt.xlabel('Time (ms)') 
plt.ylabel('Angle (degrees)') 
plt.grid('on') 
plt.show() 
0
print("maximum of |Y| is: %.4g" % np.max(np.abs(Y))) 

その他の提案を:Y = Y[:n/2+1]をではなく、

はまた、これはあなたの実際の大きさをプロットする前に、単純にこれを行うので、大きさに行う必要があることを覚えておいてくださいY = Y[range(n/2)]。 n個の入力(nが偶数である)を有する実数値データセットのフーリエ変換は、n/2 + 1個の周波数成分を有する。あなたのインデックス作成は最後のポイントを失います。 nが(あなたの場合のように)奇数ならば、それはより扱いにくくなります。

サイドノート:自分のコンピュータにあるファイルに依存しない自問自答の例を質問に付けるとよいでしょう。

関連する問題