From 462fc4f5c524ab49bd71c852bd4b530f181f0b52 Mon Sep 17 00:00:00 2001
From: Aaron LI <aly@aaronly.me>
Date: Sat, 29 Jul 2017 15:59:24 +0800
Subject: galactic/freefree: Add "dust_fraction" and "halpha_abs_th" configs

Signed-off-by: Aaron LI <aly@aaronly.me>
---
 fg21sim/configs/10-galactic.conf.spec | 12 ++++++++++++
 fg21sim/galactic/freefree.py          | 16 ++++++----------
 2 files changed, 18 insertions(+), 10 deletions(-)

diff --git a/fg21sim/configs/10-galactic.conf.spec b/fg21sim/configs/10-galactic.conf.spec
index b082a5d..116586f 100644
--- a/fg21sim/configs/10-galactic.conf.spec
+++ b/fg21sim/configs/10-galactic.conf.spec
@@ -62,6 +62,18 @@
   # The unit of the above dust map (e.g., "MJy/sr")
   dustmap_unit = string(default=None)
 
+  # Effective dust fraction in the LoS actually absorbing Halpha
+  dust_fraction = float(default=0.33, min=0.1, max=1.0)
+
+  # Halpha absorption threshold:
+  # When the dust absorption goes rather large, the true Halpha
+  # absorption can not well determined.  This configuration sets the
+  # threshold below which the dust absorption can be well determined,
+  # while the sky regions with higher absorption are masked out due
+  # to unreliable absorption correction.
+  # Unit: [mag]
+  halpha_abs_th = float(default=1.0)
+
   # The electron temperature assumed for the ionized interstellar medium
   # that generating H{\alpha} emission.
   # Unit: [K]
diff --git a/fg21sim/galactic/freefree.py b/fg21sim/galactic/freefree.py
index 2ea193c..9c9424f 100644
--- a/fg21sim/galactic/freefree.py
+++ b/fg21sim/galactic/freefree.py
@@ -75,6 +75,8 @@ class FreeFree:
         self.dustmap_path = self.configs.get_path(comp+"/dustmap")
         self.dustmap_unit = au.Unit(
             self.configs.getn(comp+"/dustmap_unit"))
+        self.f_dust = self.configs.getn(comp+"/dust_fraction")
+        self.halpha_abs_th = self.configs.getn(comp+"/halpha_abs_th")  # [mag]
         self.Te = self.configs.getn(comp+"/electron_temperature")  # [K]
         self.prefix = self.configs.getn(comp+"/prefix")
         self.save = self.configs.getn(comp+"/save")
@@ -117,29 +119,23 @@ class FreeFree:
             return
         #
         logger.info("Correct H[alpha] map for dust absorption")
-        # Effective dust fraction in the LoS actually absorbing Halpha
-        f_dust = 0.33
-        logger.info("Effective dust fraction: {0}".format(f_dust))
-        # NOTE:
+        logger.info("Effective dust fraction: {0}".format(self.f_dust))
         # Mask the regions where the true Halpha absorption is uncertain.
         # When the dust absorption goes rather large, the true Halpha
         # absorption can not well determined.
-        # Therefore, the regions where the calculated Halpha absorption
-        # greater than 1.0 mag are masked out.
-        halpha_abs_th = 1.0  # Halpha absorption threshold, unit: [ mag ]
         # Corresponding dust absorption threshold, unit: [ MJy / sr ]
-        dust_abs_th = halpha_abs_th / 0.0462 / f_dust
+        dust_abs_th = self.halpha_abs_th / 0.0462 / self.f_dust
         logger.info("Dust absorption mask threshold: " +
                     "{0:.1f} MJy/sr ".format(dust_abs_th) +
                     "<-> H[alpha] absorption threshold: " +
-                    "{0:.1f} mag".format(halpha_abs_th))
+                    "{0:.1f} mag".format(self.halpha_abs_th))
         mask = (self.dustmap.data > dust_abs_th)
         self.dustmap.data[mask] = np.nan
         fp_mask = 100 * mask.sum() / self.dustmap.data.size
         logger.warning("Dust map masked fraction: {0:.1f}%".format(fp_mask))
         #
         halphamap_corr = (self.halphamap.data *
-                          10**(self.dustmap.data * 0.0185 * f_dust))
+                          10**(self.dustmap.data * 0.0185 * self.f_dust))
         self.halphamap.data = halphamap_corr
         self._dust_corrected = True
         logger.info("Done dust absorption correction")
-- 
cgit v1.2.2