Search for $C$ violation in the decay $\eta\rightarrow\pi^0+e^++e^-$ with WASA-at-COSY

We report on the investigation of the rare decay $\eta\rightarrow\pi^0+\text{e}^++\text{e}^-$ which is of interest to study both $C$ violation in the electromagnetic interaction and to search for contributions from physics beyond the Standard Model, since the allowed decay via a two-photon intermediate state is strongly suppressed. The experiment has been performed using the WASA-at-COSY installation, located at the COSY accelerator of the Forschungszentrum J\"ulich, Germany. In total $3\times10^7$ events of the reaction $\text{p}+\text{d}\rightarrow{^3\text{He}}+\eta$ have been recorded at an excess energy of $Q = 59.8\,\textrm{MeV}$. Based on this data set the $C$ parity violating decay $\eta\rightarrow\pi^0+\gamma^*\rightarrow\pi^0+\textrm{e}^++\textrm{e}^-$ via a single-photon intermediate state has been searched for, resulting in new upper limits of $\Gamma\left(\eta\rightarrow\pi^0+\textrm{e}^++\textrm{e}^-\right)/ \Gamma\left(\eta\rightarrow\pi^++\pi^-+\pi^0\right)<3.28\times10^{-5}$ and $\Gamma\left(\eta\rightarrow\pi^0+\textrm{e}^++\textrm{e}^-\right)/ \Gamma\left(\eta\rightarrow\textrm{all}\right)<7.5\times10^{-6}\ (\textrm{CL} = 90\,\%)$, respectively.


Introduction
All strong and electromagnetic decays of the η meson are either suppressed or forbidden to first order. The η meson is, in addition, a C and P eigenstate of strong and electromagnetic interaction. This makes it well suited for the study of rare processes and the search for forbidden ones. The subject of this letter is the process η → π 0 + e + + e − via the single-photon intermediate state η → π 0 + γ * that would violate C parity conservation. The background for this process would be a twophoton process with an expected branching ratio not larger than 10 −8 according to theoretical calculations [1,2,3]. The present experimental upper limit for the branching ratio of the decay η → π 0 + e + + e − is from the seventies of the last century and amounts only to 4.5 × 10 −5 (CL = 90 %) [4]. A more stringent upper limit for the decay channel η → π 0 + e + + e − has been determined in the analysis presented in this paper. The data have been collected using the WASA-at-COSY facility and also constituted the basis for studies of other η meson decay channels already published in Ref. [5].

Experiment
The WASA-at-COSY experiment was an internal experiment operated at the accelerator COSY of the Forschungszentrum Jülich, Germany from 2006 to 2014 [6]. For the measurements discussed here, a proton beam was accelerated to a kinetic beam energy of T p = 1 GeV and collided with deuterium pellets provided by the internal pellet target. The η mesons were produced in the reaction p + d → 3 He + η.
The WASA detector setup consists of two main parts: the central detector, which was used for the reconstruction of the produced mesons and their decay particles, and the forward detector used for the measurement of the four momenta of the forward scattered 3 He nuclei. A more detailed description of the WASA-at-COSY experimental setup can be found in Ref. [5,6,7].
The data for the studies presented here were obtained in two measurement periods, one of four weeks in 2008 and one of eight weeks in 2009. For data acquisition the trigger used required a large energy loss in subsequent scintillator elements of the forward detector. Since the 3 He nucleus stemming from the reaction p + d → 3 He + η is stopped in the first layer of the WASA forward range hodoscope, a veto on the signals from the second layer was used in addition. Due to the trigger relying on information from the forward detector only, the utilized trigger was unbiased with respect to a decay mode of the η meson. In total about 3 × 10 7 events containing an η meson were recorded with 1×10 7 events originating from the 2008 period and 2×10 7 events from the 2009 period [5].

Data analysis
The analysis of the decay η → π 0 + e + + e − was based on a common analysis chain for η decay studies described in Ref. [5]. Since only very few events were expected to remain in the analysis after the event selection, an optimal choice of the selection conditions is important for the best possible result. These conditions were determined with the aid of an optimization algorithm based on Monte Carlo simulations.
Preselection. Before the selection conditions for the decay η → π 0 +e + +e − were determined, the data collected in 2008 and 2009 were preselected with conditions common to all recorded reactions. For instance, conditions on time correlations of the measured particles were used, as presented in Ref. [5]. Furthermore, to reject hits from particles that were wrongly identified as secondary particles (so-called split-offs) and electronpositron pairs from conversion of photons at the COSY beam pipe, two-dimensional cuts were utilized. More details of these conditions were published in Ref. [5].
Besides these general preselection conditions, a cut on the signature of the decay η → π 0 +e + +e − was applied requesting at least one positively and one negatively charged particle detected in the central detector, as well as at least two neutral particles originating from the π 0 meson decay π 0 → γ + γ. The last condition applied for data preselection requires the maximum considered momenta of the charged decay particles to be below p = 250 MeV/c, since the momenta of the leptons of the decay η → π 0 + e + + e − are expected to be below this value.
Monte Carlo simulations. In order to determine optimal selection conditions for the search for the decay channel η → π 0 + e + + e − , 1.8 × 10 8 Monte-Carlo events of all non-signal η decays observed yet were created with respect to their relative branching ratio [8], as well as two million events for the signal decay. These simulations were generated with the pluto++ software package [9] considering the angular distribution of p + d → 3 He + η at T p = 1 GeV according to Ref. [10]. For the various η decay channels physics models as included in  Figure 1: Invariant mass of e + e − pairs for the simulated decay η → π 0 +e + +e − . Black lined: decay via η → π 0 + γ * considering VMD. Shadowed in orange: decay according to three-particle phase space.
pluto++ were used. The reader is refered to Ref. [5] for further details.
In addition to the simulations of η decays, about 4.3 × 10 9 events for the direct pion production were created, with most events for the production reactions p + d → 3 He + π 0 + π 0 and p + d → 3 He + π + + π − , as these contribute most to the non-η background at the given kinetic beam energy. For these two-pion productions the ABC effect was incorporated into the simulations according to the model discussed in Ref. [11].
The simulations for the signal decay η → π 0 + e + + e − were generated with two different model assumptions. The first one is a decay according to pure three-particle phase space. The second is based on the vector meson dominance (VMD) model for the intermediate virtual photon. The direct decay η → π 0 + γ violates both C parity and angular momentum conservation plus global gauge invariance. Thus, there is no η → π 0 + γ on-shell contribution for the decay η → π 0 + e + + e − and the transition form factor for the off-shell contribution vanishes at zero virtuality, such that the single-photon pole is completely removed [12,13,14]. In Fig. 1 the invariant mass of the e + e − pair produced in the decay is plotted according to three-particle phase space (shadowed in orange) and the decay via η → π 0 +γ * according to the discussed model. A more detailed calculation of the model can be found in Ref. [15].
To simulate the WASA detector responses, the WASA Monte Carlo package wmc was used, which is based on geant3 [16]. The settings for the spatial, timing and energy resolution in wmc were set to agree with the resolution observed in data.
Due to the high luminosities of the WASA-at-COSY experiment, it is possible that detector responses from one event can overlap with another event. This effect was considered in the simulations and the amount of event overlap was left as a free parameter for the fit of the simulations to data (see next paragraph).
All Monte Carlo simulations were preselected with conditions identical to those for data preselection.
Data description. The choice of the selection conditions with regard to the decay channel η → π 0 + e + + e − is based on Monte Carlo simulations. It is necessary to know the contributions of the various reactions to the collected data for an optimal choice. Therefore, the 2008 and 2009 data sets were fitted separately in distributions of selected quantities by template distributions of the aforementioned Monte Carlo simulations to determine the contributions of the individual reactions to the data. In detail, these distributions are: • the missing mass m X , corresponding to the invariant mass of the proton beam and the deuteron target remaining after the 3 He four momentum has been subtracted and peaks at the η mass for the reaction p + d → 3 He + η, • the invariant mass m eeγγ of an electron-positron pair candidate and two photons, which peaks at the η mass for the decay η → π 0 + e + + e − with π 0 → γ + γ, • the invariant mass m γγ of two photons, which peaks at the π 0 mass for reactions with π 0 mesons produced, • the invariant mass m ee of an electron-positron pair candidate, • the smallest invariant mass m eγ of all four possible combinations of an electron or positron candidate and a photon and • the missing mass squared m 2 Xee , which is the invariant mass squared of the proton beam and the deuteron target remaining after the 3 He four momentum and the electron-positron pair candidate momentum have been subtracted and peaks at the π 0 mass squared for the reaction of interest.
Under the assumption of a branching ratio of the decay below the current upper limit of 4.5 × 10 −5 (CL = 90 %) [4], there are less than 150 events expected from the decay η → π 0 + e + + e − in the combined data sets after preselection, considering the preselection efficiency for the signal decay. A fit by Monte Carlo simulations including the simulated decay η → π 0 + e + + e − is consistent with zero events from this signal decay channel. Therefore, the decay η → π 0 + e + + e − was excluded from the fit. While the differential distribution for the reaction p + d → 3 He + η is well known [10], the differential distributions are known only with high uncertainties or not at all for direct multipion productions. Hence, the data were divided into ten bins in angular ranges of cos ϑ cms 3 He 5 . Monte Carlo simulations were fitted to data in the eight angular bins ranging from −1 to 0.6. The angular range 0.6 < cos ϑ cms 3 He ≤ 1 was excluded because of the lower energy resolution of the forward detector for these forward scattered 3 He nuclei. Moreover, the relative amount of background from the direct pion production is larger in this angular range, whereas less than 3 % of all p + d → 3 He + η events have a cos ϑ cms 3 He > 0.6.    The fit of the Monte Carlo simulations to the data was performed simultaneously for all angular ranges and distributions with identical scaling parameters for the simulations for all distributions within one angular range. Furthermore, the ratios for the various η decays were constrained to the branching ratios according to Ref. [8] within the given uncertainties. These were set to be identical for all angular ranges. Similarly, the amount of event overlap was included as one global fit parameter. In Fig. 2, Fig. 3, Fig. 4 and Fig. 5 the resulting Monte Carlo fits to the 2008 data are plotted for m X , m eeγγ , m γγ and m ee for the angular range 0.2 < cos ϑ cms 3 He ≤ 0.4. According to this fit most events remaining after preselection originate from the η decay η → π + + π − + π 0 , the direct p + d → 3 He + π + + π − + π 0 production and the direct two-pion production reactions. A collection of all fits is available in Ref. [15].
Selection conditions. The selection conditions for the search for the decay η → π 0 + e + + e − were based on the following quantities: • the missing mass m X to identify the production reaction p + d → 3 He + η, • the invariant mass m eeγγ of an electron-positron pair candidate and two photons to select the decay η → π 0 + e + + e − → γ + γ + e + + e − , • the invariant mass m γγ of two photons to ascertain the decay π 0 → γ + γ, • the invariant mass m ee of an electron-positron pair candidate, • the χ 2 probability of a kinematic fit with the hypothesis p + d → 3 He + γ + γ + e + + e − and • the energy loss E SEC dep of the charged particles in the central detector scintillator electromagnetic calorimeter (SEC) and their momentum p to discriminate e ± and π ± (particle identification, PID). The choice of the cut conditions was performed with 40 % of the generated Monte Carlo simulations, whereas the remaining Monte Carlo data sample was used later for the selection efficiency determination. While the graphical cut for the particle identification was chosen beforehand (see Fig. 6), as it is a common cut utilized for PID independent from the analyzed reaction, the selection conditions for the other five quantities were determined by an optimization algorithm. This algorithm is based on the relative amount of simulated signal events S R = N cut S /N The cut optimization algorithm maximizes the evaluation function by varying the selection conditions for all chosen quantities. With the aid of the cut optimization algorithm the following selection conditions were determined: m ee ≥ 0.096 GeV/c 2 and (5)

Results
After applying the selection conditions to the data, three events were left, whereas two events were expected to remain from the direct two-pion production p + d → 3 He + π 0 + π 0 according to Monte Carlo simulations. All other background reaction channels were found to give no sizeable contribution after applying the cuts. The invariant mass, m eeγγ , for these events are plotted in Fig. 7 together with simulated data. Note that the generated Monte Carlo events were scaled according to the fit to data after preselection and that the sum of all Monte Carlo events remaining after all cuts is equal to two events. The overall reconstruction efficiency for the signal decay η → π 0 + e + + e − was determined to be for a decay via η → π 0 + γ * assuming VMD, whereas the assumption of a decay according to pure three-particle phase space results in ε phase S = 0.01844 (7).
The given uncertainties are purely statistical ones. In order to calculate the upper limit for the branching ratio Γ(η → π 0 + e + + e − )/Γ(η → all), the decay channel η → π + + π − + π 0 with π 0 → γ + γ was utilized for normalization. This is a reasonable choice as this decay channel has the same signature as the signal decay and, thus, possible systematic effects introduced by differences of the signature are avoided. According to the data description by Monte Carlo simulations there were N produced η→π + π − π 0 γγ = (6.509 ± 0.018) × 10 6 (9) events in data, considering already the efficiency correction determined by Monte Carlo studies. In order to determine a final upper limit for the branching ratio of η → π 0 + e + + e − , all uncertainties have to be considered and incorporated into the calculations.
Systematics. The systematic and statistical uncertainties, which need to be considered for the upper limit determination, can be separated into uncertainties by multiplicative effects and uncertainties by offset effects. The former include an uncertainty of the reconstruction efficiency of the decay η → π 0 + e + + e − and an uncertainty in the number of η → π + + π − + π 0 → γ + γ events in data. The latter ones are uncertainties of the number of background events remaining after all cuts.
To determine the systematic uncertainty for the signal reconstruction efficiency, the resolution settings for the Monte Carlo simulations were varied within the uncertainties of the individual detector resolutions observed in data. The extracted square root of the relative variance of the reconstruction efficiency was found to be Var virtual rel = 0.059 (10) for a decay via η → π 0 +γ * assuming VMD whereas for a decay according to pure three-particle phase space one finds Var phase rel = 0.057.
In the following analysis the square root of the variance was considered as the systematic uncertainty. The uncertainty for the efficiency corrected number of η → π + + π − + π 0 → γ + γ events in data was obtained by a comparison to the efficiency corrected number determined utilizing less strict preselection conditions, namely no cuts to reject conversion or split-off events, no cut on the momentum of charged decay particles and less strict cuts on the particles' energies. Hereby a systematic uncertainty of 2.3 % was determined.
The uncertainties for the number of background events remaining after all cuts can be separated into a statistical uncertainty due to the finite number of Monte Carlo simulations and systematic uncertainties introduced by uncertainties of the fit of Monte Carlo simulations to data. The latter are dominated by differences between the Monte Carlo fit parameters for the 2008 and 2009 data sets, leading to asymmetric uncertainties. Such different fit parameters for both data sets originated mainly from different experimental settings, which affected, e.g., the event overlap due to different luminosities. To determine the overall systematic uncertainty for the number of remaining background events, the probability density functions (pdf) of the individual uncertainties were folded. The resulting pdf for the nuisance parameters λ 2008 and λ 2009 corresponds to the overall relative systematic uncertainty for the 2008 and 2009 data sets and was incorporated into the upper limit calculations. In Fig. 8 the distribution of the nuisance parameter is illustrated for the 2008 data set.
In order to investigate further possible systematic effects, the selection conditions used for the analysis were varied and the expectations according to simulations were compared to the number of events seen in data. Since the expected number of events agreed with the number of events seen in data within the statistical uncertainties, no additional systematic effect needs to be considered.
A detailed description of the uncertainty investigations is available in Ref. [15].
Upper limit. The upper limit for the relative branching ratio of the decay η → π 0 + e + + e − was calculated with the formula: with the upper limit N S,up for the number of signal events, which depends on the number of observed events and the number of expected background events. For the calculation of N S,up a Bayesian approach was chosen as given in Ref. [17] with a flat prior pdf and incorporating the determined uncertainties and the pdfs for the nuisance parameters.
These values are smaller than the previous upper limit of 4.5 × 10 −5 (CL = 90 %) [4] by a factor of six and five, respectively.

Summary
We have presented new studies with the WASA-at-COSY experiment on the C parity violating η meson decay η → π 0 + e + + e − . The obtained upper limit for the branching ratio of the decay η → π 0 + e + + e − is smaller than the previously available upper limit by a factor of five to six [4]. The results of the analysis are consistent with no events seen in data, and thus give no hint on a C violation in an electromagnetic process. Similarly, no processes from physics beyond the Standard Model are required to explain the results.
In order to further decrease this value and to continue the search for a C parity violation in an electromagnetic process, additional data were collected with WASA-at-COSY utilizing the production reaction p + p → p + p + η. Over three periods in 2008, 2010 and 2012 in total about 5 × 10 8 such events were recorded and are currently being analyzed with regard to the decay η → π 0 + e + + e − .
Besides a decay via one virtual photon, the decay η → π 0 + e + + e − could possibly occur via a hypothetical C violating dark boson U with m U < 413 MeV/c 2 where the pertinent form factor is even further suppressed (i.e. the second term in its Taylor expansion vanishes) compared with the single-photon mechanism [18]. Investigations with regard to this decay process are currently ongoing for the presented p + d → 3 He + η data sets and the p + p → p + p + η data sets recorded with WASA-at-COSY.