root/DPWavesProc/trunk/DPWavesProc/test_code/nortek_specmultiplot_compare.m

function [Winfo]=nortek_specmultiplot(dateinfo1,dateinfo2,timestep); 

%NOTE:  This function is for the Nortek AWAC profiler!

%where dateinfo1 and 2 are the strings of the first and last date extension
%of the pressure,range,etc files and timestep is the time in hours between
%each burst

%This function takes the pressure, range, orbit and sysinfo files created
%in python and outputs a Winfo data structure with the sig wave height,
%peak period, direction of peak period and dominant dir.  It also generates
%polar plots and 2d freq and directional plots of a variety of different
%estimation methods, resolutions and the nortek quickwave output. 


%first set up the time vector going from dateinfo1 to dateinfo2
Winfo.time=[];

date1=datevec(dateinfo1, 'yymmddHHMM');
date2=datevec(dateinfo2,'yymmddHHMM');

if date1 == date2
    Winfo.time=[datenum(date1)];
else
    Winfo.time=[datenum(date1)];
    newdate=date1;
    while datenum(newdate) ~= datenum(date2)
        newdate=newdate+[0 0 0 timestep 0 0];
        Winfo.time=[horzcat(Winfo.time,datenum(newdate))];
    end
end
            
      
%set up time, which is Winfo.time as date str in the yymmddHHMM format
time=datestr(Winfo.time, 'yymmddHHMM');

%set up data structure for wave info
Winfo.setup={'EMEP','IMLM' 'nortek'}';
Winfo.hsig=[];
Winfo.hconf=[];
Winfo.peakP=[];
Winfo.dirP=[];
Winfo.Ddir=[];
Winfo.Spectrum.EMEP=[];
Winfo.Spectrum.IMLM=[];
Winfo.Spectrum.nortek=[];

%load the wave data from the nortek quick wave software
wavedata=load('BPWave01.wap');

%Load the data and run the script
for i=1:length(time(:,1))
    dateinfo=time(i,:);

    pressure=strcat('pressure_',dateinfo,'.txt');
    range=strcat('range_',dateinfo,'.txt');
    orbit=strcat('orbit_',dateinfo,'.txt');
    sysinfo=strcat('sysinfo_',dateinfo,'.txt');
    nortekspec=strcat('nortekspec_',dateinfo,'.txt');

    pressure=load(pressure);
    range=load(range);
    orbit=load(orbit);
    sysinfo=load(sysinfo);
    nortekspec=load(nortekspec);

    %set up data with standard AST freq at 0.01 and dir at 2
    [ID,SM,EP]=nortek_radialtouvw(pressure,range,orbit,sysinfo,1);

    %run diwasp to generate this spectrum
    [E_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0});

    %now for IMLM
    EP.method='IMLM';
    EP.iter=3;
    [I_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0});


    %Use waveplot to plot the polar, freq and dir plots of the spectra

    [E_WI,I_WI,nortek_WI,fig1,fig2]=nortek_waveplot(E_F01_D2,I_F01_D2,sysinfo,nortekspec);
    
    saveas(fig1,['polar1_' dateinfo '.fig']);
    saveas(fig2,['dirfreq1_' dateinfo '.fig']);
    
    close all;
   
    %Fill up the Winfo data structure

    hsig=[vertcat(E_WI.hsig,I_WI.hsig,nortek_WI.hsig)];
    hconf=[vertcat(E_WI.hconf,I_WI.hconf,nortek_WI.hconf)];
    peakP=[vertcat(E_WI.tp,I_WI.tp,nortek_WI.tp)];
    dirP=[vertcat(E_WI.dtp,I_WI.dtp,nortek_WI.dtp)];
    Ddir=[vertcat(E_WI.dp,I_WI.dp,nortek_WI.dp)];
    
    Winfo.hsig=[horzcat(Winfo.hsig,hsig)];
    Winfo.hconf=[horzcat(Winfo.hconf,hconf)];
    Winfo.peakP=[horzcat(Winfo.peakP,peakP)];
    Winfo.dirP=[horzcat(Winfo.dirP,dirP)];
    Winfo.Ddir=[horzcat(Winfo.Ddir,Ddir)];
    Winfo.Spectrum.EMEP=[cat(3,Winfo.Spectrum.EMEP,E_F01_D2.S)];
    Winfo.Spectrum.IMLM=[cat(3,Winfo.Spectrum.IMLM,I_F01_D2.S)];
    Winfo.Spectrum.nortek=[cat(3,Winfo.Spectrum.nortek,nortekspec)];
    
end

%now create the time series plots of the wave info
[fig_sigh,fig_pp,fig_dp,fig_dd]=nortek_Winfo_plot(Winfo);

saveas(fig_sigh,['sigh_' dateinfo '.fig']);
saveas(fig_pp,['peakp_' dateinfo '.fig']);
saveas(fig_dp,['dirpeak_' dateinfo '.fig']);
saveas(fig_dd,['domdir_' dateinfo '.fig']);

close all;