© The Authors 2024

1. Introduction

Classical Cepheid variables in the Milky Way (MW) provide the absolute calibration of their period–luminosity (PL) relation or the Leavitt law (Leavitt & Pickering 1912) based on their individual distances derived from geometric parallaxes. The absolute calibration of the Leavitt law is crucial in order to measure extragalactic distances and determine the present expansion rate of the Universe through the ‘cosmic distance ladder’ (Freedman et al. 2001; Riess et al. 2022b). The present expansion rate in the late evolutionary universe, the local Hubble constant, is currently in discord with its early universe measurement from the Planck mission (Riess et al. 2022b; Planck Collaboration VI 2020). This tension between the Hubble constant measurements from the two extreme ends of the Universe may hint at new or missing physics in the standard cosmological model (Di Valentino et al. 2021; Abdalla et al. 2022).

The accuracy and the precision of Cepheid-based distance measurements in the traditional cosmic distance ladder are now only limited by the systematic uncertainties in the calibration of their PL relation and its metallicity dependence (Riess et al. 2022b; Bhardwaj et al. 2023). Before Gaia, the absolute calibration of Cepheid PL relations was mostly based on a limited sample of Cepheids with Hubble Space Telescope parallaxes (Benedict et al. 2007; Riess et al. 2014, 2018) or using distances determined from independent methods, such as the Baade-Wesselink or infrared surface brightness methods (Gieren et al. 1998; Fouqué et al. 2007; Storm et al. 2011; Bhardwaj et al. 2016). This is now changing thanks to increasingly more accurate and precise geometric parallaxes of thousands of Cepheids in the MW from the Gaia mission (Gaia Collaboration 2016, 2023b; Ripepi et al. 2023). Several recent studies have utilized unprecedentedly precise Gaia parallaxes to calibrate PL relations at multiple wavelengths (e.g. Groenewegen 2018; Ripepi et al. 2021, 2023; Breuval et al. 2022; Riess et al. 2022a; Cruz Reyes & Anderson 2023; Narloch et al. 2023). It has been claimed that some of these calibrations of Cepheid luminosities have reached percent-level precision at specific wavelengths (Riess et al. 2022a; Cruz Reyes & Anderson 2023).

The increasing precision of Gaia parallaxes also enabled some of the above-mentioned studies to investigate the metallicity dependence of Cepheid PL relations at multiple wavelengths. However, there is currently no consensus on the metallicity dependence of Cepheid period-luminosity-metallicity (PLZ) and period-Wesenheit-metallicity (PWZ) relations at multiple wavelengths. Most of the recent empirical studies have found a negative metallicity coefficient with values between −0.2 and −0.5 mag dex⁻¹ (e.g. Fig. 19 in Trentin et al. 2024), which was also confirmed by models for the Wesenheit relations (Anderson et al. 2016; De Somma et al. 2022). Nevertheless, these coefficients are often either weakly constrained and/or differ by a factor of two among recent empirical studies (e.g. Ripepi et al. 2021; Riess et al. 2021, 2022b; Breuval et al. 2022; Molinaro et al. 2023; Bhardwaj et al. 2023; Trentin et al. 2024). For example, Trentin et al. (2024) found a metallicity coefficient of −0.458 ± 0.052 mag dex⁻¹ in K_s-band for the MW Cepheid PL relation, while Breuval et al. (2022) found a coefficient of −0.321 ± 0.068 mag dex⁻¹ using Cepheids in the Galaxy and the Magellanic Clouds. Recently, Bhardwaj et al. (2023) determined a metallicity term of −0.43 ± 0.18 mag dex⁻¹ in K_s-band using individual metallicities of MW Cepheids. This coefficient is better constrained (−0.33 ± 0.07 mag dex⁻¹) using MW and Large Magellanic Cloud (LMC) Cepheids within a very narrow range of their mean metallicities (Δ[Fe/H] = 0.46 dex). Homogeneous high-resolution spectroscopic metallicities of MW Cepheids having a wide range of metallicity and multiband photometry are needed to complement their Gaia parallaxes to better constrain the zero-points and the metallicity coefficients of Cepheid PLZ relations.

Most of the recent studies on the empirical calibration of Cepheid PL or PLZ relations in BVIJHK_s bands utilized heterogeneous optical and infrared photometry together with homogeneous Gaia data. The optical PL relations, except those in the Gaia photometric systems, primarily used one of the largest compilations of photoelectric observations of 894 Cepheids in standard Johnson-Cousin-Kron filters (UBV(RI)c, Berdnikov 2008) or the modern CCD photometry of Cepheids in BVI bands (e.g. Berdnikov et al. 2015). Time-domain variability surveys such as the Optical Gravitational Lensing Experiment (OGLE, Udalski et al. 2018) and the All-Sky Automated Survey for Supernovae (ASAS-SN, Jayasinghe et al. 2018) have also provided optical V and/or I light curves of thousands of Cepheids in the MW, but lack multiwavelength coverage. Moreover, several ongoing large-scale variability surveys are now being carried out in the Sloan (ugriz) photometric system such as the Zwicky Transient Facility (ZTF, Bellm et al. 2019, in gri), enabling the discovery and identification of hundreds of Cepheids in the MW (Chen et al. 2020).

However, the pulsation properties of Cepheid variables and their PL relations have not yet been explored in detail in the specific Sloan filters. Hoffmann & Macri (2015) used the 8.1m Gemini North telescope to observe Cepheids and derive their PL relations in the gri filters in NGC 4258, which is one of the anchor galaxies for the extragalactic distance scale (Riess et al. 2022b). Kodric et al. (2018) presented PL relations for Cepheids in the Andromeda galaxy in Sloan-like (gri) filters using the data from the Panoramic Survey Telescope And Rapid Response System survey (Tonry et al. 2012). Recently, Adair & Lee (2023) presented the largest sample of more than 1600 Cepheids in M33 using photometry in the gri bands. Narloch et al. (2023) provided light curves of 96 MW Cepheids and their PL relations in the Sloan (gri) filters. Similar calibrations of gri band PL relations for RR Lyrae, Type II Cepheids, and the anomalous Cepheids in the globular clusters have also been provided using the ZTF data (Ngeow et al. 2022a,b,c). In view of the much anticipated Vera C. Rubin Observatory’s Legacy Survey of Space and Time (LSST, Ivezić et al. 2019), which will operate at ugrizy wavelengths, it is crucial to provide an empirical calibration of PL relations for pulsating stars in Sloan filters.

Homogeneous near-infrared (NIR) JHK_s band light curves of MW Cepheids for the largest sample of 131 Cepheids was provided by Monson & Pierce (2011). The empirical calibrations of Cepheid PL relations in JHK_s bands in the literature are either based on the above-mentioned dataset combined with older photometric measurements (e.g. Welch et al. 1984; Laney & Stobie 1994; Barnes et al. 1997) or using random phase-corrected Two Micron All Sky Survey (2MASS) observations (Skrutskie et al. 2006). NIR photometry offers several advantages over optical bands, due to lower extinction and smaller temperature variations leading to less scatter in the PL relations (see the review by Bhardwaj 2020). However, smaller variability amplitudes and near-sinusoidal light curves of Cepheid variables in the NIR bands can complicate their identification and classification, which is easier in the optical domain due to asymmetric features and larger amplitudes (Bhardwaj et al. 2015, 2017; Udalski et al. 2018). Given that all upcoming large observational facilities will mostly operate at infrared wavelengths, it is imperative to obtain homogeneous optical and infrared photometric light curves of Cepheids in the MW to fully utilize current and future Gaia parallaxes for calibrating the Leavitt law at multiple wavelengths.

This paper presents simultaneous optical (griz) and NIR (JHK_s) light curves of 78 classical Cepheids for the first time. The manuscript is the fifth in the series of the Cepheid Metallicity in the Leavitt Law (C-MetaLL) survey, which aims to provide homogeneous time-series photometry and high-resolution spectroscopy of Cepheids in our Galaxy (Ripepi et al. 2021; Trentin et al. 2023). While previous C-MetaLL publications explored the impact of new spectroscopic observations of Cepheids on their PLZ relations, this paper presents the first results of the ongoing photometric programme. The details regarding the photometric observations, the data reduction and analysis, and the photometric calibration are presented in Sect. 2. Multiwavelength light curves of Cepheids and their PL relations are discussed in Sects. 3 and 4, respectively. The impact of metallicity on the Leavitt law is discussed in Sect. 5, and the results of this work are summarized in Sect. 6.

2. Data and photometry

2.1. The sample of Cepheid variables

The sample of stars was selected to obtain complementary time-series data for Cepheids that have high-resolution spectroscopic observations as part of the C-MetaLL survey (Ripepi et al. 2021). Our initial sample consists of 80 Cepheids including 49 fundamental mode (FU) and 31 first-overtone mode (FO) stars. The majority of our Cepheids were identified as variables in the Gaia data releases. However, some of our initial targets, which were included in Ripepi et al. (2021), come from the ASAS-SN or the OGLE survey (Udalski et al. 2018; Jayasinghe et al. 2018). The pulsation periods, epochs of maximum brightness, intensity-weighted mean magnitudes in G, G_BP, and G_RP bands, and G band pulsation amplitudes for 73 Cepheids were taken from Gaia DR3 Cepheid sample (Ripepi et al. 2023). For remaining seven Cepheids, only photometric mean magnitudes in Gaia filters were taken from the source catalogue (Gaia Collaboration 2023b) since they were not included as variables in Gaia. The periods and epochs of maximum brightness, and the V/I band amplitudes for these stars were taken from the OGLE (OGLE-GD-CEP-0307/1210/1227/1286) or ASAS-SN (J163656.55-322102.2, J184816.35+004903.4/V0912 Aql, J193837.90+172322.8) catalogues (Udalski et al. 2018; Jayasinghe et al. 2018). Figure 1 displays the location of Cepheid variables used in this work.

Fig. 1. Spatial distribution of classical Cepheids used in this work in Galactic coordinates. The smaller symbol sizes represent distant Cepheids. Open circles represent stars with no metallicity measurements. The colour bar represents metallicity.

Figure 2 displays the histogram of period distribution of Cepheids in our sample. The FU mode Cepheids cover a period range of 3.15–70.80 days while FO Cepheids have periods between 0.78 and 4.91 days. These Cepheids are bright and cover an apparent magnitude range between G = 8.89 and G = 14.78 mag. The distances to our Cepheids range from 0.5 kpc to 19.7 kpc with a median distance of 3.6 kpc. High-resolution spectroscopic metallicities for 59 Cepheids were already obtained as part of the C-MetaLL survey (Ripepi et al. 2021; Trentin et al. 2023). Among the remaining stars, six more Cepheids have medium-resolution spectroscopic metallicities from Gaia Radial Velocity Spectometer (RVS, Recio-Blanco et al. 2023). However, four of these six Cepheids have poor RVS quality flags and were assigned a large uncertainty of 0.5 dex (Trentin et al. 2024). The middle panel of Fig. 2 shows the histograms of available metallicities for all 65 Cepheids.

Fig. 2. The period and metallicity distributions, and period-Wesenheit relation for MW Cepheids. Top: period distribution of 80 candidate Cepheids in the sample studied in this work. The filled and lined histograms correspond respectively to FU and FO Cepheids. Middle: histogram of metallicities for 65 Cepheids, of which 59 have homogeneous high-resolution spectroscopic metallicities from the C-MetaLL survey (Ripepi et al. 2021; Trentin et al. 2023). Bottom: period-Wesenheit relation in Gaia filters for the sample of 80 stars. The circles and squares represent fundamental and first-overtone mode Cepheids, respectively. Two extreme outliers are marked with their Gaia DR3 source ID.

To confirm the identification and classification of 80 Cepheids in our sample, we derived their PW relation, W_G = 1.90(G_BP − G_RP), in the Gaia bands (Ripepi et al. 2019). The distances to these Cepheids from Bailer-Jones et al. (2021) were used to obtain absolute Wesenheit magnitudes. We noticed two outlier stars (5325503574371253120/OGLE-GD-CEP-0307, 3318591594625903104/OGLE-GD-CEP-1286 in the PW relation shown in the bottom panel of Fig. 2. OGLE-GD-CEP-1286 is classified as rotating variable in the ASAS-SN survey (J061105.05+045303.5). The parallax uncertainties of these stars are less than 3% and these were not classified as Cepheids in the Gaia data. Therefore, we excluded these two stars from our final sample, which consists of 78 stars with 49 FU and 29 FO Cepheid variables.

2.2. Observations and data reduction

Our multiband observations of Cepheid variables were obtained between February 2021 and September 2022 using the Rapid Eye Mount (REM¹) Telescope located in La Silla, Chile. REM is a 60 cm diameter fast reacting telescope with two parallel imaging instruments: ROS2, a visible imager with four simultaneous passbands (griz), and REMIR, an infrared (JHK_s) imaging camera. The pixel scale of ROS2 and REMIR is 0.58 and 1.2 arcsec pixel⁻¹, respectively, for a total field of view of about 10 × 10 arcmin. The time-series observations were obtained in a monitoring mode where each Cepheid was observed once per night, when visible during the semester. Our observations for each epoch consisted of three dithered frames (exposures) in optical and a minimum of five dithered frames in NIR. The exposure times in optical bands were optimized depending on the brightness of the sources, and range from two to 240 s. The individual exposures were relatively short in NIR (between one and ten seconds), such that the background sky variations were negligible. For fainter targets, multiple sequences of five dithered exposures were taken in NIR bands. These Cepheids were observed on average on 24 different nights with the number of observations varying between 10 and 60.

The optical and NIR science images and associated calibrations (bias and flat frames) were downloaded from the REM data Archive². For optical images, the bias subtraction and flat fielding were performed using monthly calibration frames. For NIR images, pre-processed (dark subtraction and flat fielding) images were downloaded including the sky background and co-added images. However, the co-added NIR images were not used because several bright residual patterns were present which likely resulted from over-subtraction of sky background in several pixels. These residual variations were normalized using SExtractor (Bertin 2006) in the sky background images before applying background subtraction. Therefore, we used individual dithered frames in both optical and NIR bands for photometric data reduction.

2.3. Photometry

The photometry was performed on all dithered images separately in each filter using the DAOPHOT/ALLSTAR (Stetson 1987) and ALLFRAME (Stetson 1994) routines. In a given filter, all point-like sources were first identified using SExtractor to determine a full width at half-maximum (FWHM). An aperture photometry was then performed using DAOPHOT with an aperture equal to the median FWHM. Up to ten brightest, non-saturated stars were used to construct an empirical point-spread-function (PSF) excluding sources closer to the detector edges. We selected the first-frame (images of the first-night observations) as reference frame to create a star-list. In cases where there were not enough detections (< 10) in the first-frame, we visually selected images taken on another night to create a reference star list. The PSF photometry was obtained using ALLSTAR on all images and frame-to-frame coordinate transformations between the reference image and all epoch images were derived using DAOMATCH/DAOMASTER. This procedure required the most intervention because of the lack of a sufficient number of stars in images taken in poorer weather conditions resulted in inaccurate transformations. Those frames were excluded for which no transformation could be obtained due to poorer statistics or lack of significant overlap with the reference frame. Finally, the reference-star list and the derived coordinate transformations were used as input to ALLFRAME for performing PSF fitting across all the frames, simultaneously. The photometry on different nights was internally calibrated to the reference frame using secondary standards that have small photometric uncertainties and no epoch-to-epoch variability. These thresholds on photometric uncertainties and epoch-to-epoch variability varied between 0.01 and 0.1 mag depending on the number of sources within the field of view and their photometric precision.

The photometric calibration was performed using Gaia synthetic photometry (GSP, Gaia Collaboration 2023a) in Sloan filters (griz) and 2MASS (Skrutskie et al. 2006) catalogue in the JHK_s bands. We selected homogeneous GSP catalogue for calibrating optical photometry because there were not enough common stars for several sources in the American Association of Variable Star Observers (AAVSO) Photometric All-Sky Survey (APASS) DR10 (Henden 2019). Moreover, most of these bright Cepheids were saturated in the Sloan Digital Sky Survey catalogue. For each Cepheid location, all sources from GSP and 2MASS catalogues were extracted within the field of view of 10 × 10 arcmin. Our photometric catalogues in a given filter were then cross-matched with GSP or 2MASS within an initial tolerance of 1.0″. The common sources were used to determine a zero-point offset to calibrate our instrumental magnitudes to Sloan or 2MASS filters. The number of calibrating stars varied between five and 117 for different star or filter, and there were fewer common stars in gz-bands. The lack of number of stars covering a wide colour range precluded us from solving for a colour-term in the calibration equation. The uncertainties in the magnitude zero-point offset were propagated to the instrumental magnitudes for each Cepheid variable.

Figure 3 shows a comparison of our calibrated mean magnitudes of Cepheids (see Sect. 3.2) with the magnitudes from APASS DR10 and 2MASS for griz and JHK_s filters, respectively. The offsets are small in ri bands and increase in gz bands and at NIR wavelengths, where magnitudes from the literature are based on random-epoch observations. In the case of NIR magnitudes, the offsets increase significantly at the bright end (K_s < 8.5 mag). We note that the stars with 4 < K_s < 8.5 mag saturate in the primary 1.3 s 2MASS exposure and their NIR magnitudes come from aperture photometry on 51 ms exposure images (Skrutskie et al. 2006). The largest offset occurs at the brightest end in the H-band. For H < 6.5 mag, the median errors in the REMIR J/H/K_s magnitudes of these Cepheids are 0.05/0.09/0.08 mag, hinting at possible non-linearity and saturation at the bright end in HK_s. The typical uncertainties on 2MASS magnitudes are 0.02–0.03 mag, and thus the offsets at the bright end are still within 2σ of the combined REMIR and 2MASS uncertainties. However, we also noted increased scatter in the light curves around maximum light for Cepheids with K_s < 5 mag. In these cases, the brightest source in the REMIR images is typically the target Cepheid, and the PSF is generally constructed from fainter stars. The aperture photometry may be a better option for these brightest Cepheids, and will be explored in the subsequent photometric studies of the C-MetaLL survey.

Fig. 3. Comparison of the calibrated mean magnitudes of Cepheids in this work with the magnitudes from APASS DR10 (in griz) and 2MASS (in JHKs). The dashed and dotted lines respectively represent zero and mean offset. The mean (and the standard deviation) offsets for the bright (Ks < 8.5 mag) and faint (Ks ≥ 8.5 mag) Cepheids are also mentioned at the top left and bottom right of each panel. The filled circles and squares represent FU and FO Cepheids, respectively.

3. Optical and near-IR light curves of Cepheid variables

3.1. Pulsation periods and the phase-binned light curves

The photometric light curves of Cepheids in both optical (griz) and NIR (JHK_s) filters were utilized to determine their pulsation periods using the multiband periodogram of Saha & Vivas (2017). This hybrid algorithm for period determination using sparsely sampled light curve data at multiple wavelengths was employed to search periods in the range from 0.1 to 100 days with a stepsize of 0.001 days. An average difference of 0.002 days with a scatter of 0.012 days was found between our periods and those adopted from the literature after excluding a few spuriously determined periods. The difference in periods exceeded 0.1 days for 17 Cepheids, all of which either had large phase gaps in their light curves or poor light curve quality in at least one filter. Since all the Cepheids in our sample are bright (G < 14.5 mag), we adopted literature periods primarily from Gaia DR3 (Ripepi et al. 2023) that were based on much longer temporal baseline than our observations.

Multiwavelength light curves were phased using the adopted periods and the epochs of maximum brightness. Since we obtained multiple dithered frames in both optical and NIR filters and performed photometry on individual frames, there are often several data points at a given phase albeit with larger uncertainties due to lower signal-to-noise of individual exposures. Furthermore, the scatter at a given phase is larger for stars with a small number of secondary standards within their imaged field of view since these were used to obtain frame to frame transformations and also to calibrate the photometry. Figure 4 displays example light curves of FU and FO Cepheids with all epoch observations (in grey). The scatter in the light curves is evident which comes from larger errors due to low signal-to-noise of dithered frames and uncertainties in relative zero-points of photometry from different nights. Therefore, we decided to bin the light curves in phase improving the accuracy and precision of photometric data points in a given phase bin. We computed sliding mean values with a bin width of 0.1 in phase with 15/10 steps for FU/FO Cepheids. While 10 phase averaged data points were sufficient for sinusoidal nature of FO Cepheid light curves, we used 15 steps for FU Cepheids to fully recover the saw-tooth feature of their light curves. We found these choices recovered the sharp extrema of the light curves and also yielded the least amount of scatter between consecutive points. A weighted mean of all magnitudes was obtained from dithered frames within a given phase bin and the standard deviation of the robust mean was propagated to the photometric uncertainties.

Fig. 4. Example light curves of FU (top) and FO mode (bottom) Cepheids based on all the photometric data (in grey) from the dithered frame images. The light curves are shown in two optical (gi) and two NIR (JKs) filters. The black circles represent the weighted mean values in a given phase bin and the errors represent their standard deviations. The range of y-axis is the same in each band for a given Cepheid. The Gaia source ID, pulsation mode, and the periods are listed at the top of each panel.

Figure 4 also shows the phase-binned light curves of Cepheids overplotted on the original data points. While the variability was quite evident in the original light curves, it is more pronounced in the phase binned light curves. The phase-binning procedure was particularly useful for low-amplitude FO variables for which the larger uncertainties on photometry from individual dithered frames often dominated variability amplitudes in NIR bands. The light curves constructed from the phase binning procedure are used for further analysis in this work. These light curves were visually inspected to assign a quality flag of A,B, and C to the best, good, and poor light curves depending on the scatter and the phase coverage of the light curves. There are 52/22 stars with best and good (‘A and B’) quality flags, and only four stars have poor light curves. We note that the entire sample of Cepheids was used in the subsequent analysis irrespective of their light curve quality flags.

3.2. Light curve fitting and mean magnitudes

The light curve templates for Cepheids are not available in optical Sloan (griz) filters, but are available in NIR bands (e.g. Inno et al. 2015). Since our light curves are phase-binned, we do not expect to recover small features such as bump near the maximum light, particularly at NIR wavelengths where light curves are more sinusoidal. Therefore, we decided to analyse multiband light curves with the Fourier decomposition method (e.g. Bhardwaj et al. 2015, 2017) and not rely on a particular set of templates. The phase-binned light curves were fitted with a Fourier sine series fit in the following form:

$$\begin{array}{c}\hfill m_{\lambda }=m_{0,\lambda }+{\displaystyle \sum _{k=1}^N}{A}_{k,\lambda }sin(2k\pi ·{\mathrm{\Phi }}_{\lambda }+{\varphi }_{k,\lambda }),\end{array}$$

where m_λ is the magnitude as a function of phase (ϕ_λ) at a given wavelength (λ). The mean magnitudes (m₀) and Fourier coefficients (A_k, ϕ_k) were determined using a Fourier sine series with varying order of fit (k = 1–5). A minimum order of k = 1 was chosen to obtain purely sinusoidal fits, which are ideal for FO Cepheids, particularly at NIR wavelengths. Since the phase-binned light curves do not necessarily provide enough data points around their extrema, we find that a fifth order Fourier series was sufficient to fit the main saw-tooth characteristics and Hertzsprung progression (Hertzsprung 1926) features. Furthermore, limiting higher order fit to k = 5 also avoided possible numerical ringing in case of larger phase gaps. The best order of fit was determined using the Baart’s criteria (Baart 1982).

Figure 5 displays the Fourier-fitted light curves in grizJHK_s filters for FU and FO Cepheids with varying periods. The saw-tooth light curve shape for a short-period (P = 3.86 days) FU Cepheid and flatter maxima for Cepheids with periods close to 10 days can be noted despite the phase-binning. Most Cepheids with periods longer than 20 days typically have similar near-sinusoidal light curve shapes with decreasing amplitude as a function of period. In the case of FO Cepheids, sinusoidal Fourier series with k = 1 fitted well the majority of stars at all wavelengths. The best-fitting Fourier fits were used to derive intensity-averaged magnitude and peak-to-peak amplitude for all Cepheids except for four stars in g band and one star in K_s band. These Cepheids with lack of measurements in one filter did not have enough light curve phase points for a Fourier-fit or mean-magnitude determinations. The photometric mean magnitudes of all stars are listed in Table 1.

Fig. 5. Representative phase-binned light curves of FU and FO Cepheids with varying periods at all wavelengths (grizJHKs). The magnitudes in g, r, J, H bands were offset by −0.5, −0.2, −0.1, +0.1 for visualization purposes. The best-fitting Fourier series fits are also shown as dashed lines. The Gaia Source ID, pulsation mode, and period are listed at the top of each panel. The light curve quality flag is given at the bottom left of each panel. The uncertainties in magnitudes at a given phase point also include the scatter at that phase in the original light curves.

3.3. Period-amplitude diagrams

Figure 6 displays period-amplitude diagrams for Cepheids at multiple wavelengths. For reference, the Gaia G band amplitudes, if available, are shown in the first panel. The predicted G band amplitudes for Cepheid models representative of MW variables with Z = 0.02, Y = 0.28 are also shown for different masses (M = 4 − 11 M_⊙). These models were computed adopting canonical mass-luminosity relation and a fixed mixing length parameter (α = 1.5) for both FU and FO mode Cepheids (see De Somma et al. 2022, for details). In general, there is a good agreement between predicted and observed amplitudes for FU mode Cepheids with a maxima around log P ∼ 1.3 days. The short-period FO Cepheid models exhibit larger amplitudes than observations (De Somma et al. 2022). We note that the predicted amplitudes are sensitive to the adopted composition, mass-luminosity relations, and the convective efficiency in the pulsation models.

Fig. 6. Period-amplitude diagrams for Cepheids in multiple bands. The filled circles and squares represent FU and FO Cepheids, respectively. In the case of G band, the theoretically predicted amplitudes for metallicities (Z = 0.02) representative of Cepheids in the MW are also shown for different masses (De Somma et al. 2022). The solid and dashed lines represent FU Cepheid models, while asterisks show FO models. In grizJHKs bands, overplotted Cepheids in open star symbols are discussed in Sect. 3.3.

In grizJHK_s period-amplitude diagrams, two low-amplitude FU Cepheids (open stars) with log P ∼ 1.4 days can be seen as outliers in Fig. 6. These variables were classified as FU Cepheids in the ASAS-SN survey, but their light curves show large scatter not typical of Cepheids. Given their small V-band amplitudes of less than 0.3 mag, no obvious periodic variations were seen in our photometry despite a good phase coverage. The amplitudes of phase-binned grizJHK_s light curves are smaller than 0.1 mag, and these were not classified as variables in Gaia as well. We keep their classification as FU Cepheids since their mean magnitudes are well determined, and are used to probe their position on the PL relations. In NIR bands, four Cepheids with 1.1 < log(P) < 1.4 days exhibit relatively high amplitudes. These Cepheids have good-quality light curves and their independent G-band amplitudes also range between 0.78 and 1.12 mag.

In general, the decrease in amplitudes from the shortest-period to 10 days, an increase in amplitude with periods up to 20 days, and a reverse trend for the longer period Cepheids is also seen, which is typical for FU Cepheids (Bhardwaj et al. 2015). The separation of FU and FO Cepheids on the period-amplitude plane is noted for all wavelengths, but is most distinct in K_s-band. However, two FU Cepheids (open stars) with P < 10 days and K_s-band amplitudes of less than 0.1 mag are also populating FO amplitude cluster in all other bands. These Cepheids (5718760258978008320, 2936274153063501824) were classified by Gaia as FU mode and have periods of 4.70 and 6.39 days, and G-band amplitudes of 0.35 and 0.20 mag, respectively.

4. Absolute calibration of the Leavitt law

4.1. Gaia parallaxes

The Gaia astrometric parallaxes were obtained for all Cepheids from Gaia Collaboration (2023b) to calibrate multiband PL and PW relations. A few of these Cepheids are in the Galactic disk and anti-centre direction, and therefore the parallaxes of these distant targets have large uncertainties. The median error on parallaxes is 9.4%, but 15 Cepheids have parallax uncertainties of > 25%. Two Cepheids in our sample have negative parallaxes, and therefore the subsequent analysis is carried out in the parallax space. Several studies have found that the Gaia parallaxes are systematically smaller resulting in larger distances (Arenou et al. 2018; Groenewegen 2018; Riess et al. 2021; Bhardwaj et al. 2021; Molinaro et al. 2023). Lindegren et al. (2021) provided a parallax zero-point correction recipe that derives a shift for an individual star based on its magnitude, colour, and ecliptic latitude. The parallax correction varies between −0.075 and 0.023 mas for Cepheids in our sample. However, these parallax corrections are now known to over-correct parallaxes resulting in smaller distances (Riess et al. 2021; Cruz Reyes & Anderson 2023; Molinaro et al. 2023). We adopted Gaia parallaxes after applying the suggested corrections by Lindegren et al. (2021) and include an over-correction offset of 0.014 mas from Riess et al. (2021). The parallax uncertainties were increased by 10% following Riess et al. (2021) and include parallax correction uncertainty of 0.005 mas added in quadrature (Lindegren et al. 2021).

In addition to the uncertainties on parallaxes, we also looked at Gaia astrometric quality flags: renormalized unit weight error (RUWE) and goodness of fit (GOF). There are three stars with RUWE > 1.4 and one of these also has a GOF = 13.98 exceeding a threshold of 12.5 adopted by Riess et al. (2021). The largest RUWE among these three Cepheids is 1.55, and their parallax uncertainty are < 15%. Moreover, six Cepheids also fall near the sharp inflation point in the Lindegren et al. (2021) formulae at G = 11 mag, and their parallax corrections may have larger uncertainties. Given that our sample size is modest, we do not exclude these Cepheids from our analysis but treat them cautiously when deriving PL relations.

4.2. Reddening and extinction corrections

Groenewegen (2018) provided reddening values for more than 450 Cepheids that were primarily taken from Fernie et al. (1995). However, there are only two Cepheids (X Sct and V5567 Sgr) in our sample in common with Groenewegen (2018). Ripepi et al. (2021) derived new period-colour relations to estimate their intrinsic colours, and therefore reddening values. We used Eq. (4) of Ripepi et al. (2021) to derive intrinsic (V − I)₀ colour for all Cepheids. For apparent (V − I) colours, we first used photometric transformations from Pancino et al. (2022) to derive V and I band magnitudes using Gaia photometric data. These photometric transformations based on homogeneous Gaia data have been used to provide accurate magnitudes and colours in Johnson-Kron-Cousins photometric systems (Trentin et al. 2023) and were preferred over the heterogeneous literature compilations. The apparent V and I magnitudes were obtained for Cepheids which have (BP − RP) colours within the recommended range for photometric transformations (Pancino et al. 2022). The colour-excess values in (V − I) were converted using E(V − I) = 1.28E(B − V) (Tammann et al. 2003). Therefore, reddening E(B − V) values were obtained for 76 Cepheids, which vary between 0.092 and 2.186 mag. The uncertainty in the empirical equation used to obtain intrinsic colours of Cepheids were propagated to the errors in the reddening values.

Recently, Breuval et al. (2022) used Bayestar19³ three-dimensional reddening maps from Green et al. (2019), and period-colour relations from Riess et al. (2022b) to derive reddening values for 222 Cepheids. There is only one star in common with Breuval et al. (2022) sample. We also obtained E(B − V)_G19 reddening values for 70 Cepheids from Bayestar19 maps which varied between 0.145 and 4.017 mag. A comparison of the two sets of reddening values for 67 Cepheids suggests a median difference of E(B − V)_G19 − E(B − V) = − 0.02 mag with a scatter of 0.10 mag after excluding six stars with the largest differences. If we scale the Bayestar19 reddening values by a scale factor of 0.884 as recommended by Green et al. (2019), the median difference increases to −0.11 mag. Recently, Narloch et al. (2023) found that the Bayestar19 reddening maps are inadequate for Cepheids in their sample because the reddening values resulted in unexpectedly large scatter in the resulting PL relations. Therefore, we adopted the reddening values based on period-colour relations.

To correct the magnitudes for extinction, we adopted the Fitzpatrick (1999) reddening law assuming an R_V = 3.1. The total-to-selective absorption ratios were estimated using the dust_extinction⁴ python package. The effective central wavelength corresponding to each filter was adopted from the Spanish Virtual Observatory filter profile service⁵. The absorption ratios in different filters are provided in Table 2. These absorption ratios were used together with the E(B − V) values to apply extinction corrections to the mean-magnitudes at all wavelengths.

4.3. Period–luminosity relations

The parallax uncertainties for Cepheids in our sample are large (> 25%) for 17 stars and two stars have negative parallaxes. Moreover, 12 Cepheids have distances larger than 10 kpc. We did not adopt distances from Bailer-Jones et al. (2021) for which priors based on a three-dimensional model of the Galaxy become dominant for these distant targets. Therefore, instead of converting apparent magnitudes to absolute magnitudes using geometric distances, we decided to work in the parallax space to derive PL relations. We derived astrometry-based luminosity (ABL, Feast & Catchpole 1997; Arenou & Luri 1999) defined as:

$$\begin{array}{cc}& \text{ABL}={\overline{\omega }}_{\text{(mas)}}{10}^{0.2m_{\lambda }-2}\hfill \\ \hfill & ={10}^{0.2({\alpha }_{\lambda }+{\beta }_{\lambda }(log{P}_i-log{P}_0))},\hfill \end{array}$$(1)

where m_λ is the extinction corrected magnitude at a given wavelength for a Cepheid with period (P_i). The period (P₀) at which the zero-point is determined, is adopted at log P₀ = 1.0 days for FU Cepheids and log P₀ = 0.4 days for FO Cepheids. Since our sample size of each subtype is small, we also considered a combined sample of Cepheids by fundamentalizing the periods of FO Cepheids using the equation: P_FU = P_FO/(0.716 − 0.027log P_FO) from Feast & Catchpole (1997). The absolute zero-point was also obtained at log P₀ = 1.0 days for this combined sample. While fitting these relations in the form of Eq. (1), we iteratively removed the single largest outlier in each iteration until all residuals are within ±3σ, where σ represents the root-mean square (rms) error. We created 10⁴ random realizations of ABL fits to derive coefficients and their associated uncertainties. In this procedure, we also included the parallax zero-point offset error of 6 μas (Riess et al. 2021).

Figure 7 displays the ABL of Cepheids in multiple bands and the results of the best-fitting relation in the form of Eq. (1) are listed in Table 3. A clear decrease in the scatter and PL residuals is seen moving from optical to NIR wavelengths. Among the Cepheids that have metallicity measurements, the residuals are larger for metal-poor stars ([Fe/H] < −0.3 dex). This suggests that the metal-poor stars have higher ABL and fainter absolute magnitudes for a given period. This is consistent with recent results of a negative metallicity coefficient of PL relations in the literature (e.g. Breuval et al. 2022; Molinaro et al. 2023; Bhardwaj et al. 2023). However, the metal-poor stars, most of which are more distant, also have larger parallax uncertainties. We note that the PL relation in the g band is not shown in Fig. 7, but exhibits trends similar to other optical filters. The results given in Table 3 show that the slopes of PL relations are steeper at NIR wavelengths than optical bands, a trend that is typical for both Cepheid and RR Lyrae stars (Bhardwaj 2020).

Fig. 7. Astrometry-based luminosity of FU (circles) and FO (squares) Cepheids as a function of pulsation period with magnitudes in rizJHKs bands is shown in the top part of each subpanel. The dashed and dotted lines represent the best-fitting relation for FU and FO mode Cepheids, respectively. The middle plot in each subpanel displays the residuals of the best fit as a function of period. The dashed lines show ±3σ scatter for FU and FO modes Cepheids. The bottom plot in each subpanel shows the residuals as a function of metallicity, if available. The dashed lines in the bottom two subpanels represent zero-residual variation. The partially filled symbols represent outliers, while empty symbols represent no metallicity measurement. The colour bar represents metallicity.

4.4. Period-Wesenheit relations

The lack of accurate and independent reddening values of Cepheids in our sample can impact their PL relations. The Wesenheit magnitudes, which include a colour-term, are constructed to be reddening independent (Madore 1982). Given a reddening law, the coefficient of the colour-term of a two-band Wesenheit magnitude is derived from the total-to-selective absorption parameters at those wavelengths. We adopted the Fitzpatrick (1999) reddening law and derived PW relations for different combination of bandpasses as listed in Table 2. Figure 8 displays ABL for the PW relations for Cepheids in our sample. These relations are apparently tighter than the PL relations suggesting that the reddening values based on period-colour relations are not very accurate. The trend of larger residuals for metal-poor stars can also be seen in the bottom panels similar to PL relations. Table 4 lists the coefficients of PW relations.

Fig. 8. Same as Fig. 7, but for Wesenheit magnitudes.

We tested the PL and PW relations listed in Tables 3 and 4 for different samples after excluding light curves with poor quality flags, Cepheids with RUWE > 1.4 and GOF > 12.5, and applying different outlier rejection thresholds. No statistically significant variations are seen in the coefficients of the PL/PW relations. If we restrict the sample to parallax uncertainties smaller than 15% and/or exclude stars beyond 10 kpc, the coefficients of PL relations vary within 1σ of their quoted uncertainties, which increase with lower statistics. While using reddening values from the maps of Green et al. (2019), the PL relations exhibit significantly larger scatter, as noted by Narloch et al. (2023). The dominant source of scatter in our PL and PW relations is due to large parallax uncertainties, which are expected to improve in the future Gaia data releases.

4.5. Comparison with literature

Figure 9 displays a comparison of slopes of PL relations in optical and NIR filters with previously reported values in the recent literature. The typical trend of steeper slopes at NIR wavelengths is seen. However, we did not find a distinct linear relation between the slopes and the wavelength across all filters, as noted by Trentin et al. (2024). The slopes of ri band PL relations for the combined sample of Cepheids are in excellent agreement with the slopes of empirical PL relations for FU Cepheids derived by Narloch et al. (2023). The slope of g band PL relation in this work is steeper than that derived by Narloch et al. (2023), but the difference is still within 2σ of their quoted uncertainties. Di Criscienzo et al. (2013) provided theoretical PL relations in Sloan filters based on non-linear pulsation models representative of Cepheids in our Galaxy. The authors found their predicted slopes to be mildly steeper than the empirical slopes of PL relations derived using photometric transformation from traditional Johnson-Cousins to Sloan bands. The theoretical griz band slopes in Di Criscienzo et al. (2013) are also steeper for shorter period variables (log P < 1.0 days) than those for longer period stars. The slopes of riz band PL relations in this work are in good agreement with the theoretically predicted slopes for the entire period range models, but are significantly steeper in the g-band. The slopes of our NIR PL relations for Cepheids are in good agreement with the results of Breuval et al. (2022), Trentin et al. (2024), and Bhardwaj et al. (2023). The slopes of G, BP and RP band PL relations are in agreement with those from Breuval et al. (2022). In the case of optical PW relations, the adopted filter combinations and/or the color-coefficients are different in this work as compared to Narloch et al. (2023). Nevertheless, the slopes of W_gr PW relations are in good agreement when we consider the combined sample of Cepheids. The slope of W_G P-W relation is similar to the one derived by Breuval et al. (2022), while it is consistent within 2σ for the W_JKs PW relation.

Fig. 9. Comparison of slopes of PL relations at different wavelengths. The slopes of the Sloan band PL relations are compared with empirical (gri) and theoretical (griz) PL relations from Narloch et al. (2023) and Di Criscienzo et al. (2013), respectively. The slopes of the G, BP, RP, V, I, J, H, and Ks band PL relations are from Breuval et al. (2022), Trentin et al. (2024), and Bhardwaj et al. (2023). The shaded regions display the ±1σ standard deviation around the mean value of the slopes for the wavelengths under consideration.

5. Metallicity dependence of the Leavitt law

The influence of metallicity on the absolute magnitudes of Cepheid variables needs to be properly quantified because it is a crucial parameter in improving the overall fit of the cosmic distance ladder used to measure the Hubble constant (Riess et al. 2022a). Most of the recent studies have measured a negative sign for the metallicity coefficient of PL relation (Gieren et al. 2018; Riess et al. 2021; Ripepi et al. 2021; Breuval et al. 2022; Bhardwaj et al. 2023). While our sample of stars with spectroscopic metallicities is rather small, it is based on a homogeneously collected high-resolution spectroscopic dataset within the framework of the C-MetaLL survey, thus minimizing systematics in combining different literature [Fe/H] measurements.

The residuals of PL and PW relations in Figs. 7 and 8 show a clear trend as a function of metallicity such that metal-poor stars exhibit larger residuals. Figure 2 displays metallicity distribution of 65 stars with available [Fe/H] values for our sample of Cepheids. The most metal-poor star ([Fe/H] = −1.66 dex) in our sample comes from Gaia and has poor RVS quality flags. There are four more stars with similar quality flags and an assigned error of 0.5 dex (Trentin et al. 2024), which were excluded from this analysis. Therefore, our sample consists of 61 stars including 33 FU and 28 FO Cepheids. The [Fe/H] values of these Cepheids are between −1.11 and 0.6 dex. Given small samples of both subtypes, we only derive PLZ and PWZ relations for the combined sample of Cepheid variables. A metallicity term is added to Eq. (1) in the following form:

$$\begin{array}{cc}& \text{ABL}={\overline{\omega }}_{\text{(mas)}}{10}^{0.2m_{\lambda }-2}\hfill \\ \hfill & ={10}^{0.2({\alpha }_{\lambda }+{\beta }_{\lambda }(log{P}_i-log{P}_0)+{\gamma }_{\lambda }\text{[Fe/H]})}.\hfill \end{array}$$(2)

In order to better constrain the metallicity coefficients, we assume that the slope of the PL relation does not change with metallicity. This is a basic assumption that is applied when measuring extragalactic distances using Cepheid variables. However, the metallicity can also affect the slope of the PL relations at different wavelengths (e.g. Ripepi et al. 2021; Trentin et al. 2024). Given a small sample size, we assume a fixed slope of PL/PW relation from Tables 3 and 4 in Eq. (2) at a given wavelength. Therefore, we only solve for the absolute zero-point and the metallicity coefficient as free parameters.

Table 5 lists the coefficients of PLZ and PWZ relations. The metallicity coefficient varies significantly between −0.62 ± 0.14 mag dex⁻¹ in G band and −0.30 ± 0.11 mag dex⁻¹ in K_s band. These metallicity coefficients are systematically larger than those determined by Breuval et al. (2022) and Bhardwaj et al. (2023), but are more in agreement with previous results from the C-MetaLL survey (Ripepi et al. 2021; Trentin et al. 2024). We note that the metallicity range of MW Cepheid sample was rather small in Breuval et al. (2022) and Bhardwaj et al. (2023) with a mean value of 0.09 dex (dispersion of 0.12 dex) and 0.04 dex (dispersion of 0.09 dex), respectively. In contrast, the mean metallicity of the sample used in this work is −0.18 dex (median value of −0.11 dex) and a three times larger dispersion of 0.37 dex than Breuval et al. (2022) sample. If we restrict the sample to 38 stars with [Fe/H] values larger than −0.3 dex, the metallicity coefficients decrease (in absolute sense) considerably but exhibit larger uncertainties. For example, the metallicity coefficient in r-band becomes −0.20 ± 0.24 mag dex⁻¹ from a value of −0.46 ± 0.14 mag dex⁻¹ quoted in Table 5. These coefficients are almost zero in the PLZ/PWZ relations involving the K_s band when most metal-poor stars ([Fe/H] < −0.3 dex) are excluded. While the metallicity coefficients listed in Table 5 are weakly constrained due to lower statistics, it is important to probe this dependence more in detail. These large metallicity coefficients of the order of −0.5 mag dex⁻¹ can lead to significant biases of 0.05 mag in distance modulus determinations if the mean metallicities of the calibrator and target Cepheid PL relations differ by even Δ[Fe/H] = 0.1 dex.

5.1. Comparison of residuals of PL and PLZ relations

Figure 10 displays the residuals of PL/PW and PLZ/PWZ relations as a function of metallicity. It is clearly evident that the dominant contribution to the metallicity term comes from the metal-poor stars ([Fe/H] < −0.3 dex) which also have larger parallax uncertainties. The median parallax uncertainties for stars with [Fe/H] > −0.3 dex is 6.4% while it increases to 16.7% for Cepheids with [Fe/H] < −0.3 dex. For example, the median residuals of r band PL relation are Δ(ABL) = 0.004 and Δ(ABL) = 0.048 for metallicities higher (metal-rich) and lower (metal-poor) than −0.3 dex, respectively. When the metallicity term is included, the median residual of metal-poor stars decreases to Δ(ABL) = 0.014, while it remains the same for metal-rich stars. Similarly, the median residuals of the ABL fits for W_JKs Wesenheit relations decrease from Δ(ABL) = 0.021 to Δ(ABL) = 0.012 for metal-poor stars. However, this is also expected since the median metallicities for metal-rich and metal-poor stars are 0.10 and −0.51 dex, respectively. Therefore, the contribution to average residuals due to metallicity coefficients is smaller for metal-rich stars than for metal-poor stars. Nevertheless, it is difficult to separate the contribution of metallicity effects and parallax uncertainties to the scatter in the residuals seen for metal-poor stars in Fig. 10. If we exclude stars with larger parallax uncertainties, there are not enough metal-poor Cepheids in the sample for a proper quantification of metallicity coefficient of PL/PW relations.

Fig. 10. Comparison of residuals of r, z, Ks, Wgi, Wri, WJKs PL/PW and PLZ/PWZ relations is shown in six different panels. In each panel, the top and bottom plot shows the variation of residuals as a function of period and metallicity, respectively. The bandpasses are mentioned at the top right of each panel. The open and filled circles represent residuals of the PL/PW and PLZ/PWZ relations, respectively. The colour bar represents metallicity.

5.2. Metallicity coefficient as a function of wavelength

We compared the metallicity coefficients of PLZ/PWZ relations with recent determinations in the literature. Figure 11 displays metallicity coefficient as a function of wavelengths. Breuval et al. (2022) fitted a linear regression between γ and 1/λ and concluded that the metallicity effect is uniform over a wide range of wavelength. We do not see a strong linear correlation between metallicity coefficient and wavelength. The PLZ relations at longer wavelengths seem to suggest a marginally smaller metallicity term than at shorter wavelengths, but the uncertainties on γ values are larger for our sample of stars. The metallicity coefficients derived in this work are consistent with earlier results of Ripepi et al. (2021) and Trentin et al. (2024) within the C-MetaLL survey. We note that the approach of using individual metallicities of MW Cepheids is the same in these studies. In contrast, the results of Gieren et al. (2018), Breuval et al. (2022), and Bhardwaj et al. (2023) are based on a comparison of intercepts of PL relations for Cepheids in the MW and the Magellanic Clouds. The average value of the metallicity coefficient (−0.28 ± 0.07 mag dex⁻¹) derived in these studies is systematically smaller than those listed in Table 5. In contrast, the mean value of the metallicity coefficient (−0.47 ± 0.08 mag dex⁻¹) derived from all the values in this work, and Ripepi et al. (2021), and Trentin et al. (2024) is significantly larger. However, these metallicity coefficients are also based on a sample of Cepheids covering a wide range of metallicities. Nevertheless, the metallicity coefficients of our PLZ/PWZ relations are in agreement with most of these measurements given their large uncertainties due to a small sample of Cepheid variables.

Fig. 11. Comparison of metallicity coefficients as a function of wavelength. The metallicity coefficients of the G, BP, RP, V, I, J, H, and Ks bands are taken from Gieren et al. (2018), Ripepi et al. (2021), Breuval et al. (2022), Trentin et al. (2024), Bhardwaj et al. (2023), but no previous measurements are available in the griz bands. The shaded regions display a ±1σ standard deviation around the mean value of the metallicity coefficients from (1) this work; Trentin et al. (2024), Ripepi et al. (2021) and (2) Bhardwaj et al. (2023), Breuval et al. (2022), Gieren et al. (2018).

5.3. Parallax correction and the uncertainty in the zero-point

As shown in the previous subsections, the sample presented in this work is not optimal for the accurate and precise calibration of PLZ and PWZ relations due to low statistics, large parallax uncertainties, and the lack of accurate reddening measurements. Nevertheless, we investigate the reliability of the large metallicity coefficients and the zero points of the calibrated PLZ and PWZ relations. When compared with the results of Bhardwaj et al. (2023), the zero-points of K_s band PLZ, and W_J, Ks/W_G PWZ relations are nearly ∼0.25 mag brighter in this work. We note that both studies adopted the same parallax zero-point offset correction of −14 μas (Riess et al. 2021), but Bhardwaj et al. (2023) sample is based on MW Cepheid standards that have accurate parallaxes, low reddening, and a limited metallicity range. Therefore, the readers interested in determining Cepheid-based distances are recommended to use the calibrated relations provided in Bhardwaj et al. (2023) for metallicities closer to solar value, and those in Trentin et al. (2024) for a wide-range of metallicities.

In the recent work within the C-MetaLL survey, Trentin et al. (2024) suggested that a larger parallax zero-point offset correction results in a larger distance measurement. The authors found that both the slope and the intercept of PLZ/PWZ relation increases by 1 − 4% and the metallicity coefficient decreases (in absolute sense) if a larger parallax offset correction is adopted. For the W_J, Ks Wesenheit, Trentin et al. (2024, Table 3) found the coefficients of PWZ relation: α = −6.09 ± 0.02 mag, β = −3.29 ± 0.03 mag dex⁻¹, γ = −0.45 ± 0.05 mag dex⁻¹, with no parallax zero-point offset correction. When compared with the results in Table 5 based on a parallax zero-point offset correction of −14 μas, indeed the slope and the intercept increase and the metallicity coefficient decreases (in absolute sense), in agreement with Trentin et al. (2024). However, the zero-point is significantly brighter by ∼0.25 mag, but this difference becomes smaller when no parallax zero-point offset correction is adopted. For this latter case, we find the coefficients of W_J, Ks PWZ relation: α = −6.25 ± 0.04 mag, β = −3.42 ± 0.07 mag dex⁻¹, γ = −0.43 ± 0.10 mag dex⁻¹. We note the slope is fixed in our analysis, and the metallicity coefficient is now in excellent agreement with Trentin et al. (2024). However, the zero-point is still significantly brighter by ∼0.16 mag, due to larger parallax uncertainties. In terms of the distance to the LMC, Trentin et al. (2024, Table 3) found a value of 18.41 ± 0.02 mag using W_J, Ks Wesenheit. For the LMC Cepheids, we adopt the zero-point of 12.46 ± 0.05 mag for W_J, Ks Wesenheit from Bhardwaj et al. (2016) based on the data from Macri et al. (2015), and a mean-metallicity of −0.41 ± 0.02 (Romaniello et al. 2022). We find a LMC distance of 18.53 ± 0.06 mag for the no parallax zero-point offset correction. However this value of LMC distance increases to 18.72 ± 0.06 mag for the calibrated PWZ relation listed in Table 5 due to a smaller metallicity term and a brighter absolute zero-point. With the sample presented in this work, it is not possible to separate the contribution of parallax uncertainties (and parallax zero-point offset correction) and metallicity effects on the absolute calibration of PLZ and PWZ relations.

6. Summary

We presented new homogeneously collected light curves of 78 MW Cepheid variables at optical Sloan (griz) and NIR (JHK_s) wavelengths. These observations were obtained simultaneously using the REM telescope and present the light curves in z band for FU Cepheids and in griz band for FO Cepheids for the first time. Multiband light curves of Cepheids were obtained primarily to complement the high-resolution spectroscopic metallicities within the framework of the C-MetaLL survey (Ripepi et al. 2021). The sample of 78 Cepheids includes 49 FU and 29 FO mode variables. The light curves of Cepheids were fitted with Fourier sine series to determine accurate mean magnitude and peak-to-peak amplitudes, which were used to investigate their pulsation properties. The period-amplitude diagrams were presented for the time in Sloan filters for both subtypes of Cepheids. The mean magnitudes were used to derive PL and PW relations for Cepheid variables using their astrometry-based luminosity at multiple wavelengths.

Cepheid variables in the present study are located at a wide range of distances between 0.5 and 19.7 kpc, and therefore their Gaia parallaxes exhibit varying uncertainties with larger errors for distant objects. Moreover, some of these Cepheids are located in the Galactic disk and anti-centre direction and are significantly reddened with colour-excess values exceeding 1 mag. There are no reddening values available in the literature for most of these Cepheids, and thus the colour-excess values were obtained by determining their intrinsic colours using the empirical period-colour relations for Cepheid variables (Tammann et al. 2003; Ripepi et al. 2021). The larger uncertainties in parallaxes and reddening values and a modest sample size limits the accuracy and precision of the PL relations for these variables. The reddening uncertainties can be mitigated by employing the Wesenheit magnitudes, which results in tighter PW relations, in particular, at shorter optical wavelengths where extinction is more severe.

In addition to parallax uncertainties, metallicity variations can also contribute to the scatter in the empirical PL and PW relations. Homogeneous high-resolution spectroscopic metallicities for 59 of 78 Cepheids have already been collected as part of the C-MetaLL survey. In addition, six Cepheids have medium-resolution spectra from Gaia-RVS (Recio-Blanco et al. 2023), two of which have reliable quality flags. We investigated the residuals of astrometry-based luminosity with PL/PW relation fits as a function of their metallicities, if available. The residuals exhibit a clear trend with metallicity such that metal-poor stars have higher astrometry-based luminosity and fainter absolute magnitudes. This trend becomes significant for [Fe/H] < −0.3 dex, where the parallax uncertainties are also larger. When deriving PLZ and PWZ relations, the metallicity coefficient of these relations varies between −0.30 ± 0.11 mag dex⁻¹ in K_s band to −0.55 ± 0.12 mag dex⁻¹ in z band. While these empirical metallicity coefficients are weakly constrained due to the small sample size, they are systematically larger than previous determinations at optical and NIR filters. We compared the residuals of PL/PW and PLZ/PWZ relations to further investigate the impact of including the metallicity term, and found that the metallicity contribution predominantly affects the residuals for the most metal-poor Cepheids ([Fe/H] < −0.3 dex). The metallicity effect becomes smaller if we exclude these metal-poor Cepheids from the sample, but the uncertainties on the metallicity coefficient increase due to low statistics. The small sample size prevents us from separating the contribution to the scatter in the PL/PW relations due to the metallicity term and parallax uncertainties. These large metallicity coefficients and bright zero-points should be treated cautiously, and a larger sample of metal-poor Cepheids with more accurate parallaxes is needed to confirm these relatively large metallicity coefficients of the PLZ/PWZ relations.

The aim of the ongoing C-MetaLL survey is to increase the sample of Cepheids with homogeneous photometric and spectroscopic data to a few hundred stars, which is ideal for the absolute calibration of the Leavitt law and a proper quantification of metallicity effects at multiple wavelengths. In addition, the extension of the photometric data presented in this paper will also be useful for several other scientific goals, for example investigating the light curve structure of these variables at multiple wavelengths. The simultaneous optical and NIR light curves will be useful for a quantitative comparison with the theoretically predicted light curves, thus providing strong constraints for the input parameters to the pulsation models. Multiwavelength data will also be used to constrain the reddening values for Cepheid variables in a future study. Moreover, light curves in the optical Sloan filters will serve as templates for the identification and classification of Cepheid variables in the Vera C. Rubin Observatory’s Legacy Survey of Space and Time.

Acknowledgments

We thank the anonymous referee for the quick and constructive referee report that helped improve the manuscript. This project has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No. 886298. This research was supported by the Munich Institute for Astro-, Particle and BioPhysics (MIAPbP) which is funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s Excellence Strategy – EXC-2094 – 390783311. We acknowledge funding from INAF GO-GTO grant 2023 “C-MetaLL – Cepheid metallicity in the Leavitt law” (P.I. V. Ripepi). This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

References

All Tables

All Figures

	Fig. 1. Spatial distribution of classical Cepheids used in this work in Galactic coordinates. The smaller symbol sizes represent distant Cepheids. Open circles represent stars with no metallicity measurements. The colour bar represents metallicity.
In the text

Fig. 2. The period and metallicity distributions, and period-Wesenheit relation for MW Cepheids. Top: period distribution of 80 candidate Cepheids in the sample studied in this work. The filled and lined histograms correspond respectively to FU and FO Cepheids. Middle: histogram of metallicities for 65 Cepheids, of which 59 have homogeneous high-resolution spectroscopic metallicities from the C-MetaLL survey (Ripepi et al. 2021; Trentin et al. 2023). Bottom: period-Wesenheit relation in Gaia filters for the sample of 80 stars. The circles and squares represent fundamental and first-overtone mode Cepheids, respectively. Two extreme outliers are marked with their Gaia DR3 source ID.

In the text

Fig. 3. Comparison of the calibrated mean magnitudes of Cepheids in this work with the magnitudes from APASS DR10 (in griz) and 2MASS (in JHKs). The dashed and dotted lines respectively represent zero and mean offset. The mean (and the standard deviation) offsets for the bright (Ks < 8.5 mag) and faint (Ks ≥ 8.5 mag) Cepheids are also mentioned at the top left and bottom right of each panel. The filled circles and squares represent FU and FO Cepheids, respectively.

In the text

Fig. 4. Example light curves of FU (top) and FO mode (bottom) Cepheids based on all the photometric data (in grey) from the dithered frame images. The light curves are shown in two optical (gi) and two NIR (JKs) filters. The black circles represent the weighted mean values in a given phase bin and the errors represent their standard deviations. The range of y-axis is the same in each band for a given Cepheid. The Gaia source ID, pulsation mode, and the periods are listed at the top of each panel.

In the text

Fig. 5. Representative phase-binned light curves of FU and FO Cepheids with varying periods at all wavelengths (grizJHKs). The magnitudes in g, r, J, H bands were offset by −0.5, −0.2, −0.1, +0.1 for visualization purposes. The best-fitting Fourier series fits are also shown as dashed lines. The Gaia Source ID, pulsation mode, and period are listed at the top of each panel. The light curve quality flag is given at the bottom left of each panel. The uncertainties in magnitudes at a given phase point also include the scatter at that phase in the original light curves.

In the text

Fig. 6. Period-amplitude diagrams for Cepheids in multiple bands. The filled circles and squares represent FU and FO Cepheids, respectively. In the case of G band, the theoretically predicted amplitudes for metallicities (Z = 0.02) representative of Cepheids in the MW are also shown for different masses (De Somma et al. 2022). The solid and dashed lines represent FU Cepheid models, while asterisks show FO models. In grizJHKs bands, overplotted Cepheids in open star symbols are discussed in Sect. 3.3.

In the text

Fig. 7. Astrometry-based luminosity of FU (circles) and FO (squares) Cepheids as a function of pulsation period with magnitudes in rizJHKs bands is shown in the top part of each subpanel. The dashed and dotted lines represent the best-fitting relation for FU and FO mode Cepheids, respectively. The middle plot in each subpanel displays the residuals of the best fit as a function of period. The dashed lines show ±3σ scatter for FU and FO modes Cepheids. The bottom plot in each subpanel shows the residuals as a function of metallicity, if available. The dashed lines in the bottom two subpanels represent zero-residual variation. The partially filled symbols represent outliers, while empty symbols represent no metallicity measurement. The colour bar represents metallicity.

In the text

	Fig. 8. Same as Fig. 7, but for Wesenheit magnitudes.
In the text

Fig. 9. Comparison of slopes of PL relations at different wavelengths. The slopes of the Sloan band PL relations are compared with empirical (gri) and theoretical (griz) PL relations from Narloch et al. (2023) and Di Criscienzo et al. (2013), respectively. The slopes of the G, BP, RP, V, I, J, H, and Ks band PL relations are from Breuval et al. (2022), Trentin et al. (2024), and Bhardwaj et al. (2023). The shaded regions display the ±1σ standard deviation around the mean value of the slopes for the wavelengths under consideration.

In the text

Fig. 10. Comparison of residuals of r, z, Ks, Wgi, Wri, WJKs PL/PW and PLZ/PWZ relations is shown in six different panels. In each panel, the top and bottom plot shows the variation of residuals as a function of period and metallicity, respectively. The bandpasses are mentioned at the top right of each panel. The open and filled circles represent residuals of the PL/PW and PLZ/PWZ relations, respectively. The colour bar represents metallicity.

In the text

Fig. 11. Comparison of metallicity coefficients as a function of wavelength. The metallicity coefficients of the G, BP, RP, V, I, J, H, and Ks bands are taken from Gieren et al. (2018), Ripepi et al. (2021), Breuval et al. (2022), Trentin et al. (2024), Bhardwaj et al. (2023), but no previous measurements are available in the griz bands. The shaded regions display a ±1σ standard deviation around the mean value of the metallicity coefficients from (1) this work; Trentin et al. (2024), Ripepi et al. (2021) and (2) Bhardwaj et al. (2023), Breuval et al. (2022), Gieren et al. (2018).

In the text