Backup before removing beam_joined optimization

WiscPlanPhotonkV125/matlab_frontend/NLP_optimizer_v2.m

@@ -36,6 +36,7 @@ else
     Pat_path = varargin{1};
     path2geometry = [Pat_path, '\matlab_files\Geometry.mat'];
     path2goal = varargin{2};
+    [Goal_path,Goal_file,ext] = fileparts(path2goal);
 end
 
 str = inputdlg({'N of iterations for initial calc', 'N of iterations for full calc', ...
@@ -207,11 +208,11 @@ D_full = reshape(beamlets * w_fin, size(Geometry.data));
 %% save outputs
 NLP_result.dose = D_full;
 NLP_result.weights = w_fin;
-save([Pat_path, '\matlab_files\NLP_result.mat'], 'NLP_result');
+save([Pat_path, '\matlab_files\NLP_result_' Goal_file '.mat'], 'NLP_result');
 
 
 plot_DVH(Geometry, D_full)
-colorwash(Geometry.data, D_full, [500, 1500], [0, 66]);
+colorwash(Geometry.data, D_full, [500, 1500], [0, 36]);
 
 
 end
@@ -264,6 +265,16 @@ function penalty = eval_f(x, optGoal, nrs_i, sss_i, rgs_i)
                 temp1=max(0, (optGoal{goal_i}.nrs{nrs_i}.sss{sss_i}.rgs{rgs_i}.beamlets_pruned * x)-...
                     (optGoal{goal_i}.nrs{nrs_i}.sss{sss_i}.rgs{rgs_i}.target));
                 d_penalty = 1.0e0 * sum(temp1.*temp1);
+            case 'min_exp'
+                % penalize if achieved dose is lower than target dose
+                temp1=-min(0, (optGoal{goal_i}.nrs{nrs_i}.sss{sss_i}.rgs{rgs_i}.beamlets_pruned * x)-...
+                    (optGoal{goal_i}.nrs{nrs_i}.sss{sss_i}.rgs{rgs_i}.target));
+                d_penalty = 1.0e0 * sum(exp(temp1));
+            case 'max_exp'
+                % penalize if achieved dose is higher than target dose
+                temp1=max(0, (optGoal{goal_i}.nrs{nrs_i}.sss{sss_i}.rgs{rgs_i}.beamlets_pruned * x)-...
+                    (optGoal{goal_i}.nrs{nrs_i}.sss{sss_i}.rgs{rgs_i}.target));
+                d_penalty = 1.0e0 * sum(exp(temp1));
             case 'LeastSquare'
                 % penalize with sum of squares any deviation from target
                 % dose
@@ -302,7 +313,7 @@ function optGoal = make_robust_optGoal(optGoal, RO_params, beamlets);
     % Y - Y>0 moves image down
     % Z - in/out.
     
-    shift_mag = 3; % vox of shift
+    shift_mag = 1; % vox of shift
     nrs_scene_list={[0,0,0]};
 
     

WiscPlanPhotonkV125/matlab_frontend/get_beam_lets.m

@@ -8,7 +8,14 @@ beamlet_cell_array = read_ryan_beamlets(beamlet_batch_filename, 'ryan');
 
 numVox  = numel(Geometry.data);
 numBeamlet = size(beamlet_cell_array,2);
-beamlets = get_beamlets(beamlet_cell_array, numVox);
+
+if size(Geometry.data, 1)<129
+    batchSize = 400;
+else
+    batchSize = 200;
+end
+
+beamlets = get_beamlets(beamlet_cell_array, numVox, batchSize);
 % join beamlets into beams
 load([patient_dir '\all_beams.mat'])
 beamletOrigin=[0 0 0];
@@ -37,10 +44,9 @@ close(wbar2)
 end
 
 % ---- GET BEAMLETS ----
-function beamlets = get_beamlets(beamlet_cell_array, numVox);
+function beamlets = get_beamlets(beamlet_cell_array, numVox, batchSize);
     wbar1 = waitbar(0, 'Creating beamlet array');
     numBeam = size(beamlet_cell_array,2);
-    batchSize=400;
     beamlets = sparse(0, 0);
     for beam_i=1:numBeam
         % for each beam define how much dose it delivers on each voxel

WiscPlanPhotonkV125/matlab_frontend/goal_def_UI.m

@@ -75,9 +75,11 @@ handles.uitable1.ColumnFormat{5}{1} = 'min';
 handles.uitable1.ColumnFormat{5}{2} = 'max';
 handles.uitable1.ColumnFormat{5}{3} = 'min_sq';
 handles.uitable1.ColumnFormat{5}{4} = 'max_sq';
-handles.uitable1.ColumnFormat{5}{5} = 'LeastSquare';
-handles.uitable1.ColumnFormat{5}{6} = 'min_perc_Volume';
-handles.uitable1.ColumnFormat{5}{7} = 'max_perc_Volume';
+handles.uitable1.ColumnFormat{5}{5} = 'min_exp';
+handles.uitable1.ColumnFormat{5}{6} = 'max_exp';
+handles.uitable1.ColumnFormat{5}{7} = 'LeastSquare';
+handles.uitable1.ColumnFormat{5}{8} = 'min_perc_Volume';
+handles.uitable1.ColumnFormat{5}{9} = 'max_perc_Volume';
 
 % == populate the first entry ==
 handles.uitable1.Data = handles.uitable1.Data(1,:);
@@ -188,9 +190,9 @@ for i = 1 : size(handles.uitable1.Data, 1)
     if strcmp (handles.uitable1.Data{i,3}, 'Fixed dose')
         goal_i.D_final = str2double(handles.uitable1.Data{i,4}) * ones(size(Geometry.data));
     elseif strcmp (handles.uitable1.Data{i,3}, 'Dose map')
-        warning('Works only for NRRD (for now)')
-        
+%         warning('Works only for NRRD (for now)')
         dose_in = nrrdread(handles.uitable1.Data{i,4});
+        colorwash(Geometry.data-1000, dose_in, [-500, 500], [0, 1.5*(max(dose_in(:)))], 1, goal_i.name)
         goal_i.D_final = dose_in(goal_i.ROI_idx);
         % loading of data comes here
         % matrix size checks
@@ -247,7 +249,7 @@ if eventdata.Indices(2) == 3;
     mode = handles.uitable1.Data{eventdata.Indices(1), eventdata.Indices(2)}
     switch mode
         case 'Dose map'
-            [FileName,PathName,FilterIndex] = uigetfile('F:\021_WiscPlan_data\FET_repeat_005_1\matlab_files\*.mat', 'Select dose reference file')
+            [FileName,PathName,FilterIndex] = uigetfile('F:\021_WiscPlan_data\FET_repeat_005_1\matlab_files\*.mat', 'Select dose reference file');
             warning('loading of said profile comes here')
             handles.uitable1.Data{eventdata.Indices(1), 4} = [PathName, FileName];