# FindPelican

## PURPOSE

FindPelican - A simple test case of detection using shape script.

## SYNOPSIS

This is a script file.

## DESCRIPTION

     FindPelican - A simple test case of detection using shape script.

       Detect pelicans with a simple egret/pelican shape script template.

       The shape script template is found at "result/egret/egret_template.eps" or
     "result/egret/template.png". The detection results together with S3 and S2 maps
     are found at "result/pelican".

## CROSS-REFERENCE INFORMATION

This function calls:

- GaborFilter GaborFilter - generate Gabor filter at fixed scale and orientation.
- applyfilterMultiScale applyfilterMultiScale - Filter image by a bank of multiscale filters, with local normalization.

This function is called by:

- StartFromHere StartFromHere - The entry point of this project.

## SOURCE CODE

    0001 % FindPelican - A simple test case of detection using shape script.
    0002 %
    0003 %   Detect pelicans with a simple egret/pelican shape script template.
    0004 %
    0005 %   The shape script template is found at "result/egret/egret_template.eps" or
    0006 % "result/egret/template.png". The detection results together with S3 and S2 maps
    0007 % are found at "result/pelican".
    0008
    0009 clear
    0010 close all
    0011 addpath('common');
    0012
    0013 %% load the S2 (shape motif) dictionary and S3 (shape script) template
    0014 load template_params
    0015 load result/transformedS2Template.mat
    0016 load result/egret/shapescript.mat
    0017 category = 'pelican';
    0018
    0019 %% load the testing image
    0020 Iname = {'pelican 2_big.jpg','pelican3.jpg',...
    0021     'pelican4.jpg','pelican5.jpg','pelican6.jpg',...
    0022     'pelican7.jpg','pelican8.jpg','pelican9.jpg',...
    0023     'pelican10.jpg'};
    0024 imScale = [0.25 0.4 0.4 0.4 ...
    0025     0.22 0.36 0.5 0.4 0.3];
    0026
    0027 subsampleM1 = 1;
    0028 for iImage = 1:length(Iname)
    0029     im = imread( sprintf('image/%s/%s',category,Iname{iImage}) );
    0030     im = imresize( im,imScale(iImage),'nearest' );
    0031
    0032     if size(im,3) == 3
    0033         im = rgb2gray(im);
    0034     end
    0035
    0036     im = single(im);
    0037     S1Map = applyfilterMultiScale({im},GaborFilter, nGaborOri);
    0038
    0039     % sigmoid transform
    0040     nImage = 1;
    0041     heightS1Map = zeros(length(S1Map),1,'single');
    0042     widthS1Map = zeros(length(S1Map),1,'single');
    0043     for i = 1:length(S1Map)
    0044         heightS1Map(i) = size(S1Map{i},1);
    0045         widthS1Map(i) = size(S1Map{i},2);
    0046     end
    0047     tmpNGaborOri = 1;
    0048     for i = 1:length(S1Map)
    0049         Csigmoid(nImage, heightS1Map(i), widthS1Map(i), tmpNGaborOri, GaborSaturation, S1Map(i));
    0050     end
    0051
    0052     % compute M1 map
    0053     sizeGaborPerScale = zeros(nGaborScale,1,'single');
    0054     for iS = 1:nGaborScale
    0055         sizeGaborPerScale(iS) = size( GaborFilter{(iS-1)*nGaborOri+1}, 1 );
    0056     end
    0057     locationShift = .2;
    0058     oriShift = 1;
    0059     [M1Map M1Trace M1RowShift M1ColShift M1OriShift] = ...
    0060         CcomputeMAX1( nGaborOri, S1Map, locationShift, oriShift, sizeGaborPerScale, subsampleM1 );
    0061
    0062     % compute S2 map
    0063     S2Map = CcomputeSUM2( nGaborOri, M1Map, TransformedTemplate, subsampleS2 );
    0064
    0065     % compute M2 map
    0066     locationPerturbationFraction = 0.2; subsampleM2 = 1;
    0067     M2Map = cell(nTransform, nS2);
    0068     M2LocationTrace = cell(nTransform, nS2);
    0069     M2TransformTrace = cell(nTransform, nS2);
    0070     for iS2 = 1:nS2
    0071         sizeTemplate = zeros(nTransform,1);
    0072         for iT = 1:nTransform
    0073             sizeTemplate(iT) = sqrt( ...
    0074                 size(TransformedTemplate{iT,iS2}.sym,1)...
    0075                 * size(TransformedTemplate{iT,iS2}.sym,2) ) / subsampleS2;
    0076         end
    0077         [M2Map(:,iS2) M2LocationTrace(:,iS2) M2TransformTrace(:,iS2) storedM2LocationShift] = ...
    0078               CcomputeMAX2( templateAffMat, S2Map((iS2-1)*nTransform + (1:nTransform)), locationPerturbationFraction, ...
    0079                  int32(sizeTemplate), subsampleM2 );
    0080         storedSizeS2Template{iS2} = sizeTemplate;
    0081     end
    0082     disp('done with computing M2');
    0083
    0084     %% detection
    0085     subsampleS3 = 1;
    0086     %
    0087     % Compute SUM 3 map :
    0088     %
    0089     S3Map = CcomputeSUM3_loglinear( M2Map, S3Template, subsampleS3, nTransform );
    0090     selectedTransform = S3Template{1}.selectedTransform;
    0091     selectedRow = S3Template{1}.selectedRow;
    0092     selectedCol = S3Template{1}.selectedCol;
    0093     selectedS2Ind = S3Template{1}.selectedInd;
    0094     selectedLambda = S3Template{1}.selectedLambda;
    0095     selectedLogZ = S3Template{1}.selectedLogZ;
    0096     nF = length(selectedTransform);
    0097     %
    0098     %
    0099     %
    0100     storedS3Map = S3Map;
    0101     minval = min(S3Map{1}(:));
    0102
    0103     maxval = 10;
    0104     count = 0;
    0105     maxDetection = 1; % maximum number of detections, currently we only detect the maximum one
    0106
    0107
    0108     while maxval > -10
    0109         count = count + 1;
    0110         if count > maxDetection
    0111             break;
    0112         end
    0113
    0114         % find the best position
    0115         [maxval ind] = max(S3Map{1}(:));
    0116         col = ceil(ind/size(S3Map{1},1));
    0117         row = ind - (col-1) * size(S3Map{1},1);
    0118
    0119         overlay = single(im);
    0120         mask = zeros(size(im));
    0121         for iF = 1:nF
    0122             iT = selectedTransform(iF)+1;
    0123             iS2 = selectedS2Ind(iF) + 1;
    0124
    0125             part_row = row + selectedRow(iF);
    0126             part_col = col + selectedCol(iF);
    0127             if part_row <= 0 || part_col <= 0 || part_row > size(S3Map{1},1) || part_col > size(S3Map{1},2)
    0128                 continue;
    0129             end
    0130
    0131             % trace the deformation for each part
    0132             part_loc_trace = M2LocationTrace{iT,iS2}(part_row,part_col) + 1;
    0133             part_trans_trace = M2TransformTrace{iT,iS2}(part_row,part_col) + 1; % in matlab, the index start from 1
    0134             S2Size = storedSizeS2Template{iS2}(part_trans_trace); % size of the S2 template
    0135             part_row = part_row * subsampleS2 + floor( .5 + storedM2LocationShift(part_loc_trace,1) * S2Size ) * subsampleS2;
    0136             part_col = part_col * subsampleS2 + floor( .5 + storedM2LocationShift(part_loc_trace,2) * S2Size ) * subsampleS2;
    0137
    0138             S2T = TransformedTemplate{part_trans_trace,iS2};
    0139
    0140             for iGabor = 1:length(S2T.selectedRow)
    0141                 gabor_row = part_row + S2T.selectedRow(iGabor);
    0142                 gabor_col = part_col + S2T.selectedCol(iGabor);
    0143                 gabor_ori = S2T.selectedOri(iGabor) + 1;
    0144                 gabor_scale = S2T.selectedScale(iGabor) + 1;
    0145                 if gabor_row <= 0 || gabor_col <= 0 || gabor_row > size(im,1) || gabor_col > size(im,2)
    0146                     continue;
    0147                 end
    0148                 gabor_trace = M1Trace{(gabor_scale-1)*nGaborOri+gabor_ori}(gabor_row,gabor_col) + 1;
    0149                 gabor_dr = single(M1RowShift{gabor_scale}(gabor_trace,gabor_ori));
    0150                 gabor_dc = single(M1ColShift{gabor_scale}(gabor_trace,gabor_ori));
    0151                 gabor_do = single(M1OriShift{gabor_scale}(gabor_trace,gabor_ori));
    0152                 gabor_row = gabor_row + gabor_dr;
    0153                 gabor_col = gabor_col + gabor_dc;
    0154                 gabor_ori = gabor_do + 1; % note: gabor_do is already the transformed orientation
    0155                 while gabor_ori > nGaborOri
    0156                     gabor_ori = gabor_ori - nGaborOri;
    0157                 end
    0158                 while gabor_ori <= 0
    0159                     gabor_ori = gabor_ori + nGaborOri;
    0160                 end
    0161
    0162                 % strength of this Gabor symbol (between 0 and 1)
    0163                 S1score = S1Map{(gabor_scale-1)*nGaborOri+gabor_ori}(gabor_row,gabor_col);
    0164                 S1score = S1score - 3; % These two parameters are just for display purpose
    0165                 priorS1 = 0.5;
    0166                 strength = priorS1 * exp(S1score) / ( (1-priorS1) + priorS1 * exp(S1score) );
    0167                 strength = max(0.8,strength);
    0168
    0169                 % Gabor symbol
    0170                 patch = strength * GaborSymbol{(gabor_scale-1)*nGaborOri+gabor_ori};
    0171
    0172                 % below: the bounding box for Gabor inside the image
    0173                 top2 = gabor_row - floor(size(patch,1)/2);
    0174                 left2 = gabor_col - floor(size(patch,2)/2);
    0175                 for c = 1:size(patch,2)
    0176                     for r = 1:size(patch,1)
    0177                         c2 = left2 + c - 1;
    0178                         r2 = top2 + r - 1;
    0179                         if r2 <= 0 || r2 > size(mask,1) || c2 <= 0 || c2 > size(mask,2)
    0180                             continue;
    0181                         end
    0182                         val = patch(r,c);
    0183                         if val > mask(r2,c2)
    0184                             mask(r2,c2) = val;
    0185                         end
    0186                     end
    0187                 end
    0188             end
    0189
    0190         end
    0191
    0192         h = figure;
    0193         subplot(2,3,1);imshow(im+200*mask,[]);
    0194         subplot(2,3,2);imshow(S3Map{1},[]);
    0195         subplot(2,3,4);imshow(M2Map{selectedTransform(1)+1,selectedS2Ind(1)+1},[]);
    0196         subplot(2,3,5);imshow(M2Map{selectedTransform(2)+1,selectedS2Ind(2)+1},[]);
    0197         subplot(2,3,6);imshow(M2Map{selectedTransform(3)+1,selectedS2Ind(3)+1},[]);
    0198         mkdir(sprintf('result/%s',category));
    0199         saveas(h,sprintf('result/%s/detection%d_%s',category,count,Iname{iImage}));
    0200         save(sprintf('result/%s/detectionResult%d_%s.mat',category,count,Iname{iImage}),...
    0201             'im','mask','S3Map','M2Map', 'storedSizeS2Template');
    0202         close;
    0203
    0204         % inhibit the local region for the next detection
    0205         for r = row - floor(S3size/2/subsampleS2) : row + floor(S3size/2/subsampleS2)
    0206             for c = col - floor(S3size/2/subsampleS2) : col + floor(S3size/2/subsampleS2)
    0207                 if r <= 0 || c <= 0 || r > size(S3Map{1},1) || c > size(S3Map{1},2)
    0208                     continue;
    0209                 end
    0210                 S3Map{1}(r,c) = minval;
    0211             end
    0212         end
    0213     end
    0214 end