Cleaned up code per two new sonarcloud issues, i.e. two methods had i…

…dentical code. That is because derivative Radon-Nikodym factor is same for DAW and CAW.
VirtualPhotonics · May 2, 2024 · 906b80b · 906b80b
1 parent 6529314
commit 906b80b
Show file tree

Hide file tree

Showing 2 changed files with 22 additions and 42 deletions.
diff --git a/src/Vts/MonteCarlo/Detectors/dMCdROfRhodMuaDetector.cs b/src/Vts/MonteCarlo/Detectors/dMCdROfRhodMuaDetector.cs
@@ -158,17 +158,15 @@ protected void SetAbsorbAction(AbsorptionWeightingType awt)
             // AbsorptionWeightingType.Analog cannot have derivatives so not a case
             switch (awt)
             {
-                // note: pMC is not applied to analog processing,
-                // only DAW and CAW
+                // note: pMC is not applied to analog processing, only DAW and CAW
                 case AbsorptionWeightingType.Continuous:
-                    _absorbAction = AbsorbContinuous;
-                    break;
                 case AbsorptionWeightingType.Discrete:
                 default:
-                    _absorbAction = AbsorbDiscrete;
+                    _absorbAction = AbsorbContinuousOrDiscrete;
                     break;
             }
         }
+
         /// <summary>
         /// method to tally to detector
         /// </summary>
@@ -191,7 +189,14 @@ public void Tally(Photon photon)
             SecondMoment[ir] += photon.DP.Weight * weightFactor * photon.DP.Weight * weightFactor;
         }
 
-        private double AbsorbContinuous(IList<long> numberOfCollisions, IList<double> pathLength, IList<OpticalProperties> perturbedOps)
+        /// <summary>
+        /// The following method works for both discrete or continuous absorption weighting
+        /// </summary>
+        /// <param name="numberOfCollisions">photon number of collisions in perturbed region</param>
+        /// <param name="pathLength">photon path length in perturbed region</param>
+        /// <param name="perturbedOps">perturbed optical properties of perturbed region</param>
+        /// <returns>derivative with respect to mus</returns>
+        private double AbsorbContinuousOrDiscrete(IList<long> numberOfCollisions, IList<double> pathLength, IList<OpticalProperties> perturbedOps)
         {
             // NOTE: following code only works for single perturbed region because derivative of
             // Radon-Nikodym product needs d(AB)=dA B + A dB and this does not produce that
@@ -206,21 +211,6 @@ private double AbsorbContinuous(IList<long> numberOfCollisions, IList<double> pa
                     numberOfCollisions[i]);
         }
 
-        private double AbsorbDiscrete(IList<long> numberOfCollisions, IList<double> pathLength, IList<OpticalProperties> perturbedOps)
-        {
-            // NOTE: following code only works for single perturbed region because derivative of
-            // Radon-Nikodym product needs d(AB)=dA B + A dB and this does not produce that.
-            // Check for only one perturbedRegionIndices specified by user performed in DataStructuresValidation
-            var i = _perturbedRegionsIndices[0]; 
-            // rearranged to be more numerically stable
-            return -pathLength[i] * // dMua* factor 
-                Math.Pow(
-                    _perturbedOps[i].Mus / _referenceOps[i].Mus *
-                    Math.Exp(-(_perturbedOps[i].Mus + _perturbedOps[i].Mua - _referenceOps[i].Mus - _referenceOps[i].Mua) *
-                        pathLength[i] / numberOfCollisions[i]),
-                    numberOfCollisions[i]);
-        }
-
         /// <summary>
         /// method to normalize detector results after numPhotons launched
         /// </summary>

diff --git a/src/Vts/MonteCarlo/Detectors/dMCdROfRhodMusDetector.cs b/src/Vts/MonteCarlo/Detectors/dMCdROfRhodMusDetector.cs
@@ -157,17 +157,15 @@ protected void SetAbsorbAction(AbsorptionWeightingType awt)
             // AbsorptionWeightingType.Analog cannot have derivatives so not a case
             switch (awt)
             {
-                // note: pMC is not applied to analog processing,
-                // only DAW and CAW
+                // note: pMC is not applied to analog processing, only DAW and CAW
                 case AbsorptionWeightingType.Continuous:
-                    _absorbAction = AbsorbContinuous;
-                    break;
                 case AbsorptionWeightingType.Discrete:
                 default:
-                    _absorbAction = AbsorbDiscrete;
+                    _absorbAction = AbsorbContinuousOrDiscrete;
                     break;
             }
         }
+
         /// <summary>
         /// method to tally to detector
         /// </summary>
@@ -190,7 +188,14 @@ public void Tally(Photon photon)
             SecondMoment[ir] += photon.DP.Weight * weightFactor * photon.DP.Weight * weightFactor;
         }
 
-        private double AbsorbContinuous(IList<long> numberOfCollisions, IList<double> pathLength, IList<OpticalProperties> perturbedOps)
+        /// <summary>
+        /// The following method works for both discrete or continuous absorption weighting
+        /// </summary>
+        /// <param name="numberOfCollisions">photon number of collisions in perturbed region</param>
+        /// <param name="pathLength">photon path length in perturbed region</param>
+        /// <param name="perturbedOps">perturbed optical properties of perturbed region</param>
+        /// <returns>derivative with respect to mua</returns>
+        private double AbsorbContinuousOrDiscrete(IList<long> numberOfCollisions, IList<double> pathLength, IList<OpticalProperties> perturbedOps)
         {
             // NOTE: following code only works for single perturbed region because derivative of
             // Radon-Nikodym product needs d(AB)=dA B + A dB and this does not produce that
@@ -205,21 +210,6 @@ private double AbsorbContinuous(IList<long> numberOfCollisions, IList<double> pa
                         numberOfCollisions[i]);
         }
 
-        private double AbsorbDiscrete(IList<long> numberOfCollisions, IList<double> pathLength, IList<OpticalProperties> perturbedOps)
-        {
-            // NOTE: following code only works for single perturbed region because derivative of
-            // Radon-Nikodym product needs d(AB)=dA B + A dB and this does not produce that
-            // Check for only one perturbedRegionIndices specified by user performed in DataStructuresValidation
-            var i = _perturbedRegionsIndices[0];
-            // rearranged to be more numerically stable
-            return (numberOfCollisions[i]/ _perturbedOps[i].Mus -pathLength[i] )* // dMus* factor 
-                    Math.Pow(_perturbedOps[i].Mus / _referenceOps[i].Mus *
-                        Math.Exp(-(_perturbedOps[i].Mus + _perturbedOps[i].Mua - 
-                                   _referenceOps[i].Mus - _referenceOps[i].Mua) *
-                            pathLength[i] / numberOfCollisions[i]),
-                        numberOfCollisions[i]);
-        }
-
         /// <summary>
         /// method to normalize detector results after numPhotons launched
         /// </summary>