The Optical to Mid-Infrared Extinction Law Based on the APOGEE, Gaia DR2, Pan-STARRS1, SDSS, APASS, 2MASS and WISE Surveys

A precise interstellar dust extinction law is critically important to interpret observations. There are two indicators of extinction: the color excess ratio (CER) and the relative extinction. Compared to the CER, the wavelength-dependent relative extinction is more challenging to be determined. In this work, we combine spectroscopic, astrometric, and photometric data to derive high-precision CERs and relative extinction from optical to mid-infrared (IR) bands. A group of 61,111 red clump (RC) stars are selected as tracers by stellar parameters from APOGEE survey. The multi-band photometric data are collected from Gaia, APASS, SDSS, Pan-STARRS1, 2MASS, and WISE surveys. For the first time, we calibrate the curvature of CERs in determining CERs E(lambda-GRP)/E(GBP-GRP) from color excess -- color excess diagrams. Through elaborate uncertainties analysis, we conclude that the precision of our CERs is significantly improved (sigma<0.015). With parallaxes from Gaia DR2, we calculate the relative extinction A_GBP/A_GRP for 5,051 RC stars. By combining the CERs with the A_GBP/A_GRP, the optical -- mid-IR extinction A_lambda/A_GRP have been determined in total twenty-one bands. Given no bias toward any specific environment, our extinction law represents the average extinction law with the total-to-selective extinction ratio Rv=3.16+-0.15. Our observed extinction law supports an adjustment in parameters of CCM Rv=3.1 curve together with the near-IR power law index alpha=2.07+-0.03. The relative extinction values of HST and JWST near-IR bandpasses are predicted in 2.5% precision. As the observed reddening/extinction tracks are curved, the curvature correction needs to be considered when applying extinction correction.


The Optical Extinction
In ultraviolet (UV)/optical bands of λ < 0.9 µm, the wavelength-dependent extinction law, A λ /A V , is known to vary from one sightline to another. This variation is mainly caused by the change in dust size distributions in different environments. The various extinction curves can be approximated by a one-parameter family of curves characterized by the ratio of the total extinction to the selective extinction (Cardelli et al. 1989;hereafter CCM). Theoretically, the extinction produced by Rayleigh scattering of small grains would have a steep curve with R V ∼ 1.2, while the extinction produced by very large grains would have a flat curve with R V → ∞ (Draine 2003). Observationally, the R V value can be as small as R V ∼ 2 in some diffuse sightlines (Fitzpatrick 1999;Wang et al. 2017), or as large as R V ∼ 6 in dense molecular clouds (Mathis 1990;Fitzpatrick 1999). Sightlines toward the Galactic diffuse interstellar medium (ISM) have an average value of R V ≈ 3.1 (CCM; Draine 2003;Schlafly & Finkbeiner 2011).
Both the color excess ratio (CER) E(λ − λ 1 )/E(λ 2 − λ 1 ) and the relative extinction A λ /A λ1 are indicators of the extinction law. Based on observations and the intrinsic color indices of the targets, the color excess (CE) E(λ − λ 1 ), and the CER E(λ − λ 1 )/E(λ 2 − λ 1 ) can be derived. However, the calculation of the wavelength-dependent interstellar extinction law A λ /A λ1 is more challenging. It requires an independent determination of the extinction or the distance to the target.
As a result, many measurements in the literature used a fixed total-to-selective extinction ratio in the optical bands (e.g., R V , R I ), or a fixed relative extinction in the near-infrared (IR) bands (e.g., A J /A KS , A H /A KS ) to convert reddenings into extinction law. For example, Rieke & Lebofsky (1985) calculated the extinction law for sightlines toward the Galactic center by assuming R V = 3.09. Later, Nataf et al. (2013) derived the extinction law to the Galactic bulge using a series of fixed total-to-selective extinction ratio R I . CCM computed A V and R V by adopting a standard curve given by Rieke & Lebofsky (1985). Indebetouw et al. (2005) fitted the RC locus in near-IR colormagnitude diagrams to directly extract A H /A KS by assuming a smooth, homogeneous dust distribution representing invariant extinction law. However, this assumption is not generally applicable as pointed out by Zasowski et al. (2009). Gao et al. (2009) and Zasowski et al. (2009) adopted a fixed near-IR extinction value A H /A KS to convert CERs into extinction law A λ /A KS . The extinction law determined by this method is affected by the systematic uncertainty, which is introduced by the adopted total-to-selective extinction ratio or the near-IR relative extinction value.
The distance information provides an opportunity to independently determine the interstellar extinction, A λ /A λ1 . However, only some special lines of sight, such as the Galactic center, clusters, the accurate distances to the objects can be derived. Nishiyama et al. (2006Nishiyama et al. ( , 2009 derived the relative extinction toward the Galactic center by assuming that all RC stars are at the same distance. Fitzpatrick & Massa (2007) determined interstellar extinction curves by analyzing of 328 Galactic stars with known distances. Chen et al. (2018; hereafter Chen18) adopted classical Cepheid as a diagnostic tool to estimate the relative extinction toward the Galactic center. The uncertainties of extinction law in these measurements are dominated by the accuracy of distance.
As mentioned by Fritz et al. (2011), one possible explanation is that most of the near-IR extinction measurements in the last century were based on nearby stars (∼ 3 kpc) with low extinction (A V 5 mag). While, the near-IR extinction measured in this century based on high extinction sources which are located in the Galactic inner disk, bulge, and even the Galactic center. Another possible explanation is the prevailing systematic errors in determining the index α. When converting the near-IR CER (e.g., E(J − H)/E(J − K S )) to the relative extinction A J /A KS based on the empirical formula A λ ∝ λ −α , the selection of filter wavelengths (effective or isophotal) affects the value of α and A J /A KS . As detailedly discussed and summarized in Stead & Hoare (2009), Fritz et al. (2011, and Wang & Jiang (2014), the different choice of filter wavelengths can lead to a significant discrepancy (∼ 10%) in α and A J /A KS . More specifically, adopting the effective wavelength will lead to a larger α, comparing to the adoption of the isophotal wavelength. Stead & Hoare (2009) suggested using the effective wavelength instead of the isophotal wavelength which caused the small index value in the measurement prior to 2005. However, Fritz et al. (2011) found that using the effective wavelength as presented in Stead & Hoare (2009) would slightly overestimate the index α. Therefore, Wang & Jiang (2014) suggested using CERs to represent the near-IR extinction.
In fact, differences are existed in the near-IR CERs as well. Indebetouw et al. (2005) investigated the IR extinction for two regions in the Galactic plane with different environments and derived a constant near-IR CER E(J −H)/E(H − K S ) = 1.778 ± 0.156. This value is corroborated by later measurements: Nishiyama et al. (2006), 1.72;Wang & Jiang (2014), 1.78;and Xue et al. (2016), 1.87. However, some larger values are also reported, such as E(J −H)/E(H −K S ) = 2.08 ± 0.03 (Racca et al. 2002); 1.91 ± 0.01 (Naoi et al. 2007); 2.09 ± 0.13 (Nishiyama et al. 2009); and 1.943 ± 0.019 (Schlafly et al. 2016). To complicate matters, some studies argued that the near-IR CERs are varied. Naoi et al. (2006) investigated the extinction toward the ρ Oph cloud and Chamaeleon cloud. They found that E(J − H)/E(H − K S ) changes with increasing optical depth. Later, Zasowski et al. (2009) studied the IR relative extinction by RC stars for contiguous sightlines covering ∼ 150 • of the Galactic disk. They reported that the value of IR CER is a function of the angle from the Galactic center, and this variation trend is more obvious in mid-IR bands than that in near-IR bands. The E(J − H)/E(H − K S ) ranges from 1.95 to 2.18 around the average value 2.04 ± 0.06. However, it worths to point out that E(J − H)/E(H − K S ) does not vary much for | l |< 60 • . Recently, Wang & Jiang (2014) investigated the near-IR extinction law based on a sample of spectroscopic selected K-type giants. They reported that the near-IR CERs are universal from diffuse to dense interstellar clouds. In the measurements of CERs, there are a number of problems. One is that, the accuracy of CER depends on the purity of the sample and the accuracy of the intrinsic color index. Another is that, the slope error of the color excess -color excess diagram is usually underestimated, especially when the sample size or the extinction is small. To summarize, the values of near-IR CERs still have about 15% uncertainty and need to be further investigated.
In practice, the near-IR extinction value, such as A J /A KS , A J /E(J −K S ), is commonly used to correct the interstellar extinction. So, an independent and accurate measurement of the near-IR extinction is desired. An effort was made by Fritz et al. (2011), who used hydrogen emission lines to explore the 1 to 19 µm extinction toward the Galactic center via a distance-independent method. They derived the absolute extinction and the corresponding α = 2.11 ± 0.06. Another effort was made in recent work by Chen18, who used classical Cepheids in the direction of the Galactic center to derive the extinction law between 1 and 8 µm based on three different approaches. They suggested that the values of A J /A KS = 3.005 ± 0.031 ± 0.094, A H /A KS = 1.717 ± 0.010 ± 0.033, and α = 2.05 ± 0.07 can better describe the extinction in the inner Galactic plane. However, both works were limited to the Galactic center sightlines. In this paper, with the precise parallaxes from Gaia and accurate stellar parameters from APOGEE, we will study the optical to mid-IR extinction in large scale, not just particular sightlines.

This Work
As mentioned in the previous sections, the CERs and relative extinction results in the literature still have significant differences which can not be fully explained by the uncertainties. The precise photometric data (σ 0.01 mag) are helpful to investigate these differences. Schlafly et al. (2016) have derived accurate CERs in some optical bands by using data from the Pan-STARRS1 (PS1) survey. However, as many PS1 stars are too faint to have reliable photometric magnitudes in the 2MASS, WISE, and APOGEE surveys, the IR extinction has not been well constrained in that work. Besides, due to a lack of distances information to the target sources, they measured the reddening curve rather than extinction curve.
In this work, we try to solve these problems by combining spectroscopic, astrometric, and photometric data. The latest data release of the APOGEE survey DR14 provides us the opportunity to obtain a large and homogeneous sample of RC stars. As distance and extinction are usually degenerate, accurate distance information is critical to derive accurate extinction. The Gaia DR2 provides good trigonometric parallaxes for part of RC candidates. We gather photometric data from several of survey projects including APASS, SDSS, Pan-STARRS1, 2MASS, WISE, as well as the unprecedentedly accurate photometry data from Gaia DR2. Finally, we use RC stars to re-investigate the optical to mid-IR extinction law. Both the reddening curve and the extinction curve have been determined with high accuracy. More importantly, we discuss the potential errors in these results in detail.
The description of data sets and RC stars sample are presented in Section 2. The optical to mid-IR CERs and relative extinctions are determined in Section 3. In Section 4, we analyze the uncertainties of our extinction law in detail. We compare our extinction law with those of previous results in Section 5. The estimated Gaia extinction coefficient and predicted near-IR extinction values for bandpasses of the Hubble Space Telescope Wide Field Camera 3 (HST WFC3) and the James Webb Space Telescope (JWST) NIRCAM are also presented in Section 5. We summarize our principal conclusions in Section 6.

Data
We collect stellar parameters from APOGEE survey to construct a sample of RC candidates. We gather broad-band photometric data from the APASS, SDSS, Pan-STARRS1, 2MASS, and WISE surveys. In addition, distance and photometric information from Gaia DR2 catalog are extracted. By cross-matching of these catalogs, stellar parameter, distance and twenty-one bands photometric data from optical to IR are obtained for each star.

APOGEE
The APOGEE (Apache Point Observatory Galaxy Evolution Experiment) is a large-scale, near-IR stellar spectroscopic survey (Eisenstein et al. 2011). The high-resolution spectra (R ∼22500) provide detailed stellar atmospheric parameters (e.g., effective temperature T eff , surface gravity log g, metallicity [M/H]) and chemical abundances. The primary stellar targets of APOGEE are red giant branch (RGB) stars and RC stars in the bulge as well as faint stars (Abolfathi et al. 2018). The latest data release DR14 (Abolfathi et al. 2018) contains all data from SDSS-III (APOGEE-1) as well as two years of data from SDSS-IV (APOGEE-2). As all APOGEE data, from the beginning of APOGEE-1, were reduced using the latest data reduction pipeline, the parameters provided in DR14 are slightly different from the previous data release version.

Gaia
The Gaia mission (Gaia Collaboration et al. 2016) has released the Gaia DR2, in which more than a billion sources have trigonometric parallaxes, three bands photometry (G, G BP , G RP ) and proper motions (Gaia Collaboration et al. 2018). The G band covers the whole optical wavelength ranging from 330 to 1050 nm, while G BP band and G RP band cover the wavelength ranges of 330 -680 nm and 630 -1050 nm, respectively (Evans et al. 2018). The central wavelengths of G, G BP , G RP bands are 673, 532, 797 nm, respectively (Jordi et al. 2010). Concerning the astrometric content, for the sources with five-parameter astrometric data, the median uncertainty of the parallax is ∼ 0.04 mas for G < 14 mag sources, 0.1 mas at G = 17 mag, and 0.7 mas at G = 20 mag (Lindegren et al. 2018). Concerning the photometric content, the photometric calibrations can reach a precision as low as 2 mmag on individual measurements, the systematic effects are present at the 10 mmag level (Evans et al. 2018).

APASS
The American Association of Variable Star Observers (AAVSO) Photometric All-Sky Survey (APASS) is conducted in five filters: Landolt B, V , and Sloan g ′ , r ′ , i ′ , probing stars with V band magnitude range from 7 to 17 mag . The latest DR9 catalog covers about 99% of the sky (Henden et al. 2016). Munari et al. (2014) investigated the accuracy of APASS data and confirmed that the APASS photometry did not show any offsets or trends. As we also collect SDSS photometric data, we only adopt the B and V data from the APASS DR9 catalog.

SDSS
The Sloan Digital Sky Survey (SDSS) is both imaging and spectroscopic survey (York et al. 2000). The imaging was performed simultaneously in bandpasses u, g, r, i, and z with central wavelengths of about 370, 470, 620, 750 and 890 nm, respectively (Fukugita et al. 1996;Gunn et al. 1998). We take the photometric data from the latest data release DR14, which is the second data release of the fourth phase of the SDSS (Abolfathi et al. 2018).

2MASS
The Two Micron All Sky Survey (2MASS) is a near-IR whole-sky survey (Skrutskie et al. 2006). 2MASS point-source catalog contains photometric measurements in the J, H, K S bands with isophotal wavelengths at 1.24, 1.66, 2.16 µm, respectively. As the APOGEE objects are selected from the 2MASS, the APOGEE catalog already includes the J, H, K S measurements.

WISE
The Wide-field Infrared Survey Explorer (WISE) survey is a mid-IR full-sky survey undertaken in four bands: W 1, W 2, W 3, and W 4 bands with wavelengths center at 3.35, 4.60, 11.56, and 22.09 µm, respectively (Wright et al. 2010). The WISE photometric data are taken from the AllWISE source catalog. Since few sources in our RC sample has reliable W 4 magnitudes, we only use W 1, W 2, and W 3 data.

The Red Clump Sample
The RC stars are a group of evolved stars in the core-helium burning stage. They cover the range of spectral types G8III -K2III with effective temperature 4500 K -5300 K (Girardi 2016). As the luminosities of RC stars are fairly independent of stellar composition and age, they are standard candles and widely used to estimate distances in the Galaxy and the local group. These stars appear as a narrow strip in the color -magnitude diagram (CMD), or a clumping group in the effective temperature (T eff ) -surface gravity (log g) diagram. Therefore, they can be easily selected with photometric or spectroscopic data, and become a useful probe to study the interstellar extinction (Indebetouw et al. 2005, Gao et al. 2009, Wang et al. 2017.
On the basis of the available stellar parameters from the APOGEE DR14 survey, we try to construct a homogeneous RC sample with high purity following these steps. First, we only include the sources with spectroscopic quality S/N > 50. Besides, we limit the metallicity [M/H] > −0.5 dex to reduce the potential effects on RC absolute magnitude.
Next, RC candidates are selected based on their clumping in the T eff -log g diagram with 4550 K ≤ T eff ≤ 5120 K and 2.2 ≤ log g ≤ 2.8. After these selections, our RC sample contains 61,111 sources. 97% of these RCs have K S band magnitude in the range of 7 -12.5 mag (distance less than 6 kpc). Note that there are some contaminations, such as secondary RC (SRC) stars and RGB stars, in this RC sample. We do not remove them in determining CERs (Section 3.2), while we remove them in calculating relative extinction (Section 3.4). By cross-matching this RC sample with photometric catalogs listed in Section 2.1, the multi-band photometric data for RC stars are obtained. To guarantee the photometric precision, we select stars that satisfy the following criteria for each photometric catalog.
1. For Gaia data, we select stars with photometric error ≤ 0.01 mag and magnitude ≤ 18.0 mag in G, G BP , G RP bands 1 .
2. For APASS data, we select stars with photometric error ≤ 0.05 mag in B, V bands.
3. For SDSS data, we select stars with photometric error ≤ 0.03 mag and magnitude ≥ 14.0 mag in u SDSS , g SDSS , r SDSS , i SDSS , z SDSS bands. To remove saturated stars, we adopt the criteria |SDSS/bandpass magnitude -PS1/bandpass magnitude|≤ 0.5 mag in SDSS/g, r, i, z bands.
5. For 2MASS data, we select stars with photometric error ≤ 0.03 mag and magnitude ranging from 6.0 to 14.0 mag in J, H, K S bands.
6. For WISE data, we select stars with photometric error ≤ 0.03 mag in W 1, W 2, W 3 bands.

THE OPTICAL TO MID-IR EXTINCTION VALUES
In this section, we calculate the two indicators of the wavelength-dependent extinction law: the CERs and the relative extinction. First, we adopt the color -excess method to determine the CERs E(λ − G RP )/E(G BP − G RP ). Then, we derive the relative extinction value A GBP /A GRP by two methods. Moreover, combining the CERs with the A GBP /A GRP , the optical to mid-IR bands relative extinction values A λ /A GRP are determined.

Method
We treat the sample stars as a whole to obtain the extinction by the color -excess method. The reason why we do not calculate the extinction of the individual RC star will be discussed in Section 5.1. Briefly, this method computes the ratio k λ of two CEs E(λ − λ 1 ) and E(λ 2 − λ 1 ), which can be expressed as where A λ is the extinction in the λ band of interest, A λ1 and A λ2 are extinction in the reference λ 1 band and the comparison λ 2 band, respectively. Therefore, the relative extinction A λ /A λ1 can be derived by This method is widely applied to a group of stars with homogeneous intrinsic color indices, such as RGB stars and RC stars (Indebetouw et al. 2005;Flaherty et al. 2007;Gao et al. 2009;Zasowski et al. 2009;Wang et al. 2013;Xue et al. 2016). As seen in Equation (2), the calculation of A λ /A λ1 requires the knowledge of A λ2 /A λ1 . The near-IR extinction values A J /A KS and A H /A KS are usually used to covert the CERs into the relative extinction A λ /A KS (Section 1.1). For example, the relative extinction can be described as A λ /A KS = 1 + k λ (A J /A KS − 1), where J and K S bands are treated as the comparison λ 2 band and the reference λ 1 band. However, as discussed in Section 1.2, the A J /A KS value has 20% uncertainty under a couple of reasons. Compared to the photometric accuracy of 2MASS bands, the photometry of Gaia bands are at least a factor of 3 more precise. Therefore, we take G RP as the reference λ 1 band and G BP as the comparison λ 2 band to reduce the uncertainties in the extinction determinations. The analysis of CER uncertainties caused by adopting different bases bands will be discussed in Section 4.2.
To summarize, we calculate the CER E(λ − G RP )/E(G BP − G RP ) and the A GBP /A GRP , respectively. Then, the corresponding relative extinction A λ /A GRP are determined by 3.2. Color Excess Ratios For a group of stars, the CER is the slope k λ of a linear fit to the color excess -color excess diagram. The CE is the difference between the observed color index and the intrinsic color index E(λ 2 − λ 1 ) = (λ 2 − λ 1 ) obs − (λ 2 − λ 1 ) int . The observed color index (λ 2 − λ 1 ) obs can be easily obtained from photometric data, while the knowledge of intrinsic color index (λ 2 − λ 1 ) int needs the information of the spectral type (effective temperature) or absolute magnitude. Based on the stellar parameters from the APOGEE catalog, we first determine the T eff -intrinsic color index relations via the method adopted by Wang & Jiang (2014). The idea of this method is to consider the top 5 percentile bluest stars at a given T eff as the zero reddening star. Hence, the observed color index of these bluest stars can represent the intrinsic color index (λ − G RP ) int at the given T eff . Then, a polynomial fitting is adopted to determine the T eff -(λ − G RP ) int relations. After subtracting the intrinsic colors, CEs are determined. Figure 1 illustrates the linear fit to the color excess -color excess diagram in four bands, APASS/B, SDSS/u SDSS , PS1/g PS1 , and 2MASS/K S , respectively. The color shows the number density of RC stars, the black lines are the best fits. The distributions of the residuals of the fits, ∆E(λ − G RP ), are displayed as well.
At first glance, the color excess -color excess distribution exhibits good linearity, especially in the high-precision g PS1 band. The dispersions of residuals in B and u SDSS bands are larger than those in g PS1 and K S bands. Moreover, the triangular distribution of residuals are obvious in the low CE part of B and u SDSS bands. To reduce the residual, we further analyze the relation among the intrinsic color and metallicity or surface gravity. The intrinsic color indices are estimated by stellar parameters, including T eff , [M/H], and log g, which can be expressed as: The corresponding coefficients for each color index (λ − G RP ) int are listed in Table 1. We found that the intrinsic colors of the optical bands, such as B, u SDSS , g SDSS , g PS1 , are correlated with metallicity. The intrinsic color index (u SDSS − G RP ) int is even correlated with surface gravity. While in the other bands, the dependence on [M/H] and log g are moderate or negligible. The CEs are then determined by subtracting the intrinsic color indices (Table 1) from the observed color indices.
In addition, as seen in Figure 1, the CER in g PS1 band exhibits linear correlations with only 0.05 root mean square error (RMSE) up to about E(G BP − G RP ) = 1.5 mag. After that, the observed extinction track begins to deviate from the linear line, shown as the bowl shape in the residual distribution. We further analyze this curvature of CERs in Section 3.3.

Curvature of Color Excess Ratios
The systematic curvature of CERs seen in Figure 1 is due to the assumption of a static wavelength for each filter in fitting the color excess -color excess diagrams. When star light goes through dust, the peak of spectral energy distribution (the effective wavelength at one band) shift toward the longer wavelengths gradually, as a result, the extinction degrades. The extinction at a given bandpass filter, namely evolving filter wavelength extinction, can be expressed as where F λ (λ) is the intrinsic flux of the stellar spectra, S(λ) is the filter transmission curve. We define R(λ) = 10 −0.4A λ,0 as the extinction factor, and A λ,0 as the static wavelength extinction. According to this formula, the gradually degradation of the extinction is unavoidable, unless the width of bandpass is infinitely narrow.
In this work, we use Gaia G BP and G RP bands as the basis bands because of their excellent photometric quality, but their broad bandwidth ∆λ would cause significant curvature in the color excess -color excess diagram. To analyze the curvature, we simulate the extinction of each filter band by Equation (5). We adopt the synthetic stellar spectra F λ (λ) for a RC star with T eff = 4800 K, log g = 2.5 and [Fe/H] = −0.1 (Lejeune et al. 1997), according to the average parameters of the whole RC sample, which are T eff = 4810 ± 143 K, log g = 2.5 ± 0.1 and [Fe/H] = −0.1 ± 0.2. The spectra are convolved by filter transmission curves of each photometric system. The extinction is generated using a CCM R V = 3.1 model extinction curve with V band extinction A V,0 from 0 to 20 mag in a step of 0.005 mag. For comparison, we introduced A 0 which denotes the extinction at wavelength 5500 nm with negligible bandwidth. The left-hand panels of Figure 2 present the comparison of the modeled A λ with the A 0 for each filter, whose slope is the relative extinction A λ /A 0 . The linear lines are obviously bent at A 0 > 4 mag in some bands, such as three Gaia bands, and u SDSS , g SDSS , g PS1 bands.
The right-hand panels of Figure 2 clearly show that the extinction difference between the evolving filter wavelengths and the static wavelengths ∆A λ = A λ − A λ,0 varies with the E(G BP − G RP ). If the filter wavelengths did not evolve with the progressively extinction, the extinction tracks, in either photometric system, would have remained at ∆A λ = 0. Actually, the ∆A λ clearly deviates from zero not only in broad Gaia bands, but also in APASS B, V , SDSS u SDSS , g SDSS , and PS1 g PS1 bands (marked in Figure 2). Although the deviation of ∆A λ in the IR bands is less obvious than that in the optical bands, it does exist. When E(G BP − G RP ) > 9 mag, corresponding to IR extinction E(H − K S ) > 1.2 from the reddening law of Table 2 In the previous literature, the curvature of CERs is not significant in the color -color or color excess -color excess diagrams, due to the low photometric accuracy or low extinction. One exception is Stead & Hoare (2009), who derived the NIR reddening law by considering the curvature of NIR CERs caused by filter wavelengths that evolve with the changing spectra of progressively reddened objects. As shown in Figure 2, the curvature of CERs appears at E(G BP − G RP ) ∼ 0.5 mag (A V ∼ 1.2 mag) for some bands. Therefore, the curvature correction requires special attentions as the quality of the photometry improves.
We use the extinction tracks in Figure 2 (b) and (d) to remove the observed curvature of CERs. Figures 3 and 4 are the color excess -color excess diagrams E(λ− G BP ) vs. E(G BP − G RP ) after the curvature correction, where λ include: two APASS bands B, V , five SDSS bands u SDSS , g SDSS , r SDSS , i SDSS , z SDSS , five PS1 bands g PS1 , r PS1 , i PS1 , z PS1 , y PS1 ,  Table 2. three 2MASS bands J, H, K S , and three WISE bands W 1, W 2, W 3. The Gaia G band color excess -color excess diagram is shown separately in Figure 10. The color represents the number density of RC stars. The black lines are the best linear fits and the slopes k λ are CERs. For each λ band, the number of RC stars that participate in the linear fit, the CERs E(λ − G RP )/E(G BP − G RP ), and the dispersion of the fit σ are tabulated in the second, third, and fourth columns of Table 2, respectively. The fitting error of the slope is only 0.001 which could not represent the real error of CER. Therefore, we discuss the uncertainties of CERs in Section 4.

Relative Extinction
In the conversion of the CER E(λ − G RP )/E(G BP − G RP ) into the relative extinction A λ /A GRP by Equation (3), the knowledge of A GBP /A GRP is required. Here, we independently calculate A GBP /A GRP by two methods.
Chen18 adopted the color excess -extinction method to determine the extinction along the sightlines toward the  Galactic center region by classical Cepheids. In that work, they determined the slope of CE versus the absolute extinction plus the relevant distance modulus (DM) diagram, and then converted the slope to the relative extinction. They proved that the extinction law determined by the color excess -extinction method is consistent with those derived by other methods including color -excess method. Inspired by that work, we combine the precise parallax and photometric data from the Gaia catalog to determine A GBP /A GRP by the observed color index (G BP − G RP ) versus G RP band apparent magnitude minus distance modulus (DM) G RP − DM = A GRP + M GRP diagram. We also refer to this method as the color excess -extinction method. If all the RC stars are not affected by extinction, they would appear as a clump with a small scatter in the diagram. The RC stars affected by extinction distribute along a linear line. Therefore, the relative extinction A GBP /A GRP can be derived from the slope of the linear fit. As mentioned in Section 2.2, our RC sample contains contaminations, such as SRC stars and RGB stars. SRC stars are less luminous than RC stars, which have ignited He in non-degenerate conditions (Girardi 1999). In the color -magnitude diagram, SRC stars are located in the bluer and fainter part compared to RC stars. In our sample, most SRC stars satisfy log g > 2.5 in APOGEE T eff -log g contour. Therefore, we only select RC stars within 2.35 ≤ log g ≤ 2.5 to eliminate SRC stars and RGB stars. After that, a subsample with 30,431 RC stars is obtained. We adopt distance converted from Gaia DR2 parallax with corrections by Bailer-Jones et al. (2018), who used the Bayesian inference approach to account for the nonlinearity of the transformation and the asymmetry of the resulting probability distribution. To reduce uncertainties, we also require the fractional error of the parallax to be less than 0.1, and the parallax to be larger than 0.25. After application of these selection criteria, only 5,051 RC stars remain. Figure 5 displays the distribution of these RC stars in the (G BP − G RP ) vs. G RP − DM diagram. In the range of (G BP − G RP ) < 1.5 mag, the uncertainty of DM per extinction is large, and it may affect the result of the linear fit. Therefore, we subdivided the sample by the selection cut (G BP − G RP ) > 1.0 + 0.001 * n, where n varies from 0 to 500. The linear fit is applying to each subsample. When (G BP − G RP ) > 1.341, the fit achieves the largest coefficient of determination R 2 and we accept it as the best fit. It is shown as the black line in Figure 5, and the slope is A GRP = (1.429 ± 0.015)E(G BP − G RP ), which corresponds to A GBP /A GRP = 1.700 ± 0.007. For comparison, we also use the color -excess method to determine the extinction A GBP /A GRP . In section 3.3, we have derived E(W 2 − G RP )/E(G BP − G RP ) = −1.388 by fitting the color excess -color excess diagram. The corresponding E(W 2 − G BP )/E(G BP − G RP ) equals E(W 2 − G RP )/E(G BP − G RP ) − 1 = −2.388. Then, the relative extinction A GBP /A GRP can be expressed as With the assumption of W 2 band extinction A W2 = 0, we derived the upper limit of A GBP /A GRP = 1.720. A GBP /A GRP = 1.700 ± 0.007 determined from the color -extinction method satisfies this criterion. In the literature, the relative extinction A W2 /A KS is 0.34 ± 0.10 (Xue et al. 2016;Wang et al. 2018;Chen18). By combining the A W2 /A KS with the CERs E(K S − G RP )/E(G BP − G RP ) and E(W 2 − G RP )/E(G BP − G RP ) obtained in Section 3.1, we derive A W2 /A GBP = 0.029 ± 0.013. Plugging it into Equation (6), the derived A GBP /A GRP = 1.686 ± 0.016 from the color -excess method is consistent with that of 1.700 ± 0.007 from color -extinction method. Based on the A GBP /A GRP derived here and the CER k λ listed in Table 2, the optical to mid-IR multi-band relative extinction A λ /A GRP can be estimated by Equation (3). These results are tabulated in Table 3 (column 3). Chen18 derived the IR extinction law A λ /A KS in 2MASS, WISE and Spitzer bands using classical Cepheids projected toward the Galactic center region by three methods. We converted their relative extinction values into the A λ /A GRP and listed them in column 4 of Table 3. Our near-IR measurements are in excellent agreement with their results. Our mid-IR WISE results agree with their results at 1.2 σ. Likewise, A λ /A V (column 5) and A λ /E(G BP − G RP ) (column 6) are derived and listed in Table 3. Note that the mid-IR WISE bands A λ /A V results are refined based on more

The Residuals of Color Excess -Color Excess Diagram
Although the most RC stars fall right along the fitted lines with only 0.001 fitting errors of the slopes, these errors seem to be too small to represent the real errors of CERs. In this section, we did a residual analysis to estimate the real errors of CERs. To achieve this, we plotted the residuals (CE minus fitted functions) ∆E(λ − G RP ) distribution as a function of CE E(G BP − G RP ) diagrams (Figures 6 and 7). The color shows the number density of RC stars. The black solid lines and dashed lines are the RMSE and the mean value of residuals (∆) for stars in bins of ∆E(G BP −G RP ) = 0.1 mag, respectively. The maximums of the RMSE and the ∆ are listed in the last two columns of Table 2 named as (RMSE) max and (∆) max .
At a given band, the photometric uncertainties of the observed color indices and the dispersions of the intrinsic color indices contribute to the residuals ∆E(λ − G RP ). As seen in Figures 6 and 7, the average scatter in residuals is small as a whole. Among the twenty-one bands used in this study, residuals in Gaia bands have the smallest dispersion with the (RMSE) max = 0.008 as shown in Figures 10 (b), which is a factor of 2 better than those in PS1 bands. The reason is that the data from Gaia have the lowest photometric uncertainties, and the scatter of intrinsic color of RC stars is small as well. In PS1/grizy bands and SDSS/griz bands, the (RMSE) max are around 0.01-0.03. In the IR and short-wavelength optical bands, limited by the accuracy of the photometric data, the scatters of the residuals are obviously large. For example, the (RMSE) max is around 0.05 mag in 2MASS and WISE bands, and about 0.07 mag  in B and u SDSS bands. As for ∆, all these bands have a systematic deviation around or much less than 0.02 mag, and it validates the process of the curvature correction. We investigate the effects of the RMSE and the ∆ on the CERs, and determine the statistical and systematic errors of the CERs in Section 4.2.

Simulation
The E(λ − G RP ) vs. E(G BP − G RP ) color excess -color excess diagrams are simulated to investigate the effects of x-axis and y-axis error on slopes (CERs). To obtain the simulated CEs E(G BP − G RP ) sim and E(λ − G RP ) sim , we execute a double Gaussian function to fit the CE E(G BP − G RP ) in Figures 3 and 4. The first Gaussian distribution represents the low extinction sources located in the solar neighborhood or high Galactic latitude. The second Gaussian distribution represents the high extinction sources in the disk. Hence, the simulation value E(G BP − G RP ) sim can be expressed as: where σ 1 , µ 1 , σ 2 , and µ 2 are the fitting parameters of the double Gaussian function. Then, the simulation value where k λ is the slope of the color excess -color excess diagram listed in Table 2. This process is applied to each band shown in Figures 3 and 4, First, we test the effects of x-axis errors on the slopes. The x-axis error (i.e. error of the E(G BP −G RP )) can be inferred from the dispersion of the fit for the CER E(G − G RP )/E(G BP − G RP ) (0.005 from Table 2). To avoid underestimating the uncertainty, the x-axis error is set to be 0.005. For each band, we generated both the E(G BP − G RP ) sim and the E(λ − G RP ) sim based on the distributions of the observed CEs. After that, by fitting the E(G BP − G RP ) sim + error vs. E(λ − G RP ) sim diagrams, we determined the slope (k x ) sim , the statistical error of the slope (σ x ) sim , and the deviation (∆ x ) sim between the (k x ) sim and the k λ . Thanks to high-precision Gaia photometry, the x-axis error introduces less than 0.002 deviation to the slopes.
Then, we analyzed the impact of y-axis errors on the slopes in a similar way. Two y-axis errors were considered: one is the (RMSE) max representing the local maximum scatter; the other one is the nonlinear function which is derived by a polynomial fit to the residuals (the dashed lines in Figures 6 and 7). After that, by fitting the E(G BP − G RP ) sim vs. E(λ − G RP ) sim + errors diagrams, we determined the slope (k y ) sim , the statistical error of the slope (σ y ) sim , and the deviation (∆ y ) sim between the (k y ) sim and the k λ . By the combination of the (σ x ) sim and the (σ y ) sim , we derived the statistical error, listed as the first error item in the third column of Table 2. The sum of the deviations (∆ x ) sim and (∆ y ) sim are considered as the systematic errors. They are listed as the second error item in the third column of Table 2. In conclusion, the total uncertainties of the slopes (CERs, Table 2) are mostly less than 0.02.
In addition, we analyzed the effects of various x-axis CE error σE x on the slopes. This process can explain why we adopt Gaia bands as the x-axis. We model the color excess -color excess diagram of 2MASS K S band as an example. The corresponding values of the parameters σ 1 , µ 1 , σ 2 , and µ 2 are determined. Then, we considered five cases with various parameters (n, n 1 /n 2 , and µ 2 ) to represent different initial conditions. More specially, the total number of objects n is 500 or 50000, the ratio of low to high extinction sources n 1 /n 2 is 2/3 or 3/2, and the average reddening amount of high extinction sources µ 2 is 0.25, 0.5, or 1.0 mag. Large µ 2 denotes high extinction. Lastly, a group of x-axis CE error σE x ranging from 0.0 to 0.25 was added to CE E(G BP − G RP ) sim . By fitting the E(G BP − G RP ) sim + σE x vs. E(K S − G RP ) sim diagram, we derived the slope (k KS ) sim . In comparison to the slope k KS of E(K S − G RP )/E(G BP − G RP ) determined in Section 3.2, we estimated the percentage deviation of the slope σ slope between these two slopes by | ((k KS ) sim − k KS )/k KS |. Figure 8 shows that the error of the slope σ slope varies as the CE error σE x . Generally, the σ slope increases with the increasing of σE x . For σE x < 0.02 mag, the slowly increasing errors of the slopes are small and almost similar in different cases. While for σE x > 0.06 mag, the errors of the slopes increase dramatically. As shown in the black dashed, solid and dotted lines, the increase of µ 2 could effectively reduce the error of the slope. In addition, the σ slope of the red line with n = 500 is larger than that of the solid black line with n = 50000. It means that increasing the total number of objects can reduce the error of the slope. Furthermore, we analyzed the influence caused by the ratio of low to high extinction sources n 1 /n 2 . The σ slope of the blue line with n 1 /n 2 = 3/2 is slightly higher than that of the solid black line with n1/n2 = 2/3 from beginning to end. It means that the higher proportion of high extinction sources, the smaller error of the slope. In conclusion, the error of the slope increases as the x-axis CE error increases. For the purpose of improving the precision of CER, the bands with the best photometric quality should be set as the basis x-axis bands. When the x-axis error σE x is notable, the three parameters, average reddening amount of high extinction sources µ 2 , total number of objects n, and ratio of low to high extinction sources n 1 /n 2 , will also exacerbate uncertainties of the extinction law. Figure 8 can also be used to explain the superiority of using Gaia bands rather than 2MASS bands as the basis of CER analysis. We redo a similar analysis of E(H − K S ) vs. E(J − H) based on about 60,000 2MASS-APOGEE RCs, and the result is shown as green line in Figure 8. The error of the slope is several times larger than that caused by adopting Gaia as basis bands. If we adopt the σE x = 0.03 mag for the 2MASS CE, the corresponding 1σ error of  Figure 8. The percentage deviation of the slope σ slope varies as the x-axis CE error σEx. Five cases are considered with different parameters n, n 1 /n 2 , and µ 2 , where n is the total number of objects, n 1 /n 2 is the ratio of sources with low and high extinction, and µ 2 is the average reddening amount of high extinction sources. The green line shows the similar simulation for about 60,000 2MASS-APOGEE RCs in E(H − K S ) vs. E(J − H).
the CER is around 8.6%. Besides, when the x-axis error σE x is large, the other parameters, including the average reddening value of high extinction sources and the sample number, will also have a great impact on the slope (CER).
In addition, as mentioned in Section 3.3 and Stead & Hoare (2009), for high extinction regions, such as the Galactic center and molecular clouds, the curvature of CER becomes obvious. The curvature of CERs also influences the measurement of the slope. In conclusion, all of these are the reasons why various near-IR CERs were reported in previous works (Section 1.2). Compared to the error in the E(J − H), the error in the E(G BP − G RP ) is much smaller (< 0.005 mag) ( Figure 6). No matter what the situation is, the effect on the error of the slope is less than 1% and could be ignored. Therefore, we recommend adopting high photometric precision bands (here, G BP and G RP ) as basis bands in CER analysis to reduce the error of the slope caused by the fitting method.

Extinction of individual RC sightlines
The red clump star is a standard candle as its luminosity has relatively week dependency on the stellar composition, color and age in the solar neighborhood (Paczyński & Stanek 1998;Alves 2000;Groenewegen 2008;Girardi 2016). The wavelength-dependent extinction to each RC star sightline can be estimated by the formula A λ = m λ − M λ − 5 log d + 5 with the known values of apparent magnitude m λ , absolute magnitude M λ , and distance d. For our RC sample, with the multi-band apparent magnitudes m λ from the photometric catalogs, the distance information d from the Gaia DR2 catalog, we only need the absolute magnitude M λ to estimated the extinction A λ . Compared to optical bands, the absolute magnitude of RC stars has a reduced systematic dependence on metallicity and effective temperature in IR bands. Among the RC absolute magnitudes given in the literature, the K S band absolute magnitude M KS has the most consistent value (Ruiz-Dern et al. 2018 and reference therein).
As our goal is to determine the accurate extinction laws toward RC sightlines, we evaluate the feasibility of this single star sightline method through error analysis. To achieve this, the fractional error of the extinction (A λ ) err /A λ is estimated at the given λ band. The extinction error (A λ ) err comes from the errors of photometry (m λ ) err , the absolute magnitude (M λ ) err , and the parallax d err . Two photometric bands are taken as examples, the optical SDSS g band which has larger extinction than IR bands, and the IR 2MASS K S band which has the most consistent M KS .
The typical photometric error is ∼ 0.03 mag at g and K S bands. We select RCs by the criteria of distance less than 4 kpc and the fractional error of the parallax less than 0.1. The average error of the distance to the RC stars is ∼ 0.18 mag. The absolute magnitudes of the RC stars are M g = 1.229 ± 0.172 mag ) and M KS = −1.61 ± 0.03 mag (Alves 2000;Ruiz-Dern et al. 2018). Therefore, the typical extinction errors (A λ ) err at g and K S bands are 0.19 mag and 0.25 mag. The average extinction values of individual RC sightlines at g and K S bands are A g = 2.8 mag and A KS = 0.2 mag, respectively. The extinction uncertainty of a single RC sightline is ∼ 9% at g band, and reach ∼ 100% at K S band. Compared to the accurate extinction law determined by the colorexcess method (Table 2), the extinction uncertainties of this method are significant. Therefore, the current single star sightline analysis is inferior to the statistical color -excess method and we only treat the RC targets as a whole to investigate the extinction law.   ( The CERs have been measured by a number of works in a variety of photometric bands. So, we compare our measurements with some measurements reported in the literature. In order to compare with publications in different combinations of bands E(λ 1 − λ 2 )/E(λ 2 − λ 3 ), we convert our CERs (Table 2) into the CERs in the corresponding bands (Table 4). For example, we use three CERs E(u − G RP )/E(G BP − G RP ), E(g − G RP )/E(G BP − G RP ), and E(r − G RP )/E(G BP − G RP ) to calculate the CER E(u − g)/E(g − r). The errors of these three CERs are propagated to the E(u − g)/E(g − r) as well. The CERs derived in this work and reported in previous works are tabulated in Table 4.

The Comparison of Color Excess Ratios
As a whole, the agreement between our measurements and those of the previous publications is excellent. It is worth noting that the precision of our CERs is better than the previous results. Our SDSS band results agree well with Schlafly & Finkbeiner (2011, work (a)). We compare seven CERs through optical to IR bands to the results of Yuan et al. (2013, work (b)). Four of them are in agreement with each other within 2σ. The agreement between ours and Yuan et al. (2013) in CERs E(r − i)/E(i − z) and E(J − H)/E(H − K S ) is out of 3σ. But these two agree with work (a) and Schlafly et al. (2016, work (c) Yuan et al. (2013) is nearly twice of ours. Overall, our CERs are closely consistent with those derived by Schlafly et al. (2016, work (c)), who also reported the difference between their results and Yuan et al. (2013) Schlafly et al. (2016) pointed out the measurement in Yuan et al. (2013) was uncertain because only the low reddening objects are available in their sample.
Our measurements rely on a large sample of RC stars with precise photometry and homogeneous stellar parameters. Even though each CER listed in Table 4 contains propagated errors of three CERs, our measurements are in great agreement with previous works. More importantly, the precision of our CERs has significantly improved. To compare with model extinction curves characterized by R V , we convert our relative extinction results A λ /A GRP to A λ /A V and calculate R V value as well. The extinction values A λ /A GRP and A λ /A V are list in the third and fifth columns of Table 3. According to the definition of R V , the total-to-selective extinction ratio, we derive

The Comparison with Extinction Curves
.16 ± 0.15. As this work has no bias toward any specific environment and covers all the fields surveyed by APOGEE, our measurements represent the average extinction. It is in agreement with the average value of the Galactic diffuse ISM R V = 3.1. Fitzpatrick & Massa (2009) suggested that the CER E(K − V )/E(B − V ) can be used to estimate R V with 0.12 uncertainty. Therefore, we estimated the R V value by this method as well. The determined R V = 3.19 is also consistent with the R V determined from the total-to-selective extinction ratio. Schlafly et al. (2016) used the CER E(g PS1 − W 2)/E(g PS1 − r PS1 ) to obtain the R ′ V as a proxy of R V . Based on the definition, we derive the R ′ V = 3.21 from our measurements which is consistent with the average value 3.32 from Schlafly et al. (2016). Our R V is determined in static wavelengths (reddening/extinction after the curvature correction), so a slight difference exists when comparing with the previous R V values.
As the calculation of R V value only depends on the extinction between two optical bands, or the CER of three bands, this value can not perfectly reflect the wavelength-dependent extinction law. Therefore, we compare our optical to mid-IR multi-band extinction with different CCM model extinction curves in Figure 9. For each band, we uniformly calculate the static effective wavelength λ eff,0 through As mentioned in Section 3.3, the synthetic stellar spectra F λ (λ) are based on a RC star with T eff = 4800 K, log g = 2.5 and [Fe/H] = −0.1 (Lejeune et al. 1997), and S(λ) is the filter transmission curve. The calculated effective wavelengths λ eff,0 are tabulated in the second column of Table 3. In Figure 9 (a), our extinction results are plotted as red dots with error bars. The CCM R V = 3.1 (black line), 2.5 (green line), 2.1 (blue line) model extinction curves are also shown. Figure 9 (b) is the ratio of model extinction values from CCM R V = 3.1 (black), 2.5 (green), 2.1 (blue) to our observed extinction values A λ (Model)/A λ at the given band. As shown in Figure 9, our optical extinction law conforms with the R V = 3.1 extinction curve (< 2% deviation) in the wavelength range of 300 -550 nm. While at longer wavelengths, the new determined extinction law is significantly steeper than the R V = 3.1 extinction curve. More specifically, it agrees with the R V = 2.1 extinction curve in the wavelength range of 600 -770 nm, and agrees with the R V = 2.5 in range of 770 -1,000 nm. This feature is consistent with the R V ≈ 2.5 reported by Nataf et al. (2013), who investigated the extinction law toward the Galactic bulge based on V and I bands. In near-IR bands (0.9 µm < λ < 3 µm), none of these models can perfectly explain the quite steep trend, so a larger power-law index α is needed. To better describe the observed extinction law, we made some adjustments on the CCM R V = 3.1 extinction curve to derive the new equations given below.
Optical: 0.3 µm < λ < 1.0 µm and Y = 1/λ( µm) − 1.82, Near-infrared: 1.0 µm ≤ λ < 3.33 µm, These equations obey the form of the CCM model, while the coefficients are reanalyzed, through best fitting to our extinction in nineteen optical and near-IR bands including u, B, g SDSS , g PS1 , G BP , V , r SDSS , r PS1 , G, i SDSS , i PS1 , G RP , z PS1 , z SDSS , y PS1 , J, H, K S , and W 1 bands. The adjusted extinction curve, shown as the red dashed line in Figure 9, is consistent with the observed extinction law by better than 2.5%. Note that our new extinction law (Equations (9) and (10)) is the continuous extinction curve between 0.3 µm and 3.33 µm. The absorption features, such as ice, hydrocarbons features presenting at 3 µm with strength varying between lines of sight (e.g., Draine 2003;Fritz et al. 2011, and references therein), are not contained in the extinction curve.
The Galactic center is an ideal region to investigate the IR extinction due to the well measured Galactocentric distance (de Grijs & Bono 2016; Bland-Hawthorn & Gerhard 2016) and relatively high extinction. Along the Galactic center sightlines, steep near-IR extinction laws also have been reported in the literature. Schödel et al. (2010) derived α = 2.21 by measuring extinction between H and K bands for RC stars. Fritz et al. (2011) found α = 2.11 ± 0.06 via measurement of hydrogen emission lines. Most recently, Chen18 reported α = 2.05 ± 0.07 traced by classical Cepheids. Even steeper Galactic center extinction law has been reported by Nogueras-Lara et al. (2018) with α = 2.31 ± 0.03. Generally, our index value α = 2.07 ± 0.03 is consistent with those of Fritz et al. (2011) and Chen18. At the same time, due to much better photometry and parallax of Gaia, more pure sample of RC stars from APOGEE, and robust determination method, we recommend a steep average near-IR extinction with α = 2.07 instead of the CCM extinction with α = 1.61 in future work. Figure 10 exhibits the determination of Gaia extinction coefficient. As shown in Figure 10 (a), the distribution of RC stars in the color -color diagram exhibits good linearity and can be fitted by the black straight line up to about G BP − G RP = 2.0 mag. However, at color G BP − G RP 2.0, the distribution begins to curve. It displays great amounts of curvature for heavily reddened objects. The curve feature is even more clear in the residual distribution diagram: bottom panel of Figure 10 (a). As discussed in Section 3.3, the curvature of CER is due to the assumption of the static wavelength. Because of the broad bandwidth in the Gaia bands, the curvature is particularly evident. To obtain the real extinction at the Gaia bands, we correct the observed curvature by using our model extinction tracks (Figure 2). Figure 10 (b) displays the color excess -color excess diagram after the curvature correction. The distribution of  Table 3 as well.

The Gaia Extinction Coefficient
We compare Gaia bands extinction values with some other optical bands extinction in Figure 9. Our optical continuous extinction curve smoothly varies as wavelength. The extinctions in G BP , G, G RP bands are close to those in V, r, i bands, respectively, as their static effective wavelengths are close to each other (Table 3). This agreement proves the reliability of our Gaia extinction results.
It is worth noting that our extinction law in Table 3 and Figure 9 represents the static extinction law. Because of the non-ignorable bandwidth of the filter, the existence of curvature in reddening/extinction is unavoidable. This curvature depends on the spectral type, the filter system and the amounts of extinction. Since this curvature exacerbates with the increase of CE or extinction (Figure 2), to avoid systematic error, extinction laws are only suitable to estimate and correct extinctions for low extinction objects. For objects with moderate or heavy extinction, the extinction law needs a small correction, that is R λ,c = R λ * A λ /A λ,0 . λ denotes the band of interest, R λ,c and R λ are the corrected and the static extinction law. A λ and A λ,0 are the evolving and the static band extinction estimated by the combination of the static extinction law, the synthetic stellar spectra, and the filter transmission curve (see details in Section 3.3). Similarly, for the correction of the relative extinction, the formula is A λ1,c /A λ2,c = A λ1 /A λ2 (different to A λ1,0 /A λ2,0 ), where λ 1 and λ 2 are two bands of interest.

The Predicted Extinction in HST WFC3 and JWST NIRCAM Bandpasses
Based on the determined near-IR extinction law, we can predict the relative extinction in the other near-IR bandpasses. The relative extinction values A λ /A V and A λ /A KS for the near-IR bandpasses of the HST WFC3 and the JWST NIRCAM are evaluated. The adopted effective wavelength λ eff 2 and the predicted extinction results A λ /A V and A λ /A KS are tabulated in Table 5. The accuracy of these predicted near-IR extinction is about 2.5%. Note that the IR absorption features (Draine 2003;Fritz et al. 2011;Wang et al. 2013), such as ice features at 3.1 µm (H 2 O), 4.27 µm (CO 2 ), and 4.67 µm (CO), aliphatic hydrocarbons feature at 3.4 µm, can also affect some entries in Table 5.

CONCLUSION
We have investigated the optical to mid-IR extinction law for a group of RC stars which were selected by the stellar parameters from APOGEE survey. The multi-band photometric data are collected from Gaia, APASS, SDSS, Pan-STARRS1, 2MASS, and WISE surveys. As the extinction tracers (RC stars) cover all the fields surveyed by APOGEE, hence, our measurements represent the average extinction. Thanks to the unprecedented Gaia data, not Otherwise, extinction laws are only safe to be used for low reddened objects or photometric system with extremely narrow bandwidth. 5. In the determination of CERs by fitting the color excess -color excess diagrams, the x-axis error (including photometric error and intrinsic color error), the sample number, and the average reddening value of high extinction sources have great impact on the slope (CER), especially when the x-axis error is large. These are the reasons why various near-IR CERs were reported in previous works. Compared to the E(J − K S ), the E(G BP − G RP ) is at least a factor of 3 more precise. Therefore, we recommend adopting high photometric precision bands (here, G BP and G RP ) as basis bands in the CER analysis to reduce the error of the slope caused by the fitting method.
6. Based on the determined near-IR extinction law, we predict the relative extinction values A λ /A V and A λ /A KS for the near-IR bandpasses of the HST WFC3 and the JWST NIRCAM.