© The Authors 2023

1 Introduction

Studying solar analogs is a way to understand the evolution of stars similar to our Sun. Although the physical and evolutionary stellar characteristics of solar analogs are relatively well studied, the formation and evolution of planetary systems around these stars still lack observational constraints.

The Sun being the nearest and most studied star, where the most accurate stellar properties can be obtained, provides a way to benchmark measurements of distant stars that have similar physical properties. This fact favors a specific class of stars called solar analogs, defined as those stars that have physical properties within a certain range of solar values. Soderblom & King (1998) defined a solar analog to be a main sequence star with an effective temperature between 5278 K and 6278 K, a metallicity within the range ±0.3 dex, and it must not have a close stellar companion (arbitrarily defined with an orbital period of less than ten days).

The motivation to study solar analogs lies in the fact that stellar properties can be obtained with a greater precision since the Sun can be used as a reference (e.g., Bedell et al. 2014; Ramírez et al. 2014). Consequently, the uncertainties of the physical parameters of exoplanets can be significantly reduced for the solar analogs, which makes them interesting laboratories to understand the physics of planetary systems. In addition, stellar abundances can be obtained with greater accuracy, which informs us about the chemical composition of the environment where the planets formed. As a result, chemical clocks (e.g., da Silva et al. 2012) can provide better age diagnoses that can be used to validate the isochronal ages. This can also improve, for example, the constraints in exoplanet interior structure models based on the abundance of refractory elements (Fe, Mg, and Si) of the host star (Dorn et al. 2015; Adibekyan et al. 2021). Ultimately, studying planets around solar analogs also makes our knowledge of the conditions for the development of life less uncertain.

The Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2015) has identified more than 6000 candidate exoplanets (Guerrero et al. 2021). The transit method is efficient to detect close-in planets. Long-term radial velocity (RV) programs, such as the one carried out with the SOPHIE instrument at the Observatoire Haute-Provence, can establish the nature of candidate planets identified by TESS and can also detect additional non-transiting companions. We have identified two exoplanet candidates detected by TESS around the solar analogs TOI-1736 (TIC 408618999) and TOI-2141 (TIC 287256467). These planet candidates were identified as a TESS object of interest (TOI), which presented recurrent transit-like events in their light curves (Jenkins 2002; Jenkins et al. 2010). Our RV follow-up with SOPHIE characterizes the nature of these two planets and detects at least two additional companions in the system TOI-1736. This paper presents the detection and characterization of both systems, where we performed a detailed analysis of their spectra to obtain a refined characterization of the host stars. The stellar parameters are summarized in Table 1. We combined the TESS and RV data to characterize all planets in both systems. Finally, we present an analysis of the internal structure and possible compositions of the innermost planets.

2 Observations

2.1 TESS photometry

The TESS mission observed TOI-1736 with a cadence of 2 min in Sectors 18, 19, and 25, and with a cadence of 20 seconds in Sectors 52, 58, and 59. TOI-2141 was observed with a cadence of 2 min in Sectors 25, 26, and 52. Table 2 shows the log of TESS observations for these two objects. We obtained the TESS data products from the Mikulski Archive for Space Telescopes (MAST)², where we used the Presearch Data Conditioning (PDC) flux time series (Smith et al. 2012; Stumpe et al. 2012, Stumpe et al. 2014) processed by the TESS Science Processing Operations Center (SPOC) pipeline (Jenkins et al. 2016) versions listed in Table 2.

The SPOC searches light curves for transiting planets with an adaptive, noise-compensating matched filter (Jenkins 2002; Jenkins et al. 2010, Jenkins et al. 2020). The transit signatures of each candidate are fitted with an initial limb-darkened transit model (Li et al. 2019) and subjected to a suite of diagnostic tests (Twicken et al. 2018), all available in the TESS SPOC data validation reports (DVR). The DVR for TOI-1736 reports a candidate planet, hereafter TOI-1736 b, with an estimated radius of Rp = 3.0 ± 1.1 R_⊕ and an orbital period of P = 7.07307(1) days. The difference image centroid offsets locate the source of the transit signal within 4.4 ± 4.1 arcsec of the target star. The TESS Science Office (TSO) issued an alert for TOI 1736.01 on 27 February 2020 based on the DVR associated with the combined light curve for Sectors 18–19 (Guerrero et al. 2021).

The DVR for TOI-2141 reports two candidates, but only one passed validation tests. The unvalidated events have a periodicity of 30.6 days and are likely due to contamination by scattered light. We also note that the odd and even depth test supports a single strong feature folded on top of a much weaker and much less convincing transit-like feature for these events. The validated candidate planet, hereafter TOI-2141 b, has an estimated radius of R_p = 3.2 ± 0.3 R_⊕ and an orbital period of P = 18.26159(7) days. The centroid offsets localize the source of the transit signal within 1.7 ± 6.3 arcsec. The transit signature was first identified in the SPOC search of the Sector 26 light curve, and an alert for TOI 2141.01 was issued by TSO on 7 August 2020. Therefore, both TOI-1736 and TOI-2141 systems host at least one close-in transiting planet candidate with a size consistent with a sub-Neptune. Figures 1 and 2 depict the TESS photometry data for both respective objects.

2.2 LCOGT 1 m NEB Search

The TESS pixel scale is ~21″ pixel⁻¹, and photometric apertures typically extend out to roughly 1 arcmin, which generally results in multiple stars blending in the TESS aperture. To attempt to determine the true source of the TESS detection, we conducted ground-based photometric follow-up observations of the field around TOI-1736 as part of the TESS Follow-up Observing Program³ Sub Group 1 (TFOP; Collins 2019).

We observed a full predicted transit window of TOI-1736 b in Pan-STARRS z-short band using the Las Cumbres Observatory Global Telescope (LCOGT; Brown et al. 2013) 1.0 m network node at McDonald Observatory on UTC 2020 August 28. If the event detected in the TESS data is indeed on-target, the shallow SPOC reported depth of 380 ppm would not generally be detectable in ground-based observations. Instead, we slightly saturated TOI-1736 to enable the extraction of light curves of nearby fainter stars to attempt to rule out or identify nearby eclipsing binaries (NEBs) as potential sources of the TESS detection. The 1 m telescopes are equipped with 4096 × 4096 SINISTRO cameras having an image scale of 0″.389 per pixel, resulting in a 26′ × 26′ field of view. The images were calibrated by the standard LCOGT BANZAI pipeline (McCully et al. 2018), and photometric data were extracted using AstroImageJ (Collins et al. 2017).

To account for possible contamination from the wings of neighboring star point spread functions (PSFs), we searched for NEBs in all known Gaia DR3 and TICv8 nearby stars out to 2′.5 from TOI-1736 that are possibly bright enough in the TESS band to produce the TESS detection (assuming a 100% eclipse and 100% contamination of the TESS aperture). To attempt to account for possible delta-magnitude differences between the TESS band and the follow-up filter band, we checked stars that are an extra 0.5 magnitudes fainter in TESS-band than needed. We find that the RMS of each of the light curves of the 24 stars matching our criteria is more than a factor of 5 smaller than the expected NEB depth in the respective star. We then visually inspected each neighboring star’s light curve to ensure no obvious eclipse-like signal. All of our follow-up light curves and supporting results are available on the EXOFOP-TESS website⁴. Through our process of elimination, we find that the TESS signal must be occurring in TOI-1736 relative to known Gaia DR3 and TICv8 stars.

2.3 High contrast imaging

High-angular-resolution observations can probe close companions within ~1.2 arcsec that can create a false positive transit signal if that companion is an eclipsing binary, and which dilute the transit signal and thus yield underestimated planet radii (Cia-rdi et al. 2015). TOI-1736 was observed on October 19, 2021, by the ‘Alopeke dual-channel speckle imaging instrument on Gemini-N (PI: Howell) with a pixel scale of 0.01 arcsec pixel⁻¹ and a full width at half maximum (FWHM) resolution of 0.02 arcsec. ‘Alopeke provided simultaneous speckle imaging at 562 and 832 nm. The data were processed with the speckle pipeline (Howell et al. 2011), which yielded the 5-sigma sensitivity curves shown in Fig. 3. These observations provide a contrast at an angular separation of 0.5 arcsec of 4.2 mag at 562 nm and 6.4 mag at 832 nm. Figure 3 also shows the reconstructed image at 832 nm and its Fourier transform, which shows fringes in the power spectrum, evidence of a stellar companion to TOI-1736. The best-fit binary model for the power spectrum gives a companion at an angular separation of 0.093 ± 0.002 arcsec and position angle of 26.5 ± 1.0 deg with a magnitude difference from the primary of Amag = 2.43 ± 0.15 at 832 nm.

TOI-2141 was observed on April 25, 2021, with the 4.1-m SOAR telescope in speckle imaging. These observations provide a contrast of 5 mag at an angular separation of 1.0 arcsec in the I-band (Fig. 4). These observations did not reveal any evidence of close stellar companions to TOI-2141.

2.4 SOPHIE spectroscopy

2.4.1 Observations

SOPHIE is a high-resolution fiber-fed, cross-dispersed échelle spectrograph mounted on the 1.93-m telescope at the Obser-vatoire de Haute-Provence (OHP; Perruchot et al. 2008; Bouchy et al. 2013). It covers a wavelength domain from 387.2 nm to 694.3 nm across 39 spectral orders. We observed TOI-1736 and TOI-2141 under a program dedicated to an RV follow-up of transiting candidates (e.g., König et al. 2022; Moutou et al. 2021; Hébrard et al. 2020), where we obtained 152 spectra of TOI-1736 between 2020-08-20 and 2023-03-11 with an average peak signal-to-noise ratio (S/N) per pixel at 550 nm of 70, and 90 spectra of TOI-2141 between 2021-02-25 and 2022-09-18 with an average peak S/N of 55. Both targets have been observed in the high spectral resolution mode (HR mode, R = 75 000) of SOPHIE. Tables A.1 and A.2 present more information about these observations.

Fig. 1 TESS light curve of TOI-1736. The gray points show the TESS photometry data normalized by the median flux of each sector. The black points show the same data binned by the weighted average with bin sizes of 0.1 d. The red line shows the baseline GP model that was fitted to the binned data multiplied by the best-fit transit model for TOI-1736 b. The red triangles show the predicted central time of TOI-1736 b transits.

Fig. 2 TESS light curve of TOI-2141. The blue points show the TESS photometry data around the transits of TOI-2141 b. The red line shows the baseline GP model that was fitted to the binned data multiplied by the best-fit transit model for TOI-2141 b. The red triangles show the predicted central time of TOI-2141 b transits.

Fig. 3 Contrast ratio of TOI-1736 as a function of angular separation at 562 nm (blue line) and at 832 nm (red line) obtained from the ‘Alopeke/Gemini speckle imaging observations. The small panels show the reconstructed speckle image (right subpanel) and its power spectrum in Fourier space (left subpanel). The power spectrum shows fringes that match a binary model.

Fig. 4 Contrast ratio of TOI-2141 as a function of angular separation at I-band obtained on the 4.1-m SOAR speckle imaging observations.

2.4.2 Data reduction by the DRS

Our data have been reduced by the SOPHIE Data Reduction Software (DRS; Bouchy et al. 2009), excluding the seven first (redder) spectral orders due to their low S/N. Five spectra of TOI-1736 were not used due to their low overall S/N (below 40). The DRS automatically extracts and calibrates the spectra and computes the RV using the cross-correlation function (CCF) between the spectra and a G2-type empirical weighted numerical mask. The DRS also uses the CCF data to deliver stellar activity indicators, such as the CCF FWHM and the bisector span (BIS).

In order to improve the accuracy of SOPHIE measurements, we used the optimized procedures presented by Heidari et al. (2022, and in prep.). This includes in particular: (1) CCD charge transfer inefficiency correction (Bouchy et al. 2013); (2) correction for the moonlight contamination using the simultaneous sky spectrum obtained from the second SOPHIE fiber aperture (e.g., Pollacco et al. 2008; Hébrard et al. 2008); (3) RV constant master correction for instrumental long-term drifts (Courcol et al. 2015); and (4) correction of the instrumental short-term drifts thanks to the frequently measured drifts interpolated at the precise time of each observation.

2.4.3 The CCF analysis by sophie-toolkit

In addition to the DRS reduction, we implemented an independent CCF analysis using the methodology described in Martioli et al. (2022), where we developed the Python package sophie-toolkit⁵ to obtain the RVs and other CCF quantities from a set of reduced SOPHIE spectra (both in the e2ds or sld formats provided by the DRS). Our CCF approach is fundamentally the same as that applied by the DRS and therefore it should produce equivalent results. We compared the RVs from both CCF analyses, and despite some systematic effects in a small fraction of the data, they mostly agree within the error bars. As explained in Sect. 2.4.2, the DRS implements additional procedures to optimize the RV measurements in the SOPHIE data, which ultimately provide more accurate RVs. Thus, in our analysis, we adopted the DRS RVs.

The analysis performed by sophie-toolkit was originally designed for near-infrared observations with the SPIRou spectrograph (Donati et al. 2020), where it implements some processing steps that are aimed at mitigating the strong effects of telluric contamination and detector artifacts. Such analysis is more robust to systematics and improves some results when applied to the SOPHIE optical spectrum, where the aforementioned effects are much less significant, but still present in the data. An important step that is applied in our tool is an iterative registration of the spectra, both in flux and wavelength domains, to match each observation to a high S/N template stacked spectrum (“template” hereafter). Registration is done by shifting the spectra to the same topocentric frame using the precise RV values obtained from the CCF analysis and by scaling the fluxes using an order-by-order least squares fit to a third-order polynomial (see Martioli et al. 2022). This reduces residual systematic errors between observations, providing a differential measurement of each spectrum with respect to the template. This procedure also provides the dispersion in flux for each spectral element along the time series, which allows us to estimate statistical errors and assign proper weights to the data before performing measurements such as CCFs, RVs, and spectral indices. In addition to the CCF and RV data, sophie-toolkit provides other important products, such as the template spectrum, obtained from the median stack of all spectra in the time series. This product is saved in FITS format to store both the template and all other registered spectra. The most relevant activity indicators are also calculated and saved in a time series product, as detailed in Sect. 3.6. Appendix A presents the final RV data (from the DRS) together with the corresponding spectral quantities obtained in our analysis using the sophie-toolkit package, all compiled in Tables A.1 and A.2.

2.5 TRES spectroscopy

For spectroscopic reconnaissance, we obtained four spectra of TOI-1736 and two spectra of TOI-2141 using the 1.5 m Tillinghast Reflector Echelle Spectrograph (TRES; Fűrész et al. 2008) located at the Fred Lawrence Whipple Observatory (FLWO) in Arizona. The spectra for TOI-1736 and TOI-2141 were obtained between March 3–October 4, 2020, and between August 17 and October 20, 2020. TRES is a fiber-fed echelle spectrograph with a wavelength range of 390–910 nm and a resolving power of R = 44000. Spectra were extracted and reduced as described in Buchhave et al. (2010).

3 Stellar characterization

We carried out an analysis of the spectral and photometric data of both systems to derive the host stars’ properties as is detailed in the next sections. A summary of the final star parameters is presented in Table 1.

Fig. 5 SED fit analysis for TOI-1736. The upper panel shows results assuming a single stellar component and the lower panel shows results assuming two stellar components. The red symbols represent the observed photometric measurements, whereas the horizontal bars represent the effective passband width. The blue symbols are the model fluxes from the best-fit NextGen atmosphere model (black and gray lines). The lower panel also shows the model for the primary component (cyan line) and the secondary component (red line), separately. The insets display our model (black curve) and the observed absolute flux-calibrated Gaia spectrum (gray swathe).

3.1 Spectral energy distribution analysis

We performed an analysis of the broadband spectral energy distribution (SED) of both stars including the Gaia Early Data Release 3 parallaxes (with no systematic offset applied; see, e.g., Stassun & Torres 2021) in order to determine an empirical measurement of the stellar radius, following the procedures described in Stassun & Torres (2016) and Stassun et al. (2017, Stassun et al. 2018).

We performed a fit using NextGen stellar atmosphere models (Hauschildt et al. 1999), with the effective temperature (T_eff), metallicity ([Fe/H]), and extinction (A_V) as free parameters. The surface gravity (log 𝑔) has little influence on the broadband SED.

As a test to investigate whether the companion to TOI-1736 detected by the high contrast imaging results presented in Sect. 2.3 is real, we first fit a single-component model to the SED, as shown in the top panel of Fig. 5, where we obtained a reduced χ² of 1.1, with the best-fit parameters T_eff = 5800 ± 75 K, [Fe/H] = 0.2 ± 0.1, and A_V = 0.15 ± 0.03. According to the tool Stilism⁶ (Lallement et al. 2014), the expected extinction for the distance of ~90 pc, galactic latitude and longitude from TOI-1736 is A_V = 0.003 mag, therefore an abnormal reddening (or an infrared excess flux) was detected in our model. Then we consider a two-component stellar model to fit the SED of TOI-1736, assuming that the secondary has the same extinction (A_V = 0.003) and metallicity as the primary. The bottom panel of Fig. 5 shows these results, which give a χ² of 1.0 with the same fit parameters for the primary and T_eff = 4800 ± 200 K for the secondary. Hence, the SED of TOI-1736 provides additional support to the binary hypothesis. For TOI-2141, we adopt a single-component model as illustrated in Fig. 6, which gives a reduced χ² of 1.3 and the best fit parameters T_eff = 5750 ± 75 K, [Fe/H] = −0.1 ± 0.1, and A_V = 0.06 ± 0.04.

Integrating the SED model gives the bolometric flux on Earth (F_bol), which can be combined with the T_eff and the Gaia parallax to derive the stellar radius. For TOI-1736, the single component model gives a radius of R_★ = 1.390 ± 0.037 R_⊙ and the two-component model gives a 3σ smaller radius of R_★ = 1.241 ± 0.032 R_⊙ for the primary and R_★ = 0.748 ± 0.073 R_⊙ for the secondary. For TOI-2141, we obtained a radius of R_★ = 0.938 ± 0.027 R_⊙. Using the mass-radius relations of Torres et al. (2010), we calculate the stellar masses of 1.15 ± 0.07 M_⊙ for TOI-1736 and 0.98 ± 0.06 M_⊙ for TOI-2141.

Fig. 6 SED fit analysis for TOI-2141. Red symbols represent the observed photometric measurements, whereas the horizontal bars represent the effective width of the passband. Blue symbols are the model fluxes from the best-fit NextGen atmosphere model (black). The inset displays our model (black curve) alongside the observed absolute flux-calibrated Gaia spectrum (gray swathe).

3.2 Standard analysis of SOPHIE spectra

To obtain the stellar parameters from the SOPHIE data, we constructed an average spectrum for both stars, using all SOPHIE spectra that were unaffected by the Moon’s pollution, and performed a standard spectral analysis on them. The average spectra of TOI-1736 and TOI-2141 have an S/N at 649 nm of 1100 and 850, respectively. We use the methods presented in Santos et al. (2004) and Sousa et al. (2008) to derive the effective temperatures T_eff, the surface gravities log 𝑔, and the metallicities [Fe/H]. Using these spectroscopic parameters as input, the stellar masses were derived from the Torres et al. (2010) calibration with a correction following Santos et al. (2013). The errors were calculated from 10000 random realizations of the stellar parameters within their error bars and assuming Gaussian distributions. The parameters obtained from this analysis are presented in Table 4 labeled as “standard”.

3.3 Spectroscopic parameters from TRES spectra

We use the reduced TRES spectra to obtain the stellar parameters using the Stellar Parameter Classification tool (SPC; Buchhave et al. 2012). In short, SPC correlates each observed spectrum against a grid of synthetic spectra based on Kurucz (Kurucz 1992) atmospheric models and derives the effective temperature, surface gravity, metallicity, and rotational velocity of the star. The parameters obtained from this analysis are presented in Table 4.

3.4 Strictly differential analysis of SOPHIE spectra

We take advantage of the spectroscopic similarity between our targets and the Sun to apply a strictly differential analysis to determine atmospheric parameters and elemental abundances (e.g., Ramírez et al. 2014; Bedell et al. 2014). The differential analysis in such cases can potentially mitigate the Fe I/II excitation-ionization balance internal errors caused by possible differences in the micro-physics prescriptions, in the local determination of the continuum, or in the equivalent width (EW) measurements. Our reference solar spectrum was obtained from observations of the Moon with SOPHIE on 2023-04-01 with an S/N of 245 at 649 nm, maintaining the same instrumental configuration as the one used to observe the scientific targets.

We perform a line-by-line analysis with respect to the solar spectrum, which is based on the following steps. Given a narrow spectral range around each selected line, we normalize the stellar and solar spectrum using the same local continuum regions and then fit a Gaussian function to the solar line profiles, computing all resulting parameters (e.g., line depth, centroid, FWHM, and continuum offset). Then, we repeat the line fitting procedure for the locally normalized stellar spectrum, now fixing the continuum offset and line region information obtained from the respective solar spectrum. We employ the ionization-excitation equilibrium method where we adopted a line list from Meléndez et al. (2014), which contains excitation potentials, oscillator strengths, laboratory log 𝑔f values, and hyperfine structure corrections⁷. We included iron-peak, s- and alpha-capture elements, such as Fe, Mg, Al, Si, Ti, Y, and Ca to perform this spectroscopic analysis.

In Fig. 7, we show the relationship between the solar Fe I and Fe II EWs and those of TOI-1736 and TOI-2141. The low scatter (σEW ~ 1–2%) for both stars indicates that the differential spectroscopic equilibrium analysis with respect to the Sun is sufficiently precise for the determination of the stellar atmospheric parameters. While not an exact match to our Sun, exercising caution and maintaining a conservative approach in considering uncertainties is warranted.

We use the code q2 (Ramírez et al. 2014) to model both the stellar and solar EWs of Fe I and Fe II from our line-by-line analysis. The input-relevant line data (e.g., oscillator strength, excitation potential, etc.) are also provided. The q2 was configured to employ the Kurucz ODFNEW model atmospheres (Castelli & Kurucz 2003) and the 2019 version of the local thermodynamic equilibrium (LTE) code MOOG (Sneden 1973; Sneden et al. 2012). The stellar parameters are obtained through an iterative convergence process, where q2 iterates until it finds: 1) no systematic dependence of the derived Fe I and Fe II abundances as a function of excitation potential and reduced EW, constraining the T_eff and micro-turbulence; 2) Fe I and Fe II yielding the same abundance value, constraining the log 𝑔. Once convergence is achieved, the final atmospheric parameters are estimated together with their uncertainties, which include systematic errors, as detailed in Appendix B. For TOI-1736 we obtained T_eff = 5807 ± 13 K, log 𝑔 = 4.35 ± 0.04 dex, [Fe/H] = +0.138 ± 0.014 dex, and υ_t = 1.26 ± 0.02 km s⁻¹. For TOI-2141 we obtained T_eff = 5659 ± 12 K, log 𝑔 = 4.42 ± 0.04 dex, [Fe/H] = −0.120 ± 0.013 dex, and υ_t = 0.94 ± 0.02 km s⁻¹.

Stellar elemental abundances are obtained by q2 using local thermodynamic equilibrium (LTE) model interpolations based on the previously derived atmospheric parameters. The adopted errors for each element account for the errors from atmospheric parameter measurements and the scatter from the line-to-line abundance estimates. Table 3 summarizes the abundances obtained in our analysis.

We performed an evolutionary analysis using the derived spectroscopic parameters (T_eff, log 𝑔 and abundances), where we adopted the Yonsei-Yale evolutionary tracks (Yi et al. 2001; Demarque et al. 2004) that account for the contribution of the alpha enhancement elements to compute the ages, masses, and radii of TOI-1736 and TOI-2141. The derived atmospheric parameters are compared with those predicted by stellar models using the Bayesian framework described in Grieves et al. (2018) and Yana Galarza et al. (2021). Table 4 shows the derived values of stellar mass and radius and the derived ages are shown in Table 5.

Fig. 7 Comparison between solar and stellar EW measurements. TOI-1736 and TOI-2141 are shown in the top and bottom panels, respectively. The dashed black line is the 1:1 EWs ratio, and the solid black line is a linear fit based on the Fe I lines, represented by empty black circles. The red star symbols represent the Fe II lines, followed by their linear fit represented by the solid red line. The typical standard deviation of the average linear fit is around ~0.1–0.2 nm, which means ~2% in relative terms.

3.5 Comparison between methods for determining stellar parameters

As we have employed different methods and used data from different instruments to obtain the stellar parameters, in this section we make a comparison between these results. At the time of writing, MacDougall et al. (2023) has also published independent results for the parameters of TOI-1736 (T_eff and [Fe/H]), which we also include here for comparison.

For TOI-1736, the mean values of T_eff range from 5636 to 5807 K, where the SED, standard, and solar diff. agree within 1σ, whereas TRES and MacDougall et al. (2023) agree within 3σ with respect to the other values. All surface gravity values also agree within 3σ, although the solar diff. gives a higher value than other methods. The Gaia eDR3 trigonometric value is slightly more precise than the solar diff, but it can be contaminated by the close companion in the same way as the SED value. All metallicity values agree within 1σ. For TOI-2141, all but the SED value of T_eff agree within 1σ, and all agree within 2σ. All surface gravity and metallicity values also agree within 1σ. The parameters of this star are more similar to our Sun, being poorer in metals and more active than normal for its age.

Different instruments/methods may have different offsets for derived parameters. For this reason, we performed an absolute calibration of our differential analysis of the SOPHIE data using other stars and the Sun as a reference (see Appendix B). Since the determination is more precise and since the other methods generally agree with it, we adopt the parameters of the solar differential analysis as final values. Star masses and radii are derived from spectroscopic parameters, so they carry more or less the same level of discordance as those discussed above. For consistency, we adopted the final mass and radius values as those obtained from the solar differential analysis. As pointed out by (Tayar et al. 2022), we recognize that even the differential analysis has systematic errors in the models of stellar evolution that were not taken into account in our analysis. Therefore, we use Tayar et al. (2022)’s methods⁸ to calculate the additional errors in stellar parameters due to systematic differences in the models, as shown in Table 4. We add these errors in quadrature to the uncertainties obtained by our differential analysis, giving a more realistic uncertainty in the final values.

Turbulence velocity and rotational velocity are two quantities that depend on the adopted model and can absorb some instrumental broadening. Thus, to avoid any bias toward a specific instrument/model, we adopted their mean values.

3.6 Stellar activity and rotation

The CCF analysis performed on the SOPHIE spectra provides three quantities that are sensitive to stellar activity: the FWHM, bisector span (BIS), and RVs (e.g., Queloz et al. 2001; Boisse et al. 2009), as shown in Appendix A. In addition, we measured two well-known spectral indices that are proxies for chromospheric and coronal activity, the S-index, which relies on the emissions in the cores of the H and K Ca II lines, and Hα (e.g., Cincunegui et al. 2007). In Appendix C we present the details of our measurements of the S-index and Hα from the SOPHIE spectra. In particular, we present a calibration of the S-index to the Mount-Wilson system. The temporal variability of these activity indicators can be used to identify systematic errors in the RVs caused by spurious signals of stellar activity, thus allowing to improve the determination of planetary orbits, as will be explored in Sect. 4.

We convert the Mount-Wilson S-index to the chromospheric activity index log R′_HK using the Czesla et al. (2019) recipe that includes the photospheric correction (e.g., Mittag et al. 2013), where we use the B – V values from the magnitudes listed in Table 1 and T_eff from our analysis. The result was log ${R^{\prime }}_{\text{HK}}\left(B-V\right)=-{4.94}_{-0.25}^{+0.15}$ for TOI-1736 and log ${R^{\prime }}_{\text{HK}}\left(B-V\right)=-{4.78}_{-0.17}^{+0.13}$ for TOI-2141. Using the activity-rotation empirical calibration of Mamajek & Hillenbrand (2008) we obtained P_rot = 28 ± 5 days for TOI-1736 and P_rot = 21 ± 5 days for TOI-2141.

Also using the calibration of Mamajek & Hillenbrand (2008) we obtained gyrochronological ages of ${4.6}_{-1.3}^{+1.7}$ Gyr for TOI-1736 and ${2.7}_{-1.0}^{+1.6}$ Gyr for TOI-2141. Applying the age-chromospheric activity relation for solar analogs of Lorenzo-Oliveira et al. (2016) we estimated the ages in ${5.4}_{-1.1}^{+1.4}$ Gyr for TOI-1736 and ${4.2}_{-1.1}^{+1.6}$ Gyr for TOI-2141, which are consistent with those from gyrochronology. We adopted the activity ages as the latter, since Lorenzo-Oliveira et al. (2018) take mass and metallicity into account in their empirical relationship, proving to be more accurate for the range of parameters of our targets.

3.7 Stellar ages

Table 5 shows the derived age values of the evolutionary analysis of stellar parameters. As an additional check, we calculate the stellar ages from chemical clocks (e.g., da Silva et al. 2012; Tucci Maia et al. 2015; Delgado Mena et al. 2019). We calculated the yttrium abundance ratios ([Y/X], where X=Si I, Ti I, Ti II, or Al I) from the values reported in Table 3, where we estimated the stellar ages using the Delgado Mena et al. (2019) relations that account for T_eff and [Fe/H] effects on chemical ages. The average calibration errors are about 1.4 Gyr and the impact of T_eff-[Fe/H]-[YII/X] errors are respectively 0.5 and 0.3 Gyr for TOI-1736 and TOI-2141. The average chemical ages are shown in Table 5.

Lithium (Li) is also known to be a sensitive indicator of stellar evolution (e.g., Carlos et al. 2016; Lyubimkov 2016). To obtain Li abundance, we performed a spectral synthesis analysis using the code iSpec⁹ (Blanco-Cuaresma et al. 2014; Blanco-Cuaresma 2019) and the radiative transfer code MOOG (Sneden 1973; Sneden et al. 2012), where we use the lines in the spectral range around Li from Meléndez et al. (2012). Figure 8 illustrates our results, where we show the SOPHIE spectra and the best-fit synthetic models for TOI-1736 and TOI-2141, and also for a solar spectrum obtained from observations of the Moon. We calculated a Li abundance for the solar spectrum of A(Li) = 1.1 ± 0.2, which is consistent with the expected value for the Sun. For TOI-2141, our results are compatible with a null value, giving an upper limit of A(Li) <0.5. For TOI-1736, we obtained an abundance of A(Li) = 1.56 ± 0.16. Using the relationship between A(Li) and age from Carlos et al. (2019) we obtained an age of 4.4 ± 1.0 Gyr for TOI-1736.

All the different age determinations for both systems are presented in Table 5. There is general agreement between the ages of TOI-1736, which are in the range of 4–6 Gyr. For TOI-2141, both the chemical and evolutionary ages agree around 7.5 Gyr, which is consistent with the absence of Li, but the activity-derived age appears younger. For a final determination of the ages of both stars, we adopted the average values, where we added Tayar et al. (2022)’s model uncertainties in quadrature to the measured errors, as in Sect. 3.5, giving an age of 4.9 ± 1.3 Gyr for TOI-1736 and 6.4 ± 1.8 Gyr for TOI-2141.

Fig. 8 Spectral synthesis analysis around the Li region to obtain an estimate of the Li abundance. The points show the SOPHIE spectra of TOI-1736 (in red), TOI-2141 (in green), and a solar spectrum (in blue) from observations of the Moon. The solid lines show the respective bestfit synthetic spectra obtained from our analysis.

4 Detection and characterization of the planetary systems

4.1 Analysis of TESS photometry data

To further characterize the transits detected by TESS, we analyzed the PDCSAP flux data by applying the methods described in Martioli et al. (2021), where we fit a transit model plus a baseline polynomial in selected windows around each transit in such a way that the window size is twice the duration of the transit and with a minimum of 300 data points. We selected windows for a total of 19 transits of TOI-1736 b and four transits of TOI-2141 b. Our transit model is calculated using the BATMAN toolkit (Kreidberg 2015) and the posterior distribution of transit parameters is sampled using a Bayesian Monte Carlo Markov chain (MCMC) framework with the package emcee (Foreman-Mackey et al. 2013). We use uninformative priors for the transit parameters, as presented in Table D.1, with initial values given by those reported in the TESS DVRs and assuming circular orbits.

We divided the TESS light curve by the best-fit transit model and binned the data by the weighted average with bin sizes of 0.1 day to fit a baseline flux using Gaussian processes (GP) regression with a quasi-periodic (QP) kernel as in Martioli et al. (2022). This kernel was chosen to detect a possible variability modulated by the star’s rotation period. The priors of the GP parameters are listed in Table D.1 and the posteriors are listed in Table D.2. We opted to fix the smoothing factor and decay time values at 0.1 and 10 days, respectively, as altering these values does not impact the other parameters. Unfortunately, the periodicities obtained in our analysis are not well constrained, and their values change significantly depending on the initial guess adopted for our GP model. The low level of intrinsic photometric variability in both targets, coupled with the rotation periods expected from stellar activity (16–32 days) being closely aligned with the TESS sector duration, strongly biases the detection of periodicity. This bias is primarily due to systematic variations that occur during TESS observations, particularly near the edges of each sector, where the amplitude of these variations is greater. While there may be some intrinsic stellar variability captured by TESS photometry, the GP model obtained here may account for both systematic and intrinsic stellar phenomena in unknown proportions.

Since the periodic kernel is not essential, the primary purpose of employing the GP in this analysis is for detrending. To ensure that the choice of GP does not impact the planet parameters, we conducted experiments using both a squared exponential (SE) kernel and a model without GP. The posterior distributions of key planet parameters (T_c, a/R_★, P, and R_p/R_★) are consistent within 1-sigma for all three choices of baseline model. The results of our best-fit transit models and GP analysis using a QP kernel are illustrated in Figs. 1 and 2.

4.2 Radial velocities

The TOI-1736 RVs have a median value of −25 408.95 m s⁻¹ with an rms of 158 m s⁻¹ and a median individual error of 2 m s⁻¹. These RVs present a clear signal with an amplitude of hundreds of m s⁻¹, as can be seen in Fig. 13. An initial analysis of the generalized Lomb-Scargle (GLS) periodogram (Zechmeister & Kürster 2009) of this data show that there is a dominating signal with periodicity around 570 days and a long-term linear trend was also detected by visual inspection. To detect the TESS candidate planet, TOI-1736.01, at a period of 7.073 days, we fit the SOPHIE RV data using the online tool Data & Analysis Center for Exoplanets (DACE¹⁰), where we adopted a two-planet Keplerian model plus a linear trend. As illustrated in Fig. 9, the GLS of the RVs after subtracting the linear trend and each planet’s orbit shows a strong peak at 7.0731 d and 570.07 d for the inner and outer planet, respectively. The former agrees well with the periodicity of the transits detected by TESS, showing the RV detection of the inner planet candidate, TOI-1736 b, and an outer giant planet, TOI-1736 c. Moreover, this confirms that transits detected by TESS occur in the primary star, as the dominant flux in the SOPHIE spectra comes from the primary and not from the K3 dwarf companion. As a way to test if the SOPHIE spectra are contaminated by the flux from the companion, we performed a CCF analysis using a K0, K5, or M0 masks, which should favor the detection of spectral lines of the cooler star. Those tests resulted in slightly noisier CCF with parameters in agreement as when using the G2 mask, showing that indeed, the flux contribution of the companion to the RV measurements in the SOPHIE spectra is negligible.

The TOI-2141 RVs have a median value of −19 860.8 m s⁻¹ with an rms of 7 m s⁻¹ and a median individual error of 2.4 m s⁻¹. The RVs of this object do not appear to have any long-term trends or variability. Figure 10 shows the GLS periodogram for this data, which shows the highest peak at a period of 18.259 days, again in agreement with the periodicity of the transits detected by TESS. This shows the detection of the planet candidate TOI-2141 b in our SOPHIE data. In the following sections, we present an analysis to characterize these two planetary systems.

For both systems, we found that there is no correlation between the RVs and the four stellar activity indicators: FWHM, BIS, H_α, and S -index. Therefore, the signals detected in RV are probably due to Doppler shift and not to scenarios that imply changes in the line profile, for example, blended binaries. We have also inspected possible correlations between the RV residuals (after subtracting the fit models) and activity indicators. Only BIS shows a slight correlation of 0.29 and 0.15 for TOI-1736 and TOI-2141, respectively. We tried to de-correlate the RVs as in Boisse et al. (2009), but this did not have a significant impact on the values of the fit parameters, so we preferred to keep the final analysis without applying this correction. In Appendix C we show the time series of each activity indicator, an analysis of the GLS periodogram, and its correlations with the RV data.

Fig. 9 Generalized Lomb-Scargle periodogram for the TOI-1736 SOPHIE RVs. The top and bottom panels show the GLS for the SOPHIE RVs minus the linear trend and the fit orbits of planets c (top panel) and b (bottom panel).

Fig. 10 Generalized Lomb–Scargle periodogram for the TOI-2141 SOPHIE RVs. The blue dashed line shows the peak power at 18.3 days, which agrees with the period reported by TESS for TOI-2141 b.

4.3 Joint analysis of RVs and photometry to characterize the planets

To obtain the physical parameters of the planets, we perform a joint analysis of the photometry and RV data sets as in Martioli et al. (2022). We adopt the priors listed in Table D.1 and the initial values are those obtained by the analysis performed on the SOPHIE RVs and TESS photometry data independently.

We first fit the RV orbits and transit models jointly using the scipy.optimize.leastsq optimization tool. Then, we use the same tool to fit the RV jitter, which is a term quadratically added to the RV errors. Finally, we explore the full range covered by the priors with a Bayesian MCMC that uses 50 random walkers, with 10 000 iterations, discarding the first 3000 burn-in samples. The posterior distributions of sampled parameters are illustrated in Appendix D and the resulting best-fit parameters, obtained from the mode of the posterior and 1σ uncertainties (34% on each side of the central value) are given in Table 6. Some parameters converged to values close to the edge of their physical domains. For example, it happened for the orbital inclination (i_b ~ 90 deg) and the linear limb-darkening coefficient (u₀ ~ 0) for both transiting planets, TOI-1736 b and TOI-2141 b. We also adopted circular orbits (e = 0, ω = 90 deg) for both inner planets, as noncircular solutions resulted in higher values of the Bayesian Information Criterion (BIC), suggesting that our data cannot constrain the low values of the eccentricities of these planets.

Figures 11 and 12 illustrate the TESS photometry data for all transits detected by TESS and the best-fit transit models obtained by our joint analysis. Figures 13 and 14 illustrate the TOI-1736 RV data and the best-fit models that include the orbit of the two planets, TOI-1736 b and TOI-1736 c, and the linear trend. Figures 15 and 16 illustrate the RV data and best-fit model for TOI-2141.

Fig. 11 Transits of TOI-1736 b observed by TESS. The blue points show the TESS photometry data around the nineteen transits of TOI-1736 b, with the times relative to the central time of each transit. The red line shows the best-fit transit model and the green points show the residuals plus an arbitrary offset for better visualization.

Fig. 12 Transits of TOI-2141 b observed by TESS. The blue points show the TESS photometry data around the four transits of TOI-2141 b, with the times relative to the central time of each transit. The red line shows the best-fit transit model and the green points show the residuals plus an arbitrary offset for better visualization.

5 Gaia astrometry to constrain the mass of TOI-1736 c

The expected astrometric signal due to the orbit of the star TOI-1736 around the center of mass of the system can be estimated by the following equation: $$a_{★}=p{a}_{\text{p\hspace{0.17em}}}\left(\frac{m_{\text{p}}}{M_{★}}\right),$$(1)

where p is the parallax in milli-arcseconds, a_p is the semi-major axis of the planet’s orbit in au, m_p is the planet’s mass and M_★ is the stellar mass. For TOI-1736 c, we estimate a_★ = 0.112 ± 0.005 mas. Here, we use the method Gaia Astrometric noise Simulation To derive Orbit iNclination (GASTON; Kiefer 2019; Kiefer et al. 2019) to constrain, from Gaia DR3 (Gaia Collaboration 2021) astrometric excess noise and RV-derived orbital parameters, the orbital inclination and true mass of TOI-1736 c. The Gaia DR3 astrometric excess noise of TOI-1736 is ɛ_DR3 = 0.09221 mas. Figure 17 shows the results of our GASTON simulations of astrometric excess noise as a function of the orbital inclination, assuming the RV orbit of TOI-1736 c. We obtained an orbital inclination of $i_{\text{p}}={69}_{-17}^{+15}$ deg, which gives a semi-major axis of the star photocenter of $a_{{}_{★}}={0.119}_{-0.007}^{+0.022}$ mas. These results show that the orbit of the TOI-1736 c seems to be tilted with respect to the inner planet, but is not inconsistent, within 2σ, with a coplanar orientation. Assuming coplanarity, we computed the prospective transit epochs for TOI-1736 c; however, none of the predicted events align with TESS observations. Finally, we derived the mass of TOI-1736 c by considering our measurements of M_p sin i_p combined with the astrometric constraint. Figure 18 shows the posterior distribution for the mass of TOI-1736 c, where we obtained $M_{\text{p}}={8.7}_{-0.6}^{+1.5}$ at 1σ with a maximum value of M_p < 16.54 M_Jup at 3σ.

Fig. 13 TOI-1736 SOPHIE RVs. In the top panel, the red points show the SOPHIE RV data and the green line shows the best-fit model, including the orbits of TOI-1736 b and TOI-1736 c, plus a linear trend, which is also represented separately by the red dashed line. The bottom panel shows the residuals.

Fig. 14 TOI-1736 SOPHIE RVs in phase with the orbits of the planets. In the top panel, the gray points show the RVs phase-folded in the period of planet b, where the orbit of planet c and the linear trend have been subtracted. In the bottom panel, the gray points show the RVs phase-folded in the period of planet c, where the orbit of planet b and the linear trend has been subtracted. The red points show the binned data, where we use a bin size of 0.05 in units of orbital phase. The green lines represent the best-fit orbit model for each respective planet.

6 Discussions

The SOPHIE RV data show periodic variations at the same period as the transits detected by TESS, where the orbital fit gives a semi-amplitude of K_p = 4.1 ± 0.6 m s⁻¹ for TOI-1736 b and K_p = 6.0 ± 1.0 m s⁻¹ for TOI-2141 b. The TOI-1736 RVs also revealed an outer giant planet TOI-1736 c with an RV semi-amplitude of K_p = 201.1 ± 0.7 m s⁻¹, in addition to a long term linear trend of −0.177 m s⁻¹ d⁻¹. This trend could be due to the perturbation caused by the orbital motion of an outer planetary or stellar companion. We note that this deviation could be caused by the K-dwarf stellar companion that we detected with photometry, although we did not detect its signature in Gaia’s excess astrometric noise. Current SOPHIE data are not sufficient to distinguish at what periodicity and amplitude this variability occurs. Therefore, we have not yet inferred any characterization for this possible companion, and we will continue to monitor this star in the coming years to discover the nature of this signal in future work.

The inner planets, TOI-1736 b and TOI-2141 b, with radii of 2.44 ± 0.18 R_⊕ and 3.05 ± 0.23 R_⊕ and masses of 12.8 ± 1.8 M_⊕ and 24 ± 4 M_⊕ have similar mean densities of 4.9 ± 1.3 g cm⁻³ and 4.6 ± 1.3 g cm⁻³, respectively. The giant planet TOI-1736 c has a mass of ${8.7}_{-0.6}^{+1.5}$ M_Jup and orbits the star with a period of 570.1 ± 0.7 days and an orbital eccentricity of e = 0.362 ± 0.003. With a semi-major axis of 1.381 ± 0.017 au, this giant resides in a temperate orbit around the star, which is in the habitable zone of that star as detailed in Sect. 6.1. TOI-1736 is therefore similar to the system HD 137496 (K2-364; Azevedo Silva et al. 2022), which is also a solar analog hosting a hot, dense inner planet and an outer giant with about the same mass as TOI-1736 c and with similar orbital properties. Both are nearby and bright systems, making them interesting cases to study the evolution of planetary systems around solar analogs and, in particular, studying the occupation of their habitable zones by giant planets.

Fig. 15 TOI-2141 SOPHIE RVs. In the top panel, the red points show the SOPHIE RV data and the green line shows the best-fit model for the orbit of TOI-2141 b. The bottom panel shows the residuals.

Fig. 16 TOI-2141 SOPHIE RVs in phase with the orbit of the planet. The gray points show the RVs phase-folded in the period of planet TOI-2141 b, and the red points show the binned data, where we use a bin size of 0.05 in units of orbital phase. The green lines represent the best-fit orbit model for TOI-2141 b.

Fig. 17 Astrometric constraint to the orbital inclination of TOI-1736 c. The blue points show the excess noise data obtained from GAS-TON simulations, assuming the RV-derived parameters of TOI-1736 c. The red line shows the Gaia DR3 excess noise of TOI-1736, ɛDR3 = 0.09221 mas.

Fig. 18 TOI-1736 c mass determination. The black line shows the mass posterior distribution obtained from RV-derived parameters and GAS-TON simulations of astrometric excess noise. The solid and dashed blue lines show the median and the ±1σ range of the distribution. The red line shows the minimum mass of Mp sin ip = 8.09 MJup obtained from our RV analysis and the dotted blue line shows the 3σ upper limit of <16.54 MJup.

6.1 The exoplanets’ radiation environment

We consider the best-fit parameters from our analysis to calculate the habitable zone for these solar analogs, using the equations and data from Kopparapu et al. (2014). We obtained an optimistic lower limit (recent Venus) at 0.87 au for TOI-1736 and 0.71 au for TOI-2141, and an upper limit (early Mars) at 2.05 au for TOI-1736 and 1.68 au for TOI-2141. The runaway greenhouse limits (M_p = 1 M_⊕) ranges between 1.10 au and 1.94 au for TOI-1736, and between 0.90 au and 1.59 au for TOI-2141.

To estimate the level of radiation that each planet receives from its host star, we calculate the effective insolation flux S_eff relative to the flux incident on Earth S_⊕ as $S_{e\text{ff}}/S_{\oplus }=L_{{}_{★}}/a_{\text{p}}^2$, where L_★ is the star luminosity in solar units and a_p is the semi-major axis of the planet’s orbit in au. Figure 19 shows the insolation for exoplanets as a function of star effective temperature. TOI-1736 b and TOI-2141 b receive an insolation of 246 ± 26 S_⊕ and 50 ± 4 S_⊕, respectively, showing that these planets reside in a hot environment. We estimate the equilibrium temperature for the planets as in Heng & Demory (2013), assuming a uniform heat redistribution and an arbitrary geometric albedo of 0.1, which gives T_eq = 1076 ± 39 K and T_eq = 722 ± 23 K for TOI-1736 b and TOI-2141 b, respectively.

On the other hand, TOI-1736 c resides at a larger distance from the star, in an eccentric orbit with e_c = 0.362 ± 0.003, and a semi-minor axis of 1.287 ± 0.016 au and a semi-major axis of 1.381 ± 0.017 au. At these distances, TOI-1736 c receives an insolation of 0.81 ± 0.09 S_⊕ and 0.71 ± 0.08 S_⊕ for which we estimate an equilibrium temperature of T_eq = 258 ± 9 K and T_eq = 249 ± 9 K, respectively. This planet resides within the conservative habitable zone around the star.

Fig. 19 Exoplanet insolation compared to the incident flux on Earth. The black circles show the exoplanet data from exoplanets.eu where the sizes are proportional to the planet’s radius, the red stars show the two planets in the TOI-1736 system, and the green star shows planet TOI-2141 b. Light green dashed lines show Earth values. The insolation for both transiting planets is very high compared to Earth, resulting in an equilibrium temperature of 1076 K and 722 K for TOI-1736 b and TOI-2141 b, respectively. The giant planet TOI-1736 c receives insolation very close to Earth values, placing it in the habitable zone around the star.

6.2 Mass-radius relation for the inner planets

From mass and radius measurements of the transiting planets TOI-1736 b and TOI 2141 b, we derive their mean densities of 4.9 ± 1.3 g cm⁻³ and 4.6 ± 1.3 g cm⁻³, respectively, showing that these sub-Neptune sized planets have densities similar to the rocky planets of the Solar System. To explore further the possible chemical composition and structure of these planets we inspect the mass-radius (M-R) diagram illustrated in Fig. 20, which shows the exoplanet data from exoplanets.eu, the data from the two hot sub-Neptunes studied here, and a set of models from Zeng et al. (2019), for comparison. The M-R location of TOI-1736 b matches two models: (1) a 50% Earth-like rocky core plus a 50% H₂O layer; (2) an almost pure Earth-like rocky core with a thin 0.1% H₂ envelop. Therefore, this planet probably has a large dense rocky core surrounded by either a thick water-rich atmosphere or a thin H₂ atmosphere. TOI-2141 b is compatible with three scenarios: (1) a 50% Earth-like rocky core with a 50% H₂O layer; (2) a 49.95% Earth-like rocky core + 49.95% H₂O layer and a thin 0.1% H₂ envelop; and (3) a pure H₂O. Thus, TOI-2141 b is likely a water-rich planet. A further investigation of the internal structure of these planets can be done by including simulations constrained by the star metallicity and abundances of volatiles, which is out of the scope of this paper.

Fig. 20 Mass-radius diagram. The green circles show the exoplanet data from exoplanets.eu, the black stars show the Solar System planets and the blue and orange points show the measured values for TOI-1736 b and TOI-2141 b, respectively. Several models from Zeng et al. (2019) are also plotted for comparison.

7 Conclusions

We reported the discovery of two new planetary systems around the solar analogs TOI-1736 and TOI-2141. We monitored these systems spectroscopically with the SOPHIE instrument at OHP to obtain high-precision RV time series. These data show periodic variations at the same period as the transiting candidates detected by TESS, establishing the planetary nature of the sub-Neptunes TOI-1736 b and TOI-2141 b. Our RV data also revealed an outer giant planet in the TOI-1736 system, in addition to a long-term linear trend. This drift could be due to an additional companion that will be investigated further in future observations.

We characterized the two host stars by analyzing their high-resolution spectra combined with photometry and high-contrast imaging data. We found that TOI-1736 likely has a K3 dwarf stellar companion, which is evidenced both by the power spectrum of the speckle imaging data at 832 nm and by an infrared excess in the SED. The presence of this companion has an important impact on the determination of the star’s radius, which affects the derivation of the planet’s radius and other physical properties of the system. The flux of this possible companion has a negligible contribution to our high-resolution SOPHIE spectra.

We obtained the stellar parameters from several methods, including a SED fit analysis of the photometry data, a spectral analysis of the high-resolution SOPHIE and TRES spectra, and a differential analysis of the SOPHIE spectra using a solar spectrum from observations of the Moon as a reference. The latter provides the most accurate stellar parameters, which also improves the determination of the physical parameters of the planets. This shows the importance of studying planetary systems around solar analogs, where the solar spectrum can benchmark measurements of stellar properties and therefore provide more accurate values. We also determined elemental abundances for several species, including refractory elements (Fe, Mg, and Si) that provide constraints for the planet’s interior structure models. We found that both systems exhibit low levels of stellar activity from their log R′_HK, which are consistent with their non-variable time series of activity indices, and with the mature ages of about 5 Gyr for TOI-1736 and 7 Gyr for TOI-2141, found from chemical clocks and evolutionary analyses. The inferred radii and stellar masses for these two stars have values that differ by about ±10% of solar values, showing that they are indeed solar analogs, but not solar twins.

Therefore, these two new systems represent perspectives on how planetary systems can form and evolve around solar analogs. On the one hand, we have TOI-2141, an evolved star that hosts a warm sub-Neptune with a likely dense rocky core and a possible thick, water-rich envelope. On the other hand, we have TOI-1736, a hierarchical star system with a solar analog that hosts a hot, dense sub-Neptune and an outer giant planet in an eccentric orbit in the habitable zone.

Acknowledgements

E.M. acknowledges funding from Fundação de Amparo à Pesquisa do Estado de Minas Gerais (FAPEMIG) under project number APQ-02493-22 and research productivity grant (PQ) number 309829/2022-4 awarded by the National Council for Scientific and Technological Development (CNPq), Brazil. We thank the Observatoire de Haute-Provence (CNRS) staff for their support. This work was supported by the “Programme National de Planétologie” (PNP) of CNRS/INSU and CNES. This paper includes data collected with the TESS mission, obtained from the MAST data archive at the Space Telescope Science Institute (STScI). Funding for the TESS mission is provided by NASA’s Science Mission Directorate. We acknowledge the use of public TESS data from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center. Resources supporting this work were provided by the NASA High-End Computing (HEC) Program through the NASA Advanced Supercomputing (NAS) Division at Ames Research Center for the production of the SPOC data products. This publication makes use of The Data & Analysis Center for Exoplanets (DACE), which is a facility based at the University of Geneva (CH) dedicated to extrasolar planets data visualization, exchange, and analysis. DACE is a platform of the Swiss National Centre of Competence in Research (NCCR) PlanetS, federating the Swiss expertise in Exoplanet research. The DACE platform is available at https://dace.unige.ch. We acknowledge funding from the French ANR under contract number ANR18CE310019 (SPlaSH). X.D. and T.F. acknowledge support by the French National Research Agency in the framework of the Investissement d’Avenir program (ANR-15-IDEX-02), through the funding of the “Origin of Life” project of the Grenoble-Alpes University. This work was supported by FCT – Fundação para a Ciência e a Tecnologia through national funds and by FEDER through COMPETE2020 – Programa Operacional Competitividade e Internacionalizacão by these grants: UID/FIS/04434/2019, UIDB/04434/2020, UIDP/04434/2020, PTDC/FIS-AST/32113/2017 & POCI-01-0145-FEDER-032113, PTDC/FIS-AST/28953/2017 & POCI-01-0145-FEDER-028953, PTDC/FIS-AST/28987/2017 & POCI-01-0145-FEDER-028987. N.C.S. further acknowledges funding by the European Union (ERC, FIERCE, 101052347). Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Research Council. Neither the European Union nor the granting authority can be held responsible for them. This work makes use of observations from the LCOGT network. Part of the LCOGT telescope time was granted by NOIRLab through the Mid-Scale Innovations Program (MSIP). MSIP is funded by NSF. KAC acknowledges support from the TESS mission via subaward s3449 from MIT. S.D. is funded by the UK Science and Technology Facilities Council (grant number ST/V004735/1). S.G.S. acknowledges the support from FCT through the contract no.CEECIND/00826/2018 and POPH/FSE (EC).

Appendix A SOPHIE RVs and spectroscopic quantities

This appendix presents a compilation of all quantities derived from the SOPHIE spectra. Tables A.1 and A.2 show these quantities obtained for the spectra of TOI-1736 and TOI-2141, respectively.

Appendix B Calibration of our spectral differential analysis

The resulting stellar parameters obtained by our strictly differential analysis presented in Section 3.4 have typical internal precision of σ(T_eff) = 8 K, σ(log 𝑔) = 0.01 dex, σ([Fe/H]) = 0.01 dex, and σ(υ_t) = 0.01 km s⁻¹. To obtain the accuracy of these parameters, we conducted a series of validation tests, which are presented as supplemental material in this appendix. As a first sanity check, we visually inspect the similarity between the stellar and solar spectra for our targets, which is illustrated in Figure B.1.

Fig. B.1 Visual spectral differences between stellar and solar spectra. Solid red lines show the solar spectrum used in this work as the reference spectrum for the strictly differential analysis. The black solid lines show the stellar spectrum for TOI-1736 (upper panel) and TOI-2141 (bottom panel), and the blue points show the respective spectral differences between stellar and solar normalized fluxes.

To estimate the errors introduced by instrumental instabilities, we calculated a stacked spectrum from the mean of a randomly selected subsample containing half of the available spectra, where we performed the differential analysis. We repeat this procedure 30 times. As a result, our EWs measurements are stable at 1.5 ± 0.5% level, which translates into 7 K, 0.02 dex and 0.01 dex in T_eff, log 𝑔 and [Fe/H]. These are considered the minimum uncertainties achievable by our automatic differential technique applied to the SOPHIE data. Thus, we added these errors quadratically to the internal errors of our measurements.

To estimate possible offsets in our spectroscopic parameters and validate the error analysis, we applied the same spectroscopic analysis to observations of the Moon and to a small subset of 5 previously characterized solar analogs bracketing our sample stars’ parameters, namely: HIP 7585, HIP 28066, HIP 29432, HIP 77052, and HIP 79672. The comparison stars’ literature parameters were obtained from Spina et al. (2018) (hereafter S18) which was based on the same line list of this work and used a similar spectroscopic technique. The differences between our analysis and those from S18 are the higher S/N (500-1000) and the higher spectral resolution (R=120000) from their HARPS data. In addition, the authors performed manual EW line measurements, in contrast to our automatic procedure. Therefore, in terms of internal achievable error levels, we can assume their measurements as noise-free in comparison to ours. In Table B.1 we report the comparison between SOPHIE and S18 atmospheric parameters. The typical precision achieved by S18 are σ(T_eff) = 4 K, σ(log 𝑔) = 0.012 dex, and σ([Fe/H]) = 0.004 dex. Our results agree with their estimates with an average offset of ∆(T_eff) = 5 ± 10 K, ∆(log 𝑔) = +0.00 ± 0.02 dex, and ∆([Fe/H]) = +0.00 ± 0.02 dex. Yana Galarza et al. (2019) also validated their measurements with those in the literature using the same technique as in our work and using instrumentation with similar spectral resolution and higher signal-to-noise ratio (R=60000, S/N~400) compared to our SOPHIE data. The authors reported an rms of 12 K, 0.03 dex and 0.02 dex for T_eff, log 𝑔, and [Fe/H], respectively. In a similar fashion, S18 found an rms of σ(T_eff) = 15 K, σ(log 𝑔) = 0.03 dex, and σ([Fe/H]) = 0.017 dex. These results are in line with our differential atmospheric parameters, highlighting the robustness of our methodology applied to similar Sun-like stars.

Finally, to evaluate the accuracy of our methodology, we derived the solar atmospheric parameters from another moonlight spectrum obtained with the same instrumental setup reported in this work. Solar observations were analyzed following the same steps used to determine the parameters of solar analogs. We show in Table B.2 the accuracy test using moonlight observations and final differential atmospheric parameters for TOI-1736 and TOI-2141. The agreement with our solar reference values used in q2 confirms the differential analysis as accurate and precise, assuming there is negligible unaccounted-for residual systematics from the reduction process, for example.

Appendix C Activity indices

This appendix details our S-index and H_α measurements. It also shows the time series analysis and RV correlation of the activity indices obtained from the SOPHIE spectra of TOI-1736 and TOI-2141.

C.1 S-index

We computed the SOPHIE’s instrumental S-index as follows: $$S_{\text{SOPHIE}}=\frac{N_H+N_K}{N_V+N_R},$$(C.1)

where N_H and N_K are the integrated fluxes of two triangular band-passes in a 0.109 nm wide window centered on the H (396.85 nm) and K (393.37 nm) Ca II emission lines, and N_V and N_R are two continuum regions 2 nm wide centered at 390.107 nm and 400.107 nm, as illustrated in Figures C.3 and C.4. For the computation of S_SOPHIE, we performed flux integration using spectral values generated via Monte Carlo (MC) sampling from a normal distribution, taking into account the central value and errors associated with each spectral element. Figures C.3 and C.4 depict the posterior distributions for the MC samples.

To calibrate SOPHIE’s instrumental S-index to the Mt. Wilson Observatory (MWO) system (S_MW, Wilson 1968; Egeland et al. 2017), we downloaded archival spectra of ten selected stars that are known to have S_MW measurements in the range of 0.15 to 0.36. We only used observations in HR mode and with a peak S/N>100 in the first order, and we computed the S-index for each calibrator in the same way as described above. We least-square fit a linear relationship between S_SOPHIE and the S_MW values from Duncan et al. (1991), and then we sampled the posterior distribution of the linear coefficients using a Bayesian MCMC framework with the package emcee (Foreman-Mackey et al. 2013) as shown in figures C.1 and C.2. The best solution is expressed as follows: $$S_{\text{SOPHIE}}=0.832\pm 0.104\times S_{\text{MW\hspace{0.17em}}}+0.065\pm \mathrm{0.026.}$$(C.2)

This calibration is consistent with the previous calibration obtained by Boisse et al. (2010). The S_MW values from Duncan et al. (1991) and our measurements of S_SOPHIE for our sample of calibrators are listed in Table C.1. Solving Equation C.2 to get the S_MW for our two stars, we obtained $S_{\text{\hspace{0.17em}MW}}={0.16}_{-0.03}^{+0.04}$ for TOI-1736 and $S_{\text{\hspace{0.17em}MW}}={0.21}_{-0.03}^{+0.04}$ for TOI-2141, where the S_SOPHIE values were measured in the template spectrum of each target. We have also measured the S_MW for each individual exposure, providing the S-index time series for both stars as presented in Tables A.1 and A.2.

Fig. C.1 Calibration of S-index from the SOPHIE HR spectra. Black points show the S-index for several targets (see Table measured from the SOPHIE HR spectrum as a function of the reference MW values. The orange lines show a randomly selected sample of the fit models obtained from a posterior distribution of the model parameters, and the dashed black line show the median linear model.

Fig. C.2 Pairs plot showing the MCMC samples and posterior distributions of the coefficients for the SOPHIE S-index calibration.

C.2 Hα

The Hα index is calculated from the flux at the center of the Hα line as described by Boisse et al. (2009), that is, $$H_{\alpha }=\frac{F_{H_{\alpha }}}{F_1+F_2},$$(C.3)

where $F_{H_{\alpha }}$ is the integrated flux measured within a 0.068 nm window centered at 656.2808 nm, F₁ and F₂ are the integrated fluxes measured near the edge of the Hα wings centered at 655.087 nm and 658.031 nm with rectangular band-passes in 1.075 nm and 0.875 nm wide windows, as illustrated in Figures C.3 and C.4. In a similar way to the S-index calculation, the integration of spectral regions for Hα calculation employs a Monte Carlo sampling approach. The Hα indices measured in the template spectra are H_α = 0.2376 ± 0.0016 for TOI-1736 and H_α = 0.2488 ± 0.0023 for TOI-2141, and the individual measurements per exposure are listed in Tables A.1 and A.2.

Fig. C.3 TOI-1736 measurements of S-index and Hα from SOPHIE spectra. The top panels show the SOPHIE template spectrum around the S-index and Hα regions, showing the integration windows for the core of lines and continuum regions. The bottom panels show the posterior distribution obtained from the Monte Carlo samples calculated in the measurement of each quantity.

Fig. C.4 TOI-2141 measurements of S-index and Hα from SOPHIE spectra. The top panels show the SOPHIE template spectrum around the S-index and Hα regions, showing the integration windows for the core of lines and continuum regions. The bottom panels show the posterior distribution obtained from the Monte Carlo samples calculated in the measurement of each quantity.

C.3 Time series of activity indices

As star rotation can modulate RVs and activity indices, we can analyze these quantities together and compare their time series and respective GLS periodograms in an attempt to find possible signals that can be used to infer the star rotation period. Figures C.5 and C.6 show the time series data and the GLS periodograms for both stars, where we cannot identify any significant periodicity. The residual RVs of TOI-1736 show a peak with a false alarm probability greater than 0.001 at 25 d, which may be related to star rotation, but this peak does not appear in any other indices. As such, we could not find any strong evidence for the star’s rotation period from these data.

Fig. C.5 Time series and respective GLS periodograms for the TOI-1736 SOPHIE RVs subtracted from the best-fit model and for the following activity indices: CCF FWHM, bisector span, Hα, and S-index.

Fig. C.6 Time series and respective GLS periodograms for the TOI-2141 SOPHIE RVs subtracted from the best-fit model and for the following activity indices: CCF FWHM, bisector span, Hα, and S-index.

C.4 Correlations between activity indices and RVs

Stellar activity can change the shape of line profiles and lead to deterioration in the detection of planetary signals, or even to spurious detections. One way to inspect the impact of certain types of activity on RV measurements is to look at the correlation between the RVs subtracted from the best-fit RV orbit model and the activity indices. Here, we present these correlations for the CCF FWHM, bisector span, H_α, and S-index, as illustrated in Figures C.7 and C.8. As pointed out in the main text, there is no strong correlation between these quantities, which indicates that our RVs are not strongly affected by stellar activity.

Fig. C.7 Correlations between the TOI-1736 SOPHIE RVs subtracted from the best-fit model and the activity indices CCF FWHM, bisector span, Hα, and S-index.

Fig. C.8 Correlations between the TOI-2141 SOPHIE RVs subtracted from the best-fit model and the activity indices CCF FWHM, bisector span, Hα, and S-index.

Appendix D Priors and posterior distributions of model parameters

This appendix presents in Table D.1 the prior distributions adopted for each model parameter that we considered in our analysis. Figures D.1 and D.2 illustrate the MCMC samples and the final posterior distributions of free parameters used in our analysis.

Fig. D.1 Pairs plot showing the MCMC samples and posterior distributions of the free parameters in our joint analysis of the TESS photometry and the SOPHIE RV data of TOI-1736. The contours mark the 1σ, 2σ, and 3σ regions of the distribution. The blue crosses indicate the best-fit values for each parameter and the dashed vertical lines in the projected distributions show the median values and the 1σ uncertainty (34% on each side of the median).

Fig. D.2 Same as Figure D.1 but for the analysis of TOI-2141 data.

References

All Tables

All Figures

Fig. 1 TESS light curve of TOI-1736. The gray points show the TESS photometry data normalized by the median flux of each sector. The black points show the same data binned by the weighted average with bin sizes of 0.1 d. The red line shows the baseline GP model that was fitted to the binned data multiplied by the best-fit transit model for TOI-1736 b. The red triangles show the predicted central time of TOI-1736 b transits.

In the text

	Fig. 2 TESS light curve of TOI-2141. The blue points show the TESS photometry data around the transits of TOI-2141 b. The red line shows the baseline GP model that was fitted to the binned data multiplied by the best-fit transit model for TOI-2141 b. The red triangles show the predicted central time of TOI-2141 b transits.
In the text

	Fig. 3 Contrast ratio of TOI-1736 as a function of angular separation at 562 nm (blue line) and at 832 nm (red line) obtained from the ‘Alopeke/Gemini speckle imaging observations. The small panels show the reconstructed speckle image (right subpanel) and its power spectrum in Fourier space (left subpanel). The power spectrum shows fringes that match a binary model.
In the text

	Fig. 4 Contrast ratio of TOI-2141 as a function of angular separation at I-band obtained on the 4.1-m SOAR speckle imaging observations.
In the text

Fig. 5 SED fit analysis for TOI-1736. The upper panel shows results assuming a single stellar component and the lower panel shows results assuming two stellar components. The red symbols represent the observed photometric measurements, whereas the horizontal bars represent the effective passband width. The blue symbols are the model fluxes from the best-fit NextGen atmosphere model (black and gray lines). The lower panel also shows the model for the primary component (cyan line) and the secondary component (red line), separately. The insets display our model (black curve) and the observed absolute flux-calibrated Gaia spectrum (gray swathe).

In the text

	Fig. 6 SED fit analysis for TOI-2141. Red symbols represent the observed photometric measurements, whereas the horizontal bars represent the effective width of the passband. Blue symbols are the model fluxes from the best-fit NextGen atmosphere model (black). The inset displays our model (black curve) alongside the observed absolute flux-calibrated Gaia spectrum (gray swathe).
In the text

Fig. 7 Comparison between solar and stellar EW measurements. TOI-1736 and TOI-2141 are shown in the top and bottom panels, respectively. The dashed black line is the 1:1 EWs ratio, and the solid black line is a linear fit based on the Fe I lines, represented by empty black circles. The red star symbols represent the Fe II lines, followed by their linear fit represented by the solid red line. The typical standard deviation of the average linear fit is around ~0.1–0.2 nm, which means ~2% in relative terms.

In the text

	Fig. 8 Spectral synthesis analysis around the Li region to obtain an estimate of the Li abundance. The points show the SOPHIE spectra of TOI-1736 (in red), TOI-2141 (in green), and a solar spectrum (in blue) from observations of the Moon. The solid lines show the respective bestfit synthetic spectra obtained from our analysis.
In the text

	Fig. 9 Generalized Lomb-Scargle periodogram for the TOI-1736 SOPHIE RVs. The top and bottom panels show the GLS for the SOPHIE RVs minus the linear trend and the fit orbits of planets c (top panel) and b (bottom panel).
In the text

	Fig. 10 Generalized Lomb–Scargle periodogram for the TOI-2141 SOPHIE RVs. The blue dashed line shows the peak power at 18.3 days, which agrees with the period reported by TESS for TOI-2141 b.
In the text

	Fig. 11 Transits of TOI-1736 b observed by TESS. The blue points show the TESS photometry data around the nineteen transits of TOI-1736 b, with the times relative to the central time of each transit. The red line shows the best-fit transit model and the green points show the residuals plus an arbitrary offset for better visualization.
In the text

	Fig. 12 Transits of TOI-2141 b observed by TESS. The blue points show the TESS photometry data around the four transits of TOI-2141 b, with the times relative to the central time of each transit. The red line shows the best-fit transit model and the green points show the residuals plus an arbitrary offset for better visualization.
In the text

	Fig. 13 TOI-1736 SOPHIE RVs. In the top panel, the red points show the SOPHIE RV data and the green line shows the best-fit model, including the orbits of TOI-1736 b and TOI-1736 c, plus a linear trend, which is also represented separately by the red dashed line. The bottom panel shows the residuals.
In the text

Fig. 14 TOI-1736 SOPHIE RVs in phase with the orbits of the planets. In the top panel, the gray points show the RVs phase-folded in the period of planet b, where the orbit of planet c and the linear trend have been subtracted. In the bottom panel, the gray points show the RVs phase-folded in the period of planet c, where the orbit of planet b and the linear trend has been subtracted. The red points show the binned data, where we use a bin size of 0.05 in units of orbital phase. The green lines represent the best-fit orbit model for each respective planet.

In the text

	Fig. 15 TOI-2141 SOPHIE RVs. In the top panel, the red points show the SOPHIE RV data and the green line shows the best-fit model for the orbit of TOI-2141 b. The bottom panel shows the residuals.
In the text

	Fig. 16 TOI-2141 SOPHIE RVs in phase with the orbit of the planet. The gray points show the RVs phase-folded in the period of planet TOI-2141 b, and the red points show the binned data, where we use a bin size of 0.05 in units of orbital phase. The green lines represent the best-fit orbit model for TOI-2141 b.
In the text

	Fig. 17 Astrometric constraint to the orbital inclination of TOI-1736 c. The blue points show the excess noise data obtained from GAS-TON simulations, assuming the RV-derived parameters of TOI-1736 c. The red line shows the Gaia DR3 excess noise of TOI-1736, ɛDR3 = 0.09221 mas.
In the text

Fig. 18 TOI-1736 c mass determination. The black line shows the mass posterior distribution obtained from RV-derived parameters and GAS-TON simulations of astrometric excess noise. The solid and dashed blue lines show the median and the ±1σ range of the distribution. The red line shows the minimum mass of Mp sin ip = 8.09 MJup obtained from our RV analysis and the dotted blue line shows the 3σ upper limit of <16.54 MJup.

In the text

Fig. 19 Exoplanet insolation compared to the incident flux on Earth. The black circles show the exoplanet data from exoplanets.eu where the sizes are proportional to the planet’s radius, the red stars show the two planets in the TOI-1736 system, and the green star shows planet TOI-2141 b. Light green dashed lines show Earth values. The insolation for both transiting planets is very high compared to Earth, resulting in an equilibrium temperature of 1076 K and 722 K for TOI-1736 b and TOI-2141 b, respectively. The giant planet TOI-1736 c receives insolation very close to Earth values, placing it in the habitable zone around the star.

In the text

	Fig. 20 Mass-radius diagram. The green circles show the exoplanet data from exoplanets.eu, the black stars show the Solar System planets and the blue and orange points show the measured values for TOI-1736 b and TOI-2141 b, respectively. Several models from Zeng et al. (2019) are also plotted for comparison.
In the text

Fig. B.1 Visual spectral differences between stellar and solar spectra. Solid red lines show the solar spectrum used in this work as the reference spectrum for the strictly differential analysis. The black solid lines show the stellar spectrum for TOI-1736 (upper panel) and TOI-2141 (bottom panel), and the blue points show the respective spectral differences between stellar and solar normalized fluxes.

In the text

	Fig. C.1 Calibration of S-index from the SOPHIE HR spectra. Black points show the S-index for several targets (see Table measured from the SOPHIE HR spectrum as a function of the reference MW values. The orange lines show a randomly selected sample of the fit models obtained from a posterior distribution of the model parameters, and the dashed black line show the median linear model.
In the text

	Fig. C.2 Pairs plot showing the MCMC samples and posterior distributions of the coefficients for the SOPHIE S-index calibration.
In the text

	Fig. C.3 TOI-1736 measurements of S-index and Hα from SOPHIE spectra. The top panels show the SOPHIE template spectrum around the S-index and Hα regions, showing the integration windows for the core of lines and continuum regions. The bottom panels show the posterior distribution obtained from the Monte Carlo samples calculated in the measurement of each quantity.
In the text

	Fig. C.4 TOI-2141 measurements of S-index and Hα from SOPHIE spectra. The top panels show the SOPHIE template spectrum around the S-index and Hα regions, showing the integration windows for the core of lines and continuum regions. The bottom panels show the posterior distribution obtained from the Monte Carlo samples calculated in the measurement of each quantity.
In the text

	Fig. C.5 Time series and respective GLS periodograms for the TOI-1736 SOPHIE RVs subtracted from the best-fit model and for the following activity indices: CCF FWHM, bisector span, Hα, and S-index.
In the text

	Fig. C.6 Time series and respective GLS periodograms for the TOI-2141 SOPHIE RVs subtracted from the best-fit model and for the following activity indices: CCF FWHM, bisector span, Hα, and S-index.
In the text

	Fig. C.7 Correlations between the TOI-1736 SOPHIE RVs subtracted from the best-fit model and the activity indices CCF FWHM, bisector span, Hα, and S-index.
In the text

	Fig. C.8 Correlations between the TOI-2141 SOPHIE RVs subtracted from the best-fit model and the activity indices CCF FWHM, bisector span, Hα, and S-index.
In the text

Fig. D.1 Pairs plot showing the MCMC samples and posterior distributions of the free parameters in our joint analysis of the TESS photometry and the SOPHIE RV data of TOI-1736. The contours mark the 1σ, 2σ, and 3σ regions of the distribution. The blue crosses indicate the best-fit values for each parameter and the dashed vertical lines in the projected distributions show the median values and the 1σ uncertainty (34% on each side of the median).

In the text

	Fig. D.2 Same as Figure D.1 but for the analysis of TOI-2141 data.
In the text