Precision Light Curves from TESS Full-frame Images: A Different Imaging Approach

The fully calibrated FFIs exhibited a low-frequency residual sky background, which is meant to represent the zodiacal background and faint contamination from galaxies and unresolved stars (Stassun et al. 2018). We applied a residual background subtraction following the approach of Wang et al. (2013) and Oelkers et al. (2015). The residual background model is constructed by sampling the sky background every 32 × 32 pixels over the entire CCD, and replacing pixels with values of ADU  > 25,000, with the median sky background in each box. A model sky is then generated by interpolating between all boxes to make a "thin plate spline" (Duchon 1976). We used the IDL implementation GRID_TPS and the PYTHON implementation of Rbf to make the splines that are then subtracted from every frame.

Difference imaging requires precise frame alignment in order to produce a proper subtraction. Slight variations caused by improperly aligned frames can contribute to poorly measured fluxes, which will introduce additional dispersion in the extracted light curves. While the TESS pointing is expected to be much better than 1 pixel, spacecraft jitter is likely and has been introduced into the simulated FFIs. We use the world coordinate solution (WCS) from the image headers to align the images using the IDL implementation of HASTROM and PYTHON implementation of HCONGRID with cubic-spline interpolation.⁷

A high-quality master frame, which will be convolved and subtracted from each science frame, is required to preserve the precision of the extracted photometry. Typically, this frame is generated by median combining many individual frames with high S/N and the best seeing, obtained throughout the observing campaign. We created our master frame by median combining all 1348 images. To save machine memory, images are median combined in sets of 50, producing 26 temporary master frames, which are then combined into a single final master frame. This procedure does not apply any pixel weighting or treatment for cosmic rays, but could be improved with such measures in future reductions.

Next, we identified all 106,399 stars in the TIC with TESS magnitude (T) T < 16 that were within the master frame field of view. We then calculated the curve of growth for a variety of aperture sizes, and selected the aperture radius that optimizes the flux calculation for the majority of the stars on the frame. The flux within the aperture is calculated using a polygon approximation of the intersection between a circular aperture and a square pixel, normalized by the total area of the sum of the pixels to match the circular aperture precisely. We found this aperture radius size to be 2.5 pixels. This is consistent with the measurements of the typical PSF FWHM of 1.88 pixels, which would transcribe >99.7% of a typical star's flux within 2.39 pixels.

Fluxes for all stars were then extracted from the master frame using fixed-aperture photometry. We defined the zeropoint offset between the instrumental magnitude scale and the TIC magnitude scale as the median difference between the instrumental magnitude and the TIC T magnitude for all stars. We found the zeropoint to be 4.825 mag. We note that for a small number of stars we found a different median zeropoint offset that was ∼0.75 mag brighter than the offset above. This suggests that a second star of roughly equal brightness is present in the simulated FFIs at the same location, even though only one TIC object exists at that location. This 0.75 mag different offset was found only around some, but not all, bright stars, T < 11. These "duplicate" stars are a known issue in the TIC for version 5 and earlier, and were caused by ghost artifacts in the SDSS extended source catalog that occur near bright stars. These objects have been removed from future TIC versions. Therefore, we accepted the zeropoint offset that was consistent with the majority of stars (Stassun et al. 2018).

Typically, DIA routines use an adaptive kernel, K(x, y), defined as the combination of 2 or more Gaussians. While effective at modeling well-defined, circular PSFs, this kernel has difficulty properly fitting other PSF shapes, particularly for highly distorted stars near the edge of the frame (Bramich 2008; Miller et al. 2008). The TESS PSF varies quite significantly across the frame. Therefore, we used a Dirac δ-function kernel to compensate for the noncircular, irregular PSF shape. We defined the kernel as

$$\begin{eqnarray}&&K(x,y)=\displaystyle \sum _{\alpha =-w}^{w}\displaystyle \sum _{\beta =-w}^{w}{c}_{\alpha ,\beta }(x,y){K}_{\alpha ,\beta }(u,v),\end{eqnarray} \tag{ 1 }$$

where ${K}_{\alpha ,\beta }$ is a combination of (2w + 1)² δ-function basis vectors and K_0,0 is the centered δ function (Bramich 2008; Miller et al. 2008). We defined our basis vectors to ensure a constant photometric flux ratio between images (Alard 2000; Bramich 2008; Miller et al. 2008). In the case of $\alpha \ne 0$ and $\beta \ne 0$,

$$\begin{eqnarray}&&{K}_{\alpha ,\beta }(u,v)=\delta (u-\alpha ,v-\beta )-\delta (u,v),\end{eqnarray} \tag{ 2 }$$

while for α = 0 and β = 0,

$$\begin{eqnarray}&&{K}_{\mathrm{0,0}}(u,v)=\delta (u,v).\end{eqnarray} \tag{ 3 }$$

Stamps are taken around at most 500 bright, isolated⁸ stars to solve for the coefficients ${c}_{\alpha ,\beta }(x,y)$ using the least-squares method. We allowed ${c}_{\mathrm{0,0}}(x,y)$ to be spatially variable to compensate for variations that may correlate with position on the detector: such as poor sky subtractions, blending, or imperfections in the flat-field. Stars are identified as suitable candidates for stamps provided their measured photometric error was < 50 mmag and there were no other, brighter stars within 3 pixels of the target star. We found that we could use a 5 × 5 pixel kernel across the frame without appreciably affecting the quality of the subtraction, but significantly reducing the runtime.

We set the photometric extraction aperture at a 2.5 pixel radius (525) and set the image background to be the value determined in Section 3.1. The differential flux was then combined with the flux from the reference frame, and zero-pointed using the offset described above. Finally, the photometric errors were rescaled following a methodology similar to Lupton et al. (1989) and Kaluzny et al. (1998).

Some low-level systematics persisted in the extracted light curves, even with the care taken to properly preprocess and difference the images. While the exact source of these systematics is unknown, the likely explanations include improper image alignment, imperfect stellar stamp selection, minor variations in the calibration process, or a combination of these and other factors.

We opted to remove (detrend) these systematics from the photometry using an ensemble of light curves to craft a model of nonastrophysical signals. First, for each star, we identified the 1000 closest stars of similar magnitude (∼±0.1 mag). We have found that stars of similar magnitude tend to show similar systematics for a variety of reasons, including response to similar pixels, optical effects on a particular area of the CCD, blending effects from bright stars, and local background from sources such as zodiacal light or the Galactic plane. Atmospheric effects cause similar systematics in nearby stars in ground-based surveys, but this would not affect stars observed by TESS. Then we identify a subset of stars that, when subtracted from the light curve, decrease the root-mean-square (rms) variations of the light curve by at least 10%. Stars that do not decrease the rms are discarded from the model. If no comparison stars were found to decrease the rms of the target star, then no detrending was applied. If multiple trend stars were found to decrease the rms of the target light curve, then they were median combined into a master model.

During this process, we noticed that many stars exhibited abrupt changes in magnitude, typically after 48 observations—which is 1 day given the 30 minute FFI candence. These changes in magnitude are likely a recovery of how the simulated images were generated. While a selection of stars were fully rendered throughout the entire observing season, most stars were regenerated once per day. This effect created a visual break in the light curves every 48 data points (Jenkins et al. 2018). We removed these offsets by subtracting subsequent data points in a given trend model, and identifying observations where the deviation between data points was larger than 5σ of the mean deviation between data points. Each subsection of the trend model was then scaled to match the median magnitude of the light curve during this subset. However, if subtracting the scaled trend did not improve the rms by more than 5%, the trend during this subsection of the light curve was not scaled. This was done to ensure that true astrophysical signals would not be removed in the detrending routine.

Figure 2 shows an example of this detrending and its proficiency at removing nonastrophysical signals. While the detrending routine was proficient at reducing the rms of many stars, there were still stars that failed the detrending process (see Figure 2, right). We mention these artifact signals as a caution, but note that these events are uncommon in our tests and likely affect ∼2% of stars.

Zoom In Zoom Out Reset image size

Of course, systematics typically vary greatly between systems, and the systematics of the simulated FFIs analyzed here may not fully capture the full systematics of the final, true FFIs after TESS launch. Therefore, we propose several methods that users can implement and that we plan to include with the final release of the FFI data products. These methods include:

• 1.  
  Position based detrending: While the pointing of the telescope is expected to be exceptional, it will not be immune to slight corrections for drift. One can use the correlation of time and a star's x, y position to remove trends that appear as the star's centroid shifts from pixel to pixel. While in theory, many of these correlations should be accounted for within the δ-function kernel, in practice, we have found that for very wide-field images, the δ-function kernel is not a perfect representation across the entire image and small systematics remain on millimag levels. Given the light curves have a precision of parts-per-million, we expect there to be variations that cannot be adequately fit for by the kernel, and instead we suggest one may be able to be remove these systematics with a "residual-flat-field" as shown in Wang et al. (2013) and Vanderburg & Johnson (2014).
• 2.  
  Magnitude-based detrending: This method loosely follows the plan we outline above. Stars of similar magnitude can be combined to make a model of variations that are not astrophysical in nature. Typically, two stars of similar magnitude should not vary in identical astrophysical patterns between images, even if they were the same spectral type or variable type.
• 3.  
  Fourier-based detrending: By investigating the power spectrum of the light curve, multiple stars may show similar periods, or aliases of similar periods. These stars could be combined to create a trend pattern, if the period is known a priori to be spurious.

The quality of a differenced image can be quantified using stars that cleanly subtract in an image. Specifically, the quality of a differenced frame can be described by assessing the degree to which the deviations left in the differenced frame match the expectations for the noise from the science frame and the master frame (Alard 2000).

To show that our differenced frames match the expectations for the noise, we select a frame and normalize each pixel value by the combination of the noise from the science frame and the master frame. We defined the expected noise⁹ as $\delta =\sqrt{{I}_{N}+{R}_{N}}$, where I_N is the photon counts in the science image and R_N are the photon counts in the master frame, both normalized by the number of frames used to create each image. If the subtraction generally matches the expectations from Poisson noise, the normalized pixels should show a Gaussian distribution centered near 0 with a standard deviation of 1. Figure 3 shows the histogram of pixel values for a typical differenced frame normalized by δ. The normalized residuals show a mean of 0.132 and a standard deviation of 0.98. We believe this bias is likely due to the sky background, given each pixel is 21'' in size, and we only account for sky background during the preprocessing of the images as described in Section 3.1. It is possible that the procedure could be improved with an additional sky subtraction after the differencing step, but scatter in the normalized pixels is nearly 1, which suggests only minimal levels of noise are being introduced by this bias. We accept this as sufficient evidence that the residuals in the differenced frame are closely matching the expected Poisson noise.

Zoom In Zoom Out Reset image size

Next, we compared our resulting light curves to the current models of expected TESS precision as a function of stellar magnitude. In particular, we use the original noise model from Ricker et al. (2014) and Sullivan et al. (2015) and the updated noise model from J. Pepper et al. (2018, in preparation). These models incorporate shot noise from the star, sky background (in the case of TESS, contamination from zodiacal light, unresolved stars, and background galaxies), the read noise of the detector, and a mission-specified noise floor of 60 ppm on 1 hr timescales.¹⁰ J. Pepper et al. (2018, in preparation) included an additional estimation of the contamination by nearby stars, which we have excluded from our model because this effect is present in the simulated FFI data through physical contamination, and is manifested in our analysis by stars that deviate from the expected noise floor. These noise models differ quite substantially, particularly in the sky-dominated regime. This is principally because Ricker et al. (2014) and Sullivan et al. (2015) used an initial optimal-aperture size estimate, while J. Pepper et al. (2018, in preparation) use an updated optimal-aperture model based on lab testing of the cameras. We adopt the J. Pepper et al. (2018, in preparation) model, as a more up-to-date estimation of the photometric precision,¹¹ but we include both models for comparison.

As shown in Figure 1 (bottom left), our pipeline very satisfactorily reproduces the expected noise floor of J. Pepper et al. (2018, in preparation); the noise floor model matches the lower envelope of stars, representing nonvariable objects. Stars above the lower envelope noise floor are interpreted to represent variable objects of various types (see below).

Wide-field imagers typically have a PSF that is largely dependent on the location of the stars on the detector, and the severity of the noncircular PSF shape can be exacerbated near the edges of the frame (Pepper et al. 2007). Since the shape of the PSF is largely position-dependent on the simulated FFIs, we checked 1024 × 1024 pixel subimages in the four CCD corners to see if the precision changed as a function of detector position. As shown in Figure 1 (lower right) there is no appreciable difference on any part of the frame.

The data release notes for the simulated FFIs indicate that several hundred variable-star and planet-transit signatures were injected into the FFIs to enable checks of signal recoverability. We searched for stars in the simulated data that showed stellar variability and/or transit candidates using the variability metrics of Oelkers et al. (2018), a basic Lomb–Scargle periodicity search (Lomb 1976; Scargle 1982, LS;), and the box-least-square (BLS; Kovács et al. 2002) algorithm to identify transit-like events.

The variability metrics of Oelkers et al. (2018) were shown to work well to identify stars with large amplitude variability in time-series data. We calculated the metrics (Stetson 1996; Wang et al. 2013, rms, Δ90, Welch–Stetson J and L) for every star in the frame, and each star has their specific metric values compared to stars of similar magnitude (for rms and Δ90 metrics) or to the entire field (for J and L metrics). Stars that have metric values larger than a threshold (typically +2σ for rms and Δ90; and +3sσ for J and L) are typically considered variable. These metric values are calculated for every star, and are included in our data release.

We provide metrics for periodic signals using two methods. First, we ran a basic LS search (Lomb 1976; Scargle 1982) for periods between 0.01 and 27 days. Second, we ran a BLS search (Kovács et al. 2002) for periods between 0.3 and 27 days with 10,000 frequency steps and 100 phase bins. The LS and BLS results are included in our data release.

Typically, heuristic cutoffs will be identified using the data from the entire frame to identify the most likely astrophysical events. In the future, we plan to detect events using a methodology similar to previous work (Wang et al. 2011, 2013; Rodriguez et al. 2017; Ansdell et al. 2018; Oelkers et al. 2018). We visually inspected light curves of several stars identified in the official "Ground-Truth" lists released by NASA (Jenkins et al. 2018), and we have plotted example stars in Figure 4. Additionally, we show the precision achieved for each "Ground-Truth" object observed in our analysis in Figure 5. Some objects with large rms values are clearly not "Ground-Truth" objects, but instead were identified as poorly detrended stars (<2% of all objects). In any case, we include all calculated metrics in our data release portal (see below) to allow users to experiment with the data using different significance thresholds.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size