implement anisotropy distribution that is hierarchically sampled in s…

…igma_t/sigma_r and then converted to beta (as it's pre-processed)
sibirrer · Nov 22, 2024 · b85057c · b85057c
1 parent 9a0ea02
commit b85057c
Show file tree

Hide file tree

Showing 2 changed files with 18 additions and 3 deletions.
diff --git a/hierarc/Sampling/Distributions/anisotropy_distributions.py b/hierarc/Sampling/Distributions/anisotropy_distributions.py
@@ -1,6 +1,6 @@
 import numpy as np
 
-_SUPPORTED_DISTRIBUTIONS = ["GAUSSIAN", "GAUSSIAN_SCALED", "NONE"]
+_SUPPORTED_DISTRIBUTIONS = ["GAUSSIAN", "GAUSSIAN_SCALED", "NONE", "GAUSSIAN_TAN_RAD"]
 _SUPPORTED_MODELS = ["OM", "GOM", "const", "NONE"]
 
 
@@ -80,11 +80,14 @@ def draw_anisotropy(
                     "anisotropy parameter is out of bounds of the interpolated range!"
                 )
             # we draw a linear gaussian for 'const' anisotropy and a scaled proportional one for 'OM
-            if self._distribution_function in ["GAUSSIAN", "GAUSSIAN_SCALED"]:
+            if self._distribution_function in ["GAUSSIAN", "GAUSSIAN_SCALED", "GAUSSIAN_TAN_RAD"]:
                 if self._distribution_function in ["GAUSSIAN"]:
                     a_ani_draw = np.random.normal(a_ani, a_ani_sigma)
-                else:
+                elif self._distribution_function in ["GAUSSIAN_SCALED"]:
                     a_ani_draw = np.random.normal(a_ani, a_ani_sigma * a_ani)
+                else:
+                    # convert to beta = 1 - (sigma_t/sigma_r)^2
+                    a_ani_draw = 1 - np.random.normal(a_ani, a_ani_sigma)**2
 
                 if a_ani_draw < self._a_ani_min or a_ani_draw > self._a_ani_max:
                     return self.draw_anisotropy(

diff --git a/test/test_Sampling/test_Distributions/test_anisotropy_distribution.py b/test/test_Sampling/test_Distributions/test_anisotropy_distribution.py
@@ -28,6 +28,14 @@ def setup_method(self):
             kwargs_anisotropy_max=kwargs_anisotropy_max,
         )
 
+        self._ani_dist_tan_rad = AnisotropyDistribution(
+            anisotropy_model=anisotropy_model,
+            anisotropy_sampling=True,
+            distribution_function="GAUSSIAN_TAN_RAD",
+            kwargs_anisotropy_min=kwargs_anisotropy_min,
+            kwargs_anisotropy_max=kwargs_anisotropy_max,
+        )
+
     def test_draw_anisotropy(self):
         kwargs_anisotropy = {
             "a_ani": 1,
@@ -43,6 +51,10 @@ def test_draw_anisotropy(self):
         assert "a_ani" in kwargs_drawn
         assert "beta_inf" in kwargs_drawn
 
+        kwargs_drawn = self._ani_dist_tan_rad.draw_anisotropy(**kwargs_anisotropy)
+        assert "a_ani" in kwargs_drawn
+        assert "beta_inf" in kwargs_drawn
+
         ani_dist = AnisotropyDistribution(
             anisotropy_model="NONE",
             anisotropy_sampling=False,