-
Notifications
You must be signed in to change notification settings - Fork 0
/
Example.m
71 lines (57 loc) · 2.59 KB
/
Example.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
Fs = 500; % Sampling frequency
T = 1/Fs; % Sample time
L = 1000; % Length of signal
t = (0:L-1)*T; % Time vector
fig = figure('units','normalized','outerposition',[0 0 1 1]);
subplot(3,1,1);
%% Input signal consisting of a 20 Hz, 50 Hz and 120 Hz sine wave
x = 2*sin(2*pi*20*t) + 3*sin(2*pi*50*t) + 5*sin(2*pi*120*t);
%% Plot of input signal
subplot(4,1,1);
plot(Fs*t,x,'LineWidth',2);
xlabel('Time (ms)','FontSize',22,'FontWeight','Bold');
ylabel('Magnitude','FontSize',22,'FontWeight','Bold');
title('Input Signal Consisting of a 20Hz, 50 Hz, and 120 Hz sine wave','FontSize',32,'FontWeight','Bold');
%% Input signal with some random noise and DC offset
y = x + 2*randn(size(t)) + 0.5;
subplot(4,1,2);
plot(Fs*t,y,'LineWidth',2);
xlabel('Time (ms)','FontSize',22,'FontWeight','Bold');
ylabel('Magnitude','FontSize',22,'FontWeight','Bold');
title('Input Signal with Added Noise','FontSize',32,'FontWeight','Bold');
%% Bandpass filter the input signal using symmetric FIR with min zeros
Fstop1 = 10; % First Stopband Frequency
Fpass1 = 35; % First Passband Frequency
Fpass2 = 150; % Second Passband Frequency
Fstop2 = 167; % Second Stopband Frequency
Astop1 = 60; % First Stopband Attenuation (dB)
Apass = 1; % Passband Ripple (dB)
Astop2 = 60; % Second Stopband Attenuation (dB)
filt_design = fdesign.bandpass('fst1,fp1,fp2,fst2,ast1,ap,ast2', Fstop1, Fpass1, Fpass2, Fstop2, Astop1, Apass, Astop2, Fs);
filt_mdl = design(filt_design, 'equiripple', 'FilterStructure', 'dfsymfir','MinOrder', 'any');
filt_x = filt_mdl.filter(y);
coefficients = filt_mdl.Numerator;
%% Window the filtered data using blackman window
win_mdl = window(@blackman,L);
win_x = win_mdl.*filt_x';
wvtool(win_mdl);
%% Obtain spectral energy of noisy input
NFFT = 2^nextpow2(L);
Y = fft(y,NFFT)/L;
f = Fs/2*linspace(0,1,NFFT/2+1);
subplot(4,1,3);
plot(f,2*abs(Y(1:NFFT/2+1)),'LineWidth',2);
xlim([-2 250]);
xlabel('Frequency (Hz)','FontSize',22,'FontWeight','Bold');
ylabel('|Y(f)|','FontSize',22,'FontWeight','Bold');
title('Spectral Analysis of Noisy Input Signal','FontSize',32,'FontWeight','Bold');
%% Obtain spectral energy of filtered input
NFFT = 2^nextpow2(L);
Y = fft(win_x,NFFT)/L;
f = Fs/2*linspace(0,1,NFFT/2+1);
subplot(4,1,4);
plot(f,2*abs(Y(1:NFFT/2+1)),'LineWidth',2);
xlim([-2 250]);
xlabel('Frequency (Hz)','FontSize',22,'FontWeight','Bold');
ylabel('|Y(f)|','FontSize',22,'FontWeight','Bold');
title('Spectral Energy Windowed and Bandpass Filtered 35-200Hz','FontSize',32,'FontWeight','Bold');