restructuring model_galaxy

bagpipes/__init__.py

@@ -3,9 +3,10 @@
 from .utils import *
 from .galaxy import Galaxy
 from .model_galaxy import Model_Galaxy
-from .star_formation_history import Star_Formation_History
+from .star_formation_history import star_formation_history
 from .fit import Fit
 from .compare_fits import Compare_Fits
 from .catalogue_fit import *
 from .make_cloudy_models import *
-from .igm_inoue2014 import *
+from .igm_inoue2014 import *
+#from . import plotting

bagpipes/chemical_evolution_history.py

@@ -0,0 +1,73 @@
+from __future__ import print_function, division, absolute_import
+
+import numpy as np
+
+from . import utils
+from .star_formation_history import component_types
+
+
+class chemical_evolution_history:
+
+
+    def __init__(self, model_comp):
+
+        self.model_comp = model_comp
+
+        self.zmet_vals = utils.zmet_vals[utils.model_type]
+        self.zmet_lims = utils.zmet_lims[utils.model_type]
+
+        self.zmet_weights = {}
+
+        for name in list(self.model_comp):
+            if name in component_types or name[:-1] in component_types:
+                if (("metallicity dist" in list(self.model_comp[name])) 
+                        and self.model_comp[name]["metallicity dist"]):
+                    self.zmet_weights[name] = self.exp(self.model_comp[name])
+
+                else:
+                    self.zmet_weights[name] = self.delta(self.model_comp[name])
+
+
+    """ Delta function metallicity history. """
+    def delta(self, comp):
+
+        zmet = comp["metallicity"]
+
+        weights = np.zeros(self.zmet_vals.shape[0])
+
+        high_ind = self.zmet_vals[self.zmet_vals < zmet].shape[0]
+
+        if high_ind == self.zmet_vals.shape[0]:
+            weights[-1] = 1.
+
+        elif high_ind == 0:
+            weights[0] = 1.
+
+        else:
+            low_ind = high_ind - 1
+            width = (self.zmet_vals[high_ind] - self.zmet_vals[low_ind])
+            weights[high_ind] = (zmet - self.zmet_vals[low_ind])/width
+            weights[high_ind-1] = 1 - weights[high_ind]
+
+        return weights
+
+
+    """ P(Z) = exp(-z/z_mean). Currently no age dependency! """
+    def exp(self, comp):
+
+        mean_zmet = comp["metallicity"]
+
+        weights = np.zeros(self.zmet_vals.shape[0])
+
+        vals_hr = np.arange(0., 10., 0.01) + 0.005
+
+        factors_hr = (1./mean_zmet)*np.exp(-self.vals_hr/mean_zmet)
+
+        for i in range(zmet_weights.shape[0]):
+            lowmask = (vals_hr > self.zmet_lims[i]) 
+            highmask = (vals_hr < self.zmet_lims[i+1])
+            weights[i] = np.sum(0.01*factors_hr[lowmask & highmask])
+
+        return weights
+
+