PLFILT Low-pass filters hourly data using the pl33 filter
USAGE: [ulp]=plfilt(u)
or
[ulp,jdlp]=plfilt(u,jd)
or
[ulp,jdlp]=plfilt(u,jd,nsub)
Input: u = matrix with time series in each column
jd = julian day time vector (optional)
nsub = subsampling interval for the output low-pass data (optional)
Output: ulp = low-pass filtered data with 33 hours removed at each end.
jdlp = julian time vector corresponding to the data ulp.
Example: [ulp,jdlp]=plfilt(u,jd,4); %jdlp and ulp at 4 hour intervals
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Version 1.0 (12/4/96) Rich Signell (rsignell@usgs.gov)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%0001 function [ufilt,jdf]=plfilt(u,jd,nsamp) 0002 % PLFILT Low-pass filters hourly data using the pl33 filter 0003 % 0004 % USAGE: [ulp]=plfilt(u) 0005 % or 0006 % [ulp,jdlp]=plfilt(u,jd) 0007 % or 0008 % [ulp,jdlp]=plfilt(u,jd,nsub) 0009 % 0010 % Input: u = matrix with time series in each column 0011 % jd = julian day time vector (optional) 0012 % nsub = subsampling interval for the output low-pass data (optional) 0013 % 0014 % Output: ulp = low-pass filtered data with 33 hours removed at each end. 0015 % jdlp = julian time vector corresponding to the data ulp. 0016 % 0017 % Example: [ulp,jdlp]=plfilt(u,jd,4); %jdlp and ulp at 4 hour intervals 0018 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0019 % Version 1.0 (12/4/96) Rich Signell (rsignell@usgs.gov) 0020 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0021 0022 %half the filter weights 0023 pl33=[ 0024 -0.00027 -0.00114 -0.00211 -0.00317 -0.00427 ... 0025 -0.00537 -0.00641 -0.00735 -0.00811 -0.00864 ... 0026 -0.00887 -0.00872 -0.00816 -0.00714 -0.00560 ... 0027 -0.00355 -0.00097 0.00213 0.00574 0.00980 ... 0028 0.01425 0.01902 0.02400 0.02911 0.03423 ... 0029 0.03923 0.04399 0.04842 0.05237 0.05576 ... 0030 0.05850 0.06051 0.06174 0.06215 ]'; 0031 0032 % symmetric filter 0033 pl33=[pl33(1:34); pl33(33:-1:1)]; 0034 0035 [m,n]=size(u); 0036 if(min(m,n)==1), 0037 m=length(u); 0038 n=1; 0039 u=u(:); 0040 end 0041 if (m<length(pl33)), 0042 disp('data too short for this filter!!'); 0043 return 0044 end 0045 ufilt=zeros(m-66,n); 0046 for i=1:n, 0047 uf=conv(u(:,i),pl33); 0048 ufilt(:,i)=uf(67:length(uf)-66); 0049 end 0050 [m,n]=size(ufilt); 0051 if (nargin>1), 0052 jdf=jd(34:34+m-1); 0053 end 0054 if (nargin==3) 0055 ufilt=ufilt([1:nsamp:m],:); 0056 jdf=jdf([1:nsamp:m]); 0057 end