Changed photometric uncertainty inflation in FitModel to relative to …

…the uncertainties instead of relative to the flux, fixed issue when apply_weights was set to False, updated the Exo-Rem grid and added the high-resolution grid (exo-rem-highres)
tomasstolker · Jul 23, 2023 · d517148 · d517148
1 parent c85207d
commit d517148
Show file tree

Hide file tree

Showing 4 changed files with 83 additions and 20 deletions.
diff --git a/species/analysis/fit_model.py b/species/analysis/fit_model.py
@@ -990,6 +990,18 @@ def __init__(
                         # Set weight for photometry to FWHM of filter
                         read_filt = read_filter.ReadFilter(phot_item)
                         self.weights[phot_item] = read_filt.filter_fwhm()
+                        print(f"   - {phot_item} = {self.weights[phot_item]:.2e}")
+
+            else:
+                for spec_item in inc_spec:
+                    spec_size = self.spectrum[spec_item][0].shape[0]
+                    self.weights[spec_item] = np.full(spec_size, 1.)
+                    print(f"   - {spec_item} = {self.weights[spec_item][0]:.2f}")
+
+                for phot_item in inc_phot:
+                    # Set weight to 1 if apply_weights=False
+                    self.weights[phot_item] = 1.
+                    print(f"   - {phot_item} = {self.weights[phot_item]:.2f}")
 
         else:
             self.weights = apply_weights
@@ -1399,12 +1411,22 @@ def lnlike_func(
                 instr_check = phot_filter.split(".")[0]
 
                 if phot_filter in phot_scaling:
-                    # Inflate photometric error for filter
-                    phot_var += phot_scaling[phot_filter] ** 2 * obj_item[0] ** 2
+                    # Inflate photometry uncertainty for filter
+
+                    # Scale relative to the flux
+                    # phot_var += phot_scaling[phot_filter] ** 2 * obj_item[0] ** 2
+
+                    # Scale relative to the uncertainty
+                    phot_var += phot_scaling[phot_filter] ** 2 * obj_item[1] ** 2
 
                 elif instr_check in phot_scaling:
-                    # Inflate photometric error for instrument
-                    phot_var += phot_scaling[instr_check] ** 2 * obj_item[0] ** 2
+                    # Inflate photometry uncertainty for instrument
+
+                    # Scale relative to the flux
+                    # phot_var += phot_scaling[instr_check] ** 2 * obj_item[0] ** 2
+
+                    # Scale relative to the uncertainty
+                    phot_var += phot_scaling[instr_check] ** 2 * obj_item[1] ** 2
 
                 ln_like += -0.5 * weight * (obj_item[0] - phot_flux) ** 2 / phot_var
 
@@ -1415,14 +1437,26 @@ def lnlike_func(
                 for j in range(obj_item.shape[1]):
                     phot_var = obj_item[1, j] ** 2
 
-                    if (
-                        self.model == "powerlaw"
-                        and f"{phot_filter}_error" in param_dict
-                    ):
-                        phot_var += (
-                            param_dict[f"{phot_filter}_error"] ** 2
-                            * obj_item[0, j] ** 2
-                        )
+                    # Get the telescope/instrument name
+                    instr_check = phot_filter.split(".")[0]
+
+                    if phot_filter in phot_scaling:
+                        # Inflate photometry uncertainty for filter
+
+                        # Scale relative to the flux
+                        # phot_var += phot_scaling[phot_filter] ** 2 * obj_item[0, j] ** 2
+
+                        # Scale relative to the uncertainty
+                        phot_var += phot_scaling[phot_filter] ** 2 * obj_item[1, j] ** 2
+
+                    elif instr_check in phot_scaling:
+                        # Inflate photometry uncertainty for instrument
+
+                        # Scale relative to the flux
+                        # phot_var += phot_scaling[instr_check] ** 2 * obj_item[0, j] ** 2
+
+                        # Scale relative to the uncertainty
+                        phot_var += phot_scaling[instr_check] ** 2 * obj_item[1, j] ** 2
 
                     ln_like += (
                         -0.5 * weight * (obj_item[0, j] - phot_flux) ** 2 / phot_var

diff --git a/species/data/model_data.json b/species/data/model_data.json
@@ -134,10 +134,20 @@
     "exo-rem": {
         "parameters": ["teff", "logg", "feh", "c_o_ratio"],
         "name": "Exo-REM",
-        "file size": "706 MB",
-        "wavelength range": [0.9, 5],
-        "resolution": 10000,
-        "teff range": [1000, 2000],
+        "file size": "530 MB",
+        "wavelength range": [0.35, 250.0],
+        "resolution": 500,
+        "teff range": [400, 2000],
+        "reference": "Charney et al. (2018)",
+        "url": "https://ui.adsabs.harvard.edu/abs/2018ApJ...854..172C/abstract"
+    },
+    "exo-rem-highres": {
+        "parameters": ["teff", "logg", "feh", "c_o_ratio"],
+        "name": "Exo-REM",
+        "file size": "16 GB",
+        "wavelength range": [0.67, 250.0],
+        "resolution": 20000,
+        "teff range": [400, 2000],
         "reference": "Charney et al. (2018)",
         "url": "https://ui.adsabs.harvard.edu/abs/2018ApJ...854..172C/abstract"
     },

diff --git a/species/data/model_spectra.py b/species/data/model_spectra.py
@@ -91,6 +91,15 @@ def add_model_grid(
             "add_model of Database."
         )
 
+    elif model_name == "exo-rem":
+        warnings.warn(
+            "The Exo-Rem grid has been updated to the latest version "
+            "from https://lesia.obspm.fr/exorem/YGP_grids/. Please "
+            "consider removing the grid from the 'data_folder' if "
+            "needed such that the latest version of the grid will "
+            "be downloaded and added to the HDF5 database."
+        )
+
     if not os.path.exists(input_path):
         os.makedirs(input_path)
 

diff --git a/species/util/read_util.py b/species/util/read_util.py
@@ -301,12 +301,22 @@ def update_objectbox(
             instr_name = key.split(".")[0]
 
             if f"{key}_error" in model_param:
-                # Inflate photometric error of filter
-                var_add = model_param[f"{key}_error"] ** 2 * value[0] ** 2
+                # Inflate photometry uncertainty of filter
+
+                # Scale relative to the flux
+                # var_add = model_param[f"{key}_error"] ** 2 * value[0] ** 2
+
+                # Scale relative to the uncertainty
+                var_add = model_param[f"{key}_error"] ** 2 * value[1] ** 2
 
             elif f"{instr_name}_error" in model_param:
-                # Inflate photometric error of instrument
-                var_add = model_param[f"{instr_name}_error"] ** 2 * value[0] ** 2
+                # Inflate photometry uncertainty of instrument
+
+                # Scale relative to the flux
+                # var_add = model_param[f"{instr_name}_error"] ** 2 * value[0] ** 2
+
+                # Scale relative to the uncertainty
+                var_add = model_param[f"{instr_name}_error"] ** 2 * value[1] ** 2
 
             else:
                 # No inflation required