% ***DATA FOR ANALYSIS*** %
w = wavread( 'A4_R2.wav' );
w1 = w(1:2050);

chunkSize = 16384;

% ***NUMBER OF DATA POINTS REQUIRED*** %
outputLength   = length(w1)+chunkSize-1;
pointsComputed = 0;

% ***FILTERS*** %
strtFreq = 1;
stopFreq = 5400;
freqStep = 1;
numFilters = (stopFreq-strtFreq)/freqStep;
filterBank = bpf(strtFreq,stopFreq,freqStep,chunkSize);

% ***VARIABLES FOR ANALYSIS***%
forThis( chunkSize, numFilters ) = 0;
forNext( chunkSize, numFilters ) = 0;
a_out(outputLength, numFilters)  = 0;

% ***DESA VARIABLES*** %
d_out(outputLength, numFilters)   = 0;
toDesa(chunkSize+2, numFilters)   = 0;
thisDesaA(chunkSize, numFilters)  = 0;
thisDesaF(chunkSize, numFilters)  = 0;
twoSamps( 2, numFilters ) = 0;

% ***NUMBER OF FASTCONVOLUTIONS WE NEED TO GET THERE*** %
numChunks = ceil(max(size(w1))/chunkSize);

% ***PEAK PICKING COMPONENTS*** %
thresh = 0.005;
segment = [];
clear criticalInfo;
clear hotSpots;
criticalInfo = {};
hotSpots     = {};
oldClusters = 1;
newClusters = 1;

for i=1:numChunks
    start = (i-1)*chunkSize + 1;
    fin   = chunkSize*i;
    try
        inpData = w1( start:fin );
    catch
        inpData = w1( start:end );
    end
    
    % ***FILTER BAND ANALYSIS
    for k=1:numFilters
        disp( sprintf('convolution @ freq %d', strtFreq+freqStep*k) );
        % ***DO A SINGLE FAST OVERLAP*** %
        outData = fast_overlap( inpData, forNext(:,k), filterBank(:,k) );
        
        % ***TAKE THE NEW FINISHED DATA OUT*** %
        forThis(:,k) = outData( 1:chunkSize );
        
        % ***ADD TO DESA BLOCK*** %
        toDesa(3:end,k) = forThis(:,k);
        
        % ***TAKE THE N-POINT DESA*** %    
        [thisDesaA(:,k), thisDesaF(:,k)] = desa( toDesa(:,k), twoSamps(:,k) );
        
        % ***SET THE TWO SAMPLES*** %
        twoSamps(:,k) = forThis( end-3:end-2,k );
        
        % ***GET DESA DATA READY*** %
        toDesa(1:2,k) = forThis( end-1:end,k );
        
        % ***GET DATA READY FOR NEXT OVERLAP*** %
        forNext(:,k) = outData( chunkSize+1:end );
        
        a_out(start:fin,k) = forThis(:,k);
        d_out(start:fin,k) = thisDesaA(:,k);
        % *** I should note that there is a two sample
        % *** latency present in the DESA output.
    end

    fileNumStr = [48*ones( 1,4-length( num2str(i) ) ), num2str(i) ];
    fileName   = ['data/conv_band_' fileNumStr ];
    save( fileName, 'd_out' );
    
    pointsComputed = chunkSize*i;
    start = pointsComputed+1;
    
    disp('convolution analysis complete');
    clear criticalInfo;
    
    for l=1:chunkSize
        thisSampleNo = chunkSize*(i-1)+l;
        sampsAbove = find( thisDesaA(l,:)>thresh );
        sampsAbove = transpose( sampsAbove );
        cluster = [];
        if sampsAbove
            diffAbove   = diff( sampsAbove );
            if diffAbove
                transitions = find( diffAbove ~= 1 );
                startIndex  = 1;
                if transitions
                    for k=1:length(transitions)
                        segment = sampsAbove( [(startIndex) transitions(k)] );
                        cluster( :,size( cluster, 2 ) + 1 ) = segment;
                        startIndex = transitions(k) + 1;
                    end
                end
                segment = sampsAbove( [(startIndex) length(sampsAbove) ] );
                cluster( :,size( cluster, 2 ) + 1 ) = segment;
            end
  
        end
        criticalInfo{l} = struct( 'SampleNo', thisSampleNo, 'Clusters', cluster );
        
        numClusters = size( criticalInfo{l}.Clusters,2 );
        for n=numClusters:-1:1
            oldClusters = size( criticalInfo{l}.Clusters,2 );
            for m=1:size( hotSpots, 2 )
                if ( find( abs(hotSpots{m}.SampBand-thisSampleNo)==1 ) )
                    if is_overlap( hotSpots{m}.IndexBand, criticalInfo{l}.Clusters(:,n) )
                        hotSpots{m}.IndexBand(1) = ...
                            min( hotSpots{m}.IndexBand(1), criticalInfo{l}.Clusters(1,n) );
                        hotSpots{m}.IndexBand(2) = ...
                            max( hotSpots{m}.IndexBand(2), criticalInfo{l}.Clusters(2,n) );
                        hotSpots{m}.SampBand(2,1)= criticalInfo{l}.SampleNo;
                        criticalInfo{l}.Clusters(:,n) = [];
                        break;
                    end
                end
            end
            newClusters = size( criticalInfo{l}.Clusters,2 );
            if oldClusters==newClusters
                hotSpots{ length(hotSpots) + 1 } = ...
                    struct( 'SampBand',thisSampleNo,'IndexBand',criticalInfo{l}.Clusters(:,n) );
                criticalInfo{l}.Clusters(:,n) = [];
            end    
        end
    end
    % *** okay, so we've done all the analysis and saved all of the data.
    % *** it's now time to gather out the really significant stuff, i.e.
    % *** the clusters.
    for p=1:size( hotSpots, 2 )
    end
end

%if pointsComputed<outputLength
%    forThis = forNext( 1:outputLength-pointsComputed );
%    a_out = [a_out; forThis];
%end


