How can i calculate instantaneous heart rate from an ECG signal ?

Hello ^^ I'm trying to find instantaneous heart rate from this ECG sample.I did some filtration and used a code I've found on heart rate detection , but the heart rate plot i got is not reasonable.

Can anyone help me to know what is wrong ? or maybe have another way to obtain the heart rate ? this is the code i used to filter and obtain the HR.

D = load('test2s2.mat');
t = double(D.b1(:,1));
s = double(D.b1(:,2));
figure(1)
plot(t, s)
grid
Ts = mean(diff(t));                             % Sampling Interval (sec)
Fs = 1/Ts;                                      % Sampling Frequency (Hz)
Fn = Fs/2;                                      % Nyquist Frequency (Hz)
L = length(t);
fts = fft(s)/L;                                 % Fourier Transform
Fv = linspace(0, 1, fix(L/2)+1)*Fn;             % Frequency Vector
Iv = 1:length(Fv);                              % Index Vector
figure(2)
plot(Fv, abs(fts(Iv))*2)                        % Plot FFT
grid
xlabel('Frequency')
ylabel('Amplitude')
Wp =  [1.5  20]/Fn;                             % Passband (Normalised)
Ws =  [0.5  45]/Fn;                             % Stopband (Normalised)
Rp =  1;                                        % Passband Ripple (dB)
Rs = 50;                                        % Stopband Ripple (dB)
[n,Wn] = buttord(Wp, Ws, Rp, Rs);               % Filter Order
[b,a]  = butter(n,Wn);                          % Transfer Function
[sos,g] = tf2sos(b,a);                          % Convert To Second-Order-Section For Stability
figure(3)
freqz(sos, 4096, Fs)                            % Filter Bode Plot
sf = filtfilt(sos, g, s);                       % Filter Signal
xlim5 = [ 85   100];                            % Window (Time)
figure(5)                                       % Compare Unfiltered & Filtered Signal
subplot(2,1,1)
plot(t, s)
axis([xlim5    ylim])
title('Unfiltered EKG')
grid
subplot(2,1,2)
plot(t, sf)
axis([xlim5    ylim])
grid
title('Bandpass-Filtered EKG')

% square it detsq = sf .^ 2; figure,plot(detsq); % thresholded output detthres = zeros(length(detsq),1); % let's detect the momentary heart rate in beats per min last=0; upflag=0; pulse=zeros(length(detsq),1); p=0; for i = 1:length(detsq) if (detsq(i) > 0.25) if (upflag == 0) if (last > 0) T = i - last; p = 200/T*60; end last = i; end upflag = 100; else if (upflag>0) upflag = upflag - 1; end end pulse(i)=p; end % plot it figure,plot (t,pulse);