function get_pressure_diffs
%
%   get_pressure_diffs
%
%   Function to look at the pressure data and calculate the difference in
%   pressure between high and low tide.

%Soupy Alexander, 4/2004

%This Matlab m-file was used to create portions of U.S. Geological Survey 
%Data Series 85. 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.

%Citation:  Butman, Bradford, Alexander, P. Soupy, Anderson, S.P., 
%Lightsom, F.L., Scotti, Alberto, and Beardsley, R.C., 2004, The 
%Massachusetts Bay Internal Wave Experiment, August 1998:  Data Report:  
%U.S. Geological Survey Data Series 85, 1 DVD-ROM.

presFile = netcdf('5371tcp-a.nc', 'nowrite');
presData = presFile{'P_4023'}(:);
presTime = singlejd(presFile{'time'}(:), presFile{'time2'}(:));

[theNums, theText] = xlsread('arrival_times_summary.xls');
eventNumber = theNums(:,1);
startTime = datenum2julian(datenum('30-Dec-1899') + theNums(:,5));

theBads = find(isnan(eventNumber));
eventNumber(theBads) = [];
startTime(theBads) = [];
        
fID = fopen('pres_diffs.xls', 'w');

for indexEvent = 1:length(eventNumber)
    thisTime(indexEvent) = value2Index(presTime, startTime(indexEvent));
    
    theWindow = (thisTime(indexEvent)-10):(thisTime(indexEvent)+10);
    
    while min(presData(theWindow)) ~= presData(thisTime(indexEvent))
        thisTime(indexEvent) = min(theWindow(find(presData(theWindow) == ...
            min(presData(theWindow)))));
        theWindow = (thisTime(indexEvent)-10):(thisTime(indexEvent)+10);
    end
    
    lowPres(indexEvent) = presData(thisTime(indexEvent));
    
    thisTimeHigh(indexEvent) = thisTime(indexEvent)-30;
    
    theWindow = (thisTimeHigh(indexEvent)-10):(thisTimeHigh(indexEvent)+10);
    
    while max(presData(theWindow)) ~= presData(thisTimeHigh(indexEvent))
        thisTimeHigh(indexEvent) = min(theWindow(find(presData(theWindow) == ...
            max(presData(theWindow)))));
        theWindow = (thisTimeHigh(indexEvent)-10):(thisTimeHigh(indexEvent)+10);
    end
    
    highPres(indexEvent) = presData(thisTimeHigh(indexEvent));
    theDiff(indexEvent) = highPres(indexEvent) - lowPres(indexEvent);
    theDiffTime(indexEvent) = presTime(thisTime(indexEvent)) - ...
        presTime(thisTimeHigh(indexEvent));
    
    fprintf(fID, [num2str(eventNumber(indexEvent), '%.0f') '\t' ...
            datestr(julian2datenum(presTime(thisTime(indexEvent))), 0) '\t' ...
            datestr(julian2datenum(presTime(thisTimeHigh(indexEvent))), 0) '\t' ...
            num2str(theDiff(indexEvent), '%.3f') '\n']);
end
    