%try to do spectrogram of Hs
clear all; fclose all;


%try #1 - computes Spectra of non Nan data and then interpolate to a
%constantn df
Spectro_freq=ones(31*24,62)*NaN;
Spectro_Hs=ones(31*24,62)*NaN;
doy_all=ones(31*24,1)*NaN;
for ddd=311:341
    jd=sprintf('%03i',ddd);
    for hhh=00:23
        ic=(ddd-311)*24+hhh+1;
        [ddd hhh ic]
        hr=sprintf('%02i',hhh);
        if exist(['C:\Data\CRUISES\DYNAMO_2011\data\Riegl_corrected\Riegl_heavesplineinterpolatedtowave_v012012_leg3\Wavecorrected_',jd hr,'.prn'])
            Hs_raw=load(['C:\Data\CRUISES\DYNAMO_2011\data\Riegl_corrected\Riegl_heavesplineinterpolatedtowave_v012012_leg3\Wavecorrected_',jd hr,'.prn']);
            Hs=Hs_raw(:,end);
    %         doy_all(ic)=Hs_raw(:,1)+Hs_raw(:,2)/24+Hs_raw(:,3)/3600/24;
            doy_all(ic)=ddd+hhh/24;
            [SHs,FHs]=psd2(detrend(Hs(~isnan(Hs))),length(Hs(~isnan(Hs))),10,hamming(length(Hs(~isnan(Hs)))));[SHs1,FHs1]=specsmoo(SHs,10); 
            Spectro_freq(ic,1:length(SHs1))=FHs1;
            Spectro_Hs(ic,1:length(SHs1))=SHs1;
        end    
    end
end    

Spectro_freq(isnan(doy_all),:)=[];
Spectro_Hs(isnan(doy_all),:)=[];
doy_all(isnan(doy_all))=[];


Hs_heave2wavespline_leg3=load('C:\Data\CRUISES\DYNAMO_2011\data\Riegl_corrected\Riegl_heavesplineinterpolatedtowave_v012012_leg3\Hs_hourly_leg3.prn');
doy_leg3=Hs_heave2wavespline_leg3(:,1)+Hs_heave2wavespline_leg3(:,2)/24;
ii=find(abs(Hs_heave2wavespline_leg3(:,7))<0.06  & abs(Hs_heave2wavespline_leg3(:,5))<1 & Hs_heave2wavespline_leg3(:,6) >10);

badpts=find(abs(doy_leg3-324.4)==min(abs(doy_leg3-324.4)))-1;
Hs_heave2wavespline_leg3(badpts,:)=NaN;Spectro_freq(badpts,:)=ones(1,62)*NaN;Spectro_Hs(badpts,:)=ones(1,62)*NaN;
badpts=find(abs(doy_leg3-315.3)==min(abs(doy_leg3-315.3)));
Hs_heave2wavespline_leg3(badpts,:)=NaN;Spectro_freq(badpts,:)=ones(1,62)*NaN;Spectro_Hs(badpts,:)=ones(1,62)*NaN;
figure;plot(doy_leg3(ii),Hs_heave2wavespline_leg3(ii,3),'b.--',doy_leg3(ii),Hs_heave2wavespline_leg3(ii,4),'r-','LineWidth',2);grid;

doy_allo=repmat(doy_all,1,62);


% doy=Hs_raw(:,1)+Hs_raw(:,2)/24+Hs_raw(:,3)/3600/24;
% 
% ii=find(~isnan(Hs));
% Hs=Hs_raw(:,end);
% Hsi=interp1(doy(ii),Hs(ii),doy,'spline'); Hsi(end-20:end)=(Hsi(end-21)); Hsi(1:10)=Hsi(11);
% figure;plot(doy,Hs,'.',doy,Hsi,'r--');
% 
% [SHs,FHs]=psd2(detrend(Hs(~isnan(Hs))),length(Hs(~isnan(Hs))),10,hamming(length(Hs(~isnan(Hs)))));[SHs1,FHs1]=specsmoo(SHs,10); 
% [SHsi,FHsi]=psd2(detrend(Hsi(~isnan(Hsi))),length(Hsi(~isnan(Hsi))),10,hamming(length(Hsi(~isnan(Hsi)))));[SHsi1,FHsi1]=specsmoo(SHsi,10); 
% 
% figure;clf;loglog(FHs1,FHs1.*SHs1,'b',FHs1,FHs1.*SHs1,'r')

Spectro_freqi=Spectro_freq(38,:);
Spectro_Hsi=ones(31*24,62)*NaN;
for oo=1:length(Spectro_freq)
    try
        Spectro_Hsi(oo,:)=interp1(Spectro_freq(oo,~isnan(Spectro_freq(oo,:))),Spectro_Hs(oo,~isnan(Spectro_freq(oo,:))),Spectro_freqi,'spline');
    catch
        Spectro_Hsi(oo,:)=ones(1,62)*NaN;
    end    
end    

% doy_all(isnan(Spectro_Hsi(:,end)))=[];
Spectro_Hsi(isnan(Spectro_Hsi(:,end)),:)=[];



stdt=floor(doy_all(1)+datenum(2011,1,0));
endt=ceil(doy_all(end)+datenum(2011,1,0)); 
xs=2;
% figure;imagesc(doy_all,Spectro_freqi, real(log10(Spectro_Hsi)),[-6 -1]);  %
figure;imagesc(doy_all+datenum(2011,1,0),Spectro_freqi, real(log10(Spectro_Hsi)),[-3 1]);  %
axis xy;ylim([0 .5]);
xlabel('November (2011)');ylabel('log_{10}Hs(m^2/s)')
set(gca(gcf),'XTick',stdt:xs:endt);datetick('x','mm/dd','keepticks','keeplimits');
h = colorbar;
axes(h)


lud=find(floor(doy_all)==328);
figure;loglog(Spectro_freq(1,:), Spectro_Hs(lud(1),:),Spectro_freqi,Spectro_Hsi(lud(1),:),'r')

        