optimizeChunks.m

function [bestx, bestf] = optimizeChunks(tiepts, chunkStart, chunkEnd)

% Winnow chunk points table to only include the specified range

chunks = chunkStart:chunkEnd;
tiepts(~ismember(tiepts.chunk, chunks), :) = [];

% Calculate a time range and number of chunks in which each point appears

upoints     = unique(tiepts.point);
timeRange   = zeros(size(upoints)); 
pointFreq   = zeros(size(upoints));
for i = 1:numel(upoints)
    timeRange(i) = range(tiepts.flightTime(tiepts.point==upoints(i)));
    pointFreq(i) = numel(unique(tiepts.chunk(tiepts.point==upoints(i))));    
end

clear i

% Eliminate all points appearing in only one chunk

idx             = pointFreq < 2;
tiepts(ismember(tiepts.point, upoints(idx)), :) = [];
upoints(idx)    = [];
timeRange(idx)  = [];
pointFreq(idx)  = [];

% For each chunk, select (up to) five points with the greatest time range

selection   = false(size(upoints));

[~, iTime]  = sort(timeRange, 'descend');
upoints     = upoints(iTime);
pointFreq   = pointFreq(iTime);

for i = 1:numel(chunks)
    chunkidx    = tiepts.chunk == chunks(i);
    pointidx    = find(ismember(upoints, tiepts.point(chunkidx)));
    for j = 1:min([numel(pointidx) 30])
        selection(pointidx(j)) = true;
    end
end
        
clear iTime i frameidx pointidx j

% Create a table containing only high points

hiPoints    = upoints(selection);
highpts     = tiepts(ismember(tiepts.point, hiPoints), :);

clear i dataFile data idx datapts 

% Prepare objective function to be passed to SCEUA

idxchunk    = false(size(highpts,1), numel(chunks));
idxpoint    = false(size(highpts,1), numel(hiPoints));
hiPtFreq = pointFreq(selection);

for j = 1:numel(chunks)
    idxchunk(:,j) = highpts.chunk == chunks(j); end
for j = 1:numel(hiPoints)
    idxpoint(:,j) = highpts.point == hiPoints(j); end

functn = @(nopt, b)sceuaObjective(highpts.rawIntensity, ...
         idxpoint, idxchunk, hiPtFreq, nopt, b);
     
clear j idximage idxpoint
     
% Set up and run shuffled complex algorithm

x0      = [ones(1, numel(chunks)) zeros(1, numel(chunks))];
lb      = [1 0.8*ones(1, numel(chunks)-1) 0 -6000*ones(1, numel(chunks)-1)];
ub      = [  1.2*ones(1, numel(chunks))      6000*ones(1, numel(chunks))];
ngs     = 4;
maxn    = 5e5;
kstop   = 20;
pcento  = 0.01;
peps    = 1e-3;
iseed   = 1;
iniflg  = 1;

[bestx, bestf] = sceua(x0, lb, ub, maxn, kstop, pcento, peps, ngs, ...
                 iseed, iniflg, functn);
             
end