Skip to content
Snippets Groups Projects
Commit b38d93a0 authored by Wolfgang Ludwig's avatar Wolfgang Ludwig
Browse files

enable OAR creation of extended grains 


Signed-off-by: default avatarWolfgang Ludwig <wolfgang.ludwig@esrf.fr>
parent 3e0a59dd
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
zUtil_Deformation/gt6DCreateProjDataFromExtendedGrain.m
+ 273
− 258
View file @ b38d93a0
function [refor, estim_space_bbox_pix, estim_orient_bbox_rod] = gt6DCreateProjDataFromExtendedGrain(gr_id, phase_id, varargin)
% FUNCTION [proj, refor, or] = gt6DCreateProjDataFromExtendedGrain(gr_id, phase_id, varargin)
function [refor, estim_space_bbox_pix, estim_orient_bbox_rod] = gt6DCreateProjDataFromExtendedGrain(first, last, workingdirectory, phase_id, varargin)
% FUNCTION [proj, refor, or] = gt6DCreateProjDataFromExtendedGrain(first, last, workingdirectory, phase_id, varargin)
% proj: is a grain.proj structure
% refor: is a grain structure for the average orientation in the average
% point
% or: is a set of grain structures for the extreeme corners of the
% orientation space that should be considered
if (~exist('phase_id', 'var'))
phase_id = 1;
currentDir = pwd;
cd(workingdirectory);
if (isdeployed)
global GT_DB %#ok<NUSED,TLEV>
global GT_MATLAB_HOME %#ok<NUSED,TLEV>
load('workspaceGlobal.mat');
first = str2double(first);
last = str2double(last);
phase_id = str2double(phase_id);
end
conf = struct( ...
......@@ -22,308 +29,316 @@ function [refor, estim_space_bbox_pix, estim_orient_bbox_rod] = gt6DCreateProjDa
'use_volume_mask', true, ...
'save', false, ...
'use_extended', true, ...
'psf', []);
'psf', [], ...
'grain_ids', []);
conf = parse_pv_pairs(conf, varargin);
if (~isstruct(gr_id))
fprintf('Loading grain: ')
gr = gtForwardSimulate_v2(gr_id, gr_id, pwd, phase_id, ...
'save_grain', false, 'display_figure', false);
fprintf('Done.\n')
else
gr = gr_id;
end
p = gtLoadParameters();
refgr = gr;
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));
fwdsim_sel_refl_logic(ref_included(ref_selected)) = true;
% fwdsim_inc_refl = ref_ondet(ref_included);
fwdsim_inc_refl_logic = false(size(ref_ondet));
fwdsim_inc_refl_logic(ref_included) = true;
if (conf.include_all)
used_refl = ref_ondet;
if isempty(conf.grain_ids)
grains_list = first : last;
else
used_refl = fwdsim_sel_refl;
grains_list = conf.grain_ids;
% in case the code is deployed, we have to convert to strings
if ischar(grains_list)
grains_list = sscanf(grains_list, '%d');
end
end
[gr_center_pix, bbox_size_pix, estim_space_bbox_pix] = gt6DMergeRealSpaceVolumes(gr, 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 (isempty(conf.ospace_lims))
sampler = GtOrientationSampling(p, gr, ...
'detector_index', conf.det_index, 'verbose', conf.verbose);
r_vecs = sampler.guess_ODF_BB()';
estim_orient_bbox_rod = [min(r_vecs, [], 1), max(r_vecs, [], 1)];
bbox_size_rod = estim_orient_bbox_rod(4:6) - estim_orient_bbox_rod(1:3);
for ii = 1 : numel(grains_list)
fprintf('Loading grain: ')
gr = gtLoadGrain(phase_id, grains_list(ii));
fprintf('Done.\n')
% oversizing the orienation a bit
delta_bbox_size_rod = bbox_size_rod * (conf.ospace_oversize - 1) / 2;
estim_orient_bbox_rod = estim_orient_bbox_rod + [-delta_bbox_size_rod, delta_bbox_size_rod];
refgr = gr;
conf.ospace_lims = 2 * atand([-delta_bbox_size_rod, delta_bbox_size_rod]);
else
diff_r_vecs = tand(conf.ospace_lims / 2);
estim_orient_bbox_rod = gr.R_vector([1:3 1:3]) + diff_r_vecs;
end
bbox_size_rod = estim_orient_bbox_rod(4:6) - estim_orient_bbox_rod(1:3);
ref_ondet = refgr.proj(conf.det_index).ondet;
ref_included = refgr.proj(conf.det_index).included;
ref_selected = refgr.proj(conf.det_index).selected;
gr_center_rod = (estim_orient_bbox_rod(4:6) + estim_orient_bbox_rod(1:3)) / 2;
fwdsim_sel_refl = ref_ondet(ref_included(ref_selected));
fwdsim_sel_refl_logic = false(size(ref_ondet));
fwdsim_sel_refl_logic(ref_included(ref_selected)) = true;
% fwdsim_inc_refl = ref_ondet(ref_included);
fwdsim_inc_refl_logic = false(size(ref_ondet));
fwdsim_inc_refl_logic(ref_included) = true;
if (conf.include_all)
used_refl = ref_ondet;
else
used_refl = fwdsim_sel_refl;
end
bbox_size_deg = 2 * atand(bbox_size_rod);
[gr_center_pix, bbox_size_pix, estim_space_bbox_pix] = gt6DMergeRealSpaceVolumes(gr, conf.det_index, conf.rspace_oversize);
gr_center_mm = gtGeoSam2Sam(gr_center_pix, p.recgeo, p.samgeo, false, false);
% 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_rod, bbox_size_rod, ...
refgr, p);
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);
if (isempty(conf.ospace_lims))
sampler = GtOrientationSampling(p, gr, ...
'detector_index', conf.det_index, 'verbose', conf.verbose);
r_vecs = sampler.guess_ODF_BB()';
estim_orient_bbox_rod = [min(r_vecs, [], 1), max(r_vecs, [], 1)];
bbox_size_rod = estim_orient_bbox_rod(4:6) - estim_orient_bbox_rod(1:3);
fprintf(' Estimated orientation BBox: [%3.3f, %3.3f, %3.3f] -> [%3.3f, %3.3f, %3.3f]\n', estim_orient_bbox_rod);
fprintf(' Relative BBox size (center): [%3.3f, %3.3f, %3.3f] <- * -> [%3.3f, %3.3f, %3.3f] deg\n', conf.ospace_lims);
fprintf(' BBox size: %3.3f, %3.3f, %3.3f (deg)\n', bbox_size_deg);
fprintf('\n');
end
% oversizing the orienation a bit
delta_bbox_size_rod = bbox_size_rod * (conf.ospace_oversize - 1) / 2;
estim_orient_bbox_rod = estim_orient_bbox_rod + [-delta_bbox_size_rod, delta_bbox_size_rod];
or_abs = cat(1, ors(:).allblobs);
num_ors = numel(or_abs);
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(conf.det_index).uvw(used_refl, :);
conf.ospace_lims = 2 * atand([-delta_bbox_size_rod, delta_bbox_size_rod]);
else
uvw_tab(:, ii_g, :) = or_abs(ii_g).uvw(used_refl, :);
diff_r_vecs = tand(conf.ospace_lims / 2);
estim_orient_bbox_rod = gr.R_vector([1:3 1:3]) + diff_r_vecs;
end
end
refor = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, ...
'center', gr_center_mm, 'R_vector', gr_center_rod );
refor = gtCalculateGrain(refor, p, 'ref_omind', refgr.allblobs.omind);
refor.bb_ors = ors;
bbox_size_rod = estim_orient_bbox_rod(4:6) - estim_orient_bbox_rod(1:3);
% Let's treat those blobs at the w edge 360->0
% (from the sampled orientations perspective)
refor_ws = refor.allblobs.omega(used_refl) / gtGetOmegaStepDeg(p, conf.det_index);
uvw_tab(:, :, 3) = gtGrainAnglesTabularFix360deg(uvw_tab(:, :, 3), refor_ws, p);
gr_center_rod = (estim_orient_bbox_rod(4:6) + estim_orient_bbox_rod(1:3)) / 2;
img_bboxes_initial = [
squeeze(floor(min(uvw_tab, [], 2))), ...
squeeze( ceil(max(uvw_tab, [], 2))) ];
bbox_size_deg = 2 * atand(bbox_size_rod);
refor_ns = refor.allblobs.eta(used_refl);
% 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_rod, bbox_size_rod, ...
refgr, p);
img_sizes_initial = img_bboxes_initial(:, 4:6) - img_bboxes_initial(:, 1:3) + 1;
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);
% We avoid the vertical spots for convenience
inconvenient_etas = acosd(abs(cosd(refor_ns))) < conf.min_eta;
fprintf(' Estimated orientation BBox: [%3.3f, %3.3f, %3.3f] -> [%3.3f, %3.3f, %3.3f]\n', estim_orient_bbox_rod);
fprintf(' Relative BBox size (center): [%3.3f, %3.3f, %3.3f] <- * -> [%3.3f, %3.3f, %3.3f] deg\n', conf.ospace_lims);
fprintf(' BBox size: %3.3f, %3.3f, %3.3f (deg)\n', bbox_size_deg);
fprintf('\n');
end
img_bboxes = img_bboxes_initial(~inconvenient_etas, :);
img_sizes = img_sizes_initial(~inconvenient_etas, :);
or_abs = cat(1, ors(:).allblobs);
num_ors = numel(or_abs);
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(conf.det_index).uvw(used_refl, :);
else
uvw_tab(:, ii_g, :) = or_abs(ii_g).uvw(used_refl, :);
end
end
max_img_sizes = [max(img_sizes(:, 1)), max(img_sizes(:, 2))];
stack_imgs_oversize = min(p.fsim.oversize, conf.stack_oversize);
stackUSize = round(max_img_sizes(1) * stack_imgs_oversize);
stackVSize = round(max_img_sizes(2) * stack_imgs_oversize);
refor = struct( ...
'id', refgr.id, 'phaseid', refgr.phaseid, ...
'center', gr_center_mm, 'R_vector', gr_center_rod );
refor = gtCalculateGrain(refor, p, 'ref_omind', refgr.allblobs.omind);
refor.bb_ors = ors;
% Let's treat those blobs at the w edge 360->0
% (from the sampled orientations perspective)
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))), ...
squeeze( ceil(max(uvw_tab, [], 2))) ];
refor_ns = refor.allblobs.eta(used_refl);
img_sizes_initial = img_bboxes_initial(:, 4:6) - img_bboxes_initial(:, 1:3) + 1;
% We avoid the vertical spots for convenience
inconvenient_etas = acosd(abs(cosd(refor_ns))) < conf.min_eta;
img_bboxes = img_bboxes_initial(~inconvenient_etas, :);
img_sizes = img_sizes_initial(~inconvenient_etas, :);
max_img_sizes = [max(img_sizes(:, 1)), max(img_sizes(:, 2))];
stack_imgs_oversize = min(p.fsim.oversize, conf.stack_oversize);
stackUSize = round(max_img_sizes(1) * stack_imgs_oversize);
stackVSize = round(max_img_sizes(2) * stack_imgs_oversize);
if (conf.verbose)
if (conf.include_all)
fwdsim_inc_refl_logic_disp = fwdsim_inc_refl_logic;
fwdsim_sel_refl_logic_disp = fwdsim_sel_refl_logic;
else
fwdsim_inc_refl_logic_disp = fwdsim_inc_refl_logic(ref_included(ref_selected));
fwdsim_sel_refl_logic_disp = fwdsim_sel_refl_logic(ref_included(ref_selected));
end
fprintf('%02d) du %8d, dv %8d, dw %8d, eta: %7.3f <- used: %d, u: [%4d %4d], v: [%4d %4d], w: [%4d %4d], FwdSim: i %d, s %d\n', ...
[(1:numel(refor_ns))', img_sizes_initial, refor_ns, ~inconvenient_etas, ...
img_bboxes_initial(:, [1 4 2 5 3 6]), ...
fwdsim_inc_refl_logic_disp, fwdsim_sel_refl_logic_disp ]');
fprintf('\n');
fprintf(' Maximum images size: [%3d, %3d]\n', max_img_sizes);
fprintf('Stack images size (oversize: %1.2f): [%3d, %3d]\n', stack_imgs_oversize, stackUSize, stackVSize);
fprintf('\n');
end
if (conf.verbose)
fprintf('Loading raw images: ')
if (conf.include_all)
fwdsim_inc_refl_logic_disp = fwdsim_inc_refl_logic;
fwdsim_sel_refl_logic_disp = fwdsim_sel_refl_logic;
num_blobs = numel(used_refl);
sel_reflections = 1:num_blobs;
else
fwdsim_inc_refl_logic_disp = fwdsim_inc_refl_logic(ref_included(ref_selected));
fwdsim_sel_refl_logic_disp = fwdsim_sel_refl_logic(ref_included(ref_selected));
num_blobs = numel(ref_included);
sel_reflections = find(ref_selected);
end
fprintf('%02d) du %8d, dv %8d, dw %8d, eta: %7.3f <- used: %d, u: [%4d %4d], v: [%4d %4d], w: [%4d %4d], FwdSim: i %d, s %d\n', ...
[(1:numel(refor_ns))', img_sizes_initial, refor_ns, ~inconvenient_etas, ...
img_bboxes_initial(:, [1 4 2 5 3 6]), ...
fwdsim_inc_refl_logic_disp, fwdsim_sel_refl_logic_disp ]');
fprintf('\n');
fprintf(' Maximum images size: [%3d, %3d]\n', max_img_sizes);
fprintf('Stack images size (oversize: %1.2f): [%3d, %3d]\n', stack_imgs_oversize, stackUSize, stackVSize);
fprintf('\n');
end
sel_reflections = sel_reflections(~inconvenient_etas);
fprintf('Loading raw images: ')
if (conf.include_all)
num_blobs = numel(used_refl);
sel_reflections = 1:num_blobs;
else
num_blobs = numel(ref_included);
sel_reflections = find(ref_selected);
end
sel_reflections = sel_reflections(~inconvenient_etas);
blobs = gtFwdSimBlobDefinition('blob', num_blobs);
if (conf.use_volume_mask)
verts_rec = [ ...
estim_space_bbox_pix([1 2 3]); ...
estim_space_bbox_pix([1 2 6]); ...
estim_space_bbox_pix([1 5 3]); ...
estim_space_bbox_pix([1 5 6]); ...
estim_space_bbox_pix([4 2 3]); ...
estim_space_bbox_pix([4 2 6]); ...
estim_space_bbox_pix([4 5 3]); ...
estim_space_bbox_pix([4 5 6]); ...
];
verts_rec = bsxfun(@minus, verts_rec, gr_center_pix);
proj_mask_bls = gt6DSpreadProjectVertices2Det(refor, verts_rec, used_refl(~inconvenient_etas), p, conf.det_index);
end
blobs = gtFwdSimBlobDefinition('blob', num_blobs);
num_sel_refl = numel(sel_reflections);
for ii_i = 1:num_sel_refl
num_chars = fprintf('%02d/%02d', ii_i, num_sel_refl);
if (conf.use_volume_mask)
verts_rec = [ ...
estim_space_bbox_pix([1 2 3]); ...
estim_space_bbox_pix([1 2 6]); ...
estim_space_bbox_pix([1 5 3]); ...
estim_space_bbox_pix([1 5 6]); ...
estim_space_bbox_pix([4 2 3]); ...
estim_space_bbox_pix([4 2 6]); ...
estim_space_bbox_pix([4 5 3]); ...
estim_space_bbox_pix([4 5 6]); ...
];
verts_rec = bsxfun(@minus, verts_rec, gr_center_pix);
proj_mask_bls = gt6DSpreadProjectVertices2Det(refor, verts_rec, used_refl(~inconvenient_etas), p, conf.det_index);
end
ii_b = sel_reflections(ii_i);
num_sel_refl = numel(sel_reflections);
for ii_i = 1:num_sel_refl
num_chars = fprintf('%02d/%02d', ii_i, num_sel_refl);
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, 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, :)];
ii_b = sel_reflections(ii_i);
% Transposing to keep the same convention as spots
blob_vol = permute(blob_vol, [2 1 3]);
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, 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, :)];
blobs(ii_b).mbbsize = blob_bb(4:6);
blobs(ii_b).mbbu = [blob_bb(1), blob_bb(1) + blob_bb(4) - 1];
blobs(ii_b).mbbv = [blob_bb(2), blob_bb(2) + blob_bb(5) - 1];
blobs(ii_b).mbbw = [blob_bb(3), blob_bb(3) + blob_bb(6) - 1];
% Transposing to keep the same convention as spots
blob_vol = permute(blob_vol, [2 1 3]);
blob_size_im = [stackUSize, stackVSize, blob_bb(6)+2];
blobs(ii_b).mbbsize = blob_bb(4:6);
blobs(ii_b).mbbu = [blob_bb(1), blob_bb(1) + blob_bb(4) - 1];
blobs(ii_b).mbbv = [blob_bb(2), blob_bb(2) + blob_bb(5) - 1];
blobs(ii_b).mbbw = [blob_bb(3), blob_bb(3) + blob_bb(6) - 1];
shifts_blob = gtFwdSimGetStackShifts(stackUSize, stackVSize, blob_bb, false);
shifts_blob = [shifts_blob.u, shifts_blob.v, 1];
blob_size_im = [stackUSize, stackVSize, blob_bb(6)+2];
blobs(ii_b).intm = gtPlaceSubVolume( ...
zeros(blob_size_im, 'single'), single(blob_vol), shifts_blob);
shifts_blob = gtFwdSimGetStackShifts(stackUSize, stackVSize, blob_bb, false);
shifts_blob = [shifts_blob.u, shifts_blob.v, 1];
blobs(ii_b).bbsize = blob_size_im;
blobs(ii_b).intm = gtPlaceSubVolume( ...
zeros(blob_size_im, 'single'), single(blob_vol), shifts_blob);
blob_bb_im = [blob_bb(1:3) - shifts_blob, blob_size_im];
blobs(ii_b).bbsize = blob_size_im;
blobs(ii_b).bbuim = [blob_bb_im(1), blob_bb_im(1) + blob_bb_im(4) - 1];
blobs(ii_b).bbvim = [blob_bb_im(2), blob_bb_im(2) + blob_bb_im(5) - 1];
blobs(ii_b).bbwim = [blob_bb_im(3), blob_bb_im(3) + blob_bb_im(6) - 1];
blob_bb_im = [blob_bb(1:3) - shifts_blob, blob_size_im];
if (conf.use_volume_mask)
shift_mask = [ ...
(proj_mask_bls(ii_i).bbuim(:, 1) - blobs(ii_b).bbuim(:, 1)), ...
(proj_mask_bls(ii_i).bbvim(:, 1) - blobs(ii_b).bbvim(:, 1)), ...
1, ...
];
blobs(ii_b).bbuim = [blob_bb_im(1), blob_bb_im(1) + blob_bb_im(4) - 1];
blobs(ii_b).bbvim = [blob_bb_im(2), blob_bb_im(2) + blob_bb_im(5) - 1];
blobs(ii_b).bbwim = [blob_bb_im(3), blob_bb_im(3) + blob_bb_im(6) - 1];
mask = proj_mask_bls(ii_i).mask(:, :, ones(blobs(ii_b).bbsize(3)-2, 1));
if (conf.use_volume_mask)
shift_mask = [ ...
(proj_mask_bls(ii_i).bbuim(:, 1) - blobs(ii_b).bbuim(:, 1)), ...
(proj_mask_bls(ii_i).bbvim(:, 1) - blobs(ii_b).bbvim(:, 1)), ...
1, ...
];
blobs(ii_b).mask = gtPlaceSubVolume( ...
false(blob_size_im), mask, shift_mask);
else
blobs(ii_b).mask = gtPlaceSubVolume( ...
false(blob_size_im), true(size(blob_vol)), shifts_blob);
end
mask = proj_mask_bls(ii_i).mask(:, :, ones(blobs(ii_b).bbsize(3)-2, 1));
blob_int = sum(blobs(ii_b).intm(blobs(ii_b).mask));
blobs(ii_b).intm = blobs(ii_b).intm / blob_int;
blobs(ii_b).mask = gtPlaceSubVolume( ...
false(blob_size_im), mask, shift_mask);
else
blobs(ii_b).mask = gtPlaceSubVolume( ...
false(blob_size_im), true(size(blob_vol)), shifts_blob);
end
blobs(ii_b).intensity = blob_int;
blob_int = sum(blobs(ii_b).intm(blobs(ii_b).mask));
blobs(ii_b).intm = blobs(ii_b).intm / blob_int;
fprintf(repmat('\b', [1 num_chars]));
fprintf(repmat(' ', [1 num_chars]));
fprintf(repmat('\b', [1 num_chars]));
end
fprintf('Done.\n')
blobs(ii_b).intensity = blob_int;
proj = gtFwdSimProjDefinition();
fprintf(repmat('\b', [1 num_chars]));
fprintf(repmat(' ', [1 num_chars]));
fprintf(repmat('\b', [1 num_chars]));
end
fprintf('Done.\n')
proj.ondet = ref_ondet;
if (~conf.include_all)
proj.included = ref_included;
else
proj.included = 1:numel(used_refl);
end
proj = gtFwdSimProjDefinition();
% Now filling the remaining blobs with the blobs from refgr
% (even if they won't be used)
% ! this might change in the future ! to maybe include fwd projected
% regions on the detector
not_sel_reflections = true(size(proj.included));
not_sel_reflections(sel_reflections) = false;
not_sel_reflections = find(not_sel_reflections);
if (~conf.include_all)
blobs(not_sel_reflections) = refgr.proj.bl(not_sel_reflections);
for ii_i = 1:numel(not_sel_reflections)
ii_b = not_sel_reflections(ii_i);
blobs(ii_b).intm = [];
blobs(ii_b).mask = [];
proj.ondet = ref_ondet;
if (~conf.include_all)
proj.included = ref_included;
else
proj.included = 1:numel(used_refl);
end
else
for ii_i = 1:numel(not_sel_reflections)
ii_b = not_sel_reflections(ii_i);
blobs(ii_b).bbsize = zeros(1, 3);
blobs(ii_b).bbuim = zeros(1, 2);
blobs(ii_b).bbvim = zeros(1, 2);
blobs(ii_b).bbwim = zeros(1, 2);
blobs(ii_b).intm = zeros(stackUSize, stackVSize);
% Now filling the remaining blobs with the blobs from refgr
% (even if they won't be used)
% ! this might change in the future ! to maybe include fwd projected
% regions on the detector
not_sel_reflections = true(size(proj.included));
not_sel_reflections(sel_reflections) = false;
not_sel_reflections = find(not_sel_reflections);
if (~conf.include_all)
blobs(not_sel_reflections) = refgr.proj.bl(not_sel_reflections);
for ii_i = 1:numel(not_sel_reflections)
ii_b = not_sel_reflections(ii_i);
blobs(ii_b).intm = [];
blobs(ii_b).mask = [];
end
else
for ii_i = 1:numel(not_sel_reflections)
ii_b = not_sel_reflections(ii_i);
blobs(ii_b).bbsize = zeros(1, 3);
blobs(ii_b).bbuim = zeros(1, 2);
blobs(ii_b).bbvim = zeros(1, 2);
blobs(ii_b).bbwim = zeros(1, 2);
blobs(ii_b).intm = zeros(stackUSize, stackVSize);
end
end
end
finally_selected = false(size(proj.included));
finally_selected(sel_reflections) = true;
proj.selected = finally_selected;
spots = arrayfun(@(x){sum(x.intm, 3)}, blobs);
spots = permute(cat(3, spots{:}), [1 3 2]);
% These reconstructions usually need some padding
proj.centerpix = gr_center_pix;
proj.bl = blobs;
proj.stack = spots;
% We apply a simple oversize
vol_size = bbox_size_pix;
% vol_size = ceil(bbox_size_pix + mean(bbox_size_pix) * 0.3);
proj.vol_size_x = vol_size(2);
proj.vol_size_y = vol_size(1);
proj.vol_size_z = vol_size(3);
if (~isempty(conf.psf))
proj.psf = conf.psf;
end
refor.proj(conf.det_index) = proj;
refor.conf = conf;
if (conf.save)
fprintf('Saving the extended grain file..')
grain_filename = fullfile(p.acq.dir, '4_grains', ...
sprintf('phase_%02d', phase_id), ...
sprintf('grain_extended_%04d.mat', gr.id));
save(grain_filename, '-struct', 'refor', '-v7.3');
fprintf('\b\b: Done.\n')
fprintf('Saving to sample.mat..')
sample = GtSample.loadFromFile();
sample.phases{phase_id}.extended_params(gr.id) = ...
GtPhase.makeExtendedField(estim_space_bbox_pix, estim_orient_bbox_rod);
sample.phases{phase_id}.setUseExtended(gr.id, conf.use_extended);
sample.saveToFile();
fprintf('\b\b: Done.\n')
finally_selected = false(size(proj.included));
finally_selected(sel_reflections) = true;
proj.selected = finally_selected;
spots = arrayfun(@(x){sum(x.intm, 3)}, blobs);
spots = permute(cat(3, spots{:}), [1 3 2]);
% These reconstructions usually need some padding
proj.centerpix = gr_center_pix;
proj.bl = blobs;
proj.stack = spots;
% We apply a simple oversize
vol_size = bbox_size_pix;
% vol_size = ceil(bbox_size_pix + mean(bbox_size_pix) * 0.3);
proj.vol_size_x = vol_size(2);
proj.vol_size_y = vol_size(1);
proj.vol_size_z = vol_size(3);
if (~isempty(conf.psf))
proj.psf = conf.psf;
end
refor.proj(conf.det_index) = proj;
refor.conf = conf;
if (conf.save)
fprintf('Saving the extended grain file..')
grain_filename = fullfile(p.acq.dir, '4_grains', ...
sprintf('phase_%02d', phase_id), ...
sprintf('grain_extended_%04d.mat', gr.id));
save(grain_filename, '-struct', 'refor', '-v7.3');
fprintf('\b\b: Done.\n')
fprintf('Saving to sample.mat..')
sample = GtSample.loadFromFile();
sample.phases{phase_id}.extended_params(gr.id) = ...
GtPhase.makeExtendedField(estim_space_bbox_pix, estim_orient_bbox_rod);
sample.phases{phase_id}.setUseExtended(gr.id, conf.use_extended);
sample.saveToFile();
fprintf('\b\b: Done.\n')
end
end
cd(currentDir);
end
zUtil_Deformation/gt6DCreateProjDataFromExtendedGrains.m 0 → 100644
+ 47
− 0
View file @ b38d93a0
function gt6DCreateProjDataFromExtendedGrains(first, last, workingdirectory, phaseID, parameters, varargin)
% FUNCTION [proj, refor, or] = gt6DCreateProjDataFromExtendedGrainss(first, last, workingdirectory, phase_id, varargin)
% proj: is a grain.proj structure
% refor: is a grain structure for the average orientation in the average
% point
% or: is a set of grain structures for the extreeme corners of the
% orientation space that should be considered
cd (workingdirectory);
if (~exist('parameters','var') || isempty(parameters))
parameters = gtLoadParameters();
end
conf = struct( ...
'verbose', false, ...
'min_eta', 20, ...
'det_index', 1, ...
'ospace_oversize', 1.1, ...
'rspace_oversize', 1.1, ...
'stack_oversize', 1.2, ...
'ospace_lims', [], ...
'include_all', false, ...
'use_volume_mask', true, ...
'save', false, ...
'use_extended', true, ...
'psf', [], ...
'grain_ids', []);
conf = parse_pv_pairs(conf, varargin);
if isempty(conf.grain_ids)
grains_list = first : last;
else
grains_list = conf.grain_ids;
end
check = inputwdefault('Launch Creation of Extended Grains via OAR? [y/n]', 'y');
launch_on_OAR = strcmpi(check, 'y');
if (launch_on_OAR)
gtLaunchOnOAR('gt6DCreateProjDataFromExtendedGrain', first, last, phaseID_str, ...
parameters, 'njobs', 8, 'walltime', 3600 * 16, ...
'gpu', false, 'distribute', true, 'list', conf.grain_ids);
else
gt6DCreateProjDataFromExtendedGrain(first, last, pwd, phaseID, parameters.acq.dir, ...
parameters, 'grain_ids', conf.grain_ids);
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