0001 function qi = matRad_calcQualityIndicators(cst,pln,doseCube,refGy,refVol)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040 matRad_cfg = MatRad_Config.instance();
0041
0042 if ~exist('refVol', 'var') || isempty(refVol)
0043 refVol = [2 5 50 95 98];
0044 end
0045
0046 if ~exist('refGy', 'var') || isempty(refGy)
0047 refGy = floor(linspace(0,max(doseCube(:)),6)*10)/10;
0048 end
0049
0050
0051
0052 qi = struct;
0053 for runVoi = 1:size(cst,1)
0054
0055 indices = cst{runVoi,4}{1};
0056 numOfVoxels = numel(indices);
0057 voiPrint = sprintf('%3d %20s',cst{runVoi,1},cst{runVoi,2});
0058
0059
0060 doseInVoi = sort(doseCube(indices));
0061
0062 if ~isempty(doseInVoi)
0063
0064 qi(runVoi).name = cst{runVoi,2};
0065
0066
0067 qi(runVoi).mean = mean(doseInVoi);
0068 qi(runVoi).std = std(doseInVoi);
0069 qi(runVoi).max = doseInVoi(end);
0070 qi(runVoi).min = doseInVoi(1);
0071
0072 voiPrint = sprintf('%s - Mean dose = %5.2f Gy +/- %5.2f Gy (Max dose = %5.2f Gy, Min dose = %5.2f Gy)\n%27s', ...
0073 voiPrint,qi(runVoi).mean,qi(runVoi).std,qi(runVoi).max,qi(runVoi).min,' ');
0074
0075 DX = @(x) matRad_interp1(linspace(0,1,numOfVoxels),doseInVoi,(100-x)*0.01);
0076 VX = @(x) numel(doseInVoi(doseInVoi >= x)) / numOfVoxels;
0077
0078
0079 for runDX = 1:numel(refVol)
0080 qi(runVoi).(strcat('D_',num2str(refVol(runDX)))) = DX(refVol(runDX));
0081 voiPrint = sprintf('%sD%d%% = %5.2f Gy, ',voiPrint,refVol(runDX),DX(refVol(runDX)));
0082 end
0083 voiPrint = sprintf('%s\n%27s',voiPrint,' ');
0084 for runVX = 1:numel(refGy)
0085 sRefGy = num2str(refGy(runVX),3);
0086 qi(runVoi).(['V_' strrep(sRefGy,'.','_') 'Gy']) = VX(refGy(runVX));
0087 voiPrint = sprintf(['%sV' sRefGy 'Gy = %6.2f%%, '],voiPrint,VX(refGy(runVX))*100);
0088 end
0089 voiPrint = sprintf('%s\n%27s',voiPrint,' ');
0090
0091
0092 if strcmp(cst{runVoi,3},'TARGET') > 0
0093
0094
0095 referenceDose = inf;
0096
0097 if isstruct(cst{runVoi,6})
0098 cst{runVoi,6} = num2cell(arrayfun(@matRad_DoseOptimizationFunction.convertOldOptimizationStruct,cst{runVoi,6}));
0099 end
0100
0101 for runObjective = 1:numel(cst{runVoi,6})
0102
0103 obj = cst{runVoi,6}{runObjective};
0104 if ~isa(obj,'matRad_DoseOptimizationFunction')
0105 try
0106 obj = matRad_DoseOptimizationFunction.createInstanceFromStruct(obj);
0107 catch ME
0108 matRad_cfg.dispWarning('Objective/Constraint not valid!\n%s',ME.message)
0109 continue;
0110 end
0111 end
0112
0113
0114 if isa(obj,'DoseObjectives.matRad_SquaredDeviation') || isa(obj,'DoseObjectives.matRad_SquaredUnderdosing')
0115 referenceDose = (min(obj.getDoseParameters(),referenceDose))/pln.numOfFractions;
0116 end
0117 end
0118
0119 if referenceDose == inf
0120 voiPrint = sprintf('%s%s',voiPrint,'Warning: target has no objective that penalizes underdosage, ');
0121 else
0122
0123 StringReferenceDose = regexprep(num2str(round(referenceDose*100)/100),'\D','_');
0124
0125 VTarget95 = sum(doseInVoi >= 0.95*referenceDose);
0126 VTreated95 = sum(doseCube(:) >= 0.95*referenceDose);
0127 qi(runVoi).(['CI_' StringReferenceDose 'Gy']) = VTarget95^2/(numOfVoxels * VTreated95);
0128
0129
0130 qi(runVoi).(['HI_' StringReferenceDose 'Gy']) = (DX(5) - DX(95))/referenceDose * 100;
0131
0132 voiPrint = sprintf('%sCI = %6.4f, HI = %5.2f for reference dose of %3.1f Gy\n',voiPrint,...
0133 qi(runVoi).(['CI_' StringReferenceDose 'Gy']),qi(runVoi).(['HI_' StringReferenceDose 'Gy']),referenceDose);
0134 end
0135 end
0136
0137 matRad_cfg.dispInfo('%s\n',voiPrint);
0138 else
0139 matRad_cfg.dispInfo('%d %s - No dose information.',cst{runVoi,1},cst{runVoi,2});
0140 end
0141 end
0142
0143
0144 listOfFields = fieldnames(qi);
0145 for i = 1:size(cst,1)
0146 indices = cst{i,4}{1};
0147 doseInVoi = sort(doseCube(indices));
0148 if isempty(doseInVoi)
0149 for j = 1:numel(listOfFields)
0150 qi(i).(listOfFields{j}) = NaN;
0151 end
0152 qi(i).name = cst{i,2};
0153 end
0154 end
0155
0156 end
0157