@ USING A KALMAN FILTER TECHNIQUE TO DERIVE THE H-P FILTER
  THE INPUT IS "DAT"==MATRIX OF SERIES THAT ARE DETRENDED
                     - REMEMBER TO TAKE LOGS
  THE OUTPUT IS "UU"==MATRIX H-P DEVIATIONS FROM TREND   @

local u1, u, l, ff, t, v11, v22, v12, x, z, v, tt, d, i, m1, m2, de,
      i1, ib, e1, e2, b11, b12, b22, jj, hu;
 

l=1600;  @ SMOOTHING PARAMETER, IT CAN BE CHANGED @

u1 = dat; ff=rows(dat'); t = rows(dat);
v11 = 1; v22 = 1; v12 = 0;
i=3; v=zeros(t,3); tt =zeros(t,ff); d = zeros(t,ff);

do until i>t;
 x=v11;  z=v12;
 v11=(1/L)+4*(x-z)+v22; v12=2*x-z; v22=x; de=v11*v22-v12*v12;
 v[i,1] = v22/de; v[i,2] = -v12/de; v[i,3] = v11/de;
 x = v11+1; z=v11;
 v11 = v11-(v11*v11)/x; v22 = v22-(v12*v12)/x; v12 = v12-(z*v12)/x;
 i=i+1;
endo;

u=u1; m1=u[2,.]; m2 = u[1,.];

i=3;
do until i>t;
 x = m1; m1=2*m1-m2; m2=x;
 tt[i-1,.] = v[i,1]*m1+v[i,2]*m2;
 d[i-1,.] = v[i,2]*m1+v[i,3]*m2;
 de = v[i,1]*v[i,3]-v[i,2]*v[i,2];
 v11 = v[i,3]/de; v12 = -v[i,2]/de;
 z = (u[i,.]-m1)/(v11+1);
 m1=m1+v11*z; m2=m2+v12*z;
 i=i+1;
endo;

tt[t,.] = m1; tt[t-1,.] = m2;
m1 = u[t-1,.]; m2 = u[t,.];

i=t-2;
do until i<1;
 i1=i+1; ib=t-i+1;
 x=m1;
 m1=2*m1-m2; m2=x;
  if i>2;
   e1 = v[ib,3]*m2+v[ib,2]*m1+tt[i,.]; e2 = v[ib,2]*m2+v[ib,1]*m1+d[i,.];
   b11 = v[ib,3]+v[i1,1]; b12 = v[ib,2]+v[i1,2]; b22 = v[ib,1]+v[i1,3];
   de = b11*b22-b12*b12;
   tt[i,.] = (-b12*e1+b11*e2)/de;
  endif;
 de = v[ib,1]*v[ib,3]-v[ib,2]*v[ib,2];
 v11 = v[ib,3]/de; v12 = -v[ib,2]/de;
 z = (u[i,.]-m1)/(v11+1);
 m1 = m1 +v11*z; m2 = m2 +v12*z;
 i=i-1;
endo;

tt[1,.] = m1;
tt[2,.] = m2;
i=1;
do until i>t;
 d[i,.] = u[i,.]-tt[i,.];
 i=i+1;
endo;
hu=tt;
UU = d-ones(t,ff)*MEANC(D);
CLEAR U1;
retp(uu);
endp;