Subversion Repositories seema-scanner

Rev

Rev 219 | Rev 237 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 219 Rev 236
Line 34... Line 34...
34 
% 3D coordinates of markers in global initial alignment
 34 
% 3D coordinates of markers in global initial alignment
35 
Eg = cell(size(E));
 35 
Eg = cell(size(E));
36 
 
 36 
 
37 
% find 3D markers coordinates
37 
% find 3D markers coordinates
38 
for i=1:nSubScans
 38 
for i=1:nSubScans
39 
 
 39 
%for i=5:5
40 
    % load point cloud
 40 
    % load point cloud
41 
    pc = pcread(initialAlign(i).FileName);
 41 
    pc = pcread(initialAlign(i).FileName);
42 
    Q = pc.Location;
 42 
    Q = pc.Location;
43 
    idString = strsplit(initialAlign(i).FileName, {'.ply', '_'});
 43 
    idString = strsplit(initialAlign(i).FileName, {'.ply', '_'});
44 
    idString = idString{end-1};
 44 
    idString = idString{end-1};
Line 53... Line 53...
53 
    e1Coords = autoDetectMarkers(frame1, P1, Q);
 53 
    e1Coords = autoDetectMarkers(frame1, P1, Q);
54
54
55 
    %e0Coords = manuallyDetectMarkers(frame0, P0, Q);
 55 
    %e0Coords = manuallyDetectMarkers(frame0, P0, Q);
56 
    %e1Coords = manuallyDetectMarkers(frame1, P1, Q);
 56 
    %e1Coords = manuallyDetectMarkers(frame1, P1, Q);
57
57
58 
    %e0Coords = undistortPoints(e0Coords, camParams0);
 58 
    e0Coords = undistortPoints(e0Coords, camParams0);
59 
    %e1Coords = undistortPoints(e1Coords, camParams1);
 59 
    e1Coords = undistortPoints(e1Coords, camParams1);
60 
 
 60 
 
61 
    % match ellipse candidates between cameras based on projection
 61 
    % match ellipse candidates between cameras based on projection
62 
    E0 = projectOntoPointCloud(e0Coords, P0, Q);
 62 
    E0 = projectOntoPointCloud(e0Coords, P0, Q);
63 
    E1 = projectOntoPointCloud(e1Coords, P1, Q);
 63 
    E1 = projectOntoPointCloud(e1Coords, P1, Q);
64 
 
 64 
 
Line 182... Line 182...
182 
end
 182 
end
183 
 
 183 
 
184 
function e = autoDetectMarkers(frame, P, pointCloud)
 184 
function e = autoDetectMarkers(frame, P, pointCloud)
185 
 
 185 
 
186 
    % create mask based on morphology
 186 
    % create mask based on morphology
187 
    bw = imbinarize(rgb2gray(frame));
 187 
    g = rgb2gray(frame);
- 
 
 188 
    g(g>254) = 0;
- 
 
 189 
    bw = imbinarize(g, 'adaptive', 'Sensitivity', 10^(-50));
188 
    cc = bwconncomp(bw);
 190 
    cc = bwconncomp(bw);
189 
    rp = regionprops(cc, 'Area', 'Solidity', 'Eccentricity', 'Centroid');
 191 
    rp = regionprops(cc, 'Area', 'Solidity', 'Eccentricity', 'Centroid');
190 
    idx = ([rp.Area] > 100 & [rp.Area] < 1000 & [rp.Solidity] > 0.9);
 192 
    idx = ([rp.Area] > 100 & [rp.Area] < 1000 & [rp.Solidity] > 0.9);
191
193
192 
    initialGuesses = cat(1, rp(idx).Centroid);
 194 
    initialGuesses = cat(1, rp(idx).Centroid);
Line 232... Line 234...
232 
end
 234 
end
233 
 
 235 
 
234 
function [e, conf] = detectMarkersSubpix(frame, initGuesses, P, Q)
 236 
function [e, conf] = detectMarkersSubpix(frame, initGuesses, P, Q)
235 
 
 237 
 
236 
    % create mask based on morphology
 238 
    % create mask based on morphology
237 
    bw = imbinarize(rgb2gray(frame));
 239 
    g = rgb2gray(frame);
- 
 
 240 
    g(g>254) = 0;
- 
 
 241 
    bw = imbinarize(g);
238 
    cc = bwconncomp(bw);
 242 
    cc = bwconncomp(bw);
239 
    labels = labelmatrix(cc);
 243 
    labels = labelmatrix(cc);
240 
 
 244 
 
241 
    % project point cloud into image
 245 
    % project point cloud into image
242 
    q = [Q ones(size(Q,1),1)]*P;
 246 
    q = [Q ones(size(Q,1),1)]*P;
Line 283... Line 287...
283 
        [x, y] = meshgrid(1:151, 1:151);
 287 
        [x, y] = meshgrid(1:151, 1:151);
284 
        m((x(:)-76).^2 + (y(:)-76).^2 <= 50^2) = 1.0;
 288 
        m((x(:)-76).^2 + (y(:)-76).^2 <= 50^2) = 1.0;
285
289
286 
        % relation between marker space (px) and true marker/local plane(mm)
 290 
        % relation between marker space (px) and true marker/local plane(mm)
287 
        % true marker diameter is 1.75mm
 291 
        % true marker diameter is 1.75mm
288 
        mScale = 101/1.8; %px/mm
 292 
        mScale = 101/1.4; %px/mm
289 
        mShift = 76; %px
 293 
        mShift = 76; %px
290
294
291 
        % build homography from image to marker space
 295 
        % build homography from image to marker space
292 
        %H = fitgeotrans(q(pointMask, 1:2), mScale*Qlocal+mShift,  'projective');
 296 
        H = fitgeotrans(q(pointMask, 1:2), mScale*Qlocal+mShift,  'projective');
293 
        Hdlt = Hest_DLT([mScale*Qlocal+mShift, ones(size(Qlocal, 1), 1)]', q(pointMask,:)');
 297 
        %Hdlt = Hest_DLT([mScale*Qlocal+mShift, ones(size(Qlocal, 1), 1)]', q(pointMask,:)');
294 
        H = projective2d(Hdlt');
 298 
        %H = projective2d(Hdlt');
295
299
296 
        % bring image of marker into marker space
 300 
        % bring image of marker into marker space
297 
        imMarkerSpace = imwarp(frame, H, 'OutputView', imref2d(size(m)));
 301 
        imMarkerSpace = imwarp(frame, H, 'OutputView', imref2d(size(m)));
298 
        imMarkerSpace = rgb2gray(im2double(imMarkerSpace));
 302 
        imMarkerSpace = rgb2gray(im2double(imMarkerSpace));
299
303
Line 307... Line 311...
307 
        rege = imwarp(m, T, 'OutputView', imref2d(size(m)));
 311 
        rege = imwarp(m, T, 'OutputView', imref2d(size(m)));
308 
        %figure; imshowpair(imMarkerSpace, rege);
 312 
        %figure; imshowpair(imMarkerSpace, rege);
309
313
310 
        % measure of correlation
 314 
        % measure of correlation
311 
        confI = sum(sum(imMarkerSpace .* rege))/sqrt(sum(sum(imMarkerSpace) * sum(sum(rege))));
 315 
        confI = sum(sum(imMarkerSpace .* rege))/sqrt(sum(sum(imMarkerSpace) * sum(sum(rege))));
312 
        confI = 1.0;
 316 
        %confI = 1.0;
313
317
314 
        if confI<0.4
 318 
        if confI<0.4
315 
            continue;
 319 
            continue;
316 
        end
 320 
        end
317
321