低通滤波器实现正确还是错误?

Low pass filter implementation Correct or wrong?

我一直在研究 2 个测量旋转轴振动的传感器信号。由于信号中存在残余噪声。我试图通过去除趋势、零填充和应用低通滤波器来过滤它。下面我附上了滤波前后的信号图。过滤后的信号变化很大,这让我思考我是否真的以正确的方式进行操作。

我的 Matlab 代码是

X = xlsread(filename,'F:F');
Y = xlsread(filename,'G:G');

%Calculate frequency axis
fs = 1e6  ; % Sampling frequency (Hz)
NFFT = 2^nextpow2(length(X));  % Zero padding to nearest N power 2
df = fs/NFFT;
dt = 1/df;
%Frequency Axis defintion
f = (-(fs-df)/2:df:(fs-df)/2)';


X(2^ceil(log2(length(X))))=0;
Y(2^ceil(log2(length(Y))))=0;

%calculate time axis
T = (dt:dt:(length(X)*dt))';

subplot(2,2,1)
plot(T,X);
xlabel('Time(s)')
ylabel('X amplitude')
title('X signal before filtering')
subplot(2,2,2)
plot(T,Y);
xlabel('Time(s)')
ylabel('Y amplitude')
title('Y signal before filtering')

 X = detrend(X,0); % Removing DC Offset
 Y = detrend(Y,0); % Removing DC Offset


 % Filter parameters:
 M = length(X);    % signal length
 L = M;          % filter length   
 fc = 2*(38000/60);   % cutoff frequency


 % Design the filter using the window method:
 hsupp = (-(L-1)/2:(L-1)/2);
 hideal = (2*fc/fs)*sinc(2*fc*hsupp/fs);
 h = hamming(L)' .* hideal; % h is our filter



 % Zero pad the signal and impulse response:
 X(2^ceil(log2(M)))=0;
 xzp = X; 
 hzp = [ h zeros(1,NFFT-L) ];

 % Transform the signal and the filter:
 X = fft(xzp);
 H = fft(hzp)';
 X = X .* H;
 X = ifft(X);
 relrmserrX = norm(imag(X))/norm(X); % checked... this for zero
 X = real(X)';

 % Zero pad the signal and impulse response:
 Y(2^ceil(log2(M)))=0;
 xzp = Y; 
 hzp = [ h zeros(1,NFFT-L) ];

 % Transform the signal and the filter:
 Y = fft(xzp);
 H = fft(hzp)';
 Y = Y .* H;
 Y = ifft(Y);
 relrmserrY = norm(imag(Y))/norm(Y); % check... should be zero
 Y = real(Y)';

我画了过滤后的图,你可以在图片中看到有明显的偏差。我只想过滤噪声,但信号似乎会丢失其他成分,如果那是正确的过滤方法,我会感到有点困惑。 任何建议、提示或想法都会有所帮助。

最后我想绘制 X 与 Y 的关系图以给出轴振动的轨道。另请在下面找到另一张未过滤和过滤轨道的图片。正如您在图片中看到的那样,原始轨道也发生了变化(左图有很多噪音)。

P.S.: 我没有DSP工具箱

你的FFT和IFFT没有问题

您可以使用频率非常低的简单正弦波来测试您的代码。最终输出应该是(几乎)相同的正弦波,因为你正在对它进行低通滤波。

您已将 X(和 Y)定义为从 0 到某个值。但是,您已将 H 从 - (L-1)/2 定义为某个正值。从数学上讲,这很好,但是您只是简单地采用 H 的 fft。当您将 ffts 相乘时,Matlab 认为这与 X 具有相同的时间尺度!

因此,实际上,您已经采用了 XHfft(delta(t-d)) 的 fft,其中 d 是产生的时移。您可以在频域或时域中撤消此时移。