前提・実現したいこと
現在、勉強として信号処理プログラムを作成しています。
別のプログラム(試したことのプログラム参照)を参考に、自分なりにコードを書いています。
全く同じグラフが表示されるはずなのですが、表示できないので困っています。
現在の出力結果⤴
軸ラベル等はわざと変更しています。
発生している問題・エラーメッセージ
エラーメッセージ
該当のソースコード
Python
1import numpy as np 2import pandas as pd 3import matplotlib.pyplot as plt 4 5fs=25600 6x=np.arange(0,12800)/fs 7data=np.sin(2.0 * np.pi * 100 * x) + 0.2 * np.random.randn(len(x)) 8 9Fs = 4096 10overlap = 90 11 12#オーバーラップ処理 13#データ長[s]=データ数/サンプリング周波数 14Ts = N/fs 15 16#フレーム周期=フレームサイズ/サンプリング周波数 17Fc = Fs/fs 18 19#オーバーラップ位置=フレーム周期*(1-(オーバーラップ率/100)) 20x_ol = Fs*(1-(overlap/100)) 21 22#平均化回数(分割数) 23N_ave = int((Ts - (Fc * (overlap/100))) / (Fc * (1-(overlap/100)))) 24 25array=[] 26 27for i in range(N_ave): 28 ps=int(x_ol*i) 29 array.append(data[ps:ps+Fs:1]) 30 31#窓関数 32window=np.hamming(Fs) 33acf=1/(sum(window)/Fs) 34 35data_array=[] 36for i in range(N_ave): 37 data_array[i]=data_array[i]*window 38 39#FFT 40Datas=[] 41for i in range(N_ave): 42 43 fft_data=np.fft.fft(data_array[i]) 44 fft_data_abs=abs(fft_data) 45 fft_data_abs_half=(fft_data_abs/(Fs/2)) 46 Data=fft_data_abs_half*acf 47 48Datas.append(Data) 49t=np.arange(0,Fs)/fs 50fq = np.linspace(0,fs,Fs) 51Datas_array=np.array(Datas) 52Datas_mean=np.sqrt(np.mean(Datas_array**2,axis=0)) 53 54fig=plt.figure() 55ax1=fig.add_subplot(211) 56ax2=fig.add_subplot(212) 57#軸ラベル 58ax1.set_xlabel("time[s]",fontsize=14) 59ax1.set_ylabel("Signal_Intensity[dB]",fontsize=14) 60ax2.set_xlabel("frequency[s]",fontsize=14) 61ax2.set_ylabel("Siganl_intensity[dB]",fontsize=14) 62 63#データの範囲 64ax1.set_xlim(0,0.16) 65ax1.set_ylim(-3,3) 66ax2.set_xlim(0,200) 67ax2.set_ylim(0,1) 68 69for i in range(N_ave): 70 ax1.plot(t,data_array[i],label="signal") 71 72ax2.plot(fq[:int(Fs/2)],Datas_mean,label="signal") 73 74plt.tight_layout() 75 76plt.show() 77plt.close()
試したこと
長いですが、こちらのプログラムを参考にしています。
Python
1import numpy as np 2from scipy import signal 3from scipy import fftpack 4 5#オーバーラップ処理 6def ov(data, samplerate, Fs, overlap): 7 Ts = len(data) / samplerate #全データ長 8 Fc = Fs / samplerate #フレーム周期 9 x_ol = Fs * (1 - (overlap/100)) #オーバーラップ時のフレームずらし幅 10 N_ave = int((Ts - (Fc * (overlap/100))) / (Fc * (1-(overlap/100)))) #抽出するフレーム数(平均化に使うデータ個数) 11 12 array = [] #抽出したデータを入れる空配列の定義 13 14 #forループでデータを抽出 15 for i in range(N_ave): 16 ps = int(x_ol * i) #切り出し位置をループ毎に更新 17 array.append(data[ps:ps+Fs:1]) #切り出し位置psからフレームサイズ分抽出して配列に追加 18 return array, N_ave #オーバーラップ抽出されたデータ配列とデータ個数を戻り値にする 19 20#窓関数処理(ハニング窓) 21def hanning(data_array, Fs, N_ave): 22 han = signal.hann(Fs) #ハニング窓作成 23 acf = 1 / (sum(han) / Fs) #振幅補正係数(Amplitude Correction Factor) 24 25 #オーバーラップされた複数時間波形全てに窓関数をかける 26 for i in range(N_ave): 27 data_array[i] = data_array[i] * han #窓関数をかける 28 29 return data_array, acf 30 31#FFT処理 32def fft_ave(data_array,samplerate, Fs, N_ave, acf): 33 fft_array = [] 34 for i in range(N_ave): 35 fft_array.append(acf*np.abs(fftpack.fft(data_array[i])/(Fs/2))) #FFTをして配列に追加、窓関数補正値をかけ、(Fs/2)の正規化を実施。 36 37 fft_axis = np.linspace(0, samplerate, Fs) #周波数軸を作成 38 fft_array = np.array(fft_array) #型をndarrayに変換 39 fft_mean = np.sqrt(np.mean(fft_array ** 2, axis=0)) #全てのFFT波形のパワー平均を計算してから振幅値とする 40return fft_array, fft_mean, fft_axis 41 42import fft_function 43import numpy as np 44from matplotlib import pyplot as plt 45 46samplerate = 25600 47x = np.arange(0, 12800)/samplerate #波形生成のための間軸の作成 48data = np.sin(2.0 * np.pi * 100 * x) + 0.2 * np.random.randn(len(x)) #サイン波にランダム成分を重畳 49 50Fs = 4096 #フレームサイズ 51overlap = 90 #オーバーラップ率 52 53#作成した関数を実行:オーバーラップ抽出された時間波形配列 54time_array, N_ave = fft_function.ov(data, samplerate, Fs, overlap) 55 56#作成した関数を実行:ハニング窓関数をかける 57time_array, acf = fft_function.hanning(time_array, Fs, N_ave) 58 59#作成した関数を実行:FFTをかける 60fft_array, fft_mean, fft_axis = fft_function.fft_ave(time_array, samplerate, Fs, N_ave, acf) 61 62t = np.arange(0, Fs)/samplerate #グラフ描画のためのフレーム時間軸作成 63 64#ここからグラフ描画 65# フォントの種類とサイズを設定する。 66plt.rcParams['font.size'] = 14 67plt.rcParams['font.family'] = 'Times New Roman' 68 69# 目盛を内側にする。 70plt.rcParams['xtick.direction'] = 'in' 71plt.rcParams['ytick.direction'] = 'in' 72 73# グラフの上下左右に目盛線を付ける。 74fig = plt.figure() 75ax1 = fig.add_subplot(211) 76ax1.yaxis.set_ticks_position('both') 77ax1.xaxis.set_ticks_position('both') 78ax2 = fig.add_subplot(212) 79ax2.yaxis.set_ticks_position('both') 80ax2.xaxis.set_ticks_position('both') 81 82# 軸のラベルを設定する。 83ax1.set_xlabel('Time [s]') 84ax1.set_ylabel('Signal [V]') 85ax2.set_xlabel('Frequency [Hz]') 86ax2.set_ylabel('Signal [V]') 87 88#データの範囲と刻み目盛を明示する。 89ax1.set_xticks(np.arange(0, 2, 0.04)) 90ax1.set_yticks(np.arange(-5, 5, 1)) 91ax1.set_xlim(0, 0.16) 92ax1.set_ylim(-3, 3) 93ax2.set_xticks(np.arange(0, samplerate, 50)) 94ax2.set_yticks(np.arange(0, 3, 0.5)) 95ax2.set_xlim(0,200) 96ax2.set_ylim(0, 1) 97 98# データプロットの準備とともに、ラベルと線の太さ、凡例の設置を行う。 99for i in range(N_ave): 100 ax1.plot(t, time_array[i], label='signal', lw=1) 101 102ax2.plot(fft_axis, fft_mean, label='signal', lw=1) 103 104fig.tight_layout() 105 106# グラフを表示する。 107plt.show() 108plt.close()
補足情報(FW/ツールのバージョンなど)
プロットするデータを比較してみてはいかがですか?
「該当のソースコード」を実行したら、「N」が定義されてないので「Ts = N/fs」でエラーが出ます
自分が実行してるのと全く同じコードを提示してください
上のグラフ(ax1)で
ax1.plot(t,data_array[i],label="signal")
でプロットしてる「data_array」は、
data_array=[]
で中身空で作成して、
data_array[i]=data_array[i]*window
で掛け算してるだけだから、中身空なので、グラフに何も表示されなくて当たり前
中身空だから「data_array[i]」とかムリだし (実行したら、実際エラー出る)
コードの先頭の方で、
data=np.sin(...
で三角関数のSinで計算した結果を格納した「data」は、その後forループ内の
array.append(data[ps:ps+Fs:1])
で「array」に入れられてますが、その「array」も「data」も、その後全く使われてません
その代わりに、先に指摘したように、中身が空の「data_array」が登場して、それには何も格納してないのに、それに掛け算して、その後いろいろしてます
中身が空の「data_array」ではなく、「array」に入ってるデータを使って、いろいろしないといけないのではないのですか??