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

6D-reconstruction: enabling ospace-super-sampling with the use of W-shape-functions 


Signed-off-by: default avatarNicola Vigano <nicola.vigano@esrf.fr>
parent 86199636
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
zUtil_Deformation/Gt6DBlobReconstructor.m
+ 22
− 6
View file @ 8d42d06e
......@@ -289,10 +289,18 @@ classdef Gt6DBlobReconstructor < Gt6DVolumeToBlobProjector
num_bytes = 0;
float_bytes = GtBenchmarks.getSizeVariable(zeros(1, 1, 'single'));
num_geoms = self.get_number_geometries();
for n = 1:num_geoms
sinogram_size = [ ...
self.proj_size(1), size(self.geometries{n}, 1), self.proj_size(2)];
num_bytes = num_bytes + float_bytes * prod(sinogram_size);
if (isempty(self.ss_geometries))
for n = 1:num_geoms
sinogram_size = [ ...
self.proj_size(1), size(self.geometries{n}, 1), self.proj_size(2)];
num_bytes = num_bytes + float_bytes * prod(sinogram_size);
end
else
for n = 1:num_geoms
sinogram_size = [ ...
self.proj_size(1), size(self.ss_geometries{n}, 1), self.proj_size(2)];
num_bytes = num_bytes + float_bytes * prod(sinogram_size);
end
end
end
end
......@@ -514,8 +522,16 @@ classdef Gt6DBlobReconstructor < Gt6DVolumeToBlobProjector
fprintf(' + Backprojecting ones sinos (fake bproj)..');
c = tic();
num_geoms = self.get_number_geometries();
for ii = 1:num_geoms
self.bwd_weights{ii} = size(self.geometries{ii}, 1);
if (~isempty(self.geometries))
for ii = 1:num_geoms
self.bwd_weights{ii} = size(self.geometries{ii}, 1);
end
else
% If this happens we are using oversampling with
% shape-functions
for ii = 1:num_geoms
self.bwd_weights{ii} = numel(self.blobs);
end
end
fprintf('\b\b (%2.1f s)\n', toc(c));
end
......
zUtil_Deformation/Gt6DReconstructionAlgorithmFactory.m
+ 38
− 6
View file @ 8d42d06e
......@@ -447,11 +447,18 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
sel_allb_idxs = sampler.ondet(sampler.included(sel_bls));
omegas = ref_gr.allblobs.omega(sel_allb_idxs);
orients = sampler.get_orientations();
orients = [orients{:}];
[orients, orients_ss] = sampler.get_orientations();
num_ospace_oversampling = prod(sampler.ss_factor);
if (num_ospace_oversampling == 1)
ors_struct = [orients{:}];
else
ors_struct = cat(4, orients_ss{:});
ors_struct = [ors_struct{:}];
end
% Array of structures that contain all the info relative to
% each orientation, for each blob
grains_props = self.get_selected_props(orients, sel_bls);
grains_props = self.get_selected_props(ors_struct, sel_bls);
[w_tab, ref_ws] = self.get_shape_funcs_w_tab(shape_funcs, omegas);
[lims_blobs_w, lims_projs_w, paddings, ext_blobs_w_orig] = self.get_blob_lims(w_tab, bl, ref_ws);
......@@ -494,14 +501,39 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
psf = {};
end
tot_orient = sampler.get_total_num_orientations();
if (num_ospace_oversampling == 1)
geometries_ss = {};
offsets_ss = {};
else
geometries = reshape(geometries, num_ospace_oversampling, []);
offsets = reshape(offsets, num_ospace_oversampling, []);
geometries_ss = cell(tot_orient, 1);
offsets_ss = cell(tot_orient, 1);
for ii = 1:tot_orient
geometries_ss{ii} = geometries(:, ii);
% merging all the orientation-space super-sampling
% projection coefficients into one structure per
% orientation. This is needed in order to fwd/bwd-project
% all the different sub-orientations together.
offsets_ss{ii} = self.merge_offsets(offsets(:, ii));
end
geometries = {};
offsets = {};
end
algo_params = struct( ...
'blobs', {blobs}, ...
'geometries', {geometries}, ...
'offsets', {offsets}, ...
'geometries_ss', {{}}, ...
'offsets_ss', {{}}, ...
'geometries_ss', {geometries_ss}, ...
'offsets_ss', {offsets_ss}, ...
'psf', {psf}, ...
'good_orients', {true(1, numel(geometries))} );
'blobs_w_lims', lims_blobs_w, ...
'good_orients', {true(1, tot_orient)} );
end
function [sub_blob_slices, proj_coeffs] = get_sub_blobs(self, blobs, slices_interp, padding)
......
zUtil_Deformation/GtOrientationSampling.m
+ 20
− 8
View file @ 8d42d06e
......@@ -197,15 +197,24 @@ classdef GtOrientationSampling < handle
gr = self.lattice(ii_o).gr{ii_g};
beg_pos = gr.R_vector - half_voxel_size + half_steps;
end_pos = gr.R_vector + half_voxel_size - half_steps;
gr_ss = {};
for ii_x = beg_pos(1):steps(1):end_pos(1)
for ii_y = beg_pos(2):steps(2):end_pos(2)
for ii_z = beg_pos(3):steps(3):end_pos(3)
gr_ss{end+1} = struct( ...
'R_vector', {[ii_x, ii_y, ii_z]}, ...
range_x = beg_pos(1):steps(1):end_pos(1);
range_y = beg_pos(2):steps(2):end_pos(2);
range_z = beg_pos(3):steps(3):end_pos(3);
gr_ss = cell(self.ss_factor);
for ii_z = 1:self.ss_factor(3)
for ii_y = 1:self.ss_factor(2)
r_vecs = [ ...
range_x', ...
range_y(ii_y) * ones(self.ss_factor(1), 1), ...
range_z(ii_z) * ones(self.ss_factor(1), 1)];
for ii_x = 1:self.ss_factor(1)
gr_ss{ii_x, ii_y, ii_z} = struct( ...
'R_vector', {r_vecs(ii_x, :)}, ...
'center', {gr.center}, ...
'phaseid', {gr.phaseid}, ...
'allblobs', {[]} ); %#ok<AGROW>
'allblobs', {[]} );
end
end
end
......@@ -248,7 +257,10 @@ classdef GtOrientationSampling < handle
num_lattices = numel(self.lattice);
or_size = cell(num_lattices, 1);
for ii = 1:num_lattices
or_size{ii} = size(self.lattice(ii).gr);
or_size{ii} = [ ...
size(self.lattice(ii).gr, 1), ...
size(self.lattice(ii).gr, 2), ...
size(self.lattice(ii).gr, 3) ];
end
end
......
zUtil_Deformation/gtDefShapeFunctionsCreateW.m
+ 53
− 17
View file @ 8d42d06e
......@@ -7,30 +7,49 @@ function sf = gtDefShapeFunctionsCreateW(sampler)
omstep = gtGetOmegaStepDeg(p, sampler.detector_index);
omstep_rod = tand(omstep / 2);
sub_space_size = sampler.stats.sampling.gaps;
% Let's first compute the W extent
half_rspace_sizes = sub_space_size / 2;
ref_incs = sampler.ondet(sampler.included);
ref_inds = sampler.selected;
num_blobs = numel(find(ref_inds));
ref_gr = sampler.get_reference_grain();
pl_samd = ref_gr.allblobs.plorig(ref_incs(ref_inds), :);
if (isfield(ref_gr.allblobs, 'plcry'))
pl_cry = ref_gr.allblobs.plcry(ref_incs(ref_inds), :);
else
pl_samd = ref_gr.allblobs.plorig(ref_incs(ref_inds), :);
g = gtMathsRod2OriMat(ref_gr.R_vector');
pl_cry = gtVectorLab2Cryst(pl_samd, g)';
end
% pl_labd = gtGeoSam2Lab(pl_samd, eye(3), self.parameters.labgeo, self.parameters.samgeo, true, false)';
sub_space_size = sampler.stats.sampling.gaps;
ospace_samp_size = sampler.get_orientation_sampling_size();
ospace_samp_size = ospace_samp_size{1};
num_ospace_oversampling = prod(sampler.ss_factor);
if (num_ospace_oversampling == 1)
ospace_samp_lattice = sampler.lattice.gr;
else
ospace_samp_lattice = cat(4, sampler.lattice_ss.gr{:});
ospace_samp_lattice = reshape(ospace_samp_lattice, [sampler.ss_factor, ospace_samp_size]);
ospace_samp_lattice = permute(ospace_samp_lattice, [1 4 2 5 3 6]);
ospace_samp_size = ospace_samp_size .* sampler.ss_factor;
sub_space_size = sub_space_size ./ sampler.ss_factor;
ospace_samp_lattice = reshape(ospace_samp_lattice, ospace_samp_size);
end
g = gtMathsRod2OriMat(ref_gr.R_vector');
pl_cry = gtVectorLab2Cryst(pl_samd, g)';
% pl_labd = gtGeoSam2Lab(pl_samd, eye(3), self.parameters.labgeo, self.parameters.samgeo, true, false)';
% Let's first compute the W extent
half_rspace_sizes = sub_space_size / 2;
fprintf('Computing bounds of W shape functions: ')
c = tic();
tot_orientations = numel(sampler.lattice.gr);
tot_orientations = numel(ospace_samp_lattice);
oris_bounds_size = [size(sampler.lattice.gr, 1), size(sampler.lattice.gr, 2), size(sampler.lattice.gr, 3)] + 1;
oris_lims_min = sampler.lattice.gr{1, 1, 1}.R_vector - half_rspace_sizes;
oris_lims_max = sampler.lattice.gr{end, end, end}.R_vector + half_rspace_sizes;
oris_bounds_size = ospace_samp_size + 1;
oris_lims_min = ospace_samp_lattice{1, 1, 1}.R_vector - half_rspace_sizes;
oris_lims_max = ospace_samp_lattice{end, end, end}.R_vector + half_rspace_sizes;
oris_steps_x = linspace(oris_lims_min(1), oris_lims_max(1), oris_bounds_size(1));
oris_steps_y = linspace(oris_lims_min(2), oris_lims_max(2), oris_bounds_size(2));
......@@ -131,7 +150,7 @@ function sf = gtDefShapeFunctionsCreateW(sampler)
num_chars = fprintf('%05d/%05d', ii_g, tot_orientations);
end
gr = sampler.lattice.gr{ii_g};
or = ospace_samp_lattice{ii_g};
% Extracting ospace boundaries for the given voxel
w_tab_int = w_tabs_int(:, o_x(ii_g):o_x(ii_g)+1, o_y(ii_g):o_y(ii_g)+1, o_z(ii_g):o_z(ii_g)+1);
......@@ -145,7 +164,7 @@ function sf = gtDefShapeFunctionsCreateW(sampler)
% that we want to integrate.
r_vecs = [pos_x_exp; pos_y_exp; pos_z_exp];
r_vecs = reshape(r_vecs, 3, []);
r_vecs = gr.R_vector(ones_tot_sp, :)' + r_vecs;
r_vecs = or.R_vector(ones_tot_sp, :)' + r_vecs;
gs = gtMathsRod2OriMat(r_vecs);
% Unrolling and transposing the matrices at the same time
......@@ -169,10 +188,10 @@ function sf = gtDefShapeFunctionsCreateW(sampler)
D = Ac .^ 2 + As .^ 2;
% Precomputing useful values, like sinthetas
ominds = gr.allblobs.omind(ref_inds);
ominds = or.allblobs.omind(ref_inds);
ssp = ((ominds == 1) | (ominds == 2));
ss = ssp - ~ssp;
sinths = ss .* gr.allblobs.sintheta(ref_inds);
sinths = ss .* or.allblobs.sintheta(ref_inds);
sinths = sinths(:, ones(tot_sampling_points, 1));
CC = Acc + sinths;
......@@ -235,4 +254,21 @@ function sf = gtDefShapeFunctionsCreateW(sampler)
end
end
fprintf('Done in %g seconds.\n', toc(c))
if (num_ospace_oversampling > 1)
fprintf('Re-grouping oversampled W shape functions: ')
c = tic();
ospace_samp_size = sampler.get_orientation_sampling_size();
ospace_samp_size = ospace_samp_size{1};
tot_orient = sampler.get_total_num_orientations();
regroup = [ ...
sampler.ss_factor(1), ospace_samp_size(1), ...
sampler.ss_factor(2), ospace_samp_size(2), ...
sampler.ss_factor(3), ospace_samp_size(3) ];
sf = reshape(sf, regroup);
sf = permute(sf, [1 3 5 2 4 6]);
sf = reshape(sf, num_ospace_oversampling, tot_orient);
fprintf('Done in %g seconds.\n', toc(c))
end
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