Exclude SpeX spectra that can not be flux calibrated

tomasstolker · Jun 27, 2024 · e1d8a6c · e1d8a6c
1 parent 16206ae
commit e1d8a6c
Show file tree

Hide file tree

Showing 3 changed files with 14 additions and 6 deletions.
diff --git a/species/data/database.py b/species/data/database.py
@@ -2131,7 +2131,7 @@ def add_spectra(
 
         with h5py.File(self.database, "a") as hdf5_file:
             if f"spectra/{spec_library}" in hdf5_file:
-                del hdf5_file["spectra/" + spec_library]
+                del hdf5_file[f"spectra/{spec_library}"]
 
             add_spec_library(self.data_folder, hdf5_file, spec_library, sptypes)
 

diff --git a/species/data/spec_data/spec_spex.py b/species/data/spec_data/spec_spex.py
@@ -116,7 +116,7 @@ def add_spex(input_path: str, database: h5py._hl.files.File) -> None:
                 name = name[0].decode("utf-8")
 
             input_file = f"spex_{name}.xml"
-            data_file = Path(input_path) / input_file
+            data_file = Path(data_folder) / input_file
 
             if not data_file.exists():
                 print()
@@ -144,12 +144,12 @@ def add_spex(input_path: str, database: h5py._hl.files.File) -> None:
         if file_item.stem.startswith("spex_") and file_item.suffix == ".xml":
             table = parse_single_table(file_item)
 
-            wavelength = table.array["wavelength"]  # (Angstrom)
+            wavelength = 1e-4 * table.array["wavelength"]  # (A) -> (um)
             flux = table.array["flux"]  # Normalized units
             spec_res = table.get_field_by_id("res").value
 
-            wavelength = np.array(wavelength * 1e-4)  # (um)
-            flux = np.array(flux)  # (a.u.)
+            wavelength = np.array(wavelength)
+            flux = np.array(flux)
             error = np.full(flux.size, np.nan)
 
             # 2MASS magnitudes
@@ -217,6 +217,9 @@ def add_spex(input_path: str, database: h5py._hl.files.File) -> None:
             except KeyError:
                 sptype_nir = None
 
+            if np.isnan(h_mag):
+                continue
+
             h_flux, _ = h_twomass.magnitude_to_flux(h_mag, error=None)
             phot = h_twomass.spectrum_to_flux(wavelength, flux)  # Normalized units
 

diff --git a/species/fit/compare_spectra.py b/species/fit/compare_spectra.py
@@ -274,7 +274,6 @@ def spectral_type(
                                 chi_sq = (
                                     spec_item[indices, 1][idx_select] - c_k * flux_resample[idx_select]
                                 ) / spec_item[indices, 2][idx_select]
-                                print(chi_sq)
 
                                 g_k += np.sum(w_i * chi_sq**2)
 
@@ -285,6 +284,12 @@ def spectral_type(
                                 #     np.dot(obj_inv_cov_crop,
                                 #            spec_item[indices, 1]-c_k*flux_resample))
 
+                            if np.isnan(c_k_spec):
+                                # This can happen if the spectrum only contains NaNs
+                                # because there is an issue with the flux calibration
+                                # For example: ULAS J141623.94+134836.3 (SpeX)
+                                g_k = np.inf
+
                             # Append to the lists of results
 
                             name_list.append(item)