1
0

get_beam_lets.m 2.7 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889
  1. function [beamlets, beamlets_joined, numBeamlet, numBeam, beam_i_list] = get_beam_lets(Geometry, patient_dir)
  2. % this function loads and returns beamlets and joined beams
  3. beamlet_batch_filename = [patient_dir '\' 'batch_dose.bin'];
  4. beamlet_cell_array = read_ryan_beamlets(beamlet_batch_filename, 'ryan');
  5. numVox = numel(Geometry.data);
  6. numBeamlet = size(beamlet_cell_array,2);
  7. if size(Geometry.data, 1)<129
  8. batchSize = 400;
  9. else
  10. batchSize = 200;
  11. end
  12. beamlets = get_beamlets(beamlet_cell_array, numVox, batchSize);
  13. % join beamlets into beams
  14. load([patient_dir '\all_beams.mat'])
  15. beamletOrigin=[0 0 0];
  16. beam_i=0;
  17. beam_i_list=[];
  18. for beamlet_i = 1:numel(all_beams)
  19. newBeamletOrigin = all_beams{beamlet_i}.ip;
  20. if any(newBeamletOrigin ~= beamletOrigin)
  21. beam_i = beam_i+1;
  22. beamletOrigin = newBeamletOrigin;
  23. end
  24. beam_i_list=[beam_i_list, beam_i];
  25. end
  26. wbar2 = waitbar(0, 'merging beamlets into beams');
  27. numBeam=numel(unique(beam_i_list));
  28. beamlets_joined=sparse(size(beamlets,1), numBeam);
  29. for beam_i = 1:numel(unique(beam_i_list))
  30. beam_full = sum(beamlets(:,beam_i_list == beam_i), 2);
  31. beamlets_joined(:,beam_i) = beam_full;
  32. waitbar(beam_i/numBeam, wbar2)
  33. end
  34. close(wbar2)
  35. end
  36. % ---- GET BEAMLETS ----
  37. function beamlets = get_beamlets(beamlet_cell_array, numVox, batchSize);
  38. wbar1 = waitbar(0, 'Creating beamlet array');
  39. numBeam = size(beamlet_cell_array,2);
  40. beamlets = sparse(0, 0);
  41. for beam_i=1:numBeam
  42. % for each beam define how much dose it delivers on each voxel
  43. idx=beamlet_cell_array{1, beam_i}.non_zero_indices;
  44. % break the beamlets into multiple batches
  45. if rem(beam_i, batchSize)==1;
  46. beamlet_batch = sparse(numVox, batchSize);
  47. beam_i_temp=1;
  48. end
  49. beamlet_batch(idx, beam_i_temp) = 1000*beamlet_cell_array{1, beam_i}.non_zero_values;
  50. waitbar(beam_i/numBeam, wbar1, ['Adding beamlet array: #', num2str(beam_i)])
  51. % add the batch to full set when filled
  52. if rem(beam_i, batchSize)==0;
  53. beamlets =[beamlets, beamlet_batch];
  54. end
  55. % crop and add the batch to full set when completed
  56. if beam_i==numBeam;
  57. beamlet_batch=beamlet_batch(:, 1:beam_i_temp);
  58. beamlets =[beamlets, beamlet_batch];
  59. end
  60. beam_i_temp=beam_i_temp+1;
  61. end
  62. close(wbar1)
  63. end
  64. function show_joint_beamlets(beamlets, IMGsize, Beam_list)
  65. % this function overlays and plots multiple beamlets. The goal is to
  66. % check whether some beamlets are part of the same beam manually.
  67. beam=sum(beamlets(:, Beam_list), 2);
  68. beamImg=reshape(full(beam), IMGsize);
  69. orthoslice(beamImg)
  70. end