Partial-wave analysis of proton Compton scattering data below the pion-production threshold

Low-energy Compton scattering off the proton is used for determination of the proton polarizabilities. However, the present empirical determinations rely heavily on the theoretical description(s) of the experimental cross sections in terms of polarizabilities. The most recent determinations are based on either the fixed-$t$ dispersion relations (DR) or chiral perturbation theory in the single-baryon sector ($\chi$PT). The two approaches obtain rather different results for proton polarizabilities, most notably for $\beta_{M1}$ (magnetic dipole polarizability). We attempt to resolve this discrepancy by performing a partial-wave analysis of the world data on proton Compton scattering below threshold. We find a large sensitivity of the extraction to a few"outliers", leading us to conclude that the difference between DR and $\chi$PT extraction is a problem of the experimental database rather than of"model-dependence". We have specific suggestions for new experiments needed for an efficient improvement of the database. With the present database, the difference of proton scalar polarizabilities is constrainedto a rather broad interval: $\alpha_{E1}-\beta_{M1} = (6.8\, \ldots \, 10.9)\times 10^{-4}$ fm$^3$, with their sum fixed much more precisely [to $14.0(2)$] by the Baldin sum rule.


I. INTRODUCTION
The low-energy Compton scattering (CS) off the proton and light nuclei is the standard tool for probing the nucleon polarizabilities, see [1-5] for reviews. However, the relation between the experimental observables and polarizabilities is simple only when neglecting the higher-order terms in the low-energy expansion (LEX) of Compton amplitudes. In practice, the higher-order terms play an important role, and, for a quantitative extraction of polarizabilities from Compton scattering data, more sophisticated theoretical frameworks are being used. In the case of the proton, there are two types of "state-of-the-art" extractions: (i) based on the fixed-t dispersion relations (DR) [1, 2, 6-9] and (ii) based on chiral perturbation theory (χPT) with explicit nucleons and Delta's. The latter calculations can be divided into two types: heavy-baryon (HBχPT) [10][11][12] or manifestlycovariant (BχPT) [13,14]. The problem is that, although both DR and χPT give comparably good description of the experimental cross sections, the extracted values of polarizabilities differ, sometimes by a few standard deviations.
A notable example is provided by the magnetic dipole polarizability β M1 of the proton, which ranges from 1. 6(4) [in units of 10 −4 fm 3 , omitted in what follows] obtained in the state-of-art DR fits of the data [2,9,15] to 3.2 (5) in the χPT fits [3,12,16]. Furthermore, without using the Compton data, the BχPT at NNLO yields for the proton [13,14]: β M1 = 3.9 (7), making the discrepancy with DR more acute. Incidentally, the current PDG average [17] is basically combining the DR and HBχPT values, resulting in 2.5(4) for the proton β M1 . Their central value may serve as a compromise, but the uncertainty does not seem to reflect the spread between the DR and χPT values.
The present work is an attempt to resolve this tension in a model-independent manner, by making the partial-wave analysis (PWA) of the CS data below the pion photoproduction threshold ( 150 MeV). To our surprise, we find that the above discrepancy between DR and χPT fits is a problem of the experimental database, rather than of theoretical descriptions. As such, it calls for new experiments. We argue that new precision data for the proton CS angular distribution at backward angles and beam energy around 100 MeV are highly desirable.
The PWA is of course a good old method to study the hadronic processes at low energy. Yet, it has barely been used in proton CS (γ p → γ p). To date, the only comprehensive PWA of proton CS data remains to be the 1974 analysis of Pfeil et al. [18], in the ∆(1232) region. The region below the pion threshold has not been analyzed until now -the present PWA is the first.
Of course, this is not a first study of the Compton multipoles below threshold in general, cf. Refs. [14,[19][20][21] for calculations using DR and χPT frameworks. However, a modelindependent Compton PWA has not done until now, mainly due to the lack of accurate data. The latter problem is compensated in our analysis by the recent empirical determination of the forward Compton amplitudes through the sum rules involving the photoabsorption cross sections [22,23]. The sum rules yield two independent linear relations between the multipole amplitudes at each energy. Note that the linear relations among the partial-wave amplitudes are very rare in PWAs. They usually operate with the bi-linear relations between the amplitudes and experimental observables alone. Furthermore, we make use of the fact that the partial-wave (or multipole) amplitudes are assumed real below the threshold, that is if one neglects the radiative corrections. It is also important to build in the correct low-energy limit and treat exactly the Born contributions.
In the following we present the low-energy parametrization of the pertinent multipole amplitudes in terms of the static polarizabilities (Sec. II), and the corresponding fits of the experimental data for proton CS (Sec. III). The results are discussed in Sec. IV, and conclusions are given in Sec. V.

II. MULTIPOLE EXPANSION AND (BI-)LINEAR EMPIRICAL CONSTRAINTS
The general formalism of the multipole expansion for nucleon CS is given in [18,24], and concisely summarized in [5,Sec. 2]. The idea is that, using the rotational and discrete symmetries, the CS helicity amplitudes T σ λ ,σ λ (s,t), with σ (σ ) the helicity of initial (final) photon and λ (λ ) for the helicity arXiv:1712.05349v2 [nucl-th] 23 Apr 2018 initial (final) nucleon, admits a partial-wave expansion: with J the total angular momentum of the photon-proton system, T J σ λ ,σ λ (ω) the partial-wave amplitudes, d J λ , λ (θ ) the Wigner d-function, ω and θ the photon energy and scattering angle in the center-of-mass frame; s, t, u are Mandelstam invariants.
The partial-wave amplitudes T J (ω) are then linearly related to the amplitudes with definite parity and angular momentum l, i.e., multipoles f l± ρρ (ω), with ρ, ρ = E(lectric), or M(agnetic). The infinite sum over half-integer J is then replaced by the sum over integer l = J ∓ 1/2. Note that f 0+ ρρ = 0, by definition; hence the summation starts at l = 1.
In this work, we first write the amplitude as the sum of the Born, T Born , and the rest (non-Born)T , as illustrated in Fig. 1 (note that here the π 0 -pole contribution is part of the Born term). The same decomposition holds for the multipoles: f = f Born +f . We then truncate the multipole expansion of the non-Born amplitude at J = 3/2, whereas the Born amplitude is treated exactly. We thus retain the ten lowest non-Born multipoles, the rest are neglected. This approximation is well justified at energies below the pion production threshold (ω m π ), as the leading low-energy behavior of the non-Born multipoles is [25] Furthermore, the existing χPT calculations [13,14,19] show that the four l = 2 non-Born multipoles,f 2+ give tiny contributions below the pion threshold. In what follows we will either neglect them, or fix them to the values given by the latest BχPT calculation [14]. We shall therefore fit only the six l = 1 non-Born multipoles.
In order to build in the low-energy behavior of the non-Born multipoles [cf. Eq. (3)], we assume the following parametrization of the l = 1 multipoles in terms of static polarizabilities: where we changed the photon energy from the center-of-mass ω to the lab frame E γ = ω √ s/M. The first term in each of the square brackets of Eq. (4) is given by one of the six static polarizabilities, four of which, denoted by γ's, are spindependent. The 2nd terms are the recoil corrections (see, e.g., Ref. [14]). The 3rd terms are given by the "residual functions" f R i . The parametrization of Eq. (4) ensures the correct low-energy behavior of these multipoles. It does not imply any approximation: the six static polarizabilities as well as the residual functions are free parameters, which will next be determined from experimental data.

A. Bilinear relations: observables
Any CS observable provides bi-linear relations on CS multipoles. This is of course the usual situation for any PWA, an experimental observable, such as cross section or asymmetry, involves the amplitude squared.
Take for instance the unpolarized angular distribution dσ /dΩ, given in terms of the helicity amplitudes by dσ dΩ = 1 Substituting in here the multipole expansion of T we obtain (for J < 5/2): where c n are bilinear combinations of the multipoles. In principle, each c n can be extracted from the fit to the data, and hence one obtains 5 bilinear relations from this observable. Similarly, for the beam asymmetry, defined as where σ || and σ ⊥ are the CS cross sections with linear photonbeam polarization (parallel and perpendicular to the scattering plane), we have: thus providing 3 more bilinear relations (generally c n and d n are different).
The bilinear relations provide a system of quadratic equations for the multipoles. In reality, the angular coverage and precision of the data do not allow for unique solution of these equations, at least not at the present time. Fortunately, as discussed in what follows, the sum rules for the forward CS provide accurate linear relations, which simplify things a lot.

B. Linear relations: Sum rules
The general properties of forward CS, derived from unitarity, causality and crossing [26], allow for it to be expressed entirely in terms of integrals of total photoabsorption cross sections. In case of a spin-1/2 target such as the proton, the forward CS is characterized by two scalar amplitudes, f (ν) and g(ν), functions of the invariant ν = (s − u)/4M, which in the forward kinematics is identical to the photon lab energy E γ . The helicity amplitudes are written in terms of the scalar amplitudes as follows: where ε σ and χ λ are the photon polarization vector and the nucleon spinor, with the subscripts showing the corresponding helicities. These forward amplitudes are given by the sum rules on one hand and by the multipole expansion on the other: where σ abs Λ is the photoabsorption cross section corresponding to the total helicity Λ of the initial γ p state, 0 + is an infinitesimal positive number. The summation of the multipoles is done over L = J − 1 /2 = (l ∓ 1 /2) − 1 /2. Note also that the first term in the sum-rule expressions (due to the proton charge and anomalous magnetic moment κ) are precisely the Born contributions, whereas the integrals yield the non-Born contributions.
The empirical evaluation of the forward amplitudes f (ν) and g(ν) for proton CS has recently been performed in [22] and [23], respectively. Thus, the sum rules provide two linear relations on the multipole amplitudes. We use these relations to eliminate the residual functions f R 2 (E γ ) and f R 3 (E γ ) in Eq. (4).
The low-energy expansion of the integrals in Eq. (10a) yield sum rules for the forward combinations of static polarizabilities, such as the Baldin sum rule [27] (see, e.g., [5,Sec. 5] for more details). Based on the empirical evaluation of these sum rules, we use [22,23]: to eliminate two out of six global parameters in Eq. (4), our choice being α E1 + β M1 (so that α E1 − β M1 is fitted) and γ M1M1 . We thus are left with four global parameters and four energy-dependent functions.

III. CS DATABASE AND FITTING STRATEGY
The world database on the unpolarized angular distribution of proton CS, below the pion-production threshold, is summarized in Table I, cf. [3,28]. The number of data points con-tributed by each experiment is indicated in the column N data . The database is split into N bins = 11 energy bins, with the central values at 1 59, 69, 79, 89, 99, 109, 117, 127, 135, 143, 150 MeV. (12) We fit all these data simultaneously, hence the number of parameters is 4 + 4 N bins = 48. This is quite a large number, and we perform the fitting in two stages: 1) a Monte-Carlo swipe fitting both the static polarizabilities and the residual functions, by finding the least χ 2 for a large ensemble of parameters taken from the normal (Gaussian) distribution 2 ; 2) the χ 2 is further minimized by varying the static polarizabilities using a standard minimization routine, whilst the residual functions are kept fixed to the values determined in the Monte-Carlo swipe.
Our fit to the database of Table I results in Fit 0 of Table II. The results of Fit 1 correspond to the fit where the small (according to many existing analyses, cf. the last 3 rows of Table II) spin polarizability γ E1M2 is set to zero. The results of the two fits are consistent with each other, albeit Fit 1 provides a much better accuracy. We take it as a sign of insufficient data quality for the accurate determination of the small value of γ E1M2 , and keep γ E1M2 = 0 in our subsequent fits. The three subsequent fits in Table II are done upon various "refinements" of the database involving deletion of "outliers". Namely, 1 We have tried to optimize the number of energy-bins to minimize the number of fitting parameters. Thus, we omitted the data from the very lowenergy region (below 50 MeV), such as those of Federspiel et al. [35], which would have relatively low number of points per bin. We have also omitted two data points from the same source taken at 65.8 MeV, in order to avoid having a separate bin. 2 In the interest of full disclosure, the parameters of the normal distribution used in all of our fits are as follows. The residual functions are centered at zero with the width for f R 1,4 (E γ ) given by 10 × 10 −4 fm 3 , and the width for f R 5,6 (E γ ) given by 10 −4 fm 4 . These choices for the widths are motivated by the "natural size" argument based on the known values of the static polarizabilities, even though, the use of normal distribution allows for these functions to take any values in principle. The parameters of the normal distributions for α E1 − β M1 , γ E1E1 , γ M1E2 , and γ E1M2 are taken as  • in Fit 1 3σ , the outliers are identified according to the simple 3σ rule [41,42], i.e., as those that deviate more than 3σ from Fit 1, see Fig. 2.
• In Fit 1 , the 4 deleted outliers are: Ref. [ Hence, two of the deleted points are the same as in the previous fit and two are different. The latter two point are selected by hand such as to drive the fit closer to the BχPT-predicted cross section.
• In Fit 1 , we purge the database in accordance to what is done in χPT fits as described in [3], i.e.: omit the data of Oxley et al. [29] entirely, Bernardini et al. [32] entirely, Baranov et al. [34] at θ lab = 150 • , and Olmos de León et al. [15] at {109 MeV, 133 • }. Furthermore, as in [3], we add 5% systematic uncertainty (point-to-point in quadrature with the statistical error) to the points of Ref. [15]. Unlike [3], we do not include the floating normalization factors. Also, the points of Federspiel et al. [35] are treated as described in the footnote (i.e., omitting them below 50 and at 65.8 MeV), rather than  The proton scalar and spin polarizabilities in units 10 −4 fm 3 (scalar) and 10 −4 fm 4 (spin), obtained in the various fits described in the text, compared with the BχPT predictions [14], DR calculations [15,39] (note that only α E1 + β M1 is calculated in DR, with their difference fitted to CS data), and an experimental extraction of spin polarizabilities at MAMI [40] (performed using subtracted DRs [9]). We do not include the data on beam asymmetry in our fits, since the only data (below pion-production threshold), coming from the pilot experiment at MAMI [43], are of significantly poorer quality compared to the unpolarized data. Hopefully, the currently running A2/MAMI experiment will improve the accuracy for this observable, and thus play a crucial role in an accurate determination the magnetic polarizability β M1 [44]. At present we only verify that all our fits are in agreement with the pilot data [43], see Fig. 3. Table II presents the static polarizability values resulting from the 5 fits described in the previous section. The last column shows χ 2 /point, a measure of the fit quality. The obtained polarizabilities can be compared with the last 3 rows showing respectively the BχPT prediction, DR extraction, and the first experimental extraction of the spin polarizabilities (MAMI 2015).

IV. RESULTS AND DISCUSSION
The striking result here is that the polarizability values are fairly sensitive to the slight refinements of the database. For example, for β M1 we obtain the values ranging from 1.8(3) using the original database in Fit 1 to 3.8(4) using an improved one in Fit 1 . The latter modification of the database is similar to the one used in the χPT fits of McGovern et al. [3,12,16], which could explain why the χPT fits are significantly different from the DR fits, which in particular yield a low value of β M1 .
Let us emphasize that the BχPT row in the Table is not an extraction from CS data, but is rather a prediction, albeit of a low order [13,14]. Nonetheless, if we are to take the claimed uncertainties seriously, we must conclude that the refined databases agree somewhat better with χPT.
Besides the static polarizabilities, our fits yield the multipole amplitudes at the considered energy bins. The non-Born multipoles can equivalently be represented by the so-called dynamical polarizabilities (see, e.g., [5,Sec. 2] for definition). The blue (red) points in Fig. 4 show the dynamical polarizabilities resulting from Fit 1 (1 ). Note that the point at zero energy corresponds with the static polarizability. The error bars result from the uncertainties on the fit parameters. The results are compared with the BχPT (cyan bands) and DR (dashed lines) results. Again, we see that our solution based on the raw database (Fit 1) agrees well with DR calculation, whereas the one based on the refined database (Fit 1 ) agrees better with BχPT. Therefore, the differences between the χPT and DR results for polarizabilities are likely to be caused by deficiencies in the experimental database. How to resolve those? We first of all need to find the place where the differences among the different fits are largest. For the unpolarized cross section, the "sweet spot" is apparently at E γ 109 MeV and backward angles, see Fig. 5. At both higher and lower energies the difference among the fits quickly diminishes, cf. Fig. 6. Hence the best hope for resolving this "database consistency problem" is to obtain new precise cross-section data at energies close to 109 MeV.
Let us consider this energy region in more detail. In Fig. 5, besides the data and the results of 3 fits, we show the Born contribution (dash-dotted curve) and the BχPT prediction [13,14] (dotted curve). The deviation from the Born contribution is the effect of (dynamical) polarizabilities we are after. The polarizability contribution is at low-energy dominated by the scalar dipole polarizabilities, α E1 and β M1 , but already at 109 MeV the spin polarizabilities start to play a crucial role.
To see this, consider Table III, where the forward and backward combinations of scalar and spin polarizabilities are presented. In the fits the forward combinations are fixed by the sum rules, Eq. (11), whereas the backward combinations, α E1 − β M1 and γ π = −γ E1E1 + γ M1M1 − γ E1M2 + γ M1E2 are different from fit to fit. Fit 1 has the highest value of α E1 − β M1 and hence has the biggest deviation from the Born term at 59 MeV; the γ π value is not important at these energies. Fits 1 and 1 have α E1 − β M1 close to BχPT and as the result the three curves practically coincide at 59 MeV.
However, at 109 MeV, Fit 1 converges to Fit 1 precisely because of the different γ π value. The similar effect for Fit 1 is diminished by the difference in the value of α E1 − β M1 . Thus, at these energies the scalar and spin polarizabilities are rather entangled and cannot be extracted independently from this observable. The present PWA, on the other hand, provides a basis for a simultaneous extractions of α E1 − β M1 and the backward spin polarizability γ π .  [14] and DR [15,39] calculations.

V. SUMMARY AND CONCLUSION
We presented a first partial-wave analysis of proton Compton scattering data below the pion-production threshold (E γ 150 MeV). The only approximations, or model-dependent assumptions, we made concern the truncation of the partialwave expansion: • we account for the lowest l = 1 and 2 terms, neglecting J ≥ 5/2 contributions; • for the l = 2 multipoles, we assume the values given by the NNLO BχPT calculation [13,14], and check that the results do not change qualitatively if we put them to 0 (cf. Fit 2 variety in Table II).
The proper low-energy behavior of the (non-Born piece of) multipoles is ensured through the parameterization in terms of lowest static polarizabilities, see Eq. (4). The sum rules for the forward amplitudes impose two linear relations on the multipoles, leaving us with only four of the six amplitudes to be determined from the Compton angular-distribution data. The accuracy of the resulting solutions is significantly improved by setting (the small spin polarizability) γ E1M2 = 0 by hand. The extracted multipoles depend significantly upon very mild refinements of the world database of proton Compton scattering. The characteristic difference between the stateof-art DR and χPT analyses is likely to be explained by the database inconsistencies, rather than by differences in the theoretical framework. We claim that these inconsistencies are best to be addressed by a new precise measurement of the angular distribution at E γ ≈ 109 MeV and backward angles (cf. Fig. 5). Accurate data on polarized observables, such as the beam asymmetry, could be helpful too. The ongoing Compton scattering experiment by the A2 Collaboration at MAMI may soon provide a considerable improvement of the database, in both the angular distribution and beam asymmetry. Until then, the static polarizabilities may continue to be extracted in a rather wide range of values, manifested by our fit results in Table II.