conditional pdf $f(y\mid r)$. Assuming that the data are Gaussian distributed
around the linear model with a constant standard deviation $\sigma$383838We
observe the standard deviation to be approximately constant for the range of
interest; $r\in[100,200]$, we can write the following expression:

$f(y\mid r)\propto\text{exp}\left(-\frac{(y-g(r))^{2}}{2\sigma^{2}}\right)$
(31)

Figure 16 schematically shows the integrating process where the conditional
pdfs $f(y\mid r)$ are shown for $r=100$ and $r=200$.

Figure 16: Scaling factor determination via the empirical distance-scaling
relationship. The scaling factor distribution for groups at a distance between
$100-200$ pc depends on the conditional distribution of scaling factors at a
given distance $f(c_{x}/c_{v}\mid r)$.

By substituting Eq. (31) into Eq. (30) and solving the integral we obtain:

$f(y\mid r_{0}\leq r\leq
r_{1})\propto\text{erf}\left(\frac{y-g(r_{0})}{\sqrt{2}\sigma}\right)-\text{erf}\left(\frac{y-g(r_{1})}{\sqrt{2}\sigma}\right)$
(32)

where erf($x$) is the error function. To normalize the probability density in
Eq. (32), we compute its integral. Since both summands are of the same type,
we only need to solve the following integral:

$I(c)=\int_{-\infty}^{+\infty}\text{erf}\left(\frac{y-c}{\sqrt{2}\sigma}\right)\,dy$
(33)

The variable $c$ represents the constants $g(r_{0})$ and $g(r_{1})$. The
integral in Eq. (33) evaluates to:

$I(c)=\text{exp}\left(\frac{-(y-c)^{2}}{2\sigma^{2}}\right)\sqrt{\frac{2\sigma^{2}}{\pi}}+(y-c)\leavevmode\nobreak\
\text{erf}\left(\frac{y-c}{\sqrt{2}\sigma}\right)\leavevmode\nobreak\
\Biggr{|}_{-\infty}^{+\infty}$ (34)

Thus, the integral of Eq. (32) can be expressed in the following form:

$\int_{-\infty}^{+\infty}f(y\mid r_{0}\leq r\leq
r_{1})=n\times[I(g(r_{0}))-I(g(r_{1}))]\overset{!}{=}1$ (35)

The factor $n$ represents the normalization factor. Rearranging the resulting
terms by function type, we get the following:

$\displaystyle\begin{aligned}
I(g(r_{0}))-I(g(r_{1}))=h(y)+l(y)\leavevmode\nobreak\
\Biggr{|}_{-\infty}^{+\infty}\end{aligned}$ (36)

The functions $h(y)$ and $l(y)$ describe a sum of exponential and error
functions, respectively. The functions are defined in the following:

$\displaystyle\begin{aligned}
h(y)&=\sqrt{\frac{2\sigma^{2}}{\pi}}\left[\text{exp}\left(\frac{-(y-g(r_{0}))^{2}}{2\sigma^{2}}\right)-\text{exp}\left(\frac{-(y-g(r_{1}))^{2}}{2\sigma^{2}}\right)\right]\\\
l(y)&=(y-g(r_{0}))\leavevmode\nobreak\
\text{erf}\left(\frac{y-g(r_{0})}{\sqrt{2}\sigma}\right)-(y-g(r_{1}))\leavevmode\nobreak\
\text{erf}\left(\frac{y-g(r_{1})}{\sqrt{2}\sigma}\right)\end{aligned}$ (37)

We evaluate the summands of this primitive integral at the border
individually. First, sum Gaussians in $h(y)$ goes to zero as $y$ approaches
negative and positive infinity:

$\lim_{y\to\pm\infty}h(y)=0$ (38)

The sum of error functions can be rearranged into the following form:

$\displaystyle\begin{aligned}
l(y)&=y\left[\text{erf}\left(\frac{y-g(r_{0})}{\sqrt{2}\sigma}\right)-\text{erf}\left(\frac{y-g(r_{1})}{\sqrt{2}\sigma}\right)\right]\\\
&+g(r_{1})\leavevmode\nobreak\
\text{erf}\left(\frac{y-g(r_{1})}{\sqrt{2}\sigma}\right)-g(r_{0})\leavevmode\nobreak\
\text{erf}\left(\frac{y-g(r_{0})}{\sqrt{2}\sigma}\right)\end{aligned}$ (39)

Evaluating $l(y)$ at the borders results in the following:

$\lim_{y\to\pm\infty}y\left[\text{erf}\left(\frac{y-g(r_{0})}{\sqrt{2}\sigma}\right)-\text{erf}\left(\frac{y-g(r_{1})}{\sqrt{2}\sigma}\right)\right]=0$
(40)

and

$\displaystyle\begin{aligned}
\lim_{y\to\pm\infty}&\left[g(r_{1})\leavevmode\nobreak\
\text{erf}\left(\frac{y-g(r_{1})}{\sqrt{2}\sigma}\right)-g(r_{0})\leavevmode\nobreak\
\text{erf}\left(\frac{y-g(r_{0})}{\sqrt{2}\sigma}\right)\right]=\\\
&=\pm(g(r_{1})-g(r_{0}))\end{aligned}$ (41)

This last term is the only non-zero contribution to the integral. Its
evaluation at the lower edge results in the same but negative value to the
upper edge. Thus, the area under the curve is twice that value. The
normalization factor $n$ then becomes the following:

$\displaystyle\begin{aligned} I(g(r_{0}))-I(g(r_{1}))&=2[g(r_{1})-g(r_{0})]\\\
n&=\frac{1}{2[g(r_{1})-g(r_{0})]}=:\frac{1}{2\Delta g}\end{aligned}$ (42)

The function value difference, $\Delta g$, is always positive since $g(r)$ is
a strictly monotonically increasing function and $r_{1}>r_{0}$, see Fig. 16.
Thus, $n$ is a proper normalization factor that is non-zero and positive for
all pairs $r_{0}$ and $r_{1}$. Using this normalization constant, the
conditional pdf can be written as:

$\displaystyle f(y\mid r_{0}\leq r\leq r_{1})=\frac{1}{2\Delta
g}\left[\text{erf}\left(\frac{y-g(r_{0})}{\sqrt{2}\sigma}\right)-\text{erf}\left(\frac{y-g(r_{1})}{\sqrt{2}\sigma}\right)\right]$
(43)

Figure 17: Pdf and cdf of scaling factor conditioned on a given range of
distances. The ten red scatter points indicate samples drawn from the
10-quantile splitting procedure. We separate the pdf into 10 continuous
intervals with equal probabilities from which we derive samples as the mean
position of these intervals.

The top part of Fig. 17 shows the resulting pdf when applying Eq. (32) to
sources in Sco-Cen, where we assume a distance range of $r\in[100,200]$. Here
we can see an immediate caveat of our simple constant model uncertainty
assumption; the resulting distribution has infinite support, thus, non-zero
probability density for $f(y<0\mid r)$. Although physically meaningless, the
total probability of such events is small and, as seen below, does not
drastically influence the final sample set.

We consider sampling strategies to obtain a set of scaling factors to use in
the clustering process. Random sampling can generate almost identical
realizations, so the possible solution space might not be covered evenly.
Since we need to perform a separate clustering, run for each sample drawn,
keeping the number as small as possible is essential. To cover the space
evenly while considering the underlying probability distribution, we select a
set of 10 samples that represent 10 quantiles of the pdf. We separate the pdf
into 10 continuous intervals with equal probabilities from which we derive
samples as the mean position of these intervals.

To compute the quantiles, we determine the cdf by solving the integral over
the conditional pdf in Eq. (43); using functions $h(y)$ and $l(y)$ defined in
Eq. (37) the cdf becomes:

$\displaystyle\begin{aligned} F(y\mid r_{0}\leq r\leq r_{1})&=\int f(y\mid
r_{0}\leq r\leq r_{1})\,dy\\\ &=\frac{1}{2\Delta g}[h(y)+l(y)]+C\end{aligned}$
(44)

To obtain a proper cdf from Eq. (44), we set the constant of integration $C$
to $1/2$. Thus, the cdf becomes:

$F(y\mid r_{0}\leq r\leq r_{1})=\frac{1}{2}\left(1+\frac{1}{\Delta
g}[h(y)+l(y)]\right)$ (45)

The cdf defined in Eq. (44) for $r\in[100,200]$ is shown in the bottom part of
Fig. 17. The ten red scatter points393939The velocity scaling values are:
$c_{v}=\\{1.43,4.26,5.95,7.3,8.51,9.68,10.89,12.24,13.93,16.77\\}$. indicate
samples drawn from the 10-quantile splitting procedure where horizontal lines
indicate equal probability intervals. To invert the cdf and obtain scaling
fraction samples from $F^{-1}(y\mid r_{0}\leq r\leq r_{1})$ we used a
numerical approximation404040We made use of the open source library pynverse
v0.1.4.4 to calculate the numerical inverse of the cdf..

##  Appendix E Consensus clustering

In this section we discuss algorithmic means to identify stable cluster
solutions from an ensemble of clustering results. To identify robust
solutions, we represent all clusters across the clustering ensemble in a
graph. Each cluster of a single solution from the ensemble is represented by
one node. We connect two nodes via an edge if two respective clusters share at
least one common point. Figure 18 highlightes this step in the first two
frames. The ensemble is comprised of three clustering solutions. Its
individual runs $A$, $B$, and $C$ contain three, two, and four clusters,
respectively. Each source is classified into a single group for a given run,
resulting in disjoint sets. Thus, edges connect clusters from different runs.

Edges in the graph are weighted by the corresponding Jaccard similarity, which
measures their common overlap. The Jaccard similarity is defined as the ratio
of the intersection of two sets over their union. Typically, a value of $0.5$
or greater indicates high similarity between two sets. Since this linkage is
fundamental for determining the consensus result, it should avoid connecting
dissimilar clusters. Thus, we remove edges with a weight below $0.5$ (see Fig.
18, step 3). This threshold is quite conservative, as it can separate even
similar cluster solutions, e.g., if one cluster is a subset of the other. To
avoid over-pruning the graph, we relax the cut criterion in the following way.
Two clusters $a$ and $b$ with respective sizes $n_{a}$ and $n_{b}$, where
$n_{a}>n_{b}$, are linked if the $n_{b}$ densest points of $a$ amount to a
Jaccard similarity of $0.5$ or greater with points from cluster $b$. This
criterion guarantees connectivity between cluster extractions at different
isosurface thresholds while separating ties to clusters that randomly fragment
into multiple subclusters.

Robust clusters that exist throughout the ensemble will have strong
connections to their counterparts from different runs while having none or
very weak connections to other nodes in the graph. Thus, a robust cluster
solution builds a strong clique in this graph414141A clique in a graph is a
subset of nodes that are all connected with each other. Thus, every pair of
nodes in a clique is joined by an edge.. In contrast, unstable clusters will
have many weak connections to many other clusters from different runs but none
or very few strong ones; see Fig. 18, panel 4.

To extract all stable cluster solutions, we aim to identify all strongly
connected cliques in the graph, i.e., the consensus result. Since individual
sources and clusters can be part of several cliques we employ a voting
strategy to determine the final data partitioning (Vega-Pons & Ruiz-
Shulcloper, 2011). Each clique is represented as a vector of length $N$, where
$N$ is the number of sources in the data set. The $N$ values correspond to the
sum of the individual clusters represented as 0–1 vectors, where all entries
are 1 for sources inside a cluster and 0 elsewhere. Each source is then
associated to a single clique by maximizing the respective entry at the
source’s position (in the vector). Thus, the larger a clique, the more likely
it wins a vote. To favor robust cluster solutions, we multiply each vector by
the median of its connection strengths. This number is maximized if the
cluster is unchanged throughout different runs. This step is summarized in
panel 5 in Fig. 18.

Figure 18: Consensus clustering pipeline for a simple example using an
ensemble of three clustering solutions $A$, $B$, and $C$. (2) Clusters from
the ensemble are linked in a graph based on overlapping points. (3) Edges
between clusters are removed if the overlap between their members is
insufficient to assume a common cluster solution. (4) Cliques are extracted
representing stable cluster solutions, i.e., consensus clusters. (5) A voting
strategy determines the assignment of individual sources to cliques.

##  Appendix F Testing the independence assumption of resampled Gaia data

When testing for multimodality, SigMA simulates multiple realizations of the
Gaia data set to see how robust density dips between pairs of neighboring
peaks are. In Sect. 3.4, we discuss how these realizations are used in a
global hypothesis test that determines if samples of pairwise adjacent density
enhancements are part of a single underlying mode or not. This global
hypothesis test combines individual tests from each realization and evaluates
the global null hypothesis that no $p$-value is significant.

However, the choice of combining different $p$-values is affected by the
correlation structure between distinct tests. Typically, statistical tests
require independence across tests to guarantee proper levels of specificity
and sensitivity. In the case of the commonly used Fisher’s combination test
(Fisher, 1934), positively correlated $p$-values increase the chance of type I
errors (rejects a null hypothesis that is actually true). To investigate the
independence assumption of $p$-values across resampled data sets we consider
the effects that resampling has on the modality test procedure.

The modality test is fully described by the dimensionality of the data and
most importantly by the $k$-distances of modal and saddle points. Resampling
the data set changes the latter. The sampling of new data points is done with
Gaussians centered on mean astrometric observables with heteroscedastic error
covariance matrices available in the Gaia database. We aim to use the
covariance matrix in relation to the nearest neighbor distance to estimate the
interdependence of $p$-values.

As the entries of the covariance matrix shrink to zero, the resampled data
points converge to the original data points until they eventually coincide. At
the same time, the $p$-values across different data sets approach the same
value leading to a perfect correlation between them. On the other hand, if the
standard deviation along its principle axis (or any other direction for that
matter) extends far beyond a point nearest neighbor or even its $k$-distance,
the resampled data differs substantially from the original one, decorrelating
the $p$-values.

Figure 19 illustrates this relationship for data in the Sco-Cen box. The
histogram shows the relative uncertainty of data points, i.e., the ratio of
positional uncertainty given by the error covariance matrix over the nearest
neighbor distance. The distances are computed in the space of observed
astrometric quantities where the uncertainties are assumed Gaussian. The
majority of data points are far below unity, i.e., errors are small relative
to their absolute values. The majority (68%, i.e., $1\sigma$) has a relative
uncertainty below $0.08$ (see the 68–95–99.7% percentiles in Fig. 19 marked as
$1\sigma,\,2\sigma,\,3\sigma$, respectively). The concentration of relative
uncertainty values at zero indicates that $k$-distances across resampled data
set are strongly correlated. Thus, we cannot assume independence of $p$-values
across different samples.

Figure 19: Log-log histogram showing the relative uncertainty of data points,
i.e., the ratio of positional uncertainty given by the error covariance matrix
over the nearest neighbor distance. The majority of data points are far below
unity. The concentration of relative uncertainty values at zero indicates that
$k$-distances across resampled data set are strongly correlated. The
$1\sigma,\,2\sigma,\,3\sigma$ lines indicate the percentiles containing 68,
95, and 99.7% of the respective distribution.

##  Appendix G Distribution of point-wise densities

Instead of directly recovering clusters from phase space data, SigMA extracts
modal regions which are a mixture of field stars and cluster members. To
separate signal from background, we employ a density-based classifier that
selects cluster members as an over-density over background. To automatically
determine a density threshold we consider the distribution of field stars and
cluster members in univariate density.

As discussed in Sect. 3.5, we approximately can the treat the field star
content as a uniform distribution locally around a cluster in phase space. A
uniform distribution in phase space translates to a Gaussian distribution in
1D density space (see Sect. 3.5 for a more detailed discussion).

Cluster members are commonly modeled as multivariate Gaussians (e.g. Gagné et
al., 2014; Sarro et al., 2014; Crundall et al., 2019; Riedel et al., 2017).
Although observational findings point to more complex morphologies (e.g.
Meingast & Alves, 2019; Röser et al., 2019; Jerabkova et al., 2021; Meingast
et al., 2019; Kounkel & Covey, 2019; Cantat-Gaudin et al., 2019a; Wang & Ge,
2021; Coronado et al., 2022), the Gaussianity assumption provides a good
starting point to consider the point-wise univariate density distribution.
Figure 20 shows multiple distributions of point-wise density estimates of
$100,000$ samples drawn from an $N$-dimensional Gaussian424242The Gaussian is
chosen to have zero mean and identity covariance matrix. The left column shows
the likelihood of individual samples. It provides an assessment of the local
density under the true model. The number of samples with a relatively high
likelihood decreases exponentially as the dimension $N$ increases from one
(top row) to six (bottom row).

Neighborhood queries in high dimension are plagued by the curse of
dimensionality. As the dimensionality grows, points are increasingly isolated
making neighborhoods no longer local. This effect can already be seen for
moderate dimensions in the right column of Fig. 20. It shows point-wise
density estimates obtained via the $k$-NN technique. Around the fourth
dimension, the distribution of $k$-distances starts to converge to a normal
distribution, which incorrectly suggests an underlying uniform distribution in
$N$-dimensional input space.

We compute point-wise density estimations in five and six dimensions. Although
we cannot specifically write down a generative model for stellar groups in
phase space, we can assume that neighborhood queries are strongly affected by
the given dimensionality. Thus, we model both background and signal
contributions as Gaussians in univariate density space.

Figure 20: Distribution of point-wise density estimates of samples drawn from
an $N$-dimensional Gaussian. The left column shows the likelihood of
individual samples. It provides an assessment of the local density under the
true model. The right column shows estimated point-wise densities by employing
the $k$-NN technique. By increasing the dimensionality $N$ from one to six
(top to bottom row) the distribution of point-wise densities approaches a
Gaussian.

##  Appendix H Parameter optimization

The parameter choices of our proposed SigMA analysis pipeline are tuned to
Gaia data (see Sect. 3.3). In contrast, DBSCAN and HDBSCAN, which we use to
compare and test our clustering technique, are general clustering techniques
whose parameters we have to set. Instead of using subjective, error-prone
prior knowledge to determine suitable parametrizations, we search the space of
possible clustering results for the best result. This strategy provides a
measure of the peak performance (in case of comprehensive/absolute prior
knowledge) a clustering algorithm can achieve. Thus, a comparison against the
best results allows for a discussion on methodological advantages and
disadvantages rather than reflecting poor parameter selection.

Compared to the original data set (see Sect. 4.2.1), we generate 10 different
realizations of the “compact cluster sample” to estimate an average
performance across random configurations (see. Sect. 4.2.2). Optimal
parametrization can fluctuate between different data realizations in the
compact cluster sample case. In reality however, the compact cluster sample is
a single cluster region in which a unique parametrization should be applied
to. Thus, after the parameter search on the 10 individual data samples, we
take the median of the resulting parameters and re-run the algorithm on every
sample to obtain a clustering score to report in Table 2.

To search the space of possible model configurations we employ a grid search
in which we evaluate the clustering algorithm on a regular grid in parameters
space. Compared to other parameter tuning methods such as random search or
Bayesian optimization, a grid search has significant benefits in this
scenario. First, a large portion of the parameter space is discrete. Thus,
there is a finite grid step size that can cover the entire parameter space in
these sub-spaces. Second, the parameter spaces are low dimensional (maximally
four-dimensional) allowing densely spaced samples in each parameter axis.
Third, compared to random search and Bayesian optimization, a grid search is
deterministic and, thus, provides reproducibility. Further, as the grid is
predetermined (compared to Bayesian optimization) it allows to fully
parallelize the computation to guarantee dense sampling in reasonable times.

In the following subsections, we describe algorithm-specific parameter choices
of the employed grid search. For each algorithm, we also search optimal values
for the velocity scaling factor. We adopt values discussed in Appendix D.

### DBSCAN

DBSCAN has two main parameters, epsilon and min_samples. The parameter epsilon
describes a neighborhood radius; in particular, it is the maximum distance
within which two samples are considered to be neighbors. The parameter
min_samples denotes the minimum number of samples required in an epsilon-
neighborhood to be considered a cluster (or specifically a core point which
seeds a cluster).

Since we normalize the velocities to the spatial sub-space XYZ we can treat
epsilon as a distance in parsecs. Together with the given minimum number of
points it defines a minimum density that is needed to be considered a cluster.
We search for optimal results within a range of epsilon $\in[2,25]$ pc with a
step size of $0.5$ pc. At the same time we vary minimum number of samples in
the following range: min_samples $\in[4,40]$ with a step size of 2.

We find optimal parameters for the open cluster sample to be (epsilon,
min_samples, $c_{v}$) = (9.5, 6, 5.95). For the compact cluster sample we find
an optimal solution with (epsilon, min_samples, $c_{v}$) = (9.0, 20, 5.95)

### HDBSCAN

HDBSCAN has three main parameters, min_cluster_size, min_samples, and
cluster_selection_method. Intuitively, min_cluster_size determines the
smallest cluster sizes that HDBSCAN considers. Which points are still
associated to a cluster is determined by min_samples. By increasing
min_samples clusters are progressively more forced into denser areas leaving
more points to be declared as noise. The parameter cluster_selection_method
determines how clusters are selected from the cluster tree hierarchy, see
Appendix B for a detailed discussion.

We search for optimal results within a range of min_cluster_size $\in[20,100]$
with a step size of $2$. At the same time, we vary the minimum number of
samples min_samples in the same range and step size while requiring
min_cluster_size $\geq$ min_samples434343This is an intrinsic requirement of
the algorithm.. The selection method parameter can take the following two
values: cluster_selection_method $\in$ [“leaf”, “eom”].

We find optimal parameters for the open cluster sample to be
(min_cluster_size, min_samples, cluster_selection_method, $c_{v}$) = (60, 60,
“eom”, 8.51). For the compact cluster sample, we find an optimal solution with
(min_cluster_size, min_samples, cluster_selection_method, $c_{v}$) = (24, 18,
“eom”, 5.95).

##  Appendix I Projected velocities

Figure 21: Tangential velocities in the $v_{\alpha}$/$v_{\delta}$ plane,
showing the theoretical locations of sources with circular Galactic orbits and
LSR velocities. Shown are six different cases, while each of the lines
represents sources at all $l$ positions. The six cases are for two different
distances (100 pc, dashed lines; 200 pc, dash-dotted lines), and for three
different $b$ positions ($b=-20^{\circ}$, green; $b=0^{\circ}$, blue;
$b=25^{\circ}$, magenta). The indicated longitude positions at $l=0^{\circ}$
(box symbols) and $l=290^{\circ}$ (diamond symbols) roughly mark the eastern
and western borders of Sco-Cen. The SigMA selected Sco-Cen members are shown
as gray dots (without stability cut). See also Fig. 13 for a separation of the
clusters and for the $v_{\alpha,\mathrm{LSR}}$/$v_{\delta,\mathrm{LSR}}$
plane.

The reflex motion of the Sun influences how the observed tangential velocities
are distributed in $v_{\alpha}$/$v_{\delta}$ space. Figure 21 shows the
theoretical positions of objects if they follow a circular orbit around the
Galactic center at given positions within the Galactic potential. The orbits
are estimated within a Milky Way potential including a disk, bulge, and halo
component, using the python package galpy by Bovy (2015)
(galpy.potential.MWPotential2014; galpy.potential.vcirc) and assuming the
local standard of rest (LSR) velocity from Schönrich et al. (2010). The
projected motions are given for all Galactic longitude ($l$) positions at two
distances ($d$) of 100 pc and 200 pc and at three Galactic latitudes ($b$) of
$-20^{\circ}$, $0^{\circ}$, and $25^{\circ}$. These $d$ and $b$ ranges
encompass the Sco-Cen region, which reaches from about $l=0^{\circ}$ to
$290^{\circ}$. The members of Sco-Cen within the selected SigMA clusters are
plotted as gray dots in Fig. 21.

Overall, the young stellar groups in Sco-Cen seem to roughly follow expected
motions in our Galaxy assuming they follow the LSR velocity. The figure
additionally highlights the issues that come with the projected tangential
velocity plane $v_{\alpha}$/$v_{\delta}$, which is a function of position in
the sky and the Sun’s motions. Very nearby stellar groups, like in nearby
young local associations, could cover large areas in the sky and consequently
also in the tangential velocity plane, while being at the same time confined
in 3D velocity space (UVW). This effect can be reduced if computing the
tangential velocities relative to the LSR
$v_{\alpha,\mathrm{LSR}}$/$v_{\delta,\mathrm{LSR}}$, which eliminates the
influence of the reflex motion of the Sun (see also Sect. 2 and Appendix A). A
comparison of the two velocity spaces ($v_{\alpha}$/$v_{\delta}$ and
$v_{\alpha,\mathrm{LSR}}$/$v_{\delta,\mathrm{LSR}}$) is shown in Fig. 13,
displaying the SigMA selected members of the Sco-Cen clusters.

##  Appendix J The Gaia DR3 HRD

In this section we first describe our procedure to estimate the fraction of
possible contaminants form older stellar populations (or field stars) in the
SigMA selected Sco-Cen sample using the Gaia HRD, and second we estimate the
fraction of sub-stellar objects (brown dwarf candidates).

###  J.1 Estimating the contamination from older sources

The HRD in Fig. 14 in Sect. 5 shows a color-absolute-magnitude diagram
(equivalent to an Hertzsprung-Russel diagram, HRD) using the magnitudes from
the Gaia DR3 passbands $G$ versus $G_{\mathrm{BP}}-G_{\mathrm{RP}}$. Not all
sources have detections in all three passbands; within the SigMA selected Sco-
Cen members there are 12,724 (97%) sources which have an entry in all three
passbands. The absolute magnitude $G_{\mathrm{abs}}$ is calculated with the
distance modulus using the inverse of the parallax as distance. To estimate
how many SigMA selected Sco-Cen members are consistent with the expected ages
below $\sim$20 Myr, and which sources could be contaminants from older
populations, we first need to apply photometric quality criteria. Inferior
photometric quality mostly affects the fainter low-mass sources and shifts
them to the left in the Gaia HRD. We use the magnitude errors and the
photometric flux excess factor $C$, which are defined as follows:

$\displaystyle
G_{\mathrm{err}}=1.0857/\mathtt{phot\\_g\\_mean\\_flux\\_over\\_error}$ (46)
$\displaystyle
G_{\mathrm{RP,err}}=1.0857/\mathtt{phot\\_rp\\_mean\\_flux\\_over\\_error}$
$\displaystyle
G_{\mathrm{BP,err}}=1.0857/\mathtt{phot\\_bp\\_mean\\_flux\\_over\\_error}$
$\displaystyle
C=\mathtt{phot\\_bp\\_rp\\_excess\\_factor}=(I_{\mathrm{BP}}+I_{\mathrm{RP}})/I_{G}$
$\displaystyle C^{*}=\mathrm{corrected}\,\,\,C$

The flux excess factor $C$ (Evans et al., 2018; Riello et al., 2021) gives the
flux excess in the $G_{\mathrm{BP}}-G_{\mathrm{RP}}$ color relative to the $G$
band flux. It is recommended to use the corrected $C$ (denoted as $C^{*}$) as
given in Riello et al. (2021), which reduces the color dependence. Using these
parameters, we apply the following quality criteria to the photometry.

$\displaystyle
G_{\mathrm{err}}<0.006,\,\,G_{\mathrm{RP,err}}<0.02,\,\,G_{\mathrm{BP,err}}<0.1,\,\,C^{*}<0.3$
(47)

These cuts reduce the number of Sco-Cen members from 12,724 with complete
phototmetric information to 11,132 (leaving 87%), mainly reducing the number
of the fainter low-mass stars. If not applying any quality criteria, a
significant number of sources would be shifted toward older ages, only because
of unreliable Gaia colors.

Figure 14 shows two isochrones. The dashed line shows a 25 Myr isochrone from
Baraffe et al. (2015) (BHAC15444444http://perso.ens-
lyon.fr/isabelle.baraffe/BHAC15dir/) for Gaia DR3/EDR3 passbands.
Additionally, the solid line shows a 25 Myr isochrone from
PARSEC454545http://stev.oapd.inaf.it/cgi-bin/cmd; assuming solar metallicity
(metal fraction $z=0.0152$) and no extinction. for Gaia DR3/EDR3 passbands
(e.g., Bressan et al., 2012; Chen et al., 2014, 2015; Marigo et al., 2017;
Riello et al., 2021), including the upper–main-sequence (UMS), which is
missing in the BHAC15 models. We use both models, since BHAC15 models deliver
a better representation for low-mass stars.

To get a measure for the contamination from older sources (older than the
expected $\sim 20$ Myr), we select sources to the left of the two 25 Myr
isochrones, allowing for random scatter around the 20 Myr isochrone (in
particular young stars often show higher variability than main-sequence
stars). Additionally, we do not consider sources at the UMS, since there the
trend reverses (younger sources are to the left of the UMS). Hence, we apply a
cut at $G_{\mathrm{abs}}>3$ mag, only including fainter sources to test the
contamination from older sources.

The combined conditions deliver 737 candidate contaminants (and 10,395 young
Sco-Cen candidates) out of 11,132 sources with applied photometric quality
criteria, hence about 6–7%, which are possible contaminants from older
populations within the SigMA clusters. Considering the chosen borders, we like
to stress that this separation can only be seen as a rough estimate. In
particular, we did not consider any possible contaminants in the UMS regime,
where it is more difficult to distinguish young members from older stellar
populations. Moreover, the chosen isochrone models have intrinsic
uncertainties, and any change in metallicity or extinction is ignored in our
test. Finally, we only examine sources in the HRD in Fig. 14 with the applied
photometric quality criteria from Equ. (47). Hence, we can not make any
statement for sources with inferior photometric quality, which often also
suffer from higher astrometric uncertainties. Consequently, such sources could
also have a higher probability of having a wrong cluster membership, solely
due to the generally larger scatter in various dimensions.

To better understand the influence of the quality criteria, we repeat the
selection of old star contamination with different photometric quality cuts.
First, we consider the case of applying no cuts, using the 12,724 source with
entries in all three Gaia bands, delivering a contamination fraction from
older sources of about 12%. Next we apply somewhat looser and also stricter
quality criteria than given in Equ. (47) as follows, first showing the looser
cuts (48), and next the stricter cuts (49):

$\displaystyle
G_{\mathrm{err}}<0.01,\,G_{\mathrm{RP,err}}<0.045,\,G_{\mathrm{BP,err}}<0.25,\,C^{*}<0.5$
(48) $\displaystyle
G_{\mathrm{err}}<0.004,\,G_{\mathrm{RP,err}}<0.015,\,G_{\mathrm{BP,err}}<0.05,\,C^{*}<0.3$
(49)

With the looser cuts we get a contamination fraction of about 10% and with the
stricter cuts about 4%. It gets clear that the fraction of sources to the left
of the chosen isochrones decreases significantly (from about 12% to 4%) when
using superior photometry, which indicates that many sources indeed get
erroneously shifted to older ages when not considering the influence of
photometric errors. In conclusion, we estimate the contamination fraction from
older populations or field stars to be between 4–10%.

Additionally, we could investigate the influence of the membership stability
as delivered by the SigMA algorithm (Sect. 3.6). Figure 8 shows the influence
of different stability cuts on the “old star contamination fraction” when
using the quality criteria from Equ. (47). It can be seen that for low
stability ($\lesssim$11%) there is also a significant increase in the
contamination fraction. When only using sources with stability ¿ 11%, we would
get a “old star contamination fraction” in the range of about 2–9% for the
cases of no photometric cuts to strict photometric cuts. Hence, 2–4% can be
considered as the lower limit for contamination from older stellar populations
(or field stars), while about 10% is likely the upper limit. In conclusion,
the majority of the SigMA selected Sco-Cen members (likely more than 90%) are
sources with young Sco-Cen ages and therefore likely not contaminants from
interloping older populations or field stars.

###  J.2 Estimating the fraction of sub-stellar sources

To get an estimate of sub-stellar sources in our sample (brown dwarf
candidates), we use a 0.08 M⊙ iso-mass line in Fig. 14 (right panel), which is
extracted from BHAC15 models using ages from $10^{4}$ to $10^{10}$ yr, to get
a wide range. We select sources below 0.08 M⊙, which is defined as the
approximate hydrogen-burning limit (e.g., Baraffe et al., 1998; Burrows et
al., 2001; Freytag et al., 2010, 2012; Dieterich et al., 2014). The
uncertainties at the low-mass regime make this selection only a rough
estimate, in combination with the uncertainties of the stellar models (e.g.,
unknown metallicity, neglected extinction, different models give different
results). Additionally, all the uncertainties mentioned above in Sect. J.1
(e.g., chosen error or stability cuts) should be considered.

Using a cut 0.08 M⊙ and the quality criteria from Equ. (47) we find that there
are 1944 out of 10,395 (18.7%) potential substellar sources, when considering
only the younger sources from the middle panel in Fig. 14. This selection
indicates a fraction of sub-stellar objects of about 19% within the SigMA
clusters. If applying less strict error cuts (no cuts or Equ. (48)), the
fraction stays at about 19%, and if applying more strict error cuts (Equ.
(49)) the fraction decreases to about 12%. This is expected since more
conservative photometric error cuts affect in particular faint sources.
Changing the stability criteria does not influence these different fractions
significantly, since sources both in the stellar and substellar regime seem to
be affected almost equally. Concluding, there are about 19% of brown dwarf
candidates in the SigMA selected Sco-Cen sample.

##  Appendix K Auxiliary tables and figures

In Table 8 we give an overview of the contents of the source catalog
containing all Sco-Cen members as selected in this work, including labels for
cluster membership. The full version of the table is available online as a
machine-readable version.

We provide three additional tables, giving an overview of the literature
comparisons between the SigMA clusters and other Sco-Cen samples. In Table 9
we compare to Damiani et al. (2019) and Žerjal et al. (2021), in Table 11 we
compare to Squicciarini et al. (2021) and Miret-Roig et al. (2022b), and in
Table 10 to Kerr et al. (2021). More details on the comparisons can be found
in Sect. 5.2.

Table 8: Catalog of the 13,103 Sco-Cen members labeled for cluster membership
as identified with SigMA.

Column name | Unit | Column description
---|---|---
dr3_source_id |  | The source ID from Gaia DR3.
SigMA |  |  SigMA membership label for each source, as defined in Table 3.
stability |  | Membership stability of each source between 0–100%
distance | pc | Distance derived from inverse of the parallax
e_d_upper | pc |  Upper 1 $\sigma$ distance uncertainty determined from the 68.3 percentile of the sampled $d$ distribution
e_d_lower | pc |  Lower 1 $\sigma$ distance uncertainty determined from the 68.3 percentile of the sampled $d$ distribution
X | pc | Heliocentric Galactic Cartesian coordinate, X-axis grows positive toward the Galactic center
Y | pc | Heliocentric Galactic Cartesian coordinate, Y-axis grows positive in direction of Galactic rotation
Z | pc | Heliocentric Galactic Cartesian coordinate, Z-axis grows positive toward the Galactic North-pole
e_X_upper | pc |  Upper 1 $\sigma$ X uncertainty determined from the 68.3 percentile of the sampled X distribution
e_Y_upper | pc |  Upper 1 $\sigma$ Y uncertainty determined from the 68.3 percentile of the sampled Y distribution
e_Z_upper | pc |  Upper 1 $\sigma$ Z uncertainty determined from the 68.3 percentile of the sampled Z distribution
e_X_lower | pc |  Lower 1 $\sigma$ X uncertainty determined from the 68.3 percentile of the sampled X distribution
e_Y_lower | pc |  Lower 1 $\sigma$ Y uncertainty determined from the 68.3 percentile of the sampled Y distribution
e_Z_lower | pc |  Lower 1 $\sigma$ Z uncertainty determined from the 68.3 percentile of the sampled Z distribution
v_alpha |  km s-1 |  Tangential velocity in the direction of $\alpha$
v_delta |  km s-1 |  Tangential velocity in the direction of $\delta$
e_v_alpha_upper |  km s-1 |  Upper 1 $\sigma$ $v_{\alpha}$ uncertainty determined from the 68.3 percentile of the sampled $v_{\alpha}$ distribution
e_v_delta_upper |  km s-1 |  Upper 1 $\sigma$ $v_{\delta}$ uncertainty determined from the 68.3 percentile of the sampled $v_{\delta}$ distribution
e_v_alpha_lower |  km s-1 |  Lower 1 $\sigma$ $v_{\alpha}$ uncertainty determined from the 68.3 percentile of the sampled $v_{\alpha}$ distribution
e_v_delta_lower |  km s-1 |  Lower 1 $\sigma$ $v_{\delta}$ uncertainty determined from the 68.3 percentile of the sampled $v_{\delta}$ distribution
v_alpha_LSR |  km s-1 |  Tangential velocity in the direction of $\alpha$ and relative to the LSR
v_delta_LSR |  km s-1 |  Tangential velocity in the direction of $\delta$ and relative to the LSR
v_ERV |  km s-1 |  Estimated radial velocity, given as $\hat{v}_{r}$ in Sect. 3.5.2

464646The full machine readable version of the catalog is given online, while
a column overview is given here. We list all relevant derived parameters.
Original Gaia DR3 parameters can be queried from the Gaia Archive by using the
Gaia dr3_source_id.

Table 9: Comparing the SigMA clusters with stellar group selections from
Damiani et al. (2019) and from Žerjal et al. (2021). Only those SigMA clusters
which have cross-matches with either of the two literature samples are given
here.

SigMA |  Name (SigMA) |  Nraa$a$Number of sources from this work, for a direct comparison with the number of cross-matches given in brackets in Cols. 4–5. |  Matches with DDP19bb$b$The DPP19 group shortcuts are given for eight compact clusterings (UCL-1, UCL-2, UCL-3, Lupus 3, LCC-1, US-far, US-near), four diffuse populations (D1, D2a, D2b, US-D2), and sources that have not been assigned to any of these groups (N). The number in brackets gives the number of matches with the respective SigMA cluster. See details in Sect. 5.2.1. |  Matches with ZIC21cc$c$ZIC21 report eight sub-groups in Sco-Cen (C, E, D, F, G, T, A, U) and two additional older groups (H, I; IC 2602 and Platais 8), while there are no matches of H or I with the 37 SigMA clusters. Again, the number of matches is given in brackets. The groups contain the following numbers of sources in ZIC21: C-USco 1432, E-USco-multi 1483, D-UCL-V1062-Sco 506, F-UCL-V1062-Sco 273, G-UCL-East 1713, T-UCL-West 1057, A-LCC-North 1234, and U-LCC-South 487. See further details in Sect. 5.2.5.
---|---|---|---|---
1 | rho Oph | 535 | US-f(3)US-n(316)US-D2(68)N(20) | C-USco(405)E-USco-multi(13)
2 | nu Sco | 150 | US-n(62)US-D2(65)N(1) | C-USco(127)E-USco-multi(1)
3 | delta Sco | 691 | US-f(33)US-n(74)D1(88)D2a(8)D2b(2)US-D2(387)N(6) | G-UCL-East(3)C-USco(489)E-USco-multi(24)
4 | beta Sco | 285 | US-f(58)US-n(17)US-D2(152)N(8) | C-USco(189)E-USco-multi(19)
5 | sigma Sco | 544 | US-f(180)US-n(3)D1(14)D2a(3)US-D2(227)N(7) | G-UCL-East(2)C-USco(105)E-USco-multi(184)
6 | Antares | 502 | US-f(67)US-n(40)D1(19)D2a(5)US-D2(249)N(29) | C-USco(78)E-USco-multi(252)
7 | rho Sco | 240 | US-f(14)US-n(2)D1(159)D2a(3)US-D2(10)N(15) | G-UCL-East(48)C-USco(1)E-USco-multi(7)
8 | Sco-Body | 373 | D1(2)D2a(221)US-D2(34)N(7) | E-USco-multi(291)
9 | US-fg | 276 | D1(188)D2b(1)US-D2(7)N(12) | G-UCL-East(16)E-USco-multi(29)
10 | V1062-Sco | 1029 | UCL-1(554)D1(11)D2a(222)D2b(4)N(10) | D-UCL-V1062-Sco(499)F-UCL-V1062-Sco(228)G-UCL-East
11 | mu Sco | 54 | UCL-1(36)D1(2)D2a(5) | D-UCL-V1062-Sco(7)F-UCL-V1062-Sco(30)
12 | Libra-S | 71 | D1(1)D2a(8)D2b(13)US-D2(32) | E-USco-multi(4)
13 | Lup 1-4 | 226 | LupIII(67)D2a(47)D2b(65)N(4) | G-UCL-East(6)T-UCL-West(1)E-USco-multi(109)
14 | eta Lup | 769 | UCL-3(3)D1(549)D2a(43)D2b(10)US-D2(1)N(14) | G-UCL-East(419)E-USco-multi(15)
15 | phi Lup | 1114 | UCL-3(48)D1(627)D2a(62)D2b(148)N(38) | G-UCL-East(652)T-UCL-West(28)E-USco-multi(17)
16 | Norma-N | 42 | D1(1)N(6) |
17 | e Lup | 516 | D1(319)D2a(18)D2b(80)N(5) | G-UCL-East(349)T-UCL-West(15)
18 | UPK 606 | 131 | UCL-2(50)D1(2)D2b(57)N(1) | G-UCL-East(54)T-UCL-West(9)
19 | rho Lup | 246 | D1(17)D2b(189)N(2) | A-LCC-North (45)G-UCL-East(10)T-UCL-West(110)
20 | nu Cen | 1737 | UCL-2(2)D1(54)D2a(12)D2b(1270)US-D2(3)N(50) | A-LCC-North (70)G-UCL-East(116)T-UCL-West(790)E-US
21 | Cen-Far | 99 | D2b(41)N(1) |
22 | sig Cen | 1805 | LCC-1(1)D1(45)D2b(1384)N(43) | A-LCC-North (1077)U-LCC-South(56)T-UCL-West(66)
23 | Acrux | 394 | LCC-1(89)D1(11)D2b(242)N(4) | A-LCC-North (25)U-LCC-South(316)
24 | Musca-fg | 95 | D2b(35)N(2) | U-LCC-South(76)
25 | eps Cham | 39 |  | U-LCC-South(25)
26 | eta Cham | 30 |  | U-LCC-South(3)
27 | B59 | 32 | D2a(1)N(20) | E-USco-multi(2)
28 | Pipe-North | 42 |  | E-USco-multi(38)
29 | tet Oph | 98 | D2a(37)US-D2(6)N(2) | E-USco-multi(87)
30 | CrA-Main | 96 |  | E-USco-multi(2)
32 | Sco-Sting | 132 | D1(22)D2a(1) |
35 | L134/L183 | 24 |  | E-USco-multi(16)

474747

Table 10: Comparing the SigMA clusters with Kerr et al. (2021) clusters toward
Sco-Cen. Only those SigMA groups which have matches with KRK21 are given here,
while Oph-NorthFar is the only SigMA cluster without matches.

SigMA |  Name (SigMA) |  Nraa$a$Number of sources from this work, for a direct comparison with the number of cross-matches as given in brackets in Cols. 4–6. |  TLCbb$b$Col. 4 lists the TLC group labels if there are matches with SigMA, with the number of cross-matches in brackets. There have been only matches with the TLC groups 4, 21, and 22. |  EOMcc$c$Col. 5 lists the EOM group labels if there are matches with SigMA, with the number of cross-matches in brackets, while each EOM represents a sub-clustering within a lower level TLC group. |  LEAFdd$d$The letters in Col. 6 correspond to LEAF sub-groups, with the number of cross-matches in brackets, wile a LEAF group represents a sub-clustering within the lower level EOM group. |  Name (KRK21)ee$e$Group names from KRK21, if there is a significant overlap with SigMA. Only the (sub)group with the most significant number of cross-matches is given (in few cases more than one), as apparent from the numbers in brackets in Cols. 5 & 6.
---|---|---|---|---|---|---
1 | rho Oph/L1688 | 535 | 22(308) | 17(272) | I(109) | US-I-rho Oph
2 | nu Sco | 150 | 22(91) | 17(90) | E(54) | US-E
3 | delta Sco | 691 | 22(414) | 17(378) | F(1)H(102)I(1) | US-H
4 | beta Sco | 285 | 22(167) | 17(137) | G(29) | US-G
5 | sigma Sco | 544 | 22(296) | 17(248) | A(1)C(17)D(22) | US-C/D
6 | Antares | 502 | 22(292) | 17(239) | A(1)B(11)C(1)F(24) | US-B/F
7 | rho Sco | 240 | 22(128) | 17(64) | A(9) | US-A
8 | Scorpio-Body | 373 | 22(193) | 16(12)17(45) |  | EOM-16/US
9 | US-foreground | 276 | 22(111) | 13(30) |  | EOM-13
10 | V1062-Sco | 1029 | 22(503) | 14(20)15(347) |  | LowerSco/EOM-14
11 | mu Sco | 54 | 22(28) | 15(23) |  | LowerSco
12 | Libra-South | 71 | 22(23) |  |  |
13 | Lupus-1-4 | 226 | 22(143) | 12(102) | A(46)B(14) | Lupus-IV/III
14 | eta Lup | 769 | 22(411) | 9(15)22(102)23(6) |  | EOM-9/22/23
15 | phi Lup | 1114 | 22(391) | 17(4)19(10)23(6) |  | EOM-19/23
16 | Norma-North | 42 | 22(4) |  |  |
17 | e Lup | 516 | 22(257) | 11(1)20(76)21(8) |  | EOM-20
18 | UPK606 | 131 | 21(1)22(53) | 11(32) |  | UPK606
19 | rho Lup | 246 | 22(123) | 21(2)25(10) |  | EOM-25
20 | nu Cen | 1737 | 22(573) | 11(3)21(2)24(108)26(14)27(2) |  | EOM-24
21 | Centaurus-Far | 99 | 21(36) | 3(30) |  | Cen-South
22 | sig Cen | 1805 | 22(987) | 25(1)26(26)27(421) | C(4)D(12)E(48) | EOM-26/LCC-D/E
23 | Acrux | 394 | 22(258) | 27(208) | B(1)C(96) | LCC-C-Crux S
24 | Musca-foreground | 95 | 22(65) | 27(46) | B(16) | LCC-B
25 | eps Cham | 39 | 22(23) | 27(20) | A(17) | LCC-A-eps Cha
26 | eta Cham | 30 | 22(18) | 18(17) |  | eta Cha
27 | B59 | 32 | 22(14) | 6(13) |  | Pipe
28 | Pipe-North | 42 | 22(19) |  |  |
29 | tet Oph | 98 | 22(49) | 10(28) |  | Theia67
30 | Corona Australis | 96 | 22(53) | 8(52) |  | CrA
31 | CrA-North | 351 | 22(207) | 7(1)8(195) |  | CrA
32 | Scorpio-Sting | 132 | 22(62) | 7(11) |  | EOM-7
33 | Chamaeleon-1 | 192 | 21(101) | 1(101) |  | Cha-1
34 | Chamaeleon-2 | 54 | 21(30) | 2(30) |  | Cha-2
35 | L134/L183 | 24 | 22(6) |  |  |
36 | Oph Southeast | 61 | 4(20) |  |  | Oph Southeast

484848 See more details in Sect. 5.2.2.

Table 11: Comparing the SigMA clusters with stellar group selections from
Squicciarini et al. (2021) and from Miret-Roig et al. (2022b). Only those
SigMA groups which have cross-matches with either of the two literature
samples are given here.

SigMA |  Name (SigMA) |  Nraa$a$Number of sources from this work, for a direct comparison with the number of matches as given in brackets in Cols. 4–5. |  Matches with SGB21bb$b$The SGB21 groups (G) are numbered from 1 to 8, and their diffuse population is given with D. The number of cross-matches is given in brackets. The eight groups in SGB21 are associated with the brightest star in each group as follows: G1–$i$ Sco; G2–$\nu$ Sco B; G3–$b$ Sco; G4–HD 144273; G5–HIP 77900; G6–HIP 78968; G7–HIP 79910; G8–HD 146467. See details in Sect. 5.2.6. |  Matches with MR22cc$c$Comparison to the seven groups in US as identified by MR22. Again, the number of cross-matches is given in brackets. See details in Sect. 5.2.7.
---|---|---|---|---
1 | rho Oph | 535 | G1(428)G2(2)G3(1)G4(1)G8(2)D(51) |  $\alpha$Sco(3)$\delta$Sco(9)$\nu$Sco(1)$\sigma$Sco(67)$\rho$Oph(370)
2 | nu Sco | 150 | G1(1)G2(110)G6(10)D(23) |  $\beta$Sco(2)$\delta$Sco(2)$\nu$Sco(110)$\sigma$Sco(22)
3 | delta Sco | 691 | G1(19)G2(2)G3(390)G4(9)G5(50)G6(2)G7(1)G8(2)D(136) |  $\alpha$Sco(27)$\beta$Sco(8)$\delta$Sco(410)$\nu$Sco(29)$\pi$Sco(38)$\sigma$Sco(90)
4 | beta Sco | 285 | G1(1)G4(141)G5(11)G6(46)D(57) |  $\alpha$Sco(2)$\beta$Sco(169)$\sigma$Sco(70)
5 | sigma Sco | 544 | G1(2)G4(2)G5(104)G8(3)D(377) |  $\alpha$Sco(268)$\delta$Sco(3)$\pi$Sco(7)$\sigma$Sco(163)
6 | Antares | 502 | G1(13)G3(5)G7(44)G8(33)D(313) |  $\alpha$Sco(290)$\pi$Sco(29)$\sigma$Sco(52)$\rho$Oph(37)
7 | rho Sco | 240 | D(168) |  $\alpha$Sco(3)$\pi$Sco(180)
8 | Sco-Body | 373 |  |  $\alpha$Sco(2)$\sigma$Sco(2)
9 | US-fg | 276 | D(13) |  $\alpha$Sco(1)$\pi$Sco(194)
12 | Libra-S | 71 |  |  $\sigma$Sco(1)
14 | eta Lup | 769 |  |  $\pi$Sco(12)
15 | phi Lup | 1114 | D(2) |  $\pi$Sco(3)
20 | nu Cen | 1737 |  |  $\sigma$Sco(1)

494949

Finally, we provide additional figures of the position and velocity spaces of
the individual 37 SigMA clusters. This allows a better appreciation of the
individual cluster source distribution in each parameter space. The Figs.
22–26 are constructed as follows. The cluster names are given in the left
panel ($l,b$ panel) of each row. Each column shows one of the six different
parameter spaces for one cluster. The parameter spaces are the same as in
Figs. 10–13, namely $l$ vs $b$ (deg), $X$ vs $Y$ (pc), $X$ vs $Z$ (pc), $Y$ vs
$Z$ (pc), $v_{\alpha}$ vs $v_{\delta}$ (km s-1), and $v_{\alpha,LSR}$ vs
$v_{\delta,LSR}$ (km s-1). These axes labels are given at the top of each
column. We like to highlight that Col. 5 (tangential velocities) shows a
larger velocity range than Col. 6 (tangential velocities relative to the LSR),
where the clusters actually occupy a smaller velocity space, hence show a
smaller velocity dispersion. All SigMA selected Sco-Cen members are plotted in
gray in all panels and the given cluster is over-plotted with red dots. See
also Figs. 10–13 for an alternative view of the 37 SigMA clusters, and the
interactive 2D and 3D versions.

Figure 22: Six parameter spaces, highlighting the individual clusters in red.
Shown are the clusters SigMA 1–8 (part of US). The gray background sources are
all SigMA-selected Sco-Cen members. Cluster names are given in the left panels
of each row. The used xy-axes are given as title at the top of each column,
and tick labels are only given in the bottom row. See also Figs. 10–13 and
text for more details.  Figure 23: Same as Fig.22, but for the clusters SigMA
9–16 (part of US and UCL). Figure 24: Same as Fig.22, but for the clusters
SigMA 17–24 (part of UCL and LCC). Figure 25: Same as Fig.22, but for the
clusters SigMA 25–32 (part of LCC, Pipe, and CrA). Figure 26: Same as Fig.22,
but for the clusters SigMA 33–37 (part of Cham and NE).