function coef=tidefit(uu,tt,tstart,tstop,period)
%function coefs=tidefit(uu,tt,tstart,tstop,periods)
%        U and T are the matrices containing the variable to be fit in U, and
%        the associated times in T.  Each row of U should be a time series.
%        TSTART and TSTOP indicate the time range to be plotted.
%        RANGE is an optional argument to specify the number of bins to plot.
%        It should be a vector, for example, RANGE=[1:10] will fit bins
%        through 10.  If no range is specified, all the bins are fit.
%
%        period=[ 12.42 6.21 4.14]      default
if nargin==4
       period=[12.42 6.21 4.14]   
end
tinc=(tstop-tstart)/50;
tp=[tstart:tinc:tstop];
tp=tp(:);
ntp=length(tp);
freq=(2.*pi)*ones(period)./period;
nfreq=length(freq);
   u=uu(:);
   t=tt(:);
   ind=find(u~=NaN);
   nind=length(ind);
   if nind>=9
      u=u(ind);
      t=t(ind);
      u=u(:);
      t=t(:);
keyboard
%------ set up A -------
       A=zeros(nind,nfreq*2+1);
       A(:,1)=ones(nind,1);
       for j=[1:nfreq]
           A(:,2*j)=cos(freq(j)*t);
           A(:,2*j+1)=sin(freq(j)*t);
       end
%-------solve [A coeff = u] -----------------
       coef=A\u;
%-------generate solution components-----
       up=zeros(ntp,nfreq*2+1);
       up(:,1)=coef(1)*ones(ntp,1);
       for j=[1:nfreq]
          up(:,j+1)=coef(j*2)*cos(freq(j)*tp)+coef(j*2+1)*sin(freq(j)*tp);
       end
       up(:,nfreq+2)=sum(up(:,[1:nfreq+1])')';  % sum of all comps
       
       error=std(A*coef-u);

   plot(t,u,'o',tp,up,tp,zeros(tp),'w')
   title(sprintf('Least Squares Fit:  std dev error= %5.2g cm/s',error));
   xlabel('Time, hours');
   ylabel('cm/second')
end
pause
