matRad IPOPT callback: get jacobian structure for direct aperture optimization call jacobStruct = matRad_daoGetJacobStruct(optiProb,apertureInfoVec,dij,cst) input optiProb: option struct defining the type of optimization apertureInfoVect: aperture weights and shapes parameterized as vector dij: dose influence matrix cst: matRad cst struct output jacobStruct: jacobian of constraint function References %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0001 function jacobStruct = matRad_getJacobianStructure(optiProb,apertureInfoVec,dij,cst) 0002 % matRad IPOPT callback: get jacobian structure for direct aperture optimization 0003 % 0004 % call 0005 % jacobStruct = matRad_daoGetJacobStruct(optiProb,apertureInfoVec,dij,cst) 0006 % 0007 % input 0008 % optiProb: option struct defining the type of optimization 0009 % apertureInfoVect: aperture weights and shapes parameterized as vector 0010 % dij: dose influence matrix 0011 % cst: matRad cst struct 0012 % 0013 % output 0014 % jacobStruct: jacobian of constraint function 0015 % 0016 % References 0017 % 0018 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0019 0020 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0021 % 0022 % Copyright 2015 the matRad development team. 0023 % 0024 % This file is part of the matRad project. It is subject to the license 0025 % terms in the LICENSE file found in the top-level directory of this 0026 % distribution and at https://github.com/e0404/matRad/LICENSES.txt. No part 0027 % of the matRad project, including this file, may be copied, modified, 0028 % propagated, or distributed except according to the terms contained in the 0029 % LICENSE file. 0030 % 0031 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0032 0033 apertureInfo = optiProb.apertureInfo; 0034 0035 % jacobian structure of the dao constraints 0036 % row indices 0037 i = repmat(1:apertureInfo.totalNumOfLeafPairs,1,2); 0038 % column indices 0039 j = [apertureInfo.totalNumOfShapes+1:apertureInfo.totalNumOfShapes+apertureInfo.totalNumOfLeafPairs ... 0040 apertureInfo.totalNumOfShapes+apertureInfo.totalNumOfLeafPairs+1:apertureInfo.totalNumOfShapes+2*apertureInfo.totalNumOfLeafPairs]; 0041 0042 % -1 for left leaves, 1 for right leaves 0043 s = ones(1,2*apertureInfo.totalNumOfLeafPairs); 0044 0045 jacobStruct_dao = sparse(i,j,s, ... 0046 apertureInfo.totalNumOfLeafPairs, ... 0047 apertureInfo.totalNumOfShapes+2*apertureInfo.totalNumOfLeafPairs, ... 0048 2*apertureInfo.totalNumOfLeafPairs); 0049 0050 jacobStruct_dos_bixel = matRad_getJacobianStructure@matRad_OptimizationProblem(optiProb,apertureInfo.bixelWeights,dij,cst); 0051 % --> gives me a matrix with number of rows = num of constraints and tells 0052 % me in th columns if a beamlet has an influence on this constraint 0053 0054 % for apertures I need to check if the very beam orientation of the aperture has a bixel 0055 % that potentially influences the constraint 0056 0057 % for leaves I need to check if that particular leaf row has bixels that 0058 % potentially influence the objective which works via apertureInfo.beam(i).bixelIndMap 0059 0060 % all stuff can be done per beam direction and then I use repmat to build 0061 % up the big matrix 0062 0063 % allocate 0064 jacobStruct_dos = sparse(size(jacobStruct_dos_bixel,1),size(jacobStruct_dao,2)); 0065 0066 if ~isempty(jacobStruct_dos) 0067 0068 % first aperture weights 0069 for i = 1:apertureInfo.totalNumOfShapes 0070 currBeam = apertureInfo.mappingMx(i,1); 0071 currBixelIxInBeam = dij.beamNum == currBeam; 0072 jacobStruct_dos(:,i) = spones(sum(jacobStruct_dos_bixel(:,currBixelIxInBeam),2)); 0073 end 0074 0075 % second leaves 0076 counter = 0; 0077 for i = 1:size(apertureInfo.beam,2) 0078 for j = 1:apertureInfo.beam(i).numOfShapes 0079 for k = 1:apertureInfo.beam(i).numOfActiveLeafPairs 0080 counter = counter + 1; 0081 bixelIxInCurrRow = ~isnan(apertureInfo.beam(i).bixelIndMap(k,:)); 0082 jacobStruct_dos(:,counter+[0 apertureInfo.totalNumOfLeafPairs]) = ... 0083 repmat(spones(sum(jacobStruct_dos_bixel(:,bixelIxInCurrRow),2)),1,2); 0084 end 0085 end 0086 end 0087 0088 end 0089 0090 % concatenate 0091 jacobStruct = [jacobStruct_dao; jacobStruct_dos];