M31 refactor

frb/halos/models.py

@@ -1,7 +1,5 @@
 """ Module for DM Halo calculations
 """
-from __future__ import print_function, absolute_import, division, unicode_literals
-
 import numpy as np
 import pdb
 import warnings
@@ -544,11 +542,14 @@ class ModifiedNFW(object):
 
     Parameters:
         log_Mhalo: float, optional
-          log10 of the Halo mass (solar masses)
+          log10 of the Halo mass [dark matter + baryons] (solar masses)
+        log_MCGM: float, optional
+            log10 of the CGM mass (solar masses)
+            If provided, this sets f_hot [this is its only use]
         c: float, optional
           concentration of the halo
         f_hot: float, optional
-          Fraction of the baryons in this hot phase
+          Fraction of the total "expected" baryons in this hot (aka CGM) phase
           Will likely use this for all diffuse gas
         alpha: float, optional
           Parameter to modify NFW profile power-law
@@ -572,9 +573,21 @@ class ModifiedNFW(object):
 
 
     """
-    def __init__(self, log_Mhalo=12.2, c=7.67, f_hot=0.75, alpha=0.,
-                 y0=1., z=0., cosmo=cosmo, **kwargs):
+    def __init__(self, 
+                 log_Mhalo=12.2, 
+                 c=7.67, 
+                 f_hot:float=None, #0.75, 
+                 log_MCGM:float=None, #1e12,
+                 alpha=0., y0=1., z=0., cosmo=cosmo,fb = cosmo.Ob0/cosmo.Om0,
+                 **kwargs):
         # Init
+        if f_hot is None: 
+            if log_MCGM is None:
+                f_hot = 0.75
+                log_MCGM = np.log10(f_hot * 10**log_Mhalo*fb)
+            else:
+                f_hot = 10**log_MCGM / (10**log_Mhalo*fb)
+
         # Param
         self.log_Mhalo = log_Mhalo
         self.M_halo = 10.**self.log_Mhalo * constants.M_sun.cgs
@@ -585,22 +598,22 @@ def __init__(self, log_Mhalo=12.2, c=7.67, f_hot=0.75, alpha=0.,
         self.f_hot = f_hot
         self.zero_inner_ne = 0. # kpc
         self.cosmo = cosmo
+        self.fb = fb
 
         # Init more
-        self.setup_param(cosmo=self.cosmo)
+        self.setup_param(cosmo=self.cosmo,)
 
-    def setup_param(self,cosmo):
+    def setup_param(self,cosmo,):
         """ Setup key parameters of the model
+
+        Args:
+            cosmo: astropy cosmology
+              Cosmology of the universe
         """
         # Cosmology
-        if cosmo is None:
-            self.rhoc = 9.2e-30 * units.g / units.cm**3
-            self.fb = 0.16       # Baryon fraction
-            self.H0 = 70. *units.km/units.s/ units.Mpc
-        else:
-            self.rhoc = self.cosmo.critical_density(self.z)
-            self.fb = cosmo.Ob0/cosmo.Om0
-            self.H0 = cosmo.H0
+        self.rhoc = self.cosmo.critical_density(self.z)
+        self.H0 = cosmo.H0
+
         # Dark Matter
         self.q = self.cosmo.Ode0/(self.cosmo.Ode0+self.cosmo.Om0*(1+self.z)**3) 
         #r200 = (((3*Mlow*constants.M_sun.cgs) / (4*np.pi*200*rhoc))**(1/3)).to('kpc')
@@ -1088,7 +1101,9 @@ class M31(ModifiedNFW):
     Taking mass from van der Marel 2012
 
     """
-    def __init__(self, log_Mhalo=12.18, c=7.67, f_hot=0.75, alpha=2, y0=2, **kwargs):
+    def __init__(self, 
+                 log_Mhalo=12.18, c=7.67, 
+                 alpha=2, y0=2,f_hot=0.75, **kwargs):
 
         # Init ModifiedNFW
         ModifiedNFW.__init__(self, log_Mhalo=log_Mhalo, c=c, f_hot=f_hot,

frb/tests/test_frbhalos.py

@@ -87,6 +87,11 @@ def test_modified_NFW():
     ne = mNFW.ne(xyz)
     assert np.all(ne > nH)
 
+def test_log_MGCM():
+    # Init
+    mNFW = halos.ModifiedNFW(log_MCGM=11.0)
+    assert np.isclose(mNFW.f_hot, 0.3989835640019133)
+
 def test_milky_way():
     Galaxy = halos.MilkyWay()
     assert np.isclose(Galaxy.M_halo.to('M_sun').value, 1.51356125e+12)