0001 function resultGUI = matRad_calcCubes(w,dij,scenNum)
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 if nargin < 3
0034 scenNum = 1;
0035 end
0036
0037 resultGUI.w = w;
0038
0039
0040 for i = 1:dij.numOfBeams
0041 beamInfo(i).suffix = ['_beam', num2str(i)];
0042 beamInfo(i).logIx = (dij.beamNum == i);
0043 end
0044 beamInfo(dij.numOfBeams+1).suffix = '';
0045 beamInfo(dij.numOfBeams+1).logIx = true(size(w));
0046
0047
0048
0049 for i = 1:length(beamInfo)
0050 resultGUI.(['physicalDose', beamInfo(i).suffix]) = reshape(full(dij.physicalDose{scenNum} * (resultGUI.w .* beamInfo(i).logIx)),dij.doseGrid.dimensions);
0051 end
0052
0053
0054 if isfield(dij,'RBE')
0055 fprintf(['matRad: applying a constant RBE of ' num2str(dij.RBE) ' \n']);
0056 for i = 1:length(beamInfo)
0057 resultGUI.(['RBExDose', beamInfo(i).suffix]) = resultGUI.(['physicalDose', beamInfo(i).suffix]) * dij.RBE;
0058 end
0059 end
0060
0061
0062 if isfield(dij,'mLETDose')
0063 for i = 1:length(beamInfo)
0064 LETDoseCube = dij.mLETDose{scenNum} * (resultGUI.w .* beamInfo(i).logIx);
0065 resultGUI.(['LET', beamInfo(i).suffix]) = zeros(dij.doseGrid.dimensions);
0066 ix = resultGUI.(['physicalDose', beamInfo(i).suffix]) > 0;
0067 resultGUI.(['LET', beamInfo(i).suffix])(ix) = LETDoseCube(ix)./resultGUI.(['physicalDose', beamInfo(i).suffix])(ix);
0068 end
0069 end
0070
0071 if isfield(dij,'physicalDose_MCvar')
0072 resultGUI.physicalDose_MCvar = reshape(full(dij.physicalDose_MCvar{scenNum} * (resultGUI.w .* beamInfo(i).logIx)),dij.doseGrid.dimensions);
0073 resultGUI.physicalDose_MCstd = sqrt(resultGUI.physicalDose_MCvar);
0074 resultGUI.physicalDose_MCstdRel = resultGUI.physicalDose_MCstd ./ resultGUI.physicalDose;
0075 end
0076
0077
0078 if isfield(dij,'mAlphaDose') && isfield(dij,'mSqrtBetaDose')
0079
0080 ix = dij.bx~=0;
0081
0082 for i = 1:length(beamInfo)
0083 wBeam = (resultGUI.w .* beamInfo(i).logIx);
0084 resultGUI.(['effect', beamInfo(i).suffix]) = full(dij.mAlphaDose{scenNum} * wBeam + (dij.mSqrtBetaDose{scenNum} * wBeam).^2);
0085 resultGUI.(['effect', beamInfo(i).suffix]) = reshape(resultGUI.(['effect', beamInfo(i).suffix]),dij.doseGrid.dimensions);
0086
0087 resultGUI.(['RBExDose', beamInfo(i).suffix]) = zeros(size(resultGUI.(['effect', beamInfo(i).suffix])));
0088 resultGUI.(['RBExDose', beamInfo(i).suffix])(ix) = (sqrt(dij.ax(ix).^2 + 4 .* dij.bx(ix) .* resultGUI.(['effect', beamInfo(i).suffix])(ix)) - dij.ax(ix))./(2.*dij.bx(ix));
0089
0090 resultGUI.(['RBE', beamInfo(i).suffix]) = resultGUI.(['RBExDose', beamInfo(i).suffix])./resultGUI.(['physicalDose', beamInfo(i).suffix]);
0091
0092 resultGUI.(['alpha', beamInfo(i).suffix]) = zeros(dij.doseGrid.dimensions);
0093 resultGUI.(['beta', beamInfo(i).suffix]) = zeros(dij.doseGrid.dimensions);
0094
0095 AlphaDoseCube = full(dij.mAlphaDose{scenNum} * wBeam);
0096 resultGUI.(['alpha', beamInfo(i).suffix])(ix) = AlphaDoseCube(ix)./resultGUI.(['physicalDose', beamInfo(i).suffix])(ix);
0097
0098 SqrtBetaDoseCube = full(dij.mSqrtBetaDose{scenNum} * wBeam);
0099 resultGUI.(['beta', beamInfo(i).suffix])(ix) = (SqrtBetaDoseCube(ix)./resultGUI.(['physicalDose', beamInfo(i).suffix])(ix)).^2;
0100 end
0101 end
0102
0103
0104 resultGUI = orderfields(resultGUI);
0105
0106
0107 if any(dij.ctGrid.dimensions~=dij.doseGrid.dimensions)
0108 myFields = fieldnames(resultGUI);
0109 for i = 1:numel(myFields)
0110
0111 if numel(resultGUI.(myFields{i})) == dij.doseGrid.numOfVoxels
0112
0113
0114 resultGUI.(myFields{i}) = matRad_interp3(dij.doseGrid.x,dij.doseGrid.y',dij.doseGrid.z, ...
0115 resultGUI.(myFields{i}), ...
0116 dij.ctGrid.x,dij.ctGrid.y',dij.ctGrid.z,'linear',0);
0117
0118 end
0119
0120 end
0121 end
0122
0123 end
0124
0125
0126