Parcourir la source

Full RO is kinda working! Still need to verify/debug some thigns, but at least the function runs through.

MEDPHYSICS\pferjancic il y a 3 ans
Parent
commit
685fd94ba8
1 fichiers modifiés avec 30 ajouts et 28 suppressions
  1. 30 28
      WiscPlanPhotonkV125/matlab_frontend/helicalDosecalcSetup7_fullRO.m

+ 30 - 28
WiscPlanPhotonkV125/matlab_frontend/helicalDosecalcSetup7_fullRO.m

@@ -93,7 +93,7 @@ fieldWidth = ypmax - ypmin;
 Nrot = ceil(abs(zBow - zStern)/(pitch*fieldWidth));
 
 % Nphi = Nrot*51;  % number of angles used in the calculation
-Nphi = Nrot * N_angles;  % Grozomah
+Nphi = Nrot * N_angles;
 
 % define the limits of the angles that will be used for the calculation
 % ##phimin = 0;  % starting angle in radians
@@ -101,13 +101,16 @@ Nphi = Nrot * N_angles;  % Grozomah
 
 phi = [0:Nphi-1]/Nphi *2*pi*Nrot;
 
+Geometry.start_nominal = Geometry.start;
+
 %% account for beamlet shift
 for scenario_i = 1:numel(OptGoals.sss_scene_list)
     patient_dir = [beamlet_dir '\scenario' num2str(scenario_i)];
     mkdir(patient_dir)
     condor_folder = patient_dir;
     winxp_folder = 'winxp';
-
+    
+    Geometry.start = Geometry.start_nominal;
     % create names for condor input and output folders
     input_folder = '.';
     output_folder = '.';
@@ -185,30 +188,31 @@ for scenario_i = 1:numel(OptGoals.sss_scene_list)
     tumorMask(Geometry.ROIS{ptvInd}.ind) = 1;
 
     BW = bwdist(tumorMask);
-    tumorMaskExp = tumorMask;
+    tumorMaskExp = tumorMask;  
     tumorMaskExp(BW <= 4) = 1;
 
-    P = zeros(Mxp,Nphi);
-
-    fprintf('Checking beam''s eye view ...\n');
-    for p=1:Nphi
-        % ir and jr form the beam's eye view (BEV)
-        ir = [-sin(phi(p)); cos(phi(p)); 0];
-        jr = [0 0 1]';
-        % kr denotes the beam direction
-        kr = [cos(phi(p)); sin(phi(p)); 0];
-
-        for m=1:Mxp
-            point1 = single(-kr*SAD + [0 0 zBow + pitch*fieldWidth*phi(p)/(2*pi)]');  % source point
-            point2 = single(point1 + (SAD*kr + ir*xp(m))*10);
-            [indVisited,deffVisited] = singleRaytraceClean(tumorMaskExp,START,INC,point1,point2);
-            if ~isempty(indVisited)
-                P(m,p) = max(deffVisited);
+    % only do this for the nominal scenario, leave same beams for others.
+    if scenario_i == 1
+        P = zeros(Mxp,Nphi);
+        fprintf('Checking beam''s eye view ...\n');
+        for p=1:Nphi
+            % ir and jr form the beam's eye view (BEV)
+            ir = [-sin(phi(p)); cos(phi(p)); 0];
+            jr = [0 0 1]';
+            % kr denotes the beam direction
+            kr = [cos(phi(p)); sin(phi(p)); 0];
+            for m=1:Mxp
+                point1 = single(-kr*SAD + [0 0 zBow + pitch*fieldWidth*phi(p)/(2*pi)]');  % source point
+                point2 = single(point1 + (SAD*kr + ir*xp(m))*10);
+                [indVisited,deffVisited] = singleRaytraceClean(tumorMaskExp,START,INC,point1,point2);
+                if ~isempty(indVisited)
+                    P(m,p) = max(deffVisited);
+                end
             end
         end
-    end
-    fprintf('Finished checking BEV\n');
+        fprintf('Finished checking BEV\n');
 
+    end
     % load data required for the dose calculator
     load(kernel_file);
 
@@ -231,18 +235,16 @@ for scenario_i = 1:numel(OptGoals.sss_scene_list)
     end
 
     % account for isocenter
-    Geometry.start_nominal = single(Geometry.start - iso);
+    Geometry.start = single(Geometry.start - iso);
 
     % change Condor folder names as appropriate
     
-    
-    
     % do the isocenter shift
     shift = OptGoals.sss_scene_list{scenario_i}; % Y X Z
-    iso = [iso(1)+Geometry.voxel_size(1)*shift(1) ...
-        iso(2)+Geometry.voxel_size(2)*shift(2) ...
-        iso(3)+Geometry.voxel_size(3)*shift(3)];
-    Geometry.start = Geometry.start_nominal- iso;
+    shiftVec = [Geometry.voxel_size(1)*shift(1) ...
+        Geometry.voxel_size(2)*shift(2) ...
+        Geometry.voxel_size(3)*shift(3)];
+    Geometry.start = Geometry.start- shiftVec;
 
     % find the total number of beams
     Nbeam = Nphi*Mxp*Nyp;