SED dependent reddening

Galametz et al. (2017) investigate the impact of SED dependent reddening. The Milky Way extinction is a function of the source SED. Because we don’t know this ahead of time traditionally this was done by ‘dereddening’ the measurements assuming a flat SED. Here we investigate the impact of SED dependent reddening on the low dust COSMOS field. This impact will be larger on areas with larger reddening in wide surveys.

Equation 2 in Galametz et al. (2017) provides a means to compute reddening dependent on the SED. The essential idea is that SEDs with a blue slope across the filter will be more reddened than SEDs with a red slope or flat spectra. In the paper they recommend computing the reddening for a fixed ebv value and assuming a linear relation between the reddening and the \(E(B-V)\) value on that position of the sky.

In this notebook we use LePHARE classes to compute the SED dependent reddening on a library of SEDs and to investigate how it differs from the generic filter based approach.

[1]:
import lephare as lp
import numpy as np
from astropy.table import Table
from scipy.interpolate import interp1d
from matplotlib import pylab as plt
import struct
import timeit
import os

%matplotlib inline
User defined LEPHAREDIR is set. Code runs depend on all required
auxiliary data being present at /home/docs/.cache/lephare/data.
User defined LEPHAREWORK is set. All intermediate files will be written to:
 /home/docs/.cache/lephare/work
OMP_NUM_THREADS set to 1 (max available CPUs: 2)

1 Run prepare to give us access to library of SEDs

As is typical we set the config values and run prepare to create the required model and magnitude libraries. We update some key config values in order to switch on the SED dependent reddening. These are

  • EXT_MW_CURVE (CARDELLI[def] or NONE) Extinction curve for the Milky Way extinction (used for the energy balance when adding dust emission). Should be in $LEPHAREDIR/ext if relative.

  • EXT_ATMOS_CURVE (NONE[def] or e.g. SB_calzetti.dat) Extinction curve for the atmospheric extinction (used for the energy balance when adding dust emission). Should be in $LEPHAREDIR/ext if relative.

  • APPLY_MW_EXTINCTION (NONE[DEF], CLASSIC, GALAMETZ) Method to apply the Milky Way extinction to the templates. If CLASSIC, the extinction is applied using the E(B-V) value and the extinction curve. If GALAMETZ, the extinction is applied using the E(B-V) value and the extinction curve, but also taking into account the SED dependence of the extinction (see Galametz et al. 2017). In this case, the MW_REFERENCE_MODEL keyword must be set to define the reference SED for which the E(B-V) value is defined.

  • MW_REFERENCE_MODEL (sed/STAR/PICKLES/b5i.sed[DEF]) Reference SED for which the E(B-V) value is defined when using the GALAMETZ method to apply the Milky Way extinction. Should be in $LEPHAREDIR if relative.

  • MW_GLOBAL_EBV (NONE[def] or float) Global E(B-V) value to apply to all templates when using the GALAMETZ method to apply the Milky Way extinction. If 0, the E(B-V) value is read from the file defined by MW_EBV_FILE.

  • MW_EBV_FILE (NONE[def] or string) Name of the file containing the E(B-V) values to apply to each template when using the GALAMETZ method to apply the Milky Way extinction. The file should have two columns: source ID and E(B-V) value. Should be absolute path.

[2]:
config = lp.default_cosmos_config.copy()
bands = "ugrizy"

config.update(
    {
        "CAT_IN": os.path.join(lp.LEPHAREDIR, "examples/COSMOS_no_mw_correction.in"),
        # We look at a reduced filter list to demonstrate with ugrizy data
        "FILTER_LIST": f"cosmos/u_new.pb,hsc/gHSC.pb,hsc/rHSC.pb,hsc/iHSC.pb,hsc/zHSC.pb,hsc/yHSC.pb",
        # We use single values to deal with reduction in FILTER_LIST
        "ERR_SCALE": "0.02",  # Value from page 12 first paragraph in desprez 2023 diff fro phosphours
        "ERR_FACTOR": "1.5",  # Again from paper - diff for phosphorous - 2 lowered from 1.5 to 1. following chisquared dist
        "FILTER_CALIB": "0",
        "GLB_CONTEXT": np.sum(2 ** np.arange(len(bands))),
        "MABS_CONTEXT": np.sum(2 ** np.arange(len(bands))),
        # Reduced z grid for speed in demonstration
        "Z_STEP": ".1,0.,6.",
        # Reduced EM line dispersion for SPEED to demonstrate
        "EM_DISPERSION": "1.",
        ## THESE ARE THE VALUES THAT RELATE TO THE NEW REDDENING MODEL
        "EXT_ATMOS_CURVE": "NONE",
        "APPLY_MW_EXTINCTION": "GALAMETZ",  # "NO[DEF], CLASSIC, GALAMETZ"
        "EXT_MW_CURVE": "LMC_Fitzpatrick.dat",
        "MW_REFERENCE_MODEL": "sed/STAR/PICKLES/b5i.sed",  # The default - a B5 star
        # "MW_GLOBAL_EBV" : "0.016", # This can also be a file with id and ebv value corresponding to each source:
        "MW_EBV_FILE": os.path.join(
            lp.LEPHAREDIR, "examples/EBV_MW.in"
        ),  # FILE with IDs corresponding to the input and EBV values.
    }
)
[3]:
# Download the required data
lp.data_retrieval.get_auxiliary_data(
    keymap=config,
    # We need to have the additional files required for the Galametz method
    additional_files=[
        "examples/COSMOS_no_mw_correction.in",
        "examples/EBV_MW.in",
        "ext/LMC_Fitzpatrick.dat",
        "sed/STAR/PICKLES/b5i.sed",
    ],
)
Registry file downloaded and saved as data_registry.txt.
/home/docs/checkouts/readthedocs.org/user_builds/lephare/envs/latest/lib/python3.10/site-packages/lephare/data_retrieval.py:431: UserWarning: /home/docs/.cache/lephare/data/sed/STAR/STAR_MOD_ALL.list is present locally and will not be overwritten.
  warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/lephare/envs/latest/lib/python3.10/site-packages/lephare/data_retrieval.py:431: UserWarning: /home/docs/.cache/lephare/data/sed/GAL/COSMOS_SED/COSMOS_MOD.list is present locally and will not be overwritten.
  warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/lephare/envs/latest/lib/python3.10/site-packages/lephare/data_retrieval.py:431: UserWarning: /home/docs/.cache/lephare/data/sed/QSO/SALVATO09/AGN_MOD.list is present locally and will not be overwritten.
  warnings.warn(
Checking/downloading 418 files...
418 completed.
All files downloaded successfully and are non-empty.
Downloading file 'examples/COSMOS_no_mw_correction.in' from 'https://raw.githubusercontent.com/lephare-photoz/lephare-data/main/examples/COSMOS_no_mw_correction.in' to '/home/docs/.cache/lephare/data'.
Downloading file 'examples/EBV_MW.in' from 'https://raw.githubusercontent.com/lephare-photoz/lephare-data/main/examples/EBV_MW.in' to '/home/docs/.cache/lephare/data'.
Downloading file 'ext/LMC_Fitzpatrick.dat' from 'https://raw.githubusercontent.com/lephare-photoz/lephare-data/main/ext/LMC_Fitzpatrick.dat' to '/home/docs/.cache/lephare/data'.
Downloading file 'sed/STAR/PICKLES/b5i.sed' from 'https://raw.githubusercontent.com/lephare-photoz/lephare-data/main/sed/STAR/PICKLES/b5i.sed' to '/home/docs/.cache/lephare/data'.
Checking/downloading 4 files...
4 completed.
All files downloaded successfully and are non-empty.
[4]:
# Check the config values
config
[4]:
{'ADAPT_BAND': '5',
 'ADAPT_CONTEXT': '-1',
 'ADAPT_LIM': '1.5,23.0',
 'ADAPT_MODBIN': '1,1000',
 'ADAPT_ZBIN': '0.01,6',
 'ADDITIONAL_MAG': 'none',
 'ADD_DUSTEM': 'NO',
 'ADD_EMLINES': '0,10000',
 'AGE_RANGE': '0.,15.e9',
 'AUTO_ADAPT': 'NO',
 'CAT_FMT': 'MEME',
 'CAT_IN': '/home/docs/.cache/lephare/data/examples/COSMOS_no_mw_correction.in',
 'CAT_LINES': '0,1000000000',
 'CAT_MAG': 'AB',
 'CAT_OUT': 'zphot.out',
 'CAT_TYPE': 'LONG',
 'CHI2_OUT': 'NO',
 'COSMOLOGY': '70,0.3,0.7',
 'DZ_WIN': '1.0',
 'EBV_RANGE': '0,9',
 'EB_V': '0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.5',
 'EM_DISPERSION': '1.',
 'EM_LINES': 'EMP_UV',
 'ERR_FACTOR': '1.5',
 'ERR_SCALE': '0.02',
 'EXTERNALZ_FILE': 'NONE',
 'EXTINC_LAW': 'SMC_prevot.dat,SB_calzetti.dat,SB_calzetti_bump1.dat,SB_calzetti_bump2.dat',
 'FILTER_CALIB': '0',
 'FILTER_FILE': 'filter_cosmos',
 'FILTER_LIST': 'cosmos/u_new.pb,hsc/gHSC.pb,hsc/rHSC.pb,hsc/iHSC.pb,hsc/zHSC.pb,hsc/yHSC.pb',
 'FILTER_REP': '/home/docs/.cache/lephare/data/filt',
 'FIR_CONT': '-1',
 'FIR_FREESCALE': 'YES',
 'FIR_LIB': 'NONE',
 'FIR_LMIN': '7.0',
 'FIR_SCALE': '-1',
 'FIR_SUBSTELLAR': 'NO',
 'FORB_CONTEXT': '-1',
 'GAL_FSCALE': '1.',
 'GAL_LIB': 'LIB_CE',
 'GAL_LIB_IN': 'LIB_CE',
 'GAL_LIB_OUT': 'CE_COSMOS',
 'GAL_SED': '/home/docs/.cache/lephare/data/sed/GAL/COSMOS_SED/COSMOS_MOD.list',
 'GLB_CONTEXT': np.int64(63),
 'INP_TYPE': 'F',
 'LIB_ASCII': 'NO',
 'LIMITS_MAPP_CUT': '90',
 'LIMITS_MAPP_REF': '1',
 'LIMITS_MAPP_SEL': '1',
 'LIMITS_ZBIN': '0,99',
 'MABS_CONTEXT': np.int64(63),
 'MABS_FILT': '1,2,3,4',
 'MABS_METHOD': '1',
 'MABS_REF': '11',
 'MABS_ZBIN': '0,0.5,2,4,6',
 'MAGTYPE': 'AB',
 'MAG_ABS': '-24,-5',
 'MAG_ABS_QSO': '-30,-10',
 'MAG_REF': '3',
 'MIN_THRES': '0.02',
 'MOD_EXTINC': '0,0',
 'PARA_OUT': '/home/docs/.cache/lephare/data/examples/output.para',
 'PDZ_OUT': 'test',
 'PDZ_TYPE': 'BAY_ZG',
 'QSO_FSCALE': '1.',
 'QSO_LIB': 'LIB_QSO',
 'QSO_LIB_IN': 'LIB_QSO',
 'QSO_LIB_OUT': 'QSO_COSMOS',
 'QSO_SED': '/home/docs/.cache/lephare/data/sed/QSO/SALVATO09/AGN_MOD.list',
 'RF_COLORS': '32,4,4,13',
 'RM_DISCREPANT_BD': '500',
 'SPEC_OUT': 'NO',
 'STAR_FSCALE': '3.432E-09',
 'STAR_LIB': 'LIB_STAR',
 'STAR_LIB_IN': 'LIB_STAR',
 'STAR_LIB_OUT': 'STAR_COSMOS',
 'STAR_SED': '/home/docs/.cache/lephare/data/sed/STAR/STAR_MOD_ALL.list',
 'TRANS_TYPE': '1',
 'VERBOSE': 'NO',
 'ZFIX': 'NO',
 'ZGRID_TYPE': '0',
 'ZPHOTLIB': 'CE_COSMOS,STAR_COSMOS,QSO_COSMOS',
 'Z_INTERP': 'YES',
 'Z_METHOD': 'BEST',
 'Z_RANGE': '0.,99.99',
 'Z_STEP': '.1,0.,6.',
 'EXT_ATMOS_CURVE': 'NONE',
 'APPLY_MW_EXTINCTION': 'GALAMETZ',
 'EXT_MW_CURVE': 'LMC_Fitzpatrick.dat',
 'MW_REFERENCE_MODEL': 'sed/STAR/PICKLES/b5i.sed',
 'MW_EBV_FILE': '/home/docs/.cache/lephare/data/examples/EBV_MW.in'}
[5]:
# Because we are dealing with some low level functionality we also need the keymap in the native LePHARE format
keymap = lp.all_types_to_keymap(config)

Make the Model libraries

As for a typical run we run filter, sedtolib, and mag_mal to build the binary files. During those runs we now compute the extinction per band in addition to the band_pass_correction during mag_gal. When you initialize the PhotoZ object the band pass corrections are normalized to the reference model which is typically a B5 star.

[6]:
lp.prepare(keymap)
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor    TG-cor      VEGA M_sun(AB)   CALIB      Lb_eff    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195   -0.2745  -20.8527    6.3135       0      0.3668      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860   -0.2458  -20.7272    5.0764       0      0.4780      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466    0.2580  -21.5143    4.6523       0      0.6178      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942    0.6138  -22.2286    4.5323       0      0.7666      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169    0.7625  -22.6733    4.5147       0      0.8903      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534    0.7763  -22.9145    4.5081       0      0.9782      1.0000
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : STAR
# STAR_SED    :/home/docs/.cache/lephare/data/sed/STAR/STAR_MOD_ALL.list
# STAR_LIB    :LIB_STAR
# STAR_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_STAR.doc
# STAR_FSCALE :0.0000
#######################################
Number of templates in the list 254
Number of SED in the list (including different ages) 254
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# COSMOLOGY   :70.0000,0.3000,0.7000
# STAR_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_STAR(.doc & .bin)
# STAR_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS(.doc & .bin)
# LIB_ASCII   NO
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# CREATION_DATE Wed Oct 28 09:40:54 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : GAL
# GAL_SED    :/home/docs/.cache/lephare/data/sed/GAL/COSMOS_SED/COSMOS_MOD.list
# GAL_LIB    :LIB_CE
# GAL_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.doc
# GAL_LIB phys:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.phys
# SEL_AGE    :none
# GAL_FSCALE :1.0000
# AGE_RANGE   0.0000 15000000000.0000
#######################################
Number of templates in the list 31
Number of SED in the list (including different ages) 31
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# GAL_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_CE(.doc & .bin)
# GAL_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/CE_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# EM_LINES   EMP_UV
# EM_DISPERSION   1.0000,
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:40:54 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : QSO
# QSO_SED    :/home/docs/.cache/lephare/data/sed/QSO/SALVATO09/AGN_MOD.list
# QSO_LIB    :LIB_QSO
# QSO_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_QSO.doc
# QSO_FSCALE :1.0000
#######################################
Number of templates in the list 30
Number of SED in the list (including different ages) 30
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# QSO_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_QSO(.doc & .bin)
# QSO_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:40:54 2026
#############################################
[7]:
# This loads the library of SEDs giving us access to the reddening per model
photz = lp.PhotoZ(keymap)
#######################################
# PHOTOMETRIC REDSHIFT with OPTIONS   #
# Config file            :
# CAT_IN                 : /home/docs/.cache/lephare/data/examples/COSMOS_no_mw_correction.in
# CAT_OUT                : zphot.out
# CAT_LINES              : 0 1000000000
# PARA_OUT               : /home/docs/.cache/lephare/data/examples/output.para
# INP_TYPE               : F
# CAT_FMT[0:MEME 1:MMEE] : 0
# CAT_MAG                : AB
# ZPHOTLIB               : CE_COSMOS STAR_COSMOS QSO_COSMOS
# APPLY_MW_EXTINCTION    : GALAMETZ
# MW_GLOBAL_EBV          : -1.000000
# FIR_LIB                :
# FIR_LMIN               : 7.000000
# FIR_CONT               : -1.000000
# FIR_SCALE              : -1.000000
# FIR_FREESCALE          : YES
# FIR_SUBSTELLAR         : NO
# ERR_SCALE              : 0.020000
# ERR_FACTOR             : 1.500000
# GLB_CONTEXT            : 63
# FORB_CONTEXT           : -1
# DZ_WIN                 : 1.000000
# MIN_THRES              : 0.020000
# MAG_ABS                : -24.000000 -5.000000
# MAG_ABS_AGN            : -30.000000 -10.000000
# MAG_REF                : 3
# NZ_PRIOR               : -1 -1
# Z_INTERP               : YES
# Z_METHOD               : BEST
# MABS_METHOD            : 1
# MABS_CONTEXT           : 63
# MABS_REF               : 11
# AUTO_ADAPT             : NO
# ADAPT_BAND             : 5
# ADAPT_LIM              : 1.500000 23.000000
# ADAPT_ZBIN             : 0.010000 6.000000
# ZFIX                   : NO
# SPEC_OUT               : NO
# CHI_OUT                : NO
# PDZ_OUT                : test
#######################################
Reading input librairies ...
Read lib
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.bin
 Done with the library reading with 1830 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.bin
 Done with the library reading with 2084 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.bin
 Done with the library reading with 3975 SED read.
Read lib out
Done reading the FIR lib
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor      VEGA   CALIB    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195  -20.8527       0      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860  -20.7272       0      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466  -21.5143       0      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942  -22.2286       0      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169  -22.6733       0      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534  -22.9145       0      1.0000

At this stage the EBV values per source have already been set as the EBV file was set at the input time. We can therefore inspect each value per source.

[8]:
sources = photz.read_photoz_sources()
photz.read_mw_ebv(sources)
sources[0].mw_ebv
mw_ebv_vals = np.array([s.mw_ebv for s in sources])
plt.hist(mw_ebv_vals, bins=30)
plt.xlabel("$E(B-V)$ [mag]")
WARNING: No scaling --> No z 255
WARNING: No scaling --> No z 261
WARNING: No scaling --> No z 358
WARNING: No scaling --> No z 387
WARNING: No scaling --> No z 408
WARNING: No scaling --> No z 3369
WARNING: No scaling --> No z 5398
WARNING: No scaling --> No z 5400
WARNING: No scaling --> No z 5420
WARNING: No scaling --> No z 5432
WARNING: No scaling --> No z 5481
WARNING: No scaling --> No z 5524
WARNING: No scaling --> No z 6790
WARNING: No scaling --> No z 13469
WARNING: No scaling --> No z 13665
WARNING: No scaling --> No z 13681
WARNING: No scaling --> No z 17624
WARNING: No scaling --> No z 17688
WARNING: No scaling --> No z 17724
WARNING: No scaling --> No z 18408
WARNING: No scaling --> No z 18412
[8]:
Text(0.5, 0, '$E(B-V)$ [mag]')
../_images/notebooks_Galametz_reddening_11_2.png

2. Compute reddening

We want to compute the baseline reddening as well as the per SED values for comparison. We do this using the original method. This function simply returns a global value for each filter.

Start with traditional values

Here we just get one value per band as a baseline “CLASSIC” run. These could be applied to the input catalogue as was done or applied to the models by setting APPLY_MW_EXTINCTION to CLASSIC.

[9]:
all_filters, _, _, baseline_albd = lp.classic_extinction_values(config)
baseline_albd
[9]:
array([4.84019494, 3.66662273, 2.61482635, 1.90912797, 1.46583258,
       1.23635852])

Look at an example SED

We want to see how they vary over a filter and how the product differs from the pure filter and therefore impacts the reddening.

[10]:
# Pick one model to focus on
n_model = 20
# ex_sed=observe(photz, 5, mag)
keymap["t"] = lp.keyword("t", "Q")
qso_mag = lp.GalMag(keymap)
ex_sed = photz.fullLib[n_model]
print(f"Model {n_model} at redshift z={photz.zLib[n_model]} has its z=0 model at index {ex_sed.index_z0}\n")
ex_sed.lamb_flux = photz.fullLib[ex_sed.index_z0].lamb_flux
ex_sed.generate_spectra(photz.zLib[n_model], 1)
x, y = ex_sed.data()
plt.loglog(x, y)
Model 20 at redshift z=2.0 has its z=0 model at index 0

WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
[10]:
[<matplotlib.lines.Line2D at 0x7ab3db8c0160>]
../_images/notebooks_Galametz_reddening_15_2.png
[11]:
# The first example is a QSO
ex_sed.is_qso()
[11]:
True

We can then inspect the model dependent extinction values for any given model and compare them to the general values

[12]:
# After running the "prepare" function above, every SED in the library has the milky way extinction values in each band
ex_sed.milky_way_extinction
[12]:
[4.851102499362176,
 3.7172989513889116,
 2.646508627412381,
 1.9173622941726334,
 1.4643069933409132,
 1.2369800508024276]
[13]:
bands
[13]:
'ugrizy'
[14]:
fig, ax = plt.subplots()
filter_means = np.array([f.lambdaEff() for f in photz.allFilters])
for n, b in enumerate(bands):
    y = np.array([s.milky_way_extinction[n] - baseline_albd[n] for s in photz.fullLib])
    y = y[np.isfinite(y)]
    if len(y) < 5:
        continue

    ax.violinplot([y], positions=[filter_means[n]], widths=150, bw_method=0.1)

# FIX LIMITS EXPLICITLY

xmin = filter_means.min()
xmax = filter_means.max()
pad = 0.05 * (xmax - xmin)
ax.set_xlim(xmin - pad, xmax + pad)

ax.set_xlabel("Filter mean wavelength [$\\AA$]")
ax.set_ylabel("Difference in reddening")
plt.title("Difference between standard value and SED dependent value for all models")
# TOP AXIS
fig.canvas.draw()
ax_top = ax.twiny()
ax_top.set_xlim(ax.get_xlim())
ax_top.set_xticks(filter_means)
ax_top.set_xticklabels(bands)
[14]:
[Text(3749.546551929171, 1, 'u'),
 Text(4760.2211316989415, 1, 'g'),
 Text(6141.9732413654265, 1, 'r'),
 Text(7636.664755689224, 1, 'i'),
 Text(8906.740869095343, 1, 'z'),
 Text(9770.57514285193, 1, 'y')]
../_images/notebooks_Galametz_reddening_20_1.png

2.1 Use the python function to get the extinction values

This runs prepare and loads the values per band. These are already computed by lephare.prepare and can also be accessed by interacting directly with the lephare.PhotoZ library of SEDs. We also have a helper function Python side to perform the full calculation and return the results:

[15]:
lp.compute_model_reddening?
[16]:
import time

start = time.time()
# Get the array of albd values for each model for each band
albd_lib = lp.compute_model_reddening(keymap)
end = time.time()
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor    TG-cor      VEGA M_sun(AB)   CALIB      Lb_eff    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195   -0.2745  -20.8527    6.3135       0      0.3668      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860   -0.2458  -20.7272    5.0764       0      0.4780      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466    0.2580  -21.5143    4.6523       0      0.6178      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942    0.6138  -22.2286    4.5323       0      0.7666      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169    0.7625  -22.6733    4.5147       0      0.8903      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534    0.7763  -22.9145    4.5081       0      0.9782      1.0000
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : STAR
# STAR_SED    :/home/docs/.cache/lephare/data/sed/STAR/STAR_MOD_ALL.list
# STAR_LIB    :LIB_STAR
# STAR_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_STAR.doc
# STAR_FSCALE :0.0000
#######################################
Number of templates in the list 254
Number of SED in the list (including different ages) 254
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# COSMOLOGY   :70.0000,0.3000,0.7000
# STAR_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_STAR(.doc & .bin)
# STAR_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS(.doc & .bin)
# LIB_ASCII   NO
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# CREATION_DATE Wed Oct 28 09:40:57 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : GAL
# GAL_SED    :/home/docs/.cache/lephare/data/sed/GAL/COSMOS_SED/COSMOS_MOD.list
# GAL_LIB    :LIB_CE
# GAL_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.doc
# GAL_LIB phys:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.phys
# SEL_AGE    :none
# GAL_FSCALE :1.0000
# AGE_RANGE   0.0000 15000000000.0000
#######################################
Number of templates in the list 31
Number of SED in the list (including different ages) 31
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# GAL_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_CE(.doc & .bin)
# GAL_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/CE_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# EM_LINES   EMP_UV
# EM_DISPERSION   1.0000,
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:40:58 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : QSO
# QSO_SED    :/home/docs/.cache/lephare/data/sed/QSO/SALVATO09/AGN_MOD.list
# QSO_LIB    :LIB_QSO
# QSO_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_QSO.doc
# QSO_FSCALE :1.0000
#######################################
Number of templates in the list 30
Number of SED in the list (including different ages) 30
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# QSO_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_QSO(.doc & .bin)
# QSO_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:40:58 2026
#############################################
#######################################
# PHOTOMETRIC REDSHIFT with OPTIONS   #
# Config file            :
# CAT_IN                 : /home/docs/.cache/lephare/data/examples/COSMOS_no_mw_correction.in
# CAT_OUT                : zphot.out
# CAT_LINES              : 0 1000000000
# PARA_OUT               : /home/docs/.cache/lephare/data/examples/output.para
# INP_TYPE               : F
# CAT_FMT[0:MEME 1:MMEE] : 0
# CAT_MAG                : AB
# ZPHOTLIB               : CE_COSMOS STAR_COSMOS QSO_COSMOS
# APPLY_MW_EXTINCTION    : GALAMETZ
# MW_GLOBAL_EBV          : -1.000000
# FIR_LIB                :
# FIR_LMIN               : 7.000000
# FIR_CONT               : -1.000000
# FIR_SCALE              : -1.000000
# FIR_FREESCALE          : YES
# FIR_SUBSTELLAR         : NO
# ERR_SCALE              : 0.020000
# ERR_FACTOR             : 1.500000
# GLB_CONTEXT            : 63
# FORB_CONTEXT           : -1
# DZ_WIN                 : 1.000000
# MIN_THRES              : 0.020000
# MAG_ABS                : -24.000000 -5.000000
# MAG_ABS_AGN            : -30.000000 -10.000000
# MAG_REF                : 3
# NZ_PRIOR               : -1 -1
# Z_INTERP               : YES
# Z_METHOD               : BEST
# MABS_METHOD            : 1
# MABS_CONTEXT           : 63
# MABS_REF               : 11
# AUTO_ADAPT             : NO
# ADAPT_BAND             : 5
# ADAPT_LIM              : 1.500000 23.000000
# ADAPT_ZBIN             : 0.010000 6.000000
# ZFIX                   : NO
# SPEC_OUT               : NO
# CHI_OUT                : NO
# PDZ_OUT                : test
#######################################
Reading input librairies ...
Read lib
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.bin
 Done with the library reading with 1830 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.bin
 Done with the library reading with 2084 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.bin
 Done with the library reading with 3975 SED read.
Read lib out
Done reading the FIR lib
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor      VEGA   CALIB    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195  -20.8527       0      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860  -20.7272       0      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466  -21.5143       0      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942  -22.2286       0      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169  -22.6733       0      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534  -22.9145       0      1.0000
[17]:
print(
    f"Reddening calculation took {end-start:.2f} seconds or {(end-start)/albd_lib.shape[0]*1000:.3f} milliseconds per model"
)
Reddening calculation took 1.47 seconds or 0.369 milliseconds per model
[18]:
albd_lib.shape  # models by filters
[18]:
(3975, 6)

Look at distribution of extinction values

Note that the bluer bands have a peak at the low end where the SED has dropped out of the band due to redshifting and the code uses the reddest possible SED in such a case.

[19]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
range = [0, 5]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)
for n, b in enumerate("ugrizy"):
    data = albd_lib.T[n][~np.isnan(albd_lib.T[n]) & (albd_lib.T[n] > 0)]
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data) / counts.max(),
        label=f"${b}$",
        color=colors[n],
        linewidth=2,
        alpha=0.5,
    )
    plt.axvline(baseline_albd[n], color=colors[n])
plt.legend(loc="upper left")
plt.ylabel("Relative density", fontsize=14)
plt.xlabel("Reddening $A_{\lambda}/E_{B-V}$ [mag]", fontsize=14)
plt.xticks(fontsize=12)
plt.yticks(fontsize=12)
[19]:
(array([0. , 0.2, 0.4, 0.6, 0.8, 1. , 1.2]),
 [Text(0, 0.0, '0.0'),
  Text(0, 0.2, '0.2'),
  Text(0, 0.4, '0.4'),
  Text(0, 0.6000000000000001, '0.6'),
  Text(0, 0.8, '0.8'),
  Text(0, 1.0, '1.0'),
  Text(0, 1.2000000000000002, '1.2')])
../_images/notebooks_Galametz_reddening_27_1.png
[20]:
# Get masks to access the object types
qsos = np.array([s.is_qso() for s in photz.fullLib])
gals = np.array([s.is_gal() for s in photz.fullLib])
stars = np.array([s.is_star() for s in photz.fullLib])
ebvs = np.array([s.ebv for s in photz.fullLib])
[21]:
np.sum(np.isclose(ebvs, 0.0)), len(ebvs)
[21]:
(np.int64(3975), 3975)
[22]:
model_numbers = np.array([s.nummod for s in photz.fullLib])
[23]:
# This is the b5 star that determines the BPC normalisation
np.sum(stars & (model_numbers == 15))
[23]:
np.int64(1)
[24]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
range = [0, 5]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)
for n, b in enumerate("ugrizy"):
    data = albd_lib[qsos].T[n][~np.isnan(albd_lib[qsos].T[n]) & (albd_lib[qsos].T[n] > 0)]
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data) / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n], color=colors[n])
plt.legend(loc="upper right")
plt.ylabel("relative density")
plt.xlabel("Reddening $A_{\lambda}/E_{B-V}$ [mag]")
plt.title("QSOs")
[24]:
Text(0.5, 1.0, 'QSOs')
../_images/notebooks_Galametz_reddening_32_1.png
[25]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
range = [0, 5]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)
for n, b in enumerate("ugrizy"):
    data = albd_lib[gals].T[n][(~np.isnan(albd_lib[gals].T[n])) & (albd_lib[gals].T[n] > 0)]
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data) / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n], color=colors[n], linestyle="--")
plt.legend(loc="upper right")
plt.ylabel("relative density")
plt.xlabel("Reddening $A_{\lambda}/E_{B-V}$ [mag]")
plt.title("Galaxies")
[25]:
Text(0.5, 1.0, 'Galaxies')
../_images/notebooks_Galametz_reddening_33_1.png
[26]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
range = [0, 5]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)
for n, b in enumerate("ugrizy"):
    data = albd_lib[stars].T[n][(~np.isnan(albd_lib[stars].T[n])) & (albd_lib[stars].T[n] > 0)]
    # counts, bins = np.histogram(data, bins=bins)
    plt.hist(data, bins=bins, alpha=0.5, label=b, color=colors[n], linewidth=2)
    plt.axvline(baseline_albd[n], color=colors[n])
plt.legend(loc="upper right")
plt.ylabel("relative density")
plt.xlabel("Reddening $A_{\lambda}/E_{B-V}$ [mag]")
plt.title("Stars")
[26]:
Text(0.5, 1.0, 'Stars')
../_images/notebooks_Galametz_reddening_34_1.png

2.2 Band pass correction

We also have a new function for getting the band pass correction based on hard coded B and V bands which are run on every object cpp side. This determines the scaling of the reddening since reddening models are normalised by the \(E_{SED}(B-V)\) value of the model with respect to a reference model. There is one value per model which scales the \(E(B-V)\) on at the position on the sky.

[27]:
lp.compute_band_pass_correction?
[28]:
t1 = time.time()
band_pass_correction = lp.compute_band_pass_correction(config)  #
t2 = time.time()
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor    TG-cor      VEGA M_sun(AB)   CALIB      Lb_eff    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195   -0.2745  -20.8527    6.3135       0      0.3668      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860   -0.2458  -20.7272    5.0764       0      0.4780      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466    0.2580  -21.5143    4.6523       0      0.6178      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942    0.6138  -22.2286    4.5323       0      0.7666      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169    0.7625  -22.6733    4.5147       0      0.8903      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534    0.7763  -22.9145    4.5081       0      0.9782      1.0000
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : STAR
# STAR_SED    :/home/docs/.cache/lephare/data/sed/STAR/STAR_MOD_ALL.list
# STAR_LIB    :LIB_STAR
# STAR_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_STAR.doc
# STAR_FSCALE :0.0000
#######################################
Number of templates in the list 254
Number of SED in the list (including different ages) 254
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# COSMOLOGY   :70.0000,0.3000,0.7000
# STAR_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_STAR(.doc & .bin)
# STAR_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS(.doc & .bin)
# LIB_ASCII   NO
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# CREATION_DATE Wed Oct 28 09:41:02 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : GAL
# GAL_SED    :/home/docs/.cache/lephare/data/sed/GAL/COSMOS_SED/COSMOS_MOD.list
# GAL_LIB    :LIB_CE
# GAL_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.doc
# GAL_LIB phys:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.phys
# SEL_AGE    :none
# GAL_FSCALE :1.0000
# AGE_RANGE   0.0000 15000000000.0000
#######################################
Number of templates in the list 31
Number of SED in the list (including different ages) 31
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# GAL_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_CE(.doc & .bin)
# GAL_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/CE_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# EM_LINES   EMP_UV
# EM_DISPERSION   1.0000,
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:41:03 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : QSO
# QSO_SED    :/home/docs/.cache/lephare/data/sed/QSO/SALVATO09/AGN_MOD.list
# QSO_LIB    :LIB_QSO
# QSO_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_QSO.doc
# QSO_FSCALE :1.0000
#######################################
Number of templates in the list 30
Number of SED in the list (including different ages) 30
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# QSO_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_QSO(.doc & .bin)
# QSO_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:41:03 2026
#############################################
#######################################
# PHOTOMETRIC REDSHIFT with OPTIONS   #
# Config file            :
# CAT_IN                 : /home/docs/.cache/lephare/data/examples/COSMOS_no_mw_correction.in
# CAT_OUT                : zphot.out
# CAT_LINES              : 0 1000000000
# PARA_OUT               : /home/docs/.cache/lephare/data/examples/output.para
# INP_TYPE               : F
# CAT_FMT[0:MEME 1:MMEE] : 0
# CAT_MAG                : AB
# ZPHOTLIB               : CE_COSMOS STAR_COSMOS QSO_COSMOS
# APPLY_MW_EXTINCTION    : GALAMETZ
# MW_GLOBAL_EBV          : -1.000000
# FIR_LIB                :
# FIR_LMIN               : 7.000000
# FIR_CONT               : -1.000000
# FIR_SCALE              : -1.000000
# FIR_FREESCALE          : YES
# FIR_SUBSTELLAR         : NO
# ERR_SCALE              : 0.020000
# ERR_FACTOR             : 1.500000
# GLB_CONTEXT            : 63
# FORB_CONTEXT           : -1
# DZ_WIN                 : 1.000000
# MIN_THRES              : 0.020000
# MAG_ABS                : -24.000000 -5.000000
# MAG_ABS_AGN            : -30.000000 -10.000000
# MAG_REF                : 3
# NZ_PRIOR               : -1 -1
# Z_INTERP               : YES
# Z_METHOD               : BEST
# MABS_METHOD            : 1
# MABS_CONTEXT           : 63
# MABS_REF               : 11
# AUTO_ADAPT             : NO
# ADAPT_BAND             : 5
# ADAPT_LIM              : 1.500000 23.000000
# ADAPT_ZBIN             : 0.010000 6.000000
# ZFIX                   : NO
# SPEC_OUT               : NO
# CHI_OUT                : NO
# PDZ_OUT                : test
#######################################
Reading input librairies ...
Read lib
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.bin
 Done with the library reading with 1830 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.bin
 Done with the library reading with 2084 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.bin
 Done with the library reading with 3975 SED read.
Read lib out
Done reading the FIR lib
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor      VEGA   CALIB    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195  -20.8527       0      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860  -20.7272       0      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466  -21.5143       0      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942  -22.2286       0      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169  -22.6733       0      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534  -22.9145       0      1.0000
[29]:
print(
    f"BPC calculation took {t2-t1:.0f} seconds or {(t2-t1)/band_pass_correction.shape[0]*1000:.3f} milliseconds per model"
)
BPC calculation took 1 seconds or 0.366 milliseconds per model

Note that that function reran prepare

[30]:
band_pass_correction
[30]:
array([0.88856797, 0.93477042, 0.9319388 , ..., 0.70191983, 0.759616  ,
       0.80292392], shape=(3975,))
[31]:
plt.hist(band_pass_correction, bins=30)
plt.xlabel("bpc$_{SED}$")
[31]:
Text(0.5, 0, 'bpc$_{SED}$')
../_images/notebooks_Galametz_reddening_41_1.png
[32]:
elliptical = gals & (model_numbers < 8)

The following figures correspond to figure 3 in Galametz et al.

[33]:
plt.scatter(np.array(photz.zLib)[elliptical], band_pass_correction[elliptical])
plt.plot([0, 2], [1, 1])
plt.xlim([0, 2])
plt.ylim([0.8, 1.2])
plt.xlabel("z")
plt.ylabel("bpc$_{SED}$")
plt.title("Ellipticals")
[33]:
Text(0.5, 1.0, 'Ellipticals')
../_images/notebooks_Galametz_reddening_44_1.png
[34]:
plt.scatter(np.array(photz.zLib)[elliptical], band_pass_correction[elliptical], s=0.2)
plt.plot([0, 6], [1, 1])
plt.xlim([0, 6])
plt.ylim([0.2, 2])
plt.xlabel("z")
plt.ylabel("bpc$_{SED}$")
plt.title("Ellipticals")
[34]:
Text(0.5, 1.0, 'Ellipticals')
../_images/notebooks_Galametz_reddening_45_1.png
[35]:
albd_lib.shape, band_pass_correction.shape
[35]:
((3975, 6), (3975,))

Check the redshift dependence of the values

[36]:
# We will look at redshifts to 2, 3 and the limit for the run at 6
red_2 = np.array(photz.zLib) < 2
red_3 = np.array(photz.zLib) < 3

We adopt an example E(B-V) of 0.1 to show the impact of the band pass correction as was done in Figure 8 of Galametz et al. (2017)

[37]:
galebv_ex = 0.1
[38]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
range = [0, 0.6]
bin_width = 0.005
bins = np.arange(range[0], range[1] + bin_width, bin_width)

m = gals & red_2
for n, b in enumerate("ugrizy"):
    data = albd_lib[m].T[n]
    data = data * galebv_ex
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data),  # / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n] * galebv_ex, color=colors[n], linestyle="-.")
plt.legend(loc="upper right")
plt.ylabel("# models per 0.005 mag")
plt.xlabel("Extinction $A_X$ [mag]")
plt.title("Galaxies at z<2 without bpc")
[38]:
Text(0.5, 1.0, 'Galaxies at z<2 without bpc')
../_images/notebooks_Galametz_reddening_51_1.png
[39]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
range = [0, 0.6]
bin_width = 0.005
bins = np.arange(range[0], range[1] + bin_width, bin_width)
galebv_ex = 0.1
m = gals & red_2
for n, b in enumerate("ugrizy"):
    data = albd_lib[m].T[n]
    data = data * galebv_ex / band_pass_correction[m]
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data),  # / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n] * galebv_ex, color=colors[n], linestyle="-.")
plt.legend(loc="upper right")
plt.ylabel("# models per 0.005 mag")
plt.xlabel("Extinction $A_X$ [mag]")
plt.title("Galaxies at z<2")
[39]:
Text(0.5, 1.0, 'Galaxies at z<2')
../_images/notebooks_Galametz_reddening_52_1.png

Or alternatively use an artificial E(B-V) of 1 to plot the values scaled by the E(B-V) of the dust map

[40]:
galebv_ex = 1.0
[41]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
plt.rcParams.update({"font.size": 14})
range = [0, 6]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)
galebv_ex = 1.0
m = gals & red_3
for n, b in enumerate("ugrizy"):
    data = albd_lib[m].T[n]
    data = data * galebv_ex / band_pass_correction[m]  # *1.018
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data),  # / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n] * galebv_ex, color=colors[n], linestyle="--")
print("Galaxies to z=3")
plt.legend(loc="upper right", fontsize=8)
plt.ylabel("# models per 0.05 mag")
plt.xlabel("Extinction $A_X/E(B-V)$ [mag]")
plt.title("Galaxy templates $z<3$")
plt.xlim(range)
Galaxies to z=3
[41]:
(0.0, 6.0)
../_images/notebooks_Galametz_reddening_55_2.png
[42]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
plt.rcParams.update({"font.size": 14})
range = [0, 6]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)

m = qsos & red_3
for n, b in enumerate("ugrizy"):
    data = albd_lib[m].T[n]
    data = data * galebv_ex / band_pass_correction[m]  # *1.018
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data),  # / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n] * galebv_ex, color=colors[n], linestyle="--")
plt.legend(loc="upper right", fontsize=8)
plt.ylabel("# models per 0.05 mag")
plt.xlabel("Extinction $A_{X}/E$(B-V) [mag]")
plt.title("QSO templates $z<3$")
plt.xlim(range)
print("QSO")
QSO
../_images/notebooks_Galametz_reddening_56_1.png
[43]:
albd_lib
[43]:
array([[4.81022476, 3.66544492, 2.60017753, 1.92243948, 1.46794617,
        1.23895349],
       [4.84567432, 3.63216288, 2.68476487, 1.95107156, 1.46888891,
        1.23828616],
       [4.83604783, 3.6288114 , 2.66334591, 1.90651399, 1.46737265,
        1.23913176],
       ...,
       [4.47485   , 3.08622   , 2.29484432, 1.84253334, 1.47031246,
        1.24108416],
       [4.47485   , 3.08622   , 2.27006822, 1.88519229, 1.46963274,
        1.2415572 ],
       [4.47485   , 3.08622   , 2.24487882, 1.92556222, 1.46423814,
        1.24222454]], shape=(3975, 6))
[44]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
plt.rcParams.update({"font.size": 14})
range = [0, 6]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)
galebv_ex = 1.0
m = qsos
for n, b in enumerate("ugrizy"):
    data = albd_lib[m].T[n]
    data = data * galebv_ex / band_pass_correction[m]  # *1.018
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data),  # / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n] * galebv_ex, color=colors[n], linestyle="--")
plt.legend(loc="upper right", fontsize=8)
plt.ylabel("# models per 0.05 mag")
plt.xlabel("Extinction $A_{X}/E$(B-V) [mag]")
plt.title("QSO templates $z<6$")
plt.xlim(range)
print("QSO")
QSO
../_images/notebooks_Galametz_reddening_58_1.png
[45]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
prop_cycle = plt.rcParams["axes.prop_cycle"]
colors = prop_cycle.by_key()["color"]
plt.rcParams.update({"font.size": 14})
range = [0, 6]
bin_width = 0.05
bins = np.arange(range[0], range[1] + bin_width, bin_width)
galebv_ex = 1.0
m = gals
for n, b in enumerate("ugrizy"):
    data = albd_lib[m].T[n]
    data = data * galebv_ex / band_pass_correction[m]  # *1.018
    counts, bins = np.histogram(data, bins=bins)
    plt.hist(
        data,
        bins=bins,
        weights=np.ones_like(data),  # / counts.max(),
        alpha=0.5,
        label=b,
        color=colors[n],
        linewidth=2,
    )
    plt.axvline(baseline_albd[n] * galebv_ex, color=colors[n], linestyle="--")
plt.legend(loc="upper right")
plt.ylabel("# models per 0.05 mag")
plt.xlabel("Extinction $A_{X}/E$(B-V) [mag]")
plt.title("Galaxy templates $z<6$")
plt.xlim(range)
[45]:
(0.0, 6.0)
../_images/notebooks_Galametz_reddening_59_1.png

Look at some individual models and how the BPC changes with redshift

This is equivalent to figure 3 in Galametz et al. (2017)

[46]:
fig_size = 5
fig, ax = plt.subplots(figsize=(fig_size, fig_size))
plt.rcParams.update({"font.size": 14})
gal_models = [1, 13]
gal_names = ["Ell", "Sc"]

qso_models = [7, 28]
qso_names = ["Seyfert 1.8", "QSO"]
galebv_ex = 0.1
band_number = 1
# Lets just look at objects with intrinsic ebv =0 for clarity
model_ebvs = np.array([m.ebv for m in photz.fullLib])
range = [0, 6]
for n, g in enumerate(gal_models):
    # print(g)
    mod_mask = np.array([s.nummod == g for s in photz.fullLib])
    mod_mask &= np.array([s.is_gal() for s in photz.fullLib])
    mod_mask &= model_ebvs == 0.0
    mod_reds = np.array(photz.zLib)[mod_mask]
    # data=albd_lib[mod_mask].T[band_number]
    # data=data*galebv_ex*band_pass_correction[mod_mask]
    data = band_pass_correction[mod_mask]
    plt.plot(mod_reds, data, label=gal_names[n])

for n, q in enumerate(qso_models):
    # print(q)
    mod_mask = np.array([s.nummod == q for s in photz.fullLib])
    mod_mask &= np.array([s.is_qso() for s in photz.fullLib])
    mod_mask &= model_ebvs == 0.0
    mod_reds = np.array(photz.zLib)[mod_mask]
    # data=albd_lib[mod_mask].T[band_number]
    # data=data*galebv_ex*band_pass_correction[mod_mask]
    data = band_pass_correction[mod_mask]
    plt.plot(mod_reds, data, label=qso_names[n])

# plt.plot(range,[baseline_albd[band_number],baseline_albd[band_number]],c='k',linestyle='--')
plt.plot(range, [1, 1], c="k", linestyle="--")
plt.xlim(range)
plt.ylim([0.0, 1.2])
plt.xlabel("Redshift")
plt.ylabel("band pass correction")
plt.legend(fontsize=8)
[46]:
<matplotlib.legend.Legend at 0x7ab3bfa76470>
../_images/notebooks_Galametz_reddening_61_1.png

3. Compute photoz and compare

We want to look at some actual outputs redshifts to see how it impacts results

[47]:
# Most of the other example use "dereddened" fluxes to account for redening using the traditional method.
# Here we need to use the raw observations
input_table = Table.read(os.path.join(lp.LEPHAREDIR, "examples/COSMOS_no_mw_correction.in"), format="ascii")[
    :500
]
[48]:
cols = [input_table.colnames[0]] + input_table.colnames[3:15] + input_table.colnames[-3:]
input_table = input_table[cols]
[49]:
input_table[:5]
[49]:
Table length=5
col1col4col5col6col7col8col9col10col11col12col13col14col15col62col63col64
int64float64float64float64float64float64float64float64float64float64float64float64float64int64float64float64
11.2264058340024464e-319.171458480906563e-321.5994573122457901e-315.550314043862333e-324.273789009872182e-317.3059942026841e-322.199372687959132e-309.92560863662011e-321.938781577105002e-301.4730192590803836e-311.8850586640859814e-302.9508178430923404e-3104.8393004.0
23.907439142472311e-303.1700856255376637e-315.811632033403237e-301.30412822743503e-311.160616015380427e-291.7439220414787968e-312.042351356282506e-292.0894468509595e-313.838527337807897e-293.328404966288332e-314.928168245607709e-295.0920327599654434e-3101.16611004.0
31.72725269651295e-283.351730770889439e-319.60022921610572e-282.409428450315479e-302.458391796111975e-273.7469482404777074e-303.734917422551279e-274.1404168565905965e-304.9349166078363446e-275.5209072659251e-305.705118718570144e-276.651405949947064e-3000.164911004.0
46.161004518012826e-321.0139706124712776e-318.253574157877361e-336.741867998505554e-322.990173864155611e-318.590213552496993e-321.3406607110676255e-301.0944768116254787e-312.1480934866361488e-301.7493965537694464e-312.1680763263628393e-303.647394401392443e-3105.7423004.0
54.479272945559763e-318.644693601427378e-326.541184076450121e-316.719572984900683e-327.514292154298114e-317.904105124039268e-328.203433710829266e-319.261041729888379e-321.1335069151250403e-301.4611714515184457e-315.462789790785417e-313.068707231140827e-3105.7453003.0
[50]:
len(input_table)
[50]:
500

Run process with the reddening per model

We do this with the reddened measurements in addition to without and using the reddened models instead

[51]:
lp.prepare(config)
t1 = time.time()
out_galametz, _ = lp.process(config, input_table)
t2 = time.time()
# Run again using the classic method for comparison
lp.prepare({**config, "APPLY_MW_EXTINCTION": "CLASSIC"})
t3 = time.time()
out_classic, _ = lp.process({**config, "APPLY_MW_EXTINCTION": "CLASSIC"}, input_table)
t4 = time.time()
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor    TG-cor      VEGA M_sun(AB)   CALIB      Lb_eff    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195   -0.2745  -20.8527    6.3135       0      0.3668      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860   -0.2458  -20.7272    5.0764       0      0.4780      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466    0.2580  -21.5143    4.6523       0      0.6178      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942    0.6138  -22.2286    4.5323       0      0.7666      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169    0.7625  -22.6733    4.5147       0      0.8903      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534    0.7763  -22.9145    4.5081       0      0.9782      1.0000
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : STAR
# STAR_SED    :/home/docs/.cache/lephare/data/sed/STAR/STAR_MOD_ALL.list
# STAR_LIB    :LIB_STAR
# STAR_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_STAR.doc
# STAR_FSCALE :0.0000
#######################################
Number of templates in the list 254
Number of SED in the list (including different ages) 254
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# COSMOLOGY   :70.0000,0.3000,0.7000
# STAR_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_STAR(.doc & .bin)
# STAR_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS(.doc & .bin)
# LIB_ASCII   NO
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# CREATION_DATE Wed Oct 28 09:41:12 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : GAL
# GAL_SED    :/home/docs/.cache/lephare/data/sed/GAL/COSMOS_SED/COSMOS_MOD.list
# GAL_LIB    :LIB_CE
# GAL_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.doc
# GAL_LIB phys:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.phys
# SEL_AGE    :none
# GAL_FSCALE :1.0000
# AGE_RANGE   0.0000 15000000000.0000
#######################################
Number of templates in the list 31
Number of SED in the list (including different ages) 31
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# GAL_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_CE(.doc & .bin)
# GAL_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/CE_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# EM_LINES   EMP_UV
# EM_DISPERSION   1.0000,
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:41:12 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : QSO
# QSO_SED    :/home/docs/.cache/lephare/data/sed/QSO/SALVATO09/AGN_MOD.list
# QSO_LIB    :LIB_QSO
# QSO_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_QSO.doc
# QSO_FSCALE :1.0000
#######################################
Number of templates in the list 30
Number of SED in the list (including different ages) 30
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# QSO_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_QSO(.doc & .bin)
# QSO_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# MW_GALAMETZ   YES
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:41:13 2026
#############################################
Using user columns from input table assuming they are in the standard MEME order.
Processing 500 objects with 6 bands
#######################################
# PHOTOMETRIC REDSHIFT with OPTIONS   #
# Config file            :
# CAT_IN                 : /home/docs/.cache/lephare/data/examples/COSMOS_no_mw_correction.in
# CAT_OUT                : zphot.out
# CAT_LINES              : 0 1000000000
# PARA_OUT               : /home/docs/.cache/lephare/data/examples/output.para
# INP_TYPE               : F
# CAT_FMT[0:MEME 1:MMEE] : 0
# CAT_MAG                : AB
# ZPHOTLIB               : CE_COSMOS STAR_COSMOS QSO_COSMOS
# APPLY_MW_EXTINCTION    : GALAMETZ
# MW_GLOBAL_EBV          : -1.000000
# FIR_LIB                :
# FIR_LMIN               : 7.000000
# FIR_CONT               : -1.000000
# FIR_SCALE              : -1.000000
# FIR_FREESCALE          : YES
# FIR_SUBSTELLAR         : NO
# ERR_SCALE              : 0.020000
# ERR_FACTOR             : 1.500000
# GLB_CONTEXT            : 63
# FORB_CONTEXT           : -1
# DZ_WIN                 : 1.000000
# MIN_THRES              : 0.020000
# MAG_ABS                : -24.000000 -5.000000
# MAG_ABS_AGN            : -30.000000 -10.000000
# MAG_REF                : 3
# NZ_PRIOR               : -1 -1
# Z_INTERP               : YES
# Z_METHOD               : BEST
# MABS_METHOD            : 1
# MABS_CONTEXT           : 63
# MABS_REF               : 11
# AUTO_ADAPT             : NO
# ADAPT_BAND             : 5
# ADAPT_LIM              : 1.500000 23.000000
# ADAPT_ZBIN             : 0.010000 6.000000
# ZFIX                   : NO
# SPEC_OUT               : NO
# CHI_OUT                : NO
# PDZ_OUT                : test
#######################################
Reading input librairies ...
Read lib
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.bin
 Done with the library reading with 1830 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.bin
 Done with the library reading with 2084 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.bin
 Done with the library reading with 3975 SED read.
Read lib out
Done reading the FIR lib
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor      VEGA   CALIB    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195  -20.8527       0      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860  -20.7272       0      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466  -21.5143       0      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942  -22.2286       0      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169  -22.6733       0      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534  -22.9145       0      1.0000
WARNING: No scaling --> No z 255
WARNING: No scaling --> No z 256
WARNING: No scaling --> No z 257
WARNING: No scaling --> No z 261
WARNING: No scaling --> No z 358
WARNING: No scaling --> No z 387
WARNING: No scaling --> No z 408
# Offsets added to the modeled magnitudes (or substracted to the observed): 0.0,0.0,0.0,0.0,0.0,0.0

# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor    TG-cor      VEGA M_sun(AB)   CALIB      Lb_eff    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195   -0.2745  -20.8527    6.3135       0      0.3668      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860   -0.2458  -20.7272    5.0764       0      0.4780      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466    0.2580  -21.5143    4.6523       0      0.6178      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942    0.6138  -22.2286    4.5323       0      0.7666      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169    0.7625  -22.6733    4.5147       0      0.8903      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534    0.7763  -22.9145    4.5081       0      0.9782      1.0000
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : STAR
# STAR_SED    :/home/docs/.cache/lephare/data/sed/STAR/STAR_MOD_ALL.list
# STAR_LIB    :LIB_STAR
# STAR_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_STAR.doc
# STAR_FSCALE :0.0000
#######################################
Number of templates in the list 254
Number of SED in the list (including different ages) 254
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# COSMOLOGY   :70.0000,0.3000,0.7000
# STAR_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_STAR(.doc & .bin)
# STAR_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS(.doc & .bin)
# LIB_ASCII   NO
# MW_GALAMETZ   NO
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# CREATION_DATE Wed Oct 28 09:41:14 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : GAL
# GAL_SED    :/home/docs/.cache/lephare/data/sed/GAL/COSMOS_SED/COSMOS_MOD.list
# GAL_LIB    :LIB_CE
# GAL_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.doc
# GAL_LIB phys:/home/docs/.cache/lephare/work/lib_bin/LIB_CE.phys
# SEL_AGE    :none
# GAL_FSCALE :1.0000
# AGE_RANGE   0.0000 15000000000.0000
#######################################
Number of templates in the list 31
Number of SED in the list (including different ages) 31
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# GAL_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_CE(.doc & .bin)
# GAL_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/CE_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# EM_LINES   EMP_UV
# EM_DISPERSION   1.0000,
# MW_GALAMETZ   NO
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:41:14 2026
#############################################
#######################################
# It s translating SEDs to binary library #
# with the following options :
# Config file     :
# Library type     : QSO
# QSO_SED    :/home/docs/.cache/lephare/data/sed/QSO/SALVATO09/AGN_MOD.list
# QSO_LIB    :LIB_QSO
# QSO_LIB doc:/home/docs/.cache/lephare/work/lib_bin/LIB_QSO.doc
# QSO_FSCALE :1.0000
#######################################
Number of templates in the list 30
Number of SED in the list (including different ages) 30
WRONG NUMBER OF ARGUMENTS FOR OPTION MOD_EXTINC
We have 2 instead of 8
Use default value 0,0 for all filters
#######################################
# It s computing the SYNTHETIC MAGNITUDES #
# For Gal/QSO libraries with these OPTIONS #
# with the following options :
# Config file     :
# Filter file     : filter_cosmos
# Magnitude type     : AB
# QSO_LIB_IN    :/home/docs/.cache/lephare/work/lib_bin/LIB_QSO(.doc & .bin)
# QSO_LIB_OUT   :/home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS(.doc & .bin)
# Z_STEP   :0.1000 0.0000 6.0000
# COSMOLOGY   :70.0000,0.3000,0.7000
# EXTINC_LAW   :SMC_prevot.dat SB_calzetti.dat SB_calzetti_bump1.dat SB_calzetti_bump2.dat
# MOD_EXTINC   :0 0 0 0 0 0 0 0
# EB_V   :0.0000 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.5000
# MW_GALAMETZ   NO
# EXT_MW_CURVE   LMC_Fitzpatrick.dat
# LIB_ASCII   NO
# CREATION_DATE Wed Oct 28 09:41:14 2026
#############################################
Using user columns from input table assuming they are in the standard MEME order.
Processing 500 objects with 6 bands
#######################################
# PHOTOMETRIC REDSHIFT with OPTIONS   #
# Config file            :
# CAT_IN                 : /home/docs/.cache/lephare/data/examples/COSMOS_no_mw_correction.in
# CAT_OUT                : zphot.out
# CAT_LINES              : 0 1000000000
# PARA_OUT               : /home/docs/.cache/lephare/data/examples/output.para
# INP_TYPE               : F
# CAT_FMT[0:MEME 1:MMEE] : 0
# CAT_MAG                : AB
# ZPHOTLIB               : CE_COSMOS STAR_COSMOS QSO_COSMOS
# APPLY_MW_EXTINCTION    : CLASSIC
# MW_GLOBAL_EBV          : -1.000000
# FIR_LIB                :
# FIR_LMIN               : 7.000000
# FIR_CONT               : -1.000000
# FIR_SCALE              : -1.000000
# FIR_FREESCALE          : YES
# FIR_SUBSTELLAR         : NO
# ERR_SCALE              : 0.020000
# ERR_FACTOR             : 1.500000
# GLB_CONTEXT            : 63
# FORB_CONTEXT           : -1
# DZ_WIN                 : 1.000000
# MIN_THRES              : 0.020000
# MAG_ABS                : -24.000000 -5.000000
# MAG_ABS_AGN            : -30.000000 -10.000000
# MAG_REF                : 3
# NZ_PRIOR               : -1 -1
# Z_INTERP               : YES
# Z_METHOD               : BEST
# MABS_METHOD            : 1
# MABS_CONTEXT           : 63
# MABS_REF               : 11
# AUTO_ADAPT             : NO
# ADAPT_BAND             : 5
# ADAPT_LIM              : 1.500000 23.000000
# ADAPT_ZBIN             : 0.010000 6.000000
# ZFIX                   : NO
# SPEC_OUT               : NO
# CHI_OUT                : NO
# PDZ_OUT                : test
#######################################
Reading input librairies ...
Read lib
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/QSO_COSMOS.bin
 Done with the library reading with 1830 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/STAR_COSMOS.bin
 Done with the library reading with 2084 SED read.
Number of keywords to be read in the doc: 15
Number of keywords read at the command line (excluding -c config): 0
Reading keywords from /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.doc
Number of keywords read in the config file: 18
Keyword NUMBER_ROWS not provided
Keyword NUMBER_SED not provided
Keyword Z_FORM not provided
Reading library: /home/docs/.cache/lephare/work/lib_mag/CE_COSMOS.bin
 Done with the library reading with 3975 SED read.
Read lib out
Done reading the FIR lib
# NAME                        IDENT      Lbda_mean Lbeff(Vega)       FWHM     AB-cor      VEGA   CALIB    Fac_corr
u_new.pb                      1             3689.8      3749.5       456.0    0.6195  -20.8527       0      1.0000
gHSC.pb                       2             4850.6      4760.2      1194.2   -0.0860  -20.7272       0      1.0000
rHSC.pb                       3             6241.4      6142.0      1538.5    0.1466  -21.5143       0      1.0000
iHSC.pb                       4             7716.5      7636.7      1476.1    0.3942  -22.2286       0      1.0000
zHSC.pb                       5             8915.2      8906.7       768.0    0.5169  -22.6733       0      1.0000
yHSC.pb                       6             9801.2      9770.6       796.8    0.5534  -22.9145       0      1.0000
WARNING: No scaling --> No z 255
WARNING: No scaling --> No z 256
WARNING: No scaling --> No z 257
WARNING: No scaling --> No z 261
WARNING: No scaling --> No z 358
WARNING: No scaling --> No z 387
WARNING: No scaling --> No z 408
# Offsets added to the modeled magnitudes (or substracted to the observed): 0.0,0.0,0.0,0.0,0.0,0.0

[52]:
print(f"Time per object with GALAMETZ reddening is a factor  {(t4-t3)/(t2-t1):.2f} larger")
Time per object with GALAMETZ reddening is a factor  0.75 larger
[53]:
# How many results changed by more than 0.01?
zspec = out_classic["ZSPEC"]
z1 = out_classic["Z_BEST"]
z2 = out_galametz["Z_BEST"]
np.sum(np.abs(z1 - z2) > 0.01) / len(z1)
[53]:
np.float64(0.01)
[54]:
m = np.abs(z1 - z2) > 0.001
plt.plot([0, 2], [0, 2], c="r")
plt.scatter(z1[m], z2[m], s=1)
plt.xlabel("z (Traditional method)")
plt.ylabel("z (Galametz)")
plt.xlim([0, 2])
plt.ylim([0, 2])
[54]:
(0.0, 2.0)
../_images/notebooks_Galametz_reddening_71_1.png

How does the difference compare to general comparison with spectroscopic redshift

At these values of \(E(B-V)\) the differences due to the method are smaller than the intrinsic scatter.

[55]:
plt.plot([0, 6], [0, 6], c="r")
plt.scatter(zspec[m], z1[m], s=1)
# plt.plot([0,2],[0,2])
plt.ylabel("z (Traditional method)")
plt.xlabel("z (Spec)")
plt.xlim([0, 6])
plt.ylim([0, 6])
[55]:
(0.0, 6.0)
../_images/notebooks_Galametz_reddening_73_1.png
[56]:
plt.plot([0, 6], [0, 6], c="r")
plt.scatter(zspec[m], z2[m], s=1)
# plt.plot([0,2],[0,2])
plt.ylabel("z (Galametz)")
plt.xlabel("z (Spec)")
plt.xlim([0, 6])
plt.ylim([0, 6])
[56]:
(0.0, 6.0)
../_images/notebooks_Galametz_reddening_74_1.png

How many models changed as a result of the new method?

[57]:
same_gal = out_classic["MOD_BEST"] == out_galametz["MOD_BEST"]
same_qso = out_classic["MOD_QSO"] == out_galametz["MOD_QSO"]
[58]:
np.sum(~same_gal) / len(same_gal)
[58]:
np.float64(0.022)
[59]:
np.sum(~same_qso) / len(same_gal)
[59]:
np.float64(0.014)
[60]:
plt.scatter(z1, z2 - z1, s=5)
plt
plt.xlim([0, 6])
plt.ylim([-0.05, 0.05])
plt.plot([0, 6], [0, 0], c="r")
plt.xlabel("Z_BEST (Traditional)")
plt.ylabel("$\Delta$ Z_BEST (Galametz- Traditional)")

plt.title("Galaxies")
[60]:
Text(0.5, 1.0, 'Galaxies')
../_images/notebooks_Galametz_reddening_79_1.png

Compare stats for the traditional method and the galaxy sample

[61]:
def sigma_nmad(z1, z2):
    maskPos = (z1 > 0.02) & (z2 > 0)  # & (z1 < 1)
    delz = (z2 - z1) / (1 + z1)
    bias = np.nanmedian(delz[maskPos])
    med_delz = np.nanmedian(z2[maskPos] - z1[maskPos])
    delz_unbiased = (z2 - z1 - med_delz) / (1 + z1)
    outlier_frac = np.sum(np.abs(delz[maskPos]) > 0.15) / np.sum(maskPos)
    sigma_nmad_unbiased = 1.48 * np.nanmedian(np.abs(delz_unbiased[maskPos]))

    return {"sigma_nmad_unbiased": sigma_nmad_unbiased, "outlier_frac": outlier_frac, "bias": bias}
[62]:
# Classic
sigma_nmad(zspec, z1)
[62]:
{'sigma_nmad_unbiased': np.float64(0.05884940409665664),
 'outlier_frac': np.float64(0.11648351648351649),
 'bias': np.float64(-0.04752071949495798)}
[63]:
# Galametz
sigma_nmad(zspec, z2)
[63]:
{'sigma_nmad_unbiased': np.float64(0.058581674184708936),
 'outlier_frac': np.float64(0.11428571428571428),
 'bias': np.float64(-0.0473367373417308)}

Stats for the Galametz method and the galaxy sample

How long does the model reddening take?

[64]:
id, flux, flux_err, context, zspec, string_data = lp.table_to_data(config, input_table)
i = 0
one_obj = lp.onesource(i, photz.gridz)
one_obj.readsource(str(id[i]), flux[i], flux_err[i], context[i], zspec[i], str(string_data[i]))
photz.prep_data(one_obj)
Using user columns from input table assuming they are in the standard MEME order.
[65]:
%timeit redened_flux=one_obj.redden_flux(photz.flux,albd_lib)
3.05 ms ± 6.11 μs per loop (mean ± std. dev. of 7 runs, 100 loops each)

How does that compare to muliplying using numpy?

[66]:
%timeit red_num=photz.flux/10**(albd_lib/2.5)
1.89 ms ± 4.32 μs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)