eb_dwelltimes_traceInd_MN.m

function [ fractionHighAll, thigh, tlow, thigh_Ind, tlow_Ind, histTlow, histThigh, cumTlow, cumThigh ] = eb_dwelltimes_traceInd_MN(path, secPerFrame)

includeSingleSpikes=false;
binSizeLow=1;
binSizeHigh=5;
numTracesNotFluct=0; % I want to count the number of traces that vbFRET found no fluctuations in.
%now I go through and count Tlow and Thigh so I can make a histogram
tlow(1)=0;
thigh(1)=0; %these are the durations of the low fret and high fret states (in frames)
nTlow=0;
nThigh=0;%counters for the number of high and low states
filter=1; %increase this to set a lower limit for trace duration in frames
%%%These variables are for comparing traces
tracesLow=zeros(10,2); %these store the average low time (frames) of traces i in (i,1)
                %and the error in the mean of the low time in (i,2)
tracesHigh=zeros(10,2); %same as above with high times.

%This is for storing the fraction of frames that are high fret
%and the weights which are the lengths of the traces in frames.  
fractionHigh=zeros(1,size(path,2));
weightsFractionHigh=zeros(1,size(path,2));

for traceInd=1:size(path,1)
    %These are for getting average fraction high values
    totalPoints=0;
    totalHighPoints=0;
    nTlowFromThisTrace=0; %used to figure out if this trace fluctuated at all
    maxFRET=max(path{traceInd});
    minFRET=min(path{traceInd});
    
    currentState=0;%fret states are 1(low) or 2(high).  0 is what I'm calling the initial state (not included in statistics)
    lastChange=0;%index of frame where there was a state change 111222 -> lastChange=4
    
    tmpTraceLow=zeros(1,1);%these are arrays that hold all the low times for this trace. 
                %It is used to build tracesLow and is re-written for each
                %trace during the loop.
    tmpTraceHigh=zeros(1,1);
    nTmpTraceLow=0; %Same purpose as nTLow and nTHigh except for tmpTraceLow
    nTmpTraceHigh=0;
    
    for timeInd=2:size(path{traceInd},1)
        totalPoints=totalPoints+1;
        if path{traceInd}(timeInd)== maxFRET
            totalHighPoints=totalHighPoints+1;
        end
        if path{traceInd}(timeInd)~=path{traceInd}(timeInd-1)&&... %This is a change index
                    ((includeSingleSpikes)||...  %and we're including single spikes
                     (timeInd-lastChange>1))  %or if the change is greater than one
            if currentState==minFRET %change from a low state.
                if timeInd-lastChange < filter
                    currentState=path{traceInd}(timeInd);
                    lastChange=timeInd;
                    continue
                end
                nTlow=nTlow+1;
                nTlowFromThisTrace=nTlowFromThisTrace+1;
                tlow(nTlow)=timeInd-lastChange;
                tlow_Ind(nTlow) = traceInd;
                 
                nTmpTraceLow=nTmpTraceLow+1;
                tmpTraceLow(nTmpTraceLow)=timeInd-lastChange;
           
            elseif currentState==maxFRET % changing from a high state
               
                nThigh=nThigh+1;
                
                thigh(nThigh)=timeInd-lastChange;
                thigh_Ind(nThigh) = traceInd;
                
                nTmpTraceHigh=nTmpTraceHigh+1;
                
                tmpTraceHigh(nTmpTraceHigh)=timeInd-lastChange;
            end
           currentState=path{traceInd}(timeInd);
           lastChange=timeInd;
        end
    end
    %now adding this trace's statistics to tracesLow and tacesHigh
    tracesLow(traceInd,1)=mean(tmpTraceLow);
    tracesLow(traceInd,2)=std(tmpTraceLow)/sqrt(nTmpTraceLow);
    tracesLow(traceInd,3)=nTmpTraceLow;
    tracesHigh(traceInd,1)=mean(tmpTraceHigh);
    tracesHigh(traceInd,2)=std(tmpTraceHigh)/sqrt(nTmpTraceHigh);
    tracesHigh(traceInd,3)=nTmpTraceHigh;
    
    %now collecting statistics to find fraction high for this trace
    fractionHigh(traceInd)=double(totalHighPoints)/double(totalPoints);
    weightsFractionHigh(traceInd)=totalPoints;
    
    %if this was not a fluctuating trace add to 
    if nTlowFromThisTrace==0;
        numTracesNotFluct=numTracesNotFluct+1;
    end
end

[tmp,sortedInd]=sort(tracesLow(:,3));
tracesLow=tracesLow(sortedInd,:);
[tmp,sortedInd]=sort(tracesHigh(:,3));
tracesHigh=tracesHigh(sortedInd,:);

%calculating fraction of points that are high for all traces
fractionHighAll=sum(fractionHigh.*weightsFractionHigh)/sum(weightsFractionHigh);
errFractionHighAll=sum(...
    ((fractionHigh-fractionHighAll).*weightsFractionHigh).^2)...
    /sum(weightsFractionHigh.^2);
tlow = tlow * secPerFrame;
thigh = thigh * secPerFrame;
histTlow=hist(tlow,max(tlow)/secPerFrame);
histThigh=hist(thigh,max(thigh)/secPerFrame);

cumTlow=zeros(size(histTlow));
cumTlow(1)=histTlow(1);
indextrackerlow = [];
indextrackerlow(1) = 1;
uniquecumTlow =[];
uniquecumTlow(1) =histTlow(1);
uniquehistTlow = [];
uniquehistTlow(1) =histTlow(1);
for i=2:size(histTlow,2)
    if histTlow(i) ~=0
        indextrackerlow = [indextrackerlow , i];
        uniquecumTlow = [uniquecumTlow, cumTlow(i-1)+histTlow(i)];
        uniquehistTlow = [uniquehistTlow, histTlow(i)];
    end
    cumTlow(i)=cumTlow(i-1)+histTlow(i);
end
%and again for high
cumThigh=zeros(size(histThigh));
cumThigh(1)=histThigh(1);
indextrackerhigh = [];
indextrackerhigh (1) = 1;
uniquecumThigh =[];
uniquecumThigh(1)=histThigh(1);
uniquehistThigh = [];
uniquehistThigh(1) = histThigh(1);
for i=2:size(histThigh,2)
    if histThigh(i) ~=0
        indextrackerhigh = [indextrackerhigh , i];
        uniquecumThigh = [uniquecumThigh, cumThigh(i-1)+histThigh(i)];
        uniquehistThigh = [uniquehistThigh, histThigh(i)];
    end
    cumThigh(i)=cumThigh(i-1)+histThigh(i);
end