The Extended-Track Event Reconstruction for MiniBooNE

The extended-track reconstruction for MiniBooNE


The Booster Neutrino Experiment (MiniBooNE) searches for oscillations using the neutrino beam produced by the FNAL Booster synchrotron. The array of photomultiplier tubes (PMTs) lining the MiniBooNE detector records Cherenkov and scintillation photons from the charged particles produced in neutrino interactions. We describe a maximum likelihood fitting algorithm used to reconstruct the basic properties (position, direction, energy) of these particles from the charges and times measured by the PMTs. The likelihoods returned from fitting an event to different particle hypotheses are used to categorize it as a signal event or as one of the background processes, in particular charged current quasi-elastic scattering and neutral current production. The reconstruction and event selection techniques described here can be applied to current and future neutrino experiments using similar Cherenkov-based detection.


, , , , , 2

1 Introduction

Cherenkov detectors of have played a key role in the establishment of the phenomenon of neutrino oscillations [?, ?, ?]. These detectors typically consist of a large volume of a homogeneous transparent medium (water or mineral oil) with a high index of refraction surrounded by an array of photomultiplier tubes (PMTs). Cherenkov photons produced by charged particles emerging from the neutrino interactions and traversing the medium with (where is the index of refraction and ) are detected by the PMTs. The photons are emitted at an angle relative to the track direction, where . The radiation is azimuthally symmetric about the track direction, resulting in a ring-like pattern that can be identified on the PMT array. The quantity, spatial distribution, and arrival times of these photons provide information on the location, direction, and energy of the particle. We refer to the extraction of such information from the charge and time measurements of the PMTs as “reconstruction.” An analysis of the ring profile can also provide information on the identity of the particle and, if multiple rings are detected, the number of particles.

In this paper, we discuss the event reconstruction algorithms used in the Booster Neutrino Experiment (MiniBooNE), which searches for an excess of interactions indicative of oscillations using the Fermilab neutrino beam. The MiniBooNE detector is a sphere of radius 610.6 cm filled with Marcol 7 mineral oil () and divided into two optically isolated regions by an opaque shell of radius 575 cm, concentric with the sphere. The surface of the inner “main” region is instrumented with an array of 1280 inward-facing 8 in. PMTs which detect the light produced by the neutrino interactions. The outer “veto” region is instrumented with 240 PMTs and is used to tag charged particles that enter the detector from outside (e.g., cosmic muons) or that exit the main region. Though no scintillator was added to the Marcol 7, it scintillates weakly, resulting in the production of delayed isotropic light for particles with sub- and super-Cherenkov velocities. A more detailed description of the experiment can be found in Refs. [12] and [13].

2 Neutrino Interactions at MiniBooNE

The neutrino beam at MiniBooNE results in interactions with relatively low outgoing particle multiplicity. The largest interaction channel is the charged current quasi-elastic (CCQE) process:


which accounts for of all neutrino interactions in the MiniBooNE detector. Since the recoil proton of Eq. (1) is typically below Cherenkov threshold, only the outgoing lepton produces significant light. In modeling such events within a reconstruction algorithm, one can therefore consider only the Cherenkov and scintillation light produced by the outgoing lepton. While the recoiling nucleon can produce scintillation light, this additional source of light is not considered in the reconstruction.

A muon in the MiniBooNE detector exhibits minimum-ionizing behavior through most of its path with little chance of radiative energy loss. In contrast, electrons typically induce electromagnetic showers, resulting in additional electrons and positrons that emit Cherenkov light. Thus, muons and electrons have significantly different Cherenkov ring patterns. This difference is the basis for distinguishing a muon track from an electron track, and the reconstruction algorithm has a model for each.

The next most abundant process in MiniBooNE is single pion production, which occurs primarily via resonance or coherent scattering:


where denotes a nucleon in the case of resonance production and a nucleus in the case of coherent production. Of particular interest is neutral current (NC) production. The two photons from decay induce electromagnetic showers indistinguishable in MiniBooNE from those induced by electrons. Since the recoil nucleon/nucleus is typically below Cherenkov threshold, these NC events are well-described in the reconstruction algorithm by two electron-like tracks pointing back to a common vertex. The presence of two distinct Cherenkov rings allows for the separation of events from electron ( CCQE) events. However, misidentification (for example, due to a large energy asymmetry between the two photons or due to a small photon opening angle which leaves the Cherenkov rings overlapping) accounts for most of the reducible background in the oscillation search. Thus, successful reconstruction and identification of events is critical.

While other channels can be accommodated by the reconstruction algorithm, the appearance oscillation analysis uses only four event models: single electron track, single muon track, two tracks, and two tracks with a invariant mass.

3 The Detector Response

The Marcol 7 mineral oil used in MiniBooNE exhibits a rich array of optical phenomena despite its  m extinction length near the peak of PMT sensitivity ( nm). Cherenkov and scintillation light production is accompanied by photon absorption, fluorescence (with several possible excitation/emission spectra and lifetimes), and Rayleigh and Raman scattering. The left plot in Figure 1 summarizes the rates of various processes as a function of wavelength. The cumulative extinction rate is shown as the black line. In the near-ultraviolet region (), fluorescence processes dominate the extinction rate. The fluorescence lifetimes range from . In the visible region (), the dominant processes are Rayleigh scattering and absorption.

Figure 1: Left: Rates of optical processes in Marcol 7 as a function of wavelength. The solid black line is the overall extinction rate obtained from laboratory measurements. The dashed black line is the extrapolated rate based on in situ data. The curves labeled “Fluor ” are the excitation rates for the four identified fluorescence processes. The light blue points represent the measured rates of Rayleigh scattering, and the dashed light blue and gray lines represent extrapolated rates. Right: Reconstructed photon arrival times for R1408 (top) and R5912 (bottom) PMTs for 397 nm light flashed from the center of the detector. The black histogram is the distribution from data and the blue is the complete Monte Carlo simulation. The green (red) shows the simulation with reflections (and scattering) turned off.

The 1520 8 in. PMTs in MiniBooNE are of two types: 322 model R5912 PMTs and 1198 model R1408 PMTs, both from Hamamatsu. All of the R5912 PMTs are located in the main PMT array. The PMT time response (particularly the late-pulsing behavior) and the variation of PMT efficiency with incident angle have been characterized in external measurements using a pulsed LED [14]. The PMTs in situ have been studied using a laser calibration system that produces sub-nanosecond bursts of nearly isotropic light in the detector.

The right panel of Figure 1 shows the times of PMT hits recorded in laser calibration events relative to their expected times. The colored histograms come from Monte Carlo simulation, with the red curve showing PMT and electronics timing effects (including pre- and late-pulsing), the green curve including photon scattering, and the blue curve adding photon reflections from the PMTs and the surface of the main detector region. The resulting time structure matches well with that seen in data (black points). We note that only the earliest photoelectron’s time is reported when multiple photoelectrons are present in a hit. (A fresh hit can be seen after an electronics dead time of  ns.) The reported charge, however, does account for additional photoelectrons.

Absorption, scattering, reflections, and fluorescence processes influence both the topology and the time structure of hits. Additionally, scintillation light, emitted with a lifetime of , adds intrinsically delayed light to the prompt Cherenkov photons. The event models within the reconstruction algorithm account for all these phenomena. In the discussion, we refer to Cherenkov and scintillation light that travels undisturbed from the particle track to the PMT as “direct” light. Light that undergoes fluorescence, scattering, or reflection we term “indirect” light.

A detailed GEANT3-based [15] Monte Carlo simulation of the MiniBooNE detector serves as the central tool for developing the reconstruction algorithm’s predictive models. The accuracy of the models influences performance through resolution and biases in extracted track parameters and through the particle separation capabilities of the likelihood-based hypothesis testing. To balance the need for likelihoods that are as accurate as possible with the desire to contain the computational demands, detector and track model information extracted from the Monte Carlo simulation is tabulated, approximated, or parametrized wherever possible while minimizing any sacrifice in performance. Also, thanks to the excellent stability of the detector’s optical properties, a single configuration of the reconstruction is sufficient for use on all MiniBooNE physics data.

4 The Single Track Model

The reconstruction of a single track in the detector, whether an electron or muon, is based on a model with seven parameters:

  • the starting point:

  • the starting time:

  • the direction: ,

  • the kinetic energy: .

The starting point and time are referred to as the “vertex” of the event. The direction is described in terms of the polar angle defined with respect to the direction of the neutrino beam and the azimuthal angle about this beam axis. For the electron model, the stochastic variations in the electromagnetic showers are modeled in average – no attempt is made to account for event-by-event variations in the shower profile. This is also true of the much smaller variations in the muon propagation (e.g., straggling). This simplification allows the track to be fully specified from its initial properties summarized in the seven parameters above. We refer to this vector of parameters as .

The observables of an event are the measurements from the 1280 PMTs in the main region of the detector. For each PMT, we have:

  • a bit indicating whether the tube registered a hit

  • if the tube was hit, the measured charge of the hit

  • if the tube was hit, the measured time of the hit.

Assuming that the PMTs behave independently, the likelihood for an observed set of PMT measurements given track parameters can be expressed as:


where the products are taken over all unhit and hit PMTs and:

  • is the probability that the th tube is hit, given parameters .

  • , are the measured charge and time on the th PMT.

  • is the probability distribution function (PDF) for the measured charge on the th PMT, given , evaluated at .

  • is the PDF for the measured time given , evaluated at .

It is convenient to work with the negative logarithm of , and since the charge- and time-related portions of the likelihood are distinct, it is helpful to define:




For brevity, we refer to and as the charge and time likelihoods, respectively, although they are the negative logarithms of the likelihoods.

4.1 The Charge Likelihood

If the number of photoelectrons produced in a PMT is known, then for the observed charge is fully specified in terms of the response properties of the PMT (assumed fixed) without any reference to any other property of the detector. Further, is a Poisson variable whose mean is a function of the track parameters . These considerations motivate a two-step approach to calculating the charge likelihood.

In the first step, using the known optical photon and particle propagation properties of the detector, one determines for a given particle type (electron/muon) and set of track parameters the average number of photoelectrons that a particular PMT should observe. This quantity is referred to as , the “predicted charge” for the th PMT. In calculating the predicted charge, one considers the quantity and angular distribution of light produced by the track, the absorption, scattering, and fluorescence of light as it traverses the mineral oil, the acceptance of the PMT, and anything else that influences the mapping .

For the second step, given in Eq. (5) can be rewritten in terms of the predicted mean charges :


where and are now functions of the , with the mapping from the first step left implicit. These two functions depend only on the properties of the PMT response and are, in particular, independent of the happenings within the detector volume. This stepwise approach decouples the track and optical photon model from the PMT response model. The latter is established once using laser calibration data. The laser system provides a source of approximately isotropic and prompt light with controllable intensity. The value of for a particular laser intensity can be determined from the occupancy of the PMTs. The function for each PMT type can be determined by obtaining the distribution of measured charges seen at various intensities . These distributions incorporate all processes associated with the charge response, such as the cascade of charge down the dynode chain, the amplification and digitization of the anode signal, and the conversion of the digital readouts into calibrated charge values. Examples of obtained in this way are shown in Figure 2.

Figure 2: Examples of for photoelectrons. Note that the distribution only includes cases where the PMT registers a hit.

4.2 Calculating the Predicted Charge

We begin our description of the predicted charge calculation with a simple scenario involving a fictitious point source of isotropic direct light. This is then generalized to a line source of isotropic light, appropriate for the scintillation production of a track. We then consider a line source in which the emission has a non-trivial (but azimuthally symmetric) angular distribution, representative of Cherenkov radiation. Finally, we account for the contribution of indirect light (i.e., scattering, reflections, etc.) from both the scintillation and the Cherenkov emission. In order to simplify notation and without loss of generality, the discussion involves fixed but arbitrary PMT location and track parameters .

One needs only two parameters to describe the geometric relationship between an isotropic point source and a PMT. We choose (1) the distance between the source and the PMT and (2) the angle of incidence at the PMT, with corresponding to incidence parallel to the PMT axis. The predicted charge on the PMT can be written


where is an event-energy-dependent scintillation light yield (i.e.,  ), is a -dependent solid angle factor, is the transmission of the oil and PMT glass as a function of light propagation distance, and is the acceptance of the PMT as a function of the angle of incidence. The transmission depends on the wavelength spectrum of the light source and, thus, is specific in this case to scintillation. The wavelength-dependent PMT photocathode efficiency is also included in . The solid angle factor and PMT angular acceptance are purely geometric, independent of the particulars of the light production.

To generalize to an extended source of scintillation light, as depicted in Figure 3(left), the latter three factors in Eq. (7) become functions of , the distance along the track from its origin. We must also include an emission profile describing the distribution of scintillation production along the track. This profile satisfies . The predicted charge can now be obtained by integrating along the track’s axis:

Figure 3: Geometry of a line source of isotropic scintillation light (left) and directional Cherenkov light (right). In each case is the distance from the origin of the track, is the distance to the PMT from the point along the track, and is the angle of PMT incidence for light emitted from this point. The presence of Cherenkov light requires the additional parameter describing the light emission angle.

The dependences of , , and on and have been recast as dependences on , as the relationship is fixed for given track parameters and PMT location. Figures 5 and 5 show for 300 muons and electrons, obtained via Monte Carlo simulation. The distribution is relatively flat for muons until the end of the track, where the ionization rate per unit track length becomes larger. (Saturation effects are taken into account in the simulation.) For electrons, the distribution reflects the showering behavior of the track, with rising and falling with the number of shower particles.

Figure 4: Scintillation emission profile for muons as a function of , the distance along the track.
Figure 5: Scintillation emission profile for electrons as a function of , the distance along the track.

We now consider an extended track emitting light with a non-trivial angular distribution, as is the case for Cherenkov radiation. We introduce , the angle to the PMT from the track, as shown in Figure 3(right), and express the predicted charge due to Cherenkov radiation as:


This expression differs from its scintillation counterpart by the presence of an angular emission profile . Note that depends on in two distinct ways: the angular profile of the light changes as the track propagates and loses energy, and the angle to the PMT changes depending on which part of the track we are considering.

Figures 7 and 7 show and for simulated 300 muons and electrons. As the muon propagates, loses energy, and approaches the Cherenkov threshold, the rate of Cherenkov radiation per unit track length decreases and the Cherenkov angle becomes smaller. The scattering of the muon, which causes deviations of the track from its original direction, has been included in the Monte Carlo simulation and results in the spread of the angular distribution about the nominal . For electrons, follows the shower profile, like . The presence of the shower particles is readily apparent in the distribution, which becomes substantially wider as the “track” propagates.

Figure 6: Top: Cherenkov emission profile for muons as a function of , the distance along the track. Bottom: Two-dimensional emission profile showing the angular distribution of Cherenkov light as a function of for a muon.
Figure 7: Top: Cherenkov emission profile for electrons as a function of , the distance along the track. Bottom: Two-dimensional emission profile showing the angular distribution of Cherenkov light as a function of for a electron.

In writing Eq. (9), a few simplifications are made. A more rigorous treatment would involve , the differential probability of sending an emitted Cherenkov photon in the direction , which would be integrated over the solid angle subtended by the PMT. By assuming that is constant over the PMT face, the integration can be effected by simply multiplying the differential probability by the solid angle . Also, the azimuthal symmetry of the emission reduces the two-dimensional PDF to a one-dimensional expression, namely . We take to satisfy


for all values of , with the result that a factor of is absorbed into the definition .

4.3 Indirect Light (scattering, fluorescence, reflections)

The above formalism determines the predicted charge for light arriving at the PMTs directly from the track without any redirection. However, the detector has sources of indirect light from scattering, fluorescence, and reflection as discussed in Section 3.

The geometry for indirect light given scintillation (isotropic) emission is shown in Figure 8(left), where an infinitesimal element along the track is situated at radius from the center of the detector and at angle relative to the position of the PMT. The direct light from this track element is simply the integrand of Eq. (8):


An analytic expression for the indirect light would involve an elaborate integral over emission angles and scattering points throughout the tank. Rather than attempt this, we observe that the value of such an integral must be proportional to the source strength and must otherwise depend only on the topological variables and . The source strength dependence can be eliminated by forming a ratio of the indirect and direct light predictions:


, which we refer to as the scattering table and which we build via the detector simulation, is a property only of the detector optics and the of the track element. With this table, the indirect light contribution can be immediately incorporated into the expression for predicted charge:


where the dependence of and on has been made explicit.

Figure 8: Geometry of indirect light from an isotropic (left) and directional (right) source of light.

For Cherenkov light, the situation is more complex since the light emission is anisotropic. Two additional variables are needed to specify the direction of a vector source relative to the PMT position and the tank center. We make the following non-unique choice for these two variables, as shown in Figure 8(right):

  • : the angle between (1) the source direction vector and (2) the source-to-PMT ray (the same as defined elsewhere).

  • : The angle between (1) the plane containing the tank center, the PMT, and the source, and (2) the plane containing the track and the tank center.

The intensity of the indirect Cherenkov light is normalized with respect to a fictitious isotropic Cherenkov source with predicted charge , so that:


The total Cherenkov contribution to the mean predicted charge, including both direct and indirect components is


where is given by Eq. (9) and is given by


Implicit in this expression is the assumption that the Cherenkov angular emission profile does not change as the track propagates (and loses energy). This simplification has negligible impact, as the processes that result in indirect light tend to destroy the initial emission pattern anyway.

We now have a complete charge prediction, where the total predicted charge for the PMT is given by:


4.4 Computation of Integrals

The predicted charge integrals above would be impractical to evaluate with sufficient spatial granularity within a maximum likelihood fit, where they must be evaluated multiple times for every PMT in every event. To mitigate this, all integrations are performed and tabulated beforehand as follows.

The integrand of Eq. (8), which pertains to direct scintillation light, is the product of a “source” factor and an “acceptance” factor


As in Eq. (8), the dependence on and is recast through . If varies gradually enough, it can be approximated with a parabolic form:


An example of this is shown in Figure 9. The predicted charge then becomes


The first integral is identically unity. The remaining two depend on the energy via , but they depend on no other track parameters. This allows one to tabulate the integrals beforehand as a function of , eliminating their evaluation from the actual likelihood calculation.

Figure 9: Example of for the track geometry shown in the top figure, where a muon 3.48 meters from a PMT moves perpendicular to the PMT axis. is evaluated at three points located along the track, corresponding to the three points in the figure. In the bottom panel, the parabolic approximation for based on the three evaluated points (dashed red) is compared to the exact form (solid black).

The parabolic coefficients are obtained by evaluating Eq. (18) at three points along the track: the start of the track (), the midpoint of the Cherenkov emission profile () and twice this distance from the start (). The choice of where to evaluate is somewhat arbitrary; one could choose any three points that sample a large fraction of the range over which light is produced. The extension to indirect isotropic light is straightforward with a redefinition of that incorporates the scattering table:


An analogous three-integral expression can be formed for the Cherenkov predicted charge:


As in the scintillation case, the integrals depend on the energy of the track. However, they also depend on two parameters which define the PMT-track geometry. We choose these parameters to be the vertex-to-PMT distance, ; and the cosine of the angle-to-PMT viewed from the vertex, . Labeling the integrals , we obtain


Although the integrals depend on three parameters, it is still feasible to tabulate their values ahead of time, eliminating their explicit evaluation within the likelihood calculation.

For indirect Cherenkov light, recall that the scattering table is normalized to a fictitious isotropic Cherenkov source. Thus, the parabolic method used for indirect scintillation can be applied here to give


with two energy-dependent integrals to be tabulated.

4.5 The Time Likelihood

The time portion of the likelihood involves the PDF for the measured PMT hit time given track parameters . Some of the dependence on the track parameters can be eliminated by using the “corrected time” :


The expression removes from three terms, namely the starting time of the track , the expected time for light to propagate from the track midpoint to the PMT (with the speed of light in mineral oil given by ), and the time for the particle to propagate from its starting point to the midpoint. The dependence of on the track energy is made explicit. Defining with respect to the track midpoint (as opposed to, say, the start of the track) improves the validity of the simplifications we make below.

Given the substitution, we seek the PDF for the corrected time given arbitrary PMT-track configurations. The space of configurations is five dimensional. However, producing tables of as a function of five parameters is impractical. To reduce the task, we make the assumption that the corrected time PDF depends only on the track energy, the predicted “prompt” charge, and the predicted “late” charge, where “prompt” corresponds roughly to light arriving at the PMT directly from the track without delays induced by emission lifetimes, scattering, etc. Loosely, the assumption is that the shape of the corrected time spectrum is dominated by the physical extent of the track (characterized by its energy), and by the amount of prompt and late light reaching the PMT. The extent of the track affects the spread of possible hit times, since there is a spread of photon production times, while the amounts of prompt and late light affect the proportion of peaked “prompt” response to the tail of “late” response.

This assumption reduces the configuration to three dimensions. We make one further simplification by assuming that for a given energy , can be modeled by separate primitive prompt and late distributions, indexed by and , respectively, from which we can construct the full PDF. As a result, the PDF is indexed not by the fundamentally two-dimensional space of , but by and separately, with the full PDF calculated on-the-fly as described below. In practice, Cherenkov and scintillation primitive distributions are created and used as proxies for the desired prompt and late distributions, respectively.

Figure 10: Distributions of for direct Cherenkov light for R5912 PMTs. The primary difference between the two PMT types is that the R1480 PMTs show wider prompt timing peaks than do the newer R5912 PMTs.
Figure 11: For 300 MeV muons, the mean and width parameters from Gaussian fits like those in Figure 10 versus . The dependence is parameterized as a sixth-order polynomial in .
Figure 12: Third-level fits parametrizing the energy dependence of the second-level parameters describing the Gaussian mean and width parameters.
Figure 13: Check of primitive distributions. The parametrized likelihood distributions for Cherenkov light are compared against the actual distributions for 250 (top), 600 (middle) and 1500 (bottom) muons as a function of predicted charge. (Note that the high-, low- hits, which are less well modeled, are rare.)

The Cherenkov (prompt) primitive distributions are created by simulating particles throughout the detector with isotropically chosen directions and fixed energy . The particles are created with direct Cherenkov light only; all other sources of light (scintillation, scattering, etc.) are turned off. For each event, the true track parameters are used to evaluate the direct Cherenkov predicted charge for each hit. Histograms of the corrected time of the hits are produced for various ranges of predicted charge. Figure 10 shows three such distributions for muons. Since only direct Cherenkov light is present, the histograms show no late-time features. The shapes of the time spectra depend on the values involved, with the distributions becoming narrower and earlier with increasing predicted charge due to the increasing probability that an early photon will be recorded at the PMT. The distribution in each predicted charge range is fit to a Gaussian distribution, and the resulting Gaussian parameters (mean and width) are subsequently parametrized across using a sixth-order polynomial in . Figure 11 shows an example of these “second-level” fits for muons. The procedure is repeated at many values of , with the two second-level fits providing seven parameters for the Gaussian mean and seven parameters for the Gaussian width at each energy. The energy dependences of these 14 parameters are then fit as a function of to fourth-order polynomials in a third-level parametrization, as shown in Figure 12 for the first seven of the second-level parameters. In addition to conveniently summarizing the dependence of the time response on and , the parametrizations provide a smooth likelihood surface, as required by the minimization algorithm.

This completes the prompt primitive distributions . They are calculated for a given predicted charge at a PMT for a track with energy by evaluating the 14 third-level curves which parametrize the dependence of the second-level functions. These 14 values then determine the two functions which parametrize the dependence of the mean and width of the Gaussian parameters on at . These are evaluated at the particular value of to determine the appropriate mean and width of the Gaussian distribution which describes the Cherenkov primitive distribution appropriate for the particular values of and . The validity of the assumptions and simplifications can be tested directly by comparing the calculated primitive distributions with those actually observed in the Monte Carlo simulated events. Figure 13 shows such a comparison.

Figure 14: Check of scintillation primitive distributions. The parametrized likelihood distributions for scintillation light are compared against the actual distributions for 250 (top), 600 (middle) and 1500 (bottom) muons as a function of predicted charge.

For the scintillation (late) primitive distributions, events are generated with scintillation light only, and the histograms are created in bins of . Each histogram is fit to a sum of two exponentials convolved with a Gaussian distribution. The exponential decay constants are fixed at and , leaving three free parameters: the time origin, the Gaussian resolution, and the relative weight of the two exponentials. At a given , the three parameters are extracted in each histogram to obtain their dependence on , which is parametrized as a sixth-order polynomial in . The dependence of the three sets of seven parameters from the polynomials are parametrized as a function of by fourth-order polynomials, in analogy with the Cherenkov case. We now have , the scintillation primitive distribution. Figure 14 shows the primitive distributions compared with actual distributions from the Monte Carlo simulation for muons at , 600 and .

To obtain the complete PDF for a given PMT hit, we first divide the predicted charge into “prompt” and “late” components based on the sources:


These definitions of and encapsulate the following assumptions. First, all prompt light is due to direct Cherenkov light; all other light (including indirect Cherenkov light) follows the late scintillation-based time distribution. The of direct Cherenkov light which is included in accounts for the late pulsing of PMTs which causes promptly arriving light to appear at late times. The late time distribution is used to model indirect Cherenkov light since the fluorescence, scattering, and reflection processes which lead to indirect light have a time structure similar to that of scintillation. Further, late timing is less critical to the reconstruction, as prompt light provides essentially all the useful timing information.

The and values are used to combine the prompt and late primitive distributions to form the total PDF . In the process, one must account for the fact that a PMT can register only one hit for a given track due to the latency period of the electronics. Since the time reported by the PMT reflects the time of the first photoelectron, we assume that a PMT hit that includes a prompt photoelectron obeys the prompt primitive distribution regardless of how many late photoelectrons may follow. The Poisson distribution gives us the probabilities that no prompt or late photoelectrons are produced based on their expected mean values:


With these probabilities, one can calculate the probability that a given hit contains at least one prompt photoelectron:


We assign this as the weight for the prompt primitive distribution within the total PDF, while is used as the weight for the late primitive distribution, yielding:


where and are the prompt Cherenkov-based and late scintillation-based timing distributions, respectively. Note that Eq. (28) removes the overall probability of a hit occurring. Thus, even if the absolute probability of a prompt photoelectron is small, if .

Figure 15: Internal fit parameters for (a) a single muon or electron track and (b) two photon tracks. Each photon track includes a conversion distance parameter . The two-track parameters can be constrained such that the invariant mass of the two photons is always .

5 The Two-track Reconstruction

The single track model with seven parameters is sufficient for reconstructing and CCQE events as well as cosmic muons or their decay electrons. As discussed in Section 2, NC events require a two-track model. Figure 15 shows the 12 parameters needed to describe two photon tracks originating from a common vertex. The electron track model serves double-duty as the photon track model, under the assumption that a showering electron is indistinguishable from a showering photon apart from the conversion distance. (The mean conversion length in Marcol 7 is 67 cm.) While the two-track model is conceptually a straightforward extension of the one-track case, complexities of the likelihood space require special treatment in the likelihood maximization algorithm.

5.1 Two-track Charge Likelihood

Since the charge likelihood depends only on the measured charge and the total predicted charge at each PMT, the two-track charge likelihood can be obtained by adding together the predicted charges from the two tracks to obtain the total predicted charge at each PMT:


5.2 Two-track time likelihood

The one-track time likelihoods account for the “first photoelectron only” nature of the electronics by calculating the probability that a prompt photoelectron is present and by then weighting the prompt primitive distribution’s influence on the full time PDF accordingly. This scheme is extended to handle two tracks as follows.

The two single-track primitive distributions, and are formed for each track, where labels the track number. In anticipation of aggregating all late light later, the late “scintillation” primitive distributions are averaged to form:


Of the two tracks, one will have a midpoint that is nearer to the target PMT than the other. This track’s quantities are labeled with “n” (near); the other’s, “f” (far). In analogy with Eq. (26), we define:


The probabilities of not obtaining a photoelectron from these sources are:


from which the following weights are obtained:


where is the probability that a prompt photoelectron from the near track exists given that any photoelectron exists, etc. The weights are used to combine the two prompt primitive distributions and the averaged late primitive distribution:


resulting in the complete two-track corrected time PDF.

6 Minimization of

6.1 One-track

With the likelihood defined, the parameter set that minimizes its negative logarithm must be found. In the single track minimization process, two issues arise.

First, the energy parameter is tied to the geometry of the event via the track profiles . If the spatial parameters (, etc) are varied simultaneously with , the minimization algorithm can get confused by this correlation. The solution is to iterate the minimization process, as described below.

A second issue arises from the discrete PMT lattice, which imprints a small fluctuating signal on despite the smooth parameterizations of the input tables. Minimization algorithms that rely on gradients or that do not have controllable step sizes are thus poorly suited to the problem. Therefore, the SIMPLEX method in Minuit[19] is chosen over the MIGRAD method.

To initiate the minimization process, a vector of seed parameters is obtained from a fast fitter. In the first iteration of the minimization, the energy is held at its seeded value while the six remaining parameters (, , , , , ) are varied via Minuit/SIMPLEX to find a temporary minimum of . These six parameters are then fixed, while is freed for a second iteration. Finally, is fixed once more with the six other parameters varied (as in the initial iteration), to find the final minimum of F. The parameters from the final SIMPLEX call are returned as the best-fit set .

Figure 16: The two fitter configurations discussed in the text.

6.2 Two-track

Figure 16 shows two situations that may arise in fitting a event. The event has two photon tracks whose Cherenkov rings are represented in black, while green and red represent two ring configurations corresponding to different starting values for the parameters. In case (a), the 12 parameters are near the correct configuration. In case (b), the parameters are such that both tracks are directed toward the dominant ring. The latter case corresponds to a local minimum which may offer a worse likelihood than case (a); however, the minimization algorithm may be trapped by the fact that any small changes to the parameters result in an increased value of . For example, sweeping one track over to the smaller black ring involves passing through a region with little detected Cherenkov light. The intermediate states form a barrier of disfavored values of that the minimization algorithm code will have difficulty overcoming, especially if several parameters must be adjusted simultaneously even to make the distant configuration an improvement.

The two-track fits require a minimization approach that avoids these traps in the likelihood surface. Monte Carlo events were used to identify both common and rare trapping scenarios, and the minimization algorithm addresses each of these through the following collection of seed configurations and minimization sequences:n

  • The conversion lengths and are seeded with either 50 cm or 250 cm, leading to four possible pairs.

  • and are seeded with the results from the electron single-track fit or with one of the eight perturbations (see below).

  • and are seeded with the “best” directions from a full grid of tested directions (see below).

  • The four-vertex of the event is seeded with the four-vertex from the electron single-track fit shifted according to , and .

  • is seeded with approximately from the electron single-track fit, while the seed for is, when possible, the energy needed to give an invariant mass equal to .

For the last item, the energy seeds are based on an empirical expression that accounts for the second ring’s energy contribution to the one-track fit energy as a function of the angle between the two tracks. The function is obtained by simulating NC events and reconstructing them with the single track electron fit. Using the true photon energies () and opening angle , the relationship between and is parametrized by a linear function. This function is used to determine the seed and values for a given and .

The nine () possible seed directions for track 1 are created as follows. The spatial distribution of PMT charge is projected onto a plane perpendicular to the electron fit direction. The major axis of the covariance ellipse of the resulting two-dimensional distribution is found. The nine possible seed directions are rotations of the best electron-fit direction by , , and radians parallel to the major axis and radians perpendicular to the major axis.

The seed directions for track 2 come from either a “coarse” grid or a “fine” grid. The grid type is specified at run time.

Of the or possible combinations of track 1 and track 2 seed directions, only six are used to seed actual minimization sequences. They are the six combinations with (1) the best total likelihood (time and charge), (2) the best charge likelihood when the estimated track energies3 are similar (), (3) the best charge likelihood when the estimated track energies are dissimilar, and (4)(6) the second-best likelihood in each of these three cases. A small exception: combinations involving one of the eight perturbed track 1 directions and having estimated track energies that differ by more than a factor of five are not considered. In these energetically asymmetric cases, the electron fit direction would be a good representation of track 1 since the influence of track 2 would be small, so there is no reason to consider alternate track 1 directions.

The six best and second-best track 1 and track 2 direction and energy combinations are chosen anew for each permutation of , leaving 24 complete parameter sets in total. Each set seeds two Minuit sequences. In sequence 1, the SIMPLEX minimization routine is first run with the energy and angle parameters held fixed; once a minimum is found, SIMPLEX is called again with all parameters free. In sequence 2, only the energy parameters are held fixed in the initial SIMPLEX call. This results in 48 fit sequences, each of which ends with a completely free SIMPLEX minimization. The best parameter set (i.e., the one that results in the lowest ) from the forty-eight final SIMPLEX calls is reported as the best parameter set for the two-track fit. During all minimizations, the photon conversion points are constrained to be within the main detector volume. This prevents the fit from removing the influence of a track through an inflated conversion length. The event vertex itself, however, is not constrained.

The SIMPLEX minimizations can be performed with a constraint on the invariant mass. This is accomplished by removing as a free parameter and setting it via the relation:


where is the angle between the photon tracks. This fixed-mass mode is the actual hypothesis, while the free-mass mode allows for mass reconstruction. Both are used for identification, with the former lending its maximized likelihood and the latter providing a reconstructed mass .

7 Performance

The performance of the reconstruction algorithm on simulated neutrino interactions ( CCQE, CCQE and NC events) generated with NUANCE [16] and propagated through the detector Monte Carlo simulation is shown in Figures 1721.4 The electron and muon one-track fits are applied to the and CCQE events, respectively, while the two-track fits, with and without the mass constraint, are applied to the NC events. A minimal selection is applied to ensure that the events are contained within the main detector region and, for plots not directly related to vertex reconstruction, within the fiducial volume of the analysis.

A few observations about the figures, particularly Figure 18: (1) The resolutions for CCQE events degrade as the energy approaches the Cherenkov threshold from above, as most of the information in the event comes from the (rapidly decreasing) prompt Cherenkov light. (2) A low-energy offers little information on its direction, as all that is observed is the isotropically produced final state. Thus, the angular resolution (but not the vertex resolution) becomes poor. (3) The vertex resolution for events is worse than CCQE events since the former has spatially displaced light production thanks to the average photon conversion length. (4) Typical performance numbers for CCQE events: 20 cm vertex resolution, kinetic energy resolution, and angular resolution. (5) Typical performance numbers for CCQE events: 10 cm vertex resolution, kinetic energy resolution and angular resolution. (6) Figure 19 shows a mass resolution of .

Figure 17: Reconstructed radial event position (left) and kinetic energy (right) plotted against their true values. The top, middle, and bottom panels show CCQE, CCQE, and NC events.
Figure 18: Radius, kinetic energy, and direction resolutions. The jitter is due to limited test sample statistics in some regions. The  curves cut off at 100 MeV due to an unrelated upstream cut on the Monte Carlo sample.
Figure 19: Fitted invariant mass for Monte Carlo-simulated NC events in three low energy regions (left) and three high energy regions (right).
Figure 20: Reconstructed energies in Monte Carlo-simulated NC events plotted against true values. The left (right) panel shows the higher (lower) energy from each event. The association of each fitted track to the underlying true ’s is, in general, ambiguous. For these plots, we choose the assignment that gives the smaller combined energy and direction discrepancy.
Figure 21: Energy (left) and direction (right) resolutions for ’s in simulated NC events. For a given true energy, a dominant has expectedly better resolution than does a subordinate .

8 Particle Identification

The maximum likelihoods returned from each fit can be used for hypothesis testing. The two quantities




are used to determine for a given event whether the electron hypothesis is preferred over the muon and hypotheses. In these expressions, , and are the maximized likelihoods returned by the electron, muon, and (fixed-mass) two-track fits, respectively.

Figure 22: Distribution of ) for Monte Carlo simulated CCQE events (top) and CCQE events (bottom) as a function of reconstructed energy (from the electron hypothesis fit).
Figure 23: Left: Distribution of for Monte Carlo simulated CCQE events (top) and NC events (bottom). Right: same for . The black lines on each distribution show the selection used in the MiniBooNE oscillation analysis [12].

8.1 Electron/Muon separation

Figure 22 shows for simulated CCQE events and CCQE events the distribution of as a function of the energy reconstructed by the electron fit. For each sample, the events have been subject to a set of preselection criteria used in the appearance analysis. These require that there is only one time-cluster of PMT hits in the event, eliminating obvious CC events that produce decay electrons. Cosmic backgrounds are eliminated by requiring the minimum number of PMT hits in the main region to be and the number of veto hits to be . Furthermore, the average time of hits is required to lie within the expected beam delivery window from the Booster.

For the CCQE events, tends to take positive values, indicating that the fit to the event with the electron hypothesis is favored over the muon hypothesis. Likewise, tends to be negative for the CCQE events, indicating that the muon hypothesis fits these events better than the electron hypothesis. At high energies (), the electron/muon separation is aided by the fact that the muon pathlength grows approximately linearly with energy, while the electron shower grows more slowly.

8.2 Electron/ separation

The left plots in Figure 23 show the distributions of as a function of energy reconstructed by the single-track electron reconstruction for simulated CCQE events and NC events. The quantity tests whether a given event fits better as a single electron track or as a . Overall, electron/ separation becomes more difficult at high energies as the energies of the decay photons in events become more asymmetric and the shower fluctuations in a single electron event get large enough to mimic the presence of a second photon. Also, events in which one of the two photons leaves the detector unconverted present an essentially irreducible background at low energy.

Another quantity that can be used for electron/ separation is , the invariant mass of the two photons returned by the free-mass two-track fit. The NC events peak at the mass, whereas the CCQE events peak toward zero. This is seen in the right plots in Figure 23, where the distribution is shown as a function of reconstructed energy. The use of the true mass in the fit’s seeding procedure (Section 6.2) induces no appreciable high-mass peak in the CCQE events while helping considerably in the correct identification of NC events.

The appearance analysis at MiniBooNE utilizes both and to separate CCQE events from NC .

8.3 Particle Identification Performance

The background rejection levels reachable with these fitter-derived quantities are shown for Monte Carlo events in Figures 24 and 25. Figure 24 shows the efficiency, after fiducial volume and cosmic rejection cuts, for selecting signal CCQE events alongside the efficiency (or misidentification rate) for CCQE events. Three prototype selections are used, each chosen to give an approximately fixed signal efficiency at all energies. Note that the misidentification rate is shown with a logarithmic vertical scale. Figure 25 shows similar plots for rejection, one for the cut and one for the cut.

Figure 24: Performance of three prototype cuts for the selection. Left: Signal ( CCQE) efficiency versus the electron reconstruction’s energy. Right: Misidentification rate of CCQE events versus the electron reconstruction’s energy.
Figure 25: Performance of three prototype (top) and (bottom) cuts for the selection. Left: Signal ( CCQE) efficiency versus the electron reconstruction’s energy. Right: Misidentification rate of NC events versus the electron reconstruction’s energy.

Misidentification rates of a few percent for CCQE events are seen for signal efficiencies from . (Note that CCQE events are further reduced by about an order of magnitude by recognizing the delayed muon decay.) The two variables indicate that the misidentification rate increases with energy; this is expected, since the two decay photons become highly asymmetric energetically and since the Cherenkov rings can achieve a greater degree of overlap. The rise in the misidentification rate at low energies is due to events with an escaping photon (corresponding to the cluster of events with positive at 50200  in Figure 23). Between 200  and 600 , where most of the NC events lie, a misidentification rate of less than 1 percent can be achieved with efficiency using only the cut, while rates of a few percent can be achieved with only the cuts. In the MiniBooNE oscillation analysis, 50% signal efficiency and 1% misidentification rate are obtained by requiring both and , where the “max” functions are quadratic in the energy obtained from the single-track electron fit. These quadratic selection boundaries are indicated by the solid black lines in Figure 23. These cuts, along with an analogous cut, were tuned to optimize sensitivity to LSND-like oscillations.

9 Summary

A maximum likelihood reconstruction algorithm is used at MiniBooNE to determine track parameters under electron, muon, and hypotheses. Underlying the likelihood is a track and detector model that calculates for a given parameter set the charge and time PDFs expected for each PMT, accounting for the spatially extended production of Cherenkov and scintillation light as well as the effects of indirect light from subsequent optical processes.

The maximized likelihoods for a given event obtained under different event hypotheses are used for event selection. In particular, the ratio of the likelihoods under the electron and muon models is used to suppress CC events in the MiniBooNE oscillation search; the ratio of the likelihoods under the electron and two-track fit (with fixed invariant mass) and the invariant mass obtained from the free-mass two-track fit are used to suppress NC events. While the reconstruction has been developed within the context of the MiniBooNE oscillation search, other experiments employing similar Cherenkov detection techniques should be able to use the techniques discussed here.


  1. thanks: Now at the California Institute of Technology, Physics Department 103-33, Pasadena, CA 91125, USA.
  2. thanks: Now at the University of British Columbia, Department of Physics and Astronomy, Vancouver, BC V6T 1Z1, Canada.
  3. Recall the energy seed procedure discussed above.
  4. For comparisons of Monte Carlo to data and for discussion of the development and validation of the Monte Carlo, see Refs. [17] and [18].


  1. kamioka
  2. K. S. Hirata et al. [Kamiokande-II Collaboration], Phys. Rev. Lett. 65, 1301 (1990).
  3. Y. Fukuda et al. [Kamiokande Collaboration], Phys. Lett. B 335, 237 (1994).
  4. superk
  5. S. Fukuda et al. [Super-Kamiokande Collaboration], Phys. Lett. B 539, 179 (2002) [arXiv:hep-ex/0205075].
  6. M. H. Ahn et al. [K2K Collaboration], Phys. Rev. Lett. 90, 041801 (2003) [arXiv:hep-ex/0212007].
  7. S. Yamamoto et al. [K2K Collaboration], Phys. Rev. Lett. 96, 181801 (2006) [arXiv:hep-ex/0603004].
  8. sno
  9. Q. R. Ahmad et al. [SNO Collaboration], Phys. Rev. Lett. 87, 071301 (2001) [arXiv:nucl-ex/0106015].
  10. Q. R. Ahmad et al. [SNO Collaboration], Phys. Rev. Lett. 89, 011301 (2002) [arXiv:nucl-ex/0204008].
  11. S. N. Ahmed et al. [SNO Collaboration], Phys. Rev. Lett. 92, 181301 (2004) [arXiv:nucl-ex/0309004].
  12. A. A. Aguilar-Arevalo et al. [The MiniBooNE Collaboration], Phys. Rev. Lett. 98, 231801 (2007) [arXiv:0704.1500 [hep-ex]].
  13. A. A. Aguilar-Arevalo et al. [The MiniBooNE Collaboration], Nucl. Instrum. Meth. A 599, 28 (2009).
  14. S. J. Brice et al., Nucl. Instrum. Meth. A 562, 97 (2006).
  15. CERN Program Library Long Writeup W5013 (1993).
  16. D. Casper, Nucl. Phys. Proc. Suppl. 112, 161 (2002) [arXiv:hep-ph/0208030].
  17. R. B. Patterson, Ph.D. thesis, Princeton University, FERMILAB-THESIS-2007-19 (2007)
  18. A. A. Aguilar-Arevalo et al. [MiniBooNE Collaboration], Phys. Lett. B 664, 41 (2008) [arXiv:0803.3423 [hep-ex]].
  19. CERN Program Library Long Writeup D506
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description