% wind_ts.m  An mfile to plot the timeseries of North Carolina wind data
%            from NDBC moored buoys, NDBC C-MAN stations, NOS water level
%            stations, and the COAMPS wind model.

%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Use of this program is described in:
%
% Armstrong, B.N., Warner, J.C., Voulgaris, G., List, J.H., Thieler, E.R., 
% Martini, M.M., Montgomery, E.T., 2011, Carolinas Coastal Change Processes 
% Project Data Report for Observations near Diamond Shoals, NC from January-April 
% 2009: U.S. Geological Survey Open-File Report 2011-1156.
%
% Program written in Matlab v7.1.0 SP3
% Program ran on PC with Windows XP Professional OS.
%
% "Although this program has been used by the USGS, no warranty, 
% expressed or implied, is made by the USGS or the United States 
% Government as to the accuracy and functioning of the program 
% and related program material nor shall the fact of distribution 
% constitute any such warranty, and no responsibility is assumed 
% by the USGS in connection therewith."
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%                     
% Dependencies:
%       julian.m

clear
more off

% define data and output directories
dataDirectory = 'D:\BNA\data\Carolinas\webpage\dataDVD\Hatteras2010DVD\files';
plotDirectory = 'D:\BNA\data\Carolinas\webpage\dataDVD\Hatteras2010DVD\PLOTS';

% define data file names
dataFiles={'41025h2010.nc'};
nFiles=length(dataFiles);
name={'NDBC 41025'};
name2={'NDBC'};
%P-NAM-CONUS_12km-noaaport_best.ncd
% define parameters for plots
startTime=julian(2010,2,1,00);
stopTime=julian(2010,2,31,00);
col=get(0,'defaultaxescolororder');     %plot color order
ylab={ {'Temperature (C)'},...   %plot ylabels
       {'Pressure (hPa)'},...
       'Wave Height (m)'};
timediff=stopTime-startTime;                                %plot xlimits extend 5 percent
xlims=[(startTime-timediff*0.05) (stopTime+timediff*0.05)]; %past start time and stop time
ylims=[ 0 30;...   %plot ylimits
        990  1030;...
        0 10];
xt=[];
for yy=2010:2010
    for mm=2:2
        for dd=1:4:31
            tic=julian(yy,mm,dd,00);
            if tic>=xlims(1) & tic<=xlims(2)
                xt=[xt; tic];   % place xticks on the first of every month
            end
        end
    end
end

% initialize variables
jd=cell(nFiles,1);
temp=cell(nFiles,1);
press=cell(nFiles,1);
wvht=cell(nFiles,1);
for f=1:length(dataFiles)
    station_id=name{f};
    wind=(fullfile(dataDirectory,dataFiles{f}));
    jd=julian([ncread(wind,'Year'),ncread(wind,'Month'),...
            ncread(wind,'Day'),ncread(wind,'Hour'),ncread(wind,'Minute'),zeros(length(ncread(wind,'Minute')))]);
    windstart=find((jd(:) >= startTime),1,'first');
    windend=find((jd(:) <= stopTime),1,'last');
%     jd=jd(windstart:windend);
    temp=ncread(wind,'WaterTemp');
    temp(temp>=99)=NaN;
    press=ncread(wind,'Pressure');
    press(press>=9999)=NaN;     
    wvht=ncread(wind,'WaveHeight');
    wvht(wvht>=99)=NaN;     
        
        % create plot
        
        subplot(3,1,1)
        plot(jd, temp)
        ylabel(ylab{1});
        set(gca,'ylim',ylims(1,:))
        set(gca,'xlim',xlims,'xtick',xt,'xticklabel',[])
        for tic=1:length(xt)
            xtg=gregorian(xt(tic));
            xtd=datenum(xtg);
            xtl={datestr(xtd,6)}; %2 line label
            text(xt(tic),-2,xtl,'horizontalalignment','center',...
                                 'verticalalignment','top')
        end
        title(['NDBC Water Temperature, Pressure, and Wave Height' name(f)]);

        subplot(3,1,2)
        plot(jd,press)
        ylabel(ylab{2});
        set(gca,'ylim',ylims(2,:))
        set(gca,'xlim',xlims,'xtick',xt,'xticklabel',[])
        for tic=1:length(xt)
            xtg=gregorian(xt(tic));
            xtd=datenum(xtg);
            xtl={datestr(xtd,6)}; %2 line label
            text(xt(tic),988,xtl,'horizontalalignment','center',...
                'verticalalignment','top')
        end
        
        subplot(3,1,3)
        plot(jd,wvht)
        ylabel(ylab{3});
        set(gca,'ylim',ylims(3,:))
        set(gca,'xlim',xlims,'xtick',xt,'xticklabel',[])
        for tic=1:length(xt)
            xtg=gregorian(xt(tic));
            xtd=datenum(xtg);
            xtl={datestr(xtd,6)}; %2 line label
            text(xt(tic),-1,xtl,'horizontalalignment','center',...
                'verticalalignment','top')
        end
        % give it a title and export to .pdf
            print(fullfile(plotDirectory,[name2{f},'_tpw']),'-dpdf');            
end
