diff --git a/4_grains/GtThreshold.m b/4_grains/GtThreshold.m
index f051a1a9c84d7f610135f0acbc97aedad30660fe..3cc170a36a6a3fe482d02512d50c4e739f71829c 100644
--- a/4_grains/GtThreshold.m
+++ b/4_grains/GtThreshold.m
@@ -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;
diff --git a/5_reconstruction/GtAssembleVol3D.m b/5_reconstruction/GtAssembleVol3D.m
index af7e861ee55ed9e749eff9d6776fa127207aa0fb..2657a034f17d8b0b7ce41829ecb5ec886952536f 100644
--- a/5_reconstruction/GtAssembleVol3D.m
+++ b/5_reconstruction/GtAssembleVol3D.m
@@ -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