Merge remote-tracking branch 'refs/remotes/origin/master'

WiscPlanPhotonkV125/matlab_frontend/get_beamlets.m

@@ -12,18 +12,15 @@ ROI_idxList = unique(ROI_idxList);
 %% add idx lists for each goal!
 
 beamlet_batch_filename = [patient_dir '\' 'batch_dose.bin'];
-beamlet_cell_array = read_ryan_beamlets(beamlet_batch_filename, 'ryan');
+% beamlet_cell_array = read_ryan_beamlets(beamlet_batch_filename, 'ryan');
+beamlet_cell_array = read_ryan_beamlets(beamlet_batch_filename, 'peter', ROI_idxList);
 
 numVox  = numel(Geometry.data);
 numBeamlet = size(beamlet_cell_array,2);
-
-if size(Geometry.data, 1)<129
-    batchSize = 400;
-else
-    batchSize = 200;
-end
-
+batchSize = 200;
 beamlets = get_beamlets2(beamlet_cell_array, numVox, batchSize, ROI_idxList);
+% beamlets = beamlet_cell_array;
+
 % --- join beamlets into beams
 % load([patient_dir '\all_beams.mat'])
 % beamletOrigin=[0 0 0];
@@ -58,9 +55,13 @@ function beamlets = get_beamlets2(beamlet_cell_array, numVox, batchSize, ROI_idx
     beamlets = sparse(0, 0);
     for beam_i=1:numBeam
         % for each beam define how much dose it delivers on each voxel
-        idx=beamlet_cell_array{1, beam_i}.non_zero_indices;
+%         idx=beamlet_cell_array{1, beam_i}.non_zero_indices;
+%         [ROI_idxIntersect, ia, ~] = intersect (idx, ROI_idxList);
+%         if ~isequal(ROI_idxIntersect,idx)
+%             warning('this')
+%         end
         
-        [ROI_idxIntersect, ia, ~] = intersect (idx, ROI_idxList);
+        ROI_idxIntersect = beamlet_cell_array{1, beam_i}.non_zero_indices;
         
         if isempty(ROI_idxIntersect)
             warning(['Beamlet ' num2str(beam_i) ' is empty!'])
@@ -73,7 +74,7 @@ function beamlets = get_beamlets2(beamlet_cell_array, numVox, batchSize, ROI_idx
         end
 
 %         beamlet_batch(idx, beam_i_temp) = 1000*beamlet_cell_array{1, beam_i}.non_zero_values;
-        beamlet_batch(ROI_idxIntersect, beam_i_temp) = beamlet_cell_array{1, beam_i}.non_zero_values(ia);
+        beamlet_batch(ROI_idxIntersect, beam_i_temp) = beamlet_cell_array{1, beam_i}.non_zero_values;
         waitbar(beam_i/numBeam, wbar1, ['Adding beamlet array: #', num2str(beam_i)])
 
         % add the batch to full set when filled

WiscPlanPhotonkV125/matlab_frontend/read_ryan_beamlets.m

@@ -94,6 +94,40 @@ switch lower(mode)
                 results = 0;
             end
         end
+    case 'peter'
+        beamlets = cell(1,Nbeamlets);
+        ROI_idxList = weights;
+        wbar1 = waitbar(0, 'Creating beamlet array');
+        for k=1:Nbeamlets
+            waitbar(k/Nbeamlets, wbar1, ['Adding beamlet array: #', num2str(k)])
+            beamlet.num = fread(fid,1,'int');
+            beamlet.x_count = fread(fid,1,'int');
+            beamlet.y_count = fread(fid,1,'int');
+            beamlet.z_count = fread(fid,1,'int');
+            Nind = fread(fid,1,'int');  % number of non-zero values
+            non_zero_indices = int32(fread(fid,Nind,'int'));
+            non_zero_values = single(fread(fid,Nind,'float'));
+            % add unity to the indices to convert them from C to Matlab style indexing:
+            non_zero_indices = non_zero_indices + 1;
+            beamlet.num = beamlet.num + 1;
+            if any(non_zero_indices < 0)
+                error('index negative');
+            end
+            if any(non_zero_indices > beamlet.x_count*beamlet.y_count*beamlet.z_count)
+                error('index out of range');
+            end
+            
+            % load only the voxels within provided ROI
+            [ROI_idxIntersect, ia, ~] = intersect (non_zero_indices, ROI_idxList);
+            beamlet.non_zero_indices = ROI_idxIntersect;
+            beamlet.non_zero_values = non_zero_values(ia);
+                
+            beamlets{k} = beamlet;
+            
+        end
+        results = beamlets;
+        close(wbar1)
+        
     case 'sparse'
         nzmax = 0;
         nBmax = 500;