Study of the process $e^+e^- \to \eta\pi^0\gamma$ in the energy range $\sqrt{s} = \mbox{1.05-2.00}$ GeV with the SND detector

The process $e^+e^-\to\eta\pi^0\gamma$ is studied in the center-of-mass energy range 1.05-2.00 GeV using data with an integrated luminosity of 94.5 pb$^{-1}$ collected by the SND detector at the VEPP-2000 $e^+e^-$ collider. The $e^+e^-\to\eta\pi^0\gamma$ cross section is measured for the first time. It is shown that the dominant mechanism of this reaction is the transition through the $\omega\eta$ intermediate state. The measured cross section of the subprocess $e^+e^-\to\omega\eta\to\eta\pi^0\gamma$ is consistent with previous measurements in the $e^+e^-\to\pi^+\pi^-\pi^0\eta$ mode. It is found, with a significance of 5.8$\sigma$, that the process $e^+e^-\to\eta\pi^0\gamma$ is not completely described by hadronic vector-pseudoscalar intermediate states. The cross section of this missing contribution, which can originate from radiation processes, e. g. $e^+e^-\to a_{0}(1450)\gamma$, is measured. It is found to be 15-20 pb in the wide energy range from 1.3 to 1.9 GeV.


Introduction
This work is devoted to study of the process e + e − → ηπ 0 γ (1) in the center-of-mass energy range √ s = 1.05-2.00 GeV at the experiment with the SND detector at the VEPP-2000 e + e − collider. Previously, this process was studied near the φ-meson resonance by the SND at VEPP-2M [1], CMD-2 [2] and KLOE [3]. The dominant intermediate mechanism in this energy region is the decay φ → a 0 (980)γ. Below ( √ s = 0.920-1.004 GeV) and above ( √ s = 1.03-1.38 GeV) the φ-meson resonance the process (1) was studied by CMD-2 in Ref. [4], where a 90% confidencelevel upper limit of about 0.1 nb was set on the cross section. At higher energy, there is only the BESIII measurement of the J/ψ → ηπ 0 γ decay [5]. The ηπ 0 mass spectrum in this decay is well described by an uniform phase-space distribution. No significant signal from the J/ψ decays to a 0 (980)γ and a 2 (1320)γ was observed.
The dominant contribution to the e + e − → ηπ 0 γ cross section in the energy region under study comes from the process e + e − → ωη with the decay ω → π 0 γ. This process was measured in the BABAR [6], CMD-3 [7] and SND [8] experiments with the decay mode ω → π + π − π 0 .
In these works, the measurement was performed neglecting the interference between the ωη and other intermediate mechanisms (a 0 (980)ρ, ρ(1450)π and φη) contributing to the e + e − → π + π − π 0 η reaction. The measurements obtained under this assumption need additional verification, especially in the region √ s = 1.85-2.00 GeV, where the e + e − → ωη cross section is almost zero (< 50 pb) compared to the significant (∼ 2 nb) contribution from other mechanisms.
In this work, the most interesting is the search for radiation decays of excited vector mesons of the ρ, ω and φ families to a 0 (980)γ, a 2 (1320)γ, and a 0 (1450)γ. The measurement of these decays is important for understanding the quark structure of excited vector mesons. In particular, there are indications that the excited states of the ρ and ω mesons may contain an admixture of a vector hybrid state [9]. The widths of the radiation decays are sensitive to the hybrid admixture [10].

Detector and experiment
SND is a general-purpose non-magnetic detector [11] collecting data at the VEPP-2000 e + e − collider [12]. Its main part is a three-layer spherical electromagnetic calorimeter arXiv:2006.05465v1 [hep-ex] 9 Jun 2020 √ s = 1.05-2.00 GeV consisting of 1630 NaI(Tl) crystals. The calorimeter covers a solid angle of 95% of 4π. The energy resolution of the calorimeter for photons is σ E /E = 4.2%/ 4 E(). The angular resolution is about 1.5 • . Directions of charged particles are measured using a nine-layer drift chamber and one-layer proportional chamber in a common gas volume. The solid angle of the tracking system is 94% of 4π. A system of threshold aerogel Cherenkov counters located between the tracking system and the calorimeter is used for charged kaon identification. Outside the calorimeter, a muon detector consisting of proportional tubes and scintillation counters is placed. Simulation of the signal and the background processes takes into account radiative corrections [13]. The angular distribution of hard photon emitted from the initial state is generated according to Ref. [14]. Interactions of the particles produced in e + e − annihilation with the detector materials are simulated using the GEANT4 software [15]. The simulation takes into account variation of experimental conditions during data taking, in particular dead detector channels and beam-induced background. To take into account the effect of superimposing the beam background on the e + e − annihilation events, simulation uses special background events recorded during data taking with a random trigger. These events are superimposed on simulated events, leading to the appearance of additional tracks and photons in events.
The analysis presented in this work is based on data with an integrated luminosity of 94.5 pb −1 recorded in 2010, 2011, 2012 and 2017. These data were collected at 101 energy points in the energy region √ s = 1.05-2.00 GeV. Since the cross section of the process under study is small and relatively slowly changes with energy, the data are combined into 13 energy intervals shown in Table 2.
In this work, the process e + e − → ηπ 0 γ is studied in the five-photon final state. Therefore, it is viable to use the process e + e − → γγ for normalization. As a result of the normalization a part of systematic uncertainties associated with the hardware event selection and the beam background are canceled out. Accuracy of the luminosity measurement using the process e + e − → γγ is estimated to be 2%.

Event selection
Selection of e + e − → ηπ 0 γ → 5γ events is performed in two stages. At the first stage, we select events with exactly 5 photons with energies above 20 MeV and no charged tracks. The latter condition is ensured by requiring that the number of hits in the drift chamber is less than four. The conditions on the total energy deposition in the calorimeter (E EMC ) and the total event momentum (P EMC ) calculated using energy depositions in calorimeter crystals are imposed: To suppress cosmic-ray background, absence of a signal from the muon system is required.
The main background process is A noticeable contribution to the background comes from the QED processes e + e − → 3γ, 4γ, 5γ.
We also study background from the following reactions with multiphoton final states: Additional photons in events with three and four photons in the final state arise from splitting of the electromagnetic showers, initial state radiation, and beam background. To suppress background from the processes listed above, a kinematic fits are performed to the hypotheses e + e − → 3γ, e + e − → 5γ, e + e − → π 0 π 0 γ, and e + e − → ηπ 0 γ with the requirement of energy and momentum balance in an event. For the two latter hypotheses, the additional constraints are imposed that the invariant masses of photon pairs are equal to the masses of the π 0 and η mesons. As a result of the kinematic fit, the energies and angles of photons are refined, and the χ 2 of the proposed kinematic hypothesis is calculated. In the kinematics fits, all possible combinations of photons are tested, and the combination with the smallest χ 2 value is retained. The following conditions are applied on the obtained χ 2 values About 500 events are selected using these conditions. To estimate the background, along with the signal region determined by the conditions (6), a control region is analyzed, for which the modified condition on the χ 2 difference 10 < χ 2 ηπ 0 γ − χ 2 5γ < 60 is used.
The distributions for the four classes of events used in the fit are obtained by simulation. For the ωη class, to take into account the imperfect simulation of the ω meson line shape, a mass shift ∆m and Gaussian smearing with the dispersion ∆σ 2 are introduced into the π 0 γ mass spectrum. These parameters are determined from comparison of the ω meson peak position and width in data and simulation from the range √ s = 1.41-1.80 GeV. They are found to be ∆m = 4.0 ± 1.6 MeV and ∆σ 2 = −48 ± 48. Since ∆σ 2 agrees with zero, the π 0 γ mass resolution is not corrected.
To describe the res-ηπγ contribution, the simulation of the processes (8) is used. The distributions obtained from simulation are normalized to expected number of events and summed up. To calculate the expected number of events, we used the cross sections for the processes (8) measured in Refs. [16,17,18,19] and the ρ → π 0 γ, φ → π 0 γ, and ηγ branching ratios [20]. The contribution of the processes (9) estimated using Ref. [21,22] is found to be negligibly small. The total res-ηπγ contribution is calculated to be N res = 19.6 ± 1.3.
To describe the rad-ηπγ contribution, the simulation of the processes e + e − → a 0 (1450)γ and e + e − → a 2 (1270)γ is used. The distributions for these processes are summed with the weights 1 − α and α, respectively, where the coefficient α = 0.22 ± 0.21. The choice of the value of this coefficient is discussed below in Sec. 5. It should be noted that the fit results very weakly depend on the value of α.
The m πγ and χ 2 ηπ 0 γ −χ 2 5γ distributions for background are calculated using simulation of the processes (3)(4)(5). The distributions obtained for each process are normalized to the expected number of events and summed up. Only the shape of the distribution is used in the fit. The total number of background events is a free fit parameter.
The following relations between the numbers of signal and background events in the signal and control regions are used in the fit: N bkg = k bkg N C bkg , N C ηπγ = k sig N ηπγ , where N ηπγ and N bkg are the numbers of e + e − → ηπ 0 γ and background events in the signal region, respectively, and N C ηπγ and N C bkg are the same numbers in the control region. The coefficients k sig and k bkg are calculated using simulation in each of the 13 energy intervals. The uncertainty of k bkg is estimated by varying the cross sections of the background processes (3-5) within their errors. The value of this coefficient averaged over the energy range √ s = 1.05-2.00 GeV is k bkg = 0.53 ± 0.01. The coefficient k sig obtained from simulation is corrected to take into account the difference in the signal χ 2 ηπ 0 γ −χ 2 5γ distributions in the data and simulation. To do this, data from the energy range √ s = 1.41-1.80 GeV are used, where the cross section of the process e + e − → ωη is maximal. The mass spectrum of π 0 γ is analyzed and the number of data events in the ω meson peak is determined for the signal (N data S ) and control (N data C ) regions. The numbers of simulated e + e − → ωη events in the signal (N MC S ) and control (N MC C ) regions are also determined. The double ratio R = (N data 3 is used to correct the coefficient k sig . With this correction, this coefficient averaged over the full energy range is k sig = 0.22 ± 0.04. The simulation shows that the difference between the χ 2 ηπ 0 γ − χ 2 5γ distributions for different intermediate mechanisms of the process e + e − → ηπ 0 γ is very small. To evaluate the systematic uncertainty associated with the procedure of determining the number of events in each class, the nuisance parameters corresponding to N res , ∆m, ∆σ 2 , and the coefficients k sig and k bkg are introduced into the likelihood function with Gaussian constraints. The fit results are represented by the solid histogram in Fig. 1 (left) for the m πγ distributions, and by the dashed histogram for the χ 2 ηπ 0 γ − χ 2 5γ distribution in Fig. 1 (right). The latter take into account the scale factor R for the signal events described above, while the solid histogram in Fig. 1 shows the uncorrected distribution normalized to the number of events in the signal region. The following numbers of events of three classes are obtained in the fit: The total contribution from the background processes (3-5) calculated using simulation is N calc bkg = 123 and agrees well with the fit result.
The obtained numbers of ωη and rad-ηπγ events in 13 energy intervals are listed in Tables 2 and 4, respectively.  Table 4 also shows the distribution over the energy intervals for background events and events of the class res-ηπγ.

Model selection for the rad-ηπγ event class
As it is mentioned in the previous section, the following three processes may contribute to the rad-ηπγ class in the energy range 1.05-2.00 GeV:   η ω γ π η radγ π η resbkg Fig. 2. The ηπ 0 invariant mass distribution for data (points with error bars) from the interval √ s = 1.05-2.00 GeV. Events from the ω meson peak (700 < mπγ < 900 MeV) are rejected. Left panel: The shaded histogram shows the total distribution for the background and the processes (7) and (8). The solid, dotted and dashed histograms show the total distribution for events of all four classes, in which the distribution for the class rad-ηπγ are calculated in the models a0(1450)γ, a2(1270)γ, and a0(980)γ, respectively. Right panel: The contributions of the four classes of events are shown cumulatively. The distribution for the class res-ηπγ is obtained from the fit to the two-dimensional distribution of mπγ versus mηπ as described in the text. total distribution for the background and the processes (7) and (8). The solid, dotted and dashed histograms show the total distribution for events of all four classes, in which the distribution for the class rad-ηπγ is calculated in the models a 0 (1450)γ, a 2 (1270)γ, and a 0 (980)γ, respectively. The distributions are normalized to the numbers of events (10) found in the previous section. It is seen the observed ηπ 0 mass spectrum is best described by the model e + e − → a 0 (1450)γ.
The following model is tested to describe distributions for rad-ηπγ events: where P is, for example, a two-dimensional distribution of m πγ versus m ηπ . This simple model does not take into account the interference between the three intermediate mechanisms, but can be used to estimate the model uncertainty of the efficiency for rad-ηπγ events.
The data from the range √ s = 1.05-2.00 GeV are fitted as described in the previous section. To increase the sensitivity to the model of the intermediate states (12), the m πγ distribution for the signal region is replaced by the two-dimensional distribution of m πγ versus m ηπ . The parameters α and β are determined from the fit with the constraints that α, β and α+β vary from 0 to 1. The following values of these parameters are obtained: α = 0.22 ± 0.21 and β = 0.00 +0.08 . The m ηπ distribution obtained using simulation with these parameters is shown in Fig. 2(right).
To evaluate the significance of the rad-ηπγ signal, we compare the values of likelihood function for the fit described above (L 1 ) and the fit with the N rad ≡ 0 (L 0 ). Taking into account that the numbers of parameters in these two fit differ by three, from the difference −2(ln L 0 − ln L 1 ) = 42.3, we obtain that the significance of the observed rad-ηπγ signal (including the systematic uncertainty) is 5.8 σ.

Detection efficiency and radiative corrections
The simulation of the process e + e − → ηπ 0 γ includes initial state radiation (ISR). The detection efficiency for the process under study is calculated as a function of two parameters, √ s and normalized energy of the ISR photon, x = 2E γ / √ s, and is parametrized as follows ε r (s, x) = ε(s)g(s, x), where ε(s) ≡ ε r (s, 0). The function g(s, x) weakly depends on √ s and on the intermediate mechanism of the process e + e − → ηπ 0 γ. Its dependence on x is shown in Fig. 3. The shape of this dependence is determined by two effects. The sharp decrease of the efficiency at x corresponding to E γ = 20-30 MeV is due to the condition that an event contain exactly five photon. This condition rejects events with the ISR photon emitted at a large angle. Photons are reconstructed if their energy deposition in the calorimeter exceeds 20 MeV. A further decrease of the efficiency is due to the requirement of energy-momentum balance in an event (χ 2 5γ < 30). The visible cross section of the process e + e − → ηπ 0 γ, which is defined as σ vis = N/L, where N is the number of selected events of the process under study and L is the integrated luminosity, is related to the Born cross section σ(s) as follows: where F (x, E) is a so-called radiator function [13] describing the probability to emit extra photons with the total energy x √ s/2 from the initial state, and x max = 1 − (m η + m π 0 ) 2 /s. The formula (13) can be rewritten in the conventional form: where δ(s) is the radiative correction. Inaccuracy in simulation of the detector response for photons leads to systematic uncertainty in the detection efficiency determined using the simulation. To evaluate the efficiency corrections associated with the selection criteria, we change the boundaries of the conditions: on χ 2 5γ from 30 to 60, in χ 2 ηπ 0 γ −χ 2 5γ from 10 to 60, and on χ 2 π 0 π 0 γ −χ 2 5γ from 80 to 10, and remove the condition on χ 2 3γ . The relative change in the e + e − → ωη cross section after loosening the selection condition is taken as a correction for the detection efficiency, while the uncertainty of this correction is added to its systematic uncertainty. The total efficiency correction due to these conditions is −(5.7 ± 6.1)%. The minus sign means that the efficiency in the data is less than that in the simulation. The efficiency correction for the condition N γ = 5 is determined using five-photon events of the process e + e − → ωπ 0 , the cross section for √ s = 1.05-2.00 GeV which can be measured with the condition N γ ≥ 5, and found to be −(0.4 ± 0.2)%. In SND, photons converting into a e + e − pair in the material before the drift chamber, produce a charged track. Such events do not pass the selection criteria. Since the process under study and the process used for normalization contain different numbers of photons in the final state, improper simulation of the photon conversion lead to a shift in the measured cross section. The photon conversion probability is measured using events of the process e + e − → γγ. The corresponding efficiency correction is found to be (−0.79 ± 0.02)%. Thus, the total correction for the detection efficiency is (−6.9 ± 6.1)%.
The detection efficiency for the class rad-ηπγ is calculated in the model (12) with the coefficients α = 0.22±0.21 and β = 0.00 +0.08 . The model uncertainty of the efficiency is determined by varying the coefficients α and β within their errors. It does not exceed 3%.

Fitting the measured cross sections
To calculate the radiative correction and determine the Born cross section, the energy dependence of the visible cross section is fitted by Eq. (13). The uncertainty on the radiative correction is estimated by varying the fitted parameters within their errors. The energy dependence of the Born cross section is parametrized by a sum of contributions of two vector resonances: of the as follows where m V and Γ V are mass and width of the resonance, is the product of the branching fractions of V to e + e − and to the final state f , P f (s) is the phasespace factor.
In the fit to the e + e − → ωη → ηπ 0 γ data, the first term in Eq. (15) is assumed to be the contribution of the ω(1420) resonance, while the second effectively describes the contributions of two resonances, ω(1650) and φ(1680), which have close masses. The phase-space factor is calculated in the approximation of narrow ω meson width: P f (s) = q 3 ω (s), where q ω (s) is the ω meson momentum in the reaction e + e − → ωη. The phase ϕ between the first and second terms is chosen equal to π. The fit is performed in two variants. In the first, the V mass is fixed at the Particle Data Group (PDG) value [20] m V = 1420 MeV, while Γ V , m V , Γ V , B V and B V are free parameters. In the second, we follow the works [8,7] and fix also the ω(1420) width at Γ V = 220 MeV [20]. The fitted curves are shown in Fig. 4. The obtained values of the fit parameters are listed in Table 1. They are in good agreement with the results of the previous works [8,7].
The first variant of the fit has a better χ 2 . It is used to calculate the radiative corrections. The obtained values  Fig. 4. The energy dependence of the e + e − → ωη → ηπ 0 γ Born cross section measured in this work (filled circles). For comparison, the SND [8] (open circles), CMD-3 [7] (squares), and BABAR [6] (triangles) measurements of the e + e − → ωη cross section in the decay mode ω → π + π − π 0 are shown. These data are multiplied by the branching fraction B(ω → π 0 γ) [23]. The curves show the results of the fit described in the text.
of the radiative correction and the Born cross section are listed in Table 2. In the first column of the table we list the boundaries of the energy interval and the weighted average energy, which is calculated as where the sum is taken over the energy points entering in the energy interval, and the visible cross section is calculated using Eqs. (13) and (15) with the parameters obtained in the fit. For the cross section, the statistical and systematic errors are quoted. The main sources of the systematic uncertainty are listed in Table 3. A significant contribution to the systematic uncertainty arises from the interference of the e + e − → ωη amplitude with e + e − → ρη amplitude. The π 0 γ mass distribution for the interference term has a peak in the ω meson mass region. Therefore, in this analysis, the interference will increase or decrease the number of selected events of the process e + e − → ωη. The absolute value of the amplitude of the process e + e − → ρη is extracted from the measured e + e − → π + π − π 0 cross section [16], but the relative phase between the two amplitudes is unknown. We calculated the total cross section for the processes e + e − → ωη and e + e − → ρη with and without interference. The maximum difference between the two cross sections when varying the phase is taken as an estimate of the systematic uncertainty of the e + e − → ωη cross section. In the energy range under study it varies from 5 to 37%. The other sources of the systematic uncertainty are discussed in the previous sections.
To describe the cross section for the class rad-ηπγ, we use a model with one vector resonance decaying into a 0 (1420)γ. The energy dependence of the phase space is calculated using a Monte-Carlo event generator for the process e + e − → a 0 (1420)γ. The free parameters of the fit to the e + e − → rad-ηπ 0 γ visible cross section are m V , Γ V , B V . Their fitted values are listed in Table 1. The values of the radiative corrections determined from the fit and the values of the Born cross section for the process e + e − → Table 1. The parameters obtained in the fits to the e + e − → ωη and e + e − → rad-ηπγ cross sections. BV = B(V → e + e − )B(V → ωη) for ωη, and BV = B(V → e + e − )B(V → ηπγ) for rad-ηπγ. For the ωη final state the results for the two variants of the fit are listed in comparison with the results of the SND [8] and CMD-3 [7] experiments obtained in the ω → π + π − π 0 decay mode. The weighted average energy ( √ s) for the interval indicated in parentheses, integrated luminosity (L), detection efficiency (ε ωη ), number of selected events (N ωη ), radiative correction (1 + δ), and Born cross section (σ ωη ) for the process e + e − → ωη → ηπ 0 γ. For the cross section, the first error is statistical, and the second is systematic.  Table 3. The main sources of the systematic uncertainty on the measured e + e − → ωη and e + e − → rad-ηπ 0 γ cross sections, and the total e + e − → ηπ 0 γ cross section.  Table 4. The weighted average energy ( √ s) for the interval indicated in parentheses, integrated luminosity (L), detection efficiency (ε rad ), number of selected events (N rad ) for the process e + e − → rad-ηπ 0 γ, number of events from the processes e + e − → ρη, φη, φπ (N res ), number of background events (N bkg ), radiative correction (1 + δ), and Born cross section for the process e + e − → rad-ηπ 0 γ. For the cross section, the first error is statistical, the second is systematic.  Total η ω FIT Fig. 6. The measured total cross section for the process e + e − → ηπ 0 γ (circle) and the cross section for the process e + e − → ωη → ηπ 0 γ (squares). The curve is the result of the fit described in the text. rad-ηπγ with statistical and systematic errors are listed in Table 4. The main sources of systematic uncertainty are listed in Table 3. The fitted Born cross section is shown in Fig. 5. It is seen that events of the radiative processes are distributed over a wide energy range, from 1.3 to 1.9 GeV.
The total visible cross section for the process e + e − → ηπ 0 γ is calculated as and then fitted by Eq. (13) with the efficiency obtained for the process e + e − → ωη. The model for the Born cross section used in the fit is a sum of the models for the e + e − → ωη and e + e − → rad-ηπγ cross section. The obtained values of the Born cross section are listed in Table 5 and shown in Fig. 6 in comparison with the cross section for the intermediate state ωη. The main sources of the systematic uncertainty on the cross section are listed in Table 3. The total cross section for the process e + e − → ηπ 0 γ measured by the described method includes the contribution from the interference of the e + e − → ωη and e + e − → ρη amplitudes. Therefore, the model uncertainty associated with the interference is absent.

Summary
In the experiment with the SND detector at the VEPP-2000 collider the cross section for the process e + e − → ηπ 0 γ has been measured for the first time in the energy range from 1.05 to 2.00 GeV. The main contribution to the cross section arises from the intermediate mechanism ωη. The measured cross section for the subprocess e + e − → ωη → ηπ 0 γ agrees well with previous measurements of this cross section by SND and CMD-3 in the decay mode ω → π + π − π 0 . The significantly smaller contribution to the e + e − → ηπ 0 γ cross section from other hadronic intermediate states ρη, φη, φπ 0 , ωπ 0 , and ρπ 0 has been calculated using existing data on their production cross section. It has been found, with a significance of 5.8σ, that the process e + e − → ηπ 0 γ is not completely described by the hadronic intermediate states. We assume that the missing contribution (rad-ηπγ) arises from radiation processes, for example, e + e − → a 0 (980)γ, a 0 (1450)γ, and a 2 (1270)γ. The cross section for the process e + e − → rad-ηπ 0 γ has been measured. It is 15-20 pb in a wide energy range, from 1.3 to 1.9 GeV. The spectrum of ηπ 0 invariant masses for the events rad-ηπγ is consistent with the dominance of the intermediate mechanism a 0 (1450)γ. √ s = 1.05-2.00 GeV Table 5. Weighted average energy ( √ s) for the interval indicated in parentheses, and the total e + e − → ηπ 0 γ cross section (σ ηπγ ). For the cross section, the first error is statistical, the second is systematic.