Issue 
A&A
Volume 643, November 2020



Article Number  A18  
Number of page(s)  17  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/202037447  
Published online  27 October 2020 
Constraining the orientation of the spin axes of extrasolar minor bodies 1I/2017 U1 (‘Oumuamua) and 2I/Borisov
^{1}
Universidad Complutense de Madrid, Ciudad Universitaria,
28040
Madrid,
Spain
email: nbplanet@ucm.es
^{2}
AEGORA Research Group, Facultad de Ciencias Matemáticas, Universidad Complutense de Madrid, Ciudad Universitaria,
28040
Madrid,
Spain
Received:
5
January
2020
Accepted:
16
September
2020
Context. The orientation of the spin axis of a comet is defined by the values of its equatorial obliquity and its cometocentric longitude of the Sun at perihelion. These parameters can be computed from the components of the nongravitational force caused by outgassing if the cometary activity is well characterized. The trajectories of known interstellar bodies passing through the Solar System show nongravitational accelerations.
Aims. The spinaxis orientation of 1I/2017 U1 (‘Oumuamua) remains to be determined; for 2I/Borisov, the already released results are mutually exclusive. In both cases, the values of the components of the nongravitational force are relatively well constrained. Here, we investigate – within the framework of the forced precession model of a nonspherical cometary nucleus – the orientation of the spin axes of ‘Oumuamua and 2I/Borisov using public orbit determinations that consider outgassing.
Methods. We applied a Monte Carlo simulation using the covariance matrix method together with Monte Carlo random search techniques to compute the distributions of equatorial obliquities and cometocentric longitudes of the Sun at perihelion of ‘Oumuamua and 2I/Borisov from the values of the nongravitational parameters.
Results. We find that the equatorial obliquity of ‘Oumuamua could be about 93°, if it has a very prolate (fusiform) shape, or close to 16°, if it is very oblate (disklike). Different orbit determinations of 2I/Borisov gave obliquity values of 59° and 90°. The distributions of cometocentric longitudes were in general multimodal.
Conclusions. Our calculations suggest that the most probable spinaxis direction of ‘Oumuamua in equatorial coordinates is (280°, +46°) if very prolate or (312°, −50°) if very oblate. Our analysis favors a prolate shape. For the orbit determinations of 2I/Borisov used here, we find most probable poles pointing near (275°, +65°) and (231°, +30°), respectively. Although our analysis favors an oblate shape for 2I/Borisov, a prolate one cannot be ruled out.
Key words: methods: data analysis / methods: numerical / celestial mechanics / comets: individual: 1I/2017 U1 (‘Oumuamua) / comets: general / comets: individual: 2I/Borisov
© ESO 2020
1 Introduction
Outgassing causes nongravitational accelerations that affect both the trajectories and rotational properties of minor bodies (see for example Sekanina 1981). The spinaxis orientation of a comet is defined by two angles (see for example Fig. 1 in Yeomans et al. 2004), the equatorial obliquity (ɛ_{c}) and the cometocentric longitude or argument of the subsolar meridian at perihelion (ϕ_{c}). The subsolar meridian is the great circle intersecting the object’s rotation axis and the subsolar point, where the Sun is directly overhead as seen from the minor body.
The value of ɛ_{c} gives the angle between the orbital and equatorial planes of the comet. Retrograde comets have ɛ_{c} ∈ (90°, 180°); if the direction of rotation is consistent with that of its orbital motion, ɛ_{c} < 90° and the comet is prograde (Sekanina 1981). The value of ϕ_{c} is measured from the vernal equinox of the comet (projection of the ascending node of the orbital plane of the comet on its equator) to the subsolar meridian at perihelion in the direction of increasing true anomaly, f. Therefore, it gives the hour angle of the Sun as seen from the comet at perihelion. When ϕ_{c} ∈ (180°, 360°), the sunlit pole is the southern one; if ϕ_{c} ∈ (0°, 180°), the northern pole faces the Sun at perihelion (Sekanina 1981).
In addition to the spinaxis orientation defined by (ɛ_{c}, ϕ_{c}), the thermal lag angle (η_{c}) measures how the region that dominates outgassing is shifted behind the subsolar meridian (Yeomans et al. 2004), and its value depends on the thermophysical properties of the sublimation phenomena produced in the comet and its rotational period (Sekanina 1981). Maximum insolation takes place at comet’s noon. If η_{c} ∈ (−90°, 0°), sublimation takes place mainly when the Sun is rising as seen from the surface of the comet; when the peak of activity is shifted toward the afternoon and evening hours, as seen from the comet, η_{c} ∈ (0°, 90°). The interval η_{c} ∈ (−90°, 90°) corresponds to the comet’s dayside; outside this range, we have the comet’s nightside. Sublimation of volatile materials, such as H_{2}O or CO, generates nongravitational accelerations that can be detected by analyzing astrometric data.
The nongravitational force is the result of directional mass loss and, therefore, its detectability decreases as the rotation rate increases; the orbital evolution of fast cometary nuclei with rotation periods of a few hours is nearly unaffected by outgassing (see for example Samarasinha et al. 2004). Marsden et al. (1973) discussed a mathematical formalism that developed the model originally proposed by Whipple (1950) to quantify the nongravitational effects on the orbital evolution of comets. Vaporization of volatiles produces a force with radial, transverse, and normal components of the form A g(r) with , where r is the heliocentric distance of the comet (Marsden et al. 1973). For subsurface water ice, r_{0} = 2.808 AU, k = 4.6142, m = 2.15, n = 5.093, and α = 0.1112620426 (Marsden et al. 1973); for carbon monoxide ice, r_{0} = 5.0 AU, k = 2.6, m = 2.0, n = 3.0, and α = 0.0408373333128795 (Yabushita 1996; Micheli et al. 2018). The values of the nongravitational parameters – A_{1} (radial), A_{2} (transverse), and A_{3} (normal) – depend on those of the angles ɛ_{c}, ϕ_{c}, and η_{c}.
Assuming a spherically symmetric nucleus and applying a linear precession model, the nongravitational parameters can be computed from astrometric observations of the comet, and ɛ_{c}, ϕ_{c}, and η_{c} can be obtained from the components of the nongravitational force by using an iterative least squares algorithm, as shown by Sitarski (1990). However, many observed comets and active asteroids are known to be nonspherical; in particular, interstellar object 1I/2017 U1 (‘Oumuamua) could be a very elongated (see for example Knight et al. 2017; Meech et al. 2017; Bolin et al. 2018; Drahus et al. 2018; Fraser et al. 2018) or flattened (Mashchenko 2019) body.
Sekanina (1984) developed the forced precession model of a nonspherical cometary nucleus to reproduce the effects of the torque that anisotropic outgassing exerts. This model includes additional parameters, and the direction of the spin pole (ɛ_{c}, ϕ_{c}) now depends on the oblateness of the nucleus: p = 1 − R_{b}∕R_{a}, where R_{a} is the equatorial radius of the object and R_{b} is its polar radius. When the nonspherical body rotates on its shorter axis, p > 0, we speak of an oblate nucleus; if the rotation of the nucleus is along its longer axis, p < 0, we have a prolate one.
In this work, we use the formulae of the forced precession model of a nonspherical nucleus described by Krolikowska et al. (1998). This model assumes that the nongravitational effects due to outgassing come from one dominant active area located on the surface of an ellipsoidal object of oblateness p that is spinning about a precessing axis at a given epoch with direction defined by (ɛ_{c}, ϕ_{c}); the outgassing area reaches its maximum activity at an angle η_{c} from the subsolar meridian as seen from the surface of the object.
The relevant part of the set of Eqs. (3) in Krolikowska et al. (1998) is: (1)
where S = p (2 − p), S_{1} =S (2 − S), Λ = ϕ_{c} + f, and sin Ψ = sinɛ_{c} sinΛ. For a spherical nucleus, p = 0, the set of Eqs. (1) collapses into the one shown on page 140 of Yeomans et al. (2004) and originally derived by Sekanina (1981). In both cases, the direction cosines are also given by the expressions C_{1} = A_{1}∕A, C_{2} = A_{2}∕A, and C_{3} = A_{3}∕A, where A = . The nongravitational parameters (A_{1}, A_{2}, and A_{3}) can be computed from astrometric observations of the comet, and ɛ_{c}, ϕ_{c}, η_{c}, and p can be obtained from the components of the nongravitational force by using an iterative least squares algorithm, as shown in Krolikowska et al. (1998).
The orbit determinations of extrasolar minor bodies ‘Oumuamua and 2I/Borisov (see Tables 1–3, and 5) include the usual orbital elements – eccentricity, e, perihelion distance, q, inclination, i, longitude of the ascending node, Ω, argument of perihelion, ω, and time of perihelionpassage, τ – and the nongravitational parameters (see above). If an orbit determination is meant to reproduce actual observations accurately and to make reliable ephemeris predictions into the past and the future, one has to consider how the orbital parameters (nongravitational ones included) affect one another using the covariance matrix whose elements indicate the level to which two given parameters vary together. In this context, the mean values of the parameters and the covariance matrix define a hyperellipsoid in multidimensional space. Such data can be used to derive distributions of the angles (ɛ_{c}, ϕ_{c}, and η_{c}), and the oblateness parameter, p, compatible with the observations via Monte Carlo random search starting with input data generated as described by de la Fuente Marcos & de la Fuente Marcos (2015).
Here, we outline a procedure to derive the spinaxis orientation and thermal lag angle of a comet from its nongravitational orbit determination that includes the mutual uncertainties between the various parameters. This approach is applied to study the orientation of the spin axes of extrasolar minor bodies ‘Oumuamua and 2I/Borisov. This paper is organized as follows. In Sect. 2, we present the data used in our analyses, discussing and validating our techniques. In Sect. 3, we applyour methodology to ‘Oumuamua’s data. The spinaxis orientation of 2I/Borisov is studied in Sect. 4. Our results are placed within the context of previous work in Sect. 5 and our conclusions are summarized in Sect. 6.
2 Data description and methods
In this work, we find two types of assumptions that affect the results obtained. On the one hand, we have those linked to the model used to constrain the three angles (ɛ_{c}, ϕ_{c}, and η_{c}), and the value of the oblateness; on the other, we find those associated with the public orbit determination, that has not been computed by us, such as the one species dominating the outgassing.
The assumptions behind the forced precession model of a nonspherical cometary nucleus are described in detail in Sekanina (1984) and Krolikowska et al. (1998). The ones affecting the orbit determination are validated numerically by their ability to reproduce the available observations and predict the future positions of the object. If the orbit determination is capable of providing the actual position of the object within a reasonable level of uncertainty (often less than a few arcseconds), then one has to assume that the starting hypotheses (used to compute the orbit determination) are sufficiently justified. However, it is certainly possible that some starting hypotheses could be formally inconsistent with evidence from other sources and still produce numerically correct results, in the sense of being consistent with the observations.
2.1 Data sources
In this paper, we use publicly available data from Jet Propulsion Laboratory’s (JPL) SmallBody Database (SBDB, Giorgini 2015)^{1} to investigate the spinaxis orientations of 1I/2017 U1 (‘Oumuamua) and 2I/Borisov. The data include the orbital elements and nongravitational parameters shown in Tables 1–3, and their associated covariance matrices. The data were used to construct the distributions of the angles ɛ_{c}, ϕ_{c}, and η_{c}, and in some cases that of the oblateness. Additional data have been retrieved from JPL’s SBDB or the Minor Planet Center (MPC, Rudenko 2016) using the tools provided by the Python package Astroquery (Ginsburg et al. 2019) and the Python module sbpy (Mommert et al. 2019).
JPL’s SBDB explicitly states that the orbit determination of ‘Oumuamua in Table 1 is identical to solution 7c in Micheli et al. (2018) that assumes a nongravitational acceleration that uses normalization factors appropriate for CO driven sublimation (see above). However, we have to acknowledge that outgassing of carbon based molecules, specifically CO and CO_{2} has been definitively ruled out as accelerants by NASA’s Spitzer Space Telescope nondetection of ‘Oumuamua (Trilling et al. 2018). The authors of this analysis suggest that outgassing of another gas species, perhaps H_{2}O, may be responsible for the observed nongravitational acceleration; however, the producer of the orbit determination in Table 1, D. Farnocchia, is listed as a coauthor of both Micheli et al. (2018) and Trilling et al. (2018). In their Table 1, ‘Oumuamua ISSI Team et al. (2019) show that the upper limits for CO from different sources are clearly orders of magnitude lower than those of other species. In any case, this inconsistency may not have any major impact on our results.
The primary objective of this work is not in contesting the orbit determinations computed by others, but in exploring their associated effects on the possible rotational properties of the objects. As we use the numerical mean values and their associated covariance matrices, and they are able to reproduce observational results (the astrometry) within reasonable accuracy levels, we consider their numerical values as essentially valid even if the original hypotheses used to compute them could be invalid. The subject of the actual driver of the detected outgassing remains controversial in the case of ‘Oumuamua(see for example Flekkøy et al. 2019; Hui & Knight 2019; Seligman & Laughlin 2020; Hoang & Loeb 2020).
As for the orbit determinations of 2I/Borisov in Tables 2 and 3, JPL’s SBDB states that the standard model – in other words, assuming a nongravitational acceleration that uses normalization factors appropriate for H_{2}O driven sublimation (see above) – has been applied to compute the nongravitational parameters. A water ice driver is supported by Sekanina (2019). However, it has already been confirmed that this comet exhibits an unusually high CO abundance (Cordiner et al. 2020). Although the orbital elements in Tables 2 and 3 are, in statistical terms, mutually consistent, we cannot say the same about the nongravitational parameters. The orbit determination in Table 2 – for which the last observation used was acquired on January 6, 2020 – shows a statistically significant (above the 6.8σ level) value of the nongravitational normal acceleration parameter (but the others are statistically compatible with zero). This is often the case when cometary outgassing comes from a single active region located at arotation pole and when the obliquity of the orbit plane to the equatorial plane is near 50° or 130° (see Sect. 3.1 in Yeomans et al. 2004).
In sharp contrast, the orbit determination in Table 3 has no value of the normal nongravitational parameter; therefore, we have to assume that it is compatible with zero. In addition, the value of A_{2} is consistent with that in Table 2 but now the dominant nongravitational acceleration is radial (statistically significant above the 17σ level), as it is in the case of ‘Oumuamua (see Table 1). The orbit determination in Table 3 includes observations acquired in January and February. Assuming that both orbit determinations are correct, they may represent two snapshots in the dynamical life of 2I/Borisov that may reflect the rapidly evolving nature of the outgassing processes taken place at the surface of the comet. In fact, there is robust observational evidence to support this interpretation: a sequence of outbursts was observed by multiple observers in FebruaryMarch (Bolin et al. 2020a; Drahus et al. 2020; Jewitt et al. 2020b,c; Zhang et al. 2020).
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar object 1I/2017 U1 (‘Oumuamua).
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar comet 2I/Borisov (I).
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar comet 2I/Borisov (II).
2.2 Methodology
The first step of the approach used in this paper is an implementation of the Monte Carlo using the Covariance Matrix (MCCM, Bordovitsyna et al. 2001; Avdyushev & Banshchikova 2007) methodology in which a Monte Carlo process generates control or clone orbits based on the nominal orbit, but adding random noise on each orbital parameter by making use of the covariance matrix. Considering a covariance matrix as computed by JPL’s Solar System Dynamics Group, Horizons OnLine Ephemeris System, the vector including the mean values of the orbital parameters at a given epoch t_{0} is of the form v = (e, q, τ, Ω, ω, i, A_{1}, A_{2}, A_{3}). If C is the covariance matrix at the same epoch associated with the nominal orbital solution that is symmetric and positivesemidefinite, then C = A A^{T} where A is a lower triangular matrix with real and positive diagonal elements, A^{T} is the transpose of A. In the case studied here, these matrices were 9 × 9 and the expressions of the individual elements a_{ij} of matrix A are shown in Appendix A.
If r is a vector made of univariate Gaussian random numbers – components r_{i} with i = 1, 9 generated using the BoxMuller method (Box & Muller 1958) – the sets of initial orbital elements of the control orbits can be computed using the expressions (assuming the structure provided by JPL’s Horizons OnLine Ephemeris System pointed out above), v_{c} =v + A r: (2)
This approach and set of equations have already been used to study the dynamics of 2I/Borisov as presented in de León et al. (2020).
In order to reconstruct the distributions of ɛ_{c}, ϕ_{c}, and η_{c}, and the oblateness, we generated 10^{5} instances of a given orbit (unless stated otherwise) using the set of Eqs. (2). For each one of them, we estimated the associated values of the angles (and p, if not imposed by additional evidence) by applying a Monte Carlo random search algorithm. The random search uses sets of values of the three angles (and p in relevant cases) that are generated randomly (uniformly) within the relevant intervals (see Sect. 1) until a stopping criterion is met. For example, if we consider the direction cosine for the radial component of the force (see equations in Sect. 1): (3)
that is also given by the expression C_{1} = A_{1}∕A, where A = . We can construct C_{1c} using data from MCCM and then try values of ɛ_{c}, ϕ_{c}, η_{c}, and p (to compute S and S_{1}) on Eq. (3) until C_{1c} − C_{1} < Δ, where Δ is a certain tolerance value that in our calculations and for practical reasons was set to 0.0001 (our results do not depend significantly on the tolerance value). A set of ɛ_{c}, ϕ_{c}, η_{c}, and p is regarded as a valid solution if C_{ic} − C_{i} < Δ, with i = 1, 3 (the radial, transverse, and normal direction cosines). We applied this technique and stopping criterion to obtain 10^{5} sets of angles and p (one per orbital instance) that we used to generate the distributions of angles and oblateness. This approach is robust when working with multimodal distributions that may signal the existence of multiple solutions. An iterative procedure like the one discussed by Sitarski (1990) may not be able to obtain all the statistically viable solutions.
Our technique delivers the distribution of spinaxis orientations defined by (ɛ_{c}, ϕ_{c}), but other alternative approaches discussed in the literature provide the ecliptic coordinates, (λ, β), of the spinaxis orientation. We can explore the associated distribution of pole directions by taking into account the dot product of the unit vector in the direction of the rotation pole – u_{r} = (cosβ_{p} cosλ_{p}, cosβ_{p} sinλ_{p}, sinβ_{p}) – with that of the orbital pole – u_{o} = (sini_{c} sinΩ_{c}, −sini_{c} cosΩ_{c}, cosi_{c}) – that gives u_{r} ⋅u_{o} = cosɛ_{c}. Therefore, the expression that can be used to apply a random search technique with a stopping criterion to find (λ_{p}, β_{p}) is: cos β_{p} cosλ_{p} sini_{c} sinΩ_{c} − cosβ_{p} sinλ_{p} sini_{c} cosΩ_{c} + sinβ_{p} cosi_{c} = cosɛ_{c}. From here, we can compute the equivalent distribution in equatorial coordinates, (α, δ), using the conversion tools provided by Astropy (Astropy Collaboration 2013, 2018).
In order to analyze the results, we produced histograms using the Matplotlib library (Hunter 2007) with two different sets of bins computed using Astropy (Astropy Collaboration 2013, 2018): one by applying the Freedman and Diaconis rule (Freedman & Diaconis 1981) and a second one using the Bayesian Blocks technique (Scargle et al. 2013). While the Freedman and Diaconis rule produces a constant bin size, the Bayesian Blocks technique allows varying bin widths that work well for complicated distributions like the ones studied here. We consider two different statistically motivated histograms to show that our conclusions do not depend on the choice of bin width. On the other hand and instead of using frequencybased histograms, we considered counts to form a probability density so the area under the histogram will sum to one. In addition, we construct diagrams to visualize better the statistical significance of the spinaxis orientations – as (ɛ_{c}, ϕ_{c}), (λ_{p}, β_{p}), or (α_{p}, δ_{p}) maps – representing kerneldensity estimates using Gaussian kernels (see for example Scott 1992) implemented by SciPy (Virtanen et al. 2020).
The spinaxis orientations obtained using the methodology described above are referred to the same epoch as the associated orbit determination. We have not studied the time evolution of the rotational properties of minor bodies as described by, for example, Krolikowska et al. (1998). However, for an object following a hyperbolic path and in absence of outbursts, fragmentation events, or close encounters with massive planets, the overall spinaxis orientation is expected to remain relatively stable as it is often the case for each apparition (not for multiple consecutive apparitions)of wellstudied comets, see Sects. 3.5–3.8 in Yeomans et al. (2004). In this context, we may interpret the results obtained for a given epoch as the most probable values of the rotational properties around that epoch.
Fig. 1 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 51P/Harrington. Top left panel: histogram of equatorial obliquity, ɛ_{51P/Harrington}, with bins computed using the Freedman and Diaconis rule, while top right panel: uses the Bayesian Blocks technique. Similarly, the cometocentric longitude, ϕ_{51P/Harrington}, histograms displayed in the lower panels also use the Freedman and Diaconis rule (left) and Bayesian Blocks (right). Histograms are based on data from Table 4 (see the text for details). 
2.3 Validation: comet 51P/Harrington
The nongravitational motion of comet 51P/Harrington has been studied by Sitarski (1996). This comet has experienced multiple fragmentation events since its discovery on August 14, 1953 (Sekanina 1997; Kidger & Manteca 2002). The work by Sitarski (1996) is referred to the epoch February 12, 1995 and the dataarc spans from discovery time to January 28, 1995. The values of the relevant properties are ɛ_{51P/Harrington} = 77. °6 ± 1. °1, ϕ_{51P/Harrington} = 92° ± 4°, η_{51P/Harrington} = 18. °4 ± 1. °4, and p_{51P/Harrington} = −0.32 ± 0.08 (Sitarski 1996), and they have been computed using Sekanina’s forced precession model of a nonspherical cometary nucleus. Sitarski (1996) considered asymmetric cometary activity with respect to perihelion in his calculations and A = 7.9192 × 10^{−9} ± 1.987 × 10^{−10} AU d^{−2}, but the values of the individual nongravitational acceleration parameters were not provided.
Table 4 shows the latest orbit determination available for comet 51P/Harrington. JPL’s SBDB states that the orbit determination has been computed using the standard model (see above). When comparing the data in Table 4 with those in Tables 1–3, we realize that the statistical relevance of the nongravitational parameters is substantially higher in the case of comet 51P/Harrington. It is also important to emphasize that 51P/Harrington is a periodic comet that has been studied for multiple decades (its orbital period is 7.16 yr); interstellar objects can only be studied once and they have comparatively short visibility windows.
Applying the procedure outlined in Sect. 2 using as input data those in Table 4 and the associated covariance matrix, we generated 10^{4} instances of a given orbit using the set of Eqs. (2). The results for the spinaxis orientation (ɛ_{51P/Harrington}, ϕ_{51P/Harrington}) are shown in Fig. 1. The top panels display the distribution of the equatorial obliquity that is unimodal with median and 16th and 84th percentiles of ; the bottom panels show a bimodal distribution for the cometocentric longitude with average maxima at about 55° and 235° (therefore, separated by 180°).
These results indicate that the effective equatorial plane of 51P/Harrington during its most recent perihelion was nearly perpendicular to its orbital plane; on the other hand, if ϕ_{51P/Harrington} ~ 55° the northern pole of the comet was facing the Sun at perihelion, or the southern one if ϕ_{51P/Harrington} ~ 235°. Our result for the value of the equatorial obliquity is fully consistent with that in Sitarski (1996), in particular with the trend shown in his Fig. 1. The value ofϕ_{51P/Harrington} in Fig. 1, bottom panels, is however inconsistent with the one in Sitarski (1996) and the trend in his Fig. 2; the median and 16th and 84th percentiles of ϕ_{51P/Harrington} are . It can be argued that the differences in ϕ_{51P/Harrington} could be theresult of the fragmentation events pointed out above.
Figure 2 shows an asymmetric distribution for η_{51P/Harrington}, with no clear maximum. The median and 16th and 84th percentiles of the distribution are . This impliesthat the maximum outgassing was taking place when the Sun was rising as seen from the surface of the comet. Again, this result is inconsistent with the one in Sitarski (1996), but it might be a byproduct of the fragmentation events pointed out above if the fresh surface has a different value of the thermal inertia.
Figure 3 shows a bimodal distribution for the oblateness of the nucleus of 51P/Harrington. The median and 16th and 84th percentiles of the distribution are , with maxima at about 0.06 and −0.6. Our results favor an oblate nucleus although a prolate one cannot be discarded. Considering the large dispersion, our median value is statistically consistent with the results in Sitarski (1996) and our secondary maximum is marginally consistent with his result. Again, the fact that the surface of the comet has experienced multiple alterations since the data used by Sitarski (1996) were acquired and that some of these episodes may have even altered the shape of 51P/Harrington make it difficult to validate our results against those in Sitarski (1996). This is particularly obvious when we consider that the value of A in Table 4 is nearly five times larger than the one computed by Sitarski (1996). In addition, the model used to derive the values of the parameters in Table 4 considered symmetric – not asymmetric as in Sitarski (1996) – cometary activity with respect to perihelion.
Taking into account the dot product of the unit vector in the direction of the rotation pole with that of the orbital pole and its relationship with ɛ_{51P/Harrington} as pointed out above, we have computed the distributions of (λ_{p}, β_{p}) and (α_{p}, δ_{p}). Figure 4 shows the resulting distributions in ecliptic coordinates and Fig. 5 shows the maps of statistical significance computed as described in the second to last paragraph of Sect. 2.2. The most probable spinaxis direction of 51P/Harrington in equatorial coordinates could be (260°, −13°) or (120°, +20°). On the other hand, the equatorial coordinates of the Sun when the comet reached perihelion on August 12, 2015 were approximately (141. °8, +15. °0), somewhat consistent with the orientation of one of the poles.
Heliocentric orbital elements and 1σ uncertainties of comet 51P/Harrington.
Fig. 2 Distribution of thermal lag angle for 51P/Harrington. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 4. 
Fig. 3 Distribution of values of the oblateness of the nucleus for 51P/Harrington. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 4. 
Fig. 4 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 51P/Harrington. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 4. 
3 Interstellar object 1I/2017 U1 (‘Oumuamua)
We applied the procedure outlined in Sect. 2 to produce the distributions of ɛ_{c}, ϕ_{c}, η_{c}, and p_{c} for interstellar object 1I/2017 U1 (‘Oumuamua) considering the orbit determination in Table 1 that is discussed by Micheli et al. (2018). The model favored by these authors assumes that ‘Oumuamua’s cometlike outgassing is driven by CO, not H_{2}O (see Sect. 1 and the discussion in Sect. 2). The results for the spinaxis orientation (ɛ_{‘Oumuamua}, ϕ_{‘Oumuamua}) are displayed in Fig. 6, which summarizes the results of 10^{5} Monte Carlo experiments.
Figure 6, top panels, shows the distribution of the equatorial obliquity. The median and 16th and 84th percentiles of the distribution are . This result indicates that the effective equatorial plane of ‘Oumuamua was nearly perpendicular to its orbital plane. Figure 6, bottom panels, shows a bimodal distribution with average maxima at about 140° and 320° (therefore, separated by 180°). If ϕ_{‘Oumuamua} ~ 140° the northern pole of ‘Oumuamua was facing the Sun at perihelion, if ϕ_{‘Oumuamua} ~ 320° was the southern one. In general, our methodology cannot always distinguish which pole was facing the Sun when ‘Oumuamua reached perihelion.
Figure 7 also shows a unimodal distribution for η_{‘Oumuamua}. The median and 16th and 84th percentiles of the distribution are . This impliesthat the maximum outgassing was taking place mainly at ‘Oumuamua’s noon, when the Sun was at its highest as seen from the surface of the interstellar visitor.
Figure 8 shows a unimodal distribution for the oblateness of ‘Oumuamua with median and 16th and 84th percentiles of . This result slightly favors a fusiform shape for ‘Oumuamua over a disklike one. The actual shape of ‘Oumuamua remains a controversial subject (see for example the discussion in Vazan & Sari 2020). Although multiple studies have pointed out that the available observational data strongly favor that ‘Oumuamua has a very elongated shape (see for example Knight et al. 2017; Bolin et al. 2018; Drahus et al. 2018; Fraser et al. 2018), Belton et al. (2018) stated that its shape might be fusiform, but a disklike appearance could not be ruled out.
Mashchenko (2019) presented a very detailed modeling of the light curve of ‘Oumuamua to show that his bestfitting model, with a probability of 91%, was a thin disk with p = 0.836; his second choice, with a probability of 16%, was a thin spindle with p = −6.7. The thindisk model was significantly more successful at reproducing the observed variations in apparent magnitude over time. However, the results presented by Mashchenko (2019) have been contested by Zhang & Lin (2020) by arguing that multiple lines of evidence favor the very prolate or cigarlike shape against the very oblate or disklike shape. In particular, they pointed out that the fusiform case is energetically more stable. In the following, we test the impact of both models on the values of the rotational parameters of ‘Oumuamua. We first assume that p = −6.7 together with the data in Table 1, then we consider p = 0.836.
Figures 9 and 10, show the distributions of ɛ_{‘Oumuamua}, ϕ_{‘Oumuamua} and η_{‘Oumuamua} under the assumption of p = −6.7; we generated 10^{4} instances of a given orbit. The results for the orientation of the pole are similar to those in Fig. 6: the equatorial obliquity becomes 93° ± 3° and ϕ_{‘Oumuamua} ~ 140° or ~ 320°. However, Fig. 10 shows a relatively flat distribution for η_{‘Oumuamua}; in other words, outgassing takes place evenly throughout a rotation period. The actual location of the pole is studied in Figs. 11 and 12: under the thinspindle scenario, the most probable spinaxis direction of ‘Oumuamua in equatorial coordinates could be (280°, +46°), with approximate Galacticcoordinates l = 75°, b = +21°, that point slightly above the Galactic disk, toward the constellation of Lyra.
Figures 13 and 14, show the distributions of ɛ_{‘Oumuamua}, ϕ_{‘Oumuamua}, and η_{‘Oumuamua} for 10^{4} instances of the orbit, under the assumption of p = 0.836. The results for the orientation of the pole are very different from those in Fig. 6: the equatorial obliquity is now and ϕ_{‘Oumuamua} ~ 67° or ~247°. However, the distribution of η_{‘Oumuamua} in Fig. 14 resembles that in Fig. 7; maximum outgassing would take place at ‘Oumuamua’s noon with η_{‘Oumuamua} = 3° ± 2°. The location of the pole under the thindisk scenario is shown in Figs. 15 and 16; the most probable spinaxis direction of ‘Oumuamua in equatorial coordinates could be (312°, −50°). This is consistent with the estimate of the orbital pole of ‘Oumuamua in de la Fuente Marcos & de la Fuente Marcos (2017a,b): compare Fig. 16, middle panel, with the middle panel in Fig. 1 of de la Fuente Marcos & de la Fuente Marcos (2017b), both in ecliptic coordinates. The approximate Galactic coordinates of the rotational pole are l = 349°, b = −39°, that point slightly below the Galactic bulge, toward the constellation Indus.
The rotational properties provided by our modeling are very different depending on how extreme the shape of ‘Oumuamua is; this suggests that the spinaxis orientation of ‘Oumuamua is tightly controlled by its shape. In the Solar System, most cometary nuclei are closer to a prolate configuration, not oblate (Samarasinha et al. 2004). Although the work by Mashchenko (2019) is widely regarded as the best available at the moment, the arguments put forward by Zhang & Lin (2020) may indicate that a final answer on ‘Oumuamua’s shape still remains elusive.
Fig. 5 Gaussian kernel density estimation of the spinaxis orientations computed using the SciPy library (Virtanen et al. 2020). Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of comet 51P/Harrington available from the MPC. 
Fig. 6 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 1I/2017 U1 (‘Oumuamua). Top left panel: histogram of equatorial obliquity, ɛ_{‘Oumuamua}, with bins computed using the Freedman and Diaconis rule, while top right panel: uses the Bayesian Blocks technique. Similarly, the cometocentric longitude, ϕ_{‘Oumuamua}, histograms displayed in the lower panels also use the Freedman and Diaconis rule (left) and Bayesian Blocks (right). Histograms are based on data from Table 1 (see the text for details). 
Fig. 7 Distribution of thermal lag angle for 1I/2017 U1 (‘Oumuamua). Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 
Fig. 8 Distribution of values of the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua). Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 
Fig. 9 Same as Fig. 6 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to −6.7 as suggestedby Mashchenko (2019). 
Fig. 10 Same as Fig. 7 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to −6.7 as suggestedby Mashchenko (2019). 
Fig. 11 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 1I/2017 U1 (‘Oumuamua) assuming a very prolate shape, with p = −6.7. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 
Fig. 12 Gaussian kernel density estimation of the spinaxis orientations computed using the SciPy library (Virtanen et al. 2020) for a very prolate case (p = −6.7). Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of interstellar object 1I/2017 U1 (‘Oumuamua) available from the MPC. 
4 Interstellar comet 2I/Borisov
It is unlikelythat additional data on 1I/2017 U1 (‘Oumuamua) may eventually emerge; in sharp contrast, interstellar comet 2I/Borisov is still being actively observed. The two most recent orbit determinations of 2I/Borisov are shown in Tables 2 and 3 and both consider nongravitational effects driven by sublimation of water ice (see Sect. 1). However, Ye et al. (2020) argued that a model based on sublimation of CO could be more consistent with the observations compared to a H_{2}O model. In addition, Bodewits et al. (2020) have pointed out that the coma of 2I/Borisov contains significantly more CO than H_{2}O gas, with abundances of at least 173%, more than three times higher than previously measured for any comet in the inner (<2.5 AU) Solar System. However, the water production rate prior to perihelion was found to increase faster than for most known dynamicallynew comets (Xing et al. 2020); a hyperactive nucleus was also favored by McKay et al. (2020). Again and as considered above, an orbit determination may produce numerically correct results (in the sense of correct ephemeris predictions) even if the hypotheses used to carry out the calculations are physically questionable or even invalid.
The orbit determination in Table 2 corresponds to an epoch nearly a month past perihelion and shows that only the normal acceleration may be playing a significant role on the dynamical evolution of 2I/Borisov as the other components have values statistically compatible with zero. Using data from Table 2 as input, weobtained the distributions in Figs. 17 and 18. The distribution of ɛ_{2I/Borisov} is unimodal but wide (see Fig. 17, top panels); there is one statistically significant maximum and the median and 16th and 84th percentiles are . Figure 17, bottom panels, shows a multimodal distribution for ϕ_{2I/Borisov}; the absolute maxima are found at about 90° and 270°, with secondary maxima at 40° and 220° (therefore, each pair separated by 180°), and the associated dispersions are relatively small. Figure 18 displays a multimodal distribution for η_{2I/Borisov}, with a dominant maximum at about −130°, which is on the nightside of the comet.
The most probable value of the equatorial obliquity of 2I/Borisov is sometimes found for periodic comets in the Solar System, when outgassing comes from a single active region located at a rotation pole (see Sect. 3.1 in Yeomans et al. 2004). Surprisingly, outgassing appears to be restricted to the nightside of the comet. Most comets tend to exhibit afternoon activity (Sekanina 1981). The lack of dayside activity seems problematic when considering that ‘Oumuamua’s activity was restricted to the dayside for the very oblate case; however, outbursts from the nightside of comets 9P/Tempel 1 (Farnham et al. 2007, 2013), 103P/Hartley 2 (A’Hearn et al. 2011; Bruck Syal et al. 2013), and 67P/ChuryumovGerasimenko (Knollenberg et al. 2016; Rinaldi et al. 2019) have been observed.
Figure 19 shows that the most probable value for the oblateness parameter for 2I/Borisov corresponds to that of a very oblate body with median and 16th and 84th percentiles of ; the peak of the distribution is found at p = 0.95. This is even more extreme than the value of p = 0.836 favored by Mashchenko (2019) for 1I/2017 U1 (‘Oumuamua).
The location of the pole is shown in Figs. 20 and 21. The most probable spinaxis direction of 2I/Borisov in equatorial coordinates could be (275°, +65°), see Fig. 21, bottom panel, with approximate Galactic coordinates l = 95°, b = +28° that point slightly above the Galactic disk, toward the constellation Draco.
Our results are quite different when considering the most recent orbit determination of 2I/Borisov in Table 3. The distribution of ɛ_{2I/Borisov} in Fig. 22, top panels, is still unimodal but now the median and 16th and 84th percentiles are 90° ± 52°. The distribution of ϕ_{2I/Borisov} in Fig. 22, bottom panels, is multimodal as the one in Fig. 17, bottom panels, but now the dominant maxima are found at about 100° and 280°. The distribution of η_{2I/Borisov} in Fig. 23 is very different from that in Fig. 18 and shows peaks at −20° and 20°, with median and 16th and 84th percentiles of 0° ± 29°; the maximum of activity takes place around noon as seen from the comet.
Figure 24 shows that the distribution of the values for the oblateness parameter of 2I/Borisov is bimodal with an absolute maximum of about 0.34 and a secondary peak at −0.36. Although our analysis favors an oblate shape for 2I/Borisov, a prolate one cannot be ruled out.
The location of the pole is shown in Figs. 25 and 26. The most probable spinaxis direction of 2I/Borisov in equatorial coordinates could be (231°, +30°), see Fig. 26, bottom panel, with approximate Galactic coordinates l = 47°, b = +56° that point well above the Galactic disk, toward the constellation Corona Borealis.
Fig. 13 Same as Fig. 6 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to 0.836 as suggested by Mashchenko (2019). 
Fig. 14 Same as Fig. 7 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to 0.836 as suggested by Mashchenko (2019). 
Fig. 15 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 1I/2017 U1 (‘Oumuamua) assuming a very oblate shape, with p = 0.836. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 
Fig. 16 Gaussian kernel density estimation of the spinaxis orientations computed using the SciPy library (Virtanen et al. 2020) for a very oblate case (p = 0.836). Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of interstellar object 1I/2017 U1 (‘Oumuamua) available from the MPC. 
Fig. 17 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 2I/Borisov. These results correspond to the orbit determination of comet 2I/Borisov in Table 2. 
Fig. 18 Distribution of thermal lag angle for 2I/Borisov. These results correspond to the orbit determination of comet 2I/Borisov in Table 2. 
Fig. 19 Distribution of the oblateness parameter for 2I/Borisov. These results correspond to the orbit determination of comet 2I/Borisov in Table 2. 
Fig. 20 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 2I/Borisov. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 2. 
5 Discussion
The main objective of the research presented here is not to argue against or in favor of the various orbit determinations (and their underlying outgassing models) used to generate input data for our analysis, but to point out the consequences and implications that a given nongravitational solution may have on certain comet properties such as the spinaxis orientation, the thermal lag angle, or the value of the oblateness. On the other hand, both 1I/2017 U1 (‘Oumuamua) and 2I/Borisov are singleapparition minor bodies and it is unclear how many of the trends observed forwellstudied periodic comets and discussed in the literature can be expected to be present for these extrasolar visitors. Sekanina (1981) argues that most periodic comets tend to have values of the obliquity close to 90° and the spinaxis tends to point toward the Sun at perihelion, which results in the most extreme insolation possible. In addition, the forced precession model of a nonspherical cometary nucleus provides better results when applied to data from several consecutive apparitions of the body under study (Krolikowska et al. 1998), but this is clearly not possible when considering interstellar objects.
Fig. 21 Gaussian kernel density estimation of the spinaxis orientations of 2I/Borisov computed using the SciPylibrary (Virtanen et al. 2020) and based on data from Table 2. Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of 2I/Borisov available from the MPC. Already published spinaxis directions are plotted as filled orange symbols with white error bars (when available): triangles (Ye et al. 2020), square (Manzini et al. 2020), star (Bolin et al. 2020b), and circle (Kim et al. 2020). 
5.1 1I/2017 U1 (‘Oumuamua)
Although noestimates of the orientation of ‘Oumuamua’s rotation pole have been made public yet, Fitzsimmons et al. (2018) assumed an obliquity of 0° in their analyses and Micheli et al. (2018) considered 45°, but additional modeling by Zhou (2019) indicated that ‘Oumuamua’s effective obliquity when it entered the Solar System may have been 93. °5 (this result has been removed from the final, published version of this work, Zhou 2020). The value assumed by Fitzsimmons et al. (2018) is marginally consistent with the one computed here for a very oblate shape (see Fig. 13, top panel); the value obtained by Zhou (2019) coincides with our estimate in the case of a very prolate shape (see Fig. 9, top panel).
Seligman et al. (2019) argued that ‘Oumuamua’s outgassing is made of H_{2}O molecules and that the venting of this material is expected to be collimated toward the Sun. We found a small lag angle between the subsolar meridian and the direction of maximum mass ejection for the case of a very oblate shape, but the very prolate case shows no statistically significant value of the lag angle, which can be interpreted as lack of correlation between level of activity and elevation of the Sun over the local horizon. On the other hand, the equatorial coordinates of the Sun when ‘Oumuamua reached perihelion on September 9, 2017 were approximately (167. °8, +5. °2), therefore none of our predictions regarding the spinaxis direction are close to the position favored (see above) by Sekanina (1981).
Rafikov (2018) showed that under standard conditions, if ‘Oumuamua is an elongated body experiencing outgassing, it should have followed a rapid and violent rotational evolution leading to rotational fission and probably catastrophic disruption prior to its discovery. The fact is that ‘Oumuamua was discovered and Seligman et al. (2019) explained that being longterm structurally stable, outgassing, and having an extremely oblong appearance could be simultaneously possible if the active region tracks the subsolar point (maximal venting occurs when both are aligned). The modeling presented by Mashchenko (2019) is consistent with the discussion in Seligman et al. (2019).
Our analysis of the rotational properties of ‘Oumuamua in Sect. 3 shows that, for the same orbit determination in Table 1, a very prolate or cigarlike shape with p = −6.7 leads to outgassing taking place evenly throughout the rotation period (see Fig. 10); in sharp contrast, a very oblate or disklike shape with p = 0.836 leads to a strong correlation between time of maximum activity or outgassing and alignment between active region and subsolar point (see Fig. 14) as demanded by Seligman et al. (2019). On the other hand, the Solar System is no strange to very oblate or flattened bodies, 486958 Arrokoth (2014 MU_{69}), a transNeptunian contact binary, has been found to have an average oblateness of p = 0.695, with 0.808 for the large lobe and 0.582 for the small lobe (Amarante & Winter 2020).
It can be argued that the values of the parameters obtained when no value of p is imposed (, ) lead to a trivial solution to the problem under study here. The effective equatorial plane of ‘Oumuamua was nearly perpendicular to its orbital plane and the maximum outgassing was taking place at ‘Oumuamua’s noon, when the Sun was at its highest as seen from the surface of the interstellar visitor. In this case and if ‘Oumuamua is fusiform, its body would be parallel to its orbital plane and the thermal lag could be effectively zero, then the outgassing is confined only to the radial direction and the two poles are indistinct. This is caused by the fact that both the transverse and normal components of the force are statistically compatible with zero. The inclusion of information on p leads to more complex solutions.
The orbit determination in Table 1 corresponds to 74 d past perihelion and shows that only the radial acceleration may have played a significant role on the dynamical evolution of ‘Oumuamua as it visited the inner Solar System because the other components have rather uncertain values, statistically compatible with zero (Micheli et al. 2018). ‘Oumuamuais in a tumbling rotational state (Belton et al. 2018; Drahus et al. 2018; Fraser et al. 2018). Belton et al. (2018) argued for two periods – a short one of 8.7 h and a long one (slow precession) of 54.5 h that are the result of tumbling – and a long, narrow rodlike shape. Drahus et al. (2018) found a variable light curve with a periodicity of 7.56 h consistent with a very elongated shape resulting from a catastrophic collision in the distant past. Extensive reviews on ‘Oumuamua’s properties can be found in Hainaut et al. (2018) and ‘Oumuamua ISSI Team et al. (2019).
Fig. 22 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 2I/Borisov. Same as Fig. 17 but for the orbit determination of comet 2I/Borisov in Table 3. 
Fig. 23 Distribution of thermal lag angle for 2I/Borisov. Same as Fig. 18 but for the orbit determination of comet 2I/Borisov in Table 3. 
Fig. 24 Distribution of the oblateness parameter for 2I/Borisov. Same as Fig. 19 but for the orbit determination of comet 2I/Borisov in Table 3. 
Fig. 25 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 2I/Borisov. Same as Fig. 20 but based on data from Table 3. 
5.2 2I/Borisov
Although no estimates on the spinaxis orientation of ‘Oumuamua have been presented in the literature, there are a number of released results for 2I/Borisov although they are all mutually exclusive. Ye et al. (2020) found two equally preferred locations (see Fig. 7 in their paper) for the pole’s right ascension and declination: (340°, +30°) and (205°, −55°). Kim et al. (2020) discussed a persistent asymmetry in the dust coma that could be best explained by a thermal lag on the rotating nucleus causing peak mass loss to take place during the comet nucleus afternoon; these authors calculated a value for the obliquity of 2I/Borisov of 30° and estimated a pole direction of (α, δ) = (205°, +52°). Manzini et al. (2020) computed a spinaxis direction (α, δ) = (260° ± 15°, −35° ± 10°). Bolin & Lisse (2020) estimated a value of the oblateness ≤−0.5 and argued that a jet is close to the rotation pole that points toward (α, δ) = (322° ± 10°, +37° ± 10°), mentioning that the rotation of the nucleus of comet 2I/Borisov seems to occur on a single axis and that it appears to be not chaotic.
Together with our own estimates in the form of maps of statistical significance, the various estimates for the spinaxis orientation of 2I/Borisov are plotted in Figs. 21 and 26, bottom panels, as filled orange symbols with white error bars (when available): triangles (Ye et al. 2020), square (Manzini et al. 2020), star (Bolin & Lisse 2020), and circle (Kim et al. 2020). The published solutions for the spinaxis orientation of 2I/Borisov are quite different and they may signal actual shortterm changes in the values of the rotational properties of this object or problems with the interpretation of some of the observed features. Our estimate based on Table 2 (see Fig. 21, bottom panel) is inconsistent with any of the other four, but the one based on the orbit determination in Table 3 (see Fig. 26, bottom panel) is marginally consistent with the result obtained by Kim et al. (2020) although the associated values of the obliquity are quite different, 30° versus 90°. Our value of the oblateness for the solution in Table 3 is inconsistent with the one computed by Bolin & Lisse (2020). On the other hand, the equatorial coordinates of the Sun when 2I/Borisov reached perihelion on December 8, 2019 were approximately (254. °6, 22. °7), so none of our predictions are close to the position favored by Sekanina (1981). This location is also inconsistent with the estimates provided by Ye et al. (2020), Bolin & Lisse (2020), and Kim et al. (2020); it is, however, somewhat consistent with the spinaxis direction computed by Manzini et al. (2020).
Most of the observed properties of interstellar comet 2I/Borisov have been found to be remarkably similar to those of known Solar System comets (see for example de León et al. 2019, 2020; Fitzsimmons et al. 2019; Opitom et al. 2019; Cremonese et al. 2020; Guzik et al. 2020; Jewitt et al. 2020a; Kareta et al. 2020; Yang et al. 2020). Recent estimates, placed the rotation period of 2I/Borisov at about 5.3 h or perhaps 10.6 h (Bolin & Lisse 2020). Gladman et al. (2019) had found a very slow rotation rate with a period of about 13 d. If confirmed, this would be the longest rotation period ever observed in a comet. Slowly rotating Solar System comets, such as 1P/Halley and 109P/SwiftTuttle, have rotation periods of about 2.8 d (see for example Table 1 in Samarasinha et al. 2004); comet 29P/SchwassmannWachmann 1 – well known for its enormous and random outbursts – was reported to have a rotation period of about 5 d (see for example Sekanina 1981), but recent determinations place its spin rate close to 1.4 d (see for example Table 1 in Samarasinha et al. 2004). JPL’s SBDB indicates that, among known comets, P/2006 HR30 (Siding Spring) has the longest rotation period with 2.9 d.
The thermal lag angle is relatively small for our results based on Table 3 and this suggests that outgassing for 2I/Borisov appears to have been active only when the Sun was rather high over the local horizon. As for the large lag angle obtained in the case of the results based on Table 2, it suggests that an insulating crust of nonvolatile material might have been present over most of the surface of this comet prior to the sequence of outbursts pointed out above; no significant change in the slightly reddish color of the comet was observed from late September 2019 through late January 2020 (Hui et al. 2020). Alternatively, this could be the result of its putative long period as the effects of the outgassing from the comet’s dayside may cancel out and only discrete nightside outbursts may result in a net contribution to the value of nongravitational normal acceleration parameter. Yet another scenario may involve an active vent close to the south pole of the comet that, due to the viewing geometry, remained on the nightside of the comet during the time interval covered by the observations. Both Manzini et al. (2020) and Bolin & Lisse (2020) argued that such a scenario may have been possible.
An improved orbit determination of 2I/Borisov was published on August 21, 2020 (see Table 5), when this work was under review, that also considers nongravitational effects driven by sublimation of water ice (see Sect. 1). This solution is close to that in Table 3 although it assummes asymmetric nongravitational forces including an additional parameter, dτ, to account for an asymmetry of the comet’s outgassing (see for example Yeomans et al. 2004). We have repeated the calculations discussed in Sect. 4 using the data in Table 5 as input and found results that are similar (see Fig. 27) to those obtained for the orbit determination in Table 3. The inclusion of the nongravitational perihelion offset has no effects on the mathematical expressions used to perform our calculations (see Sect. 1).
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar comet 2I/Borisov (III).
6 Conclusions
In this paper, we have shown that the most probable values of the equatorial obliquity and the cometocentric longitude of the Sun at perihelion as well as the value of the thermal lag angle and the oblateness of a comet can be computed from an orbit determination that includes nongravitational parameters using a Monte Carlo random search approach within the framework of the forced precession model of a nonspherical cometary nucleus. The algorithm receives input from a MCCM process thatincludes the uncertainties of the orbit determination in a selfconsistent way. The resulting combined algorithm is applied to investigate the spinaxis orientations of extrasolar minor bodies 1I/2017 U1 (‘Oumuamua) and 2I/Borisov.
For ‘Oumuamua and using data from Table 1, we found that its equatorial obliquity could be about 93°, if it has a very prolate (fusiform) shape, or close to 16°, if it is very oblate (disklike). Outgassing for ‘Oumuamua appears to have been active only when the Sun was highest over the local horizon, as the thermal lag angle is small, if the body has a thindisk appearance; for a cigarlike shape, outgassing may have been evenly distributed in time. Our calculations suggest that the most probable spinaxis direction of ‘Oumuamua in equatorial coordinates could be (280°, +46°) if very prolate or (312°, −50°) if very oblate. Our analysis favors a prolate shape. No previous estimates of the direction of the polar axis of ‘Oumuamua are available.
For 2I/Borisov, three different orbit determinations gave values of 59° (using data from Table 2) and 90° (using data from Tables 3 and 5) for its obliquity, with most probable poles pointing near (275°, +65°) and (231°, +30°), respectively. Although our analysis favors an oblate shape (p = 0.34) for 2I/Borisov, a prolate one cannot be ruled out. Our estimates for the spinaxis direction of 2I/Borisov are inconsistent with published determinations (except perhaps with the one in Kim et al. 2020), which already were mutually exclusive.
No further data may be collected from ‘Oumuamua, but new data on 2I/Borisov will be made public and they may lead to better, morereliable estimates of its rotational properties. The new data may improve what is already known and/or perhaps provide new snapshots of the evolution of its rotation.
Acknowledgements
We thank the referee for her/his constructive, detailed, and insightful report that included very helpful suggestions regarding the presentation of this paper and the interpretation of our results, J. de León, J. Licandro, and M. SerraRicart for comments on 2I/Borisov’s ongoing cometary activity, and A. I. Gómez de Castro for providing access to computing facilities. Part of the calculations and the data analysis were completed on the Brigit HPC server of the “Universidad Complutense de Madrid” (UCM), and we thank S. Cano Alsúa for his help during this stage. This research was partially supported by the Spanish “Ministerio de Economía y Competitividad” (MINECO) under grant ESP201787813R. In preparation of this paper, we made use of the NASA Astrophysics Data System, the ASTROPH eprint server, and the MPC data server. This research made use of Astropy^{2}, a communitydeveloped core Python package for Astronomy (Astropy Collaboration 2013, 2018).
Appendix A Elements of the matrix A
If the elements of the covariance matrix C are written as c_{ij} and those of A as a_{ij} (so C = A A^{T}), where those are the entries in the ith row and jth column, they are related by the following expressions:
References
 A’Hearn, M. F., Belton, M. J. S., Delamere, W. A., et al. 2011, Science, 332, 1396 [NASA ADS] [CrossRef] [Google Scholar]
 Amarante, A., & Winter, O. C. 2020, MNRAS, 496, 4154 [CrossRef] [Google Scholar]
 Astropy Collaboration (Robitaille, T. P., et al.) 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Astropy Collaboration (PriceWhelan, A. M., et al.) 2018, AJ, 156, 123 [Google Scholar]
 Avdyushev, V. A., & Banshchikova, M. A. 2007, Sol. Syst. Res., 41, 413 [NASA ADS] [CrossRef] [Google Scholar]
 Belton, M. J. S., Hainaut, O. R., Meech, K. J., et al. 2018, ApJ, 856, L21 [NASA ADS] [CrossRef] [Google Scholar]
 Bodewits, D., Noonan, J. W., Feldman, P. D., et al. 2020, Nat. Astron., 4, 867 [CrossRef] [Google Scholar]
 Bolin, B. T., & Lisse, C. M. 2020, MNRAS, 497, 4031 [CrossRef] [Google Scholar]
 Bolin, B. T., Weaver, H. A., Fernandez, Y. R., et al. 2018, ApJ, 852, L2 [CrossRef] [Google Scholar]
 Bolin, B. T., Bodewits, D., Lisse, C. M., et al. 2020a, ATel 13613, 1 [Google Scholar]
 Bolin, B. T., Lisse, C. M., Kasliwal, M. M., et al. 2020b, AJ, 160, 26 [CrossRef] [Google Scholar]
 Bordovitsyna, T., Avdyushev, V., & Chernitsov, A. 2001, Celest. Mech. Dyn. Astron., 80, 227 [NASA ADS] [CrossRef] [Google Scholar]
 Box, G. E. P., & Muller, M. E. 1958, Ann. Math. Stat., 29, 610 [CrossRef] [Google Scholar]
 Bruck Syal, M., Schultz, P. H., Sunshine, J. M., et al. 2013, Icarus, 222, 610 [NASA ADS] [CrossRef] [Google Scholar]
 Cordiner, M. A., Milam, S. N., Biver, N., et al. 2020, Nat. Astron., 4, 861 [CrossRef] [Google Scholar]
 Cremonese, G., Fulle, M., Cambianica, P., et al. 2020, ApJ, 893, L12 [CrossRef] [Google Scholar]
 de la Fuente Marcos, C., & de la Fuente Marcos, R. 2015, MNRAS, 453, 1288 [NASA ADS] [CrossRef] [Google Scholar]
 de la Fuente Marcos, C., & de la Fuente Marcos, R. 2017a, Res. Notes Amer. Astron. Soc., 1, 5 [CrossRef] [Google Scholar]
 de la Fuente Marcos, C. & de la Fuente Marcos, R. 2017b, Res. Notes Amer. Astron. Soc., 1, 9 [CrossRef] [Google Scholar]
 de León, J., Licandro, J., SerraRicart, M., et al. 2019, Res. Notes Amer. Astron. Soc., 3, 131 [NASA ADS] [CrossRef] [Google Scholar]
 de León, J., Licandro, J., de la Fuente Marcos, C., et al. 2020, MNRAS, 495, 2053 [CrossRef] [Google Scholar]
 Drahus, M., Guzik, P., Waniak, W., et al. 2018, Nat. Astron., 2, 407 [NASA ADS] [CrossRef] [Google Scholar]
 Drahus, M., Guzik, P., Udalski, A., et al. 2020, ATel 13549, 1 [Google Scholar]
 Farnham, T. L., Wellnitz, D. D., Hampton, D. L., et al. 2007, Icarus, 187, 26 [NASA ADS] [CrossRef] [Google Scholar]
 Farnham, T. L., Bodewits, D., Li, J.Y., et al. 2013, Icarus, 222, 540 [NASA ADS] [CrossRef] [Google Scholar]
 Fitzsimmons, A., Snodgrass, C., Rozitis, B., et al. 2018, Nat. Astron., 2, 133 [NASA ADS] [CrossRef] [Google Scholar]
 Fitzsimmons, A., Hainaut, O., Meech, K. J., et al. 2019, ApJ, 885, L9 [NASA ADS] [CrossRef] [Google Scholar]
 Flekkøy, E. G., Luu, J., & Toussaint, R. 2019, ApJ, 885, L41 [CrossRef] [Google Scholar]
 Fraser, W. C., Pravec, P., Fitzsimmons, A., et al. 2018, Nat. Astron., 2, 383 [NASA ADS] [CrossRef] [Google Scholar]
 Freedman, D., & Diaconis, P. 1981, Prob. Theory Rel. Fields, 57, 453 [Google Scholar]
 Ginsburg, A., Sipőcz, B. M., Brasseur, C. E., et al. 2019, AJ, 157, 98 [NASA ADS] [CrossRef] [Google Scholar]
 Giorgini, J. D. 2015, IAUGA, 22, 2256293 [Google Scholar]
 Gladman, B., Boley, A., & Balam, D. 2019, Res. Notes Amer. Astron. Soc., 3, 187 [CrossRef] [Google Scholar]
 Guzik, P., Drahus, M., Rusek, K., et al. 2020, Nat. Astron., 4, 53 [CrossRef] [Google Scholar]
 Hainaut, O. R., Meech, K. J., Micheli, M., et al. 2018, The Messenger, 173, 13 [Google Scholar]
 Hoang, T., & Loeb, A. 2020, ApJ, 899, L23 [CrossRef] [Google Scholar]
 Hui, M.T., & Knight, M. M. 2019, AJ, 158, 256 [CrossRef] [Google Scholar]
 Hui, M.T., Ye, Q.Z., Föhring, D., et al. 2020, AJ, 160, 92 [CrossRef] [Google Scholar]
 Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [Google Scholar]
 Jewitt, D., Hui, M.T., Kim, Y., et al. 2020a, ApJ, 888, L23 [CrossRef] [Google Scholar]
 Jewitt, D., Mutchler, M., Kim, Y., et al. 2020b, ATel 13611, 1 [Google Scholar]
 Jewitt, D., Kim, Y., Mutchler, M., et al. 2020c, ApJ, 896, L39 [CrossRef] [Google Scholar]
 Kareta, T., Andrews, J., Noonan, J. W., et al. 2020, ApJ, 889, L38 [CrossRef] [Google Scholar]
 Kidger, M. R., & Manteca, J. 2002, Earth Moon Planets, 90, 153 [CrossRef] [Google Scholar]
 Kim, Y., Jewitt, D., Mutchler, M., et al. 2020, ApJ, 895, L34 [CrossRef] [Google Scholar]
 Knight, M. M., Protopapa, S., Kelley, M. S. P., et al. 2017, ApJ, 851, L31 [CrossRef] [Google Scholar]
 Knollenberg, J., Lin, Z. Y., Hviid, S. F., et al. 2016, A&A, 596, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Krolikowska, M., Sitarski, G., & Szutowicz, S. 1998, A&A, 335, 757 [Google Scholar]
 Lin, H. W., Lee, C.H., Gerdes, D. W., et al. 2020, ApJ, 889, L30 [CrossRef] [Google Scholar]
 Manzini, F., Oldani, V., Ochner, P., et al. 2020, MNRAS, 495, L92 [CrossRef] [Google Scholar]
 Marsden, B. G., Sekanina, Z., & Yeomans, D. K. 1973, AJ, 78, 211 [NASA ADS] [CrossRef] [Google Scholar]
 Mashchenko, S. 2019, MNRAS, 489, 3003 [NASA ADS] [CrossRef] [Google Scholar]
 McKay, A. J., Cochran, A. L., Dello Russo, N., et al. 2020, ApJ, 889, L10 [CrossRef] [Google Scholar]
 Meech, K. J., Weryk, R., Micheli, M., et al. 2017, Nature, 552, 378 [NASA ADS] [CrossRef] [Google Scholar]
 Micheli, M., Farnocchia, D., Meech, K. J., et al. 2018, Nature, 559, 223 [NASA ADS] [CrossRef] [Google Scholar]
 Mommert, M., Kelley, M., de ValBorro, M., et al. 2019, J. Open Source Softw., 4, 1426 [CrossRef] [Google Scholar]
 Opitom, C., Fitzsimmons, A., Jehin, E., et al. 2019, A&A, 631, L8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 ‘Oumuamua ISSI Team, Bannister, M. T., Bhandare, A., et al. 2019, Nat. Astron., 3, 594 [CrossRef] [Google Scholar]
 Rafikov, R. R. 2018, ApJ, 867, L17 [CrossRef] [Google Scholar]
 Rinaldi, G., Formisano, M., Kappel, D., et al. 2019, A&A, 630, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rudenko, M. 2016, IAU Symp., 318, 265 [Google Scholar]
 Samarasinha, N. H., Mueller, B. E. A., Belton, M. J. S., et al. 2004, Comets II (Cambridge, UK: Cambridge University Press), 281 [Google Scholar]
 Scargle, J. D., Norris, J. P., Jackson, B., et al. 2013, ApJ, 764, 167 [NASA ADS] [CrossRef] [Google Scholar]
 Scott, D. W. 1992, Multivariate Density Estimation: Theory, Practice, and Visualization (New York: John Wiley & Sons) [Google Scholar]
 Sekanina, Z. 1981, ARA&A, 9, 113 [Google Scholar]
 Sekanina, Z. 1984, AJ, 89, 1573 [NASA ADS] [CrossRef] [Google Scholar]
 Sekanina, Z. 1997, A&A, 318, L5 [Google Scholar]
 Sekanina, Z. 2019, ArXiv eprints [arXiv:1911.06271] [Google Scholar]
 Seligman, D., & Laughlin, G. 2020, ApJ, 896, L8 [CrossRef] [Google Scholar]
 Seligman, D., Laughlin, G., & Batygin, K. 2019, ApJ, 876, L26 [CrossRef] [Google Scholar]
 Sitarski, G. 1990, Acta Astron., 40, 405 [NASA ADS] [Google Scholar]
 Sitarski, G. 1996, Acta Astron., 46, 29 [NASA ADS] [Google Scholar]
 Trilling, D. E., Mommert, M., Hora, J. L., et al. 2018, AJ, 156, 261 [NASA ADS] [CrossRef] [Google Scholar]
 Vazan, A., & Sari, R. 2020, MNRAS, 493, 1546 [CrossRef] [Google Scholar]
 Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Methods, 17, 261 [Google Scholar]
 Whipple, F. L. 1950, ApJ, 111, 375 [NASA ADS] [CrossRef] [Google Scholar]
 Xing, Z., Bodewits, D., Noonan, J., et al. 2020, ApJ, 893, L48 [CrossRef] [Google Scholar]
 Yabushita, S. 1996, MNRAS, 283, 347 [NASA ADS] [CrossRef] [Google Scholar]
 Yang, B., Kelley, M. S. P., Meech, K. J., et al. 2020, A&A, 634, L6 [CrossRef] [EDP Sciences] [Google Scholar]
 Ye, Q., Kelley, M. S. P., Bolin, B. T., et al. 2020, AJ, 159, 77 [CrossRef] [Google Scholar]
 Yeomans, D. K., Chodas, P. W., Sitarski, G., et al. 2004, Comets II (Cambridge, UK: Cambridge University Press), 137 [Google Scholar]
 Zhang, Y., & Lin, D. N. C. 2020, Nat. Astron., 4, 852 [CrossRef] [Google Scholar]
 Zhang, Q., Ye,Q., & Kolokolova, L. 2020, ATel 13618, 1 [Google Scholar]
 Zhou, W. H. 2019, ArXiv eprints [arXiv:1911.12228v2] [Google Scholar]
 Zhou, W. H. 2020, ApJ, 899, 42 [CrossRef] [Google Scholar]
All Tables
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar object 1I/2017 U1 (‘Oumuamua).
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar comet 2I/Borisov (I).
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar comet 2I/Borisov (II).
Heliocentric and barycentric orbital elements, and 1σ uncertainties of interstellar comet 2I/Borisov (III).
All Figures
Fig. 1 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 51P/Harrington. Top left panel: histogram of equatorial obliquity, ɛ_{51P/Harrington}, with bins computed using the Freedman and Diaconis rule, while top right panel: uses the Bayesian Blocks technique. Similarly, the cometocentric longitude, ϕ_{51P/Harrington}, histograms displayed in the lower panels also use the Freedman and Diaconis rule (left) and Bayesian Blocks (right). Histograms are based on data from Table 4 (see the text for details). 

In the text 
Fig. 2 Distribution of thermal lag angle for 51P/Harrington. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 4. 

In the text 
Fig. 3 Distribution of values of the oblateness of the nucleus for 51P/Harrington. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 4. 

In the text 
Fig. 4 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 51P/Harrington. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 4. 

In the text 
Fig. 5 Gaussian kernel density estimation of the spinaxis orientations computed using the SciPy library (Virtanen et al. 2020). Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of comet 51P/Harrington available from the MPC. 

In the text 
Fig. 6 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 1I/2017 U1 (‘Oumuamua). Top left panel: histogram of equatorial obliquity, ɛ_{‘Oumuamua}, with bins computed using the Freedman and Diaconis rule, while top right panel: uses the Bayesian Blocks technique. Similarly, the cometocentric longitude, ϕ_{‘Oumuamua}, histograms displayed in the lower panels also use the Freedman and Diaconis rule (left) and Bayesian Blocks (right). Histograms are based on data from Table 1 (see the text for details). 

In the text 
Fig. 7 Distribution of thermal lag angle for 1I/2017 U1 (‘Oumuamua). Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 

In the text 
Fig. 8 Distribution of values of the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua). Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 

In the text 
Fig. 9 Same as Fig. 6 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to −6.7 as suggestedby Mashchenko (2019). 

In the text 
Fig. 10 Same as Fig. 7 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to −6.7 as suggestedby Mashchenko (2019). 

In the text 
Fig. 11 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 1I/2017 U1 (‘Oumuamua) assuming a very prolate shape, with p = −6.7. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 

In the text 
Fig. 12 Gaussian kernel density estimation of the spinaxis orientations computed using the SciPy library (Virtanen et al. 2020) for a very prolate case (p = −6.7). Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of interstellar object 1I/2017 U1 (‘Oumuamua) available from the MPC. 

In the text 
Fig. 13 Same as Fig. 6 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to 0.836 as suggested by Mashchenko (2019). 

In the text 
Fig. 14 Same as Fig. 7 but assuming a value for the oblateness of interstellar object 1I/2017 U1 (‘Oumuamua) equal to 0.836 as suggested by Mashchenko (2019). 

In the text 
Fig. 15 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 1I/2017 U1 (‘Oumuamua) assuming a very oblate shape, with p = 0.836. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 1. 

In the text 
Fig. 16 Gaussian kernel density estimation of the spinaxis orientations computed using the SciPy library (Virtanen et al. 2020) for a very oblate case (p = 0.836). Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of interstellar object 1I/2017 U1 (‘Oumuamua) available from the MPC. 

In the text 
Fig. 17 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 2I/Borisov. These results correspond to the orbit determination of comet 2I/Borisov in Table 2. 

In the text 
Fig. 18 Distribution of thermal lag angle for 2I/Borisov. These results correspond to the orbit determination of comet 2I/Borisov in Table 2. 

In the text 
Fig. 19 Distribution of the oblateness parameter for 2I/Borisov. These results correspond to the orbit determination of comet 2I/Borisov in Table 2. 

In the text 
Fig. 20 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 2I/Borisov. Left panel: histogram with bins computed using the Freedman and Diaconis rule, while right panel: uses the Bayesian Blocks technique. Histograms are based on data from Table 2. 

In the text 
Fig. 21 Gaussian kernel density estimation of the spinaxis orientations of 2I/Borisov computed using the SciPylibrary (Virtanen et al. 2020) and based on data from Table 2. Top panel: (ɛ_{c}, ϕ_{c}) map. The results in ecliptic coordinates, (λ_{p}, β_{p}), are shown in the middle panel and those in equatorial coordinates, (α_{p}, δ_{p}), are displayed in the bottom panel that includes the actual observations (in red) of 2I/Borisov available from the MPC. Already published spinaxis directions are plotted as filled orange symbols with white error bars (when available): triangles (Ye et al. 2020), square (Manzini et al. 2020), star (Bolin et al. 2020b), and circle (Kim et al. 2020). 

In the text 
Fig. 22 Distributions of equatorial obliquity and cometocentric longitude of the Sun at perihelion for 2I/Borisov. Same as Fig. 17 but for the orbit determination of comet 2I/Borisov in Table 3. 

In the text 
Fig. 23 Distribution of thermal lag angle for 2I/Borisov. Same as Fig. 18 but for the orbit determination of comet 2I/Borisov in Table 3. 

In the text 
Fig. 24 Distribution of the oblateness parameter for 2I/Borisov. Same as Fig. 19 but for the orbit determination of comet 2I/Borisov in Table 3. 

In the text 
Fig. 25 Distributions of spinaxis orientations, (λ_{p}, β_{p}), for 2I/Borisov. Same as Fig. 20 but based on data from Table 3. 

In the text 
Fig. 27 Same as Fig. 26, bottom panel, but based on data from Table 5. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.