© P. A. B. Galli et al. 2020

1. Introduction

The Lupus molecular cloud complex is one of the closest and largest low-mass star-forming regions of the southern sky. It consists of several subgroups associated with different molecular clouds labelled as Lupus 1–9 that exhibit distinct morphologies as revealed by CO surveys and extinction maps (Hara et al. 1999; Cambrésy 1999; Tachihara et al. 2001; Dobashi et al. 2005). Although they are part of the same complex, Lupus clouds differ significantly regarding their star formation activity. Some clouds have dense concentrations of T Tauri stars (e.g. Lupus 3) while others (e.g. Lupus 7–9) display no signs of ongoing star formation (Comerón 2008).

Early surveys leading to the identification of Lupus stars mostly focused on the molecular clouds with active star formation (i.e. Lupus 1–4), but more recently a few young stellar objects (YSOs) have also been discovered in Lupus 5 and 6 (Manara et al. 2018; Melton 2020). Most of the hitherto known classical T Tauri stars (CTTSs) in Lupus were identified by Schwartz (1977) and Hughes et al. (1994) based on their strong Hα and infrared excess emissions. Succeeding studies identified a more dispersed and older population of weak-emission line T Tauri stars (WTTSs) from ROSAT X-ray pointed observations surrounding the Lupus clouds, which greatly exceed the number of CTTSs located in the immediate vicinity of the molecular clouds (Krautter et al. 1997; Wichmann et al. 1997a,b). The census of Lupus stars was later expanded by Merín et al. (2008) based on infrared observations collected with the Spitzer Space Telescope as part of the cores to disk (c2d) legacy project. In a subsequent study, Comerón et al. (2009) revealed a new population of cool stars and brown dwarfs potentially associated with Lupus based on their spectral energy distributions (SEDs) derived from optical and infrared photometry. López Martí et al. (2011) confirmed most of these sources to be genuine members of the Lupus region based on their proper motions and showed the importance of using kinematic information (e.g. proper motions, parallaxes, and radial velocities) to assess membership.

The distance to Lupus has undergone several revisions and has become a matter of debate in the recent decades although it has always been recognised as one of the closest star-forming regions to the Sun. Distance determinations in the literature range from 100 pc (Knude & Hog 1998) to 360 pc (Knude & Nielsen 2001). Franco (1990) and Hughes et al. (1993) estimated the distances of 165 ± 15 pc and 140 ± 20 pc, respectively, based on the interstellar reddening of field stars. Lombardi et al. (2008) used a more robust method based on 2MASS wide field extinction maps and derived the distance of 155 ± 8 pc with a depth of ${51}_{-35}^{+61}$ pc, which the authors explained as the Lupus clouds being at different distances. It is important to mention that these distances are only average estimates based on indirect methods because distances for individual stars derived from trigonometric parallaxes existed (until recently) for only a few stars in this region. For example, Bertout et al. (1999) used trigonometric parallaxes of five stars from the HIPPARCOS catalogue (ESA 1997) to estimate the distance to Lupus. Three of them, which are associated with Lupus 1, 2 and 4 have a mean parallax of ϖ = 6.79 ± 1.50 mas that roughly defines a distance of ${147}_{-27}^{+42}$ pc. The other two stars, located in Lupus 3, have a mean parallax of ϖ = 4.38 ± 0.67 mas, which puts this cloud much farther away at a distance of ${228}_{-30}^{+41}$ pc. Alternatively, kinematic distances to individual stars in this region were derived in the past from the convergent point method under the assumption that the Lupus stars are comoving (Makarov 2007; Galli et al. 2013). The resulting distances confirmed the important depth effects (i.e. distance variations along the line of sight) reported in previous studies.

Lupus is located between the Upper Scorpius (US) and Upper Centaurus-Lupus (UCL) subgroups of the Scorpius-Centaurus (Sco-Cen) association, but despite the close proximity in the sky it represents a more recent star formation episode in this region (Preibisch et al. 2008). In this context, the existence of a more dispersed and older population of WTTSs near the Lupus clouds as reported in previous studies is not surprising, but their association with the younger population concentrated in the immediate vicinity of the molecular clouds needs to be confirmed. The second data release of the Gaia space mission (Gaia-DR2, Gaia Collaboration 2018) provides the best astrometry available to date to further assess the membership of the historically known members. When combined with radial velocity information these data can be used to constrain not only the 3D positions of the stars, but also their 3D space motions. The stellar content of the Sco-Cen association was recently revised by Damiani et al. (2019) based on Gaia-DR2 data, but a dedicated study of the Lupus clouds is still missing in the literature. The proper motions and parallaxes in the Gaia-DR2 catalogue are more precise by one or two orders of magnitude compared to the ground-based astrometry used in previous studies of the Lupus region (see e.g. Makarov 2007; López Martí et al. 2011; Galli et al. 2013), and they will therefore allow us to obtain a clean sample of members and a more accurate picture of the overall 6D structure of the complex.

This paper is one in a series as part of the Dynamical Analysis of Nearby Clusters project (DANCe, Bouy et al. 2013). Here, we investigate the census of stars and kinematic properties of the Lupus association of young stars based on Gaia-DR2 data. It is organised as follows. In Sect. 2 we compile the list of Lupus stars from previous studies in the literature that is used in Sect. 3 to perform a new membership analysis based on a novel methodology developed by our team (Sarro et al. 2014; Olivares et al. 2019). In Sect. 4 we revisit several properties of the Lupus subgroups (e.g. distance, spatial velocity, age, spatial distribution, and initial mass function) based on our new sample of cluster members. Finally, we summarise our conclusions in Sect. 5.

2. Sample of Lupus stars from previous studies

We construct our initial sample of stars based on the lists of members and candidate members that have been associated with the Lupus star-forming region in previous studies. We take the samples of (i) 73 stars from Hughes et al. (1994), (ii) 136 stars from Krautter et al. (1997), (iii) 92 stars from Wichmann et al. (1997a), (iv) 48 stars from Wichmann et al. (1997b), (v) 159 stars from Merín et al. (2008), (vi) 248 stars from Tables 4 and 9 of Comerón et al. (2009), (vii) 82 stars from Tables 2 and 4 of López Martí et al. (2011), and (viii) 69 stars from Damiani et al. (2019). This results in a sample of 508 stars after removing the sources in common among the several lists of Lupus stars. We find proper motions and parallaxes for 441 stars of this sample in the Gaia-DR2 catalogue. Figure 1 shows the location of the stars in this sample with respect to the main star-forming clouds of the Lupus complex. We use the boundaries defined in Fig. 3 of Hara et al. (1999) to assign the sources in the Lupus sample to the corresponding molecular cloud of the complex based on their sky positions. This procedure does not confirm membership to the corresponding clouds, but allows us to compare the stellar proper motions of the stars projected against the various clouds. As discussed in Sect. 1, it is apparent that the stars identified by Krautter et al. (1997) and Wichmann et al. (1997a,b) constitute a more dispersed “off-cloud” population that surrounds the molecular clouds.

Fig. 1. Location of the Lupus stars identified in previous studies overlaid on the extinction map of Dobashi et al. (2005) in Galactic coordinates. Colours and symbols denote the samples of stars from each study. The red dashed lines indicate the position of the main molecular clouds defined by Hara et al. (1999).

Similarly, we compile a list of known members of the adjacent Sco-Cen association that will be used in the forthcoming discussion. We restrict this sample to the to stars in the UCL subgroup that are located in the general region of the main star-forming clouds of the Lupus complex (334° < l <  342° and 5° < b <  25°). Our UCL sample is based on the lists of sources previously identified in this region by de Zeeuw et al. (1999), Preibisch et al. (2008), Pecaut & Mamajek (2016) and Damiani et al. (2019). It includes1352 stars with measured proper motions and parallaxes in the Gaia-DR2 catalogue.

Figure 2 shows the distribution of proper motions for the Lupus and UCL stars compiled from the literature. This preliminary analysis reveals the existence of (at least) three populations of stars in the Lupus sample from the literature. First, there is a background population of sources projected towards Lupus 1, 3, and 4 with proper motions typically smaller than 10 mas yr⁻¹ which are more likely to be unrelated to the Lupus star-forming region as anticipated in past studies (López Martí et al. 2011; Galli et al. 2013). Second, there is a more diffuse population of stars with proper motions that overlap with the UCL stars in this region of the sky. Indeed, many of them have been identified as Sco-Cen members in the recent study conducted by Damiani et al. (2019). Finally, we note an overdensity of stars around (μ_α cos δ, μ_δ)≃(−11, −23) mas yr⁻¹ (see highlighted region in Fig. 2) including stars projected towards the four molecular clouds of the complex (Lupus 1–4) and with proper motions consistent with membership in Lupus (see e.g. Galli et al. 2013). The latter will be used in our forthcoming analysis (see Sect. 3) to search for additional members to the Lupus association and investigate the kinematic properties of this star-forming region.

Fig. 2. Proper motion of Lupus and UCL stars identified in previous studies. The different colours denote the subgroups of the Lupus sample. The cyan rectangle indicates the locus of Lupus candidate stars in the space of proper motions (see discussion in Sect. 2).

3. Membership analysis

The membership analysis performed in this study for the Lupus star-forming region is based on the methodology developed by Sarro et al. (2014) and Olivares et al. (2019). The main steps of our analysis are briefly summarised in this section and we refer the reader to the original papers for further details.

First, we define the representation space of our membership analysis (i.e. the set of observables) that we use to classify the sources as cluster members or field stars. The representation space includes both astrometric and photometric features of the stars given in the Gaia-DR2 catalogue. It does not include sky positions since the Lupus association is spread over a relatively large sky region without any clear over-density with respect to the field population in terms of stellar positions (see Fig. 1). In our case the inclusion of stellar positions in the representation dilutes the discriminant power of parallaxes, proper motions and photometry in the membership analysis, thus resulting in higher contamination rates. Furthermore, we do not use the blue photometry (BP) from Gaia-DR2 due to some inconsistencies reported in the literature in the BP system (see e.g. Maíz Apellániz & Weiler 2018). We measured the relative importance of the photometric features using a random-forest classifier and concluded that G_RP is a more discriminant magnitude than G for identifying members. So, our representation space for the membership analysis is defined by μ_α cos δ, μ_δ, ϖ, G_RP, and G − G_RP.

The region of the sky downloaded from the Gaia-DR2 catalogue to perform our membership analysis is defined by 334° ≤ l ≤ 342° and 5° ≤ b ≤ 25° to include the molecular clouds of the complex with active star formation (Lupus 1–4). This input catalogue contains 14 471 847 sources (10 470 632 of them with complete data in the chosen representation space) whose membership status with respect to the Lupus region will be investigated in our analysis.

We modelled the field population using Gaussian Mixture Models (GMM)¹. The field model was computed only once at the beginning of our membership analysis and held fixed during the process. We tested field models with 20, 40, 60, 80, 100, 120, 140 and 160 components, and inferred their parameters using a random sample of 10⁶ sources from our input catalogue. We compute the Bayesian information criteria (BIC) for each one of the models and take the one (with 80 components) that returns the minimum BIC value as the optimum model for our analysis. The main objective our membership analysis is to distinguish the Lupus association from the rest of the sources in the region surveyed by our study which we collectively denote as field population. It includes background sources, foreground stars and stellar groups at similar distances of Lupus, but with different properties (e.g. the Sco-Cen subgroups). The number of GMM components used to model the field are valid in a statistical sense and no attempt has been made to assign a physical interpretation or characterise the various components of the field (as done e.g for the cluster itself in the forthcoming sections).

The cluster is modelled with GMM in the astrometric space and a principal curve in the photometric space. The principal curve defines the cluster isochrone with a spread at any point that is given by a multivariate Gaussian. Both the principal curve and its spread are initialised with an input list of cluster members. The cluster model is built iteratively based on the initial list of cluster members that is provided to the algorithm in the first iteration of the process. The initial list of cluster members that we use consists of a sample of 88 stars compiled from the literature and located within the proper motion locus indicated in Fig. 2. This sample contains sources projected against the four main molecular clouds of the complex (Lupus 1–4) that are more likely to be associated with the Lupus region than to the Sco-Cen association (see Sect. 2). It contains 13 stars in Lupus 1, three stars in Lupus 2, 47 stars in Lupus 3, nine stars in Lupus 4, and 16 off-cloud stars. This first list of “candidate” members is indeed needed to start the membership analysis, but it can be incomplete and somewhat contaminated since it will be updated in the following iterations. Its main purpose is to locate the cluster locus in the space of astrometric observables and to define the empirical isochrone of the cluster in the photometric space to start the membership analysis. The algorithm uses the input list of members to infer the astrometric and photometric parameters of the cluster model. Then, the method assigns membership probabilities to the sources and classifies them into cluster and field stars based on a probability threshold p_in previously defined by the user. At each step of the algorithm the marginal (i.e. data independent) class probabilities are estimated as the fraction of sources in each category obtained in the previous iteration. The resulting list of cluster members is used as input for the next iteration and this procedure is repeated until the list of cluster members remains fixed after successive iterations.

Once that our solution has converged we generate synthetic data (based on the properties of the cluster and field models) and evaluate the performance of our classifier to define an optimum probability threshold p_opt (see Olivares et al. 2019). Finally, we classify the sources in our catalogue used for the membership analysis as members (p ≥ p_opt) and non-members.

Table 1 shows the results of our membership analysis using different values for the user-defined probability threshold p_in. The true positive rate (TPR, i.e. fraction of cluster members generated from synthetic data that are recovered by the algorithm) and contamination rate (CR, i.e. fraction of field stars generated from synthetic data that are classified as cluster members by the algorithm) given in this table allow us to evaluate the quality of our membership analysis. However, we caution the reader that these figures derived from synthetic data are not absolute numbers but only estimates that are valid in the absence of the true distributions (and given the assumed cluster model). The solution obtained with p_in = 0.5 exhibits the highest CR and contains fewer cluster members compared to the solutions derived with p_in = 0.6 and p_in = 0.7. The latter is explained from the more conservative p_opt threshold that results due to the low initial threshold. On the other hand, the solutions derived from p_in = 0.8 and p_in = 0.9 have lower CRs (and higher TPRs), but a significantly reduced sample of cluster members. The solution obtained with p_in = 0.6 offers the best compromise between the number of clusters members in the sample and the performance of our classifier (we note that the TPR and CR are compatible with the solutions obtained with the highest p_in values). We have therefore adopted this solution (with 137 stars) as our final list of cluster members for the present study of the Lupus region. Table A.1 lists the 137 stars identified in our membership analysis together with other results that will be discussed in the forthcoming sections of this paper. We also provide in Table A.2 the membership probabilities for all sources in the field using the p_in values investigated in our analysis, so that the reader can select cluster members based on other criteria that are more specific to his/her scientific objectives.

We tested the dependency of our membership analysis on the initial list of members that is provided to the algorithm in the first iteration. We selected the stars for the initial list using a somewhat narrower range of proper motions that concentrates most of the Lupus 3 candidate stars compiled from previous studies (−13.0 <  μ_α cos δ <  −8.0 and −25.0 <  μ_δ <  −21.0) and repeated the membership analysis as described before. We concluded based on similar arguments presented for the first run that the solution obtained with p_in = 0.6 returns the best compromise between the sample size (141 members) and performance of our classifier (TPR = 0.91 ± 0.02 and CR = 0.13 ± 0.04) as compared to the other solutions with different p_in thresholds. These numbers are consistent with the values obtained in the first run (see Table 1) and small variations are expected to occur given that our classifier exhibits a non-zero CR. We recovered 95% of the stars obtained in the first run and confirm that our result for the membership analysis is not dependent on the initial list of stars that we use (as long as they define the general properties of the Lupus association).

Figures 3 and 4 show the proper motions, parallaxes and the colour-magnitude diagram of the Lupus stars selected in our membership analysis. The empirical isochrone that we derive from our membership analysis is given in Table A.3. Our sample of cluster members includes stars in the magnitude range from G ≃ 8 mag to G ≃ 18 mag. Most of the cluster members identified in our study are located in Lupus 3 and 4 as illustrated in Fig. 5. We note that Lupus 1 and 2 host only a minor fraction of cluster members. Interestingly, our sample includes three and one stars projected towards Lupus 5 and 6, respectively. We confirm one star, namely Gaia DR2 6021420630046381440, from Manara et al. (2018) and another three stars (Gaia DR2 5995157042469914752, Gaia DR2 5997091667544824960 and Gaia DR2 5997086204346520448) from Melton (2020) which are located in Lupus 5 and 6 as bona fide members of the region. We also confirm the existence of a more dispersed population of stars spread over the entire complex that is less numerous than the off-cloud population reported in pre-Gaia studies (see e.g. Krautter et al. 1997; Wichmann et al. 1997a,b; Galli et al. 2013).

Fig. 3. Proper motions and parallaxes of the 137 Lupus stars identified in our membership analysis. The stars are colour-coded based on their membership probabilities which are scaled from zero to one. Triangles indicate the sources with RUWE ≥ 1.4 (see Sect. 4.1).

Fig. 4. Colour-magnitude diagram of the sample of Lupus stars identified in our membership analysis. The black solid line denotes the empirical isochrone of the Lupus association derived in this study (see Table A.3). The stars are colour-coded based on their membership probabilities which are scaled from zero to one. Triangles indicate the sources with RUWE ≥ 1.4 (see Sect. 4.1). The arrow indicates the extinction vector of AV = 1 mag converted to the Gaia bands using the relative extinction values computed by Wang & Chen (2019).

Fig. 5. Location of the Lupus stars overlaid on the extinction map of Dobashi et al. (2005) in Galactic coordinates. Colours and symbols denote the initial list of candidate stars compiled from the literature (88 stars) and the final members identified from our membership analysis (137 stars). The rectangles indicate the position of the main molecular clouds defined by Hara et al. (1999).

The Venn diagram represented in Fig. 6 shows that we have confirmed only 90 stars from previous studies as members of the Lupus association with our methodology. One important fraction of the literature sample with available data in Gaia-DR2, namely 351 stars, is rejected in our membership analysis. Damiani et al. (2019) classified 105 of them as Sco-Cen members. Among the remaining 246 sources we verify that most of them (i.e. 148 stars) have small parallaxes (ϖ <  3 mas) confirming that they are indeed background stars unrelated to the Lupus region (as already anticipated in Sect. 2, see Fig. 2). We find more background sources at lower Galactic latitudes as also observed by Manara et al. (2018). The remaining sample of 98 sources rejected from the literature includes 78 stars that are mostly spread well beyond the location of the Lupus clouds in a region not covered by our survey (l <  334°) and exhibit different proper motions and parallaxes compared to the cluster members identified in this paper. The other 20 stars are indeed projected towards the molecular clouds of the complex (mostly in Lupus 3), but only one of them, namely Gaia DR2 5997013327325390592, has proper motion and parallax consistent with membership in Lupus. This star has G_RP = 11.77 mag and G − G_RP = 1.65 mag which puts it significantly above the empirical isochrone defined by the other cluster members (see e.g. Fig. 4). It was first proposed to be a member of the Lupus association by Damiani et al. (2019) and we found no more information about this star in the literature. We suspect that it is a binary or a high-order multiple system to explain its position in the colour-magnitude diagram.

Fig. 6. Venn diagram comparing the number of Lupus stars compiled from the literature (with Gaia-DR2 data), Sco-Cen stars in the UCL subgroup, and the sample of confirmed Lupus members identified in this paper from our membership analysis.

A recent study conducted by Teixeira et al. (2020) and published during the revision process of this paper reports on the discovery of new YSOs with circumstellar discs in the Lupus complex. Two points are worth mentioning regarding the comparison of our sample of cluster members with their results. First, the survey conducted by Teixeira et al. (2020) covers a larger area in the sky than our study and extends clearly beyond the location of the Lupus molecular clouds. They identified 60 new YSOs that are potential members of the Lupus complex (these sources are labelled as “Lupus and/or UCL” in Table 4 of that study), but only 39 of them are included in the field selected for our membership analysis in this paper. Second, the methodology used in the two studies to select cluster members differs significantly. Teixeira et al. (2020) selected their sample of cluster members from proper motion cuts based on Lupus stars previously identified in the literature. As shown in Fig. 6, many of the Lupus stars identified in the literature are probable field stars (unrelated to the Lupus clouds) and the adopted proper motion constraints in that study to select new members (i.e. −27 <  μ_α cos δ <  3 mas yr⁻¹ and −33 <  μ_δ <  −13 mas yr⁻¹) extend to a region in the proper motion vector diagram that is mostly populated by UCL stars (see Fig. 2). As a result, we found only 3 stars (out of the sample of 39 sources) that are in common with the study of Teixeira et al. (2020). The remaining 36 sources listed by Teixeira et al. (2020) and included in our analysis have mostly small membership probabilities suggesting that they are unrelated to Lupus. Indeed, the authors themselves concluded in their study that the older stars in the sample are most probably UCL members.

In this study, we report on the discovery of 47 new members of the Lupus association, which represents an increase of more than 50% with respect to the number of confirmed members already known from previous studies. Most of them (i.e. 43 stars) had previously been assigned to the Sco-Cen association by Damiani et al. (2019). One reason to explain this is that the authors included only Lupus 3 in their membership analysis so that the Lupus stars belonging to the other subgroups were either assigned to Sco-Cen or to the field population. The 43 stars populate all subgroups of the Lupus complex (except for Lupus 3) and our membership analysis focused on the properties of the whole region was able to identify them as members. In the following we use our new sample of 137 members to revisit some properties of the Lupus star-forming region.

4. Properties of the Lupus subgroups

4.1. Refining the sample of Lupus stars

We refine our list of cluster members by selecting the stars with the best astrometry in the sample to accurately determine the properties of the Lupus subgroups. We remove 24 stars with poor astrometric solution based the re-normalised unit weight error (RUWE) criterion² (i.e. with RUWE ≥ 1.4). These sources exhibit the largest uncertainties in proper motions and parallaxes (see Fig. 3) so that their membership status will require further investigation with future data releases of the Gaia space mission. This step reduces the sample of cluster members to 113 stars.

4.2. Proper motions and parallaxes

Figure 7 shows the distribution of proper motions and parallaxes for the several subgroups of the Lupus complex. The mean values of each parameter are listed in Table 2. A first inspection of the sample reveals that the stars projected towards Lupus 3 exhibit different proper motions (in right ascension) that are shifted by about 1–2 mas yr⁻¹ with respect to the stars in Lupus 1, 2 and 4. This proper motion shift translates into 0.8–1.5 km s⁻¹ in the tangential velocities of the subgroups. On the other hand, we see no important differences between the proper motion components in declination and parallaxes among the various subgroups of the Lupus complex.

Fig. 7. Proper motions and parallaxes for the Lupus subgroups in our sample of cluster members.

We performed a two sample Kolmogorov-Smirnov test to quantitatively confirm our findings. We compare the distribution of proper motions and parallaxes from Lupus 3 stars with respect to the stars in the other molecular clouds (Lupus 1, 2 and 4). We find the p-values of 1.97 × 10⁻⁶, 0.64, and 0.23, for the analyses based on μ_α cos δ, μ_δ and ϖ, respectively. Thus, we conclude that the proper motion component in right ascension is the most distinctive astrometric feature that distinguishes the stars projected towards Lupus 3 from the other subgroups in the region. The few stars in Lupus 5 and 6 that have been retained in our sample have proper motions that are more consistent with Lupus 3 stars.

This small difference in proper motions among the subgroups is seen now for the first time with the more precise Gaia-DR2 data that is available nowadays. Previous studies about the kinematic properties of the Lupus region (see e.g. Makarov 2007; López Martí et al. 2011; Galli et al. 2013) used proper motions measured from the ground with a typical precision of about 2 mas yr⁻¹ that was not enough to detect subtle differences of the stars located in different subgroups. The median uncertainties in proper motions and parallaxes for the Lupus stars in our sample based on Gaia-DR2 data are σ_{μα cos δ} = 0.15 mas yr⁻¹, σ_μδ = 0.09 mas yr⁻¹ and σ_ϖ = 0.08 mas, respectively. This represents an improvement in precision of about one order of magnitude with respect to pre-Gaia studies.

4.3. Radial velocities

We search for radial velocity information for the stars in our sample using the CDS databases and data mining tools (Wenger et al. 2000). Our search in the literature is based on Wichmann et al. (1999), Gontcharov (2006), Torres et al. (2006), James et al. (2006), Guenther et al. (2007), Galli et al. (2013), Frasca et al. (2017), and the Gaia-DR2 catalogue. We find radial velocities for 52 stars in the sample of 113 sources with good astrometry (see Sect. 4.1). For a few stars we find more than one single radial velocity measurement and in such cases we use the most precise result.

The radial velocity distribution of the Lupus stars in our sample is shown in Fig. 8. We use the interquartile range (IQR) criterium to identify the outliers that lie over 1.5 × IQR below the first quartile or above the third quartile of the distribution. Only the radial velocity of Gaia DR2 5997015079686716672 satisfies this condition. We discard this radial velocity measurement (V_r = 15.9 ± 0.7 km s⁻¹, Frasca et al. 2017), but we still retain the star in the sample because its proper motion and parallax are consistent with membership in Lupus.

Fig. 8. Kernel density estimate of the distribution of radial velocity measurements for the sample of 52 members with available data. We used a kernel bandwidth of 2 km s−1. The tick marks in the horizontal axis mark the individual radial velocity values of each star.

We find that 20 candidate members in our sample of 113 stars have been identified as binaries or multiple systems in previous studies (Ghez et al. 1997; Merín et al. 2008). Nineteen of them have radial velocity information in the literature, but only EX Lup (Gaia DR2 5996902860784332800) shows periodic radial velocity variations most probably due to a substellar companion (Kóspál et al. 2014). Analogously, we retain this star in the sample, but we do not use its radial velocity in the forthcoming analysis. Altogether, this reduces the sample of sources with radial velocity information to 50 stars.

Table 3 lists the mean radial velocity values for the Lupus subgroups (except for Lupus 5 and 6 which have no published data). We find no significant difference between the mean radial velocities of the subgroups within the uncertainties. The admittedly large scatter of radial velocities that is reported here is most probably related to the precision of the individual radial velocity measurements. The median uncertainty of the radial velocities in our sample is 2.2 km s⁻¹. Removing the binaries from the sample has negligible impact in reducing the radial velocity scatter. On the other hand, we note that the median uncertainty on the tangential velocities derived from Gaia-DR2 data is only 0.3 km s⁻¹ and 0.4 km s⁻¹ in the right ascension and declination components, respectively. Thus, the precision of the radial velocity measurements reported in the literature is currently the main limitation to investigate the dynamics (e.g. internal motions, expansion and rotation effects) of the Lupus complex.

4.4. Distance and spatial velocity

In this section we re-visit the distance and kinematic properties of the Lupus region based on our new sample of cluster members, the more precise Gaia-DR2 astrometry and the radial velocity measurements collected from the literature. The first step consists in correcting the Gaia-DR2 parallaxes by the zero-point shift of −0.030 mas that is present in the published data (see e.g. Lindegren et al. 2018). There are different estimates of the zero-point correction for the Gaia-DR2 parallaxes in the literature that range from −0.031 ± 0.011 mas (Graczyk et al. 2019) to −0.082 ± 0.033 mas (Stassun & Torres 2018), but we verified that the impact of this correction in our solution is negligible given the close proximity of the Lupus clouds. In addition, we also add 0.1 mas yr⁻¹ in quadrature to the formal uncertainties on the proper motions of the stars to take the systematics errors of the Gaia-DR2 catalogue into account (see e.g. Lindegren et al. 2018). This procedure is likely to overestimate the proper motion uncertainty for some stars in the sample given that the systematic errors in the Gaia-DR2 catalogue depend on the position, magnitude and colour of each source (see e.g. Luri et al. 2018). These corrections are necessary for the purpose of estimating distances and spatial velocities with the corresponding uncertainties taking into account systematic errors, but they do not affect our results of the membership analysis.

Then, we used Bayesian inference to convert the observed parallaxes and proper motions into distances and 2D tangential velocities. The prior that we used for the angular velocity (i.e. 2D tangential velocity) is a beta function following the online tutorials available in the Gaia archive (see e.g. Luri et al. 2018). We investigated two types of prior families for the distance. The first one is based on purely statistical probability density distributions (hereafter, statistical priors) and the second type is inspired on astrophysical assumptions (hereafter, astrophysical priors) e.g. the luminosity and density profiles of stellar clusters. The statistical priors investigated in this study are based on the Uniform and Gaussian distributions. The astrophysical priors include the surface brightness profile derived from star counts in the Large Magellanic Cloud by Elson et al. (1987, hereafter, the EFF prior) and the King’s profile distribution derived from globular clusters (King 1962, hereafter, the King prior). We used the Kalkayotl code (Olivares 2019; Olivares et al. 2020)³ to compute the distances with these different prior families. We refer the reader to that paper for more details on the properties, validation and performance of these priors. We are aware about the existence of another prior in the literature, namely the exponentially decreasing space density prior (hereafter, the EDSD prior, Bailer-Jones 2015; Astraatmadja & Bailer-Jones 2016), but we do not use it in our study. The EDSD prior is related to the distribution of stars in the Galaxy and it was proposed to be used in the context of large samples with wide distributions of parallaxes and uncertainties which is not the case of the Lupus sample.

We find no significant differences in the resulting distances obtained with different priors (Uniform, Gaussian, EFF and King) which agree within 1σ of the corresponding distance uncertainties. Thus, we confirm that our results are not affected by our choice of the prior. We have therefore decided to use the distances obtained with the Uniform prior, which is the simplest one, and report them as our final results. The distance and spatial velocity for individual stars in our sample are given in Table A.1 together with other results derived in this paper. We use the resulting distances to compute the 3D position XYZ of the stars in a right-handed reference system that has its origin at the Sun, where X points to the Galactic centre, Y points in the direction of Galactic rotation, and Z points to the Galactic north pole. Then, we convert the 2D tangential velocity and radial velocity of the stars into the UVW Galactic velocity components in this same reference system following the transformation outlined by Johnson & Soderblom (1987).

Table 4 lists the distance and spatial velocity of the Lupus subgroups in our sample. The large uncertainties of our distance results observed for some subgroups can be explained by the small number of stars used in each solution. Our results reveal that the Lupus subgroups are located at the same distance (of about 160 pc) within the uncertainties derived for each solution. The only exception is Lupus 6 whose distance estimate given in Table 4 is based on the parallax of one single star. Despite the common distance to the Sun, the subgroups (Lupus 1–4) are separated by about 8–20 pc among themselves in the space of 3D positions as illustrated in Fig. 9.

Fig. 9. 3D spatial distribution of the 113 members of the Lupus association. The different colours denote the various subgroups of the Lupus star-forming region. The dashed lines connect each point in the plot to its projection on the XY plane (with Z = 0). The reference system used to represent the 3D spatial distribution of the stars is the same as described in the text of Sect. 4.4.

The distances derived in this paper are still consistent with the results obtained by Dzib et al. (2018) despite the different samples of stars that have been used in each study. Our results do not confirm the large spread in distances of tens of parsec along the line of sight reported in previous studies based on kinematic parallaxes (Makarov 2007; Galli et al. 2013) and HIPPARCOS data (see e.g. Bertout et al. 1999; Comerón 2008). One reason to explain the discrepancy with previous findings in the pre-Gaia era is the more precise and accurate data that is available nowadays to investigate the Lupus region. Another reason is the different samples of Lupus stars used in each study. For example, Makarov (2007) and Galli et al. (2013) included the more dispersed stars identified by Krautter et al. (1997) and Wichmann et al. (1997a,b) in their solutions (see Fig. 1) which are mostly rejected in our membership analysis. In addition, many of the stars projected towards the molecular clouds and included in their analyses appear now to be background sources unrelated to Lupus in light of the Gaia-DR2 data (see Sect. 3).

The mean spatial velocities of the Lupus subgroups given in Table 4 confirm that they move at a common speed. We therefore detect no significant relative motion between the subgroups at the level of a few km s⁻¹ which is currently determined by the precision of radial velocity data. The off-cloud population is located at the same distance and moves with the same velocity of the stars projected towards Lupus 1–4. This confirms that they are comoving with the on-cloud population and therefore belong to the same association.

We now compare the distance and kinematics of Lupus with the two adjacent subgroups of the Sco-Cen association (US and UCL). The relative motion of Lupus with respect to the space motion of US and UCL derived by Wright & Mamajek (2018) is (ΔU, ΔV, ΔW) = (2.4, −0.7, −0.4)±(0.4, 0.2, 0.1) km s⁻¹ and (ΔU, ΔV, ΔW) = (2.1, 2.4, −1.6)±(0.4, 0.2, 0.1) km s⁻¹, respectively (in the sense Lupus “minus” US or UCL). This implies a relative bulk motion of 2.5 ± 0.4 km s⁻¹ and 3.6 ± 0.3 km s⁻¹ with respect to US and UCL, respectively. Damiani et al. (2019) revealed that Sco-Cen is made up of several substructures and the distances of the stars range from about 100 to 200 pc. The closest substructure of the Sco-Cen association to the Lupus molecular cloud complex is a compact group of 551 stars near V1062 Sco and named as UCL-1 (Röser et al. 2018; Damiani et al. 2019). The mean parallax of ϖ = 5.668 ± 0.011 mas and distance of $174.7_{-0.6}^{+0.5}$ pc put UCL-1 further away from the Sun as compared to Lupus (see Table 4). Thus, despite the common location in the plane of the sky Lupus and Sco-Cen stars are located at different distances and exhibit distinct spatial velocities. Interestingly, we note that the observed difference between the spatial velocity of the Lupus clouds, UCL and US is comparable to the difference in spatial velocity of the Sco-Cen subgroups among themselves. The Lupus clouds could therefore be regarded as another subgroup of the Sco-Cen complex that resulted from a more recent star formation episode.

4.5. Relative ages of the subgroups

In this section we compare the relative age of the subgroups based on the location of the stars in the HR-diagram compared to evolutionary models and the fraction of disc-bearing stars in each population. First, we use the Virtual Observatory SED Analyzer (VOSA, Bayo et al. 2008) to fit the SED, derive the effective temperatures and bolometric luminosities of the stars. We use the BT-Settl (Allard et al. 2012) grid of models to fit the SED and left the extinction A_V as a free parameter (in the range of 0–10 mag). We use the optical photometry from Gaia-DR2 and infrared photometry from the 2MASS (Cutri et al. 2003), AllWISE (Wright et al. 2010) and Spitzer c2d (Merín et al. 2008) surveys. We provide ourselves the photometric data of the stars to the VOSA service when possible to avoid erroneous cross-matches when querying these catalogues with the system interface. We cross-matched our sample of stars with the 2MASS and AllWISE catalogues using the source identifiers given in the auxiliary tables TMASS_BEST_NEIGHBOUR and ALLWISE_BEST_NEIGHBOUR available in the Gaia archive. We derived effective temperatures and bolometric luminosities for 110 stars in the sample with available photometric data. The two stars of Lupus 2 in our sample are among the three rejected sources with insufficient photometric data to fit the SED with VOSA.

We compare the effective temperatures obtained from the SED fits with those derived from the spectral type of the stars compiled from the literature (Hughes et al. 1994; Krautter et al. 1997; Comerón et al. 2009, 2013; Alcalá et al. 2017). We find spectral classification for only 69 stars (i.e. 61% of the sample) and convert them to effective temperatures using the tables provided by Pecaut & Mamajek (2013). We observe a good agreement between the two effective temperature estimates with a mean deviation of −32 ± 31 K (in the sense, result from spectral type “minus” result from SED fit) and rms of 260 K. We have therefore decided to use the effective temperature estimates derived from the SED fits in the following analysis to work with a homogeneous data set and use this information for a more significant number of stars in the sample.

The resulting HR-diagram is shown in Fig. 10. We compare the location of the stars in the HR-diagram with the Baraffe et al. (2015) and Siess et al. (2000) models that combined together cover the entire mass range of our sample which ranges from about 0.03 to 2.4 M_⊙. Gaia DR2 5997082081177906048 (HR 5999) is a binary early-type star and the most massive member in our sample (m = 2.432 ± 0.122 M_⊙, see Vioque et al. 2018). The remaining sources in the sample are less massive than 1.4 M_⊙ and fall into the mass domain covered by the Baraffe et al. (2015) models. Our results suggest that the Lupus stars in the sample are mostly younger than 5 Myr. Many sources also lie above the 1 Myr isochrone and some of them are likely to be binaries or high-order multiple systems that will demand further investigation in future studies. Although age determination at these early stages of stellar evolution is rather uncertain, we detect no evident age gradient or segregation among the Lupus subgroups from the HR-diagram analysis. Table A.4 provides the age estimates inferred from the Baraffe et al. (2015) and Siess et al. (2000) models together with the results obtained from the SED fit with VOSA.

Fig. 10. HR diagram of the Lupus star-forming region with the Baraffe et al. (2015) and Siess et al. (2000) evolutionary models. The solid and dashed lines represent isochrones and tracks, respectively, with the corresponding ages (in Myr) and masses (in M⊙) indicated in each panel. The most massive star in our sample (HR 5999) is not shown to improve the visibility of the low-mass range that includes the other stars in our sample. The error bars are smaller than the size of the symbols.

Alternatively, we compare the age of the subgroups based on the fraction of disc-bearing stars in each population. We compute the spectral index α (Lada et al. 1987) of each source in the sample with available photometric data in the infrared between 2 and 20 μm. This restricts the sample to 104 stars with AllWISE and Spitzer data. Then, we classify the stars into Class I (embedded source, α >  0), Class II (disc-bearing star, −2 <  α <  0) or Class III (optically thin or no disc, α <  −2). We compare our results derived from the spectral index with the classification scheme proposed by Koenig & Leisawitz (2014) based on infrared colours and we obtain the same SED subclass for all stars in the sample.

Figure 11 shows that the distribution of the spectral index α (see Table A.1) is bimodal suggesting the existence of two populations (disc-bearing and discless stars) with approximately the same number of sources. We verified that the bimodality of the distribution of spectral indices is not an artefact caused by the different number of points and wavelength range of the photometric data available for each star to compute the spectral indices. A similar result was also observed for the Chamaeleon I star-forming region (see Fig. 11 of Luhman et al. 2008). Previous studies suggested that the early disappearence of circumstellar discs could be related to environmental effects imposed by the presence of massive stars that produce strong UV radiation, stellar winds, and supernova explosions (see e.g. Walter et al. 1994; Martin et al. 1998). However, we see no dependency of the spectral index on the position of the stars in our sample and any nearby OB star surrounding the Lupus clouds to support this scenario. As shown in Fig. 10 we note the existence of only a few stars older than 10 Myr in our sample which could be potential contaminants (as expected based on the performance of our classifier, see Table 1). Thus, the hypothesis of contamination by older field stars does not explain the bimodality of spectral indices in our sample given that the two populations have the same number of stars, similar ages and are younger than the potential contaminants from the Sco-Cen association. Alternatively, we investigated the dependency of the spectral indices on the age and colour of the stars. Figure 12 shows that the age distribution of the Class II and Class III stars in Lupus overlap. The median age of Class III stars inferred from the Baraffe et al. (2015) and Siess et al. (2000) stellar models is about 3 Myr which yields a rough estimate of the typical disc lifetime in the Lupus association. However, it should also be noted from Fig. 12 that we observe an excess of Class III stars at cooler temperatures (red colours) suggesting that the survival time of circumstellar discs may also depend on other stellar parameters. For example, Galli et al. (2015) used an empirical disc evolution model to determine the lifetime of circumstellar discs in Lupus in terms of the mass of the star. According to their model the average lifetime of a circumstellar disc around a star with 0.1 M_⊙ is of the order of 1 Myr which could explain the early disapperance of circumstellar discs for some stars in our sample.

Fig. 11. Histogram of the spectral index α computed for 104 stars in the sample from infrared photometry.

Fig. 12. Dependency of the spectral index on the age inferred from the Baraffe et al. (2015) models (left panel) and colour of the stars using Gaia photometry (right panel).

In Table 5 we compare the fraction of disc-bearing stars in the Lupus subgroups with their isochronal ages inferred from the Baraffe et al. (2015) and Siess et al. (2000) models. The following discussion about isochronal ages is obviously restricted to the stars covered by the evolutionary models used in this study (see Fig. 10). This reduces the sample of stars with age estimates to only a few members in some subgroups, but the current analysis is still useful to provide a relative dating among them. For example, we note that Lupus 3 has similar fractions of disc-bearing and discless stars. As shown in Fig. 13, the populations of YSOs in Lupus 3 overlap in the 3D space of positions with no evident segregation between the two subclasses making it the most extended subgroup of the complex with about 8 pc in the X direction (see also Fig. 9). It is interesting to note that the median ages of the Lupus subgroups inferred from isochrones range from about 1 to 3 Myr. The only exception shown in Table 5 is Lupus 5, but our sample has only two stars. One of them, namely Gaia DR2 6021420630046381440, appears to be older than 10 Myr. It was first identified as a Lupus 5 member by Manara et al. (2018), but the authors also raised the possibility that the targets selected in their study could be dispersed members of the Lupus 3 cloud. This would explain the discrepant age estimate that we observe for this star compared to the other source in our Lupus 5 sample, namely Gaia DR2 6021414479652887296, that has an age of about 1–2 Myr (see Fig. 10) and is more consistent with the age reported for the other subgroups. We can therefore conclude that the Lupus subgroups are coeval. Interestingly, we also note that the off-cloud population of Lupus has similar ages compared to the on-cloud stars. The common age reinforces our previous conclusion in Sect. 4.4 based on spatial velocities that the off-cloud stars belong to the same association of the stars projected towards the Lupus clouds.

Fig. 13. Spatial distribution of the Lupus stars in our sample of cluster members. Colours indicate the different subgroups of the Lupus complex. Filled and open symbols denote Class II and Class III stars, respectively.

4.6. Initial mass function

In the following we use our new sample of members to discuss the initial mass function (IMF) of the Lupus association. As illustrated in Fig. 10 many stars in our sample fall into an age and mass domain that are not covered by the pre-main sequence stars evolutionary models published in the literature. We have therefore decided to use the distribution of spectral types in the sample as a proxy for the IMF. This procedure also allows us to compare our results in Lupus with other star-forming regions investigated in previous studies that adopted a similar approach.

We convert the effective temperatures derived from the SED fits (see Table A.4) to spectral types using the tables given by Pecaut & Mamajek (2013). We adopt the conversion from effective temperature to spectral type for 5–30 Myr old stars given in Table 6 of that study for most stars in our sample (with temperatures ranging from 2880 to 7280 K) and use the tabulated values for dwarf stars (Table 5 of Pecaut & Mamajek 2013) for the remaining stars as an approximation. We estimate the resulting spectral types to be accurate by about 2 subclasses based on the rms of 260 K that we find when comparing the effective temperatures derived from the SED fit with the spectral classification of the stars with available information in our sample (see Sect. 4.5).

We estimate the completeness limits of our new sample of Lupus stars based on the completeness limits of the Gaia-DR2 catalogue which is the only source of data used for the membership analysis. The Gaia-DR2 catalogue is complete between G = 12 and G = 17 mag (see e.g. Gaia Collaboration 2018). The extinction is variable in the region covered by our survey (see Fig. 5) and the measurements given in the Gaia-DR2 catalogue for individual stars in our sample range from A_G = 0.8 to A_G = 2.5 mag with a median value of 1.5 mag. The Gaia-DR2 completeness limit of G = 17 mag corrected for the maximum observed extinction of the stars in our sample at the distance of the Lupus clouds (see Table 4) yields the absolute magnitude of G_abs = 8.5 mag. This translates into 0.2 M_⊙ using the 3 Myr isochrone from the BT-Settl models (Allard et al. 2012) and we therefore estimate our sample to be complete to this mass limit.

Figure 14 shows the distribution of spectral types in the Lupus association compared to the Taurus and US star-forming regions using the spectral classifications derived by Esplin & Luhman (2019) and Luhman (2018), respectively. It is apparent that the sample of US exhibits a surplus of early-type stars which appear in small number in Taurus and are restricted to only one member (namely, HR 5999) in Lupus. Otherwise, we note that the distribution of spectral types in Lupus resembles the distribution of Taurus and US. We therefore argue that the IMF shows little variation for late-type stars among the three regions despite the different number of stars in the various samples. The samples of stars in Taurus and US show an important number of members with substellar masses which have been identified using ancillary data from optical and infrared surveys (in addition to the Gaia-DR2 catalogue). When we restrict the Taurus and US samples to the members that have been observed by the Gaia satellite we find a similar shape for the distribution of spectral types at the faint end of the IMF that we obtain in this study for the Lupus association.

Fig. 14. Distribution of spectral types of the stars in Lupus (this paper), Taurus (Esplin & Luhman 2019) and US (Luhman 2018) as a proxy for the IMF. The dashed line in the upper panel indicates the completeness limit (with the corresponding mass) of the Lupus sample investigated in this study.

Our study based on Gaia-DR2 data addressed many of the past uncertainties regarding membership, distance and the kinematic properties of the Lupus association returning the most complete census of stars down to about 0.2 M_⊙. This step was made necessary to further investigate the population of Lupus stars. Future studies using ancillary data to complement the Gaia catalogue will search for the least-massive members and significantly increase the number of substellar objects. This will extend the IMF of the Lupus association down to the planetary-mass regime as done in other star-forming regions.

5. Conclusions

In this study, we have revised the census of stars and properties of the various subgroups in the young stellar association located in the Lupus star-forming region based on Gaia-DR2 data. We applied a probabilistic method to infer membership probabilities of more than 10⁷ sources over a field of 160 deg² encompassing the main star-forming clouds of the Lupus complex. We identified 137 stars (spread over the molecular clouds Lupus 1–6) that are probable members of the Lupus association. Our analysis confirms 90 stars from the literature that were previously associated with the Lupus region and adds 47 new members to list (an increase of more than 50% with respect to the number of known members). We confirm that about 80% of the historically known members of the region with available astrometry in the Gaia-DR2 catalogue are more likely to be background sources unrelated to Lupus, field stars or members of the adjacent Sco-Cen association in light of our new membership analysis.

Our results on the 6D structure of the Lupus region show that the different subgroups of the complex (defined by the sample of stars projected towards the various molecular clouds) are located at the same distance (of about 160 pc) and move with the same spatial velocity. This confirms that the Lupus subgroups are comoving and belong to the same association. The HR-diagram shows that most stars in our sample are younger than 5 Myr and they cover the mass range from about 0.03 to 2.4 M_⊙. Our results show that there is a superposition of ages between disc-bearing and discless stars in the Lupus association. The median age of the Lupus subgroups ranges from about 1 to 3 Myr and we therefore conclude that they are coeval. Both age indicators (fraction of disc-bearing stars and isochronal ages) reveal that the Lupus association is younger than the population of YSOs in the Corona-Australis star-forming region recently investigated by our team with the same methodology (Galli et al. 2020). Our new sample of Lupus stars is complete down to 0.2 M_⊙ and the IMF shows little variation compared to other star-forming regions.

The main limitation to investigate the dynamics of the Lupus complex (e.g. relative motion of the subgroups, expansion and rotation effects) is currently hampered by the precision on the radial velocities available in the literature (about 2 km s⁻¹). We therefore encourage astronomers to perform high-resolution spectroscopy of Lupus stars to allow for more precise studies on the kinematic properties of this region in conjunction with the future data releases of the Gaia space mission.

Acknowledgments

We would like to thank the referee for carefully reading our manuscript and providing constructive comments. This research has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No 682903, P. I. H. Bouy), and from the French State in the framework of the “Investments for the future” Program, IdEx Bordeaux, reference ANR-10-IDEX-03-02. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, funded by the National Aeronautics and Space Administration. This publication makes use of VOSA, developed under the Spanish Virtual Observatory project supported by the Spanish MINECO through Grant AyA2017-84089. VOSA has been partially updated by using funding from the European Union’s Horizon 2020 Research and Innovation Programme, under Grant Agreement n^o 776403 (EXOPLANETS-A).

References

Appendix A: Additional tables (online material)

All Tables

All Figures

	Fig. 1. Location of the Lupus stars identified in previous studies overlaid on the extinction map of Dobashi et al. (2005) in Galactic coordinates. Colours and symbols denote the samples of stars from each study. The red dashed lines indicate the position of the main molecular clouds defined by Hara et al. (1999).
In the text

	Fig. 2. Proper motion of Lupus and UCL stars identified in previous studies. The different colours denote the subgroups of the Lupus sample. The cyan rectangle indicates the locus of Lupus candidate stars in the space of proper motions (see discussion in Sect. 2).
In the text

	Fig. 3. Proper motions and parallaxes of the 137 Lupus stars identified in our membership analysis. The stars are colour-coded based on their membership probabilities which are scaled from zero to one. Triangles indicate the sources with RUWE ≥ 1.4 (see Sect. 4.1).
In the text

Fig. 4. Colour-magnitude diagram of the sample of Lupus stars identified in our membership analysis. The black solid line denotes the empirical isochrone of the Lupus association derived in this study (see Table A.3). The stars are colour-coded based on their membership probabilities which are scaled from zero to one. Triangles indicate the sources with RUWE ≥ 1.4 (see Sect. 4.1). The arrow indicates the extinction vector of AV = 1 mag converted to the Gaia bands using the relative extinction values computed by Wang & Chen (2019).

In the text

	Fig. 5. Location of the Lupus stars overlaid on the extinction map of Dobashi et al. (2005) in Galactic coordinates. Colours and symbols denote the initial list of candidate stars compiled from the literature (88 stars) and the final members identified from our membership analysis (137 stars). The rectangles indicate the position of the main molecular clouds defined by Hara et al. (1999).
In the text

	Fig. 6. Venn diagram comparing the number of Lupus stars compiled from the literature (with Gaia-DR2 data), Sco-Cen stars in the UCL subgroup, and the sample of confirmed Lupus members identified in this paper from our membership analysis.
In the text

	Fig. 7. Proper motions and parallaxes for the Lupus subgroups in our sample of cluster members.
In the text

	Fig. 8. Kernel density estimate of the distribution of radial velocity measurements for the sample of 52 members with available data. We used a kernel bandwidth of 2 km s−1. The tick marks in the horizontal axis mark the individual radial velocity values of each star.
In the text

	Fig. 9. 3D spatial distribution of the 113 members of the Lupus association. The different colours denote the various subgroups of the Lupus star-forming region. The dashed lines connect each point in the plot to its projection on the XY plane (with Z = 0). The reference system used to represent the 3D spatial distribution of the stars is the same as described in the text of Sect. 4.4.
In the text

Fig. 10. HR diagram of the Lupus star-forming region with the Baraffe et al. (2015) and Siess et al. (2000) evolutionary models. The solid and dashed lines represent isochrones and tracks, respectively, with the corresponding ages (in Myr) and masses (in M⊙) indicated in each panel. The most massive star in our sample (HR 5999) is not shown to improve the visibility of the low-mass range that includes the other stars in our sample. The error bars are smaller than the size of the symbols.

In the text

	Fig. 11. Histogram of the spectral index α computed for 104 stars in the sample from infrared photometry.
In the text

	Fig. 12. Dependency of the spectral index on the age inferred from the Baraffe et al. (2015) models (left panel) and colour of the stars using Gaia photometry (right panel).
In the text

	Fig. 13. Spatial distribution of the Lupus stars in our sample of cluster members. Colours indicate the different subgroups of the Lupus complex. Filled and open symbols denote Class II and Class III stars, respectively.
In the text

	Fig. 14. Distribution of spectral types of the stars in Lupus (this paper), Taurus (Esplin & Luhman 2019) and US (Luhman 2018) as a proxy for the IMF. The dashed line in the upper panel indicates the completeness limit (with the corresponding mass) of the Lupus sample investigated in this study.
In the text