FEDReD II : 3D Extinction Map with 2MASS and Gaia DR2 data

Aims. We aim to map the 3D distribution of the interstellar extinction of the Milky Way disk up to distances larger than those probed with the Gaia parallax alone. Methods. We apply the FEDReD (Field Extinction-Distance Relation Deconvolver) algorithm to the 2MASS near-infrared photometry together with the Gaia DR2 astrometry and photometry. This algorithm uses a Bayesian deconvolution approach, based on an empirical HR-diagram representative of the local thin disk, in order to map the extinction as a function of distance of various fields of view. Results. We analysed more than 5.6 million stars to obtain an extinction map of the entire Galactic disk within $|b|<0.24^{\circ}$. This map provides information up to $5~\mathrm{kpc}$ in the direction of the Galactic centre and at more than $7~\mathrm{kpc}$ in the direction of the anticentre. This map reveals the complete shape of structures known locally, such as the Vela complex or the split of the local arm. Furthermore our extinction map shows many large"clean bubbles"especially one in the Sagittarius -- Carina complex, and four others which define a structure that we nickname the butterfly.


Introduction
The interstellar extinction attenuates the light coming from background objects. Moreover this attenuation is a function of the wavelength of radiation, which causes the reddening phenomenon. The study of extinction is mandatory to recover the absolute magnitude and intrinsic colours of Galactic or extragalactic stars.
Mapping the extinction also provides access to the distribution of dust. This component of the Galactic disk is itself an interesting key to understand the evolution of the Milky Way. In fact, high dust density areas are expected to be associated with high star formation regions. Mapping the extinction of the Galactic disk is thus a way to study the spatial structure of spiral arms.
In the past decades, several methods have been developed in order to draw 3D extinction maps of the Galactic disk. With limited resolution and distances, the first results were obtained by Fitzgerald (1968); Neckel & Klare (1980); Arenou et al. (1992); Hakkila et al. (1997). Marshall et al. (2006) compared the 2MASS photometry to the stellar population synthesis of the Besançon Model (Robin et al. 2012) to obtain extinction density. This approach was also used by Chen et al. (2013) and Schultheis et al. (2014) which added Glimpse and VVV catalogue to their dataset. Sale et al. (2014) applied a hierarchical Bayesian method to infer stellar properties and extinction profile with the IPHAS survey. Also using a Bayesian method, Green et al. (2014) derived an extinction map from Pan-STARRS data, later combined with 2MASS (Green et al. 2015, 2018  inversion on APOGEE data to obtain density map of the solar neighbourhood. Chen et al. (2019b) used a random forest algorithm on the crossmatch of Gaia DR2, 2MASS and WISE, trained on spectroscopic data, to infer the local extinction density.
In this work we use the Field Extinction-Distance Relation Deconvolver (FEDReD) algorithm described in Babusiaux et al. (2020). This is a Bayesian deconvolution algorithm which uses an empirical H-R diagram to study photometry and parallax, taking into account the completeness of the field of view under study, to derive both the extinction as a function of distance and the stellar distance distribution. We apply the algorithm to data of the Two Micron All Sky Survey (2MASS, Skrutskie et al. 2006) crossmatched with the second Data Release of Gaia (DR2, Gaia Collaboration et al. 2018). This paper is organised as follows : in Section 2 we present the data and filters used to select stars. Section 3 is a quick sum up of the FEDReD method. In Section 4 we detail our method to merge every field of view results into a self-consistent map. Section 5 presents our extinction map, with the details of the different visible features and their relation with the other components of the Galactic disk.

Data
This work uses photometry and astrometry based on two surveys, Gaia and 2MASS. Gaia DR2 data provides high precision parallax and photometry in G, G BP and G RP bands 2MASS provides photometry in three near -infrared bands: J, H and K s . As Gaia DR2 and 2MASS are full sky surveys, we are able to draw an extinction map of the entire Galactic disk.
The 2MASS catalogue is used as a basis for our analysis, as its completeness is easier to model than the Gaia one, the latter being more dependent on the crowding as well as on the Gaia's scanning law. Every star that we study has therefore 2MASS 1 arXiv:2007.03734v1 [astro-ph.GA] 7 Jul 2020 C. Hottier et al.: FEDReD II : 3D Extinction Map with 2MASS and Gaia DR2 data photometry which can be completed with Gaia-DR2 parallax, photometry or both.
To select stars in the 2MASS catalogue, we used the ph qual flag and kept every star with at least a ph qual=D for each photometric band.
We use the Gaia-2MASS crossmatch provided by Marrese et al. (2019). When a 2MASS source had more than one Gaia best neighbour, we did not associate any Gaia information to the 2MASS source.
For the Gaia photometry (Evans et al. 2018), we do not use G BP and G RP photometry for stars affected by crowding issues using the filter phot bp rp excess factor > 1.3 + 0.06 × (G BP − G RP ) 2 and we do not use G BP information for the faint stars with G BP > 18 which are affected by background underestimation (see Evans et al. 2018;Arenou et al. 2018). Moreover, we add quadratically 10 mmag to the uncertainties to take into account the systematics. We also take into account the 3 mmag/mag drift on the G band (Arenou et al. 2018;Weiler 2018).

FEDReD in a nutshell
In this work, we use the Field Extinction-Distance Relation Deconvolver (FEDReD) method, presented in Babusiaux et al. (2020). We just sum up here the general process of FEDReD to analyse a line of sight (LoS) and infer the relation between extinction and distance.
The algorithm works in two separated steps. The first one deals with the individual analysis of each star contained in the LoS. It looks for P O j | A 0 , D , i.e. the likelihood of an observed star O j (considering its apparent magnitudes and possibly its parallax) to be at a distance D with extinction A 0 (absorption at 550 nm). To compute this probability, FEDReD compare the apparent photometry of the star to an empirical Hertzprung-Russel diagram based on Gaia DR2 representative of local thin disk stars. We take into account the colour and extinction dependant extinction coefficients by using Danielski et al. (2018) models with the same coefficients as used in Lallement et al. (2019).
Once we have computed every stars' density P O j | A 0 , D in the LoS, we merge them to obtain an estimate of the joint distribution of extinction and distance of the entire field of view, P (A 0 , D), with a Bayesian iterative Richardson-Lucy deconvolution (Richardson 1972;Lucy 1974). In other words, the prior of P k (A 0 , D) of the k th iteration is the posterior of the k − 1 th iteration. To initiate the process, we build a simple prior P 0 (A 0 , D) by multiplying two prior distributions, the distance distribution of stars, P 0 (D) and the distribution of extinction given the distance, P 0 (A 0 | D). The prior P 0 (D) just follows a square law of the distance to take into account the cone effect. The prior P 0 (A 0 | D) is null where A 0 > 10×D and flat elsewhere, this condition being experimentally verified using the map of Lallement et al. As we take into account the completeness, the result of the deconvolution is actually the estimate of the probability distri-butionP (A 0 , D | S ) where S is the completeness of the LoS, which is here the completeness of the 2MASS photometry. To model the completeness we estimate the probability distribution of P (S | A 0 , D) using the empirical HR diagram and a rough completeness model. This distribution is used both during the deconvolution process to obtainP (A 0 , D | S ) and from it to derive the finalP (A 0 , D) distribution.
To obtain the relation A 0 (D) from the previous distribution, we use a Monte-Carlo process to draw monotonic increasing relations in the distribution ofP (A 0 , D). We assess the probability of each Monte-Carlo Solution (MCS) and we keep the 1 000 best MCSs, which correspond to the 1 000 best relations A 0 (D) of the LoS. Finally, we fit a constrained median cubic spline through the MCSs using cobs R library (Ng & Maechler 2007), to get the best fit for the given LoS.
We discretise the extinction -distance space that we probe. We choose here for the extinction a sampling from 0 to 30 mag with a step of 0.05 mag. Concerning the distance space it covers 0.1 to 30 kpc but the step is linear in the distance modulus space with a width of 0.05 mag.

Merging line of sight results
To obtain an extinction map of the Galactic disk, we split it into small fields of view and we analyse them separately with FEDReD. To ensure the continuity between LoS, half of each LoS is shared with its 2 neighbour fields. So in practice it means that each star is contributing to the information of two contiguous LoS. This allows an efficient post-processing. Indeed, the best fit output of FEDReD can be polluted by two major effects. On one hand, the deconvolution can lead to very noisy P (A 0 , D), in particular in crowded fields where outliers can be numerous. On the other hand, several dust clumps could be present in the field of view and create an angular differential extinction which is difficult to handle with a single median spline fit. Thus we use the MCS information coming from neighbours fields of view to remove outliers solutions and smooth the results.
To do so, we remove every MCS which is not included in the envelope drawn by the maximum and minimum values of the two neighbour MCS envelops. We repeat the operation until no more MCS is removed, 30 to 60 iterations being usually needed to clean every LoS. This step converges thanks to the fact that two neighbour fields overlap, ensuring the consistency of two consecutive fields of view.
Once our sample is cleaned we can draw a map by picking randomly one of the remaining MCS for each LoS. We smooth this map by averaging each field with its two direct neighbours applying weights 1/2, 1/4 and 1/4 for centre field and the sides respectively, as each LoS is shared with its 2 neighbour fields. We draw 1 000 of these smoothed maps. We then fit a median constrained cubic spline on solutions using cobs to get the A 0 (D) relation of each LoS.
To obtain the map of extinction density a 0 we decumulate the A 0 (D) relations by processing the difference between subsequent distance bins and normalise by the width of each distance bin.
Finally, we determine minimum and maximum valid distance for each field of view. This distance interval is defined by the observability of a red clump star in the infrared photometric band, considering the completeness of the field of view. We use red clump stars because their small intrinsic dispersion of absolute magnitude and colour leads to individual distributions P O j | A 0 , D being more peaked than other stellar types, so they bring more constraints (Babusiaux et al. 2020). We process the theoretical apparent magnitude of a red clump star by using the absolute magnitudes inferred by Ruiz-Dern et al. the best MCS for the A 0 (D) relation. Valid distances correspond to the ones where the red clump apparent magnitude is between the saturation and the completeness magnitudes of the field.

Results
The method described above has been applied to 3 764 LoS. Each LoS is centred at the Galactic latitude b = 0 • with a latitude width of 0.48 • . The longitude width and position are defined so as to obtain fields of view which contain 3 000 stars and share 1 500 stars with each neighbours. This overlap with adjacent fields is necessary to properly clean MCS by the neighbours minimum, maximum envelope. The average width of a LoS in first and fourth quadrant is 0.12 • whereas the average width is 0.43 • in the second and third quadrant. This difference is due to the larger number of observed stars in the central region of the Galaxy. We represent on Figure 1 the resulting extinction density map in the Galactic plane. The white area centred on the Sun corresponds to the too close distances where the red clump is saturated (see previous section).

Uncertainty on the extinction and on the extinction density
By using the clean sample of MCSs we can obtain the minimum and maximum cumulated extinction at each distances (A 0min (D) and A 0max (D)).
We also compute an error map of the extinction density. Basically, we use our sample of MCSs cleaned of outliers and we apply a bootstrap technique. To do so, we randomly keep a subsample of MCSs for each line of sight and we apply the algorithm previously described to obtain a bootstrapped map. We compute 100 of these bootstrapped maps and then process the standard deviation at each (l, D) location to obtain the extinction density uncertainty. Figure 2 presents the relative uncertainty.
This error map mostly represents the sampling error, so it underestimates the true uncertainty of our results. The few locations where the relative uncertainty is really high correspond to areas with almost no extinction or areas with a high-increasing rate of extinction density (i.e. a high da 0 /dD) so the algorithm can not locate exactly the beginning of the cloud. We also see that the uncertainty increase at large distances towards the anticentre.

Main Structures
On Figure 3 we label the main features visible on our map. We also overplot some other ISM tracers on the map, in order to obtain a better view of the ISM structures. In Figure 4  The biggest structure in our map is the Sagittarius-Carinae complex, which extends in a cone between = 300 • and = 30 • . Inside this high extinction area we can notice a cavity centred on = 335 • , D = 2 kpc. This cavity presents an artefact which seems to be a Finger-of-God. This complex also contains a high density of molecular clouds and Hii regions except in the cavity, which confirms the location and the size of this clean area.
Several high extinction structures related to spiral arms can be seen in the first quadrant. The Local arm, at = 80 • , presents a very strong extinction consistent with the local arm masers of Reid et al. (2019). The split (Lallement et al. 2019), between = 30 • and = 40 • , is almost parallel to the Local arm. We can also notice the extinction overdensity at = 60 • between D = 1.5 kpc and D = 3.5 kpc which corresponds to the spur (Xu et al. 2018).
In the second quadrant, two main extinction overdensities appear. The first one is at = 111 • and D = 3 kpc and is associated with the Cassiopeia -Cepheus complex (Ungerechts et al. 2000). It coincides with the Reid et al. (2019) masers associated with the Perseus arm, and it is also well-marked by Hii regions. The second overdensity is the Cameliopardalis-Cassiopeia cloud (Chen et al. 2014) which is at = 145 • and 1 ≤ D ≤ 3 kpc.
The third quadrant contains many small high extinction areas associated with the local arm according to the masers or to the Radcliffe wave (Alves et al. 2020). There is also the Gemini molecular cloud (Carpenter et al. 1995) at ( = 190 • , D = 2 kpc) related to masers of the Perseus arm.
The Vela cloud overlaps the third and the fourth quadrants. It exhibits a strong overdensity in the foreground and is prolonged by an empty bubble surrounded by a thin extinction edge. The distant boundary of this bubble is marked by molecular clouds and an HII region.
Near the Vela cloud, at = 282 • and D = 1.6 kpc, the Carina complex (Zhang et al. 2001) also presents a strong foreground structure and is prolonged up to D = 6 kpc. This large elongation is also drawn by molecular clouds but the foreground structure only appears in extinction and is also well visible in Lallement et al. (2019).
In addition to all of these structures, the Sun appears surrounded by four bubbles without extinction, without taking into account the local bubble. The first one is delimited by the Sagittarius-Carinae complex and the split. The second is centred on = 165 • .The third one lies between the Local arm masers's and the Vela cloud. Finally, the longest one is surrounded by the Vela cloud and Carina-complex on one side and by the Sagittarius-Carinae on the other side. Those four extinction bubbles draw a structure that we will nickname "the Butterfly". Those four bubbles also present a lack of every other ISM tracers. While most of the masers and Hii regions used to create these models have a footprint on the dust map ( Fig.6 and 5), the dust behaviour is very patchy and does not lead to obvious continuous spiral arm footprints.
The local arm is very well defined by the extinction in the first quadrant, however, its path in the third/fourth quadrant is not obvious. According to masers, it seems that the Local arm does not really reach the Sun and miss the Vela complex to follow the Radcliffe wave (Alves et al. 2020). On the other hand, following the Hou & Han (2014) model, the Vela complex is part of the Local arm.
The Sagittarius arm shows a strong extinction but it is cut in two parts by the first clean bubble. Moreover, it appears that the split of the Local arm, even if it seems to be connected to the local arm at close distance, is more related to the Sagittarius arm at larger distance. The spur looks like an extinction bridge between the Local and the Sagittarius arms.
The Scutum arm crosses a high extinction area, however, we cannot distinguish a clear separation between the Scutum and the Sagittarius extinction areas. They are merged in the Sagittarius-Carina complex.
The Perseus arm is mostly visible at masers's locations, it presents a very patchy structure (Baba et al. 2018) and is only noticeable thanks to the visual guide of Reid et al. (2014). Furthermore, the two masers of this arm at the bottom of the third quadrant are not even visible in our extinction map because their extinction signature is too small for our spatial resolution.
The outer arm is not visible in our extinction map. This is due to its high distance from the Galactic plane.

Comparison with other work
We compared our extinction map to one of the most complete available extinction map of the solar neighbourhood, by Lallement et al. (2019). They basically use the same data as our study (2MASS and Gaia-DR2) limited to stars with a relative parallax error better than 20%. Their results are presented in the reference extinction at 550nm, as we do in this work. We present in Figure 8    The reader can notice the very good agreement between the two maps within the confidence limit of Lallement et al. (2019) results. The presence of Finger-of-God on both maps is also apparent, though at different locations on the two maps, allowing us to really identify them as spurious.
The foreground part of the Vela cloud is roughly identical on both works, which confirms this crescent shape. This particular    However, we detect structures and prolongation of structures at larger distance which did not appear in previous results in the

Conclusion
We have used the Gaia DR2 data and 2MASS photometry to estimate the extinction within the Galactic plane. We used data of about 5.6 million stars to produce an extinction map, reaching about 4 kpc in the direction of the Galactic centre and more than 6 kpc in other directions. Thanks to this result we are able to confirm general structures of the solar neighbourhood, already revealed by recent previous works (Lallement et al. 2019;Green et al. 2019;Chen et al. 2019b). The Local arm and the Split represent two structures, close to the Sun, with a high extinction separated by a clean corridor. Nevertheless our map reveals a possible relation between the Split and the Sagittarius arm at larger distance. The Galactic centre direction is dominated by the Sagittarius -Carina complex. Because of this high extinction, as well as crowding, we are not able to go beyond this complex. We do not distinguish the Sagittarius and the Scutum arm within this complex. The Perseus arm extinction component seems very fragmented, even if the location of masers which trace this arm are visible in the dust. The Vela and Carina complex are the two strongest extinction areas of the forth quadrant. Furthermore, FEDReD's map reveals the extinction prolongation of the Carina Complex and the bubble structure behind the Vela cloud. Finally, we also observe four empty bubbles close to the Sun, and we reveal their ends for two of them.
Data corresponding to the maps presented here are available in Tables 1 and 3. We also provide the values of the cumulated extinction A 0 and corresponding asymmetric uncertainties in Tables 2, 4 and 5. The entire version will be available as a machine-readable form at the CDS In future works, we will also explore the Galactic disk at higher and lower latitudes. This will allow a better study of some structures, for example the Outer arm which is known to be warped.
To explore larger distances, and in particular the Galactic centre and the start of spiral arms, we will use deeper nearinfrared surveys, such as UKIDSS (Lucas et al. 2008) and VISTA (Minniti et al. 2010). Moreover we look forward for the future third data release of Gaia mission which will provide better parallax and photometry constraints. Table 1. Some columns and rows of the a 0 extinction density (in mag/kpc) of the map presented in Figure 1. Each column corresponds to a LoS, the longitude is given in the first row (in degree). The first column corresponds to the heliocentric distance in kiloparsecs. Blank spaces represent spatial locations outside of the distance confidence intervals (see section 4). The full table is available in electronic form at the CDS.  Table 3. Some columns and rows of the a 0 density uncertainty σ a 0 (in mag/kpc) presented in Figure 2 (see section 5.1). Each column corresponds to a LoS, the longitude is the first row (in degree). The first column corresponds to the heliocentric distance in kiloparsecs. Blank spaces represent spatial locations outside of the distance confidence intervals. The full table is available in electronic form at the CDS.  Table 4. Some columns and rows of the maximum values of the extinction A 0 (see section 5.1). Each column corresponds to a LoS, the longitude is the first row (in degree). The first column corresponds to the heliocentric distance in kiloparsecs. Blank spaces represent spatial locations outside of the distance confidence intervals. The full table is available in electronic form at the CDS.  Table 5. Some columns and rows of the minimum values of extinction A 0 (see section 5.1). Each column corresponds to a LoS, the longitude is the first row (in degree). The first column corresponds to the heliocentric distance in kiloparsecs. Blank spaces represent spatial locations outside of the distance confidence intervals. The full table is available in electronic form at the CDS.