© F. Murgas et al. 2022

1 Introduction

Ultra-short-period (USP) planets are defined by the community as planetary-mass objects with orbital periods shorter than one day (Sahu et al. 2006; Sanchis-Ojeda et al. 2013). Of the ~5000 exoplanets discovered to date¹, close to 120 are USPs, and only ~30 of them have both mass and radius measurements. Some of the first discovered USP planets are Corot-7b (Léger et al. 2009, P = 0.85 days), 55 Cnc e (Dawson & Fabrycky 2010; Winn et al. 2011, P = 0.74 days), Kepler-10b (Batalha et al. 2011, P = 0.84 days), and Kepler-78b (Sanchis-Ojeda et al. 2013, P = 0.36 days). The USP planets with the shortest orbital period discovered so far are the planet candidate KOI-1843.03 (Ofir & Dreizler 2013) and K2-137b (Smith et al. 2018); these two planets revolve around M dwarfs in just 0.18 days (~4.3 h).

The first systematic search for transiting USPs was conducted by Sanchis-Ojeda et al. (2014). Using Kepler data, they found 106 USP candidates, 6 of which have P < 6 h. The authors noted that the planetary radii of these objects were rarely larger than 2.0 R_⊕. Winn et al. (2018) also found that the radius distribution of USPs declines sharply around 2 R_⊕, and they proposed that this sharp decline is attributed to photoevaporation of the planet atmosphere that is caused by radiation from the host star.

Although the period cutoff (P < 1 day) is arbitrary, this type of objects appears to have some distinct characteristics that may indicate that they are a distinct population. For example, the occurrence rates of USP planets seem to depend on spectral type. Sanchis-Ojeda et al. (2014) found that the USP occurrence rate falls from 1.1 ± 0.4% for M dwarfs to 0.15 ± 0.05% for K dwarfs. In the same work, they establish that the occurrence rate for G dwarfs is 0.51 ± 0.07%, meaning that about 1 in 200 solar-type stars hosts a USP planet (for planets larger than 0.84 R_⊕).

There are other characteristics that separate USP planets from other types of short-period exoplanets such as hot Jupiters. Winn et al. (2017) found that the metallicity distribution for stars that host hot Jupiters is different from that of stars that ho planets: Hot-Jupiter host stars tend to be more metal rich than the stars that have USP planets. Winn et al. (2018) also pointed out that, unlike hot Jupiters, USP planets can be found in multiple planet systems. When USP planets are found in multiple planet systems, the ratio of the orbital periods between the inner planet and its nearest neighbor is typically higher than 4 (P₂/P₁ ≳ 4) (Steffen & Farr 2013; Winn et al. 2018; Pu & Lai 2019). This period ratio is higher than what is seen in Kepler multiplanet systems (Fabrycky et al. 2014). Additionally, for systems of multiple transiting planets, the dispersion of the orbital inclination of the different transiting planets has been found to be higher when a USP planet is part of the system (Dai et al. 2018).

The origin of USP planets is still an open issue. Their close-in orbit usually places them inside the dust-sublimation region around their host star, meaning that they are unlikely to have formed in their present-day orbits. Hence, in all the proposed scenarios, these objects formed in a wider orbit around their host star and migrated to their current position. Petrovich et al. (2019) proposed that secular interactions in multiplanet systems (e.g., N_planets > 3) can affect the inner planet (with P in the 5–10 days range) and push it into a highly eccentric orbit. It is eventually tidally captured by the star in a short-period orbit that is circularized over time as a result of planetary tides. Pu & Lai (2019) explored another formation mechanism in which a small rocky planet is born with a period of a few days and moderate eccentricity (e ≳ 0.1) in a multiplanetary system; the outer planets tidally interact with each other and with the innermost planet, damping the eccentricity to a value close to zero and shrinking the semi-major axis in a quasi-equilibrium state. At the end of this process, the innermost planet becomes a USP planet, while the second planet, an Earth- or super-Earth-sized planet, stabilizes itself on a 10 day orbit. TOI-500 is the first four-planet system for which this mechanism has been proven to work (Serrano et al. 2022). Millholland & Spalding (2020) proposed an obliquity-driven tidal migration mechanism. There, in a system of planets with strong mutual inclinations, planetary obliquities and tides become excited in a positive-feedback loop that forces inward migration, until a condition is reached in which the high obliquities are tidally destabilized and migration stalls.

Because USP planets are rare, those with well-determined parameters are especially useful for comparison with the predictions made by theoretical models. The NASA-sponsored Transiting Exoplanet Survey Satellite (TESS; Ricker et al. 2014) is a space telescope equipped with four cameras observing an area of sky of 24° × 96° degrees. Its main goal is to search for transiting exoplanets around bright stars (5 < I_C < 13 mag) by observing the same part of the sky almost uninterrupted during ~27 days. Launched in April 2018, TESS is past its original 2-yr mission and is currently in its extended mission, which was approved to continue until the end of September 2022. Although USPs are rare, the TESS ~27-day observation cycles offer a unique opportunity to find this type of objects, and because the mission focused on bright stars, it allows Doppler mass measurements through radial velocity follow-up observations.

Here we report the discovery of HD 20329b, a USP planet around a bright (V = 8.74 mag, J = 7.5 mag) G-type star, which was discovered using TESS data. This paper is organized as follows: in Sect. 2 we describe the TESS data, spectroscopic follow-up, and high resolution imaging of HD 20329. In Sect. 3, we describe the methods we used to determine the stellar parameters, the light curve, and the radial velocity fitting procedure. In Sect. 4, we present the parameters of HD 20329b and place this planet in the context of known USP planets. Section 5 presents the conclusions of this work.

2 Observations

2.1 TESS photometry

HD 20329 (TIC 333657795; Stassun et al. 2018) was observed by TESS from 20 August 2021 until 16 September 2021 (sector 42) and from 16 September 2021 until 12 October 2021 (sector 43). For the stellar coordinates and magnitudes, see Table 1. For sector 42, the target was observed on camera 4 CCD 3, while in sector 43, the star was placed on camera 2 CCD 1. For each TESS sector, the star was observed for ~25 days, and the images were stacked using a 2-min cadence mode.

TESS observes a field continuously for ~27 days. At every orbit perigee (~13 days), science operations are interrupted, and the data are sent to Earth for processing. The raw images were processed by the Science Processing Operations Center (SPOC) at NASA Ames Research Center. The SPOC pipeline (Jenkins et al. 2016) performs image calibration and data-quality control, extracts photometry for all the TESS target stars in the field of view, and searches the extracted light curves for transit signatures.

The data reduction process of TESS time series starts by using simple aperture photometry (SAP: Morris et al. 2020) to generate an initial light curve. Then the Presearch Data Conditioning (PDC) pipeline module (Smith et al. 2012; Stumpe et al. 2014) removes some instrumental systematic effects from the time series. Transit events are searched with the wavelet-based matched filter described in Jenkins (2002) and Jenkins et al. (2020), and are then fit to transit models, including the contribution made by stellar limb-darkening effects (Li et al. 2019). Finally, a set of diagnostic tests are applied to the light curves to establish whether the detected transit events have a planetary origin (Twicken et al. 2018). After this process, the TESS science office reviews the transit signature and promotes the candidate as a TESS object of interest (TOI) if it has a likely planetary origin. In the case of HD 20329, the planetary candidate was assigned the TOI identification TOI-4524.01, and the community was alerted to it in October 2021. The SPOC transit depth and orbital period of TOI-4524.01 were 210 ± 0.60 ppm and P = 0.926014 ± 0.00005 days, making this TOI a USP planet candidate.

The TESS light curves were analyzed independently using the Détection Spécialisée de Transits (DST; Cabrera et al. 2012) pipeline. Variability in the PDC-SAP light curve was first removed using a Savitzky-Golay filter (Savitzky & Golay 1964; Press et al. 1992), and transit searches were performed. A transit signal with an orbital period of 0.92653 ± 0.00011 day and a transit depth of 204 ± 17 ppm was detected, consistent with the signal detected by SPOC.

We used the TESS PDC-SAP light curves for the transit analysis presented in this work. The PDC-SAP curves were corrected for instrumental systematic effects and include some correction for flux contamination from nearby stars. The TESS sector 42 and 43 observations are publicly available at the Barbara A. Mikulski Archive for Space Telescopes (MAST²). Figure 1 shows the TESS target pixel file (TPF) images around HD 20329 for sectors 42 and 43, and the photometric aperture used by TESS is highlighted with red squares.

2.2 Ground-based seeing-limited photometry with MuSCAT2

TESS has a relatively large pixel scale (21″ pixel⁻¹), hence the detected transit event might originate from another star located inside the TESS photometric aperture. Centroid analysis results from the TESS pipeline using images from sectors 42 and 43 combined eliminate the possibility that nearby catalog stars caused the transit signature. In particular, they exclude the 17th mag star 37″ southwest of HD 20329, the 12th mag star 49’″ northwest of it, as well as some very dim stars slightly over 25’″ to the south-southwest.

When TESS announced HD 20329 as an object of interest, the transit depth found by the pipeline was ~0.230 mmag (210 ppm), meaning that the transit event probably cannot be detected on HD 20329 with the medium-sized ground-based telescopes that are typically used to rule out false positives. Nonetheless, it is possible to rule out other stars in the field as causing the transit signal based on seeing-limited photometry.

HD 20329 was observed on the nights of 4 and 5 December 2021 UT with the simultaneous multicolor imager MuSCAT2 (Narita et al. 2019) mounted on the 1.5 m Telescopio Carlos Sánchez (TCS) at Teide Observatory, Spain. MuSCAT2 has four CCDs with 1024 × 1024 pixels, and each camera has a field of view of 7.4′ × 7.4′ (pixel scale of 0.44″ pixel⁻¹). The instrument is capable of taking images simultaneously in g′, r′, i′, and z_s bands with short read-out times. For both nights, the raw data were reduced by the MuSCAT2 pipeline (Parviainen et al. 2019), which performs standard image calibration and aperture photometry and is capable of modeling the instrumental systematics present in the data while simultaneously fitting a transit model to the light curve.

On the night of 4 December 2021 UT, the telescope was defocused to avoid saturation of HD 20329. On the night of 5 December 2021 UT, we also defocused the telescope, but we let the target star saturate in order to detect fainter stars in the field. For both nights, the exposure times were set to 5 s for all the MuSCAT2 bands. For the first night, we were unable to detect the transit on target, likely due to the small transit depth of the event; on the second night HD 20329 was saturated on purpose so that no transit detection could be performed. Inside a field of view with a radius of 2.5′, four stars identified by Gaia (DR3; Gaia Collaboration 2022) can produce a transit signal with the depth detected by TESS. Of these stars, we were only able to rule out the nearby star (TIC 333657797, V = 12.5 mag) as an eclipsing binary causing the transit events. The other three Gaia sources were too faint and are not detected in the MuSCAT2 images.

Fig. 1 TESS target pixel file image of HD 20329 observed in sector 42 (left) and sector 43 (right), made with tpf plotter (Aller et al. 2020). The pixels highlighted in red show the aperture used by TESS to create the light curves. The position and sizes of the red circles represent the position and TESS magnitudes of nearby stars, respectively.

2.3 Ground-based high-resolution imaging

Part of the validation process of transiting exoplanets is the assessment of possible flux contamination by nearby companions (bound or unbound to the target star) and its effect on the derived planetary radius (Ciardi et al. 2015). For this reason, we observed HD 20329 with a combination of high-resolution imaging resources including optical speckle and near-infrared (NIR) adaptive optics (AO). Astrometric data from Gaia EDR3 (Gaia Collaboration 2021) was also used to provide additional constraints on the presence of undetected stellar companions as well as wide companions.

2.3.1 SOAR optical speckle imaging

We observed HD 20329 using speckle imaging on 20 November 2021 UT with the 4.1 m Southern Astrophysical Research (SOAR) telescope (Tokovinin 2018). The observations were made as part of the SOAR TESS survey (see Ziegler et al. 2020 for details). The data were acquired with the HRCam instrument (Cousins I band). The observation reached a sensitive threshold to ensure the detection of stars 5.8 mag. fainter than the target at an angular distance of 1″. No nearby stars around HD 20329 were detected using a search radius of 3″. Figure 2 shows the 5σ detection sensitivity and speckle autocorrelation functions (ACF) from the SOAR observations.

Fig. 2 Fig. 2. Speckle sensitivity curve and ACF (inset) from SOAR HRCam observations of HD 20329. No nearby stars were detected within 3″of HD 20329 with HRCam.

2.3.2 Palomar NIR AO imaging

HD 20329 was observed by the adaptive-optics instrument PHARO (Hayward et al. 2001; field of view of ~25″, pixel scale of 0.025″ per pixel) mounted on the 5 m Hale telescope on 11 November 2021. The observations were made using the narrow-band NIR Br − γ filter and applying a five-point dither pattern. The raw data were reduced using standard procedures (flat, dark, and sky calibration; calibrated science images were coadded). The 5σ sensitive curve presented in Fig. 3 was determined using injection of artificial sources following the procedure outlined in (Furlan et al. 2017). No stellar companions were detected in PHARO observations.

Fig. 3 Palomar NIR AO imaging and sensitivity curves for HD 20329 taken in the Bry filter. The images were taken in good seeing conditions, and we reach a contrast of 7 mag. fainter than the host star within 0.5″. Inset: image of the central portion of the data, centered on the star.

2.3.3 Gaia assessment

We used Gaia EDR3 measurements to search for stellar companions around HD 20329 following Luque et al. (2022; see Sect. 3.2.3 of Luque et al. 2022 and references therein). There are no nearby stars around the position HD 20329 with similar astrometric properties (i.e., parallaxes and proper motions; Mugrauer & Michel 2020, 2021) that could indicate that they are bound to the host star. Additionally, the Gaia renormalized unit weight error (RUWE), a metric used to measure the level of astrometric noise caused by a gravitationally bound unseen companion (e.g., Ziegler et al. 2020), is consistent with a single-star model for the case of HD 20329 (EDR3 RUWE = 0.92, below the 1.4 RUWE threshold for a multiple star system).

2.4 Spectroscopic observations

Between 29 November 2021 (UT) and 30 January 2022 (UT), we collected 120 spectra with the High Accuracy Radial velocity Planet Searcher for the Northern hemisphere (HARPS-N: λ ∈ (378−691) nm, R ≈ 115 000, Cosentino et al. 2012) mounted at the 3.58 m Telescopio Nazionale Galileo (TNG) of Roque de los Muchachos Observatory in La Palma, Spain, under the observing program CAT21A_119. The exposure time was set to 259−1800 s, based on weather conditions and scheduling constraints, leading to a S/N per pixel of 20−136 at 5500 Å. The spectra were extracted using the off-line version of the HARPS-N data reduction software (DRS) pipeline (Cosentino et al. 2014), version 3.7. Doppler measurements and spectral activity indicators (CCF_BIS, CCF_FWHM, CCF_CTR and Mont-Wilson S-index, S_MW) were measured using an online version of the DRS, the YABI tool³, by cross-correlating the extracted spectra with a G2 mask (Baranne et al. 1996). We also used serval code (Zechmeister et al. 2018) to measure relative RVs by the template-matching, chromatic index (CRX), differential line width (dLW), and Ha, sodium Na D1 & Na D2 indexes. The uncertainties in the relative RVs measured with serval are in the range 0.5−4.4 m s⁻¹, with a mean value of 1.2 m s⁻¹. The uncertainties in the absolute RVs measured with the online version of DRS (YABI) are in a the range 0.6–6.1 m s⁻¹, with a mean value of 1.3 m s⁻¹. Table Al gives the time stamps of the spectra in BJD_TDB, absolute RVs, and spectral activity indicators (CCF_BIS, CCF_FWHM, CCF_CTR, and Mont-Wilson S-index, S_MW) measured with YABI. Table A.2 gives the relative RVs and spectral activity measurements (CRX, dLW, and Ha, and sodium Na Dl & Na D2 indexes) measured with serval. In the joint RV and transit analysis presented in Sect. 3.2, we used relative RVs measured from HARPS-N spectra with serval by the template-matching technique.

Figure 4 shows the generalized Lomb-Scargle (GLS; Zechmeister & Kürster 2009) periodogram for the radial velocities and activity indices. The false-alarm probability (FAP) levels were computed using the bootstrap method (Murdoch et al. 1993; Hatzes 2016) with 10 000 iterations. The peak of the RVs GLS periodogram is located at P = 0.9263 days, in agreement with the orbital period of the transiting candidate announced by TESS.

Fig. 4 GLS periodogram (Zechmeister & Kürster 2009) for HARPS-N radial velocity measurements, (α) Radial velocity measurement residuals after fitting the transiting-planet signal, (b) Stellar activity indices of HD 20329 ((c) to (k)). The peak with the maximum power in the periodogram for the RV measurements corresponds to a period of P = 0.9263 days (vertical green line), in agreement with the period of the transiting-planet candidate reported by TESS. The dash-dotted blue and segmented red line represent the 0.01% and 0.1% FAP levels, respectively.

3 Analysis

3.1 Stellar parameters

The stellar parameters for HD 20329 are presented in Table 1. The analysis of the coadded HARPS-N stellar spectrum was carried out by using the BACCHUS code (Masseron et al. 2016), which relies on the MARCS model atmospheres (Gustafsson et al. 2008) and atomic and molecular line lists from Heiter et al. (2021). In brief, the surface gravity was determined by requiring ionization balance of the Fe I lines and the Fe II line. A microturbulence velocity was also derived by requiring no trend of Fe line abundances against their equivalent widths. The output metallicity is represented by the average abundance of the Fe I lines. An effective temperature of 5596 ± 50 K was derived by requiring no trend of the Fe I lines abundances against their respective excitation potential.

We compared our stellar parameters to those found in the literature in the Hypatia Catalog (Hinkel et al. 2014). HD 20329 was spectroscopically observed by Ramírez et al. (2013) and Bedell et al. (2017). The effective temperatures and iron-content determined here are in excellent agreement with both Ramírez et al. (2013) and Bedell et al. (2017), especially when [Fe/H] is normalized to the same solar scale. While the surface gravities of Ramírez et al. (2013) and Bedell et al. (2017) disagree within their own errors, the value that we determined (4.40 ± 0.07) overlaps with the two literature determinations within errors. The [Fe/H] of HD 20329, in addition to other iron-peak elements (Cr, Mn, and Ni), is somewhat subsolar. However, many elements, especially the α-elements (C, O, Si, Ca, Ti), within HD 20329 are supersolar (Ramírez et al. 2013; Bedell et al. 2017). In combination with these abundance trends, the stellar kinematics indicate that HD 20329 likely originated from the thick disk according to a conservative kinematic prescription by Bensby et al. (2003) and as noted within Bedell et al. (2017).

The stellar rotation (v sin i) was estimated by measuring the average of the Fe line broadening after subtracting the instrument and the thermal, collisional, and microturbulence broadening. This velocity leads to an estimate of the rotational period of 15 ± 3 days (assuming sin i = 1). Based on the Noyes et al. (1984) and Mamajek & Hillenbrand (2008) activity-rotation relations and using (B − V) of 0.670 and the log measured with YABI (−5.03 ± 0.03), we estimated a rotation period of HD 20329 for 28.6 ± 5.8 days and 30.7 ± 3.3 days, respectively. The discrepancy between the two indicators may be explained by the fact that for relatively long rotation periods, its impact on the line broadening becomes very weak and its disentanglement from other sources of broadening becomes more uncertain. Although the empirical calibrations of chromospheric activity and rotational periods have their own issues, we nevertheless favor a period of ~30 days. We tried to establish the rotation period of HD 20329b using long-term photometry data from several surveys. We were only able to access the photometry from the ASAS-SN public light-curve archive⁴ (Shappee et al. 2014; Kochanek et al. 2017), but found no conclusive signs of photometric modulation attributable to the rotation period of the star (see Appendix B).

In a second step, we used the Bayesian tool PARAM (Rodrigues et al. 2014, 2017) to derive the stellar mass, radius, and age using the spectroscopic parameters and the updated Gaia EDR3 luminosity along with our spectroscopic temperature. The resulting error radius is shown in Table 1 and appears to be particularly small (1.8%). The very precise luminosity provided by Gaia allows us to constrain stellar parameters such as the radius very well. However, these Bayesian tools underestimate the error budget because they do not take the systematic errors between one set of isochrones to the next into account because of the various underlying assumptions in the respective stellar evolutionary codes. In order to attempt to take these systematic errors into account, we combined the results of the two sets of isochrones provided by PARAM (i.e., MESA and Par-sec) and added the difference between the two sets of results to the error budget provided by PARAM. However, although using two set of isochrones may mitigate underlying systematic errors, our formal error budget for radius and luminosity maybe still be underestimated, as demonstrated by Tayar et al. (2022). For solar-type stars such as HD 20329, absolute errors may rather be up to 4%, 2%, 5%, and 20% for radius, luminosity, mass, and age, respectively.

Despite its nearly solar metallicity, the derived age from the isochrones indicates that the star is old (11 Gyr). In parallel, the HARPS-N spectrum allows us to check other indices of age, notably the chromospheric activity of the Ca H&K lines and the Li abundance. We did not observe any sign of chromospheric activity in the cores of the Ca H and K lines (hence leading to a relatively high log ). Using the activity-age relation of Mamajek & Hillenbrand (2008), we also found the age of HD 20329 to be in a range of 4–8 Gyr, which is quite consistent with the age of 9–13 Gyr as determined by isochrone fitting. In addition, we did not observe any lithium line, which supports the idea that the old age of the stars permits a slow depletion of all the lithium in its shallow convective external envelope.

3.2 Transit light curve and radial velocity model fit

To obtain the planetary mass and radius, we fit the light curves and radial velocity measurements simultaneously following the procedure described in Murgas et al. (2021). In summary, the transits were modeled with PyTransit⁵ (Parviainen 2015) adopting a quadratic limb-darkening (LD) law. The fitted LD coefficients (u₁ and u₂) were weighted against the predicted values computed by LDTK⁶ (Parviainen & Aigrain 2015) while using the Kipping (2013) coefficient parameterization. The HARPS-N radial velocity measurements were fit using RadVel⁷ (Fulton et al. 2018). For the transits, we set as free parameters the planet-to-star radius ratio R_p/R_*, the quadratic LD coefficients, the central time of the transit T_c, the planetary orbital period P, the stellar density ρ_*, the transit impact parameter b, and in the case of noncircular-orbit models, the eccentricity e and argument of the periastron ω. The free parameters for the radial velocity fit were the radial velocity semi-amplitude (K_RV), the host star systemic velocity (γ₀), and the instrumental jitter (σ_RVjitter). For the RV modeling, the orbital period, central transit time, eccentricity, and argument of the periastron were also set free, but were taken to be global parameters in common with the transit model. We also introduced a term to model the slope seen in the RV measurements.

As a result of their short orbital periods, it is expected that the orbits of USPs become circularized over time. Nevertheless, we decided to fit a noncircular orbit to the RVs, setting as global free parameters the eccentricity and the argument of the periastron using the parameterization sin(ω) and cos(ω) (parameter limits [−1, 1]). With this parameterization, we sample values of e ∈ [0, 1] and ω ∈ [0, 2π].

We modeled the residual correlated noise in TESS and HARPS-N data using Gaussian processes (GPs; e.g., Rasmussen & Williams 2010; Gibson et al. 2012; Ambikasaran et al. 2015). For TESS light curves, we chose a Matern 3/2 kernel, (1)

where ∣t_i − t_j∣ is the time between epochs in the series, and the hyperparameters c₁ and τ₁ were set free. Because of its flexibility, the Matern 3/2 kernel is a commonly used kernel for modeling TESS systematics.

The radial velocity package RadVel allows modeling correlated noise using GPs. For the radial velocity time series, we chose an exponential squared kernel (i.e., a Gaussian kernel), (2)

where t_i − t_j is the time between epochs in the series, and the hyperparameters c₂ and τ₂ were set free. We chose a Gaussian kernel to model the red noise in the RVs since it is a simple model and we saw no signs of RV signals induced by the rotation period of the star in our RV measurements periodogram, nor evident correlations between the RVs and the measured activity indices. This is probably caused by the fact that we likely covered only two rotation periods of the star (RV baseline of 62 days) if the stellar rotation is around 30 days based on activity relations.

For the fitting procedure, we constrained the prior values for the orbital period and central time of the transit based on the results of applying Transit Least Squares (TLS, Hippke & Heller 2019) to the TESS time series. Then we maximized a posterior function for the joint data set (TESS plus HARPS-N) using PyDE⁸, and started a Markov chain Monte Carlo (MCMC) using Emcee (Foreman-Mackey et al. 2013) to explore the parameter space. We used 160 chains to fit 19 free parameters, 2000 iterations as a burn-in phase, and 8000 iterations for the main MCMC phase. After the MCMC was complete, we computed the percentiles corresponding to the median and 1σ limits (from the median) of the distribution for each variable. We adopted these values as the final parameter estimates and their corresponding uncertainties.

We repeated this procedure to test whether the data supported a circular or an eccentric model and the use of GPs for modeling the RVs. We compared four models: (1) a circular-orbit model without GPs for the RV data, (2) a circular-orbit model including GPs (for TESS and the RVs), (3) a noncircular-orbit model using GPs only for TESS photometry, and (4) a noncircular-orbit fit including GPs (for TESS and HARPS-N RVs). Then we computed the model comparison metric Bayesian information criterion (BIC; Schwarz 1978) to determine the best approach to fitting our data. We found that the preferred model based on these criteria was the circular-orbit model including GPs for both TESS and RV measurements, with a BIC difference between models of ΔBIC = −7.6 when compared to the second best model, that is, the circular model without GPs for the RVs (see Table C.1 for these fit results). Hence, all the results presented in the following sections refer to the circular-orbit model including GPs for modeling the red noise in TESS and HARPS-N measurements. Figure C.1 shows the posterior distribution for the final fitted orbital parameters.

Our joint fit analysis found a linear trend in the residuals of the RV measurements after subtracting the velocity change induced by HD 20329b (see Fig. 6 mid panel). This trend may be caused by another planet in the system. Based on the TESS SAP light curves, HD 20329 appears to be a star with low photometric modulation. In addition to the low value of the log activity index, HD 20329 is potentially a good candidate for long-term radial velocity monitoring to search for other planets in this system.

We tested for any noticeable curvature in the residuals of the RV measurements. For this purpose, we added a new parameter to the RV modeling (circular orbit including GPs) to take this change into account. The results of the joint fit with this parameterization are presented in Table C.2. From this fit, we find m s−³, that is, consistent with 0, meaning that no significant curvature is found in the residuals of the fit. This was also supported by the BIC model selection criteria values, which supported an RV model with a linear trend as the best fit.

Fig. 5 HD 20329 TESS light curves and best model fit. Top and bottom left panels: TESS sector 42 and 43 photometry. The red line shows the best model fit with and without systematic effects, and the blue and black points are TESS unbinned and binned data points (with a bin size of ~1 h). Individual transit events of HD 20329b are marked with vertical black lines. Right panel: TESS phase-folded light curves and best-fitting model (black line) after subtracting the photometric variability in the two TESS sectors. The red points are binned data points.

Fig. 6 Radial velocity measurements of HD 20329 taken with HARPS-N. (a) RV time series and best-fitting model (blue line) including red noise. The best model shown here was computed using the median values of the posterior distribution for each fitted parameter. The orange shaded area around the blue line represent the 1σ uncertainty levels of the fitted model. (b) Residuals of the fit after subtracting the single-planet Keplerian function and systematic noise component from the RVs. (c) Radial velocity measurements in phase after subtracting the red noise and best-fitting model (blue line).

4 Results and discussion

The final parameter values and 1σ uncertainties are presented in Table 2. Figures 5 and 6 show the TESS light curves and the radial velocity measurements made by HARPS-N, respectively, and the best model found by the joint fit.

We find that HD 20329b has an orbital period of P = 0.926118 ± 0.00005 days, and according to our model selection, an orbital eccentricity consistent with 0, as expected for this type of planets. Using the stellar parameter values derived from a coadded HARPS-N spectrum, we determined that HD 20329b has a radius of R_p = 1.72 ± 0.07 R_⊕, a mass of M_p = 7.42 ± 1.09 M_⊕, and a bulk density of ρ_p = 8.06 ± 1.53 g cm⁻³. Assuming a Bond albedo of 0.3, we find that this planet has an equilibrium temperature of T_eq = 1958 ± 25 K.

Figure 7 shows the mass versus radius distribution for known transiting planets (parameters taken from TEPcat; Southworth 2011). USP planets around stars with T_eff ≤ 4000 K and T_eff > 4000 K are represented by orange circles and blue triangles, respectively. The figure also includes the rocky planet models of Zeng et al. (2016, 2019) with an equilibrium temperature of 2000 K. The composition models shown in Fig. 7 are planets with pure iron cores (100% Fe), Earth-like rocky compositions (32.5% Fe plus 67.5% MgSiO₃), and Earth-like composition planet cores with 0.1% H₂ gaseous envelopes. The position of HD 20329b in the mass-radius diagram agrees with other known USP planets with radii smaller than 2.0 R_⊕, and it most likely indicates a rocky composition with little to no atmosphere.

Fig. 7 Mass-radius diagram for transiting planets with mass determinations with a precision better than 30% (parameters taken from the TEPcat database; Southworth 2011). The position of HD 20329b is shown by the red star, USP planets orbiting M-type stars (Teff ≤ 4000 K) are marked with the orange circles, USP planets around stars with Teff > 4000 K are represented by blue triangles, and non-USP planets around other types of stars are marked by gray squares. The lines in the massradius diagram represent the composition models of Zeng et al. (2016, 2019): planets with pure iron cores (100% Fe, brown line), Earth-like rocky compositions (32.5% Fe plus 67.5% MgSiO3, dashed green line), and a 99.9% Earthlike rocky core (32.5% Fe plus 67.5% MgSiO3) with a 0.1% H2 envelope with a temperature of 2000 K (dash-dotted pink line).

Fig. 8 TSM vs. ESM (Kempton et al. 2018) values for known transiting planets. The USP planets (P < 1 day) with radii smaller than 2.0 R⊙ are marked by circles. The color and size of the circles represent the equilibrium temperature and radius of the planet, respectively. USPs with radii larger than 2.0 R⊙ are represented with purple pentagons, and small planets (Rp < 2.0 R⊙) with orbital periods longer than one day are marked with gray squares. The position of HD 20329b is marked by the red star and arrow. We did not label all USP planets in the plot for easy viewing.

4.1 Prospects for atmospheric characterization

We used data from NASA Exoplanet Archive⁹ to compute the transmission spectroscopy metric (TSM) and emission spectroscopy metric (ESM) as defined by Kempton et al. (2018) for known transiting planets. These two metrics are used as indicators for the feasibility of detecting atmospheric signals during a transit (TSM) or detecting the secondary eclipse of a transiting planet (ESM) with the James Webb Space Telescope (JWST; Gardner et al. 2006). We find a TSM value of 45.7 for HD 20329b. However, Kempton et al. (2018) recommend a threshold for the TSM value higher than 90 for planets in the radius range of 1.5 < R_p < 10 R_⊕; this would mean that HD 20329b is not an ideal candidate for transmission spectroscopy follow-up if the planet holds an atmosphere. On the other hand, we find an ESM = 10.2. This is a favorable indicator that the secondary eclipse of this planet might be detected with JWST given that Kempton et al. (2018) established an ESM threshold for rocky worlds of 7.5. We explore the possibility of detecting the secondary transit and phase variations of HD 20329b using TESS data in Sect. 4.2.

Figure 8 shows the TSM and ESM values for known transiting planets. We focused on USP period planets and planets with radii smaller than 2 R_⊕. The planets with the most similar TSM and ESM values to HD 20329b are GJ 367b and K2-141b. GJ 367b was discovered by Lam et al. (2021). The planet is a sub-Earth that orbits an M dwarf with an orbital period of P = 0.321 days (7.7 h) and has a mass and radius of 0.546 ± 0.078 M_⊕ and 0.718 ± 0.054 R_⊕. The mean bulk density of this planet is similar to that of HD 20329b with ρ_p = 8.1 ± 2.2 g cm⁻³. Because GJ 367b orbits an M dwarf, this planet receives less flux than HD 20329b, and their equilibrium temperatures differ by ~360 K (T_eq = 1597 ± 39 K for GJ 367b assuming A_B = 0.3). K2-141b was discovered independently by Malavolta et al. (2018) and Barragán et al. (2018) using data from the Kepler extended mission K2 (Howell et al. 2014). This planet is a rocky super-Earth orbiting a K dwarf with a period of P = 0.280 days (~6.7 h), and with a mass and radius of 5.08 ± 0.041 M_⊕ and 1.51 ± 0.05 R_⊕. K2-141b has a similar mean bulk density (ρ_v = 8.1 ± 1.1 g cm⁻³) and equilibrium temperature ( K, Barragán et al. 2018) as HD 20329b.

Fig. 9 TESS observations from sectors 42 and 43 folded over the orbital phase (in days) and binned (for visualization) with the median posterior phase-curve model (black line) and its 16th and 84th percentile limits (blue shading). The eclipse is modeled assuming a uniform stellar disk with a depth scaled from 0 (eclipse) to 1 (out of eclipse).

4.2 Secondary eclipse in TESS data

In this section, we describe our efforts to explore the detection of the secondary transit and phase variations of HD 20329b using TESS data. We modeled the TESS light curves from sectors 42 and 43 using the PhaseCurveLPF phase-curve model implemented in PyTransit (see Parviainen et al. 2021). The phase-curve model includes the effects from thermal emission, reflection, ellipsoidal variation, and Doppler boosting, but we forced the ellipsoidal variation and Doppler boosting amplitudes to zero with narrow normal priors because the planet-to-star mass ratio is too low for these two effects to have any practical importance. Furthermore, we set the geometric albedo to zero (again, using a narrow normal prior) because the reflection and emission components are strongly degenerate, and we are primarily interested in determining whether an eclipse signal exists in the light curves.

With these exceptions, the phase curve model includes 12 free parameters (transit center, orbital period, stellar density, impact parameter, cos ω, sin ω, planet-star area ratio, logio dayside planet-star flux ratio from emission, logio night-side planet-star flux ratio from emission, emission peak offset, and two limb-darkening coefficients). We set weakly informative normal priors on the transit center and orbital period based on the information in ExoFOP, and an informative normal prior, N(0.88, 0.07), on the stellar density based on the stellar characterisation described in Sect. 3.1. We also forced a circular orbit by setting tight zero-centered normal priors on cos ω and sin ω and constrained the limb-darkening coefficients using priors calculated with LDTK. The baseline flux variations were modeled as a time-dependent GP using Celerite (Foreman-Mackey et al. 2017). Each TESS sector was assigned a separate GP with its own hyperparameters, and the hyperparameters were kept free and marginalized over during the posterior sampling. The day- and nightside logio flux ratios had uniform priors from −3 to 0.

The phase-curve model parameter posteriors were obtained as usual. We first determined the global posterior maximum using a global optimization method, and then sampled the posterior using MCMC sampling.

The results from the phase-curve modeling (Fig. 9) support the existence of an eclipse signal in the TESS data. The logio dayside flux ratio posterior has its median at −0.98, and the day-side flux ratio estimate¹⁰ is 11 ± 8%, which corresponds to an eclipse depth of 20 ± 14 ppm, as shown in Fig. 10. The nightside flux ratio is not constrained. The dayside posterior does not exclude a no-eclipse scenario (i.e., with a planet-star flux ratio of zero), but does support a relatively significant planetary emission signal. We used logio flux ratio as a sampling parameter, and a uniform prior on it sets a reciprocal (log-uniform) prior on the flux ratio itself. We would expect a uniform logio flux ratio posterior distribution in the absence of an eclipse signal, which would allow us to give only an upper limit for the flux ratio. However, the logio dayside flux ratio posterior has a strong mode with a tail toward lower values, which indicates that models with little to no eclipse signal can also explain the observations.

Curious about the possible eclipse detection, we mapped the brightness temperatures and geometric albedos that might explain the dayside flux ratio posterior. We modeled the flux ratio as the sum of emission and reflection components. The emission component was calculated using the BT-Settl spectra¹¹, where we assumed the host star to have an effective temperature of 5640 K and kept the planet temperature as a free parameter. The reflection component had the geometric albedo as a free parameter, and the remaining parameters (semi-major axis and planet-star radius ratio) were fixed to the posterior median values from the phase-curve modeling.

We estimated the brightness temperature and geometric albedo posterior via MCMC sampling using a Gaussian kernel density estimate (kde) of the logio dayside flux ratio posterior as a target distribution (i.e., the log likelihood is based on the kde log pdf evaluated for the logio flux ratio model values). The final geometric albedo versus brightness distribution is shown in Fig. 11. The dayside flux ratio posterior from the phase-curve modeling suggests a very high brightness temperature for the planet that can be reduced only by high geometric albedos. The median posterior brightness temperature for A_g < 0.25 is ≈3500 K, while for 0.25 < A_g < 0.5, it is ≈3300 K. However, the median posterior value does not characterize the marginalized brightness temperature posterior well because it changes from a distribution with a clear mode closer to a uniform distribution as the geometric albedo increases. The upper limit for the brightness temperature decreases with increasing geometric albedo (as it should), but the 99th posterior percentile is ~4000 K over the A_g range from 0 to 1. To support this conclusion, we also modeled the light curve with the transit and light-curve modeler code (Csizmadia 2020), which obtained similar results (see Appendix D.1).

Fig. 10 Posterior probability densities for the dayside logio planet-to-star flux ratio used as a sampling parameter and the flux ratio.

4.3 Additional planets in the system

As pointed out in Sanchis-Ojeda et al. (2014) and Winn et al. (2018), USP planets are often found in multiple-planet systems. We searched for additional transit events in TESS data using Transit Least Squares (TLS; Hippke & Heller 2019). We used the results of the joint model fit to mask the transit events of planet b and to remove the photometric noise in the TESS time series modeled by the GPs. We also tested median filters with different window sizes to remove the photometric variability in the light curves and ran TLS multiple times. We found no evidence of other significant transit events in TESS data. The transits that TLS detected are likely spurious signals, with at least one transit occurring at the beginning or end of a block of TESS observations where systematic effects due to the pointing of the satellite are stronger. In addition, the detected events presented a low signal detection efficiency (SDE) in comparison to HD 20329b, which TLS detected.

The time baseline of our RV measurements is 62 days because we do not detect a significant curvature in the RV residuals (see end of Sect. 3.2). This means that the object that causes the RV slope must have an orbital period longer than at least four times the baseline. Using Kepler’s third law and assuming M_* ≫ M_p, an object with P = 120 days has an orbital semi-major axis of 0.46 AU. Based on the maximum ΔRV difference of the residuals and setting this value as the semi-amplitude of an object with orbital period P = 120 days, we find a minimum mass of M_p sin(i) = 0.17 M_Jup.

If the radial velocity change seen in the RV residuals is caused by a stellar mass object, this might affect the projected motion of the star in the sky. Brandt (2021) used the HlPPARCOS (ESA 1997) and Gaia EDR3 (Lindegren et al. 2021) data to create a catalog of proper motions with a time baseline of ~24 yr (1991.25, 2015.5) to identify astrometrically accelerating stars. For HD 20329, the Gaia proper motions are µ₀ = (111.77 ± 0.030, −202.41 ± 0.025) mas yr⁻¹ and the derived HlPPARCOS-Gaia proper motions are µ_HG = (111.73 ± 0.035, −202.43 ± 0.026) mas yr⁻¹, meaning that the star presented a change in proper motion of Δµ = ∣µ_G − µ_G∣ = 0.047 ± 0.044 mas yr⁻¹ over 24 yr. This null acceleration is consistent with the Gaia EDR3 assessment discussed in Sect. 2.3, in which Gaia measurements favor a single-star model. The slope in the RV residuals is therefore likely caused by a substellar object.

Fig. 11 Posterior probability densities for the geometric albedo and brightness temperature.

5 Conclusions

We reported the discovery of HD 20329b, an ultra-short-period transiting planet around a solar-type star. Observations made by the NASA TESS mission led to the initial detection of the transits. Follow-up radial velocity measurements taken with the HARPS-N spectrograph made it possible to confirm the planetary nature of the transiting object and to establish its mass and radius.

From a coadded HARPS-N mean stellar spectrum, we estimate that the host star has an effective temperature of T_eff = 5596 ± 50 K, a metallicity of [Fe/H] = −0.07 ± 0.06 dex, and a derived stellar mass and radius of = 0.90 ± 0.05 M_⊙ and R_* = 1.13 ± 0.02 R_⊙.

We analyzed TESS photometric data from sectors 42 and 43 and radial velocity measurements taken with the HARPS-N spectrograph. The transit observations made by TESS and the radial velocity measurements were fit simultaneously using an MCMC procedure that included Gaussian processes to model the systematic effects present in the light curves and RV measurements. We find that HD 20329b has a radius of R_p = 1.72 ± 0.07 R_⊕ and a mass of M_p = 7.42 ± 1.09 M_⊕, and that it orbits its star with a period of 0.926118 ± 0.50 × 10⁻⁴ days. We derived a mean bulk density of ρ_p = 8.06 ± 1.53 g cm⁻³, indicating a likely rocky composition. We note that after subtracting the planetary signal of the USP planet from the RV measurements, the RV residuals present a slope that might indicate the presence of an additional low-mass object orbiting HD 20329.

The ESM for HD 20329b indicates that this planet presents favorable conditions for secondary transit follow-up with the JWST. We used a simple phase-curve model including reflected light and thermal emission to search for the secondary transit and phase variations in TESS light curves. Our results support the existence of a significant, although not conclusive, eclipse signal in the TESS data, with a dayside flux ratio of 11% and a relatively strong planetary emission signal. Our modeling indicates a brightness temperature of ~3500 K for low geometric albedo values (A_g < 0.25) and an upper limit on the brightness temperature of ~4000 K over the range A_g ∈ [0, 1]. Precise observations, preferable in the IR, are needed to confirm these results.

Overall, HD 20329b is a new addition to the ~ 120 currently known USP planets. It presents very favorable metrics for a secondary transit detection with the JWST and for radial velocity follow-up to search for additional planets in the system.

Acknowledgements

T.M. acknowledges financial support from the Spanish Ministry of Science and Innovation (MICINN) through the Spanish State Research Agency, under the Severo Ochoa Program 2020-2023 (CEX2019-000920-S) as well as support from the ACIISI, Consejería de Economía, Conocimiento y Empleo del Gobierno de Canarias and the European Regional Development Fund (ERDF) under grant with reference PROID2021010128. R.L. acknowledges funding from University of La Laguna through the Margarita Salas Fellowship from the Spanish Ministry of Universities ref. UNI/551/2021-May 26 under the EU Next Generation funds and financial support from the Spanish Ministerio de Ciencia e Innovación, through project PID2019-109522GB-C52, and the Centre of Excellence “Severo Ochoa” award to the Instituto de Astrofísica de Andalucía (SEV-2017-0709). PK acknolwedges the support from grant LTT-20015. K.W.F.L. was supported by Deutsche Forschungsgemeinschaft grants RA714/14-1 within the DFG Schwerpunkt SPP 1992, Exploring the Diversity of Extrasolar Planets. HJD acknowledges support from the Spanish Research Agency of the Ministry of Science and Innovation (AEI-MICINN) under the grant PID2019-107061GB-C66, DOI: 10.13039/501100011033. C.M.P. and J.K. gratefully acknowledge the support of the Swedish National Space Agency (DNR 65/19 2020-00104). This paper includes data collected by the TESS mission, which are publicly available from the Mikulski Archive for Space Telescopes (MAST). Funding for the TESS mission is provided by NASA’s Science Mission directorate. 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. We acknowledge the use of public TESS data from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center. This research has made use of the Exoplanet Follow-up Observation Program website, which is operated by the California Institute of Technology, under contract with the National Aeronautics and Space Administration under the Exoplanet Exploration Program. This publication makes use of data products from the AAVSO Photometric All Sky Survey (APASS). Funded by the Robert Martin Ayers Sciences Fund and the National Science Foundation. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. This work made use of tpfplotter by J. Lillo-Box (publicly available in www.github.com/jlillo/tpfplotter), which also made use of the python packages astropy, lightkurve, matplotlib and numpy. This work is partly supported by JSPS KAKENHI Grant Numbers JP17H04574, JP18H05442, JST CREST Grant Number JPMJCR1761, and the Astrobiology Center of National Institutes of Natural Sciences (NINS) (Grant Number AB031010). This paper is based on observations made with the MuSCAT2 instrument, developed by ABC, at Telescopio Carlos Sánchez operated on the island of Tenerife by the IAC in the Spanish Observatorio del Teide.

Appendix A Radial velocity data

Radial velocity and activity indices measurements are available electronically at CDS.

Appendix B Stellar rotation from ASAS-SN photometry

Inhomogeneities in the stellar surface (e.g., spots or plages) that appear and disappear out of sight as the star rotates can imprint flux variations that can be detected by photometric measurements. We obtained V-band photometry for HD 20329b from ASAS-SN public light-curve archive¹² (Shappee et al. 2014; Kochanek et al. 2017). The ASAS-SN archive has 251 observations spanning a time baseline of 1832 days (t_start = 2456618.931 HJD, t_end = 2458451.844 HJD). The ASAS-SN photometry is shown in Figure B.1. The GLS periodogram of the photometry presents several peaks with about equal power. The strongest peak corresponds to a period of 423 days. We fit the ASAS-SN photometry using the package for Gaussian Processes Celerite (Foreman-Mackey et al. 2017) and fit the photometry with the kernel (B.1)

where ∣t_i − t_j∣ is the difference between two epochs of observations, B, C, L are positive constants, and P_rot is the stellar rotational period (see Foreman-Mackey et al. 2017 for details). For the fit, we considered each ASAS-SN camera (bd and bh) as an independent instrument. Each data set had as free parameters the constants B, C, L, but shared a common rotational period parameter.

The fitting procedure consisted of a global optimization of a log likelihood function. Then we sampled the posterior distribution of the kernel we used to fit the photometry with Emcee (Foreman-Mackey et al. 2013) using 120 chains and 10000 iterations for the main MCMC. The final parameter values (median and 1σ uncertainties) were estimated from the posterior distribution.

Figure B.1 presents the photometry (top panel) and the posterior distribution of the fitted rotation period parameter (bottom right panel). The distribution of the rotation period parameter presents a long a tail of possible values, but the median value is days.

Fig. B.1 Ground-based long-term photometry of HD 20329. Top panel: ASAS-SN V-band measurements from cameras bd and bh. Bottom left panel: GLS periodogram (Zechmeister & Kürster 2009) of the V-band photometry. Bottom right panel: Posterior distribution of the rotation period parameter from the fit.

Appendix C Light curve and radial velocity joint fit

In this section, we present the tests and results from the joint fit of TESS photometric data and the HARPS-N radial velocity measurements. Figure C.1 shows the parameter distributions for the best-fit model, table C.1 presents the results of the fit models we tested in our model selection analysis, and table C.2 presents the results of a fit including a second-order term to model the linear trend seen in the RV residuals of the best-fit model.

Fig. C.1 Correlation plot for the fit orbital parameters. Limb-darkening coefficients and systematic effect parameters were intentionally left out for easy viewing. The blue lines give the median values of the parameters.

Appendix D Light-curve fit with TLCM

We also fit the light curve with the TRANSIT AND LIGHT CURVE MODELLER code described in detail in Csizmadia (2020). We briefly summarize its main features here. This code is able to perform a joint radial velocity and light-curve fit. The transit and occultation events can be modeled with different limb-darkening laws. The out-of-transit variations can be modeled by including the ellipsoidal, reflection, and beaming effects (for details, see Csizmadia et al. 2021). Contamination effects and eccentric orbit are included. The instrumental noise and stellar variability in the light curve are modeled by wavelets that are simultaneously fit with the light-curve model (for a more detailed description, see Csizmadia et al. 2021).

In the fit, the free parameters were the scaled semi-major axis (a/R_star), the planet-to-star radius ratio (R_planet/R_star), the impact parameter (b), the planet-to-star surface brightness ratio, the epoch, the period, the RV amplitude for the beaming effect, the planetary geometric albedo, the mass ratio for the ellipsoidal effect, the wavelet noise parameters σ_w and σ_r, and the two limb-darkening coefficient combinations for the quadratic limb-darkening law u₊ = u_a + u_b and u₋ = u_a − u_b. Since the derived eccentricity seems to be compatible with zero (see Table 2), we assumed a circular orbit for the planet.

The solution was optimized with the genetic algorithm (Geem et al. 2001) first, then a differential evolution MCMC analysis (Nelson et al. 2014; Sherri et al. 2017) was performed to estimate the median of the posterior distribution. The uncertainties were obtained using the usual 16-84% rule. The results are reported in Table D.1.

The combined model + noise model fit is presented graphically in the left panel of Figure D.1. For visualization purposes, we subtracted the wavelet-based red-noise component from the light curve, then we phase-folded and plotted it in the right panel of Figure D.1. We binned the phase-folded, red-noise-corrected light curve into 50 bins (bin size ~27 minutes) to show the occultation event (secondary transit) at phase 0.5. Although there is a flux loss at this phase, we did not detect the occultation of the planet at a 3σ significance level. This detection is just tentative. We argue for this in the following way: First, a similar flux loss is visible in this binned light curve at phase 0.25, where we do not expect such a flux drop (Figure D.2). Second, the observed geometric albedo is , meaning that it cannot be derived from this light-curve set. We also measured the surface brightness ratio of the planet and the star to be , which is compatible with zero. When blackbody radiation is asumed for the star and the planet, the median value would mean a temperature of 3300 K for the planet. This agrees well with the analysis presented in Section 4.2. The occultation depth is 30 ± 54 ppm, which does not exclude a nondetection of the secondary transit.

We conclude that we can give only a 1σ upper limit of 84 ppm for the occultation depth in the system of HD 20329b with TLCM. The error bars given by TLCM are slightly larger than the values in Table 2. The derived planet-to-star radius ratio agrees within 1σ of the error bar with the values reported in Table 2, and all other parameters also agree reasonably well with the finally accepted values.

Fig. D.1 HD 20329 TESS photometry and TLCM fit. Left panel: TESS normalized SAP-fluxes vs time. Black dots are the observations. The red curve is the result of the simultaneous transit + occultation + beaming + reflection + ellipsoidal + wavelet-based noise model fit. Top right panel: Primary transit of HD 20329b (black dots) and best model fit (red line) after subtracting the red-noise component from the light curve (top). Bottom right panel: Residuals of the fit.

Fig. D.2 TESS binned light curve of HD 20329b (50 bins with a ~27 min. bin size), zoomed at the out-of-transit part (phase curve). See details in Appendix D. TLCM normalized the light curve to phase 0.25.

References

All Tables

All Figures

	Fig. 1 TESS target pixel file image of HD 20329 observed in sector 42 (left) and sector 43 (right), made with tpf plotter (Aller et al. 2020). The pixels highlighted in red show the aperture used by TESS to create the light curves. The position and sizes of the red circles represent the position and TESS magnitudes of nearby stars, respectively.
In the text

	Fig. 2 Fig. 2. Speckle sensitivity curve and ACF (inset) from SOAR HRCam observations of HD 20329. No nearby stars were detected within 3″of HD 20329 with HRCam.
In the text

	Fig. 3 Palomar NIR AO imaging and sensitivity curves for HD 20329 taken in the Bry filter. The images were taken in good seeing conditions, and we reach a contrast of 7 mag. fainter than the host star within 0.5″. Inset: image of the central portion of the data, centered on the star.
In the text

Fig. 4 GLS periodogram (Zechmeister & Kürster 2009) for HARPS-N radial velocity measurements, (α) Radial velocity measurement residuals after fitting the transiting-planet signal, (b) Stellar activity indices of HD 20329 ((c) to (k)). The peak with the maximum power in the periodogram for the RV measurements corresponds to a period of P = 0.9263 days (vertical green line), in agreement with the period of the transiting-planet candidate reported by TESS. The dash-dotted blue and segmented red line represent the 0.01% and 0.1% FAP levels, respectively.

In the text

Fig. 5 HD 20329 TESS light curves and best model fit. Top and bottom left panels: TESS sector 42 and 43 photometry. The red line shows the best model fit with and without systematic effects, and the blue and black points are TESS unbinned and binned data points (with a bin size of ~1 h). Individual transit events of HD 20329b are marked with vertical black lines. Right panel: TESS phase-folded light curves and best-fitting model (black line) after subtracting the photometric variability in the two TESS sectors. The red points are binned data points.

In the text

Fig. 6 Radial velocity measurements of HD 20329 taken with HARPS-N. (a) RV time series and best-fitting model (blue line) including red noise. The best model shown here was computed using the median values of the posterior distribution for each fitted parameter. The orange shaded area around the blue line represent the 1σ uncertainty levels of the fitted model. (b) Residuals of the fit after subtracting the single-planet Keplerian function and systematic noise component from the RVs. (c) Radial velocity measurements in phase after subtracting the red noise and best-fitting model (blue line).

In the text

Fig. 7 Mass-radius diagram for transiting planets with mass determinations with a precision better than 30% (parameters taken from the TEPcat database; Southworth 2011). The position of HD 20329b is shown by the red star, USP planets orbiting M-type stars (Teff ≤ 4000 K) are marked with the orange circles, USP planets around stars with Teff > 4000 K are represented by blue triangles, and non-USP planets around other types of stars are marked by gray squares. The lines in the massradius diagram represent the composition models of Zeng et al. (2016, 2019): planets with pure iron cores (100% Fe, brown line), Earth-like rocky compositions (32.5% Fe plus 67.5% MgSiO3, dashed green line), and a 99.9% Earthlike rocky core (32.5% Fe plus 67.5% MgSiO3) with a 0.1% H2 envelope with a temperature of 2000 K (dash-dotted pink line).

In the text

Fig. 8 TSM vs. ESM (Kempton et al. 2018) values for known transiting planets. The USP planets (P < 1 day) with radii smaller than 2.0 R⊙ are marked by circles. The color and size of the circles represent the equilibrium temperature and radius of the planet, respectively. USPs with radii larger than 2.0 R⊙ are represented with purple pentagons, and small planets (Rp < 2.0 R⊙) with orbital periods longer than one day are marked with gray squares. The position of HD 20329b is marked by the red star and arrow. We did not label all USP planets in the plot for easy viewing.

In the text

	Fig. 9 TESS observations from sectors 42 and 43 folded over the orbital phase (in days) and binned (for visualization) with the median posterior phase-curve model (black line) and its 16th and 84th percentile limits (blue shading). The eclipse is modeled assuming a uniform stellar disk with a depth scaled from 0 (eclipse) to 1 (out of eclipse).
In the text

	Fig. 10 Posterior probability densities for the dayside logio planet-to-star flux ratio used as a sampling parameter and the flux ratio.
In the text

	Fig. 11 Posterior probability densities for the geometric albedo and brightness temperature.
In the text

	Fig. B.1 Ground-based long-term photometry of HD 20329. Top panel: ASAS-SN V-band measurements from cameras bd and bh. Bottom left panel: GLS periodogram (Zechmeister & Kürster 2009) of the V-band photometry. Bottom right panel: Posterior distribution of the rotation period parameter from the fit.
In the text

	Fig. C.1 Correlation plot for the fit orbital parameters. Limb-darkening coefficients and systematic effect parameters were intentionally left out for easy viewing. The blue lines give the median values of the parameters.
In the text

Fig. D.1 HD 20329 TESS photometry and TLCM fit. Left panel: TESS normalized SAP-fluxes vs time. Black dots are the observations. The red curve is the result of the simultaneous transit + occultation + beaming + reflection + ellipsoidal + wavelet-based noise model fit. Top right panel: Primary transit of HD 20329b (black dots) and best model fit (red line) after subtracting the red-noise component from the light curve (top). Bottom right panel: Residuals of the fit.

In the text

	Fig. D.2 TESS binned light curve of HD 20329b (50 bins with a ~27 min. bin size), zoomed at the out-of-transit part (phase curve). See details in Appendix D. TLCM normalized the light curve to phase 0.25.
In the text