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

Assembling: added the possibility to assemble DMVOLs from 6D reconstructions 


Signed-off-by: default avatarNicola Vigano <nicola.vigano@esrf.fr>
parent d43cde17
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
4_grains/GtThreshold.m
+ 1
− 1
View file @ a0bfe686
......@@ -449,7 +449,7 @@ classdef GtThreshold < handle
% for bb start values
grain.seg = gtCrop(segvol, ceil(segbb));
grain.Area = regprop.Area;
segbb(1:3) = grain.proj.shift + floor(segbb(1:3));
segbb(1:3) = grain.proj.shift + ceil(segbb(1:3)) - 1;
grain.segbb = segbb;
end
grain.threshold = threshold;
......
5_reconstruction/GtAssembleVol3D.m
+ 102
− 70
View file @ a0bfe686
......@@ -125,6 +125,10 @@ classdef GtAssembleVol3D < handle
grainDirsModel = fullfile('4_grains', 'phase_%02d', 'grain_%04d.mat');
obj.localPar.cache.createModel(grainDirsModel, 'grain');
end
if (~obj.localPar.cache.isModel('grain_details'))
grainDirsModel = fullfile('4_grains', 'phase_%02d', 'grain_details_%04d.mat');
obj.localPar.cache.createModel(grainDirsModel, 'grain_details');
end
if (obj.localPar.deal_with_twins)
disp('dealing with twins')
......@@ -153,19 +157,28 @@ classdef GtAssembleVol3D < handle
end
% Let's initialize the volume
volume = zeros(obj.parameters.acq.bb(3), ...
obj.parameters.acq.bb(3), ...
obj.parameters.acq.bb(4));
grain_ids = zeros(obj.parameters.acq.bb([3, 3, 4]), 'int16');
if (obj.localPar.deal_with_twins)
volume = obj.calculatePhaseVolumeTwins(volume, sampleInfo, phaseID, list);
grain_ids = obj.calculatePhaseVolumeTwins(grain_ids, sampleInfo, phaseID, list);
obj.savePhaseTwinOutput(phaseID);
else
volume = obj.doAssemblePhaseVolume(volume, sampleInfo, phaseID, list);
grain_ids = obj.doAssemblePhaseVolume(grain_ids, sampleInfo, phaseID, list);
end
% store volume
obj.volumes.phases{phaseID} = volume;
obj.volumes.phases{phaseID}.grains = grain_ids;
if (strcmpi(obj.parameters.rec.grains.algorithm(1:2), '6D'))
dmvol = zeros([obj.parameters.acq.bb([3, 3, 4]), 3]);
dmvol = obj.doAssemblePhaseDmVolume(dmvol, sampleInfo, phaseID, list);
if (ismember(obj.localPar.overlaps, {'zero', 'conflict'}))
% set overlapping regions to zero - this could be another label colour
dmvol(grain_ids == -1 | grain_ids == 0) = 0;
end
obj.volumes.phases{phaseID}.dmvol = dmvol;
end
if (obj.localPar.autosave)
obj.savePhaseVolume(phaseID);
......@@ -184,54 +197,6 @@ classdef GtAssembleVol3D < handle
end
end % end function assembleAllPhases
function assembleInExistingVolume(obj, phaseID, filename)
% GTASSEMBLEVOL3D/ASSEMBLEINEXISTINGVOLUME Discouraged function to write
% a new grain map, to a previously existing grain map.
if (~exist('phaseID', 'var'))
phaseID = 1;
end
[sampleInfo, list] = obj.initialiseAssembling(phaseID);
% look for the file
if (~exist('filename', 'var'))
filename = obj.localPar.add_to_filename;
end
absolute_file_path = fullfile(obj.parameters.acq.dir, ...
'5_reconstruction', filename);
if exist(filename, 'file')
volume = edf_read(filename);
% check if it's in the 5_reconstruction directory
elseif exist( absolute_file_path, 'file')
volume = edf_read(absolute_file_path);
% otherwise, complain
else
error('ASSEMBLE:wrong_argument', ...
['Cannot find the volume to add to: %s\ncheck the' ...
' filename / path\n'], ...
filename)
end
% if we want to add grains to an existing volume, this is loaded
% in. Then labels in the list are set to zero, before the grains
% are "re-added" to the volume. I don't think this will deal
% with overlapping grains correctly.
if (max(volume(:)) > 0)
for ii = 1:length(list)
volume(obj.volume == list(ii)) = 0;
end
end
volume = obj.doAssemblePhaseVolume(volume, sampleInfo, phaseID, list);
% Dumb saving.
name = fullfile('5_reconstruction', filename);
fprintf('writing volume %s\n', name);
edf_write(volume, name, 'uint16');
end % end function assembleInExistingVolume
function assembleCompleteVolume(obj)
% GTASSEMBLEVOL3D/ASSEMBLEFINALVOLUME
......@@ -239,20 +204,30 @@ classdef GtAssembleVol3D < handle
sampleInfo = obj.localPar.cache.get('sample', {}, 'sample');
% Let's initialize the phase IDs volumes
volumePIDs = zeros( obj.parameters.acq.bb(3), ...
obj.parameters.acq.bb(3), ...
obj.parameters.acq.bb(4), 'int16' );
volumePIDs = zeros(obj.parameters.acq.bb([3, 3, 4]), 'int16');
% Let's initialize the grains IDs volumes
volumeGIDs = zeros( obj.parameters.acq.bb(3), ...
obj.parameters.acq.bb(3), ...
obj.parameters.acq.bb(4), 'int16' );
volumeGIDs = zeros(obj.parameters.acq.bb([3, 3, 4]), 'int16');
assemble_6D = strcmpi(obj.parameters.rec.grains.algorithm(1:2), '6D');
if (assemble_6D)
volumeDM = zeros([obj.parameters.acq.bb([3, 3, 4]), 3]);
end
phasesNum = length(sampleInfo.phases);
gauge = GtGauge(phasesNum, 'Processing phases: ');
for ii = 1:phasesNum
gauge.incrementAndDisplay();
% no real need to cache phase volumes
phaseVol = gtMATVolReader(sampleInfo.phases{ii}.volumeFile);
phase_file_name = sampleInfo.phases{ii}.volumeFile;
phase_file_name = regexp(phase_file_name, ':', 'split');
phase_file = matfile(phase_file_name{1});
phaseVol = phase_file.vol;
assemble_phase_6D = assemble_6D ...
&& ismember('dmvol', fieldnames(phase_file));
if (assemble_phase_6D)
phaseDmVol = phase_file.dmvol;
end
grainMap = phaseVol;
phaseVol = int16(phaseVol ~= 0) * ii;
......@@ -262,6 +237,9 @@ classdef GtAssembleVol3D < handle
maskVol = gtMATVolReader([maskFile ':mask']);
phaseVol = phaseVol .* int16(maskVol);
grainMap = grainMap .* int16(maskVol);
if (assemble_phase_6D)
phaseDmVol(~maskVol) = 0;
end
end
switch (obj.localPar.overlaps)
......@@ -272,10 +250,18 @@ classdef GtAssembleVol3D < handle
volumeGIDs = volumeGIDs + grainMap;
volumeGIDs(overlappingPoints) = -1;
if (assemble_phase_6D)
volumeDM = volumeDM + phaseDmVol;
volumeDM(overlappingPoints) = -1;
end
case {'summed'}
volumePIDs = volumePIDs + phaseVol;
volumeGIDs = volumeGIDs + grainMap;
if (assemble_phase_6D)
volumeDM = volumeDM + phaseDmVol;
end
end
end
gauge.delete();
......@@ -295,6 +281,7 @@ classdef GtAssembleVol3D < handle
obj.volumes.complete.phases = volumePIDs;
obj.volumes.complete.grains = volumeGIDs;
obj.volumes.complete.dmvol = volumeDM;
% Autosaving
if (obj.localPar.autosave)
......@@ -314,9 +301,11 @@ classdef GtAssembleVol3D < handle
name = obj.checkAndSaveExistingFile(fileName);
% Small trick to keep consistency in field name in the written file.
outStruct = [];
volField = 'vol';
outStruct.(volField) = int16(obj.volumes.phases{phaseID});
outStruct = struct( ...
'vol', int16(obj.volumes.phases{phaseID}.grains) );
if (isfield(obj.volumes.phases{phaseID}, 'dmvol'))
outStruct.dmvol = obj.volumes.phases{phaseID}.dmvol; %#ok<STRNU>
end
% We lock, to avoid other processes to write to this file:
% Actually this is commented out because it should happen
......@@ -330,7 +319,7 @@ classdef GtAssembleVol3D < handle
obj.checkSampleInfo(sampleInfo, phaseID);
% Store info about volume in the metadata
sampleInfo.phases{phaseID}.volumeFile = [name ':' volField];
sampleInfo.phases{phaseID}.volumeFile = [name ':vol'];
obj.localPar.cache.set('sample', {}, sampleInfo, 'sample');
% Flush metadata
fprintf('Flushing sample.mat\n');
......@@ -341,7 +330,7 @@ classdef GtAssembleVol3D < handle
% Saving the data
fprintf('Writing phase %02d volume to: %s\n', phaseID, name);
save(name, '-struct', 'outStruct', volField, '-v7.3');
save(name, '-struct', 'outStruct', '-v7.3');
end % end function savePhaseVolume
......@@ -373,8 +362,8 @@ classdef GtAssembleVol3D < handle
% Saving the data
fprintf('Writing complete volume to: %s\n', name);
outStruct = obj.volumes.complete;
save(name, '-struct', 'outStruct', 'phases', 'grains', '-v7.3');
outStruct = obj.volumes.complete; %#ok<NASGU>
save(name, '-struct', 'outStruct', '-v7.3');
end % end function saveCompleteVolume
......@@ -495,7 +484,7 @@ classdef GtAssembleVol3D < handle
function volume = doAssemblePhaseVolume(obj, volume, sampleInfo, phaseID, list)
% GTASSEMBLEVOL3D/DOASSEMBLEPHASEVOLUME Function that effectively assembles phase volume
grainsNum = length(list);
gauge = GtGauge(grainsNum, sprintf('Phase %02d: ', phaseID));
......@@ -530,6 +519,49 @@ classdef GtAssembleVol3D < handle
gauge.delete();
end % end function doAssemblePhaseVolume
function dmvol = doAssemblePhaseDmVolume(obj, dmvol, sampleInfo, phaseID, list)
% GTASSEMBLEVOL3D/DOASSEMBLEPHASEVOLUME Function that effectively assembles phase volume
grainsNum = length(list);
gauge = GtGauge(grainsNum, sprintf('Phase %02d: ', phaseID));
for ii = list;
gauge.setExtraText(sprintf('(grain %04d)', ii));
gauge.incrementAndDisplay();
if (~sampleInfo.phases{phaseID}.getSelected(ii))
fprintf('\nSkipping deselected grain: %04d\n', ii);
gauge.rePrint();
continue;
end
try
volBBox = sampleInfo.phases{phaseID}.getBoundingBox(ii);
seg_vol = obj.localPar.cache.get('grain', {phaseID, ii}, 'seg');
odf6d = obj.localPar.cache.get('grain_details', {phaseID, ii}, 'ODF6D');
segbb = [volBBox(1:3) - odf6d.shift + 1, volBBox(4:6)];
for ii_dim = 1:3
single_dim_dmvol = odf6d.voxels_avg_R_vectors(:, :, :, ii_dim);
single_dim_dmvol = gtCrop(single_dim_dmvol, segbb);
single_dim_dmvol(~seg_vol) = 0;
dmvol(:, :, :, ii_dim) = gtPlaceSubVolume( ...
dmvol(:, :, :, ii_dim), single_dim_dmvol, ...
volBBox(1:3), [], obj.localPar.overlaps, false);
end
catch mexc
gtPrintException(mexc, sprintf('\nGrain %d failed!!\n', ii));
gauge.rePrint();
obj.redo = [obj.redo ii];
end % end try catch loop
% Trigger callback (usually used in graphical applications)
gtTriggerUpdateCallback(obj.localPar.update_status_cb);
end % for loop through all grains
gauge.delete();
end % end function doAssemblePhaseVolume
function volumeUpdate = calculatePhaseVolumeTwins(obj, volume, sampleInfo, phaseID, list)
% GTASSEMBLEVOL3D/CALCULATEPHASEVOLUMETWINS
% obj has the sigma as a property - not sure if useful like this
......
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