FEDReD III : Unraveling the 3D structure of Vela

Context. The Vela complex is a region of the sky that gathers several stellar and interstellar structures in a few hundred square degrees. Aims. Gaia data now allows us to obtain a 3D view of the Vela interstellar structures through the dust extinction. Methods. We used the FEDReD (Field Extinction-Distance Relation Deconvolver) algorithm on near-infrared 2MASS data, cross-matched with the Gaia DR2 catalogue, to obtain a 3D cube of extinction density. We applied the FellWalker algorithm on this cube to locate clumps and dense structures. Results. We analysed 18 million stars on $450~\mathrm{deg}^2$ to obtain the extinction density of the Vela complex from 0.5 to 8~kpc at $\ell\in[250\degr,280\degr]$ and $b\in[$-10$\degr,5\degr]$. This cube reveals the complete morphology of known structures and relations between them. In particular, we show that the Vela Molecular Ridge is more likely composed by three substructures instead of four, as suggested by the 2D densities. These substructures form the shell of a large cavity. This cavity is visually aligned with the Vela supernova remnant but located at a larger distance. We provide a catalogue of location, distance, size and total dust content of ISM clumps that we extracted from the extinction density cube.


Introduction
A large number of structures are known towards the Vela constellation. A view from the Improved Reprocessing of the IRAS Survey (Miville-Deschênes & Lagache 2005, IRIS) of this sky area is displayed in Fig. 1. This is the Vela complex. All of its structures reside in the region defined by 240 • ≤ ≤ 280 • and 6 • ≤ b ≤ −20 • (Pettersson 2008). As stellar structures, we count several OB association as well as R associations (Pettersson 2008). However, the Vela complex also presents several large interstellar structures.
The Gum nebula (Gum 1952) is the largest structure of the area (Pettersson 2008). It is a ring shell structure with an angular diameter of 34 • centred at ( , b) = (262 • , −3 • ) (Reynoso & Dubner 1997) at a distance of 450pc. The nebula is expanding and shares this expansion with the cometary globules (Woermann et al. 2001) which are a set of dark clouds with tails pointing to the centre of the Gum nebula (Pettersson 2008).
The IRAS Vela shell (hereafter IVS) is a ring shell first noticed in IRAS data by Sahu (1992). It is centred at ( , b) = (263 • , −7 • ) at the distance of 450pc (Pettersson 2008). This shell is related to the Vela OB2 association and may have a common progenitor and may share their history (Cantat-Gaudin et al. 2019, and references therein).
It is canonically accepted that all these structures are quite close to the Sun. The Vela complex nevertheless also presents a large ridge in the background. It was first noticed in CO Survey (Dame et al. 1987;May et al. 1988) and named Vela Molecular Ridge (VMR) by Murphy & May (1991). This structure can be split in four parts, A, B, C and D, according to CO peaks. Liseau et al. (1992) estimated the distances to VMR A, C and D at 0.7 ± 0.2kpc and VMR B at 2kpc, using infrared photometry which could imply that VMR B is not related to the A, C and D parts. More recently, Massi et al. (2019) used stellar parallax from Gaia DR2 (Gaia Collaboration et al. 2018) to estimate the VMR C distances farther at D = 0.95 ± 0.05 kpc.
As shown by the previous paragraph, the Vela region has been studied through several markers and techniques, but the extinction component, especially in 3D, is poorly exploited compared to other tracers. Franco (2012) analyses the interstellar reddening in six distinct areas using photometry to study the largest structures of the Vela complex except the densest parts of the complex.
On the other hand, there are some 3D extinction maps, but some of them are focused on other regions, such as Chen et al. (2013) or Schultheis et al. (2014), which compare photometric surveys to the Besançon model (Robin et al. 2012) to study the Galactic bulge. Sale et al. (2014) and Green et al. (2018) both used a Bayesian approach, respectively on IPHAS survey and on Pan-STARRS and 2MASS, to derive extinctions. Rezaei Kh. et al. (2018) analysed APOGEE survey with a non-parametric 3D inversion to obtain the dust density. But these three techniques used northern sky surveys, so they are not able to reach Vela.
There are also some 3D extinction maps which cover the Vela direction. Capitanio et al. (2017) and Lallement et al. (2019) applied the 3D inversion technique of Vergely et al. (2001) on composite dust proxies using respectively the Gaia DR1 parallaxes and the 2MASS crossmatch with Gaia DR2. Chen et al. (2019) used a random forest algorithm on a dataset built with Gaia DR2, WISE and 2MASS to compute the extinction density. Hottier et al. (2020)  with the FEDReD algorithm ) and derived the extinction density in the Galactic plane.
However, as far as we know, there is no study focussing on the 3D extinction structure of the Vela complex. In a previous work , we noticed that the Vela complex presents large structures that reach distances up to 5kpc. We also discussed in Hottier et al. (2020) the possibility that the Vela complex could belong to the local arm. This belonging seems to be confirmed by the spiral arms locations found by Khoperskov et al. (2020). Studying the Vela complex could therefore give the opportunity to study a spiral arm viewed from the inside.
In this work we used FEDReD to probe the 3D extinction density distribution towards Vela. In section 2 we present the analysed dataset. Section 3 sums up the FEDReD algorithm and develops the differences from Hottier et al. (2020). In Sec 4 we present the method of clump extraction. Section 5 describes and analyses the three dimensional density map and the clumps and cavities extracted from it.

Data
In this study, we use the infrared photometry in bands J, H and K from the 2MASS survey (Skrutskie et al. 2006) and we combine them with the photometry in bands G, G BP and G RP and the astrometry both coming from the Gaia DR2 (Gaia Collaboration et al. 2018). We use the same methodology as Hottier et al. (2020) to filter and merge data of these surveys, so we will just briefly review them here.
We use the 2MASS near-infrared photometry as principal data, that is our dataset completeness will be the 2MASS one. In practice, every star included in our dataset has ph_qual≥D in all three J, H and K bands. Once we have selected the stars in 2MASS, we use the Marrese et al. (2019) crossmatch to potentially add Gaia DR2 parallaxes and photometry.
To filter the Gaia photometry (Evans et al. 2018), we do not use G BP and G RP when phot_bp_rp_excess_factor > 1.3 + 0.06 × (G BP − G RP ) 2 and G BP > 18 according to Evans et al. (2018); Arenou et al. (2018). For the astrometric information (Lindegren et al. 2018), we use the equation 1 of Arenou et al. (2018), we correct the parallax zero point of −0.03 mas and we do not use the parallax when + 3 × σ < 0 to remove spurious astrometric solutions.

Extinction map with FEDReD
To analyse this dataset, we used the FEDReD algorithm. The entire algorithm description as well as tests on mock and observed data are presented in Babusiaux et al. (2020). Thus we will just briefly explain the main steps of the algorithm and differences from Hottier et al. (2020).

Field of view analysis
FEDReD analyses photometric and astrometric data, field of view by field of view, and infers both the evolution of the extinction as a function of distance and the stellar density distribution. To do so, it works in two steps.
The first step is the processing for each star of the likelihood of this observed star being at the distance D with extinction A 0 (extinction at 550nm) P (O | A 0 , D). This distribution is processed by comparing the apparent photometry of the stars to an empirical HR diagram built from 2MASS and Gaia DR2 (see Babusiaux et al. 2020, section 3.1 for details), it also uses parallax information when it is available.
Once the likelihood of each star computed, FEDReD applies a Bayesian deconvolution to obtain the joint distribution of extinction and distances P (A 0 , D). This deconvolution also takes into account the completeness of the field of view, which is estimated from the observed near-infrared photometry distribution.
To initialise this deconvolution, we use two simple priors. The prior on the distance distribution is a square law of the distance, corresponding to the cone effect. Concerning the extinction given the distances, P (A 0 | D), we use an uniform distribution (contrary to Hottier et al. 2020).
From the joint distribution of extinction and distance, FE-DReD generates Monte-Carlo Solutions (MCSs) of the increasing relation A 0 (D) drawn following the probability distribution P (A 0 , D). As Red Clump stars are the ones providing the strongest constrains on the distance/extinction, we restricted our results to the distance interval [D min , D max ] where red clump stars are observed by 2MASS, that is the distance at which a red clump star saturates (D min ) or is fainter than the completeness limit (D max ). Unlike Hottier et al. (2020), this distance interval restriction is also used within the algorithm determining the A 0 (D) and not just in the post-processing. At the end of this process we obtain 1000 MCSs by field of view.

Merging fields of view into extinction cubes
The fields of view are 0.38 • wide in longitude and latitude, and the centers are spaced by 0.38 • in latitude and longitude. Therefore each field of view overlaps its first neighbour (top, bottom, left and right) by half of its angular surface and overlaps its second neighbour (four diagonals) by the quarter of its angular surface. This means that each star is inside three different field of view, which allows us a good continuity between fields of view.
To merge results from each field of view into a consistent extinction cube, we use the exact same algorithm as in Hottier et al. (2020). Firstly, we iteratively clean the MCS samples of each field of view by using neighbour fields' MCS envelopes as upper and lower limits, the convergence of this process being provided by the overlapping of fields. Secondly, we randomly draw one MCS per field of view (inside their respective clean pool) to obtain a relation between the extinction A 0 and the distance D for each field of view. Then, at each distance bin we smooth the extinction value using the eight neighbour fields of view, to obtain the extinction cube. We randomly draw 100 of these cubes. Finally, we use a constrained cubic spline fit (Ng & Maechler 2007) to obtain the median relation of extinction as a function of distance of each field of view. These relations are then decumulated and normalised by the distance width of bins to obtain the extinction density a 0 in each voxel of the cube.

Extinction uncertainty
As in Hottier et al. (2020), we also compute the uncertainty of our extinction and extinction density cubes. Concerning the extinction A 0 we use the sample of MCSs after the cleaning process by the neighbour envelope to get the maximum and minimum values of extinction at each distance bin (A 0min (D) and A 0max (D)).
To estimate the extinction density uncertainty, we use the exact same algorithm explained in section 3.2 but on a bootstrapped MCS sample instead. We build 100 bootstrap merged cubes and we compute the standard deviation of extinction density a 0 for each voxel to get the uncertainty on the extinction density. As discussed Hottier et al. (2020) (section 5.1), this un-certainty map mostly represents the sampling error, and it underestimates the true uncertainty of our results.

Clump Extraction
To obtain distance, shape and size of extinction density clumps, we interpolate our data cube on a regular Galactic Cartesian grid 1 with a voxel size of 10 pc.
Using this new interpolated cube, we first tried to process the iso-surface density to locate clumps. This allowed us to spatially constrain some clumps but it required very sensitive settings for each clump we looked for. Moreover, this technique has difficulty resolving distinct clumps, and local peaks can hide the true shape and size of a clump.
Dendograms have already been used to extract clumps and molecular clouds (Goodman et al. 2009;Chen et al. 2019). This algorithm uses hierarchical tree to segregate structures, showing relations between each clump. To do so it only uses the local value of pixels and links voxels with their neighbours if values are compatibles. But as we work with a density cube processed by line of sight, the extinction density can leak to larger distances and create some fingers of god structures.
To avoid the extraction of structures only dominated by fingers of god effect, we use the FellWalker algorithm 2 (Berry 2015) which identifies clumps by analysing the local gradient. In a few words, this algorithm "walks" through the entire volume, voxel by voxel. To choose the next voxel it always chooses the steepest path. When it reaches a local maximum, it checks nearby for a higher point (in order to avoid extrema due to noise). If this higher point exists the algorithm "jumps" and so on. If a higher spot does not exist, the algorithm has just reached the top of a new clump. When the entire volume is mapped the algorithm check if clumps can be merged. Finally, it can reassign a voxel to the most common clump around it.
As the finger of god effect also impact the gradient of extinction density, the FellWalker algorithm is also affected. While a minimum value based algorithm (such as dendogram) just wipe out low extinction structures to avoid elongations, the FellWalker algorithm ability of setting up both a minimal value and a minimum gradient can mitigate the extraction of finger of god effect when looking for the edge of clumps, without missing the ones with low values of extinction. Nevertheless, the clumps segregation is depending on the parameter settings of the algorithm.
With this in mind, we manually tune the parameters of the algorithm in order to avoid elongations and extract parts of the VMR. The minimum extinction density to be a part of a clump is 3 mag kpc −1 . The size of the maximum jump is set to 4 voxels which corresponds to 0.04 kpc. The minimum initial density variation is 0.5 mag kpc −1 from one voxel to another. Two clumps may be merged if the depth of the col between them is smaller than 3 mag kpc −1 . The clump outline process is done with a cube of 4 voxel. Finally, we removed clumps which contain only one or two voxels as those are just high extinction density spots induced by the decumulation process.

Processing clump parameters
The result of the FellWalker algorithm is a 3D mask for each clump. From this mask we directly obtain angular and distance bound of clumps as well as the volumes of the clumps. To obtain the center of a clump, we compute the barycentric position using extinction density as weight.
In order to estimate the total dust content A c of a clump we have to integrate the extinction density a 0 over the volume of the clump: if we discretise the equation we get: with i the i th voxel in the clumps, l i , b i and d i the position of the voxel, a 0i the extinction density in this voxel and δd i the distance width of the voxel. As these voxels follow the cube we obtained by merging the fields of view, δl i and δb i are constant over the cube δl = δb = 0.18 • . In practice, as we looked for clumps in a regular cartesian cube, we used the Eq 1 in cartesian referential frame with dV = (0.01kpc) 3 and a 0 the value of extinction density in each voxel. Note that these A c can be used to estimate the masses of clumps, following the method described in Chen et al. (2020). Some of the clumps are located on the edge of our extinction cube. We mark them as "not full" and their estimated extinction is only a lower estimate while their distance only corresponds to what we observe in our extinction cube.

Uncertainties on the parameters
To estimate the uncertainty on each parameter, we used the random cubes that we have built to compute the density uncertainty. We apply on them the same FellWalker algorithm with the same parameters. We cross identify clumps from their barycentric position and we estimate each parameter uncertainty using their standard deviation.
We kept, in the final catalog, clumps which are detected in at least 99% of the bootstraps. Indeed, some large clumps are identified as several small ones in some bootstraps, preventing the usage of a 100% threshold.

Extracting cavities
We also used the FellWalker algorithm to spatially constrain density cavities. As this algorithm is designed to search peaks rather than valleys, we inverted the values of the density cube following this equation : max(a 0 (l, b, d))−a 0 (l, b, d). The cavities in the "inverted" data cube appear more as a plateau than as a mountain so we did not use the minimum steepness parameter. The minimum value parameter is set in order to map the area of the original cube with extinction less than 0.5 mag. As for the clump extraction we set the maximum jump to 4 voxels.
Inverting the data cube implies that every area without extinction is detected as a cavity, so areas above and below the complex are considered as many little cavities which merge into a single one during the merging process. Moreover, the shell of cavities that we can see by eye are porous, which leads them to be merged with the clean empty part below and above the complex. To avoid this phenomenon, we clean our cavity sample before the merging step. To do so we remove from the sample the cavities with a median voxel values inferior to 0.1 mag.kpc −1 . Then we merge the cavities with a col height threshold at 0.03 mag.kpc −1 .
Finally, we remove cavities on the edge of our field of view. Indeed, in the case of clumps, if a structure is detected on the edge it is for sure a clump. But in the case of cavity detection, we are looking for area with low extinction density. In the case of cavity inside our data cube we can ensure that this area is really surrounded by an extinction shell, while a cavity on the edge of the cube can be completely open.
As for clump, we also processed uncertainty on cavities parameters, but due to the large volume of the cavities and the thin shell around them, the tuning of the parameters is very sensitive. This implies that we do not recover cavities in bootstrap as often as we did for clumps. So, in order to provide uncertainties, we lower the threshold at 50% to keep cavities in final catalog.

Extinction column
We present in Fig. 2 a view of the cumulative extinction A 0 of the Vela complex up to 5 kpc, roughly the distance limit of Vela complex noticed in Hottier et al. (2020), with the same boundaries as in Fig. 1. We also overplotted the main structures of Vela.
The VMR appears as the strongest extinction structure of the entire complex and is approximately located at 260 • ≤ ≤ 272 • and −4 • ≤ b ≤ 3 • . While the overall shape of the VMR is quite similar to the IRAS view shown in Fig. 1, it differs in the details, as the strongest emissions at 60 µm do not coincide with the areas with the strongest cumulative extinction. As in Fig. 1, the IVS is very faint but using a visual guide, it is still visible as a diffuse filament. On the other hand, due to the relatively small area covered by our study, the Gum nebula footprint cannot be distinguished.

Extinction density towards Vela
In Fig. 3 we present the extinction density a 0 cut along constant Galactic latitudes. Fig. 4 shows views of the extinction density at several distances. We also add on those two figures, the contours of clumps that we segmented with the FellWalker algorithm.
On both figures, we do not provide the extinction density close to the Sun because of the red clump visibility criterion (see Sect 3.1). To see the closeby Gum Nebula and the VSNR, we used instead the Lallement et al. (2019) extinction cube. They used the same data as this work (2MASS and Gaia DR2) restricted to stars with relative error on the parallax less than 20%, and they apply the same extinction law (see Appendix. A for a more detailed comparison). In Fig. 5 we plot the extinction density that they obtain at the canonical distance of the Gum Nebula and VSNR. Indeed, at D = 0.25 kpc we saw a shell which has the size of the VSNR. We also see a large shell which seems to correspond to the Gum Nebula at D = 0.25 kpc whereas we cannot distinguish any footprint at D = 0.45 kpc. This could mean that the Gum Nebula is closer than usually indicated in the literature.

VMR and Clumps
The FellWalker algorithm found 14 clumps. Table 1 presents their parameters. We identify them by a simple integer id sorted by their volumes. We visually identify clumps to the VMR parts (column Notes of Table 1). The VMR C (clump 3) goes beyond the edge of our extinction cube but we can see that we detect the main part of it. It begins at D = 0.5 kpc and its barycenter is at D = 0.9 ± 0.09 kpc, which is consistent with the Liseau et al. (1992) estimation at 0.7 ± 0.2 kpc. This distance is also inside the confidence interval of Massi et al. (2019)(0.950 ± 0.050 kpc) and our distance interval contains every molecular distances of Zucker et al. (2020) labeled VMRC, which range between 0.86 and 0.97 kpc. However VMR C is very large, its longitude footprint goes from = 256 • to more than = 280 • and its maximum distance reaches D = 1.6 kpc.
VMR D (clump 1), the biggest of this study, is farther than VMR C. Its barycenter is at D = 1.9 ± 0.4 kpc and its front part is at D = 0.75 kpc, while Liseau et al. (1992) only see the front part of it and estimate VMR C,D and A to be roughtly at the same distance of 0.7 ± 0.2 kpc.
Concerning VMR A and VMR B, they are identified as a unique clump (clump 2) 3 . The reason why this clump is split into two pieces in the literature, is the overlapping with VMR C in the foreground which induces a projection effect. Nevertheless, VMR A is one of the biggest clump of the Vela region and it seems to be related to VMR C and VMR D at high distances (D > 1.2kpc) where they form a shell (see Sect 5.4). We can also note than this clump contains the molecular cloud RCW 38 which is located at D = 1.6 kpc by Zucker et al. (2020).
As far as we know, the other clumps have not yet been identified and named. Clump 8 is a large clump, and it is bigger than what we process as it goes beyond the edge of our extinction cube. Clump 6 presents a relatively high extinction density and reasonable apparent angle but it is barely visible in the column view (Fig. 2) because it is overlapped by VMR D.

Cavities
The FellWalker algorithm found 9 cavities. As for clumps, we identify them by a simple integer ID sorted by volume. We represent the contour of this cavities in distance views on Figure 6 and their parameters are presented in Table 2.
The biggest one (cavity 1) is centered on ( , b) = (266.1 • , −3.1 • ). It is very visible on IRAS data (Fig. 1) and on the cumulated extinction view (Fig. 2). Its angular location fits the VSNR angular location but the center of this cavity is at Fig. 3. Extinction density a 0 at different Galactic latitudes b. x and y represent galactic coordinates, the Sun is at (0, 0) and the Galactic centre direction is to the right. We also add the contour of the clumps, their parameters being presented in Table 1. D = 3.1 ± 0.35 kpc and its closest "wall" is at D = 1.2 kpc, behind the foreground part of VMR C (clump 3) and its shell is composed by all the VMR parts. On the other hand, VSNR is at 0.25 kpc with an upper limit at 0.49 kpc (Cha et al. 1999), in the foreground of the VMR. We saw its local extinction fingerprint in Fig 5. So, it appears that this shell, described in the 2D views as VSNR, visible in infrared and in extinction column, is actually composed of two distinct structures physically separated by VMR C.
As Cavity 1, Cavities 2 and 3 are under the main extinction density of the galactic plane and their low latitudes shells are very thin. While the others are located inside the galactic plane. So far, we have not found the origins of these cavities.

Conclusion
We used the FEDReD algorithm to analyse photometry and parallax from 2MASS and Gaia DR2 of stars in the Vela complex direction. This provided us with the distribution of the extinction density in 3D. Consequently we were able to obtain the distance and the 3D shape of known structures of this area.
For this purpose we applied the FellWalker algorithm, which produced very good results when applied to extracting clumps. However it turns out not to be the best designed tool to extract cavities as it requires a very specific tuning. For further inves-tigation on cavities, we recommend using or developing other approaches, although we have not found a better one so far.
Nevertheless, we managed to spatially measure 14dense extinction clouds (Table 1) and 9 cavities (Table 2). We then visually identified these clouds using the projections of the known structures. Doing this, we determined the distance of Vela Molecular Ridge components. It appears that the split of VMR into 4 parts described in the literature is more a construction due to a face-on view of very large structures, where the foreground is composed by VMR C and VMR D. At high longitude (i.e. the western side) the so called VMR A and VMR B are actually one clump only and prolongate a large shell completed with the background part of VMR C and VMR D. Moreover, it came out that what seems to be the infrared shell of the VSNR is actually mainly due to a background cavity surrounded by the components of the VMR.
This study of the Vela complex in extinction reveals that this area contains structures of different types and in interaction. In fact this complex is probably a part of the local arm (Hou & Han 2014;Hottier et al. 2020;Khoperskov et al. 2020). Therefore Vela is a perfect lab to study a spiral arm from the inside.