Astrometry and Occultation Predictions to Trans-Neptunian and Centaur Objects Observed within the Dark Energy Survey

The very first step consists of obtaining the necessary information—pointing, observing date, location in the DES database, among others—on all CCD images acquired during the first three years of observations within the DES. This was done through easyaccess (Carrasco Kind et al. 2018), a friendly structured query language (SQL)-based tool to query the DES database. The result from such a query was a file containing the metadata from 4,292,847 CCD images. This file then feeds into the Sky Body Tracker (SkyBoT; Berthier et al. 2006).

SkyBoT is a project aimed at providing a virtual observatory tool useful to prepare and analyze observations of solar system objects. In addition to the web-interface service it offers, queries are also possible from the command line. The basic inputs to a cone search,⁴¹ for instance, are IAU identification of the observatory, J2000 pointing coordinates of a given CCD image, observation date, and a region centered on the pointing coordinates. All of these data come from the metadata previously mentioned. The output is a text or VOTable file format with pieces of information on all of the known small solar system bodies inside the given region, such as their J2000 astrometric right ascensions and declinations, V magnitudes, names and numbers (when they are numbered), and dynamical classes, among others. Table 1 lists the total number of TNOs and Centaurs found in the DES images as well as the expected number of objects for which positions can be determined from them. As we will see later in the text, these expected numbers (column 4 in particular) were surpassed.

The result of the search with the SkyBoT was a file having 1,708,335 entries, most of them of around 140,000 main-belt asteroid objects in more than 1.5 million CCD images. These objects, in addition to a few thousand members of other dynamical classes also found in the images, are being treated separately.

Note that the detection of a TNO or Centaur is not expected for all of the selected CCD images. Objects that are faint (V ≳ 24.0) in the DES images, or images taken under non-transparent sky, may not provide a detectable signal of the target. The most frequent exposure time of the DES frames presented here is 90 s (see Morganson et al. 2018).

Our astrometric tool is the Platform for Reduction of Astronomical Images Automatically (PRAIA; Assafin et al. 2011) package. PRAIA was conceived to determine photometry and accurate positions from large numbers of CCD images as unsupervised as possible. Its use and performance have been reported by various works (see, for instance, Assafin et al. 2013; Thuillot et al. 2015; Gomes-Júnior et al. 2016) from reference frame to solar system studies. The reference catalog used here for astrometry is the Gaia Data Release 2 (Lindegren et al. 2018). All differences in R.A. as well as all uncertainties related to measurements along R.A. are multiplied by the cosine of the decl.

A Intel(R) Xeon(R) CPU E5-2650 v42.20 GHz configuration, using 40 cores, reduces 1000 CCD images in 20 minutes from a parallelized run of PRAIA. A total of 12,561 CCD images were treated here.

The presence of distortion effects, also known as the field distortion pattern (FDP), are expected in detectors with large FOVs such as that of the DECam. Common solutions are, e.g., the use of a high-degree polynomial (not always recommended) to relate CCD and gnomonic coordinates of reference stars, the brute-force determination of a distortion mask (e.g., Assafin et al. 2010), and the construction of an empirical model that takes into consideration effects due to the atmosphere and the instrument. This last one was the solution adopted here to correct for the FDP.

Such a solution (hereafter C0) is based on the model developed by Bernstein et al. (2017) and was the first step toward the determination of positions. C0 provides corrections for the instrumental distortion effects including color terms from the optics, delivering an astrometric solution for the DECam with rms errors below 10 mas. This astrometric solution is obtained from a parametric model that considers the celestial coordinates of an object and its respective pixel coordinates along with a set of observing circumstances (e.g., object's color, exposure time, filter), profiting from internal comparisons of around 40 million high signal-to-noise ratio measurements of stellar images. A first degree polynomial can be subsequently used to relate CCD and gnomonic coordinates of reference stars, providing reliable solutions from fields with low star densities. Observed positions will be sent to the MPC.

The refinement of orbits is obtained with the code numerical integration of the motion of an asteroid (NIMA; Desmars et al. 2015). NIMA starts from existing orbital parameters and then iteratively corrects the state vector from the differences between observations and computed positions through least squares. NIMA adopts a specific weighing scheme that takes into account the estimated precision of each position (σ_i), depending on the observatory and stellar catalog used as reference to determine the observed positions and the number of observations obtained during the same night in the same observatory (N_i) as well as a possible bias due to the observatory (b_i). The final variance of observation i is given by ${\omega }_{i}^{2}={N}_{i}{b}_{i}^{2}+{\sigma }_{i}^{2}$. As a consequence, the weight is given by $1/{\omega }_{i}^{2}$. This weighing scheme is particularly relevant when we consider old epoch positions that do not use the Gaia catalog as a reference.

The values used in the NIMA weighing scheme are described in Desmars et al. (2015) and were consolidated before the release of the astrometric data from the Gaia mission. Therefore, the code was improved to profit from the DES observations and from the Gaia releases. In this way, we have adopted ${\sigma }_{i}={b}_{i}=0\buildrel{\prime\prime}\over{.} 125$ for observations reduced with the Gaia DR1 and σ_i = b_i = 01 for observations reduced with the Gaia DR2. We emphasize that the latter is the case of DES observations presented here.

It is possible to run NIMA, with the help of few scripts, in an unsupervised way so that it is suitable for a pipeline. One of its outputs is the object ephemeris in a format (bsp—binary Spacecraft and Planet Kernel) that can be readily used by the SPICE/NAIF tools (Acton 1996; Acton et al. 2018) to derive the state vector of a given body at any time.

The prediction of an occultation event is given by prediction maps that show where and when, on the Earth, such an event can be observed. This involves the knowledge of the Earth's position in space, the geocentric ephemeris of the occulting body, and a set of stellar positions in the neighborhoods of the sky path of the occulting object as seen by a geocentric observer (see details in Assafin et al. 2010). Note that, with the astrometry from Gaia, the uncertainties in predictions rest completely upon the accuracy of the ephemerides.

A dedicated website, as presented in the next section, provides these occultations maps where many events occurring during daylight are also shown. This is done so that we are aware of even those ones that can be observed near the Earth terminator.

The determination of positions of TNOs and Centaurs from the DES images was subject to at least three constraints. The first one is that the ephemeris position of the target falls inside a box size of 4'' × 4'' centered on its observational counterpart. The second is an iterative 3σ filtering on the offsets, as obtained from the differences between observations and a reference ephemeris, to eliminate outliers. The third constraint is based on a brief inspection of the magnitudes as obtained from the DES database for each filter. Differences larger than Δ = 0.9 mag between the brightest and faintest values in each filter, when multiple measurements were available, were investigated and eventually eliminated. This value of Δ takes into account a maximum variation of σ_S = 0.15 (absolute value) in the magnitude due to the object's rotation, a maximum uncertainty of σ_M = 0.1 in the observed magnitude, and a maximum variation of σ_P = 0.25 (absolute value) in the observed magnitude due to the phase angle. In other words, ${\rm{\Delta }}\sim 3\times \sqrt{{{\sigma }_{S}}^{2}+{{\sigma }_{M}}^{2}+{{\sigma }_{P}}^{2}}$.

These constraints were expected to provide a reliable identification of the solar system objects in the images with minimum elimination of good data. However, a preliminary orbit fitting of some objects still showed the presence of real outliers (misidentifications). To solve this, a fourth filter was applied to our data and affected mostly those sources whose ephemerides presented large uncertainties (extension and doubtful sources; see Section 4.2). This filter has as an input the offsets that remained from the application of the previous filters and works as follows.

First, a mean (m₀) and a standard deviation (s₀) are obtained from a sigma-clipping iterative process, where σ is a low value (1.5 in the present case). The adopted standard deviation is the largest value between 10 mas and s₀ as given by the sigma-clipping iterations. Then, any offset within N times the adopted standard deviation from the mean was kept. Most frequently, N = 5 was used.

As a result from this process, misidentifications of TNOs and Centaurs from the images were reduced to a minimum, although real outliers can still be found mostly among the doubtful sources.

Our results in astrometry are organized in Tables 5–7 (Appendix), and the respective source distribution in the sky can be seen in Figure 1.

Zoom In Zoom Out Reset image size

Table 5 (main) considers those sources for which the 1σ ephemeris uncertainty (σ_E) in both R.A./decl. is smaller than or equal to 2'' for TNOs and Centaurs and the number of observations (N) is greater than or equal to 3. Table 6 (extension) considers those sources for which the ephemeris uncertainty is 2'' < σ_E ≤ 12'' and ${\text{}}N\geqslant 5$. Table 7 (doubtful) considers the remaining sources. All of the ephemeris uncertainties used in these tables were obtained from JPL on 2018 April 27 and are referred to 2014 January 1 at 0 hr UTC. Note that these uncertainties are given as they appear in their respective ephemerides, that is, 3σ values.

Note that the choice of the 4'' square box, although somewhat arbitrary, is a good compromise within the organization of our results to keep reliable source identifications in Tables 5 and 6, most of them in Table 5. Few objects would have moved from Tables 6 to 5 if we had opted, for instance, for a 5'' or 6'' square box. This is so because objects in Table 6 frequently have at least one coordinate (R.A./decl.) with a large ephemeris uncertainty when compared to the respective columns in Table 5. In any case, as shown later, objects in Table 6 are also a contribution to orbit refinement.

4.2.1. The Extension Table: Rationale

Most (90%) of the CCD images treated here have less than 1100 sources. Knowing that the size of one CCD in the DECam is ∼9' × 18', we can consider that there is one field object,⁴² on average, inside a box of 24'' × 24''. In this way, it is expected that a box of this size centered on the ephemeris (calculated) position of an object in Table 6 contains the respective observed position and a field star. If any of them fall inside a box of 4'' × 4'' around the ephemeris position, then this observed position is flagged as an eligible target. If not eliminated by the other steps of the filtering process, then this observed position is selected to refine the respective orbit.

We adopted the number five as the minimum number of filtered (see Section 4.1) selected positions that an object with an ephemeris uncertainty of 2'' < σ_E ≤ 12'' must have to appear in the extension table. Orbits for the objects in this table do not have the same quality as those for objects in Table 5. However, as illustrated by Figure 2 (compare it to Figure 6 panel (a), shown later in the text), the five or more positions of each object in that table are a relevant contribution to the refinement of their respective orbits.

Zoom In Zoom Out Reset image size

In the astrometric analysis of these images, it is interesting to introduce here the concept of limiting magnitude, as presented by Neilsen et al. (2015) and also discussed by Morganson et al.(2018).

The limiting magnitude is that at which the magnitude of a star is measured with an uncertainty of 0.1 mag. It can be shown to be related to a quantity τ by

$$\begin{eqnarray}&&{m}_{\mathrm{lim}}={m}_{0}+1.25\mathrm{log}\tau ,\end{eqnarray} \tag{ 1 }$$

where τ is a scaling factor to the actual exposure time (given by the image header). As a consequence, an effective exposure time can be defined as τ × nominal exposure time. The τ quantity and the limiting magnitude, therefore, can be used as a quality parameter for a given image. In order to determine the limiting magnitude in the r-band shown in Figures 3 and 4, the value m₀ = 23.1 was taken from Neilsen et al. (2015) and the values of τ were obtained directly from the DES database for each CCD (Morganson et al. 2018).

Zoom In Zoom Out Reset image size

The accuracy of the observations for the objects presented in Tables 5–6 (columns 5 and 6) is illustrated by Figure 3, where the average limiting magnitude (22.9) in the r-band (dashed line) sets a rough limit in the upper panels from which the uncertainties become larger, mainly when the number of observations is low. It also shows that the sources with a large number (hundreds) of observations have magnitudes that are close to or fainter than this limiting magnitude.

Two relevant features are shown by Figure 3. First, the lower panels show that, even in frames with the shortest exposure time (90 s), we detect sources with r as faint as ∼24.0 with a quality that is comparable to those from frames with an exposure time of 400 s thanks to the excellent quality of the images. It is worth mentioning that the faintest objects are more than 1 mag fainter than the average limiting magnitude in the r-band. Second, it is also possible to note that the range of uncertainties in R.A. is wider than that in decl. This feature most probably results from the fact that the ephemeris uncertainties (columns 3 and 4, Tables 5–7) are, on average, larger in R.A. than in decl., since we do not verify such a large difference between our measurements in R.A. and decl. as discussed below.

The standard deviations in Tables 5–7 (columns 5 and 6), obtained from the differences between the observed positions and those from the respective JPL ephemeris, is a common way to express the positional accuracy of solar system targets. These differences vary as a function of time so that, in the present study, the standard deviations provided by these columns numerically overestimate the internal accuracy (or repeatability) of the astrometric measurements.

A second astrometric empirical model (hereafter C1), also developed by the DES collaboration and based on Bernstein et al. (2017), provides improved astrometric solutions for all of the good-quality wide-survey DES exposures for years one through four of the survey. From C1, instrumental solutions are believed accurate to smaller than 3 mas rms per coordinate (see Bernstein et al. 2017). As a consequence, every DES astrometric measurement will be limited by the stochastic atmospheric distortions, typically ∼10 mas rms in a single exposure within this solution. Note that, as compared to C0, C1 is available to a smaller set of DES exposures.

We compared the positions we determined for TNOs and Centaurs to all those ones resulting from C1. This comparison is summarized in Table 2, where all of the differences we found between our results and those from C1 were kept. It is important to note, however, that C1 does not provide a solution for all CCDs. We stress that C1 is only used to provide a more realistic estimate of the internal accuracy of our measurements as well as a comparison between our positions and those from the most recent astrometric empirical model developed by the DES collaboration. C1 does not participate in any of the astrometric determinations provided here.

The standard deviations shown in Table 2 (columns 4 and 5) are a more reliable estimate of the internal accuracy of our measurements, as compared to those obtained in Tables 5–6. This internal accuracy is given by the standard deviation of the measurements, not of the mean. Therefore, the small systematics between both solutions (columns 2 and 3) cannot be considered negligible. Part of them, at least, may be explained by the fact that the empirical model is based on the Gaia Data Release 1 (Gaia DR1; Lindegren et al. 2016). It is also worth mentioning that, when our positions are referred to the Gaia DR1 (that is, the Gaia DR1 is used as reference for astrometry), the values of these standard deviations in R.A. and decl. are more similar to each other.

On the other hand, a realistic estimate of the final positional accuracy of the targets (or how accurate their equatorial coordinates are given in the International Celestial Reference Frame (Ma et al. 1998)) can be obtained from the root mean square (rms) of the reference stars, as given by the differences between their observed and catalog positions, and the precision in the determination of the object's centroid. The latter, as well as the rms of the reference stars for different filters and magnitude ranges, are provided by Table 3. In this context, this final accuracy to both equatorial coordinates is obtained, at the 1σ level, from the quantity

$$\begin{eqnarray}&&{\sigma }_{F}=\sqrt{{{\sigma }_{C}}^{2}+{{\sigma }_{R}}^{2}},\end{eqnarray} \tag{ 2 }$$

where σ_C is the uncertainty in the determination of the objects' centroid and σ_R is the rms of the reference stars. For the r filter, for instance, 12 mas < σ_F < 20 mas.

When dealing with solar system objects, the mid-exposure time (time of the shutter opening plus half of the exposure time) is of particular importance. DECam has a shutter that takes a while (about 1 s) to cross the focal plane, so the actual mean of the exposed time depends on the position in the focal plane. To compensate for this feature, the mid-exposure time was obtained by adding

$$\begin{eqnarray}&&0.5\times (\mathrm{exposure}\,\mathrm{time}+1.05\,{\rm{s}})\end{eqnarray} \tag{ 3 }$$

to the value of the Modified Julian Date (MJD) as read from the image headers (see Flaugher et al. 2015). This becomes particularly relevant when dealing with objects in the inner solar system.

In Figure 4 we show the detection efficiency as measured by the number of observed positions divided by the number of images for a given object. This figure has contributions from all of the images matched to objects in Tables 5 and 6, including those taken under non-photometric sky. This efficiency justifies the more favorable detection statistics shown in Table 4 (column 3) as compared to the initial estimates given by Table 1. It is true that this latter, as opposed to Table 4, considers only those objects for which the uncertainty in the ephemeris is ≤2''. However, Table 5 alone, with 114 entries, corroborates this better performance.

Zoom In Zoom Out Reset image size

Orbit refinement is a straightforward process with NIMA, once positions are obtained. One ephemeris (bsp format) file is provided for each of the 177 TNOs and each of the 25 Centaurs (see Table 4), from which the J2000 equatorial heliocentric state vector of each body at any time⁴³ can be obtained with the help of the SPICE/NAIF tools.

As far as stellar occultations are concerned, it is enough to be aware of an occultation event one or two years in advance so that the object's ephemeris can be more intensively refined, if necessary, and the respective observation missions for the occultation can be organized. In this way, these ephemerides should be sufficiently accurate for 1–2 yr after the most recent observations and constant updates must be provided. Ideally, we consider an ephemeris to be sufficiently accurate when its 1σ uncertainty is smaller than the angular size of the respective occulting body and very few objects—(10199) Chariklo and Pluto among them—profit from such ephemerides. Observations like those from the DECam are invaluable to change this scenario.

One disadvantage of the bsp files is that they do not carry information on uncertainties. Our dedicated website provides an orbit quality table in which uncertainties are given in steps of six months to each target. These uncertainties vary from few to hundreds of milliarcseconds, depending mainly on the astrometric quality of the current epoch of observations.

The result of an ephemeris refinement is illustrated by Figures 2 (object from Table 6) and 6 panel (a) (object from Table 5). They compare the refined orbit with its counterpart from JPL and show the uncertainty of the refined orbit along with the recently observed positions of the respective solar system body. Among others, it helps to have a first idea of the work still needed to reach suitable uncertainties for successful predictions.

The waving pattern seen in Figure 6 panel (a) is a common feature. It is a consequence of the different heliocentric distances of the solar system bodies as determined from NIMA and JPL combined with the Earth's motion around the Sun. Deep sky surveys like the DES also play a relevant role to improve the determination of these distances by providing observations at different phase angles.

Orbits determined in this work can be found from  http://lesia.obspm.fr/lucky-star/des/nima. For each object, a text file lists the positions determined here as well as the respective observational history from AstDys⁴⁴ (MPC, if the object is not found in the AstDys) that were used to determine the orbit. The 1σ orbit uncertainty (${\sigma }_{\alpha }\cos \delta$ and σ_δ) is given for a period of two years in steps of six months from the last observation. Orbits themselves are available in the bsp format. Details on the pages content are provided in a README file.

One important feature of surveys like DES is the possibility to provide a better insight on dynamical theories as the number of objects on which such theories may be employable increase through new discoveries. This is illustrated with the help of Figure 5.

Zoom In Zoom Out Reset image size

Considering explicitly the osculating elements, it is interesting to note that the MPC lists, to date, 48 objects with q > 40 au and a > 50 au. They constitute a conspicuous population of detached objects, for which mechanisms capable of increasing their perihelia is a subject of interest. Three of these—2013 VD24, 2014 QR441, and 2005 TB190—were observed by the DES, the first two being discovered by the survey. All of them are shown in Figure 5.

Gomes (2011) showed that there is a path between a scattering particle, induced by the migration of the giant planets, and the stable orbit similar to that of 2004 XR190 (black square in Figure 5, object not observed by the DES). This path results from a combination of Neptune's migration and mean motion resonance (MMR) plus Kozai resonance. One of the features of this dynamical path is that the new stable orbits escape the MMR of Neptune. The discovery of more objects by deep sky surveys with q > 40 au and a > 50 au may help to confirm this dynamical path.

2013 VD24 (close to the 5:2 resonance) and 2014 QR441 (close to the 7:2 resonance) are potentially among these objects. Numerical integrations of the equations of motion are necessary to check if they are not trapped in the resonances indicated in Figure 5. A more detailed study is ongoing.

A dedicated website also provides access to occultation prediction maps for the TNOs and Centaurs in this work.

These maps can be found at  http://lesia.obspm.fr/lucky-star/des/predictions along with a link to specific ongoing campaigns where intense astrometric efforts are done to orbit improvement. These specific campaigns are those for which worldwide alerts are sent. The basic pieces of information given by the maps are as illustrated by Figure 6(b).

Zoom In Zoom Out Reset image size

Prediction maps, plots with ephemeris uncertainties, as well as the respective ephemerides (bsp files) are available and are constantly updated at the websites mentioned earlier in the text.