%computes the tide level before and after a decrease in temperature of at
%least rateTc.  Returns allOut={indices,risingOut,falseOut}.  The indices
%are for all critical temperature decreases (TD*s). Rising out is the tide level
%for TD*'s when the tide is increasing. FalseOut is for TD*'s when the tide
%isn't incoming.  ESL is called from averageESL, or you can use it alone.
%e.g. ESL(tempData, tideData, rateTc), and you can add an extra column if you want to
%include say times, or waveheights.




%tempData=array of sequential temperatures at different sites: one site per column.
%tideData=array of tide Data corresponding to the same times as tempData. to get tide data
%rateTc=critical rate in degrees per sample interval.
%you can interpolate a tide record that coincides with your temperature data by using ftinterp or interp1 
%in matlab and tide data from NOAA/NOS http://co-ops.nos.noaa.gov/data_res.html 



function allOut=ESL(varargin)

if nargin==3,
    tempData=varargin{1}; 
    tideData=varargin{2};
    rateTc=varargin{3};
    otherData=varargin{1}; %if no other data, then use tideData as otherData
elseif nargin==4,
    tempData=varargin{1}; 
    tideData=varargin{2};
    rateTc=varargin{3};
    otherData=varargin{4};
else
    disp('wrong number of arguments for ESL');
    help ESL
    return 
end
    
    
    
    
indices={}; %cell array for cooling indices.
indicesRising={}; %subset of indices where there is a risinig tide
%computing temperature difference matrix
tempDataOffset=tempData(1:size(tempData,1)-1,:); %temp matrix without last row
dTemp=tempData(2:size(tempData),:)-tempDataOffset; 
      %tempDataShifted - [temp matrix without first row]
      %e.g., if tempData=[1,2;3,4;5,6;7,8], dTemp=[3,4;5,6;7,8]-[1,2;3,4;5,6]

%get indices if cooling faster than rateTc, put in cell array indices      
for i=1:size(dTemp,2)
    indices{i}=find(dTemp(:,i)<-rateTc);
        %2 columns, first is the high temp, and second 
        %is the low temp, so
        %tempData(indices{i}(:,2),1)-tempData(indices{i}(:,1),1)is more negative than the critical rate
        %if the change is due to the tide coming in,
        %tideData(indices{1}(:,2),1)>tideData(indices{1}(:,1),1), and ESL
        %lies between the two.
        tideOut{i}=[tideData(indices{i}),tideData(indices{i}+1)]; %the tide levels bracketing the temperature change.    
        
        otherOut{i}=otherData(indices{i}); %
        
end

for j=1:size(dTemp,2);
    if size(tideOut{j},1)>0
        diff=tideOut{j}(:,2)-tideOut{j}(:,1);
        ok=find(diff>=0); %only ESL if it's on a rising tide.  
        notOK=find(diff<0);
        risingOut{j}=tideOut{j}(ok,:);
        indicesRising{j}=indices{j}(ok);
        risingOther{j}=otherOut{j}(ok,:);
        falseOut{j}=tideOut{j}(notOK,:);
    else 
        risingOut{j}=[];
        falseOut{j}=[];
    end
end
allOut={indicesRising,risingOut,falseOut,risingOther};