S Udry

X Dumusque

C Lovis

D Segransan

R F Diaz

W Benz

F Bouchy

A Coffinet

G Lo Curto

M Mayor

D Ségransan

C Mordasini

F Motalebi

F Pepe

D Queloz

N C Santos

A Wyttenbach

R Alonso

A Collier Cameron

M Deleuil

P Figueira

M Gillon

C Moutou

D Pollacco

E Pompei

The HARPS search for southern extra-solar planets. XLII. Eight HARPS multi-planet systems hosting 20 super-Earth and Neptune-mass companions

Keywords: Planetary systems -Techniques: RVs -Techniques: spectroscopy -Methods: data analysis -Stars: individual: HD 20003, HD 20781, HD 20782

    The HARPS search for southern extra-solar planets. XLII. Eight HARPS multi-planet systems hosting 20 super-Earth and Neptune-mass companions

Introduction

The radial velocity (RV) planet search programs with the HARPS spectrograph on the ESO 3.6-m telescope [START_REF] Pepe | Optical and IR Telescope Instrumentation and Detectors[END_REF]Mayor et al. 2003) have contributed in a tremendous way to our The analysis of the radial-velocity measurements were performed using the Data and Analysis Center for Exoplanets (DACE) developed in the frame of the Swiss NCCR PlanetS and available for the community at the following address: https://dace.unige.ch/

The HARPS RV measurements discussed in this paper are available in electronic form at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsweb.u-strasbg.fr/cgi-bin/qcat?J/A+A/.

knowledge of the population of small-mass planets around solartype stars. The HARPS planet-search program on Guaranteed Time Observations (GTO, PI: M. Mayor) was on-going for 6 years between autumn 2003 and spring 2009. The high-precision part of this HARPS GTO survey aimed at the detection of very low-mass planets in a sample of quiet solar-type stars already screened for giant planets at a lower precision with the CORALIE Echelle spectrograph mounted on the 1.2-m Swiss telescope on the same site (Udry et al. 2000). The GTO was then continued within the ESO Large Programs 183.C-0972, 183.C-1005 and 192.C-0852 (PI: S. Udry), from 2009 to 2016. Within these programs, HARPS has allowed for the detection (or has contributed to the detection) of more than 100 extra-solar planet candidates (see detections in Díaz et al. 2016a;Moutou et al. 2015;Lo Curto et al. 2013;Dumusque et al. 2011a;Mayor et al. 2011;Pepe et al. 2011;Moutou et al. 2011;Lovis et al. 2011b). In particular, HARPS has unveiled the existence of a large population of low-mass planets including super-Earths and hot Neptunes previous to the launch of the Kepler satellite which provided us with an overwhelming sample of thousands of smallsize transiting candidates (Coughlin et al. 2016;Mullally et al. 2015;Borucki et al. 2011;Batalha et al. 2011). A preliminary analysis of the HARPS data showed at that time that at least 30 % of solar-type stars were hosting low-mass planets on short-period orbits with less than 50 days [START_REF] Lovis | IAU Symposium[END_REF]. A comprehensive analysis of our high-precision sample combined with 18 years of data from CORALIE allowed us later to precise this occurrence rate: about 50% of the stars surveyed have planets with masses below 50 M ⊕ on short to moderate period orbits (Mayor et al. 2011). Furthermore, a large fraction of those planets are in multi-planetary systems. This preliminary statistics of hot super-Earth and Neptune frequency is now beautifully confirmed by the impressive results of the Kepler mission.

With the RV technique, the variation of the velocity of the central star due to the perturbing effect of small-mass planets becomes very small, of the order or even smaller than the uncertainties of individual measurements. The problems to solve and characterize the individual systems are then multi-fold, requiring to disentangle the planetary from the stellar, instrumental, and statistical noise effects. Efficient statistical techniques, mainly based on a Bayesian approach, have been developed to optimise the process and thus the outcome of ongoing RV surveys (see e.g., Dumusque et al. 2016;Díaz et al. 2016b, and references therein). A large number of observations is however paramount for a complete probe of the planetary content of the system, and to take full advantage of the developed technique of analysis. In this context, focusing on the most observed, closest and brightest stars, an on-going HARPS is continuing the original observing efforts and providing us with an unprecedented sample of well observed stars. In this paper we describe 8 planetary systems hosting 20 planets. The detection of these planets has been announced in Mayor et al. (2011) 1 studying statistical properties of the systems discovered with HARPS. Most of them are super-Earths and Neptunian planets on relatively short periods, and member of a multi-planet system. We present here the orbital solutions for one 4-planet system around HD 20781, three 3-planet systems around HD 20003, HD 31527 and HD 136352, and four 2-planet systems around HD 21693, HD 45184, HD 51608, and HD 134060. In addition, we give updated parameters for the very eccentric planet orbiting HD 20782 (Jones et al. 2006), derived by combining HARPS, CORALIE and the published UCLES data. The paper is organised as follows. In Sec. 2 we discuss the primary star properties. Radial-velocity measurements and orbital solutions of each system are presented in Sects. 3 and 5, while Sec. 4 describes the framework used for analysing the data. We provide concluding remarks in Sec. 6.

Stellar characteristics

This section provides basic information about the stars hosting the planets presented in this paper. Effective temperature, gravity and metallicity are derived from the spectroscopic analysis of HARPS spectra and provided in Sousa et al. (2008). We used the improved Hipparcos astrometric parallaxes re-derived by van Leeuwen (2007) to determine the absolute V-band mag-nitude using the apparent visual magnitude from Hipparcos (ESA 1997). Metallicities, together with the effective temperatures and M V are then used to estimate basic stellar parameters (ages, masses) using theoretical isochrones from the grid of Geneva stellar evolution models, including a Bayesian estimation method (Mowlavi et al. 2012).

Individual spectra were also used to derive measurements of the chromospheric activity S-index, log(R HK ) and Hα, following a similar approach as used by Santos et al. (2000) and Gomes da [START_REF] Gomes Da Silva | [END_REF]. Using Mamajek & Hillenbrand (2008), we estimate rotational periods for our stars from the empirical correlation of stellar rotation and chromospheric activity index (Noyes et al. 1984). We also derived the v sin (i) from a calibration of the Full Width at Half Maximum (FWHM) of the HARPS Cross-Correlation Function (CCF) following a standard approach (e.g. Santos et al. 2002). All those extra indicators are used to disentangle small-amplitude planetary signals from stellar "noise" (see e.g. the cases of CoRoT-7 and α Centauri B; Queloz et al. 2009;Dumusque et al. 2012). The derived stellar parameters are summarized in Table 1.

HARPS RV measurements

Radial velocities presented here have been obtained with the HARPS high-resolution spectrograph installed on the 3.6m ESO telescope at La Silla Observatory (Mayor et al. 2003). The longterm m s -1 RV precision is ensured by nightly ThAr calibrations [START_REF] Lovis | Precision Spectroscopy in Astrophysics[END_REF]). On the short timescale of a night, the high precision is obtained using simultaneous ThAr (from 2003 to 2013) or Fabry-Perot étalon (since 2013) reference calibrations. The data reduction was performed with the latest version of the HARPS pipeline. In addition to the barycentric radial velocities with internal error bars, the reduction provides the Bisector Inverse Slope (BIS) of the HARPS CCF (Pepe et al. 2002;Baranne et al. 1996), as defined by Queloz et al. (2000), the FWHM and the Contrast of the CCF, as well as the Ca II activity indexes S and log(R HK ).

By construction of the HARPS high-precision sample from low-activity stars in the CORALIE volume-limited planet-search sample, the eight stars discussed here present low activity levels. The average values of the log(R HK ) activity index, estimated from the chromospheric re-emission in the Ca II H and K lines at λ = 3933.66 Å and 3968.47 Å, are low, ranging from -5.01 to -4.87. Activity-induced RV jitter on stellar rotation time scales, due to spots and faculae on the stellar surface, is thus expected to remain at a low level. The potential influence of stellar activity on RVs is nevertheless scrutinised closely when long-term variation of activity indexes are observed (magnetic cycle) or when the planetary and stellar rotation periods (or harmonics) are of similar values.

Low-mass planets are very often found in multi-planet systems (e.g. [START_REF] Lovis | [END_REF]Lovis et al. , 2011b;;Udry & Santos 2007;Mayor et al. 2009;Vogt et al. 2010;Díaz et al. 2016b;Motalebi et al. 2015;Latham et al. 2011;Lissauer et al. 2011Lissauer et al. , 2014;;[START_REF] Fabrycky | [END_REF], for examples of RV results and Kepler findings). The several components in the system give rise to often complex, low-amplitude RV signals, not easy to solve for. The optimal observing strategy for a given star is a priori unknown, the relevant planetary periods possibly extending over three orders of magnitude. An adequate and efficient observing strategy has been developed coping with the need to accumulate a large number of measurements and to probe different timescales of variation. We typically follow our targets every night during a few initial Notes. a Star already known as a planet host for which we give here an updated orbit including HARPS observations. observing runs, and after gathering a couple of dozens of observations we compare the observed dispersion with the expected jitter for the considered stellar spectral type. For significant variations we then continue monitoring at the same cadence when high-frequency variations are seen, or we adapt the frequency to possible variations at longer timescales. The majority of the stars in our HARPS high-precision program have been followed since 2003, gathering observations spanning more than 11 years. We are not considering here observations obtained after May 2015 when a major upgrade of the instrument was implemented (change of the optical fibres)2 . We also removed measurements obtained on the nights JD=2455115 and 2455399 as an unexplained instrumental systematic produced RV measurements that were off by more than 10 sigma on several stars observed those nights. Each RV measurement corresponds normally to a 15-minute HARPS exposure. This long exposure time allows to mitigate the short-timescale variations induced by stellar oscillations and therefore to improve RV precision (Santos et al. 2004;Dumusque et al. 2011c). For bright stars (V ≤ 6.5), the exposures are split in several sub-exposures within the 15 minutes in order to avoid CCD saturation while keeping the same strategy to mitigate stellar oscillations. The different time series presented in this paper are composed of binned points calculated through a weighted average of all points taken within an hour, so that all the observations taken within 15 minutes are binned together. The stars presented here have between 178 and 245 observations typically spread over ∼4000 days and with a sampling allowing the detection of planets with periods from below 1 day to the full span of the measurements. The obtained SNR at 550 nm typically ranges from 100 to 250 (up to 400 for HD 134060), depending on the star and weather conditions. The corresponding quantified uncertainties on the RVs range then from 0.33 to 0.76 m s -1 , including photon noise, calibration errors and instrumental drift uncertainty (see Table 2). This does not include other instrumental systematics like telescope centering and guiding errors, which are expected to be small but difficult to estimate. Possible additional uncertainties not included in this estimate might originate from the RV intrinsic variability of the star (jitter) due to stellar oscillations, granulation (Dravins 1982;Dumusque et al. 2011c) and magnetic activity (Dumusque et al. 2014b(Dumusque et al. , 2011b;;Meunier et al. 2010;Desort et al. 2007;Saar & Donahue 1997). Those extra errors induced by instrumental systematics or stellar signals will be taken into account when modeling the RVs (see Sec. 4).

The analysis of long series of HARPS radial velocities have recently revealed that, at the level of precision of HARPS, stitching of the detector (i.e. difference between the inter-pixel size and the gap between imprinted patches of pixels, see Sect. 4.3) may introduce a residual signal in the radial velocities at periods close to 6 months or 1 year when a strong stellar line crosses A&A proofs: manuscript no. udry2018 (Dumusque et al. 2015). To avoid this effect, new sets of radial velocities for the stars presented in this paper have been obtained by removing the corresponding zone of the spectra in the template used for the cross correlation.

The final one-hour binned HARPS rsadial velocities and log(R HK ) measurements are displayed in the left column of Fig. 1 for the 8 planet-host stars discussed in the paper. These velocities and the parameters inferred from the spectra and cross correlation functions are provided in electronic form at CDS. A sample of these data is provided in Table 3. The statistics of the RV series are listed in Table 2.

General approach of the data analysis

The data analysis presented in this paper is performed using a set of online tools available from the DACE platform3 . The RV tool on this platform allows to upload any RV measurements and then perform a multiple Keplerian adjustment to the data using the approach described in Delisle et al. (2016). Once preliminary Keplerian plus drift parameters are found using this iterative approach, it is possible to run a full Bayesian MCMC analysis using an efficient algorithm, allowing various choices of models for signals (and noise) of different origins (Díaz et al. 2016a(Díaz et al. , 2014)).

For each system, we performed a full MCMC analysis, probing the following set of variables for planetary signals: log P, log K, √ e cos ω, √ e sin ω and λ 0 , each one corresponding to the period, the RV semi-amplitude, the eccentricity, the argument of periastron and the mean longitude at a given reference epoch. We used √ e cos ω and √ e sin ω as free parameters rather than the eccentricity and the argument of periastron because they translate into a uniform prior in eccentricity (Anderson et al. 2011). The mean longitude (λ 0 ) is also preferred as a free parameter (instead of the mean anomaly or the date of passage through periastron) since this quantity is not degenerated at low eccentricities.

When no hint of magnetic activity can be seen in the different activity observables (S-index, log(R HK ), Hα-index), we fitted a RV model composed of Keplerians plus a polynomial up to the second order. The MCMC analysis is performed with uniform priors for all variables, with the exception of the stellar mass for which a Gaussian prior is used based on the information given in Table 1. We chose an error in stellar mass of 0.1 M for all systems to propagate the stellar mass error to the estimation of the planet masses. To take into account uncertainty due to instrumental systematics and/or stellar signals not estimated by the reduction pipeline, we included in the MCMC analysis a white-noise jitter parameter, σ JIT , that is quadratically added to the individual RV error bars.

When a magnetic cycle is detected in the different activity observables (S-index, log(R HK ), Hα-index), we decided, in addition to the model described above, to include two extra components in our RV model to account for the RV variation induced by this magnetic cycle. As explained in Meunier et al. (2016) and Dumusque et al. (2011a), the variation of the total filling factor of spots and faculae along a magnetic cycle change the total amount of stellar convective blueshift, as it is reduced inside spots and faculae due to strong magnetic fields, which thus change the absolute RV of the star. A positive correlation between the different activity observables and the RV is thus expected. In this paper, we decided to use the method proposed by Meunier & Lagrange (2013) to mitigate the impact of magnetic cycles, i.e. to adjust a linear correlation between RV and one of the activity index. Here we chose to use the log(R HK ) (Vaughan et al. 1978;Wilson 1968;Noyes et al. 1984). A larger number of spots and faculae are present on the stellar surface during high-activity phases of the magnetic cycle, which implies a stronger stellar jitter due to those surface features coming in and out of view. To account for this stellar jitter that changes in amplitude along the magnetic cycle, we included in the MCMC analysis two extra white-noise jitter parameters, σ JIT LOW and σ JIT HIGH , that correspond to the jitter during the minimum and the maximum phases of the magnetic cycle. We therefore replace σ JIT described in the precedent paragraph by σ JIT LOW + (σ JIT HIGH -σ JIT LOW ).norm( log(R HK )), where norm(log(R HK )) corresponds to log(R HK ) normalized from 0 to 1. For each RV measurement, a new jitter parameter is derived according to the activity level and is quadratically added to the corresponding RV error bar. This is similar to the approach adopted in Díaz et al. (2016a). The list of all the parameters probed by our MCMC is given in Table 4.

Before running a full MCMC analysis to obtain reliable posteriors for the orbital parameters of the planets present in the RV measurements, we first iteratively look for significant signals in the data using the approach of Delisle et al. (2016), which gives us a good approximation of the orbital parameters that are used as initial conditions for the MCMC analysis. The iterative approach follows several independent but complementary steps.

Removing long-term trends

Long-term trends in the data, due either to long-period companions (stellar or planetary) or magnetic cycles (Dumusque et al. 2011a, andLovis et al. (2011a) that provides a broad view of magnetic cycle behavior for the entire HARPS high-precision RV survey), might perturb the detection of planets on shorter-period orbits due to the spread of the energy of the signal at high frequencies in the Generalized Lomb-Scargle periodogram (GLS, Zechmeister & Kürster 2009;Scargle 1982) and aliases of these signals appearing at shorter periods. We fit a polynomial up to the second order to account for a non-resolved long-period companion. In addition, if a significant long-period signal is seen in -Stellar origin. The signal is due to stellar intrinsic phenomena. The most common one is the variation of the measured RVs due to spectral line asymmetries induced by spots and faculae coming in and out of view on the surface of the star, modulating the signal over a star rotational period (Haywood et al. 2016(Haywood et al. , 2014;;Dumusque et al. 2014a;Meunier et al. 2010;Desort et al. 2007;Santos et al. 2000;Saar & Donahue 1997). In order to recognize RV variations of intrinsic stellar origin, we compare the periods of the derived orbital solutions with the rotational periods of the stars and their harmonics estimated using the log(R HK ) average activity level (see Table 1, Mamajek & Hillenbrand 2008;Noyes et al. 1984). We also compare the orbital periods with the variation timescales of spectroscopic activity indicators derived from the CCF, such as the BIS and the FWHM, and derived from spectral lines sensitive to activity, such as the Ca II H and K lines (S-index and log(R HK )) or the Hα line (Hα-index). -Instrumental origin. As already mentioned earlier, spurious signals of small amplitudes with periods close to one year or one of its harmonics (half a year, a third of a year) can be created by a discontinuity in the wavelength calibration introduced by the stitching of the detector. The 4k×4k HARPS CCD is composed of 32 blocks of 512×1024 pixels. When each block were imprinted to form the detector mosaic, the technology at the time was not precise enough to ensure that pixels between blocks had the same size than intra-pixels. Therefore, every 512 pixels in the spectral direction, the CCD presents pixels that differs in size (Wilken et al. 2010).

[ m s -1 ] U instrumental jitter of instrument X ]-∞,∞[ σ JIT [ m s -1 ] U stellar jitter (If only one instrument is used, this parameter includes σ X ) ]-∞,∞[ γ X [km/s] U Constant velocity offset of instrument X ]-∞,∞[ lin [ m s -1 yr -1 ] U Linear drift ]-∞,∞[ quad [m 2 s -2 yr -2 ] U Quadratic drift ]-∞,∞[ log (P) log([days]) U Logarithm of the period [0,∞[ log (K) log([ m s -1 ]) U Logarithm of the RV semi-amplitude [0,∞[ √ e cos ω - U ]-1,1[ √ e sin ω - U ]-1,1[ λ 0 = M 0 + ω [deg] U Mean longitude (M 0 = mean anomaly) [0,
Block stitching may introduce a residual signal in the RVs at periods close to 6 months or 1 year when strong stellar lines crosse block stitchings due to the yearly motion of the Earth around the Sun. As mentioned in Sect. 3, in order to mitigate this effect new sets of RVs for the stars presented in this paper have been obtained by removing from the correlation masks used the spectral lines potentially affected. For each star, the systemic velocity will shift the stellar spectrum on the CCD, therefore an optimization of the correlation mask is done on a star-by-star basis (for more information, see Dumusque et al. 2015). -Observational limitations. Aliases have to be taken into account. The ones due to one year or one day sampling effects are well known (e.g., Dawson & Fabrycky 2010). They apply on all signals and not only on the planetary ones.

Description and analysis of individual systems

In this section, we present eight planetary systems including Neptune-mass and super-Earth planets. Presence of planets around these stars and preliminary system characterizations have been announced at the "Extreme Solar System II" conference held at Grand Teton USA in September 2011 and published in Mayor et al. (2011). We present here the detailed analysis of each system with updated data, describing and discussing the derived planetary system characteristics, as well as other (i.e. non-planetary) signals in the data. The typical photon-noise plus calibration uncertainty of the observations is 0.74 ms -1 , well below the observed dispersion of the RVs at 5.35 m s -1 . The RVs with their GLS periodogram and the log(R HK ) time-series are displayed in Fig. 1. A long-period variation is observed in the RV data, as well as in log(R HK ), indicative of a magnetic cycle effect on the velocities. To correct the velocities for this effect, we modeled the long-term variation by a Keplerian with parameters fixed and determined from the log(R HK ), except for the amplitude that was left free to vary.

After correction for the long-period variation due to the magnetic cycle plus fitting a second-order polynomial drift, two peaks very clearly emerge in the periodogram well above the 0.1 % p-value limit, at periods around 11.9 and 33.9 days (Fig. 2). Once those two signals are fitted for along with the magnetic cycle and the second-order polynomial drift, a significant signal at 184 days appears in the GLS periodogram of the RV residuals, with strong aliases at 127 and 359 days. A global fit including a second-order polynomial drift, a Keplerian for the magnetic cycle and three Keplerians for signals at 11.9, 33.9 and 184 days allows us to model all the variations seen in the RVs. No significant signal (with p-values smaller than 5 %) are then present in the RV residuals (Fig. 2, bottom-right plot).

After this preliminary phase looking for significant signals in the data, we searched for the best-fit parameters with an MCMC, using a model composed of a second-order polynomial drift, a linear correlation of the RVs with log(R HK ) to adjust the magnetic cycle effect, three Keplerians to fit for the signals at 11.9, 33.9 and 184 days, and two jitters that correspond to the instrumental plus stellar noise at the minimum and maximum of the magnetic cycle (see Sec. 4). This model converged to a stable solution. However, a significant signal in the residuals was still present near 3000 days (Fig. 3), a period similar to the one of the magnetic cycle, an indication that our modeling of the long period variation was not satisfactory, very probably because the magnetic cycle of HD 20003 is not fully covered, as seen in Fig. 1. We finally decided then to adjust a model composed of 3 Keplerians at 11.9, 33.9 and 184 days, an extra Keplerian at long-period to absorb the effect of the magnetic signal plus any long-term trend that could be present in the data, and one jitter corresponding to the instrumental plus stellar noise. In this case, the GLS periodogram of the residuals is extremely similar to the last plot in Fig. 2 and therefore no more significant signal is present in those RV residuals. Note that the period of the signal adjusting the magnetic cycle changed however from 4038 days (log(R HK )) to 3298 days (RVs). Although we do not understand fully the long-term signal present in the RVs of HD 20003, adjusting it with different models gives similar orbital parameters for the signals at 11.9, 33.9 and 184 days. We adopted as final solution the model with the magnetic cycle modeled by a long-period Keplerian with free parameters.

The best-fit for each inner planet, the signal at 184 days, and the RV residuals are displayed in Fig. 4. The best-fit parameters can be found in Table 5. A careful look at the activity indicators, log(R HK ), BIS SPAN and FWHM, after removing the longperiod signal induced by the stellar magnetic cycle, reveals no significant peaks that match the planetary signals found in our analysis (see Fig. A.1). We are therefore confident that those signals are not induced by stellar activity.

The correlation between the activity index log(R HK ) and the RV residuals when removing all the detected signals except the magnetic cycle effect can be seen in the left panel of Fig. 5. Note that this dependence in activity is included in our model when adjusting the long-period Keplerian.

A&A proofs: manuscript no. udry2018 Fig. 2. GLS Periodogram of the RV residuals of HD 20003 at each step after removing, from top-left to bottom-right, the magnetic cycle effect plus a second-order polynomial drift, then the two planets at 11 and 33 days one after the other, and finally the signal at 184 days. The GLS periodogram of the raw RVs is shown in Fig. 1.

Table 5. Best-derived solution for the planetary system orbiting HD20003. For each parameter, the median of the posterior distribution is considered, with error bars computed from the MCMC chains using a 68.3% confidence interval. σ O-C corresponds to the weighted standard deviation of the residuals around this best solution. All the parameters probed by the MCMC can be found in Annex, in Table B

A 3rd planet with a period of 184 days?

Although the signals at 11.9 and 33.9 days are clearly induced by planets orbiting around HD 20003, the signal at 184 days is more difficult to interpret as it is very close to half a year, with strong aliases at 127 days and 359 days (close to 1/3 and 1 year).

As explained in Sec. 4.3, signals at a year or harmonics of it can be induced by a discontinuity in the wavelength calibration introduced by tiny gaps between the different quadrants of the detector. However, the RV data that we analyze here have been corrected for this effect following the work of Dumusque et al. (2015). This correction has been already applied to the RVs of several stars, and has always been successful in removing the spurious signal. The fact that the RV residuals after removing the two planets, and the long-period effect of the magnetic cycle do not correlate with the barycentric Earth velocity (BERV, see right panel of Fig. 5) disfavors the hypothesis that the 184day signal is due to a discontinuity in the wavelength calibration introduced by tiny gaps between the different quadrants of the detector. In conclusion, we would be inclined to claim that this signal at 184 days is a real planetary signal. However, although 5. intensively tested, there is always the possibility that the data reduction system does not correct for all the instrumental effects in this particular case. We leave the signal as a potential one, and encourage other teams to confirm it using other pipelines or other facilities than HARPS.

Discussion of the planetary system

A surprising aspect of this system is the high eccentricity of the inner planet, e=0.38, while the planet at 33.9 days has an eccen-tricity close to 0. From a planet formation or a dynamical point of view, this is not straightforward to explain. There might be actually several possibilities to explain the observed configuration.

It is out of the scope if this paper to check all of them in details, we just mention here the basic ideas.

1) The first and simplest possibility is that the high eccentricity of the planet at 11.9 days actually hides a planet at half the orbital period (Anglada-Escudé et al. 2010), in (or close to) a 2:1 resonance configuration. We therefore tried to fit a model with an extra planet with initial period at 5.95 days, setting the eccentricities of the 5.95 and 11.9-day planets to zero or letting them free to vary. Model comparison shows that, in both cases, the more complex solution is disfavoured with a ∆BIC of 8.3 and 24.3, respectively. It seems therefore that the significant eccentricity of planet b is real.

2) The two planets are close to the 3:1 resonance. An increase of the eccentricity could have happened if, in the past, the two planets have been migrating trapped in the resonance. Then a specific event made the two planets go slightly out of resonance on the nowadays observed orbits. Possible scenarios for this event could be instability when the gas disappeared or the presence of an additional planet (e.g. resonant chain) that was ejected when the eccentricity increased (scattering).

3) Invoking additional bodies in the system raises the possibility of an external companion that could induce a Kozai-lidov effect on the inner planets in case of a large mutual inclination (Kozai 1962;Lidov 1961). Such a companion could actually be hinted by our difficulty of modelling the long-period variation observed in the radial velocities. It is however not very clear to understand the differential effect on the two inner planets.

4) Finally, a very appealing potential explanation might rely on the spin-orbit evolution of planets feeling the tidal effect of the central star and the gravitational perturbation of another planetary companion. Correia et al. (2012) have shown that "under some particular initial conditions, orbital and spin evolution cannot be dissociated and counter-intuitive behaviors can be observed, such as the secular increase of the eccentricity. This effect can last over long timescales and may explain the high eccentricities observed for moderate close-in planets".

HD 20781: A packed system with 2 super-Earths and 2

Neptune-mass planets HD 20781 is part of a visual binary system including another star, HD 20782, known to host a planet on a 595-day very eccentric orbit (Jones et al. 2006). We take the opportunity of the discovery of a compact system of small planets around HD 20781 to provide here as well an updated solution for the planet around HD 20782 using HARPS data.

Analysis of the HD 20781 system

HD 20781 was part of the original high-precision HARPS GTO survey and the star has been then followed for more than 11 years (4093 days). Over this time span, we gathered a total of 226 high signal-to-noise spectra (<S/N> of 112 at 550 nm) corresponding in the end to 216 RV measurements binned over 1 hour. As reported in Table 2, the typical precision of individual measurements is 0.76 m s -1 including photon noise and calibration uncertainties. The raw RV rms is significantly higher, at 3.41 m s -1 , pointing towards additional variations in the data of potentially stellar or planetary origin, assuming the instrumental effects are kept below the photon-noise level of the observations. As a first approach we looked at the RV and activity index time series shown in Fig. 1. No significant long-term variation is observed in log(R HK ) data and no long-term variation is visible in the GLS periodogram of the velocity time series. We conclude that there is no noticeable sign of a magnetic activity cycle for this star. The average value of log(R HK ) at -5.03 ± 0.01 is also very low with a very small dispersion, similar to the Sun at minimum activity.

Due to the small activity level and the large number of observations, the GLS periodogram of the velocity series is very clean, with significant peaks clearly coming out of the noise background. So, even if the final characterization of the planetary system parameters is performed through a Bayesian-based MCMC approach, a step-by-step analysis of the system, characterizing and then removing one planet after the other from the data, will provide an excellent illustration of the significance of the planet detection in this system. The most prominent peak in the GLS periodogram is at a period around 20 days. Deriving a Keplerian solution for this planet and then removing the corresponding signal from the raw RVs makes a clear signal at 86 days appear in the GLS periodogram of the residuals (Fig. 6). Keeping on with the same approach, we can clearly identify, sequentially, significant signals first at 13.9 and then at 5.3 days. It has to be noted here that none of the significant periods in the data are close to the stellar rotational period, 46.8 days, estimated from the activity index (Table 1).

As described in Sec. 4, the final determination of the planetary system parameters is performed using a MCMC sampler and a model composed of four Keplerians representing the planetary signals, and an extra white-noise jitter to consider potential stellar or instrumental noise not included in the RV error bars (the σ JIT parameter, see Table 4). Phase-folded planetary solutions, as well as the RV residuals around the best solution, are displayed in Fig. 7. The best-fit parameters are reported in Table 6.

Looking at the periodograms of the log(R HK ), the BIS SPAN and the FWHM in Fig. A.2, to check if any announced planet matches any signal in the activity indicators, we see that the log(R HK ) times series presents signals at 115, 81 and 68 days, and the FWHM time series at 380 days. This later signal is probably due to interaction with the window functions, creating a signal near a year. The former signals are more difficult to interpret as they are not compatible with the estimated rotational period of the star, 46.8 days (Table 1). In particular, we can note that the 86-day Neptune-mass planet has a period close to the detected peaks in the log(R HK ) time-series at 81 days, and furthermore is close to the 1-year alias of the 115-d signal (days). The amplitude of these signals however very low, at the level of ∼0.01 dex, 20 times smaller than the solar magnetic cycle variation. Such signals can therefore probably not be responsible for the 2.6 m s -1 periodic variation detected in the RVs at 86 days.

Our conclusion is thus that the system HD 20781 hosts two inner super-Earths with periods of 5.3 and 13.9 days and two outer Neptune-mass planets with periods of 29 and 86 days. With small eccentricities, stability is not a concern for this system of small-mass planets in a moderately compact configuration.

HD 20782: More planets in this visual binary system

The star HD 20782 is the brightest companion of the HD 20781-HD 20782 binary system. It is known to harbor a 595-day very eccentric Jupiter planet (Jones et al. 2006). A total of 71 high signal-to-noise spectra (<S/N> of 181 at 550 nm) were obtained with HARPS on this target, which translates into 68 RV measurements after binning the data over 1 hour. Because of the very large semi-amplitude and eccentricity of the RV signal induced by the planet, we decided to include in addition to the HARPS very precise measurements the lower precision RV data obtained with UCLES (published in Jones et al. 2006) and CORALIE.

We searched for the best-fit parameters using a MCMC sampler. The solution converges to an extremely eccentric Jupitermass planet with a period of 597 days. The phase-folded planetary solution, as well as the RV residuals around the best solution and their GLS periodogram, are displayed in Fig. 8. The best-fit parameters are reported in Table 6. Although the eccentricity, the amplitude and the argument of periastron are compatible within one sigma with the values presented in Jones et al. (2006), the period and argument of periastron are not compatible, even if different by less than 2%. This can be explained by the addition of the CORALIE and HARPS data that allows to sample much better the periastron passage.

The high eccentricity of HD 20782 b has commonly been explained via the Kozai-Lidov mechanism (Kozai 1962;Lidov 1961) induced by HD 20781. 3). Those signals are also present in the FWHM and the bisector span of the CCF, although less significant. The signal at ∼740 days, close to two years is probably due to the sampling of the data, and the 33.5-day signal is likely the stellar rotation period. This value is compatible with the rotation period derived from the mean activity log(R HK ) level, i.e. 36 days (Table 1).

Looking at the raw RVs and their GLS periodogram in Fig. 1, it is clear that the observed magnetic cycle has an impact on the measured RV measurements. To remove the RV contribution of the magnetic cycle, we remove from the RVs a Keplerian model that has the same parameters as the Keplerian fitted to log(R HK ) but with an amplitude free to vary. The GLS periodogram of the RV residuals after correcting for the magnetic cycle effect are displayed in the top panel of Fig. 9. A highly significant signal at 54 days is present in the data. When removing this signal by fitting a Keplerian, an extra signal at 23 days is seen in the residuals (middle panel of Fig. 9). No signal with p-value smaller than 1 % appears in the residuals of a two-Keplerian model; we therefore stop looking for extra signals. Note however that the most significant signal left in the GLS periodogram corresponds to a period of 16 days, likely the first harmonic of the stellar rotation period, which is expected from stellar activity (Boisse et al. 2011).

After this preliminary phase looking for significant signals in the data, we search for the best-fit parameters with an MCMC, using a model composed of a linear correlation of the RVs with log(R HK ) to adjust the magnetic cycle effect, two Keplerian functions for the signals at 23 and 54 days, and two jitter values correspond to the instrumental plus stellar noise at the minimum and maximum of the magnetic cycle (see Sec. 4). The best-fit for each planet and the RV residuals are displayed in Fig. 10, and the best-fit parameters can be found in Table 7. None of the two detected planet corresponds to signals found in the activity indicators (see Fig. A.3).

When looking at the RV residuals, the scatter is still rather high even after removing all the significant signal detected in the data. Although in Hipparcos the star is catalogued as a G8 dwarf, the spectroscopic survey of nearby stars NSTAR finds that that HD 21693 is a G9IV-V therefore a slightly evolved star (Gray et al. 2006). Evolved stars presents higher photometric and RV jitter associated with more significant granulation, which might explain this significant residual jitter (Bastien et al. 2014(Bastien et al. , 2013;;Dumusque et al. 2011c).

The correlation between the activity index log(R HK ) and the RV residuals removing only the two-planet solution can be seen in Fig. 11. Note that this correlation is considered in the MCMC model we used to fit the RV data.

Our analysis of the HARPS RV measurements of HD 21693 provides strong evidence that 2 Neptune-mass planets orbit the star, with periods of 22.7 and 53.7 days. With such periods, this planetary system is close to a 5:2 resonance but with a period ratio smaller than the 5:2 commensurability. If the system formed through a convergent migration scenario trapping the planets in the 5:2 resonance, some process has to be invoked to push the system out of the resonance configuration after the disappearance of the protoplanetary disk. 7.

HD 31527: A 3-Neptune system

In total, 257 high signal-to-noise spectra (<S/N> of 180 at 550 nm) of HD 31527 were gathered over a time span of 11 years (4135 days). This results in 245 observations of the star when data are binned over 1 hour. The 3.19 m s -1 observed dispersion of the RVs is much larger than the typical RV precision of 0.64 m s -1 , pointing towards the existence of extra signals in the data. The raw RVs, their GLS periodogram and the calcium activity index of HD 31527 are shown in Fig. 1. The calcium activity index log(R HK ) does not show any significant variation as a function of time, therefore we do not expect the Table 7. Best-fitted solution for the planetary system orbiting HD 21693. For each parameter, the median of the posterior is considered, with error bars computed from the MCMC chains using a 68.3 % confidence interval. σ O-C corresponds to the weighted standard deviation of the residuals around this best solutions. All the parameters probe by the MCMC can be found in Annex, in Table B RVs to be affected by long-period signals generally induced by magnetic cycles. The mean activity level of the star is equal to log(R HK )= -4.96, close to solar minimum. The RVs should therefore be exempt as well of activity signal at the rotational period of the star and its harmonics due to active regions on the stellar surface (Boisse et al. 2011). This is confirmed when looking at the periodograms of the different activity indicators in Fig. A.4. Only a signal at 400 days in the BIS SPAN is sig-nificant. This signal, not too far from a year, might be due the interplay between the time series and the window function. As we can see in the GLS periodogram of the raw RVs (Fig. 1), two extremely significant signals appear at 17 and 52 days. After fitting these signals with a two-Keplerian model, a third significant signal at 271 days can be seen in the RV residuals (middle panel of Fig. 12). Finally, the residuals of a three-Keplerian model do not show any signal with p-value smaller than 1 % (bottom panel of Fig. 12). The fact that no signal is present at the estimated rotation period of the star (19 days, Table 1) or its harmonics confirms that the RVs of HD 31527 are not affected by significant activity signal.

After this first search for significant signals in the data, we fitted, using a MCMC sampler, a three-Keplerian model to the data including a white-noise jitter component to account for stellar and instrumental uncertainties not included in the RV error bars. The best-fit for each planet, as well as the RV residuals, can be seen in Fig. 13. We report in addition the best-fit parameters in Table 8. Fig. A.4).We therefore conclude that HD 31527 hosts 3 Neptune-mass planets, with periods of 16.6, 51.2 and 272 days. The star is a G2 dwarf like the Sun, therefore the outer planet in this system lies on an orbit between the ones of Venus and the Earth, therefore in the habitable zone of its host star (Selsis et al. 2007). This planet, with a minimum mass of 13 Earth-masses, is however most likely composed of a large gas envelope (Rogers 2015;Wolfgang & Lopez 2015;Weiss & Marcy 2014), except if it is similar to Kepler-10 c in composition (Dumusque et al. 2014b). 8.

None of the signals announced here matches signals in the different activity indicators (see

Table 8. Best-fitted solution for the planetary system orbiting HD 31527. For each parameter, the median of the posterior distribution is considered, with error bars computed from the MCMC chains using a 68.3 % confidence interval. σ O-C corresponds to the weighted standard deviation of the residuals around the best solution. All the parameters probeb by the MCMC can be found in Annex, in Table B

HD 45184: A system of two close-in Neptunes

We gathered a total of 309 high signal-to-noise spectra (<S/N> of 221 at 550 nm) of the G1.5 dwarf HD 45184 during a time span of 11 years (4160 days). This results in 178 RV measurements, when the data are binned over 1 hour. The average precision of the RVs is 0.41 m s -1 considering only photon noise and calibration uncertainties. The raw RV rms is much higher, 4.72 m s -1 , which implies that significant signals are present in the data. In Fig. 1, we display the raw RVs and their GLS pe-riodogram, and the log(R HK ) time series. We see that a significant magnetic cycle affects log(R HK ), with values ranging from -5.00 to -4.86 with a periodicity of 5 years. This magnetic cycle is therefore smaller in amplitude and period than the solar cycle, despite the similarities of the derived stellar parameters with those of the Sun. To see if significant signals were present in the calcium activity index despite the long-period magnetic cycle, we fitted the log(R HK ) time series with a Keplerian. In the residuals, a strong signal at 20 days is present, likely corresponding A&A proofs: manuscript no. udry2018

Fig. 14. From top to bottom: GLS periodogram of the RVs after removing the effect induced by the stellar magnetic cycles and then one planet after the other in the analysis of HD 45184. The GLS periodogram of the raw RVs is shown in Fig. 1.

to the stellar rotation period (see Fig. A.5). This value is fully compatible with the rotation estimated using the log(R HK ) average level (21.5 days, see Table 1; Mamajek & Hillenbrand 2008;Noyes et al. 1984).

Looking at the raw RVs of HD 45184 and their GLS periodogram in Fig. 1, we see that the magnetic cycle observed in log(R HK ) has an influence on the RVs. To remove the RV contribution of the magnetic cycle, we fitted the log(R HK ) with a Keplerian, and removed the same Keplerian from the RVs leaving the amplitude as a free parameter. In the residuals, displayed in the top panel of Fig. 14, we see a significant signal at 6 days. Once this signal is removed by fitting a Keplerian with the corresponding period, another signal at 13 days appears (see middle panel of Fig. 14). Finally after removing a two-Keplerian model to the RVs corrected for the magnetic cycle effect, no signal with p-value smaller than 10 % remains. We therefore stop looking for extra signals in the data. Note that although not significant, the highest peak in the GLS periodogram of the RV residuals is at 18.6 days, likely the imprint of the stellar rotation period. The RVs are therefore slightly affected by stellar activity, however at a level that is not perturbing the detection and characterization of the two planets at 6 and 13 days.

After the preliminary phase of looking for significant signals, we fitted the RVs using a MCMC sampler and a model composed of a linear correlation of the RVs with log(R HK ) to adjust the magnetic cycle effect, two Keplerians to fit for the signals at 5.9 and 13.1 days, and two jitters that correspond to the instrumental plus stellar noise at the minimum and maximum of the magnetic cycle (see Sec. 4). Each planet with its best fit can be seen in Fig. 15, as well as the RV residuals after the best solution has been removed. The best-fit parameters are reported in Table 9. None of the signals at 5.9 and 13. 9.

In Fig. 16, we show the RV residuals after removing the bestfit solution for planets b and c as a function of the log(R HK ). The observed strong correlation indicates that most of the residuals are due to activity-related effects and motivates the use of our model that includes a correlation between log(R HK ) and RVs to mitigate the effect of long-term activity.

HD 51608: 2 Neptune-mass planets

Over a time span of 11 years (4158 days), 218 high signalto-noise spectra (<S/N> of 133 at 550 nm) of HD 51608 were gathered with HARPS, resulting in a total of 216 measurements binned over 1 hour with a typical photon-noise and calibration uncertainty of 0.62 m s -1 . This value is significantly below the observed 4.07 m s -1 dispersion of the RVs, pointing towards the existence of additional signals in the data. The raw RVs, their GLS periodogram and the calcium activity index of HD 51608 are shown in Fig. 1. A small, albeit significant, long-term variation can be seen in log(R HK ), with values in log(R HK ) ranging from -5.04 to -4.96 and with a period of 11 years. Although the period of this magnetic cycle is very similar to that of the Sun, its amplitude is much lower. After fitting this long-period sig-Fig. 16. RV residuals when removing all the detected signals except the magnetic cycle effect plotted as a function of the activity index log(R HK ) for HD 45184. The observed correlation indicates that most of the RV residual variation is due to activity-related effects.

Table 9. Best-fitted solution for the planetary system orbiting HD 45184. For each parameter, the median of the posterior distribution is considered, with error bars computed from the MCMC chains using a 68.3% confidence interval. σ O-C corresponds to the weighted standard deviation of the residuals around the best solution. All the parameters probed by the MCMC can be found in Annex, in Table B nal in log(R HK ), a signal with a p-value of ∼5 % and a period of 37 days is detected in the log(R HK ) residuals, and can also be seen in the BIS SPAN and the FWHM of the CCF (see Fig. A.6). This is likely a signature of stellar activity as the mean log(R HK ) level gives an estimated rotation period of 40±4 days (Table 1; Mamajek & Hillenbrand 2008;Noyes et al. 1984).

In the raw RVs, very strong signals at 14 and 96 days are present (see Fig. 1). Once fitting a two Keplerian model to account for those signals, a long-period signal with a p-value Fig. 17. From top to bottom: GLS periodogram of the RV residuals of HD 51608 after removing the effect of the magnetic cycle and, at each step, after removing one planet after the other. The GLS periodogram of the raw RVs is shown in Fig. 1. smaller than 0.1 % appear in the GLS periodogram (see middle panel of Fig. 17). This signal is induced by the stellar magnetic cycle, and we remove it as in the precedent cases by fitting a Keplerian to the log(R HK ), and removing the same Keplerian from the RVs leaving the amplitude as a free parameter. After fitting the effect of the two planets plus the magnetic cycle, no significant signal with p-value smaller than 5 % is left in the residuals.

After this preliminary stage of checking significant signals in the data, we searched for the best-fit parameters using a MCMC sampler, and selecting a model composed of a linear correlation between the RVs and log(R HK ) two Keplerians to fit for the planetary signals at 14.1 and 96.0 days, and two jitters that correspond to the instrumental plus stellar noise at the minimum and maximum of the magnetic cycle (see Sec. 4). The best-fit solution for the two planets are shown in Table 10 and The correlation between the activity index log(R HK ) and the RV residuals removing only the two-planet solution can be clearly seen in Fig. 19. Note that this correlation is considered in the model we used to fit the RV data.

With masses of 14.3 and 12.8 Earth-masses, the two planets orbiting HD 51608 are likely Neptune-type planets. There is however still the possibility to have (one of) them similar in composition to Kepler-10 c (Dumusque et al. 2014b).

HD 134060: A short-period Neptune on an eccentric

orbit with a long-period more massive companion

Data analysis

A total of 335 high signal-to-noise spectra (<S/N> of 199 at 550 nm) of HD 134060 have been gathered over a time span One very significant signal at 3.3 days can be seen in the GLS periodogram of the raw RVs (see Fig. 1). Once this signal is fitted with a Keplerian, another significant peak appears at 1292 days, as can be seen in Fig. 20. After fitting simultaneously those two signals, nothing is left in the RV residuals with p-values smaller than 10 %.

To get the best possible orbital parameters for those two planets with reliable error bars, we perform an MCMC analysis with a model composed of two Keplerians plus a white-noise jitter to account for stellar and instrumental uncertainties not included in the RV error bars. The best-fit parameters can be found in Fig. 19. RV residuals around the best derived solution without considering the magnetic cycle effect plotted as a function of the activity index log(R HK ) for HD 51608. The observed correlation indicates that most of the RV residual variation is due to activity-related effects.

Table 10. Best-fitted solution for the planetary system orbiting HD51608. For each parameter, the median of the posterior distribution is considered, with error bars computed from the MCMC chains using a 68.3 % confidence interval. σ O-C corresponds to the weighted standard deviation of the residuals around the best solution. All the parameters probed by the MCMC can be found in Annex, in Table B

σ (O-C) [m s -1 ] 1.60 log (Post) -409.1398 +2.4423 -3.1572
Table 11. The planetary signals, folded in phase, can be seen in Fig. 21 along with the RV residuals shown in the bottom plot.

Origin of the high eccentricity of the inner planet?

The MCMC converges to a solution with an inner planet of minimum mass 10.1 M ⊕ and period 3.27 days, on a relatively high eccentricity orbit with e = 0.45. An inner eccentric planets in a system with an additional outer planet is reminiscent of the system HD 20003 above. The main differences are that in HD 134060 the outer companion is farther out and more massive, and the inner planet is experiencing stronger tides from the star making the circularization timescale shorter. Nevertheless some of the Table 11. Best-fitted solution for the planetary system orbiting HD 134060. For each parameter, the median of the posterior distribution is considered, with error bars computed from the MCMC posteriors using a 68.3 % confidence interval. The value σ (O-C) corresponds to the weighted standard deviation of the residuals around the best solution. All the parameters probed by the MCMC can be found in Annex, in potential mechanisms proposed for explaining the high eccentricity might be considered as well.

1) The high eccentricity can hide a planetary system in 2:1 resonance, thus the existence of another planet at half its orbital period (Anglada-Escudé et al. 2010). We therefore tried to fit a model with an extra planet at 1.65 days, fixing the eccentricities at zero or leaving them free to vary. In both cases, the more complex solution is disfavoured by Bayesian model comparison with a ∆BIC of 3.1 and 23.5, respectively. We therefore keep the simplest solution with the relatively high eccentricity of the inner planet.

2) This planet has a long-period companion that has a minimum mass three times larger. In the case of the two planetary orbits being strongly non-coplanar it is likely that the long-period planet is perturbing its inner companion through a Kozai-Lidov mechanism causing libration of its orbit (Kozai 1962;Lidov 1961). During this process, the eccentricity of the inner planet can reach very high values. The inner planet starts then to interact with its host star during close approaches, inducing a circularization of the inner orbit. Because of the conservation of the total angular momentum, eccentricity can increase only if the inclination of the orbit changes. Inner planets under the influence of a Kozai-Lidov mechanism are therefore likely to be on inclined orbits relative to the stellar rotational plane. This can be measured using the Rossiter-McLaughlin effect if by chance the planet transits its host star. This has been checked for HD 134060 by [START_REF] Gillon | [END_REF] using the Spitzer space telescope with unfortunately a null result. Another difficulty here is that for a planet at 3.3 days, the circularization timescale is normally very short for typical values of the tidal dissipation parameters, therefore a priori preventing us from observing the system in such a configuration.

3) Scattering is always a valid possibility when high eccentricities are involved. As for the previous point, preventing the system to then circularize is also the challenge.

A&A proofs: manuscript no. udry2018 Fig. 22. From top to bottom: GLS periodogram of the residuals at each step, after removing one planet after the other in the analysis of HD 136352. The GLS periodogram of the raw RVs is shown in Fig. 1.

4) The explanation invoking spin-orbit coupling during the evolution of the planet orbit under the influence of tides from the central star and gravitational perturbation from the outer companion (Correia et al. 2012) might finally be the most appealing process. It has however still to be demonstrated that the process works for the present set of parameters. As often, the challenge is in the balance of timescales. To help, HD 134060 seems to be a relatively young sub-giant with an age estimate of 1.75 Gyr (Delgado Mena et al. 2015).

In any case, with the high eccentricity of its inner planet HD 134060 provides a very interesting system to study planetary formation and evolution.

HD 136352: A 3-planet system

HD 136352 was part of the original high-precision HARPS GTO survey and the star has been followed for nearly 11 years (3993 days). Over this time span, we gathered a total of 649 high signal-to-noise spectra (<S/N> of 231 at 550 nm) corresponding in the end to 240 RV measurements binned over 1 hour. As reported in Table 2, the typical precision of individual measurements is 0.33 m s -1 including photon noise and calibration uncertainties, an order of magnitude smaller than the observed raw RV rms of 2.74 m s -1 . This suggests that significant signals, of stellar or planetary origin, are present in the data. As a first approach we looked at the log(R HK ) activity index time series in Fig. 1. No significant long-term variation is observed in the log(R HK ) data and no long-term variation is visible in the GLS periodogram of the velocity time series. We conclude that there is no noticeable sign of a magnetic activity cycle for this star. The average value of log(R HK ) at -4.95 is low with a small dispersion of ∼ 0.01, close to the Sun at minimum activity. No significant effect of stellar activity in the RV measurements is thus expected for this star. This is confirmed by the fact that no significant signals appears in the periodogram of the different activity indicators in Fig. A.8. Due to the small activity level and the large number of observations, the GLS periodogram of the velocity series is actually very clean, with peaks at 27.6, 11.6 and 108 days, in order of decreasing significance (see Fig. 22). After fitting those three signals with Keplerians, a study of the GLS periodogram of the RV residuals shows a peaks at 123 days with a p-value between 1 and 0.1 %, thus an interesting signal that we will consider in the MCMC analysis.

The best orbital parameters for the three planets orbiting HD 136352 are searched for using a MCMC sampler using a model composed of three Keplerians and an extra white-noise jitter to account for instrumental and stellar uncertainties not included in the RV error bars. Phase-folded planetary solutions are displayed in Fig. 23, as well as the RV residuals around the best solution. The best-fit parameters are reported in Table 12.

Our analysis of HD 136352 converges to the detection of three planet orbiting this G4V star. With minimum masses of 4.8, 8.6 and 10.8 M ⊕ , HD 136352 host three super-Earth on orbits ranging from 11 to 108 days. The two inner planets, with periods of 11.6 and 27.6 days, are close to a 5:2 resonance and, contrary to the two planets in HD 21693, with a period ratio larger than 2.5.

As mentioned above, an interesting signal is present in the residuals at 123 days. Despite the period being suspicious (1/3 of a year), we ran another MCMC trial including this fourth signal as an additional planet (Keplerian). The fit converged towards a non-eccentric signal with an amplitude of 0.65 m s -1 and a period of 122.6 days. However, when comparing the threeand four-planet model solutions, the case with three planets is strongly favored, with a ∆BIC of 39.6. We also tried to add to the three-planet model a polynomial of the first or the second order to check if this 123-day signal could be due to a very-long period companion, whose orbit is not covered by the data, creating an alias at a period of one year or one of its harmonics. This approach, however, did not reduce the signal at 123 days and was disfavoured by a model comparison using the BIC. We also looked for other possibilities of the three planet scenario, fitting the aliases of the 11 and 27.6-day signals, 0.91 and 0.96 day, respectively. These other possibilities are also ruled out with difference in BIC>17. Finally we tested the sensitivity of the GLS periodogram to outliers. We found that by simply removing two of them, for example JD 2455411 and 24556168, the amplitude of the peak found at 123 days goes above a p-value of 1 %. Therefore, several arguments points in the direction of this interesting, albeit not conclusively significant, signal at 123 days being more likely an artefact induced by noise in the data or interaction with the window function rather than a bona-fide planet. We therefore keep for HD 136352 the three planet solution, with periods of 11.6, 27.6 and 108 days.

Conclusion

We have reported the discovery and detailed characterization of eight planetary systems hosting twenty low-mass planets, plus a possible additional planet on a 6-month orbit, discovered with the HARPS Echelle spectrograph mounted on the 3.6-m ESO telescope located at La Silla Observatory. We also improved the characterization of the extremely eccentric Jupiter orbiting HD20782. 12.

Table 12. Best-fitted solution for the planetary system orbiting HD 136352. For each parameter, the median of the posterior distribution is considered, with error bars computed from the MCMC posteriors using a 68.3 % confidence interval. The value σ (O-C) corresponds to the weighted standard deviation of the residuals around the best solution. All the parameters probed by the MCMC can be found in Annex, in Table B.9. See Table 4 -Very small-mass planets with minimum masses below 6 M ⊕ , that are found orbiting HD 20781 and HD 136352, on short periods, less than 15 days. -More massive planets in the super-Earth to Neptune transition regime, spanning a range in minimum masses from 8 to 17 M ⊕ , and a range in periods from a few days to nearly a year. They represent the bulk of our detections. -Massive planets found on long-period orbits.

The RV technique sensitivity goes down when moving towards small-mass and long-period planets. It is therefore not surprising that most of our detections are in the second mass category. The lack of massive objects on short-period orbits is already well established, and the detection of extremely smallmass planets, although numerous from Kepler statistics (e.g. Coughlin et al. 2016) and previous RV surveys (e.g. Butler et al. 2017;Mayor et al. 2011), is challenging for RV surveys because of the perturbing effect of additional stellar signals (Dumusque et al. 2016) and the preponderance of appearance of these planets in multiple systems.

The detection of these planets had been announced in Mayor et al. (2011) studying statistical properties of the systems discovered with HARPS. However, with more data in hand, it was possible to discover new bona-fide planets. This is the case for the two inner super-Earth orbiting HD20781 and the 13.1-day period super-Earth orbiting HD45184. There is also the 184-day signal found in the timeseries of HD20003 even if, as explained in Sec. 5.1, we cannot exclude an instrumental origin to this signal. These new detections show that gathering more data helps in detecting small-mass planets on short-period orbits as well as long-period signals. Although expensive, characterizing those signals is therefore a necessary task to get an as complete as possible view of planetary systems.

The most important targets in today's exoplanet studies are the transiting small-mass planets orbiting bright stars. They allow for exquisite determination of their orbital and physical properties, bringing priceless constrains to our understanding of planet formation, of the establishment of planetary system architecture, and of the diversity of planet internal compositions. Those candidates would also be excellent candidates for further atmospheric characterization. A dedicated Spitzer survey looking for potential transit events induced by the shortest-period planets in our sample was conducted by [START_REF] Gillon | [END_REF]. Unfortunately, after searching for transits of the innermost planets orbiting HD 20003, HD 20781, HD 31527, HD 45184, HD 51608 and HD 134060, no detection was reported in this study. The search will continue, in particular with the Swiss-ESA led CHEOPS satellite to be launched end of 2018 [START_REF] Fortier | Space Telescopes and Instrumentation 2014: Optical, Infrared, and Millimeter Wave[END_REF].

Some systems presented here are interesting in terms of formation process and architecture. HD 20003 hosts two Neptunemass planets with periods of 11.9 and 33.9 days, thus close to a 3:1 commensurability. This configuration might explain the relatively high eccentricity measured on the innermost planet in this system, especially if some resonance is involved in the process. We note also that the two Neptune-mass planets orbiting HD 21693 and the two innermost super-Earth found around HD 136352 are close to a 5:2 resonance. Finally, HD 134060 is also an interesting dynamical system as it harbors a small-mass planet on a fairly eccentric 3-day orbit, accompanied by a more massive long-period planet at ∼1300 days. Several possibilities for the formation of such configurations have been mentioned in the text. The list is probably not exhaustive and further theoretical work is need to analyse these systems in a deeper way.

Here we just checked that the stability of the proposed solutions is verified over a few tens of thousands of years, avoiding the systems to quickly fall apart.

  @unige.ch Based on observations made with the HARPS instrument on the ESO 3.6 m telescope at La Silla Observatory under the GTO program 072.C-0488 and Large program 193.C-0972/193.C-1005/.

Fig. 1 .

 1 Fig. 1. Left column: From top to bottom, HARPS RV measurements as a function of barycentric Julian Date obtained for HD 20003, HD 20781, HD 21693, HD 31527, HD 45184, HD 51608, HD 134060 and HD 136352. Error bars only include photon and calibration noise. Middle column: log(R HK ) activity indicator as a function of time for the same stars. Right column: GLS periodogram of the corresponding RV measurements.

5. 1 .

 1 HD 20003: Eccentric, close-in Neptune-mass planets close to the 3:1 commensurability 5.1.1. Radial velocity analysis An ensemble of 184 high signal-to-noise observations (<S/N> of 110 at 550 nm) of HD 20003 have been gathered covering about 11 years (4063 days).

Fig. 3 .

 3 Fig. 3. GLS Periodogram of the RV residuals of HD 20003 obtained from the MCMC solution including a linear correlation between RV and log(R HK ) to account for the magnetic cycle, a second-order polynomial drift and three Keplerian.

Fig. 4 .

 4 Fig. 4. Phase-folded RV measurements of HD 20003 with the best Keplerian solution represented as a black curve for each of the signals in the data. From top-left to bottom-right: planet b, planet c and the signal at 184 days. The residuals around the solution are displayed in the lower-right panel. Contrarily to Fig. 1 error bars include now the instrumental + stellar jitter from the MCMC modeling. Corresponding planetary orbital elements are listed in Table5.

Fig. 5 .

 5 Fig. 5. Left: Activity index log(R HK ) as a function of the RV residuals when removing all the detected signals except the magnetic cycle effect for HD 20003. The observed correlation indicates that most of the RV residual variation is due to activity-related effects. Right: Barycentric Earth velocity projected on the line of sight (BERV) as a function of the RV residuals around the best derived solution without considering the 184-day signal.The fact that no correlation can be observed disfavours the hypothesis that the 184-day signal is due to a discontinuity in the wavelength calibration introduced by tiny gaps between the different quadrants of the detector.

Fig. 6 .

 6 Fig. 6. GLS periodogram of the residuals at each step, after removing one planet after the other in the analysis of HD 20781 (from left to right and top to bottom). The GLS periodogram of the raw RVs is shown in Fig. 1.

Fig. 7 .

 7 Fig. 7. Phase-folded RV measurements of HD 20781 with the best-fit solution represented as a black curve for each of the signal found in the data (from left to right and top to bottom: planet b, c, d and e). Error bars include photon and calibration noise, as well as instrumental + stellar jitters derived from the MCMC modeling. The residuals around the solution are displayed in the lower panel. Corresponding orbital elements are listed in Table6.

Fig. 8 .

 8 Fig. 8. Best keplerian solution for the eccentric planet orbiting HD 20782. Top: Phase-folded RV measurements with the best solution represented as a black curve. Bottom: RV residuals around the best Keplerian solution.

  Fig. 9. Top: GLS periodogram of the residual RVs of HD 21693 after removing the RV contribution of the magnetic cycle. Middle and bottom: GLS periodogram of the residuals at each step, after removing one planet after the other in the analysis. The GLS periodogram of the raw RVs is shown in Fig. 1.

Fig. 10 .

 10 Fig. 10. Phase-folded RV measurements of HD 21693 with the best planet solution represented as a black curve for each of the signal in the data (top to bottom: planet b and planet c). The residuals around the solution are displayed in the lower panel. Corresponding orbital elements are listed in Table7.

Fig. 11 .

 11 Fig. 11. RV residuals when removing all the detected signals except the magnetic cycle effect plotted as a function of the activity index log(R HK ) for HD 21693. The observed correlation indicates that most of the RV residual variation is due to activity-related effects.

Fig. 12 .

 12 Fig.12. Top: GLS periodogram of the residual RVs of HD 31527 after removing the best-fit Keplerian to account for the significant signal at 17 days seen in the raw RVs (bottom right panel of Fig.1). Middle and bottom: GLS periodogram of the residuals at each step, after removing the second and third planet from the RVs.

Fig. 13 .

 13 Fig. 13. Phase-folded RV measurements of HD 31527 with the best Keplerian solution for planet b, c, d and e represented as black curves (from left right and top to bottom). Error bars include photon and calibration noise, as well as a jitter estimated from the MCMC analysis. The residuals around the solution are displayed in the bottom right panel. Corresponding orbital elements are listed in Table8.

  1 matched signals in the different activity indicators (see Fig. A.5), and therefore those signal are associated with bona-fide planets.

Fig. 15 .

 15 Fig. 15. Phase-folded RV measurements of HD 45184 with, from top to bottom, the best Keplerian solution for planets b and c, and the corresponding residuals. Error bars include photon and calibration noise, as well as a jitter effect (stellar + instrumental) determined in the analysis. Corresponding orbital elements are listed in Table9.

  illustrated in Fig. 18, along with the RV residuals. The two signals detected in RVs are not matching any significant signal in the different activity indicators (see Fig. A.6) and are therefore associated to bona-fide planets.

Fig. 18 .

 18 Fig. 18. Phase-folded RV measurements of HD 51608 with, from top to bottom the best Keplerian solution for planet b, c and the residuals around the best-fitted solution. Corresponding orbital elements are listed in Table

  Fig. 20. From top to bottom: GLS periodogram of the residuals at each step, after removing one planet after the other in the analysis of HD 134060. The GLS periodogram of the raw RVs is shown in Fig. 1.

Fig. 21 .

 21 Fig. 21. Phase-folded RV measurements of HD 134060 with the best fitted solution for planet b and c represented as a black curve. The RV around the best-fit solution are displayed in the lower panel. Corresponding orbital elements are listed in Table11.

  Fig. 23. Phase-folded RV measurements of HD 136352 with, from left to right and top to bottom the best fitted solution for planet b, c and d. The residuals around the solution are displayed in the lower right panel. Corresponding orbital elements are listed in Table12.

Fig. 24 .

 24 Fig. 24. Summary of the confirmed planetary systems with indication of their periods (x-axis), masses (symbol size) and eccentricities (colour code).

Fig

  Fig. A.2. Same as Fig. A.1 but for HD20781

Table 1 .

 1 Observed and inferred stellar parameters for the stars hosting planetary systems described in this paper.

	Parameters	HD 20003	HD 20781 a	HD 21693	HD 31527	HD 45184	HD 51608	HD 134060	HD 136352
	Sp. Type (1)	G8V	K0V	G8V	G2V	G1.5V	G7V	G3IV	G4V
	V (1)	8.39	8.48	7.95	7.49	6.37	8.17	6.29	5.65
	B -V (1)	0.77	0.82	0.76	0.61	0.62	0.77	0.62	0.63
	π [mas] (1)	22.83±0.65	28.27±1.08	30.88±0.49	25.93±0.60	45.70±0.40	28.71±0.51	41.32±0.45	67.51±0.39
	M V	(1)	5.22	5.70	5.40	4.87	4.65	5.51	4.38	4.83
	T eff [K] (2)	5494±27	5256±29	5430±26	5898 ±13	5869±14	5358±22	5966±14	5664±14
	[Fe/H] (2)	0.04 ± 0.02 -0.11 ± 0.02	0.0 ± 0.02	-0.17 ± 0.01 0.04 ± 0.01 -0.07 ± 0.01 0.14 ± 0.01 -0.34 ± 0.01
	log (g) (2)	4.41±0.05	4.37±0.05	4.37 ±0.04	4.45±0.02	4.47±0.02	4.36±0.05	4.43 ±0.03	4.39±0.02
	M [M ] (2)	0.875	0.70	0.80	0.96	1.03	0.80	1.095	0.81
	L [L ] (2)	0.72±0.03	0.49±0.04	0.62±0.02	1.20±0.03	1.13±0.01	0.57±0.02	1.44±0.02	0.99±0.01
	log(R HK ) (3)	-4.97±0.05	-5.03±0.01	-4.91±0.05	-4.96±0.01	-4.91±0.01	-4.98±0.02	-5.00±0.01	-4.95±0.01
	v sin (i) [km s -1 ] (3)	1.9	1.1	1.6	1.9	2.1	1.3	2.6	< 1
	P rot [days] (4)	38.9±4.0	46.8±4.4	35.2±4.0	20.3±2.9	21.5±3.0	40.0±4.0	23.0±2.9	23.8±3.1

Notes. (1) Astrometric and visual photometric data from the Hipparcos Catalogs

(ESA 1997;van Leeuwen 2007)

. (

2

) From

Sousa et al. (2008)

spectroscopic analysis. (3) Parameter derived using HARPS spectra or CCF. (4) From the calibration of

Mamajek & Hillenbrand (2008)

. a HD 20782, the stellar visual companion of HD 20781, hosts a planet as well. The corresponding stellar parameters are given in

Jones et al. (2006)

.

Table 2 .

 2 General

		N spectr	< SNR >	N meas	∆T	σ RV	< ε RV >
			550 nm	1-hour bin	[days]	[m/s]	[m/s]
	HD 20003	184	110	184	4063	5.35	0.74
	HD 20781	225	112	216	4093	3.41	0.76
	HD 20782 a	71	181	68	4111	36.86	0.56
	HD 21693	210	141	210	4106	4.72	0.60
	HD 31527	256	180	245	4135	3.19	0.64
	HD 45184	308	221	178	4160	4.72	0.41
	HD 51608	218	133	216	4158	4.07	0.62
	HD 134060	335	199	155	4083	3.68	0.40
	HD 136352	649	231	240	3993	2.74	0.33

statistics of the HARPS observations of the planethost stars presented in this paper, with the number of individual spectra observed, <S/N> the average signal-to-noise ratio at 550 nm of those spectra, the number of measurements obtained after binning the data over 1 hour, ∆T the time span of the observations, σ RV the rms value of the RV set and < ε RV > the mean RV photon-noise + calibration error.

Table 3 .

 3 HARPS RVs and parameters inferred from the spectra and cross correlation functions for the 9 planet-host stars discussed in the paper.

	JDB [-2400000 days]	RV [km/s]	ε RV [m/s]	FWHM Contrast [km/s] [%]	BIS [m/s]	log(R HK ) SN50
	HD 20003 53728.576344	-16.09631 0.00054 6.88881	52.452	-0.03595	-4.8864	141.40
	HD 20003 53763.545747	-16.09548 0.00048 6.88610	52.481	-0.03718	-4.8729	165.50
	HD 20003 53788.524407	-16.09975 0.00052 6.88607	52.510	-0.03588	-4.8968	151.20
	...
	the gaps between pixel patches due to the yearly motion of the
	Earth around the Sun

Table 4 .

 4 Parameters probed by the MCMC. The symbols U and N used for the priors definition stands for uniform and normal distributions, respectively.

	Parameters	Units	Priors	Description

  360[Parameters probed by MCMC with magnetic cycle (in addition to the previous ones, except σ JIT LOW and σ JIT HIGH that replace σ JIT

			Udry et al.: 20 super-Earths and hot Neptunes detected with HARPS
	σ JIT LOW	[ m s -1 ]	U	stellar jitter at the minimum of the magnetic cycle (If only one
				instrument is used, this parameter includes σ X )	]-∞,∞[
	σ JIT HIGH	[ m s -1 ]	U	stellar jitter at the maximum of the magnetic cycle (If only one
				instrument is used, this parameter includes σ X )	]-∞,∞[
	log(R HK ) lin	[ m s -1 ]	U	slope of the correlation between the RV and log(R HK ) (The log(R HK ) variation is normalized between 0 and 1, thus this parameter is in [ m s -1 ])	]-∞,∞[
	Physical Parameters derived from the MCMC posteriors (not probed)
	P	[d]	-	Period
	K	[m s -1 ]	-	RV semi-amplitude
	e	-	-	Orbit eccentricity
	ω	[deg]	-	Argument of periastron
	T P	[d]	-	Time of passage at periastron
	T C	[d]	-	Time of transit
	Ar	[AU]	-	Semi-major axis of the relative orbit
	M.sin i	[M Jup ]	-	Mass relative to Jupiter
	M.sin i	[M Earth ]	-	Mass relative the Earth
	the GLS periodogram of the calcium activity index log(R HK ), we fit the RVs with a polynomial up to the second order, plus	value is smaller than 1 % 4 . If the p-values of some signals in the GLS periodogram are smaller than 1%, we adjust a Ke-
	a Keplerian model reproducing the log(R HK ) variation, leaving only the amplitude as a free parameter. Note that when correct-	plerian orbital solution to the signal presenting the smallest p-value following the approach of Delisle et al. (2016) with
	ing magnetic cycle effect in RVs when adopting a step-by-step	the detected period as a guess value.
	approach, we do not consider a simple linear correlation between RV and log(R HK ) as explained above because significant plane-tary signals not yet removed from the data might destroy any existing correlation. There is no a priori reason why the mag-netic cycle should look like a Keplerian, however a Keplerien has more degrees of freedom than a simple sinusoid and can there-fore better estimate the long-term variation seen in log(R HK ). Regarding signal significance, a signal is considered worth look-ing into when its p-value, which gives the probability that the	-Additional planets in the systems are then considered if other significant signals are present in the GLS periodogram of the RV residuals. Note that each time a Keplerian signal is added to the model, a global fit including all the previously detected signals is performed. -We stopped when no more signals in the residuals present a p-value smaller than 1 % and we finally visually inspected the solutions.
	signal appears just by chance, is smaller or equal to a chosen
	threshold. Experience has shown us that a p-value of 1 % pro-	4.3. Origin of the different signals found
	vides a good empirical threshold. Note that in this paper, p-values are estimated in a Monte Carlo approach by bootstrapping 10'000 times over the dates of the observations.	After the main periodic signals have been recognized in the data, it is important to identify the ones that are very likely not from planetary origin. They may be of different natures: stellar, in-
					strumental or observational. Here is a brief description of some
	4.2. Detection of periodic signals in the data	cases encountered:
	-Planetary signals are searched for in the GLS periodogram
	of the RV residuals after correcting long-term trends. As for
	magnetic cycles, we consider a signal as significant if its p-

Table 6 .

 6 Best-fitted solution for the planetary system orbiting HD 20781 and HD 20782. For each parameter, the median of the posterior distribution is considered, with error bars computed from the MCMC posteriors using a 68.3 % confidence interval. The value σ (O-C) corresponds to the weighted standard deviation of the residuals around the best solution. All the parameters probed by the MCMC can be found in Annex, in Tables B.2 and B.3. See Table4for definition of the parameters.

	Param.	Units	HD20781b	HD20781c	HD20781d	HD20781e	HD20782b
	P	[d]	5.3135 +0.0010 -0.0010	13.8905 +0.0033 -0.0034	29.1580 +0.0102 -0.0100	85.5073 +0.0983 -0.0947	597.0643 +0.0256 -0.0256
	K	[m s -1 ]	0.91 +0.15 -0.15	1.81 +0.16 -0.16	2.82 +0.17 -0.16	2.60 +0.14 -0.14	118.43 +1.78 -1.78
	e		0.10 +0.11 -0.07	0.09 +0.09 -0.06	0.11 +0.05 -0.06	0.06 +0.06 -0.04	0.95 +0.001 -0.001
	ω	[deg]	84.04 +141.70 -108.41	7.44 +53.87 -70.21	60.99 +30.79 -30.03	70.59 +61.40 -67.58	143.58 +0.56 -0.66
	T P	[d]	55503.2027 +2.0979 -1.5934	55503.5204 +2.0425 -2.6815	55511.3258 +2.4394 -2.4382	55513.3912 +14.3623 -16.0090	55247.0150 +0.0770 -0.0837
	T C	[d]	55503.2888 +0.2361 -0.2384	55506.4386 +0.3546 -0.4385	55513.2432 +0.4250 -0.4410	55517.4714 +1.3596 -1.4838	55246.1712 +0.0706 -0.0706
	Ar	[AU]	0.0529 +0.0024 -0.0027	0.1004 +0.0046 -0.0051	0.1647 +0.0076 -0.0083	0.3374 +0.0155 -0.0170	1.3649 +0.0466 -0.0495
	M.sin i	[M Jup ]	0.0061 +0.0012 -0.0011	0.0168 +0.0022 -0.0021	0.0334 +0.0038 -0.0037	0.0442 +0.0049 -0.0049	1.4878 +0.1045 -0.1066
	M.sin i	[M Earth ]	1.93 +0.39 -0.36	5.33 +0.70 -0.67	10.61 +1.20 -1.19	14.03 +1.56 -1.56	472.83 +33.22 -33.88
	offset UCLES [m s -1 ]					5.4431 +0.9302 -0.9370
	γ COR98	[m s -1 ]					39928.1039 +1.8747 -1.9160
	γ COR07	[m s -1 ]					39930.6146 +2.0828 -2.1696
	γ COR14	[m s -1 ]					39956.5569 +1.5850 -1.6638
	γ HARPS	[m s -1 ]		40369.2080 +0.1147 -0.1104		39964.8070 +0.2193 -0.2275
	σ (O-C)	[m s -1 ]		1.45		2.34
	log (Post)			-397.1395 +3.0943 -3.8525		-379.6596 +2.6521 -3.2805

  .4.

	Param.	Units	HD21693b	HD21693c
	P	[d]	22.6786 +0.0085 -0.0087	53.7357 +0.0312 -0.0309
	K	[m s -1 ]	2.20 +0.22 -0.22	3.44 +0.20 -0.20
	e		0.12 +0.09 -0.08	0.07 +0.06 -0.05
	ω	[deg]	-91.04 +50.57 -50.35	-17.34 +47.50 -55.35
	T P	[d]	55492.0549 +3.0968 -3.1310	55528.4064 +7.0762 -8.1820
	T C	[d]	55480.7554 +0.7454 -0.6973	55543.5822 +1.0676 -1.2491
	Ar	[AU]	0.1455 +0.0058 -0.0063	0.2586 +0.0103 -0.0113
	M.sin i	[M Jup ]	0.0259 +0.0034 -0.0033	0.0547 +0.0056 -0.0056
	M.sin i	[M Earth ]	8.23 +1.08 -1.05	17.37 +1.77 -1.79
	γ HARPS	[m s -1 ]	39768.8113 +0.1425 -0.1471
	σ (O-C)	[m s -1 ]	2.05
	log (Post)		-440.9160 +2.2958 -3.1188

  .5.

	Param.	Units	HD31527b	HD31527c	HD31527d
	P	[d]	16.5535 +0.0034 -0.0035	51.2053 +0.0373 -0.0368	271.6737 +2.1135 -2.2471
	K	[m s -1 ]	2.72 +0.13 -0.13	2.51 +0.14 -0.14	1.25 +0.17 -0.16
	e		0.10 +0.05 -0.05	0.04 +0.05 -0.03	0.24 +0.13 -0.13
	ω	[deg]	41.12 +29.46 -34.41	-23.25 +68.94 -152.43	179.00 +31.11 -26.20
	T P	[d]	55499.5453 +1.3130 -1.5818	55526.3434 +9.7876 -21.5943	55718.8091 +21.3880 -17.3759
	T C	[d]	55501.4585 +0.2640 -0.2667	55542.0635 +0.7471 -0.9373	55670.8324 +11.7715 -13.7509
	Ar	[AU]	0.1254 +0.0041 -0.0045	0.2663 +0.0088 -0.0095	0.8098 +0.0273 -0.0293
	M.sin i	[M Jup ]	0.0329 +0.0028 -0.0028	0.0445 +0.0040 -0.0039	0.0372 +0.0053 -0.0052
	M.sin i	[M Earth ]	10.47 +0.89 -0.87	14.16 +1.28 -1.23	11.82 +1.70 -1.64
	γ HARPS	[m s -1 ]		25739.7025 +0.0952 -0.0971
	σ (O-C)	[m s -1 ]		1.41
	log (Post)			-439.4138 +2.7350 -3.3929

  .6.

	Param.	Units	HD45184b	HD45184c
	P	[d]	5.8854 +0.0003 -0.0003	13.1354 +0.0026 -0.0025
	K	[m s -1 ]	4.26 +0.23 -0.23	2.36 +0.23 -0.23
	e		0.07 +0.05 -0.05	0.07 +0.07 -0.05
	ω	[deg]	145.80 +49.43 -47.76	-197.97	-80.87
	T P	[d]	55500.2509 +0.7996 -0.7790	55497.4412 +4.4543 -2.9188
	T C	[d]	55499.4150 +0.1050 -0.0903	55494.8821 +0.3364 -0.3065
	Ar	[AU]	0.0644 +0.0020 -0.0021	0.1100 +0.0034 -0.0036
	M.sin i	[M Jup ]	0.0384 +0.0033 -0.0032	0.0277 +0.0034 -0.0032
	M.sin i	[M Earth ]	12.19 +1.06 -1.03	8.81 +1.09 -1.02
	γ HARPS	[m s -1 ]	-3757.6506 +0.1562 -0.1595
	σ (O-C)	[m s -1 ]	2.15
	log (Post)		-382.5338 +2.3828 -2.9637

  .7.

	Param.	Units	HD51608b	HD51608c
	P	[d]	14.0726 +0.0016 -0.0016	95.9446 +0.1555 -0.1366
	K	[m s -1 ]	3.95 +0.16 -0.16	2.36 +0.17 -0.17
	e		0.09 +0.04 -0.04	0.14 +0.07 -0.07
	ω	[deg]	117.45 +30.27 -28.91	-165.07 +29.34 -33.93
	T P	[d]	55494.5239 +1.1795 -1.1073	55498.3725 +7.6262 -9.0495
	T C	[d]	55493.6039 +0.1684 -0.1646	55474.0824 +2.5542 -2.5650
	Ar	[AU]	0.1059 +0.0043 -0.0046	0.3809 +0.0153 -0.0164
	M.sin i	[M Jup ]	0.0402 +0.0038 -0.0037	0.0450 +0.0051 -0.0048
	M.sin i	[M Earth ]	12.77 +1.20 -1.19	14.31 +1.63 -1.53
	γ HARPS	[m s -1 ]	39977.2351 +0.1159 -0.1147

  Table B.8. See Table 4 for definition of the parameters.

	Param.		HD134060b	HD134060c
	P	[d]	3.2696 +0.0001 -0.0001	1291.5646 +48.0333 -44.2197
	K	[m s -1 ]	4.61 +0.22 -0.22	1.65 +0.24 -0.23
	e		0.45 +0.04 -0.04	0.11 +0.13 -0.07
	ω	[deg]	-98.23 +6.61 -6.75	-132.73 +121.28 -55.51
	T P	[d]	55499.6542 +0.0420 -0.0407	55232.8160 +444.3209 -186.5492
	T C	[d]	55498.1943 +0.1065 -0.1037	56057.2704 +61.2663 -48.5676
	Ar	[AU]	0.0444 +0.0013 -0.0014	2.3928 +0.0929 -0.0951
	M.sin i	[M Jup ]	0.0318 +0.0025 -0.0024	0.0922 +0.0139 -0.0133
	M.sin i	[M Earth ]	10.10 +0.79 -0.75	29.29 +4.43 -4.24
	γ HARPS	[m s -1 ]	37987.9484 +0.1512 -0.1489
	σ (O-C)	[m s -1 ]		1.64
	log (Post)		-304.2166 +2.1293 -2.9181

  for definition of the parameters.

	Param.	Units	HD136352b	HD136352c	HD136352d
	P	[d]	11.5824 +0.0024 -0.0025	27.5821 +0.0089 -0.0086	107.5983 +0.2796 -0.2669
	K	[m s -1 ]	1.59 +0.13 -0.13	2.65 +0.14 -0.14	1.35 +0.15 -0.15
	e		0.14 +0.08 -0.08	0.04 +0.05 -0.03	0.09 +0.10 -0.07
	ω	[deg]	-185.64 +35.97 -37.67	-110.22 +74.58 -75.76	170.89 +63.23 -76.56
	T P	[d]	55496.7009 +1.1393 -1.2084	55490.0572 +5.7741 -5.8036	55494.9775 +18.6996 -22.9117
	T C	[d]	55494.3913 +0.3534 -0.3259	55477.8619 +0.4105 -0.3730	55472.4695 +4.6515 -3.2834
	Ar	[AU]	0.0934 +0.0037 -0.0040	0.1666 +0.0066 -0.0072	0.4128 +0.0163 -0.0175
	M.sin i	[M Jup ]	0.0151 +0.0018 -0.0018	0.0340 +0.0034 -0.0033	0.0270 +0.0037 -0.0036
	M.sin i	[M Earth ]	4.81 +0.59 -0.56	10.80 +1.08 -1.05	8.58 +1.18 -1.15
	γ HARPS	[m s -1 ]		-68709.0314 +0.0937 -0.0964
	σ (O-C)	[m s -1 ]		1.35
	log (Post)			-420.5739 +2.5027 -3.1788

Table B .

 B 1. Parameters probed by the MCMC used to fit the RV measurements of HD20003. The maximum likelihood solution (Max(Like)), the median (Med), mode (Mod) and standard deviation (Std) of the posterior distribution for each parameter is provided, as well as the 68.3 % (CI(15.85),CI(84.15)) and 95.45 % (CI(2.275),CI(97.725)) confidence intervals. The prior for each parameter can be of type: U: uniform, N: normal, SN:split normal, T N:truncated normal Article number, page 28 of 31 Udry et al.: 20 super-Earths and hot Neptunes detected with HARPSTable B.4. Same as Table B.1 for HD 21693 Table B.5. Same as Table B.1 for HD 31527 Article number, page 29 of 31 A&A proofs: manuscript no. udry2018 Table B.6. Same as Table B.1 for HD 45184 Table B.7. Same as Table B.1 for HD 51608 Article number, page 30 of 31 Udry et al.: 20 super-Earths and hot Neptunes detected with HARPS Table B.8. Same as Table B.1 for HD 134060 Table B.9. Same as Table B.1 for HD 136352

													A&A proofs: manuscript no. udry2018
	Table B.2. Same as Table B.1 for HD 20781
		Param.	Units			Max(Like)	Med	Mod	Std	CI(15.85)	CI(84.15)	CI(2.275)	CI(97.725)	Prior
				Param. Param. Param.				Units Units Units		Max(Like) Max(Like) Max(Like)	Med Med Med	Mod Mod Mod	Std Std Std Likelihood	CI(15.85) CI(15.85) CI(15.85)	CI(84.15) CI(84.15) CI(84.15)	CI(2.275) CI(2.275) CI(2.275)	CI(97.725) CI(97.725) CI(97.725)	Prior Prior Prior
	.	log (Post)					-298.850697	-304.216566	-304.307005	Likelihood Likelihood Likelihood 2.241852 -307.134697	-302.087246	-310.611617	-300.536513
		log (Post) log (Post) log (Post) log (Like)				-388.542397 -435.045817 -375.948188 -297.620659 -302.893646 -397.139540 -440.915985 -382.533843 -302.472853 -396.801832 -440.662807 -381.953342 2.232698 3.086207 2.436011 2.376383 -305.807282 -400.992042 -444.034771 -385.497577 -300.789692 -394.045214 -438.620178 -380.151077 -309.374131 -405.822507 -447.941122 -389.208051 -299.237149 -391.632479 -436.845354 -378.372729
		Param. log (Like) log (Like) log (Like) log (Prior)	Units	Max(Like) -388.379068 -434.994431 -375.907261 -1.230038 -1.253327 Med -396.816561 -440.747390 -382.434308 -1.227765 Mod -396.641259 -440.523392 -381.860465 0.254368 Std 3.045776 2.426114 2.358405 -1.556962 CI(15.85) -400.618909 -443.837350 -385.390539 -1.023282 CI(84.15) -393.744957 -438.463925 -380.066476 -2.140093 CI(2.275) -405.246050 -447.707065 -389.064124 -0.820972 CI(97.725) -391.364866 -436.680607 -378.293195	Prior
		M	log (Prior) log (Prior) log (Prior) [M ]			-0.163328 -0.051386 -0.040927 1.099460 1.094932 -0.271536 -0.133026 -0.076679	Likelihood 0.188786 0.126005 0.074946 0.088166 -0.192668 -0.077867 -0.029834 1.086339 0.993140 -0.532565 -0.308420 -0.179107 1.195814 -0.133491 -0.041395 -0.022112 0.898980 -0.959781 -0.592687 -0.352708 1.298334 -0.053132 -0.007046 -0.003342 N(1.095, 0.1)
		σ JIT	log (Post) M M M [m s -1 ]	[M ] [M ] [M ]	-355.440573 0.798309 0.627847 1.058757 1.58 1.65	-364.186739 0.700500 0.799037 1.029104	-365.130519 0.675260 0.801958 1.010914 1.63	3.064118 0.088466 0.087416 0.087356 0.22	-367.902792 0.599877 0.699111 0.930545 1.43	-360.869032 0.801716 0.898724 1.127460 1.35	-372.251620 0.500287 0.600636 0.830397 1.31	-358.366379 0.905060 0.995794 1.231321 1.04	N(0.7, 0.1) N(0.8, 0.1) N(1.03, 0.1) U
		log (Like) σ JIT σ JIT LOW σ JIT LOW γ HARPS [m s -1 ] [m s -1 ] [m s -1 ] [m s -1 ] 37987.94 -354.388255 1.14 1.27 1.38 37987.95 -363.200339 1.27 1.26 1.76	-364.552124 1.21 1.25 1.90 37987.91	3.075716 0.19 0.44 0.64 0.13	-366.886986 1.09 0.73 0.73 37987.80	-359.786081 1.01 0.79 1.27 37988.10	-371.228449 1.00 0.14 0.11 37987.65	-357.313993 0.69 0.38 0.83 37988.25	U	U U U
		log (Prior) γ HARPS σ JIT HIGH σ JIT HIGH log (P)	[d]	-1.052318 40369.22 2.69 2.33 0.514500 [m s -1 ] [m s -1 ] [m s -1 ] 0.514500 -1.000952 40369.21 2.56 2.74	-1.023733 40369.17 2.70 2.77 0.514499 0.000011 0.258084 0.10 0.56 0.68	-1.317939 40369.10 1.86 1.70 0.514488	-0.740211 40369.32 2.06 2.17 0.514512	-1.736133 40368.98 0.85 0.93 0.514478	-0.528578 40369.44 1.60 1.62 0.514526	U	U U U
		M log (P) γ HARPS γ HARPS log (K)	[M ] [d] [m s -1 ] [m s -1 ] [m s -1 ]	0.947719 0.725379 39768.81 -3757.68 0.67 0.66	0.877722 0.725377 39768.81 -3757.65	0.888341 0.725315 39768.82 -3757.69 0.66	0.086110 0.000076 0.13 0.14 0.02	0.777570 0.725292 39768.66 -3757.81 0.64	0.972946 0.725459 39768.95 -3757.49 0.68	0.677797 0.725196 39768.52 -3757.98 0.62	1.077890 0.725567 39769.09 -3757.34 0.70	N(0.875, 0.1) U U U U
		σ JIT γ HARPS log (P) log (K) √ e. cos ω log (K) RHK index lin RHK index lin e. cos ω √ e. cos ω log (P) log (P) e. sin ω √ e. sin ω λ 0 log (P) log (K) √ e. cos ω √ e. sin ω log (K) λ 0 √ √ [deg] [m s -1 ] [m s -1 ] [d] [m s -1 ] [m s -1 ] [ m s -1 ] [ m s -1 ] -0.125031 1.45 -16104.48 1.073639 0.58 -0.021425 0.05 11.066005 10.723380 -0.096070 1.59 -16104.11 1.073654 0.58 -0.028442 -0.04 10.860734 10.634926 -0.087114 1.63 -16104.25 1.073651 0.58 -0.023510 -0.04 10.732702 10.630724 0.067809 0.34 0.31 0.000052 0.02 0.077999 0.07 0.440301 0.607624 -0.173750 1.17 -16104.39 1.073593 0.56 -0.117180 -0.12 10.361899 9.941268 -0.018918 1.24 -16103.74 1.073713 0.61 0.059704 0.03 11.369837 11.330657 -0.246029 0.77 -16104.65 1.073536 0.53 -0.209618 -0.22 9.877045 9.258677 0.063097 0.93 -16103.07 1.073771 0.63 0.142977 0.09 11.871009 11.996205 -0.061140 0.043862 0.039181 0.230069 -0.243492 0.299379 -0.460848 0.499573 [d] 1.355767 1.355617 1.355612 0.000146 1.355449 1.355780 1.355285 1.355943 [d] 0.769779 0.769779 0.769776 0.000020 0.769758 0.769802 0.769736 0.769826 -0.673348 -0.660554 -0.664313 0.027638 -0.690593 -0.628249 -0.720202 -0.591182 0.216991 -0.013862 -0.028798 0.215428 -0.264618 0.241274 -0.455895 0.436077 [deg] -135.345360 -134.251473 -137.352873 9.387901 -144.854378 -123.506806 -156.133446 -112.221519 [d] 1.142791 1.142717 1.142690 0.000092 1.142611 1.142819 1.142490 1.142923 [m s -1 ] 0.35 0.34 0.34 0.04 0.30 0.38 0.24 0.42 -0.033137 -0.003440 -0.067300 0.179167 -0.215702 0.205008 -0.392233 0.358369 -0.247432 -0.267481 -0.370194 0.162410 -0.411532 -0.042976 -0.522567 0.169768 [m s -1 ] 0.63 0.63 0.63 0.02 0.61 0.65 0.58 0.67 -60.971221 -60.219046 -60.857964 2.817351 -63.378221 -57.046679 -66.686867 -53.740675 √ e. cos ω -0.220553 -0.170838 -0.275286 0.126595 -0.294336 -0.000876 -0.378620 0.156389 log (P) [d] 3.115614 3.111116 3.110321 0.013896 3.095987 3.126974 3.080422 3.146266 √ e. sin ω 0.182613 0.117241 0.149405 0.131397 -0.053453 0.252324 -0.200167 0.349677 log (K) [m s -1 ] 0.19 0.22 0.22 0.05 0.15 0.27 0.08 0.33	U U U U U	U U U U U U U U U U U U U U U U U U U U
		√ λ 0 e. sin ω log (P) log (K) √ e. cos ω log (P) √ e. sin ω e. sin ω √ log (K) λ 0 λ 0 e. cos ω √ e. cos ω log (P) log (P) e. sin ω √ e. sin ω λ 0 log (K) √ e. cos ω log (K) e. cos ω √ λ 0 √ √ [deg] [deg] [d] [m s -1 ] [d] [m s -1 ] [deg] [deg] 0.018089 -0.596953 225.658516 1.530369 0.49 0.089332 1.464692 -0.025692 -0.059901 0.29 34.107432 -230.080798 -0.131288 -0.593533 225.988820 1.530505 0.50 0.163526 1.464758 -0.063953 0.065831 0.26 35.450889 -229.352706 -0.033631 -0.597847 225.701948 1.530577 0.50 0.213761 1.464755 -0.075276 0.087988 0.24 33.905458 -229.816814 0.237542 0.037382 3.096407 0.000282 0.03 0.146216 0.000134 0.137633 0.173593 0.03 4.533112 2.694368 -0.408103 -0.633941 222.316287 1.530165 0.47 -0.035613 1.464609 -0.218240 -0.146692 0.22 30.349272 -232.391963 0.151476 -0.548640 229.436249 1.530811 0.53 0.306443 1.464910 0.105252 0.260735 0.29 40.588204 -226.300698 -0.601802 -0.670370 218.679062 1.529813 0.44 -0.203105 1.464452 -0.335603 -0.312367 0.18 25.065784 -235.479746 0.365222 -0.497859 232.672381 1.531065 0.55 0.412917 1.465071 0.237976 0.404779 0.33 46.053381 -223.032360 0.327288 0.203741 0.321124 0.172122 -0.029464 0.370760 -0.219811 0.484086 [d] 1.730230 1.730263 1.730205 0.000221 1.730013 1.730515 1.729784 1.730800 [d] 1.118470 1.118443 1.118431 0.000075 1.118359 1.118527 1.118271 1.118619 -0.193191 -0.102419 -0.157397 0.182951 -0.296909 0.130699 -0.452167 0.321316 0.021983 0.026817 0.016131 0.159055 -0.169199 0.211170 -0.318653 0.343092 [deg] -84.722610 -83.737346 -84.982562 4.380473 -88.530458 -78.631060 -93.607328 -73.479892 [m s -1 ] 0.54 0.54 0.54 0.02 0.51 0.56 0.48 0.59 0.269425 0.196029 0.243736 0.136666 0.005735 0.324520 -0.157091 0.416385 [m s -1 ] 0.40 0.37 0.37 0.04 0.33 0.41 0.28 0.45 -0.147627 -0.073911 -0.085546 0.176247 -0.277303 0.140983 -0.420074 0.304579 -57.886308 -59.605407 -62.436000 6.753625 -66.940199 -51.615752 -74.094300 -42.949848	U U U	U U U U U U U U U U U U U U U U U U U U
		√ λ 0 e. sin ω log (K) λ 0 λ 0 √ e. cos ω Param. Units [deg] [m s -1 ] [deg] [deg] Max(Like) 0.339711 120.340507 0.45 -206.679521 -129.227327 0.139677 Med log (P) [d] 2.262470 √ e. sin ω 0.246274	0.219084 124.671746 0.45 -207.980760 -127.324650 0.145327 2.263926 0.263747	0.272743 124.146671 0.45 -207.807573 -126.147264 0.162099 Mod Std 0.133488 3.478923 0.02 2.957609 4.833219 0.120006 CI(15.85) 0.030981 120.635953 0.42 -211.285700 -132.789170 -0.010957 CI(84.15) 0.334583 128.547476 0.48 -204.615563 -121.833139 0.265346 CI(2.275) -0.157897 116.645239 0.40 -214.868823 -138.259570 -0.162692 CI(97.725) 0.421322 132.285391 0.50 -201.248502 -116.150675 0.358965 2.263963 0.002232 2.261303 2.266297 2.258300 2.268859 0.306374 0.113409 0.114332 0.358597 -0.077325 0.434823 Likelihood	Prior	U U U U U U U U
		log (K) λ 0 log (Post)		[m s -1 ] [deg] -413.840756 0.23 -77.213730 -420.573939 0.21 -78.248282 -420.519231 0.20 -78.600903 2.525979 0.05 2.996830 -423.752736 0.15 -81.624509 -418.071226 0.26 -74.896708 -427.592767 0.08 -84.983023 -416.139887 0.30 -71.294902	U U
		√ log (P) e. cos ω log (Like)			[d]	0.191246 1.931848 -413.418243 -420.301368 0.167121 1.932003	0.209789 1.932093 -420.080839 2.517934 0.189206 0.000430 -423.494074 -0.078520 1.931522 -417.796941 0.364574 1.932502 -427.287283 -0.282064 1.931041 -415.883659 0.506428 1.932995	U U
	√ λ 0 e. sin ω log (P) log (K) log (K) Param. Param. log (Prior) √ e. cos ω √ e. sin ω λ 0 log (Post) log (Like) M [M ] [deg] [d] [m s -1 ] [m s -1 ] Units Units -0.422514 0.419103 197.675128 3.582659 0.74 0.43 Max(Like) Max(Like) -0.225282 0.213797 196.751843 3.606144 0.76 0.41 Med Med 0.148239 0.050362 0.268142 0.155008 [deg] 12.584852 15.024637 -431.457393 -439.413776 -0.181747 0.320208 196.887430 3.595842 0.75 0.41 Mod Mod 0.168419 0.195134 6.258785 0.028188 0.02 0.02 Std Std 0.090292 0.146393 -0.456156 -0.049677 189.645165 3.580559 0.73 0.39 CI(15.85) CI(15.85) -0.087916 0.402924 203.799985 3.640255 0.78 0.44 CI(84.15) CI(84.15) -0.795628 -0.274331 182.342310 3.557643 0.71 0.37 CI(2.275) CI(2.275) -0.025262 0.537274 211.192757 3.702132 0.80 0.46 CI(97.725) CI(97.725) -0.126620 0.219662 -0.260249 0.340229 0.226421 0.139655 -0.035773 0.293236 -0.188818 0.387469 14.339843 3.408448 11.079349 18.883891 7.424984 22.814245 Likelihood -439.404479 2.693424 -442.806691 -436.678805 -446.790781 -434.536716 -430.869912 -438.920421 -438.066741 2.744421 -442.319982 -436.104557 -446.392727 -433.960853 0.836645 0.811051 0.810685 0.088144 0.710931 0.910337 0.609205 1.008754 N(0.81, 0.1) U U U U U Prior Prior U U U Likelihood log (Post) -402.833168 -409.139845 -408.668989 2.465262 -412.297050 -406.697560 -415.940872 σ JIT [m s -1 ] 1.43 1.41 1.41 0.30 1.02 1.08 0.76 0.81 U -404.885802 log (Like) -402.629727 -408.960474 -408.454977 2.467315 -412.097485 -406.520548 -415.757993 γ HARPS [m s -1 ] -68709.12 -68709.03 -68709.05 0.08 -68709.13 -68708.94 -68709.22 -68708.84 U -404.678963 √ e. cos ω 0.028669 -0.069460 -0.042689 0.161723 -0.257669 0.126252 -0.420089 0.262372 U √ e. sin ω -0.205147 -0.146710 -0.208581 0.122078 -0.262944 0.019978 -0.347443 0.162734 U λ 0 [deg] 137.392067 139.338913 137.305526 3.530496 135.946353 143.669909 132.569853 149.602355 U Table B.3. Same as Table B.1 for HD 20782 Param. Units Max(Like) Med Mod Std CI(15.85) CI(84.15) CI(2.275) CI(97.725) Prior Likelihood log (Post) -372.933997 -379.659609 -379.184909 2.652630 -382.940155 -377.007478 -387.071255 -375.024009 log (Like) -367.930148 -374.651836 -374.180899 2.652279 -377.934595 -371.997235 -382.071125 -370.018358 log (Prior) -5.003849 -5.008809 -5.011567 0.009534 -5.019713 -4.997733 -5.030582 -4.987338 M [M ] 0.985227 0.950073 0.916830 0.087670 0.850317 1.050787 0.755244 1.148009 N(0.95, 0.1) σ COR07-DRS -3-4 [m s -1 ] 5.33 6.45 5.20 1.77 4.84 8.72 3.61 12.03 U σ COR14-DRS -3-8 [m s -1 ] 4.32 4.66 4.01 1.59 3.08 6.58 1.49 9.29 U σ COR98-DRS -3-3 [m s -1 ] 5.26 5.88 5.33 1.74 4.14 7.98 2.33 10.72 U σ HARPS [m s -1 ] 1.42 1.49 1.47 0.21 1.25 1.71 1.00 1.95 U σ UCLES -PU B-2006 [m s -1 ] 0.06 1.60 0.00 1.05 0.49 2.95 0.07 4.41 U σ JIT [m s -1 ] 1.29 0.94 1.00 0.46 0.37 1.47 0.06 1.93 U log (Prior) -0.587481 -0.421477 -0.228174 0.322560 -0.872877 -0.155258 -1.460968 -0.044428 M [M ] 0.964548 0.960466 0.935623 0.087688 0.861241 1.058519 0.761776 1.159099 N(0.96, 0.1) σ JIT [m s -1 ] 1.24 1.29 1.29 0.30 0.95 0.93 0.77 0.62 U γ HARPS [m s -1 ] 25739.68 25739.70 25739.69 0.09 25739.61 25739.80 25739.51 25739.90 U log (P) [d] 1.218857 1.218889 1.218885 0.000080 1.218798 1.218979 1.218706 1.219071 U log (K) [m s -1 ] 0.43 0.43 0.43 0.02 0.41 0.46 0.39 0.48 U √ e. cos ω 0.242808 0.210416 0.235029 0.104713 0.075031 0.311766 -0.082660 0.383529 log (Prior) -0.203442 -0.163596 -0.147311 0.100964 -0.296183 -0.073908 -0.485685 log (P) [d] 1.063786 1.063799 1.063829 0.000081 1.063703 1.063890 1.063608 1.063979 U -0.023113 M [M ] 0.980358 0.800314 0.788707 0.088192 0.701466 0.900776 0.600351 0.999712 log (K) [m s -1 ] 0.23 0.20 0.19 0.03 0.16 0.24 0.12 0.27 U N(0.8, 0.1) σ JIT LOW [m s -1 ] 1.20 1.24 1.21 0.42 0.76 0.78 0.16 0.39 √ e. cos ω -0.342695 -0.314967 -0.410810 0.145521 -0.437347 -0.116554 -0.524238 0.125811 U U σ JIT HIGH [m s -1 ] 2.04 1.94 1.99 0.52 1.30 1.49 0.31 U √ e. sin ω -0.090658 0.033773 0.091648 0.156779 -0.152919 0.214589 -0.306061 0.352353 U 1.06 [m s -1 ] 39977.30 39977.24 39977.20 0.10 39977.12 39977.35 39977.00 39977.47 U λ 0 [deg] 277.740369 277.024300 275.232882 4.443846 272.023532 282.083700 266.558482 287.198688 U γ HARPS RHK index lin [ m s -1 ] 3.902459 4.129760 4.057943 0.501237 3.567302 4.713965 2.972341 5.281426 log (P) [d] 1.440604 1.440627 1.440596 0.000121 1.440492 1.440768 1.440359 U 1.440909 U log (K) [m s -1 ] 0.43 0.42 0.42 0.02 0.40 0.45 0.37 U 0.47 U √ e. sin ω 0.236848 0.188237 0.211870 0.120112 0.021413 0.298912 -0.126407 0.384773 U λ 0 [deg] 52.547369 51.262624 51.323989 2.622456 48.258036 54.220352 45.382590 57.230513 U log (P) [d] 1.709426 1.709315 1.709215 0.000282 1.709003 1.709631 1.708664 1.709963 U log (K) [m s -1 ] 0.40 0.40 0.40 0.02 0.37 0.42 0.35 0.45 U √ e. cos ω 0.023117 0.081468 0.086880 0.144497 -0.097212 0.242734 -0.237589 0.358056 log (P) [d] 1.148361 1.148375 1.148387 0.000043 1.148326 1.148424 1.148275 1.148474 U log (K) √ e. cos ω -0.070645 -0.042094 -0.084481 0.126899 -0.182267 0.113998 -0.296099 0.237968 U [m s -1 ] 0.61 0.60 0.60 0.02 0.58 0.61 0.56 0.63 U √ √ e. sin ω -0.212716 -0.113587 -0.193076 0.134373 -0.259032 0.060329 -0.348573 U 0.197605 e. cos ω -0.118480 -0.123787 -0.183326 0.102911 -0.229128 0.007749 -0.308200 0.133668 U √ e. sin ω 0.335597 0.239601 0.269015 0.095793 0.113091 0.322123 -0.059036 U [deg] 19.705439 19.919436 19.941697 2.590743 16.987116 22.886544 14.008324 25.935025 U λ 0 0.387245 λ 0 [deg] 257.734884 257.274281 256.797681 2.077620 254.933378 259.655465 252.551371 262.119813 log (P) [d] 2.030837 2.031806 2.031583 0.000968 2.030727 2.032933 2.029701 U 2.034155 U log (K) [m s -1 ] 0.12 0.13 0.12 0.04 0.08 0.17 0.03 U 0.22 U √ e. sin ω -0.164323 0.014834 -0.010018 0.126248 -0.134951 0.164944 -0.248677 0.274415 U λ 0 [deg] -208.406874 -207.583870 -207.914089 2.785788 -210.745542 -204.465797 -213.962526 -201.198653 U log (P) [d] 2.434889 2.434048 2.433252 0.003107 2.430440 2.437413 2.426290 2.440704 U log (P) [d] 1.981979 1.982021 1.981746 0.000578 1.981402 1.982724 1.980805 1.983426 U log (K) √ e. cos ω -0.489331 -0.199704 -0.273372 0.203550 -0.406230 0.080327 -0.542331 0.280257 U [m s -1 ] 0.37 0.37 0.37 0.03 0.34 0.40 0.30 0.43 U √ √ e. sin ω 0.029048 0.037665 0.048956 0.162724 -0.157033 0.224967 -0.309643 U 0.370264 e. cos ω -0.370075 -0.326983 -0.359983 0.111428 -0.420448 -0.181879 -0.492776 0.046683 U √ e. sin ω -0.080548 -0.087455 -0.136843 0.149558 -0.250395 0.099871 -0.366481 0.262974 U λ 0 [deg] 185.975741 188.581732 187.252622 5.246122 182.652564 194.619420 176.751015 200.688274 U
		γ COR07-DRS -3-4 γ COR14-DRS -3-8 γ COR98-DRS -3-3 e. sin ω √ e. cos ω √ log (K) λ 0	[m s -1 ] [m s -1 ] [m s -1 ] [m s -1 ] [deg]	39931.00 39958.04 39929.00 0.097884 -0.541593 0.11 -160.261307	39930.61 39956.56 39928.10 0.006756 -0.450052 0.10 -158.719759	39930.34 39956.34 39927.80 0.058079 -0.536450 0.09 -159.674199	1.92 1.46 1.68 0.173485 0.169210 0.05 3.735633	39928.44 39954.89 39926.19 -0.200976 -0.588081 0.04 -162.883613	39932.70 39958.14 39929.98 0.195484 -0.235463 0.15 -154.431466	39925.88 39952.82 39924.15 -0.403356 -0.692747 -0.03 -167.054572	39935.32 39959.96 39932.04 0.363079 0.104560 0.21 -149.921315	U U U U	U U U
		γ HARPS λ 0			[m s -1 ] [deg] -111.177767 39964.88	39964.81 -110.997045	39964.82 -113.124819	0.20 5.834754	39964.58 -117.548307	39965.03 -104.311167	39964.35 -124.016867	39965.25 -97.254697	U	U
		offset UCLES -PU B-2006	[m s -1 ]		5.30	5.44	5.26	0.83	4.51	6.37	3.53	7.36	U
		log (P)					[d]	2.776014	2.776021	2.776014	0.000016	2.776003	2.776040	2.775987	2.776056	U
		log (K)				[m s -1 ]		2.07	2.07	2.07	0.01	2.07	2.08	2.06	2.09	U
		√	e. cos ω							-0.780714	-0.784083	-0.787192	0.005216	-0.789895	-0.777061	-0.793328	-0.773612	U
		√	e. sin ω							0.582695	0.578516	0.573850	0.006592	0.571141	0.587284	0.566660	0.591712	U
		λ 0						[deg]	295.804494	296.111799	296.327571	0.488659	295.460223	296.661315	295.142950	296.986011	U
													Article number, page 27 of 31

This paper presents the characterization of half of the planets announced inMayor et al. (2011). The latter paper was then waiting for the detailed analysis to be published before being resubmitted in parallel to the present paper.

The change of the fibres induces a small RV offset, different for each target. This change also corresponds to a change of scientific focus (and contributors) of our HARPS programs. After the change, only a few new measurements were obtained for our targets, not enough to precisely characterize the instrument offset and sample long-term signals.

The DACE platform is available at http://dace.unige.ch. The online tools used to analyse radial-velocity data can be found in the section Observations => RVs

The goal of our program is to derive reliable statistical distributions of orbital and planet properties that can be used to constrain planet formation models. Signals of less significance are of course of great interest but require more observations to be confirmed as bona fide planets.

Article number, page 26 of 31 Udry et al.: 20 super-Earths and hot Neptunes detected with HARPS

Acknowledgements. We thank the Swiss National Science Foundation (SNSF) and the Geneva University for their continuous support to our planet search programs. This work has been in particular carried out in the frame of the National Centre for Competence in Research 'PlanetS' supported by SNSF. XD is grateful to the Society in Science-The Branco Weiss Fellowship for its financial support. N.C.S. was supported by Fundação para a Ciência e a Tecnologia (FCT, Portugal) through the research grant through national funds and by FEDER through COM-PETE2020 by grants UID/FIS/04434/2013 & POCI-01-0145-FEDER-007672 and PTDC/FIS-AST/1526/2014 & POCI-01-0145-FEDER-016886, as well as through Investigador FCT contract nr. IF/00169/2012/CP0150/CT0002 funded by FCT and POPH/FSE (EC). PF acknowledges support by Fundação para a Ciência e a Tecnologia (FCT) through Investigador FCT contracts of reference IF/00169/2012/CP0150/CT0002, respectively, and POPH/FSE (EC) by FEDER funding through the program "Programa Operacional de Factores de Competitividade -COMPETE". PF further acknowledges support from Fundação para a Ciência e a Tecnologia (FCT) in the form of an exploratory project of reference IF/01037/2013/CP1191/CT0001. ACC acknowledges support from STFC consolidated grant ST/M001296/1. This research has made use of the SIMBAD database and of the VizieR catalogue access tool operated at CDS, France, and used the DACE platform developed in the frame of PlanetS (https://dace.unige.ch).

Appendix A: Periodogram of the activity and CCF indicators

To check if any signal detected in RVs corresponds to signals measured in activity indicator, we show here the Lomb-Scargle periodograms of the log(R HK ), the BIS SPAN, and the FWHM.

Those indicators have been shown to be sensitive to activity, and therefore any signal appearing both in the RVs and at least one of those indicators might be induced by stellar activity. Because magnetic cycles will be seen as a long-period significant signals in all those indicators, we removed any long-term signal either by fitting a Keplerian to adjust at best the observed magnetic cycle, as in Fig. 1, or by adjusting a second order polynomial.

Appendix B: Parameters probed by MCMC

We present here the parameter estimates issued from the MCMC modeling of the planetary systems.