close all
NWCUSP=load('NWC_USP.txt');
VTXUSP=load('VTX_USP.txt');
PAMT=load('Pam_track_NOAA_mat.txt');
files = gunzip('gshhs_c.b.gz','./');
filename = files{1};
indexfile = gshhs(filename, 'createindex');
latM=[70 -70];
lonM=[70 210];
load coast
if (1==1)
    T=[];
    LatL=[];
    LonL=[];
    AL=[];
    for d=9:15
        filename=sprintf(['A201503/A201503%02d' '.loc'],d);
        L=load(filename)';
        T=[T d*24*3600+L(1,:)*3600+L(2,:)*60+L(3,:)];
        LatL=[LatL L(4,:)];
        LonL=[LonL L(5,:)];
        LonL(LonL<-150)=LonL(LonL<-150)+360;
        AL=[AL L(6,:)];
    end
end
k=0;
for d=9:15
    for h=0:23
        k=k+1;
        t0=d*24*3600+(h-0.5)*3600;
        II=and(T>=t0,T<t0+3600);
        %II=and(II,abs(AL)>=25);
        %II=and(II,and(LatL>=-40,LatL<10));
        II=and(II,and(LatL>=-70,LatL<70));
        %II=and(II,and(LonL>=140,LonL<200));
        II=and(II,and(LonL>=70,LonL<210));
        II=find(II);
            filename=sprintf('MTS/MTS21503%02d%02dIR1.pgm',d,h);
            filename
            JL=floor(LonL(II)-70)+1;
            IL=floor(LatL(II)+70)+1;
            %C=double(imread(filename));
            C=zeros(141,141);
            for i=1:length(IL)
                C(IL(i),JL(i))=C(IL(i),JL(i))+AL(II(i));
                C(IL(i)+1,JL(i))=C(IL(i)+1,JL(i))+AL(II(i));
                C(IL(i),JL(i)+1)=C(IL(i),JL(i)+1)+AL(II(i));
                C(1+IL(i),1+JL(i))=C(1+IL(i),1+JL(i))+AL(II(i));
            end
            %C=log10(C);
            H=imagesc(lonM,-latM,C)
            set(gca,'YDir','normal');
            cm=min(min(C));
            cM=max(max(C));
            cM=5e3;
            set(gca,'CLim',[cm cM])
            %axis([140 200 -40 10])
            axis([70 210 -70 70])
            %colormap gray
            hold
            plot(long,lat,'Color','y','Linewidth',2)
            plot(PAMT(:,2),PAMT(:,1),'Color','k','Linewidth',3 )
            plot(NWCUSP(:,2),NWCUSP(:,1),'Color','r','Linewidth',3)
            plot(VTXUSP(:,2),VTXUSP(:,1),'Color','r','Linewidth',3)
            set(gca,'DataAspectRatio',[1 1 0.1]);
            set(gca,'FontSize',20);
            %for j=II
            %    if (AL(j)>0)
            %        plot(LonL(j),LatL(j),'xb','Linewidth',3,'MarkerSize',log10(floor(abs(AL(j))+4)*10));
            %    else
            %        plot(LonL(j),LatL(j),'xb','Linewidth',3,'MarkerSize',log10(floor(abs(AL(j))+4)*10));
            %    end
            %end
            hold 
            
            title(sprintf('MTS2 IR1 channel 2015/03/%02d %02d:00 UTC',d,h));
            set(gcf,'Position',[63           1        1858        1133]);
            A=getframe(gcf);
            imwrite(A.cdata,sprintf('SAM_%04d.jpg',k));
            %pause
    end
end