Transiting exocomets detected in broadband light by TESS in the $β$ Pictoris system
Sebastian Zieba Konstanze Zwintz Matthew A. Kenworthy Grant M. Kennedy
 

   
      Transiting exoplanet around 

   

   Universität Innsbruck, Institut für Astro- und Teilchenphysik,
             Technikerstraße 25, A-6020 Innsbruck

              konstanze.zwintz@uibk.ac.at
         Leiden Observatory, Leiden University, P.O. Box 9513, 2300 RA Leiden, The Netherlands
        Department of Physics, University of Warwick, Gibbet Hill Road, Coventry CV4 7AL, UK
             

   

  
  
   
  
   We search for signs of falling evaporating bodies (FEBs, also known as exocomets) in photometric time series obtained for  after fitting and removing its δ Scuti-type pulsation frequencies. 
  
   Using photometric data obtained by the TESS satellite we determined the pulsational properties of the exoplanet host star  through frequency analysis. We then pre-whitened the 54 identified δ Scuti p-modes and investigated the residual photometric time series for the presence of FEBs.
  
   We identify three distinct dipping events in the light curve of  over a 105-day period. These dips have depths from 0.5 to 2 millimagnitudes and durations of up to 2 days for the largest dip. These dips are asymmetric in nature and are consistent with a model of an evaporating comet with an extended tail crossing the disc of the star.
  
   We present the first broadband detections of exocomets crossing the disc of β Pictoris, complementing the predictions made 20 years earlier by <cit.>. No periodic transits are seen in this time series. These observations confirm the spectroscopic detection of exocomets in calcium H and K lines that have been seen in high resolution spectroscopy.

   

   Transiting exocomets detected in broadband light by TESS 

   in the β Pictoris system
    S. Zieba 1
         K. Zwintz1
          M. A. Kenworthy2
          G. M. Kennedy 3
          
    Received 26 March 2019; accepted 
==================================================================================================






§ INTRODUCTION


Our voyage of extrasolar planetary system discovery and exploration arguably began with the image of the circumstellar disc of β Pictoris <cit.>, which showed that the infrared excesses discovered by the Infrared Astronomical Satellite  (IRAS) <cit.> were not “shells”, but discs on scales similar to the extent of our solar system.

β Pictoris (HD 39060) has a V magnitude of 3.86 and lies at a distance of 19.76 pc <cit.>. It is a member of the β Pictoris moving group with an age of ∼23 Myr <cit.>. 
β Pictoris harbours a warped debris disc composed of dust and gas, which is seen nearly edge-on and has an outer extent ranging from 1450 AU to 1835 AU <cit.>. In 2010, the planet β Pictoris b was detected using observations with the NaCo instrument at the Very Large Telescope (VLT)  by <cit.>. 

Because of its close proximity and circumstellar disc, the β Pictoris system can be considered an ideal test bed to study the formation and evolution of planetary systems, including minor bodies such as exocomets and exomoons. 
Indeed, the first reports of comets around stars other than our Sun came from spectroscopic observations of the β Pictoris system, which identified transient absorption features that change in time <cit.>. Similar detections have also been made towards other stars, for example HD 172555 <cit.>. Models based on these detections have included cometary comae and have placed the exocomets at distances of a few to several tens of stellar radii, which was confirmed by their observed acceleration <cit.>. Subsequent work theoretically predicted that occultations of the stellar light by these exocomets could be detected in broadband photometric time series. <cit.> and <cit.> estimated that the depths of such signals would be of the order of tenths of a percent in flux and suggested the detection of an exocomet passage in front of β Pictoris in observations from 1981 <cit.>.

To date just a few stars show photometric transits consistent with the predictions. The first of these were discovered by <cit.> around the two stars KIC 3542116 and KIC 11084727. The dips seen for these stars show a remarkable similarity to the models of <cit.>, which strengthens the exocomet interpretation. A more recent study by <cit.> describes an algorithm to identify transit features with an ingress that is more rapid than the egress, i.e. such as expected for exocometary bodies with a long tail. This work recovered the events discovered by <cit.>, found one more around HD 182952, and argued that the host stars are likely young (∼100Myr old). For completeness, we note that the deep and irregular dimming events seen towards KIC 8462852 have also been interpreted as having a cometary origin <cit.>, although the shape and depth of these events are very different to those discovered by <cit.> and <cit.>.


In this letter we present observations of  from TESS and show that the duration and depth of a dipping event in Sector 6 are consistent with the transit of an exocomet and agree with the predictions made by <cit.>.






§ OBSERVATIONS



We used data collected by the Transiting Exoplanet Survey Satellite <cit.> from 19 October 2018 to 1 February 2019 in Sectors 4 through 7.
This leads to a baseline of about 105 days with a duty cycle of 85% for the four sectors.
 (TIC 270577175, T = 3.82 mag) is one of the preselected targets for which “short cadence” 2-minute data are provided.
We investigated the instrumental behaviour of TESS during the β Pictoris observations in detail to exclude any systematic effects or instrumental artefacts. A detailed description is given in the Appendix.





§ RESULTS




 §.§ Fitting the stellar pulsations of Beta Pictoris
 

 is known to show δ Scuti-like pulsations <cit.>, which we needed to fit and remove in order to search for signs of exocomets. 
Therefore, we conducted a frequency analysis using the software package Period04 <cit.> and extracted all pulsation frequencies down to a signal-to-noise ratio of 4 following <cit.>.

From this, we identified 54 significant p-modes in the frequency range from 23 and 76. A one-day zoom into the light curve illustrating the pulsational behaviour is shown in the upper panel of Figure <ref>.

A detailed description, the full list of frequencies (see Table <ref>) and the corresponding amplitude spectrum (see Figure <ref>) are provided in the Appendix.

We subtracted the multi-sine fit using the 54 frequencies, amplitudes, and phases from the original TESS light curve and investigated the residuals (lower panel in Figure <ref>) for the presence of FEBs.






 §.§ Exocomet transit modelling


The pulsation removed light curve shows three notable transit features, as highlighted in Figure <ref>.

All three transits show a triangular or sawtooth shape that has a steep increase of absorption followed by an exponential decay back to the full flux level of the star.

The three transits have different shapes and are not periodic. Assuming a circular orbit with period P>105 days and a mass of 1.80 M_⊙  for β Pictoris <cit.> leads to a minimum orbital radius of about 0.5 AU. It is almost certain that the comet orbits are highly eccentric, so this non-periodicity is not particularly surprising.

<cit.> carried out a comprehensive modelling study of exocomet light curves that includes modelling of the tail ejecta and its particle size distribution and the subsequent acceleration due to radiation pressure from the star for different distributions of the sizes of the ejected particles.

Not knowing the geometry of the exocomet orbit leads to degeneracies in the possible range of particle distributions and orbital parameters, so in this Letter we adopt a simpler 1-D model as used in <cit.>. While that model was designed with a transiting planet in mind, the only practical difference in this work is that the orbital period and eccentricity are largely unconstrained. We therefore assume a circular orbit with period P and orbital radius a for the modelling, but report the results in terms of the transverse velocity v_trans of the comet across the face of the star.

The angle between the observer, centre of the star, and the comet θ is related to the time t by θ= 2π(t-t_mid)/P, where t_mid is the time of the transit following the nomenclature of <cit.>. The star is modelled as a limb-darkened circular disc with radius R_*=1.53 R_⊙ <cit.> and linear limb darkening coefficient u=0.275 for β Pictoris <cit.>.


The exocomet is assumed to be optically thin, and its extinction cross-section ρ is expressed in units of stellar area. 

The model has a hard front edge and the cross section then drops exponentially  behind the comet with normalisation constant c_e and an exponential factor λ so that ρ(Δθ)=c_e e^-λ(Δθ), where Δθ=(θ - θ') is the angle between the position of the comet and an arbitrary point on the circular orbit of the comet.

The characteristic length of the tail is therefore 1/λ.


The comet crosses the limb-darkened disc of the star with an impact parameter b and chord length r_c, related by b=[1-(r_c/2R_*)]^1/2.

The total extinction of the star is the convolution of ρ with the stellar disc, and so the resultant intensity I of the star is


    I(θ) = 1 - ∫^2π_0 ρ(θ-θ') i (θ',r̂_c) dθ'.


The function i(θ',r̂_c) is the limb-darkened intensity of the star along a given chord defined by

    i(θ',r̂_c)=1-u [ 1- a/R_*√(sin^2(r̂_c/2)-sin^2 θ') ],


which is set to zero for r̂_c/2<θ ' < 2π - r̂_c/2, where
θ' is the angle between the centre of the star and a point along the orbit of the comet and the angular size of the chord as seen from the comet r̂_c= arcsin(r_c/2a).


 With five parameters (t_mid, b, c_e, λ, P)  the model is underconstrained in that transits with more rapid ingress at large b (shorter transit chords) are compensated for by longer orbital periods and shorter comet tails, thus yielding the same ingress or egress slope in the model light curve. Fixing b=0 constrains P and λ, although if the comet were not assumed to have a hard front edge these would again be degenerate.

We use the Python software package emcee <cit.> to perform a Markov Chain Monte Carlo (MCMC) fit to the binned light curve around the largest transit feature.

We average the photometry over 30-minute bins, and calculate the error bars on the binned photometry by determining the root mean square of each set of photometric points and dividing by the square root of the number of points per bin.

The photometric light curve and the best-fit model are shown for b=0 in Figure <ref>, along with the residuals of the fit in the lower panel. The resultant triangle plot showing the parameter distributions and correlations between the variables is shown in Figure <ref>.

The ±1σ errors of the parameters are computed by measuring the interval between 16% and 84% of the merged distribution, and the best-fit values and their uncertainties are listed in Table <ref>.











§ DISCUSSION


As noted above, the comet orbit is not constrained by our model because the chord across the star is unknown, yielding a near-total degeneracy between b, λ and v_trans, which is derived from P. By setting b=0 the transverse velocity is constrained, for our assumption of a hard front edge for the comet, and has a value of v_trans=19.6 km s^-1 (i.e. equivalent to the circular velocity at about 3 AU). This transverse velocity is surprisingly slow given that the radial velocities of the comets detected in calcium absorption are typically tens to hundreds of km s^-1 and imply distances at transit of only a few tens of stellar radii. This discrepancy therefore implies a high eccentricity (i.e. near-radially plunging) orbit if the spectroscopic and photometric events were from objects on similar orbits. While high eccentricity orbits are expected, a further issue is that the characteristic transit shape seen in this work and by <cit.> is not expected for comets with arbitrary pericentre orientations <cit.>. Further orbital constraints may therefore be derived from the TESS data, which we leave for more detailed future modelling work. 

The characteristic tail length 1/λ is small. However, as with the velocity this is also a transverse length, and in addition the physical size depends on the distance to the star. If we assume 1 AU, then this length is ∼2 × 10^8 km, or about a tenth of the stellar diameter, indicating that most of the dust is concentrated near the comet nucleus. The close match of the data with the models of <cit.> suggests that those physical models are a reasonable starting point to further interpret these data.

We chose not to include forward scattering for our exocomet model, even though <cit.> and others include it in their models.

Our hesitation is that the long-term variation in the flux of  over Sector 6 seems to not extend or be seen in the adjacent sectors.

An initial examination suggests that there is a long-lived forward scattering halo associated with transit A (illustrated in the rightmost panel of Figure <ref>), but it is not centred on the transit itself.

Forward scattering peaks are seen in transit B and C. However because these transits are at much smaller significance than transit A, there is a much larger degeneracy in the model fitting for these two dips, and exploring whether these dips have a common orbital connection with transit A is left to a future paper.



§ CONCLUSIONS


We report the first photometric broadband detection of an exocomet transit in the  system, consistent with the prediction of <cit.>.

We fit a simple cometary tail model to the largest of these events, seen at BTJD 1486.4, and find that it is consistent with an exponentially decaying optical depth tail convolved with the limb-darkened disc of the star.

Future work includes a more comprehensive modelling of all the transits seen in the TESS light curve, and a  study of whether high resolution spectral observations were taken simultaneously with this transit to confirm the cometary hypothesis.

A comprehensive photometric campaign was conducted during the Hill sphere transit of β Pictoris b <cit.>, which includes bRing <cit.>, ASTEP   <cit.>, and BRITE-Constellation <cit.> photometry.

Searching this data may yield more candidate broadband transits.

A dedicated cubesat mission to monitor , similar to the PICSAT project  <cit.> would uncover more of these events, and photometry at different wavebands would help determine the size distribution of the tail particles <cit.>.



This paper includes data collected by the TESS mission, which are publicly available from the Mikulski Archive for Space Telescopes (MAST).

Funding for the TESS mission is provided by the NASA Explorer Program.

This work has made use of data from the European Space Agency (ESA) mission Gaia (<https://www.cosmos.esa.int/gaia>), processed by the Gaia Data Processing and Analysis Consortium (DPAC; <https://www.cosmos.esa.int/web/gaia/dpac/consortium>).

Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

GMK is supported by the Royal Society as a Royal Society University Research Fellow.

We made use of the software package Period04 <cit.>, the Python programming language <cit.>, and the open-source Python packages numpy <cit.>, matplotlib <cit.>, astropy <cit.>, lightkurve <cit.>, and emcee <cit.>.

An on-line repository with materials used in this work is available at <https://github.com/sebastian-zieba/betaPic_comet>.



aa







§ INSTRUMENTAL BEHAVIOUR






 §.§ TESS photometric time series

 (TIC 270577175, T = 3.82 mag) is one of the preselected targets for which “short cadence” 2-minute data is provided. This Candidate Target List (CTL) is a subset of the TESS Input Catalogue[<https://mast.stsci.edu/portal/Mashup/Clients/Mast/Portal.html>] <cit.> with about 200 000 targets for TESS 2-minute cadence observations, which were primarily chosen in order to maximise the yield of transiting exoplanets.


Figure <ref> shows the complete TESS  light curve.

The beginnings of each  sector are delineated with a vertical red line.

The visible gaps are related to the data downlink at the perigee of the orbit of TESS, which occurs every 13.7 days and during which TESS halts its observations for about one day.

 was observed with CCD 1 of Camera 4 during Sector 4, with CCD 4 of Camera 3 during Sector 5 and with CCD 3 of Camera 3 during Sectors 6 and 7.
As the angle between the boresight of these cameras and both the Earth and the Moon was never smaller than 37^∘, we do not expect strong features created by scattered light from these two bodies during these observations. 




We used the 2-minute Presearch Data Conditioning <cit.> light curve from the Science Processing Operations Center (SPOC) pipeline <cit.>, which was originally developed for the Kepler mission <cit.>.

These light curves were corrected for systematics by the SPOC pipeline.

During the “momentum dumps”, i.e. thruster firings that reduce the speed of the reaction wheels and decrease the pointing jitter, TESS loses the fine attitude control mode for about 15 minutes (see the TESS Instrument Handbook[<https://archive.stsci.edu/missions/tess/doc/TESS_Instrument_Handbook_v0.1.pdf>]).

Consequently, the pointing is less stable during these particular times, resulting in potential changes in the photometric fluxes of the target stars.

Momentum dumps took place every ≈2.5 days in Sector 4, every three days in Sector 5, and every 3.125 days in Sector 6 and 7 (red triangles in Figure <ref>).

We omitted data collected in the vicinity of momentum dumps in the subsequent analysis.

The data products were accessed and modified with the python package lightkurve <cit.>, which retrieves the data from the MAST archive[<https://archive.stsci.edu/tess/>].

We removed every measurement with a non-zero “quality” flag (see §9 in the TESS Science Data Products Description Document[<https://archive.stsci.edu/missions/tess/doc/EXP-TESS-ARC-ICD-TM-0014.pdf>]), which marks anomalies like cosmic ray events or instrumental issues. 

Three days between BTJD 1421 and BTJD 1424 were removed owing to higher than normal rates of spacecraft jitter.

This occurred just after an instrument anomaly (see the Data Release Note of Sector 4 for more information[<https://archive.stsci.edu/tess/tess_drn.html>]) and is visible in Figure <ref> as the first gap, which is noticeably bigger than the others. 

Finally, we normalised the four sectors by dividing each of the sectors by their respective median flux and we combined these into one light curve.





 §.§ Identification of possible transit events



Sector 5 saw an automated creation of a data validation (DV) report, where the pipeline identified two decreases in brightness at about BTJD 1442 and BTJD 1459.

These dips are indicated by the leftmost two downward-pointing black arrows in Figure <ref>.

The third arrow represents a clearly visible change in brightness for about two days at BTJD 1486.
This arrow is asymmetric in shape and very similar to the predicted light curve generated by a comet with an extended tail transiting the disc of the star, as hypothesised by <cit.> and <cit.>.

We used the apertures which were created by the pipeline visible as the white boxes in Figure <ref>.

The apertures have an elongated shape due to the blooming effect with a typical size of about 30 rows by 4 columns.

To determine the robustness of these measurements, custom apertures of differing geometries were used to extract the photometry of  - all three dips appear in these reductions.

As a quick check, we used the TESScut function of lightkurve, which allows us to create a light curve out of the full frame images (FFIs).

Those products are generated every 30 minutes by TESS.

However, we have to consider that those - in contrast to the light curves provided by the pipeline - do not have an aperture mask and are not background corrected.

After trying out different aperture masks and background masks for , we report that the dimming features are also distinctly visible in the FFIs. 

We then inspected the postage stamp images, which have a cadence of 2 minutes.

As the background is subtracted from the flux of the star, a bright object passing by, such as an asteroid,  could temporarily increase the background flux and thus decrease the flux of the star.

However, nothing like this is visible in the postage stamps during the dimming events.

The times of the event do not coincide with momentum dumps (as can be seen in Figure <ref>) and the jitter of TESS also did not increase.

We checked for similar dimmings in the stars closest to  that have short cadence data.

TIC 270574544 and TIC 270626557 have a separation of 810 and 1000 arcseconds relative to , which corresponds to 39 and 48 pixels, respectively.

Those two stars were observed in the same sectors as  and on the same CCDs and no similar decrease in brightness is visible.

Different background apertures in the vicinity of  show no dips either.

These last two checks also rule out a “rolling shutter” effect.

The brightest star besides  in the pipeline aperture is Gaia DR2 4792773113918130560 in Sector 5 and Gaia DR2 4792772873399766528 in Sector 6 <cit.>. 

Those stars are too faint to create those dimming events. 

Stars in the postage stamps of  with a luminosity high enough to create those dips are indicated with red circles in Figure <ref>.

Finally we checked the possibility of Intra Camera Crosstalk as the cause of the dimmings. 

A bright object can create a decrease in flux in other parts of the same CCD (see §6.8.6 of the TESS Instrument Handbook[<https://archive.stsci.edu/missions/tess/doc/TESS_Instrument_Handbook_v0.1.pdf>]).

Those echoes can be found on the same row and 1024 columns higher than the bright source and they were only seen in CCD 3 of each camera.

The big dip in Sector 6 was observed with CCD 3, however 1024 columns above  lies outside the imaging area of this CCD.
We conclude that  is the true host of the dimming events.









§ DETAILS ON THE PULSATION FREQUENCY ANALYSIS


For the pulsation frequency analysis we removed all Sector 4 data as the PDCSAP flux exhibits strong variations.

As these are not visible in the SAP light curve, they are most likely artefacts of the reduction.

Consequently, this part of the light curve would introduce artificial frequencies and, hence, was omitted for the pulsational analysis.

The total time baseline of the TESS light curve used for the pulsational analysis is 78 days with a duty cycle of 80.5%.

The original amplitude spectrum of β Pictoris based on TESS data is shown in Figure <ref>, where the detected pulsation frequencies are represented in red.

The 54 δ Scuti-type frequencies, amplitudes and phases detected in the TESS data for β Pictoris are listed in Table <ref>. Frequency, amplitude and phase errors are calculated using the relations by <cit.>.
The zero point for the phase calculation is at BTJD = 1410.9035 d.








