function x = bpf(T1,T2,sp,x)
%BPF BandPass filter between T1 < period < T2
% BPF(T1,T2,SP,X) filters the input vector X, at the sampling
% rate SP, using a 2nd order bandpass digital Butterworth
% filter, constructed to retain information between
% the input periods T1<T2.  All arguments are REQUIRED.
% SP must be the sample rate in hours.  The MATLAB routine 
% FILTFILT (zero-phase digital filtering) is used to perform 
% the actual filtering.
%
%  INPUT: T1,T2 -  Period Band to pass (hours, T1<T2)
%         sp    -  sampling rate (hours)
%         x     -  signal
%
%  OUTPUT: y   -  band passed version of input x
%
% CALL AS: y=bpf(10,40,sp,x); % for band pass between 10 and 40 hours

% Brian Blanton, Fall 1998
%

if nargin==0 & nargout==0
   disp('Call as: y=bpf(T1,T2,sample_rate,x)')  
   return
end

w1=cutoff(sp,T1);
w2=cutoff(sp,T2);
[b,a] = butter(2,[w2 w1]);
x = filtfilt(b, a, x);