The HARPS search for southern extra-solar planets XLVII. Five Jupiter-mass planets in long-period orbits, one highly irradiated Neptune, one brown dwarf, and ﬁve stellar binaries

Context. The long-term ongoing HARPS radial velocity survey of extra-solar planets initiated in 2003 provides a unique data set with a 19 year baseline that allows the detection of long-period exoplanets, brown dwarfs


Introduction
Of the currently more than 5300 known exoplanets, only 6% have periods of longer than 3 years, and less than 3% have peri-0100.C-0097, 0101.C-0379, 0102.C-0558, 0103.C-0432 106.21R4.001, 108.222V. 001, 183.C-0972, 192.C-0852 and 196.C-1006. Article number, page 1 of 22 arXiv:2306.04420v1 [astro-ph.EP] 7 Jun 2023 A&A proofs: manuscript no. 46203 ods longer than 10 years 1 . 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 year 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 2 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, Demangeon et al. 2021, Sreenivas et al. 2022. The HARPS volume-limited survey (up to 57.5pc) 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 years, 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.3M J , a > 1AU), 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 32 +24 −16 %; 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.03M J to ∼ 4M 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 1 See The Extrasolar Planets Encyclopaedia, http://exoplanet.eu 2 High Accuracy Radial velocity Planet Searcher. 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 longperiod 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 Section 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.

Observations
The observations used in the present study were carried out using the HARPS spectrograph at the ESO 3.6 meter 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 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 R HK , 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 m s −1 , 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 HD61383.

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 et al. 2022 (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 R HK 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 (log R HK > −5.0). 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 R HK and P rot as an estimate of their uncertainties.

Radial velocity data and orbital solutions
The RV variations are examined using tools provided by the Data and Analysis Center for Exoplanets (DACE) 3 . We use the data obtained before the fibre upgrade in 2015 as the reference systemic velocity, which we denote as γ 03 and introduce an offset constant for the data obtained after the fibre upgrade, which we denote as γ 15 . 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 detrending 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 periodograms, 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. (2014Díaz et al. ( , 2016. A magnetic cycle is discussed in subsection 4.1, as an example of how to distinguish stellar activity from long-period exoplanets. In subsection 4.2 we present the uncertain solution of HD3964, 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 subsection 4.3, while those of the brown dwarf candidates and binaries are presented in subsections 4.4 and 4.5. For all solutions presented, we fixed the HARPS instrumental error to 0.75 m s −1 . 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 section 5 classify some of the brown dwarf companions as stellar binaries; we therefore introduce these as candidates.

Magnetic cycle of HD11608
HD11608 was part of our analysis because the RVs show a prominent peak above 0.1% FAP at ∼ 3950 days in the periodogram. 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 R HK (see Fig. 1), which all have peaks in the periodogram above 0.1% FAP at ∼ 3300 days. After detrending for log R HK 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 Figure 2. The corresponding Keplerian solution has either an unlikely offset between the two datasets (21 m s −1 ) or, when fixing the offset to 11 m s −1 -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 HD11608 are most likely caused by stellar activity. If there is a companion present, its RV signal is heavily entangled with its magnetic cycle.

Uncertain orbital solution of HD3964
HD3964 is a relatively quiet star (log R HK = −4.87). The log R HK 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 Figure 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 HD3964, but we cannot ex-  clude 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 Figure 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 HD3964 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 HD3964 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 m s −1 before, 9.8 m s −1 after), which is 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.

HD94771
The observations of the quiet relatively evolved (1.9 R ) star HD94771 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 year 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 m s −1 stellar jitter is towards the lower limit for a 1.2 M star with log(g) ∼4 cm s −2 (Luhn et al. 2020).
We conclude that neither stellar activity nor the rotational period of HD94771 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 HD94771 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 HD94771 b, and so we favour the latter explanation. The observations cover three periods of the orbit of HD94471b, as visible in Figure 6.

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 Sindex, Hα-index, log R HK , 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 periodogram, and there are only some minor changes in the model (∆χ 2 red = 0.01). 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 Figure 8. There are no strong indications of an additional companion, a drift does not improve the model, and 2.4 m s −1 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.

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.4M 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 (-14 m s −1 ) between the HARPS data before and after the fibre upgrade Lo Curto et al. (2015) and a stellar jitter of 10 m s −1 . When including BD-210397 b, the stellar jitter reduces to 8.9 m s −1 , the offset is more plausible at 19 m s −1 , and the χ 2 red 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 1periodogram 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 Figure 10. The period of BD-210397 c is not well constrained; more observations are required. As the log R HK value of BD-210397 is undefined, we cannot directly conclude as to whether or not the high stellar jitter (8.9 m s −1 ) 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 −1 ) 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(g) = 4.67 cm s −1 , the stellar jitter is expected to be within the range of 3-9 m s −1 (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 m s −1 to ∼7 m s −1 . 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 Figure 10. Where for HARPS the average (absolute) O-C is 8 m s −1 , HIRES on average varies from the model by 10 m s −1 . The HIRES data do agree with the model found but cover a very limited time span.

HD74698
The RVs of the quiet star HD74698 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 −1 . 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 m s −1 ) and P c = 3449 days (K c = 5.3 m s −1 ). 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.5 m s −1 ); it does not correspond to any activity indicator. Apart from DACE, we used two independent programs: KIMA (Faria et al. 2018) and the 1 -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 1 -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. 1 -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, Figure 12 the phase-folded solution for HD74698 b, and Table 3 the corresponding orbital parameters. Neither signal originates from stellar activity or the rotational period. We conclude that HD74698 b and c are exoplanetary 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.

Brown dwarf candidates
The mass limits of brown dwarfs, are subject to debate; between ∼ 13M 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. (2010Sahlmann et al. ( , 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, longterm RV surveys like ours provide a means for their detection.

HD62364
The observations of the low-metallicity star HD62364 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 HD62364 b is at the edge of the range of brown dwarfs, but since its mass is higher with the inclination found (see section 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 HD62364 b, as is visible in Figure 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.0 m s −1 stellar jitter is within the expected range for a 1.2 M star with log(g) ∼4 (Luhn et al. 2020).

HD56380
The inactive star HD56380 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 Figure  14. We conclude that the RV variation is not caused by activity or the star's rotational period but by a companion. The 1.2 m s −1 stellar jitter of the Keplerian model is low for a 0.8M star with log(g) ∼4.5 and a drift does not enhance the fit.

HD221638
The relatively quiet star HD221638 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 detrending 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 ∆χ 2 red = 0.08, im-Article number, page 9 of 22   plying 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 Figure 15.

HD33473A
The companion of the evolved star 4 (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 (∆χ 2 red = 2.25). 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 −1 . 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 (∆ log R HK ∼ −0.4), 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 years 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 HD33473B, the stellar companion at a separation of 10 arcsecs (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 Figure 16; more observations could further improve the orbital solution. We refer to HD33473C as the companion to HD33473A, which was previously designated as HD33473A b by (Moutou et al. 2011).

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  Consequently, the secondary is expected to be at V 0 . We inspected various combinations of V 0 and q; by searching for the deepest peak at the varying location of V 0 , we found which combination best matches the CCFs. The range of V 0 is defined by what is feasible according to the RV model. In this section, we present the finding of one stellar binary (HD11938) confirmed using this method, and constrain its inclination i. The other presented stellar binary (HD61383) cannot be confirmed this way due to the blending of the primary and secondary peaks in the CCF.

HD11938
There are 25 observations of the active star HD11938 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 HD11938B we fixed the offset between the HARPS data before and after the fibre upgrade to 15 m s −1 as suggested by Lo Curto et al. (2015) for K4/5V spectral types. There is a clear signal in the periodogram of 35 years. The stellar binary is also visible in the CCFs, where the best solution is at a V 0 of 40.7 ± 0.2 km s −1 (= γ 03 ). The mass ratio q in that case is 0.44 ± 0.03, corresponding to a secondary mass of 0.33 ± 0.03 Fig. 17: Averaged CCF residuals of HD11938 for different radial velocity shifts. The deepest peak (orange curve) occurs for a shift corresponding to a mass ratio of 0.44.
M . Figure 17 shows the average CCF residual for this V 0 and q combination in red. HD11938 is an active star (logR HK ∼ −4.5) 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 −1 , is therefore not caused by activity. We fix γ 03 to the found systemic velocity and its error margins V 0 = 40.7 ± 0.2 km s −1 in three iterations (i.e. 40.5, 40.7, and 40.9 km s −1 ). The model presented in Table 6 corresponds to 40.7 km s −1 , and the presented errors correspond to the models found for 40.5 and 40.9 km s −1 . 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 Figure 18. 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.

HD61383
There are 66 HARPS observations of the metal-poor quiet star HD61383 that pass the DRS QC. The average photon noise is 3.1 m s −1 . To increase phase coverage, the analysis also includes 20 observations from CORALIE (Queloz et al. 2001) with an average photon noise of 6.3 m s −1 . The CORALIE instrumental error is fixed at 5.0 m s −1 . The phase of the data is not fully covered; this would require at least 20 more years of observations, as is visible in Figure 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 −1 , as HD61383 is a G3V star. There is a strong peak visible in the periodogram at 39.1 years 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 longterm 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 Figure 19. ments from the Gaia Observation Forecast Tool 5 . All presented targets have an expected S/N of above 6.2, the minimum value needed to retrieve an astrometric orbit from combined astrometry 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).

True mass and inclination from astrometry
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 , 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, 6 and fix the offset between the radial velocity data to the offsets presented in section 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, HD74698 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 HD74698 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, HD74698 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 (HD56380B, HD221638B, and HD33473C). HD61383 has an edge-on orbit (i = 86.2 • ), and its true mass is very close to the minimum mass. For HD11938, 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 HD11938, 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 years 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 renormalised 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 HD56380 (1.8348), which corresponds to one of our most massive companions (0.36 M ). The other massive companion HD11938B does not present a large RUWE, probably due to its very long orbital period (35 years).
Considering that HD56380B is our second most massive companion, we decided to reinspect the CCF following the method described in section 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.

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 HD33473C discovered earlier by Moutou et al. (2011). An overview of the detections is presented in Figure 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 (HD11608) and a possible exoplanet signal with strong entanglement with its magnetic cycle (HD3964). 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 HD3964 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.
HD56380 and HD61383 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 HD56380 and HD61383 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 HD74698 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. HD3964 has a remaining RV signal of Fig. 20: Overview of the discovered companions, including the uncertain planetary signal HD3964 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 1M . The giant planets (0.3M J < M < 13M J ) 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. 4 m s −1 , 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. HD94771 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 HD74698. As BD-210397 includes a high jitter and HD74698 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 7 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 exoplanets. Our results are in agreement. Only one of the proposed exoplanets (HD94771 b) has an eccentricity of larger than 0.3, while HD62364 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 Figure 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. 7 Defining 'confirmed' as companions with well-constrained masses and periods (> 3σ).
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 4M 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 2 sin i = 2.01M 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 years allows us to detect Saturn-mass companions with periods of 10 years and brown dwarfs and stellar binaries with periods of up to 50 years. 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 HD11938B using the secondary peak in the combined CCFs, and of HD62364 b, HD5680B, HD221638B, HD33473C, HD11938B, and HD61383B via absolute astrometry. As byproducts of this survey, the brown dwarfs and binaries demon-strate the effectiveness of a long baseline. The impact of continuing with observations is made clear by the difference between the orbital solution of HD33473C 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.