© The Authors 2025

1 Introduction

Trans-Neptunian Objects (TNOs) are remnants of the original protoplanetary disk that orbit beyond Neptune. Although some of the largest TNOs (such as Pluto or Eris) might have experienced significant evolutionary changes, many smaller objects are believed to have remained largely unaltered, offering valuable information about the composition and conditions of the early Solar System (Nesvorný 2018).

Studying TNO shapes gives us insight into their bulk composition. It is generally believed that there is a material-dependent critical size above which TNOs behave as fluids and adopt hydrostatic equilibrium shapes (Tancredi & Favre 2008). This has not been thoroughly tested for TNOs. Moreover, a precise determination of the shape can reveal the internal structure and determine ice-rock ratios (Rambaux et al. 2017).

Surface properties can also be used as a proxy for TNO composition (Barucci & Merlin 2020). Two relevant material-dependent parameters are the geometric albedo p_V and the beaming parameter η. The former measures the surface reflecting properties (Karttunen et al. 2017); the latter is known to be dependent on surface roughness and thermal inertia effects (Lellouch et al. 2013). However, these properties are hard to constrain if the size of the object is unknown.

A stellar occultation occurs when a TNO passes in front of a distant star, thus blocking its light temporarily. With observers in different locations, a stellar occultation can put strong constraints on the size and shape of the occulting TNO. The accuracy of this technique can be of the order of kilometres (Ortiz et al. 2020b), only paralleled by in situ space missions.

Occultations enable the computation of the geometric albedo p_V when combined with independent measurements of the absolute magnitude H_V, as in Santos-Sanz et al. (2021); Kretlow et al. (2024a). Occultations also allow us to determine the apparent (projected) axis ratio. This technique has revealed that some TNOs exhibit shapes that deviate from the three-dimensional hydrostatic equilibrium expected for objects with a uniform composition (Ortiz et al. 2017; Vara-Lubiano et al. 2022). The projected axis ratio has only been tabulated for a handful of TNOs (see detailed information in Sect. 4.2). Acquiring a statistically significant sample of projected axis ratios is a first step if we want to infer what the three-dimensional shape distribution is in the trans-Neptunian Belt.

Stellar occultations have also revealed the presence of rings (Braga-Ribas et al. 2014; Ortiz et al. 2017, 2023), topographic features (Dias-Oliveira et al. 2017; Rommel et al. 2023), and binarity (Buie et al. 2020; Leiva et al. 2020), and have provided accurate astrometric positions for the object down to a few milliarcseconds (Rommel et al. 2020). In the case of double stars, occultations can constrain the relative geometry and relative brightness of the stars (Leiva et al. 2020).

The object (470316) 2007 OC10 is a TNO discovered at Palomar Observatory on 2007 July 22. It presents an orbital eccentricity of e = 0.29 and an inclination of i = 21.7°. The perihelion lies at 35.53 AU (data from the IAU Minor Planet Center). According to the classification criteria presented in Gladman et al. (2008) and the computations in Volk & Laerhoven (2024), it belongs to the dynamical class of detached TNOs: non-resonant and non-scattering objects with orbital eccentricities e > 0.24. We highlight that (229762) G!kún||'hòmdímà (provisional designation 2007 UK126) is the only other detached TNO with published shape analysis through stellar occultations, with an area-equivalent diameter of D_eq = ${638}_{-14}^{+28}$ km (Benedetti-Rossi et al. 2016).

The TNO 2007 OC10 was a target of the TNOs are Cool programme (Müller et al. 2010), conceived to derive physical and thermal properties of TNOs with Herschel and Spitzer data. The fit of radiometric models to 2007 OC10 data provided an area-equivalent diameter (sometimes referred to as effective diameter) of D_eq = 309 ± 37 km and a geometric albedo of p_V = ${12.7}_{-2.8}^{+0.4}\%$, adopting an absolute magnitude of H_V = 5.43 ± 0.10 (Santos-Sanz et al. 2012). Comparing occultation data with thermal measurements can provide further insights into the physical properties of the object under study. Discrepancies between thermal and occultation projected size have been proposed as a way to discover unresolved satellites (Ortiz et al. 2020a; Rommel et al. 2023). In addition, occultation results can be used to improve our knowledge of the beaming parameter (Harris 1998).

In this paper, we present the results from the occultation produced by 2007 OC10 on 2022 August 22. The occultation was observed from several sites on the Mediterranean Basin (Sect. 2). We acquired complementary photometry to constrain the absolute magnitude (Sect. 2.4). We adopt a parametric model (Sect. 3.1) and a Markov chain Monte Carlo (MCMC) technique to derive a probability distribution of position, size, and shape parameters (Sect. 3.2). Finally, we compare our size measurement with thermal data to refine thermal parameters (Sect. 4.3). Given the interest and importance of finding close satellites around TNOs that cannot be resolved with ALMA or space-based telescopes, we constrain any putative undetected satellites (Sect. 5).

2 Observations

2.1 Prediction of the stellar occultation

The occultation event was predicted by the Collaborative Occultation Resources and Archive (CORA) service (Kretlow et al. 2024b). According to Jet Propulsion Laboratory¹ (JPL) ephemeris (solution # 12), 2007 OC10 would produce an occultation of the star identified as 2727866328215869952 in Gaia DR3 on 2022 August 22 at 01:32:37 UTC ± 0.56 min. The upper half of Table 1 summarises the occulted star astrometric information from Gaia DR3 catalogue. The lower half shows the derived astrometric position and uncertainties, propagated at the occultation epoch, used in the occultation analysis.

The original predicted shadow path, adopting a single star model, crossed continental Europe (see Fig. 1). It would be accessible by small-class telescopes due to the considerable brightness of the star, V = 15.23. In order to further refine this already promising occultation prediction, astrometry of 2007 OC10 was acquired in the days prior to the event. The aim was to update the ephemeris of the TNO, reducing uncertainties in astrometric position. Eighteen images of 2007 OC10 were taken from the Liverpool 2 m telescope in Roque de los Muchachos Observatory (La Palma, Spain) between 2022 August 5 and 2022 August 6. The 4k × 4k IO:O instrument was used, with SDSS-r filter. Sixty more images were taken from the 1.23 m telescope in Calar Alto Obervatory (Almería, Spain, abbreviated as CAHA) from 2022 August 4 to 2022 August 8. We used the 4k × 4k DLR-MKIII CCD instrument, with no filter. The exposure times ranged from 300 s to 400 s in both cases.

The average seeing on the observations from Calar Alto was 1.2 arcsec, whereas in Roque de los Muchachos it was 1.0 arcsec. The signal-to-noise ratio was, on average, of S/N=16.4 for the images of 2007 OC10 taken from Calar Alto and S/N=16.6 for the case of Liverpool images.

The obtained offsets in the cross-track direction (perpendicular to the motion of 2007 OC10) from JPL #12 ephemeris were Δ_Liverpool = −54 ± 8 mas and Δ_CAHA = −50 ± 7 mas respectively. The negative signs indicate a shift southward in the celestial sphere. The updated predictions are those shown in Fig. 1.

The 1σ region for the original prediction in Fig. 1 is derived from the (RA, Dec) uncertainties available in JPL. For the updated predictions, they represent the formal uncertainties derived from the standard deviation of the observed astrometric points. Both updated predictions agree within these 1σ values. It should be noted that there is a ~60 mas difference in the observed offset in the direction of motion. This is most likely due to low accuracy in the time-stamping of the images, which is not uncommon for professional telescopes. Hence, we did not update the expected closest approach time and use the CORA reported value (extracted from JPL #12 ephemeris).

Fig. 1 Comparison of occultation predictions. At this time, the double nature of the star was unknown. The continuous lines represent the predicted path of the centre of the shadow. The 1σ regions along the path are represented by the broken lines. The red diamonds represent the stations that collected data. The grey squares represent the sites that suffered from bad weather or whose data were not used in the analysis. Some stations have been slightly shifted in position for improved visualisation, but the exact coordinates used in the analysis are those in Table C.1.

2.2 Observation of the stellar occultation

On the 2022 August 22, 18 stations attempted an observation of the event as reported in the Tubitak Occultation Portal website (Kilic et al. 2022). Their locations are represented in Fig. 1. On Table C.1 there is detailed information about each observing site. Cancelada, Forcarei and Sabadell suffered from bad weather and no data could be collected. The rest of the stations obtained FITS images series, with the exception of Guirguillano Observatory, which recorded a video, and Las Cumbres Observatory (LCO, Tenerife), which used the drift scan technique. This technique consists in exposing for a long time without tracking, so that the stars in the field leave a trail on the image.

From the 15 data sets obtained, 13 were used in the occultation analysis. Data from Observatorio de Aras de los Olmos (OAO) was excluded because it consisted on 300 s exposure time images, whereas the estimated maximum occultation duration was 12 s assuming a 300 km diameter. The drift scan from LCO was also discarded because of its low S/N.

Botorrita experienced synchronisation problems between the GPS and the acquisition software. The system time was stamped on the images and showed a drift of about 3 seconds with respect to the GPS. Çanakkale data also shows problems with the time stamps because these are truncated to the nearest second. However, after data reduction, it is clear that both sites did not detect the occultation (see Fig. 2). Therefore, shifts in the timings of these data sets do not alter the result of the analysis.

For the rest of the stations, the synchronisation of the computer time to GPS or internet NTP servers worked correctly. Their timings can be trusted to 0.1 s at least (see Ortiz et al. 2017, 2020a).

The data from Guirguillano station was captured at 10.24 s integrations and stored as a 25 fps video file. This is a usual technique achieved by repeating the integration into consecutive video frames. Each video frame contains timestamps from a GPS receiver, allowing the start of the integration and the exposure time during the analysis to be reconstructed. The standard processing consists of extracting the individual video frames, grouping and averaging all the frames corresponding to the same integration. The result is then saved into FITS files, storing the start of the exposure and the restored exposure time as metadata. Additional characteristics of the video data and procedure details are described in Benedetti-Rossi et al. (2016) and Buie & Keller (2016).

Fig. 2 Light curves of the stellar occultation (normalised flux vs time at mid-exposure). The origin of time is set to be the time of closest approach: 2022-08-22T01:32:37 UTC. They are ordered by distance perpendicular to the predicted shadow path (as predicted by CAHA data). The error bars are too small to be seen in the plot in most light curves, but they can be checked in the supplementary data on CDS.

2.3 Data reduction: Occultation light curves

To obtain light curves for the final 13 data sets, we used our own Interactive Data Language (IDL) based routines for photometry. These routines are run on the FITS files, which have previously been dark- and flat-field calibrated, whenever calibration images were available.

These routines perform aperture photometry on the target star and on several comparison stars. Firstly, a file containing the time information is generated, with the mid-exposure time for each frame. Then, for every frame the integrated flux of the target is divided by the sum of the integrated fluxes of all comparison stars, providing a differential light curve. These fluxes are integrated using a circular aperture and are background corrected (background values computed integrating in an annulus). Finally, the resulting light curve is normalised so that the baseline level (outside any potential flux drop) is equal to one. The selection of the aperture radius was made manually for each case to maximise the light curve S/N. The same principle was applied to the selection of comparison stars.

The obtained light curves are plotted in Fig. 2. They can be accessed on the CDS archive. There are clear occultation features in the curves from GAL Hassin, Serra la Nave (SLN) and Artemis, without secondary drops. In a classic approach, they would be considered positive detections. However, TAR2 light curve is much lower in S/N (Artemis aperture is 100 cm, TAR2 is 42 cm) and the flux drop is not apparent. But this telescope is only 300 m away from Artemis (both at the Teide Observatory). Therefore, it should also be a positive detection.

A stellar flux drop of 75% was observed at GAL Hassin and SLN stations, whereas the drop of only 25% was measured at Artemis. This result was unexpected, as the predicted magnitude drop, based on the known magnitudes of the star (Table 1) and the TNO (from JPL ephemeris used in CORA), was estimated to be 5.7 magnitudes, equivalent to a 99.5% drop in flux. Liverpool and Oukaimeden observed no flux drop, and they are located in between Sicily (GAL Hassin and SLN stations) and Tenerife (Artemis and TAR2 stations), according to the direction of the shadow path (see Fig. 1). These are clear indications that the target star is an optical double or a binary: The telescopes in Sicily observed the occultation of the brighter member of the pair and those in Tenerife observed the occultation of the fainter member. As seen from Liverpool and Oukaimeden, 2007 OC10 passed in between the two stars, not occulting any of them. Therefore, we implemented a double star model in the analysis of the event (see Sect. 3.1).

The fact that the star is double is coherent with the Gaia DR3 Re-normalized Unit Weight Error value (RUWE) being >1.4 (see Table 1). The RUWE is expected to be close to 1.0 when a single-star model fits the astrometric observations adequately. A value noticeably higher than 1.0 (e.g. >1.4) may suggest that the source is either not a single star or presents challenges for the astrometric solution (Castro-Ginard et al. 2024).

2.4 Photometric measurements

Photometric observations of 2007 OC10 were made with the IO:O instrument of the 2 m Liverpool Telescope at La Palma Observatory, the 2k×2k CCD camera of the 1.5 m telescope at Sierra Nevada Observatory and the 4k×4k DLR CCD camera of the 1.23 m telescope at Calar Alto observatory. Exposure times ranged from 360 s to 500 s. The observations were performed with SDSS r, Johnson Cousins R, Johnson V, and with no filters. The data reduction and analysis were the same as those described in Ortiz et al. (2020a).

For the absolute calibration we used the Gaia DR3 stars in the field of view of each image and the Riello et al. (2021) transformation equations to the Johnson–Cousins system as explained in Morales et al. (2022). The uncertainty in the absolute calibrations changed slightly from night to night with the mean value at 0.05 mag.

The observed magnitudes were translated into reduced magnitudes by subtracting 5 log₁₀(R_hΔ), where R_h and Δ were the heliocentric and geocentric distances of 2007 OC10. For the absolute magnitude and phase slope determination, we used only the most recent data, those since JD 2459641.5 up to JD 2460559.5, to avoid any possible small change due to the slight variation of the aspect angle in several years. From a linear fit of the reduced magnitudes versus phase angle, the H_V determined is 5.40 ± 0.02 and the phase slope parameter is 0.18 ± 0.02 mag/deg. The fit was performed by removing outliers above and below 3σ. The H_R and phase slope derived from the data are 4.85 ± 0.02 and 0.204 ± 0.014 mag/deg respectively. The V − R is thus 0.55 ± 0.03.

The H_V magnitude in the JPL Horizons system is 5.24, but this is assuming a fixed phase slope. Given the very heterogeneous sources of the data, this value is highly uncertain. Nevertheless, it is not too far from the value determined here.

In Perna et al. (2013), a value of H_V = 5.36 ± 0.13 is reported in Table 3, based on their own observations, but this assumes a phase slope of 0.12 ± 0.06 mag/degree according to the table. However, when we correct their M_V(1, 1, 1.1) values reported in their Table 2² and take the mean, the obtained value is 5.59 ± 0.02. Given that the observations were obtained at 1.1 degrees of phase angle, the final value is (using the reported phase slope) H_V = 5.46 ± 0.07, not 5.36. If we use our phase slope parameter of 0.18 mag/deg we obtain H_V = 5.40 ± 0.03, perfectly compatible with our own derived value.

In Alvarez-Candal et al. (2019) photometric data of 2007 OC10 was also part of their data set. From an internal communication, they obtained values H_V = 5.71 ± 0.14 and H_R = 4.99 ± 0.14, which are also compatible with our results.

From our own photometric measurements, the peak-to-valley amplitude in the rotational light curve seems to be ≲0.1 mag nitudes. Due to this small photometric variation, no rotational period could be confidently recovered from our data.

3 Methodology

3.1 The model

During the stellar occultation, 2007 OC10 was at a geocentric distance of Δ = 35.532 AU. Hence, the Fresnel scale $F=\sqrt{\lambda \mathrm{\Delta }/2}$ at a wavelength of λ = 600 nm is F = 1.3 km. The stellar diameter can be estimated from the formulae in van Belle (1999). Using the magnitudes in Table 1, the angular diameter at the distance of 2007 OC10 is 0.14 km assuming a main sequence star, 0.27 km assuming a variable star and 0.17 km assuming a super giant. The data from the occultation revealed a double star, so these values represent an upper limit for the individual diameters.

The shortest exposure time in our data set is 0.3 s, from the Artemis Observatory. In that time interval, the sky plane covered distance was 7.1 km, according to the mean shadow velocity of 23.6 km/s. This value is considerably larger than that of the Fresnel scale or any possible stellar diameter. As a result, our light curves are mainly dominated by exposure times. We consequently modelled the occultation with a purely geometric approach; considering the stars as point-like sources and with no diffraction effects.

The canvas for our modelling is the sky plane. This plane is perpendicular, at any instant, to the radius vector connecting the centre of 2007 OC10 and the photocentre of the target stars (for a detailed explanation, see Elliot et al. 1978 and references therein). We set the origin of this plane to be the predicted barycentre of 2007 OC10 according to some reference ephemeris. In our case, we choose as reference ephemeris the updated version from CAHA data. The Cartesian axes are oriented so that X corresponds to the eastwards direction and Y corresponds to the northwards direction (see Fig. 3).

The projected shape of 2007 OC10 will be modelled as an ellipse. In light of the observed data, we need to include the presence of two stars in the occultation (see discussion in Sect. 2.3). Consequently, our model depends on the following 8 free parameters:

• a is the semi-major axis of the ellipse.

• b/a is the axis ratio (b represents the semi-minor axis).

• ϕ is the tilt angle of the ellipse: from the East direction to the semi-major axis. This angle is restricted to the [−90°, 90°] range, being positive when measured towards the North.

• (x_c, y_c) represents the position of the ellipse on the sky plane. These two parameters represent the offset of 2007 OC10 with respect to the reference ephemeris. Given that the offset is approximately constant during the short duration of the occultation, its value provides a new astrometric measurement of the 2007 OC10 for further refining its orbit.

• f₁ represents the normalised flux of the brightest star (primary star). As discussed in Sect. 2.3, the expected flux drop was 99.5% considering a single star. Since our light curve with the highest S/N is the one from Artemis with S/N=30, this residual brightness from 2007 OC10 will not be distinguishable from noise and is neglected. Therefore, the normalisation of the light curves implies 1 = f₁ + f₂, where f₂ = 1 − f₁ is the normalised flux of the fainter star (secondary star).

• (x₁, y₁) define the position of the primary star on the sky plane, taking the photocentre as origin: (x₁, y₁) = r_{star 1} − r_photocentre. The position relative to the photocentre for the secondary star (x₂, y₂) is now completely determined by $$(x_2,y_2)=-\frac{f_1}{1-f_1}(x_1,y_1).$$(1)

Fig. 3 Schematic representation of the model. The dashed line represents the path followed by the target star photocentre (αast, δast in Table 1) on the sky plane, as seen from a certain site. The stars are represented by the yellow dots; the size of the dot is proportional to the brightness. The origin on the sky plane (black dot) is determined by the theoretical astrometric position of 2007 OC10 according to our reference ephemeris. In this particular configuration of parameters, the site for which the chord is represented would observe the occultation, with a residual flux in the drop of f1 (secondary star occulted, primary star visible).

3.2 The Bayesian approach

The traditional analysis of a stellar occultation is usually performed as follows. Firstly, the obtained light curves are identified as positive (flux drop detected) or negative (no detection). In most cases, diffraction or extensive source effects are negligible (see Sect. 3.1), so a square-well is fitted to all the positive light curves. This determines the ingress and egress times (beginning and end of the drop). Then, they are compared to synthetic times produced by certain shape models. The differences between synthetic and measured times are minimised following an iterative χ² approach. Many TNOs have been successfully modelled with this procedure (Santos-Sanz et al. 2022; Pereira et al. 2023).

However, this methodology can be impractical for moderate to low signal-to-noise ratios in some of the light curves, as is the case in the present work. This situation can be aggravated in the scenario of a double star being occulted. If one of the stars is significantly fainter, stations detecting only the occultation of this star might not be able to distinguish the small flux drop from the noise. Therefore, we adopted an approach where each data point from each light curve is treated individually, comparing them to synthetic light curves produced by a certain model. In this way, a flux drop smaller than noise level due to the occultation of the faint star is naturally covered, without the need to decide if the occultation was actually detected or not. This is the case of TAR2 data. Ingress and egress times should have been computed in a classic approach. Clearly, this is not an easy task when the data quality is not optimal (see Fig. 2). Worse still, in the absence of Artemis data, the TAR2 light curve could have been mistakenly identified as negative, yielding to misleading conclusions. In this aspect, our methodological approach is advantageous, as no a priori assumptions need to be made about the nature of the light curves.

Although this is not the present case, the situation in which a single station observes the blended occultation of both stars would be treated the same way with this approach. In contrast, determining ingress and egress times for each star would imply fitting increasingly complex stepped functions, which might be cumbersome in noisy data sets.

To fit the synthetic data curves produced by our model to the observed ones, we follow a methodology similar to that of Leiva et al. (2020) or Strauss et al. (2021). We use Bayesian statistics to determine a posterior probability density function (posterior pdf or simply posterior) for the parameters described before under our adopted model. The aim is to characterise this posterior. This approach is advantageous because it provides a better understanding of parameter correlations, which are likely to appear in double star occultations. By characterising the full posterior distributions, we are also able to comprehensively quantify uncertainties (Gregory 2005).

The posterior pdf P(θ|D) is constructed from the likelihood function ℒ(|θ), a prior distribution of the model parameters P(θ) and a normalisation factor P(D): $$P(\theta \mid D)=\frac{\mathcal{L}(D\mid \theta )P(\theta )}{P(D)}$$(2)

Here D represents the obtained data set, comprising each normalised flux value f_i for each light curve. θ is our 8-dimensional parameter vector. The likelihood measures the plausibility of observing the data set D given the parameter values θ. The prior distribution reflects our a priori knowledge about the model parameters. In this work, we consider the simplest case where P(θ) is the product of the marginal prior distributions for each parameter (no correlations).

Modelling the normalised flux errors σ_i as normally distributed and uncorrelated, the likelihood function is given by (Gregory 2005) $$\mathcal{L}(D\mid \theta )=(2\pi {)}^{-N/2}\left\{\prod _{i=1}^N\frac{1}{{\sigma }_i}\right\}\times \{\exp [-\sum _{i=1}^N\frac{{(f_i-m(t_i\mid \theta ))}^2}{2{\sigma }_i^2}]\},$$(3)

where N is the total number of data points and m(t_i|θ) is the modelled light curve value under the parameters θ at time t_i, corresponding to the data point f_i.

To characterise the posterior pdf, we applied a Markov chain Monte Carlo (MCMC) method to obtain a representative sample from P(θ|D). We adopt the use of Python package emcee 3.1.6 (Foreman-Mackey et al. 2013), which implements the affine invariant sampling algorithm described in Goodman & Weare (2010).

In our samplings, we used n_w = 512 random walkers, which are initially distributed according to the selected prior. The sampler is initially run for n_burn steps to ensure convergence of the sampling. This is the burn-in phase, where the walkers explore the parameter space. We checked that the number of steps in this phase is at least 10 times larger than the estimated autocorrelation time in each parameter, to ensure that we are not sampling far from convergence. From the last position of the burn-in phase, the sampler is run for n_iter additional steps to obtain a representative sample. The final sample size is n_w × n_iter (we discard the burn-in phase).

Our MCMC samplings were performed with n_burn = 3000 and n_iter = 1000. Since P(θ|D) does not show signs of being multi-modal (just one set of parameters satisfactorily explains the data), it is sensible to present the median of the corresponding marginal distributions as the nominal values of each parameter. The provided errorbars are taken from 68% maximum-density credible intervals. That is, the shortest interval containing 68% of the sample elements, corresponding to the 1σ region.

4 Occultation results

4.1 Bayesian analysis: Prior selection and results

In Bayesian analysis, prior selection is a sensitive matter as it can significantly influence the posterior in cases of low S/N data. The only a priori information we have about 2007 OC10 is the published result from thermal data. Population-wise, we also have certain information from published TNO occultations about the shape distribution. We report the results of two different samplings taking these into account. We aim to assess the relative influence of the priors and occultation data in the posterior.

Our preferred sampling (sampling 1) will have uniform priors in all parameters. This is the most uninformative choice and serves only a practical purpose of limiting the parameter-space region the walkers must explore. We selected this as our preferred sampling, since we intend to compare the results on the object size from thermal data reported in Santos-Sanz et al. (2012). Including these as priors would bias our sampling beforehand, invalidating the comparison.

Our sampling 2 includes as priors the area-equivalent diameter from thermal data and shape distribution of TNOs from stellar occultations. Since D_eq is not directly a parameter of our model, we need to introduce prior distributions in a and b/a so that the resulting prior is compatible with a normal centred in the nominal D_eq value and σ equal to the reported uncertainty. To do so, we considered that the best practical choice for the b/a prior would be a normal fit to the observed b/a values in other TNO occultations. Then, the prior in semi-major axis can be chosen accordingly.

In Table 2, we collect all the published values for projected axis ratio in TNO stellar occultations. We acknowledge that the sample of TNOs with published occultations is biased towards the largest objects. This highlights the need to accumulate more data if we intend to perform statistically significant analyses.

The fitted truncated normal between 0 and 1 for these axis ratio values is N(μ = 0.889, σ = 0.143) (see Appendix A for details). If a normal prior N(μ = 169 km, σ = 17 km) is chosen for a, the resulting a priori distribution for D_eq resembles the desired distribution from thermal data N(μ = 309 km, σ = 37 km): The mean is μ = 308.1 km, the median is M = 308.2 km and the standard deviation is σ = 37.1 km. To further confirm the compatibility, we performed a Kolmogorov-Smirnov test comparing a sample of this a priori D_eq distribution against N(μ = 309 km, σ = 37 km). With a sample size of 1000, the test yields a p-value=0.61. Therefore, it is unjustified to reject the null hypothesis: The sample comes from the desired distribution.

For sampling 1 and sampling 2, we selected uniform priors in the position parameters (x_c, y_c). We know that the actual position (x_c, y_c) of 2007 OC10 on the sky plane has to be between the two flux-drop areas. And it should be closer to the drop from SLN and GAL Hassin, sites where the occultation by the primary star was detected. Since f₁ ≃ 0.75, from Eq. (1) we can infer that (x_c, y_c) should be three times closer to this flux drop area than the one from Artemis and TAR2 (secondary star). With this in consideration and visually inspecting the light curves plotted on the sky plane (see Fig. 4), we choose an uniform prior in the region [−900, −300] km × [−500, 400] km, confidently covering the actual value of (x_c, y_c). This limits the futile exploration of regions where the ellipse positions on the sky plane are far away from the actual flux drop areas.

The situation for (x₁, y₁) is similar to that of (x_c, y_c), since the double nature of the star was unknown prior to the occultation. Hence, we base the prior selection in practical purposes. From visual inspection, the primary star position relative to the photocentre should be in the x₁, y₁ < 0 region, to cause an occultation visible from SLN and GAL Hassin (see Fig. 3). Consequently, a conservative uniform prior is chosen in the region [−150, −500] km × [0, −350] km to limit the walkers from exploring unnecessary parameter space volume.

For the rest of parameters that have not been discussed, the priors are the same for sampling 1 and sampling 2:

• The possible orientation of the projected shape was unknown a priori, hence we set a completely uninformative prior for ϕ: U[−90°, 90°].

• As the star was an unknown double, there is no a priori information about the possible value of f₁. Hence, we set an uninformative prior U[0.5, 1].

We represent the resulting samples in a pairwise marginal plot in Fig. D.1. From this representation it is clear that the samplings show no signs of multi-modality, so there is no degeneracy in the solutions. The clearest correlations appear between the positional parameters x_c, x₁ and f₁, which is not unexpected. Slight variations of the position of 2007 OC10 and the primary star can produce the same shadow, as long as f₁ is adjusted so that the secondary star shadow also stays in place. However, since the value f₁ is accurately determined by Artemis light curve, this variation is very well constrained. Another (less marked) correlation can be found between b/a and y_c: if the position of 2007 OC10 on the sky plane was further south, a moderately rounder projected shape would be admissible.

The nominal parameter values with the corresponding uncertainties for both samplings are given in Table 3. Figure 4 offers a qualitative visualisation of the result of our preferred sampling 1. A sub-sample of 100 parameter vectors is extracted randomly from the full posterior sample. Those belonging to the 68% credible region of each parameter space are shown by drawing the flux drop area (i.e. the projected elliptical limb) on the sky plane expected for this parameter configuration. The expected shadows for the nominal parameters in Table 3 are also represented for reference.

Along both samplings, the obtained values for the model parameters are compatible within the provided (1σ) credible intervals. This is an indication that the sampling is heavily influenced by the obtained occultation data, and not by the chosen priors.

The axis ratio parameter is the one suffering from greater variation between sampling 1 and sampling 2. This is expected since the prior in b/a for sampling 2 gives preference to rounder objects, favouring a higher median in the posterior marginalised distribution for b/a. However, this axis ratio prior is based on a scarce sample, without significant representation of objects of the same size or dynamical class as 2007 OC10. Therefore, we prefer to report sampling 1 result as more realistic and unbiased.

The differences in axis ratio due to different priors are translated in a difference of the obtained area-equivalent diameters: D_eq = ${330}_{-55}^{+56}$ km in sampling 1 and D_eq = ${343}_{-34}^{+40}$ km in sampling 2. In Fig. 5, it is shown how the different prior distributions in D_eq affect the resulting sample. The combination of a ~ U[100, 300] km and b/a ~ U[0.0, 1.0] (sampling 1) results in a D_eq prior with a peak at ~200 km, but very long tails compared to the N(309, 37) km prior distribution in sampling 2. Despite the notable differences in these priors, the results obtained are compatible within the credible interval 1σ. Again, this supports our claim that the posterior density is primarily influenced by the likelihood (determined by the occultation data) and not the priors.

Fig. 4 Occultation data plotted on the sky plane with the results of sampling 1. The X − and Y-axes are pointing towards east and north, respectively. The arrow indicates the relative direction of motion of the occulted stars. Each line segment represents an exposure. The transparency of each segment corresponds to the normalised flux measured in that integration (see Fig. 2). No transparency means flux = 1. In the light curve from Çanakkale there is an abnormally long dead-time between two concrete exposures. The corresponding blank space should not be mistaken for a flux drop. The pair of black ellipses represent the nominal solution reported in Table 3 for sampling 1. Their centres are marked with red dots. The thinner grey ellipses represent the 1σ uncertainties from the extracted sample, as discussed in the text.

Fig. 5 Histogram comparison of the posterior samples in Deq. The priors are plotted for reference (dashed lines). Sampling 1 is our preferred sampling, with uniform priors in all our model parameters. However, as Deq is not directly a model parameter, the prior in this value is not uniform. In sampling 2, the prior in Deq is equivalent to the value reported in Santos-Sanz et al. (2012) from thermal data. This prior pdf peaks at a higher Deq value than the prior pdf of sampling 1. Consequently, the posterior pdf in sampling 2 also peaks at a higher value than in sampling 1. However, both distributions yield compatible results within 1σ (68% credible intervals). Moreover, the result from sampling 1 is also compatible with the result from thermal data (sampling 2 prior).

4.2 Projected axis ratio of the object

The determination of the projected axis ratio of the object b/a = 0.58 ± 0.16 (or equivalently, b = ${141}_{-23}^{+24}\text{km}$) is of great relevance. This adds another data point to the scarce data set on projected TNO shapes (see Table 2). The importance of this value is twofold.

Firstly, the scarcity of the data is further emphasised when dividing the sample according to dynamical classes. This division might be of interest if we seek to identify differences in the shape distribution according to origin, since we know that the evolutional history varies between dynamical classes (Nesvorný 2018). As stated previously, this is only the second detached TNO to be measured with stellar occultations.

Secondly, 2007 OC10 is only the third-smallest TNO with a measured projected shape through occultations. This is of remarkable interest, as shapes are expected to depend notably on size (Tancredi & Favre 2008). Theoretically, there should be a certain critical diameter above which shapes should be governed by hydrostatic equilibrium (Chandrasekhar 1987). However, there is no empirical evidence of where exactly this limit is located. In fact, even the existence of a critical diameter is to be confirmed, since some large (D_eq > 500 km) TNOs depart from hydrostatic shapes of equilibrium (Ortiz et al. 2017; Vara-Lubiano et al. 2022). Therefore, it is important to explore the lower end of the size parameter space in order to acquire sufficient data to test this hypothesis with some statistical significance.

4.3 Surface properties: Geometric albedo and beaming parameter

In this section, we combined our occultation results from the elliptical model with photometric and thermal measurements to refine our knowledge about the geometric albedo p_V and beaming parameter η. We contrast them with the values derived in the literature from the TNOs are Cool survey (Santos-Sanz et al. 2012; Lellouch et al. 2013).

We acknowledge that the projected area (and consequently D_eq) at the instant of the occultation could differ from the one at the epoch of thermal measurements. Due to the absence of information on the three-dimensional shape of 2007 OC10 and, in particular, its rotation axis, we can only statistically estimate the variation in projected area caused by rotation and changes in aspect angle. The statistical analysis is performed in Appendix B. While the difference in projected area between the two epochs may be significant, for the purposes of this study, we assume that the projection observed during the occultation and thermal measurements remains the same.

The geometric albedo and beaming parameter values obtained in Santos-Sanz et al. (2012); Lellouch et al. (2013) from thermal data come from a modified NEATM model (Harris 1998). The only difference with the original NEATM model is that zero phase angle is assumed. The NEATM model assumes a spherical body in instantaneous equilibrium with the solar radiation, in the limit of zero thermal inertia. A socalled beaming parameter η is introduced to account for surface roughness and the fact that the thermal inertia is non-zero. The two unknown parameters that completely determine the thermal emission in the NEATM model are the geometric albedo p_V and the beaming parameter η. The absolute magnitude H_V of the object is also needed to determine the diameter of the sphere when combined with p_V. The absolute magnitude is obtained from independent photometric measurements.

Within the TNOs are Cool survey, the thermal flux from (470316) 2007 OC10 was observed with Herschel/PACS at three different wavelengths: 70 μm, 100 μm and 160 μm. As discussed in Santos-Sanz et al. (2012), only these three data points were considered insufficient to fit p_V and η simultaneously. Hence, the authors adopted an approach in which the beaming parameter is fixed to η = 1.20 ± 0.35 (the mean value known at the time for TNOs; Stansberry et al. 2008), so they fit p_V as the only free parameter.

Using the absolute magnitude derived from our own observations of H_V = 5.40 ± 0.02 (see Sect. 2.4), the obtained value for the geometric albedo is p_V = ${11.2}_{-5.0}^{+2.1}\%$. This is consistent with p_V = ${12.7}_{-2.8}^{+0.4}\%$ from Santos-Sanz et al. (2012). A H_V value of 5.43 ± 0.10 was used on their work, similar to that derived in the present paper. If we adopt their value of absolute magnitude, the variation in our result is negligible (see Table 4).

Lellouch et al. (2013) reanalysed the same thermal data for 2007 OC10 reported in Santos-Sanz et al. (2012) but fitting both parameters p_V and η, instead of using a fixed η value (see Table 4). The uncertainties are considerably larger than in Santos-Sanz et al. (2012), and η is barely constrained. The beaming parameter for 2007 OC10 is effectively unknown. Occultation data can help constrain this value. Because we have now independently determined p_V, we can fit η as the only free parameter. For comparability, we adopt the same modified NEATM model, phase integral approximation, and fitting method used in Santos-Sanz et al. (2012) and Lellouch et al. (2013).

Adopting the most recent and precise absolute magnitude value H_V = 5.40 ± 0.02 derived in this work, the resulting beaming parameter is η = ${1.42}_{-0.58}^{+0.75}$. The value obtained for H_V = 5.43 ± 0.10 is also reported in Table 4 for completeness. There is a considerable reduction in the uncertainty of η with respect to the value reported in Lellouch et al. (2013). Our result shows that 2007 OC10 has a beaming parameter compatible with the mean of the population. However, the relative error is still ~50% in η, preventing us from drawing any further conclusions.

4.4 Double star characterisation

As a by-product of the shape analysis of 2007 OC10, we have constrained the relative position of the occulted stars at the epoch of the occultation: 2022-08-22T01:32:37 UTC. The fitting of (x₁, y₁) and f₁ completely determines the relative position of both stars. From our preferred sampling, we computed an angular separation of ${57}_{-11}^{+4}$ mas and a position angle of the secondary star with respect to the primary star (measured from north towards the East) of θ = ${56}_{-17}^{+3}$ degrees.

In addition, the fitting of f₁ allows us to conclude that the instantaneous magnitude difference between stars was ${1.18}_{-0.07}^{+0.07}$ magnitudes.

4.5 Astrometric constraints

Occultations are a source of mas-level astrometry for TNOs (Rommel et al. 2020). In this particular case, the occultation also enables us to constrain the astrometric position of both stars. All these astrometric constraints rely on the precise determination of the double star photocentre at the occultation epoch (see Table 1). To compute this position (α_ast, δ_ast) the parallax and proper motion needs to be taken into account. However, the duplicity nature of the star was unknown when computing these values, which might be a source of systematic error. As discussed in Leiva et al. (2020), in case both stars are gravitationally bound, then the Gaia DR3 proper motion reflects a combination of the star system barycentric motion and the movement of the photocentre relative to the barycentre. On the contrary, if the stars are just aligned by chance, then Gaia DR3 proper motion represents a combined effect of the individual proper motion and parallax of each star. Trying to discern between both cases is neither possible with the occultation data nor in the scope of this work.

We report the instantaneous astrometric positions of 2007 OC10 and both stars in Table 5. We note that our model assumes a fixed position for the photocentre with no uncertainty. Therefore, the uncertainties in Table 5 combine the uncertainties from the occultation sampling result and those from the Gaia astrometric catalogue (α_ast, δ_ast in Table 1). The latter ((0.1 mas)) are assumed to be normally distributed and independent, which is a fair approximation since the dominating source of error comes from the sampling ((1–10 mas)).

5 Comparison with thermal data. Upper limit on the size of a putative satellite

There are no secondary flux drops in the occultation light curves that suggest a secondary object, but Fig. 4 shows that the sky plane area covered by the data has about 1500 km in width, with gaps between sites as large as 500 km. A secondary object in the gaps between observing sites or outside the band covered by the occultation data would remain undetected. The thermal data, taken with Herschel/PACS, have an angular resolution of 5 arcsec (Poglitsch et al. 2010). The angular resolution corresponds to about 130 000 km at the object’s distance during the thermal measurements (Santos-Sanz et al. 2012), a resolution too coarse to detect any of the known trans-Neptunian binaries (Noll et al. 2020). There are no published imaging data from the Hubble Space Telescope (HST), which could reach resolutions of about 1000 km at the distance of the Kuiper belt. Considering the limited angular resolution of the occultation and thermal imaging data, we cannot directly rule out the presence of a yet undetected secondary object.

Fitting a single-object model to a hypothetical double-object to the thermal measurements will result in a systematically larger thermal diameter to account for the additional flux from the secondary. Neglecting differences due to changing geometry (see discussion in Sect. 4.3 and Appendix B), an unresolved secondary object will show up as a larger thermal diameter with respect to the occultation data. This has been previously proposed for the TNOs 2002 MS4 and 2002 TC302 (Rommel et al. 2023; Ortiz et al. 2020a). Table 4 shows that the area-equivalent diameter from our elliptical model is compatible, within error bars, with the thermally derived diameter values. The relatively large uncertainties in the diameter from both independent measurements do not discard cases of a larger thermal diameter than the occultation diameter. Thus, we cannot categorically rule out a secondary, although the similarity in both thermal and occultation values allows us to constrain its maximum size.

In this section, we adopt a novel statistical approach to derive the upper limit on the size of a hypothetical secondary object, using the occultation-derived diameter for the primary and thermal data to constrain the diameter of the secondary. We adopted a simple model consisting of primary and secondary objects modelled with the same modified NEATM from Santos-Sanz et al. (2012). Both objects have the same geometric albedo and beaming parameter, differing only on their diameters, with the constrain that the combined brightness of both objects must yield an absolute magnitude of the system of H_V = 5.4. The model has three independent parameters: The diameter of the primary D_eq, the diameter of the secondary d_eq, and the beaming parameter η.

We implemented a Bayesian approach to derive posterior distributions similar to the occultation analysis. The likelihood has the same functional form as Eq. (3), where f_i and σ_i are the thermal fluxes and their corresponding uncertainties, and m(t_i|θ) is the modelled flux (sum of the two bodies) at each of the N = 3 wavelengths t_i for the model parameters θ = (D_eq, d_eq, η). As the equivalent diameter derived from the occultation is a direct measurement, we used the posterior probability distribution in Fig. 5 (sampling 1) as the diameter of the primary object, leaving only d_eq and η as free parameters to explore.

We used the TNO sample of beaming parameters from the TNOs are Cool programme (Lellouch et al. 2013, Table 3) to derive a prior distribution for the beaming parameter η, excluding centaurs to be consistent with the previous analysis (see Appendix A for details). Since we have no a priori information about the putative satellite, we chose a uniform prior U[0, 400] km for d_eq, excluding exploring the unphysical region of negative diameters.

Figure 6 shows the pairwise marginal plot obtained after performing a MCMC sampling of the posterior distribution using the emcee package with n_w = 256 random walkers, n_iter = 25 000 iterations steps, and n_burn = 10 000 burn in steps. We verified that the burn in phase is at least 10 times the estimated autocorrelation time to ensure convergence of the posterior distribution.

The posterior pdf reaches its peak at d_eq = 0 km indicating that the occultation and thermal data are better explained by a single object in our simple model. The presence of a sizeable satellite cannot be completely ruled out with the available data as indicated by the posterior of d_eq having an extended wing. In Fig. 6, we plot the 80%, 90% and 95% percentiles of the marginalised distribution in D_eq.. They represent quantitative upper limits for the size of a putative satellite. The interpretation is that, according to our sampling and model, if 2007 OC10 has a satellite there is a 95% probability that the satellite is smaller than 261.8 km in diameter, a 90% probability that it is smaller than 225.3 km and a 80% that it is smaller than 178.2 km.

The lower-left panel in Fig. 6 shows a strong correlation between η and D_eq.. A higher beaming parameter implies a reduced thermal emission, and hence larger objects can be allocated. The diameters of the primary and secondary are mildly correlated, as smaller primary diameters allow the presence of a larger satellite. We also observe that the posterior distribution of the beaming parameter shows a secondary peak at η ~ 2.3, resulting in a higher median η than the single-body model. The extended wing in the joint marginal distribution of η and D_eq (Fig. 6, lower-left panel) indicates that the feature in the posterior of η is produced by the distribution of D_eq and the right-skewed prior in η (see Appendix A). The larger values in D_eq allow the larger values of η in the prior to show up in the posterior.

Fig. 6 Pairwise marginal plot for the sample of our thermal binaryobject model MCMC run. The marginalised distribution in Deq is the same as the posterior from the occultation (sampling 1 in Fig. 5). We show it in red to highlight that this parameter is not being updated in this sampling. The regions depicted in the 2D histograms encircle the 39.4% (1σ) and 78.8% (2σ) of the samples. The coloured lines represent quantiles of the marginalised deq distribution. These are interpreted as upper limits for the diameter of the satellite.

6 Summary and conclusions

On 2022 August 22, the detached trans-Neptunian object (470316) 2007 OC10 caused an occultation of the Gaia DR3 star 2727866328215869952. The event was observed from 18 sites, resulting in 13 light curves (3 lost due to bad weather, 2 discarded due to data quality). Partial flux drops of ~25% and ~75% were observed from different locations. This is a distinctive sign that two different stars have been occulted, confirming the discovery of a double star. No secondary features indicating rings, debris, or satellites were detected in the light curves.

This is the third-smallest TNO with an occultation-derived shape to date, third only to Arrokoth (Buie et al. 2020) and 2003 UY117 (Kretlow et al. 2024a). It is clear that this sample should be further enlarged if we intend to determine the critical diameter (material-dependent) that sets a limit above which TNOs adopt hydrostatic equilibrium shapes, as proposed in Tancredi & Favre (2008). This is not an easy task since the bulk of TNOs similar in size to 2007 OC10 have ~300 mas uncertainty in their ephemerides (Ortiz et al. 2020b), which translates to ~8000 km in uncertainty for the shadow path. Consequently, we must carry out dedicated astrometric campaigns, such as those conducted in the present work, to correctly predict stellar occultations.

G!kún||'hòmdímà is the only other detached TNO with a published physical characterisation through occultations. This population is believed to have suffered a dynamical evolution similar to the scattered TNOs and substantially different to the Cold Classicals (Benedetti-Rossi et al. 2019). Whether the difference in evolutionary paths is somehow reflected in the shape is still impossible to determine because of the scarcity of the sample.

An eight-parameter model was applied to the occultation light curves using a Bayesian approach, employing the emcee Python package to derive the posterior probability distribution. Our samplings were run for enough time to ensure convergence. A prior sensitivity analysis shows that the parameter posteriors are strongly dependent on the collected occultation data and are robust against variation in prior distributions.

From our preferred sampling we have determined a projected axis ratio of b/a = 0.58 ± 0.16, which is the smallest value to date for a TNO (see Table 2). The obtained area-equivalent diameter is D_eq = ${330}_{-55}^{+56}$ km and the geometric albedo p_V = ${11.2}_{-5.0}^{+2.1}\%$. These are compatible with the thermal data from Herschel/PACS, available in Santos-Sanz et al. (2012).

The present work demonstrates that stellar occultations represent a great opportunity to complement, and sometimes reinterpret, the thermal data from the TNOs are Cool programme. We were able to reduce the uncertainty on the value of the beaming parameter for 2007 OC10 from η = ${1.17}_{-0.64}^{+1.11}$ (Lellouch et al. 2013) to ${1.42}_{-0.58}^{+0.75}$. Proceeding similarly for other TNOs are Cool objects only observed with Herschel/PACS (poorly constrained η) would be of great interest. In this way, a more statistically significant analysis could be performed on the surface properties of TNOs from the distribution of beaming parameters and the correlation with other parameters.

The results from this paper also evidence that stellar occultations enable milliarcsecond-level astrometry for distant Solar System objects (as put forward in Rommel et al. 2020) and the stellar components in the case of double sources (see Table 5). As shown, this level of accuracy can be achieved in the case of moderate S/N light curves, even if they only scan a limited region of the body’s limb (Fig. 4).

The fact that the occultation-derived diameter is compatible with that from thermal data is, in principle, evidence against the existence of any putative large satellite. Satellites would increase the measured thermal flux and the estimated thermal diameter would be systematically enlarged. The discussion usually ends here for past works analysing stellar occultation results for TNOs are Cool targets. We implemented a novel methodology to quantitatively set an upper limit on the size of a putative satellite for objects with thermal and occultation data.

Although our data are well explained by a single object model, the presence of a satellite of size comparable to the main body cannot be completely ruled out due to the uncertainties in the occultation size and thermal flux. We derive an upper limit of 225.3 km for a putative satellite. A better constrain in the diameter from occultations, thermal measurements with a better wavelength coverage or smaller thermal flux uncertainties would put tighter constrains in a putative satellite. Although our analysis involves a simple model, we emphasise that simple conclusions from the compatibility of thermal and occultation diameters must be drawn with caution.

Imaging from HST could detect a satellite if the projected separation is greater than ~1000 km at 2007 OC10 distance. A close-in orbit is also possible since the occultation sites have several >250 km gaps where the satellite could be allocated. Determining the fraction of tight binaries is of great interest since it can give us insights into key parameters of Solar System evolution, such as the lifespan of the massive disk (Nesvorný & Vokrouhlický 2019). The methodology outlined in this paper offers a valuable approach to addressing this issue, particularly in cases involving well-constrained occultations and thermal data with good wavelength coverage and S/N.

To further refine the physical characterisation of 2007 OC10, it would be ideal to determine the rotation period and establish the rotational phase at the occultation epoch. This would help us correct for some systematic effects since the analysis of the beaming parameter assumes that the occultation and photometric measurements were taken at the same phase as the occultation. However, this assumption is reasonable since the rotational light curve amplitude of 2007 OC10 is less than 0.1 magnitudes. Therefore, the rotational variability of the area-equivalent diameter is ≲5%.

The determination of a rotational period would also help to constrain the three-dimensional shape, allowing us to reanalyse the thermal data with a thermophysical model to better constrain the object surface properties.

Further stellar occultations would also greatly contribute to improving 2007 OC10 models. Inconveniently, 2007 OC10 is currently moving in a region very sparsely populated by stars, reducing the chances of detecting future events, but still some opportunities are available for medium and large telescopes in the 2025–2030 time frame.

Data availability

Occultation light curve is available at the CDS via anonymous ftp to cdsarc.cds.unistra.fr (130.79.128.5) or via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/697/A157

Acknowledgements

The authors thank the referees and the editor for their valuable feedback. This publication is funded by the Spanish Ministry of Universities through the university training programme FPU2022/00492. The authors of IAA-CSIC acknowledge financial support from the Severo Ochoa grant CEX2021-001131-S funded by MCIN/AEI/10.13039/501100011033. Part of this work was supported by the Spanish projects PID2020-112789GB-I00 from AEI and Proyecto de Excelencia de la Junta de Andalucía PY20-01309. P.S-S. acknowledges financial support by the Spanish grants PID2022-139555NB-I00 and PDC2022-133985-I00 funded by MCIN/AEI/10.13039/501100011033 and by the European Union “NextGenerationEU”/PRTR. This work is based on observations collected at the Centro Astronómico Hispano en Andalucía (CAHA) at Calar Alto, operated jointly by Junta de Andalucía and Consejo Superior de Investigaciones Científicas (IAA-CSIC). This work is based on observations made with the Liverpool Telescope operated on the island of La Palma by the Instituto de Astrofísica de Canarias in the Spanish Roque de los Muchachos Observatory. This work uses data from GAL Hassin – Centro Internazionale per le Scienze Astronomiche, Via della Fontana Mitri, 90010, Isnello (PA), Italy. IST60 telescope and its equipment are supported by the Republic of Turkey Ministry of Development (2016K12137) and Istanbul University with the project numbers BAP-3685, FBG-2017- 23943. Based on observations collected at the Observatorio de Sierra Nevada, operated by the Instituto de Astrofísica de Andalucía (IAA-CSIC). RIS acknowledges financial support from the Spanish Ministry of Science and Innovation through the project PID2022-138896NA-C54. TSR acknowledges funding from Ministerio de Ciencia e Innovación (Spanish Government), PGC2021, PID2021-125883NB-C21. FBR acknowledges CNPq grant 316604/2023-2. This work was (partially) supported by the Spanish MICIN/AEI/10.13039/501100011033 and by “ERDF A way of making Europe” by the “European Union” through grant PID2021-122842OB-C21, and the Institute of Cosmos Sciences University of Barcelona (ICCUB, Unidad de Excelencia ‘María de Maeztu’) through grant CEX2019-000918-M. This work makes use of observations from the Las Cumbres Observatory global telescope network. Partially based on observations made with the Tx40 telescope at the Observatorio Astrofísico de Javalambre in Teruel, a Spanish Infraestructura Cientifico-Técnica Singular (ICTS) owned, managed, and operated by the Centro de Estudios de Física del Cosmos de Aragón (CEFCA). Tx40 is funded with the Fondos de Inversiones de Teruel (FITE). J.d.W. and M.I.T. gratefully acknowledge financial support from the Heising-Simons Foundation, Dr. and Mrs. Colin Masson and Dr. Peter A. Gilman for Artemis, the first telescope of the SPECULOOS network situated in Tenerife, Spain.

Appendix A Prior distributions

Appendix A.1 Occultation sampling: b/a

To fit a truncated normal distribution to the b/a available data points in Table 2, we performed a maximum likelihood estimation. To find the point (μ,σ) where this maximum is achieved, we implemented the Broyden, Fletcher, Goldfarb, and Shannon (BFGS) method (Mikosch et al. 2006) with the mean and standard deviation of the sample as initial points. The resulting optimal truncated normal has parameters μ = 0.889 and σ = 0.143. It is plotted in Fig. A.1 together with a histogram of the sample of b/a values.

Fig. A.1 Histogram of the available sample of projected axis ratios in published TNO occultations together with our fitted truncated normal. The latter is the used prior distribution in b/a in sampling 2.

Then, to obtain a sample of D_eq values from this b/a prior combined with our N(μ = 169 km, σ = 17 km) prior in a, we extracted a 1000 size sample from both prior distributions and combine them through D_eq = 2a²(b/a). The result is a 1000 size sample of diameters which is used in the K-S test against the distribution from the thermal data.

Appendix A.2 Thermal data sampling: η

In this case, as it is clear form the histogram in Fig. A.2, a normal distribution is not appropriate. Therefore, we implemented a Gaussian Kernel Density Estimation with the gaussian_kde class of scipy.stats Python package. For the bandwidth selection, we implemented Scott’s rule (Scott 2015). The resulting prior pdf for η is shown in Fig. A.2

Fig. A.2 Histogram of the available sample of beaming parameters from TNOs are Cool (Centaurs excluded; Lellouch et al. 2013) together with our Gaussian kernel density estimated pdf. The latter is the used prior distribution in η in the sampling in Sect. 5.

Appendix B Projected area change between occultation and thermal measurement epoch

To statistically estimate the change in projected area between the mean epoch of the radiometric measurements with Herschel/PACS, 2010 October 17, and the occultation data, 2022 August 22, we considered an ellipsoidal three-dimensional shape for 2007 OC10. We extracted 100000 random ellipsoids and rotation poles. The rotation poles are extracted from an isotropic distribution on the sphere. The axis ratios of the ellipsoid a/b and c/a are randomly selected from uniform distributions of U[0.5, 1] and U[0.3,b/a] respectively, excluding only very elongated shapes from our analysis. The sampled shapes and rotation poles were then filtered. We exclude those projecting an ellipse incompatible with an axis ratio 0.58 ± 0.16 on 2022 August 22, as observed in the occultation. We also exclude those that, on this same epoch, show a rotational light curve amplitude larger than 0.1 mag (considering flux proportional to the projected area). For those ellipsoids that are compatible with our observations, we then extracted a random rotational phase to compute the projected area on 2010 October 17, A_thermal. This area can be compared with that projected in the occultation A_occ.

In Fig. B.1, the resulting distribution of area differences (relative to A_occ) is shown. According to this analysis, there is a 50% probability that the area difference is less than 10%, and a 86% probability that the area difference is less than 25%. However, large area differences are plausible.

Fig. B.1 Histogram of the relative area difference (|Athermal − Aocc|/Aocc) for our random sample of ellipsoids and rotation poles. The epoch for Athermal is 2010 October 17 and the rotational phase is random. The epoch for Aocc is 2022 August 22 and the rotational phase is such that a projected axis ratio compatible with the occultation is achieved.

Appendix C Observing sites data

Appendix D Pairwise marginal plot for the obtained posterior sample

Fig. D.1 Pairwise marginal plot for our two MCMC samples of the posterior density. There are no signs of other modes representing alternative solutions. In the bidimensional plots, the represented region encircles 39.4% of the samples (representing the 1σ region). The red shading in these regions is a bidimensional histogram for our preferred sampling 1.

References

All Tables

All Figures

Fig. 1 Comparison of occultation predictions. At this time, the double nature of the star was unknown. The continuous lines represent the predicted path of the centre of the shadow. The 1σ regions along the path are represented by the broken lines. The red diamonds represent the stations that collected data. The grey squares represent the sites that suffered from bad weather or whose data were not used in the analysis. Some stations have been slightly shifted in position for improved visualisation, but the exact coordinates used in the analysis are those in Table C.1.

In the text

Fig. 2 Light curves of the stellar occultation (normalised flux vs time at mid-exposure). The origin of time is set to be the time of closest approach: 2022-08-22T01:32:37 UTC. They are ordered by distance perpendicular to the predicted shadow path (as predicted by CAHA data). The error bars are too small to be seen in the plot in most light curves, but they can be checked in the supplementary data on CDS.

In the text

Fig. 3 Schematic representation of the model. The dashed line represents the path followed by the target star photocentre (αast, δast in Table 1) on the sky plane, as seen from a certain site. The stars are represented by the yellow dots; the size of the dot is proportional to the brightness. The origin on the sky plane (black dot) is determined by the theoretical astrometric position of 2007 OC10 according to our reference ephemeris. In this particular configuration of parameters, the site for which the chord is represented would observe the occultation, with a residual flux in the drop of f1 (secondary star occulted, primary star visible).

In the text

Fig. 4 Occultation data plotted on the sky plane with the results of sampling 1. The X − and Y-axes are pointing towards east and north, respectively. The arrow indicates the relative direction of motion of the occulted stars. Each line segment represents an exposure. The transparency of each segment corresponds to the normalised flux measured in that integration (see Fig. 2). No transparency means flux = 1. In the light curve from Çanakkale there is an abnormally long dead-time between two concrete exposures. The corresponding blank space should not be mistaken for a flux drop. The pair of black ellipses represent the nominal solution reported in Table 3 for sampling 1. Their centres are marked with red dots. The thinner grey ellipses represent the 1σ uncertainties from the extracted sample, as discussed in the text.

In the text

Fig. 5 Histogram comparison of the posterior samples in Deq. The priors are plotted for reference (dashed lines). Sampling 1 is our preferred sampling, with uniform priors in all our model parameters. However, as Deq is not directly a model parameter, the prior in this value is not uniform. In sampling 2, the prior in Deq is equivalent to the value reported in Santos-Sanz et al. (2012) from thermal data. This prior pdf peaks at a higher Deq value than the prior pdf of sampling 1. Consequently, the posterior pdf in sampling 2 also peaks at a higher value than in sampling 1. However, both distributions yield compatible results within 1σ (68% credible intervals). Moreover, the result from sampling 1 is also compatible with the result from thermal data (sampling 2 prior).

In the text

Fig. 6 Pairwise marginal plot for the sample of our thermal binaryobject model MCMC run. The marginalised distribution in Deq is the same as the posterior from the occultation (sampling 1 in Fig. 5). We show it in red to highlight that this parameter is not being updated in this sampling. The regions depicted in the 2D histograms encircle the 39.4% (1σ) and 78.8% (2σ) of the samples. The coloured lines represent quantiles of the marginalised deq distribution. These are interpreted as upper limits for the diameter of the satellite.

In the text

	Fig. A.1 Histogram of the available sample of projected axis ratios in published TNO occultations together with our fitted truncated normal. The latter is the used prior distribution in b/a in sampling 2.
In the text

	Fig. A.2 Histogram of the available sample of beaming parameters from TNOs are Cool (Centaurs excluded; Lellouch et al. 2013) together with our Gaussian kernel density estimated pdf. The latter is the used prior distribution in η in the sampling in Sect. 5.
In the text

	Fig. B.1 Histogram of the relative area difference (|Athermal − Aocc|/Aocc) for our random sample of ellipsoids and rotation poles. The epoch for Athermal is 2010 October 17 and the rotational phase is random. The epoch for Aocc is 2022 August 22 and the rotational phase is such that a projected axis ratio compatible with the occultation is achieved.
In the text

	Fig. D.1 Pairwise marginal plot for our two MCMC samples of the posterior density. There are no signs of other modes representing alternative solutions. In the bidimensional plots, the represented region encircles 39.4% of the samples (representing the 1σ region). The red shading in these regions is a bidimensional histogram for our preferred sampling 1.
In the text