The outer disc in shambles: blind detection of Monoceros and ACS with Gaia's astrometric sample

The astrometric sample of Gaia allows us to study the outermost Galactic disc, the halo and their interface. It is precisely at the very edge of the disc where the effects of external perturbations are expected to be the most noticeable. Our goal is to detect the kinematic substructure present in the halo and at the edge of the Milky Way (MW) disc, and provide observational constraints on their phase-space distribution. We download, one HEALpix at a time, the proper motion histogram of distant stars, to which we apply a Wavelet Transformation to reveal the significant overdensities. We then analyse the large coherent structures that appear in the sky. We reveal a sharp yet complex anticentre dominated by Monoceros (MNC) and the Anticentre Stream (ACS) in the north, which we find with an intensity comparable to the Magellanic clouds and the Sagittarius stream, and by MNC south and TriAnd at negative latitudes. Our method allows us to perform a morphological analysis of MNC and ACS, both spanning more than 100$^\circ$ in longitude, and to provide a high purity sample of giants with which we track MNC down to latitudes as low as $\sim$5$^\circ$. Their colour-magnitude diagram is consistent with extended structures at a distance of $\sim$10-11 kpc originated in the disc, with a very low ratio of RR Lyrae over M giants, and kinematics compatible with the rotation curve at those distances or only slightly slower. We present a precise characterisation of MNC and ACS, two previously known structures that our method reveals naturally, allowing us to detect them without limiting ourselves to a particular stellar type and, for the first time, using only kinematics. Our results allow future studies to model their chemo-dynamics and evolution, thus constraining some of the most influential processes that shaped the MW.


Introduction
Most of the studies that discovered new substructures within the second data release (DR2, Gaia Collaboration et al. 2018b) of the Gaia mission (Gaia Collaboration et al. 2016) used the full 6D phase-space sample, and their impact on our current understanding of the Milky Way (MW) and its history is undeniable. Good examples of that are the work of Belokurov et al. (2018), Gaia Collaboration et al. (2018a, Haywood et al. (2018) and Helmi et al. (2018), which identified a large group of stars accreted in the last major merger event of the MW that took place ∼10 Gyr ago. Another example is the advance in the study of the moving groups and the possibility to now visualise the kinematic substructure directly in the plane of Galactocentric radii against rotational velocity with the ridges Kawata et al. 2018;Ramos et al. 2018;Laporte et al. 2019b;Fragkoudi et al. 2019;Khanna et al. 2019). Nevertheless, this sample is limited to G< 13 mag approximately (above that, the completeness drops significantly), restricting the explotation of the kinematic data to a volume of ∼3 kpc radius from the Sun. Despite some attempts to extend the kinematic maps to further distances either by usemail: p.ramos@unistra.fr ing statistical corrections to the parallax (López-Corredoira & Sylos Labini 2019; López-Corredoira et al. 2020), or by adding photometric (Anders et al. 2019) or spectroscopic information (Liu et al. 2017;Wang et al. 2019), these only have a significant amount of stars up to Galactocentric radii of ∼16 kpc.
In contrast, the 5D sample (only astrometry and no radial velocity) is more than two orders of magnitude larger than the 6D one. Its power is exemplified by the work of, for instance, Castro-Ginard et al. (2018 by discovering hundreds of new open clusters throughout the disc, Malhan et al. (2018) and Ibata et al. (2019) by revealing several new tidal streams in the halo, or by the large sample of halo stars selected using a combination of photometry and proper motions provided in Koppelman & Helmi (2020). Another good example is the detection of the Sagittarius (Sgr, Ibata et al. 1994) stream using mostly its kinematic signature Ibata et al. 2020;Ramos et al. 2020).
The astrometric sample reaches down to G ∼21 mag, expanding the volume probed significantly 1 , meaning that we can use it to trace kinematic structures well into the halo. Among the different stellar systems that we expect to find within the 5D sample we count globular clusters (e.g., Baumgardt et al. 2019), streams (e.g., Belokurov et al. 2006), dwarf galaxies and even ultra faint dwarf galaxies (e.g., Willman et al. 2005;Belokurov et al. 2007;Koposov et al. 2015). Moreover, this sample also covers the outermost regions of the MW disc, where Newberg et al. (2002) reported, almost two decades ago, the presence of a peculiar population above the mid-plane of the Galaxy, bluer than the thick disc and clearly appreciable as an overdensity of Main Sequence Turn-off (MSTO) stars at a distance of ∼10 kpc from the Sun. Known as Monoceros (MNC, also referred to as Galactic anticentre stellar structure or GASS), this structure was observed to span more than a hundred degrees in longitude (100 • < l < 270 • , see, e.g., Rocha-Pinto et al. 2004;Morganson et al. 2016) both in the north and south hemispheres. During the past two decades, there has been an intense debate over its origin, in part due to the difficulties of confronting the data with the different models available (Slater et al. 2014). Out of the many possible mechanisms proposed by Ibata et al. (2003), there have been mainly two leading hypotheses: accretion and disc perturbation.
The idea that MNC is the tidal debris of an accreted satellite was based on its morphology (it looks like a stream) and on its metallicity and kinematics (e.g., Yanny et al. 2003;Crane et al. 2003;Wilhelm et al. 2005;Conn et al. 2005), partially supported by the simulations of Helmi et al. (2003) and Peñarrubia et al. (2005). This lead to the hunt for its progenitor and, after discarding the Canis Major (Martin et al. 2004) over-density as a candidate (e.g., Momany et al. 2006;Rocha-Pinto et al. 2006;Carballo-Bello et al. 2020), none has yet been found, despite the attempts to detect the continuation of the hypothetical tidal stream at other Galactic latitudes (Conn et al. 2007(Conn et al. , 2008. Nevertheless, upper limits for the total mass of the progenitor have been calculated (e.g., Guglielmo et al. 2018).
On the other hand, several works have shown that the close passage of a satellite can induce significant substructure in the outer parts of the disc (e.g., Younger et al. 2008;Purcell et al. 2011;Gómez et al. 2016). The interaction with a dwarf galaxy as massive as Sgr could cause some of the disc material to move to more inclined and eccentric orbits, and produce a stream of stars consistent with the observations (see also Kazantzidis et al. 2008, where, instead of a single satellite, the perturbers are 6 dark matter subhalos of masses ∼10 10 M ). More recently, the simulations by Laporte et al. (2018Laporte et al. ( , 2019a have shown that it is possible to create extended structures similar to MNC as well as feather-like structures during a satellite encounter, while at the same time reproducing qualitatively part of the phase-space substructure observed at the solar neighbourhood. The detection of a vertical wave-like pattern in the disc (Widrow et al. 2012) that propagates almost radially Schönrich & Dehnen 2018), in agreement with the simulations of Gómez et al. (2013), and the discovery of the phase-space spiral and consequent confirmation that our Galaxy is undergoing phase-mixing ) further supports the perturbative scenario.
Other structures in the outer disc are the Anticentre Stream (ACS) and Eastern Band Structure (EBS, Grillmair 2006), or the Triangulum-Andromeda (TriAnd1 and TriAnd2) overdensities Rocha-Pinto et al. 2004;Martin et al. 2007). The connection between all of these and MNC has also been subject to scrutiny for many years and is still not entirely clear (but see the models of Xia et al. 2015;Sheffield et al. 2018). For instance EBS, which was described as an independent structure by Grillmair (2011), has now been suggested to be part of the MNC ring by Deason et al. (2018) using a combination of Gaia and SDSS (York et al. 2000) data. de Boer et al. (2018 re-analysed the kinematics of MNC and ACS with SDSS astrometry calibrated with Gaia DR1 to provide accurate kinematic maps, showing that they have similar yet clearly distinct kinematic trends that can be used to establish the processes that form them. Very recently, Laporte et al. (2020) studied the [Mg/Fe]-[Fe/H] distribution of ACS and MNC with a combination of Gaia DR2 data and LAMOST-SEGUE-APOGEE to, guided also by their previous simulations of an isolated MW interacting with Sgr, intercede in favour of a disc origin for the two structures. In their work, they also show with colour-magnitude diagrams (CMDs) that both have a conspicuous red clump (RC), in contrast with previous studies that mainly focused on the Main Sequence (MS), the MSTO or the 2MASS M-giants. Here we aim to provide an independent detection and characterisation of these structures that can help us clarify their true extent and 3D morphology, as well as their nature. A deeper understanding of the events that lead to the observed stellar distribution in the anticentre could be used to constrain the orbit and mass of Sgr, as well as its effect on the gas and stars of our Galaxy.
In this work, we search for substructure following the strategy devised by Antoja et al. (2015b). The original goal of this method was to detect ultra faint dwarf galaxies in the halo using the fact that these should create, simultaneously, an overdensity in proper motion space and in the sky. Here we use the first half of the methodology, that is, its application in proper motion space only, to find the kinematic substructure at large heliocentric distances. This approach allows us to scan the whole Celestial sphere systematically, homogeneously, and using a statically robust technique that can distinguish small but significant overdensities in proper motion space, and track their changes as we move with Galactic longitude and latitude. As a result, our all-sky maps are dominated by three large structures: the Magellanic clouds, the Sgr stream (as reported in Antoja et al. 2020), and MNC-ACS. Our goal is to map and study the morphology and kinematics of the structures in the Galactic anticentre, taking advantage of our methodology which allows us detect them and obtain a large set of members with (almost) no prior information.
This paper is organised as follows. In Sect. 2 we describe the strategy used to process the large amount of data available with Gaia. Section 3 then enumerates the different systems detected with our method and, in Sect. 4, we focus on characterising the complex kinematics of the anticentre, specially in the north where we observe MNC and ACS. We discuss the implications of our findings in Sect. 5, and finally present our conclusions in Sect. 6.

Data and methods
In this work we use the same sample and methodology presented in Antoja et al. (2020, hereafter A20), which we reproduce here for convenience. We exploit the full Gaia DR2 (Gaia Collaboration et al. 2018b), not restricting ourselves to any magnitude limit other than the one intrinsic to the instruments, and applying only the following two filters: where LVL is the level of the HEALpix grid (here, 5), HPNUM is the HEALpix to be processed, and BINSIZE is the size of the histogram binning (here 0.24 mas yr −1 ). Then, we apply the Wavelet Transformation (WT, Starck & Murtagh 2002), followed by a peak detection algorithm, to each of the 12 288 histograms we have downloaded to detect the significant kinematic structures (assuming Poisson noise). The resulting wavelet coefficient of each peak is then, by construction, proportional to its density in proper motion space. To simplify the analysis, we only keep one structure at each HEALpix, the one with the highest relative intensity (WT /N hp ): 2 Hosted at: https://gea.esac.esa.int/archive/ where N hp is the total number of sources in the HEALpix and is used to normalise the wavelet coefficient (for more details on the method, see A20). Figure 1 shows the Mollweide projection of the sky in Galactic coordinates coloured by the relative intensity (Eq. 2) of the highest peak in the proper motion histogram. By selecting only the overdensity with the largest intensity present at each proper motion histogram, we can focus on the dominant kinematic structure of the HEALpix 3 . The normalisation used in Eq. 2 compensates the density gradient of the Galaxy and gives more contrast to the structures at higher latitudes. This figure reveals a wealth of substructure that cannot be seen with a simple density map of our sample. For instance, Fig. 2 contains the number of sources that pass our filters at each HEALpix (top panel) where we can only identify the Magellanic clouds, some globular clusters and the imprints of the extinction (which is shown in the middle panel for comparison, Schlegel et al. 1998) or the Gaia scanning law (bottom panel). With Fig. 1 we have been able to detect:

Global map of the substructures
-The Magellanic clouds: Their angular size in the kinematic maps, in contrast with their apparent angular size in the star counts map of Fig. 2 (top panel), is larger and shows the true extent of these systems as already noted in, e.g. Gaia Collaboration et al. (2018c). We also detect substructure within them (not shown here), and we are able to recover some of the globular clusters 4 that orbit the Large Magellanic cloud, like the recently detected Gaia 3 (Torrealba et al. 2019). -The Sgr stream: The core of this dwarf galaxy and its stream are also clearly visible in our map. Interestingly enough, we do not observe it in the top panel of Fig. 2, which highlights the difficulty, even in the Gaia era, to detect this structure with just stellar counts. -Almost vertically mirrored to the Sgr stream, we note a feature with a stream-like shape. The CMDs of the sources that produce these peaks do not present any coherent isochronelike shape. Instead, they appear clumped around the faint and blue corner of the diagram. Added to the fact that their proper motions are always normally distributed around the origin, regardless of the position in the sky, this leads us to conclude that these sources are actually quasars. Although quasars are ubiquitous and should not produce a band in the sky, the scanning law of Gaia favours certain regions of the sky as can be seen in the bottom panel of Fig. 2 with the astrometric_gof_al 5 that quantifies the quality of the astrometric solution. In these parts, the astrometric uncertainties of the quasars are low enough so that they become more relevant than the diffuse halo population of stars. Since we did not expect to obtain such a clear signal from the quasars, we did not remove them beforehand. Nevertheless, our results are not affected by their presence. With the list of extended objects that will be published in Gaia DR3, we will be able to remove these objects up-front already within the queries. -Nearby galaxies: Apart from the Magellanic clouds, we also detect M31 and M33, the latter clearly visible in Fig. 1 at (l, b) ∼ (133 • , -31 • ). -Dwarf spheroidals: Our method is able to detect several dwarf spheroidals, like Fornax, Sextans or Sculptor, as well as fainter ones like Draco. -Globular clusters: We recover 51 globular clusters from the Bica et al. (2019) catalogue. From their 200 objects classified as globular clusters, more than half are in the bulge where we do not detect any either because the contrast with the foreground is too low or due to the cut in parallax applied. -Ultra Faint Dwarf galaxies: We do not recover any of the known Ultra Faint Dwarfs galaxies based solely on their kinematic signature. Their proper motion uncertainties are too large and there are too few members (see Massari & Helmi 2018) to produce a significant overdensity. Nevertheless, we note that when we apply the full methodology described in Antoja et al. (2015b), that includes the search of peaks in the sky and not only in proper motion as done here, we can effectively recover most of them. -Anti-centre: Apart from all the substructure we find in the halo, our methodology reveals complex kinematic substructures towards the anticentre of the MW, dominated by two arch-like features in the north Galactic hemisphere. After comparing with the extinction map shown in the middle panel of Fig. 2, we confirm that these features are not aligned with regions of high absorption. Also, we have checked how the astrometric_gof_al map (bottom panel of Fig. 2) superposes to the WT intensity map, from which we conclude that the shape of the bottom arch is artificially enhanced by the scanning law. The cavity at ∼180 • , b ∼ 20 • coincides with a region poorly sampled by Gaia and therefore the intensity (proportional to stellar counts) is lower.
We note that some of these structures also appear in Fig. 3, where we colour the sky according to the proper motion of the highest peak (the peak used to colour Fig. 1). The most clear one is the Sgr stream, which we have analysed in detail in A20. Also, the structure at latitude b∼35 • (140 • < l < 200 • ) appears as a conspicuous arch in the proper motion map. We devote the following sections to the analysis and characterisation of the kinematic substructure present at the outer disc, focusing mostly in the north where we observe these two conspicuous arches already mentioned.

Monoceros and ACS
In Fig. 4 we present a zoom-in of Fig. 1 towards the anticentre and show our selection of the two structures that appear after colouring the sky according to the relative intensity of the dominant structure in proper motion. To build this selection in an objective manner, we first apply a Gaussian softening (two sigmas) of the 2D image to erase the HEALpix limits, and then apply a bi-directional Sobel filter 6 to reveal edges. By doing so, the two arches are cleanly separated at all longitudes. The final step is to select only the HEALpix whose Sobel intensity is above a certain threshold (0.0035 for the bottom arch, 0.0040 for the top one). However, if we applied this selection blindly we would obtain a long list of HEALpix that comprises several structures. Instead, we first draw a rectangle around each arch in an appropriate coordinate system. This coordinate system, different for the bottom and top arch, is obtained by rotating the Celestial sphere with respect to the Galactic reference frame until the structure lies roughly flat at zero latitude in the new reference frame. The resulting final selections are the contours shown in 6 Included in the Python package Scikit-image (van der Walt S et al. 2014) the bottom plot of Fig. 4. The structures can be seen to continue beyond the contours defined, especially for the feature at lower latitudes, but we focus on the regions where they are the most intense. We add three horizontal lines that represent an approximated latitude limit of each structure at ell ∼180 • . By comparing the shape and location in the sky of these structures we note that they match with the MNC ring (bottom) and ACS (top) (e.g., Newberg et al. 2002;Grillmair 2006;Slater et al. 2014;Morganson et al. 2016). The patches we obtain are also in good agreement with the regions delineated by Laporte et al. (2020) but are much more concise. In contrast to previous works, though, since we are not relying on counts but instead we detect these structures in relative intensity (Eq. 2), their morphology appears sharper and well defined. For instance, we observe a MNC structure that has a clear arch like shape 7 extending from ∼120 • to ∼230 • in longitude, where it meets the disc at a latitude of ∼10 • . Nevertheless, we stress again that in the case of MNC, although the structure is physical, the sobel filter is enhancing the edge caused by the scanning law of Gaia. In contrast, the ACS is thinner, stays above MNC for the ranges of longitudes where we detect it and has its strongest signal at l ∼140 • , where MNC has already almost merged with the disc. This is the most precise picture of the anticentre available to date, thanks to the introduction of kinematic information in the detection of these structures.
Once we have identified the MNC and ACS regions, we explore their kinematics and CMDs in more detail in Fig. 5. Each row contains the results for different regions: the first row is MNC, the second corresponds to the list of HEALpix that fall between MNC and ACS (hereafter, the bridge), then the third is ACS and, finally, the fourth is a region above ACS. The last row is an example of what we would expect from a Galaxy with no substructure, obtained from a mock catalogue (Appendix B). With this exercise we can evaluate the continuity of these structures, compare their characteristics with nearby regions in the sky where we do not see an enhancement in relative intensity (c.f. Fig. 4) and contrast them with the predictions of a MW model. In the first column of this plot we have aggregated all the stars that, in their respective HEALpix, fall within the highest intensity proper motion peak. We refer to such stars as peak stars. To show what we would see if we had not done this kinematic selection, the grey contours on top represent the CMD of all the stars of the region. The second and third columns contain, respectively, the trends of µ l and µ b with Galactic longitude for all the stars within the region. Here the black line encircles the stars selected kinematically, that is, the outer-contour of the volume that the peak stars occupy in this space. To provide some contrast with a fiduciary galaxy, in the bottom row we repeat the same process for the the particles in the mock catalogue that fall within the ACS footprint. In this case, the peak stars are selected according to the position and size of the peaks detected in the data, which is why the contours of panels h (i) and n (o) are so similar.
The first thing that we note is the presence of a Giant branch all the way from MNC to ACS, that disappears once we explore latitudes larger than ∼40 • . The fact that we see a well defined RC means that these stars share a similar distance which, based on their magnitudes (∼15.5 mag) and Galactic latitudes (b > 15 • ), puts them well above the mid-plane of the Galaxy (z > 2 kpc) at a height larger than the scale height of even the thick disc.  If we compare the observed CMDs with the CMD of the mock catalogue, we note that we do not expect many stars in this region of the diagram as the nearby giants have been already removed with the cut in parallax (these are bright enough to have a reliable parallax) and the farther ones are not in the model since there are not many stars at such heights/distances.
We also note, accompanying the Giant branch, a dense clump of stars appearing in panels a-d-g which is bluer than the MS of the disc seen in the mock. Newberg et al. (2002) already reported that the main sequence turn off of MNC was bluer than the thick disc and we detect the same behaviour for the stars in the peaks found within MNC and ACS. This group of stars is consistent with being the MS of an isochrone containing the RC discussed above. The rest of the stars that fall outside said isocrhone seem to follow the contours of the CMD obtained with all the stars (no kinematic selection, grey contours in those panels), and are most likely nearby, dwarf field stars that overlap with these structures in the proper motion plane. In contrast to previous works (e.g. Newberg et al. 2002;Ivezić et al. 2008;Xu et al. 2015;Thomas et al. 2019), where they detect an overdensity in counts for a given population, MS, MSTO, blue stragglers or M-giants, here we have unveiled the whole sequence in the CMD by performing a blind kinematic selection of the stars instead.
In the bottom row of Fig. 5 we have included a few curves that represent the proper motion that Gaia would measure for a star at given distance if it only had azimuthal velocity 8 . In orange (red), this distance is 10 kpc (4 kpc) and in solid (dashed) the rotation velocity is 200km s −1 (220km s −1 ). A structure that is too near, like the brown curve, does not match the contour delineated by the peak stars (black contours), while a structure that does not rotate, as could be the halo, would have to be at a distance larger than 50 kpc in average to fall within the black lines. And even then, its shape would not be compatible with the data. While we note that other combinations of distances and velocities could produce a similar shape (even if the result is not physically supported), the dashed red line shows a good agreement with our observations and corresponds to a structure at ∼10 kpc P. Ramos et al.: Monoceros, ACS and other outer disc structures with Gaia astrometry The proper motion maps contain all the stars in the region, and the black line is the zero-contour of the peak stars, that is, the stars that fall within the highest intensity proper motion peak of their HEALpix. First row: MNC. Second row: Bridge between MNC and ACS. Third row: ACS. Forth row: above ACS. Fifth row: Same region as ACS but for the mock particles, selecting the stars for the CMD according to the contours of panels h and i. In the bottom panel we also include the proper motions expected from a structure 30 • above the plane at a given distance, 4 (red) or 10 (orange) kpc, rotating at a given velocity, 200 (solid) or 220 (dashed) km s −1 , but with no radial or vertical velocity.
rotating slightly slower than the disc. In other words, the peak stars in ACS (the same applies to MNC) have proper motions that change with Galactic longitude in a way that is compatible with a structure at a distance of ∼10 kpc rotating at a speed similar to the disc or slower, in agreement with the analysis by de Boer et al. (2018). By comparing the CMDs inside and outside the patches defined in Fig. 4, it is clear that, even though the MS of these structures is the dominant fraction, by focusing only on the giants we can gain contrast with the MW foreground. Therefore, we  Fig 4). A sudden increase of the ratio can be seen in the part where ACS is the more intense.
introduce another tag, apart from the one that we have already been using to separate stars inside and outside the proper motion peaks. The stars will be called giants whenever they are redder than G BP − G RP > 1 mag and their apparent magnitude smaller (brighter) than the line: where we have used the slope calculated in Romero- Gómez et al. (2019) to follow the extinction vector, and the zero-point is adjusted by eye to reduce the contamination from the disc while preserving the RC as much as possible. We note, however, that we are still selecting some faint, red dwarfs at all latitudes, the great majority of which are not classified as peak stars (for sources redder than 2 mag in G BP − G RP and G>14 mag, only ∼100 out of ∼8000 in MNC and ∼30 out of ∼4000 in ACS). This means that the giants tagged also as peak stars are more likely to be true giants, whereas field stars tagged as giants have a larger probability of being nearby red dwarfs. Figure 6 shows the fraction of giants inside the peaks (i.e., stars tagged as giants and peak stars simultaneously) with respect to all the stars tagged as giants as a function of latitude. Given that our classification is rather rough, we should treat these as simple estimates and focus on the trends instead. What we observe is that the parts where MNC and ACS have the strongest signal in relative intensity (c.f. Fig. 4) coincide with the regions where this ratio is the highest. We have already seen that the relative intensity of ACS decreases with with longitude, and here we note that the ratio of giants also diminishes moving from one curve to the other. We also observe that the bridge keeps a constant ratio, showing that it is just the region where the tail of the two structures overlap. Finally, we note that our patch around ACS is too broad, as the ratio drops abruptly at b ∼37 • , coinciding with the place where we observe a discontinuity in the kinematics 9 .
Tables 1 and 2, available online, contain the list of sources classified simultaneously as peak stars and as giants for, respectively, MNC (10 079 sources) and ACS (2 104 sources).
A&A proofs: manuscript no. aanda Table 1. MNC stars classified both as peak and giants (top 2 rows). The first column contains the source_id followed, in columns two and three, by the right ascension and declination of the star. Columns four to seven contain the proper motions (ICRS) and the corresponding uncertainties. The eight and ninth are the apparent magnitude in the G band and the Gaia colour G BP − G RP , respectively. Then, the Galactic coordinates, and b are given in columns twelve and thirteen. Finally, in the last column, we give the absorption in the G band (see text).

Anticentre region: North vs South
Above, we have proven that our methodology can detect kinematic substructures and isolate them effectively from the rest of the disc. Also, that the giants inside the proper motion peaks are good tracers of MNC and ACS since the contamination in that region of the CMD is expected to be very low once we have filtered by parallax and kinematics. Therefore, we can use the location of the RC to trace the structures in physical space. To do so, we now explore the changes in apparent magnitude of the stars that we have tagged as peak stars, for different ranges of longitude, both in the north and south hemispheres. Figure 7 shows the distribution of apparent magnitudes of peak stars with respect to Galactic latitude using a histogram normalised by bins of b, for the range 130 • < l < 150 • . This figure confronts the data (top) with the expectations from the mock (bottom). In this case, the particles selected in the mock correspond to the peaks that are detected in the mock itself (in contrast to panel m of Fig. 5 where we used the peaks detected in the data). We find an overdensity of stars in the north at a magnitude ∼16 that corresponds to the RC seen in Fig. 5. It is most intense above b >30 • and corresponds to the ACS. We see it extending rather continuously down to b ∼10 • where it merges with the disc, following an arch that is compatible with the increase in extinction.
In the south, we observe an excess of bright stars (G <17 mag) at latitudes between 15 • and 25 • with respect to the mock. Interestingly enough, the intensity maps ( Fig. 1) do not show an enhancement as is the case for the north, not even compared to the mock map (Fig. B.1). Based on their apparent magnitudes and location in the sky, it is very likely that these stars form the diffuse stellar population detected by Ibata et al. (2003) and that is sometimes called MNC South. This is also the region where the TriAnd overdensities have been reported Rocha-Pinto et al. 2004;Martin et al. 2007) and, as we show below, we do detect it with our method. A detailed study of these structures and a comparison with previous studies (e.g., Fig. 4 from Perottoni et al. 2018) is out of the scope of this work, but we will revisit it once the eDR3 Gaia data release (Brown 2019) is made public and we have more and better astrometric data.
In Fig. 8 we now present the difference between the data and the mock in the same plane of apparent magnitude against Galactic latitude, for different bins in longitude. From right to left, these are: 130 • < l < 150 • , 150 • < l < 170 • , 170 • < l < 190 • , and 190 • < l < 210 • . In panels (d) and (h) we see the subtraction of the left panels of Fig. 7 from the right panels (i.e., data minus mock after properly normalising the histograms). As already mentioned, we see two distinct overdensities in the south, a bright (14< G <15.5 mag) one corresponding to MNC south and a fainter one (G∼17 mag) corresponding to TriAnd, at a magnitude consistent with the most recent determinations of its heliocentric distance (e.g., Bergemann et al. 2018). These features cannot be seen so clearly at other longitudes except for panel (e) where the diffuse overdensity appears again (latitudes between ∼20 • and ∼30 • , brighter than ∼17 mag). While we cannot dis- card that the southern structures are indeed discontinuous, taking into account the corrugations in the disc reported by Xu et al. (2015) the most likely scenario is that we can only detect them with our method where the extinction is low enough. With the next Gaia releases we will be able to assess better their continuity.
In the north (top panels), we note that ACS decreases its intensity and shifts to lower latitudes when we move towards the third quadrant of the Galaxy, as we see also in Fig. 4. In the intermediate panels (b and c) we observe two concentrations at different latitudes, but similar apparent magnitudes, corresponding to MNC and ACS, whose tails overlap forming the bridge that we mentioned above. More importantly, we see MNC extending more and more towards lower latitudes, keeping roughly the same apparent magnitude throughout. This is interesting as MNC is usually hard to trace so deep into the disc due to the foreground stars. And yet, in panel (a), using our kinematic selection its RC can be traced down to a latitude ∼5 • .
By measuring the median G for the Giant stars only, selected according to Eq. 3, we can investigate the relative distance of these structures. To do so, however, since the latitudes probed by MNC and ACS change with longitude, we focus only in the range 130 • < l < 170 • where they remain rather flat and the bridge is quite wide (see Fig. 4). The median of the G magnitude and the associated one sigma interval of uncertainity at different bins in latitude is shown in Fig. 9 for the Giant stars inside the proper motion peaks. To compensate for the effects of extinction, we have first corrected the apparent magnitude using the G BP − G RP colour of each source individually and the prescription detailed in Appendix A of Ramos et al. (2020). What we observe is that, below b ∼28 • where we identify MNC, the RC is brighter than above b ∼31 • where ACS begins. Since we have used the integrated extinction up to infinity (Schlegel et al. 1998), the separation that we observe is an upper limit: if we assume that the extinction applied to ACS is correct, since it is at a  Fig. 9. Apparent magnitude of giants peak stars in the anticentre region (130 • < l < 170 • ) as a function of Galactic latitude, after correcting for extinction (see text). The error bars denote the 1σ uncertainty on the median computed as σ π 2N , where σ is the standard deviation of the apparent magnitude in the bin. Vertical lines represent the approximate limits of each structure in that range of Galactic longitudes (see Fig. 4), and the right axis represents the distance to a RC star with the apparent magnitude shown in the left axis. The horizontal lines correspond to the median G magnitude for the giants peak stars within MNC (cyan) and ACS (orange). The shaded areas contain the ±3σ interval of uncertainty on the median and they extend from the minimum to the maximum latitude of the peak stars within each patch (the vertical dashed lines serve only as an orientation). As can be seen, the ACS is fainter than MNC and this translates to a difference in distance of ∼1 kpc. higher latitude, then MNC could actually be less extincted than assumed, and therefore be intrinsically fainter than the value we are recovering. Nevertheless, since both structures are quite far, this effect should be small and we can safely conclude that ACS is farther away than MNC.
To be more quantitative, we measure now the median G, corrected for extinction, for all the peak giant stars in this longitude range. The result is the shaded areas shown in Fig. 9 where we can clearly see that ACS is, once we convert the differ-ence in magnitude to distance, roughly 1 kpc farther away than MNC, with a discrepancy of more than 3σ. Also, we have estimated the median distance to each of the two structures. In doing so, we assume that the median apparent magnitude measured (horizontal lines in Fig. 9) corresponds to the magnitude of the RC. By imposing that the absolute magnitude of a RC stars is M G =0.495 mag (Ruiz-Dern et al. 2018), we obtain the following median distances and their statistical uncertainties 10 : D MNC ∼ 10.6±0.1 kpc, and D ACS ∼ 11.7±0.2 kpc.
Nevertheless, without a precise calibration of each individual star, and its extinction, we cannot investigate the changes in distance with longitude and latitude which is key to reveal the 3D shape of these structures. We make a first attempt to study the distribution of the structures we detect along the line of sight by cross-matching the peak stars with StarHorse (Anders et al. 2019), a catalogue of Bayesian derived astrophysical parameters obtained from the photometry of Gaia , Pan-STARRS1, 2MASS, and AllWISE combined. We download all the stars in the anticentre (100 • < l < 260 • and -60 • < b < 60 • ) with SH_OUTFLAG equal to "00000", as recommended in (Anders et al. 2019), and with a distance (50th percentile) less than 20 kpc. From the 13 098 038 peak stars in the north, we find 429 565 in StarHorse. In the south, the cross-match returns 514 167 stars out of the 13 669 647 peak stars. Most of them, however, are faint dwarfs found at low latitudes, closer than 10 kpc, whose parallax quality is not good enough to discard them with the filter presented in Sect. 2. If instead we restrict ourselves to the giant peak stars, then we find 378 955 (out of 1 449 250) in the north and 458 403 (out of 1 286 132) in the south. Figure 10 shows the distribution of StarHorse distances as a function of Galactic latitude for the four ranges of longitude explored above. The first thing we note is the effect of our selection function as the nearby giants are missing and a wall of stars at a distance of ∼6 kpc is formed. The tails extend up to ∼15 kpc, point beyond which StarHorse distance uncertainties become too large. Compared with the corresponding figure for the mock (Fig. B.3), where we see the disc extending much farther away, we note a clear excess of stars in the data at latitudes larger than 20 • and at a distance of >7 kpc. We associate these to MNC and ACS. In the right column, in the range of longitudes where ACS is more intense, we see it clearly separated from MNC and slightly farther away. As we shift our view towards the third Galactic quadrant, these structures recede, becoming less prominent and shifting to lower latitudes (as we showed above). MNC covers a large range of latitudes and connects smoothly with the disc, but the lack of stars and the uncertainties prevent us from determining if there is a distance gradient with latitude or not.
The south does not show the same structures as the north but we note that, at least for panels (f) and (g), the extinction is higher than in the north, which could block our line of sight. We note that we do not recover so clearly the structures detected in panel (h) of Fig. 7, probably due to low statistics. However, we do observe an increase of stars in panel e in the form of a diffuse distribution of distant stars at latitudes between -20 • and -30 • , coinciding with the location of the structure S200-24-19.8 10 Here, the dominant source of uncertainty is the systematic errors, which are not included in the error bars given. The more important ones are i) the assumption that the median magnitude is the magnitude of the RC, and ii) not using 3D extinction maps but instead correcting with the integrated extinction to infinity. Other sources of systematic uncertainty are the error on the absolute magnitude of the RC, contamination from stars that are not giants, or errors in the extinction map. reported in Newberg et al. (2002) and also the detection by Xu et al. (2015) that they associate with the aforementioned TriAnd.

RR Lyrae to M giant ratio
With a kinematically selected sample of stars for both MNC and ACS, we can now check if their population is consistent with having been born in an extragalatic system or not. All known MW dSphs have a large fraction of RR Lyrae stars (Vivas & Zinn 2006) and a low one of M giants (Price-Whelan et al. 2015), whereas the opposite happens with the Galactic discs. Hence, we will follow Price-Whelan et al. (2015), where they used the ratio between the number of RR Lyrae and M giant stars in TriAnd Martin et al. 2007) to argue that these structures were probably disc stars kicked-out by an external perturbation. For that, we need to estimate the number of RR Lyrae and M giants within MNC and ACS. Based on previous studies, we expect a low number of RR Lyrae in these structures (Kinman et al. 2004) but a large number of M-giants (e.g., Rocha-Pinto et al. 2004). Now the question is how low is the ratio between the two populations.
Starting with the RR Lyrae, we use the catalogue described in Mateu et al. (2020) which combines the VariClassifier and Specific Objects Studies (SOS) catalogues from Gaia DR2 (Holl et al. 2018;Clementini et al. 2019) with the ASAS-SN-II catalogue Jayasinghe et al. (2019), providing optimal completeness at the bright end (G < 15). We select only those sources that fall within the sky patches defined in Fig. 4, for a total of 900 ( 800) RR Lyrae in MNC (ACS), of which 253 stars fall within the one standard deviation range around the median distance in the case of MNC (6.7 to 16.7 kpc) and another 253 for ACS (7.8 to 17.7 kpc). Of these, only 12 (6) are also consistent with the kinematic signature we have detected, in the space of l-b-pmrapmdec. Based on these results, and after correcting for the completeness of the RR Lyrae catalogue (estimated at 80% at these magnitudes in Mateu et al. 2020), we conclude that, at most, MNC has 15 RR Lyrae and ACS, no more than 8.In parallel, we cross-match the list of RR Lyrae with our sample to see the effect that the cut in parallax 11 has and we observe that, from the 253 (253) stars that we had, only 123 (123) remain in MNC (ACS). Finally, if we now keep only those classified as peak stars (i.e., probable MNC/ACS members) we find only 1 and 2 RR Lyrae for, respectively, MNC and ACS. These figures are lower and, even if we correct by the 48.6% reduction in completeness caused by the cut in parallax, the maximum amount of RR Lyrae in MNC (ACS) would be 3 (5).
On the other hand, the M giants are much more numerous. We use the official Gaia cross-match with 2MASS and the selection proposed by Majewski et al. (2003) to obtain the corresponding G BP − G RP colour cut necessary for MNC and ACS, finding that we can select confidently M giants in our sample of candidates with the following limits: G < 15.5mag and G BP − G RP > 1.5mag. These cuts result in a total of 959 M giants for MNC and 155 for ACS. Of course, these values are not corrected for the completeness of the sample, as opposed to those corresponding to the RR Lyrae. Therefore, if we take the maximum number of RR Lyrae that these structures can have, the ratio f RR:MG that we provide becomes a strict upper-limit: f RR:MG < 1.5% for MNC, f RR:MG < 5.2% for ACS.
The fractions obtained are consistent with previous independent estimates (Sheffield et al. 2018)  much lower than the expected values for extragalactic systems like Sagittarius or the LMC, for which we expect a fraction ∼50% (Price- Whelan et al. 2015). In fact, it is hard to reconcile these values with the hypothesis that MNC and ACS are tidal tails of an accreted satellite. It would require either a young system that managed to reach high metallicities (perhaps formed from already metal-enriched gas), or the effect of a dynamical mechanism that could segregate RR Lyrae from giants within the stream, something that we do not observe in other streams like the Sagittarius stream (e.g., Antoja et al. 2020;Ibata et al. 2020;Ramos et al. 2020).

Discussion
Although recent work favours a disc origin for these structures, the debate over its origin (the alternative being that these are the tidal debris of an accreted MW satellite) is still on going. Based solely on the morphology that our method allows to observe, both MNC and ACS could very well be different wraps of the same tidal stream. If that was the case then we should be able to see at least a hint of continuity in the south, unless the tails only emerge from behind the disc at Galactic longitudes where the disc is already too dense for our method to detect them. Peñarrubia et al. (2005) presented an N-body model fitted to the observations of MNC available at the time, which was later used by Slater et al. (2014) to show that there is a broad agreement with the PanSTARRS-1 (Chambers et al. 2016) data. We find that the arch described by the debris generated with their model is too wide to explain MNC as we detected it 12 but we note the presence of a tail (top right panel of Fig. 5 from Slater et al. 2014) that resembles ACS. Comparing the morphology now with the N-body simulation of Laporte et al. (2019a) 13 , where these 12 Yet, the shape that we recover could be affected by the Gaia scanning law, as we mentioned in 4.1 13 We tried to compare also with the simulations of Kazantzidis et al. (2008, Fig. 5) and Gómez et al. (2016, Fig. 4) but there are too few "feathers" appear as a result of the interaction with Sgr, we note that the overall agreement is good (top left panel of their Fig. 1). For instance, the difference in latitude between the two structures and the fact that the one on top stops abruptly at a given longitude. Nonetheless, the "feather" corresponding to ACS obtained with their N-body simulations is thicker than the observed one and does not present a higher density of stars close to the turning point of the vertical oscillation (the highest point in latitude) as we see in the data. However, this could simply be due to resolution limitations and the fact that these simulations were not meant to be an exact match to the MW.
To explore a little bit deeper the simulations of Laporte et al. (2018), in Fig. 11 we show the plane of Galactic latitude against Heliocentric distance of the particles that are at Galactocentric radius > 18 kpc. To do that, we locate the Sun at (x,y)=(-8, 0) kpc and look at the particles with 120 • < < 240 • . The snapshots are chosen such that we can see the outer disc at the beginning of the simulation (first row), right after the first pericentric passage of Sgr (second), at the time of first apocentre (third), right after the second pericentre (fourth) and, finally, at the second apocentre (fifth). The first thing we note is that, while the first pericentre passage seems to have a similar effect in both simulations, causing some material to be ejected from the system (third row), the difference becomes more noticeable with the second passage. In particular, we can see how in the H2 simulation a couple of "feathers" appear in the fourth row, with roughly constant declination (∼20 • and ∼40 • ) and a large extension in distance. These correspond to the MNC-and ACS-like structures reported in Laporte et al. (2019a). The shape that we observe for MNC in the Gaia data seems more concentrated in distance while being more extended in latitude (see Fig. 7). Nevertheless, a more detailed analysis of the distances in our data, as mentioned in previous sections, is still necessary to produce a more precise 3D characterisation of MNC and ACS. particles to make a good assessment given the level of detail that the data is providing. Secondly, we note that, even at the timescales of ∼100 Myr, the large-scale distribution of stars in the whole outer disc of the simulated MW changes significantly with time in both the radial and vertical directions. In turn, this would suggest that the phase-space configuration of the outer disc and of its feathers at present time can in principle pose strong constraints on the timeevolution of the perturbation (assuming that they were caused by a single perturber). The data that we have obtained in this work, being precise, continuous and covering a large range of longitudes and latitudes, can be used for describing MNC and ACS with analytical/semi-analytical models like, for instance, with the analytic method presented in Weinberg (1998). This efficient way of exploring the parameter space could be used to obtain a quantitative measure of the goodness of fit for each model, beyond the small set of simulations analysed here. Such an exercise is crucial as it would allow us to predict where the continuation of these structures should appear, both in Galactic coordinates as well as in distance, providing a way to search for them actively. More importantly, we could quantify the mass of the perturber and the history of its orbit. Or even, as mentioned in Laporte et al. (2019a), use MNC and ACS to constrain the Galactic potential: the rotation curve at that distance, its slope, the shape of the Dark Matter halo, etc.
Another option to compare models objectively is to add particles generated with N-body models, one that could reproduce the observed MNC and ACS, to a mock catalogue of the Galaxy. In this work we have used the Rybizki et al. (2018) catalogue as an example of Galaxy without substructure, but we now know that it was not so representative of the MW and, also, that it underestimated the observational errors (Rybizki et al. 2020). The recent Rybizki et al. (2020) catalogue has fixed much of these issues and provides a mock to use with the up-coming Gaia EDR3 (Brown 2019). With these approach, and taking into account the nuances of our methodology, we can attempt a quantitative comparison with the models. One of the key parameters to generate the N-body models would be the stellar mass contained inside these structures, which is currently poorly measured (Morganson et al. 2016). In this work, we have attempted to quantify the fraction of the disc that is within MNC and ACS in Fig. 6, but it is just a rough estimate based solely on the number of giants.
We have also studied the kinematic information obtained with the proper motion of the peaks, and we see that MNC and ACS rotate slightly slower than the disc at the solar position, in agreement with de Boer et al. (2018). Nevertheless, this alone does not prove that these structures were once part of the disc. Since the two mechanisms proposed (extra-galactic/internal) have distinct formation time-scales and chemical properties, by exploring also their CMDs we can more easily distinguish between both. In this sense, the ratio of RR Lyrae to M giants that we find (<5%) is unlikely for structures composed of tidal debris from an accreted satellite. This adds to the recent studies of Bergemann et al. (2018) and Laporte et al. (2020), where it is shown that the abundances, the distribution in the [Mg/Fe] vs [Fe/H] plane and the mean metallicity of MNC and ACS are inconsistent with the extragalatic scenario.
If indeed MNC and ACS were once part of the disc, then we can now use them as chemical fossils, an idea already exposed in Laporte et al. (2020). After they were kicked out of the disc, the stellar formation of these structures most likely came to a halt as any gas that initially accompanied the stars must have quickly settled back to the disc thanks to its efficient energy dissipation mechanisms. As a result, their current population is a frozen relic of the outskirts of the MW at the time when the perturbation occurred. With enough spectroscopic abundances we could learn about the gas that dwelled at the edge of our Galaxy some Gigayears ago and use that information to constrain the chemo-dynamical models of the MW.

Conclusions
The application of the WT to the proper motion space has proven extremely useful to reveal the kinematic substructure of the halo and outer disc. By removing most of the foreground with a simple yet effective cut in parallax, our method is able to efficiently detect kinematic substructure in the halo and even external galaxies like M33 or the Magellanic clouds, several dwarf spheroidals and dozens of globular clusters, as well as the Sgr stream. It has also revealed the sharpest picture of the anticentre, with MNC and ACS appearing as the third most prominent structures in the distant sky (only after the Magellanic clouds and Sgr).
We have been able to blindly detect the whole MNC north as well as the ACS from l ∼ 120 • to l ∼ 230 • . Our findings are in good agreement with previous studies like Laporte et al. (2020), who also used DR2 data to investigate these structures. Nevertheless, we have been able to characterise their morphology with great detail, which is crucial for obtaining the orbital parameters of these groups of stars. We observe MNC with an arch-like shape, broader at small longitudes and becoming thinner towards larger longitudes. Nevertheless, the RC stars that we have selected can be seen to span a wider range of latitudes, therefore a detailed study of the selection function of Gaia and the extinction is needed to confirm how much of this shape is cause by the scanning law. ACS can be seen at larger latitudes than MNC throughout the whole longitude range where we detect them, and has a maximum of relative intensity when it reaches the highest latitude at l∼140 • (consistent with a pile up of stars at the maximum height in the orbit), and stops abruptly at a longitudes of ∼110 • . This behaviour, added to the fact that we do not observe a clear continuity in the south, favours the perturbative scenario proposed by Ibata et al. (2003), later supported by the simulations of many other authors (e.g., Gómez et al. 2016;Laporte et al. 2019a). Moreover, the kinematics of these features, which differ from the bulk motion of the disc stars that lay in front, are compatible with a low eccentricity orbit at ∼10 kpc that rotates similarly to the disc.
By analysing the apparent magnitude of the RC stars selected by proper motion we have been able to trace MNC down to a latitude of ∼5 • , closer to the disc than ever before. Also, by measuring the median apparent magnitude of the RC stars of each structure and converting to Heliocentric distance, we have determined that ACS (∼11.7 kpc) is roughly 1 kpc farther away from the Sun than MNC (∼10.6 kpc). This actually means that both structures are at roughly the same Galactocentric radius (but at heights above the disc of, respectively, ∼6.5 kpc and ∼4.5 kpc). Also, we have shown that MNC and ACS, despite being different structures, are extended in distance and in the sky, and their tails overlap both in the 3D physical space as well as in kinematic space.
In the south, we have found a diffuse population of giants at 130 • < l < 150 • and 190 • < l < 210 • , coinciding with the regions of low extinction, that we do not observe in the mock catalogue nor in the north. Their apparent magnitudes span the range 14 < G < 15.5 mag which implies distances for a RC star not affected by extinction between 5 and 10 kpc. These could be related with the vertical wave described in Xu et al. (2015), and are most likely the so-called MNC south first reported by Ibata et al. (2003). On the other hand, in the longitude range 130 • < l < 150 • we have observed a faint trace of RC at G ∼17 mag (heliocentric distance ∼16 kpc) that most likely corresponds to the TriAnd overdensity. Nevertheless, due to the contamination of nearby stars close to the disc and large the distance uncertainties, we have not been able to explore the morphological connection between these structure and MNC north.
Studies like this will benefit the most from the next Gaia release (EDR3) which is expected to contain proper motions twice as precise (on average), and increase the number of stars at the faintest magnitudes (Brown 2019). As a result, the structures that we detect will become more concentrated in the proper motion space, and in turn produce stronger signals in our maps of relative intensity. Also, the effects of the scanning law should diminish as a result of the additional year of observations. With it, we should be able to remove the foreground contamination more efficiently and detect the anticentre structures continuously at all latitudes, providing a direct observation of their 3D morphology. Moreover, the WEAVE spectroscopic survey will also observe this region, so we could potentially obtain radial velocities and abundances for a large fraction of the giants in our sample. This means that we will be able to trace MNC and ACS, probably even TriAnd, more clearly and deeper, and obtain a less contaminated sample of members.
The challenge now is to find the way to use our data to form a coherent and unified picture of the outer disc, constraining the properties of the different agents involved (MW disc, dark matter halo, Sgr), as well as to prepare these methods to work with the future samples that will soon become available. In some of the figures, e.g. Figs 5 or B.3, we select the particles not according to the peak obtained with the mock but with coordinates of the peaks detected in the data. In doing so we can check what is distribution in the CMD or in distance of the particles that have the observed kinematics.