© The Authors 2023

1 Introduction

Of the currently more than 5300 known exoplanets, only 6% have periods of longer than 3 yr, and less than 3% have periods longer than 10 yr¹. The radial velocity (RV) method, one of the leading methods in providing long-period exoplanet detections, was responsible for most of them. Long-term RV surveys using precise spectrographs made this accomplishment possible, including but not limited to: the over 12-yr historical ELODIE program initiated in 1994 (Mayor & Queloz 1995; Naef et al. 2005), the SOPHIE large program started in 2006 (Bouchy et al. 2009), the Anglo-Australian Planet Search launched in 1998 (Tinney et al. 2001), and the HARPS² volume-limited survey, from which the observations in this paper originate. From these different RV surveys, Jupiter analogues and long-period companions were discovered and their properties scrutinised (e.g. Boisse et al. 2012; Moutou et al. 2014; Wittenmyer et al. 2016; Rey et al. 2017; Bryan et al. 2019; Feng et al. 2019; Kiefer et al. 2019; Dalal et al. 2021; Demangeon et al. 2021; Sreenivas et al. 2022). The HARPS volume-limited survey (up to 57.5 pc) started as part of the HARPS Guaranteed Time Observations (GTO) in 2003 (Mayor et al. 2003) and is targeting low-activity, solar-type dwarf stars with spectral types from late F to early M (Lo Curto et al. 2010). The baseline of the survey, which is over 19 yr, allows the detection of long-period signals, including gas giants beyond the ice line. The detection of lighter planets and other long-period companions such as brown dwarfs is also possible, although these targets are not the direct aim of the program.

The detection of long-period planets, and in particular of ‘cold-Jupiters’ (M > 0.3 M_J, a > 1 AU), strongly correlates with the presence of inner super-Earths; more specifically, Zhu & Wu (2018) quote a conditional probability of super-Earths of 90% ± 20%. Rosenthal et al. (2022) are probing a different class of outer planets with semi-major axis a in the range 0.3 to 3 AU and in their recent publication present a value of ; as they are looking at different samples, this is effectively not in contrast.

In any case, the host stars to our detected long-period planets are potential sources of yet-undetected super-Earths. Our survey is not designed to detect these super-Earths, and our precision and sampling are not optimal. Outer giant planets might stabilise planetary systems, as for example shown by the Nice model (Morbidelli 2018) in which the gas giants are believed to shield the inner planets from collisions; Jupiter and Saturn for instance could be responsible for the existence of Earth and the other inner planets in their current orbits.

Hot giant planets tend to orbit metal-rich host stars, a correlation known as the giant planet–metallicity correlation (Fischer & Valenti 2005; Osborn & Bayliss 2020). As the metallicity of the host star is representative of the metallicity content in the protoplanetary disc, this result is an important piece of observational evidence in support of the core-accretion scenario (Pollack et al. 1996). In this formation scenario, a large metal core, which is expected to be more easily formed in a metal-rich disc (i.e. more planetesimals), efficiently accretes gas to form a gas giant. While hot Jupiters are more frequent around metal-rich stars, for longer-period giant planets this is not necessarily the case. Adibekyan et al. (2013) showed that planets (from ∼0.03 M_J to ∼4 M_J) formed in metal-poor systems generally have longer periods than those formed in metal-rich systems. These authors suggest a metal-poor disc may form the giant planets further out and/or formation may start later, which would mean the planets undergo less migration. Short- and long-period giant planets may have the same formation but different migration histories (Osborn & Bayliss 2020). More data are needed to further constrain the formation and evolution of long-period giant planets. With this research, we are probing this still largely unknown population of exoplanets.

The RV method observes the stellar wobble induced by the gravitational pull of the companion orbiting its host star. Stellar activity can mimic this wobble and create false-positive detections. As our program focuses on finding Jupiter-mass long-period companions, it does not require high precision and is not hampered by the short-period (∼minutes) p-mode oscillations. For long-period planets, the greatest difficulties stem from the entanglement of the RV with the stellar magnetic cycle and the RV-induced signal from binary stars (or multiple-star systems). Although the magnetic cycle can inconveniently dominate the RV variation, the signal induced by the magnetic cycle correlates well with the stellar activity indices. The correlation can therefore be used to determine whether or not there is entanglement with the magnetic cycle and then be fitted and removed if the long-term activity index evolves smoothly (Dumusque et al. 2011).

In search of gas giants, we discuss 12 stars in this paper. Each signal is thoroughly examined to ensure the variation is not caused by stellar activity. Section 2 contains information on the observations. In Sect. 3 we present the stellar parameters of the host stars. In Sect. 4, we discuss one magnetic cycle and derive the orbital solutions for each star from the RV observations. In Sect. 5, we present the combination of astrometric measurements from HIPPARCOS and Gaia with RVs in order to constrain the mass and inclination of the heaviest companions. Our results are discussed in Sect. 6.

2 Observations

The observations used in the present study were carried out using the HARPS spectrograph at the ESO 3.6-m telescope at La Silla Observatory, Chile (Mayor et al. 2003). The presented signals are part of a long-term ongoing exoplanet survey that started in 2003 as part of the HARPS GTO program (Mayor 2003) and is still continuing today (Lo Curto 2010). In May 2015, the fibre link of HARPS was upgraded (Lo Curto et al. 2015). With the upgrade, the instrumental profile was changed significantly. As this might induce an offset, we considered the data before and after the fibre upgrade as coming from two independent instruments. The magnitude of the offset varies for different stars and affects the CCF-FWHM as well.

The RVs result from version 3.5 of the HARPS data reduction software (DRS), which cross-correlates each spectrum with a numerical stellar template, matching its spectral type. The pipeline derives several cross-correlation function (CCF) parameters: the RVs, the full width at half maximum (FWHM), the bisector span, and the contrast. Other results are the Mount Wilson S-index, the chromospheric emission ratio log , and the activity indices based on the Na, Ca, and Hα lines.

The number of measurements and the basic characteristics of the measurements are summarised in Table 1. Our program is targeting mainly Jupiter-mass planets and therefore requires only a moderate RV precision of ∼3 ms⁻¹, which corresponds to a signal-to-noise ratio (S/N) of about 40. As this goal is not always obtained, we include observations with an S/N of larger than 25. This limit is used as additional quality control (QC); because the exposure times are defined for an S/N of 40, values below 25 indicate a poor observation (bad seeing, bad weather, etc.). The total number of measurements N_meas excludes the observations with a S/N below 25 and observations that did not pass the DRS QC. The DRS QC checks for the reliability and accuracy of the data obtained by verifying the instrumental stability, data-reduction process, and consistency with known calibration sources (ESO 2019). The public HIRES (Vogt et al. 1994) RV data are added to the analysis of BD-210397, and CORALIE observations (Queloz et al. 2001) are added to the study of HD 61383.

3 Stellar properties

All stars presented in this paper are main sequence stars, like our Sun. Table 2 summarises the parameters of the stars analysed in the present paper. The spectral type, V and B – V values originate from the HIPPARCOS catalogue (ESA 1997), where the magnitudes are dereddened by Gomes da Silva et al. (2021). The parallax with the derived distance is taken from Gaia DR3 (Gaia Collaboration 2023). The spectroscopic analysis of Delgado Mena et al. (2017) provides log ɡ (corrected by Gomes da Silva et al. 2021), T_eff, and [Fe/H] for all stars, except BD-210397 (Yee et al. 2017) and HIP54597 (Sousa et al. 2011), which are too cold for the parameter estimation by Delgado Mena et al. (2017). The parameters M_V, L, and R_⋆ are determined from the above values using the bolometric correction from Flower (1996). The age and mass estimates are obtained by Delgado Mena et al. (2019), apart from HIP54597 and BD-210397, which are calculated via theoretical isochrones (Bressan et al. 2012).

The average FWHM is obtained from the HARPS data from before 2015 and cross-correlated with G2 masks, as the v sin i approximation adapted from Santos et al. (2002) is calibrated on the same data period and spectral type. The FWHM uncertainty is the standard deviation and is therefore an indicator of activity in the RV signal. The average chromosphere emission ratio log results from the Noyes et al. (1984) method. The rotational period is averaged and estimated from the activity-Rossby relations described by Mamajek & Hillenbrand (2008), implementing the convective overturn time from Noyes et al. (1984). These relations are defined for stars with (). Here, we use the same method for relatively quiet stars to get an approximation of the rotational period. Long-term companions have periods that do not coincide with the rotational period of the stars and we therefore do not require a very precise value. We use the standard deviations of log and P_rot as an estimate of their uncertainties.

4 Radial velocity data and orbital solutions

The RV variations are examined using tools provided by the Data and Analysis Center for Exoplanets (DACE)³. We use the data obtained before the fibre upgrade in 2015 as the reference systemic velocity, which we denote as γ₀₃ and introduce an offset constant for the data obtained after the fibre upgrade, which we denote as γ₁₅. To ensure the observed RV variation is caused by a companion, its correlations with the stellar activity indicators and the bisector span are inspected. Delisle et al. (2018) describe the computation of activity cycles as part of a detrend-ing process to remove systematic effects from RV data. A model is fitted to the systematic effects – including activity cycles –and is subtracted from the RV data to improve the sensitivity with which planets are detected. As our main focus is on long-period companions, inspection of the correlations is sufficient. For periods shorter than a few months, the RV variation is also compared to the rotational period of the star and the bisector span to infer whether the signal might originate from starspots or plages. After fitting away the activity signals detected in the peri-odograms, further periodic RV variations in the data are searched for using the false alarm probability (FAP) of the periodogram as an indicator of the statistical significance of the candidate signal. The Keplerian model is computed as described in Delisle et al. (2016). The values of the FAP are calculated analytically following (Baluev 2008). When fitting multiple Keplerian solutions, the stellar activity indicators are again examined before each new fit. An MCMC is used to specify the orbital solution and its uncertainties, the algorithm for which is defined in Díaz et al. (2014, 2016).

A magnetic cycle is discussed in Sect. 4.1, as an example of how to distinguish stellar activity from long-period exoplanets. In Sect. 4.2 we present the uncertain solution of HD 3964, where the magnetic cycle hinders the accurate determination of the RV amplitude. The orbital solutions of the exoplanets without strong entanglement with the star’s magnetic cycle are presented in Sect. 4.3, while those of the brown dwarf candidates and binaries are presented in Sects. 4.4 and 4.5. For all solutions presented, we fixed the HARPS instrumental error to 0.75 m s⁻¹. This value is close to some of the lowest residual RMSs we obtain from the literature; see e.g. Lovis et al. (2006). The true mass and inclination derived from astrometry in Sect. 5 classify some of the brown dwarf companions as stellar binaries; we therefore introduce these as candidates.

4.1 Magnetic cycle of HD 11608

HD 11608 was part of our analysis because the RVs show a prominent peak above 0.1% FAP at -3950 days in the peri-odogram. However, after inspecting the stellar activity correlation values, we find that it is more likely activity induced.

The observations show high (≳0.5) correlations with many activity indicators, such as the S-index, Ca-index, and log (see Fig. 1), which all have peaks in the periodogram above 0.1% FAP at ~3300 days. After detrending for log with a Gaussian low-pass filter (timescale 0.5 yr) the peak reduces to below 10% FAP. When using a longer timescale of 1 yr, the peak remains above 1% FAP, but the residuals correlate with the CCF-Contrast, with a coefficient of −0.21. When subsequently detrending the CCF-Bissector, all correlation coefficients fall below 0.1 and the RV peak stays above 10% FAP, as visible in Fig. 2. The corresponding Keplerian solution has either an unlikely offset between the two datasets (21 m s⁻¹) or, when fixing the offset to 11 m s⁻¹ − which is a more plausible offset for K1/K2V (Lo Curto et al. 2015) – does not converge and finds a highly eccentric orbit with a period of decades. As neither solution is convincing, we conclude that the RV variations of HD 11608 are most likely caused by stellar activity. If there is a companion present, its RV signal is heavily entangled with its magnetic cycle.

4.2 Uncertain orbital solution of HD 3964

HD 3964 is a relatively quiet star (). The periodogram has a peak at 23.5 days, which corresponds to the stellar rotation period (see Table 2). The RV variation shows a significant correlation with Hα, as is visible in the correlation values in Fig. 3. The Hα-index has a peak in the periodogram around 1150 days, which is similar to the detected long-period RV variation at 1086 days. There is a clear magnetic cycle visible in the RV variation of HD 3964, but we cannot exclude a distant companion as well. After detrending for the Hα-index, the RV variation peak in the periodogram remains above 0.1% FAP level at a period of 1086 days, as is visible in Fig. 4, where the first periodogram is before and the second is after detrending the Hα-index. The two other strong peaks (both above 10% FAP before detrending) are the aliases of the 1086 day period. As the periods of the magnetic cycle and possible companion are similar, there is no way to accurately determine the RV amplitude. Therefore, we present the orbital solution of HD 3964 b with caution.

If the magnetic cycle is not intertwined with the companion signal, the detected long-period variation at 1086 days can be attributed to the presence of a distant Jupiter-like planet of 0.58 M_J minimum mass. The period of the RV variation does not correspond to the rotational period of the star. Figure 5 shows the RV variation of HD 3964 with the proposed orbital solution versus time; the residuals (O–C) are included.

The residuals show a correlation with the CCF-Bissector with a correlation coefficient of ∼0.3. However, the CCF-Bissector does not show any peaks above 10% FAP in its periodogram. The average variation in RV signal changes after the 2015 fibre upgrade (6.2 ms⁻¹ before, 9.8 ms⁻¹ after), which is also visible in Fig. 5. The stellar jitter of the model is equal to 3.1ms⁻¹ and the σ(O–C) is 3.49 ms⁻¹.

A drift does not improve the model; therefore, if there is a second companion, it should be at a short period. More extensive RV follow-up observations could provide more insight into the origin of the remaining fluctuations. Complementary RV measurements in the near-infrared (NIR) may help to disentangle stellar activity using the chromaticity of stellar active regions.

Fig. 1 Correlation between log and the RV values of HD 11608.

Fig. 2 Periodograms of the RV data of HD 11608 before (top) and after (bottom) detrending for and the CCF-Bissector with a Gaussian low-pass filter (timescale 1 yr). The horizontal lines are the 10%, 1%, and 0.1% FAP levels.

Fig. 3 Correlation between Hα and the RV values of HD 3964.

Fig. 4 Periodograms of HD 3964 for the initial RV data (up) and the RV data detrended for Hα (down). The horizontal lines are the 10%, 1%, and 0.1% FAP levels.

4.3 Planetary systems

4.3.1 HD94771

The observations of the quiet relatively evolved (1.9 R_⊙) star HD 94771 reveal no strong correlation with stellar activity indicators. However, there is a strong periodic RV variation above FAP 0.1% for 2164 days plus its 1 day and 1 yr alias. The Keplerian solution corresponds to a companion with a minimum mass of 0.53 M_J. There is no strong suggestion of a second companion in the current data: adding a drift does not improve the model and the 2.2 ms⁻¹ stellar jitter is towards the lower limit for a 1.2 M_⊙ star with log(ɡ) ~ 4 cm s⁻² (Luhn et al. 2020). We conclude that neither stellar activity nor the rotational period of HD 94771 induces the observed RV variation. The signal is caused by a giant exoplanet. The low stellar jitter might imply the rotation of the star is seen pole on; if the whole system is misaligned, the mass is much higher. In combination with the eccentric orbit, this suggests HD 94771 b is potentially a brown dwarf. However, it is also possible that the star is in a quiet part of its cycle. As mentioned in Sect. 5, there is no strong astrometric signal visible for HD 94771 b, and so we favour the latter explanation. The observations cover three periods of the orbit of HD 94471b, as visible in Fig. 6.

Fig. 5 Radial velocity observations of HD 3964 versus time and the residuals to the best-fit model after detrending. HARPS data from before the fibre upgrade in 2015 are shown in blue and those from after the upgrade are shown in orange.

4.3.2 HIP54597

HIP54597 is a relatively quiet star. Only one observation shows excessive activity in comparison to the others. Therefore, after the pre-selection of the S/N > 25 and the DRS QC, this observation was excluded from the analysis. Stellar activity indicators S-index, Hα-index, , and P_rot suggest a modulation above 10% FAP level at a slightly shorter period (~2825 days) than the RV variation (3250 days), and have correlation coefficients of approximately −0.25.

However, after detrending for the S-index and the CCF-FWHM, the RV signal stays above FAP 0.1% (see Fig. 7); moreover, the detrending does not significantly affect the peri-odogram, and there are only some minor changes in the model (). The correlation is small enough and the difference in periods large enough to conclude that the observed RV variation is not caused by stellar activity or the rotational period of the star, but by a companion with a minimum mass of 2.01 M_J. The Keplerian solution is shown in Fig. 8. There are no strong indications of an additional companion, a drift does not improve the model, and 2.4 m s⁻¹ stellar jitter is relatively low for a K5V star (Luhn et al. 2020). The correlations suggest that there may be a magnetic cycle at play.

Fig. 6 Radial velocity observations of HD 94771 versus time and the residuals to the best-fit model. Line colours are as in Fig. 5.

Fig. 7 Radial velocity periodogram of HIP54597 after detrending for the S-index and the CCF-FWHM. The horizontal lines are the 10%, 1%, and 0.1% FAP levels.

4.3.3 bd-210397

There are 83 HARPS observations of the late K star BD-210397 that pass the DRS QC. One data point was re-reduced with a K5 mask to match the rest of the observations. The HARPS observations show a negative correlation with the Hα-, Ca-, Na-,and S-index, but their periodicity (respectively ~5800, ~5500, ~5400, and ~5250 days) is not similar to the RV variation (1891 and 6360 days). After detrending for the S-index, the correlation coefficients are reduced to below 0.1 and the two RV variations remain visible in the periodogram (see Fig. 9). We conclude that the mentioned correlations are not the origin of the RV variation. We include 18 public HIRES (Vogt et al. 1994) observations to the analysis for better convergence of the fit.

There are two long-period RV variations visible in the observations, which are attributed to a companion with a minimum mass of 0.7 M_J and a period of 1891 days (BD-210397 b) and another companion (BD-210397 c) of 2.4 M_J minimum mass and a period of 6360 days. When fitting BD-210397 c first, the peak of BD-210397 b reduces from above 0.1% to below 10% FAP. However, the model including only BD-210397 c comes with an unrealistic offset (−14ms⁻¹) between the HARPS data before and after the fibre upgrade Lo Curto et al. (2015) and a stellar jitter of 10ms⁻¹. When including BD-210397 b, the stellar jitter reduces to 8.9 ms⁻¹, the offset is more plausible at 19 ms⁻¹, and the improves from 25 to 17 (without including stellar jitter). The ℓ1 -periodogram (Hara et al. 2017), also finds the presence of the 1891 day signal after fitting the 6360-day signal. The ℓ1-periodogram is a variant of the Lomb-Scargle periodogram, but instead of minimising the ℓ2 norm (sum of squares), it minimises the ℓ1 norm (sum of absolute values). This allows the ℓ1 periodogram to be less sensitive to outliers and non-Gaussian noise.

The model is presented in Table 3, and the solution is visible in Fig. 10. The period of BD-210397 c is not well constrained; more observations are required. As the value of BD-210397 is undefined, we cannot directly conclude as to whether or not the high stellar jitter (8.9ms⁻¹) is caused by activity; although heavily uncertain, the age (6.2 ± 4.7 Gyr) suggests it is not. However, the standard deviation of the FWHM (±0.022 km s⁻¹) is large in comparison to the other targets presented in this paper, which indicates that BD-210397 might be on the more active side. Combining the S -index amplitude with the stellar mass 0.679 M_⊙ and log(ɡ) = 4.67 cm s⁻¹, the stellar jitter is expected to be within the range of 3–9 m s⁻¹ (Luhn et al. 2020). The stellar jitter found is high but within the expected range.

Another explanation for the stellar jitter could be a very short-period companion. The periodogram indeed shows a ∼0.5-day signal below 10% FAP level, with its aliases also present. When adding this signal, the jitter is reduced from ∼8.9 ms⁻¹ to ∼7 ms⁻¹. A short-period companion helps to reduce the jitter, but there are not enough data points to be sure. To determine whether the high stellar jitter is caused by another companion, by sampling or by aliasing effects, follow-up measurements with short time intervals are required.

The HIRES and HARPS O–C values are visible in Fig. 10. Where for HARPS the average (absolute) O–C is 8 m s⁻¹, HIRES on average varies from the model by 10 m s⁻¹ . The HIRES data do agree with the model found but cover a very limited time span.

Fig. 8 Radial velocity observations of HIP54597. Line colours are as in Fig. 5.

Fig. 9 Radial velocity data periodogram of BD-210397 after detrending for the S-index. The horizontal lines are the 10%, 1%, and 0.1% FAP levels.

4.3.4 HD74698

The RVs of the quiet star HD 74698 do not show a high correlation with stellar activity indicators apart from the CCF-FWHM (0.26), which has a period ∼8500 days, and disappears when changing the offset between the two datasets. This is even more evident when binning the data every 60 days. As it is dependent on offset, the CCF-FWHM is not detrended and the offset is fixed to 15 m s⁻¹. This value is expected for a G5V star (Lo Curto et al. 2015), reduces the correlation between the RV and CCF-FWHM to almost zero, and is the best-fit offset found by the binned data. There are two signals above 0.1% FAP level visible in the RV periodogram, P_b = 15.017 days (K_b = 5.8 ms⁻¹) and P_c = 3449 days (K_c = 5.3ms⁻¹). After adding the first Keplerian model, the correlations with the stellar activity indicators remain low, indicating that both periods do not originate from activity. There is a third signal visible in the periodogram at ∼1000 days (K ∼ 6.5ms⁻¹); it does not correspond to any activity indicator. Apart from DACE, we used two independent programs: KIMA (Faria et al. 2018) and the ℓ₁-periodogram (Hara et al. 2017), to verify the model. Both programs also suggest the presence of the 1000 day signal. KIMA finds the three-planet solution to be the most likely, and ℓ₁-periodogram finds all three periods but does not suggest the 1000-day signal as a prominent one. As the 1000 day signal is below 10% FAP level in the DACE periodogram and strongly depends on the offset between the two datasets, it does not appear sufficiently robust to be included in the analysis.

ℓ₁ -periodogram is also used to calculate the periodograms of the stellar activity indicators, none of the stellar activities correspond to P_b, P_c, or the possible 1000 day signal. Figure 11 shows the best-fit model, Fig. 12 the phase-folded solution for HD 74698 b, and Table 4 the corresponding orbital parameters. Neither signal originates from stellar activity or the rotational period. We conclude that HD 74698 b and c are exoplane-tary companions. There are large variations still visible in the residuals, which are potentially caused by a 1000-day period companion. Additional observations can provide more insight into this compound system, for which there is substantial evidence of a third companion. As the minimum masses of the two found signals are in the range from Neptune to Saturn, this star is a very interesting target for follow-up observations.

4.4 Brown dwarf candidates

The mass limits of brown dwarfs, are subject to debate; between ∼13 M_J (limit of deuterium fusion) and 80 M_J (limit of hydrogen fusion) is the generally applied range (Spiegel et al. 2011; Chabrier & Baraffe 2000). However, Sahlmann et al. (2010, 2011b) suggest that planet masses can go up to 25 M_J Here we apply the commonly used limit of 13 M_J in order to differentiate brown dwarfs from massive planets. The eccentricities of brown dwarfs are usually larger than those of exoplanets. As they are faint, they are difficult to detect by direct imaging; here, long-term RV surveys like ours provide a means for their detection. The orbital solutions are presented in Table 5.

Fig. 10 Radial velocity curves for BD-210397, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively (right). The public HIRES data are labelled in yellow in the left panel.

Fig. 11 Radial velocity observations of HD 74698, including two signals: the 3530-day period planet with a minimum mass of 0.38 MJ and a 15.017-day period planet with a minimum mass of 0.08 MJ. HARPS data from before and after the fibre upgrade in 2015 are marked in blue and orange, respectively.

4.4.1 HD62364

The observations of the low-metallicity star HD 62364 reveal a strong high-eccentric RV variation. With a 5138 day period, this signal corresponds to a companion with a minimum mass of 12.7 M_J The minimum mass of HD 62364 b is at the edge of the range of brown dwarfs, but since its mass is higher with the inclination found (see Sect. 5) and high eccentricity is more common among brown dwarfs, we conclude that this companion is probably a brown dwarf. The HARPS spectra cover the entire phase of HD 62364 b, as is visible in Fig. 13. There is no strong correlation with any stellar activity indicator and the rotational period is too short to create the observed fluctuation. There is no drift present and the remaining 3.0ms⁻¹ stellar jitter is within the expected range for a 1.2 M_⋆ star with log(ɡ) ∼ 4 (Luhn et al. 2020).

Fig. 12 Phase-folded RV curves for HD 74698 b after subtracting HD 74698 c. HARPS data before and after the fibre upgrade in 2015 are marked in blue and orange.

4.4.2 HD 56380

The inactive star HD 56380 exhibits a strong RV variation above 0.1% FAP level at 3254 days, corresponding to a 33.2 M_J minimum mass companion, without strong stellar activity indicator correlations. The phase is well covered, as is visible in Fig. 14. We conclude that the RV variation is not caused by activity or the star’s rotational period but by a companion. The 1.2ms⁻¹ stellar jitter of the Keplerian model is low for a 0.8 M_⊙ star with log(ɡ) ∼ 4.5 and a drift does not enhance the fit.

Fig. 13 Radial velocity curves for HD 62364, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively.

Fig. 14 Radial velocity curves for HD 56380, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively.

4.4.3 HD221638

The relatively quiet star HD 221638 shows an RV variation with a period of ∼2 days that correlates with the CCF-FWHM with a coefficient of 0.43 and with the S-index with a coefficient of 0.20. None of the periods come close to the RV variation at 6064 days, which stays above 0.1% FAP level even after detrend-ing the CCF-FWHM. The long-period RV variation corresponds to a companion with a minimum mass of 53 M_J. Adding a linear drift to the orbital parameters improves the fit by , implying the presence of a companion with an even longer period. After detrending the CCF-FWHM and adding a drift, there is no remaining stellar jitter. The phase is fully covered and shown in Fig. 15.

Fig. 15 Radial velocity curves for HD 221638, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively.

4.4.4 HD33473a

The companion of the evolved star⁴ (2.3R_⊙) HD33473A was discovered by Moutou et al. (2011). However, the Keplerian solution was incomplete as the observations covered a small fraction of the period and assumed an additional long-term drift.

The increased number of observations allows us to present a significantly improved orbital solution (). As the phase is not fully covered and the model has difficulties converging when the offset is allowed to vary, we fix the offset between the two datasets to 14 m s⁻¹. This offset corresponds to an approximation for its spectral type G3V, derived from the values mentioned for G2V and G4V in Lo Curto et al. (2015). When excluding one observation with excessive activity (), the RV variation correlates with the CCF-Contrast with a coefficient of 0.44 (no period above 10% FAP) and with the CCF-FWHM with a coefficient of −0.49. The CCF-FWHM shows a long period below 10% FAP, which appears to be increasing outside of the range of the span of the observations.

After detrending the CCF-FWHM, the strong signal of 50 yr remains visible in the periodogram. This period corresponds to a companion with a minimum mass of 38.3 M_J, which is significantly heavier than initially thought (7.2 ± 0.3 M_J). The rotational period of the star (32 days) cannot account for the signal, nor can HD 33473B, the stellar companion at a separation of 10 arcsec (Chanamé & Ramírez 2012). Though a stellar companion justifies the inclusion of a long-term drift in the Keplerian model, we decided not to include one, as a linear drift does not improve the model. The phase is still not entirely covered, as shown in Fig. 16; more observations could further improve the orbital solution. We refer to HD 33473C as the companion to HD 33473A, which was previously designated as HD 33473A b by Moutou et al. (2011).

Fig. 16 Radial velocity curves for HD 33473A, which was originally discovered by Moutou et al. (2011). Blue and orange colours indicate HARPS data before and after the fibre upgrade in 2015, the light blue line is the orbital solution as found by Moutou et al. (2011), and the black line is the orbital solution as proposed in this paper.

Fig. 17 Averaged CCF residuals of HD 11938 for different radial velocity shifts. The deepest peak (orange curve) occurs for a shift corresponding to a mass ratio of 0.44.

4.5 Binaries

Long-term RV monitoring allows the detection of binaries not yet resolved by direct imaging because they are too faint and their separation is too small. However, as their periods are on the order of decades, the phases of the binaries are not always fully covered. To confirm the origin of the RV variation, we searched for a second component in the CCFs, and resolve binaries as SB2; see Bouchy et al. (2016) for an in-depth explanation of the applied method. Briefly, the CCFs are recomputed with an M mask to optimise the detectability of the low-mass companion. The CCFs are shifted to the systemic velocity V0 and their average is subtracted to remove the first component. The expected radial velocity of the second component is defined by the radial velocity of the main component V1 , the systemic velocity V₀, and the mass ratio q = M_pl/M_⋆. The residual CCFs are shifted to this expected value and again averaged. Consequently, the secondary is expected to be at V₀. We inspected various combinations of V₀ and q; by searching for the deepest peak at the varying location of V₀, we found which combination best matches the CCFs. The range of V₀ is defined by what is feasible according to the RV model. In this section, we present the finding of one stellar binary (HD 11938) confirmed using this method, and constrain its inclination i. The other presented stellar binary (HD 61383) cannot be confirmed this way due to the blending of the primary and secondary peaks in the CCF.

4.5.1 HD 11938

There are 25 observations of the active star HD 11938 that pass the DRS quality control and have S/N > 25. One measurement was reduced with a G2 mask instead of a K5 mask and was therefore reprocessed. As there are only 25 observations that do not fully cover the phase of HD 11938B we fixed the offset between the HARPS data before and after the fibre upgrade to 15 m s⁻¹ as suggested by Lo Curto et al. (2015) for K4/5V spectral types. There is a clear signal in the periodogram of 35 yr. The stellar binary is also visible in the CCFs, where the best solution is at a V₀ of 40.7 ± 0.2 km s⁻¹ (= γ₀₃). The mass ratio q in that case is 0.44 ± 0.03, corresponding to a secondary mass of 0.33 ± 0.03 M_⊙. Figure 17 shows the average CCF residual for this V₀ and q combination in red.

HD 11938 is an active star () with correlations with all stellar activity indicators, but none have periodical signals above 10% FAP level in the periodogram. This strong RV variation signal, with an amplitude of 2.41 km s⁻¹, is therefore not caused by activity. We fix γ₀₃ to the found systemic velocity and its error margins V₀ = 40.7 ± 0.2 km s⁻¹ in three iterations (i.e. 40.5, 40.7, and 40.9 km s⁻¹). The model presented in Table 6 corresponds to 40.7 km s⁻¹, and the presented errors correspond to the models found for 40.5 and 40.9km s⁻¹. This is done because the MCMC model drifts away from the found systemic velocity when using a prior. We note that the period is heavily dependent on the systemic velocity. The solution is visible in Fig. 18. The found minimum mass msini = 0.21 M_⊙ corresponds to an inclination of 39.33° ± 0.07°. The model does not include detrending or a drift, as both do not improve the model. The relatively large 9 ms⁻¹ stellar jitter agrees with an active star (without detrending). We conclude that the rotation period cannot account for the signal nor can the stellar activity. The RV variation belongs to a stellar binary companion, separated by 262 mas from its host star. The Gaia archive does not show a secondary within this distance range.

4.5.2 HD61383

There are 66 HARPS observations of the metal-poor quiet star HD 61383 that pass the DRS QC. The average photon noise is 3.1 m s⁻¹. To increase phase coverage, the analysis also includes 20 observations from CORALIE (Queloz et al. 2001) with an average photon noise of 6.3 ms⁻¹. The CORALIE instrumental error is fixed at 5.0 m s⁻¹. The phase of the data is not fully covered; this would require at least 20 more years of observations, as is visible in Fig. 19. Due to the incomplete coverage of the phase, the model has difficulties converging. We fix the offset between the HARPS data to 14 m s⁻¹, as HD 61383 is a G3V star. There is a strong peak visible in the periodogram at 39.1 yr corresponding to a minimum mass of 0.18 M_⊙.

The HARPS observations weakly correlate with almost all stellar activity indicators and though there are some with long-term trends, the amplitude of the RV variation is too large to be caused by stellar activity. The binary is barely visible in the CCFs as a secondary peak, but the amplitude is small in comparison to the noise, as it is in most CCFs blended with the primary peak. With the currently available data, we present the orbit, as shown in Table 6 and Fig. 19.

Fig. 18 Radial velocity curves for HD 11938, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015.

5 True mass and inclination from astrometry

Astrometry is a promising technique for measuring the inclinations of our systems, and the single measurement precision of Gaia in particular should allow such measurements (10 µas, based on Gaia EDR3; Gaia Collaboration 2021). The S/N reached by Gaia can be estimated by combining the astrometric signal (approximately for M_pl ≪ M_⋆) with the single measurement precision and the number of Gaia measurements from the Gaia Observation Forecast Tool⁵. All presented targets have an expected S/N of above 6.2, the minimum value needed to retrieve an astrometric orbit from combined astrom-etry and radial velocity data (Sahlmann et al. 2011a). However, it is not certain that Gaia covers the full orbit, and as the presented periods are on the order of several years, depending on the companion mass, the signal might get included in the proper motion. The Gaia excess noise values and their significance are an indication of the astrometric signal. For all presented targets, the excess noise significance is larger than 2, meaning there is a significant astrometric signal seen by Gaia (Lindegren et al. 2012).

As the individual Gaia observations have not yet been released, we are limited to the currently available proper motion and RA/Dec positions. To derive the true mass and inclination, we use the Python package orvara (Brandt et al. 2021), which can fit Keplerian orbits by utilising a combination of radial velocity, relative astrometry, and available absolute astrometry data. Here we combine the radial velocity observations presented in this paper and the absolute astrometry data that comes from the HIPPARCOS-Gaia Catalog of Accelerations (HGCA; Brandt 2021). In addition, we use an extension of orvara that allows priors on orbital periods and semi-major axes⁶, and fix the offset between the radial velocity data to the offsets presented in Sect. 4, as orvara otherwise finds values outside of the expected offset range. See Appendix A for the resulting corner plots, and Appendix B for the proper motion plots.

For most of the exoplanets, the amplitude and period are too low to properly fit the HIPPARCOS-Gaia proper motions. Apart from HIP54597 b, HD 74698 c, and BD-210397 c, the mass and inclination of the exoplanets cannot be constrained. For HIP54597 b, the constraint is weak and varies between an inclination of 120° and 60°. This might increase the mass of the planet by 10%–30% but not significantly more. The constraint of BD-210397 c is slightly better, with peaks at 50° and 130°, favouring the former, and a mass of 2.7 M_J. The inclination of HD 74698 c is 90° ± 33°; though this is not well-constrained, it shows the minimum mass is likely equal to its true mass. For HIP54597 b, HD 74698 c, and BD-210397 c, the inclination is not well-constrained. Their orvara results are therefore not included in Table 3. However, we do conclude that their companions have masses in the planetary range. The convergence is better for our more massive companions, and with the found models, three of the four brown dwarf candidates are found to be stellar binaries (HD 56380B, HD 221638B, and HD 33473C). HD 61383 has an edge-on orbit (i = 86.2°), and its true mass is very close to the minimum mass. For HD 11938, orvara finds an inclination of 64.0°, corresponding to a mass of 0.386 M_⊙, which is slightly higher than the mass found by the secondary peak in the CCF. Orvara does not find the same RV orbital parameters as the simple RV analysis, which causes discrepancies between all presented minimum masses and true masses. For the vast majority of our stars, the difference is negligible. For HD 11938, on the other hand, orvara finds a different systemic velocity. This stems from the fact that the orbit is incomplete and in the RV model the systemic velocity is fixed to the value found by the second component in the CCF. When the systemic velocity is let free, the period can increase up to 41 yr and the minimum mass increases to 0.30 M_⊙, which corresponds better to the model found by orvara.

The true masses found are in agreement with the Gaia renor-malised unit weight error (RUWE) values, which are expected to be 1 for well-behaved single-star solutions. Values significantly greater than 1, with a limit generally placed at 1.4, either indicate the star is not a single star or is otherwise problematic. All presented stars have RUWE values below 1.4, except HD 56380 (1.8348), which corresponds to one of our most massive companions (0.36 M_⊙). The other massive companion HD 11938B does not present a large RUWE, probably due to its very long orbital period (35 yr).

Considering that HD 56380B is our second most massive companion, we decided to reinspect the CCF following the method described in Sect. 4.5 and identified a second component compatible with a massive companion. However, considering that both components are always blended (within 2 km s⁻¹ ), the method does not allow us to derive a precise mass ratio.

Fig. 19 Radial velocity curves for HD 61383, with purple indicating CORALIE data, and blue and orange showing HARPS data before and after the fibre upgrade in 2015.

Fig. 20 Overview of the discovered companions, including the uncertain planetary signal HD 3964 b (white star). Star markers indicate companions discussed in this paper, and colours differentiate between short-period (black) and long-period giants (green), brown dwarfs (orange), and stellar binaries (yellow). Arrows indicate the true masses as found by orvara. Round markers are previously detected (using the RV method) exoplanets and brown dwarfs with well-constrained (>3σ) masses and periods from the exoplanet.eu archive. The grey shaded regions indicate the RV precision detection limits 50 cm s−1 and 5 m s−1 for 1 M⊙. The giant planets (0.3 MJ < M < 13 MJ) are highlighted where the difference between the short-period (blue-shaded) and long-period (green-shaded) objects is defined at 100 days. The brown dwarf region is shaded orange.

6 Discussion and conclusions

In this paper, we discuss the RV variations of 12 stars. We present the discovery of six exoplanets, one uncertain exoplanet candidate, one brown dwarf, four stellar binaries and the improved orbital solution of the stellar binary HD 33473C discovered earlier by Moutou et al. (2011). An overview of the detections is presented in Fig. 20.

The biggest difficulty for long-period companions arises from the RV signal entanglement with the star’s magnetic cycle. Therefore, for all proposed companions, we discuss the possibility of the RV variation originating from activity. For comparison, we also present an activity-induced signal (HD 11608) and a possible exoplanet signal with strong entanglement with its magnetic cycle (HD 3964). These two stars show the strength of investigating the correlation with stellar activity indicators. Detrending reduces the activity-induced component but makes it difficult to constrain the orbital period. As HD 3964 has an RV period (1086 days) similar to the Hα-index period (1150), there is no way to accurately determine the RV signature of the possible companion.

HD 56380 and HD 61383 are relatively old (>12 Gyr). In principle, this is not inconsistent with their average lifetime, but it is uncommon to get so close to the upper age limit for these spectral types. A possible explanation could be that the light of the secondary was blended in the observations, causing the age of the primary star to be overestimated. The ages originate from Delgado Mena et al. (2019), who use isochrones based on T_eff, [Fe/H], and V magnitudes. The latter is sensitive to blending. As HD 56380 and HD 61383 are binaries, and given that we find that the secondary peak in the RV CCF is blended with the primary peak for both sources, the age could indeed be influenced by the secondary.

After determining the orbital solutions, overlap with transit data was examined for all companions as expected by the combination of the period and time of periastron. There are no CHEOPS light curves for the presented stars, and the TESS light curves have no superposition with the potential transits. The only companion to have superposition with the TESS sectors is HD 74698 b, which does not show a transit or periodicity of 15 days: the orbit is not well-aligned with our line of sight. As the potential transits have uncertainties on the order of months, and for some of our targets even on the order of years, we also used a more general approach and looked at the Lomb Scargle periodograms of the TESS light curves. Again, no periodicities were found, meaning there are no observed transits for our targets, but also that there are no activity-induced signals with the same periods as our companions.

Long-period companions provide interesting systems for follow-up observations. HD 3964 has a remaining RV signal of 4 m s⁻¹ , which is likely caused by activity given that the CCF-Bissector is correlated to the residuals, but this could be better constrained by follow-up observations. HD 94771 does not show any strong signs of multiple companions and has a relatively high eccentricity, making the stability of a potential multi-planet system less likely. For HIP54597, there is no significant evidence of additional planets, but there is for both BD-210397 and HD 74698. As BD-210397 includes a high jitter and HD 74698 shows a 1000 day period in the periodogram, both targets are very interesting targets for follow-up observations. Apart from BD-210397, none of the exoplanet’s Keplerian models show a high jitter. For BD-210397, the activity is unknown, but according to the S -index, stellar mass, and surface gravity, the jitter falls within the expected range (Luhn et al. 2020). If the jitter is caused by a companion, it is a very short period signal (∼0.5 days); this can be verified by observations at short time intervals.

It is common for brown dwarfs to have higher eccentricities than exoplanets. According to the NASA exoplanet archive, 94% of the confirmed exoplanets⁷ have eccentricities e ≤ 0.5, versus 69% of the brown dwarfs. For e ≤ 0.3, this value decreases to 84% for the exoplanets and to 56% for the brown dwarfs. Brown dwarfs are less likely to have eccentricities below 0.3 than exo-planets. Our results are in agreement. Only one of the proposed exoplanets (HD 94771 b) has an eccentricity of larger than 0.3, while HD 62364 has a high eccentricity (0.607). When comparing long-period giant planets to short-period giant planets, placing the limit at 100 days (see also Fig. 20), a similar trend occurs: 90% of the short-period giant planets have eccentricities e ≤ 0.3, versus 72% of the long-period giant planets.

According to the planet-metallicity correlation, host stars with higher metallicities are more likely to host a giant planet (Fischer & Valenti 2005), which is particularly true for hot Jupiters (Osborn & Bayliss 2020). The short-period giant planets (<100 days) have parent stars with an average metallicity of 0.11 dex. This value is 0.04 dex on average for stars hosting long-period giant planets, a similar value to our small sample of presented stars hosting (long-period) giant planets (0.05 dex). This is in agreement with the findings of Adibekyan et al. (2013): giant planets orbiting metal-poor stars have longer periods than those in metal-rich systems. The findings of (Santos et al. 2017) highlight the possibility that stars with massive giant planets (M ≳ 4 M_J) do not follow the same metallicity trend as stars with lower mass giant planets. The stars hosting massive giant planets are on average more metal-poor. Though all planets discussed in the present paper have minimal masses of below 4 M_J, this may indicate that HIP54597 b ([Fe/H] = −0.22, m₂ sin i = 2.01 M_J) falls in the massive giant category; however, the results from orvara show this is unlikely.

The discovery of very long-period companions requires years of observations. With the present research, we probe a region of a largely unknown population. The current baseline of 19 yr allows us to detect Saturn-mass companions with periods of 10 yr and brown dwarfs and stellar binaries with periods of up to 50 yr. Though the RV observations do not cover the full phase of the found stellar binaries, we are able to constrain the orbital parameters and confirm the origin of HD 11938B using the secondary peak in the combined CCFs, and of HD 62364 b, HD 5680B, HD221638B, HD33473C, HD 11938B, and HD61383B via absolute astrometry. As by-products of this survey, the brown dwarfs andbinaries demonstrate the effectiveness ofa long baseline. The impact of continuing with observations is made clear by the difference between the orbital solution of HD 33473C found by Moutou et al. (2011) and the solution for this object presented here. Incomplete orbits are prone to errors. Following indications by Rosenthal et al. (2022) and Zhu & Wu (2018), the long-period Jupiters detected in this work are prime targets for the detection of low-mass inner planets.

Acknowledgements

The authors thank the ESO staff at La Silla for their diligent and competent help during the observations and for the effort to maintain the instrument operating and stable for so many years. We thank Emanuela Pompei for providing helpful comments. The HARPS spectrograph was built by the contributions of the Swiss FNRS, the Geneva University, the French Institut National des Sciences de l’Univers (INSU) and ESO. This research made use of the Simbad database, operated at the CDS, Strasbourg, France. This work has been carried out within the framework of the NCCR PlanetS supported by the Swiss National Science Foundation. NCS acknowledges the support from the European Research Council through grant agreement 101052347 (FIERCE). This work was supported by FCT Fundação para a Ciência e a Tecnologia through national funds and by FEDER through COMPETE2020 Programa Operacional Competitividade e Internacionalização by these grants: UIDB/04434/2020; UIDP/04434/2020. This publication makes use of The Data & Analysis Center for Exoplanets (DACE), which is a facility based at the University of Geneva (CH) dedicated to extrasolar planets data visualisation, exchange and analysis. DACE is a platform of the Swiss National Centre of Competence in Research (NCCR) PlanetS, federating the Swiss expertise in Exoplanet research. The DACE platform is available at https://dace.unige.ch.

Appendix A Orvara corner plots

Fig. A.1 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 3964 b. There are no strong constraints on the secondary mass or inclination.

Fig. A.2 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 94771 b. There are again no strong constraints on the secondary mass or inclination.

Fig. A3 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HIP54597 b. Though the mass and inclination are not well-constrained, the inclination varies between 120° and 60°. The mass of the planet might increase by ~10%–30% but not significantly more.

Fig. A.4 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of BD-210397 c. The inclination peaks at 50° and 130°, favouring 50°. The found mass is 2.7 MJ. This corner plot is for the long-period (7063 d) planet; the fit also includes BD-210397 b.

Fig. A.5 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 74698 с. The inclination is not well-constrained at 90° ± 33°; the true mass is equal to the found minimal mass. The fit also includes HD 74698 b.

Fig. A.6 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 62364 b. The mass and inclination are well-constrained.

Fig. A.7 Comer plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 5638OB. The mass and inclination are well-constrained.

Fig. A.8 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 221638B. The mass and inclination are well-constrained.

Fig. A.9 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 33473C. The mass and inclination are well-constrained.

Fig. A.10 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 11938B. The mass and inclination are well-constrained.

Fig. A.11 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 61383B. The mass and inclination are well-constrained.

Appendix B Orvara proper motion plots

Fig. B.1 Proper motion measurements from HIPPARCOS and Gaia (grey markers) and fitted proper motions in right ascension and declination. The black line is the orbit with the highest likelihood and the 50 coloured lines are randomly drawn orbits from the posterior distribution, where the colour corresponds to the mass of the companion. The bottom figure shows the difference between observed and calculated astrometry.

References

All Tables

All Figures

	Fig. 1 Correlation between log and the RV values of HD 11608.
In the text

	Fig. 2 Periodograms of the RV data of HD 11608 before (top) and after (bottom) detrending for and the CCF-Bissector with a Gaussian low-pass filter (timescale 1 yr). The horizontal lines are the 10%, 1%, and 0.1% FAP levels.
In the text

	Fig. 3 Correlation between Hα and the RV values of HD 3964.
In the text

	Fig. 4 Periodograms of HD 3964 for the initial RV data (up) and the RV data detrended for Hα (down). The horizontal lines are the 10%, 1%, and 0.1% FAP levels.
In the text

	Fig. 5 Radial velocity observations of HD 3964 versus time and the residuals to the best-fit model after detrending. HARPS data from before the fibre upgrade in 2015 are shown in blue and those from after the upgrade are shown in orange.
In the text

	Fig. 6 Radial velocity observations of HD 94771 versus time and the residuals to the best-fit model. Line colours are as in Fig. 5.
In the text

	Fig. 7 Radial velocity periodogram of HIP54597 after detrending for the S-index and the CCF-FWHM. The horizontal lines are the 10%, 1%, and 0.1% FAP levels.
In the text

	Fig. 8 Radial velocity observations of HIP54597. Line colours are as in Fig. 5.
In the text

	Fig. 9 Radial velocity data periodogram of BD-210397 after detrending for the S-index. The horizontal lines are the 10%, 1%, and 0.1% FAP levels.
In the text

	Fig. 10 Radial velocity curves for BD-210397, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively (right). The public HIRES data are labelled in yellow in the left panel.
In the text

	Fig. 11 Radial velocity observations of HD 74698, including two signals: the 3530-day period planet with a minimum mass of 0.38 MJ and a 15.017-day period planet with a minimum mass of 0.08 MJ. HARPS data from before and after the fibre upgrade in 2015 are marked in blue and orange, respectively.
In the text

	Fig. 12 Phase-folded RV curves for HD 74698 b after subtracting HD 74698 c. HARPS data before and after the fibre upgrade in 2015 are marked in blue and orange.
In the text

	Fig. 13 Radial velocity curves for HD 62364, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively.
In the text

	Fig. 14 Radial velocity curves for HD 56380, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively.
In the text

	Fig. 15 Radial velocity curves for HD 221638, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015, respectively.
In the text

	Fig. 16 Radial velocity curves for HD 33473A, which was originally discovered by Moutou et al. (2011). Blue and orange colours indicate HARPS data before and after the fibre upgrade in 2015, the light blue line is the orbital solution as found by Moutou et al. (2011), and the black line is the orbital solution as proposed in this paper.
In the text

	Fig. 17 Averaged CCF residuals of HD 11938 for different radial velocity shifts. The deepest peak (orange curve) occurs for a shift corresponding to a mass ratio of 0.44.
In the text

	Fig. 18 Radial velocity curves for HD 11938, with blue and orange indicating HARPS data before and after the fibre upgrade in 2015.
In the text

	Fig. 19 Radial velocity curves for HD 61383, with purple indicating CORALIE data, and blue and orange showing HARPS data before and after the fibre upgrade in 2015.
In the text

Fig. 20 Overview of the discovered companions, including the uncertain planetary signal HD 3964 b (white star). Star markers indicate companions discussed in this paper, and colours differentiate between short-period (black) and long-period giants (green), brown dwarfs (orange), and stellar binaries (yellow). Arrows indicate the true masses as found by orvara. Round markers are previously detected (using the RV method) exoplanets and brown dwarfs with well-constrained (>3σ) masses and periods from the exoplanet.eu archive. The grey shaded regions indicate the RV precision detection limits 50 cm s−1 and 5 m s−1 for 1 M⊙. The giant planets (0.3 MJ < M < 13 MJ) are highlighted where the difference between the short-period (blue-shaded) and long-period (green-shaded) objects is defined at 100 days. The brown dwarf region is shaded orange.

In the text

	Fig. A.1 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 3964 b. There are no strong constraints on the secondary mass or inclination.
In the text

	Fig. A.2 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 94771 b. There are again no strong constraints on the secondary mass or inclination.
In the text

	Fig. A3 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HIP54597 b. Though the mass and inclination are not well-constrained, the inclination varies between 120° and 60°. The mass of the planet might increase by ~10%–30% but not significantly more.
In the text

	Fig. A.4 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of BD-210397 c. The inclination peaks at 50° and 130°, favouring 50°. The found mass is 2.7 MJ. This corner plot is for the long-period (7063 d) planet; the fit also includes BD-210397 b.
In the text

	Fig. A.5 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 74698 с. The inclination is not well-constrained at 90° ± 33°; the true mass is equal to the found minimal mass. The fit also includes HD 74698 b.
In the text

	Fig. A.6 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 62364 b. The mass and inclination are well-constrained.
In the text

	Fig. A.7 Comer plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 5638OB. The mass and inclination are well-constrained.
In the text

	Fig. A.8 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 221638B. The mass and inclination are well-constrained.
In the text

	Fig. A.9 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 33473C. The mass and inclination are well-constrained.
In the text

	Fig. A.10 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 11938B. The mass and inclination are well-constrained.
In the text

	Fig. A.11 Corner plot with the orvara-derived primary star mass, companion mass, semimajor axis, eccentricity, and inclination of HD 61383B. The mass and inclination are well-constrained.
In the text

Fig. B.1 Proper motion measurements from HIPPARCOS and Gaia (grey markers) and fitted proper motions in right ascension and declination. The black line is the orbit with the highest likelihood and the 50 coloured lines are randomly drawn orbits from the posterior distribution, where the colour corresponds to the mass of the companion. The bottom figure shows the difference between observed and calculated astrometry.

In the text