Calibration of the charge and energy loss per unit length of the MicroBooNE liquid argon time projection chamber using muons and protons
MicroBooNE collaboration C. Adams M. Alrashed R. An J. Anthony J. Asaadi A. Ashkenazi S. Balasubramanian B. Baller C. Barnes G. Barr V. Basque M. Bass F. Bay S. Berkman A. Bhanderi A. Bhat M. Bishai A. Blake T. Bolton L. Camilleri D. Caratelli I. Caro Terrazas R. Carr R. Castillo Fernandez F. Cavanna G. Cerati Y. Chen E. Church D. Cianci E. O. Cohen J. M. Conrad M. Convery L. Cooper-Troendle J. I. Crespo-Anadon M. Del Tutto D. Devitt A. Diaz L. Domine K. Duffy S. Dytman B. Eberly A. Ereditato L. Escudero Sanchez J. Esquivel J. J. Evans R. S. Fitzpatrick B. T. Fleming N. Foppiani D. Franco A. P. Furmanski D. Garcia-Gamez S. Gardiner V. Genty D. Goeldi S. Gollapinni O. Goodwin E. Gramellini P. Green H. Greenlee R. Grosso L. Gu W. Gu R. Guenette P. Guzowski P. Hamilton O. Hen C. Hill G. A. Horton-Smith A. Hourlier E. C. Huang R. Itay C. James J. Jan de Vries X. Ji L. Jiang J. H. Jo R. A. Johnson J. Joshi Y. J. Jwa G. Karagiorgi W. Ketchum B. Kirby M. Kirby T. Kobilarcik I. Kreslo I. Lepetic Y. Li A. Lister B. R. Littlejohn S. Lockwitz D. Lorca W. C. Louis M. Luethi B. Lundberg X. Luo A. Marchionni S. Marcocci C. Mariani J. Marshall J. Martin-Albo D. A. Martinez Caicedo K. Mason A. Mastbaum N. McConkey V. Meddage T. Mettler K. Miller J. Mills K. Mistry T. Mohayai A. Mogan J. Moon M. Mooney C. D. Moore J. Mousseau M. Murphy R. Murrells D. Naples R. K. Neely P. Nienaber J. Nowak O. Palamara V. Pandey V. Paolone A. Papadopoulou V. Papavassiliou S. F. Pate A. Paudel Z. Pavlovic E. Piasetzky D. Porzio S. Prince G. Pulliam X. Qian J. L. Raaf A. Rafique L. Ren L. Rochester H. E. Rogers M. Ross-Lonergan C. Rudolf von Rohr B. Russell G. Scanavini D. W. Schmitz A. Schukraft W. Seligman M. H. Shaevitz R. Sharankova J. Sinclair A. Smith E. L. Snider M. Soderberg S. Soldner-Rembold S. R. Soleti P. Spentzouris J. Spitz M. Stancari J. St. John T. Strauss K. Sutton S. Sword-Fehlberg A. M. Szelc N. Tagg W. Tang K. Terao R. T. Thornton M. Toups Y. -T. Tsai S. Tufanli T. Usher W. Van De Pontseele R. G. Van de Water B. Viren M. Weber H. Wei D. A. Wickremasinghe Z. Williams S. Wolbers T. Wongjirad K. Woodruff M. Wospakrik W. Wu T. Yang G. Yarbrough L. E. Yates G. P. Zeller J. Zennamo C. Zhang
physics.ins-det hep-ex

=1
  
















       
   
         
 
 


conditions* 
  
>l<@= >X
MicroBooNE Collaborationj]C. Adamsl]M. Alrashedk]R. Anc]J. Anthonycc]J. Asaadip]A. Ashkenazihh]S. Balasubramaniani]B. Ballerq]C. Barnest]G. Barro]V. Basqueb]M. Bassdd]F. Bayi]S. Berkmano]A. Bhanderiz]A. Bhatb]M. Bishaim]A. Blakel]T. Boltong]L. Camillerii]D. Caratellif]I. Caro Terrazasp]R. Carri]R. Castillo Fernandezi]F. Cavannai]G. Ceratia]Y. Chenu]E. Churchg]D. Cianciaa]E. O. Cohenp]J. M. Conradx]M. Converyhh]L. Cooper-Troendleg]J. I. Crespo-Anadónt]M. Del Tuttom]D. Devittp]A. Diazx]L. Dominei]K. Duffyv]S. Dytmanh]B. Eberlya]A. Ereditatoc]L. Escudero Sanchezz]J. Esquivelo]J. J. Evansq]R. S. Fitzpatrickhh]B. T. Flemingj]N. Foppianihh]D. Francoo]A. P. Furmanskio]D. Garcia-Gamezi]S. Gardinerg]V. Gentya]D. Goeldibb]S. Gollapinnio]O. Goodwini]E. Gramellinio]P. Greeni]H. Greenleee]R. Grossoff]L. Gub]W. Guj]R. Guenetteo]P. Guzowskiz]P. Hamiltonp]O. Heno]C. Hilll]G. A. Horton-Smithp]A. Hourliern]E.-C. Huangx]R. Itayi]C. Jamesc]J. Jan de Vriesb]X. Jiv]L. Jianghh]J. H. Joe]R. A. Johnsonb]J. Joshig]Y.-J. Jwag]G. Karagiorgii]W. Ketchumb]B. Kirbyi]M. Kirbyi]T. Kobilarcika]I. Kreslok]I. Lepeticb]Y. Lim]A. Listerk]B. R. Littlejohni]S. Lockwitza]D. Lorcan]W. C. Louisa]M. Luethii]B. Lundberghh]X. Luoi]A. Marchionnii]S. Marcocciff]C. Marianigg]J. Marshallj]J. Martin-Alboy]D. A. Martinez Caicedoee]K. Masond]A. Mastbaumo]N. McConkeyl]V. Meddagea]T. Mettlerd]K. Milleree]J. Millso]K. Mistrybb]A. Mogani]T. Mohayaip]J. Moonf]M. Mooneyi]C. D. Mooreq]J. Mousseauff]M. Murphyo]R. Murrellsv]D. Naplesl]R. K. Neelyw]P. Nienaberm]J. Nowaki]O. Palamaraff]V. Pandeyv]V. Paolonep]A. Papadopoulour]V. Papavassiliour]S. F. Patel]A. Paudeli]Z. Pavlovicaa]E. Piasetzkyo]D. Porzioj]S. Princez]G. Pulliamb]X. Qiani]J. L. Raafl]A. Rafiquer]L. Renx]L. Rochesterf]H.E. Rogersg]M. Ross-Lonergana]C. Rudolf von Rohrhh]B. Russellhh]G. Scanavinid]D. W. Schmitzi]A. Schukraftg]W. Seligmang]M. H. Shaevitzee]R. Sharankovaa]J. Sinclairc]A. Smithi]E. L. Sniderz]M. Soderbergo]S. Söldner-Remboldt,j]S. R. Soletii]P. Spentzourisq]J. Spitzi]M. Stancarii]J. St. Johni]T. Straussg]K. Suttonr]S. Sword-Fehlbergo]A. M. Szelcs]N. Taggbb]W. Tangx]K. Teraon]R. T. Thorntoni]M. Toupsx]Y.-T. Tsaihh]S. Tufanlix]T. Ushert,j]W. Van De Pontseelen]R. G. Van de Waterb]B. Virena]M. Weberb]H. Weiv]D. A. Wickremasinghecc]Z. Williamsi]S. Wolbersee]T. Wongjiradr]K. Woodruffi]M. Wospakriki]W. Wui]T. Yangbb]G. Yarbroughp]L. E. Yatesi]G. P. Zelleri]J. Zennamob]C. Zhang[a]Universität Bern, Bern CH-3012, Switzerland[b]Brookhaven National Laboratory (BNL), Upton, NY, 11973, USA[c]University of Cambridge, Cambridge CB3 0HE, United Kingdom[d]University of Chicago, Chicago, IL, 60637, USA[e]University of Cincinnati, Cincinnati, OH, 45221, USA[f]Colorado State University, Fort Collins, CO, 80523, USA[g]Columbia University, New York, NY, 10027, USA[h]Davidson College, Davidson, NC, 28035, USA[i]Fermi National Accelerator Laboratory (FNAL), Batavia, IL 60510, USA[j]Harvard University, Cambridge, MA 02138, USA[k]Illinois Institute of Technology (IIT), Chicago, IL 60616, USA[l]Kansas State University (KSU), Manhattan, KS, 66506, USA[m]Lancaster University, Lancaster LA1 4YW, United Kingdom[n]Los Alamos National Laboratory (LANL), Los Alamos, NM, 87545, USA[o]The University of Manchester, Manchester M13 9PL, United Kingdom[p]Massachusetts Institute of Technology (MIT), Cambridge, MA, 02139, USA[q]University of Michigan, Ann Arbor, MI, 48109, USA[r]New Mexico State University (NMSU), Las Cruces, NM, 88003, USA[s]Otterbein University, Westerville, OH, 43081, USA[t]University of Oxford, Oxford OX1 3RH, United Kingdom[u]Pacific Northwest National Laboratory (PNNL), Richland, WA, 99352, USA[v]University of Pittsburgh, Pittsburgh, PA, 15260, USA[w]Saint Mary's University of Minnesota, Winona, MN, 55987, USA[x]SLAC National Accelerator Laboratory, Menlo Park, CA, 94025, USA[y]South Dakota School of Mines and Technology (SDSMT), Rapid City, SD, 57701, USA[z]Syracuse University, Syracuse, NY, 13244, USA[aa]Tel Aviv University, Tel Aviv, Israel, 69978[bb]University of Tennessee, Knoxville, TN, 37996, USA[cc]University of Texas, Arlington, TX, 76019, USA[dd]TUBITAK Space Technologies Research Institute, METU Campus, TR-06800, Ankara, Turkey[ee]Tufts University, Medford, MA, 02155, USA[ff]Center for Neutrino Physics, Virginia Tech, Blacksburg, VA, 24061, USA[gg]University of Warwick, Coventry CV4 7AL, United Kingdom[hh]Wright Laboratory, Department of Physics, Yale University, New Haven, CT, 06520, USAmicroboone_info@fnal.gov We describe a method used to calibrate the position- and time-dependent response of the MicroBooNE liquid argon time projection chamber anode wires to ionization particle energy loss. The method makes use of crossing cosmic-ray muons to partially correct anode wire signals for multiple effects as a function of time and position, including cross-connected TPC wires, space charge effects, electron attachment to impurities, diffusion, and recombination. The overall energy scale is then determined using fully-contained beam-induced muons originating and stopping in the active region of the detector. Using this method, we obtain an absolute energy scale uncertainty of 2% in data. We use stopping protons to further refine the relation between the measured charge and the energy loss for highly-ionizing particles. This data-driven detector calibration improves both the measurement of total deposited energy and particle identification based on energy loss per unit length as a function of residual range. As an example, the proton selection efficiency is increased by 2% after detector calibration. Calibration of the charge and energy loss per unit length of the MicroBooNE liquid argon time projection chamber using muons and protons
    [
    December 30, 2023
========================================================================================================================================



§ INTRODUCTION


Liquid argon time projection chamber (LArTPC) technology provides both particle tracking and energy loss per unit length reconstruction with high resolution. The basic working principle of a LArTPC neutrino detector is that neutrinos first interact with an argon nucleus and produce charged and neutral secondary particles.  The charged secondary particles travel through liquid argon and mainly lose their energy by ionizing and exciting argon atoms.  The ionization electrons travel under an applied electric field to a set of anode wire planes. The charge is measured on the anode wire planes in order to reconstruct particle trajectories and energies. Excited argon atoms also produce scintillation light, which is detected by the photon detectors. In this paper, we address the calibration of the charge collected at the anode wires for the MicroBooNE LArTPC <cit.>.

MicroBooNE is a 170-ton LArTPC with dimensions 2.56×2.33×10.36 m^3 (horizontal drift dimension×height×length). The nominal electric field inside the TPC is 0.273 kV/cm, which leads to a nominal electron drift velocity of 0.11 cm/μs. The drifted charge from particle interactions is read out in three planes with a plane spacing and a wire pitch of 3 mm. The 3456 collection plane wires are vertical and the angle between the induction and collection plane wires is 60 degrees. The readout window size is 4.8 ms and the ADC sampling rate is 2 MHz. Scintillation photons are observed by 32 photo-multipliers (PMTs) <cit.>. MicroBooNE started collecting neutrino events in Booster Neutrino Beam (BNB) at Fermilab in October 2015.

The first step of energy loss per unit length (dE/dx) reconstruction in LArTPC detectors involves the extraction of charge information from the signals (waveforms) on the anode plane wires. The MicroBooNE experiment uses several signal processing techniques <cit.> including noise filtering and signal deconvolution for charge extraction. However, the total charge extracted in this way normally does not equal the total charge produced from ionization for a number of reasons. The reasons are distortions in detector response due to cross-connected TPC channels <cit.>, space charge effects (SCE) <cit.>, electron attachment to impurities <cit.>, diffusion <cit.> and recombination <cit.>. To trace back  the exact amount of charge released from the original interaction, we have to correct for each of these effects starting from the ionization charge per unit length in the direction of the particle's travel, dQ/dx, as reconstructed from the signal
collected on the TPC wires. Throughout this paper dx refers to the length of a short segment of the particle track, which should not be confused with an infinitesimal distance in x (drift) direction of the TPC.




The data-driven detector calibration consists of several steps. In the first step we correct the position- and time-dependence of the detector response to ionization charge using data from cosmic ray muons (CR) which enter the detector from either the anode or cathode side and exit through the opposite face (crossing CR). This correction process is known as the dQ/dx calibration or equalization calibration of the detector. Once the detector response is corrected to be uniform throughout the TPC and in time, we determine a calibration for converting dQ/dx to dE/dx, by using stopping muons from neutrino interactions. Here stopping muons from neutrino interactions refers to the muons produced in the charged current interactions of neutrinos which ultimately decay inside the detector with a Bragg peak in their dE/dx profile. This process is known as the overall dE/dx calibration of the detector. In the last step, stopping protons from neutrino interactions are used to further refine the relation between the measured charge and the energy loss for highly ionizing particles.

The MicroBooNE coordinate system is shown in figure 1 (left). MicroBooNE uses a right-handed coordinate system in which the y axis is vertical, the x axis is horizontal, perpendicular to the anode and cathode planes, and the z=x× y is also horizontal and along the beam direction. The TPC signal formation in MicroBooNE is illustrated in figure 1 (right), where the induction wire planes are referred to as the “U” and “V” planes and the collection wire plane is referred to as the “Y” plane <cit.>. 
In this analysis we focus on calibrating only the collection wire plane of the detector, as this is the wire plane predominantly used for calorimetry in a LArTPC. The calibration of induction planes is generally more difficult because the response of the induction wires is highly dependent on the angle of the tracks relative to the wires. This is caused by the cancellation of overlapping induction signals for large-angle tracks due to the bipolar signal shape. New techniques are being developed to improve the reconstruction of the induction plane signals <cit.>, which will allow reliable dE/dx measurements on induction plane wires in the future. 

The two-step calibration procedure presented in this paper is similar to calibration techniques developed for other calorimeters such as the MINOS detectors <cit.>. 



§ DQ/DX CALIBRATION


 §.§ Introduction

The goal of the dQ/dx calibration procedure is to make the detector response uniform in space and over time. There are many effects that can produce a nonuniform detector response. The dominant effects are described in the following sections.



  §.§.§ Cross-connected TPC channels




Cross-connected wires, which can distort the electric field between wire planes <cit.>, affect the dQ/dx response of roughly 20% of the detector volume.
Figure <ref> shows the distribution of dQ/dx in the collection plane as a function of Y and Z coordinates in data. The highlighted diagonal region represents a “shorted-U” region where multiple U plane channels are shorted to one or more V plane channels. The highlighted vertical region represents a “shorted-Y” region where multiple Y plane channels are shorted to one or more V plane channels. 
By calibrating out these effects we recover the affected areas for useful physics.



  §.§.§ Space charge effects (SCE)

Since MicroBooNE is a surface-based detector, there is a significant flux of CR tracks in the detector volume. Because of this, there is a significant accumulation of slow-moving positive argon ions inside the detector, which is enough to distort the uniformity of the drift electric field. These distortions in the electric field inside the TPC have two significant effects:

  *  distortions in the magnitude of the drift electric field
compared to a uniform electric field, and

  *  distortions in the electric field direction compared to the nominal
direction parallel to the x-axis.

When the magnitude of the electric field is nonuniform, the recombination of electrons and ions is affected. This recombination effect <cit.> is sensitive to changes in the electric field. When the drift electric field is relatively low, recombination becomes dominant compared to the recombination at higher electric fields. Moreover, space charge effects can lead to spatial distortions in the trajectories of reconstructed particle tracks and electromagnetic showers. The positive argon ions built up in the detector tend to drag ionization electrons closer to the middle of the detector.
The space charge effects are expected to be stronger at the detector edges transverse to the drift direction. The cumulative effect of the space charge in track reconstruction leads to squeezing of the reconstructed track in transverse directions and bending towards the cathode. See Figure <ref>.



Figure <ref> shows how space charge effects <cit.> implemented in the MicroBooNE simulation change dQ/dx values over the entire drift distance of the MicroBooNE LArTPC. As seen in the right plot of figure <ref>, the dQ/dx values closer to the cathode are higher compared to that at the anode. The accumulation of positive ions causes the electric field magnitude closer to the cathode to be approximately 10% higher than at the anode. The higher field suppresses electron-ion recombination near the cathode. In addition, due to spatial distortions, tracks reconstructed closer to the cathode are bent and squeezed making reconstructed dx values smaller. Thus we observe a higher collected charge per unit track length closer to the cathode.



  §.§.§ Electron attachment to impurities

When a cloud of ionization electrons drifts to the anode, electronegative contaminants such as H2O and O2 can capture some of the drifting electrons, reducing the dQ/dx measured at the wire planes.
The capture rate is inversely dependent on the strength of the applied electric field, as at higher electric field magnitude, ionization electrons have a greater chance of making it to the anode plane before attaching to electronegative impurities.

Equation <ref> describes the depletion of a cloud of ionization
electrons due to capture by electronegative contaminants as the electrons drift towards the anode: 

    n_e(t_collected)/n_e(t_0) = exp(-(t_collected-t_0)/τ),

where t_0 is the start time, n_e(t_0) is the initial number of electrons at time t_0, n_e(t_collected) is the number of electrons collected by anode plane wires after a time t_collected and τ is the electron lifetime.

The electron lifetime depends on the amount of electronegative contaminants present in the medium, where higher electron lifetime corresponds to lower argon contamination levels. Figure <ref> shows the effect of electron attachment to impurities on dQ/dx values under different purity conditions. In high argon purity conditions space charge effects become dominant, causing dQ/dx to increase from anode to cathode. In low argon purity conditions dQ/dx value drops from the anode to the cathode due to electron attachment to impurities.


  §.§.§ Diffusion


Both longitudinal and transverse diffusion can be studied.
The cloud of ionization electrons tends to get smeared out in the direction of the drift because of longitudinal diffusion. This widens and lowers the pulse height of the signal at longer drift distance, which can lead to loss of signal if the pulse-height is below  reconstruction threshold. The effect of longitudinal diffusion is small in practice because the electronics shaping time is selected to be comparable to the diffusion smearing. In addition, transverse diffusion charge can allow charge to be detected on
multiple wires, which can smear the detected signal and reduce the resolution of charge reconstruction. This smearing is not corrected in the signal processing procedure used for this analysis. 
 



  §.§.§ Temporal variations

The detector response can change over time because of effects such as drift of the electronics gains, changes in temperature, different running conditions, etc. Note the time referred to here is the calendar time, not to be confused with drift time in the TPC. The most significant time-dependent change affecting MicroBooNE calibrations are changes in argon purity.



  §.§.§ dQ/dx equalization strategy


The general strategy of the dQ/dx calibration is to separate detector nonuniformities into yz plane, x (drift) direction, and calendar time variations, and calibrate them in sequence using CR. More details will be discussed in section <ref>. 

CR are the standard candle for uniform energy deposition throughout the detector. These CR have typical momenta in the range of 4 - 5 GeV, which results in a peak dE/dx of ∼ 1.7 MeV/cm . 


It should be noted that the calibration scheme described here is an approximation because contributions such as space charge effects cannot be completely factorized in separate yz plane and drift directions: for a given x value, the variations in the yz plane are different. Ideally we should carry out this calibration by voxelizing the detector into small three-dimensional cubes and derive a calibration constant for each cube. But limited statistics make the voxelation approach impractical. Moreover, spatial distortions introduced by space charge effects that impact track reconstruction would not be completely addressed by such a voxel-based calibration scheme.



 §.§ Data sample


The dQ/dx calibration of the detector is carried out using both data and Monte Carlo (MC) simulation. For MC simulation we use samples of simulated CORSIKA CR events <cit.> overlaid with neutrino interactions simulated with the  genie generator <cit.>. The simulation of particle propagation in the MicroBooNE detector is based on Geant4 <cit.>. The drift of ionization electrons to the wire planes and propagation of scintillation light to the PMTs are modeled in LArSoft <cit.>, which induces effects such as recombination, diffusion, space charge effects and electronics readout including cross-connected TPC channels. The electron lifetime is set to infinity in the simulation so there is no electron attachment to impurities. For data, we use data collected in the MicroBooNE detector with a trigger coincident with the beam and with a random trigger anti-coincident with the beam
from February to October of 2016, and a second period from September 2017 to March 2018.

In both MC and data we use the Pandora <cit.> pattern recognition program, combined with a Kalman filter implemented in the LArSoft framework to fit the CR tracks. 



 §.§ Event selection




For the dQ/dx calibration of the detector, we choose to use the anode-cathode crossing CR. The main reason is that the anode-cathode crossing CR span the entire drift distance, which make them a very valuable sample to study any effects that depend on the drift distance. 
Each CR in the selection is “tagged” with an initial time t_0, corresponding to the time that the electrons produced nearest the anode are detected. Crossing CR have a wide spatial but limited angular coverage. Figure <ref> shows an example of anode-cathode crossing CR in MicroBooNE data. Figure <ref> shows the x projected track length distribution of all CR. It has been estimated that ∼ 0.13% of the tracks are anode-cathode crossing CR. Out of a total of 69429 tracks from a single run 5709 taken on April 1, 2016, 87 are crossing CR and 52 tracks survived after angular cuts. 



The tracks used in the dQ/dx calibration satisfy the following selection criteria:

  *  250 cm < Track projected x length < 270 cm: Any track that satisfies this selection cut is considered to be an anode-cathode crossing CR. A 20 cm wide window is selected in order to account for imperfect reconstruction of track start and end positions.


  *  The absolute value of track angle θ_XZ (see figure <ref>) should not be in the range of 75^∘ to 105^∘: This selection ensures that we are excluding tracks which are nearly orthogonal to the wire planes, which can be mis-reconstructed <cit.>. See figure <ref>.


  *  The absolute value of track angle θ_YZ (see figure <ref>) should not be in the range of 80^∘ to 100^∘: This selection helps to remove tracks that are nearly parallel to the collection plane wires. See figure <ref>. The difference between MC and data in these figures is due to the fact that the MC simulation of detector effects such as space charge effects and the wire field response is not in a perfect agreement with data.


 §.§ Analysis method

The dQ/dx calibration of the collection plane response is carried out in three separate steps.

  *  Detector calibration in the yz plane. This step aims to remove the effects of space charge, cross-connected TPC channels, and transverse diffusion.
This step is performed using all crossing tracks occurring over a period of several months.


  *  Detector calibration in the drift (x) direction, which may be time-dependent. This step removes effects of electron attachment to impurities, space charge, and longitudinal diffusion. A separate correction is derived for each day of data.


  *  Detector calibration in time. This step removes any temporal variations in the overall detector response, and is only performed for data. There is currently no time dependence in the MC simulation.


  §.§.§ Method for detector calibration in the yz plane

We segment the yz plane into 5 cm by 5 cm cells.
Each three-dimensional (3D) hit of each crossing CR track 
is assigned to a cell based on the y and z coordinates of the hit as
reconstructed by Pandora, without correction for SCE-induced distortion.
For each cell containing more than five hits, we calculate a median 
charge deposit per unit length, [(dQ/dx)(y_i,z_i)]_Local, where y_i and z_i are
the central coordinates of each cell. Figure <ref> shows the distribution of number of hits for each yz plane cell in the collection plane for both MC and data respectively. The statistics is low in the peripheral cells because of the spatial distortion caused by space charge. However, this is not a problem since 1) yz correction factors are still quite uniform near the TPC edges as shown in figures <ref> and <ref>; 2) all the analyses will define a fiducial volume to exclude activities in the peripheral region. The use of the median dQ/dx mitigates effects of delta ray contamination, misreconstructed calorimetric information, and long tails of the dQ/dx distribution.
We also calculate a global median charge deposit per unit length, [dQ/dx]_Global using
all 3D hits of all crossing CR tracks. The correction factor for each cell is
defined as


    C(y_i,z_i) = [dQ/dx]_Global/[(dQ/dx)(y_i,z_i)]_Local .


If a cell has five or fewer space points, the correction factor for that cell is calculated from the C(y_i,z_i) values of nearby cells.



  §.§.§ Method for detector calibration in the drift direction

We segment the drift direction into bins of 10 cm in data and bins of 5 cm in MC.
These bin sizes are chosen to minimize statistical fluctuations in the correction factors.
Each 3D hit of each crossing CR track on a given day
is assigned to a bin based on the reconstructed x coordinate of the hit,
with x position corrected only for the t_0 of the CR as described in section
<ref>, and charge deposit corrected only by the C(y_i,z_i) factor.
For each bin containing more than five hits, we calculate a median 
charge deposit per unit length, [(dQ/dx)'(x_i)]_Local, where x_i is the
central position of the bin, and the prime on (dQ/dx)' denotes that the yz correction has been applied.
The correction factor for each bin in the drift direction is defined as


    C(x_i) = [(dQ/dx)']_Global/[(dQ/dx)'(x_i)]_Local ,


where [(dQ/dx)']_Global is the global median dQ/dx value after applying the yz correction to each 3D hit.

In data, we derive separate drift direction correction factors for each day as space charge effects and electron lifetime can change over time. To get a reliable set of correction factors, we need a data sample which contains enough statistics. For this purpose, we derive correction factors in the drift direction only for the days where there are more than 40 crossing tracks after all the angular cuts. In the case of lack of statistics for a given day, we use the drift direction correction factors derived for the immediate neighboring day. For MC we derive only a single set of drift direction correction factors for the collection plane as there is no time variation of space charge effects and the electron lifetime is set to infinity. 



  §.§.§ Time dependent calibration of the detector

After applying the yz and drift direction factors to each 3D hit, the 
daily global median value of dQ/dx of all hits on crossing CR tracks varies over
time, as shown in figure <ref>. 
Dips lasting several days reflect times when the argon purity was low, e.g., during recirculation pump maintenance.
We choose a reference value 
[dQ/dx]_Ref≡ 210 ADC/cm for the global median dQ/dx in the collection plane, since all global median dQ/dx values calculated for each day fluctuate around this value after yz- and drift-dependent correction. We then define the overall time dependent correction as


    C(t) = [dQ/dx]_Ref/[(dQ/dx)”(t)]_Global ,


where [dQ/dx]_Ref is the reference dQ/dx value in the collection plane, and [(dQ/dx)”(t)]_Global is the global median dQ/dx value in the collection plane after correcting for yz plane and drift direction irregularities for each day. 



Values of dQ/dx fully calibrated to the reference value are obtained by multiplying the reconstructed values by the factor C(y_i,z_i) C(x_i) C(t). For MC simulated data without temporal variation we define C(t)=1.































 §.§ Results


  §.§.§ dQ/dx calibration (Simulation)

The dQ/dx calibration of the detector in MC is performed using a single large dataset described in section <ref>. Figure <ref> shows the extracted yz correction factors C(y,z) and drift direction x correction factors C(x). The variation of the correction factors with respect to the drift distance could be attributed to electron attachment to impurities as well as space charge effects as shown in Figure <ref>. The bump in the last bin of the drift direction correction distribution is due to low statistics in that bin.



  §.§.§ dQ/dx calibration (data)

The yz and drift direction x correction factors are derived for two separate datasets as described in section <ref>.
In the first dataset, we have data combined from February 2016 to October 2016; and in the second dataset, we have data combined from September 2017 to March 2018. See figure <ref> 
 for the variation of the yz correction factors for data taken from February to October in 2016 and drift direction x correction factors for data taken on February 25, 2016 in the collection plane. The distributions for data taken from September 2017 to March 2018 are similar.













The drift direction corrections factors and time corrections are derived for the time period from February 2016 to October 2016 and the time period from September 2017 to March 2018. C(t) is calculated for each day. Due to intermittent periods of detector downtime, some days are skipped in this measurement. See Figure <ref> for the variation of the time-dependent correction factor over time. Figure <ref> shows the impact of the dQ/dx calibration on the dQ/dx distribution in both MC and data. 

In both data and MC, we see that the width of the central peak of
the dQ/dx distribution has narrowed by 2%, and the significant
excess of low dQ/dx values in the data before calibration has
been corrected.  This represents a significant improvement in the
measurement of dE/dx, particularly in regions requiring the
largest yz corrections, which is important for the measurement of
total deposited energy and particle identification based on
energy loss per unit length as a function of residual range.



 §.§ Angular dependence study

The calibration procedure described in this paper does not make any angular-dependent corrections because of the limited angular coverage of anode-cathode crossing CR. The angular distributions of the anode-cathode crossing CR tracks are shown in figure <ref>. In order to more closely study the angular dependence of the dQ/dx calibration, we divide the crossing CR into four angular subsamples:

  * 0^∘<|θ_XZ|<50^∘,

  * 50^∘<|θ_XZ|<75^∘,

  * 105^∘<|θ_XZ|<130^∘ and

  * 130^∘<|θ_XZ|<180^∘.

and derive the yz and x correction factors for each subsample following the same procedure described above. All of the other selection criteria remain the same.

For each yz cell in each subsample, we compute the relative difference of the yz correction factor for the subsample C_subsample(y_i,z_i) from the correction factor found using all crossing CR C(y_i,z_i).
Figure <ref> shows the distribution of the relative differences of the yz correction factors relative to the complete sample, [C_subsample(y_i, z_i)  - C(y_i,z_i)] / C(y_i,z_i),
for each angular subsample.

The shift in the mean of the distribution away from zero is indicative of the amount of bias observed in a given angular bin. 
The peaks of the four distributions all agree within 2%.



Figure <ref> shows the relative differences in the x correction factors between the different angular subsamples and the combined sample. The maximum variation is 1.5%. 



Based on the above studies, we assign a 2% systematic error on the yz correction and a 1.5% systematic error on the x correction for the collection plane. This study shows very little angular dependence for the sample selected for the dQ/dx calibration. However, because of the limited angular coverage of this sample, more studies are underway to further refine the angular dependence of the detector response. 




§ DE/DX CALIBRATION


 §.§ Introduction

After the dQ/dx correction is complete we can determine the absolute energy loss per unit length scale. To move from the calibrated dQ/dx to dE/dx, we use the modified box model <cit.> for recombination:


    (dE/dx)_calibrated=exp((dQ/dx)calibrated/Ccalβ^' Wion/ρℰ)-α/β^'/ρℰ,

withCcal is a calibration constant used to convert ADC values to number of electrons,
 
Wion    23.6 x 10-6 MeV/electron (work function of                                argon), 
ℰ    0.273 kV/cm (MicroBooNE drift electric field), 
ρ    1.38 g/cm3 (liquid argon density at a pressure 124.106 kPa), 
β'    0.212 (kV/cm)(g/cm2)/MeV, and 
α    0.93. 


The last two parameters were measured by the ArgoNeuT experiment <cit.> at an electric field of 0.481 kV/cm.
The modified box model is applied at MicroBooNE's electric field of 0.273 kV/cm. 

According to the above equation, precise determination of the calibration constant Ccal which translates ADC/cm to (number of electrons)/cm, is important in determining the absolute energy scale.  The charge Q is measured by the integral of the pulse of the deconvoluted anode response signal where the simulated electronics and field responses are removed.  Ccal is normalized so that the unit (ADC) corresponds to 200 electrons. In the case where the detector response is perfectly modeled ( e.g. in the simulation), the calibration constant Ccal should be exactly 1/200. Here the goal of the dE/dx calibration is to determine the calibration constant Ccal using stopping muons as the standard candle, because they have a well-understood energy loss profile. A sample of stopping muons can either be identified from cosmic data or neutrino interactions. The method we use for dE/dx calibration is discussed in section <ref>.



 §.§ Data sample

The dE/dx calibration of the detector is performed on both MC and data samples. For MC we use a CORSIKA CR simulation overlaid with neutrino interactions to obtain stopping muons induced by the neutrino interactions. In data we use datasets collected from February 2016 to October 2016 and from September 2017 to March 2018, with the final states of the inclusive charged current ν_μ interactions <cit.>,
to obtain candidate stopping muons from charged current interactions. In both MC and data we use the Pandora <cit.> pattern recognition program to reconstruct the neutrino and cosmic interactions. 




 §.§ Event selection

In the dE/dx calibration of the detector, we employ stopping muons induced by neutrino interactions. The motivation is to get a sample of muons from neutrino interactions with known muon energy loss per unit length so that a  portion of the track can be selected that corresponds closely to a minimum ionizing particle (MIP), for which the dE/dx is very well understood theoretically to better than 1%.
This sample contains the particles with the energy regime and the angular profile well representing the signatures of neutrino interactions.
Further, most events have forward-going muons, providing us with plenty hits in the collection plane.
The well known starting time of the stopping muons, determined by the arrival of the neutrinos produced in BNB, allows us to apply the dQ/dx calibration corrections.





We select the candidate muon tracks with the reconstructed start and end points 10 cm away from the TPC border in each dimension.
The tracks with the reconstructed angles, 75^∘ < θ_XZ < 105^∘ and 80^∘ < θ_YZ < 100^∘, are discarded to exclude tracks pointing into the wire planes and tracks nearly parallel to the collection plane wires, respectively.
To further ensure the quality of the reconstruction and the feature of minimum ionization, the tracks are required to have a minimum length of 150 cm.
In the MC sample, the candidate muon track has to match the true stopping muon coming from a neutrino interaction, and the distance between the start (end) point of the reconstructed track and that of the true muon has to be shorter than 5 cm.
In the data sample, the ratio of the median dQ/dx value in the last 5 cm segment to that of the first 5 cm segment of the candidate track is required to be greater than 1.5, ensuring the selected tracks to have the profile of a stopping muon.
The selection criterion based on the dQ/dx profile is determined from a MC simulation sample, as shown in figure <ref>. The non-muon contamination in the selected data sample is estimated to be less than 1% using MC.


































 §.§ Analysis method

The following explains the procedure to perform the dE/dx calibration both in MC and data:


  *  Start with a set of tracks that satisfy the selection in section <ref>.




  *  Segment the last 200 cm residual ranges (residual range is defined as the path of the particle before its stopping point) of tracks into 5 cm bins, which leads to a total of 40 residual range bins. If the track length is less than 200 cm, the entire track is used, otherwise, the last 200 cm segment with respect to the stopping point is used. 


  *  Loop over all the 3D hits of each selected track and fill the residual range
bins with dE/dx values derived using equation <ref> by setting the calibration constant Ccal to an arbitrary value. Moreover, we consider only the 3D hits which are separated by 0.3 to 0.4 cm; as the minimum spacing between 3D hits is set by the wire spacing of 0.3 cm, this selects forward going muons.











  *  After looping over all the tracks, fit each dE/dx distribution to a Landau-convoluted Gaussian function <cit.> to extract the most probable dE/dx value (MPV) representing that particular residual range bin.



  *  Plot the MPV dE/dx values against the kinetic energy of the particle; for each residual range bin, we take the middle bin value as the representative residual range value of that bin and transform that to kinetic energy. In this transformation we use a cubic spline fit to the tabulated values of  CSDA (continuous slowing down approximation) residual range  vs. kinetic energy for stopping muons in liquid argon <cit.>.


  *  Compare the curve generated in the previous step with the prediction made by the Landau-Vavilov function <cit.> in the kinetic energy range of the muons from 250 MeV to 450 MeV, which is in the region expected for MIPs, and a χ2 value is calculated using equation <ref>. The Landau-Vavilov function describes the energy loss probability distribution for a particle in a given medium. The most probable energy loss of a particle is dependent on the thickness of the energy absorber. To get the predicted MPV dE/dx, we set the absorber thickness to be 0.35 cm, the average value of dx.


    χ^2 = ∑((MPV(dE/dx)Predicted-MPV(dE/dx)Measured)^2/σ^2) ,

where we sum all the data points in the kinetic energy region of 250 MeV to 450 MeV.  

For MC, σ^2=δ_fit^2, where δfit is the uncertainty associated with the MPV dE/dx extracted by fitting a Landau-convoluted Gaussian function to the energy distribution. 

For data, σ^2=δ_fit^2 + δ_recombination^2, where δrecombination is the systematic error associated with recombination model uncertainties (here we take the uncertainty to be 1.5% of the measured dE/dx. This uncertainty is coming from the 1σ level uncertainties of the two free parameters β^' and α in the modified box model for recombination in  equation <ref>).

Note the MC sample was simulated and reconstructed with the same recombination model. Therefore, there is no systematic uncertainty associated with the recombination model for MC.


  *  Iterate through the steps iii to vi described above several times to generate tabulated set of data between a given calibration constant (Ccal) and χ2 value (here we change the parameter Ccal in each iteration and calculate the χ2 value. The total number of iterations in this procedure is the total number of data points in the plots shown in figure <ref>).


  *  In the final step, plot the χ2 values generated against calibration constants and fit that distribution with a second order polynomial to get the calibration constant which corresponds to the lowest χ2 value. 



  *  With the newly derived calibration constant (Ccal), calculate the dE/dx values using equation <ref> with calibrated dQ/dx values as the input. 


After the dE/dx calibration, we compare the newly derived dE/dx values with uncorrected dE/dx values to see the effects of applying the algorithm. (See section <ref>).















 §.§ Results

The absolute calibration of the detector is performed on both MC and data samples. Figure <ref> shows χ^2 - χ_Min^2 vs. the calibration constant Ccal for the collection plane. The best fit value is the one that gives the minimal χ^2. The uncertainty on the extracted value of Ccal is given by Δχ^2 = χ^2 - χ_Min^2 = 1. 



The calibration constants derived for the collection plane are shown in table <ref>. The MC calibration constant uncertainty is statistical only. The data calibration constant uncertainty includes both the statistical uncertainty and the recombination uncertainty. The difference between data and MC calibration constants is due to the fact that the MC simulation of detector response is not in a perfect agreement with data.



Figure <ref> shows the comparison between the predicted most probable energy loss <cit.> with the fitted most probable energy loss using the calibration constants shown in table <ref> for stopping muons both in MC and data. Here the disagreement between the fitted distribution and theory at lower kinetic energies in data is mainly due to the recombination model uncertainty. Results can be applied to all time ranges, since temporal variations are taken into account by the relative dQ/dx corrections. The first fitted point in the distributions corresponds to the end of the tracks. Most of the time the end of the track can happen in between two wires. Therefore calorimetric information for the end points of tracks is not properly reconstructed. The significant deviation between prediction and our fitted distribution for MC in the first bin is due to this reason.



  
The absolute energy-scale calibration can be validated in a data-driven way by comparing the range-based energy to that obtained via calorimetry for selected stopping muon candidates in data and MC. For each selected stopping muon track, the kinetic energy calculated by range and by calorimetry is computed, and the relative difference between these methods is shown in figure <ref>. Here, in calculating the kinetic energies, the entire track is used. The agreement between the two, approximately 2% for data (1% in MC), gives confidence in the proper absolute energy-scale calibration to within this level of agreement, and provides a closure test for this calibration procedure. The tails on the positive side are caused by contributions of other particles to the muon track near the neutrino interaction vertex. A MC study shows that after the calibration procedure, the bias in (E_calo-E_range)/E_range reduces from -3.2% to -1.2% and the resolution improves from 8.1% to 7.7%.



§ ENERGY CORRECTION WITH PROTONS


 §.§ Introduction


After the nonuniformities in detector response have been corrected for and the absolute energy scale has been determined as described in the previous sections, we perform studies using a proton sample to investigate the effect of electron-ion recombination.

The calibration of dQ/dx and dE/dx presented in the previous sections sets the correct scale for energy loss per unit length for MIPs, for which dE/dx remains nearly a constant. As particles approach the end of their range, they become increasingly ionizing. Electron-ion recombination introduces a non-linear relationship between dQ/dx and dE/dx, which has been extensively studied <cit.>. In this study, 
we use stopping protons to measure the non-linear relationship between dQ/dx and dE/dx. 
The variation of dE/dx with the proton energy allows us to measure the recombination effect for a wide range of stopping powers up to about 20 MeV/cm. 











The detector simulation used by MicroBooNE employs the modified box model by default to simulate the recombination effect, with the parameters obtained by the ArgoNeuT experiment <cit.>. The charge deposition dQ/dx at each Geant4 step is calculated using the energy loss dE/dx following:


    dQ/dx(e/cm) =ln(dE/dxβ^'/ρℰ+α)/β^'/ρℰ W_ion ,

with the parameters Wion, ℰ, ρ, β' and α, already described in the inverse of this equation, eq. <ref>.




In this analysis, we compare the proton dQ/dx vs.dE/dx measurements with two recombination models: the modified box model as described in eq. <ref> and Birks' law, following eq. <ref> with A_B=0.8 and k=0.0486 (kV/cm)(g/cm^2)/MeV, as measured by the ICARUS collaboration <cit.>: 

    dQ/dx(e/cm)=A_B/W_ion(dE/dx/1+k/ρℰdE/dx).

The two models are parameterized by ArgoNeuT and ICARUS at different electric fields from the operating electric field at MicroBooNE. Therefore, it is interesting to compare both models with the MicroBooNE measurement to study the dependence of recombination on electric field. 


To perform these studies, we use the calibrated dQ/dx information and calculate dE/dx using an empirical fit to Bethe-Bloch energy loss using the range information,  (dE/dx)_range. The (dE/dx)_range values for protons are calculated using the residual range information per hit, as described in Ref. <cit.> according to eq. <ref>
    (dE/dx)_range=A· R^b

where R is the residual range and A and b are empirical constants (the constant A in eq. <ref> should not be confused with A_B from Birks' law in eq. <ref>). We use the values obtained by ArgoNeuT, A = 17 MeV/cm^1-b and b = -0.42 in the calculation. The dE/dx values calculated in the way are independent of the charge measurement, and therefore can be compared to the calibrated dQ/dx values to derive recombination model parameters.



 §.§ Event selection


In this study we use protons emerging from neutrino interactions with argon. We choose a sample of events with one muon and two protons. This sample is referred to as 1μ2p. Such events are produced almost exclusively by the neutrino interactions, rather than the cosmic ray interactions, as verified using data taken when the BNB beam was off. Having three tracks sharing the same vertex improves the resolution of vertex and track reconstruction. This sample provides a highly pure proton sample with well reconstructed track trajectories, which is ideal for the recombination measurement. 




The event selection is based on a set of cuts designed to select charged-current ν_μ interactions producing 1μ2p final state in the TPC fiducial volume (FV). We take advantage of the very specific number of charged particles and the absence of e^±/π^0 in the final state, which yields a sample of well reconstructed tracks.

For the 1μ2p events, the signature is three tracks originating from the same vertex, one consistent with a muon and the other two consistent with protons. The selection criteria are described as follows:


  *  Beam window: 32 PMTs are used in MicroBooNE to record scintillation light information from particles traveling inside the TPC. In order to select neutrino interactions, we require there should be a scintillation light signal detected in the PMTs coincident with the neutrino beam spill.  


  *  Vertex: The reconstructed neutrino vertex is required to be in the FV. The FV is defined as a rectangular box shaped as follows:  the boundary from the anode plane and  the  cathode  plane is 10 cm,  the boundary from the upstream and downstream ends is 10 cm, and the boundary from the top and bottom of the TPC is 20 cm.

It is the same FV as used by the MicroBooNE physics analysis reported in <cit.>.

The more stringent requirement in the vertical direction ensures a higher CR background removal.





  *  Three tracks cut: Only three tracks are connected to the vertex within a 5 cm radius. 
If more than one neutrino candidate is present in an event, the one where the sum of the distance between the start position of each track to each other is smallest is selected. This suppresses cosmic ray background mimicking neutrino interactions.


  *  Proton containment: The leading particle, defined as the longest track, is taken as the muon candidate; the remaining two tracks are taken as proton candidates. These two proton candidates are required to be contained in the TPC FV to aid in particle identification (PID).


  *  Minimum number of hits: For each track we require a minimum number of 10 hits in total, summing all three planes. In order to have the most reliable calorimetry information, we also require at least 5 hits in the collection plane. This requirement is applied to all tracks in a 1μ2p event.


  *  PID cut: We use the χ^2 test with respect to a proton hypothesis in the dE/dx vs. residual range curves, from Bethe-Bloch predictions, to discriminate MIP with respect to non-MIP particles. The predicted proton dE/dx vs. residual range curve is calculated using the Geant4 simulation, where the residual range is sampled with a step size of 0.08 cm.


For each selected track (either in MC or data), we reconstruct dE/dx per hit in the collection plane as described in section <ref> and calculate the χ^2 between the reconstructed and the Geant4-determined mean dE/dx at each given residual range bin. These values are then summed for all hits on the selected track. 

We exclude the first and last hits on the track from this calculation because of the uncertainty in the dx calculation. This is due to the fact that the dx calculation for the first and last hits may be wrong since the exact start and end positions of the track between the wires are unknown. 
It also mitigates the complications caused by energy deposition overlap at the interaction vertex from other particles ( e.g. short protons) produced by the neutrino interaction.


The value obtained is then normalized by the number of degrees of freedom (d.o.f.), which corresponds to the number of hits in the collection plane. This calculation can be determined for several particle assumptions, but in this analysis we exclusively use the proton hypothesis. 


    PID = χ^2_proton/d.o.f. = ∑_Track hits((dE/dxmeasured- dE/dxtheory)/σ_dE/dx_)^2/d.o.f.,



where σ_dE/dx is the estimated resolution in measuring a fixed value of dE/dx, the subscript  measured and  theory corresponds to the dE/dx measured at the wires using the charge deposition information and the predicted by Geant4 respectively. Values used for resolution in dE/dx are taken from studies performed by the ArgoNeuT Collaboration. This does not take into account the 3 mm wire spacing in MicroBooNE (compared to the 4 mm wire spacing in ArgoNeuT), which leads to an overestimated resolution in dE/dx when compared to the resolution we obtain from the MicroBooNE MC studies. We prefer to use the values from ArgoNeuT as a conservative approach since the recombination model is not tuned based on MicroBooNE data at this stage yet. 

In order to select a sample of highly ionizing particles, an initial guess of recombination correction is needed, for which we employ the modified box model with the parameters measured by ArgoNeuT. This model is sufficiently close to the final measurement to ensure a good separation of protons and muons.

Taking into account both the particle discrimination (muon  vs. proton) and the data/MC shape agreement, we apply a PID requirement that identifies protons as the tracks with χ^2proton/d.o.f.<88, and muons as tracks with χ^2proton/d.o.f>88. The values are optimized by maximizing the efficiency and purity for the proton identification.


 



Table <ref> shows the purity and efficiency of the muon and proton selections, estimated from the MC sample. Purity is defined as the quantity of particles from an specific species which has been correctly identified, divided by the total number of particles of that species in the sample. For comparison, the same values are shown before applying the equalization procedure described in section <ref>. Note that in the MC simulation, we reconstruct dE/dx using the recombination parameters obtained by ArgoNeuT, as they are the values used in the simulation. 
The equalization calibration improves both the purity and efficiency of the muon and proton selections. 

  *  Shower veto: We remove events with a shower reconstructed within 15 cm of the vertex, both to reduce π^0 contamination and to reduce background caused by tracks misreconstructed as showers.


In this analysis, a total of 226 proton candidates are selected from the neutrino data.



 §.§ Analysis method


With the proton sample described we study two different recombination models:
the modified box model, eq. <ref>, and Birks' law, eq. <ref>. The comparison of these two models with the MicroBooNE proton dataset (described in section <ref>) is shown in figure <ref>. In this figure, the measured dQ/dx, after calibration, is shown with respect to the theoretically predicted energy deposit (dE/dx)_range, computed from the proton residual range following eq. <ref>, for both data (left) and MC (right). Both model predictions are calculated at the nominal electric field of 0.273 kV/cm for MicroBooNE.





We observe a good agreement of data with Birks' law and the modified box model at low dE/dx in figure *fig:protons2D_defaultBirksBox_data, as expected. Low dE/dx corresponds to the MIP region, which has been calibrated with muons previously. At higher dE/dx, data shows a lower recombination effect than the one predicted by the modified box model. On the other hand, Birks' law agrees much better with data at higher ionization. This is not surprising since Birks' law was tuned with ICARUS data taken at drift fields of 0.2, 0.35 and 0.5 kV/cm while the modified box model was tuned with ArgoNeuT data taken at a single drift field of 0.481 kV/cm. As expected, the MC distribution agrees well with the modified box model using the ArgoNeuT parameters, shown as the black curve in figure *fig:protons2D_defaultBirksBox_mc, since it is the same model used for the simulation.

Because the parameters β^' in the modified box model and k in Birks' law are coupled to the electric field ℰ in eqs <ref> and <ref>, it is the quantity β^'/ℰ or k/ℰ that we are able to measure. The accumulation of positive ions produced by CR in the liquid argon leads to a nonuniform electric field. 
Even though the nonuniformities are largely reduced by the equalization procedure described in section <ref>, there could still be residual nonuniformities in the electric field distribution.
The discrepancies shown in figure *fig:protons2D_defaultBirksBox_data could be due to the fact that the recombination model parameters derived at other experiments are not applicable at the MicroBooNE electric field, or that the average electric field at MicroBooNE is different from the nominal electric field of 0.273 kV/cm.
For the recombination measurements reported in this paper, we adopt an effective approach by fixing the electric field at 0.273 kV/cm and tune the other model parameters to improve the agreement with the MicroBooNE proton data.




We obtain the  effective recombination parameters through a χ^2 fit. We construct a χ^2 function based on  the mean dQ/dx value per each dE/dx bin in data and MC. We then minimize the χ^2 function to obtain the best fit parameters (α and β^' for the modified box model, and A_B and k for Birks' law). For data, dQ/dx in (e/cm) is calculated by dividing the calibrated dQ/dx in (ADC/cm), as described in section <ref>, by the calibration constant C_cal, as described in section <ref>;  i.e.
    dQ/dx(e/cm)^data = dQ/dx(ADC/cm)^data/C_cal,


and dE/dx is calculated using the residual range following eq. <ref>. For MC, dE/dx is also calculated using residual range and dQ/dx is converted from (dE/dx)_range using the recombination model (the modified box model eq. <ref> or Birks' law eq. <ref>). The χ^2 is defined as:









    χ^2 = ∑_dE/dx bins(dQ/dx(e/cm)^data-dQ/dx(α,β^',(dE/dx)range)^MC/√((σ_(dQ/dx)^data^2+σ_(dQ/dx)^MC^2)))^2,

where σ_(dQ/dx)^data and σ_(dQ/dx)^MC are the statistical errors from data and MC respectively.



 §.§ Results


We perform the χ^2 fit using both the modified box model and Birks' law. The best fit results are summarized in table <ref>, including the statistical uncertainties from the fit. To obtain these results, we assume an electric field of 0.273 kV/cm, which is the nominal electric field in MicroBooNE. Under this assumption, the parameters that are not coupled to the electric field, α and A_B, are consistent with the measurements from ArgoNeuT and ICARUS, while the parameters that are coupled to the electric field, β^' and k, are in a disagreement with the measurements from ArgoNeuT and ICARUS. The effect of the obtained result can be clearly observed in figure <ref>, which is the same as figure *fig:protons2D_defaultBirksBox_data now showing also the model curves using  effective recombination parameters. We observe a better agreement of the data to the two models with the new parameters.






The  effective recombination parameters are validated using the MicroBooNE ν_μ CC inclusive data sample. The calculated dE/dx values using the  effective recombination parameters are compared with the Bethe-Bloch prediction with respect to residual range in figure *fig:allprotons_dEdxvsrr_OnBeam_ub. As a comparison, we show in figure *fig:allprotons_dEdxvsrr_OnBeam_argo the reconstructed dE/dx using the modified box model with the ArgoNeuT parametrization. Improvement between the calculated dE/dx distribution and the Bethe-Bloch prediction for protons can be seen when using the  effective recombination parameters (right), with respect to the same model using the parameters obtained by ArgoNeuT (left). 










The impact of the full calibration procedure on proton dE/dx reconstruction is shown in figure <ref>. For figures *fig:comparison_protondEdx_a, *fig:comparison_protondEdx_b and *fig:comparison_protondEdx_c, the distributions show the dE/dx per hit of the proton candidates in data compared to those from the simulation. Simulated dE/dx is calculated using calibrated dQ/dx and the modified box model with ArgoNeuT parameters. Data dE/dx in figures *fig:comparison_protondEdx_a is calculated using uncalibrated dQ/dx and the modified box model with ArgoNeuT parameters. Data dE/dx in figures *fig:comparison_protondEdx_b is calculated using calibrated dQ/dx and the modified box model with ArgoNeuT parameters. The equalization calibration clearly improves the dE/dx reconstruction, but there is a discrepancy at high dE/dx values. Data dE/dx in figures *fig:comparison_protondEdx_c plot is calculated using calibrated dQ/dx and the modified box model with the MicroBooNE  effective recombination parameters as reported in table <ref>. The agreement is improved for high dE/dx values.  The agreement for low dE/dx values gets slightly worse. This will be improved by retuning the dE/dx scale using stopping muons, as described in section <ref>, after applying the  effective recombination parameters in the future. 

We also show the dE/dx per hit distributions for selected protons and muons from data in figures *fig:comparison_protondEdx_d. In this plot, the equalization calibration and the  effective recombination correction are applied. We can clearly see from the figure that protons and muons can be well identified and have different widths in dE/dx as expected. 



§ CONCLUSIONS

In this paper, we describe for the first time a method that calibrates the calorimetric response of the ionization signal in a LArTPC using muons. We use crossing muons to remove spatial and temporal variations in the detector response. We use stopping muons to determine the absolute energy loss per unit length scale for minimum ionizing particles.  We use stopping protons to further refine the relation between the measured charge and energy loss for highly-ionizing particles. The most probable value of calibrated energy loss agrees with expectation. The calibrated energy loss per unit length is used in particle identification and energy reconstruction. This method can be followed by other LArTPC experiments to calibrate their detectors. 



§ ACKNOWLEDGEMENTS

This document was prepared by the MicroBooNE collaboration using the
resources of the Fermi National Accelerator Laboratory (Fermilab), a
U.S. Department of Energy, Office of Science, HEP User Facility.
Fermilab is managed by Fermi Research Alliance, LLC (FRA), acting
under Contract No. DE-AC02-07CH11359.  MicroBooNE is supported by the
following: the U.S. Department of Energy, Office of Science, Offices
of High Energy Physics and Nuclear Physics; the U.S. National Science
Foundation; the Swiss National Science Foundation; the Science and
Technology Facilities Council (STFC), part of the United Kingdom Research and Innovation; and The Royal Society (United Kingdom).  Additional support for the laser calibration system and cosmic ray tagger was provided by the Albert
Einstein Center for Fundamental Physics, Bern, Switzerland.

99uBDetJINST
R. Acciarri et al. (MicroBooNE Collaboration), “Design and construction of the MicroBooNE detector”, JINST 12 (2017), P02017 [arXiv:1612.05824].


Conrad:2015xta 
  J. Conrad  et al. [MicroBooNE Collaboration],
  “The Photomultiplier Tube Calibration System of the MicroBooNE Experiment,”
  JINST  10, no. 06, T06001 (2015)
  doi:10.1088/1748-0221/10/06/T06001
  [arXiv:1502.04159 [physics.ins-det]].
  
  

uBsignalJINST
R. Acciarri et al. (MicroBooNE Collaboration), “Ionization Electron Signal Processing in Single Phase LAr TPCs I: Algorithm Description and Quantitative Evaluation with MicroBooNE Simulation”, JINST 13 (2018), P07006 [arXiv:1802.08709].

uBsignalJINST_2
R. Acciarri et al. (MicroBooNE Collaboration), “Ionization Electron Signal Processing in Single Phase LAr TPCs II: Data/Simulation Comparison and Performance in MicroBooNE”, JINST 13 (2018), P07007 [arXiv:1804.02583].

uBnoiseJINSTpre
R. Acciarri et al. (MicroBooNE Collaboration), “Noise Characterization and Filtering in the MicroBooNE Liquid Argon TPC”, JINST 12 (2017) P08003 [arXiv:1705.07341].

sc2
M. Mooney, “The MicroBooNE Experiment and the impact of Space Charge Effects”, Proceeding of DPF 2015 Meeting of APS, [arXiv:1511.01563].

ICARUS1O2
E. Buckley et al.,“A study of ionization electrons drifting over large distances in liquid argon”, Nucl. Instrum. and Methods in Physics Research A 275, 364-372 (1989). 

ICARUS2O2
A. Bettini et al.,“A study of the factors affecting the electron lifetime in ultra-pure liquid argon”, Nucl. Instrum. and Methods in Physics Research A 305, 177-186 (1991).

ICARUSpm
S. Amoruso et al. (ICARUS Collaboration), “Analysis of the liquid argon purity in the ICARUS T600 TPC”, Nucl. Instrum. and Methods in Physics Research A 516, 68-79 (2014).

longbo
C. Bromberg et al. (LongBo Collaboration), “Design and operation of LongBo: a 2 m long drift liquid argon TPC”, JINST 10 (2015) P07015.

argoneut
C. Anderson et al. (ArgoNeuT Collaboration), “The ArgoNeut detector in the Numi low-energy beam at Fermilab", JINST 7 (2012) P10019. 

diffusion2015Li
Y. Li et al., Measurement of Longitudinal Electron Diffusion in Liquid Argon, arXiv:1508.07059, Nucl. Instrum. Meth. A816 P160-170 (2016)

recomb
R. Acciarri et al. (ArgoNeuT Collaboration), “A study of electron recombination using highly ionizing particles in the ArgoNeuT Liquid Argon TPC”, JINST 8 (2013) P08005. 

recombICARUS
S. Amoruso et al.(ICARUS Collaboration), “Study of electron recombination in liquid argon with the ICARUS TPC", Nucl. Instrum. Meth. A523 P275-286 (2004)

Hartnell:2005uq 
  J. J. Hartnell,
  “Measurement of the Calorimetric Energy Scale in MINOS,”
  FERMILAB-THESIS-2005-51,
  doi:10.2172/875528.
 
asic_IEEE_2011
G. D. Germonimo, et al., “Front-End ASIC for a Liquid Argon TPC", IEEE Trans. Nucl. Sci 58 (2011) 1376-85.  

Heck:1998vt 
  D. Heck, G. Schatz, T. Thouw, J. Knapp and J. N. Capdevielle,
  FZKA-6019.
  

Andreopoulos:2009rq 
  C. Andreopoulos  et al.,
  
  Nucl. Instrum. Meth. A  614, 87 (2010).
  
  

Agostinelli:2002hh 
  S. Agostinelli  et al. [GEANT4 Collaboration],
  “GEANT4: A Simulation toolkit,”
  Nucl. Instrum. Meth. A  506, 250 (2003).
  doi:10.1016/S0168-9002(03)01368-8
  
  

larsoft
  E. L. Snider and G. Petrillo,
  “LArSoft: Toolkit for Simulation, Reconstruction and Analysis of Liquid Argon TPC Neutrino Detectors,”
  J. Phys. Conf. Ser.  898, no. 4, 042057 (2017).
  doi:10.1088/1742-6596/898/4/042057
  
  

  
pandora
R. Acciarri et al. (MicroBooNE Collaboration), “The Pandora Multi-Algorithm Approach to Automated Pattern Recognition of Cosmic Ray Muon and Neutrino Events in the MicroBooNE Detector”, arXiv:1708.03135, Eur. Phys. J. C78, 1, 82 (2018).

numu_cc_filt
R. Acciarri et al. (MicroBooNE Collaboration),“First Measurement of Muon Neutrino Charged Current Neutral Pion Production on Argon with the MicroBooNE LAr TPC", Phys. Rev. D99 (2019) 091102(R) [arXiv:1811.02700].

numu_cc_Marco
R. Acciarri et al. (MicroBooNE Collaboration),“First Measurement of Inclusive Muon Neutrino Charged Current Differential Cross Sections on Argon at Enu ∼0.8 GeV with the MicroBooNE Detector", Phys. Rev. Lett. 123, 131801 (2019) [arXiv:1905.09694].

lgfit
Convoluted Landau and Gaussian Fitting Function,
<https://root.cern.ch/root/html/tutorials/fit/langaus.C.html>csda_res_range
Donald E. Groom, Nikolai V. Mokhov, Sergei I. Striganov, Muon Stopping Power and Range Tables 10 MeV - 100 TeV, Atomic Data and Nuclear Data Tables, Vol. 76, No. 2, July 2001.

landau_vavilov


H. Bichsel,
  “Straggling in Thin Silicon Detectors,”
  Rev. Mod. Phys.  60, 663 (1988).
  doi:10.1103/RevModPhys.60.663
  
  

Thomas:1987zz
J. Thomas and D. A. Imel, “Recombination of electron-ion pairs in liquid argon and liquid xenon”, Phys. Rev. A36, 614-616 (1987).
  
Birks:1964zz
  John B. Birks, “The Theory and practice of scintillation counting”,
  Pergamon Press (1964).