function [sp,df] = tsmooth(p,kw)
%
% s = tsmooth(p,kw) computes a spectral estimate by smoothing periodogram 
% p, using a TRIANGULAR window of full (odd) width = kw points 
%
% see pdgm.m for details of suitable periodogram calculation.
% see also rsmooth.m, csmooth.m for box & cosine-bell smoothing windows
%
% MIM Feb/90
%     May-92 a call [sp,df] = tsmooth(p,kw) now gives an extra o/p vector, df,
%     the first element of which is the no. deg. of freedom in the Chi-square
%     approx. to the dist. of the spectral estimates, sp.
%

if kw < 1, kw=1;msg='warning, zero length window call',end;
kworig=kw;
%
% generate triangular window of full (odd) width = original
% kw (if it wasn't odd go to next higher odd no.). window ->0
% at lags +/- (kw/2+1) 
kw=fix(kw/2);k0=kw+1;
w=ones(2*kw+1,1);
for k=1:kw,
 w(k0+k)=1-k/(kw+1);
 w(k0-k)=1-k/(kw+1);
end;
w=w/sum(w);
%%msg='check window',keyboard;
%
%
%
% reflect periodogram about end points .....
% original periodogram now starts with zero freq. pt at index = kw+1
% Note use of conjugate in reflected bits so we can handle cross-periodograms
% This won't affect auto-periodograms.
np=length(p);
pp=ones(np+(2*kw),1);
pp(kw+1:kw+np)=p(1:np);
pp(kw:-1:1)=conj(pp(kw+1+1:kw+1+kw));
pp(kw+np+1:kw+np+kw)=conj(pp(kw+np-1:-1:kw+np-kw));
pp(kw+1)=0;%should already be so if mean correction done at pdgm stage;
%%msg='check reflection in pp',keyboard;
%
%
% run window along pp...
sp=conv(pp,w);
%
% select desired bit.....
sp = sp(2*kw:2*kw+np-1);
%
% adjust first kw+1 points for zero in periodogram at f=0
for k=0:kw
  fixx=1/(1-w(k0+k));
  sp(k+1)=sp(k+1)*fixx;
end

df=(1+0*sp)*(2/sum(w.*w));
