Skip to content
Snippets Groups Projects
Commit 81452f18 authored by Nicola Vigano's avatar Nicola Vigano
Browse files

6D-Twin-reconstruction: fixed selection of twin blobs from parent's shared...

6D-Twin-reconstruction: fixed selection of twin blobs from parent's shared selection, and matching of the two

Signed-off-by: default avatarNicola Vigano <nicola.vigano@esrf.fr>
parent 18a290d3
No related branches found
No related tags found
No related merge requests found
...@@ -81,8 +81,8 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -81,8 +81,8 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% ones that were selected among the included, because the % ones that were selected among the included, because the
% information we have in the twins about the spots that % information we have in the twins about the spots that
% coincide are about the included and not the selected % coincide are about the included and not the selected
selected_pos_in_included = zeros(size(ref_inc)); ref_inc_to_sel = zeros(size(ref_inc));
selected_pos_in_included(ref_sel) = 1:numel(find(ref_sel)); ref_inc_to_sel(ref_sel) = 1:numel(find(ref_sel));
blobs_w_interp = cell(num_regions, 1); blobs_w_interp = cell(num_regions, 1);
...@@ -97,7 +97,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -97,7 +97,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% the shared reflections % the shared reflections
if (ii_o_reg > 1) if (ii_o_reg > 1)
[shared, parent_pos_shared_bls] = ... [shared, parent_pos_shared_bls] = ...
self.get_shared_bls(selected_pos_in_included, sampler, ii_o_reg); self.get_shared_bls(ref_inc_to_sel, sampler, ii_o_reg);
blobs_w_interp{ii_o_reg}(shared) = blobs_w_interp{1}(parent_pos_shared_bls); blobs_w_interp{ii_o_reg}(shared) = blobs_w_interp{1}(parent_pos_shared_bls);
...@@ -141,21 +141,21 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -141,21 +141,21 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% We should first go through all the shared blobs in each % We should first go through all the shared blobs in each
% orientation and extend them if needed, and then extend the % orientation and extend them if needed, and then extend the
% positions pointed by the orientations in the shared blobs % positions pointed by the orientations in the shared blobs
all_shared_blobs = find(ref_gr.proj(self.det_index).shared_bls_twins(ref_sel)); all_shared_blobs_sel = find(ref_gr.proj(self.det_index).shared_bls_with_twins(ref_sel));
% The reference num_interps % The reference num_interps
shared_blobs_w_interp = zeros(num_parent_blobs, 1); shared_blobs_w_interp = zeros(num_parent_blobs, 1);
shared_blobs_w_interp(all_shared_blobs) = blobs_w_interp{1}(all_shared_blobs); shared_blobs_w_interp(all_shared_blobs_sel) = blobs_w_interp{1}(all_shared_blobs_sel);
% A record for the W limits of shared blobs between parent and % A record for the W limits of shared blobs between parent and
% twins % twins
shared_bls_w_lims = zeros(num_parent_blobs, 2, num_regions); shared_bls_w_lims = zeros(num_parent_blobs, 2, num_regions);
shared_bls_w_lims(all_shared_blobs, :, 1) = algo_params(1).blobs_w_lims(all_shared_blobs, :) ./ shared_blobs_w_interp(all_shared_blobs, [1 1]); shared_bls_w_lims(all_shared_blobs_sel, :, 1) = algo_params(1).blobs_w_lims(all_shared_blobs_sel, :) ./ shared_blobs_w_interp(all_shared_blobs_sel, [1 1]);
rec_shared_bls_w_lims = shared_bls_w_lims(:, :, 1); rec_shared_bls_w_lims = shared_bls_w_lims(:, :, 1);
for ii_o_reg = 2:num_regions for ii_o_reg = 2:num_regions
[shared, parent_pos_shared_bls] = ... [shared, parent_pos_shared_bls] = ...
self.get_shared_bls(selected_pos_in_included, sampler, ii_o_reg); self.get_shared_bls(ref_inc_to_sel, sampler, ii_o_reg);
shared_bls_w_lims(parent_pos_shared_bls, :, ii_o_reg) = algo_params(ii_o_reg).blobs_w_lims(shared, :) ... shared_bls_w_lims(parent_pos_shared_bls, :, ii_o_reg) = algo_params(ii_o_reg).blobs_w_lims(shared, :) ...
./ shared_blobs_w_interp(parent_pos_shared_bls, [1 1]); ./ shared_blobs_w_interp(parent_pos_shared_bls, [1 1]);
...@@ -165,35 +165,52 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -165,35 +165,52 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
max(rec_shared_bls_w_lims(parent_pos_shared_bls, 2), shared_bls_w_lims(parent_pos_shared_bls, 2, ii_o_reg)) ]; max(rec_shared_bls_w_lims(parent_pos_shared_bls, 2), shared_bls_w_lims(parent_pos_shared_bls, 2, ii_o_reg)) ];
end end
to_be_changed = any(rec_shared_bls_w_lims(all_shared_blobs, :) ~= shared_bls_w_lims(all_shared_blobs, :, 1), 2); to_be_changed = any(rec_shared_bls_w_lims(all_shared_blobs_sel, :) ~= shared_bls_w_lims(all_shared_blobs_sel, :, 1), 2);
rec_rescaled_shared_blobs_sizes = rec_shared_bls_w_lims(all_shared_blobs, 2) - rec_shared_bls_w_lims(all_shared_blobs, 1) + 1; rec_rescaled_shared_blobs_sizes = rec_shared_bls_w_lims(all_shared_blobs_sel, 2) - rec_shared_bls_w_lims(all_shared_blobs_sel, 1) + 1;
fprintf('Shared blobs:\n')
% And now let's extend those blobs that need it! % And now let's extend those blobs that need it!
for ii_b = 1:numel(all_shared_blobs) for ii_b = 1:numel(all_shared_blobs_sel)
shared_blob_pos_inc = all_shared_blobs_sel(ii_b);
% Collecting IDs of other regions using the said blobs
reg_ids = zeros(num_regions, 1);
for ii_o_reg = 2:num_regions
[~, parent_pos_shared_bls] = ...
self.get_shared_bls(ref_inc_to_sel, sampler, ii_o_reg);
if (~isempty(find(parent_pos_shared_bls == shared_blob_pos_inc, 1)))
curr_or = sampler(ii_o_reg).get_reference_grain();
reg_ids(ii_o_reg) = curr_or.id;
end
end
reg_ids = reg_ids(reg_ids > 0);
twins_string = sprintf(' %d', reg_ids);
fprintf('%02d) Blob in parent: %d, shared with twins:%s\n', ii_b, shared_blob_pos_inc, twins_string)
% if they changed, let's enlarge them % if they changed, let's enlarge them
if (to_be_changed(ii_b)) if (to_be_changed(ii_b))
shared_blob_pos = all_shared_blobs(ii_b); proj_shift = shared_bls_w_lims(shared_blob_pos_inc, 1, 1) ...
- rec_shared_bls_w_lims(shared_blob_pos_inc, 1);
proj_shift = shared_bls_w_lims(shared_blob_pos, 1, 1) ... curr_blob_size = size(blobs{shared_blob_pos_inc});
- rec_shared_bls_w_lims(shared_blob_pos, 1);
curr_blob_size = size(blobs{shared_blob_pos});
new_size = [curr_blob_size(1), rec_rescaled_shared_blobs_sizes(ii_b), curr_blob_size(3)]; new_size = [curr_blob_size(1), rec_rescaled_shared_blobs_sizes(ii_b), curr_blob_size(3)];
fprintf('* Enlarging blob: %d (shared %d), depth: %d -> %d, with shift: %d\n', ... fprintf(' * Enlarging blob, depth: %d -> %d, with shift: %d\n', ...
shared_blob_pos, ii_b, curr_blob_size(2), new_size(2), proj_shift) curr_blob_size(2), new_size(2), proj_shift)
% Creating the new blob % Creating the new blob
new_blob = zeros(new_size, 'like', blobs{shared_blob_pos}); new_blob = zeros(new_size, 'like', blobs{shared_blob_pos_inc});
new_blob = gtPlaceSubVolume(new_blob, blobs{shared_blob_pos}, [0, proj_shift, 0]); new_blob = gtPlaceSubVolume(new_blob, blobs{shared_blob_pos_inc}, [0, proj_shift, 0]);
blobs{shared_blob_pos} = new_blob; blobs{shared_blob_pos_inc} = new_blob;
% Fixing the shift in parent's projection % Fixing the shift in parent's projection
if (~isempty(geometries_ss)) if (~isempty(geometries_ss))
for ii_o = 1:numel(offsets_ss) for ii_o = 1:numel(offsets_ss)
for ii_ss = 1:numel(offsets_ss{ii_o}) for ii_ss = 1:numel(offsets_ss{ii_o})
inds_of_shared_blobs = offsets_ss{ii_o}{ii_ss}.blob_offsets == shared_blob_pos; inds_of_shared_blobs = offsets_ss{ii_o}{ii_ss}.blob_offsets == shared_blob_pos_inc;
offsets_ss{ii_o}{ii_ss}.proj_offsets(inds_of_shared_blobs) ... offsets_ss{ii_o}{ii_ss}.proj_offsets(inds_of_shared_blobs) ...
= offsets_ss{ii_o}{ii_ss}.proj_offsets(inds_of_shared_blobs) + proj_shift; = offsets_ss{ii_o}{ii_ss}.proj_offsets(inds_of_shared_blobs) + proj_shift;
...@@ -201,7 +218,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -201,7 +218,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
end end
else else
for ii_o = 1:numel(offsets) for ii_o = 1:numel(offsets)
inds_of_shared_blobs = offsets{ii_o}.blob_offsets == shared_blob_pos; inds_of_shared_blobs = offsets{ii_o}.blob_offsets == shared_blob_pos_inc;
offsets{ii_o}.proj_offsets(inds_of_shared_blobs) ... offsets{ii_o}.proj_offsets(inds_of_shared_blobs) ...
= offsets{ii_o}.proj_offsets(inds_of_shared_blobs) + proj_shift; = offsets{ii_o}.proj_offsets(inds_of_shared_blobs) + proj_shift;
...@@ -209,14 +226,12 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -209,14 +226,12 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
end end
end end
end end
if (any(to_be_changed)) fprintf('\n')
fprintf('\n')
end
base_offset_shift = num_parent_blobs; base_offset_shift = num_parent_blobs;
for ii_o_reg = 2:num_regions for ii_o_reg = 2:num_regions
[shared, parent_pos_shared_bls, shared_pos] = ... [shared, parent_pos_shared_bls, shared_pos] = ...
self.get_shared_bls(selected_pos_in_included, sampler, ii_o_reg); self.get_shared_bls(ref_inc_to_sel, sampler, ii_o_reg);
% This shift is related to the new size of the shared % This shift is related to the new size of the shared
% blobs, to accomodate every other sampled region in % blobs, to accomodate every other sampled region in
...@@ -249,7 +264,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -249,7 +264,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% Treating the shared blobs % Treating the shared blobs
inds_of_blobs(~to_be_shifted) ... inds_of_blobs(~to_be_shifted) ...
= selected_pos_in_included(shared_pos(inds_of_blobs(~to_be_shifted))); = ref_inc_to_sel(shared_pos(inds_of_blobs(~to_be_shifted)));
proj_coeffs_ss{ii_o}{ii_ss}.blob_offsets = inds_of_blobs; proj_coeffs_ss{ii_o}{ii_ss}.blob_offsets = inds_of_blobs;
...@@ -285,7 +300,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -285,7 +300,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% Treating the shared blobs % Treating the shared blobs
inds_of_blobs(~to_be_shifted) ... inds_of_blobs(~to_be_shifted) ...
= selected_pos_in_included(shared_pos(inds_of_blobs(~to_be_shifted))); = ref_inc_to_sel(shared_pos(inds_of_blobs(~to_be_shifted)));
proj_coeffs{ii_o}.blob_offsets = inds_of_blobs; proj_coeffs{ii_o}.blob_offsets = inds_of_blobs;
...@@ -353,7 +368,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle ...@@ -353,7 +368,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
curr_sel = sampler(ii_o_reg).selected; curr_sel = sampler(ii_o_reg).selected;
% This contains the position of the shared blobs in % This contains the position of the shared blobs in
% parent's numbering % parent's numbering
shared_pos = curr_or.proj(self.det_index).shared_bls_parent(curr_sel); shared_pos = curr_or.proj(self.det_index).shared_bls_with_parent(curr_sel);
% This indicates the ones that are actually shared among % This indicates the ones that are actually shared among
% the twin's blobs % the twin's blobs
shared = shared_pos > 0; shared = shared_pos > 0;
......
...@@ -128,14 +128,6 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -128,14 +128,6 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
grs(ii_g) = gtGrainAllblobsFilterOrder(grs(ii_g), refgr_omind); grs(ii_g) = gtGrainAllblobsFilterOrder(grs(ii_g), refgr_omind);
end end
% At the moment we only support the configuration: parent + twins
% refl_matches links the positions (i.e the hklsp type) of shared reflections
% i.e. refl_matches= [2, 0; 0, 1] indicates that the reflection in
% twin(1).allblobs(2) shares the spot with parent.allblobs(1) and
% twin(2).allblobs(1) shares the spot with parent.allblobs(1)
[refl_matches_ref, refl_matches_twin] = gt6DGetMatchingReflections(refgr, grs(2:end), conf.angular_toll);
refgr_proj = refgr.proj(conf.det_index); refgr_proj = refgr.proj(conf.det_index);
if (numel(conf.use_raw_images) == 1) if (numel(conf.use_raw_images) == 1)
...@@ -152,10 +144,6 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -152,10 +144,6 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
ref_selected = refgr_proj.selected; ref_selected = refgr_proj.selected;
end end
if (conf.verbose > 1)
produce_matching_reflections_table(grs, conf, refl_matches_ref)
end
if (conf.verbose > 1 && false) if (conf.verbose > 1 && false)
produce_beam_intensity_renorm_figure(refgr, p, ref_ondet, ref_included, ref_selected) produce_beam_intensity_renorm_figure(refgr, p, ref_ondet, ref_included, ref_selected)
end end
...@@ -185,6 +173,8 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -185,6 +173,8 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
samp_ors = cell(num_grains, 1); samp_ors = cell(num_grains, 1);
size_refl_list = size(grs(1).allblobs.omega); size_refl_list = size(grs(1).allblobs.omega);
estim_orient_bbox = zeros(num_grains, 6);
% Determining orientation-space bounding boxes. They are assumed to be % Determining orientation-space bounding boxes. They are assumed to be
% detached from each others. % detached from each others.
for ii_g = 1:num_grains for ii_g = 1:num_grains
...@@ -194,17 +184,17 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -194,17 +184,17 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
'detector_index', conf.det_index, 'verbose', conf.verbose); 'detector_index', conf.det_index, 'verbose', conf.verbose);
r_vecs = sampler.guess_ODF_BB()'; r_vecs = sampler.guess_ODF_BB()';
estim_orient_bbox = [min(r_vecs, [], 1), max(r_vecs, [], 1)]; estim_orient_bbox(ii_g, :) = [min(r_vecs, [], 1), max(r_vecs, [], 1)];
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3); bbox_size_rod = estim_orient_bbox(ii_g, 4:6) - estim_orient_bbox(ii_g, 1:3);
% oversizing the orienation a bit % oversizing the orienation a bit
delta_bbox_size_rod = bbox_size_rod * (conf.ospace_oversize - 1) / 2; delta_bbox_size_rod = bbox_size_rod * (conf.ospace_oversize - 1) / 2;
estim_orient_bbox = estim_orient_bbox + [-delta_bbox_size_rod, delta_bbox_size_rod]; estim_orient_bbox(ii_g, :) = estim_orient_bbox(ii_g, :) + [-delta_bbox_size_rod, delta_bbox_size_rod];
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3); bbox_size_rod = estim_orient_bbox(ii_g, 4:6) - estim_orient_bbox(ii_g, 1:3);
bbox_size_deg = 2 * atand(bbox_size_rod); bbox_size_deg = 2 * atand(bbox_size_rod);
fprintf(' - Estimated extent: [%g, %g, %g] -> [%g, %g, %g]\n', estim_orient_bbox); fprintf(' - Estimated extent: [%g, %g, %g] -> [%g, %g, %g]\n', estim_orient_bbox(ii_g, :));
fprintf(' BBox size: %3.3f, %3.3f, %3.3f (deg)\n', bbox_size_deg); fprintf(' BBox size: %3.3f, %3.3f, %3.3f (deg)\n', bbox_size_deg);
refor = struct( ... refor = struct( ...
...@@ -213,23 +203,14 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -213,23 +203,14 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
refor = gtCalculateGrain(refor, p, 'ref_omind', refgr_omind); refor = gtCalculateGrain(refor, p, 'ref_omind', refgr_omind);
refor.bb_ors = gt6DCalculateFullSpaceCorners(estim_space_bbox_mm, ... refor.bb_ors = gt6DCalculateFullSpaceCorners(estim_space_bbox_mm, ...
bbox_size_mm, estim_orient_bbox, bbox_size_rod, refgr, p); bbox_size_mm, estim_orient_bbox(ii_g, :), bbox_size_rod, refgr, p);
samp_ors{ii_g} = refor; samp_ors{ii_g} = refor;
% Populating the extended projection data-structures % Populating the extended projection data-structures
if (ii_g > 1) if (ii_g == 1)
extended_projs(ii_g) = find_matches_with_refor(extended_projs(1), ... % We first deal with the parent grain, and we will create the
grs(ii_g), refl_matches_ref(:, ii_g-1), refl_matches_twin(:, ii_g-1), conf); % datastructures: refl_matches_ref and refl_matches_twin
% (2) twin variants might share reflections which are not
% shared with the parent - we currently do not attempt to
% detect and handle these
extended_projs(1).selected = extended_projs(1).selected ...
| extended_projs(ii_g).shared_to_be_selected_in_parent;
extended_projs(1).allreflbool_selected(ref_ondet(ref_included)) = extended_projs(1).selected;
else
extended_projs(ii_g).ondet = ref_ondet; extended_projs(ii_g).ondet = ref_ondet;
extended_projs(ii_g).included = ref_included; extended_projs(ii_g).included = ref_included;
extended_projs(ii_g).selected = ref_selected; extended_projs(ii_g).selected = ref_selected;
...@@ -244,6 +225,36 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -244,6 +225,36 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
extended_projs(ii_g).allreflbool_ondet = bool_ref_ondet; extended_projs(ii_g).allreflbool_ondet = bool_ref_ondet;
extended_projs(ii_g).allreflbool_included = bool_ref_inc; extended_projs(ii_g).allreflbool_included = bool_ref_inc;
extended_projs(ii_g).allreflbool_selected = bool_ref_sel; extended_projs(ii_g).allreflbool_selected = bool_ref_sel;
% At the moment we only support the configuration: parent + twins
% refl_matches links the positions (i.e the hklsp type) of shared reflections
% i.e. refl_matches= [2, 0; 0, 1] indicates that the reflection in
% twin(1).allblobs(2) shares the spot with parent.allblobs(1) and
% twin(2).allblobs(1) shares the spot with parent.allblobs(1)
fprintf(' - Determining matching reflections with other regions..')
dis_angle = zeros(1, num_grains-1);
for ii_g_2 = 2:num_grains
dis_angle(ii_g_2-1) = gtDisorientation(grs(1).R_vector', grs(ii_g_2).R_vector', symm);
end
angular_tollerance = max(conf.angular_toll, dis_angle * 1.05);
[matches_ref_to_twin, matches_twin_to_ref] = gt6DGetMatchingReflections(refgr, grs(2:end), angular_tollerance);
fprintf('\b\b: Done.\n')
if (conf.verbose > 1)
produce_matching_reflections_table(grs, conf, matches_ref_to_twin)
end
else
extended_projs(ii_g) = find_matches_with_refor(extended_projs(1), ...
grs(ii_g), matches_ref_to_twin(:, ii_g-1), matches_twin_to_ref(:, ii_g-1), conf);
% (2) twin variants might share reflections which are not
% shared with the parent - we currently do not attempt to
% detect and handle these
extended_projs(1).selected = extended_projs(1).selected ...
| extended_projs(ii_g).shared_to_be_selected_in_parent;
extended_projs(1).allreflbool_selected(ref_ondet(ref_included)) = extended_projs(1).selected;
end end
local_ondet = extended_projs(ii_g).ondet; local_ondet = extended_projs(ii_g).ondet;
...@@ -259,7 +270,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -259,7 +270,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
% one of them - so it would lead to inconsistent projections % one of them - so it would lead to inconsistent projections
for ii_g = 2 : num_grains for ii_g = 2 : num_grains
extended_projs(ii_g) = find_matches_with_refor(extended_projs(1), ... extended_projs(ii_g) = find_matches_with_refor(extended_projs(1), ...
grs(ii_g), refl_matches_ref(:, ii_g-1), refl_matches_twin(:, ii_g-1), conf); grs(ii_g), matches_ref_to_twin(:, ii_g-1), matches_twin_to_ref(:, ii_g-1), conf);
end end
samp_ors = [samp_ors{:}]; samp_ors = [samp_ors{:}];
...@@ -286,8 +297,8 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -286,8 +297,8 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
% Let's treat those blobs at the w edge 360->0 % Let's treat those blobs at the w edge 360->0
% (from the sampled orientations perspective) % (from the sampled orientations perspective)
refor_ws = refor.allblobs.omega(local_inc_refl) / gtGetOmegaStepDeg(p, conf.det_index); or_ws = samp_ors(ii_g).allblobs.omega(local_inc_refl) / gtGetOmegaStepDeg(p, conf.det_index);
uvw_tab{ii_g}(:, :, 3) = gtGrainAnglesTabularFix360deg(uvw_tab{ii_g}(:, :, 3), refor_ws, p); uvw_tab{ii_g}(:, :, 3) = gtGrainAnglesTabularFix360deg(uvw_tab{ii_g}(:, :, 3), or_ws, p);
img_bboxes{ii_g} = [ ... img_bboxes{ii_g} = [ ...
reshape(floor(min(uvw_tab{ii_g}, [], 2)), [], 3), ... reshape(floor(min(uvw_tab{ii_g}, [], 2)), [], 3), ...
...@@ -297,12 +308,27 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -297,12 +308,27 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
max_img_sizes(ii_g, :) = max(img_sizes{ii_g}(:, 1:2), [], 1); max_img_sizes(ii_g, :) = max(img_sizes{ii_g}(:, 1:2), [], 1);
if (conf.verbose) if (conf.verbose || true)
refor_ns = refor.allblobs.eta(local_inc_refl); or_ns = samp_ors(ii_g).allblobs.eta(local_inc_refl);
tbl_out_format = '%02d)du %8d, dv %8d, dw %8d, eta: %7.3f, omega: %7.2f <- used: %d, u: [%4d %4d], v: [%4d %4d], w: [%4d %4d], sel: %d';
tbl_output = [(1:numel(or_ns))', ...
img_sizes{ii_g}, or_ns, or_ws, ...
~inconvenient_etas{ii_g}, ...
img_bboxes{ii_g}(:, [1 4 2 5 3 6]), ...
extended_projs(ii_g).selected ];
fprintf('\nImages for orientation %d:\n', ii_g) fprintf('\nImages for orientation %d:\n', ii_g)
fprintf('%02d)du %8d, dv %8d, dw %8d, eta: %7.3f <- used: %d, u: [%4d %4d], v: [%4d %4d], w: [%4d %4d]\n', ... if (ii_g == 1)
[(1:numel(refor_ns))', img_sizes{ii_g}, refor_ns, ~inconvenient_etas{ii_g}, img_bboxes{ii_g}(:, [1 4 2 5 3 6]) ]'); fprintf([tbl_out_format '\n'], ...
tbl_output');
else
pos_of_shared_blobs = zeros(size(extended_projs(ii_g).shared_refl_inc));
pos_of_shared_blobs(extended_projs(ii_g).shared_refl_inc) = extended_projs(ii_g).shared_refl_pos_in_parent_inc;
fprintf([tbl_out_format ', shared: %d\n'], ...
[tbl_output, pos_of_shared_blobs]');
end
end end
else else
max_img_sizes(ii_g, :) = [ ... max_img_sizes(ii_g, :) = [ ...
...@@ -328,7 +354,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -328,7 +354,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
for ii_b = 1:num_blobs for ii_b = 1:num_blobs
num_chars = fprintf('%03d/%03d', ii_b, num_blobs); num_chars = fprintf('%03d/%03d', ii_b, num_blobs);
if (ii_g == 1 || ~extended_projs(ii_g).shared_refl(ii_b)) if (ii_g == 1 || ~extended_projs(ii_g).shared_refl_inc(ii_b))
blobs(ii_b) = load_blob( ... blobs(ii_b) = load_blob( ...
img_bboxes{ii_g}(ii_b, :), img_sizes{ii_g}(ii_b, :), ... img_bboxes{ii_g}(ii_b, :), img_sizes{ii_g}(ii_b, :), ...
stackUSize, stackVSize, p, conf); stackUSize, stackVSize, p, conf);
...@@ -338,7 +364,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -338,7 +364,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
end end
% Using the same because it gives the same size/uvw-limits % Using the same because it gives the same size/uvw-limits
if (ii_g > 1) 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); blobs(extended_projs(ii_g).shared_refl_inc) = samp_ors(1).proj(conf.det_index).bl(extended_projs(ii_g).shared_refl_pos_in_parent_inc);
end end
else else
fprintf(' - Loading segmented blobs..') fprintf(' - Loading segmented blobs..')
...@@ -347,16 +373,16 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -347,16 +373,16 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
blobs = grs(ii_g).proj.bl; blobs = grs(ii_g).proj.bl;
else else
blobs = gtFwdSimBlobDefinition('blob', num_blobs); blobs = gtFwdSimBlobDefinition('blob', num_blobs);
blobs(~extended_projs(ii_g).shared_refl) = grs(ii_g).proj(conf.det_index).bl(~extended_projs(ii_g).shared_refl_included_pos); blobs(~extended_projs(ii_g).shared_refl_inc) = grs(ii_g).proj(conf.det_index).bl(~extended_projs(ii_g).shared_refl_inc);
% Using the same because it gives the same size/uvw-limits % Using the same because it gives the same size/uvw-limits
blobs(extended_projs(ii_g).shared_refl) = samp_ors(1).proj(conf.det_index).bl(extended_projs(ii_g).shared_refl_pos_in_parent); blobs(extended_projs(ii_g).shared_refl_inc) = samp_ors(1).proj(conf.det_index).bl(extended_projs(ii_g).shared_refl_pos_in_parent_inc);
if (conf.verbose) if (conf.verbose)
parent_ids = grs(1).fwd_sim(conf.det_index).spotid(ref_included(extended_projs(ii_g).shared_refl_pos_in_parent))'; parent_ids = grs(1).fwd_sim(conf.det_index).spotid(ref_included(extended_projs(ii_g).shared_refl_pos_in_parent_inc))';
twin_ids = grs(ii_g).fwd_sim(conf.det_index).spotid(extended_projs(ii_g).included(extended_projs(ii_g).shared_refl))'; twin_ids = grs(ii_g).fwd_sim(conf.det_index).spotid(extended_projs(ii_g).included(extended_projs(ii_g).shared_refl_inc))';
fprintf('Shared included blobs:\n') fprintf('Shared included blobs:\n')
fprintf(' - %d: parent blob_id %d, twin blob_id %d\n', ... fprintf(' - %d: parent blob_id %d, twin blob_id %d\n', ...
[1:numel(find(extended_projs(ii_g).shared_refl)); parent_ids; twin_ids]) [1:numel(find(extended_projs(ii_g).shared_refl_inc)); parent_ids; twin_ids])
end end
end end
end end
...@@ -366,7 +392,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -366,7 +392,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
c = tic(); c = tic();
proj_bl_masks = project_masks(samp_ors(ii_g), extended_projs(ii_g), grs(1).fwd_sim(conf.det_index).gv_verts, p, conf); 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); blobs = assign_masks(blobs, proj_bl_masks, ~extended_projs(ii_g).shared_refl_inc);
end end
fprintf('\b\b: Done in %g seconds.\n - Equalizing blob sizes..', toc(c)) fprintf('\b\b: Done in %g seconds.\n - Equalizing blob sizes..', toc(c))
...@@ -407,12 +433,12 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -407,12 +433,12 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
% work, but it should be replaced by a more advanced global % work, but it should be replaced by a more advanced global
% solution for handling shared blobs by different orientation boxes % solution for handling shared blobs by different orientation boxes
if (ii_g > 1) if (ii_g > 1)
proj.shared_bls_parent = zeros(size(extended_projs(ii_g).shared_refl)); proj.shared_bls_with_parent = zeros(size(extended_projs(ii_g).shared_refl_inc));
proj.shared_bls_parent(extended_projs(ii_g).shared_refl) = extended_projs(ii_g).shared_refl_pos_in_parent; proj.shared_bls_with_parent(extended_projs(ii_g).shared_refl_inc) = extended_projs(ii_g).shared_refl_pos_in_parent_inc;
samp_ors(1).proj(conf.det_index).shared_bls_twins(extended_projs(ii_g).shared_refl_pos_in_parent) = true; samp_ors(1).proj(conf.det_index).shared_bls_with_twins(extended_projs(ii_g).shared_refl_pos_in_parent_inc) = true;
else else
proj.shared_bls_twins = false(num_blobs, 1); proj.shared_bls_with_twins = false(num_blobs, 1);
end end
samp_ors(ii_g).proj(conf.det_index) = proj; samp_ors(ii_g).proj(conf.det_index) = proj;
...@@ -449,7 +475,8 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat ...@@ -449,7 +475,8 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
fprintf('Saving to sample.mat..') fprintf('Saving to sample.mat..')
sample = GtSample.loadFromFile(); sample = GtSample.loadFromFile();
cl_info = GtPhase.makeCluster(grs_list, ... cl_info = GtPhase.makeCluster(grs_list, ...
cat(1, samp_ors(:).R_vector), estim_space_bbox_pix, [], 1); cat(1, samp_ors(:).R_vector), ...
estim_space_bbox_pix, estim_orient_bbox, 1);
sample.phases{phase_id}.setClusterInfo(grs_list, cl_info); sample.phases{phase_id}.setClusterInfo(grs_list, cl_info);
sample.saveToFile(); sample.saveToFile();
fprintf('\b\b: Done.\n') fprintf('\b\b: Done.\n')
...@@ -548,7 +575,6 @@ function produce_matching_reflections_table(grs, conf, refl_matches) ...@@ -548,7 +575,6 @@ function produce_matching_reflections_table(grs, conf, refl_matches)
acosd(abs(cosd(grs(1).allblobs.eta(vis_sel)))) < conf.min_eta, ... acosd(abs(cosd(grs(1).allblobs.eta(vis_sel)))) < conf.min_eta, ...
]'); ]');
end end
end end
function produce_beam_intensity_renorm_figure(refgr, p, ref_ondet, ref_included, ref_selected) function produce_beam_intensity_renorm_figure(refgr, p, ref_ondet, ref_included, ref_selected)
...@@ -611,9 +637,8 @@ function eproj = get6DExtendedProjDefinition(num_structs) ...@@ -611,9 +637,8 @@ function eproj = get6DExtendedProjDefinition(num_structs)
'allreflbool_ondet', repl_cell, ... 'allreflbool_ondet', repl_cell, ...
'allreflbool_included', repl_cell, ... 'allreflbool_included', repl_cell, ...
'allreflbool_selected', repl_cell, ... 'allreflbool_selected', repl_cell, ...
'shared_refl', repl_cell, ... 'shared_refl_inc', repl_cell, ...
'shared_refl_pos_in_parent', repl_cell, ... 'shared_refl_pos_in_parent_inc', repl_cell, ...
'shared_refl_included_pos', repl_cell, ...
'shared_to_be_selected_in_parent', repl_cell); 'shared_to_be_selected_in_parent', repl_cell);
end end
...@@ -624,7 +649,7 @@ function eproj = find_matches_with_refor(ref_eproj, twingr, matches_ref_to_twin, ...@@ -624,7 +649,7 @@ function eproj = find_matches_with_refor(ref_eproj, twingr, matches_ref_to_twin,
eproj = get6DExtendedProjDefinition(); eproj = get6DExtendedProjDefinition();
size_refl_list = size(ref_eproj.allreflbool_ondet); size_allbls_refl_list = size(ref_eproj.allreflbool_ondet);
ref_ondet = ref_eproj.ondet; ref_ondet = ref_eproj.ondet;
ref_included = ref_eproj.included; ref_included = ref_eproj.included;
...@@ -643,7 +668,7 @@ function eproj = find_matches_with_refor(ref_eproj, twingr, matches_ref_to_twin, ...@@ -643,7 +668,7 @@ function eproj = find_matches_with_refor(ref_eproj, twingr, matches_ref_to_twin,
% bool_ref_shared(find(refl_matches_ref > 0)) = true; % bool_ref_shared(find(refl_matches_ref > 0)) = true;
% list of twin reflections for which a shared spot has been % list of twin reflections for which a shared spot has been
% predicted % predicted
bool_twin_shared = false(size_refl_list); bool_twin_shared = false(size_allbls_refl_list);
bool_twin_shared(matches_twin_to_ref > 0) = true; bool_twin_shared(matches_twin_to_ref > 0) = true;
twingr_proj = twingr.proj(conf.det_index); twingr_proj = twingr.proj(conf.det_index);
...@@ -652,50 +677,55 @@ function eproj = find_matches_with_refor(ref_eproj, twingr, matches_ref_to_twin, ...@@ -652,50 +677,55 @@ function eproj = find_matches_with_refor(ref_eproj, twingr, matches_ref_to_twin,
twin_inc = twingr_proj.included; twin_inc = twingr_proj.included;
twin_sel = twingr_proj.selected; twin_sel = twingr_proj.selected;
bool_twin_ondet = false(size_refl_list); bool_twin_ondet = false(size_allbls_refl_list);
bool_twin_ondet(twin_ondet) = true; bool_twin_ondet(twin_ondet) = true;
bool_twin_inc = false(size_refl_list); bool_twin_inc = false(size_allbls_refl_list);
bool_twin_inc(twin_ondet(twin_inc)) = true; bool_twin_inc(twin_ondet(twin_inc)) = true;
bool_twin_sel = false(size_refl_list); bool_twin_sel = false(size_allbls_refl_list);
bool_twin_sel(twin_ondet(twin_inc(twin_sel))) = true; bool_twin_sel(twin_ondet(twin_inc(twin_sel))) = true;
% parent indices of shared reflections that are ondet for the twin % parent indices of shared reflections that are ondet for the twin
ref_shared_ondet = false(size_refl_list); ref_shared_ondet = false(size_allbls_refl_list);
ref_shared_ondet(matches_twin_to_ref(bool_twin_shared & bool_twin_ondet)) = true; ref_shared_ondet(matches_twin_to_ref(bool_twin_shared & bool_twin_ondet)) = true;
% ...and which were included for the parent % ...and which were included for the parent
ref_shared_ondet_inc = bool_ref_inc & ref_shared_ondet; ref_shared_ondet_inc = bool_ref_inc & ref_shared_ondet;
ref_shared_selected = ref_shared_ondet_inc & bool_ref_sel; ref_shared_selected = ref_shared_ondet_inc & bool_ref_sel;
% translate this back into twin indices % translate this back into twin indices
shared_refl = false(size_refl_list); twin_shared_refl = false(size_allbls_refl_list);
shared_refl(matches_ref_to_twin(ref_shared_ondet_inc)) = true; twin_shared_refl(matches_ref_to_twin(ref_shared_ondet_inc)) = true;
bool_twin_inc_redundant = bool_twin_inc & shared_refl; bool_twin_inc_redundant = bool_twin_inc & twin_shared_refl;
% If a shared & twin included spot is selected in parent, select % If a shared & twin included spot is selected in parent, select
% it for the twin, as well % it for the twin, as well
bool_twin_sel(bool_twin_inc(matches_ref_to_twin(ref_shared_selected))) = true; % XXX - NOTE: next step would be to include these spots in the twin,
% even if they were not included.. but we have to take care that we
% handle correctly the copy from the previous data-structure, in case
% that we use segmented blobs
bool_twin_sel(matches_ref_to_twin(ref_shared_selected)) = bool_twin_inc(matches_ref_to_twin(ref_shared_selected));
twin_sel = bool_twin_sel(twin_ondet(twin_inc)); twin_sel = bool_twin_sel(twin_ondet(twin_inc));
bool_ref_matched_twin_inc = false(size_refl_list); bool_ref_matched_twin_inc = false(size_allbls_refl_list);
bool_ref_matched_twin_inc(matches_twin_to_ref(bool_twin_inc_redundant)) = true; bool_ref_matched_twin_inc(matches_twin_to_ref(bool_twin_inc_redundant)) = true;
shared_to_be_selected_in_parent = bool_ref_matched_twin_inc & ~bool_ref_sel; shared_to_be_selected_in_parent = bool_ref_matched_twin_inc & ~bool_ref_sel;
eproj.shared_to_be_selected_in_parent = shared_to_be_selected_in_parent(...
ref_ondet(ref_included)); eproj.shared_to_be_selected_in_parent ...
= shared_to_be_selected_in_parent(ref_ondet(ref_included));
% We now find the blobs in the old twin bl structure that are % We now find the blobs in the old twin bl structure that are
% shared with the parent % shared with the parent
eproj.shared_refl_included_pos = bool_twin_inc_redundant(twin_ondet(twin_inc)); eproj.shared_refl_inc = bool_twin_inc_redundant(twin_ondet(twin_inc));
% We cut shared_refl down to the size of the included % We cut shared_refl down to the size of the included
eproj.shared_refl = shared_refl(twin_ondet(twin_inc)); eproj.shared_refl_inc = twin_shared_refl(twin_ondet(twin_inc));
% We have to find the included IDs of the corresponding shared % We have to find the included IDs of the corresponding shared
% reflections % reflections
refl_id_in_parent = zeros(size_refl_list); refl_id_in_parent_inc = zeros(size_allbls_refl_list);
refl_id_in_parent(ref_ondet(ref_included)) = (1:numel(ref_included)); refl_id_in_parent_inc(ref_ondet(ref_included)) = (1:numel(ref_included));
eproj.shared_refl_pos_in_parent = refl_id_in_parent(matches_twin_to_ref(bool_twin_inc_redundant)); eproj.shared_refl_pos_in_parent_inc = refl_id_in_parent_inc(matches_twin_to_ref(bool_twin_inc_redundant));
eproj.ondet = twin_ondet; eproj.ondet = twin_ondet;
eproj.included = twin_inc; eproj.included = twin_inc;
......
function [ref_match, twin_match] = gt6DGetMatchingReflections(ref_or, ors, angular_toll, cryst) function [matches_ref_to_twin, matches_twin_to_ref] = gt6DGetMatchingReflections(ref_or, ors, angular_tolls)
% ref_match contains the indices of shared reflections in the % ref_match contains the indices of shared reflections in the
% twin_gr.allblobs list: i.e. [0 10; ...] would tell us that % twin_gr.allblobs list: i.e. [0 10; ...] would tell us that
% ref.allblobs(1) shares a spot with twin_or(2).allblobs(10) % ref.allblobs(1) shares a spot with twin_or(2).allblobs(10)
...@@ -7,8 +7,15 @@ function [ref_match, twin_match] = gt6DGetMatchingReflections(ref_or, ors, angul ...@@ -7,8 +7,15 @@ function [ref_match, twin_match] = gt6DGetMatchingReflections(ref_or, ors, angul
% ref_gr.allblobs list: i.e. [10 0; ...] would tell us that % ref_gr.allblobs list: i.e. [10 0; ...] would tell us that
% twin_or(1).allblobs(1) shares a spot with ref_gr.allblobs(10) % twin_or(1).allblobs(1) shares a spot with ref_gr.allblobs(10)
if (~exist('angular_toll', 'var') || isempty(angular_toll)) angular_tolls = reshape(angular_tolls, 1, []);
angular_toll = 1; % degrees if (~exist('angular_toll', 'var') || isempty(angular_tolls))
angular_tolls = 1; % degrees
end
num_ors = numel(ors);
if (numel(angular_tolls) == 1 && num_ors > 1)
angular_tolls = angular_tolls(1, ones(num_ors, 1));
end end
ref_ws = ref_or.allblobs.omega; ref_ws = ref_or.allblobs.omega;
...@@ -16,36 +23,33 @@ function [ref_match, twin_match] = gt6DGetMatchingReflections(ref_or, ors, angul ...@@ -16,36 +23,33 @@ function [ref_match, twin_match] = gt6DGetMatchingReflections(ref_or, ors, angul
ref_tt = ref_or.allblobs.thetatype; ref_tt = ref_or.allblobs.thetatype;
ref_Lfactors = ref_or.allblobs.lorentzfac; ref_Lfactors = ref_or.allblobs.lorentzfac;
w_toll = angular_toll * ref_Lfactors; num_spots = numel(ref_ws);
n_toll = angular_toll; angular_tolls = angular_tolls(ones(num_spots, 1), :);
w_toll = angular_tolls .* ref_Lfactors(:, ones(num_ors, 1));
n_toll = angular_tolls;
ors_allblobs = [ors(:).allblobs]; ors_allblobs = [ors(:).allblobs];
ors_ws = cat(2, ors_allblobs(:).omega); ors_ws = cat(2, ors_allblobs(:).omega);
ors_ns = cat(2, ors_allblobs(:).eta); ors_ns = cat(2, ors_allblobs(:).eta);
ors_tt = cat(2, ors_allblobs(:).thetatype); ors_tt = cat(2, ors_allblobs(:).thetatype);
num_ors = numel(ors);
num_spots = numel(ref_ws);
ones_array = ones(num_spots, num_ors); ones_array = ones(num_spots, num_ors);
ref_match = zeros(num_spots, num_ors); matches_ref_to_twin = zeros(num_spots, num_ors);
twin_match = zeros(num_spots, num_ors); matches_twin_to_ref = zeros(num_spots, num_ors);
row_ind = [];
for ii = 1 : num_spots for ii_refl = 1:num_spots
ref_ws_array = ref_ws(ii) * ones_array; ref_ws_array = ref_ws(ii_refl) * ones_array;
ref_ns_array = ref_ns(ii) * ones_array; ref_ns_array = ref_ns(ii_refl) * ones_array;
ref_tt_array = ref_tt(ii) * ones_array; ref_tt_array = ref_tt(ii_refl) * ones_array;
tt_equal = ref_tt_array == ors_tt; tt_equal = ref_tt_array == ors_tt;
w_diff = abs(mod(ref_ws_array - ors_ws + 180, 360) - 180); w_diff = abs(mod(ref_ws_array - ors_ws + 180, 360) - 180);
n_diff = abs(mod(ref_ns_array - ors_ns + 180, 360) - 180); n_diff = abs(mod(ref_ns_array - ors_ns + 180, 360) - 180);
[row_ind, col_ind] = find((w_diff < w_toll(:,ones(num_ors,1))) & (n_diff < n_toll) & tt_equal); [twin_refl_inds, twin_inds] = find((w_diff < w_toll) & (n_diff < n_toll) & tt_equal);
if ~isempty(row_ind)
for jj = 1 : numel(row_ind) for jj = 1:numel(twin_refl_inds)
ref_match(ii, col_ind(jj)) = row_ind(jj); matches_twin_to_ref(twin_refl_inds(jj), twin_inds(jj)) = ii_refl;
twin_match(row_ind(jj), col_ind(jj)) = ii; matches_ref_to_twin(ii_refl, twin_inds(jj)) = twin_refl_inds(jj);
end
end end
end end
end end
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment