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

6D-reconstruction: added NW shape function reconstruction for near-field data

parent bf1d0e64
No related branches found
No related tags found
No related merge requests found
......@@ -15,6 +15,7 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% Different types of shape functions:
% 'none' : no shape functions are used
% 'w' : Only omega is considered
% 'nw' : (eta, omega) shape functions
% 'nw2uvw' : (eta, omega) shape functions -> UVW shape functions
% 'uvw' : Using FED's projector for UVW shape functions
shape_functions_type = 'none';
......@@ -52,6 +53,10 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
shape_funcs = gtDefShapeFunctionsCreateW(sampler);
algo_params = self.get_algorithm_params_shape_funcs(sampler, shape_funcs, self.shape_functions_type);
shape_funcs = {};
case 'nw'
shape_funcs = gtDefShapeFunctionsCreateNW(sampler);
algo_params = self.get_algorithm_params_shape_funcs(sampler, shape_funcs, self.shape_functions_type);
shape_funcs = {};
case 'nw2uvw'
shape_funcs = gtDefShapeFunctionsCreateNW(sampler, 2, 'data_type', 'single');
shape_funcs = gtDefShapeFunctionsNW2UVW(sampler, shape_funcs);
......@@ -458,6 +463,9 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
case 'w'
[geometries, offsets] = self.compute_proj_geom_shapefuncs_w( ...
bl, ref_gr, ors_allblobs, shape_funcs, lims_blobs_w);
case 'nw'
[geometries, offsets] = self.compute_proj_geom_shapefuncs_nw( ...
bl, ref_gr, ors_allblobs, shape_funcs, lims_blobs_w);
case {'nw2uvw', 'uvw'}
geometries = self.compute_proj_geom_numinterp(bl, ref_gr, ors_allblobs);
offsets = self.make_proj_geom_shape_funcs_uvw(...
......@@ -763,11 +771,11 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
num_orients = numel(ors_allblobs);
bb_bls_uv_oi = repmat(bb_bls_uv, [num_orients, 1]);
dvecsam = cat(1, ors_allblobs(:).dvecsam);
dvec_sam = cat(1, ors_allblobs(:).dvecsam);
rot_w_l2s = cat(3, ors_allblobs(:).srot);
proj_geom = gtGeoProjForReconstruction(...
dvecsam, rot_w_l2s, ref_gr.center, bb_bls_uv_oi, [], ...
dvec_sam, rot_w_l2s, ref_gr.center, bb_bls_uv_oi, [], ...
self.parameters.detgeo(self.det_index), ...
self.parameters.labgeo, ...
self.parameters.samgeo, ...
......@@ -883,16 +891,24 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% Assign offsets
blob_offsets = cell(1, num_blobs);
proj_coeffs = cell(1, num_blobs);
ws = cell(1, num_blobs);
for ii_b = 1:num_blobs
proj_coeffs{ii_b} = reshape(sf(ii_b).intm, 1, []);
valid_coeffs = proj_coeffs{ii_b} > eps('single');
proj_coeffs{ii_b} = proj_coeffs{ii_b}(valid_coeffs);
ws{ii_b} = sf(ii_b).bbwim(1):sf(ii_b).bbwim(2);
blob_offsets{ii_b} = ii_b * ones(1, sf(ii_b).bbsize);
ws{ii_b} = ws{ii_b}(valid_coeffs);
blob_offsets{ii_b} = ii_b * ones(1, sum(valid_coeffs));
end
ws = [ws{:}];
sino_offsets = 1:numel(ws);
blob_offsets = [blob_offsets{:}];
proj_offsets = ws - lims_blobs_w(blob_offsets, 1)' + 1;
proj_coeffs = cat(1, sf.intm)';
proj_coeffs = [proj_coeffs{:}];
offsets{ii_o}.('sino_offsets') = sino_offsets;
offsets{ii_o}.('blob_offsets') = blob_offsets;
......@@ -906,12 +922,99 @@ classdef Gt6DReconstructionAlgorithmFactory < handle
% Computing geometries
rot_w_l2s = gtMathsRotationTensor(ws' .* om_step, rot_comp_w);
dveclab = ors_allblobs(ii_o).dvec(blob_offsets, :);
dvecsam = gtGeoLab2Sam(dveclab, rot_w_l2s, ...
dvec_lab = ors_allblobs(ii_o).dvec(blob_offsets, :);
dvec_sam = gtGeoLab2Sam(dvec_lab, rot_w_l2s, ...
self.parameters.labgeo, self.parameters.samgeo, true);
proj_geom = gtGeoProjForReconstruction(...
dvec_sam, rot_w_l2s, ref_gr.center, ...
bb_bls_uv(blob_offsets, :), [], ...
self.parameters.detgeo(self.det_index), ...
self.parameters.labgeo, ...
self.parameters.samgeo, ...
self.parameters.recgeo(self.det_index), ...
'ASTRA_grain');
if (self.volume_downscaling > 1)
proj_geom(:, 4:12) = proj_geom(:, 4:12) / self.volume_downscaling;
end
geometries{ii_o} = double(proj_geom);
end
end
function [geometries, offsets] = compute_proj_geom_shapefuncs_nw(self, bl, ref_gr, ors_allblobs, sf_nw, lims_blobs_w)
bls_bbuim = cat(1, bl.bbuim);
bls_bbvim = cat(1, bl.bbvim);
bls_bbsize = cat(1, bl.bbsize);
bb_bls_uv = [bls_bbuim(:, 1), bls_bbvim(:, 1), bls_bbsize(:, [1 2])];
num_blobs = numel(bl);
num_orients = numel(ors_allblobs);
rot_comp_w = gtMathsRotationMatrixComp(self.parameters.labgeo.rotdir', 'col');
om_step = gtGetOmegaStepDeg(self.parameters);
geometries = cell(num_orients, 1);
offsets = cell(num_orients, 1);
for ii_o = 1:num_orients
sf = sf_nw{ii_o};
bl_bbsizes = cat(1, sf.bbsize);
% Assign offsets
blob_offsets = cell(1, num_blobs);
proj_coeffs = cell(1, num_blobs);
ns = cell(1, num_blobs);
ws = cell(1, num_blobs);
for ii_b = 1:num_blobs
ns{ii_b} = linspace(sf(ii_b).bbnim(1), sf(ii_b).bbnim(2), bl_bbsizes(ii_b, 1));
ws{ii_b} = sf(ii_b).bbwim(1):sf(ii_b).bbwim(2);
ns{ii_b} = reshape(ns{ii_b}, [], 1);
ns{ii_b} = ns{ii_b}(:, ones(sf(ii_b).bbsize(2), 1));
ns{ii_b} = reshape(ns{ii_b}, [], 1);
ws{ii_b} = ws{ii_b}(ones(sf(ii_b).bbsize(1), 1), :);
ws{ii_b} = reshape(ws{ii_b}, [], 1);
proj_coeffs{ii_b} = reshape(sf(ii_b).intm, 1, []);
valid_coeffs = proj_coeffs{ii_b} > eps('single');
ns{ii_b} = ns{ii_b}(valid_coeffs);
ws{ii_b} = ws{ii_b}(valid_coeffs);
proj_coeffs{ii_b} = proj_coeffs{ii_b}(valid_coeffs);
blob_offsets{ii_b} = ii_b * ones(1, sum(valid_coeffs));
end
ws = cat(1, ws{:});
ns = cat(1, ns{:});
sino_offsets = 1:numel(ws);
blob_offsets = [blob_offsets{:}];
proj_offsets = (ws - lims_blobs_w(blob_offsets, 1))' + 1;
proj_coeffs = [proj_coeffs{:}];
offsets{ii_o}.('sino_offsets') = sino_offsets;
offsets{ii_o}.('blob_offsets') = blob_offsets;
offsets{ii_o}.('proj_offsets') = proj_offsets;
offsets{ii_o}.('proj_coeffs') = proj_coeffs;
offsets{ii_o}.proj = struct( ...
'sino_offset', sino_offsets, 'blob_offset', blob_offsets, ...
'proj_offsets', proj_offsets, 'proj_coeffs', proj_coeffs );
% Computing geometries
rot_w_l2s = gtMathsRotationTensor(ws .* om_step, rot_comp_w);
t = ors_allblobs(ii_o).theta(blob_offsets);
pl_lab = gtGeoPlLabFromThetaEta(t, ns, self.parameters.labgeo);
dvec_lab = gtFedPredictDiffVecMultiple(pl_lab', self.parameters.labgeo.beamdir')';
dvec_sam = gtGeoLab2Sam(dvec_lab, rot_w_l2s, ...
self.parameters.labgeo, self.parameters.samgeo, true);
proj_geom = gtGeoProjForReconstruction(...
dvecsam, rot_w_l2s, ref_gr.center, ...
dvec_sam, rot_w_l2s, ref_gr.center, ...
bb_bls_uv(blob_offsets, :), [], ...
self.parameters.detgeo(self.det_index), ...
self.parameters.labgeo, ...
......
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