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

6D-Reconstruction: updated Cluster/Twins/Extended creation functions to share... 

6D-Reconstruction: updated Cluster/Twins/Extended creation functions to share more code and support multi-detector

Signed-off-by: default avatarNicola Vigano <nicola.vigano@esrf.fr>
parent 1fea31c7
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 203 additions and 256 deletions
zUtil_Deformation/gt6DCalculateFullSpaceCorners.m 0 → 100644
+ 30
− 0
View file @ 873df95b
function orients = gt6DCalculateFullSpaceCorners(estim_space_bbox_mm, bbox_size_mm, estim_orient_bbox, bbox_size_rod, refgr, p)
orients = cell(2, 2, 2, 2, 2, 2);
for ii_x1 = 0:1
base_x1 = estim_space_bbox_mm(1:3) + [ii_x1 * bbox_size_mm(1), 0, 0];
for ii_x2 = 0:1
base_x2 = base_x1 + [0, ii_x2 * bbox_size_mm(2), 0];
for ii_x3 = 0:1
base_x3 = base_x2 + [0, 0, ii_x3 * bbox_size_mm(3)];
for ii_r1 = 0:1
base_r1 = estim_orient_bbox(1:3) + [ii_r1 * bbox_size_rod(1), 0, 0];
for ii_r2 = 0:1
base_r2 = base_r1 + [0, ii_r2 * bbox_size_rod(2), 0];
for ii_r3 = 0:1
base_r3 = base_r2 + [0, 0, ii_r3 * bbox_size_rod(3)];
% fprintf(' center: %f, %f, %f, r_vec: %f, %f, %f\n', base_x3, base_r3)
orients{ii_r3+1, ii_r2+1, ii_r1+1, ii_x3+1, ii_x2+1, ii_x1+1} = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, ...
'center', base_x3, 'R_vector', base_r3);
end
end
end
end
end
end
orients = gtCalculateGrain_p(orients, p, 'ref_omind', refgr.allblobs.omind);
orients = [orients{:}];
end
\ No newline at end of file
zUtil_Deformation/gt6DCreateProjDataFromExtendedGrain.m
+ 28
− 78
View file @ 873df95b
......@@ -13,6 +13,7 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
conf = struct( ...
'verbose', false, ...
'min_eta', 20, ...
'det_index', 1, ...
'ospace_oversize', 1.1, ...
'rspace_oversize', 1.1, ...
'ospace_lims', [], ...
......@@ -31,20 +32,10 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
p = gtLoadParameters();
refgr = gr;
refgr_omind = refgr.allblobs.omind;
refgr_omind_ind = [ ...
find(refgr_omind == 1); find(refgr_omind == 2); ...
find(refgr_omind == 3); find(refgr_omind == 4) ];
if (isfield(refgr.proj, 'ondet'))
ref_ondet = refgr.proj.ondet;
ref_included = refgr.proj.included;
ref_selected = refgr.proj.selected;
else
ref_ondet = refgr.ondet;
ref_included = refgr.included;
ref_selected = refgr.selected;
end
ref_ondet = refgr.proj(conf.det_index).ondet;
ref_included = refgr.proj(conf.det_index).included;
ref_selected = refgr.proj(conf.det_index).selected;
fwdsim_sel_refl = ref_ondet(ref_included(ref_selected));
fwdsim_sel_refl_logic = false(size(ref_ondet));
......@@ -58,12 +49,8 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
used_refl = fwdsim_sel_refl;
end
vol_half_size = [gr.proj.vol_size_y, gr.proj.vol_size_x, gr.proj.vol_size_z] / 2;
space_bbox = [ ...
gr.proj.centerpix - vol_half_size * conf.rspace_oversize, ...
gr.proj.centerpix + vol_half_size * conf.rspace_oversize ];
estim_space_bbox_pix = [floor(space_bbox(1:3)), ceil(space_bbox(4:6))];
bbox_size_pix = estim_space_bbox_pix(4:6) - estim_space_bbox_pix(1:3);
[gr_center_pix, bbox_size_pix, estim_space_bbox_pix] = gt6DMergeRealSpaceVolumes(gr, conf.det_index);
gr_center_mm = gtGeoSam2Sam(gr_center_pix, p.recgeo, p.samgeo, false, false);
estim_space_bbox_mm = [ ...
gtGeoSam2Sam(estim_space_bbox_pix(1:3), p.recgeo, p.samgeo, false, false), ...
......@@ -71,7 +58,8 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
bbox_size_mm = estim_space_bbox_mm(4:6) - estim_space_bbox_mm(1:3);
if (isempty(conf.ospace_lims))
sampler = GtOrientationSampling(p, gr, 'verbose', conf.verbose);
sampler = GtOrientationSampling(p, gr, ...
'detector_index', conf.det_index, 'verbose', conf.verbose);
r_vecs = sampler.guess_ODF_BB()';
estim_orient_bbox = [min(r_vecs, [], 1), max(r_vecs, [], 1)];
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3);
......@@ -85,44 +73,17 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
end
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3);
gr_center_rod = (estim_orient_bbox(4:6) + estim_orient_bbox(1:3)) / 2;
bbox_size_deg = 2 * atand(bbox_size_rod);
% Let's now compute the bb on the images, by computing for each corner
% of the space bb, the position on the detector of each corner of the
% orientation bb.
or = {};
for ii_x1 = 0:1
base_x1 = estim_space_bbox_mm(1:3) + [ii_x1 * bbox_size_mm(1), 0, 0];
for ii_x2 = 0:1
base_x2 = base_x1 + [0, ii_x2 * bbox_size_mm(2), 0];
for ii_x3 = 0:1
base_x3 = base_x2 + [0, 0, ii_x3 * bbox_size_mm(3)];
for ii_r1 = 0:1
base_r1 = estim_orient_bbox(1:3) + [ii_r1 * bbox_size_rod(1), 0, 0];
for ii_r2 = 0:1
base_r2 = base_r1 + [0, ii_r2 * bbox_size_rod(2), 0];
for ii_r3 = 0:1
base_r3 = base_r2 + [0, 0, ii_r3 * bbox_size_rod(3)];
% fprintf(' center: %f, %f, %f, r_vec: %f, %f, %f\n', base_x3, base_r3)
or{end+1} = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, ...
'center', base_x3, 'R_vector', base_r3); %#ok<AGROW>
end
end
end
end
end
end
or = gtCalculateGrain_p(or, p);
or = [or{:}];
num_ors = numel(or);
for ii_g = num_ors:-1:1
or(ii_g) = gtGrainAllblobsFilterOrder(or(ii_g), refgr_omind_ind);
end
ors = gt6DCalculateFullSpaceCorners( ...
estim_space_bbox_mm, bbox_size_mm, ...
estim_orient_bbox, bbox_size_rod, ...
refgr, p);
if (conf.verbose)
fprintf('\n');
......@@ -134,11 +95,11 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
fprintf('\n');
end
or_abs = cat(1, or(:).allblobs);
or_abs = cat(1, ors(:).allblobs);
uvw_tab = zeros(numel(used_refl), num_ors, 3);
for ii_g = 1:num_ors
if (isfield(or_abs(1), 'detector'))
uvw_tab(:, ii_g, :) = or_abs(ii_g).detector(1).uvw(used_refl, :);
uvw_tab(:, ii_g, :) = or_abs(ii_g).detector(conf.det_index).uvw(used_refl, :);
else
uvw_tab(:, ii_g, :) = or_abs(ii_g).uvw(used_refl, :);
end
......@@ -146,25 +107,13 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
refor = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, ...
'center', estim_space_bbox_mm(1:3) + bbox_size_mm / 2, ...
'R_vector', estim_orient_bbox(1:3) + bbox_size_rod / 2 );
refor = gtCalculateGrain(refor, p);
refor = gtGrainAllblobsFilterOrder(refor, refgr_omind_ind);
num_images = gtGetTotNumberOfImages(p);
refor_ws = refor.allblobs.omega(used_refl) / p.labgeo.omstep;
'center', gr_center_mm, 'R_vector', gr_center_rod );
refor = gtCalculateGrain(refor, p, 'ref_omind', refgr.allblobs.omind);
% Let's treat those blobs at the w edge 360->0
% (from the sampled orientations perspective)
opposite_ws = mod(refor_ws + num_images / 2, num_images);
gt_ref_int = opposite_ws > num_images / 2;
lt_ref_int = ~gt_ref_int;
opposite_ws_repmat = opposite_ws(:, ones(1, size(uvw_tab(:, :, 3), 2)));
uvw_tab(gt_ref_int, :, 3) = uvw_tab(gt_ref_int, :, 3) ...
- num_images .* (uvw_tab(gt_ref_int, :, 3) > opposite_ws_repmat(gt_ref_int, :));
uvw_tab(lt_ref_int, :, 3) = uvw_tab(lt_ref_int, :, 3) ...
+ num_images .* (uvw_tab(lt_ref_int, :, 3) < opposite_ws_repmat(lt_ref_int, :));
refor_ws = refor.allblobs.omega(used_refl) / gtGetOmegaStepDeg(p, conf.det_index);
uvw_tab(:, :, 3) = gtGrainAnglesTabularFix360deg(uvw_tab(:, :, 3), refor_ws, p);
img_bboxes_initial = [
squeeze(floor(min(uvw_tab, [], 2))), ...
......@@ -214,7 +163,7 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
end
sel_reflections = sel_reflections(~inconvenient_etas);
blobs(1:num_blobs) = gtFwdSimBlobDefinition();
blobs = gtFwdSimBlobDefinition('blob', num_blobs);
num_sel_refl = numel(sel_reflections);
for ii_i = 1:num_sel_refl
......@@ -223,7 +172,7 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
ii_b = sel_reflections(ii_i);
bb = [img_bboxes(ii_i, 1:2), img_sizes(ii_i, 1:2)];
blob_vol = gtGetRawRoi(img_bboxes(ii_i, 3), img_bboxes(ii_i, 6), p.acq, bb);
blob_vol = gtGetRawRoi(img_bboxes(ii_i, 3), img_bboxes(ii_i, 6), p.acq, bb, conf.det_index);
blob_vol(blob_vol < 0) = 0;
blob_vol(isnan(blob_vol)) = 0;
blob_bb = [img_bboxes(ii_i, 1:3), img_sizes(ii_i, :)];
......@@ -298,6 +247,8 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
end
fprintf('Done.\n')
proj = gtFwdSimProjDefinition();
proj.ondet = ref_ondet;
if (~conf.include_all)
proj.included = ref_included;
......@@ -338,8 +289,7 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
spots = permute(cat(3, spots{:}), [1 3 2]);
% These reconstructions usually need some padding
proj.centerpix = estim_space_bbox_pix(1:3) + bbox_size_pix / 2;
% proj.centerpix = gtGeoSam2Sam(refor.center, p.samgeo, p.recgeo, false, false);
proj.centerpix = gr_center_pix;
proj.bl = blobs;
proj.stack = spots;
% We apply a simple oversize
......@@ -348,8 +298,8 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
proj.vol_size_y = vol_size(1);
proj.vol_size_z = vol_size(3);
refor.proj = proj;
refor.bb_ors = or;
refor.proj(conf.det_index) = proj;
refor.bb_ors = ors;
if (conf.save)
fprintf('Saving the extended grain file..')
......
zUtil_Deformation/gt6DCreateProjDataFromGrainCluster.m
+ 40
− 90
View file @ 873df95b
......@@ -13,6 +13,8 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
conf = struct( ...
'verbose', false, ...
'min_eta', 15, ...
'det_index', 1, ...
'ospace_oversize', 1.1, ...
'save', false );
conf = parse_pv_pairs(conf, varargin);
......@@ -51,96 +53,49 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
refgr = grs(1);
refgr_omind = refgr.allblobs.omind;
refgr_omind_ind = [ ...
find(refgr_omind == 1); find(refgr_omind == 2); ...
find(refgr_omind == 3); find(refgr_omind == 4) ];
for ii_g = num_grains:-1:2
grs(ii_g) = gtGrainAllblobsFilterOrder(grs(ii_g), refgr_omind_ind);
grs(ii_g) = gtGrainAllblobsFilterOrder(grs(ii_g), refgr_omind);
end
if (isfield(refgr.proj, 'ondet'))
ref_ondet = refgr.proj.ondet;
ref_included = refgr.proj.included;
ref_selected = refgr.proj.selected;
else
ref_ondet = refgr.ondet;
ref_included = refgr.included;
ref_selected = refgr.selected;
end
ref_ondet = refgr.proj(conf.det_index).ondet;
ref_included = refgr.proj(conf.det_index).included;
ref_selected = refgr.proj(conf.det_index).selected;
fwdsim_sel_refl = ref_ondet(ref_included(ref_selected));
space_bboxes = zeros(num_grains, 6);
for ii_g = 1:num_grains
vol_half_size = [grs(ii_g).proj.vol_size_y, grs(ii_g).proj.vol_size_x, grs(ii_g).proj.vol_size_z] / 2;
space_bboxes(ii_g, :) = [ ...
grs(ii_g).proj.centerpix - vol_half_size, ...
grs(ii_g).proj.centerpix + vol_half_size ];
end
estim_space_bbox_pix = [floor(min(space_bboxes(:, 1:3), [], 1)), ceil(max(space_bboxes(:, 4:6), [], 1))];
bbox_size_pix = estim_space_bbox_pix(4:6) - estim_space_bbox_pix(1:3);
% Real-space volume size
[gr_center_pix, bbox_size_pix, estim_space_bbox_pix] = gt6DMergeRealSpaceVolumes(grs, conf.det_index);
gr_center_mm = gtGeoSam2Sam(gr_center_pix, p.recgeo, p.samgeo, false, false);
estim_space_bbox_mm = [ ...
gtGeoSam2Sam(estim_space_bbox_pix(1:3), p.recgeo, p.samgeo, false, false), ...
gtGeoSam2Sam(estim_space_bbox_pix(4:6), p.recgeo, p.samgeo, false, false) ];
bbox_size_mm = estim_space_bbox_mm(4:6) - estim_space_bbox_mm(1:3);
% Orientation-space volume size
for ii_g = num_grains:-1:1
sampler = GtOrientationSampling(p, grs(ii_g), 'verbose', conf.verbose);
sampler = GtOrientationSampling(p, grs(ii_g), ...
'detector_index', conf.det_index, 'verbose', conf.verbose);
odfw_R_vectors{ii_g} = sampler.guess_ODF_BB()';
end
r_vecs = cat(1, odfw_R_vectors{:});
estim_orient_bbox = [min(r_vecs, [], 1), max(r_vecs, [], 1)];
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3);
% oversizing the orienation a bit
delta_bbox_size_rod = bbox_size_rod * 0.05;
% oversizing the orienation BB a bit
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];
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3);
bbox_size_deg = 2 * atand(bbox_size_rod);
% Let's now compute the bb on the images, by computing for each corner
% of the space bb, the position on the detector of each corner of the
% orientation bb.
or = {};
for ii_x1 = 0:1
base_x1 = estim_space_bbox_mm(1:3) + [ii_x1 * bbox_size_mm(1), 0, 0];
for ii_x2 = 0:1
base_x2 = base_x1 + [0, ii_x2 * bbox_size_mm(2), 0];
for ii_x3 = 0:1
base_x3 = base_x2 + [0, 0, ii_x3 * bbox_size_mm(3)];
for ii_r1 = 0:1
base_r1 = estim_orient_bbox(1:3) + [ii_r1 * bbox_size_rod(1), 0, 0];
for ii_r2 = 0:1
base_r2 = base_r1 + [0, ii_r2 * bbox_size_rod(2), 0];
for ii_r3 = 0:1
base_r3 = base_r2 + [0, 0, ii_r3 * bbox_size_rod(3)];
% fprintf(' center: %f, %f, %f, r_vec: %f, %f, %f\n', base_x3, base_r3)
or{end+1} = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, ...
'center', base_x3, 'R_vector', base_r3); %#ok<AGROW>
end
end
end
end
end
end
or = gtCalculateGrain_p(or, p);
or = [or{:}];
num_ors = numel(or);
gr_center_rod = (estim_orient_bbox(4:6) + estim_orient_bbox(1:3)) / 2;
for ii_g = num_ors:-1:1
or(ii_g) = gtGrainAllblobsFilterOrder(or(ii_g), refgr_omind_ind);
end
bbox_size_deg = 2 * atand(bbox_size_rod);
if (conf.verbose)
f = figure();
ax = axes('parent', f);
hold(ax, 'on');
gt6DPlotOrientationBBox(ax, cat(1, or(:).R_vector));
gt6DPlotOrientationBBox(ax, [estim_orient_bbox(1:3); estim_orient_bbox(4:6)]);
for ii_g = 1:num_grains
scatter3(ax, grs(ii_g).R_vector(1), grs(ii_g).R_vector(2), grs(ii_g).R_vector(3), 30);
end
......@@ -157,37 +112,30 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
fprintf('\n');
end
or_abs = cat(1, or(:).allblobs);
% Let's now compute the bb on the images, by computing for each corner
% of the space bb, the position on the detector of each corner of the
% orientation bb.
ors = gt6DCalculateFullSpaceCorners( ...
estim_space_bbox_mm, bbox_size_mm, ...
estim_orient_bbox, bbox_size_rod, ...
refgr, p);
or_abs = cat(1, ors(:).allblobs);
num_ors = numel(ors);
uvw_tab = zeros(numel(fwdsim_sel_refl), num_ors, 3);
for ii_g = 1:num_ors
if (isfield(or_abs(1), 'detector'))
uvw_tab(:, ii_g, :) = or_abs(ii_g).detector(1).uvw(fwdsim_sel_refl, :);
else
uvw_tab(:, ii_g, :) = or_abs(ii_g).uvw(fwdsim_sel_refl, :);
end
uvw_tab(:, ii_g, :) = or_abs(ii_g).detector(conf.det_index).uvw(fwdsim_sel_refl, :);
end
refor = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, 'gids', grs_list, ...
'center', estim_space_bbox_mm(1:3) + bbox_size_mm / 2, ...
'R_vector', estim_orient_bbox(1:3) + bbox_size_rod / 2 );
refor = gtCalculateGrain(refor, p);
refor = gtGrainAllblobsFilterOrder(refor, refgr_omind_ind);
num_images = gtGetTotNumberOfImages(p);
refor_ws = refor.allblobs.omega(fwdsim_sel_refl) / p.labgeo.omstep;
'center', gr_center_mm, 'R_vector', gr_center_rod );
refor = gtCalculateGrain(refor, p, 'ref_omind', refgr_omind);
% Let's treat those blobs at the w edge 360->0
% (from the sampled orientations perspective)
opposite_ws = mod(refor_ws + num_images / 2, num_images);
gt_ref_int = opposite_ws > num_images / 2;
lt_ref_int = ~gt_ref_int;
opposite_ws_repmat = opposite_ws(:, ones(1, size(uvw_tab(:, :, 3), 2)));
uvw_tab(gt_ref_int, :, 3) = uvw_tab(gt_ref_int, :, 3) ...
- num_images .* (uvw_tab(gt_ref_int, :, 3) > opposite_ws_repmat(gt_ref_int, :));
uvw_tab(lt_ref_int, :, 3) = uvw_tab(lt_ref_int, :, 3) ...
+ num_images .* (uvw_tab(lt_ref_int, :, 3) < opposite_ws_repmat(lt_ref_int, :));
refor_ws = refor.allblobs.omega(fwdsim_sel_refl) / gtGetOmegaStepDeg(p, conf.det_index);
uvw_tab(:, :, 3) = gtGrainAnglesTabularFix360deg(uvw_tab(:, :, 3), refor_ws, p);
img_bboxes = [
squeeze(floor(min(uvw_tab, [], 2))), ...
......@@ -222,7 +170,7 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
fprintf('Loading raw images: ')
num_blobs = numel(ref_included);
blobs(1:num_blobs) = gtFwdSimBlobDefinition();
blobs = gtFwdSimBlobDefinition('blob', num_blobs);
sel_reflections = find(ref_selected);
sel_reflections = sel_reflections(~inconvenient_etas);
......@@ -234,7 +182,7 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
ii_b = sel_reflections(ii_i);
bb = [img_bboxes(ii_i, 1:2), img_sizes(ii_i, 1:2)];
blob_vol = gtGetRawRoi(img_bboxes(ii_i, 3), img_bboxes(ii_i, 6), p.acq, bb);
blob_vol = gtGetRawRoi(img_bboxes(ii_i, 3), img_bboxes(ii_i, 6), p.acq, bb, conf.det_index);
blob_vol(blob_vol < 0) = 0;
blob_bb = [img_bboxes(ii_i, 1:3), img_sizes(ii_i, :)];
......@@ -318,7 +266,9 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
spots = arrayfun(@(x){sum(x.intm, 3)}, blobs);
spots = permute(cat(3, spots{:}), [1 3 2]);
proj.centerpix = estim_space_bbox_pix(1:3) + bbox_size_pix / 2;
proj = gtFwdSimProjDefinition();
proj.centerpix = gr_center_pix;
proj.bl = blobs;
proj.stack = spots;
vol_size = bbox_size_pix + mean(bbox_size_pix) * 0.3;
......@@ -333,8 +283,8 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDataFr
finally_selected(sel_reflections) = true;
proj.selected = finally_selected;
refor.proj = proj;
refor.bb_ors = or;
refor.proj(conf.det_index) = proj;
refor.bb_ors = ors;
if (conf.save)
fprintf('Saving the cluster file..')
......
zUtil_Deformation/gt6DCreateProjDataFromTwinnedGrainFwdProj.m
+ 38
− 88
View file @ 873df95b
......@@ -14,9 +14,10 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
conf = struct( ...
'verbose', false, ...
'det_ind', 1, ...
'det_index', 1, ...
'min_eta', 15, ...
'rspace_oversize', 1.1, ...
'ospace_oversize', 1.1, ...
'check_spots', false, ...
'stack_oversize', 1.4, ...
'save', false );
......@@ -106,9 +107,13 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
end
refgr = grs(1);
grs(2:end) = reorder_reflections(grs(2:end), refgr);
refgr_proj = refgr.proj(conf.det_ind);
refgr_omind = refgr.allblobs.omind;
for ii_g = 2:numel(grs)
grs(ii_g) = gtGrainAllblobsFilterOrder(grs(ii_g), refgr_omind); %#ok<AGROW>
end
refgr_proj = refgr.proj(conf.det_index);
ref_ondet = refgr_proj.ondet;
ref_included = refgr_proj.included;
......@@ -162,7 +167,6 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
if (conf.verbose)
vis_sel = ref_ondet(ref_included);
% vis_sel = ref_ondet;
fprintf('%02d) w1 %6.2f, w2 %6.2f, n1 %6.2f, n2 %6.2f, w_inds (%d<->%d), hkls [%2d %2d %2d]<->[%2d %2d %2d] <- diff: %7.2f (!! %d, L %g) [bad n1 %d n2 %d]\n', ...
[(1:numel(vis_sel))', ...
grs(1).allblobs.omega(vis_sel), ...
......@@ -181,62 +185,57 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
]');
end
vol_size = [ ...
refgr_proj.vol_size_y, ...
refgr_proj.vol_size_x, ...
refgr_proj.vol_size_z ];
vol_half_size = vol_size / 2;
space_bbox = [ ...
refgr_proj.centerpix - vol_half_size * conf.rspace_oversize, ...
refgr_proj.centerpix + vol_half_size * conf.rspace_oversize ];
estim_space_bbox_pix = [floor(space_bbox(1:3)), ceil(space_bbox(4:6))];
bbox_size_pix = estim_space_bbox_pix(4:6) - estim_space_bbox_pix(1:3);
[gr_center_pix, bbox_size_pix, estim_space_bbox_pix] = gt6DMergeRealSpaceVolumes(grs, conf.det_index, conf.rspace_oversize);
gr_center_mm = gtGeoSam2Sam(gr_center_pix, p.recgeo, p.samgeo, false, false);
estim_space_bbox_mm = [ ...
gtGeoSam2Sam(estim_space_bbox_pix(1:3), p.recgeo, p.samgeo, false, false), ...
gtGeoSam2Sam(estim_space_bbox_pix(4:6), p.recgeo, p.samgeo, false, false) ];
bbox_size_mm = estim_space_bbox_mm(4:6) - estim_space_bbox_mm(1:3);
if (conf.verbose)
fprintf('\n');
fprintf('Estimated spatial voxel BBox: [%3d, %3d, %3d] -> [%3d, %3d, %3d]\n', estim_space_bbox_pix);
fprintf(' BBox size: %3d, %3d, %3d (%f, %f, %f mm)\n', bbox_size_pix, bbox_size_mm);
end
img_sizes = zeros(num_grains, 2);
for ii_g = num_grains:-1:1
sampler = GtOrientationSampling(p, grs(ii_g), 'verbose', conf.verbose);
sampler = GtOrientationSampling(p, grs(ii_g), ...
'detector_index', conf.det_index, 'verbose', conf.verbose);
r_vecs = sampler.guess_ODF_BB()';
estim_orient_bbox = [min(r_vecs, [], 1), max(r_vecs, [], 1)];
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3);
% oversizing the orienation a bit
delta_bbox_size_rod = bbox_size_rod * 0.05;
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];
bbox_size_rod = estim_orient_bbox(4:6) - estim_orient_bbox(1:3);
bbox_size_deg = 2 * atand(bbox_size_rod);
if (conf.verbose)
fprintf('\n');
fprintf('Estimated spatial voxel BBox: [%3d, %3d, %3d] -> [%3d, %3d, %3d]\n', estim_space_bbox_pix);
fprintf(' BBox size: %3d, %3d, %3d (%f, %f, %f mm)\n', bbox_size_pix, bbox_size_mm);
fprintf(' Estimated orientation BBox: [%3.3f, %3.3f, %3.3f] -> [%3.3f, %3.3f, %3.3f]\n', estim_orient_bbox);
fprintf(' BBox size: %3.3f, %3.3f, %3.3f (deg)\n', bbox_size_deg);
fprintf('%d) Estimated orientation BBox: [%3.3f, %3.3f, %3.3f] -> [%3.3f, %3.3f, %3.3f]\n', ii_g, estim_orient_bbox);
fprintf(' BBox size: %3.3f, %3.3f, %3.3f (deg)\n', bbox_size_deg);
fprintf('\n');
end
refor = struct( ...
'id', grs(ii_g).id, 'phaseid', grs(ii_g).phaseid, ...
'center', estim_space_bbox_mm(1:3) + bbox_size_mm / 2, ...
'R_vector', grs(ii_g).R_vector );
refor = gtCalculateGrain(refor, p);
'center', gr_center_mm, 'R_vector', grs(ii_g).R_vector );
refor = gtCalculateGrain(refor, p, 'ref_omind', refgr_omind);
refor.bb_ors = get_6D_space_corners(estim_space_bbox_mm, ...
refor.bb_ors = gt6DCalculateFullSpaceCorners(estim_space_bbox_mm, ...
bbox_size_mm, estim_orient_bbox, bbox_size_rod, refgr, p);
img_sizes(ii_g, :) = [size(grs(ii_g).proj.stack, 1), size(grs(ii_g).proj.stack, 3), ];
img_sizes(ii_g, :) = [ ...
size(grs(ii_g).proj(conf.det_index).stack, 1), ...
size(grs(ii_g).proj(conf.det_index).stack, 3), ];
samp_ors(ii_g) = refor;
end
samp_ors = reorder_reflections(samp_ors, refgr);
[stackUSize, stackVSize] = get_stack_UV_size(cat(1, img_sizes), p, conf);
......@@ -245,13 +244,16 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
bool_ref_sel = false(size(refgr.allblobs.omega));
bool_ref_sel(ref_ondet(ref_included(ref_selected))) = true;
% At the moment we only support the configuration: parent + twins
refl_matches = gt6DGetMatchingReflections(refgr, grs(2:end), 2);
for ii_g = 1:num_grains
fprintf('Loading raw images for grain %d: ', ii_g)
if (ii_g > 1)
twingr = grs(ii_g);
twingr_proj = twingr.proj(conf.det_ind);
twingr_proj = twingr.proj(conf.det_index);
twin_ondet = twingr_proj.ondet;
twin_inc = twingr_proj.included;
......@@ -268,17 +270,7 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
% twin, and were included for the parent
bool_test_inc = bool_ref_inc & bool_twin_ondet;
inc_ref_omegas = samp_ors(1).allblobs.omega(bool_test_inc);
inc_ref_Lfac = samp_ors(1).allblobs.lorentzfac(bool_test_inc);
inc_ref_etas = samp_ors(1).allblobs.eta(bool_test_inc);
inc_twin_omegas = samp_ors(ii_g).allblobs.omega(bool_test_inc);
inc_twin_etas = samp_ors(ii_g).allblobs.eta(bool_test_inc);
w_diff = abs(mod(inc_ref_omegas - inc_twin_omegas + 180, 360) - 180);
n_diff = abs(mod(inc_ref_etas - inc_twin_etas + 180, 360) - 180);
test_shared = (w_diff < inc_ref_Lfac) & (n_diff < 5);
test_shared = refl_matches(bool_test_inc, ii_g-1);
% let's add the shared as included
shared_parent = false(size(refgr.allblobs.omega));
......@@ -338,7 +330,9 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
spots = arrayfun(@(x){sum(x.intm, 3)}, blobs);
spots = permute(cat(3, spots{:}), [1 3 2]);
proj.centerpix = estim_space_bbox_pix(1:3) + bbox_size_pix / 2;
proj = gtFwdSimProjDefinition();
proj.centerpix = gr_center_pix;
proj.bl = blobs;
proj.stack = spots;
vol_size = bbox_size_pix + mean(bbox_size_pix) * 0.3;
......@@ -354,15 +348,15 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
proj.shared_bls_parent = zeros(size(shared_parent));
proj.shared_bls_parent(shared_parent) = find(shared_parent_pos);
samp_ors(ii_g).proj(conf.det_ind).shared_bls_twins(shared_parent_pos) = true;
samp_ors(1).proj(conf.det_index).shared_bls_twins(shared_parent_pos) = true;
else
proj.shared_bls_twins = false(num_blobs, 1);
end
samp_ors(ii_g).proj(conf.det_ind) = proj;
samp_ors(ii_g).proj(conf.det_index) = proj;
if (conf.check_spots && ii_g > 1)
samp_ors(ii_g).proj.selected = gtGuiGrainMontage(samp_ors(ii_g), [], true);
samp_ors(ii_g).proj.selected = gtGuiGrainMontage(samp_ors(ii_g), [], true, conf.det_index);
end
end
......@@ -397,49 +391,6 @@ function [samp_ors, estim_space_bbox_pix, estim_orient_bbox] = gt6DCreateProjDat
end
end
function orients = get_6D_space_corners(estim_space_bbox_mm, bbox_size_mm, estim_orient_bbox, bbox_size_rod, refgr, p)
orients = cell(2, 2, 2, 2, 2, 2);
for ii_x1 = 0:1
base_x1 = estim_space_bbox_mm(1:3) + [ii_x1 * bbox_size_mm(1), 0, 0];
for ii_x2 = 0:1
base_x2 = base_x1 + [0, ii_x2 * bbox_size_mm(2), 0];
for ii_x3 = 0:1
base_x3 = base_x2 + [0, 0, ii_x3 * bbox_size_mm(3)];
for ii_r1 = 0:1
base_r1 = estim_orient_bbox(1:3) + [ii_r1 * bbox_size_rod(1), 0, 0];
for ii_r2 = 0:1
base_r2 = base_r1 + [0, ii_r2 * bbox_size_rod(2), 0];
for ii_r3 = 0:1
base_r3 = base_r2 + [0, 0, ii_r3 * bbox_size_rod(3)];
% fprintf(' center: %f, %f, %f, r_vec: %f, %f, %f\n', base_x3, base_r3)
orients{ii_r3+1, ii_r2+1, ii_r1+1, ii_x3+1, ii_x2+1, ii_x1+1} = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, ...
'center', base_x3, 'R_vector', base_r3);
end
end
end
end
end
end
orients = gtCalculateGrain_p(orients, p);
orients = [orients{:}];
orients = reorder_reflections(orients, refgr);
end
function orients = reorder_reflections(orients, refgr)
refgr_omind = refgr.allblobs.omind;
refgr_omind_ind = [ ...
find(refgr_omind == 1); find(refgr_omind == 2); ...
find(refgr_omind == 3); find(refgr_omind == 4) ];
for ii_g = 1:numel(orients)
orients(ii_g) = gtGrainAllblobsFilterOrder(orients(ii_g), refgr_omind_ind);
end
end
function [u_size, v_size] = get_stack_UV_size(img_sizes, p, conf)
max_img_sizes = [max(img_sizes(:, 1)), max(img_sizes(:, 2))];
......@@ -496,4 +447,3 @@ end
zUtil_Deformation/gt6DGetMatchingReflections.m 0 → 100644
+ 28
− 0
View file @ 873df95b
function match = gt6DGetMatchingReflections(ref_or, ors, angular_toll)
if (~exist('angular_toll', 'var') || isempty(angular_toll))
angular_toll = 1; % degrees
end
ref_ws = ref_or.allblobs.omega;
ref_ns = ref_or.allblobs.eta;
ref_Lfactors = ref_or.allblobs.lorentzfac;
w_toll = angular_toll * ref_Lfactors;
n_toll = angular_toll;
ors_allblobs = [ors(:).allblobs];
ors_ws = cat(2, ors_allblobs(:).omega);
ors_ns = cat(2, ors_allblobs(:).eta);
ones_num_ors = ones(numel(ors), 1);
ref_ws = ref_ws(:, ones_num_ors);
ref_ns = ref_ns(:, ones_num_ors);
w_toll = w_toll(:, ones_num_ors);
w_diff = abs(mod(ref_ws - ors_ws + 180, 360) - 180);
n_diff = abs(mod(ref_ns - ors_ns + 180, 360) - 180);
match = (w_diff < w_toll) & (n_diff < n_toll);
end
zUtil_Deformation/gt6DMergeRealSpaceVolumes.m 0 → 100644
+ 39
− 0
View file @ 873df95b
function [gr_center_pix, vol_size_pix, bbox_pix] = gt6DMergeRealSpaceVolumes(grs, det_index, oversize)
% num_ors = numel(grs);
% vol_sizes = zeros(num_ors, 3);
% vol_centers = zeros(num_ors, 3);
% for ii = 1:num_ors
% vol_sizes(ii, :) = [ ...
% grs(ii).proj(det_index).vol_size_y, ...
% grs(ii).proj(det_index).vol_size_x, ...
% grs(ii).proj(det_index).vol_size_z ];
% vol_centers(ii, :) = grs(ii).proj(det_index).centerpix;
% end
projs = arrayfun(@(x){x.proj(det_index)}, grs);
projs = [projs{:}];
vol_sizes = [cat(1, projs.vol_size_y), cat(1, projs.vol_size_x), cat(1, projs.vol_size_z)];
vol_centers = cat(1, projs.centerpix);
vol_half_sizes = vol_sizes / 2;
min_corners = vol_centers - vol_half_sizes;
max_corners = vol_centers + vol_half_sizes;
full_vol_min_corners = min(min_corners, [], 1);
full_vol_max_corners = max(max_corners, [], 1);
gr_center_pix = (full_vol_max_corners + full_vol_min_corners) / 2;
vol_size_pix = full_vol_max_corners - full_vol_min_corners;
if (exist('oversize', 'var'))
vol_size_pix = vol_size_pix * oversize;
end
vol_size_pix = ceil(vol_size_pix);
vol_half_size_pix = vol_size_pix / 2;
bbox_pix = [gr_center_pix - vol_half_size_pix, gr_center_pix + vol_half_size_pix];
end
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment