The code below generates a $3-$beat ECG signal and then adds powerline interference at $60$ Hz along with its harmonics at $180$ Hz and $300$ Hz$.$ Run the code in MATLAB and explain the purpose of the chosen filter.

Then, add an additional harmonic at $400$ Hz to the noise signal and modify the filter by adding a zero at this frequency to remove the $400$ Hz harmonic from the ECG signal. Include all $5$ resulting figures with a brief explanation for each.

$$

MATLAB CODE:

clear all, close all

$\%$ Step $1$: Define parameters

fs $=1000$; $\%$ Sampling frequency $(500$ Hz$)$

T $=3$; $\%$ Duration in seconds for $3$ beats

t $=0$:$1/$fs:T$-(1/$fs$)$; $\%$ Time vector

$\%$ Step $2$: Create sine signal $($pretend ECG$)$

single$\_$beat $=$ sin$(2*$ pi $*1*$ t$(1$:end$/3))$; $\%$ Replace with a sine wave for testing

single$\_$beat $=$ single$\_$beat $/$ max$($abs$($single$\_$beat$))$; $\%$ Normalize signal

num$\_$beats $=3$;

ecg$\_$signal $=$ repmat$($single$\_$beat, $1,$ num$\_$beats$)$;

t$\_$ecg $=0$:$1/$fs:$($length$($ecg$\_$signal$)-1)/$fs;

$\%$ Step $3$: Add power line interference

f$\_$noise$1=60$; f$\_$noise$2=180$; f$\_$noise$3=300$;

noise $=0\mathrm{.}2*$ sin$(2*$ pi $*$ f$\_$noise$1*$ t$\_$ecg$)+...$

$0\mathrm{.}1*$ sin$(2*$ pi $*$ f$\_$noise$2*$ t$\_$ecg$)+...$

$0\mathrm{.}05*$ sin$(2*$ pi $*$ f$\_$noise$3*$ t$\_$ecg$)$;

noisy$\_$ecg $=$ ecg$\_$signal $+$ noise;

$\%$ Step $4$: Plot noisy ECG signal

figure; plot$($t$\_$ecg, noisy$\_$ecg, $\text{'}$r$\text{'})$;

title$(\text{'}$Noisy Sine Signal'$)$;

xlabel$(\text{'}$Time $($s$)\text{'})$; ylabel$(\text{'}$Amplitude$\text{'})$; grid on;

$\%$ Step $5$: Power Spectral Density

$[$psd$\_$values, f$]=$ pwelch$($noisy$\_$ecg, $[],[],[],$ fs$)$;

figure; plot$($f$,10*$log$10($psd$\_$values$))$;

title$(\text{'}$Power Spectral Density'$)$; xlabel$(\text{'}$Frequency $($Hz$)\text{'})$;

ylabel$(\text{'}$Power$/$Frequency $($dB$/$Hz$)\text{'})$;

grid on;

$\%$ Step $6$: Design comb filter

f$0=[$f$\_$noise$1,$ f$\_$noise$2,$ f$\_$noise$3]$;

z $=[]$;

for i $=1$:length$($f$0)$

theta $=2*$ pi $*$ f$0($i$)/$ fs;

z $=[$z$,$ exp$(1$j$*$theta$),$ exp$(-1$j$*$theta$)]$;

end

b $=$ poly$($z$)$; a $=1$; $\%$ Filter coefficients

$\%$ Step $7$: Apply comb filter

filtered$\_$ecg $=$ filter$($b$,$ a$,$ noisy$\_$ecg$)$;

$\%$ Step $8$: Plot filtered ECG

figure; plot$($t$\_$ecg, filtered$\_$ecg, $\text{'}$b$\text{'})$;

title$(\text{'}$Filtered Sine Signal'$)$;

xlabel$(\text{'}$Time $($s$)\text{'})$; ylabel$(\text{'}$Amplitude$\text{'})$;

grid on;

$\%$ Step $9$: Frequency Response

figure; freqz$($b$,$ a$,1024,$ fs$)$;

title$(\text{'}$Frequency Response of Comb Filter'$)$;

grid on;