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

6D-Reconstruction: added initial twin reconstruction 


Signed-off-by: default avatarNicola Vigano <nicola.vigano@esrf.fr>
parent 2648a93b
No related branches found
No related tags found
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
5_reconstruction/gtReconstruct6DLaunchAlgorithm.m
+ 1
− 1
View file @ d2422b5a
......@@ -5,7 +5,7 @@ function [algo, good_or] = gtReconstruct6DLaunchAlgorithm(sampler, rec_opts, par
rec_factory = Gt6DReconstructionAlgorithmFactory(parameters, ...
'volume_downscaling', rec_opts.volume_downscaling, ...
'rspace_oversize', rec_opts.rspace_oversize );
[algo, good_or] = rec_factory.getSubBlobReconstructionAlgo(sampler, rec_opts.num_interp);
[algo, good_or] = rec_factory.getReconstructionAlgo(sampler, rec_opts.num_interp);
mem_cons = [ ...
algo.get_peak_memory_consumption(rec_opts.algorithm), ...
......
5_reconstruction/gtReconstructGrainTwinCluster.m 0 → 100644
+ 111
− 0
View file @ d2422b5a
function gtReconstructGrainTwinCluster(grain_ids, phase_id, parameters, varargin)
% gtReconstructGrainOrientation 6D reconstructions on a GPU machine
% gtAstraReconstructGrain(grainIDs, phaseID, [parameters])
% -------------------------------------------------------
if (~exist('parameters', 'var') || isempty(parameters))
parameters = gtLoadParameters();
end
parameters.fsim.check_spot = false;
rec_opts = gtReconstruct6DGetParamenters(parameters);
phase_dir = fullfile(parameters.acq.dir, '4_grains', ...
sprintf('phase_%02d', phase_id));
grains_str_ids = sprintf('_%04d', grain_ids);
grain_cluster_file = fullfile(phase_dir, ...
sprintf('grain_twin_cluster%s.mat', grains_str_ids));
fprintf('Loading the cluster file..')
grs = load(grain_cluster_file);
grs = grs.samp_ors;
fprintf('\b\b: Done.\n')
num_grains = numel(grs);
for ii_g = 1:num_grains
bb_gvdm = cat(1, grs(ii_g).bb_ors(:).R_vector)';
sampler(ii_g) = GtOrientationSampling(parameters, grs(ii_g)); %#ok<AGROW>
sampler(ii_g).make_even_simple_grid('cubic', rec_opts.grid_edge, bb_gvdm, 1)
if (rec_opts.super_sampling > 1)
sampler(ii_g).make_supersampling_simple_grid([1 2 3], rec_opts.super_sampling);
end
end
[algo, good_or] = gtReconstruct6DLaunchAlgorithm(sampler, rec_opts, parameters);
vols = algo.getCurrentSolution();
% Restoring initial volume size (depending on the rounding)
if (rec_opts.volume_downscaling > 1)
fprintf('Expanding volumes..')
c = tic();
vols = gtMathsUpsampleVolume(vols, rec_opts.volume_downscaling);
fprintf('\b\b: Done (%f seconds).\n', toc(c))
end
ranges = zeros(num_grains, 2);
ranges(1, 1) = 1;
for ii_g = 2:num_grains
ranges(ii_g, 1) = ranges(ii_g-1) + sampler(ii_g-1).get_total_num_orientations();
end
ranges(1:end-1, 2) = ranges(2:end, 1) - 1;
ranges(end, 2) = numel(vols);
for ii_g = num_grains:-1:1
or_sizes = sampler(ii_g).get_orientation_sampling_size();
gr_vols = vols(ranges(ii_g, 1):ranges(ii_g, 2));
gr_good_or = good_or(ranges(ii_g, 1):ranges(ii_g, 2));
[avg_R_vecs, avg_R_vecs_int] = sampler(ii_g).getAverageOrientations(gr_vols, gr_good_or);
avg_R_vec = sampler(ii_g).getAverageOrientation(gr_vols, gr_good_or);
s_g_odf = reshape(sampler(ii_g).getODF(gr_vols, gr_good_or), or_sizes{1});
vol_size = size(avg_R_vecs_int);
shift = gtFwdSimComputeVolumeShifts(grs(ii_g).proj, parameters, vol_size);
[kam, gam(ii_g)] = gtDefComputeKernelAverageMisorientation(avg_R_vecs, avg_R_vecs_int);
[igm, gos(ii_g)] = gtDefComputeIntraGranularMisorientation(avg_R_vecs, avg_R_vecs_int, 'R_vector', grs(ii_g).R_vector);
ODF6D(ii_g) = struct( ...
'voxels_avg_R_vectors', {avg_R_vecs}, ...
'intensity', {avg_R_vecs_int}, ...
'shift', {shift}, ...
'R_vectors', {sampler(ii_g).get_R_vectors()}, ...
'single_grain_ODF', {s_g_odf}, ...
'single_grain_avg_R_vector', {avg_R_vec}, ...
'kernel_average_misorientation', {kam}, ...
'intra_granular_misorientation', {igm} );
end
cluster_rec_file = fullfile(phase_dir, ...
sprintf('grain_twin_cluster_details%s.mat', grains_str_ids));
% Saving and cleaning at the same time
if (exist(cluster_rec_file, 'file'))
cl_rec = load(cluster_rec_file);
cl_rec.ODF6D = ODF6D;
else
cl_rec = struct('ODF6D', ODF6D); %#ok<NASGU>
end
save(cluster_rec_file, '-struct', 'cl_rec', '-v7.3');
if (algo.verbose)
clear('ODF6D')
[proj_blobs, proj_spots] = algo.getProjectionOfCurrentSolution();
for ii_g = num_grains:-1:1
gr_vols = vols(ranges(ii_g, 1):ranges(ii_g, 2));
ODF6D(ii_g) = struct( ...
'orientation_volumes', {gr_vols}, ...
'fwd_projected_blobs', {proj_blobs}, ...
'fwd_projected_spots', {proj_spots}, ...
'grain_average_misorientation', {gam(ii_g)}, ...
'grain_orientation_spread', {gos(ii_g)}, ...
'compute_statistics', {algo.get_statistics()} );
end
grain_full_details_file = fullfile(phase_dir, ...
sprintf('grain_twin_cluster_full_details%s.mat', grains_str_ids));
save(grain_full_details_file, 'ODF6D', '-v7.3');
end
end
zUtil_Deformation/Gt6DReconstructionAlgorithmFactory.m
+ 100
− 16
View file @ d2422b5a
......@@ -21,7 +21,15 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
self.parameters = parameters;
end
function [algo, good_orients, blobs] = getSubBlobReconstructionAlgo(self, sampler, num_interp, varargin)
function [algo, good_orients, blobs] = getReconstructionAlgo(self, sampler, num_interp, varargin)
if (numel(sampler) == 1)
[algo, good_orients, blobs] = self.getGrainReconstructionAlgo(sampler, num_interp, varargin{:});
else
[algo, good_orients, blobs] = self.getTwinReconstructionAlgo(sampler, num_interp, varargin{:});
end
end
function [algo, good_orients, blobs] = getGrainReconstructionAlgo(self, sampler, num_interp, varargin)
algo_params = self.get_algorithm_params(sampler, num_interp);
% Build Empty volumes
......@@ -31,19 +39,6 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
if (self.volume_downscaling > 1)
volume_size = ceil(volume_size / self.volume_downscaling);
num_geoms = numel(algo_params.geometries);
for ii = 1:num_geoms
algo_params.geometries{ii}(:, 4:12) = algo_params.geometries{ii}(:, 4:12) / self.volume_downscaling;
end
if (~isempty(algo_params.geometries_ss))
for ii = 1:num_geoms
geom_ss = algo_params.geometries_ss{ii};
for ii_ss = 1:numel(geom_ss)
geom_ss{ii_ss}(:, 4:12) = geom_ss{ii_ss}(:, 4:12) / self.volume_downscaling;
end
algo_params.geometries_ss{ii} = geom_ss;
end
end
end
good_orients = algo_params.good_orients;
......@@ -60,6 +55,78 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
'psf', algo_params.psf, ...
varargin{:} );
end
function [algo, good_orients, blobs] = getTwinReconstructionAlgo(self, sampler, num_interp, varargin)
num_grains = numel(sampler);
for ii_g = 1:num_grains
algo_params(ii_g) = self.get_algorithm_params(sampler(ii_g), num_interp); %#ok<AGROW>
end
% Build Empty volumes
ref_gr = sampler(1).get_reference_grain();
spacing = mean([ref_gr.proj.vol_size_y, ref_gr.proj.vol_size_x, ref_gr.proj.vol_size_z]) * (self.rspace_oversize - 1);
volume_size = round([ref_gr.proj.vol_size_y, ref_gr.proj.vol_size_x, ref_gr.proj.vol_size_z] + spacing);
if (self.volume_downscaling > 1)
volume_size = ceil(volume_size / self.volume_downscaling);
end
blobs = cat(1, algo_params(:).blobs);
psf = cat(1, algo_params(:).psf);
good_orients = cat(1, algo_params(:).good_orients);
geometries = cat(1, algo_params(:).geometries);
geometries_ss = cat(1, algo_params(:).geometries_ss);
offsets = algo_params(1).offsets;
offsets_ss = algo_params(1).offsets_ss;
base_offset_shift = numel(algo_params(1).blobs);
for ii_g = 2:num_grains
temp_gr = sampler(ii_g).get_reference_grain();
shared = temp_gr.proj.shared_parent(sampler(ii_g).selected);
orient_offsets = algo_params(ii_g).offsets;
for ii_o = 1:numel(orient_offsets)
to_be_summed = ~shared(orient_offsets{ii_o}.blob_offsets);
orient_offsets{ii_o}.blob_offsets(to_be_summed) ...
= orient_offsets{ii_o}.blob_offsets(to_be_summed) + base_offset_shift;
for ii_p = 1:numel(orient_offsets{ii_o}.proj)
ii_b = orient_offsets{ii_o}.proj(ii_p).blob_offset;
if (~temp_gr.proj.shared_parent(ii_b))
orient_offsets{ii_o}.proj(ii_p).blob_offset ...
= ii_b + base_offset_shift;
end
end
end
offsets = cat(1, offsets, orient_offsets);
if (~isempty(geometries_ss))
orient_offsets_ss = algo_params(ii_g).offsets_ss;
for ii_o = 1:numel(orient_offsets_ss)
for ii_ss = 1:numel(orient_offsets_ss{ii_o})
to_be_summed = ~shared(orient_offsets{ii_o}{ii_ss}.blob_offsets);
orient_offsets_ss{ii_o}{ii_ss}.blob_offsets(to_be_summed) ...
= orient_offsets_ss{ii_o}{ii_ss}.blob_offsets(to_be_summed) + base_offset_shift;
end
end
offsets_ss = cat(1, offsets_ss, orient_offsets_ss);
end
base_offset_shift = base_offset_shift + numel(algo_params(ii_g).blobs);
end
algo = Gt6DBlobReconstructor(volume_size, blobs, ...
'data_type', self.data_type, ...
'geometries', geometries, ...
'offsets', offsets, ...
'ss_geometries', geometries_ss, ...
'ss_offsets', offsets_ss, ...
'use_astra_projectors', true, ...
'volume_ss', self.volume_downscaling, ...
'psf', psf, ...
varargin{:} );
end
end
methods (Access = protected)
......@@ -246,6 +313,21 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
offsets_ss = {};
end
if (self.volume_downscaling > 1)
for ii = 1:num_geoms
geometries{ii}(:, 4:12) = geometries{ii}(:, 4:12) / self.volume_downscaling;
end
if (~isempty(geometries_ss))
for ii = 1:num_geoms
geom_ss = geometries_ss{ii};
for ii_ss = 1:numel(geom_ss)
geom_ss{ii_ss}(:, 4:12) = geom_ss{ii_ss}(:, 4:12) / self.volume_downscaling;
end
geometries_ss{ii} = geom_ss;
end
end
end
if (isfield(bl, 'psf'))
psf = arrayfun(@(x){ permute(x.psf, [1 3 2]) }, bl);
elseif (isfield(sampler.ref_gr.proj, 'psf'))
......@@ -261,7 +343,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
'geometries_ss', {geometries_ss}, ...
'offsets_ss', {offsets_ss}, ...
'psf', {psf}, ...
'good_orients', {good_orients} );
'good_orients', {good_orients'} );
end
function [sub_blob_slices, proj_coeffs] = get_sub_blobs(self, blobs, slices_interp, padding)
......@@ -441,6 +523,8 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
else
geometries = {grains_props(:).proj_geom};
end
geometries = reshape(geometries, [], 1);
num_geoms = numel(geometries);
w_int_tab(:, :, 1) = floor(w_tab);
......@@ -451,7 +535,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
slice_pos_in_blob = w_int_tab - repmat(extreemes_blobs_w(:, 1), [1 num_geoms 2]) + 1;
% Assign offsets
offsets = cell(1, num_geoms);
offsets = cell(num_geoms, 1);
for ii = 1:num_geoms
valid_proj_dirs = ...
(w_tab(:, ii) >= abs_lims(:, 1)) ...
......
zUtil_Deformation/gt6DCreateProjDataFromTwinnedGrain.m 0 → 100644
+ 606
− 0
View file @ d2422b5a
This diff is collapsed.
zUtil_Deformation/test_generation/do6D.m
+ 2
− 2
View file @ d2422b5a
......@@ -186,7 +186,7 @@ switch(mode)
recon_algo_factory = Gt6DReconstructionAlgorithmFactory();
% [reconstructor, blobs] = recon_algo_factory.getBlobReconstructionAlgo( ...
% test_data.bl, blobs_on_det, sampler, '2D');
[reconstructor, blobs] = recon_algo_factory.getSubBlobReconstructionAlgo( ...
[reconstructor, blobs] = recon_algo_factory.getReconstructionAlgo( ...
test_data.bl, blobs_on_det, sampler, '2D', num_interp);
case 'testPeterSubBlob3D'
cd('~/data/testdata/')
......@@ -230,7 +230,7 @@ switch(mode)
end
recon_algo_factory = Gt6DReconstructionAlgorithmFactory();
[reconstructor, blobs] = recon_algo_factory.getSubBlobReconstructionAlgo( ...
[reconstructor, blobs] = recon_algo_factory.getReconstructionAlgo( ...
test_data.bl, blobs_on_det, sampler, '3D', num_interp);
case 'peterDataset'
cd('/mntdirect/_data_id19_degisher/mi1026/id18/AlLi2b_vert_ ')
......
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