First concurrent extraction of the leading-order scalar and spin proton polarizabilities

We performed the first simultaneous extraction of the six leading-order proton polarizabilities. We reached this milestone thanks to both new high-quality experimental data and an innovative bootstrap-base fitting method. These new results provide a self-consistent and fundamental benchmark for all future theoretical and experimental polarizability estimates.


INTRODUCTION
Understanding the hadron structure in the nonperturbative regime of quantum chromodynamics (QCD) is one of the major challenge of modern physics. We can classify hadrons in terms of their global properties, such as mass and spin, but we can not fully explain how these properties emerge from the underlying dynamics of the hadron's interior. A clean probe to investigate the internal structure of hadrons is the Compton scattering process that gives access to observables with a clear interpretation in terms of structure-dependent properties of the hadrons. In particular, real Compton scattering (RCS) at low energies is parametrized by polarizabilities that describe the response of the charge and magnetization distributions inside the nucleon to an applied quasistatic electromagnetic field. These structure constants are fundamental properties of the nucleon and their determination has driven a relevant experimental effort in the last few years [1][2][3][4][5].
The effective multipole interactions for the coupling of the electric ( E) and magnetic ( H) fields of the photon with the internal structure of the nucleon is described at leading order in terms of the electric (α E1 ) and magnetic (β M 1 ) scalar polarizabilities [6,7]: while the four spin polarizabilities (γ ... ) show up in the subleading terms: where σ are the proton's Pauli spin matrices,˙ E = ∂ t E and E ij = 1 2 (∇ i E j + ∇ j E i ) are partial derivatives with respect to time and space, respectively.
In addition to being fundamental properties of the nucleon, polarizabilities play a profound role in precision atomic physics, in the evaluation of the nuclear corrections to atomic energy levels [8][9][10][11][12] and in astrophysics, influencing neutron star properties [13].
Despite their evident importance in a broad range of physics topics, up to now a self-consistent experimental extraction of all the different polarizability values has not been possible, due to the scarce quality of the available database (see, for instance, Ref. [14]). In all existing fits of the RCS data, some of the polarizabilities have been fixed either using theoretical calculations [15][16][17][18] or empirical evaluations from other reactions [19], or, at most, have been constrained to vary within certain intervals [14]. The situation has recently improved with the first measurements of the double-polarization observables Σ 2x [1] and Σ 2z [3] and new data for the unpolarized differential cross sections and the single-polarized Σ 3 asymmetry [2,4,5]. The beam asymmetry is defined as [6]: where dσ (⊥) is the polarized cross section obtained with a photon beam polarized parallel (or perpendicularly) to the scattering plane and an unpolarized target. In a similar way, the double-polarization asymmetries can be defined as: where dσ R(L) ±x is the polarized cross section obtained with circular right-handed (left-handed) photon polarization and target spin aligned transversely (±x) with respect to the incident beam direction, while dσ R(L) ±z is obtained with the target spin aligned longitudinally (±z) with respect to the incident beam direction.
In particular, the work from Mornacchi et al. [4] provides the highest statistics proton RCS single dataset ever obtained below the pion photoproduction threshold, with 60 unpolarized differential cross section points and 36 beam asymmetry points measured over a large angular range and with small statistical and systematic errors. Therefore, it represents a significant improvement for a more accurate extractions of all the different polarizability values.
In this paper we present the first consistent and simultaneous fit of the six leading-order static proton polariz-arXiv:2204.13491v2 [hep-ph] 15 Oct 2022 abilities. It has been obtained thanks to both the new experimental data and an innovative bootstrap-base fitting method [20]. This algorithm has already been deployed successfully for the extraction of the proton scalar dynamical and static polarizabilities from low-energy RCS data [14,21]. The theoretical framework used for this extraction is based on fixed-t subtracted dispersion relations (DRs) [7,22,23]. The theoretical uncertainties associated to the model dependence of our results are also evaluated by using, as input, pion photoproduction amplitudes obtained from three different partial wave analyses (PWAs) of the available experimental data: BnGa-2019 [24], MAID-2021 [25,26], and SAID-MA19 [27]. This is the very first time that such a comprehensive and self-consistent study on the simultaneous extraction of all the six leading-order proton polarizabilities from RCS data is performed.

DATABASE SELECTION AND FIT PROCEDURE
The proton RCS database used for this work consists of two main sets: the unpolarized differential cross section data, and the (single and double) polarization asymmetries. The former can be further divided into low-and high-energy data, namely data for which the incoming photon energy E γ is below or above the pion photoproduction threshold (∼ 150 MeV), respectively.
For the low-energy set, in addition to the new data from Refs. [4,5], we used the same selection extensively discussed in a previous work by Pasquini et al. [14], which includes datasets from Refs. [28][29][30][31][32][33][34][35][36][37][38][39]. For the high-energy set, thanks to the DR model used for the theoretical framework, we were able to consider data measured up to E γ = 300 MeV (corresponding to the 2π photoproduction threshold), thus extending the range used in the fits of Refs. [14,16]. Furthermore, only for these high-energy data, we decided to narrow the selection to the newgeneration experiments, namely the measurements performed using tagged photon facilities (see, for instance, the review of Ref. [40]). Their main advantage, in addition to a more reliable photon flux determination, is that the incoming photon energy is known with a resolution of a few MeV. This is an essential ingredient to reject the overwhelming background coming from the single π 0 photoproduction channel that has a cross section by two orders of magnitude higher than Compton scattering process, and can mimic the Compton signature when one of the two photons coming from the π 0 decay escapes the particle detection. The available datasets for the unpolarized cross section come from two different facilities: MAMI [41][42][43][44][45][46] (with also few data above threshold from Ref. [39]) and LEGS [47,48].
The polarization observables have enhanced sensitivity to the spin polarizabilities [23,49], hence are crucial for the extraction of these structure constants. The adopted datasets include three different polarization observables: [3], and Σ 3 both below [2,4,5] and above [48] pion photoproduction threshold.
The fit to extract the leading-order scalar and spin polarizabilities was performed using a bootstrap-based method [20] that consists of randomly generating N Monte Carlo replicas of the fitted experimental database, where each data point e The indices i, j, and b run over the number of data points in each dataset, the number of datasets, and the bootstrap replica, respectively; r i,j,b is a random number extracted from the normal distribution N (0, 1), and δ j,b is a random variable that accounts for the effect of the common systematic errors, independently for each dataset. From each of these simulated databases, a set of fitted parameters is extracted. The mean and standard deviation of the obtained distributions give then the best value and the error for each of the fitted parameters. This technique offers several advantages compared to other fitting procedures, especially when different datasets are used together, as in the present work: i) a straightforward inclusion of common systematic uncertainties without any a priori assumption on their distributions and without introducing any additional fit parameter; ii) the probability distribution of the fit parameters (often non-gaussian) is obtained directly by the procedure itself; iii) the uncertainties on possible nuisance model parameters are easily and directly taken into account in the sampling procedure; iv) the correct fit p-value is always provided when systematic uncertainties are present and in all the other cases when the goodness-of-fit distribution is not given by the χ 2 -distribution.
In the first step of our analysis, we checked the consistency of the selected database by looking at the distribution of the normalized residual for each of the largest datasets (e.g., with more than 40 data points): the unpolarized cross section from the A2 [4], the TAPS [39], and the LARA [44] Collaborations; and the unpolarised cross section and beam asymmetry from the LEGS [48] Collaboration.
Since the cross section data from LARA and LEGS are known to be in significant disagreement between each other (see, for instance, Ref. [15]), we performed a preliminary test by alternatively including LARA or LEGS data in the fit database and simultaneously fitting all six polarizabilities using the MAID-2021 [25,26] multipole solution.
For each of the two configurations, we took the output polarizability best-values, calculated the residual distribution for each dataset, and produced a probability plot [50] for assessing the residual normal distribution. In all cases, the residuals were found to follow fairly well the expected normal distribution for all the selected datasets (see Fig.1 of the Supplemental Material [51], where the probability plots refer to the test without the LEGS data), except for the unpolarized cross section from both the LEGS and LARA Collaborations, as shown in Fig. 1.
We repeated the same test using the SAID-MA19 [27] and BnGa-2019 [24] PWAs and obtained very similar results. Given also the small fit p-values, i.e., 5·10 −4 and 1 · 10 −2 for the LEGS and LARA data, respectively, with all the PWA inputs, we excluded both sets from the database of the present work.
Although the data for the photon asymmetry Σ 3 above threshold are from the same dataset of the LEGS unpolarized cross section, they give a consistent residual plot (see Fig. 1 of the Supplemental Material [51]). Such an agreement is not surprising since, as noted, e.g., in Ref. [15], most of the possible systematic biases cancel out in an asymmetry observable.
Taking all these considerations together, 25 datasets for a total of 388 points were included in the fit. The angular and energy coverage of all the used datasets are listed in the Supplemental Material [51].

RESULTS AND DISCUSSION
A total of N = 10 4 bootstrapped samples of the database was generated and the minimization was performed at the end of each iteration with all six proton polarizabilities treated as free parameters. For convenience, we used as actual fit parameters some linear combina-tions of the scalar and spin polarizabilities: The last term, γ π , is the sum of the dispersive contribution γ disp π , to be fitted to the data, and the pion-pole contribution, fixed to γ π 0 −pole π = −46.7 [40]. The choice of the fit parameters allows for a direct comparison, as a consistency check, of the fit results for α E1 + β M 1 and γ 0 to the available experimental predictions of the Baldin [39,[55][56][57] and Gellmann-Goldberger-Thirring (GGT) [55][56][57][58] sum rule values, respectively, obtained using data for the total photoabsorption cross section. When present, point-topoint systematic errors were added in quadrature to the statistical errors, while common systematic scale factors were treated as in the previous bootstrap extractions by Pasquini et al. [14,21], namely they are assumed to follow a uniform distribution, unless otherwise specified in the original publication. Moreover, when multiple systematic sources are given, the final error is the product of all the generated random uniform variables.
We performed the minimizations by using the nonlinear least-squares fitting routines of the GSL library [59]. As a consistency check of this procedure, both the gradient and the simplex methods were used as best-fit algorithms and identical results were obtained. The entire procedure was performed for the three PWAs. The obtained distributions for each parameter are reported in the Supplemental Material [51]. The parameter distributions obtained using the three different PWA inputs have the same shape and differ only for a shift in their central values. For this reason, we evaluated the central polarizability values as the mathematical average of the three different sets of fit values. Additionally, the largest of the differences between each set of fit values and the average was used to estimate an additional model error (conservatively considered as a standard deviation) due to dependence on the PWA used as input in the DRs. The resulting best-fit values are: The quoted fit errors are the 68% confidence level (CL) and include the contribution of both the statistical and systematic uncertainties of the experimental data. The additional model-dependent systematic uncertainties, evaluated as explained above, are given in rms units, and have to be added in quadrature to the previous ones to get the overall values of the estimated systematic uncertainties. All different fits gave, within rounding errors, a mini-  [14], HBChPT [16], and BChPT [17] refer to three extractions within DR [7,22,23], HBChPT [18], and BChPT [52] frameworks, respectively. In all cases, most of the spin polarizabilities were fixed to a given value. For the four spin polarizabilities in panels c) to f), some of the existing extractions and theoretical predictions are shown below and above the dotted line, respectively. In particular, A2 refers to the most recent experimental extraction [3], where both scalar polarizabilities were fixed to the PDG values [53]. DR [7], HBChPT [16,18], and BChPT [54] above the dotted line refer to theoretical predictions using different approaches.
mum value of the fit function equal toχ 2 red = 1.13. As mentioned before, the expected goodness-of-fit distribution is not given by the χ 2 function, because of the correlations between points of a same dataset introduced by the systematic uncertainties. This density was then evaluated in the framework of the bootstrap technique [20] and a p-value = 24% was estimated from its cumulative distribution. A plot of this function and the fit correlation matrix are reported in Fig. 3 and Table 3 of the Supplemental Material [51], respectively.
The obtained value of both α E1 + β M 1 and γ 0 are in agreement, within the quoted errors, with the available estimates of the Baldin and GGT sum rule values listed in Refs. [56,57].
The final residual distribution, the average χ 2 -perpoint for each subset, and the comparison between the results of the fit with a sample of experimental data from different observables, both below and above pion threshold, are collected in Figs. 3-5 of the Supplemental Material [51]. Taken together, all these results confirm the validity of our overall fit procedure.
From the fit results reported in Eq. (6), we obtained the 68% CL intervals of the six proton static polarizabil- where the meaning of the fit and model errors is the same as in Eq. (6).
The proton polarizability values reported in Eq. (7) are shown as blue points in Fig. 2 (bottom row in each panel). The blue horizontal bars represent the errors given by the fit procedure. The increase of the overall systematic uncertainties due to the inclusion of the model errors is shown in red. In Fig. 2a and Fig. 2b, the vertical gray bands give the average 68% CL interval on α E1 and β M 1 as evaluated by the Particle Data Group (PDG) [53]. In the same figures, the new results are compared with some of the existing global extractions of α E1 and β M 1 using DRs [14], HBChPT [16], and BChPT [17], respectively, where at least three spin polarizabilities were kept fixed. In Fig. 2c, Fig. 2d, Fig. 2f, and Fig. 2e our new results are compared with the last experimental extraction [3], where the two scalar polarizabilities were fixed, and three theoretical calculations within DRs [7], HBChPT [16,18], and ChPT [17]. The high quality, as well as the importance, of these new results is highlighted in all the plots. In fact, they provide a self-consistent extraction of the six leading-order static proton polarizabilities without any fitting approximation, and with errors that are competitive with those of all the existing evaluations, that were all performed by constraining some of the polarizabilities to reduce the uncertainty on the ones of interest. In particular, the new results from the A2 Collaboration on the unpolarized cross section [4] were fundamental in reducing the correlations between α E1 − β M 1 and γ π , and between α E1 − β M 1 and γ M 1M 1 .
These new results then are a significant benchmark for all future theoretical and experimental polarizability estimates. However, the fit parameters have still relevant uncertainties and, as it can be seen from the values reported in the Supplemental Material [51], there is still a slightly high correlation between γ M 1M 1 and γ π (ρ (γ M 1M 1 −γπ) 0.80). These are clear indications of the need for new dedicated measurements, that will be discussed in detail in a dedicated forthcoming publication.
In summary, we presented the first simultaneous and self-consistent extraction of the six leading-order static proton polarizabilities. The fit was performed using a bootstrap-based technique combined with a fixed-t subtracted Dispersion Relation model for the theoretical calculation, using three different PWA solutions as input.
The obtained values have an error that is competitive with the existing extractions, that were all obtained with the inclusion of constraints on some of the fit parameters. These results provide new important information to our understanding of the internal electromagnetic proton structure, and should be used as input for further experimental and theoretical extractions.
We are grateful to Igor Strakovsky and Viktor Kashevarov, who provided us with the SAID-