© ESO 2020

1. Introduction

The study of the physical properties and evolution of massive stars (M > 8−9 M_⊙) is crucial for many aspects of our understanding of the Universe. They play an important role in the chemodynamical evolution of the galaxies (Matteucci 2012) and were key players in the epoch of reionization of the Universe (Bromm et al. 1999; Abel et al. 2000). They are the precursors of hyperenergetic supernovae, long-duration γ-ray burst (see Langer 2012, and references therein), and the recently detected gravitational wave events (e.g., Abbott et al. 2016, 2017; Ackley et al. 2020). Their high luminosities make them observable individually at large distances, and they are thus optimal tools for access to invaluable information about abundances and distances in galaxies at up to several megaparsec (e.g., Urbaneja et al. 2003; Castro et al. 2008; Kudritzki et al. 2013). Moreover, through their feedback into the interstellar medium in the form of ultraviolet radiation and stellar winds, massive stars critically affect the star formation process by both triggering the formation of new generations of stars and stopping mass accretion in the surrounding forming stars.

Most massive stars are found within or are linked to young open clusters and the so-called OB associations (Lada & Lada 2003; Portegies Zwart et al. 2010). These stellar groupings are therefore perfect laboratories to study them.

Humphreys (1978) compiled a catalog of all known supergiants (Sgs) and O stars in associations and clusters of the Milky Way, including over 1000 objects of this type. Among the list of associations quoted in that paper, Per OB1, which also includes the famous h and χ Persei double cluster, clearly stands out as one of the richest. In particular, it is one of the few Galactic OB associations in which, given its age (∼13−14 Myr, Slesnick et al. 2002; Li et al. 2019), a massive star population covering a wide range of evolutionary stages can be found (e.g., it harbors 23 red Sgs and several dozen blue Sgs). In addition, it is also relatively close to us (d ∼ 2.2−2.4 kpc, Babusiaux 2018a; Davies & Beasor 2019) and is characterized by a moderate extinction (E(B − V) ∼ 0.6, Slesnick et al. 2002). This unique combination of characteristics makes Per OB1 a very interesting testbed for the study of a large interrelated population of evolved massive stars from an evolutionary point of view.

Per OB1 has attracted the attention of the astrophysical community for many years and has been the subject of studies from many different fronts. We highlight the investigation of how the association could have been formed (Lee & Lim 2008); the membership of stars to the association (Humphreys 1970, 1978; Garmany & Stencel 1992; Lee & Lim 2008; Mel’Nik & Dambis 2009) and, in particular, to h and χ Persei (Uribe et al. 2002; Currie et al. 2010); the characterization of the kinematics of the region (Mel’nik & Dambis 2017; Zhong et al. 2019); the identification of blue Sg binaries (Abt & Levy 1973); or the spectroscopic characterization of different samples of blue stars in the region (including the determination of rotational velocities, stellar parameters and surface abundances; Slettebak 1968; Lennon et al. 1988; Kendall et al. 1995, 1996; Strom et al. 2005; Li et al. 2019), also reaching the red Sg domain (Gazak et al. 2014).

Despite all the information compiled about the Per OB1 association, and particularly, h and χ Persei, we still lack a complete homogeneous empirical characterization (that also takes environmental and kinematical information into account) of the physical and evolutionary properties of its massive star population. This is the main objective of this series of papers, which is based on a set of high-quality observations including high-resolution, multi-epoch spectroscopy (mostly gathered in the framework of the IACOB project, see Simón-Díaz et al. 2015 and references therein), and astrometric information delivered by the Gaia mission (Brown 2018b; Lindegren et al. 2018). The compiled empirical information resulting from the analysis of this observational dataset will allow us to proceed in our understanding of massive star evolution, and also investigate some long-standing and new open questions in this important field of stellar astrophysics. These questions include the evolutionary status of the blue supergiants, or the effect that binarity and rotation have on the evolution of massive stars.

In this first paper, we carry out a membership analysis of a sample of 88 blue and red Sgs located within 4.5 deg from the center of the Per OB1 association, and we also investigate some of its kinematical properties. In Sect. 2 we present the sample of stars and the main characteristics of the compiled observations. In Sect. 3 we describe the strategy we followed to derive reliable radial velocities (RVs). Section 4 presents the results extracted from the analysis of the observations, mainly referring to parallaxes, proper motions, RV measurements, and the identification of spectroscopic signatures of binarity and other types of spectroscopic variability phenomena. In Sect. 5 we use all these results to establish and apply our membership criteria to all stars in the sample, and we also identify outliers for each of the considered quantities, in particular, binary and runaway stars. We also analyze some global features of Per OB1, and discuss some individual cases of interest. The main conclusions of this work and some future prospects are provided in Sect. 6.

2. Sample definition and observations

In this section, we describe the process we have followed to build the sample under study, and to compile the associated observations. The latter mainly refers to high-quality spectroscopy obtained with the FIbre-fed Echelle Spectrograph (FIES) (Telting et al. 2014) and the High Efficiency and Resolution Mercator Echelle Spectrograph (HERMES) (Raskin et al. 2011) high-resolution spectrographs attached to the 2.56 m Nordic Optical Telescope (NOT) and the 1.2 m Mercator telescope, respectively, and astrometric and photometric data delivered by the Gaia mission in the second data release (DR2, Brown 2018b; Lindegren et al. 2018; Evans et al. 2018).

2.1. Sample definition

The final sample of targets considered for this work comprises 88 blue and red massive stars located within 4.5 deg from the center of the Per OB1 association (as defined in Mel’nik & Dambis 2017). To restrict the sample to the most massive stars, the luminosity classes (LCs) were limited to bright giants (Gs) and Sgs (LC II and I, respectively) in the case of the O- and B-type stars, and to Sgs when we refer to A- and later-type stars. In addition, the sample includes a few O and early-B Gs for which we already had available observations in the IACOB spectroscopic database (see Sect. 2.2).

Table A.1 summarizes the list of targets, separated and ordered by spectral type (SpT). We note that the quoted spectral classifications were carefully revised using the spectra with the best signal-to-noise ratio (S/N) of our own spectroscopic observations (see Sect. 2.2) following the criteria explained in Negueruela (in prep.) and Dorda et al. (2018), for the case of the blue and red Sg samples, respectively. In addition, Fig. 1 shows their location on the sky. We also indicate as a large green circle the search area of 4.5 deg around the center of Per OB1, marked as a green cross. Most stars, including those from the h and χ Persei double cluster (indicated as two small green circles), are concentrated along the diagonal of the image. In addition, our sample includes four stars lying within one degree from the center of IC 1805 (the Heart nebula, located in the top left corner of the figure). The top panel in Fig. 2 depicts the histogram of SpT of the sample, which shows that the majority of stars are B Sgs.

Fig. 1. Sky map with the complete sample of stars. Purple, blue, cyan, golden, and red symbols represent the O-, B-, A-, F-, and K- and M-type stars, respectively. This color code is the same in all the plots unless otherwise specified. The central green cross denotes the center of the Per OB1 association taken from Mel’nik & Dambis (2017). The large green circle indicates a 4.5-degree circle around the center. The small green circles show the positions of h and χ Persei. The background image, used for reference, was taken from DSS-red.

To assemble this sample of stars, we considered several bibliographic sources, including the works by Humphreys (1978), Garmany & Stencel (1992), Currie et al. (2010), and Gazak et al. (2014). In a first step, we used the Topcat¹ Virtual Observatory tool to cross-match all the stars that are quoted in these four papers and fulfilled the criteria indicated above and the list of targets with spectra available in the IACOB spectroscopic database. In a second step, we tried to obtain new spectra of as many of the missing stars as possible using the NOT or Mercator telescopes (see Sect. 2.2).

From the original lists of luminous stars in Galactic OB associations quoted in Humphreys (1978) and Garmany & Stencel (1992), we found 207 targets that are located within 4.5 deg of the center of Per OB1. Only 109 of these fulfill our luminosity class criteria; the rest are either dwarfs, (sub)giants, or do not have a defined luminosity class. Our sample includes 82 of these, but we miss spectra for another 12 (7 B and 5 M Sgs).

We also used the list of targets quoted in the extensive study of the stellar population of h and χ Persei by Currie et al. (2010) to find suitable candidates. From the complete list of several ten thousand stars, only 23 were found to have luminosity classes I or II. We currently have spectra for 17 of them. Of the remaining 6 (all of them B Sgs), one was identified previously when we cross-matched our observations with the list of targets in Humphreys (1978) and Garmany & Stencel (1992). This means that we lack spectra for another 5 blue Sgs at the time of writing.

Last, our sample includes all the red Sgs of those listed in Gazak et al. (2014). In summary, the sample of stars we discuss here comprises all the blue and red Sgs (except for 12 B and five M Sgs, listed at the end of Table A.1 for future reference) that are quoted in the abovementioned papers and are located within 4.5 deg around the center of Per OB1. Further notes on the actual completeness of our sample can be found in Sect. 4.3.

2.2. Spectroscopic observations

The spectroscopic observations of the stars in the sample come from different observing runs performed between November 2010 and December 2019 using either the FIES (NOT) or the HERMES (Mercator) instruments.

The first observations, comprising an initial sample of B, A, and M Sgs in Per OB1 selected from Humphreys (1978), were obtained in 2010 during an observing run of three nights with Mercator (PI. M.A. Urbaneja). The O stars in the sample were targeted by the IACOB project (P.I. S. Simón-Díaz) as part of a more general objective of observing all O stars in the Northern Hemisphere up to V_mag = 9. These observations, obtained with both HERMES and FIES, include a minimum of three epochs per target (see more details in Holgado et al. 2018, 2020; Holgado 2019). We also benefit from the multi-epoch observations available for a subsample of O and B Sgs as gathered by the IACOB project as part of a subproject aimed at investigating line-profile variability phenomena in the OBA Sg domain and its relation with pulsational-type phenomena (see, e.g., Simón-Díaz et al. 2010, 2017, 2018; Aerts et al. 2017, 2018). The time span of these observations covers several years. We also count on multi-epoch observations of red Sgs obtained during several of our observing runs with HERMES. Last, all these observations have more recently been complemented by FIES spectroscopy obtained as part of the time granted to A. de Burgos in 2018 by the Spanish time-allocation committee, and through internal service observations performed in 2019 and 2020 by A. de Burgos. In addition, we were able to obtain a new epoch for a large fraction of stars in the sample during an observing run with Mercator in December 2019.

FIES is a cross-dispersed high-resolution échelle spectrograph mounted at the 2.56 m Nordic Optical Telescope (NOT), located at the Observatorio del Roque de los Muchachos on La Palma, Canary Islands, Spain. The observations made with FIES were taken with different fibers/resolutions from R ∼ 25 000 to R ∼ 67 000, and with a wavelength coverage of 370−830 nm.

HERMES is a fibre-fed prism cross-dispersed échelle spectrograph mounted at the 1.2 m Mercator Telescope, also located at the Observatorio del Roque de los Muchachos. It provides a spectral resolution of R ∼ 85 000 and wavelength coverage of 377−900 nm, similar to FIES.

The FIES and HERMES spectrographs provide good mechanical and thermal stability that allows for a good precision in RV measurements. For FIES, the RV y accuracy² for the high-resolution fiber has been proved to be 5−10 m s⁻¹, regardless of the atmospheric conditions. For the medium-resolution fiber under poor conditions, the precision reaches 150 m s⁻¹. In the case of HERMES, the precision obtained for the low- and high-resolution fibers is 2.5 and 2 m s⁻¹, respectively (Raskin et al. 2011). In both cases this precision is well above the precision required for this work, as we expect variations of several km s⁻¹ for the blue Gs/Sgs, and a few km s⁻¹ for the red Sgs.

All the spectra were reduced using the FIESTool (Stempels & Telting 2017) and HermesDRS³ dedicated pipelines. Both pipelines provide merged wavelength-calibrated spectra. In addition, we used our own programs, implemented in IDL, to normalize the spectra and compute the heliocentric velocity to be applied to each spectrum before the associated RV was measured (see Sect. 3).

As indicated above, we have multi-epoch spectroscopy for a large fraction of the stars in our sample. The bottom panel in Fig. 2 summarizes this characteristic of our observations, showing the histogram of the collected number of spectra per star. In addition, Table A.2 quotes all those stars for which we have five or more spectra. This table includes the time span covered by the spectra, together with the total number of spectra for each of these stars, separated by SpT. It is important to remark that the cadence of the spectra taken for each star is very inhomogeneous, as they were gathered during different observing runs, as described at the beginning of this section.

Fig. 2. Histograms by SpT separated with colors (top) and number of spectra separated by SpT and stacked (bottom).

2.3. Photometric and astrometric data

For all the stars in the sample, Table A.1 quotes the Gaia G_mag and BP_mag − RP_mag, parallaxes (ϖ) and proper motions (μ_α,μ_δ), as well as associated errors, retrieved from Gaia DR2. Sources in the Gaia catalog were identified using Topcat, defining a radius threshold of 2 arcsec.

We adopted a parallax zero-point offset of −0.03 mas (see Lindegren et al. 2018), which is already applied to all values quoted in Table A.1 and used to generate the various figures in the paper. We note, however, that some other authors push this value up to −0.08 mas (see Stassun & Torres 2018; Davies & Beasor 2019).

The Gaia DR2 renormalized unit weight error (RUWE) is also included in the last column of Table A.1. The value of this quantity is used to estimate the goodness of the Gaia astrometric solution for each individual target. Following recommendations by the Gaia team for the known issues⁴, we decided to adopt a RUWE = 1.4 to distinguish between good and bad solutions.

Seven stars (or 8% of the sample) have an associated RUWE higher than this value. Their parallaxes and proper motions are indicated in parentheses in Table A.1. Hereafter, we call them stars with “unreliable astrometry” or “unreliable astrometric solution”. For all the stars with a RUWE < 1.4, the top panel of Fig. 3 shows the G_mag against the Gaia error in parallax, and the bottom panel shows the Gaia error in total proper motion against the Gaia error in parallax.

Fig. 3. Top: Gmag against the Gaia error in parallax. Bottom: Gaia error in total proper motion against the error in parallax. Both panels include all stars in our sample except for the seven targets with Gaia RUWE > 1.4 (see Sect. 2.3).

The G_mag of the stars in our working sample ranges between 5.1 and 9.7 mag. It has been shown that bright sources (G_mag < 6) also result in unreliable astrometric solutions because of uncalibrated CCD saturation (Lindegren et al. 2018). In the sample, four stars have magnitudes lower than 6, and they are discussed in detail in Sect. 5.1. In order to verify the Gaia DR2 parallaxes and proper motions for the brightest stars in the sample, we also retrieved the values provided in the HIPPARCOS (van Leeuwen 2007), and TGAS (Michalik et al. 2015) catalogs; however, the results were not better.

The Gaia errors in parallax range between 0.032 and 0.121 mas, while the errors in total proper motion range between 0.046 and 0.308 mas yr⁻¹. Six stars have uncertainties in parallax ∼0.08 mas or larger. The same have uncertainties in total proper motion larger than 0.18 mas yr⁻¹. They are all red Sgs except for HD 14489 (the A Sg in the upper right corner). The explanation for their large errors lies in the combined effect of large size and variability for the red Sgs, and the high brightness for HD 14489 (G_mag = 5.1). In both cases, the Gaia astrometric solution is affected (see Sect. 5.1). Of the red Sgs, HD 14528 (in the upper right corner) has the largest errors and also a relatively high RUWE value (1.25), followed by HD 14489, which in comparison has a RUWE = 0.81. In particular, for HD 14528, we adopt the results from Asaki et al. (2010) from this point on, who used the very long baseline interferometry (VLBI) technique to derive the astrometric parameters.

For the stars with RUWE > 1.4, errors in parallax range between 0.086 and 0.384 mas with a mean of 0.170 mas, and errors in total proper motion range between 0.171 and 0.501 mas yr⁻¹, with a mean of 0.347 mas yr⁻¹. We note that as expected, all these stars have larger errors than those associated with the main concentration of stars in the bottom panel of Fig. 3.

3. Radial velocity measurements

We first generated various suitable lists of spectral lines, optimized for the different SpT, using information available in the Atomic Line List interface⁵ (van Hoof 2018), and the SpectroWeb⁶ database. Each line list comprises a few to several dozen strong (log(gf) > −0.5), unblended lines covering the full 390−650 nm spectral window (or 510−870 nm in the case of the red Sgs).

For early and mid O-type stars, a few lines of N III-V and O III were used. In addition, we also included some He I lines to compensate for the lower number of available metal lines. For the late O-type stars, we added some lines of Si IV and O II. The situation improves for the B and A Sgs, were a much larger sample of lines is available, including lines from Si II-IV, N II-III, O II-III, S II-III, C II, Mg II, and Fe II. Last, in the case of red Sgs, we mostly used lines from Mg I, Ti I, Fe I, Ca I, Cr I, Ni I, and V I.

We then used our own tool (developed in Python 3.6) to perform a RV analysis. For each star, the corresponding list of lines was selected based on its SpT. For each line, an iterative normalization of the surrounding local continuum was made. Then, each line was fit to either a Gaussian or a Gaussian plus a rotational profile, depending on the first estimate of the full width at half maximum (FWHM) of the line. The measured central wavelength was then used to calculate the RV of each individual line in the initial line list (see above). From all the identified lines we removed those with equivalent widths lower than 25 mÅ directly before we carried out an iterative sigma clipping (using a threshold of 2σ) to remove potential poorly fit lines or incorrect identifications. The RVs of the surviving lines were then averaged, and we calculated the standard deviation of the final RV. This process was repeated for each spectrum and for each star in the sample.

The measurements of the individual RVs, together with the number of lines used for each spectrum, are listed in Table A.5. For O-type stars, the average number of lines is 12, the final average number of lines after sigma clipping is 6, and the typical uncertainties associated with the dispersion of RV measurements obtained after sigma clipping is ∼3.9 km s⁻¹. For the B-type stars, these values are 37 and 22 lines and ∼0.9 km s⁻¹, respectively. For A/F-type stars, they are 42 and 32 lines, and ∼0.26 km s⁻¹. Finally, for the K/M-type stars, they are 31 and 24 lines, and ∼0.17 km s⁻¹. This error is larger for the O-type stars for two main reasons: the first is that fewer lines are available, and the second reason is related to the broadening of the diagnostic lines, which is much larger for the O-type stars than in the cooler B, A, and red Sgs.

The RV results for the spectra with the best S/N are shown in the last column of Table A.1. For each star, we also searched for double-lined spectroscopic binaries (SB2) by looking at different key diagnostic lines (e.g., He Iλ5875, Si IIIλ4552, O IIIλ5592, C IIλ4267, and Mg IIλ4481).

We were able to measure individual RVs for the two components in three of the five SB2. We used the spectrum of maximum separation between them. Their values are listed in Table A.1.

For each star with four or more spectra, an average RV was calculated as the mean of the RVs obtained for each individual spectrum. In addition to the associated standard deviation, the peak-to-peak amplitude of variability in RV (RV_PP) was calculated as the difference between the highest and lowest individual RVs, and its error was calculated as the square root of the sum of the squares of the individual uncertainties. The results for the stars for which multi-epoch spectroscopy is available are listed in Table A.2.

Last, we also visually inspected the line-profile variability in each star with available multi-epoch spectroscopy. By doing this we were able to identify those cases in which any detected variability is more likely due to stellar oscillations than to (single-line) spectroscopic binarity (see Sect. 4.4.2).

4. Results

Figure 4 summarizes all the compiled information on astrometry and RVs (except for the information we extracted from the multi-epoch spectroscopy). The top panel of the figure shows the position of the stars in the sky, and the corresponding proper motions are indicated with arrows. For reference, we also indicate the location of the h and χ Persei double cluster (green circles at the center of the image) and the Galactic plane (dashed yellow line).

Fig. 4. Top panel: sky map of all the stars in the sample overplotted over a DSS-red image of the region. Dashed yellow line indicates the galactic plane, and the central green cross marks the center of Per OB1 (as defined in Mel’nik & Dambis 2017). Green circles indicated the location of the h and χ Persei double cluster. Colored vectors indicate the individual proper motion of each star. Middle and bottom panels: parallax and RV of the spectrum with the highest S/N, respectively, for each star in the sample against their position in right ascension. Open circles and square symbols indicate stars that deviate more than 2σ from the mean of the distribution of parallaxes and RVs, respectively (see Sects. 4.1 and 4.4.1). Stars with bad astrometry (see Sect. 2.3) are indicated with a plus, and no proper motion vectors are overplotted.

This image is complemented with another two panels, in which the distribution of parallaxes and RVs (as derived from the best S/N spectrum of each star) is plotted against the right ascension (middle and bottom panels, respectively). These two panels allow us to better identify the location in the sky of the outliers of both distributions, and to easily connect the information of the three investigated quantities.

From a first visual inspection of this summary figure, it becomes clear that generally speaking, the stars in our sample (including those located in the h and χ double cluster) belong to a connected population in terms of proper motions, parallaxes, and RVs. In addition, there is a non-negligible number of outliers that we discuss in detail in the next sections. They are potential nonmembers of the Per OB1 association, and/or runaway stars and binary systems.

4.1. Parallaxes and proper motions

Figure 5 shows again the results for proper motions and parallaxes (ϖ) from a different perspective. The central panel of the figure depicts the combined distribution of these two quantities, this time using the modulus of the proper motion (μ), defined as the square root of the sum of the squares of the proper motion in right ascension and declination. Stars labeled “unreliable astrometric solution” (see Sect. 2.3) are excluded from this figure.

Fig. 5. Total proper motions against parallax for the sample of stars except for those labeled “unreliable astrometry” (see Sect. 2.3). The 2σ boundaries of the distribution are shown as a rectangle. Empty colored circles show outliers of the distribution of any of the two quantities, and the associated uncertainties are overplotted. The mean and standard deviation obtained from the stars within the 2σ box are shown in the top right corner.

Most of the stars are grouped together around ϖ ≈ 0.4 mas and μ ≈ 1.2 mas yr⁻¹. This is also shown in the left and bottom panels of Fig. 5, where histograms of both parallax and total proper motion are shown.

An iterative 2σ clipping of these distributions results in ϖ = 0.40 ± 0.07 mas, and μ = 1.22 ± 0.26 mas yr⁻¹, and the identification of a total of 18 outliers (i.e., deviating more than 2σ from the mean of the distribution). The 2σ boundaries of the distribution (0.265 < ϖ < 0.540 mas, and 0.706 < μ < 1.740 mas yr⁻¹, respectively) and the outliers are highlighted in Fig. 5. The latter are also indicated in the second and third columns of Table A.4 and are discussed in Sect. 5.1.

These results assume that no different local substructures exist in the region, especially in terms of parallax. To investigate this statement further, we show again in Fig. 6 an image of the region with the proper motions overplotted, but this time using the mean proper motion obtained by considering the 16 stars located within 15 arcmin from the center of h and χ Persei, respectively, and having reliable astrometry (see the black arrow in the bottom right corner of the figure, corresponding to μ_α cos δ = −0.47 and μ_δ = −0.99 mas yr⁻¹).

Fig. 6. Same as the top panel of Fig. 4, but this time, the individual proper motion of each star is referred to the mean proper motion of 16 stars in h and χ Persei with good astrometric solution (black arrow in the bottom right corner of the top panel).

This figure is complemented with the information provided in Table 1, where we summarize the resulting means and standard deviations of parallaxes and proper motions when the sample is divided into circular regions around the center of h and χ Persei. The first region only includes the double cluster. The other regions extend outward by one degree each, starting at a distance of 30 arcmin from the center of the double cluster.

Based on the results presented in this section, we conclude the following: There is some empirical evidence of the existence local substructures in the spacial distribution of proper motions (see further discussion in Sect. 5.2.1). These subgroups of stars have a compatible distribution of parallaxes and proper motions. As a result, this justifies the decision to use the whole sample of stars to obtain the mean values and standard deviations of these two quantities to characterize this population of stars, as well as to identify potential outliers in parallax (i.e., nonmembers) and proper motion (i.e., runaway stars).

4.2. Comparison with previous works

4.2.1. Distance

We have obtained an average value for the parallax of ϖ = 0.398 ± 0.066 mas (adopting a zero offset of −0.03 mas). This value represents the mean of all stars in the sample with good astrometric solution that are not outliers in parallax and total proper motion.

Based on the corrected computed distances to these stars from Bailer-Jones et al. (2018), we obtain an average of d = 2566 ± 432 pc. This is compatible with the distance obtained using the inverse of our derived parallax: d = 2510 ± 415 pc. When we assume this distance, the projected distance extends up to ∼180 pc for the furthest stars in the association. In particular for the stars in the double cluster used in Table 1, we obtain a distance of d = 2340 ± 328 pc using the inverse of the derived parallax.

This result agrees well with previous estimates for h and χ Persei using different approaches. In addition, it also indicates that the parallax zero-point offset correction proposed by Lindegren et al. (2018) is adequate.

To give some examples, Uribe et al. (2002) obtained an average double cluster distance of d = 2014 ± 46 pc using the ZAMS fitting approach. Currie et al. (2010) used main-sequence stars with a very large sample and obtained a distance to each cluster of d_h = 2290${}_{-82}^{+87}$ pc and d_χ = 2344${}_{-85}^{+88}$ pc. The previously mentioned work by Asaki et al. (2010) estimated a distance to HD 14528 of d_h = 2420${}_{-90}^{+110}$ pc using high-precision interferometric observations. More recently, Babusiaux (2018a) published mean parallaxes for a broad selection of open clusters using Gaia DR2 including h and χ Persei. By applying a −0.03 mas zero-point offset, they obtained d_h = 2239 pc and d_χ = 2357${}_{-85}^{+88}$ pc. Finally, Davies & Beasor (2019) estimated the distance to h Persei in d_h = 2250${}_{-140}^{+160}$ pc, adopting an offset of −0.05 mas for the Gaia parallaxes.

The aim of this work is not to provide a better estimate, but to ensure that the stars selected here based on their parallax belong to the association. Only a few works provide distances to the Per OB1 association. For instance, Shull & Danforth (2019) used the photometric distance and Gaia parallaxes for a selection of O-type stars to derive a distance to the association of d = 2.99 ± 0.14 kpc and d = 2.47 ± 0.57 kpc, respectively.

4.2.2. Proper motions

For the stars that are not outliers in proper motion and parallax, we obtain mean values and standard deviation for the individual components of the proper motion of μ_α cos δ = −0.51 ± 0.48 mas yr⁻¹, μ_δ = −1.00 ± 0.31 mas yr⁻¹. This result agrees quite well with previous results obtained in the literature by other authors and different samples of stars. For example, Zhong et al. (2019) investigated a sample of more than 2100 stars (covering a much wider range in mass than our study) located within 7.5 degrees around the h and χ Persei double cluster. They found for each cluster μ_α cos δ = −0.71 ± 0.18 mas yr⁻¹ and μ_δ = −1.12 ± 0.17 mas yr⁻¹, respectively. Similar results were also obtained by Mel’nik & Dambis (2017) and Li et al. (2019).

4.3. Completeness of the sample

As indicated in Sect. 2.1, our sample of 88 stars includes almost all blue and red Sgs (LC I and II) quoted in Humphreys (1978), Garmany & Stencel (1992), Currie et al. (2010), and Gazak et al. (2014), plus a few LC III objects (Gs) with late-O and early-B spectral types. In particular, from a total of 107 targets quoted in these four papers that meet our selection criteria, we only lack spectra for 12 B and 5 M Sgs.

To further evaluate the completeness of our sample, we benefit from photometry provided by Gaia and the results about parallaxes and proper motions described in Sect. 4.1. To this aim, we retrieved all the stars in the Gaia DR2 catalog with G_mag brighter than 10.5 whose parallaxes and total proper motions lie within 2σ of the distributions depicted in Fig. 5. We then removed all stars with RUWE larger than 1.4, and those classified by the SIMBAD Astronomical Database⁷ as dwarfs or subgiants (luminosity classes V and IV).

The results are presented in a color-magnitude diagram (CMD) in Fig. 7, where we use the same color-code as in previous figures for the stars in our sample, but this time, we also highlight the 17 stars that are classified as LC III or II-III stars in light green. For reference purposes, we also include a A_v = 1.7 mag reddening vector and two reddened 14 Myr isochrones⁸ (solid lines) shifted to a distance of 2.5 kpc (or, equivalently, a distance modulus of 12 mag). The values of reddening for the isochrones (A_v = 1.0 and 2.7, respectively) were selected to embrace the main-sequence band, corresponding to the region of the CMD with higher density of gray points in the bottom left corner.

Fig. 7. Color-magnitude diagram (using Gaia photometry) of stars located within 4.5 degrees from the center of the Per OB1 association. Colored stars shows the stars in our sample, gray circles represent the remaining stars from Gaia, green diamonds show 17 blue and red supergiants quoted in the literature for which we lack spectra (see the last part of Table A.1). Two isochrones and a reddening vector are also included for reference purposes. See Sect. 4.3 for explanation.

From inspection of this figure we can conclude that the level of completeness in our sample is very high, specially when we concentrate on the region of the CMD where the blue and red Sgs are located (purple, dark blue, cyan, and red stars). Interestingly, we also find that a high percentage of the 12 B Sgs quoted in Humphreys (1978), Garmany & Stencel (1992), and Currie et al. (2010) are likely B Gs, instead of B Sgs. These refer to all green diamonds with G_mag < 9, most of them classified as B Sgs in Currie et al. (2010) (see the last rows of Table A.1).

In addition, Fig. 7 allows us to conclude that the blue and red Sg population of Per OB1 is affected by a variable reddening that ranges from A_v ∼ 1.0 to 2.7 mag (in agreement with previous findings by Li et al. 2019), and that the age associated with the blue and red Sg population is not compatible because the higher mass present in the 14 Myr isochrone is ∼14 M_⊙, while all O, B, and A Gs/Sgs included in our sample are expected to have masses higher than 20 M_⊙. This latter result will be further investigated in the next paper of this series, after information about the stellar parameters of the full working sample is included.

4.4. Radial velocities

By following the strategy described in Sect. 3, we obtained RV estimates for all the available spectra in our sample of stars. These measurements are used (1) to investigate the RV distributions resulting from the best S/N spectra, (2) to provide empirical constraints on intrinsic spectroscopic variability typically associated with the various types of stars, and (3) to identify spectroscopic binaries and runaway candidates.

4.4.1. Best S/N spectra

The bottom panel of Fig. 4 shows the RVs of all stars in the sample, obtained from the best S/N spectra, as a function of the position of the stars in right ascension. The associated distributions, this time separated by SpTs, are depicted in the form of histograms in Fig. 8, with the mean and standard deviation associated with each RV distribution (after performing an iterative 2σ clipping) indicated at the top of the various panels. The corresponding outliers in each distribution are indicated as open squares in the bottom panel of Fig. 4 and listed in the fourth column of Table A.4.

Fig. 8. RV distributions associated with the different SpT groups resulting from the analysis of the best S/N spectra. The orange bin in the second panel from the bottom is HD 12842, the F Sg.

From a visual inspection of Fig. 8 we can conclude that except for the case of O-type stars, which has a flatter and more scattered distribution, the other three distributions are quite similar (when the outliers are eliminated), following a more or less clear Gaussian shape. (For the A/F-type stars, only the two situated on the right-most side of Fig. 8 are outliers. The consequence of having fewer stars than for the B- and K/M-type stars results in a poorer Gaussian shape). The mean values of these three distributions are compatible within the uncertainties, with a difference smaller than 2 km s⁻¹. Interestingly, the standard deviation of the distributions significantly drops from O- to B- and A-type stars, and continues to decrease to the K/M-type stars (see further notes in Sects. 4.4.2 and 4.4.3).

4.4.2. Multi-epoch spectra: intrinsic variability

As indicated in the bottom panel of Fig. 2, we have more than one spectrum for 73 of the stars in the sample. These observations can be used to identify binaries; however, as extensively discussed in Simón-Díaz et al. (in prep.), the effect of intrinsic variability also needs to be taken into account to minimize the spurious detection of single-line spectroscopic binaries (SB1) in the blue supergiant domain (see also further notes regarding the red supergiant domain in Patrick et al. 2019, 2020).

Some examples of the type of spectroscopic variability phenomena produced by stellar oscillations or the effect of a variable stellar wind in the OBA Sg domain can be found in Fullerton et al. (1996), Prinja et al. (2004, 2006), Kaufer et al. (2006), Kraus et al. (2015), Simón-Díaz et al. (2017, 2018), Aerts et al. (2017, 2018), for example. This effect is also illustrated in Fig. 9 using a subsample of 15 stars in PerOB1 for which five or more spectra are available, and whose detected variability in RV is more likely produced by intrinsic variability than by the orbital motion in a binary system (see Table A.2 and further notes in Sect. 4.4.3).

Fig. 9. Measured RVs (subtracted from their mean) for a sample of 15 stars (ordered by SpT) for which we have five or more spectra, and whose detected variability in RV is more likely produced by intrinsic variability than by the orbital motion in a binary system.

These results warn us about the dangers of using a single snapshot observation to associate the outliers detected in the RV distributions shown in Fig. 8 with potential runaway stars or spectroscopic binaries. Some of these cases might even correspond to a single measurement in a specific phase of the intrinsic variability of the star instead of being associated with the orbital motion in a binary system or with a single star with an anomalous RV due to an ejection event. They also partially explain why the standard deviation of the RV distributions presented in Fig. 8 becomes smaller when moving from the blue to the red Sgs. This is just a consequence of the behavior of the characteristic amplitude of spectroscopic variability with SpT (see Table 2 and Simón-Díaz et al., in prep.). Last, it also affects the fraction of detected SB1 stars using multi-epoch observations, or the final sample of outliers in RV (see further notes in Sects. 4.4.3 and 5.2.3, respectively).

To evaluate the effect that including information about multi-epoch spectroscopy has on the identification of outliers in the RV distribution, we have repeated the same exercise as in the case of the single-snapshot observations (Sect. 4.4.1), but modifying the individual measurements (obtained from the analysis of the best S/N spectra) of stars for which four or more spectra are available by the mean of the multi-epoch RV measurements. Results of this exercise are presented in the “RV multi” column of Table A.4. Although the number of stars with a modified outlier status in RV is small in this specific example (only HD 13402 and HD 12953), the results presented in Fig. 9 indicate that it could have been larger if other epochs of the time series had been selected as single-snapshot observations.

4.4.3. Multi-epoch spectra: spectroscopic binaries

Given the stability of the FIES and HERMES instruments, and the accuracy reached in the RV measurements for most of the stars in the sample with multi-epoch spectroscopy, it might be tempted to assign the SB1 status to all stars showing a RV_PP above a few km s⁻¹. However, as indicated in Sect. 4.4.2, the intrinsic variability in single blue supergiants can reach amplitudes of a few dozen km s⁻¹; hence, many of these identification may lead to spurious results.

To avoid this situation as much as possible, and in order to identify the most secure candidates to be SB1, we performed a careful inspection of the type of line-profile variability detected in each of the stars with more than one spectrum. To this aim, we mainly considered the following diagnostic lines, whenever available: He Iλ5875, Si IIIλ4552, O IIIλ5592, C IIλ4267, and Mg IIλ4481. For the case of the two red supergiants with multi-epoch observations, we found that the measured RV_PP is lower than 2 km s⁻¹, which we directly attribute to intrinsic variability.

The list of clearly identified SB1 is presented in Table A.3. In addition to the four SB1 stars quoted there, we found five SB2 systems (some of them directly detected from a single-snapshot observation) and labeled “LPV/SB1?” another five cases in which we are not entirely sure if the detected variability is due to binarity or intrinsic variability. All this information is also added to Table A.4 (column “Spec. variability”).

We also performed a bibliographic search for previously identified binaries in our sample of blue and red supergiants. We mainly concentrated in the works by Holgado et al. (2018, 2020) and Maíz Apellániz et al. (2019) for the case of O-type stars, and Abt & Levy (1973), Zakirov & Shaidullin (1985), and Laur et al. (2017) for the B supergiant sample. In addition, we made use of The International Variable Star Index (VSI)⁹.

In total, we found that six out of our sample of ten detected SB1 or SB2 systems from this work were previously identified in any of these references as such (HD 16429 is actually a triple system McSwain 2003). This implies four newly detected binaries: HD 13969, HD 14476, HD 17378 (all SB1), and HD 13402 (SB2). We also found three binaries in the literature that were not detected from our available spectra because of short time-coverage: BD +56578, an eclipsing binary (Tarasov et al. 2016; Laur et al. 2017), plus HD 17603 and HD 14956, identified as SB1 by Holgado et al. (2018) and Abt & Levy (1973), respectively. All of them are labeled “(lit.)” in the corresponding column of Table A.4.

The stars classified as “LPV/SB1?” are HD 13036, HD 13854, HD 13267, HD 14542, HD 12953, and HD 17378. HD 13036 (B0.2 III), HD 13267 (B6 Iab) and HD 14542 (B8 Iab) have RV_PP = 10−14 km s⁻¹; although this value is at the boundary of the expected variability due to pulsations, which may indicate an SB1 classification, we cannot conclude after visual inspection of their line-profile variability. HD 13854 (B1 Ia-Iab), mostly looks like a pulsational variable, but we do not discard the possibility entirely that this star might be a SB1. We note that Abt & Levy (1973) provide RV_PP = 24.8 km s⁻¹; however, they did not consider it as a binary. For HD 12953 (A1 Iae), we measured RV_PP = 10.4 km s⁻¹, the largest variability in the A supergiant sample; however, after visual inspection of its line profile variability, we cannot conclude whether this is a SB1 system. Abt & Levy (1973) found RV_PP = 15.8 km s⁻¹ for this star, which would favor that it is an SB1. HD 17378 (A6 Ia) has RV_PP = 8 km s⁻¹, which is large enough to consider it as potential binary. However, we only have three spectra.

Last, we found that although HD 14956 (B2 Ia) was classified as an SB1 with a period of P = 175 days, and RV_PP = 27.0 km s⁻¹ (Abt & Levy 1973), and Laur et al. (2017) classified this star as α Cygni variable, we do not see such signs of SB1 variations, as we measure RV_PP = 5.5 km s⁻¹. However, we do not have enough spectra (three) to discard this possibility.

These results about detected spectroscopic binaries, along with the RV distributions obtained from the analysis of the best S/N spectra (i.e., obtained from a single-snapshot observation), allow us to evaluate the extent to which these distributions can be used to identify spectroscopic binaries among the outliers. We find that only four out of all the SB1/SB2 systems detected by means of multi-epoch spectroscopy are outliers in the abovementioned distributions. In addition, some outliers have not been detected as spectroscopic binaries although more than four spectra are available for them (e.g., HD 13268, O8.5 IIIn, RV = −106.2 km s⁻¹). These results can be explained when we take into account (1) that the best S/N of some of the spectroscopic binaries correspond to an orbital phase in which the RV is close to the systemic velocity, and (2) that some outliers in RV might be runaways and not necessarily binaries. The latter situation is the case of HD 13268, a well-known runaway star (see also Sect. 5.1). This means that if a given star is an outlier in RV, it is useful to first investigate its runaway nature (by means of its proper motion) before marking it as a potential spectroscopic binary, and vice versa; for example, although the measured RV of the B1 Ib-II star HD 14052 (−90.8 km s⁻¹) deviates more than 3σ from the mean, this star is not an outlier in proper motion, and so we may conclude that it is more likely a spectroscopic binary than a runaway. This is confirmed through access to multi-epoch spectroscopy. We further discuss the percentage of spectroscopic binaries in our sample of stars in Sect. 5.2.3.

5. Discussion

Table A.4 compiles and summarizes some information of interest for the discussion about membership and final identification of spectroscopic binaries and runaway stars. Columns “ϖ” and “μ” indicate if a given star is part of the bulk distribution of parallaxes and proper motions, respectively, or if it is detected as an outlier of these distributions (Sect. 4.1). Columns “RV best” to “RV final” provide similar information for the case of RV estimates obtained from the best S/N spectra (Sect. 4.4.1) for stars with four or more spectra (Sect. 4.4.2), or the final distribution of RVs (Sect. 5.2.1), respectively. In all these cases, different symbols are used to identify secure or doubtful cases.

For completeness, we also add to Table A.4 information about confirmed spectroscopic binaries, our final decision on cluster membership status (columns “Spec. variability” and “Member”), as well as some other comments of interest for the final interpretation of results (column “Comments”).

5.1. Cluster membership

As discussed in Sect. 4.1, most of the stars in the sample with good astrometry (81 stars) are grouped together in the proper motion versus parallax diagram. The mean and standard deviation of the distribution of these two quantities are ϖ = 0.40 ± 0.07 mas, and μ = 1.22 ± 0.26 mas yr⁻¹, respectively. All the stars that are located within the 2σ boundaries of the distribution (64 in total) are directly considered as members and labeled with filled circles in columns “ϖ” and “μ” of Table A.4. The remaining 17 stars are marked with an open circle or a cross in Table A.4 depending on whether they deviate by 2−3σ or more than 3σ, respectively. We note that in this case, columns “ϖ” and “μ” include information about the remaining seven stars that were not included in Fig. 5: those labeled “unreliable astrometry” (or RUWE > 1.4). Because the information about parallaxes and proper motions is uncertain for them, we exclude these stars for the moment and mark them using brackets surrounding the corresponding symbols in columns “ϖ” and “μ” of Table A.4.

These are not the only stars with unreliable parallaxes. Figures 3 and 5 include a small sample of 6 K/M-type supergiants that despite a RUWE value well below 1.4 have larger errors than the rest of stars in the sample, and interestingly, all of them are systematically shifted to larger parallaxes (although except for one, all have total proper motions within the 2σ boundaries and proper motion vectors compatible with the bulk of member stars, Fig. 4). They are also all marked with brackets in Table A.4.

This is likely connected to an already known problem that affects the reliability of the Gaia DR2 astrometric solution. In brief, as pointed out by Pasquato et al. (2011) and Chiavassa et al. (2018), the position of the centroid changes on timescales of several months or a few years because of the large size and strong intrinsic photocentric variability of red supergiants. This effect leads to unreliable parallaxes and errors.

A particular example of interest regarding this issue with the astrometric solution of Gaia for the case of red supergiants is the highly variable star HD 14528 (S Per, ϖ_Gaia = 0.25 ± 0.12 mas, μ_total Gaia = 2.57 ± 0.31 mas yr⁻¹). This star has an average angular size of 6.6 mas (Richards et al. 2012). It was monitored for six years by Asaki et al. (2010) with VLBI. The authors obtained an independent parallax of 0.413 ± 0.017 mas, which is just at the center of the distribution. We therefore cannot discard completely that these six K/M-type supergiants, which are outliers in parallax using data from Gaia DR2, are members of Per OB1.

The last star that we place in brackets is the A-type supergiant HD 14489. This is the brightest stars in our sample, with G_mag = 5.1. As shown in Fig. 3, this star also has much larger errors in parallax and proper motions than the bulk of stars in the sample. This may be related to the current limitation of Gaia DR2 regarding the reliability of the astrometric solutions for stars brighter than G_mag ≲ 6 (Lindegren et al. 2018). Another three stars share this issue, but their associated astrometric errors are much smaller and their magnitudes are close to G_mag = 6; therefore we decided to consider their astrometric solutions reliable.

Taking all this information into account, we decided to following strategy below to evaluate the membership to Per OB1 of each star in our sample. Stars with reliable values on parallax and proper motion (i.e., not marked with parentheses in Cols. 2 and 3 of Table A.4) are considered as confirmed members if they do not deviate more than 2σ from the mean of the distribution of parallaxes. For stars with unreliable values of parallax and proper motion (i.e., highlighted with brackets in columns “ϖ” and “μ” of Table A.4), we adopted the following: if they are not outliers in parallax, they are considered likely members; if they are outliers in parallax, we consider them candidate members, except for the K/M-type stars, which remain likely members because of the arguments provided above. Last, stars with reliable astrometry that are outliers in parallax (as well as those stars in IC 1805, see below) are considered nonmembers.

Most of the stars are properly classified using these criteria. However, a few cases deserve further attention.

HD 13022 (O9.7 III) and HD 12842 (F3 Ib). These two stars are classified as members following the guidelines above, but they are outliers in proper motion (Fig. 5). Interestingly, they have a very small proper motion compared to the rest of the stars in the sample (see in Fig. 4 the two stars with very small vectors located at (RA, Dec) ∼ (32, 58.5) deg). Awaiting a more detailed study of these two stars, we continue considering them members for the moment.

HD 16691 (O4 If), HD 15642 (O9.5 II-IIIn), HD 13745 (O9.7 II(n)), and HD 13268 (ON8.5 IIIn). These four O-type stars are clear outliers in proper motion (see Figs. 4 and 5). We consider the first three runaway members because their parallaxes lie within the 2σ boundaries. The fourth (HD 13268) is an interesting case; although this star has a somewhat larger parallax, it has a RV of ∼105 km s⁻¹. Therefore, given its spectral classification and this high RV pointing to us, it can still be considered a runaway member of Per OB1. This star is a well-known fast-rotating nitrogen-rich O-type runaway star (e.g., Abt et al. 1972; Mathys 1989; Simón-Díaz & Herrero 2014; Martins et al. 2015; Cazorla et al. 2017a,b).

HD 14322 (B8 Iab). This star is an outlier in parallax with a value of ϖ = 0.21 ± 0.04 mas. Although the TGAS catalog provides a value for it of ϖ = 0.44 ± 0.38 mas (within the boundaries of ϖ), the error is much larger. This inconsistency caused us to modify its status from nonmember to member candidate while awaiting Gaia DR3.

HD 14489 (A1 Ia). This is a bright A-type star (G_mag = 5.1), outlier in parallax, and with the largest parallax error. Although it has a RUWE = 0.81, we do not trust its Gaia astrometry, as explained before, because of its brightness. The result from TGAS provides a parallax of ϖ = 0.45 ± 0.94 mas, and although it is within the adopted boundaries of Per OB1, the error is very large. This star is also an outlier in RV and close to the 2σ boundary in proper motion. Therefore we decide to label it a runaway member candidate.

BD +56724 (M4-M5 Ia-Iab). This star has the largest parallax in Fig. 5, and a RUWE = 0.93. Although the reliability of Gaia DR2 parallaxes for the K/M-type stars may be low, its large deviation from the mean of the distribution could mean that this star is not a member. It is also an outlier in proper motion, but its magnitude and RV are similar to other red supergiants in the sample. We therefore retain this star as member candidate for the moment.

HD 15570 (O4 If), HD 15558 (O4.5 III(f)), HD16429 (O9 II(n)), and BD +60493 (B0.5 Ia). All these stars are located within or in the surroundings of IC 1805. Interestingly, all of them but one are O-type stars. Although they are located within the 2σ boundaries of the parallax and proper motion distribution (except for HD 16429, but this is a triple system with a RUWE = 8.8), we decided to mark them nonmembers based on their separated location in the sky and their direct connection with the surrounding H II region. They seem to be linked to a younger star-forming region located at higher galactic latitudes (but at the same distance). Most of them are also outliers in RV (see Table A.4), but this is likely due to their binary nature.

The final result of this classification, also taking into account the comments on some individual stars presented above, is summarized in column “Member” of Table A.4. In total, we have 70 confirmed members, 9 likely members, 5 member candidates, and 4 nonmembers. Interestingly, only stars in IC 1805 are finally classified nonmembers. The remaining 84 stars likely belong to the Per OB1 association (although some of them are identified as runways, see Sect. 5.2.3).

5.2. Kinematics

In Sects. 4.1 and 4.4 we provided a global overview of the results about proper motions and RVs for the complete sample of stars, also including some information about identified spectroscopic binaries. In this section we discuss these results more in detail. We also refer to Mel’nik & Dambis (2017), Zhong et al. (2019), Melnik & Dambis (2020) for complementary (and in some cases more detailed) information about the global and internal kinematical properties of stars in the Per OB1 association.

5.2.1. Proper motions

Figure 5 and the top panel of Fig. 4 provide a global overview of the distribution of proper motions in the whole star sample. These figures show that (except for a few outliers) most of the stars in our sample that are located below the Galactic plane (among them, those in h and χ Persei) can be considered a dynamically connected population of stars. This result perfectly agrees with previous findings by Lee & Lim (2008). Using proper motions from the HIPPARCOS mission, these authors showed that the luminous members of the Per OB1 association exhibit a bulk motion away from the Galactic plane, such that their average velocity increases with height above the Galactic plane.

Furthermore, inspection of the results for proper motions (relative to h and χ Persei) and parallaxes presented in Fig. 6 and Table 1 allows us to conclude that the distributions of parallaxes associated with stars located at increasing distances from the center of the double cluster are all compatible (at least we find no clear subgrouping in terms of parallax at least given the accuracy of Gaia DR2 astrometry – except maybe the stars in h and χ Persei because their parallax is somewhat larger or their distance is somewhat closer). We also conclude that the mean and standard deviation of the distribution of proper motions in h and χ Persei is much smaller than in the more extended population.

Moreover, the spatial distribution of proper motions in the extended population of blue and red supergiants in Per OB1 does not follow an expanding structure centered in the h and χ Persei double cluster. Instead, the local proper motions of most of the stars located north of these clusters seem to point outwards from an imagined center located at about 1 degree north of the double cluster (see also Mel’nik & Dambis 2018; Melnik & Dambis 2020). These results from the proper motion, linked with the results by Lee & Lim (2008) mentioned above, are compatible with a scenario in which the halo population of blue and red supergiants around the double cluster has been formed from a more diffuse region of interstellar material compared to the denser region associated with the clusters themselves.

In addition, four O-type stars south of the region can be clearly considered runaways based on the size and direction of their proper motions (see also Fig. 5). Interestingly, their proper motion vectors do not point outward from h and χ Persei, but to a far more extended region of the Galactic plane (see also the discussion in Sect. 5.2.3). Finally, as also indicated in Sect. 5.1, the stars located within or near IC 1805 likely belong to a younger population of stars that is not necessarily connected with the remaining stars in Per OB1.

5.2.2. Radial velocities

Figure 8 and the bottom panel of Fig. 4 summarize the RV results obtained with the best S/N spectra. The analysis of these spectra has allowed us to characterize the RV distributions for the different SpT groups and to identify potential spectroscopic binary systems and runaway stars among the outliers of the distributions (see columns “RV best” and “RV multi” in Table A.4). We then illustrated in Sect. 4.4.2 (see also Fig. 9) the importance of incorporating information from the analysis of multi-epoch observations for the correct interpretation of the RV distribution, and in particular, to avoid the spurious identification of spectroscopic binaries (either from a single epoch or from multi-epoch observations) due to the effect on the measured RVs of the intrinsic variability caused by stellar oscillation and/or wind variability in the blue supergiant domain.

Last, we learned that after eliminating outliers associated with confirmed spectroscopic binaries (via multi-epoch spectroscopy) and runways (via proper motions), the RV distributions for the B, A/F, and K/M Sgs are fairly compatible in terms of mean values and standard deviations. In addition, we found that most of the O-type stars in the sample are either (1) runaways, as detected from the proper motions, (2) spectroscopic binaries, or (3) are considered nonmembers because they are located nearby IC 1805, far away from the main distribution of stars in Per OB1. As a result, the RV distribution of the O-type sample is remarkably broader than for those associated with the other SpT.

We now take all these results into account to provide final information about RVs in Fig. 10 and in the column “RV final” of Table A.4. To do this, we first replaced the list of measurements obtained from the best S/N spectra by the mean value resulting from the analysis of the multi-epoch observations for those cases for which we have more than one spectrum. Then we used this list of values, except for all the SB2 binaries, and the stars identified as nonmembers (see Sect. 5.1), to obtain the mean and standard deviation by performing an iterative 2σ clipping.

Fig. 10. Radial distance from the center of Per OB1 against the mean RV for the stars with more than one spectrum, or RV for the stars with only one spectrum, excluding the SB2 binaries. The stars within 1 degree from IC 1805 were also excluded. The filled gray rectangle shows 2σ of all the RVs, excluding those from stars identified as binaries (SB1 or SB2). The filled colored circles denote stars identified as SB1. The filled colored crosses denote starts outside 2σ from the mean that are not identified as binaries. The colored stars show the remaining stars. The colored error bars present RVPP for stars with multi-epoch data, except for SB1. The open triangles show stars that are outliers in proper motion and were therefore identified as runaways.

The results of this process are presented in Fig. 10, where the RVs of all stars that were not excluded from the list are presented against the radial distance from the center of Per OB1. The obtained mean and standard deviation are shown in the top right corner (−42.9 ± 3.5 km s⁻¹) of the figure, and the horizontal gray band indicates the 2σ boundaries.

Although most of the stars in the sample are concentrated within the central 100 pc, we observe that except for a few cases, the remaining stars also lie within the 2σ boundaries. Therefore, once more, and as was suggested by Currie et al. (2010) and Zhong et al. (2019), the extended population of blue and red Sgs in Per OB1 (up to 200 pc, i.e., relatively far away from h and χ Persei) seems to have a common origin in terms of kinematics. No global gradient (as a function of distance to the center of the association) or local substructures are observed in the distribution of RVs.

As indicated above, we obtain $\overline{\mathrm{RV}}$ = −42.9 ± 3.5 km s⁻¹ using the whole sample of stars that are not excluded from the list. Regarding h and χ Persei, we obtained average values of RV_χ Per = −44.4 ± 1.4 km s⁻¹ and RV_h Per = −41.1 ± 2.6 km s⁻¹, respectively. These results agree well with those previously obtained by other authors. For the association as a whole, Mel’nik & Dambis (2017) provided a mean value of $\overline{\mathrm{RV}}$ = −43.2 ± 7.0 km s⁻¹ using available information of member stars from the TGAS catalog. For the individual clusters, Liu et al. (1991) provided RV_χ Per = −44.4 ± 0.7 km s⁻¹ and RV_h Per = −46.8 ± 1.7 km s⁻¹, respectively, using a sample of cluster stars (mainly early type, more specifically, B- and A-type stars). In the particular case of h Per, the fact that we have only three suitable stars to compute the mean may explain the poorer agreement.

As in the case of the analysis of the best S/N spectra and the multi-epoch observations, we provide in column “RV final” of Table A.4 a list of identifiers to separate the outliers of the distribution of final values of RV from stars within the 2σ boundaries. This information is used in the next section to determine additional SB1 stars that have not previously been identified based on the available multi-epoch spectra.

5.2.3. Runaway and binary stars

The last two columns of Table 3 summarize the final number and percentage of identified runaways in each of the four SpT groups. We proceed as follows to obtain this. First, we assigned the runaway status to all outliers in proper motion (meaning that the magnitude of any of their individual components, or the total proper motion, deviates by more than 2σ from the mean of the corresponding distribution). Then we considered the possibility of identifying additional runaway stars through their RVs. In this case, we decided to only label them (if not detected as SB2) clear runaways if their RV deviates by more than 15 km s⁻¹ from the mean of the final RV distribution presented in Fig. 10. We made this decision based on two arguments: the first refers to the result presented in the fifth column of Table 3; namely, the number of identified runaways when all outliers in RV are considered that deviate by more than 2σ is too large when compared with those detected through PMs. The second argument is based on the results of RV_PP that is expected to be produced by intrinsic variability, which can be as high as 10−20 km s⁻¹ in some cases (see Table 2).

These arguments are supported by the fact that it is very unlikely to find a runaway star that is an outlier in RV, but not in at least one of the components of the proper motion. In contrast, as described above, intrinsic variability can lead to single-snapshot RV measurement that can easily deviate by up to 10−20 km s⁻¹ (or, equivalently, about 4−5σ in this specific sample of stars). Furthermore, this situation can be even more dramatic for large-amplitude SB1 systems for which only a low number of spectra is available. It is therefore more likely that a star that is not an outlier in proper motion but is an outlier in RV is an SB1 than a runaway. Alternatively, if the deviation in RV is smaller than the typical intrinsic variability corresponding to the associated SpT, it might not even be a binary star.

A practical example of the latter situation is the BN2 II-III star BD +56578, for which we only have three spectra that cover a very short time-span (one day). This star is not an outlier in proper motion, but is an outlier in RV (it deviates by 13σ from the mean of the RV distribution). Based on what we have described above, this star should be labeled as a potential binary, a suspicion that is confirmed from the literature (Laur et al. 2017).

Following these arguments, all stars in Fig. 10 whose RV measurements deviate by more than 2σ (i.e., which lie outside the gray band) and up to 10−20 km s⁻¹ and that have not previously been detected as runaways through proper motions (open triangles) or as spectroscopic binaries through multi-epoch spectroscopy (filled circles) are quite likely single pulsating stars.

Overall, we identify a total of 11 runaway stars. The group of stars with a larger number of runaways (45%) are the O-type stars. This is followed by the B and K/M Sgs, with 5% each. (We note, however, that the runaway status of the M Sg BD +56724 can be a spurious result because it is based on the Gaia DR2 proper motion, which may not be as reliable as for the other stars because of the problems regarding size and variability of the red Sgs). Last, the lower percentage of runaways is found for the A Sgs, with only 1 or 2%, depending on whether we trust the Gaia DR2 proper motion of the bright star HD 14489, which also has a much larger parallax than the remaining stars in Per OB1.

It thus becomes clear again that a high percentage of the O-type stars in the Per OB1 region can be considered a dynamically distinct group. However, in contrast to previous assumptions (see, e.g., Walborn 2002), the fact that all of them are found within the 2σ boundaries of the parallax distribution indicates that they belong to the same grouping as the remaining blue and red supergiants in Per OB1, and not to a more distant, dispersed association. Although further confirmation is needed, the most likely origin of the O-star runaways is a dynamical kick by a supernova explosion in a previously bounded binary system. This hypothesis is reinforced by the fact that none of the detected runaways are identified as binary systems (and the other way round).

Table 4 summarizes the results for the detected binaries, again separated by SpT group, and this time differentiating the B Sgs from the B Gs because they represent the evolutionary descendants of main-sequence stars in two different mass domains. We refer to Sect. 4.4.2 for a description of how the SB1 and SB2 stars where identified. The targets labeled “SB1?” include targets fulfilling any of the two following criteria. On the one hand, stars with five or more spectra for which we cannot clearly decide whether the detected line-profile variability is due to intrinsic variability or orbital motion. On the other hand, following the arguments above, we identified stars as “SB1?” whose RV_PP is larger than the typical intrinsic variability expected for their SpT (see Fig. 9).

The main conclusions from inspection of the results presented in Table 4 (and Table 3) are summarized as follows. First, the percentage of detected spectroscopic binaries decreases toward later SpT, or equivalently as the massive star evolution proceeds. This result agrees with recent findings by Barbá et al. (2017), Patrick et al. (2019, 2020), Simón-Díaz et al. (in prep.). When we assume that the detected runaways indicate a past binary evolution, the total percentage of clear binaries (excluding those labeled “SB1?”) would decrease from ∼60% to ∼15% when the O star and B Sg samples are compared, and further below ∼5% when the cooler Sgs are considered. Second, while the decreasing tendency remains in both cases, the exact behavior of the percentage of detected spectroscopic binaries is different depending on whether we also include the stars labeled “SB1?” stars. Therefore it is critical to confirm or dismiss our suspicion that most of the stars with RV_PP below 10−15 km s⁻¹ are actually single pulsating stars and not spectroscopic binaries. Access to multi-epoch data for the whole sample of star is therefore crucial to obtain reliable empirical information about the relative percentage of binaries throughout the massive star evolution. Finally, as an aside, the percentage of spectroscopic binaries is much higher among the B Gs than in the B Sgs.

This clearly shows that any further attempt to interpret the empirical properties of this sample of massive stars in an evolutionary context must take into account that a large fraction of the O stars is or likely has been part of a binary or multiple system. In addition, some of the other more evolved targets may also have been affected by binary evolution.

6. Summary and future prospects

Our study has provided all the necessary environmental information that will be used in a forthcoming paper, in which we will also incorporate results obtained from a quantitative spectroscopic analysis of the whole sample (including stellar parameters and surface abundances) to perform a complete homogeneous characterization of the physical and evolutionary properties of the massive star population of the Per OB1 association.

In this paper, we have studied a sample of 88 massive stars located within 4.5 deg from the center of the Per OB1 association using high-resolution multi-epoch spectroscopy, and astrometric information from the Gaia second data release (DR2).

We have investigated membership of all star in the sample to the Per OB1 association, resulting in 70 members, 9 likely members, and another 5 candidates that require further investigation, while the other 4 were considered nonmembers as they belong to IC 1805.

We have found eight clear and two likely runaway stars, most of them O-type stars. We also identified 5 SB1 and five SB2 stars (these include three and one new binary systems, respectively), plus another 11 potential SB1 stars that we propose are single pulsating stars.

To obtain these results, we took their parallaxes and proper motions (as compiled from Gaia DR2) into account, and the RV estimates obtained from the available multi-epoch and/or single snapshot spectra. In addition, we also considered the reliability of the astrometry provided by Gaia through the RUWE value, the potential decrease in reliability of Gaia astrometry in the case of the red Sgs because of their large size and photocentric variability, and the expected amplitude of spectroscopic variability produced by stellar pulsations and/or wind variability when spectroscopic binaries are identified based on their RV measurements.

We have also analyzed some global properties of the sample and obtained averages in parallax, total proper motion, and RV of ϖ = 0.40 ± 0.07 mas, μ = 1.22 ± 0.26 mas yr⁻¹ (μ_α cos δ = −0.50 ± 0.48, μ_δ = −0.99 ± 0.31), and −42.9 ± 3.5 km s⁻¹. All these results agree relatively well with previous studies based on different stellar samples comprising the Per OB1 association (some of them focused on the h and χ Persei clusters).

Generally speaking, no important differences are detected in the distribution of parallaxes, proper motions, and RVs when stars in h and χ Persei or the full sample are considered, which suggests a very extended dynamically interrelated population. However, a few clear outliers in the proper motion and RV distributions are also found. A large fraction of these are O-type stars (almost 50%). The further analysis of their proper motions and RVs indicates that they are runaway stars, probably resulting from the kick of a supernova explosion in a previously bounded binary system.

Finally, we have found that the percentage of secure binaries decreases from the hotter to the cooler Sgs. In particular, this percentage decreases from 15% to 10% when the O star and B Sg samples are compared (or alternatively, from 60% to 15% when we consider the runaway stars as previous binaries), and it practically vanishes in the A/F and K/M Sgs. Further investigation of the potential connection between this result and merging processes that occur during the evolution of massive stars is an interesting direction of future work.

Acknowledgments

Based on observations made with the Nordic Optical Telescope, operated by NOTSA, and the Mercator Telescope, operated by the Flemish Community, both at the Observatorio de El Roque de los Muchachos (La Palma, Spain) of the Instituto de Astrofísica de Canarias. We acknowledge funding from the Spanish Government Ministerio de Ciencia e Innovación through grants PGC-2018-091 3741-B-C211/C22, SEV 2015-0548, and CEX2019-000920-S and from the Canarian Agency for Research, Innovation and Information Society (ACIISI), of the Canary Islands Government, and the European Regional Development Fund (ERDF), under grant with reference ProID2017010115. This research made use of the SIMBAD, operated at Centre de Données astronomiques de Strasbourg, France, and NASA’s Astrophysics Data System. The background images were taken from The STScI Digitized Sky Survey (Copyright link).

References

All Tables

All Figures

Fig. 1. Sky map with the complete sample of stars. Purple, blue, cyan, golden, and red symbols represent the O-, B-, A-, F-, and K- and M-type stars, respectively. This color code is the same in all the plots unless otherwise specified. The central green cross denotes the center of the Per OB1 association taken from Mel’nik & Dambis (2017). The large green circle indicates a 4.5-degree circle around the center. The small green circles show the positions of h and χ Persei. The background image, used for reference, was taken from DSS-red.

In the text

	Fig. 2. Histograms by SpT separated with colors (top) and number of spectra separated by SpT and stacked (bottom).
In the text

	Fig. 3. Top: Gmag against the Gaia error in parallax. Bottom: Gaia error in total proper motion against the error in parallax. Both panels include all stars in our sample except for the seven targets with Gaia RUWE > 1.4 (see Sect. 2.3).
In the text

Fig. 4. Top panel: sky map of all the stars in the sample overplotted over a DSS-red image of the region. Dashed yellow line indicates the galactic plane, and the central green cross marks the center of Per OB1 (as defined in Mel’nik & Dambis 2017). Green circles indicated the location of the h and χ Persei double cluster. Colored vectors indicate the individual proper motion of each star. Middle and bottom panels: parallax and RV of the spectrum with the highest S/N, respectively, for each star in the sample against their position in right ascension. Open circles and square symbols indicate stars that deviate more than 2σ from the mean of the distribution of parallaxes and RVs, respectively (see Sects. 4.1 and 4.4.1). Stars with bad astrometry (see Sect. 2.3) are indicated with a plus, and no proper motion vectors are overplotted.

In the text

Fig. 5. Total proper motions against parallax for the sample of stars except for those labeled “unreliable astrometry” (see Sect. 2.3). The 2σ boundaries of the distribution are shown as a rectangle. Empty colored circles show outliers of the distribution of any of the two quantities, and the associated uncertainties are overplotted. The mean and standard deviation obtained from the stars within the 2σ box are shown in the top right corner.

In the text

	Fig. 6. Same as the top panel of Fig. 4, but this time, the individual proper motion of each star is referred to the mean proper motion of 16 stars in h and χ Persei with good astrometric solution (black arrow in the bottom right corner of the top panel).
In the text

Fig. 7. Color-magnitude diagram (using Gaia photometry) of stars located within 4.5 degrees from the center of the Per OB1 association. Colored stars shows the stars in our sample, gray circles represent the remaining stars from Gaia, green diamonds show 17 blue and red supergiants quoted in the literature for which we lack spectra (see the last part of Table A.1). Two isochrones and a reddening vector are also included for reference purposes. See Sect. 4.3 for explanation.

In the text

	Fig. 8. RV distributions associated with the different SpT groups resulting from the analysis of the best S/N spectra. The orange bin in the second panel from the bottom is HD 12842, the F Sg.
In the text

	Fig. 9. Measured RVs (subtracted from their mean) for a sample of 15 stars (ordered by SpT) for which we have five or more spectra, and whose detected variability in RV is more likely produced by intrinsic variability than by the orbital motion in a binary system.
In the text

Fig. 10. Radial distance from the center of Per OB1 against the mean RV for the stars with more than one spectrum, or RV for the stars with only one spectrum, excluding the SB2 binaries. The stars within 1 degree from IC 1805 were also excluded. The filled gray rectangle shows 2σ of all the RVs, excluding those from stars identified as binaries (SB1 or SB2). The filled colored circles denote stars identified as SB1. The filled colored crosses denote starts outside 2σ from the mean that are not identified as binaries. The colored stars show the remaining stars. The colored error bars present RVPP for stars with multi-epoch data, except for SB1. The open triangles show stars that are outliers in proper motion and were therefore identified as runaways.

In the text