Constraints on millicharged particles with low threshold germanium detectors at Kuo-Sheng Reactor Neutrino Laboratory

Relativistic millicharged particles ($\chi_q$) have been proposed in various extensions to the Standard Model of particle physics. We consider the scenarios where they are produced at nuclear reactor core and via interactions of cosmic-rays with the earth's atmosphere. Millicharged particles could also be candidates for dark matter, and become relativistic through acceleration by supernova explosion shock waves. The atomic ionization cross section of $\chi_q$ with matter are derived with the equivalent photon approximation. Smoking-gun signatures with significant enhancement in the differential cross section are identified. New limits on the mass and charge of $\chi_q$ are derived, using data taken with a point-contact germanium detector with 500g mass functioning at an energy threshold of 300~eV at the Kuo-Sheng Reactor Neutrino Laboratory.


I. INTRODUCTION
The origin of electric charge quantization is one of the profound intriguing mysteries of the nature [1]. Studies are made to explore the theoretical foundation of electric charge quantization in Kaluza-Klein higher dimensional theory [2], the existence of magnetic monopoles [3], and grand unified theories [4]. Despite consistency of charge quantization with all experimental data, the absence of evidence for magnetic monopoles and grand unification continues to motivate searches for the existence of nonquantized charged particles commonly known as millicharged particles and denoted by χ q with mass m χ q in this article [5]. Such particles can be naturally obtained via including an extra Abelian gauge U HS ð1Þ (the subscript denotes "hidden sector") into the standard model (SM) gauge groups [6,7]. The SM particles are not charged under this new gauge group, while the χ q under U HS ð1Þ acquire small electric charge (δe 0 ) due to the kinetic mixing of SM photon and HS dark photon, where δ is the charge fraction of χ q and e 0 is the standard electron charge. The parameter space defined by ðm χ q ; δÞ is strongly constrained by astrophysical observations [1,5,8,9], cosmic microwave background [10,11], big bang nucleosynthesis [12,13], and direct laboratory experiments [14][15][16].
The theme of this article is to report new direct laboratory limits on χ q using data acquired by an n-type point-contact germanium (PCGe) detector with sub-keV sensitivity and excellent energy resolution. The PCGe technology [17] is being used in the studies of coherent neutrino nucleus elastic scattering, as well as searches of "light" WIMPs [18][19][20][21] and neutrinoless double beta decay.
This article is organized as follows. Three cases of χ q are selected for investigations, whose origins and flux estimations are given in Sec. II. The interaction of χ q with matter and the derivation of cross sections are discussed in Sec. III. Physics constraints on χ q derived from data taken with the PCGe at the Kuo-Sheng Reactor Neutrino Laboratory (KSNL) [17,[22][23][24][25] are presented in Sec. IV.

II. SOURCES OF MILLICHARGED PARTICLES
A. Nuclear reactor χ q Nuclear reactor cores are powerful sources of γ rays and therefore have been used for searches of axions [26,27], dark photons [28,29], and millicharged neutrinos [30]. Approximately half of the γ-radioactivity are of "prompt" origins directly from highly excited fission fragments, while the rest are from radiative de-excitation of the daughter nuclei, inelastic neutron scattering and capture of neutron by core materials [26]. The total prompt γ-ray spectrum, in the units of MeV −1 s −1 , for FRJ-1 (Merlin) reactor core is approximated (E γ ≥ 200 keV) by [31] This spectrum has been adopted in earlier reactor axion [26] and dark photon searches [28]. The light-χ q can be produced through Compton-like processes, where γ rays of OðMeVÞ energy scatter off electrons in the reactor core [32]. The differential production cross section of χ q -χ q via a Compton-like process, as depicted schematically in Fig. 1, can be estimated from well-known lepton-pair production process in the field of an electron. To adapt the lepton-pair production formula in the case of χ q production, the χ q -χ q vertex is parameterized by δe and the lepton mass is replaced by m χ q , where m e is a mass of the electron; E χ q and E γ are, respectively, energies of χ q and γ rays, while Eχ q ¼ E γ − E χ q . The differential χ q -flux (ϕ χ q ) is determined from the convolution of the reactor γ-ray spectrum and differential production cross section normalized by the total cross section (σ tot ), where R is a distance of detector from the center of reactor core. The factor 2 in Eq. (3) comes from the fact that χ q particles produce in pairs and both can interact at the detectors. Although the reactor core produces γ rays of energy up to 10 MeV, the 1-5 MeV interval of the reactor γ-ray spectrum is expected to dominate the χ q -χ q pair production process of Eq. (2) [32]. The Compton scattering cross section (σ c ) of γ rays in this energy interval is dominated over atomic photoelectric absorption and the pair production cross section, even for high Z-materials. Therefore, σ tot ≃ σ c is a reasonable approximation for the 1-5 MeV energy interval. A typical reactor core with 1 GW thermal power produces O (10 20 ) photons per second at OðMeVÞ energies. These γ rays can interact with electrons producing 2.3 × 10 18 × δ 2 number of χ q per second at the reactor core.

B. Atmospheric χ q
High-energy cosmic rays can produce relativistic-χ q when they interact with nucleus in the Earth's atmosphere [33,34]. Bremsstrahlung and multiperipheral (which include the interaction of two virtual photons) are main production mechanisms of χ q in the collision of charged particles and their cross sections are proportional to δ 2 α 4 Z 2 and δ 4 α 4 Z 2 , respectively [33]. The energetic cosmic rays may produce χ q with masses inaccessible to current particle accelerators.
The investigation of cosmic rays can provide constraints on χ q via interactions between χ q and detectors. The energy loss of χ q through excitation and ionization is proportional to δ 2 , which is much lower than the minimum ionizing particles of unit charge under the similar conditions. These exotic χ q are therefore also called lightly ionizing particles and fractional charged particles in the literature [34][35][36][37]. Since these particles deposit little energy at detectors, the high resolution and sub-keV threshold of PCGe make it an ideal detector for their searches.
The mass range of cosmogenic produced χ q is unknown due to unknown production conditions. Therefore, the experimental sensitivity is usually expressed in terms of the integral incoming flux (I χ q ) in the units of cm −2 s −1 sr −1 as a function of δ [34].

C. Dark cosmic rays-DM-χ q
If self-interaction between dark matter (DM) is mediated by massless dark photon, the gauge invariance mandates that DM would be multicomponent [38]. Two oppositely FIG. 1. The production of χ q -χ q via Compton-like mechanism based on the kinetic mixing of dark photon with the SM photon. charged massive fermion interacting through a massless Uð1Þ HS gauge boson would give rise to bound states which, implies atomic DM [39]. The atomic dark sector also requires the existence of a dark matter-antimatter asymmetry in the early universe to set the relic abundance of DM. The self-interaction of dark matter is an interesting scenario which has the potentials to explain and solve the small scale issues of the cold collisionless dark matter paradigm without modifying the physics at large scale [40,41]. The dark sector particles may have both neutral and ionized components. The latter component may be due to incomplete recombination of primordial DM gas [42] and reionization by sources such as starlight and supernova explosions [43], which can efficiently overcome the binding energy of dark atoms. The ionized components of DM could be accelerated analogous to standard cosmic rays [44].
The diffuse-shock acceleration is widely accepted as being responsible for accelerating the standard cosmic rays at the blast waves of supernova remnant [45]. Ambient charged particles captured by supernova shock wave fronts can be accelerated to high energies during the lifetime of a supernova remnant. The key feature of this theory [46], usually known as first order Fermi acceleration, is that the fractional energy gain by charged particles via each crossing of supersonic shock front and scattered back from turbulent magnetic fields associated with the shock front is proportional to the first order of shock velocity. This mechanism also naturally explains the commonly observed power law cosmic-ray energy spectrum. The maximum energy attained by χ q of charge fraction δ is simply product of rate of energy gain and time spent in shock. The upper limit of maximum energy can be obtained by assuming a minimum diffusion length equal to gyroradius of χ q A typical supernova can accelerate a proton to 10 6 GeV energy level with the following parameters such as the shock wave speed V S ∼ 10 4 km s −1 , total acceleration time t A ≈ 10 3 years, and magnetic field B ∼ 10 −10 T. If the gyroradius of the χ q as DM candidates is smaller than the length of shock waves, they can be efficiently accelerated by the Fermi mechanism. Therefore, accelerated DM-χ q may accompany ordinary cosmic rays arriving on the Earth. The precise prediction of the DM-χ q "dark cosmic-ray" flux requires knowledge of the injection spectrum. An estimation is obtained in Ref. [44] by using results of ion acceleration in shocks. The energy spectrum in the units of cm −2 s −1 sr −1 GeV −1 , follows a simple power law, where α ¼ 2.7 is the differential spectral index of the ordinary cosmic-ray flux. The maximal attainable energy of DM-χ q would be suppressed by a factor of δ from the maximum energy of the cosmic-ray proton.

III. ATOMIC IONIZATION WITH MILLICHARGED PARTICLES
The millicharged χ q 's are relativistic and can interact electromagnetically with matter via atomic ionization through the t-channel process depicted in Fig. 2(a). The cross sections for the analogous atomic ionization processes due to neutrino electromagnetic interactions were recently derived [30,47]. The case of dark matter sterile neutrino as a nonrelativistic probe was also investigated [48]. This work expands these studies to one with having a relativistic probe with finite mass. The differential cross section with respect to the energy transfer (T) can be expressed as where ⃗ k 1 and ⃗ k 2 are, respectively, the 3-momenta of the incoming and outgoing χ q , q 2 is the squared 4-momentum transfer, R L and R T are the longitudinal and transverse response functions, while V L and V T are the longitudinal and transverse kinematic factors for electromagnetic interaction. The response functions can, in principle, be calculated by many-body wave functions. In practice, for most cases, the calculations are highly nontrivial and schemes such as free electron approximation (FEA) and equivalent photon approximation (EPA) provide good estimations at certain kinematic regions.
The cross section evaluation with FEA assumes electrons occupying atomic orbitals as individual and independent. The electrons are inactive when energy transfer is below the ionization threshold for their corresponding orbitals. Sharp step functions are therefore produced in the cross section. Derivations with FEA take the specific kinematics case q 2 ¼ −2m e T (where T is the energy transfer) fixed by energy-momentum conservation. For high-energy scattering where the atomic binding energy is negligible, FEA is a good approximation. The differential cross section for χ q elastic scattering off electron is given by However, FEA largely underestimates the scattering cross section at T comparable to the atomic transition energy scale [49].
In kinematical region, where the momentum transfer by highly relativistic charged particles is small, EPA is a good approximation in the cross section calculations. In the limit of q 2 → 0, the contribution from the longitudinal polarized virtual photons vanishes and that from transversely polarized virtual photons approaches the cross section due to real photons. Furthermore, the on-shell transverse response function is directly related to the total photo-absorption cross section σ γ abs ðTÞ ¼ As a result, the differential cross section of Eq. (7) can be approximated by where the transverse part for Coulomb interaction is described by and the scattering angle θ depends on ⃗ k 1 , ⃗ k 2 and ⃗ q. Upon integration after a change of variable from dðcos θÞ to d⃗ q 2 , the EPA cross section for ultrarelativistic particle can be written as The derived differential cross sections of m χ q ¼ 1 keV at E χ q ¼ 1 MeV on Ge target under FEA and EPA schemes are depicted in Fig. 2(b). Although the EPA is a good approximation in the most interesting sub-keV region of T, it underestimates the scattering cross section above a few keV region of T where FEA works well. The EPA and FEA schemes therefore serve as conservative approximations in the region near and away from ionization thresholds, respectively. Accordingly, the scattering cross sections are combined in this analysis such as The uncertainties of the differential cross section are estimated via comparisons with an ab initio many-body calculation using multiconfiguration relativistic randomphase approximation (MCRRPA) which has been successfully applied to photoexcitation, photoionization, and neutrino-induced ionization of germanium and xenon atomic systems [30,50,51]. As illustrated in Fig. 2(b), the combined cross sections of Eq. (13) from EPA and FEA is 30% less than those due to MCRRPA calculations at energy away from the discrete excitation thresholds. This scheme will give conservative estimate in deriving limits. The "smoking gun" experimental signatures of χ q are observation of K-and L-shell peaks at the specific binding energies in the energy spectra, with known intensity ratios. This feature requires excellent energy resolution of the PCGe detectors where keV X-ray peaks can be resolved.

IV. EXPERIMENTAL CONSTRAINTS
Analysis results presented in this article are based on data taken with an n-type PCGe at KSNL. The facility [22][23][24][25] is located at a distance of 28 m from a 2.9 GW thermal power nuclear reactor core with 30 meter-water-equivalent overburden. The n-type PCGe is enclosed by a NaI(Tl) anti-Compton (AC) detector of mass 38.3 kg. This setup is placed inside the 50 ton shielding structure equipped with cosmic-ray (CR) veto scintillator panels. Data with n-type PCGe are selected for its low threshold at 300 eV, and the absence of anomalous surface events which serves to minimize complications to the analysis [17]. A total of 124.2=70.3 kg-days of reactor ON/OFF data are adopted for this analysis.
Every Ge-trigger is categorized by AC þð−Þ ⊗ CR þð−Þ , where the superscript -(+) denotes anti-coincidence (coincidence) of the Ge-signals with the AC and CR detectors [18]. The AC − ⊗ CR − events are uncorrelated with other active veto detector systems and are therefore candidates of neutrino, WIMP and other exotic events. The Reactor ON-OFF residual spectrum and the raw rates for the AC − ⊗ CR − events are depicted in Figs. 3(a) and 3(b), respectively.
The expected differential count rates due to possible χ q interaction with matter are obtained by integrating the χ q -flux formulas of Sec. II with the differential cross sections of Sec. III: where ρ A is atomic number density per unit target mass and ðE min ; E max Þ are the (minimum, maximum) energy of χ q . Constraints from each of the three discussed χ q channels are derived from the measured AC − ⊗ CR − spectra after subtraction of (i) internal radioactivity due to K=L-shell X-ray peaks from cosmogenically activated isotopes in the Ge-target, and (ii) a flat background estimated from ambient high-energy γ rays, following background understanding and analysis procedures from earlier work on similar detectors and configurations [18][19][20][21].
The attenuation of relativistic-χ q between source and detector is estimated by via the radiation length (X e 0 ) of high energy electrons pre-dominantly through bremsstrahlung in matter. The analogous radiation length (X χ 0 ) of the relativistic-χ q is given by where m e is mass of the electron. The KSNL reactor source is approximated by 10 m of water (X e 0 ¼ 36.08 cm), 10 m of concrete (X e 0 ¼ 11.55 cm) and 50 cm of lead (X e 0 ¼ 0.55 cm), while the atmospheric source is given by 20 km of air (X e 0 ¼ 29.89), 10 cm of concrete and 50 cm of lead. The attenuation effects place upper bounds to the excluded regions in the analysis.

A. Reactor-χ q
The reactor ON/OFF comparison provides a sensitive probe for laboratory searches on exotics particles. A positive signal of reactor-χ q discussed in Sec. II A would manifest itself as excess of events with 1=T 2 profile in the ON-OFF "residual" energy spectra. The kinematic threshold of χ q pair production in the field of electron is 2m χ q ðm χ q þ m e Þ=m e , therefore the search of reactor-χ q is only sensitive to m χ q < 1 MeV. The millicharged χ q produced by 1 to 5 MeV γ rays are found to be broadly distributed over the energy range 0 < E χ q < 5 MeV, with a few percent below 50 keV.
More than 50% of the γ-fluxes of nuclear reactor cores are due to prompt-γ of highly excited fission fragments, with slight dependence on the actual fuel composition [26]. This introduces an uncertainty to prompt γ-flux Eq. (1), giving a 10% systematic error to reactor-χ q fluxes.
The expected rate of χ q producing in reactor at Gedetector, in the units of kg  3. (a) Reactor ON-OFF residual spectrum of n-type PCGe with AC − ⊗ CR − selection. The best-fit 2σ region of possible χ q interactions are shown as the yellow band, with an excluded scenario of reactor-associated χ q at specified (m χ q ¼ 1 keV, δ ¼ 3 × 10 −6 ) is superimposed. (b) The total AC − ⊗ CR − spectrum showing a flat background due to ambient high-energy γ-rays and the L-shell X-rays from internal radioactivity. Excluded scenarios of atmospheric and DM-χ q at specified (δ ¼ 2.5 × 10 4 , I δ ¼ 1 × 10 −2 cm −2 s −1 sr −1 ) and (m χ q ¼ 1 keV, δ ¼ 1.06 × 10 −7 ), respectively, are superimposed. on δ as a function of m χ q , the minimum-χ 2 fit method is adopted. Systematic uncertainty in flux (AE10%) is taken into account. The input are varied within one-RMS range, and the most conservative constraints are adopted. The excluded (m χ q , δ) parameter space at 90% C.L. is thus derived and displayed in Fig. 4(a) as a solid deep-red region. Previous laboratory limits [5,[14][15][16]44] as well as bounds from cosmological and astrophysical model [52], are superimposed. The residual spectrum with best-fit for m χ q ¼ 1 keV is depicted in Fig. 3(a) and the 2σ uncertainty band is superimposed. These new results provide improved limits at light-χ q of m χ q < 4 × 10 −4 GeV=c 2 .

B. Atmospheric-χ q
The high energy cosmic rays are capable of creating massive χ q via interaction with the Earth's atmosphere, as discussed in Sec. II B. In this analysis, we focus on χ q with δ < 10 −3 and the signal would be neutrino-like, where χ q interact only once with the Ge-target, without producing signals at CR and AC veto detectors.
The upper limits on I χ q are derived using the most conservative assumption that predicted signal cannot exceed the measured background AC − ⊗ CR − rate depicted in Fig. 3(b). The differential cross section in Eq. (12) is independent of E χ q for ultrarelativistic condition. The theoretically expected spectrum of χ q -induced events is where ϕ χ q is the only free parameter used to fit the measured spectrum. The corresponding upper limits on I χ q as a function of δ are displayed in Fig. 4(b) and constraints from previous experiments [54][55][56][57][58] are also superimposed. The constraints are applicable to m χ q > 4 × 10 −4 GeV=c 2 , below which attenuation of χ q -fluxes by radiation effects before reaching the detector would limit the sensitivities. The differential cross section due to χ q -interaction has enhancement in sub-keV region of T as shown in Fig. 2(b), thereby extending the lower reach of δ from 10 −3 to new domain of 10 −6 . An excluded scenario of (δ ¼ 2.5 × 10 −4 , Fig. 3(b) as illustration.
C. Dark cosmic ray DM-χ q As discussed in Sec. II C, dark matter ion can be accelerated by supernova shock waves. In the absence of injection spectrum, the prediction of dark cosmic-ray flux is strongly model-dependent, in which a conservative systematic uncertainty of AE50% is assigned. Subsequent analysis follows the same procedures as described in the reactor-χ q case. The expected yield of the DM-χ q events in the detector is calculated by mean of Eq. (14) adopting flux parameter space for millicharged particles with a massless dark photon. Cosmological and astrophysical bounds from Ref. [52] are denoted as dotted lines. The direct laboratory limits from Accelerators [5], OPOS [14], SLAC [15], LHC [16], LSND [53], MiniBooNE [53] and SK-ν [44] experiments are represented as shaded regions. New 90% C.L. excluded regions of this work with χ q from reactor and dark primary cosmic rays are shown as red and blue shaded areas, respectively. The dotted lines correspond to the upper bounds of the exclusion regions, due to complete attenuation of χ q before reaching the detector. (b) Excluded parameter space at 90% C.L. on incoming flux of χ q from secondary dark cosmic rays versus its charge fraction δ. Results from other benchmark experiments such as LSD [54], Kamiokande-II [55], MACRO [56], CDMS II [57], and MAJORANA [58] are superimposed. The constraints are applicable to m χ q > 4 × 10 −4 GeV=c 2 , below which attenuation effects before reaching the detectors would limit the sensitivities. from Eq. (5). The minimum energy (E min ) of dark cosmic rays is adopted with the assumption of their being ultrarelativistic, or E min > γm χ q . The value of γ is assigned to be 10, where EPA and exact calculation for cross section are consistent [59]. The maximum energy of dark cosmic rays is taken as E max ¼ δ × 10 6 GeV.
The evaluation of upper limits for various mass follows the conservative assumption that the predicted signal cannot exceed the measured rates. For fixed m χ q , δ is the only free parameter derived by minimum-χ 2 fit method with a flat background assumption. The model-dependent excluded regions on δ at 90% C.L. for m χ q between 10 −6 GeV=c 2 to 1 GeV=c 2 are depicted in Fig. 4(a) as shaded blue region. We note that the bounds are strongly dependent on the DM-χ q model adopted and therefore not at comparable solid footing as those for reactor-χ q . The blue dotted-lines denote boundaries of exclusion due to the attenuation of relativistic-χ q between source and detector. The lower reach of m χ q is defined by E min > γm χ q , the minimum energy required by DM-χ q to produce recoil electron above the minimum detectable energy. Example of an excluded spectrum is shown in Fig. 3(b) for m χ q ¼ 1 keV and δ ¼ 1.06 × 10 −7 .

V. CONCLUSION AND PROSPECTS
The hidden sector with massless gauge boson allows possibilities of multicomponent dark matter. Its ionic constituents can acquire small charges (millicharge) under the hidden sector gauge group. In this paper, we have derived the direct constraints on χ q with low-threshold point-contact germanium detectors under the scenarios of χ q produced (i) at nuclear power reactors, (ii) as products of cosmic-ray interactions, and (iii) as dark matter particles accelerated by supernova shock. The sub-keV sensitivity leads to improving direct limits of δ at small m χ q and extending the lower reach of δ to 10 −6 .