% CDIP MONTHLY SEA SURFACE TEMPERATURE (SST) - MATLAB

% Plots a single month (or any specified time-period) of Sea Surface Temperature (SST) data for selected station
% See http://cdip.ucsd.edu/themes/cdip?pb=1&bl=cdip?pb=1&d2=p1&d2=p70&u3=s:093:st:1:v:temperature for example monthly SST plot.

%-----------------------------------------------------------------------
% To access NetCDF file directly from THREDDS:
%	- Download the nctoolbox from: https://github.com/nctoolbox
%		1. Download and extract into a folder
%		2. Open matlab and at the matlab prompt cd to the extracted folder where setup_nctoolbox lives. (Use pwd to see where you currently are).
%		3. Run setup_nctoolbox
%		4. Optionally edit startup.m and add lines in as described in the install instructions.
%-----------------------------------------------------------------------

% CHANGE DIRECTORY TO 'NCTOOLBOX' DIRECTORY
cd nctoolbox

% OPEN NCTOOLBOX IN MATLAB
setup_nctoolbox



% USER ENTERS STATION NUMBER AND START/END DATES FOR PLOT
stn = '093';
startdate = '01/01/2011 00:00';
enddate = '01/31/2011 23:59';



% CONNECT TO THREDDS SERVER AND OPEN NETCDF FILE
urlbase = 'http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/';  % Set 'base' THREDDS URL and pieces to concatenate with user-defined station number (set above)
p1 = 'p1/';
urlend = 'p1_historic.nc';
dsurl = strcat(urlbase,stn,p1,stn,urlend);
ds = ncdataset(dsurl);

% PRINT LIST OF VARIABLES IN NETCDF FILE
% ds.variables

% GET BUOY NAME AND TRANSPOSE TO HORIZONTAL STRING
buoyname = ds.data('metaStationName');
buoytitle = transpose(buoyname);

% EXTRACT MONTH/YEAR FROM 'STARTDATE' 
formatOut = 'mmm yyyy';
monthtitle = datestr(datenum(startdate,'mm/dd/yyyy'),formatOut);


% CONVERT START/END DATES TO MATLAB SERIAL NUMBERS
startser = datenum(startdate, 'mm/dd/yyyy HH:MM');
endser = datenum(enddate, 'mm/dd/yyyy HH:MM');


% CALL TIME VARIABLE 
timevar = ds.data('sstTime');
timeconvert = ds.time('sstTime',timevar); % Convert UNIX timestamps to Matlab serial units

% FIND INDEX NUMBERS OF CLOSEST TIME ARRAY VALUES TO START/END DATES
diffstart = abs(timeconvert-startser);
[startidx startidx] = min(diffstart); % index of closest value
closeststart = timeconvert(startidx); % closest value

diffend = abs(timeconvert-endser);
[endidx endidx] = min(diffend); % index of closest value
closestend = timeconvert(endidx); % closest value


% ASSIGN NAMES TO VARIABLES TO BE USED IN COMPENDIUM PLOT, AND SPECIFY TIMERANGE
timecut = ds.time('sstTime',ds.data('sstTime',startidx,endidx)); % Create a subset 'time' variable from startdate to enddate, and convert to Matlab serial units
timedt = datetime(timecut,'ConvertFrom','datenum'); % Convert Matlab serial units 'timecut' section to datetime objects

sst = ds.data('sstSeaSurfaceTemperature',startidx,endidx);




%% SET UP FIGURE FOR PLOTTING
   
% MAKE SUBPLOTS AND ADD DATA

xvals = timedt; % define x-axis variable (time)
yvals = sst; % define y-variable (sst)

figure % Create new figure

% SUBPLOT
subplot(1,1,1)
plot(xvals,yvals,'DatetimeTickFormat','dd')
ylim([0 25])
grid on
grid minor
title(monthtitle)
ylabel('SST, C')
xlabel('Day')

suptitle(buoytitle) % Set main plot title