THE CANADA–FRANCE ECLIPTIC PLANE SURVEY—L3 DATA RELEASE: THE ORBITAL STRUCTURE OF THE KUIPER BELT^*

Jewitt et al. (1996) provided one of the first estimates of the size and shape of the Kuiper Belt. In this early work they considered the biases against detections of high-inclination objects and determined that the intrinsic width of the Kuiper Belt is likely around ∼30° FWHM, much broader than anticipated. As surveys of the Kuiper Belt have continued, observers are becoming more aware that a correct assessment of all observational bias is critical if one is to correctly measure the intrinsic distribution of material in the distant solar system. Over the past decade, a number of surveys (Jewitt et al. 1998; Larsen et al. 2001; Trujillo et al. 2001; Gladman et al. 2001; Allen et al. 2002; Jones et al. 2006) have further refined these statistics and Elliot et al. (2005) provide the survey of the Kuiper Belt with the largest number of detected and tracked objects. A compilation of the detection rates and subpopulation estimates from these KBO surveys can be found in the Kavelaars et al. (2008). These previous surveys have provided exciting insights into this region of the solar system.

The precise determination of a KBO orbit often requires an observational arc of many years in length. The long time commitment and frequent observations often lead to observational biases that are difficult to account for (see Jones et al. 2006; Kavelaars et al. 2008, for further details). A survey which desired to provide a detailed and accurate description of the underlying Kuiper Belt population must pay very close attention to survey characterization. In the Canada-France-Ecliptic Plane Survey (CFEPS) we have taken great pains to characterize and track the survey detections so that one can make robust statements regarding the underlying orbital distribution in the Kuiper Belt. In the CFEPS project we have tracked ∼85% of our characterized detections (see Section 3.1), to-date this is the highest fraction of tracked objects for any KBO survey. Though the sample we report here is only ∼5% of the world inventory of KBOs with known orbits, we are able to make a number of interesting and new statements regarding the intrinsic orbital structure of the Kuiper Belt region because of our careful attention to the survey characteristics.

In Sections 2 and 3, we describe the observations and characterization of the CFEPS L3 KBO sample. In Section 4, we explain the use of our survey simulator and statistical methods and in Section 5 we describe our parametrized models of the KBO orbit distribution and provide estimates of the total number of objects in each subpopulation modeled.

Observations where acquired as part of a queue-scheduled survey at the CFHT. Our fields were chosen to lie within 2 degrees of the ecliptic plane and the CFHT QSO (queue service observer) imaged those fields that were at opposition whenever our priority on the facility was higher than that of competing programs. This "pressure" balancing approach to the initial discovery observations resulted in about one "block" of observations per CFHT queue observing session (nominally the 7 days leading and trailing the new moon). As a result, the CFEPS sky coverage will span nearly the entire R.A. range with breaks near the galactic plane crossings (where we did not request observations) and at locations of very high queue pressure (typically May and October). The sky coverage of our L3 release is 94 square degrees (see Table 1)

The LS-VW observations are acquired as part of the CFHT queue system. For each field observed on a photometric night the CFHT QSO provides calibrated images using their ELIXIR processing software (Magnier & Cuillandre 2004). Our photometry is reported in the AB zero-point system defined by ELIXIR combined with an average color term for the camera run of the observation. All of the CFEPS discovery observations were acquired in photometric conditions. Using differential aperture photometry we measure the total flux for each of our detected objects on all photometric nights and these fluxes are reported in Table 3.

The photometric calibration of the discovery triplets to a common reference frame and determination of our detection efficiency is required for our survey simulator analysis. To characterize our detection efficiency we inserted artificial sources, moving on the sky with rates typical of KBOs, into our CCD images and then determine our efficiency of recovering these artificial sources by running the images through MOP, the same detection pipeline as used in discovery processing. Most of our discovery observations were made using the g' filter at CFHT, with some made using the r' filter. For fields observed in r' we shifted the limits to a nominal g' value using a color of (g' − r') = 0.73, typical of the (g' − r') in our sample (see Table 3). Our detection and characterization processes are fully described in Petit et al. (2004) and Jones et al. (2006).

Those objects in each block that have a magnitude above the block's 40% detection probability are considered to be part of the L3 characterized sample. Our pipeline detection has a human-operator confirmation step for real objects; the detection efficiency, however, is determined later by implanting hundreds of artificial objects per CCD and these sources are found without the operator confirmation step. Determining the detection efficiency of the artificial objects using a human-operator was too time consuming for the ∼10⁵ artificial objects used to characterize the detection efficiency of the L3 fields. Below ∼40% the detection efficiencies determined by human operators and MOP diverge (Petit et al. 2004); since characterization is critical to the CFEPS goals, we are unable to utilize the sample faintward of the measured 40% detection efficiency level The characterized L3 sample is those 57 objects (of the 73 discovered) whose flux at detection is above the 40% detection threshold for their block; this represents roughly 80% of the discovery sample, consistent with the shape of the KBO luminosity function and typical decays in detection efficiency. These 57 KBOs make up our characterized sample (see Table 3 for a list of KBOs considered part of this sample).

Our discovery and tracking observations were made using short exposures designed to maximize the efficiency of detection and tracking of the KBOs in the field. These observations do not provide the high-precision flux measurements necessary for possible classification using broadband-colors of KBOs and we do not comment on this aspect of the L3 CFEPS sample.

Tracking during the first opposition was done using the built-in follow up of the LS-VW project. Subsequent tracking, over the next three oppositions, occurred at a variety of facilities, including CFHT. The observational efforts are summarized in Table 4. In spring 2006 the CFEPS project made an initial data release of the complete observing record for the L3 objects (before all the refinement observations for all objects were complete). The L3 release was reported to the MPC (Gladman et al. 2006; Kavelaars et al. 2006a, 2006b) and additional followup that has occurred since the 2006 release has also been reported to the MPC. Detailed data for the L3 objects can be found on the CFEPS web site at http://www.cfeps.net where updated electronic tables of all of the information in this paper can be obtained. The correspondence between CFEPS internal designations and MPC designation can be determined using Tables 5 and 6. The tracking observations provide sufficient information to allow reliable orbits to be determined such that ephemeris errors are smaller then a few tens of arc-seconds over the next five years.

Of the 73(57) KBOs in our L3 discovery (characterized) sample 55 (48) have been tracked through 3 oppositions or more (i.e., not lost) and their orbits are now known to a precision of Δa/a < 0.1% and can be reliably classified into orbital subpopulations (see below). The very high tracking fraction (fraction of characterized sample objects tracked) for the CFEPS project (84% of KBOs in our L3 characterized sample have been tracked to 3 or more oppositions) is a testament to careful planning of our survey strategy. The initial tracking of KBOs discovered by CFEPS is through blind return to the discovery fields to ensure that there is no orbital dependency in the tracked fraction. We do find, however, that the tracked fraction is a function of the magnitude of the KBO and have characterized this bias. For the L3 fields we find

$$\begin{equation*} f_t(m_g) =\left\lbrace \begin{array}{@{}ll} 1.0 & (m_g < 22.8) \\ \ms 1.0 -0.25(m_g - 22.8) & (m_g > 22.8) \end{array}\right. \end{equation*}$$

where f_t is the tracked fraction. The tracked fraction remains well above 50% down to the characterized limit of the survey blocks. We have also examined the magnitude dependence of the tracked fraction of our pre-survey discoveries (Jones et al. 2006) and find

$$\begin{equation*} f_t(m_g) =\left\lbrace \begin{array}{@{}ll} 1.0 & (m_g < 24.5)\\ \ms 0.0 & (m_g > 24.5). \end{array}\right. \end{equation*}$$

Our pre-survey discovery observations were made in the "R" band and we have transformed our pre-survey limits to g_AB, for use in our survey simulator, using a constant color offset of (g_AB − R) = 0.8 (Hainaut & Delsanti 2002).

Kuiper Belt objects are subdivided into broad dynamical classes based on their orbital elements and dynamical behavior. Using the procedure fully described in Gladman et al. (2008) (similar to that described in Chiang et al. 2003), we report the L3 sample classifications as of March 2008 (including all refinement observations to date). Briefly, based on all available astrometry we determine the distribution of elements consistent with the observation for each object. A best fit and two extremal "clones" consistent with the observations are integrated for 10⁷ years and their orbital evolution is used to provide the classification. If the three orbital evolutions provide the same answer, the classification is deemed to be "secure"; if not the classification is insecure and further refinement observations are needed. The results of this exercise are given in Tables 5 and 6. Based on current knowledge of the L3 sample seven objects remain insecure (even though these have 5 opposition observational arcs!), and all of these are due to their proximity to a resonance border where the remaining astrometric uncertainty makes it unclear if the object is actually resonant. For these "insecure" objects we list them in the category shown by two of the three clones and note the other possible classification (see Table 5 and 6).

We find that 33 (60%) of our L3 tracked sample are in the classical belt, 18 (33%) objects are in a mean-motion resonance with Neptune, 8 (15%) of which are plutinos, the remaining sample consists of 2 (4%) objects on scattering orbits and 2 (4%) in the detached population. Given these orbital classification in our observed population we now determine a process for the parametrization the intrinsic population of the Kuiper Belt implied by our observations; the wildly different detectability of these populations in a flux-limited survey implies that the "true" population ratios may be different from the observed ones.

Our goal is to determine the intrinsic population of KBOs that is consistent with our observations. A standard approach would be to use these observations to constrain some reasonably agreed upon parametrized model of the orbital distribution of KBOs. There is, unfortunately, no such model; the Kuiper Belt has arrived in the current state by a complex sequence of events related to planetesimal accretion, gravitational erosion, probably planet migration, possibly resonance sweeping, and perhaps sculpting by large rogue objects. The only models we have today to describe the KBO population are the results of long term numerical integration, with limited resolution of all but the most numerous dynamical subpopulations, and no possibility to vary parameters to fit the observations. To circumvent this limitation here, we rely on ad-hoc parametrization of the various dynamical subpopulations. Our hope is that these parametrization will guide future dynamical modeling.

We proceed by examining both the orbital distribution of the L3+Pre sample as well as the distribution of orbital elements as reported by the MPC. By examining the observed orbit distributions of each of the Kuiper Belt subpopulations we can then make some initial guesses as to likely parametrization of the intrinsic population. Then, through extensive use of the survey simulator, we pare away those intrinsic models that are inconsistent with our characterized sample. Regardless of how well we might reproduce the observations, one should remember that the parametrization used here is completely ad hoc, and should not be viewed as resulting from any specific model of the formation and evolution of the KBO population.

Although we make some progress here, a reliable description of the internal orbital structure (of say the insides of the 3:2 resonance) will require factors of many to orders of magnitude larger sample sizes than we have available from our calibrated L3+Pre survey catalog. We anticipate that the complete CFEPS sample (L3, L4, and L5), which will contain roughly 200 well-characterized orbits, will be available before the end of 2009.

Once we have decided on a model, given either by a set of parametrized distributions, or as the result of a full-fledged numerical simulation of the formation and evolution of the outer solar system, we pass the model orbit catalog into our survey simulator and obtain a set of simulated detections. These detections represent those KBOs our survey would detect if the input model is a reasonable representation of the intrinsic KBO population. Given these simulated detections, we compare their orbital-element distribution with those of the real objects.

Ideally we would compare the multidimensional distribution of all orbital elements simultaneously; given the small sample sizes being considered, a multidimensional approach is not warranted. Instead we rely on a series of one-dimensional distribution tests, using either the Anderson–Darling (AD) statistic (for elements a, e, i, q, Δ) or the Kuiper-modified Kolmogorov–Smirnov (KKS) D-statistic for cyclic variables (elements M, N, ω). The AD test statistic is chosen since this statistic is more responsive to the tails of a distribution than the more traditional Kolmogorov–Smirnov test statistic.

Using the survey simulator we generate many thousands of simulated detections for a given intrinsic model and use these large simulated populations to bootstrap the relevant AD and KKS statistics¹². We then proceed as usual for these tests by rejecting models at a given confidence level of X% when the probability of finding an AD or KKS statistic value larger than the one observed is smaller than (100 − X)/100, i.e., the observed statistic is larger than X% of our bootstrap cases for that model.

Our AD and KKS statistics are computed for one dimensional distributions, but we assess our models based on the set of orbital-element (a, e, i, q) and observation circumstance (r, m, M) distributions. For a given model, we compute the rejection confidence level for each of the measured statistics and retain the largest rejection confidence as the rejection confidence level for that model as a whole. That is, we take the model orbit-element or observation circumstance distribution that least matches the observed distribution as the indicator of quality for a particular input model, determining which orbit distribution least matches the data via the AD and KKS statistics. For example, consider only the a, e, and i distributions for some hypothetical Kuiper Belt model, $\cal H$: if the a distribution AD statistic is larger than 60% of the bootstrap cases (i.e., the model cannot be rejected), the e distribution AD statistic is larger than 95% of cases (i.e., marginally rejected) and the i distribution AD statistic is larger than 99.9% of the bootstrap values (i.e., the i distribution strongly rejects the model), then we reject model ${\cal H}$ at the 99.9% confidence level, based on the rejection of the i distribution. For a model distribution to be considered acceptable, all the simulator-observed model orbital-element distributions must be consistent with the survey-observed orbital-element distributions, i.e., none of the element distributions may reject the model.

The observed element and orbit element distribution for Kuiper Belt objects are correlated, most importantly; the distribution of distance at detection, r, and peri-center q are coupled to the eccentricity (e), semimajor axis (a), and luminosity (m) distributions. In our "worst case" rejection process we sometimes reject a model based on its r or q distribution, even though the individual a, e, m distributions may be marginally acceptable. In these cases the r and q distributions reveal that, although the distributions of the orbit elements are acceptable, their coupled distribution (which a survey is most sensitive to via the observed element distributions) is not acceptable. The use of coupled orbit element distributions raises the specter of "double testing": sometimes we reject based on a or e and other times on the basis of r or q and these are, in-fact, measuring the same quantities with different couplings to the survey circumstance. However, we use a and e rejections to specifically tune the limits of our a an e distributions, and when examining the global acceptability of the model we then revert to the q and r distributions. In this way we avoid "double testing."

Once a set of model orbital element distributions have been determined, essentially creating the parametrized Kuiper Belt, we can use these distribution functions to determine the total absolute population of the Kuiper Belt, subject to our model. We report our determination of the number of objects with absolute magnitudes H_g < 10.0 for each subpopulation. The CFEPS survey did not detect any objects with H_g>10 and thus this value represents the global limit on our size sensitivity.

For each dynamical subpopulation, we draw objects from our model population until the survey simulator detects the same number of objects as in the real L3+Pre sample for that dynamical class. We repeat this process 2000 times and for each of these realizations we count the number of objects in the simulated population needed to account for the number of objects in the L3+Pre sample (see Figure 3).

Zoom In Zoom Out Reset image size

We prefer to quote N(H_g < 10.0) as this number is independent of albedo. To convert our estimates to the number of KBOs larger than 100 km diameter one must assume an albedo and determined the measured luminosity for the Sun in a particular waveband. We denoted this albedo/luminosity dependent population estimate as N(D_p>100 km). For objects that are visible only via reflected sunlight,

$$\begin{equation} {\rm H}_x = m_{\odot }(x) + 42.38 - 2.5\log _{10}(p_x D^2) \end{equation} \tag{ 1 }$$

where p_x is the x-band albedo of the TNO and m_☉(x) is the Sun's AB x-band magnitude. For the g' filter, m_☉(g') = −26.47 (Willmer 2008). Scaling off a reference albedo of 0.05 and a reference diameter of 100 km yields:

$$\begin{equation} {\rm H}_g = 9.16 - 2.5 \log _{10}\left(\frac{p_g}{0.05}\left(\frac{D}{100\,{\rm km}}\right)^2\right). \end{equation} \tag{ 2 }$$

The CFEPS observations are only sensitive to the brightest members of the KBO population and the luminosity function of this population is well represented by a uniform power-law size distribution (Petit et al. 2006). Starting from our N(H_g < 10) estimates (which corresponds to D_p ≃ 70 km for p_g = 0.05), the shift to a larger diameter is performed by a multiplication by 10^α×ΔH where ΔH = 10 − H(D_p = 100 km) is the difference between H = 10 and the magnitude corresponding to a diameter D_p = 100 km. We adopt α = 0.72 for the slope of the magnitude distribution for all our modeling (Petit et al. 2006), leading to N(D_p>100 km) = 0.25*N(H_g < 10). Many previous estimates of N(D_p>100 km) were made assuming p ≃ 0.04, in which case N(D_p>100 km) = 0.37*N(H_g < 10.0). The magnitude of the corrections due to: albedo-dependent color (small), the slope of the magnitude distribution (moderate), and absolute albedo (large) render such comparisons difficult. Using N(H_g < 10) escapes from the albedo uncertainties and is only weakly dependent on our assumed value of α.

The main classical Kuiper Belt is well described as a population with a two component inclination distribution (see Levison & Stern 2001; Brown 2001; Elliot et al. 2005) and these two components also appear to have distinct surface properties (Tegler & Romanishin 1998, 2000; Trujillo & Brown 2002). The classical belt population also appears to have a limited radial extent, initially Dones (1997) noted that, even with the small handful of objects then known, either the luminosity function of Kuiper Belt objects is quite steep or the population does not extend beyond a ∼ 50 AU. Jewitt et al. (1998) found that a classical Kuiper Belt with a flat radial distribution beyond 50 AU was inconsistent with their observations while Gladman et al. (2001) concluded that the observed radial distribution of KBOs in their survey was consistent with a radii distribution that declines beyond 50 AU. Further observations by Allen et al. (2001) demonstrated a lack of objects on circular orbits beyond 50 AU, and if the classical belt does extended beyond this distance then the surface density must drop dramatically. Additional strong evidence for an outer-edge to the classical Kuiper Belt came from Trujillo et al. (2001) and Trujillo & Brown (2001) and the more recent analysis of Hahn & Malhotra (2005) indicates that the lack of objects on circular orbits with a>50 AU is consistent with a classical Kuiper Belt that is truncated at a ∼ 45 AU.

Jones et al. (2006) demonstrated that a filled phase-space model of the main classical belt, one that filled the available phase space down to a minimum perihelion of q_min = 38 AU, was not rejected by their small, but well-characterized, sample. We passed this same model¹³ through the survey simulator, configured for the L3+Pre sample. We now reject this model at more than 99.9% confidence, dotted line in panel A of Figure 4. (Figure 4 presents the a, e, i CDFs of the main classical belt model distributions as observed through our survey simulator.) The orbital information available within the CFEPS project has now reached the quality level where more precise modeling of the internal structure of the main classical belt can be tested.

Zoom In Zoom Out Reset image size

5.1.1. An Initial Inclination Distribution

Recognizing that the observed inclination distribution of classical Kuiper Belt objects follows the shape of a two component Gaussian, Brown suggested an intrinsic model distribution of

$$\begin{equation} f_2 (i) = \sin i \left(a_m \exp \left(\frac{ -i^2 }{2 \sigma _c^2} \right) + (1-a_m) \exp \left(\frac{ -i^2 }{2 \sigma _h^2} \right) \right) \end{equation} \tag{ 3 }$$

where σ_c is the width of the "cold-component," σ_h is that of the "hot-component," and a_m is a weighting factor, not the fraction of the population that is in the cold-component: to avoid potential confusion we will quote our results below in terms of the fraction f_h of the intrinsic population in the hot-component. We initially adopt the best-fit values for the classical population as reported in Brown (2001): σ_c = 2.2, σ_h = 17.0 and f_h = .8.

To test this inclination distribution using our survey simulator approach we require a complete description of the distribution of all six orbital elements. We combined Brown's two component inclination distribution with a uniform a/e distribution, bounded by 40 < a < 47 AU and q>38 AU, and uniform mean longitude, peri-center and node distributions (solid line in panel A of Figure 4). The Brown (2001) model provides a marginally acceptable match to the inclination distribution of the L3+Pre sample of main classical belt objects. The model a and e CDFs when observed through the survey simulator, however, do not match the L3+Pre sample at greater than 99.9% confidence level and thus we reject this model. The a, e distributions of our 33 classical belt objects are not matched by a uniform a, e distribution for the main classical belt. A more complex a, e distribution is required.

5.1.2. The Semimajor Axis Distribution

5.1.3. The Eccentricity Distribution

5.1.4. Determining the Ratio of Hot to Cold Main Classical Belt Objects

We now have a parametrized model of the main classical Kuiper Belt with a set of parameters that provides a good match to the L3+Pre sample, how well can we constrain the parameters of this model? Although there are many parameters (regarding the internal orbital structure of the classical belt) that could be explored, such as the relative weighting of different semimajor axes or eccentricities, we have chosen to concentrate on the parameters f_h and σ_h which are relevant to cosmogonical models currently being discussed in the literature.

For our input model orbit distribution we generate two separate subpopulations of the main classical Kuiper Belt using the following recipe.

• 1.  
  Select if the object will be from the hot or cold-component such that f_h fraction of drawn objects are hot and generate an inclination drawn from the corresponding distribution. Note that in principle objects of any inclination could come from either component.
• 2.  
  If this is a cold-component object, choose a random semimajor axis on the interval (42.4,45.0); if hot use the interval (40.0,47.0) and if the object falls in the ν₈ (a < 42.4 AU; i < 12°) draw another hot (i, a) pair.
• 3.  
  Draw a random e on the interval [0,0.24) (objects with e ⩾ 0.24 that are on stable orbits are considered members of the detached population and not considered in this model), uniformly distributed for cold objects but distributed with probability P(e) ∝ e for hot objects (Jones et al. 2006); if this would result in a perihelion q < 39(38) AU for the cold (hot) classical Kuiper Belt objects then return to step 2.
• 4.  
  Draw random values for Ω, ω, and M, uniformly distributed between [0, 360).
• 5.  
  Draw a random value of the H magnitude from an exponential distribution with α = 0.72.

Each such drawn object is then passed through the CFEPS survey simulator to determine if it would have been detected and tracked. The orbital-element distribution for a large number of simulated tracked objects is then compared to that of the L3+Pre tracked objects that satisfy the constraints (40 < a < 47) AU and (q>38) AU. This restriction on q is needed since our current sample is too small to allow us to accurately model the complex stable phase space in the 35 < q < 38 AU zone.

We ran a grid of 2050 models using the above prescription covering a range of (f_h, σ_h) pairs, constrained to the value of σ_c = 2.2 from Brown (2001) and the a/e distributions as described above. We take as our "best fit" those model parameters (f_h and σ_h) for which the AD statistic of the i, q and r distributions is minimized. We find our best-fit models is f_h = 0.6 and σ_h = 15°. Our "confidence boundaries" are given by those parameters which produce models with i, q and r AD statistics that are smaller than 68%, 95%, and 99% of the boot-strapped values of the statistics, Figure 5. The contours of Figure 5 explore, in essence, the range of plausible f_h and σ_h parameters for our intrinsic model when tested against the L3+Pre sample.

Zoom In Zoom Out Reset image size

An examination of Figure 5 reveals that our current sample only weakly constrains the width of the main-belt's hot component. This poor constraint is a reflection of the small number of large-inclination objects in our current sample, along with the fact that a given discovery cannot be uniquely assigned to either of these two overlapping components.

If our "thick" model of Figure 5 applies, then based on our survey simulator analysis, all objects in our L3+Pre sample with i < 10° are more likely to be cold-component objects (in the large-i tail of that population) than hot-component members. Thus, in a thick hot-component model, our current sample likely contains only two objects from the hot population and these two TNOs provide only a weak constraint on the hot component's width. If one were to (unjustifiably) assign all the i>7° TNOs to the hot population, the thick model can be rejected at >99% confidence. Thus, although external evidence (from other surveys) suggests that σ_h ≃ 15°, our current sample provides only a loose confirmation of this result.

In our modeling we have not explored the relation between the hot main classical Kuiper Belt and detached component of the belt which has e>0.24 and a>50 AU. Our current statistics on the distant belt are very limited, however the a>50 AU detached component of the Kuiper Belt may actually be the extension of the hot main classical belt.

5.1.5. Population Estimates

Following the procedure described in Section 4.3, and using the best fit parameters described in the proceeding sections, we compute a population estimate for the main classical belt, giving (see Table 5):

$$\begin{eqnarray*} &&N_{\rm classical} ({\rm H}_g < 10.0) = \big(120^{+50}_{-46}\big) \times 10^3 \end{eqnarray*}$$

where the uncertainties reflect a 95% confidence limit assuming the underlying orbital model and its parameter values are correct. This population estimate is model dependent, and varies strongly with the width σ_h of the hot-component. As seen in Figure 5, increasing σ_h requires a larger fraction of hot-component objects, which increases the population needed to match the observations. This is due to the detection bias against finding hot objects in an ecliptic survey. As a comparison, Table 7 gives the population estimates derived from two other models within the 95% contour. For a "thin disk" model we take the the parameters that minimize the fractional size of the hot population, σ_h = 10° and f_h = 0.2, which produces a population estimate, for the hot-component, that is slightly more than half our best-fit model estimate. For a "thick disk" model we take the acceptable model with f_h ∼ 0.7 which gives a population estimate 2-sigma larger than the estimate based on our best-fit parameters. In all three cases the inclination distribution of the cold-component is held constant at σ_c = 22 and we note that the total cold population (which is (1 − f_h) × N(H_g < 10)) ranges between 50,000–60,000 for these models. Unsurprisingly the total size of the better-sampled cold population is relatively well determined (∼10%), while the poorly sampled hot-component contains most of the uncertainty (factors of several).

Table 7. Model Dependent Population Estimates

	N(Hg < 10) (1000 s)	Lower	Upper
Inner Classical Belt
Two-component P(i)	2	0.3	7
One-component P(i)	8	1	27
Main Classical Belt (σh, fh)
(10,0.2)	74	48	105
(15,0.6)	126	80	177
(40,0.7)	183	112	258
Plutinos
Nominal model	25	13	50
Kozai sub-pop	30	15	55
α = 0.54	12	6	23

Notes. Estimates are given for our model for each subpopulation within the Kuiper Belt. The values in the upper and lower columns are the upper and lower bounds on our 2-σ confidence region for the model-dependent population estimate. All values are in thousands of objects. The Kozai subpopulation shows how forcing half of all plutinos with i>15° to librate in the Kozai resonance causes a small increase in the population estimate.

Download table as:  ASCIITypeset image

Figure 6 shows a representation of our underlying model relative to the real L3+Pre detections and a set of simulated detections. In the classical belt's a = 45–47 AU region the survey simulator produces few detections and those that do occur have large eccentricities due to the perihelion flux bias (see also Jones et al. 2006).

Zoom In Zoom Out Reset image size

There are eight plutinos in the L3+Pre sample. The plutino population is of historical importance as the first recognized resonant subpopulation of the Kuiper Belt (Tholen et al. 1994), and is by far the most numerous in the flux-limited observational catalogs. Plutinos are forced by their resonant argument (e.g., Malhotra 1995, 1996) to come to perihelion away from Neptune (precisely ±90° for a zero libration-amplitude plutino).

Compared to the non-resonant objects, modeling the plutino population is complicated by the need to respect the longitude relations embodied in the resonant argument. For objects in 3:2 mean motion resonance with Neptune the resonant angle ϕ₃₂ has the form:

$$\begin{equation} \phi _{32} = 3 \lambda - 2 \lambda _N - \varpi \end{equation} \tag{ 4 }$$

where λ and ϖ are the object's mean longitude and longitude of perihelion and λ_N is Neptune's mean longitude.

Plutinos librating in the resonance will have their resonant angle ϕ₃₂ librate around 180° with some amplitude L₃₂. For example, Pluto's resonant angle has a libration amplitude of L₃₂ ≃ 84° (Milani et al. 1989). At maximum libration, the angle between the pericenter direction of Pluto (i.e., λ = ϖ) and Neptune is (180° − 84°)/2 = 48° or (180° + 84°)/2 = 132°; at these extrema the librating relative perihelion longitude reverses direction. This inability to approach Neptune causes the "hole" seen in the on-ecliptic projection of Figure 7.

Zoom In Zoom Out Reset image size

The resonant argument, ϕ₃₂, oscillates sinusoidally around 180°, thus passing more time at the two extrema ϕ₃₂ ± L₃₂ than near 180° itself. The observational impact of this effect is that objects with libration amplitudes of L₃₂ have a bias to being found near (180 ± L₃₂)/2 degrees of longitude away from Neptune. Flux bias weights against the discovery of plutinos far from their current pericenter locations, and so when looking at a particular longitude from Neptune one will tend to find more objects for whom the observed longitude relative to Neptune is such that their resonant argument is at the extrema. Indeed, the L3 plutinos have libration amplitudes similar to the longitudes where they were found. In practice plutinos have L₃₂ < 130°, and so no plutino can be at pericenter in fields that are close to the direction (or antidirection) of Neptune, explaining why the Presurvey, L3q and L3y blocks yielded no plutinos (see Table 1 and Figure 7).

Because of the 1/r⁴ nature of the reflected flux and the strong bias to detection at certain longitudes, correct modeling of the plutino population requires knowledge of the ecliptic longitudes and flux limits of the surveys. The change in distance from perihelion to aphelion for an e≃0.25 plutino (30–50 AU) drops its apparent magnitude by Δm = 2 magnitudes. Because of the slope of the luminosity function, this change in flux causes a drop by a factor of 10^αΔm ∼ 30 in the plutino fraction at longitudes where they are preferentially at aphelion. The dependence of the heliocentric distances at detection on the size distribution is illustrated in Figure 8; if the size distribution is flatter a larger fraction of objects above the flux limit are found at larger distance. The median detection distance varies by more than 2 AU if the power law H-magnitude distribution spans the range of values given in the literature. Therefore, one cannot determine the size distribution of the plutinos from surveys that cover a range of solar longitudes unless one knows the ecliptic longitude of the surveys where plutinos were both detected and not detected. Our characterized survey provides this information, but with only 8 plutinos in the L3+Pre sample, we cannot uniquely constrain this large parameter space since the detected r distribution also depends on the e and ϕ₃₂ distributions.

Zoom In Zoom Out Reset image size

We proceeded to find a satisfactory model of the intrinsic plutino orbital element distribution in a two-step process. First we produce hypothetical parametric representations of the intrinsic orbital-element distributions based on available observations and dynamical modeling and then we make small adjustments to these distributions until they are not rejected by the L3+Pre sample.

The a, e limits are taken as the borders of the stable resonance (Morbidelli 1997). A Gaussian e-distribution with mean of 0.21 and standard deviation of 0.06 provides a reasonable match to the observed eccentricity distribution of the plutinos listed in the 3:2 SSBN07 classification (Gladman et al. 2008) and provides an initial model for this element (i.e., we use eccentricity distribution of the MPC reported plutino detections as a proxy for the intrinsic one, with no attempt to correct for detection biases).

Brown (2001) concluded that a single component Gaussian inclination distribution of width 10⁺³₋₂ degrees satisfactorily represented that plutino population, so we begin with this 10-degree width.

The libration amplitude, L₃₂, distribution is taken to be a match to the compilation of Lykawka & Mukai (2007); a symmetric "triangle" where the most likely value is 65°, falling to zero at 0° and 130°. The resonant angle ϕ₃₂ cannot be chosen uniformly in the range (180 − L₃₂, 180 + L₃₂) since ϕ₃₂ oscillates sinusoidally in this range; accordingly, we time-weight the ϕ₃₂ oscillation when producing our model distribution of plutinos. Longitudes of node and mean anomalies were chosen randomly. Having picked a resonant angle ϕ₃₂, Equation (6) then forces the value of the plutino's longitude of perihelion ϖ.

This straightforward model produced a/e and distance at detections distributions that are rejected at better than the 99% confidence and an inclination distribution rejected at the 95% confidence level (see Figure 9). To first order the inclination distribution is not strongly coupled to the other distributions, so we proceed by first searching for models that reproduce the observed inclination distribution and then focus on determining the (a, e) and libration amplitude distributions.

Zoom In Zoom Out Reset image size

5.2.1. The Inclination Distribution of the Plutinos

5.2.2. The a/e and Libration Amplitude Distributions

5.2.3. The Plutino Size Distribution

5.2.4. Population Estimate

Our nominal model yields an estimated plutino population

$$\begin{equation} N_{\rm plutino} ({\rm H}_g < 10.0) = 25^{+25}_{-12} \times 10^3 \end{equation} \tag{ 5 }$$

(95% confidence, model dependent) as the total plutino population (Table 7).

With only eight detected plutinos, our population estimates are only good to about a factor of two due to the large parameter space which affects detectability. For example, although none of our eight L3 plutinos librate in the Kozai resonance, other high-inclination plutinos do (Lykawka & Mukai 2007). If we assume half of the plutinos with i>15° have their arguments of pericenter librate around 90 or 270 degrees with a width of 50 degrees, then we can increase the estimated population of plutinos to ∼30 000 (Table 7) without violating the lack of a Kozai population in the L3+Pre sample. This occurs because only a few of the L3 blocks probe the longitudes where these plutinos would cross the ecliptic and those that did cross our field of view would not be at pericenter and thus harder to detect. This further reinforces the need for knowledge of the survey longitudes in order to estimate the size of resonant populations.

Changing the nominal model by flattening the size distribution to α = 0.54 causes the population estimate to drop substantially to N(H_g < 10)12, 000, since the plutinos are detectable at greater distances.

Given all the model dependencies, we conclude that the plutino population is ∼25, 000 objects with H_g < 10, to a factor of 2. The nominal estimate indicates that the intrinsic plutino fraction is about 20% of the main classical belt, whose estimate is also uncertain to a factor of 2. Since the majority of Kuiper Belt surveys have been conducted near the ecliptic plane, the fact that the plutino population lacks a cold-component results in an apparent plutino fraction (relative to the classical belt) in the MPC catalog being smaller than the intrinsic one, a fact that some readers may find counter intuitive.

The only previous published plutino population estimate, based on analysis of a KBO survey (Trujillo et al. 2001) is N(D_p>100) = 1400, accurate to a factor of two. Converting our estimate to the same size limit assuming albedo of 0.04 yields an estimate of 9,000, also to a factor of two precision (neglecting the albedo uncertainty). Clearly we estimate a significantly larger plutino population; this is result is due to the fact that the plutino "bias corrections" applied previously assumed the same inclination distribution for the plutinos as for the classical belt. Lacking a cold-component the plutino population is relatively less detectable than the classical belt. We estimate that the plutino population is 40%–50% as numerous as the cold classical main belt population.

The large fraction of plutinos, compared to the main classical belt seems to indicate that Neptune did, at some point, have a slow migration phase which allowed the 3:2 resonance to sweep up a substantial population, although quantitative estimates of the ratio of plutinos to classical objects is lacking except in Hahn & Malhotra (2005). The lack of a low-inclination component in the plutino population, however, appears to be a critical problem for the concept of slow migration into a cold disk.

The inner classical belt is defined as the non-resonant, non-scattering orbits with semimajor axes between Neptune and the 3:2 resonance (Gladman et al. 2008). Duncan et al. (1995) first showed, via long-term numerical simulation, that interior to the a = 39.4 AU semimajor axis of the 3:2, there is a region of orbital stability (from a≃ 36–39 at low e), and remarked with surprise that no objects were known in that region. Petit & Gladman (2003) suggested that the lack of such objects could be due to systematic follow-up bias in early surveys, and showed a handful of such objects derived from more systematic tracking campaigns. Kavelaars et al. (2008) and Parker et al. (2007) discuss this tracking bias in more detail.

Based on current observations and numerical-stability simulations the number of inner classical belt objects is anticipated to be low since the volume of phase space available is much smaller than for the main classical belt. In Gladman et al. (2008), 16 of the 510 TNOs (3%) are classified as inner classical belt objects, which is in rough agreement with our having two such discoveries in the 63-object L3+Pre catalog (L3w06 and L3y14PD). With only two characterized detections we currently provide only a rough population estimate based on a simple phase-space model and must leave exploring the internal orbital structure of this region for the future.

Guided by the results of Duncan et al. (1995) and unpublished integrations by one of the authors (B. Gladman), we have modeled the inner classical belt with the following algorithm. The a/q elements are uniformly distributed in the region a = 37–39 AU and q>35.3 AU. We choose orbital inclinations uniformly from the same double-Gaussian inclination model as for our nominal classical main belt model, σ_c = 22, σ_h = 15°, and f_h = 0.6 (see Section 5.1), except that we exclude the region 7° < i < 22° in order to account for the strong destabilizing presence of the ν₈ secular resonance which passes through the inner classical belt at these inclinations (Knezevic et al. 1991; Duncan et al. 1995). The angles ω, Ω, and M are uniformly distributed between (0,360) degrees.

The number of inner classical belt objects is estimated, subject to our model constraints, using the procedure described in Section 4.3. We find the inner classical belt population to be (see Table 7)

$$\begin{equation*} N_{{\rm inner}} ({\rm H}_g < 10) = 2200^{+5000}_{-1900} \end{equation*}$$

The uncertainty range is the 95% confidence region and is model dependent. The large range of uncertainty is due to the presence of only two detections in the L3+Pre sample. Our two L3+Pre inner-belt detections both have low (i < 5°) inclinations. For the two component i-distribution given above, the "hole" caused by the ν₈ resonance results in a very low fraction of detected hot-component inner classical belt TNOs in an ecliptic survey. We calculate that if there is a cold-component then the cold inner-belters should make up ∼96% of the inner-belt CFEPS detections. If, however, the inner belt has only a hot population then low-i (i < 5°) KBOs still represent ∼70% of the expected detections. Under the assumption of a hot only model, (i.e., f_h = 1) the size of the intrinsic population required to match the L3+Pre sample would quadruple (Table 7) due to the lower detectability of the intrinsic population. Our 2 low-inclination detections are, of course, consistent with either of these hypothesis.

If we use the CFEPS L3+Pre survey simulations as a proxy for generic KBO surveys conducted near the ecliptic then we might expect the MPC fraction of low and high inclination inner-belt objects discovered near the ecliptic would match that of the simulator. This simulator by proxy appraoch is not too risky, since the non-resonant objects obey no phase relationship with Neptune. Of the 11 inner classical belt objects in the MPC database (Gladman et al. 2008) that were discovered within 2° of the ecliptic, 7 have inclinations of less than 5°. This 64% ratio is close to that expected for an inner-belt containing no cold-component (70%) and in (mild) conflict with the two component model (90%), leading us to suggest that the inner classical belt consists of a purely hot inclination population.

We are cautious not to draw too strong a conclusion from this examination of the MPC database since we expect that inner-belt objects with small inclinations may have been mistakenly classified as plutinos and lost, biasing the MPC against low-i inner-belters. Even with these uncertainties, the inner-belt population is unlikely to exceed 2% of the classical belt population, and may well be smaller than 1%. Reproducing the size and inclination structure of this population will be an important constraint on outer solar system formation models.

Based on the L3 sample we have determined that the cold- and hot-components of the main classical belt have distinct eccentricity and semimajor axis distributions. Our modeling of the orbit distribution indicates that the cold main classical belt is isolated in semimajor axis and likely does not extend the full distance to the 2:1 resonance; the cold belt appears to end at ∼46 AU.

Migration models (see Hahn & Malhotra 2005; Chiang et al. 2003; Gomes et al. 2004, for example) produce outer edges to the cold and hot populations by using an initially truncated disk (usually with a pre-migration edge well interior to the current 47.4 AU location of the 2:1 resonance). Many of these models have sought primarily to generate an outer edge at a = 47.4 AU. The fact that we have found different edges (in semimajor axis space, which makes for a fuzzy edge in heliocentric distance) for the hot and cold populations is a new constraint for formation models.

In some outer solar system formation model the radial extent of the cold disk arises as a consequence of the time-scale for the damping of Neptune's eccentricity (see Levison et al. 2008, for example) while the extent of Neptune's migration is controlled by initial radial extent of the planetesimal belt of the solar system (Gomes et al. 2004). Alternatively, the cold classical main Kuiper Belt could be the in situ remains of the primordial solar system and the extent of this population may be indicating the edge of primordial disk. Regardless, formation scenarios are now left to the challenge of producing a cold core whose eccentricity and semimajor axis distributions are distinct from the hot population.

The separate eccentricity and semimajor axis distribution of the hot and cold main classical belt populations add to the growing number of compelling observational differences: different colors (Tegler & Romanishin 1998), different binary fractions (Stephens & Noll 2006), different albedos (Grundy et al. 2005), and perhaps different size distributions (Bernstein et al. 2004; Elliot et al. 2005). Together these imply that the two populations had separate formation locations or processes. Based on the dynamics of the cold-component we feel it likely that this population formed in situ and are more comfortable with models that transfer the hot-component from some other region of the solar system.

The hot main classical belt population has an eccentricity distribution considerably hotter than that of the cold-component, with a semimajor axes distribution that does not exhibit the restriction of the cold-component. Intriguingly, we believe it plausible that the hot-component has the same inclination distribution in the inner classical belt, main classical belt, plutino population, and perhaps other populations not modeled here due to lack of detections.

The hot main classical belt distribution is, in our opinion, compellingly like the high-inclination population produced in Levison et al. (2008) and also somewhat similar to that of Hahn & Malhotra (2005). Testing these models quantitatively will require passing the results of the model simulations through the CFEPS simulator and comparing the resulting orbital distributions to the CFEPS L3+Pre detections. To facilitate quantitative comparison of models with the CFEPS results, the authors provide full details of the survey simulator on the CFEPS Web site (http://www.cfeps.net).

Examining the distribution of MPC inner belters lends evidence to our suggestion that the inner belt lacks a low-inclination component. The authors are not aware of any outer solar system model that predicts a change in the inclination distribution between the inner and main classical belts. Additional evidence for the hot-component membership of the inner belt could come from a census of the colors of the members of this population.

Models of resonance capture during a smooth Neptune migration (Malhotra 1993) indicate that if the large eccentricities of current plutinos are the result of "eccentricity pumping" during resonance migration, the 3:2 resonance must have swept through the inner belt. Resonance capture alone, however, does not strongly affect the inclinations of captured objects (Gomes 2003; Hahn & Malhotra 2005). Thus, if the inner belt had a low-inclination component during the putative resonance sweeping, they should have been easily captured and we should find this cold (low-i) population in the current plutino population; it seems not to be there today. Based on this reasoning, we postulate that the inner belt never had a low-inclination component or it was removed prior to the resonance sweeping migration of Neptune.

In the so-called "Nice" model (Levison et al. 2008) the low-inclination component of the main classical belt is produced by resonance overlap and diffusion in the zone between Neptune and the 2:1 resonance. This mechanism may be very effective at creating a low-i component in the inner belt or at least appears to produce an inner belt with the same inclination distribution as the main classical belt (see Figures 6, 9, and 10 in Levison et al. 2008). More detailed modeling of this scenario and additional observations are needed before a conclusion can be drawn regarding the structure of the inner belt region predicted.

Is our estimated plutino population large enough to be the source of the Jupiter-family comets (JFCs)? Levison & Duncan (1997) examined the JFC supply problem from the Kuiper Belt as a whole and Morbidelli (1997) quantitatively explored the hypothesis of the gradual gravitational erosion of the plutino population as a supply source for the JFCs and estimated that there would need to be ∼4.5 × 10⁸ objects currently in the 3:2 resonance for this population to be a major source of JFCs.

Our estimate of the plutino population, N(H_g < 10) ∼ 3 × 10⁴, must be scaled by the size-distribution relation to compare with the estimate from Morbidelli (1997). We take, as a plausible estimate, the size of a typical JFC to be D ∼ 5 km (Fernández et al. 1999). Using Equation (4) with an albedo of $p_{{g^{\prime }}}=0.05$ gives H_g ∼ 15.7 for a typical, D ∼ 5 km, JFC resulting in ΔH_g = 5.7 between the JFC size and the size limit of our population estimate; assuming a uniform power law for the plutino luminosity function with α = 0.72 implies an intrinsic population of N(d>5 km) 4 × 10⁸ plutinos, resulting in the plutinos being a plausible source for the JFC population.

Due to ambiguities in the size distribution slope, the location at which the steep size slope for large objects rolls over to a shallower slope (Bernstein et al. 2004; Petit et al. 2006; Fraser et al. 2008) and the true size of a typical JFC, these numbers should be treated with caution. For example, if JFCs are 10 km rather than 5 km, the corresponding 3:2 population of 3 × 10⁷ objects is insufficient as a supply source.

The sensitivity of these numbers to the unknown value of the typical JFC size is especially worrying, and results in little confidence in a conclusion either for or against a resonant source; we remark that if a resonant source were thought plausible, the now-realized abundant population of objects in high-order resonances beyond the 3:2 (Hahn & Malhotra 2005; Gladman et al. 2008; Lykawka & Mukai 2007) must be seriously examined as an additional source of gradually-destabilized objects which might rival or exceed the 3:2 alone. Arguments favoring the scattering disk as the most likely source of the JFCs (Duncan & Levison 1997) continue to be appealing. Lastly, if either the size distribution of plutinos is shallower, α ∼ 0.52 (Elliot et al. 2005; Hahn & Malhotra 2005) or if the size distribution rolls over to some flat slope as expected then the number of small plutinos drops by a factor of 10 and there are many orders of magnitude too few plutinos available to act as the source of JFCs.