%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% test parameters for particle finding

% camera calibration model
calib_file				= 'F:\mixer\lpt\calibration\21-06-30\camera-model.sheet.mat';
% pattern of source images to load
src_dir					= 'F:\mixer\lpt\run29\';
src_file_pat			= '%s/frame_%04d.cam%d.tif';
src_file_loader			= @(fidx, cidx) load_tif(src_file_pat, src_dir, fidx, cidx-1);
src_range				= 0 : 179;
% pattern of background image to load
bg_file_pat				= '%s/background.cam%d.tif';
bg_file_fun				= @(cam) sprintf(bg_file_pat, src_dir, cam-1);
% camera to test
cam_idx					= 2;

% particle detection and matching options
mopts					= struct();
mopts.algorithm			= 'template';
if strcmp(mopts.algorithm, 'steerable')
	% steerable filter matching options
	mopts.steerable				= struct();
	mopts.steerable.r_filter	= 8;			% scale of ridge filter, px
	mopts.steerable.im_threshold= 0.5;			% threshold for particle detection
	mopts.steerable.b_overlap   = true;			% enable detection of overlapping particles?
	mopts.steerable.n_theta		= 16;			% angular resolution for multi-particle matching
	mopts.steerable.l_max_px	= 100;			% longest length of rod, px
elseif strcmp(mopts.algorithm, 'template')
	% template matching options
	mopts.template              = struct();
	mopts.template.kernel		= -fspecial('log', [5 5], 1);
	mopts.template.threshold    = .1;
	mopts.template.mask_fun		= @(img) img > 0.5;	% mask function
elseif strcmp(mopts.algorithm, 'circle')
	% hough transform circle matching options
	mopts.circle				= struct();
	mopts.circle.threshold		= 100;			% edge threshold for circle detection
	mopts.circle.radius			= [6 12];		% radius in px for circular search
	mopts.circle.polarity		= 'bright';       % polarity (for shadowgraph, is dark)
	mopts.circle.binning		= 1;			% factor by which to downsample image in search (helps accelerate search)
else
	error('unrecognised algorithm');
end

% build basic options structure
opts			= struct();
opts.matching	= mopts;
opts.frame_range= src_range;
opts.frame_rate = 1;

%% load camera calibration
fprintf('load camera calibration: %s\n', calib_file);
setup			= load(calib_file);
camera			= setup.camera;
i_axis			= camera(1).i_axis;
j_axis			= camera(1).j_axis;
n_i             = length(i_axis);
n_j             = length(j_axis);

%% define image preprocessing
fprintf('Loading background image ... ');
preproc			= cell(4, 3);
for c = 1 : 3
	bg_file			= bg_file_fun(c);
	im_dark			= 0*double(imread(bg_file,1))';
	im_white		= double(imread(bg_file,2))' - im_dark;	

	preproc{1,c}	= struct('type', 'ssbg', 'bg', im_dark);
	preproc{2,c}	= struct('type', 'levelize', 'white', im_white);
	preproc{3,c}	= struct('type', 'clip', 'threshold', [0 1]);
 	preproc{4,c}	= struct('type', 'invert', 'offset', 1);
	%preproc{4,c}	= struct('type', 'median', 'size', [3 3]);
end
topts.preprocessing	= preproc;
fprintf('done\n');

%% build tracker
tracker			= ShadowTrack(setup, opts);

%% create figure
figure(1);
clf;
h_im = imagesc(i_axis, j_axis, zeros(n_i,n_j)');

hold on;
h_ax = gca;
set(h_ax,'color','b');
axis equal tight xy;
colormap gray;

%% iterate over images
n_frames = length(src_range);
for f = 1 : n_frames
    %% load an image
    img				= src_file_loader(src_range(f), cam_idx);
    img_preproc		= image_preprocess(img, preproc(:,cam_idx));
    
    %% mask image
    img_mask		= mopts.template.mask_fun(img_preproc);
    
    %% find particles
    p2d				= tracker.find_particles(img_preproc);
    im_ctr			= [p2d.centroid];
    im_amp			= [p2d.intensity];
    
    %% display image
    h_im.CData      = img_preproc';
    h_im.AlphaData = ~img_mask';
    old = findobj(h_ax,'type','line','UserData',f-10);
    delete(old);
    set(findobj(h_ax,'type','line'), 'marker', '.', 'color', 'r')
    plot(h_ax, im_ctr(1,:), im_ctr(2,:), 'gx', 'UserData', f);
    pause();
end
    
    

