Updated insensitive blocking layers at detector top and bottom.

Default thickness of the blocking layers is now 0 (ie, no blocking layers) to generate the same results as the previous code without blocking layers. It is recommended to change BLOCKING_LAYER_TOP to 0.0008f cm to have a more realistic DQE(0). Small bug in the layer definition corrected too.
DIDSR · May 25, 2021 · 103114a · 103114a
1 parent 1b00b0b
commit 103114a
Showing 1 changed file with 7 additions and 4 deletions.
diff --git a/MC-GPU_kernel_v1.5b.cu b/MC-GPU_kernel_v1.5b.cu
@@ -37,6 +37,7 @@
 //
 ////////////////////////////////////////////////////////////////////////////////
 
+
 // ** This software is described in the following reference (please cite it in yuor papers):
 //      Andreu Badal, Diksha Sharma, Christian G.Graff, Rongping Zeng, and Aldo Badano, Mammography and breast  
 //      tomosynthesis simulator for virtual clinical trials, Computer Physics Communications 261, p. 107779 (2021) 
@@ -49,8 +50,9 @@
 //   be tracked and might be detected somewhere else. The insensitive top layer causes a measurable
 //   drop in DQE(0), but it does not affect MTF as implemented. Pixel values will be reduced (less  
 //   energy detected per history. [Reference: Zhou et al., Med. Phys. 34, 1098-1109 (2007)]
-#define BLOCKING_LAYER_TOP    0.0008f  // [cm] = 8 micron. Thickness layer closer to source.            !!BLOCKING_LAYER!!
-#define BLOCKING_LAYER_BOTTOM 0.0008f  // [cm] = 8 micron. Thickness layer further from source.         !!BLOCKING_LAYER!!
+#define BLOCKING_LAYER_TOP    0.0000f  // [cm] Thickness layer closer to source.    Example: 0.0008f for a 8 micron layer (0.0 == no layer).  !!BLOCKING_LAYER!!
+#define BLOCKING_LAYER_BOTTOM 0.0000f  // [cm] Thickness layer further from source. Example: 0.0008f for a 8 micron layer (0.0 == no layer).  !!BLOCKING_LAYER!!
+
 
 ////////////////////////////////////////////////////////////////////////////////
 //!  Initialize the image array, ie, set all pixels to zero
@@ -1514,7 +1516,8 @@ __device__ inline void tally_image(float* energy, float3* position, float3* dire
       //      scatter_state=0 for non-scattered, =1 for Compton, =2 for Rayleigh, and =3 for multiple scatter.
 
       // - Do not count the energy deposited inside the top or bottom blocking (dead) layers of the detector (fluorescence still generated).  !!BLOCKING_LAYER!!
-      if ((position->y > BLOCKING_LAYER_BOTTOM) && (position->y < (detector_data_SHARED->scintillator_thickness-BLOCKING_LAYER_TOP)))      // !!BLOCKING_LAYER!!
+      //   After rotation, the detector top layer starts at Y=0 and grows towards positive Y (bottom), with radiation expected from negative Y, moving towards positive Y.
+      if ((position->y > BLOCKING_LAYER_TOP) && (position->y < (detector_data_SHARED->scintillator_thickness-BLOCKING_LAYER_BOTTOM)))      // !!BLOCKING_LAYER!!
 
       atomicAdd(( image +                                                               // Pointer to beginning of image array
                 (int)(*scatter_state) * detector_data_SHARED->total_num_pixels +         // Offset to corresponding scatter image
@@ -1538,7 +1541,7 @@ __device__ inline void tally_image(float* energy, float3* position, float3* dire
         // -- Tally fluorescence energy in the corresponding pixel, unless escaped:        
         position->y = position->y + dist_detector*direction->y;
 
-        if ((position->y > BLOCKING_LAYER_BOTTOM) && (position->y < (detector_data_SHARED->scintillator_thickness-BLOCKING_LAYER_TOP)))   // !!BLOCKING_LAYER!!
+        if ((position->y > BLOCKING_LAYER_TOP) && (position->y < (detector_data_SHARED->scintillator_thickness-BLOCKING_LAYER_BOTTOM)))   // !!BLOCKING_LAYER!!
         {
           position->x = position->x + dist_detector*direction->x;
           pixel_coord_x = __float2int_rd((position->x - detector_data_SHARED->offset.x + 0.5f*detector_data_SHARED->width_X) * detector_data_SHARED->inv_pixel_size_X);  // CUDA intrinsic function converts float to integer rounding down (to minus inf)