python中的FFT谱图

Question

我有python3.4。在

我发送了一个2MHz（例如）频率，并在一段时间内接收到空化（直到我停止测量）。 我想得到一个频谱图（空化与频率），更有趣的是空化在次谐波（1MHz）频率下的频谱图。在

数据保存在sdataA（=气蚀）和t（=测量时间）中

我试着在FFTA中保存fft

FFTA = np.array([])
FFTA = np.fft.fft(dataA)
FFTA = np.append(FFTA, dataA)

我有实数和复数 然后我只取了一半（从0到1MHz），保存了真实而复杂的数据。在

^{pr2}$

我试图通过以下方法获得频率：

FFTAfreqs = np.fft.fftfreq(len(yAre))

但它完全错了（我用print (FFTAfreqs)打印了数据）

我还绘制了数据，但还是错了：

plt.plot(t, FFTA[range(n)].real, 'b-', t, FFTA[range(n)].imag, 'r--')
plt.legend(('real', 'imaginary'))
plt.show()

我怎样才能输出在次谐波（1MHz）频率下的空化频谱图？在

编辑：

数据示例：

查看“dataA”和“time”示例：

dataA = [6.08E-04,2.78E-04,3.64E-04,3.64E-04,4.37E-04,4.09E-04,4.49E-04,4.09E-04,3.52E-04,3.24E-04,3.92E-04,3.24E-04,2.67E-04,3.24E-04,2.95E-04,2.95E-04,4.94E-04,4.09E-04,3.64E-04,3.07E-04]
time = [0.00E+00,4.96E-07,9.92E-07,1.49E-06,1.98E-06,2.48E-06,2.98E-06,3.47E-06,3.97E-06,4.46E-06,4.96E-06,5.46E-06,5.95E-06,6.45E-06,6.94E-06,7.44E-06,7.94E-06,8.43E-06,8.93E-06,9.42E-06]

编辑二： 在@Martin示例中，我尝试了以下代码，请告诉我是否正确。在

如果数据和时间保存为h5文件（或我已经发布的数据）

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd

dfdata = pd.read_hdf("C:\\data_python\\DataA.h5")
dft = pd.read_hdf("C:\\data_python\\time.h5")
dft_cor = int((len(dft)-2)*4.96E-6)   # calculating the measured time

fs = 2000000 #sampling frequency 2MHz
CHUNK = 10000
signal_time = dft_cor # seconds

def sine(freq,fs,secs):
    data=dfdata
    wave = np.sin(freq*2*np.pi*data)
    return wave

a1 = sine(fs,fs,120)
a2 = sine(fs/2,fs,120)

signal = a1+a2
afft = np.abs(np.fft.fft(signal[0:CHUNK]))
freqs = np.linspace(0,fs,CHUNK)[0:int(fs/2)]
spectrogram_chunk = freqs/np.amax(freqs*1.0)

# Plot spectral analysis
plt.plot(freqs[0:1000000],afft[0:1000000]) # 0-1MHz
plt.show()

number_of_chunks = 1000

# Empty spectrogram
Spectrogram = np.zeros(shape = [CHUNK,number_of_chunks])

for i in range(number_of_chunks):
    afft = np.abs(np.fft.fft(signal[i*CHUNK:(1+i)*CHUNK]))
    freqs = np.linspace(0,fs,CHUNK)[0:int(fs/2)]
    spectrogram_chunk = afft/np.amax(afft*1.0)
    try:
        Spectrogram[:,i]=spectrogram_chunk
    except:
        break


import cv2
Spectrogram = Spectrogram[0:1000000,:]
cv2.imshow('spectrogram',np.uint8(255*Spectrogram/np.amax(Spectrogram)))
cv2.waitKey()
cv2.destroyAllWindows()

Answer 1

似乎你的问题不在Python中，而是在理解什么是谱图。在

---

谱图是对信号进行谱分析的序列。在

1）你需要把信号分成块。在

2）对这些块进行光谱分析并将其粘在一起。在

示例：

您有1秒的音频录制（44100赫兹采样）。这意味着录音将有1s*44100->；44100个样本。定义块大小=1024（例如）。在

对于每个块，您将执行FFT，并将其粘贴到2D矩阵中（X轴-块的FFT，Y轴-块编号，）。44100个样本/块~44个FFT，每个FFT覆盖1024/44100~0.023秒的信号

块越大，光谱图越精确，但“实时性”越差。在

数据块越小，频谱图的准确率越低，但是当你“更频繁”测量频率时，你会有更多的测量值。在

如果你需要1MHZ——实际上，你不能使用任何高于1MHZ的东西，你只需要得到一半的FFT数组——而哪一半并不重要，因为1MHZ只是采样频率的一半，而FFT镜像的是高于采样频率1/2的任何东西。在

---

关于FFT，你不需要复数。你想做的

快速傅里叶变换=np.abs（FFT）#Edit-我刚注意到你用了'.real'，但我会把它保留在这里

因为你想要实数。在

谱图的制备谱图示例

150HZ和300HZ的音频信号

import numpy as np
import matplotlib.pyplot as plt


fs = 44100#sampling frequency
CHUNK = 10000
signal_time = 20 # seconds

def sine(freq,fs,secs):
    data=np.arange(fs*secs)/(fs*1.0)
    wave = np.sin(freq*2*np.pi*data)
    return wave

a1 = sine(150,fs,120)
a2 = sine(300,fs,120)

signal = a1+a2
afft = np.abs(np.fft.fft(signal[0:CHUNK]))
freqs = np.linspace(0,fs,CHUNK)[0:int(fs/2)]
spectrogram_chunk = freqs/np.amax(freqs*1.0)

# Plot spectral analysis
plt.plot(freqs[0:250],afft[0:250])
plt.show()

number_of_chunks = 1000

# Empty spectrogram
Spectrogram = np.zeros(shape = [CHUNK,number_of_chunks])

for i in range(number_of_chunks):
    afft = np.abs(np.fft.fft(signal[i*CHUNK:(1+i)*CHUNK]))
    freqs = np.linspace(0,fs,CHUNK)[0:int(fs/2)]
    #plt.plot(spectrogram_chunk[0:250],afft[0:250])
    #plt.show()
    spectrogram_chunk = afft/np.amax(afft*1.0)
    #print(signal[i*CHUNK:(1+i)*CHUNK].shape)
    try:
        Spectrogram[:,i]=spectrogram_chunk
    except:
        break


import cv2
Spectrogram = Spectrogram[0:250,:]
cv2.imshow('spectrogram',np.uint8(255*Spectrogram/np.amax(Spectrogram)))
cv2.waitKey()
cv2.destroyAllWindows()

单块谱分析

谱图