Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add additional test for PaCER function #69

Open
wants to merge 10 commits into
base: develop
Choose a base branch
from
193 changes: 152 additions & 41 deletions test/suite/testPaCER.m
Original file line number Diff line number Diff line change
Expand Up @@ -79,18 +79,18 @@
assert(norm(skelSkelmms_Mask_new{k} - refData_brainMask.skelSkelmms_Mask_Ref{k}) < tol);
end

%% test different electrode type:
%% test different electrodes type:
% load the reference data
refData_electrodeType = load([refDataPath filesep 'refData_PaCER_electrodeType.mat']);

% define the input argument
niiCT_electrodesType = niiCT_PostOP_new;
xml_Plan_new;

% generate the new output with Medtronic 3387 electrode type.
[elecModels_Medtronic3387_new, elecPointCloudsStruct_Medtronic3387_new, intensityProfiles_Medtronic3387_new, skelSkelmms_Medtronic3387_new] = PaCER(niiCT_electrodesType, 'medtronicXMLPlan', xml_Plan_new, 'electrodeType', 'Medtronic 3387');

% compare the new data against the reference data using XML plan and
% providing electrode type (Medtronic 3387)
% compare the new data against the reference data using XML plan and providing electrode type (Medtronic 3387)
structureComparison(elecModels_Medtronic3387_new, refData_electrodeType.elecModels_Medtronic3387_ref)
assert(isequal(elecPointCloudsStruct_Medtronic3387_new, refData_electrodeType.elecPointCloudsStruct_Medtronic3387_ref))
for k=1:length(intensityProfiles_Medtronic3387_new)
Expand All @@ -103,8 +103,7 @@
% generate the new output with Medtronic 3389 electrode type.
[elecModels_Medtronic3389_new, elecPointCloudsStruct_Medtronic3389_new, intensityProfiles_Medtronic3389_new, skelSkelmms_Medtronic3389_new] = PaCER(niiCT_electrodesType, 'medtronicXMLPlan', xml_Plan_new, 'electrodeType', 'Medtronic 3389');

% compare the new data against the reference data using XML plan and
% providing electrode type (Medtronic 3389)
% compare the new data against the reference data using XML plan and providing electrode type (Medtronic 3389)
structureComparison(elecModels_Medtronic3389_new, refData_electrodeType.elecModels_Medtronic3389_ref)
assert(isequal(elecPointCloudsStruct_Medtronic3389_new, refData_electrodeType.elecPointCloudsStruct_Medtronic3389_ref));
for k=1:length(intensityProfiles_Medtronic3389_new)
Expand All @@ -128,18 +127,8 @@
assert(norm(skelSkelmms_Boston_new{k} - refData_electrodeType.skelSkelmms_Boston_ref{k}) < tol);
end

%% test the warning messages
% test if slice thickness is greater than 1 mm
warningMessage = 'Slice thickness is greater than 1 mm! Independent contact detection is most likly not possible. Forcing contactAreaCenter based method.';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_PostOP_new}))

% test if slice thickness is greater than 0.7 mm
warningMessage = 'Slice thickness is greater than 0.7 mm! Independet contact detection might not work reliable in this case. However, for certain electrode types with large contacts spacings you might be lucky.';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_brainMask_new}))

%% test different case
% load reference data including the new condition for
% peakWaveCenter,'displayProfiles' true, 'displayMPR' true and final degree setup to 1
% load reference data including the new condition for peakWaveCenter,'displayProfiles' true, 'displayMPR' true and final degree setup to 1
refData_case = load([refDataPath filesep 'refData_PaCER_different_case.mat']);

% if contact detection set to peakWaveCenter
Expand All @@ -155,45 +144,167 @@
assert(norm(skelSkelmms_peakWaveCenter_new{k} - refData_case.skelSkelmms_peakWaveCenter_ref{k}) < tol);
end

% if 'displayProfiles' is true
[elecModels_displayProfiles_new, elecPointCloudsStruct_displayProfiles_new, intensityProfiles_displayProfiles_new, skelSkelmms_displayProfiles_new] = PaCER(niiCT_brainMask_new,'displayProfiles', true);
% if final degree setup to 1
[elecModels_degree_new, elecPointCloudsStruct_degree_new, intensityProfiles_degree_new, skelSkelmms_degree_new] = PaCER(niiCT_brainMask_new,'finalDegree', 1);

% compare reference data and new outputs
structureComparison(elecModels_displayProfiles_new, refData_case.elecModels_displayProfiles_ref)
assert(isequal(elecPointCloudsStruct_displayProfiles_new, refData_case.elecPointCloudsStruct_displayProfiles_ref))
structureComparison(elecModels_degree_new, refData_case.elecModels_degree_ref)
assert(isequal(elecPointCloudsStruct_degree_new, refData_case.elecPointCloudsStruct_degree_ref))
for k=1:length(intensityProfiles_degree_new)
assert(norm(intensityProfiles_degree_new{k} - refData_case.intensityProfiles_degree_ref{k}) < tol);
end
for k=1:length(skelSkelmms_degree_new)
assert(norm(skelSkelmms_degree_new{k} - refData_case.skelSkelmms_degree_ref{k}) < tol);
end

%% if 'displayProfiles and 'displayMPR' are true
% load reference data (with displayMPR and displayProfiles)
refData_MPR = load([refDataPath filesep 'refData_PaCER_DisplayMPR.mat']);

% generate new output
[elecModels_displayProfiles_new, elecPointCloudsStruct_displayProfiles_new, intensityProfiles_displayProfiles_new, skelSkelmms_displayProfiles_new] = PaCER(niiCT_brainMask_new,'displayProfiles', true, 'displayMPR', true);

% compare reference data and new data
structureComparison(elecModels_displayProfiles_new, refData_MPR.elecModels_displayProfiles_ref)
assert(isequal(elecPointCloudsStruct_displayProfiles_new, refData_MPR.elecPointCloudsStruct_displayProfiles_ref))
for k=1:length(intensityProfiles_displayProfiles_new)
assert(norm(intensityProfiles_displayProfiles_new{k} - refData_case.intensityProfiles_displayProfiles_ref{k}) < tol);
assert(norm(intensityProfiles_displayProfiles_new{k} - refData_MPR.intensityProfiles_displayProfiles_ref{k}) < tol);
end
for k=1:length(skelSkelmms_displayProfiles_new)
assert(norm(skelSkelmms_displayProfiles_new{k} - refData_case.skelSkelmms_displayProfiles_ref{k}) < tol);
assert(norm(skelSkelmms_displayProfiles_new{k} - refData_MPR.skelSkelmms_displayProfiles_ref{k}) < tol);
end

% if 'displayMPR' is true
[elecModels_displayMPR_new, elecPointCloudsStruct_displayMPR_new, intensityProfiles_displayMPR_new, skelSkelmms_displayMPR_new] = PaCER(niiCT_brainMask_new,'displayMPR', true);
%% testing additional parameters
%load the reference data generated by using additional inputs in PaCER function: reverseDir,contactAreaCenter , peak, peakWaveCenter
refData = load([refDataPath filesep 'refData_PaCER_additional_parameter.mat']);

% compare reference data and new outputs
structureComparison(elecModels_displayMPR_new, refData_case.elecModels_displayMPR_ref)
assert(isequal(elecPointCloudsStruct_displayMPR_new, refData_case.elecPointCloudsStruct_displayMPR_ref))
for k=1:length(intensityProfiles_displayMPR_new)
assert(norm(intensityProfiles_displayMPR_new{k} - refData_case.intensityProfiles_displayMPR_ref{k}) < tol);
%load a new CT input
niiCT_PaCER_9 = NiftiMod([inputDataPath filesep 'PaCER_ct_post_OK_9.nii.gz']);

% generate the new outputs with 'reverseDir', 'contactAreaCenter' and 'peak'
[elecModels_new_9, elecPointCloudsStruct_new_9, intensityProfiles_new_9, skelSkelmms_new_9] = PaCER(niiCT_PaCER_9, 'reverseDir', true, 'contactAreaCenter','peak');

% compare the generated data against the ref data
structureComparison(refData.elecModels_ref_9, elecModels_new_9)
assert(isequal(elecPointCloudsStruct_new_9, refData.elecPointCloudsStruct_ref_9));
for k=1:length(intensityProfiles_new_9)
assert(norm(intensityProfiles_new_9{k} - refData.intensityProfiles_ref_9{k}) < tol);
end
for k=1:length(skelSkelmms_displayMPR_new)
assert(norm(skelSkelmms_displayMPR_new{k} - refData_case.skelSkelmms_displayMPR_ref{k}) < tol);
for k=1:length(skelSkelmms_new_9)
assert(norm(skelSkelmms_new_9{k} - refData.skelSkelmms_ref_9{k}) < tol);
end

% if final degree setup to 1
%[elecModels_degree_ref, elecPointCloudsStruct_degree_ref, intensityProfiles_degree_ref, skelSkelmms_degree_ref] = PaCER(niiCT_brainMask_new,'finalDegree', 1);
[elecModels_degree_new, elecPointCloudsStruct_degree_new, intensityProfiles_degree_new, skelSkelmms_degree_new] = PaCER(niiCT_brainMask_new,'finalDegree', 1);
% load a new CT input
niiCT_PaCER_8 = NiftiMod([inputDataPath filesep 'PaCER_ct_post_OK_8.nii.gz']);

% compare reference data and new outputs
structureComparison(elecModels_degree_new, refData_case.elecModels_degree_ref)
assert(isequal(elecPointCloudsStruct_degree_new, refData_case.elecPointCloudsStruct_degree_ref))
for k=1:length(intensityProfiles_degree_new)
assert(norm(intensityProfiles_degree_new{k} - refData_case.intensityProfiles_degree_ref{k}) < tol);
% generate the new output with 'reverseDir', true, 'contactAreaCenter' and 'peakWaveCenter' additional inputs
[elecModels_new_8, elecPointCloudsStruct_new_8, intensityProfiles_new_8, skelSkelmms_new_8] = PaCER(niiCT_PaCER_8, 'reverseDir', true, 'contactAreaCenter','peakWaveCenter');

% compare the generated data against the ref data
structureComparison(refData.elecModels_ref_8, elecModels_new_8)
assert(isequal(elecPointCloudsStruct_new_8, refData.elecPointCloudsStruct_ref_8));
for k=1:length(intensityProfiles_new_8)
assert(norm(intensityProfiles_new_8{k} - refData.intensityProfiles_ref_8{k}) < tol);
end
for k=1:length(skelSkelmms_degree_new)
assert(norm(skelSkelmms_degree_new{k} - refData_case.skelSkelmms_degree_ref{k}) < tol);
for k=1:length(skelSkelmms_new_8)
assert(norm(skelSkelmms_new_8{k} - refData.skelSkelmms_ref_8{k}) < tol);
end

% load a new CT input
niiCT_PaCER_6 = NiftiMod([inputDataPath filesep 'PaCER_ct_post_OK_6.nii.gz']);

% generate the new output with 'reverseDir', 'contactDetectionMethod' and 'peak' input arguments
[elecModels_new_6, elecPointCloudsStruct_new_6, intensityProfiles_new_6, skelSkelmms_new_6] = PaCER(niiCT_PaCER_6, 'reverseDir', true, 'contactDetectionMethod','peak');

% compare the new generated data against the ref data
structureComparison(refData.elecModels_ref_6, elecModels_new_6)
assert(isequal(elecPointCloudsStruct_new_6, refData.elecPointCloudsStruct_ref_6));
for k=1:length(intensityProfiles_new_6)
assert(norm(intensityProfiles_new_6{k} - refData.intensityProfiles_ref_6{k}) < tol);
end
for k=1:length(skelSkelmms_new_6)
assert(norm(skelSkelmms_new_6{k} - refData.skelSkelmms_ref_6{k}) < tol);
end

% load a new CT input
niiCT_PaCER_5 = NiftiMod([inputDataPath filesep 'PaCER_ct_post_OK_5.nii.gz']);

% generate the new output with 'reverseDir', 'contactDetectionMethod' and 'contactAreaCenter' input arguments
[elecModels_new_5, elecPointCloudsStruct_new_5, intensityProfiles_new_5, skelSkelmms_new_5] = PaCER(niiCT_PaCER_5, 'reverseDir', true, 'contactDetectionMethod','contactAreaCenter');

% compare the new generated data against the ref data
structureComparison(refData.elecModels_ref_5, elecModels_new_5)
assert(isequal(elecPointCloudsStruct_new_5, refData.elecPointCloudsStruct_ref_5));
for k=1:length(intensityProfiles_new_5)
assert(norm(intensityProfiles_new_5{k} - refData.intensityProfiles_ref_5{k}) < tol);
end
for k=1:length(skelSkelmms_new_5)
assert(norm(skelSkelmms_new_5{k} - refData.skelSkelmms_ref_5{k}) < tol);
end

%% Use extra dataset to test error message and warning message
% test if PaCER thrown an error when no electrode artifact are found in the CT supplied
w = warning ('off','all');
niiCT_inputs = NiftiMod([inputDataPath filesep 'PaCER_ct_post_OK_2.nii.gz']);
try
[elecModels_new, elecPointCloudsStruct_new, intensityProfiles_new, skelSkelmms_new] = PaCER(niiCT_inputs,'noMask',false);
catch ME
assert(length(ME.message) > 0)
end
w = warning ('on','all');

% catch error message if medtronicXMLPlan is used but no corresponding file is provided
w = warning ('off','all');
niiCT_error_1 = NiftiMod([inputDataPath filesep 'PaCER_postop_error_1.nii.gz']);
try
[elecModels_new, elecPointCloudsStruct_new, intensityProfiles_new, skelSkelmms_new] = PaCER(niiCT_error_1, 'medtronicXMLPlan');
catch ME
assert(length(ME.message) > 0)
end
w = warning ('on','all');

% test the warning messages
% test if PaCER thrown specific warning message depending on the dataset used
warningMessage = ' invPolyArcLength3: given arcLength is negative! Forcing t=0. This is wrong but might be approximatly okay for the use case! Check carefully!';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_error_1}))

% use additional dataset to test warning messages
% load a new CT data
niiCT_error_2 = NiftiMod([inputDataPath filesep 'PaCER_postop_error_2.nii.gz']);

% test if warning message appears when the model is loaded
warningMessage = 'checkNifti: qform_code == sform_code, however the transformation defined in the qform differes from the sform! This might indicate a serious flaw in the nifti header and lead to unexpected results as different tools/algorithms might deal differently with this situation. Fix the nifti header of your file before continuing.';
assert(verifyFunctionWarning('checkNiftiHdr', warningMessage, 'inputs', {niiCT_error_2}))

% test if CT planes in Z direction are not exactly aligned
warningMessage = 'CT planes in Z direction are not exactly aligned. Trying with 0.1 mm tolerance';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_error_2}))

% test if electrode type and contact between them are not detected
warningMessage_1 = 'determineElectrodeType: Could NOT detect electrode type! Electrode contact detection might by flawed. To low image resolution (to large slice thickness)!? Set electrode type manually if you want to continue with this data';
warningMessage_2 = 'Could NOT detect independent electrode contacts. Check image quality. ';
assert(verifyFunctionWarning('PaCER', warningMessage_1, 'inputs', {niiCT_error_2}))
assert(verifyFunctionWarning('PaCER', warningMessage_2, 'inputs', {niiCT_error_2}))

% test if slice thickness is greater than 1 mm
warningMessage = 'Slice thickness is greater than 1 mm! Independent contact detection is most likly not possible. Forcing contactAreaCenter based method.';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_PostOP_new}))

% test if Uncommon fraction of CT data in threshold range
warningMessage = 'Uncommon fraction of CT data in threshold range (15-60 HU). Trying to compensate. Make sure to use "soft tissue" reconstruction filters for the CT (e.g. J30 kernel) if this fails. ';
assert(verifyFunctionWarning('extractBrainConvHull', warningMessage, 'inputs', {niiCT_PostOP_new}))

% test if No electrode specification given
warningMessage = 'No electrode specification given! Set electrodeType option! Trying to estimate type by contactAreaWidth only which might be wrong!';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_PostOP_new}))

% test if PaCER cannot detect independent electrode contacts
warningMessage = 'Could NOT detect independent electrode contacts. Check image quality. ';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_PostOP_new}))

% test if slice thickness is greater than 0.7 mm
warningMessage = 'Slice thickness is greater than 0.7 mm! Independet contact detection might not work reliable in this case. However, for certain electrode types with large contacts spacings you might be lucky.';
assert(verifyFunctionWarning('PaCER', warningMessage, 'inputs', {niiCT_brainMask_new}))

%% change back to the current directory
cd(currentDir);
2 changes: 1 addition & 1 deletion test/testAll.m
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@
else
% on the CI, always reset the path to make absolutely sure, that we test
% the current version
restoredefaultpath()
restoredefaultpath();
launchTestSuite = true;
end

Expand Down