function disper01(X, Y)
%
%   the function disper01(X, Y) calculates and plots the dispersion
%   of a group of floats whose x,y coordinates are the rows of 
%   matrices X and Y respectively (each column being a float).
%   Reference:  
%     Davis, R. E., P. D. Killworth and J. R. Blundell,1996: Comparison
%   of Autonomous Lagrangian Circulation Explorer and fine resolution 
%   Antarctic model reults in the South Atlantic.  J. Geophys. Res.,
%   101, C1, 855-884.
%
%------------------------------------------------------------------------
%   Written by Vassilis Zervakis, in June 1996, while at
%              AINCO-Interocean, Madrid, Spain
%------------------------------------------------------------------------

%   Calculation of coordinates of center of mass
[ktotal,NFL] = size(X);
xc = zeros(ktotal,1);
yc = zeros(ktotal,1);
Nc = zeros(ktotal,1);
for jj = 1:ktotal
  for ifl = 1:NFL
     if X(jj,ifl) ~= 999.,
        Nc(jj) = Nc(jj) + 1;
        xc(jj) = xc(jj) + X(jj,ifl); 
        yc(jj) = yc(jj) + Y(jj,ifl); 
     end
  end
  if Nc(jj) == 0.,
     xc(jj) = 999.;
     yc(jj) = 999.;
    else
     xc(jj) = xc(jj) / Nc(jj);
     yc(jj) = yc(jj) / Nc(jj);
  end
end
%
%   Calculation of dispersion of all floats
disper = zeros(ktotal,1);
for jj = 1:ktotal
  if Nc(jj) <= 1,
    disper(jj) = 0.;
   else,
    for ifl = 1:NFL
     if X(jj,ifl) ~= 999.,
      disper(jj)=disper(jj)+(latlondi(X(jj,ifl),Y(jj,ifl),xc(jj),yc(jj)))^2;
     end
    end
    disper(jj) =  sqrt(disper(jj) / (Nc(jj)-1));
   end
end
figure
plot(disper);
ylabel('dispersion, km');
xlabel('days');
mbox('OK','Use the mouse to define range of','computation of diffusivity')
[x1,y1] = ginput(2);
aa = zeros(ktotal,2);
bb = zeros(ktotal,1);
ncount = 0;
for jj = 1:ktotal
  if (x1(1) <= jj) & (x1(2) >= jj),
     ncount = ncount + 1;
     aa(ncount,1) = jj*24.*3600.;
     aa(ncount,2) = 1;
     bb(ncount,1) = disper(jj)*100000.;
  end
end
rho = mean(bb(1:ncount));
fit = (aa \ bb) /2 ;
Kfit = rho*(fit(1)) /2;
mbox('OK',['Eddy diffusivity is' num2str(Kfit)],' cm^2/s') ;
