6D-Reconstruction-Twins: added masking out what is outside the parent convex shape for raw images

Signed-off-by: Nicola Vigano <nicola.vigano@esrf.fr>

zUtil_Deformation/gt6DCreateProjDataFromTwinnedGrainFwdProj.m

@@ -292,6 +292,9 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
             % Using the same because it gives the same size/uvw-limits
             if (ii_g > 1)
                 blobs(extended_projs(ii_g).shared_refl) = samp_ors(1).proj(conf.det_index).bl(extended_projs(ii_g).shared_refl_pos_in_parent);
+
+                proj_bl_masks = project_masks(samp_ors(ii_g), extended_projs(ii_g), grs(1).fwd_sim(conf.det_index).gv_verts, p, conf);
+                blobs = assign_masks(blobs, proj_bl_masks, ~extended_projs(ii_g).shared_refl);
             end
         else
             fprintf(' - Loading segmented blobs..')
@@ -689,4 +692,37 @@ function blob = load_blob(img_bboxes, img_sizes, stackUSize, stackVSize, p, conf
     blob.intensity = blob_int;
 end
 
+function proj_bl_masks = project_masks(refor, eproj, parent_verts, p, conf)
+    sel = eproj.ondet(eproj.included);
+    proj_bl_masks = gt6DSpreadProjectVertices2Det(refor, parent_verts, sel, p, conf.det_index);
+end
+
+function blobs = assign_masks(blobs, proj_mask_blobs, sel)
+    blobs_u_lims = cat(1, blobs(sel).bbuim);
+    blobs_v_lims = cat(1, blobs(sel).bbvim);
+
+    proj_masks_u_lims = cat(1, proj_mask_blobs(sel).bbuim);
+    proj_masks_v_lims = cat(1, proj_mask_blobs(sel).bbvim);
+
+    shifts = [ ...
+        (proj_masks_u_lims(:, 1) - blobs_u_lims(:, 1)), ...
+        (proj_masks_v_lims(:, 1) - blobs_v_lims(:, 1)) ];
+    shifts(:, 3) = 0;
+
+    inds = find(sel);
+    for ii = 1:numel(inds)
+        ii_b = inds(ii);
+
+        mask = uint8(proj_mask_blobs(ii_b).mask);
+        mask = mask(:, :, ones(blobs(ii_b).bbsize(3), 1));
+
+        blobs(ii_b).mask = gtPlaceSubVolume( ...
+            zeros(blobs(ii_b).bbsize, 'uint8'), mask, shifts(ii, :) );
+
+        blobs(ii_b).intm(~blobs(ii_b).mask) = 0;
+        blob_int = gtMathsSumNDVol(blobs(ii_b).intm);
+        blobs(ii_b).intensity = blob_int;
+    end
+end
+
 

zUtil_Deformation/gt6DSpreadProjectVertices2Det.m

+ function proj_bls = gt6DSpreadProjectVertices2Det(gr, verts_rec, ondet, parameters, det_ind, verbose)
+% FUNCTION proj_bls = gt6DSpreadProjectVertices2Det(gr, verts_rec, ondet, parameters, det_ind, verbose)
+%
+
+    if (~exist('verbose', 'var') || isempty(verbose))
+       verbose = false;
+    end
+
+    detgeo = parameters.detgeo(det_ind);
+    labgeo = parameters.labgeo;
+    samgeo = parameters.samgeo;
+    recgeo = parameters.recgeo(det_ind);
+
+    omstep = gtGetOmegaStepDeg(parameters, det_ind);
+
+    r_vecs = cat(1, gr.bb_ors(:).R_vector);
+    r_vecs = unique(r_vecs, 'rows');
+
+    if (size(r_vecs, 1) ~= 8)
+        error('gt6DSpreadProjectVertices2Det:wrong_argument', ...
+            'Incoming bb orientations have %d unique orietations (8 were expected)', ...
+            size(r_vecs, 1))
+    end
+
+    ors = cell(8, 1);
+    for ii_o = 1:8
+        ors{ii_o} = struct('id', gr.id, 'phaseid', gr.phaseid, ...
+            'center', gr.center, 'R_vector', r_vecs(ii_o, :));
+    end
+    ors = gtCalculateGrain_p(ors, parameters, 'ref_omind', gr.allblobs.omind);
+    ors = [ors{:}];
+
+    num_blobs = numel(ondet);
+    proj_bls = gtFwdSimBlobDefinition('blob', num_blobs);
+
+    dvecs_lab = zeros(num_blobs, 3, 8);
+    omegas = zeros(num_blobs, 8);
+    rot_l2s = zeros(3, 3, num_blobs, 8);
+
+    for ii_o = 1:8
+        dvecs_lab(:, :, ii_o) = ors(ii_o).allblobs.dvec(ondet, :);
+        omegas(:, ii_o) = ors(ii_o).allblobs.omega(ondet);
+        rot_l2s(:, :, :, ii_o) = ors(ii_o).allblobs.srot(:, :, ondet);
+    end
+    rot_s2l = permute(rot_l2s, [2 1 3 4]);
+
+    if (verbose)
+        uv = gr.allblobs.detector(det_ind).uvw(ondet, 1:2);
+    end
+
+    ones_verts = ones(size(verts_rec, 1), 1);
+    gc_sam = gr.center(ones_verts, :);
+
+    verts_sam = repmat(gtGeoSam2Sam(verts_rec, recgeo, samgeo, false) + gc_sam, [8, 1]);
+
+    % Now we project the vertices to the spots and take the convex hull
+    for ii_b = 1:num_blobs
+        rot_s2l_w = rot_s2l(:, :, ii_b .* ones_verts, :);
+        rot_s2l_w = reshape(rot_s2l_w, 3, 3, []);
+        verts_lab = gtGeoSam2Lab(verts_sam, rot_s2l_w, labgeo, samgeo, false);
+        dvec_vers = dvecs_lab(ii_b .* ones_verts, :, :);
+        dvec_vers = reshape(permute(dvec_vers, [1 3 2]), [], 3);
+
+        verts_det_pos = gtFedPredictUVMultiple([], dvec_vers', verts_lab', ...
+            detgeo.detrefpos', detgeo.detnorm', detgeo.Qdet, ...
+            [detgeo.detrefu, detgeo.detrefv]')';
+
+        if (verbose)
+            f = figure();
+            ax = axes('parent', f);
+            hold(ax, 'on')
+            scatter(ax, verts_det_pos(:, 1), verts_det_pos(:, 2));
+            scatter(ax, uv(ii_b, 1), uv(ii_b, 2), 30, 'r');
+            hold(ax, 'off')
+        end
+
+        verts_idxs = convhull(verts_det_pos);
+        verts_det_pos = verts_det_pos(verts_idxs, :);
+
+        bb_2d = [floor(min(verts_det_pos, [], 1)), ceil(max(verts_det_pos, [], 1))];
+        bb_2d = [bb_2d(1:2), (bb_2d(3:4) - bb_2d(1:2) + 1)];
+
+        mask = zeros(bb_2d(3:4));
+        verts_mask_pos = round(verts_det_pos) - bb_2d(ones(numel(verts_idxs), 1), 1:2) + 1;
+        verts_mask_indx = verts_mask_pos(:, 1) + bb_2d(3) * (verts_mask_pos(:, 2) - 1);
+        mask(verts_mask_indx) = 1;
+        mask = bwconvhull(mask);
+
+        proj_bls(ii_b).bbuim = [bb_2d(1), (bb_2d(1) + bb_2d(3) - 1)];
+        proj_bls(ii_b).bbvim = [bb_2d(2), (bb_2d(2) + bb_2d(4) - 1)];
+        proj_bls(ii_b).bbwim = [ ...
+            floor(min(omegas(ii_b, :)) / omstep), ...
+            ceil(max(omegas(ii_b, :)) / omstep) ];
+        proj_bls(ii_b).bbsize = [bb_2d(3:4), 1];
+        proj_bls(ii_b).mask = mask;
+
+        if (verbose)
+            f = figure();
+            ax = axes('parent', f);
+            imshow(mask, 'parent', ax)
+
+            pause
+        end
+    end
+ end