% met_ts.m  An mfile to plot the timeseries of South Carolina
%           air temperature, sea temperature, barometric pressure,
%           significant wave height, and average wave period from NDBC
%           moored buoys, NDBC C-MAN stations, and NOS water level stations.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 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:
%       the netcdf toolbox
%       julian.m
%       timeplt.m
%       dolandscape.m
%       fixpaper.m
%       suptitle.m
%       smart_interp.m

clear
more off

% define data and output directories
dataDirectory = 'D:\BNA\data\Carolinas\webpage\dataDVD\Hatteras2009DVD\supporting data';
plotDirectory = 'D:\BNA\data\Carolinas\webpage\dataDVD\Hatteras2009DVD\PLOTS\MET';

% define data file names
dataFiles = { 'http://dods.ndbc.noaa.gov/thredds/dodsC/data/stdmet/41025/41025h2009.nc' };   % NDBC (met) buoy station
sta_type = 'NDBC';
sta_id = '41025';
nFiles = length(dataFiles);

% define parameters for plots
ncvarnames = {'air_temperature','sea_surface_temperature','air_pressure','wave_height','average_wpd'};
names = {'air_temp','sea_temp','bp','Hs','T'};
startTime=julian(2009,1,8,00);
stopTime=julian(2009,05,15,00);
jd_int = startTime:1/24:stopTime;
ndt = length(jd_int);
col = get(0,'defaultaxescolororder');     %plot color order
ylab = { {'Average Wave','Period (sec)'},...   %plot ylabels
         {'Significant','Wave Height (m)'},...
         {'Air','Pressure (hPa)'},...
         {'Sea Surface','Temperature','(degrees C)'},...
         {'Air','Temperature','(degrees C)'} };
timediff = stopTime-startTime;                                %plot xlimits extend 5 percent
xlims = [(startTime-timediff*0.05) (stopTime+timediff*0.05)]; %past start time and stop time
ylims = [ 2 10;...   %plot ylimits
          0 6;...
       980 1040;...
          0 30;...
          -10 30];
xt=[];
for yy=2009:2009
    for mm=1:12
        tic=julian(yy,mm,01,00);
        if tic>=xlims(1) & tic<=xlims(2)
            xt=[xt; tic];
        end
    end
end

% initialize variables
air_temp_all = nan(ndt,nFiles);
sea_temp_all = nan(ndt,nFiles);
bp_all = nan(ndt,nFiles);
Hs_all = nan(ndt,nFiles);
T_all = nan(ndt,nFiles);

for f=1:nFiles
    
        nc = netcdf(dataFiles{f});
        sta_type = 'NDBC';
        sta_id = '41025';
        station = [sta_type,' ',sta_id]; %station ID
        
        % get time
        jd = (nc{'time'}(:)/3600/24+julian(1970,01,01,0));
        metidx(1) = find(jd(:) >= startTime,1,'first');
        metidx(2) = find(jd(:) <= stopTime,1,'last');
        jd = jd(metidx(1):metidx(2));
        
        % get data, nan the bad stuff, and interpolate onto 
        % standard time array
        for n = 1:length(ncvarnames);
            names{n}
            ncobj = nc{ncvarnames{n}};
            if isempty(ncobj)
                data_int = nan(size(jd_int'));
            else
                theFillValue = ncobj.FillValue_(:)
                data = ncobj(metidx(1):metidx(2));
                bads = find(data >= theFillValue);
                data(bads) = nan;
                data_int = smart_interp(jd, data, jd_int, 3)';
            end
            eval([names{n},'= data_int;'])
        end
        nc = close(nc);


        % call timeplt.m to create individual plots
        % (clearing global timeplt*, the ResizeFcn,
        % and the Tag are necessary to release
        % timeplt.m's control over your figure!)
        xdata = jd_int;
        ydata = [T Hs bp sea_temp air_temp];
        panelOrd = [1 2 3 4 5]; 
        figure
        han = timeplt2(xdata,ydata,panelOrd);
        clear global timeplt*
        set(gcf,'ResizeFcn',[])
        set(gcf,'Tag','')
        dolandscape
        fixpaper
 
        % reformat xlims, ylims, ylabels, and line color
        for h = 1:length(han)
            axes(han(h))
            grid on
            set(gca,'xlim',xlims,'xtick',xt)
            set(gca,'ylim',ylims(h,:))
            ylabhan = get(gca,'Ylabel');
            set(ylabhan,'String',ylab{h})
            pos=get(ylabhan,'position');
            linehan1 = findobj(gca,'Type','line','color','b');
            set(linehan1,'color',col(f,:))
            linehan2 = findobj(gca,'Type','line','color','k');
            if ~isempty(linehan2)
                delete(linehan2)
            end
                  
        end

        % reformat xlabels
        axes(han(1))
        oldhan = findobj(gca,'Tag','xtick labels');
        delete(oldhan);
            set(gca,'XTickLabel',{''})

        for tic = 1:length(xt)
            xtg = gregorian(xt(tic));
            xtd = datenum(xtg);
            xtl = {[datestr(xtd,3),' 1'];datestr(xtd,10)}; %2 line label
            text(xt(tic),1,xtl,'horizontalalignment','center',...
                                 'verticalalignment','top')
        end
        
        % give it a title and export to .pdf
        set(gcf,'PaperPositionMode','Auto')
        suptitle({station,' Temperature, Pressure, and Wave Observations'});
        print(fullfile(plotDirectory,['temp_press_wvs_',sta_id]),'-dpdf')

        % get data ready for combined plot
        air_temp_all(:,f) = air_temp;
        sea_temp_all(:,f) = sea_temp;
        bp_all(:,f) = bp;
        Hs_all(:,f) = Hs;
        T_all(:,f) = T;

end            
