%Powerline noise elimination
%MATLAB tutorial series (Part 2.2)

%% Step 1: plot the signal

% Load ecg signal from .mat file
load('TNS_2_3_2019_example1.mat','ecg')
% Fs is given
Fs=250;
% Ts sampling period
Ts=1/Fs;
% Length of the signal
N=length(ecg);
% Ending time
Tmax=(N-1)*Ts;
% Time vector
t=0:Ts:Tmax;
% Plot the original signal
figure(1)
plot(t,ecg)
xlabel('Time (s)')
title('ECG + noise')


%% Step 2: plot one-sided magnitude spectrum

% Compute fft
ECG=fft(ecg);
% Take abs and scale it
ECG2=abs(ECG/N);
% Pick the first half
ECG1=ECG2(1:N/2+1);
% Multiply by 2 (except the DC part), to compenseate 
% the removed side from the spectrum.
ECG1(2:end-1) = 2*ECG1(2:end-1);
% Frequency range
F = Fs*(0:(N/2))/N;
% Plot single-sided spectrum
figure(2)
plot(F,ECG1,'LineWidth',2.5) 
title('Single-Sided Amplitude Spectrum')
xlabel('f (Hz)');


%% Step 3: viusal ...


%% Step 4: I already designed my filt ...

% load the filter object
load Notch50;
% removing the noise
pure_ecg=filter(Notch50,ecg);


%% Step 5: plot the results

figure(3)
plot(t,pure_ecg,'LineWidth',2); title('without noise')             
xlabel('time (s)')
ylabel('amplitude')
% Zoom in
figure(4)
plot(t(500:length(t)/6),ecg(500:length(t)/6),'LineWidth',2.5); title('Noisy and without noise')             
xlabel('time')
ylabel('amplitude')
hold on
plot(t(500:length(t)/6),pure_ecg(500:length(t)/6),'LineWidth',2.5)


%% Filter visualization

fvtool(Notch50)







