% function kino(flg,jlong,n)
    flg=1;
    n=3;
    jlong=2150;
% program kino_tp.m     
% animation of the QG model output: predicted streamfunction (SLA) field 
% and referenced SLA field (OA output).
% $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
%                                                                      $
%           PROGRAM KINO_TP                                            $
%                                                                      $
% $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
% set up input parameters:
% origin of the coordinate system
       f=pi*64.2/180;
%      yn=9.0;
%      xn=22.95;
       yn=10.2728;
       xn=23.4591;
% BATS site coordinates
       xbats=-64.17;
       ybats=31.74;
% set julian day, model day, the beginning and the ending days of the movie  
       jul=datenum(1995,10,31); % (1992,10,27)
       t=jul-datenum(1992,1,1)
       jbeg=1
       jend=jbeg+jlong
% set colormap
       cmapname='jet';
       ncolbase=64;
       nband=10;
       
       [cmap]=gen_cmap(ncolbase,nband,cmapname);
      
       colormap(cmap);
% set caxis limits
       cval1=[-20:4:20];
       cval2=[-10:4:10];
       cval3=[1:11];
% create array geographical coordinates rotated cclockwise 64 deg     
       for k=1:1:25
       for l=1:1:31
       xi(k,l)=(xn)+(-(k-1)*sin(f)+(l-1)*cos(f))*0.366;
       yi(k,l)=(yn)+( (k-1)*cos(f)+(l-1)*sin(f))*0.366;
       end
       end
       xi=-70+(xi-.85)*10/34;
       yi=28.25+(yi-1)*8/32;      
% "diamond" perimeter
       xv=[xi(1,1);xi(1,31);xi(25,31);xi(25,1);xi(1,1)];
       yv=[yi(1,1);yi(1,31);yi(25,31);yi(25,1);yi(1,1)];  
% set Box Aspect Ratio 
       lonmax=-min(xi(:));
       lonmin=-max(xi(:));
       latmax=max(yi(:));
       latmin=min(yi(:));
       weight=111*cos(pi*(latmin+(latmax-latmin)*.5)/180);
       xc=(lonmax-lonmin)*weight;
       yc=(latmax-latmin)*111;
       relyx=yc/xc;
% load source file to draw  T/P and ERS tracks
       if(flg==1)
       load /usr/people/valery/altim/aviso/tp36ers2a13_1;            
       va1=tp36ers2a13_1;
       clear tp36ers2a13_1;
       load /usr/people/valery/altim/aviso/tp12ers1c;            
       va1=[tp12ers1c;va1];
       clear tp12ers1c;
       load /usr/people/valery/altim/aviso/tp36ers2a13_2;            
       va1=[va1;tp36ers2a13_2];
       clear tp36ers2a13_2;
       end
       load /usr/people/valery/altim/aviso/berm.dat;
% load dynheig
       load /usr/people/valery/bval16/dynh20.dat;
       dyn=dynh20(:,3);
% load source files: computed and referenced streamfunction and vorticity
       load /usr/people/valery/p2disk/qgsbl/exam/dh3a15b5f.vbc;
       a=dh3a15b5f;
       a=a/2.3162;
       clear dh3a15b5f;
       load /usr/people/valery/p2disk/qgsbl/exam/dh3a15b5f.vbr;
       b=dh3a15b5f;
       b=b/2.3162;
       clear dh3a15b5f;
% calculate rms difference (cm)   
       rsl=sqrt(mean((a(1:2151,:)-b(1:2151,:)).^2,2));
% dimensional mean
      rms=mean(rsl);
% referenced SLA (streamfunction) rms
       sf=sqrt(mean((b(1:2151,:)).^2,2));
% normalized RMS difference
       v=rsl./mean(sf)*100;
% normaliaed mean
       nmrs=mean(v)
 set(gcf,'color','k');
       hh=1;
% start drawing SLA field
       for j=[jbeg:jend];
       as=reshape(a(j,:),31,25);
       as=as';
       bs=reshape(b(j,:),31,25);
       bs=bs';
       cs=as-bs;
% draw predicted SLA
       
       subplot(2,2,1)
       contourf3(xi,yi,as,cval1);
       hold on;
       plot(berm(:,1),berm(:,2),'color','k','LineWidth',[3.]);
%      plot(xl,yl,'k-','LineWidth',[4.0]);
       if(isempty(find((j+299)==dyn)))
       plot(xbats,ybats,'r.','MarkerSize',22);
       else
       plot(xbats,ybats,'o','MarkerSize',10,'Markerfacecolor',...
       'w','Markeredgecolor','r');
       end
       text(xbats,ybats,['  BATS '],'HorizontalAlignment','left',...
       'FontWeight','bold','FontSize',10);
       hold off
%       axis square;
       set(gca,'PlotBoxAspectRatio',[1 relyx*1 1]);
       caxis([-20 20]);
              
       h1=colorbar;
   

       set(h1,'XColor','w','YColor','w');
       set(h1,'FontWeight','bold','FontSize',10);
       set(h1,'ytick',cval1);
%       set(h1,'YTickLabel',[-20 -12 -4 0 4 12 20]); 
       set(h1,'ticklength',[0.05 0.025]);
       set(gca,'XColor','w','YColor','w');
       set(gca,'FontWeight','bold','FontSize',12);
%      xlabel('Longitude (deg E)','color','w','FontWeight','bold','FontSize',8);
       ylabel('Latitude (deg N)','color','w','FontWeight','bold','FontSize',12);
       title(['Predicted SLA (cm)'],'color','g',...
       'FontAngle','italic','FontWeight','bold','FontSize',12);

       set(gca,'XGrid','on','color','k','YGrid','on','color','k');
% find TP and ERS tracks for 7 days period      
       if(flg==1)
       tmin=t-n;
       tmax=t+n+1;
       v1=va1(find(va1(:,1)>=tmin & va1(:,1)<=tmax),:);
       er1=v1(find(v1(:,5)==15),:);
       tp1=v1(find(v1(:,5)==1),:);
       end; 
% draw observed SLA
       subplot(2,2,2)
       contourf3(xi,yi,bs,cval1);
       hold on
       plot(berm(:,1),berm(:,2),'color','k','LineWidth',[3.]);
%      plot(xl,yl,'k-','LineWidth',[4.0]);
       if(flg==1)
       plot(tp1(:,3),tp1(:,2),'w.','markersize',7.0);      
       plot(er1(:,3),er1(:,2),'g.','markersize',7.0);      
       end
       if(isempty(find((j+299)==dyn)))
       plot(xbats,ybats,'r.','MarkerSize',22);
       else
       plot(xbats,ybats,'o','MarkerSize',10,'Markerfacecolor',...
       'w','Markeredgecolor','r');
       end
       text(xbats,ybats,['  BATS '],'HorizontalAlignment','left',...
       'FontWeight','bold','FontSize',10);
       hold off
%       axis square;
       set(gca,'PlotBoxAspectRatio',[1 relyx*1 1]);
       caxis([-20 20]);
       
       h2=colorbar;
%
       set(h2,'XColor','w','YColor','w');
       set(h2,'FontWeight','bold','FontSize',10);
       set(h2,'ytick',cval1);
%       set(h2,'YTickLabel',[-20 -12 -4 0 4 12 20]);
       set(h2,'ticklength',[0.05 0.025]);
       set(gca,'XColor','w','YColor','w');
       set(gca,'FontWeight','bold','FontSize',12);
%      xlabel('Longitude (deg E)','color','w','FontWeight','bold','FontSize',8);
%      ylabel('Latitude (deg N)','color','w','FontWeight','bold','FontSize',12);
       title(['Observed SLA (cm)'],'color','g',...
       'FontAngle','italic','FontWeight','bold','FontSize',12);
       set(gca,'XGrid','on','color','k','YGrid','on','color','k');
% draw predicted - observed SLA difference
       subplot(2,2,3)
       contourf3(xi,yi,cs,cval1);
       hold on;
       plot(berm(:,1),berm(:,2),'color','k','LineWidth',[3.]);
%      plot(xl,yl,'k-','LineWidth',[4.0]);
       if(isempty(find((j+299)==dyn)))
       plot(xbats,ybats,'r.','MarkerSize',22);
       else
       plot(xbats,ybats,'o','MarkerSize',10,'Markerfacecolor',...
       'w','Markeredgecolor','r');
       end
       text(xbats,ybats,['  BATS '],'HorizontalAlignment','left',...
       'FontWeight','bold','FontSize',10);
       hold off
%       axis square;
       set(gca,'PlotBoxAspectRatio',[1 relyx*1 1]);
       caxis([-20 20]);

       h3=colorbar;

       set(h3,'XColor','w','YColor','w');
       set(h3,'FontWeight','bold','FontSize',10);
       if(j==1)
       set(h3,'ytick',cval3);
       set(h3,'yticklabel',[-20 -16 -12 -8 -4 0 4 8 12  16  20]);
       else
       set(h3,'ytick',cval1);
%       set(h3,'yticklabel',[-10  -6  -2  2  6  10]);
       end
       set(h3,'ticklength',[0.05 0.025]);
       set(gca,'XColor','w','YColor','w');
       set(gca,'FontWeight','bold','FontSize',12);
       xlabel('Longitude (deg E)','color','w','FontWeight','bold','FontSize',12);
       ylabel('Latitude (deg N)','color','w','FontWeight','bold','FontSize',12);
       title(['SLA Difference (cm)'],'color','g',...
       'FontAngle','italic','FontWeight','bold','FontSize',12);
       set(gca,'XGrid','on','color','k','YGrid','on','color','k');
% draw plot
       subplot(2,2,4)
       xx=300:2450;      
       h5=area(xx,rsl);
% draw gap in ERS data
       xz=722:1:1214;
       yz(1:493)=8.5;
        hold on
       plot(xz,yz,'-','Linewidth',2.5)
       y1=0:.01:rsl(j);
       x1=j+299;      
       plot(x1,y1,'w.','MarkerSize',5);
       plot(x1,rsl(j),'*w','Markersize',6)       
       hold off
       tt=text(680,9,['NO ERS DATA'],'Fontsize',7);
%      text(370,90,['NRMS=',num2str(nmrs,'%5.1f'),'%'],'HorizontalAlignment','left','fontsize',6);
       axis([300 2450 0 10]);
       grid on
       set(gca,'XTick',[366 731 1096 1461 1827 2192])
       set(gca,'XTickLabel',[93 94 95 96 97 98])
%
       set(gca,'XColor','w','YColor','w');
       set(gca,'FontWeight','bold','FontSize',10);
       set(gca,'XGrid','on','YGrid','on');
%       set(gca,'YTick',[0 2 4 6 8 10]);
%       set(gca,'YTickLabel',[0 2 4 6 8 10]);

       xlabel('Time (Years)','color','w','FontWeight','bold','FontSize',12);
       ylabel('RMS Diff  (cm)','color','w','FontWeight','bold','FontSize',12);
       title([datestr(jul,1)]...
       ,'color','g','FontWeight','bold','FontSize',14);
       colormap(cmap);
       set(gcf,'InvertHardcopy','off');
% set path for output files
       jj=sprintf('/whoinet/deltorchio/p3/valery/qg_tp/st%i',j);
% or use a default directory
%       jj=sprintf('st%i',j);
       pp=['print ' jj '.ppm -dppm'];
       eval(pp);
       t=t+1;
       jul=jul+1;
%      input('?');
       end
% **************************************************************

