%
% This script plots an individual track on a satellite map background, with
% colors indicating intensity and dates marked every 24 hours. 
% Copyright 2009 WindRiskTech L.L.C.
%--------------------------------------------------------------------------
params
if ~exist('templot.mat','file') || ~exist('nt','var')
    setplot;
end
load templot
load temp latstore longstore
clf('reset')
%
%  Read in image
%
imagename=strcat('map_',imageres);
try
    load(imagename);
    %A=imread('eocene_map.jpg');
catch err
    errordlg('Cannot find image file; may neeed to run initialize.m')
    return
end    
%
[ay, ax, dum]=size(A);
aspect=ax/ay;
%       
image(A)
%
hold on
x(1:iq)=ax.*(longstore(nt,1:iq)+360-displacement)/360;
y(1:iq)=ay-ay.*(latstore(nt,1:iq)+90)./180;
%
if iq > 1
%
ntrack=1;
mini=iq;
np2=0;
x2=zeros(1,iq);
y2=x2;
v2=x2;pres2=x2;vpot2=x2;shear2=x2;monthm2=x2;daym2=x2;hour2=x2;
vtrans2=x2;
for i=1:iq
    if x(i) > ax
        ntrack=2;
        mini=min(i,mini);
        np2=np2+1;
        x2(np2)=x(i)-ax;
        y2(np2)=y(i);
        v2(np2)=v(i);
        pres2(np2)=pres(i);
        vpot2(np2)=vpot(i);
        shear2(np2)=shear(i);
        monthm2(np2)=monthm(i);
        daym2(np2)=daym(i);
        hour2(np2)=hour(i);
        vtrans2(np2)=vtrans(i);
    end
end
x2(np2+1:end)=[];
y2(np2+1:end)=[];
% 
if cspline == 'y'
   a=size(x);
   b=a(1,2);
   bx=1:b;
   tint=floor(6/trackres);
   tintin=1/tint;
   bxm=1:tintin:b;
   x=spline(bx,x,bxm);
   y=spline(bx,y,bxm);
   v=spline(bx,v,bxm);
   pres=spline(bx,pres,bxm);
   vpot=spline(bx,vpot,bxm);
   shear=spline(bx,shear,bxm);
   monthm=spline(bx,monthm,bxm);
   daym=spline(bx,daym,bxm);
   hour=spline(bx,hour,bxm);
   vtrans=spline(bx,vtrans,bxm);
end   
a=size(x);
b=a(1,2);
if colorcoded > 0
    xt=zeros(1,2);
    yt=zeros(1,2);
   for i=2:b
      xt(1)=x(i-1);
      xt(2)=x(i);
      yt(1)=y(i-1);
      yt(2)=y(i);
      vta=(1852/3600)*v(i-1);
      if colorcoded == 1
         if vta < 33
            vtcolor=vcolor(1,:);
         elseif vta >= 33 && vta < 45
            vtcolor=vcolor(round(0.2*qcolor),:);
         elseif vta >= 45 && vta < 54.4
            vtcolor=vcolor(round(0.47*qcolor),:);
         elseif vta >= 54.4 && vta < 62.4
            vtcolor=vcolor(round(0.6*qcolor),:);
         elseif vta >= 62.4 && vta < 69.5
            vtcolor=vcolor(round(.77*qcolor),:);
         else
            vtcolor=vcolor(round(0.95*qcolor),:);
         end
      else   
         vtr=round(0.4267.*(qcolor/64)*(3600.*vta./1852-30));
         %vtr=round(vta)-15;
         %vtr=round(vtr*(qcolor/64));
         vtr=min(vtr,qcolor);
         vtr=max(vtr,1);
         vtcolor=vcolor(vtr,:);
      end   
      plot(xt,yt,'color',vtcolor,'linewidth',thick)
      lc=vtcolor;
   end   
else
   plot(x,y,lc,'linewidth',thick)
end
bd=size(xd);
bd=bd(1,2);
if strcmp(dplot,'y')  
    offset=zeros(1,bd);
if bd > 1   
    ul=zeros(1,bd); 
    vl=zeros(1,bd);
   for kl=2:bd-1
      ul(kl)=0.5.*(xd(kl+1)-xd(kl-1));
      vl(kl)=0.5.*(yd(kl+1)-yd(kl-1));
   end
   ul(1)=4.*(x(2)-x(1));
   vl(1)=4.*(y(2)-y(1));
   ul(bd)=4.*(x(b)-x(b-1));
   vl(bd)=4.*(y(b)-y(b-1));
else
   ul(1)=x(2)-x(1);
   vl(1)=y(2)-y(1);
end   
for kl=1:bd
   offset(kl)=-0.2;
   if day(kl)>9
      offset(kl)=-0.5;
   end
end   
vabs=0.001+sqrt(ul.^2+vl.^2);
dx=15.*(ax/3200).*textdfacsat.*((vl./vabs)+offset);
dy=15.*(ay/1600).*textdfacsat.*(-ul./vabs);
xdb=xd+dx;
ydb=yd+dy;
g=text(xdb,ydb,date);
set(g,'fontname',gfont,'fontweight',gfontweight,'fontsize',textfsizesat,'color',textcolorsat)
end
ua=x(a(1,2))-x(a(1,2)-1);
va=y(a(1,2))-y(a(1,2)-1);
uabs=sqrt(ua.^2+va.^2)+0.001;
ua=10*(ax/3600)*ua/uabs;
va=10*(ay/1600)*va/uabs;
xq=x(a(1,2));
yq=y(a(1,2));
%
g1=quiver(xq,yq,ua,va,10);
set(g1,'maxheadsize',5);
set(g1,'showarrowhead','off')
set(g1,'showarrowhead','on')
set(g1,'linewidth',thick,'color',lc);
if strcmp(cplot,'y')
   plot(xd,yd,'.','color',ccol,'MarkerSize',csize)   
end   
%
end
%
if ntrack == 2
    iq=np2;
    clear x y
    x=x2;
    y=y2;
    v=v2;
    pres=pres2;
    vpot=vpot2;
    shear=shear2;
    monthm=monthm2;
    daym=daym2;
    hour=hour2;
    vtrans=vtrans2;
n=0;
m=0;
iswitch=0;
for i=1:iq
    if hour(i) <= 0.5
        m=m+1;
    end
end
xd=zeros(1,m);
yd=zeros(1,m);
day=zeros(1,m);
month=zeros(1,m);
date=cellstr(num2str(daym(1:m)));
m=0;
for i=1:iq
    if hour(i) <= 0.5
       m=m+1;
       xd(m)=x(i);
       yd(m)=y(i);
       day(m)=daym(i);
       month(m)=monthm(i);
       daystring=num2str(day(m));
       date(m)=cellstr(daystring);
    end      
end
%
if strcmp(cspline,'y')
   a=size(x);
   b=a(1,2);
   bx=1:b;
   tint=floor(6/trackres);
   tintin=1/tint;
   bxm=1:tintin:b;
   x=spline(bx,x,bxm);
   y=spline(bx,y,bxm);
   v=spline(bx,v,bxm);
   pres=spline(bx,pres,bxm);
   vpot=spline(bx,vpot,bxm);
   shear=spline(bx,shear,bxm);
   monthm=spline(bx,monthm,bxm);
   daym=spline(bx,daym,bxm);
   hour=spline(bx,hour,bxm);
   vtrans=spline(bx,vtrans,bxm);
end   
a=size(x);
b=a(1,2);
if colorcoded > 0
    xt=zeros(1,2);
    yt=zeros(1,2);
   for i=2:b
      xt(1)=x(i-1);
      xt(2)=x(i);
      yt(1)=y(i-1);
      yt(2)=y(i);
      vta=(1852/3600)*v(i-1);
      if colorcoded == 1
         if vta < 33
            vtcolor=vcolor(1,:);
         elseif vta >= 33 && vta < 45
            vtcolor=vcolor(round(0.2*qcolor),:);
         elseif vta >= 45 && vta < 54.4
            vtcolor=vcolor(round(0.47*qcolor),:);
         elseif vta >= 54.4 && vta < 62.4
            vtcolor=vcolor(round(0.6*qcolor),:);
         elseif vta >= 62.4 && vta < 69.5
            vtcolor=vcolor(round(.77*qcolor),:);
         else
            vtcolor=vcolor(round(0.95*qcolor),:);
         end
      else   
         vtr=round(0.4267.*(qcolor/64)*(3600.*vta./1852-30)); 
         %vtr=round(vta)-15;
         %vtr=round(vtr*(qcolor/64));
         vtr=min(vtr,qcolor);
         vtr=max(vtr,1);                 
         vtcolor=vcolor(vtr,:);
      end   
      plot(xt,yt,'color',vtcolor,'linewidth',thick)
      lc=vtcolor;
   end   
else
   plot(x,y,lc,'linewidth',thick)
end
bd=size(xd);
bd=bd(1,2);
if dplot == 'y'  
    offset=zeros(1,bd);
if bd > 1   
    ul=zeros(1,bd); 
    vl=zeros(1,bd);
   for kl=2:bd-1
      ul(kl)=0.5.*(xd(kl+1)-xd(kl-1));
      vl(kl)=0.5.*(yd(kl+1)-yd(kl-1));
   end
   ul(1)=4.*(x(2)-x(1));
   vl(1)=4.*(y(2)-y(1));
   ul(bd)=4.*(x(b)-x(b-1));
   vl(bd)=4.*(y(b)-y(b-1));
else
   ul(1)=x(2)-x(1);
   vl(1)=y(2)-y(1);
end   
for kl=1:bd
   offset(kl)=-0.2;
   if day(kl)>9
      offset(kl)=-0.5;
   end
end   
vabs=0.001+sqrt(ul.^2+vl.^2);
dx=15.*(ax./3200).*textdfacsat.*((vl./vabs)+offset);
dy=15.*(ay./1600).*textdfacsat.*(-ul./vabs);
xdb=xd+dx;
ydb=yd+dy;
g=text(xdb,ydb,date);
set(g,'fontname',gfont,'fontweight',gfontweight,'fontsize',textfsizesat,'color',textcolorsat)
end
%
ua=x(a(1,2))-x(a(1,2)-1);
va=y(a(1,2))-y(a(1,2)-1);
uabs=sqrt(ua^2+va^2)+0.001;
ua=10*(ax/3600)*ua/uabs;
va=10*(ay/1600)*va/uabs;
xq=x(a(1,2));
yq=y(a(1,2));
g1=quiver(xq,yq,ua,va,10);
set(g1,'maxheadsize',5);
set(g1,'showarrowhead','off')
set(g1,'showarrowhead','on')
set(g1,'linewidth',thick,'color',lc);
if cplot == 'y'
   plot(xd,yd,'.','color',ccol,'MarkerSize',csize)   
end   
%
if colorcoded > 0 && cbar == 'y'
   if colorcoded == 1
      ucolor=[vcolor(1,:);vcolor(round(0.2*qcolor),:);vcolor(round(0.47*qcolor),:); ...
          vcolor(round(0.6*qcolor),:);vcolor(round(0.77*qcolor),:);vcolor(round(0.95*qcolor),:)];
      aa={'TS';'1';'2';'3';'4';'5'};
      bb=1/6;
      tstring='Category';
   else
      ucolor=vcolor;
      aa=40:10:180;
      bb=1/14;
      tstring=wunitslabel;
   end   
   colormap(ucolor);
   rts=colorbar('horiz');
   posit=get(gca,'position');
   rr=get(rts,'position');
   xl=get(gca,'xlim');
   yl=get(gca,'ylim');
   yfrac=(yl(2)-yl(1))/(xl(2)-xl(1));
   rr2=rr;
   rr2(3)=0.75*rr(3);
   rr2(3)=max(rr2(3),0.4);
   rr2(1)=posit(1)+0.5*posit(3)-0.5*rr2(3);
   rr2(4)=0.015;
   %
   set(rts,'position',rr2)
   set(rts,'ticklength',0)
   set(rts,'fontname',gfont,'fontweight','bold')
   set(rts,'xtick',0:bb:1,'xticklabel',aa)
   set(gca,'position',posit)
   set(get(rts,'xlabel'),'String',tstring,'fontweight','bold')
end
%
end
title(tit3,'fontsize',11,'fontweight','bold','interpreter','none')
%
% Label map
%
xbeg=(incr-mod(displacement,incr))*ax/360;
xend=ax*(1-incr/360)+xbeg;
if xbeg == 0
    xend=ax;
end    
xtic=xbeg:incr.*ax/360:xend;
xtic(1)=xbeg+0.005.*ax;
set(gca,'xtick',xtic);
xinc=displacement+incr:incr:360+displacement;
xsize=max(size(xinc));
xtlab=cell(1,xsize);
for i=1:xsize
    if xinc(i) <= 180
        xtlab(i)={strcat(num2str(xinc(i)), 'E')};
    elseif xinc(i) > 180 && xinc(i) <360
        xtlab(i)={strcat(num2str(360-xinc(i)),'W')};
    else
        xtlab(i)={strcat(num2str(xinc(i)-360),'E')};
    end
end    
set(gca,'xticklabel',xtlab)
ytic=0:incr.*ay/180:ay;
ytic(1)=0.005.*ay;
set(gca,'ytick',ytic)
yinc=90:-incr:-90;
ysize=max(size(yinc));
ytlab=cell(1,ysize);
for i=1:ysize
    if yinc(i) > 0
        ytlab(i)={strcat(num2str(yinc(i)), 'N')};
    elseif yinc(i) < 0
        ytlab(i)={strcat(num2str(-yinc(i)),'S')};
    else
        ytlab(i)={num2str(yinc(i))};
    end
end    
set(gca,'yticklabel',ytlab)
set(gca,'fontname',gfont,'fontweight','bold')
set(gca,'tickdir','out')   
%
% Add grid lines to map
%
xline=zeros(max(size(ytic)),ax);
yline=zeros(max(size(ytic)),ax);
for i=1:max(size(ytic))
    xline(i,1:ax)=1:ax;
    yline(i,1:ax)=ytic(i);
end
xline=transpose(xline);
yline=transpose(yline);
hold on
line(xline,yline,'linestyle',gridlstyle,'linewidth',glinewidth,'color',gridcolor)
clear xline yline
xline=zeros(max(size(xtic)),ay);
yline=zeros(max(size(xtic)),ay);
hold off
hold on
for i=1:max(size(xtic))
    yline(i,1:ay)=1:ay;
    xline(i,1:ay)=xtic(i);
end
xline=transpose(xline);
yline=transpose(yline);
line(xline,yline,'linestyle',gridlstyle,'linewidth',glinewidth,'color',gridcolor)
hold off
%
xmin=floor(min(x));
xmax=ceil(max(x));
ymin=floor(min(y));
ymax=ceil(max(y));
delx=xmax-xmin;
dely=ymax-ymin;
xmin=xmin-bfac*delx;
xmin=max(xmin,1);
xmax=xmax+bfac*delx;
xmax=min(xmax,ax);
ymin=ymin-bfac*dely;
ymin=max(ymin,1);
ymax=ymax+bfac*dely;
ymax=min(ymax,ay);
delx=xmax-xmin;
dely=ymax-ymin;
delx2=aspect*dely;
dely2=delx/aspect;
if delx2 > delx
    delxf=delx2;
    delyf=dely;
else
    delxf=delx;
    delyf=dely2;
end    
xmaxn=0.5*(xmin+xmax+delxf);
xminn=0.5*(xmin+xmax-delxf);
xmax=xmaxn;
xmin=xminn;
ymaxn=0.5*(ymin+ymax+delyf);
yminn=0.5*(ymin+ymax-delyf);
ymax=ymaxn;
ymin=yminn;
xmin=floor(xmin);
xmax=ceil(xmax);
ymin=floor(ymin);
ymax=ceil(ymax);
xmax=min(xmax,ax);
set(gca,'xlim',[xmin xmax],'ylim', [ymin ymax]) 
%
% Add button to access documentation of script
%
if strcmp(Docdisp,'y')
    try
        delete(hButton)
    catch
    end    
    tempname='trmap';
    mname = strcat('graphics_docs/html/',tempname,'_doc.html');
    comstr=['open ',mname,';delete(hButton)'];
    labelStr = '<html><center><a href="">Doc';
    hButton = uicontrol('string',labelStr,'pos',[10,10,50,20],'callback',comstr); 
end        