function rr = myecg2hr(ecgt,ecg, threshfilt)

% determine sampling rate
fs = 10/(ecgt(11)-ecgt(1));

% generate thresholding filter
if(nargin < 3) 
    n = 5;              % order of 
    k = -10;            % sensitivity setting
    w = fs/9;           % how wide the threshold should be
    
    w = 2^ceil(log(w)/log(2));
    l = floor(w/n);
    threshfilt = real(ifft(fft(ones(2*l,1),2*w).^n));
    b = k/sum(threshfilt);
    threshfilt = b*threshfilt;
    threshfilt = max(threshfilt)-threshfilt;
end

% specify times of "good" beat
t1 = 42.23;
t2 = 42.26;

% specify pass band for ecg data
lowpasscutoff = 350;
highpasscutoff = 15;

% bandpass filter the ecg
ecgfilt = lowpass(lowpasscutoff,fs,highpass(highpasscutoff,fs,ecg));

% find time points for "good" beat
times = find((ecgt > t1) & (ecgt <= t2));
matched = fft(ecgfilt(times), length(ecgfilt));

% match filter
Fecgfilt2 = fft(ecgfilt).*matched;
ecgfilt2 = real(ifft(Fecgfilt2));

% theshold the signal
thresh = real(ifft(fft(abs(ecgfilt2),length(ecgfilt2)+length(threshfilt)).*fft(threshfilt,length(ecgfilt2)+length(threshfilt))));
n1 = floor(length(threshfilt)/2);
n2 = length(ecgfilt2);
thresh = thresh(1:(n2-n1+1));
thresht = ecgt(n1:n2);
cond = ecgfilt2(n1:n2) > thresh;
cond = diff(cond) == 1;
rrt = ecgt(cond);
rr = diff(rrt);

% display results of filtering
figure(1)
%subplot(411);
%plot(ecgt(1:length(threshfilt)),threshfilt);
subplot(311);
plot(ecgt,ecg);
subplot(312);
plot(ecgt,ecgfilt);
subplot(313);
plot(ecgt,ecgfilt2);

% display results of thesholding
figure(2);
subplot(311);
hold off;
n2
n1
size(ecgt)
size(threshfilt)
plot(ecgt,ecgfilt2,thresht,thresh,ecgt(1:length(threshfilt)),threshfilt/max(abs(threshfilt)));
hold on;
stem(rrt,mean(thresh)*ones(size(rrt)),'r');
subplot(312);
plot(rrt(1:length(rr)), rr);

% computer heart rate
[nt nv] = hrtach(0, rr*1000, ones(size(rr)), 64);
subplot(313);
plot(nt,nv);
