Axionlike particles searches in reactor experiments

Reactor neutrino experiments provide a rich environment for the study of axionlike particles (ALPs). Using the intense photon flux produced in the nuclear reactor core, these experiments have the potential to probe ALPs with masses below 10 MeV. We explore the feasibility of these searches by considering ALPs produced through Primakoff and Compton-like processes as well as nuclear transitions. These particles can subsequently interact with the material of a nearby detector via inverse Primakoff and inverse Compton-like scatterings, via axio-electric absorption, or they can decay into photon or electron-positron pairs. We demonstrate that reactor-based neutrino experiments have a high potential to test ALP-photon couplings and masses, currently probed only by cosmological and astrophysical observations, thus providing complementary laboratory-based searches. We furthermore show how reactor facilities will be able to test previously unexplored regions in the $\sim$MeV ALP mass range and ALP-electron couplings of the order of $g_{aee} \sim 10^{-8}$ as well as ALP-nucleon couplings of the order of $g_{ann}^{(1)} \sim 10^{-9}$, testing regions beyond TEXONO and Borexino limits.


Introduction
In recent years there has been an increasing interest in weakly/feebly coupled light physics. Arguably the main motivation has been given by searches at the LHC, which have provided null results of new physics at the TeV scale. Weakly coupled hidden sectors are a rather plausible explanation for the absence of such signals. They are motivated as well by the fact that the absence of dark matter direct detection signals can be justified if the new secluded sector involves sub-MeV degrees of freedom. A common feature of these scenarios is the presence of light scalars/pseudoscalars, usually called axionlike particles (ALPs), with masses and couplings to Standard Model (SM) particles which can span over many orders of magnitude. Examples of the latter are QCD axions [1][2][3][4], a widely explored extension of the SM addressing the strong CP problem. They can also arise as a consequence of spontaneously broken approximate symmetries, of which examples include familons and Majorons [5][6][7]. In these scenarios, ALPs are the pseudo-Nambu-Goldstone bosons of theories where a symmetry group G ⊃ SU (3) × SU (2) × U (1) is spontaneously broken at a high scale, leaving behind light (pseudo)scalar states.
ALPs may generate signals in a large variety of experimental environments. They can as well affect the evolution of stars and of core-collapse supernova explosions and can potentially affect the early universe dynamics. They can, for instance, contribute to the number of relativistic degrees of freedom or they can release entropy in the heat bath thus disrupting BBN predictions (see e.g. [8,9] and references therein). In the high-mass regime, and at the laboratory level, ALPs have been searched for in beam-dump experiments as well as in radiative Υ decays at CLEO and BaBar [10][11][12][13][14][15]. They have been as well investigated at accelerator experiments through different event topologies, first at LEP and more recently at the LHC [16][17][18][19]. Prospects of exploring regions not yet tested in this mass regime include Belle II, NA62, NA64, SeaQuest, MATHUSLA, CODEX-b, FASER and SHiP [20][21][22][23][24][25][26][27][28].
ALPs have been searched for as well in reactor experiments, motivated by theoretical arguments pointed out first by Weinberg and subsequently by Donnelly, Freedman, Lytel, Peccei and Schwartz [3,29] 1 . To the best of our knowledge, ALP reactor searches looking for various nuclear magnetic transitions were pioneered at the Bugey reactor [30], and followed by searches carried out by Zehnder and Lehmann et. al. in which instead ALPs were searched for in nuclear magnetic transitions of 127 Ba * and 65 Cu * [31,32]. Most recently, the TEXONO collaboration extended the analysis of Ref. [30] through a comprehensive investigation at the Kuo-Sheng reactor nuclear power plant, by also employing nuclear magnetic transitions along with inverse Compton-like and Primakoff scattering for detection [33].
At present, there is a rapidly developing physics program aiming at measuring coherent elastic neutrino-nucleus scattering (CEνNS) using reactor neutrinos, which offers exciting possibilities to explore physics beyond the SM. To this purpose, several new experiments with sub-keV threshold capabilities are in preparation in The Americas, Western and Eastern Europe and Asia, including-but not limited to-CONNIE, MINER, CONUS, ν-cleus, Ricochet, RED-100 and TEXONO [34][35][36][37][38][39][40] (for a detailed review see Ref. [41]). Apart from CEνNS measurements, this program is aiming at improving upon previous measurements of the weak mixing-angle, the root-mean-square radii of neutron distributions for different isotopes, as well as to look for new physics signatures in the form of e.g. neutrino magnetic moments, light vector and scalar mediators and neutrino generalized interactions (see e.g. [42][43][44][45][46][47][48][49][50][51][52][53][54]). Recently, it has been pointed out that these technologies also open up a new direction for ALP searches [55], further motivating our present work.
In this paper, we investigate the feasibility and extent of such searches by including leading ALP production and detection mechanisms at reactor neutrino experiments. In addition to nuclear magnetic transitions, considered by previous reactor searches, here we take into account additional ALP production mechanisms such as Compton-like and Primakoff scattering. We furthermore consider detection through their inverse processes as well as via axio-electric absorption and ALP decays to diphoton and electron pairs. Our analysis follows minimal ALP assumptions in the sense that it considers only those processes that allow production and detection through the same ALP-SM coupling, the exception being production through nuclear de-excitation. In this 1 These arguments were applicable to the QCD axion and so these searches. The experimental results, however, hold as well for ALPs. case we consider detection through either ALP-photon or ALP-electron couplings, following to a large extent the analysis carried out by TEXONO in Ref. [33]. Rather than sticking to a particular CEνNS reactor-based experiment, our analysis is formulated as generic as possible, and therefore being representative of what current or near-future technologies could achieve.
Due to kinematic constraints, reactor-based searches will cover regions up to ALP masses of order of 10 MeV. Thus they will not cover some of the regions that Belle II, NA62, NA64, SeaQuest, MATHUSLA, CODEX-b, FASER and SHiP will have access to. However they can potentially provide valuable information on regions yet unexplored by laboratory experiments, while in some cases they can actually explore regions so far not tested at all. Since CEνNS reactor experiments will run anyway, an ALP research program will serve towards widening the physics reach of these facilities in a direction complementary to other ALP search strategies.
The remainder of this paper is organized as follows. In Sec. 2 we present the interactions our analysis deals with and we briefly discuss some representative QCD axion models. In Sec. 3 we discuss the ALPs production and detection mechanisms considered in the present study, relevant at reactor experiments. In Sec. 4 we present our main results, the sensitivities achievable in ongoing or near-future experiments, along with sensitivities expected at the next-generation CEνNS reactor experiments. Finally, in Sec. 5 we summarize and present our conclusions.

Axionlike particle interactions
Axionlike particles (ALP) are pseudoscalars that feebly couple to the SM particles. They arise in a variety of new physics scenarios, including those related with dark matter in which the pseudoscalar mediates interactions between the dark and visible sectors [56,57]. In general, they can be expected in scenarios where the spontaneous breaking of a global symmetry takes place. In that regard, the Majoron could be thought as a prototypical example linked to spontaneous lepton number symmetry breaking, and ultimately with the origin of Majorana neutrino masses [58]. From a more theoretical point of view, they can be as well understood as low energy manifestations of string theories [59,60].
Strictly speaking, the SM degrees of freedom to which the ALPs couple depend on the specific model realization. However, from a purely phenomenological point of view the ALP couplings can be parameterized in terms of dimension-five effective operators. Here we consider couplings to photons, electrons and nucleons, for which the interactions can be written according to where n refers to the proton-neutron isospin doublet and g (0) ann (g (1) ann ) to the isosinglet (isotriplet) axion-nucleon couplings, while the first term in the Lagrangian involves the ALP kinetic and mass terms as well as the scalar potential.
It is worth noting that the dimension-five operator nature of the couplings to photons is evident, while for fermions this is not the case. Indeed, couplings to fermions are dictated by derivative interactions (see e.g. [8,9]) 2) with f a the ALP decay constant and C af dimensionless strength couplings. Operators in the Lagrangian (2.2) are of dimension five, but they can be cast into those described in the Lagrangian of Eq. (2.1) with the aid of the fermion equation of motionf ∂ µ γ 5 γ µ f = −im ff γ 5 f , combined with integration by parts. In terms of the decay constant and fermion masses, the couplings in the Lagrangian (2.2) are thus given by So far the discussion applies to ALPs, which are not related with the solution to the strong CP problem and for which their mass and decay constant are not related. Reactor experiments are best suited to test parameter space regions of ALPs rather than of the QCD axion (see Sec. 4). Nevertheless, in what follows we briefly provide a few relations derived in QCD axion models that will enable comparing the extent at which they could be tested (for comprehensive reviews see e.g. [8,9]). The interactions in (2.1) along with (2.3) hold as well for the QCD axion, but in this case the relation m a 5.7 10 12 GeV f a µeV (2.4) applies. This expression combined with Eq. (2.3) leads to relations which can be mapped into the axion mass-coupling plane in the form of stripes whose widths depend on the UV completion accounting for the axion effective Lagrangian, namely The dimensionless parameters in Eqs. (2.5) and (2.6) read (see e.g. Ref. [9])

9)
C a,sea = 0.038c 0 s + 0.012c 0 where E and N correspond to the QCD and electromagnetic anomaly coefficients and the second term in Eq. (2.11) arises from one-loop corrections. Barring the case E/N=0, the axion couples to (gluons) photons and electrons (g aee = 0 even for c 0 e = 0). The coefficients c 0 f are model dependent and in some instances are zero at the tree level.
The most renowned invisible axion UV completions are the KSVZ and DFSZ models [61][62][63][64], or variants of them. To compare with our sensitivities we select few of them as follows: for KSVZ-type models we consider those with c 0 u i = c 0 d i = c 0 e = 0 and E/N ⊂ [44/3, 5/3], which cover a large variety of hadronic axion models [65]. For the case of DSFZ-type models we instead consider the DFSZ-I and DFSZ-II realizations, for which the relevant parameters are given by DSFZ-I:

ALPs production and detection mechanisms at reactors
Photons are abundantly produced in nuclear power plants. Processes responsible for γ emission are fission, decay of fission products, capture processes in fuel and other materials, inelastic scattering in the fuel, and decay of capture products [66]. The prompt photon flux can be described by an approximation originally derived from an analysis of γ radiation in the FRJ-1 research reactor [67] and which has been recently used in dark photon phenomenological analyses [68,69]. Once produced, these photons can interact with the fuel material which we assume to be 235 U. This choice is motivated by the fact that most of the reactor experiments our analysis applies to, use this uranium isotope (totally or partially) as fuel material. Possible interactions fall in three category processes: (i) a first category with final-state photons, which involves Rayleigh (elastic) or Compton (inelastic) scattering; (ii) a second category that accounts for possible photon absorption by the material, and that involves photo-electric absorption and e-pair production by the nucleus or electron field, (iii) photon-ALP conversion. Specifically, categories (i) and (ii) are SM processes, while category (iii) requires a SM-ALP coupling. The relevant Feynman diagrams for processes (i) and (ii) are shown in Fig. 1. It is worth noting that processes of type (iii) should involve an initial-and final-state photon and ALP, respectively. Apart from Primakoff photon-ALP conversion (diagram (a) in Fig. 2), the other two proceed from the SM Compton and Rayleigh processes by changing the final-state photon to an ALP (see diagrams (b) and (c) in Fig. 2). The leading mechanisms are Compton-like production and nuclear de-excitation. All these processes involve independent SM-ALP couplings: Primakoff production is enabled by g aγγ , Compton-like production by g aee and nuclear de-excitation by g ann (n = neutron, proton). Primakoff, Compton-like and nuclear de-excitation can also serve as detection mechanisms. The other possible detection processes follow from the SM photoelectric absorption and e-pair production by changing γ → a, as shown in Fig. 3. In the former case, axio-electric absorption [76], an ALP is absorbed by an atom resulting in electron emission. In the latter, an ALP is absorbed and an electronpair is emitted. Depending on the ALP lifetime, detection can proceed as well through ALP decay processes: a → γγ, a → e + e − and a → nn (since for the last two processes m a > 2m e and m a > 2m n should be satisfied, decay to nucleons is not a viable possibility in nuclear reactor experiments). Table 2 shows the different production and detection mechanisms together with the relevant couplings. From  ALPs production mechanisms at a nuclear reactor plant. Process (a) corresponds to Primakoff photon-ALP conversion, process (b) to Compton-like scattering and process (c) to nuclear de-excitation. Depending on the type of SM-ALP coupling (g aXX , X = γ, e, n, where n refers to nucleons), in addition to the processes shown in Fig. 1, photons interacting with the reactor fuel can be subject to these processes.

Experiment
this table one can see that self-contained one-parameter analyses (axion mass aside) of ALP production should be done as follows: • Analyses with g aγγ should involve production through Primakoff and detection through inverse Primakoff and decay into a γ−pair final states.
• Analyses with g aee should include production through Compton-like and detection through inverse Compton-like, axio-electric, e−pair production in e field and decay into an electron-positron pair final states.  Figure 3. ALPs detection mechanisms, along with inverse Primakoff, inverse Compton-like and nuclear excitation (inverse processes in Fig. 2). Process (a) corresponds to the axioelectric process (photo-electric absorption equivalent, diagram (c) in Fig. 1), while process (b) to axion e−pair production (e−pair production in the field of the nucleus or electron equivalent, diagram (d) in Fig. 1). As in the production case, whether such detection mechanisms are operative depends on the type of coupling the analysis is subject to. See Tab. 2 for details.
• Analyses with g ann should be done considering production through nuclear deexcitation and detection through nuclear excitation.
We relax our analysis from considering ALP e−pair production in the field of N or e − (diagram (b) in Fig. 3). Although we do not calculate the ALP e−pair production cross sections, we expect that the detection through these processes should produce a rather small event yield, similar to their SM counterparts which are suppressed (see Fig. 4). We do not consider either ALP nuclear excitation in detection. Hence, in our analysis of ALP production through nuclear de-excitation (Sec. 3.3) we consider detection mechanisms which involve other ALP couplings (either g aee or g aγγ ).
Phenomenological ALP studies follow a three-step analysis: (i) production, (ii) propagation towards the detection chamber and (iii) detection. The relevant production mechanisms assumed in the present work, namely, Primakoff conversion, Compton-like scattering or nuclear de-excitation, generate an ALP flux which propagates and can potentially reach the detector depending on the ALP lifetime, determined by its survival probability. Neglecting possible interactions with the detector shielding material, ALPs reaching the detector will then interact with the target material yielding-in principle-a measurable signal. The quantities that intervene in this three-step process then combine through a convolution, that defines the differential event yield from which the number of counts at the detector is calculated.
Possible signals depend on the kind of process ALPs are subject to. For instance, in the axio-electric process ALP absorption by the material produces phonons (from the temperature rise of the absorber) and electron-hole pairs leading to heat and ionization signals. Ideally, detection through the axio-electric process should then be done using dielectric crystals, whose properties are optimized for such type of signatures. DAMA, CDMS and more recently the SuperCDMS Soudan experiment

Process
Coupling Prod Det

Decay processes
Process Coupling Prod Det n-pair final state a → n + n g ann Table 2. ALPs production and detection mechanisms. Check marks (crosses) are used to identify (exclude) production and/or detection processes. Color refers to processes that should be-in principle-included in a self-contained one-parameter analysis. Note that decay to nucleon pairs is forbidden by kinematic arguments.
have followed this strategy, i.e. DAMA using NaI crystals, while CDMS and Su-perCDMS Soudan using a cryogenic detector based on high-purity germanium and silicon crystals [77][78][79]. The different detector technologies as well as detector specifications of typical experiments to which our analysis can apply are shown in Tab. 3. These parameters motivate the generic numbers used in our analysis: germanium as material target, m det = 10 kg and L = 10 m (see Sec. 4).

Production and detection through photon-ALP coupling
Let us now focus on Primakoff-like production and detection processes. In full generality, the determination of the ALP flux involves a convolution (over ingoing photon energies) of the photon flux and the photon-ALP conversion differential cross section, properly normalized to account for the processes (i) and (ii). With the ALP flux at hand, and after assuring the ALP survival probability is sizeable enough, the detection proceeds through a second convolution (over ingoing ALP energies) of the ALP flux and the ALP-photon conversion differential cross section. If the ALP  Table 3. Detectors to which our analysis may apply, along with relevant parameters for ALP detection. With germanium as a target material, m det = 10 kg and L = 10 m are rather representative of the experiments our analysis can cover. CCD refers to chargecoupled devices.
undergoes decay within the detector chamber, photon signals can be triggered as well by its decay. In such a scenario the ALP detection is therefore a combination of inverse-Primakoff scattering and decay contributions, with the latter being mainly relevant in the "high" ALP mass region (the decay partial width is proportional to the ALP mass). Let us discuss this in more detail. Accounting for the probability that an ALP of energy E a is emitted after a γ − 235 U interaction, the ALP flux proceeds then from the following convolution where σ Tot = σ SM + σ P P , with the SM cross section being σ SM = σ R + σ C + σ PE + σ e−pair , while σ P P ≡ σ γ→a is the total Primakoff scattering cross section for ALP production. Results for the SM components are obtained from the XCOM: Photon cross section database [80] and shown in Fig. 4 for completeness. The integration  limits are determined by the kinematic relation 2 found for cos θ = +1 (cos θ = −1) and M N being the nuclear mass. Differences between the upper and lower limits are in the sub-keV ballpark, almost independent of the ALP mass. Kinematically, this implies that a photon with a given energy E γ produces ALPs within a rather narrow energy range.
The ALP survival probability, assuring that the ALP flux reaches the detector, can be cast as [55] P surv = e −LEa/| pa|τ , (3.4) where L refers to the distance between the reactor and detector and τ to the ALP lifetime in the fixed target frame, determined in turn by the ALP total decay width: Thus, for a photon with energy E γ the consequent ALP distribution can be approximated with a monochromatic one 3 This observation allows to determine the ALP flux with a simplified version of For the calculation of the total Primakoff scattering cross section we rely on the expression [81] a result that does not involve any approximation (see discussion below). Here s, t are Mandelstam variables and αZ 2 F 2 (t) is the effective electric charge to which the photon couples to. Its value is controlled by F 2 (t), a form factor that we normalize to unity. For a sufficiently small exchanged momentum |q| = √ −t ≤ 7.39m 2 e 7.0×10 −3 fm −1 , the DeBroglie wavelength of the scattering process is large enough so to partially involve the atomic electron cloud [82]. In that case the effective electric charge is diminished (screened), with the reduction depending on how small the exchanged momentum becomes. For such regime we use the following atomic form factor [83] with a 0 calculated in the Thomas-Fermi atom using the Moliere representation (3.10) Since we choose t = q 2 < 0, contrary to what Ref. [83] does, Eq. (3.9) results from shifting t → −t in the form factor provided by that reference. For a larger exchanged momentum |q| = √ −t > 7.39m 2 e 7.0 × 10 −3 fm −1 the DeBroglie wavelength becomes smaller and hence the atomic electron cloud is not relevant anymore. However as q increases the nuclear electric charge diminishes as well, 3 For definiteness, in what follows the cross sections are written in the form σ Y X , where X =(P, D) indicates the corresponding mechanism e.g. production or detection and Y =(P, C, A) indicates the nature of the process e.g. Primakoff-like, Compton-like and axio-electric, respectively. Similarly, the ALP flux is given in the form therefore the atomic form factor has to be replaced by a nuclear form factor which we choose to be of the Helm form [84] 4 . It is worth emphasizing that while the inclusion of a nuclear form factor has a negligible effect on our results, the inclusion of the atomic form factor is instead important, in particular in the low ALP mass region m a 10 4 eV. This can be understood from the nature of the Primakoff process, for which the largest amount of events are found at small t or along the forward direction.
The calculation of the Primakoff cross section in Ref. [82] is valid in the forward scattering limit as well as in the limit E γ = E a and m a → 0. In contrast, the result in Eq. (3.8) does not involve any of those approximations, even if E γ = E a is still valid to a large degree, as we have already pointed out. However, departures from forward scattering as well as a massive ALP can produce sizable deviations, particularly relevant at photon energy threshold 5 .
For the integration of the Primakoff differential cross section in Eq. (3.6) we use the following limits [85] where p 1 (k 1 ) refer to the initial-state photon (final-state axion) three-momentum given by with s = M 2 N + 2M N E γ . Then, the differential event rate consists of two terms, one from the inverse-Primakoff scattering and another one from ALP decays a → γγ Here, the first term results from a convolution of the ALP flux and the ALP-photon differential cross section 6 , namely which taking into account Eq. (3.6) allows writing the ALP yield through Primakoff scattering according to 15) where N T = N A /m molar with N A = 6.022 × 10 23 /mol and m molar the target nuclei molar mass in kg/mol (N T thus measures the amount of nuclei/kg). Moreover, ∆t stands for the data-taking time and m det refers to the detector mass in kg. Integration limits are E min a = 25 keV and E max a = 10 MeV (determined by the photon natural threshold at a reactor). The second term in Eq. (3.13) is instead given by where P decay accounts for the probability that the decay occurs within the detector, and reads where L det stands for the detector length and A = L 2 det denotes the detector transverse area.

Production and detection through electron-ALP coupling
The MeV photons produced in a reactor may also scatter off electrons in the materials of the reactor core and they can produce ALPs via the Compton-like process γe − → ae − . This process requires the individual SM-ALP coupling g aee to be non-zero. The ALPs so produced can be subsequently detected via three possible processes relying on the same g aee coupling: inverse Compton-like, axio-electric and decay into e + e − pairs. The total number of observed ALP events is computed similarly to the Primakoff case.
As a first step, we calculate the Compton-produced ALP flux via a convolution of the reactor photon flux dΦ γ /dE γ and the production cross section dσ C P /dE a , as where we assume E γ ,min = (2m e m a +m 2 a )/(2m e ), E γ ,max = 10 MeV and σ Tot = σ SM + σ C P . The differential Compton-like cross section for ALP production has threshold s > (m a + m e ) 2 and is given by [86] dσ

19) where s is the Mandelstam variable and
(3.20) The total Compton-like production cross section is obtained by integrating dσ P /dE a over the integration limits [87]: One then gets the differential number of ALP events with a second convolution, of the ALP flux dΦ a /dE a with the detection cross section. The latter will be different depending on the detection mechanism: namely inverse Compton-like, axio-electric, and e + e − decay, respectively. In this case, the total ALP decay width reads: while, the relevant detection cross sections can be cast as [88][89][90][91] dσ where σ PE denotes the SM photoelectric cross section of the detector material. Here, θ is the scattering angle, y = 2m e E a + m 2 a and β = | p a |/E a 7 . The integration limits E a,min (max) are functions of E γ : Finally, the total number of observed ALP events from Compton-like production in a reactor neutrino experiment would be for the detection mechanisms inverse Compton-like, axio-electric, and e + e − decay, respectively. The integration limits E a,1 (2) are extracted from the ALP flux range, whereas the final integration over E γ is performed using the limits E γ,min = 25 keV-which we assume to be the minimum detectable photon energy-and E γ,max that can be obtained by inverting Eq. (3.26), leading to the simple relation E γ,max = y max /(2(E a,max + m e − | p a |)).

Production through nucleon-ALP coupling
A heavy nucleus disintegrates through α or β decay or fission processes, while the resulting daughter nucleus is typically left in an excited state. It can then undergo fission, provided the state is above the excitation energy for this process to take place, otherwise it will in most cases decay by γ emission. The study of nuclear de-excitation follows from longitudinal (L), transverse (T ), transverse-longitudinal (T L) and transverse-transverse (T T ) nuclear responses. Their multipole expansion (projection in the angular momentum basis) define the Coulomb multipole operator M JM (q), the transverse electric and magnetic multipole operators T el JM (q) and T mag JM (q) and the longitudinal multipole operator L JM (q). The angular momentum quantum numbers J and M are determined by the angular momentum conservation conditions |J f − J i | ≤ J ≤ J f + J i and M = M f − M i , where J i and J f stand for the momenta of the initial and final states involved in the nuclear excitation (or de-excitation) process (for details see e.g. [92,93]).
For virtual photons -as e.g. those involved in electron-nucleus scattering processes -all responses matter, while for real photons only transverse electric and magnetic responses are relevant. These account for electric (magnetic) transitions EJ (MJ), with multi-polarities determined by J. The lowest orders correspond to radiation emission from charge (J = 0, only electric), dipole (J = 1) and quadrupole (J = 2). In addition to the angular momentum conservation, transitions are as well subject to parity conservation, for instance an EJ transition produces a photon with parity Π γ = (−1) J , while a MJ transition a photon with parity Π γ = (−1) J+1 . Only photons for which Π γ = Π i Π f is satisfied are allowed. If parity is not changed, ∆Π = 0, only even (odd) electric (magnetic) transitions are possible. If on the contrary ∆Π = 0, only odd (even) electric (magnetic) transitions are available.   If ALPs couple to nucleons, nuclear de-excitation can proceed as well through ALP emission. Since ALPs are pseudoscalars (|J Π = |0 − ) they can be emitted only in magnetic transitions [3], in contrast to photons which are produced in electric transitions as well. The possible quantum numbers they carry away are therefore |0 − , |1 + , |2 − , |3 + , ... [29]. To the best of our knowledge, searches for axions in magnetic transitions were pioneered by Zehnder in a nuclear decay experiment which used the 664 keV M4 transition from the |11/2 − excited state to the |3/2 + ground state of 127 Ba * [31]. Subsequent searches were conducted at the Bugey reactor (∼ 2.8 GW) using the neutron capture M1 isovector transition along with M4 transitions in 97 Nb, 91 Y, 137 Ba and 135 Xe (see left graph in Fig. 5) [30]. Searches focusing on the 1115 keV M1 transition from the |5/2 − excited state to the |3/2 − ground state of 65 Cu were performed as well [32], while the same approach was later followed by Avignone et. al. [76]. More recently, similar searches were pursued at the Kuo-Sheng Nuclear Power Station by the TEXONO collaboration, which in addition to the transitions used at Bugey, considered a M1 transition in 7 Li (see left graph in Fig. 5) [33].
In our analysis we employ the fluxes reported by TEXONO: below we discuss in more detail the magnetic transitions from which they are generated. Thermal neutron capture on proton in the cooling water, p + n → d + γ, produces a monochromatic line at 2230 keV (deuteron binding energy). The deuteron ground state has magnetic dipole and electric quadrupole moments, the emitted γ is therefore mainly M1. On the other hand, the 7 Li * de-excitation involves M1+E2 transitions from the |1/2 − first excited state to the |3/2 − ground state, and so again the transition is mainly M1. For 91 Y * and 97 Nb * de-excitations from the first excited state to the ground state, as shown in Fig. 5, correspond to M4 transitions. Finally, 135 Xe * and 137 Ba * de-excitations involve M4 transitions from the second excited state to the ground state. Note that for 135 Xe, de-excitation from the |1/2 + first excited state to the ground state corresponds to M1+E2 transitions. If it was available, the flux from this transition (mainly M1) would be much larger than that from the M4 transition. The corresponding fluxes are shown in the right graph of Fig. 5.
ALP production depends on the ALP-nucleon couplings and on the ALP mass, while the relevant ALP flux is determined through the competition of the axionto-photon emission rates, as first calculated in [29]. In contrast to M4 transitions, neutron capture isovector M1 transitions (pn→ dγ) depend only on kinematics and the vector ALP-nucleon coupling [94] where µ 1 refers to the isovector magnetic moment which follows from the neutron and proton magnetic moments, µ 1 = µ p − µ n = 4.71 in nuclear magneton units µ N = 1/2m p (m p is the proton mass) [76]. For MJ transitions this ratio instead involves nuclear physics input, and in the long-wavelength limit it reads [33] Γ a Γ γ MJ = 1 πα (3.32) Here, µ 0 refers to the isosinglet magnetic moment normalized to µ N and obtained as µ 0 = µ p + µ n = 0.88 [76]. The remaining δ, η and κ parameters are nuclear-structure dependent. In particular, the multipole mixing ratio δ is determined by the EJ -to-MJ (J = J + 1) transition probabilities ratio, while κ and η are defined by ratios of the orbital ( L) and spin ( σ) nuclear operators reduced matrix elements, namely [76] where J i,f refer to initial and final state angular momenta and the sum runs over the total number of nucleons. Values of δ, η and κ for the 7 Li * M1 transition have been calculated in Refs. [95] and also adopted by the TEXONO collaboration in their analysis in [33]. The transition is an admixture of M1+E2, for which the transition probability ratio is small δ 1 and therefore it can be neglected, while the estimated values for η and κ are 0.5 and 1.0. For the remaining transitions δ = 0, as expected since they are purely magnetic, while values for η and κ are not available. However, using results for 65 Cu, 57 Fe, 55 Mn and 23 Na-calculated in Refs. [76,96]-TEXONO estimated these values to be η = 1 and κ = −3 for 91 Y and 97 Nb and η = −1 and κ = 1 for 135 Xe and 137 Ba [33]. We adopt these values in our analysis.
With the values of the nuclear-structure-dependent parameters specified, we can now write the ALP flux for the i−th transition (3.34) with Φ i γ the photon flux per fission for the i−th transition specified in the right graph in Fig. 5 and R f the fission rate, which we assume that proceeds entirely from 235 U and estimate to be 1.248 × 10 20 fissions/s for a 4 GW nuclear reactor. The ALP yield at the detector depends on the detection mechanism (see Tab. 2). In the presence of only g (a) ann couplings, the detection would be possible only through ALP nuclear excitation. This will require the detector to be built from one of the daughter isotopes from which the ALP was produced, to enable ALP absorption. Since this does not seem to be a viable possibility, in this case we abandon our minimal assumption analysis and instead we rely on detection through Primakoff, Compton-like, axio-electric processes and decays. The nuclear de-excitation analysis is therefore sensitive to coupling pairs and the ALP mass, in contrast to the Primakoff and Compton-like analyses which determine sensitivities on a single coupling and the ALP mass.
The ALP detection rate for the i−th transition proceeds as always from a convolution of the ALP flux and the differential cross section of the corresponding detection process. For detection through inverse Primakoff, Eq. (3.14) applies and so we write Using the ALP flux in Eq. (3.34) the ALP yield then becomes The same applies as well for detection through inverse Compton-like and axio-electric processes, for which we then write 37) and to ALPs detected through their decays into diphoton or electron-positron pairs:  Table 4. Production and detection parameters used in our analysis. Upper values refer to current (or near-future) experiments, while those in the lower row to a next-generation experiment. These values-used for definitiveness-should be understood as representative of the actual experimental parameters shown in Tabs. 1 and 3. P refers to power, PM to production material and TM to target material. Our analysis includes two background hypotheses and we assume the detector to be 100% efficient.

Experimental sensitivities
In this section we discuss sensitivities that the experiments listed in Tab. 1 and 3 could in principle achieve. Rather than focusing on a particular case we adopt a generic analysis whose results we believe provide a general picture. To do so we fix the relevant experimental production and detection parameters as shown in Tab. 4.
Upper row values provide sensitivities for ongoing or near-future experiments, while lower row values represent what could be achieved in a next-generation experiment.
Given an ALP coupling, our analysis includes the processes shown in Tab. 2, the exception being the case of ALP production through nuclear de-excitations for the reasons we have discussed in Sec. 3.3. Calculation of sensitivities is done with the aid of the following χ 2 parameter function Here by N ALP we generically refer to the ALP event yield, obtained through any of the processes that we have discussed in Sec. 3, while N bkg refers to the background hypothesis adopted (see Tab. 4). We begin our analysis by comparing the size of the ALP flux generated by Primakoff and Compton-like scattering as well as by M1 transitions (neutron capture in proton and 7 Li * de-excitation). Figure 6 shows the result for parameter choices as displayed in the plot. One can see that a large amount of ALPs can be produced through any of the mechanisms we consider. Despite being parameter dependent, the size of the flux can become large in certain regions of the parameter space, as shown in Fig. 6. Indeed, reactor-based neutrino experiments have the potential of reaching large sensitivities, further motivating our present study. Note that M4 transitions produce a rather comparable flux to that from M1 de-excitation processes. However we found that their event yield is small (see Fig. 9), and so in the ALP-nucleon coupling analysis we consider only p(n, d)γ. With large ALP fluxes, a large ALP event yield crucially depends on whether they can reach and decay within the detector -which is accounted for by the probabilities in Eq. (3.4). While a number of ALPs can as well interact with the detector shielding thus diminishing the expected ALP flux, for the sake of simplicity we neglect those interactions and we assume the interaction probability to be suppressed (P det → 0). In order to demonstrate the importance of the different detector variables, we have calculated differential event rates for the Primakoff and decay a → γγ cases (representative of models with dominant g aγγ coupling 8 ) varying the reactordetector distance L, the detector mass m det and the detector length L det . Results are shown in Fig. 7.
One can see that the reactor-detector distance plays an important role regardless of the process, as expected given the quadratic flux dependence on that variable. Moving from 10 m to 2 m enhances the event rate by a factor 25, due to the usual isotropic reduction factor 1/(4πL 2 ). For a + N → γ + N scattering, the detector mass plays an important role as well although its effect is less pronounced than the reactor-detector distance. For a → γγ decay the detector length is instead a key variable. As shown in the right graph in Fig. 7 Figure 7. Left graph: ALP differential event rate for the case of production and detection through Primakoff scattering (applicable to ALP models with dominant coupling only to photons) for different values of the reactor-detector baseline and detector masses. Right graph: Same as in left graph, but for detection through ALP decay, a → γγ. The reactordetector baseline as well as the detector length L det vary in this case.
large impact in the event rate, resulting in an enhancement of about two orders of magnitude. From these results one can fairly conclude that, if backgrounds are well understood and under control, ideally a short reactor-detector distance combined with a large detector mass and length would maximize the ALP event rate.
As we have already stressed, for scenarios with dominant g aee coupling three possible detection mechanisms are possible: Compton-like scattering, axio-electric absorption and ALP decay to electron-positron pairs. To compare their relative importance we have calculated their differential event rates for particular parameter choices. Results are shown in Fig. 8 plotted versus the relevant kinematic variable in each case. Compton-like and axio-electric absorption have about the same energy behavior, although in the former case the relevant variable is the photon energy while in the latter is the electron recoil energy. As it will become clear when discussing g aee and g aee × g (1) ann sensitivities, due to their close energy behavior they contribute similarly to the event yield. A comparison of these processes with ALP decay to electron-positron pairs shows that at large energies and for large ALP masses the decay dominates over the absorption and scattering processes. Inverse Primakoff compared with Compton-like scattering extends all over the photon energy range, but produces less events. This results in ALP-electron sensitivities being slightly better. The same applies for decays.
For scenarios with dominant g ann coupling the event yield depends on the detection mechanism. For illustration, we have calculated for particular parameter space point values the different yields. Results are displayed in Fig. 9. First of all one can see the difference between M1 and M4 transitions. In all cases the event yields from neutron capture in proton and 7 Li * de-excitation exceeds that from 91 Y * by more than one order of magnitude, a trend that applies for the other M4 transitions as well.
Although the ALP event rate from decay to final state photons has no kinematic constraint, its size does depend on the ALP mass as expected from Eq. (3.5). The absence of a kinematic limit allows all transitions to generate such a process, with a monochromatic event yield characterized by the energy of the corresponding transition. Decay to electron-positron pairs is instead subject to the kinematic constraint m a > 2m e , and so only p(n, γ)d can induce such process. Although involving different numerical values, comparison of the left upper graph with the lower left graph in Fig. 9 allows to see that a → γγ can play an important role, actually dominant once m a 10 4 eV (this will become apparent when discussing g aγγ experimental sensitivities). Comparison of the right upper graph and the middle and right lower graphs enables establishing the same conclusion, decay to electron pair dominates as soon as the channel is kinematically open. Comparison of the middle and right lower graphs in turn demonstrates that Compton-like scattering plays a more important role than axio-electric absorption does, with the event rate exceeding about an order of magnitude.
We now turn to sensitivities that could be achieved in models with dominant g aγγ coupling. With the production and detection specifications we have assumed, the region of interest that can be explored corresponds to g aγγ 10 −6 GeV −1 and m a ⊂ [1, 10 7 ] eV. In this region there are various laboratory as well as astrophysical and cosmological constraints that apply. Relevant laboratory limits in the "high" mass region follow from beam-dump experiments, invisible Υ decays and e + e − → inv + γ searches [10,11,[13][14][15]97]. In the "low" mass region instead limits from CAST and SUMICO are the most important [98,99]. Laboratory limits shown in Fig. 10 are axio-electric Figure 9. Event yields for ALP production through magnetic transitions and detection through: decays a → γγ and a → e + e − (left and right upper graphs), inverse-Primakoff scattering a + N → γ + N (left lower graph), Compton-like scattering a + e → γ + e (middle lower graph) and axio-electric absorption a + e + Z → γ + e + Z (right lower graph). The event yield has been calculated for a 4 GW reactor power plant assuming L = 10 m. Here, g eff ann = g (0) ann κ + g (1) ann .
adapted from Ref. [100]. Although our analysis applies for ALPs, for completeness we show the region of hadronic QCD axion models defined by E/N ⊂ [5/3, 44/3] [65]. Constraints from astrophysical sources are derived from SN1987A and horizontal branch (HB) stars [101][102][103][104][105][106]. SN1987A limits follow from SN energy loss, which in the absence of ALPs is driven by neutrino emission, and from the visible signal that ALP emission and subsequent decay to photons will generate. Given the temperature of the SN environment, these bounds extend up to ALP masses of the order of 10 8 eV. Horizontal branch stars limits follow, instead, from stellar cooling arguments. The presence of the new coupling opens a channel that accelerates the consumption of helium resulting in a reduction of the HB stars lifetime [107]. It is worth emphasizing that the SN limit is subject to large uncertainties, mainly due to the lack of a detailed understanding of the SN1987A. Actually it has been recently pointed out that it might not apply at all [108]. HB stars bounds as well as CAST+SUMICO [98,[109][110][111] limits can be substantially relaxed if the photon-ALP coupling and ALP mass depend on enviromental conditions such as the matter density and temperature, see e.g. [112] and references therein.
The triangular white region in Fig. 10, the so-called cosmological triangle, is subject as well to cosmological constraints of which BBN and N eff are particularly relevant [113][114][115]. Under standard cosmological arguments (ΛCDM) with a high inflationary scale that region is therefore ruled out. Departures from standard arguments, however, open this spot in parameter space. Possible scenarios include low-reheating temperature (T R ) models, extra contributions to the radiation energy density and non-vanishing neutrino chemical potentials [114,115]. Let us take the case of low-reheating models. If T R T fo , where T fo is the Primakoff freeze-out temperature, ALPs are abundant in the heat bath. Photons produced by the decay of this thermal ALP distribution will then dilute the baryon density and disrupt N eff , the former affecting BBN. If instead T R T fo -and ALPs are not produced by e.g. inflaton decays-ALPs will be much less abundant and so the amount of photons injected in the heat bath [115], thus having a small effect on BBN and N eff .
From Fig. 10 one can see that 90%CL sensitivities will cover regions where CAST and SUMICO as well as the SN 1987A and HB stars cooling arguments have placed limits. As we have mentioned, these bounds can be evaded if environmental effectsrelevant in stellar interiors-are at work. Reactor experiment measurements therefore will be able to test whether such effects are present. Sensitivities cover as well the cosmological triangle, thus implying that this type of measurements will test the cosmological hypotheses behind the constraints in that region. These searches therefore will be complementary to those that will be carried out at BELLE-II (50 ab −1 ) [116]. Undergoing and near-future experiments will partially cover that region and nextgeneration experiments will complement those measurements, allowing to cover the full region .
In models with dominant g aee coupling there are as well limits in the region that  reactor experiments can conver, g aee 10 −8 and in the same ALP mass range as that of the g aγγ case. In the "low" mass window Red Giants cooling arguments as well as ALP searches at the EDELWEISS III experiment place the most relevant bounds [117,118]. At "high" ALP masses beam dump experiments are important [119]. These limits along with the sensitivities that could be achieved in current, near-future or next generation experiments are displayed in Fig. 11. We show as well the region covered by the DFSZ-I and DFSZ-II QCD axion models.
From Fig. 11 one can see that the region m a ⊂ [5, 20] × 10 5 eV is not subject to constraints. Inverse Compton-like scattering and axio-electric absorption will allow testing couplings down to ∼ 10 −5 in that region, with the precise value determined by how well backgrounds are kept under control. Smaller values might be achievable in next generation experiments. Decay to electron-positron pairs extends the region that can be covered in that mass window, extending sensitivities down to couplings of order g aee ∼ 10 −8 . Inverse Compton-like scattering and axio-electric absorption cover as well regions subject to limits from Red Giants cooling arguments and EDELWEISS  Figure 11. Electron-ALP coupling 90%CL sensitivities along with relevant laboratory and astrophysical limits that apply in the region of interest. Black contours correspond to 90%CL sensitivities achievable with ongoing or near-future experiments (see Tab. 2). Purple contours instead indicate 90%CL sensitivites that could be achieved in a nextgeneration xenon detector. Solid (dashed-dotted) contours refer to sensitivities under the assumption of detection via inverse Compton-like scattering (axio-electric absorption) and different background rates. See the text for more details. For completeness we show as well the region of DFSZ-I and DFSZ-II QCD axion models.
III. Since the latter searches for solar ALPs, sensitivities in the region m a 5×10 5 eV will test as well possible scenarios (hypotheses) where environmental effects could allow circumventing those limits.
Finally our results for ALP models with dominant (tree level) g (1) ann and g aγγ couplings (g (1) ann and g aee couplings) are shown in the left graph (right graph) in Fig. 12 9 . Relevant constraints follow from TEXONO [33], Borexino [120,121] and BGO [122,123], that cover ALP mass values up to 10 6 eV. Constraints from solar ALP searches using 57 Fe first excited state and detection through axio-electric absorption apply in the case of the g aee × g (1) ann combination. Those limits follow from a variety of experiments, namely EDELWEISS III [118], PANDAX-II [124], CDEX [125] and the MAJORANA DEMONSTRATOR [126]. Due to kinematic constraints they cover the ALP mass range [1, 14.4 × 10 3 ] eV. In the case of the g aγγ × g (1) ann combination we include HB stars cooling and SN1987A energy loss limits, which should be interpreted as subject to the constraints the individual couplings  Tab. 3). Purple contours instead indicate 90%CL sensitivites that could be achieved in a next-generation xenon detector. Solid (dotted) contours refer to sensitivities under the assumption different background rates, (see the legend). Right graph: Same as in left graph but for electron-ALP times nucleon-ALP couplings sensitivities. Here solid (dotdashed) curves correspond to detection via inverse Compton-like scattering and decays into electron-positron pairs (axio-electric absorption), with different background rates. The solar axion region includes limits from EDELWEISS III, MAJORANA DEMONSTRATOR, CDEX and PandaX-II [118,[124][125][126]. For completeness we include the region of DFSZ-I and DFSZ-II QCD axion models.
are subject to, g aγγ 2 × 10 −6 GeV −1 (see Fig. 10) and g (1) ann 5 × 10 −6 [127]. As we have done in the previous cases, we include as well the region covered by the DFSZ-I and DFSZ-II QCD axion models.
In the case of detection through inverse Primakoff and ALP decay to photons, 90%CL sensitivity contours in the low mass range include regions covered by TEX-ONO and Borexino searches. For ALP masses above ∼ 4×10 5 eV and up to kinematic threshold ∼ 2 × 10 6 eV, ALP decay will widen searches in parameter space. Furthermore 90%CL sensitivities go inside some of the HB stars and SN1987A regions. These measurements therefore can partially-and potentially-test the hypotheses behind these constraints. Detection through g aee −related processes also cover Borexino and TEXONO regions, but more importantly will access an ALP mass stripe within [8,20] × 10 5 eV and down to couplings of order 10 −15 . We have also shown that next-generation experiments have the capability to substantially improve those sensitivities.
There is an infinite set of individual couplings g aXX (X = γ, e) and g (1) ann that satisfy the condition g aXX ×g (1) ann = const, which is not explicit in our results in Fig. 12. The interplay of the individual couplings on sensitivities can be seen in Fig. 13, which shows the 90%CL sensitivities calculated in the two-dimensional coupling plane for  Figure 13. 90%CL sensitivities in ALP production through magnetic de-excitation and detection through g aγγ and g aee related processes. Results are displayed in the couplings plane for a fixed ALP mass. See text for more details. Bounds on the isovector axionnucleon coupling g (1) ann are adapted from [128]. For completeness we include the region of DFSZ-I and DFSZ-II QCD axion models.
fixed ALP mass along with the limits the individual couplings are subject to. For illustrative purposes, we have fixed m a = 1.2 MeV (10 keV) for the left (right) graphs in Fig. 13, depicting sensitivities for ALPs detected via decay (scattering) processes, respectively. Indeed, decays are more relevant at higher masses, while scattering processes lead to almost constant number of events for masses m a 10 4 (10 5 ) eV for Primakoff (Compton-like, axio-electric) detection mechanisms. From these plots one can see that most of the parameter space inside these sensitivities is already probed by a variety of different experiments, from laboratory to astrophysics and cosmology [128,129]. In the case of decays however, reactor searches could probe some small regions of parameter space currently uncovered, at g (1) ann 10 −9 and g aee ∼ 10 −6 , or in the "cosmological triangle" region at g (1) ann ∼ 10 −6 and g aγγ ∼ 5 × 10 −6 . In the rest of parameter space, reactor searches would still give valuable information, by providing complementary laboratory probes especially in the regions currently tested only by astrophysical observations.

Conclusions
Nuclear reactor power plants produce a large amount of photons and so they are suitable sources for ALPs production. In this paper we have studied such possibility by considering production through three possible mechanisms: Compton-like and Primakoff scattering as well as nuclear de-excitation in the fuel material. Motivated by the current CEνNS experimental program, which aims at observing CEνNS induced by reactor neutrino fluxes using various low-threshold technologies, we have analyzed the prospects for detection of ALPs in such experiments. We have assumed a 10 kg germanium detector, which we regard as representative of current (CONUS [36]) or near-future CEνNS experiments (including MINER [35], Ricochet [38], νGeN [70] and TEXONO [40]). We have also explored future capabilities of next generation experiments, assuming a ton-scale xenon detector.
Concerning ALPs detection, we have considered processes controlled by the same couplings that determine their production. Namely, we have focused on inverse-Primakoff scattering and ALP decays into photons for production through Primakoff scattering (processes controlled by g aγγ ), while we have investigated inverse Compton-like scattering, axio-electric absorption and decay into electron pairs for ALPs produced through Compton-like scattering (processes controlled by g aee ). For nuclear de-excitation we have instead considered detection through processes controlled by g aγγ and g aee couplings.
We have calculated 90%CL sensitivities to ALPs in the mass range (1 − 10 7 eV), and we have confronted them to existing constraints. Our main results are displayed in Figs. 10-13, and they can be summarized as follows. For models with dominant g aγγ coupling, these experiments will probe a region in the ALP mass range around 1 MeV, presently constrained only by cosmological observations. Moreover, since reactor searches will cover regions currently constrained by CAST+SUMICO, SN 1987A and HB stars cooling arguments, they can also potentially test ALPs environmentdependent properties in stellar media. For ALP models with dominant ALP-electron coupling, reactor searches will probe ALP masses in the range [5, 20] × 10 5 eV and ALP-electron couplings extending down to g aee ∼ 10 −8 . They will also cover regions constrained by Red Giants cooling arguments as well as by EDELWEISS III solar axions searches, thus testing also in this case the potential relevance of ALP environmental effects.
Searches for ALPs produced in nuclear de-excitation require the presence of two different ALP-SM couplings. If detection proceeds through ALP-photon coupling induced-processes, they can improve current limits set by Borexino solar axion searches and by the TEXONO reactor experiment in the ALP mass range ∼ [0.2, 2] MeV. Moreover these searches could probe regions constrained by SN1987A data and so the different hypotheses this limit comes along with. Detection through ALP-electron coupling induced-processes will cover regions already explored by other experiments such as BGO, Borexino and TEXONO. Still, a small region in parameter space (currently unexplored) with ALP masses ∼ MeV and couplings g (1) ann ∼ 10 −9 could be tested.
To summarize, we have considered the possibility that ongoing and near-future reactor-based CEνNS experiments can be used for ALP searches. Of course these searches could be run in any other reactor-based experiment. However since CEνNS experiments will run anyway, they will probably be the most suited environment for these type of searches. We have shown that sensitivities can cover interesting regions in the ALP parameter space not yet explored by laboratory experiments, thus motivating the inclusion of ALP searches in the physics program of such reactor neutrino experiments.