M = importdata('C:\Users\mingfaifong\Documents\raw traces\anchovy6_E-69_VC-1.txt','\t',6);

%% filter out line noise

fs = 1/(mode(diff(M.data(:,1)))*0.001);             %#sampling rate
f0 = 60;                %#notch frequency
fn = fs/2;              %#Nyquist frequency
freqRatio = f0/fn;      %#ratio of notch freq. to Nyquist freq.

notchWidth = 0.1;       %#width of the notch

%Compute zeros
zeros = [exp( sqrt(-1)*pi*freqRatio ), exp( -sqrt(-1)*pi*freqRatio )];

%#Compute poles
poles = (1-notchWidth) * zeros;

figure;
zplane(zeros.', poles.');

b = poly( zeros ); %# Get moving average filter coefficients
a = poly( poles ); %# Get autoregressive filter coefficients

figure;
freqz(b,a,32000,fs)

%#filter signal x
x = M.data(:,2);
y = filter(b,a,x);

%% line noise removal (try 2)

% Design a filter with a Q-factor of Q=35 to remove a 60 Hz tone from 
% system running at 300 Hz.
Wo = 60/(10000/2);  BW = Wo/35;
BW = 0.1;
[b,a] = iirnotch(Wo,BW); 
x = M.data(:,2);
y = filter(b,a,x);

%fvtool(b,a,x);

%%
hold off;
plot(M.data(:,1)-min(M.data(:,1)),M.data(:,2),'-k');
hold on;
plot(M.data(:,1)-min(M.data(:,1)),y-50,'-k');

ylim([-150 0]);
xlim([0 2000]);
