© ESO 2023

1. Introduction

Open clusters are gravitationally bound groups of a few tens of stars to over a thousand stars that formed from the same cloud of dust and gas. The similar distance, age, and metallicity of the cluster members means that a group of clusters at different ages can be viewed as snapshots of stellar evolution across time. This makes them an ideal laboratory to study the universality of the initial luminosity and mass functions (Salpeter 1955; Miller & Scalo 1979; Kroupa 1998, 2001; Scalo et al. 1998; Chabrier 2001; Bastian et al. 2010), nuclear timescales, and the evolution of the stellar multiplicity rates.

Over 1100 open clusters were already discovered in our Galaxy (Mermilliod 1995; Dias et al. 2002; Kharchenko et al. 2013) before the advent of the Gaia astrometric mission (Gaia Collaboration 2016). The arrival of the Gaia data releases (Gaia Collaboration 2017, 2018, 2021, 2023) have resulted in improved lists of cluster members as well as a discovery of new clusters (e.g., Cantat-Gaudin et al. 2018a; Castro-Ginard et al. 2019; Kounkel & Covey 2019; Sim et al. 2019; Jaehnig et al. 2021), stellar streams (Borsato et al. 2020), tails, and coronae around clusters (Meingast et al. 2021; Moranta et al. 2022).

Precise parallaxes and proper motions from Gaia have also improved the search for members in the kinematic parameter space. Different clustering techniques, in mostly five-dimensional space (i.e., missing radial velocities), have been used to evaluate membership in open clusters. Each of these techniques has their own biases and data quality cuts. For example, hierarchical density-based spatial clustering of applications with noise (HDBSCAN; McInnes et al. 2017) helped reveal a number of new clusters, associations, and comoving groups within 1 kpc (Kounkel & Covey 2019). Tarricq et al. (2022) applied the same method to a search for cluster members up to 50 pc around their centres in Gaia EDR3 data and detected vast coronae around almost all the clusters, as well as tidal tails in 71 open clusters. Castro-Ginard et al. (2018, 2019, 2020) used DBSCAN in Gaia DR2 and He et al. (2022) used it in Gaia EDR3 to discover numerous new clusters in the Milky Way. Unsupervised photometric membership assignment in stellar clusters (UPMASK; Krone-Martins & Moitinho 2014) was applied to the proper motions and parallaxes in the TGAS catalogue (Michalik et al. 2015) by Cantat-Gaudin et al. (2018b) to characterise 128 known open clusters in the solar neighbourhood. Cantat-Gaudin et al. (2018a) used the same method on a compiled list of thousands of known or putative clusters from the literature in the Gaia DR2 catalogue (Gaia Collaboration 2018) to obtain a list of members and cluster parameters for 1229 clusters. Recently, Perren et al. (2022) used the same algorithm for the analysis of the 25 most distant (> 9 kpc) catalogued open clusters. In contrast with (H)DBSCAN and UPMASK, Extreme Deconvolution Gaussian Mixture Models (Bovy et al. 2011) take astrometric uncertainties into account. With this method, Olivares et al. (2019) analysed the oldest cluster in the solar neighbourhood, Ruprecht 147; Price-Whelan et al. (2019) discovered a disrupting young open cluster in the halo of the Milky Way; and Jaehnig et al. (2021) provided membership lists for 431 open clusters. Perryman et al. (1998) developed a convergent point technique to retrieve the Hyades members in the HIPPARCOS data (Perryman et al. 1997). The membership of this cluster was later improved with the same method using the TGAS catalogue (Reino et al. 2018) and Gaia DR2 (Lodieu et al. 2019b). Lodieu et al. (2019a) improved memberships of the α Per, Pleiades, and Praesepe clusters with the same algorithm.

While stellar mass is the main parameter of the evolution of a star, its determination usually depends on models and their assumptions, except in the case of spectroscopic or eclipsing binaries. Observations have shown that the shape of the initial mass function appears universal among star-forming regions and open clusters, at least formally within the error bars (Bastian et al. 2010). Likewise, the luminosity function constructed from direct observables – stellar magnitudes that are model-free (unlike stellar masses) – reflects the cluster’s initial mass function, the history of its dynamical evolution, and its age. For example, when stars fuse their hydrogen fuel and become giants, they increase their luminosity. Similarly, the dynamical evolution causes mass segregation and the gradual evaporation of low-mass members due to the interactions with the external massive objects while moving around the Galaxy.

The main goal of this study is to produce catalogues of members for a volume-limited (distances less than 500 pc) and age-limited (younger than 1 Gyr) sample of open clusters outside the Galactic plane, taking advantage of the latest Gaia Data Release 3 (DR3; Gaia Collaboration 2023). We do not make any magnitude cuts to reach the faintest member objects, and we include bright HIPPARCOS stars that are missing parallaxes and proper motions in Gaia in order to reach maximal completeness. Our aim is to study the evolution of the luminosity functions, dividing our sample of clusters into three age groups. We will make all the catalogues and plots available publicly through an interactive website for the scientific community.

The paper is organised as follows. Section 2 introduces the main source of data and the clusters under scrutiny. Section 3 describes the procedure to infer membership of any source to a given cluster. We describe the catalogues of members in Sect. 4 and provide details about the website with interactive plots in Sect. 5. We analyse the luminosity functions of the clusters in Sect. 6 and conclude in Sect. 7.

2. Data

In this section, we describe the selection procedure to compile a list of clusters, which is the main focus of this study. Then, we continue with the summary of the data to evaluate individual stellar membership for each cluster. This includes stellar positions, distances, proper motions, and radial velocities as well as their mass estimates and supplementary photometric measurements.

2.1. Selection of clusters

Our preliminary list of clusters is based on the known clusters from Kharchenko et al. (2013), Cantat-Gaudin et al. (2018a), Tarricq et al. (2021), and Dias et al. (2021) as well as the WEBDA website¹. To avoid crowded and highly reddened regions of the Milky Way, we only selected clusters that are located outside the Galactic plane (|b|> 10 deg). We excluded clusters beyond 500 pc due to the increasing parallax uncertainties of their members. We set an upper literature age limit of 1 Gyr since our objective is to focus on young clusters, and there are not many older open clusters within 500 pc. For example, in the list of Friel (1995), we find one such cluster – NGC 752 at ∼430 pc with the age of 1.61 ± 0.03 ± 0.05 Gyr (Sandquist et al. 2023). However, we added well-known clusters such as IC 2602, IC 2391, and α Persei. These three clusters are within b = 10 deg but have been included because they are nearby and young and have been studied in detail in the literature (e.g., Randich et al. 2001; Lodieu et al. 2005, 2012, 2019a; Dobbie et al. 2010; Nisak et al. 2022). The following clusters/regions Sigma Orionis (Garrison 1967; Barrado y Navascués et al. 2001; Caballero 2018), IC 348 and NGC 1333 (Luhman et al. 2016; Olivares et al. 2023), Feigelson 1 (also known as ϵ Cha Cluster; Dickson-Vandervelde et al. 2021), Platais 2 and 5 (Platais et al. 1998), UBC 19 (Castro-Ginard et al. 2018), Melotte 227 (Epstein 1968), Collinder 65 (Yen et al. 2018), Collinder 70 (Caballero & Solano 2008), and ASCC 100 (Kharchenko et al. 2005) were removed because they are not isolated, and many of them are found in very young complex regions (e.g., Taurus, Orion, Upper Scorpius and Rho Ophiuchi) that display an extended structure that is often non-trivial to isolate from the surrounding environment. This selection left us with 49 clusters, which are presented in a list in Table 1 with their equatorial and Galactic coordinates, distances, and numbers of initial candidates with and without radial velocities. All clusters are unique, and a few of them are located relatively close to each other. In rare cases, their members might partially overlap. We explain this issue in Sect. 4.1.

2.2. Stellar positions and velocities

We downloaded the catalogue with astrometric measurements of the Milky Way stars from the Gaia DR3 database (Gaia Collaboration 2016, 2023). We limited the selection to stars with parallax uncertainties better than 20% but applied no constraints on any other parameters except for parallax to encompass a parameter space large enough to include all possible candidate members for all clusters. Due to the fact that stars with 1/parallax just beyond 500 pc have relatively large parallax uncertainties (typical parallax_error of 12% as opposed to 9% for stars between 400 and 500 pc and 4% for stars between 100 and 200 pc) and can thus still be consistent with the distance of a cluster at 500 pc, we set the parallax limit to 1.6 mas, which corresponds to distances up to 625 pc. Our input catalogue therefore includes right ascension, declination, parallax, proper motions, and associated uncertainties for 44 million stars.

To avoid the minimal risk of negative parallaxes for our most distant stars, parallax uncertainties needed to be taken into account in the distance determination (e.g., see Fig. 3 in Luri et al. 2018). Stellar distances were obtained from Bailer-Jones et al. (2021), who used parallaxes from Gaia EDR3² (Gaia Collaboration 2021) with a prior constructed from a 3D model of the Galaxy. In particular, we worked with their geometric distances r_med_geo and their lower and upper limits r_lo_geo and r_hi_geo, respectively. However, their Galactic prior with a scale length varies across the sky. If this scale distance is far away and the parallax uncertainty is high, this may result in a mode distance considerably away from the true one. For parallaxes with an error of 20%, this translates to a distance error of 10 − 20% for a cluster at 500 pc, as shown in Fig. 1 of Olivares et al. (2020). Despite the exclusion of the clusters in the prior of Bailer-Jones et al. (2021), the effect on the derived distances is negligible. The difference with the inverted parallax distance is minimal, typically less than 1 pc for stars closer than 400 pc and 4 pc for stars beyond this limit.

The latest Gaia DR3 catalogue (Katz et al. 2023) contains radial velocities for 33 million stars. The overall radial velocity completeness ratio in our sample is 14%, and 44% of the stars with radial velocity have uncertainty better than 10%. The completeness varies over spectral types, from 92% for G stars to 2% for M stars, with typical uncertainties of 1.44 km s⁻¹ for G and 4.60 km s⁻¹ for M stars, respectively.

While the majority of our clusters are composed of relatively faint stars, the nearest ones contain stars visible to the naked eye. Due to the saturation limit at the bright end, 20% of stars brighter than magnitude three do not have an entry in Gaia EDR3 (Fabricius et al. 2021). Additionally, a fraction of bright sources, for example Alcyone, the brightest star in the Pleiades cluster, lack measurements of parallaxes and proper motions. Fortunately, the saturation limit in Gaia’s predecessor, the HIPPARCOS space telescope, was much brighter. We therefore took parallaxes, proper motions, and their uncertainties from the new reduction of the HIPPARCOS measurements (van Leeuwen 2007). We used the cross-match with Gaia DR3 that is available in the Gaia archives and added the measurements to our list in the following way: If a Gaia source brighter than G = 6 mag lacked parallax, we inserted the HIPPARCOS data for this star. Furthermore, there is a small number of HIPPARCOS stars that do not have a counterpart in Gaia. We appended such stars brighter than H_p = 6 mag to our list and noted the underlying assumption of a perfect cross-match. This might not be true in case of binary stars unresolved in HIPPARCOS. We checked the offsets of the HIPPARCOS catalogue with respect to Gaia for the entire overlapping sample (100 000 stars) and provided the median values. There is an offset of 0.1087 mas in the parallax. This is below the typical parallax uncertainty of 0.3600 mas in HIPPARCOS but greater than a typical Gaia parallax error of 0.0662 mas in this magnitude range. The offsets of proper motions are 0.0972 and −0.0165 mas yr⁻¹ in right ascension and declination, respectively. These values are not only smaller than HIPPARCOS uncertainties (0.3400 and 0.2900 mas yr⁻¹, respectively) but also smaller or comparable with the Gaia uncertainties (0.0658 and 0.0606 mas yr⁻¹, respectively).

2.3. Stellar masses

One of the steps in our membership estimation algorithm described below in Sect. 3.2 is the determination of the systemic velocity for each cluster (i.e., the mean velocity of the centre of mass). This means that we needed to estimate individual stellar masses. Approximate mass computation is sufficient for our needs, assuming a spherical symmetry of the clusters and the number of members large enough to provide a robust estimate. This presumption is supported by Perryman et al. (1998) and Reino et al. (2018), who concluded that the chosen weighting scheme in mass determination does not have a large influence on the resulting barycentre. Similarly, assuming an isotropic distribution of binaries in clusters, their mass determination technique does not have an impact on the barycentre (Lodieu et al. 2019b). Unresolved binary stars in our catalogue thus have masses determined in the same way as single objects. We note that this does, however, skew the total mass of the cluster.

Mass estimation using the luminosity–mass or magnitude–mass relations requires a particular model for each cluster to meet its age and metallicity. To simplify the process, we adopted a colour–mass relation that is valid for a wide range of masses and ages below 1 Gyr. A comparison of different isochrones with ages between 50 Myr and 1 Gyr in Fig. 1 shows no large variation that would exceed our requirements of rough mass estimates. We thus chose the models at 200 Myr.

Fig. 1. Comparison of isochrones for a range of ages from 50 Myr to 1 Gyr for the BHAC15 (dashed lines) and PARSEC models (solid lines). The polynomial fit used as a colour–mass relation is shown in magenta. While we only plot the relation for G − RP > −0.2, the fit extends to bluer colours in order to cover the most massive stars in our input catalogue.

Our mass estimation is based on the colour–mass relation from two models that are joined together at the point of their best overlap at G − R_P = 0.8 mag. The vast majority of cluster member candidates have positive G − R_P colours. We took the isochrones with solar metallicity from Baraffe et al. (2015) for G − R_P > 0.8 mag and PARSEC isochrones (the PAdova and TRieste Stellar Evolution Code; Bressan et al. 2012; Chen et al. 2014, 2015; Tang et al. 2014; Marigo et al. 2017; Pastorelli et al. 2019, 2020) for stars with for G − R_P < 0.8 mag. We fit an 11th order polynomial to these two isochrones to provide a relation that can be applied to all stars in our list. To avoid the extrapolation problems outside the models, we set the lower mass limit to 0.05 M_⊙. For dwarf stars with G − R_p < 0, we followed a similar procedure but took the 50 Myr isochrones.

The masses of HIPPARCOS stars are based in the B − V–mass relation from an updated table³ originally published by Pecaut & Mamajek (2013, hereafter we cite the table as P13). We fit a 21st order polynomial to this relation and estimated the masses of HIPPARCOS stars using their B − V colours. Again, we reinforced all masses to be greater than 0.05 M_⊙. We set masses for a few individual stars above 5 M_⊙ manually from Tetzlaff et al. (2011) in mass determination with both colours.

Some of our older clusters contain white dwarfs. Typically, about 4% of objects in the input catalogue are found in the white dwarf region of the colour–magnitude diagrams. We defined all stars fainter than M_G = 10(G − R_p)+5 as white dwarfs and assigned them a typical mass of a DA white dwarf (hydrogen lines), 0.59 M_⊙ (Kepler et al. 2007), which represent 80% of all white dwarfs. This sample partially overlaps with a list of white dwarfs identified in Gaia EDR3 (Gentile Fusillo et al. 2021).

2.4. Field population

We compared the luminosity functions and mass distributions of the clusters with the field population. We thus prepared a sample of field stars in the following way: We selected stars from the Gaia DR3 catalogue within 25 pc with small parallax errors by requiring parallax > 40 mas and parallax_error < 8 mas in the Gaia query. We excluded stars with either missing or noisy photometry in the G and R_p filters. In particular, we required phot_rp_mean_flux_over_error and phot_g_mean_flux_over_error to be more than 100. These selection criteria resulted in 4760 stars that we proclaim as the field population. We computed the total mass within the given distance to the Sun. This cumulative mass distribution – the total mass within the given distance to the Sun – reflects a uniform distribution with no over-densities. We tested this claim by taking the same sample and simulating a uniform distribution of the stars within a sphere with the same volume. In practice, we kept the original locations of the stars but shuffled their masses. The cumulative mass distributions of both samples are practically the same and correspond to a mass density of 0.033 M_⊙ pc⁻³ (see Fig. 2). This sample was assumed to be complete down to L4 objects, while the completeness for the most distant clusters in our list does not exceed mid-M dwarfs. For comparison, we expected it to be complete to M9/L0 in the Hyades (see our reasoning in Sect. 2.3) and ∼M7 in the Pleiades (Lodieu et al. 2019a,b). This means that membership lists for most of the clusters in this work lack brown dwarfs. We note that the comparison with the field is not ideal for more distant clusters, since these samples might not come from the same environment and are affected by completeness and contamination issues.

Fig. 2. Tidal radius determined by a cross-section between the cumulative radial mass distribution of the cluster and M(r) from Eq. (8) (blue line). The tidal radius for single stars (black) represents a minimal value, while a synthetic sample where all binaries have an equal mass (red) indicates the upper limit. The difference between these tidal radii is about 10%. The field distribution was added for comparison (purple and green for the actual and simulated field, respectively; see Sect. 2.4).

2.5. Photometric data

We cross-matched our catalogues of cluster members (see Sect. 4.1) with the photometric surveys described in Table 2 to complement the Gaia astrometry and photometry of all the stars. We used cross-matches provided in the Gaia database for 2MASS, AllWISE, and Pan-STARRS. For SDSS, UKIDSS, and VISTA, we performed a coordinate cross-match in TopCat within the radius of 3 arcsec.

3. Membership determination

The membership determination technique used in this work exploits the fact that stars in open clusters display a common location and velocities, as opposed to the uniform distribution of field stars in a Cartesian parameter space. The method takes advantage of the full 6D space (i.e., stellar positions and velocities) and the 5D space, when radial velocities are missing. It is based on the comparison of the transverse and radial velocities of a candidate star with the values expected for a member of a cluster at the location of a star. The expected values were estimated from the cluster’s systemic velocity and the position of the candidate star with respect to the cluster. This approach was developed by Perryman et al. (1998) to retrieve the members of the Hyades cluster from the HIPPARCOS data with a later application to Gaia DR1 (Reino et al. 2018). Lodieu et al. (2019a,b) utilised the same method on the Gaia DR2 data to improve the memberships of the Hyades, α Persei, the Pleiades, and Praesepe clusters. The specifics of this technique are described in the aforementioned papers. In this section, we describe the input data for each cluster and follow it with a summary of our membership selection algorithm, which is an improved version of the one used by Lodieu et al. (2019b).

The terms we use to describe different samples of stars in the next sections are the following: ‘Preliminary members’ are used only in the membership determination method but are not part of the results; see Sect. 3.1.1. ‘Bona fide members’ are the most reliable members of the cluster, as these are stars within 1 tidal radius of the cluster centre (see Sect. 4). Finally, ‘cluster members’ include stars up to 3 tidal radii.

3.1. Input catalogues for each cluster

The input catalogues for each cluster consist of two parts, the list of preliminary candidate members and the list of objects in the cluster region. The purpose of the catalogue with the preliminary candidate members is the computation of the initial systemic velocity and barycentre of the cluster, which were later refined in the iterative loop described in Sect. 3.2. As for the list of objects in the cluster region, it encompasses a large volume around the cluster in order to include all potential candidate members. This list is a simple subset of our entire Gaia catalogue with 44 million stars, outlined in Sect. 2.2, and its purpose is to speed up the membership selection algorithm. We added Cartesian coordinates for stellar velocities and positions as required in some steps of the membership determination algorithm. In the following subsections, we describe the preparation process for both datasets and the coordinate transformations.

3.1.1. Preliminary candidate members

The catalogue of preliminary candidate members was prepared manually for each cluster to ensure a coherent catalogue. The goal was to find enough bona fide members to estimate the initial systemic velocity and barycentre of a cluster in a robust way. This means that our goal was to have a low contamination rate and that we were not concerned with the completeness level.

In the first step, we limited the data to the region of the sky centred on a cluster with coordinates reported in the literature. This region typically spans an angle in the sky that corresponds to a diameter of 60 pc at the distance of the cluster centre. Since the typical radius of a bound cluster is less than 10 pc, this diameter should include all bona fide members. We note that not all clusters appear as clear over-densities in the sky. Therefore, we identified each of them as an overdensity in the proper motion space with the help of the literature values. We noted an occasionally complex and crowded nature for the proper motion space with asymmetric and sometimes partially overlapping over-densities, due to the rather large selected area in the sky. We refer the reader to Appendix D for a description of clusters with peculiarities. The cut in the proper motion space is circular with a typical radius of 0.7 mas yr⁻¹. We manually adjusted this radius after a visual inspection to cover the obvious overdensity but stayed conservative in order to limit the contamination to the minimum.

To further help eliminate interlopers, we removed outliers in the distance versus total proper motion space. This step revealed a clear, narrow sequence in the colour–magnitude diagram that represents a coeval population of stars.

The selection of initial candidate members was an iterative process because some individual and less populated clusters appeared as somewhat sparse groups of stars in the sky with only a small number of members. We modified the parameters of the cuts as needed to optimise the selection. Finally, we manually removed the outliers in the colour–magnitude diagram, taking out all white dwarfs and (sub)giants as well as any stars with radial velocities significantly different from the average value in the sample (e.g., more than 20 km s⁻¹). This left us with a sample of stars that are very likely members of the cluster. The number of initial candidate members and the number of stars with a radial velocity in the initial sample is listed in Table 1. The fraction of stars with a radial velocity in the initial sample has increased thanks to Gaia DR3 and varies from 20 to 70%, depending on the spectral type, with a typical value of 30%. These numbers enable good statistics for determining a robust initial systemic velocity of the cluster.

3.1.2. Objects in the cluster region

To ensure the volume of each of the input catalogues was large enough to enclose all potential candidate members, we performed wide and simple cuts to our fully downloaded Gaia catalogue (Sect. 2.2). In particular, we only made cuts in the right ascension and declination in order to limit the catalogue to the cluster centre and its vicinity in the sky. We took the right ascension and declination from the literature and kept all the stars that are located within an angular span and have a radius approximately three times the estimated radius of the cluster. We visually compared the distribution of this input catalogue and the initial list of candidate members in the sky to make sure the volume is large enough. We did not make any additional cuts in the distance in order to mitigate the risk of missing potential members with large distance uncertainties (except for the distance cuts explained in Sect. 2.2). The typical number of objects in the cluster region varies from a few 10 000 to a few 100 000 stars, depending on the part of the sky and distance from the Galactic plane.

3.1.3. Coordinate transformations

We transformed the spatial observables α and δ as well as the distance d (in parsecs) from the input catalogue of 44 million stars to the Galactic Cartesian system of physical distances with the centre in the Sun. The position of the ith star is described by b_i = (b_xi,  b_yi,  b_zi). The velocity vector v_i = (v_xi,  v_yi,  v_zi) was computed by the transformation of the observed transverse and radial velocities into the equatorial coordinate frame with

$$\begin{array}{c}\hfill \left(\begin{array}{c}{v}_{\mathit{xi}}\\ v_{\mathit{yi}}\\ v_{\mathit{zi}}\end{array}\right)={\mathbf{R}}_i\left(\begin{array}{cc}& V_{\alpha \ast i}\\ & V_{\delta i}\\ & V_{\mathit{Ri}}\end{array}\right).\end{array}$$(1)

The radial velocity is denoted by V_Ri. The transverse velocities V_α * i and V_δi are defined as

$$\begin{array}{cc}\hfill V_{\alpha \ast i}& ={\mu }_{\alpha \ast i}{A}_v/{\pi }_i,\mathrm{and}\hfill \end{array}$$(2)

$$\begin{array}{cc}\hfill V_{\delta i}& ={\mu }_{\delta i}{A}_v/{\pi }_i,\hfill \end{array}$$(3)

where μ_δi and μ_α * i  =  μ_αi cos δ_i are proper motions and ${\widehat{\pi }}_i$ is a parallax computed from the distance r_med_geo. The constant A_v = 4.74047 km yr s⁻¹ was used to convert the units and express V_α * i and V_δi in kilometres per second. The rotation matrix R_i is given by

$$\begin{array}{c}\hfill {\mathbf{R}}_i=\left(\begin{array}{ccc}-sin{\alpha }_i& -sin{\delta }_{i}cos{\alpha }_i& cos{\delta }_{i}cos{\alpha }_i\\ cos{\alpha }_i& -sin{\delta }_{i}sin{\alpha }_i& cos{\delta }_{i}sin{\alpha }_i\\ 0& cos{\delta }_i& sin{\delta }_i\\ \end{array}\right).\end{array}$$(4)

3.2. Membership determination algorithm

Stellar clusters appear as over-densities in the spatial and kinematic parameter space. The process of identifying cluster members thus operates in the full 6D parameter space. The algorithm takes the kinematic properties of the stars in the cluster region and compares them with the systemic velocity of the group. However, the members of a cluster must be known in order to determine its systemic properties. To break this cyclic dependency, we used a two-step procedure that iteratively improves both the membership list and the systemic cluster properties. In the first iteration, we predefined the systemic cluster properties from the initial list of bona fide members and searched for the candidate members. In the next iteration, the new membership list was used to improve the cluster properties and in turn update the catalogue of candidate members. This iteration was repeated until converged. We proceed next with the details of each iterative step, convergence criterion, and the treatment of stars with missing radial velocities.

3.2.1. Systemic properties of the cluster

The systemic properties include the barycentre, systemic velocity, and half-mass radius of the cluster. The computation of the barycentre and velocity was done in the Cartesian coordinate system by

$$\begin{array}{c}\hfill {\mathit{b}}_{\mathrm{c}}=\frac{\mathrm{\Sigma }{m}_i{\mathit{b}}_i}{\mathrm{\Sigma }{m}_i},{\mathit{v}}_{\mathrm{c}}=\frac{\mathrm{\Sigma }{m}_i{\mathit{v}}_i}{\mathrm{\Sigma }{m}_i},\end{array}$$(5)

where b_c denotes a barycentre and v_c marks the systemic velocity vector of the cluster. The weights m_i are stellar masses (see Sect. 2.3). While the assumption of the entirely single nature of the stars keeps the total cluster mass at its lower limit, it does not affect the barycentre significantly. This is because the contributions of the binary stars cancel out due to their expected isotropic distribution.

We note that radial velocity measurements are not available for all the stars. Therefore, the list of stars used to compute the systemic velocity is a subsample of the catalogue employed in the determination of the barycentre.

The initial cluster systemic velocity and barycentre were determined from the preliminary list of cluster members (Sect. 3.1.1). In the later iterative steps, we took candidate members within the 3σ confidence interval (see Sect. 3.2.2 for definition) and the half-mass radius of the cluster to refine the systemic cluster parameters.

3.2.2. Search for cluster members

Once the systemic cluster properties were determined, we compared the kinematic properties of each object in the input catalogue with those of the cluster. In particular, we computed the expected transverse and radial velocities for a candidate member of a cluster at the location of the candidate star:

$$\begin{array}{c}\hfill \left(\begin{array}{c}{V}_{\alpha \ast i}^e\\ V_{\delta i}^e\\ V_{\mathrm{R}i}^e\end{array}\right)={\mathbf{R}}_i^{-1}\left(\begin{array}{c}{v}_x\\ v_y\\ v_z\end{array}\right).\end{array}$$(6)

The matrix R_i is defined in Eq. (4). The velocity vector in the right-hand side corresponds to the cluster systemic velocity v_c = (v_x,  v_y,  v_x).

To evaluate the membership of a star, we first computed the difference between its observed and expected transverse and radial velocities and denoted it with a vector z_i. To account for the measurement uncertainties and the correlations between the observables, we followed Perryman et al. (1998) and prepared a covariance matrix for the observed values and another matrix for the expected values of the star. The sum of both matrices Σ represents the confidence region and was used to compute the scaled distance between the observed and expected velocity vector by

$$\begin{array}{c}\hfill c={\mathit{z}}^T{\mathrm{\Sigma }}^{-1}\mathit{z}.\end{array}$$(7)

Assuming a normal distribution, c² follows an χ² distribution for a given number of degrees of freedom. For three degrees of freedom (radial velocity and 2D proper motion), a 3σ confidence interval corresponds to c = 14.16. For two degrees of freedom, when radial velocity is missing, c equals to 11.83. In other words, stars with a c value equal to or lower than the one corresponding to the selected 3σ confidence interval are considered bona fide members of the cluster. We note that this part of the algorithm does not take into account over-densities in the positions with respect to the cluster centre.

3.2.3. Iterations

The algorithm iterates between the computation of the systemic cluster properties and membership determination until the convergence criterion is achieved. This happens when the absolute difference between the distance of the cluster’s barycentre to the Sun in the previous and the new step drops below 0.01 pc, which typically occurred after the second iteration.

Our definition of cluster membership consists of two criteria. The first criterion is that a star must be within a 3σ confidence interval in the c value. The second criterion is that only stars that are gravitationally bound to the cluster can be included, that is, stars within the tidal radius of the cluster as per definition. We describe the determination of the tidal radius in Sect. 3.2.4.

The tidal radius is always larger than the half-mass radius. To make the estimate of the systemic cluster properties more robust, we recomputed the barycentre and cluster velocity with the objects within the tidal radius. Finally, we took the new cluster properties and used them to reiterate the tidal radius.

3.2.4. Tidal radius

Cluster members are subject to the gravitational fields of the cluster itself and the Galaxy. In a simplified view, stars with sufficient energy are able to escape the cluster at the point where the gravitational pulls between the cluster and the Galaxy balance out (i.e., the first Lagrangian point). We defined the distance of the Lagrangian point from the cluster centre as a tidal radius, assuming that stars are in a circular orbit. We took Eq. (3) from Röser et al. (2011), which describes the relation between the tidal radius r and the mass of the cluster within the tidal radius M_cluster

$$\begin{array}{c}\hfill M_{}\text{cluster}=\frac{4A(A-B)}{G}{r}^3.\end{array}$$(8)

Here, G is a gravitational constant, and A = 14.5 km s⁻¹ kpc⁻¹ and B = −13.0 km s⁻¹ kpc⁻¹ are Oort’s constants, from Piskunov et al. (2006), determined from 581 clusters within 2.5 kpc.

We solved the equation numerically by comparing the relation from Eq. (8) with the cumulative radial mass distribution of the cluster. The intersection of both curves corresponds to the tidal radius. The example of the Pleiades is shown in Fig. 2, where we only use stars with a high membership probability (i.e., 3σ confidence). The figure reveals that the slope of the cumulative radial mass distribution is steeper within the tidal radius, indicating a larger spatial density of the stars. For comparison, we added the distribution of the field stars in the solar neighbourhood (see Sect. 2.4) that reflects a uniform distribution of stars with a density of 0.033 M_⊙ pc⁻³.

The cumulative radial mass distribution of the 3σ candidate members gives a lower limit of the tidal radius because our mass determination method assumes all stars are single (Sect. 2.3). If we account for binaries, the total mass and thus the tidal radius would be larger. To estimate the upper limit of the tidal radius, we took the members of the Pleiades, manually selected all stars seen by eye above the single cluster sequence (Fig. 3), and proclaimed them as equal mass binaries by doubling their mass. This gave a tidal radius of 12.5 pc that is 11% larger than the value from the sample with single stars (Fig. 2). However, the true impact of the binaries is likely smaller because not all binaries are equal mass objects. Additionally, clusters appear elongated along the line of sight due to the parallax uncertainties, and this affects the determination of their tidal radius. Consequently, not all bound stars are actually included as bona fide members by our selection criteria. However, we chose to keep this procedure because we only use members within 1 tidal radius to estimate the global cluster properties (e.g., cluster centre in the physical and velocity space). In our studies of the luminosity functions, we included members up to 3 tidal radii. At such distances to the cluster centre, the stellar density drops, and the distance-to-centre cut does not have such a large dependence on the accuracy of the tidal radius.

Fig. 3. Colour–magnitude diagram for the Pleiades cluster with bona fide and candidate members within 3 tidal radii. The top axis indicates approximate spectral types from the relation in P13. Stars with high ruwe values (i.e., potential multiples) are marked with green circles. Triangles denote white dwarfs from Gentile Fusillo et al. (2021). We found one white dwarf among the bona fide members (EGGR 25).

4. Clusters and their members

Our membership determination algorithm (Sect. 3) produced catalogues of bona fide members within 1 tidal radius and within a 3σ confidence interval in c for all the clusters. We also report candidates within 3 tidal radii for all the clusters in our tables (Sect. 4.1) in order to be as inclusive as possible and not reject evaporated objects and sources belonging to the corona and tidal tails.

We divided the list of clusters is into two subsamples (bound and unbound) because tidal radius computation was not possible for the clusters in the unbound sample. Indeed, they appear too sparse and do not reach the minimum total mass that is necessary for a system to be gravitationally bound (and not be disrupted by the Galactic forces). We manually set their radii to a value that approximately includes most of the obvious members, based on the visual inspection of their overdensity in the proper motion space. We note that their radii are only rough estimates and should be refined for a more detailed analysis.

We further split the bound clusters into three age categories, each of them being represented by one of three benchmark clusters IC 2602, the Pleiades, and the Hyades. These three age classes cover different stages in the nuclear and dynamical evolution of the clusters. At the age of ${46}_{-5}^{+6}\mathrm{Myr}$ (Dobbie et al. 2010), IC 2602 still contains pre-main-sequence stars in the low-mass end. Stars in the Pleiades at the age of 112 ± 5 Myr (Dahm 2015) have already settled into the main sequence, while the low-mass members of the older Hyades (640${}_{-49}^{+67}$ Myr; Lodieu et al. 2019b) might have already evaporated. The age ranges of the groups are < 50 Myr for the IC 2602, 100–200 Myr for the Pleiades, and ≥500 Myr for the Hyades group. The age classification of the clusters into three groups is based on the comparison of their colour–magnitude diagrams with those of the three benchmark clusters. We present the results in Table B.1. In total, there are seven clusters in the IC 2602 age group (Alessi 13, IC 2391, UBC 7, Collinder 135, UBC 9, and Stephenson 1), ten clusters similar to the Pleiades (α Persei, Platais 3, Blanco 1, ASCC 41, ASCC 101, NGC 2516, Alessi 10, Alessi 24, and NGC 1039), and six clusters comparable to the Hyades (Coma Berenices, NGC 2632, Alessi 3, NGC 1901, and Mamajek 4).

In Sect. 4.1, we describe the catalogues and evaluate the membership classification algorithm by comparing our members of the three benchmark clusters with the literature. In Sect. 4.2, we outline the systemic cluster properties, provide their coordinates in the position and velocity spaces, and estimate their velocity dispersions. We comment on completeness and contamination rates in Sect. 4.3.

4.1. Catalogues of cluster members

A schema of the membership catalogue is presented in Table C.1. Full membership lists for all clusters are available as supplementary material in the Vizier database⁴ and on our website, together with the visualisation tools (Sect. 5). We present all members in one catalogue and note that a small percentage of stars (0.8%) appear in the catalogue twice. This duplication is due to their membership in close cluster pairs, resulting in their simultaneous assignment to two clusters (for more details, see Appendix D). We expanded the catalogues with the photometric data from the external databases as outlined in Sect. 2.5. Below, we briefly describe the three benchmark clusters, IC 2602, Pleiades and Hyades, and compare their members from this work with the catalogues from the literature. We note that while other works performed more strict data quality cuts – Cantat-Gaudin et al. (2018a), for example, used stars brighter than GaiaG = 18 mag – we merely applied filters on the parallax quality and thus reached the detection limit of Gaia. Additionally, with the aim of comparing the luminosity functions in mind, we made sure that all our clusters were analysed homogeneously (i.e., all with the same technique). We chose the convergent point method (Perryman et al. 1998) that was designed for astronomy, while most of the literature is based on machine learning techniques that were not originally created for astrophysics.

4.1.1. Pleiades

The Pleiades is one of the most prominent and best-studied clusters. Due to its proximity, large number of members, and high proper motion, it stands out from background stars. Lodieu et al. (2019a) provide a comprehensive list of works related to a historical debate about the distance and age of the Pleiades. The same paper used Gaia DR2 data to retrieve 1248 sources within the tidal radius of 11.6 pc. The authors estimated a distance of 135.2 ± 0.4 pc and based its proper motion on 342 stars with radial velocities. They used white dwarfs to estimate the age of the Pleiades to ${132}_{-27}^{+26}\mathrm{Myr}$, in agreement with the lithium depletion boundary age from Stauffer et al. (1998), Barrado y Navascués et al. (2004); and Dahm (2015).

In this work, we found 1355 bona fide members within a tidal radius of 11.28 ± 0.03 pc. While our tidal radius is slightly smaller than the value from Lodieu et al. (2019a), we report about 9% more members (within the respective tidal radii). This is likely due to more resolved binaries thanks to the longer Gaia baseline. Our distance of 135.0 ± 0.3 pc is in excellent agreement with the aforementioned result.

To make a fair comparison, we took stars within 3σ in the c value and with a distance to the cluster centre of less than 11.6 pc in both samples. The overlap between the samples is 1156 stars; 77 stars are found in the Lodieu et al. (2019a) list but are not members in this work. Additionally, this work increases the sample by 223 stars, which is likely due to now-resolved binary stars.

We compared our catalogue with the members from Cantat-Gaudin & Anders (2020) and found 1045 objects in common. There are 43 stars found only in their catalogue and 334 only in our list. The reason for a larger number of bona fide members in our catalogue is twofold. On one hand, Gaia DR3 resolved more binaries, and on the other hand, Cantat-Gaudin & Anders (2020) only considered stars with a five-parameter astrometric solution in the G = 5–18 mag range.

The narrow sequence in the colour–magnitude diagram in Fig. 3 – independent of the kinematics used in classification – reinforces the credibility of the result and is consistent with the coeval nature of the sample typical for clusters. The contamination rate is higher for cooler stars due to their intrinsically low luminosities, binarity, and lack of radial velocity measurements, but the overall contamination rate remains small. Low mass outliers above the main sequence in the G − R_p colour fall on the main sequence in other colour spaces.

The binary sequence in the cool part of the colour–magnitude diagram is clearly visible in Fig. 3. These stars together with other outliers typically have ruwe > 1.4, indicating either binarity or problems with the astrometric solution.

The significant increase in the number of sources with radial velocity measurements in Gaia DR3 enables a robust statistical study of clusters in the velocity space. Except for a few of the brightest stars, practically all stars earlier than ∼M4 have a radial velocity measurement. This accounts for about one third of the Pleiades members (480 objects). We inferred a systemic cluster velocity of 32.21 ± 0.27 km s⁻¹ in the Galactic reference frame.

4.1.2. Hyades

The nearest star cluster to the Sun is the Hyades, at the distance of 47.03 ± 0.20 pc (Lodieu et al. 2019b). Its age is estimated to be 640${}_{-49}^{+67}\mathrm{Myr}$ based on white dwarfs (Lodieu et al. 2019b) and 650 ± 70 Myr from its lithium depletion boundary (Lodieu et al. 2018; Martín et al. 2018). There are 385 stars within its tidal radius of 9 pc (Lodieu et al. 2019b), but many more are found in its tidal tails. Indeed, Röser et al. (2019) reports 501 members to be within two tidal radii, and 529 stars to be in its asymmetric tidal tails extending up to 170 pc. Similarly, N-body simulation of the Hyades has shown that the total extent of its tails might reach almost 1 kpc and might thus contain more mass (Jerabkova et al. 2021).

We inferred a distance of 46.47 ± 0.14 pc, which is slightly lower than Lodieu et al. (2019b) but in agreement within 2σ. Our tidal radius (8.4  ±  0.02 pc) is also slightly lower than previous estimates (panel a in Fig. 4), but the total number of members (427) is slightly larger. We note that our input candidate list is centred on α = 67.2149 deg and δ = 16.5569 deg and spans 70° in diameter (three times its tidal diameter), while the outskirts of the Hyades extend beyond this limit (Jerabkova et al. 2021). We therefore focussed only on the cluster core and its immediate vicinity. Figure 5 shows the colour–magnitude diagram of the Hyades and reveals seven white dwarfs that are known members of the cluster.

Fig. 4. Cumulative radial mass distributions for all clusters. These plots were used to determine tidal radii. The term M(r) is a theoretical curve that describes the location of the first Lagrange point. The dashed lines represent clusters whose tidal radii we could not determine with our technique, so we set their characteristic radii manually.

Fig. 5. Colour–magnitude diagram for the Hyades cluster with bona fide and candidate members within 3 tidal radii. The top axis indicates approximate spectral types from the relation in P13. Stars with high ruwe values (i.e., potential multiples) are marked with green circles. Triangles denote white dwarfs from Gentile Fusillo et al. (2021). We found seven known white dwarfs among the bona fide members.

To compare our sample with the literature, we prepared a catalogue from Lodieu et al. (2019b) by limiting the data to 3σ in c. We took stars with a distance to the cluster centre within 9 pc in both cases. There are 346 stars in common between the samples. Additionally, 14 stars were found only in the Lodieu et al. (2019b) catalogue and 108 stars only in this work. We attribute the increase of the sample in this work to a larger fraction of resolved binaries. We found a similar but slightly smaller overlap in other works. For example, membership lists from Röser et al. (2019) and (Jerabkova et al. 2021) respectively have 333 and 318 stars in common with our work. On the other hand, we found 121 and 136 stars from this work that are not in these two catalogues, respectively.

4.1.3. IC 2602

The youngest benchmark cluster in our list is IC 2602, also known as the Southern Pleiades. Its lithium depletion age is estimated to be 46${}_{-5}^{+6}\mathrm{Myr}$ (Dobbie et al. 2010) with a distance to 151.6${}_{-1.8}^{+1.9}\mathrm{pc}$ (Nisak et al. 2022). We estimated its tidal radius to 7.2 pc in Fig. 4 and its distance to 150.39 ± 0.72 pc, which are in agreement with the literature. Its location in the Galactic plane (b = −5 deg) makes membership selection difficult. Thus, our membership algorithm results in many interlopers with large photometric uncertainties (σ_G ≥ 0.15) at the faint end, which we removed. Since these interlopers have larger relative parallax uncertainties than the more reliable members, we excluded all stars fainter than GaiaG = 8 mag and with parallax_error/parallax > 0.075. Figure 6 shows the colour-magnitude diagram of the cluster.

Fig. 6. Colour–magnitude diagram for the IC 2602 cluster with bona fide and candidate members within 3 tidal radii. The top axis indicates approximate spectral types from the relation in P13. Stars with high ruwe values (i.e., potential multiples) are marked with green circles.

Our catalogue contains 308 members within 1 tidal radius and 619 stars within 3 tidal radii, while Cantat-Gaudin & Anders (2020) and Nisak et al. (2022) report 311 and 451 members, respectively; both are based on Gaia DR2. The overlap between Nisak et al. (2022) and our sample is 416 stars, while their catalogue contains 44 additional stars that are not in our list. In our comparison with Cantat-Gaudin & Anders (2020), we found 295 objects in common, 324 stars in our list that are not in theirs, and 20 stars in their list that are not in ours.

4.2. Catalogue of cluster characteristics

We present the summary of the cluster properties in Table B.1 for bound clusters and Table B.2 for unbound clusters. The tables include the cluster position in the physical and velocity field, expressed in both observable and Cartesian space. We added tidal radius, the number of bona fide cluster members N, tidal mass, and velocity dispersion in all three directions. We computed the velocity dispersions for each cluster. Table B.1 lists them separately for transversal and radial velocities due to large uncertainties of the latter. We determined a 2D dispersion σ_2D from transversal velocities.

To check if cluster properties evolve with time, we computed the median values for clusters in each age bin, but we observed no clear trends in the median tidal radius and their dispersion, the minimal and maximal tidal radius, nor in 2D dispersion in velocity and the number of bona fide members. This might be due to our age bins being too wide and the fact that while all clusters lose mass with time, the surviving clusters are those with an initially larger mass.

4.3. Completeness and contamination rate

In this section, we address the completeness and contamination rate in our catalogues of cluster members. The completeness of the Gaia DR3 catalogue is affected by the saturation limitations in the bright end and the intrinsic faintness of the low-mass stars. Bright stars missing in Gaia are assumed to be covered in our list by the insertion of the HIPPARCOS catalogue, as described in Sect. 2. The completeness ratio at GaiaG magnitude 20 is 92.2% for a five- and six-parameter solution (Lindegren et al. 2021). Additionally, our input catalogues are affected by the parallax quality cut that correlates with the magnitude and the fact that radial velocities are only available to G = 16 mag.

In Gaia, ‘relatively few sources are found with separations less than about 0.6 arcsec’ (Lindegren et al. 2021). We checked the crowded regions in the core of the most populated cluster in this work, NGC 2516, for a possible source of incompleteness. Stellar pairs with a separation of less than 5 arcsec (less than 1% of the bona fide members) show a rather uniform distribution in the sky. We concluded that stellar density does not affect the incompleteness.

The adopted membership determination algorithm in this work assumes that members have the same motion as the cluster. While the measurement uncertainties in proper motions are small, relatively large errors in parallaxes and radial velocity smear out the internal kinematics, for example, the rotation of the cluster. Our estimate is that this algorithm works well for stars within the tidal radius and those that are still found in the vicinity of the cluster and evaporated in the direction of the cluster motion. However, the technique might not be adequate to identify tidal tails in all directions and large distances from the cluster centre. Stars that escaped in a direction perpendicular to the cluster motion will not be detected. Additionally, the results become less trivial to interpret for clusters that have neighbours in both velocity and physical space, for example, Platais 4 (see Appendix D.3). In this case, the algorithm cannot distinguish between the two clusters and detects stars from both clusters as if there was only one large cluster.

Binary stars orbit their common centre of mass and cause periodic oscillations in their proper motions and radial velocities. Ideally, the velocity of the centre of their mass should be compared with the cluster in our algorithms, but we avoided this step due to the large number of stars in our sample. Instead, we tested the membership detection of binaries in the Pleiades by analysing the distribution of their c values. While Gaia resolves wide binary stars with a separation > 0.4 arcsec, ruwe > 1.4 is a good proxy for close binary objects. Stars with ruwe > 1.4 have a slightly broader distribution of the c values, with the main peak at c = 1 and a smaller peak at around c = 8. This is still well within our 3σ confidence interval (c ≈ 14 for three degrees of freedom). Spectroscopic binaries should be detected by Gaia, but if the common centre of mass is significantly displaced (e.g., in non-equal mass objects), we may miss them. We concluded that the 3σ criterion is good enough to incorporate the majority of binary stars, with a possible exception to the non-equal mass binaries.

The contamination rate remains small for bound clusters. We provide a quantitative estimate from the number of outliers in the colour–magnitude diagram. Contamination in our most distant cluster (NGC 1039 at b = −15.6deg) is up to 2% within 1 tidal radius and 9 − 15% within 3 tidal radii, depending on how conservative we are in estimating the intrinsic spread in the lower main sequence. We note that this is only an approximate estimate, as it is often non-trivial to distinguish between, for example, subdwarfs or overluminous low-mass members and true outliers that might otherwise be reliable kinematic members with problematic photometry. We automated the contamination estimation procedure and applied it to all clusters in the following way: We fit an 11th order polynomial to the cluster sequence within 1 tidal radius and between 0 < G − R_P < 1.5 and computed a root-mean-square of the sequence. We counted stars that are more than 3 RMS away from the sequence as contaminants. We note that this is only an approximate estimate, as the polynomial fit is not an ideal solution. We then computed the contamination rate as the number of outliers over the total number of stars in this colour range. The contamination rates are below 5% for most of the clusters. We found no obvious correlation with the Galactic latitude. Most of our clusters are within 20 deg from the Galactic plane. Those further away have, on average, contamination rate less than 2%. Further, we checked crowdedness in the cluster cores as a possible source of incompleteness and concluded that it is negligible since the pairs with a separation of less than 5 arcsec show a uniform distribution in the sky and are not concentrated in the cluster cores.

5. The website

Visualisation of the multidimensional parameter space is challenging, especially if the catalogues are large and the phenomena of interest are observed on various scales. To present our results in an adequate way for all the clusters, we designed a website⁵. It includes interactive and 3D visualisation tools and allows the user to select the quantities on the axes of the plot for the cluster of interest. The list of available columns includes Gaia data and the external photometric measurements described in Sect. 2.5. We added a database of PARSEC isochrones to overplot over the colour–magnitude diagrams.

One of the highlights is the interactive link between our plots and the Aladin sky atlas (Bonnarel et al. 2000; Boch & Fernique 2014)⁶. In particular, a click on a given star in our plot centres the Aladin plugin on the chosen object. This feature can be useful in many cases, for example, to inspect the astrometric aspect of close binaries and of faint or bright objects. We linked every cluster page with the catalogue of its members to the Vizier database and enabled a direct download of the data for further analysis.

The website was built to enable the professional community to have a quick look at well-known clusters. At the same time, we hope to bring the latest scientific results to the public and inspire new citizen science and exciting student projects.

6. Luminosity functions

The luminosity function represents the number of stars per magnitude bins within a given cluster. It is subject to a temporal variation due to nuclear and dynamical processes of the stellar population. For example, the low-mass stars in the young cluster IC 2602 are still found in the pre-main-sequence stage and thus appear more luminous than their older counterparts. In the later phases of a cluster’s life, its orbital motion in the Galaxy creates interactions and results in the gradual evaporation of weakly bound and low-mass members. For example, the Hyades lost a significant fraction of M dwarfs from its core. In fact, a significant number of its members very likely drifted away to its tidal tails (Röser et al. 2019; Jerabkova et al. 2021). In this work, we also studied the timescale of dynamical evolution and estimated the percentage of M dwarfs lost from clusters at different stages of their life.

To address this question of mass loss, we first re-evaluated the ages of the clusters in the bound subsample from Table B.1 whose ages span ∼50 to ∼650 Myr in order to cover the most significant period in the evolution of the clusters (Sect. 4). We discuss their constructed spectral type distributions in Sect. 6.1 and their luminosity functions in Sect. 6.2. For comparison, we added the field population from the immediate Solar neighbourhood described in Sect. 2.4 but manually removed white dwarfs because they correspond to objects with larger initial masses.

6.1. Spectral type distributions

The low-resolution version of the luminosity function is a spectral type distribution. We divided stars into spectral types based on their G − R_P colour from the Gaia DR3 catalogue. The spectral type classification is based on P13. To test the reliability of this conversion, we took the sample of field stars within 8 pc from Kirkpatrick et al. (2012) and compared their spectral types with the classification from the colour conversion. Out of 141 dwarf stars that were cross-matched with the Gaia DR3 catalogue, we found eight cases with discrepant spectral types. This is a relatively small error (6%); thus, we concluded that the classification using the G − R_P colour relation is sufficient for our purpose.

Figure 11 from Kirkpatrick et al. (2012) presents a sample of all objects within 8 pc in three different ways: the total mass per spectral type, the number of stars per spectral type, and the spatial mass density per spectral type. To reproduce and compare our results with their figure, we took our mass estimates from Sect. 2.3 and counted the number of stars per spectral type. We note that binary stars in our work are treated as single objects in the mass determination and thus represent a lower mass limit. We show the results for our field sample within 25 pc (described in Sect. 2.4) and the benchmark clusters in Fig. 7. In order to compare the clusters with the field population, we overplot the field sample scaled to the total mass and number of late-G stars. To trace the change in the distribution in greater detail, we split the spectral types into early- and late-type classes (i.e., from 0–4 and 5–9 subtype bins).

Fig. 7. Spectral type and total mass distribution of the field dwarfs within 25 pc and three benchmark clusters, IC 2602, the Pleiades, and the Hyades. Spectral types G, K, and M are split into early and late subtypes. The field population was used to compute the expected number of stars for each cluster with respect to the number of late-G dwarfs (grey bars).

As expected, the number of objects and the total mass grow towards the late spectral types, with a peak around M3–M4 dwarfs. This corresponds to the peak of the present day mass function in the solar neighbourhood and open clusters (Kroupa 1998, 2001; Chabrier 2001; Bastian et al. 2010). We quantified this observation by computing the fraction of M dwarfs in the entire sample. We did this for masses and the number of stars and present the results in Table 3.

The mass and number ratio of M dwarfs in IC 2602 and the Pleiades match the field population (∼40% of the total cluster mass and 80% in terms of a star count). In contrast, only 65% of stars in the Hyades are M dwarfs, representing 30% of the total cluster mass. Assuming the same initial mass function as for IC 2602 and the Pleiades, it is clear that the Hyades lost its M dwarfs because of evaporation during its multiple rotations around the Galaxy and the dynamical processes that lead to the formation of the tidal tails that stretch up to nearly 1 kpc (Jerabkova et al. 2021). The number ratio between the early- and late-M dwarfs in the Hyades is twice as large as in IC 2602 and the Pleiades. This means that the late-type M dwarfs are missing in the Hyades. Their absence cannot be attributed to Gaia’s incompleteness, since the Hyades is the closest of all clusters (Lodieu et al. 2019b).

The low-mass stars in the youngest clusters are in their pre-main-sequence phase. Clusters with the Pleiades age can already contain white dwarfs, while clusters in the oldest group with the Hyades age stars can produce (sub)giants. The number of A-type stars in a cluster decreases with age, as expected from stellar evolution models. On the low-mass end, we observed that the numbers of M dwarf members is similar between the field, the Pleiades, and IC 2602 within the Poisson error bars, while the Hyades shows a clear deficit of low-mass and very low mass stars. We quantified this fact by calculating the ratio of early-M dwarfs to G-type stars in the three benchmark clusters (Table 3). This substantial loss of M dwarfs is the result of interaction during the motion of the cluster in our Galaxy, which yielded the dynamical evaporation of the lower mass members. This appears to be natural in the evolution of a young cluster into an older one, as simulated by Kroupa et al. (2001).

6.2. Luminosity functions

A luminosity function is defined as the distribution of the luminosity in a sample (i.e., a distribution of stellar magnitudes). Its profile is subject to the intrinsic properties of a population, such as its initial mass function, age, and star formation history, as well as the dynamical processes that cause gradual evaporation of the low-mass stars in the population.

We determined absolute GaiaG magnitudes using the Bailer-Jones et al. (2021) distances and assumed a negligible extinction since the clusters are located close to the Sun. To check our assumption, we compiled E(B − V) values from the literature for the bound clusters in Appendix A. Clusters beyond 400 pc (e.g., Alessi 10, NGC 2516, and NGC 1039) are the most reddened clusters in our sample, with the E(B − V) reaching 0.2 mag, and we found no correlation with their age. Three clusters are missing the literature value, and two of them (UBC 7: 40 − 50 Myr Kovaleva et al. 2020, and UBC 9, classified as young in this work) are estimated to be very young. We thus provide our own estimate of the upper limit of reddening and extinction based on dustmaps (Green 2018). The large model reddening values exceed those measured in other works. The upper limit for E(B − V) is 0.21 for UBC 7, but it is based on other clusters. We deduced that the value is likely much smaller. Assuming R_V = 3.1, the typical extinction A_V is much less than one magnitude, and we thus concluded that it does not significantly affect our luminosity functions. A possible exception is Alessi 10, which has an extinction of A_V = 0.6 that is approximately half of the magnitude bin in our luminosity function, as described below. We computed the luminosity functions in the range from 0 to 15 mag with a fixed step of one magnitude and normalised them by the area under the curve after counting the number of objects per magnitude bin.

Figure 8 shows the luminosity functions of the clusters in each age group together with the benchmark clusters for comparison. In the figure, we plot all members within 3 tidal radii. Within each age group, the luminosity functions seem to be similar to each other, as they practically overlap with the one of the benchmark cluster. To quantify this comparison, we used a Kolmogorov–Smirnov test for two samples from SciPy. The test was performed pair-wise within the common magnitude completeness limits of each cluster pair. We present the p-values for all pairs in Fig. 9, where we sort clusters according to their age group. The p-values indicate the probability that the two samples, in our case the luminosity functions of the two clusters in question, come from the same distribution. We note that off the diagonal (the diagonal compares clusters to themselves), the values are generally lower than 0.5, with a few exceptions. Because it is not trivial to interpret the absolute values and the test might sometimes be affected by the low-number statistics in the bright end of the luminosity functions, we focussed on the groups as a whole. We most often found pairs with a relatively high p-value in the youngest age group. This means that clusters in this age group are more similar to each other than in the intermediate and old age groups. At the same time, it is clear that there is no significant similarity between the young and older clusters from the Pleiades and Hyades age groups, with an exception being cluster Alessi 10. On the other hand, there is no significant difference between the Pleiades and Hyades age groups. Clusters NGC 1901, NGC 1039, and Platais 3 are affected by incompleteness and low-number statistics. The average p-values for both age groups, the Pleiades and the Hyades, appear to be lower, and one is less likely to find a pair from the same age group with a high degree of similarity.

Fig. 8. Luminosity functions for clusters in each age group. These functions have been normalised such that the area under the curve equals one. We include all members within 3 tidal radii. The solid black lines show the benchmark cluster, and shaded areas represent the Poissonian uncertainty. We sort the clusters by their distance to more easily trace the completeness limit (denoted by vertical lines). The y-axes show the p-value from the Kolmogorov–Smirnov comparison, with the benchmark cluster within the common completeness limits.

Fig. 9. Pair-wise Kolmogorov–Smirnov test p-values for the luminosity functions of the bound clusters. We sort the clusters by their age group and present them with the following colour code: black for the IC 2602, red for the Pleiades, and blue for the Hyades age group.

Moreover, there are seven clusters that bear similarities to clusters from another age group. These clusters are Alessi 24, ASCC 101, and Platais 3 from the intermediate age group and Melotte 111, Mamajek 4, Melotte 25, and Alessi 3 from the old age group.

To summarise, clusters in the youngest age group are the most similar to each other, and the degree of similarity later decreases. These observations led us to draw a parallel with the mass functions in star-forming regions, open clusters, and globular clusters, as they are universal within current observational uncertainties (see review by Bastian et al. 2010). Based on our data, the luminosity function appears to be universal at least within 500 pc from the Sun. We reiterate that our cluster selection function only depends on the distance (500 pc), upper age limit of 1 Gyr, and the distance from the Galactic plane (|b|> 10 deg), and the membership selection bias was minimised by the fact that we used the same method on all of the clusters. Our results are consistent with the earlier study of Piskunov et al. (2008), who observed a strong similarity between the luminosity functions in nearby clusters and other galaxies, and agree with the findings in young regions, as Galli et al. (2020) state that their Gaia members in the Taurus, Upper Scorpius, and Lupus regions show a similar shape for the distribution of spectral types at the faint end of the initial mass function.

If the luminosity function is universal at least up to the Hyades age, then all our clusters must have been born with similar luminosity functions. However, the luminosity functions of our intermediate and old age clusters differ from the young clusters, and this might be due to many reasons. One reason is stellar evolution that, for example, causes massive stars to increase their luminosity when they turn into giants. On the other hand, since single stars of the same mass and metallicity evolve at the same rate, clusters of similar advanced age (e.g., the Hyades age) should still have similar luminosity functions. Our clusters are single (with the exception being the pair Collinder 135 – UBC 7, see Appendix D.1), and the reason for the increased dissimilarity with age might be in their internal dynamical evolution, which depends on their density. However, further analysis with models and simulations would be needed to provide more evidence to support this claim.

7. Conclusions

We selected 49 open clusters younger than 1 Gyr and closer than 500 pc and reassessed their membership using Gaia DR3 and HIPPARCOS 6D astrometric information (positions, parallaxes, proper motions, and radial velocities) over the full range of validity of both missions. We discarded some known clusters in the Galactic plane. We provide the physical parameters of each cluster (3D position and velocity, distance, proper motion, age, tidal radius and mass, number of members) and an updated list of bona fide members, up to 1 tidal radius, and candidates, up to 3 tidal radii. We will make the new list of members public in a dedicated webpage⁷ as well as on Vizier (the Centre de Données de Strasbourg) for use by professional and amateur astronomers and educators. We studied the evolution of the luminosity function as a function of age and found that luminosity functions of the youngest clusters differ from the older populations, are more likely similar to each other, and show a greater degree of similarity than older clusters. We explained this observation with the universal luminosity function within the volume of our sample (500 pc). Clusters with ages similar to the Pleiades or Hyades show a lesser degree of similarity than younger clusters. If the luminosity function is universal, at least up to the age of the Hyades, one of the reasons for increased diversity of these single clusters could be internal dynamical evolution, but future work is needed to support this claim.

Acknowledgments

We thank the anonymous reviewer and acknowledge their valuable contribution that helped to significantly improve this paper. MŽ and NL acknowledge support from the Consejería de Economía, Conocimiento y Empleo del Gobierno de Canarias and the European Regional Development Fund (ERDF) under grant with reference PROID2020010052. MŽ, NL, JO, VB and ELM acknowledge support from the Agencia Estatal de Investigación del Ministerio de Ciencia e Innovación (AEI-MCINN) under grant PID2019-109522GB-C53. APG acknowledges support from the grant PID2020-120052GB-I00 financed by MCIN/AEI/10.13039/501100011033. Co-funded by the European Union (ERC, SUBSTELLAR, project number 101054354). Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Research Council. Neither the European Union nor the granting authority can be held responsible for them. This research has made use of the Simbad and Vizier databases, and the Aladin sky atlas operated at the centre de Données Astronomiques de Strasbourg (CDS), and of NASA’s Astrophysics Data System Bibliographic Services (ADS). This research has made use of “Aladin sky atlas” developed at CDS, Strasbourg Observatory, 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 Two Micron All Sky Survey, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation. Based on observations obtained as part of the VISTA Hemisphere Survey, ESO Program, 179.A-2010 (PI: McMahon). The UKIDSS project is defined in Lawrence et al. (2007). UKIDSS uses the UKIRT Wide Field Camera (WFCAM; Casali et al. 2007). The photometric system is described in Hewett et al. (2006), and the calibration is described in Hodgkin et al. (2009). The pipeline processing and science archive are described in Irwin et al. (in prep.) and Hambly et al. (2008). 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, and NEOWISE, which is a project of the Jet Propulsion Laboratory/California Institute of Technology. WISE and NEOWISE are funded by the National Aeronautics and Space Administration. The Pan-STARRS1 Surveys (PS1) and the PS1 public science archive have been made possible through contributions by the Institute for Astronomy, the University of Hawaii, the Pan-STARRS Project Office, the Max-Planck Society and its participating institutes, the Max Planck Institute for Astronomy, Heidelberg and the Max Planck Institute for Extraterrestrial Physics, Garching, The Johns Hopkins University, Durham University, the University of Edinburgh, the Queen’s University Belfast, the Harvard-Smithsonian Center for Astrophysics, the Las Cumbres Observatory Global Telescope Network Incorporated, the National Central University of Taiwan, the Space Telescope Science Institute, the National Aeronautics and Space Administration under Grant No. NNX08AR22G issued through the Planetary Science Division of the NASA Science Mission Directorate, the National Science Foundation Grant No. AST-1238877, the University of Maryland, Eotvos Lorand University (ELTE), the Los Alamos National Laboratory, and the Gordon and Betty Moore Foundation. Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the U.S. Department of Energy Office of Science. The SDSS-III website is http://www.sdss3.org/. SDSS-III is managed by the Astrophysical Research Consortium for the Participating Institutions of the SDSS-III Collaboration including the University of Arizona, the Brazilian Participation Group, Brookhaven National Laboratory, Carnegie Mellon University, University of Florida, the French Participation Group, the German Participation Group, Harvard University, the Instituto de Astrofisica de Canarias, the Michigan State/Notre Dame/JINA Participation Group, Johns Hopkins University, Lawrence Berkeley National Laboratory, Max Planck Institute for Astrophysics, Max Planck Institute for Extraterrestrial Physics, New Mexico State University, New York University, Ohio State University, Pennsylvania State University, University of Portsmouth, Princeton University, the Spanish Participation Group, University of Tokyo, University of Utah, Vanderbilt University, University of Virginia, University of Washington, and Yale University. Software: astropy (Price-Whelan et al. 2018), NumPy (Harris et al. 2020), IPython (Pérez & Granger 2007), topcat (Taylor 2005) and Matplotlib (Hunter 2007).

References

Appendix A: Extinction

In Table A.1, we present model extinction and reddening values for each bound cluster and compare them with the values in the literature. We estimated the upper E(B-V) limit for each star using the 2D maps of dust reddening with the GNILC method from the Planck Collaboration Int. XLVIII 2016 that is integrated in dustmaps (Green 2018). To avoid the model uncertainties for very hot and cool stars, we applied a photometric filter and only kept stars with Gaia’s 0 < B_P − R_P < 1. We assumed R_V = 3.1.

Appendix B: Cluster properties

Table B.1 presents the characteristics of bound clusters, while unbound objects are listed in Table B.2. Cartesian coordinates and velocity dispersions are given in Table B.3.

Appendix C: Catalogue of cluster members

We describe a schema of the table with cluster members in C.1. We note that some source_id entries are repeated. This happens in cases where stars are simultaneously candidate members of two different clusters that overlap in their outskirts, for example, in the case of the known pair Collinder 135 and UBC 7 (see Appendices D.1 and D.2).

Appendix D: Individual clusters

D.1. A pair: Collinder 135 and UBC 7

Clusters Collinder 135 and UBC 7 are a physical pair of clusters (Kovaleva et al. 2020). They have very similar proper motions, but their distances to the Sun are 294 pc and 273 pc, respectively. Due to their proximity – they are found only 23 pc apart – some of their candidate members at 3 tidal radii overlap (52 stars).

D.2. A triplet: NGC 2516, UPK 540, and UPK 533

The Southern Beehive (NGC 2516) is known to be enveloped by a large 500 pc halo (Bouma et al. 2021). While we only list candidates up to its 3 tidal radii, we find a partial overlap between Southern Beehive, UPK 540, and UPK 533. While they lie 54 pc apart from each other, UPK 540 and UPK 533 are found 50 and 90 pc away from NGC 2516, respectively.

Out of 54 bona fide members in UPK 540, 25 are found in the halo of NGC 2516 in (Bouma et al. 2021), while there are no stars from UPK 533 in the halo. At the same time, UPK 533 is located very close to another cluster, NGC 2547. They are close but clearly separable in the proper motion space, but we found a significant overdensity of stars from NGC 2547 within the 3 tidal radii of UPK 533. Based on the colour-magnitude diagram of UPK 533, the contamination of its 44 bona fide members is likely small.

Due to their small separation, our list of members includes repeated stars that belong to two clusters. There are 158 repeated stars from the pair UPK 540–NGC 2516 and five from the pair UPK 540–UPK 533.

D.3. Platais 4

This cluster is a weak overdensity next to a massive cluster in the proper motion space. Both clusters are located in a similar region of the sky. We believe that these factors led to a false characterisation of this cluster. The results for this cluster are not convincing, especially because we do not have a 3D velocity for it.

D.4. UBC 31

This cluster is located in a very complex region of the proper motion space with no clear distinction from its close neighbours, NGC 1579 and IC 348. The bona fide members of this cluster suffer from a higher contamination rate that increases significantly for candidate members beyond 1 tidal radius.

D.5. UBC 159

One of the most distant clusters in our list, UBC 159 is located at around 500 pc from the Sun. This is a young cluster with its members in the pre-main-sequence stage. Due to its large distance, the contamination rate seems to be higher than for other clusters.

All Tables

All Figures

	Fig. 1. Comparison of isochrones for a range of ages from 50 Myr to 1 Gyr for the BHAC15 (dashed lines) and PARSEC models (solid lines). The polynomial fit used as a colour–mass relation is shown in magenta. While we only plot the relation for G − RP > −0.2, the fit extends to bluer colours in order to cover the most massive stars in our input catalogue.
In the text

Fig. 2. Tidal radius determined by a cross-section between the cumulative radial mass distribution of the cluster and M(r) from Eq. (8) (blue line). The tidal radius for single stars (black) represents a minimal value, while a synthetic sample where all binaries have an equal mass (red) indicates the upper limit. The difference between these tidal radii is about 10%. The field distribution was added for comparison (purple and green for the actual and simulated field, respectively; see Sect. 2.4).

In the text

Fig. 3. Colour–magnitude diagram for the Pleiades cluster with bona fide and candidate members within 3 tidal radii. The top axis indicates approximate spectral types from the relation in P13. Stars with high ruwe values (i.e., potential multiples) are marked with green circles. Triangles denote white dwarfs from Gentile Fusillo et al. (2021). We found one white dwarf among the bona fide members (EGGR 25).

In the text

	Fig. 4. Cumulative radial mass distributions for all clusters. These plots were used to determine tidal radii. The term M(r) is a theoretical curve that describes the location of the first Lagrange point. The dashed lines represent clusters whose tidal radii we could not determine with our technique, so we set their characteristic radii manually.
In the text

Fig. 5. Colour–magnitude diagram for the Hyades cluster with bona fide and candidate members within 3 tidal radii. The top axis indicates approximate spectral types from the relation in P13. Stars with high ruwe values (i.e., potential multiples) are marked with green circles. Triangles denote white dwarfs from Gentile Fusillo et al. (2021). We found seven known white dwarfs among the bona fide members.

In the text

	Fig. 6. Colour–magnitude diagram for the IC 2602 cluster with bona fide and candidate members within 3 tidal radii. The top axis indicates approximate spectral types from the relation in P13. Stars with high ruwe values (i.e., potential multiples) are marked with green circles.
In the text

	Fig. 7. Spectral type and total mass distribution of the field dwarfs within 25 pc and three benchmark clusters, IC 2602, the Pleiades, and the Hyades. Spectral types G, K, and M are split into early and late subtypes. The field population was used to compute the expected number of stars for each cluster with respect to the number of late-G dwarfs (grey bars).
In the text

Fig. 8. Luminosity functions for clusters in each age group. These functions have been normalised such that the area under the curve equals one. We include all members within 3 tidal radii. The solid black lines show the benchmark cluster, and shaded areas represent the Poissonian uncertainty. We sort the clusters by their distance to more easily trace the completeness limit (denoted by vertical lines). The y-axes show the p-value from the Kolmogorov–Smirnov comparison, with the benchmark cluster within the common completeness limits.

In the text

	Fig. 9. Pair-wise Kolmogorov–Smirnov test p-values for the luminosity functions of the bound clusters. We sort the clusters by their age group and present them with the following colour code: black for the IC 2602, red for the Pleiades, and blue for the Hyades age group.
In the text