python-3.x - 用python绘制两个傅立叶系数

Question

我正在使用 python 中的离散傅立叶变换。在我的 python 代码中，如何在复平面上绘制两个傅立叶系数。

我已经看到在 matlab 中他们使用以下代码来执行此操作，他们使用 dsearchn 来绘制它，如下所示：

% create the signal
srate  = 1000; % hz
time   = 0:1/srate:2; % time vector in seconds
pnts   = length(time); % number of time points
signal = 2.5 * sin( 2*pi*4*time ) ...
       + 1.5 * sin( 2*pi*6.5*time );

% prepare the Fourier transform
fourTime = (0:pnts-1)/pnts;
fCoefs   = zeros(size(signal));

% compute frequencies vector
hz = linspace(0,srate/2,floor(pnts/2)+1);

%% plot two Fourier coefficients

coefs2plot = dsearchn(hz',[4 4.5]');

% extract magnitude and angle
mag = abs(fCoefs(coefs2plot));
phs = angle(fCoefs(coefs2plot));

figure(2), clf
plot( real(fCoefs(coefs2plot)) , imag(fCoefs(coefs2plot)) ,'o','linew',2,'markersize',10,'markerfacecolor','r');

% make plot look nicer
axislims = max(mag)*1.1;
set(gca,'xlim',[-1 1]*axislims,'ylim',[-1 1]*axislims)
grid on, hold on, axis square
plot(get(gca,'xlim'),[0 0],'k','linew',2)
plot([0 0],get(gca,'ylim'),'k','linew',2)
xlabel('Real axis')
ylabel('Imaginary axis')
title('Complex plane')

有人告诉我，我可以使用 scipy: scipy.spatial.cKDTree中的包，但我不知道如何在 python 代码中以 matlab 代码中的示例来实现它。谁能帮我

提前致谢

标签： python-3.xmatlabscipy

Answer 1

下面的绘图命令可能不是很优雅，但我希望我能得到你想要的。'dsearchn' 的用法与傅里叶变换无关，而是在寻找大特征。如果你像我一样绘制整个光谱，你可以直观地找到这些特征。在图像 2-6 中，y 轴是标题，x 轴是频率，以 Hz 为单位。

%matplotlib inline
import numpy as np
import matplotlib.pyplot as p
from ipywidgets import *

from numpy.fft import * 

# create the signal
srate  = 1000 # Hz   
numsec=1 
time   = np.linspace(0,numsec,srate*numsec+1) # time vector in seconds

signal = 2.5 * np.sin( 2*np.pi*4*time ) + 1.5 * np.sin( 2*np.pi*6.5*time );

p.figure(figsize=(14,6))
p.subplot(241)
p.plot(time,signal,'.-',lw=1,ms=2)
p.title('signal')

#  compute frequencies 
fourier = np.fft.fft(signal)
n = signal.size
timestep = 1/srate
freq = np.fft.fftfreq(n, d=timestep)

p.subplot(242)
p.plot(freq,np.abs(fourier),'.-',ms=2,lw=0.5)
p.plot(freq[4],np.abs(fourier[4]),'.-',ms=15 )
p.plot(freq[6],np.abs(fourier[6]),'.-',ms=15 )
p.title('full spectrum (abs)')

p.subplot(243)
p.plot(freq[:20],np.abs(fourier)[:20],'.-')
p.plot(freq[4],np.abs(fourier[4]),'.-',ms=15 )
p.plot(freq[6],np.abs(fourier[6]),'.-',ms=15 )
p.title('zoomed in (abs)')

p.subplot(244)
p.plot(freq[:20],np.angle(fourier)[:20],'.-')
p.plot(freq[4],np.angle(fourier[4]),'.-',ms=15 )
p.plot(freq[6],np.angle(fourier[6]),'.-',ms=15 )
p.title('zoomed in (phase)')

p.subplot(245)
p.plot(freq[:20],np.real(fourier)[:20],'.-')
p.plot(freq[4],np.real(fourier[4]),'.-',ms=15 )
p.plot(freq[6],np.real(fourier[6]),'.-',ms=15 )
p.title('zoomed in (real)')

p.subplot(246)
p.plot(freq[:20],np.imag(fourier)[:20],'.-')
p.plot(freq[4],np.imag(fourier[4]),'.-',ms=15 )
p.plot(freq[6],np.imag(fourier[6]),'.-',ms=15 )
p.title('zoomed in (imag)')

p.subplot(248)
p.plot(np.real(fourier)[:50],np.imag(fourier)[:50],'.')
p.plot(np.real(fourier[4]),np.imag(fourier[4]),'.',ms=15 )
p.plot(np.real(fourier[6]),np.imag(fourier[6]),'.',ms=15 )
p.xlabel('real axis')
p.ylabel('imag.axis')
p.title('complex plane')