Skip to content

Commit bb8ec01

Browse files
JenHardtJenHardt
committed
added RBE test for regular and 4d plus a small bug fix
1 parent 01e3dd8 commit bb8ec01

File tree

2 files changed

+102
-9
lines changed

2 files changed

+102
-9
lines changed

‎matRad/doseCalc/+DoseEngines/matRad_TopasMCEngine.m‎

Lines changed: 12 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -246,13 +246,13 @@ function writeAllFiles(obj,ct,cst,stf,machine,w)
246246
end
247247

248248
% Get alpha beta parameters from bioParam struct
249-
for i = 1:length(obj.bioParameters.AvailableAlphaXBetaX)
250-
if ~isempty(strfind(lower(obj.bioParameters.AvailableAlphaXBetaX{i,2}),'default'))
251-
break
252-
end
249+
if isfield(obj.bioParameters, 'tissuseAlphaX')
250+
obj.bioParameters.AlphaX = obj.bioModel.tissueAlphaX(1);
251+
obj.bioParameters.BetaX = obj.bioModel.tissueBetaX(1);
252+
end
253+
if numel(obj.bioParameters.AlphaX)>1
254+
matRad_cfg.dispWarning('!!! Only a unique alpha/beta ratio supported at the moment. Found multiple, only the first one will be used !!!!');
253255
end
254-
obj.bioParameters.AlphaX = obj.bioParameters.AvailableAlphaXBetaX{5,1}(1);
255-
obj.bioParameters.BetaX = obj.bioParameters.AvailableAlphaXBetaX{5,1}(2);
256256

257257
end
258258
if obj.scorer.LET
@@ -294,8 +294,11 @@ function writeAllFiles(obj,ct,cst,stf,machine,w)
294294
% Save used RBE models
295295
if obj.scorer.RBE
296296
obj.MCparam.RBE_models = obj.scorer.RBE_model;
297-
[obj.MCparam.ax,obj.MCparam.bx] = matRad_getPhotonLQMParameters(cst,prod(ct.cubeDim),obj.MCparam.numOfCtScen);
298-
obj.MCparam.abx(obj.MCparam.bx>0) = obj.MCparam.ax(obj.MCparam.bx>0)./obj.MCparam.bx(obj.MCparam.bx>0);
297+
[obj.MCparam.ax,obj.MCparam.bx] = matRad_getPhotonLQMParameters(obj.cstDoseGrid,prod(ct.cubeDim),obj.VdoseGrid);
298+
obj.MCparam.abx = arrayfun(@(scen) zeros(size(obj.MCparam.bx{scen})), 1:obj.MCparam.numOfCtScen, 'UniformOutput',false);
299+
for scen=1:obj.MCparam.numOfCtScen
300+
obj.MCparam.abx{scen}(obj.MCparam.bx{scen}>0) = obj.MCparam.ax{scen}(obj.MCparam.bx{scen}>0)./obj.MCparam.bx{scen}(obj.MCparam.bx{scen}>0);
301+
end
299302
end
300303

301304
% fill in bixels, rays and beams in case of dij calculation or external calculation
@@ -1380,7 +1383,7 @@ function writeScorers(obj,fID,beamIx)
13801383
end
13811384
fprintf(fID,'\n%s\n\n',scorerTxt);
13821385
scorerInString = {'tabulatedAlpha', 'tabulatedBeta', 'RBE', 'McNamaraAlpha', 'McNamaraBeta', 'WedenbergAlpha', 'WedenbergBeta'};
1383-
for ixWrite = ixToWrite4DInterplay
1386+
for ixWrite = ixToWrite4D
13841387
fprintf(fID,'b:Sc/%s%i/Active = Tf/%s%i/Active/Value\n', scorerInString{ixWrite},PhaseNum,scorerInString{ixWrite},PhaseNum);
13851388
fprintf(fID,'s:Tf/%s%i/Active/Function = "Step"\n', scorerInString{ixWrite}, PhaseNum);
13861389
fprintf(fID,'dv:Tf/%s%i/Active/Times= %i ',scorerInString{ixWrite},PhaseNum, obj.MCparam.cutNumOfBixel(beamIx));

‎test/doseCalc/test_TopasMCEngine.m‎

Lines changed: 90 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -74,6 +74,54 @@
7474
rmdir(folderName,'s'); %clean up
7575
end
7676

77+
78+
function test_TopasMCdoseCalcBasicRBE
79+
% test if all the necessary output files are written vor a couple of cases.
80+
% i am not using the default number of histories for testing her, insted 1e6.
81+
% Because the files are just writen and not simulated so we dont care about simulation time.
82+
% To few histories may result in wierd behavior in the topas interface, i.e if a beam
83+
% recieves no histories because there are not enough to be distributed accros the spots,
84+
% it causes an error
85+
radModes = DoseEngines.matRad_TopasMCEngine.possibleRadiationModes;
86+
matRad_cfg = MatRad_Config.instance();
87+
88+
if moxunit_util_platform_is_octave
89+
confirm_recursive_rmdir(false,'local');
90+
end
91+
92+
for i = 2:numel(radModes)
93+
switch radModes{i}
94+
case 'protons'
95+
RBEmodel = {'mcn', 'wed'};
96+
case {'helium', 'carbon'}
97+
RBEmodel ={'libamtrack','lem'};
98+
end
99+
matRad_cfg = MatRad_Config.instance();
100+
load([radModes{i} '_testData.mat']);
101+
102+
pln.propDoseCalc.engine = 'TOPAS';
103+
pln.propDoseCalc.externalCalculation = 'write';
104+
pln.propDoseCalc.numHistoriesDirect = 1e6;
105+
pln.propDoseCalc.scorer.RBE = true;
106+
pln.propDoseCalc.scorer.RBE_model = RBEmodel;
107+
pln.bioModel = matRad_bioModel(radModes{i},'none');
108+
resultGUI = matRad_calcDoseForward(ct,cst,stf,pln, ones(1,sum([stf(:).totalNumOfBixels])));
109+
110+
folderName = [matRad_cfg.primaryUserFolder filesep 'TOPAS' filesep];
111+
folderName = [folderName stf(1).radiationMode,'_',stf(1).machine,'_',datestr(now, 'dd-mm-yy')];
112+
%check of outputfolder exists
113+
assertTrue(isfolder(folderName));
114+
%check if file in folder existi
115+
assertTrue(isfile([folderName filesep 'matRad_cube.dat']));
116+
assertTrue(isfile([folderName filesep 'matRad_cube.txt']));
117+
assertTrue(isfile([folderName filesep 'MCparam.mat']));
118+
for j = 1:pln.propStf.numOfBeams
119+
assertTrue(isfile([folderName filesep 'beamSetup_matRad_plan_field' num2str(j) '.txt']));
120+
assertTrue(isfile([folderName filesep 'matRad_plan_field' num2str(j) '_run1.txt']));
121+
end
122+
rmdir(folderName,'s'); %clean up
123+
end
124+
77125
function test_TopasMCdoseCalcMultRuns
78126
numOfRuns = 5;
79127
radModes = DoseEngines.matRad_TopasMCEngine.possibleRadiationModes;
@@ -138,7 +186,47 @@
138186
confirm_recursive_rmdir(false,'local');
139187
end
140188

189+
% physical Dose
190+
for i = 1:numel(radModes)
191+
if ~strcmp(radModes{i},'photons')
192+
load([radModes{i} '_testData.mat']);
193+
[ct,cst] = matRad_addMovement(ct, cst,5, numOfPhases,[0 3 0],'dvfType','pull');
194+
pln.bioModel = matRad_bioModel(radModes{i},'none');
195+
resultGUI.w = ones(1,sum([stf(:).totalNumOfBixels]))';
196+
timeSequence = matRad_makeBixelTimeSeq(stf, resultGUI);
197+
timeSequence = matRad_makePhaseMatrix(timeSequence, ct.numOfCtScen, ct.motionPeriod, 'linear');
198+
pln.propDoseCalc.engine = 'TOPAS';
199+
pln.propDoseCalc.externalCalculation = 'write';
200+
pln.propDoseCalc.calc4DInterplay = true;
201+
pln.propDoseCalc.calcTimeSequence = timeSequence;
202+
pln.propDoseCalc.numHistoriesDirect = 1e6;
203+
resultGUI = matRad_calcDoseForward(ct,cst,stf,pln, resultGUI.w);
204+
205+
folderName = [matRad_cfg.primaryUserFolder filesep 'TOPAS' filesep];
206+
folderName = [folderName stf(1).radiationMode,'_',stf(1).machine,'_',datestr(now, 'dd-mm-yy')];
207+
%check of outputfolder exists
208+
assertTrue(isfolder(folderName));
209+
%check if file in folder existi
210+
assertTrue(isfile([folderName filesep 'MCparam.mat']));
211+
for j = 1:pln.propStf.numOfBeams
212+
assertTrue(isfile([folderName filesep 'beamSetup_matRad_plan_field' num2str(j) '.txt']));
213+
assertTrue(isfile([folderName filesep 'matRad_plan_field' num2str(j) '_run1.txt']));
214+
assertTrue(isfile([folderName filesep 'matRad_cube_field' num2str(j) '.txt']));
215+
for k = 1:numOfPhases
216+
assertTrue(isfile([folderName filesep 'matRad_cube' num2str(k) '.dat']));
217+
end
218+
end
219+
rmdir(folderName,'s'); %clean up
220+
end
221+
end
222+
%RBExDose
141223
for i = 1:numel(radModes)
224+
switch radModes{i}
225+
case 'protons'
226+
RBEmodel = {'mcn', 'wed'};
227+
case {'helium', 'carbon'}
228+
RBEmodel ={'libamtrack','lem'};
229+
end
142230
if ~strcmp(radModes{i},'photons')
143231
load([radModes{i} '_testData.mat']);
144232
[ct,cst] = matRad_addMovement(ct, cst,5, numOfPhases,[0 3 0],'dvfType','pull');
@@ -151,6 +239,8 @@
151239
pln.propDoseCalc.calc4DInterplay = true;
152240
pln.propDoseCalc.calcTimeSequence = timeSequence;
153241
pln.propDoseCalc.numHistoriesDirect = 1e6;
242+
pln.propDoseCalc.scorer.RBE = true;
243+
pln.propDoseCalc.scorer.RBE_model = RBEmodel;
154244
resultGUI = matRad_calcDoseForward(ct,cst,stf,pln, resultGUI.w);
155245

156246
folderName = [matRad_cfg.primaryUserFolder filesep 'TOPAS' filesep];

0 commit comments

Comments
 (0)