

REstar = [2400 1800 1200 600 0 -600 -1200 -1800];
%tbp0 = [1020.4 1020.55 1020.75 1020.95  1021.05 1021.0 1021.25 1021.3];
nR = length(REstar);
Bo=1000;

index = [20 80 140 200 280 360];
ni = length(index);

XMAX = nan*ones(nR,ni);
YMAX = nan*ones(nR,ni);
TIME = nan*ones(nR,ni);

for j=1:nR,
    for i=1:ni,
        [REstar(j)  index(i)]
  [F,x,y,t]  = read_binary_file(9,REstar(j),Bo,index(i),'f','dump');
%[P,x,y,t]  = read_binary_file(9,REstar(i),Bo,index0,'p');
%[Z,x,y,t]  = read_binary_file(9,REstar(i),Bo,index0,'z');
  ii = find(F>0.2);  F(ii)=1;   ii=find(F<=0.21);  F(ii)=0;

  R = read_surface_extract(9,REstar(j),Bo,t);
  Xd = R(:,3);  
  Yd = R(:,4);
  ii=find(Xd>x(1) & Xd<x(end) & Yd<0.9);
  Xd=Xd(ii); Yd=Yd(ii);

  [jj,ii] = find(F);
  Y=y(jj);
  [ymx,imx]=max(Y)
  ymax = y(jj(imx));
  xmax = x(ii(imx));

  ii = find( abs(Xd-xmax)<1 &  abs(Yd-ymax)<0.4 );
  Xd2 = Xd(ii);  Yd2 = Yd(ii);
  [ymax2,imx2] = max(Yd2);
  xmax2 = Xd2(imx2);
  XMAX(j,i) = xmax2;
  YMAX(j,i) = ymax2;  
  TIME(j,i)=t;
  figure(1)
  imagesc(x,y,F)
  axis xy
  hold on
  plot(Rpc.Location(:,1),Rpc.Location(:,2),'m.')
  plot(xmax2,ymax2,'r*','markersize',12)
  hold off
  title(sprintf('REstar = %f,  index=%3d,  time = %.2f',REstar(j),t,index(i)))
  pause(0.1)
  
end;
end;

figure(2),clf
C= 1.26491106;
tt = TIME(2,:);
xt = C*(tt-1001)+15;

XMAXP = XMAX - ones(8,1)*xt;
subplot(2,1,1)
plot(tt,XMAX,'*-',tt,xt)

subplot(2,1,2)
plot(tt,XMAXP,'*-')

figure(3),clf
subplot(2,1,1)
plot(tt,YMAX,'*-')
legend('2400','1800','600','0','-600','-1200','-1800','location','northwest')
subplot(2,1,2)
plot(XMAX,YMAX,'*-')

return

%hold on
%patch(Xf,Yf,'g')

figure(2),clf
imagesc(x,y,P)
colorbar('vert')
caxis([-0.01 0.01])
axis xy

figure(3)
imagesc(x,y,BW)
axis xy

