Quasielastic axial-vector mass from experiments on neutrino-nucleus scattering

We analyze available experimental data on the total and differential charged-current cross sections for quasielastic neutrino and antineutrino scattering off nucleons, measured with a variety of nuclear targets in the accelerator experiments at ANL, BNL, FNAL, CERN, and IHEP, dating from the end of sixties to the present day. The data are used to adjust the poorly known value of the axial-vector mass of the nucleon.


Introduction
A precise knowledge of the cross sections for chargedcurrent induced quasielastic scattering (QES) of neutrinos and antineutrinos on nuclear targets is a pressing demand of the current and planning next generation experiments with accelerator and atmospheric neutrino beams, aiming at the further exploration of neutrino oscillations, probing nonstandard neutrino interactions, searches for proton decay, and related phenomena.
The quasielastic cross sections are very sensitive to the poorly known shape of the weak axial-vector form factor F A (Q 2 ) of the nucleon. Adopting the conventional dipole approximation, this form factor is determined by the axialvector coupling g A = F A (0) and the phenomenological parameter M A , the so-called axial-vector (dipole) mass related to the root-mean-square axial radius by The experimental values of M A extracted from neutrino and antineutrino scattering data and from the more involved and vastly model-dependent analyses of charged pion electroproduction off protons, show very wide spread, from roughly 0.7 to 1.2 GeV with the formal weighted averages [1,2] M A = 1.026 ± 0.021 GeV from ν µ , ν µ experiments, 1.069 ± 0.016 GeV from π electroproduction.
In the likelihood analysis, we use the most accurate phenomenological parametrizations for the vector form factors of the nucleon [72,73], we take into account all known sources of uncertainties, in particular, the systematic errors in the energy spectra of ν µ and ν µ beams. For description of nuclear effects we apply the standard RFG model. We examine possible difference between the values of M A extracted from ν µ and ν µ data, and cross-check our results with the data on Q 2 distributions measured in several experiments.

.1 Structure functions and cross section
Let us first summarize the well-known phenomenology for describing the hypercharge conserved quasielastic reactions on free nucleon targets ν ℓ (k) + n(p) → ℓ − (k ′ ) + p(p ′ ), ν ℓ (k) + p(p) → ℓ + (k ′ ) + n(p ′ ). (1) Here k, k ′ , p, and p ′ denote the four-momenta and ℓ stays for e, µ, or τ . In this paper, we will neglect the protonneutron mass difference, 2 since the resulting correction, in the ν µ /ν µ case, exclusively works near the reaction threshold and practically negligible for the energies of our current interest. The general formulas which take this effect into account, were derived in Ref. [74] (assuming T and C invariance) and in Refs. [75,76] (avoiding these assumptions). The double differential cross-section for these processes is a convolution of spin-averaged leptonic and hadronic tensors L αβ and W αβ : 1 The νe, νe, ντ , and ντ beams from past and current accelerator experiments are not appropriate for measuring the QES cross sections. 2 While our computer code operates with the most general formulas and relevant kinematics.
Here G F is the Fermi coupling, q = k − k ′ is the fourmomentum transferred from the incoming (anti)neutrino to the nucleon, Q 2 = −q 2 , M W is the mass of intermediate W -boson; E ν , E ℓ , P ℓ = E 2 ℓ − m 2 ℓ , and θ ℓ are, respectively, the incident (anti)neutrino energy, outgoing lepton energy, momentum, and scattering angle in the lab frame, m ℓ is the lepton mass. The leptonic tensor defined by the product of the weak leptonic currents, is given by where the upper (lower) sign is for ν ℓ (ν ℓ ). Assuming the isotopic invariance, the hadronic tensor is defined by the six structure functions W i (Q 2 ): where M is the mass of the "isoscalar" nucleon. Then combining Eqs. (3) and (4) yields where κ = m ℓ /2M .
In order to connect the structure functions with the nucleon form factors, we define the charged hadronic current for the QES process (see, e.g., Ref. [77]): Here V ud is the ud transition element from the Cabibbo-Kobayashi-Maskawa quark-mixing matrix and The form factors F i are in general complex functions of Q 2 . After standard calculations one finds with , The only difference between this result and that from Ref. [77] is in the relative sign of the terms in ω 6 which does not contribute to the QES cross section. 3 Inserting Eqs. (5) and (8) into Eq. (2) gives the commonly known formula for the differential cross section for reactions (1) on free nucleon targets:

Induced scalar and tensor form factors
The quoted formulas take into account the nonstandard G parity violating axial and vector second-class currents (SCC) which induce the nonzero scalar and tensor form factors F S and F T . The most robust restrictions on the SCC couplings F S,T (0) come from the studies of β decay of complex nuclei (see, e.g., Refs. [78,79] and quoted therein references). However, these studies are almost insensitive to the SCC effects at nonzero Q 2 . The latter were investigated in several (anti)neutrino experiments at BNL [22,3 According to Llewellyn Smith, the functions ω ′ 5 = ω5 − ω2 and ω6 are, respectively, the real and imaginary parts of a unique function. Our examination does not confirm this property for the general case of nonvanishing second-class current induced form factors FS and FT . 25,26,27] (Q 2 1.2 GeV 2 ) and in the IHEP-ITEP spark chamber experiment at Serpukhov [68] (Q 2 2.4 GeV 2 ), adopting the ad hoc dipole parameterizations The strongest (but yet not too telling) 90% C.L. upper limit for the axial SCC strength ξ T has been obtained at the BNL AGS ν µ experiment [27] as a function of the "tensor mass" M T , assuming conservation of vector current (CVC) (that is ξ S = 0), and simple dipole form for the vector and axial form factors with M V = 0.84 GeV and M A = 1.09 GeV. The limit ranges between 0.78 at M T = 0.5 GeV to about 0.11 at M T = 1.5 GeV. In so much as the contribution of the scalar form factor into the QES cross section is suppressed by (m µ /M ) 2 ≈ 0.01, the 90% C.L. constraint to the vector SCC strength ξ S is even less impressive: ξ S < 1.9, assuming ξ T = 0, M S = 1 GeV, and the same M V and M A as above.
Below, keeping in mind this vagueness, we will assume the time and charge invariance of the hadronic current. Under this standard assumption, all the form factors are real functions of Q 2 and to relate elastic and inelastic form factors, and imposes quark-hadron duality asymptotic constraints at high momentum transfers where the quark structure dominates. The parametrization is based on the same datasets as were used by Kelly [80], updated to include some recent experimental results. Quark-hadron duality implies that the squared ratio of neutron and proton magnetic form factors should be the same as the ratio of the corresponding inelastic structure functions F n 2 and F p 2 in the limit ξ p,n = 1: Here d and u are the partonic density functions. The authors fit the data under the two assumptions: d/u = 0 and d/u = 0.2. One more duality-motivated constraint is the equality applied for the highest Q 2 data points for the neutron electric form factor included into the BBBA(07) fit. The GKex(05) model is in fact a modification of the QCD inspired vector dominance model (VDM) by Gari and Krüempelmann (GK) [81] extended and fine-tuned by Lomon [82,83] in order to match the current and consistent earlier experimental data. The data set used by Lomon includes the polarization transfer measurements, which are directly related to the ratios of electric to magnetic form factors, and differential cross section measurements of the magnetic form factors. The electric form factors derived from the Rosenbluth separation of the differential cross section are only used for the lower range of Q 2 where the magnetic contributions are less dominant. Among several versions of the parametrization considered by Lomon, we chose the latest one "GKex(05)" described in Ref. [73]. This version incorporates the data that has become available since the publication [83] and has a bit better χ 2 . The fitted parameters agree with the known constraints and the model is consistent with VDM at low Q 2 , while approaching perturbative QCD behavior at high Q 2 . The quark-hadron duality constraint is not imposed. Figure 1 shows a comparison of the GKex(05) and BBBA(07) parametrizations for the form factors G p,n E and G p,n M divided by the standard dipole G D , against the experimental data extracted using either the Rosenbluth separation or polarization transfer techniques (including a series of double-polarization measurements of neutron knock-out from a polarized 2 H or 3 He targets). The data assemblage is borrowed from Refs. [84,85,86,87] and recent reviews [88,89]. It is seen from the figure that the models are numerically close to each other at low momentum transfers covered by experiment, but diverge at high Q 2 . The most serious disagreement between the models is in the neutron electric form factor at Q 2 2 GeV 2 . In section 4, we examine how the model differences affect the extracted value of the axial mass.

Axial-vector and induced pseudoscalar form factors
For the axial and pseudoscalar form factors we use the conventional parametrizations [77] where F A (0) = g A is the axial coupling, m π is the charged pion mass, and M A is the axial-vector mass treated as a free parameter. In fact, Eq. (10) is a conjecture inspired by the hypothesis of partial conservation of the axial current (PCAC), expectation that the form factor F P is dominated by the pion pole near Q 2 = 0, and the "technical" condition which is obviously fulfilled for the experimental lower limit of M A . Since the pseudoscalar contribution enters into the cross sections multiplied by (m ℓ /M ) 2 , the uncertainty caused by this approximation may only be important for ν τ /ν τ induced reactions (especially in the low-Q 2 range, see, e.g., Refs. [90,91]) and it is insignificant for reactions induced by electron and muon (anti)neutrinos.

Constants
The most precise determination of V ud comes from superallowed nuclear beta decays (0 + → 0 + transitions). We adopt the weighted average of the nine best measured superallowed decays V (SA) ud = 0.97377 ± 0.00027 recommended by the Particle Data Group (PDG) [92]. Note that this value is consistent with that of the PIBETA experiment at PSI [93], V (PIBETA) ud = 0.9728 ± 0.0030, obtained from the measured branching ratio for pion beta decay π + → π 0 e + ν.
For the axial-vector and Fermi coupling constants, we use the standard PDG averaged values: g A = −1.2695 ± 0.0029 and G F = 1.16637 × 10 −5 GeV 2 [92]. In several papers (see, e.g., Ref. [94] and references therein) it is suggested to use the value G ′ F = 1.1803×10 −5 GeV 2 obtained from 0 + → 0 + nuclear β decays, rather than the standard G F obtained from muon β decay. The coupling constant G ′ F subsumes the bulk of the inner radiative corrections. However, some neutrino experiments already take the radiative corrections into account (sometimes in quite different ways) in the measured cross sections. That is why, in this study, we simply add the corresponding difference (of about 2%) to the overall uncertainty of the fit. Note that using the G ′ F instead of G F would lead to a few percent decrease of the output value of M A .

Relativistic Fermi gas model
Since the main part of the experimental data on the QES cross sections for nuclear targets was not corrected for nuclear effects, we must take these into account in our calculations. In the present work, we use the RFG model by Smith and Moniz [10] incorporated as a standard tool into essentially all neutrino event generators employed in accelerator and astroparticle neutrino experiments.
According to RFG, the hadronic tensor W αβ given by Eq. (4) must be replaced with the tensor T αβ , which describes the bound nucleon. This tensor is of the same Lorentz structure as W αβ and is defined by the six invariant nuclear structure functions T i (Q 2 ). Thus, in the in the lab. frame where p lab = (M t , 0), M t is the mass of the target nucleus, and . The function n i (p) is the Fermi momentum distribution of the target nucleons, satisfying the normalization condition n i (p)dp = 1.
The factor 1 − n f (p + q) (the unoccupation probability) takes into account the Pauli blocking for the outgoing nucleon. The relative velocity v rel which represents the flux of incident particles, is given by Explicitly defining the three-momenta q, p, and p, q = (0, 0, |q|) , p = (sin θ k , 0, cos θ k ) |q|, p = (sin θ p cos φ p , sin θ p sin φ p , cos θ p ) |p|, one obtains v rel = [E p − |p| (cos θ k cos θ p + sin θ k sin θ p sin ϕ p )] /M t , where is the total energy of the bound nucleon and ǫ b is the effective binding energy. The angle θ k is defined by For determining the angle θ p , one can use the energy conservation law defined by delta-function is the total energy of the outgoing nucleon. Then the condition must be obeyed. The nuclear structure functions are the linear combination of the W i and can be straightforwardly calculated from Eqs. (4) and (11): All the fits are done with the CERN function minimization and error analysis package "MINUIT" (version 94.1) [110], taking care of getting an accurate error matrix. The errors of the output parameters quoted below correspond to the usual one-standard-deviation (1σ) errors (MINUIT default).
For the analysis, we have selected the most statistically reliable measurements of the total and differential cross sections for each nuclear target, which were not superseded or reconsidered (due to increased statistics, revised normalization, etc.) in the posterior reports of the same experimental groups. Finally, we include into the global fit the data on the total cross sections from Refs. [17,21,22,32,36,37,40,50,55,60,61,68,71] and the data for the differential cross sections from Refs. [50,57,64,65,68,71,101]. The remaining data are either obsolete, or exhibit uncontrollable systematic errors and/or fall well outside the most probable range determined through the fit of the full dataset; the value of χ 2 evaluated for each subset of the rejected data usually exceeds (3 − 4) NDF.
Since the differential cross sections dσ/dQ 2 were measured, as a rule, within rather wide ranges of the energy spectra of ν µ and ν µ beams, we use only the data from such experiments, in which the spectra were known (measured or calculated and then calibrated) with reasonably good accuracy. All the energy spectra (borrowed from Refs. [50,56,68,101,104,111,112]) necessary for numerical averaging of the calculated differential cross sections and distributions were parametrized. To avoid the loss of accuracy, the precision of these parametrizations was chosen to be at least an order of magnitude better than the experimental accuracy of the spectra themselves. For a verification, we have estimated the mean energies of the beams for different energy intervals, and have compared these against the published values.
The analyses were performed for neutrino and antineutrino data separately, and for the full set of the ν and ν data together. For each fit, we have included the data for either total or differential cross sections, as well as for the cross sections of both types together. The main results of the analysis are summarised in Tables 2 and 3 and illustrated in Figs. 2-19. Let us discuss these results in details.

Main results of the global fit
As is seen from Table 2, the differences between the values of M A extracted from the fits of each type, performed with the BBBA(07) and GKex(05) models for the vector form factors vary between 0.3% and 1.3% that is less than or of the order of one standard deviation in the M A extractions and is comparable with the accuracy of the most precise measurements of the electric and magnetic form factors. The values of χ 2 /NDF are essentially the same for BBBA(07) and GKex(05). The differences in the M A values obtained with the two versions of the BBBA(07) model corresponding to d/u = 0 and 0.2 (the latter is not shown in the table) are less than 0.2% that is practically negligible. Therefore, in the following we will solely discuss the d/u = 0 case.
The M A values obtained from the fits to the differential cross sections are systematically lower those obtained from the total cross sections. The differences amount ∼ 1.5% (∼ 5.7%) for ν µ (ν µ ) that is (especially in antineutrino case) above the statistical error of the fit and is caused mainly by uncertainties in the energy spectra of ν µ and ν µ and, in lesser extent, in the nuclear effects. Figures 2 and 3 show a compilation of the available data on the total QES cross sections for the following nuclear targets: hydrogen [21], deuterium [14,15,16,17 [37], steel [12], propane [43], freon [40,45,48,50,62,70,71,104], and also propane-freon [52,55,56] and neonhydrogen [33,34,36] mixtures. The recent MiniBooNE 2007 datapoint [5] (carbon target) estimated from the reported value of M A is also shown in Fig. 2 for comparison.
The compilation does not include obviously obsolete data (e.g., ANL 1972 [13], CERN HLBC 1965/1966 [38,39]), as well as the data identical to those reported in the posterior publications of the same experimental groups  Table 2. Values of MA (given in GeV), extracted by fitting the νµ, νµ, and νµ + νµ data on total and differential QES cross sections, using the BBBA(07) and GKex(05) models for the vector form factors of the nucleon. The χ 2 /NDF values for each fit are shown in parentheses. preliminary and are reproduced here by permission of the NOMAD Collaboration.
All the deuterium data quoted in Fig. 2 and freon data in Fig. 3 were converted to a free nucleon target by the experimenters. 4 The BNL 1981 experiment [22] had reported the E ν and Q 2 dependencies of M A extracted from a fit of the experimental Q 2 distribution rather than the cross section; we quote the BNL 1981 cross section recalculated from M A by Kitagaki et al. [32]. Similarly, the FNAL 1984 rectangle [33,34] and FNAL 1987 datapoint [36] were calculated by the experimenters (for free proton target) using the M A value extracted from the measured Q 2 distribution of ν µ events recorded in the Fermilab 15' bubble chamber filled with a heavy neon-hydrogen mixture. The data from several freon experiments (e.g., [40,45,97]) reported in the original papers in units cm 2 per nucleon of freon nucleus, were converted to the standard units.
All solid curves shown in the figures were calculated using the BBBA(07) model for vector form factors with d/u = 0 and always correspond to the best fit value obtained from the global fit of neutrino and antineutrino data on the total and differential cross sections (see Table  2). We do not show the cross sections calculated with the GKex(05) model since the difference will be practically invisible.
The dashed curves in Fig. 2 are calculated with the M A values extracted from the best fit to the (preliminary) NOMAD total cross section data alone [61]: 4 The nuclear corrections applied to the deuterium data under consideration, were treated according to Singh [113]. The nuclear effects for the freon data were modeled using a Fermi gas approach.
both agree with the global fit value (12). Note that these results were obtained with the GKex(05) vector form factors. Fitting the NOMAD data with the BBBA(07) form factors increases M ν A and M ν A by about 0.8 and 0.9%, respectively, that still remains well within the errors quoted in (13).
As is seen from the figures, the obtained result, despite the non-optimal χ 2 and large spread of the data, is not in conflict with the main part of the data excluded from the global fit. Moreover, it well agrees with the world averaged value of obtained in Ref. [3] as a result of their reanalysis of the "raw" data from ν µ d and ν According to the global fit (see Table 2), the difference between the values of M ν A and M ν A obtained by fitting the neutrino and antineutrino data separately, reaches about 3.5% for BBBA(07) and about 3.2% for GKex(05) that is above the statistical error in determination of M ν A and M ν A . However, taking into account the systematic difference between the fits of total and differential cross section data, as well as high values of χ 2 /NDF, this difference cannot be considered statistically significant. Furthermore, the fit to the antineutrino data is not stable relative to including/excluding some data subsets. In particular, as is seen from Fig. 2, the total NuTeV cross sections per nucleon bound in iron, averaged over the energy range E ν,ν = 30 ÷ 300 GeV Fig. 2 by rectangles) notably exceed the corresponding best fit curves whereby the NuTeV data [37] strongly affects the global fit values of M ν A and M ν A . To clarify this point further, we have performed additional fits, in which the datasets obtained in experiments with non-active targets have been removed. Namely, we excluded the highest energy NuTeV total cross section data (iron target) [37] and the data on differential cross sections measured with the IHEP-ITEP spark chamber detector with aluminium filters [63,65,68], since these experiments do not have an active target to measure recoil hadrons and surely remove resonance background. In order to minimize possible uncertainties in nuclear corrections, the lowest-energy CERN 1967 total cross section data (freon target) [40] were also excluded from these fits. The results of this analysis are summarized in Table 3. It is seen that the additional reduction of the dataset essentially decreases the resulting values of M A . Concurrently it improves the statistical quality of the fits to the total cross section data, while slightly increases the χ 2 /NDF for the fit to the differential cross sections. Besides that, the M A values extracted from the total and differential cross sections become bit more consistent. The differences between M ν A and M ν A [-65 MeV for BBBA(07) and -75 MeV for GKex(05)] become opposite in sign to those obtained from our "default" fit performed with the full dataset. However, both M ν A and M ν A values are still compatible, within the 1σ deviation, with the average value of M ν+ν A . So we may reckon that (i) the axial mass extraction is rather responsive to the choice of the data subsets and (ii) the current experimental data cannot definitely confirm or disconfirm possible difference between the axial masses extracted from experiments with neutrino and antineutrino beams. Similar fit performed for the differential cross section data only, from which all the ν µ d data were excluded, leads to an increase of M ν A by about 4.2% (4.4%) for BBBA(07) (GKex(05)). However, the statistical error of this fit increases too. Including into this fit the non-deuterium data on total cross sections diminish the increase of M ν A to about 1.2% for both BBBA(07) and GKex(05). Hence, the above conclusions remain essentially unchanged.

Further details on differential cross section data
As is known from the comparison with the low-energy electron-nucleus scattering data, the RFG description of the low-Q 2 region is not enough accurate especially at energies below ∼ 2 GeV (for recent discussion, see, e.g., Refs. [115,116] and references therein). Moreover, the shape of dσ/dQ 2 at Q 2 0.1 GeV 2 is slowly sensitive to variations of M A (see below). Thus, in order to minimize possible uncertainties due to nuclear effects, the points with Q 2 < 0.15 GeV 2 were rejected from the fit of the differential cross section dataset. Leaving these points in the fit would lead to a decrease of the output values of M ν A , M ν A , and M ν,ν A obtained from the dσ/dQ 2 dataset by, respectively, 1.8, 3.3, and 2.2% for BBBA(07) and 2.0, 4.0, and 2.6% for GKex(05) form factors. The corresponding decrease of M A derived from the full dataset (σ and dσ/dQ 2 ) is clearly less essential: respectively, 0.7, 1.3, and 0.9% for BBBA(07) and 0.7, 1.5, and 1.0% for GKex(05).
Of course, the mentioned uncertainty still remains in the RFG calculations of the total cross sections, since the contribution from the low-Q 2 region is essential at low energies. To illustrate this, we show in Fig. 4 the relative contribution of the region Q 2 < Q 2 1 into the total cross section, R Q 2 1 = σ Q 2 < Q 2 1 /σ, as a function of Q 2 1 , evaluated for ν µ and ν ν QE interactions with carbon at several (anti)neutrino energies using M A = 1 GeV. 5 It is seen that for neutrino-nucleus interactions R 0.25 as Q 2 1 < 0.15 GeV 2 and E ν > 0.7 GeV that is for all energies of our current interest. As a result, a few percent error expected in dσ/dQ 2 due to inaccuracy of the RFG model for the low-Q 2 region, becomes nearly negligible in the total cross section. However it is not the case for antineutrino interactions, for which the ratio R Q 2 1 = 0.15 GeV 2 becomes reasonably small (R 0.3) only for E ν 2 GeV. Therefore the lower energy antineutrino total cross section data may bias an uncontrolled (while still small) additional uncertainty. Fortunately, the major part of the data participated in the global fit satisfies the above conditions and our examination demonstrates that the related uncertainty is not weighty. Figures 5-7 (a) and 8-11 represent the spectrum-averaged differential cross sections for several nuclear targets: deuterium (Fig. 5) [57], aluminium (Figs. 6 and 7 (a)) [63, 65,68], freon (Figs. 8, 9 and 10) [50,70,71,101,104], and propane-freon mixture (Fig. 11) [55,52]. In Fig. 7 (b) we show (for illustrative purposes only) the axial-vector form factor extracted in the IHEP-ITEP spark chamber experiment [68]. All the quoted data, except those from Ref. [70] (superseded by the data from the more recent publication by the SKAT Collaboration [71]), model-dependent IHEP-ITEP data on F A (Q 2 ) [68], and a few rejected low-Q 2 datapoints, participate in the global fit. We show the cross sections calculated with M A obtained by individual fits to the data of each experiment alone and compare these against the cross sections evaluated with the global-fit value of M A . All the details are recounted in the captions and legends of the figures. The comparison demonstrates that the individual and global fits generally Table 3. The same as in Table 2 but after exclusion of the datasets from experiments with non-active targets (NuTeV 1984 [37], IHEP-ITEP 1981,82,85 [63, 65, 68]) and the lowest-energy data of CERN 1967 [40] (see text for details).

BBBA(07)
GKex(05) do not contradict each other. The differences are within the experimental errors and are not of systematic nature.
As a further test of the global fit, we show in Fig.  12 the flux-weighted differential cross sections dσ(ν µ n → µ − p)/dy and dσ(ν µ p → µ + n)/dy (divided by energy), which were measured with the Gargamelle bubble chamber filled with liquid freon and exposed to the wide-band CERN-PS ν µ and ν µ beams. Several analyses of these data samples are available from the literature (see Refs. [46,47,49,101] and also Ref. [105] for a review). Figure 12 shows two representative versions taken from Refs. [47] and [101] -the preliminary and final results of the GGM experiment, respectively. The data are shown for the five narrow instrumental ranges: 1 − 2, 2 − 3, 3 − 5, 5 − 11, and 5 − 20 GeV. The measured cross sections were converted from freon to a free nucleon target by the experimenters, after accounting for Fermi motion of the nucleons and Pauli suppression of quasielastic events.
For a qualitative comparison, we have performed individual fits to the GGM data, separately for neutrino and antineutrino differential cross sections. In order to reduce possible error introduced by RFG calculations of nuclear effects, the energy range of 1 − 2 GeV has been excluded from this likelihood analysis. As is seen from the figure, the M A value extracted from the neutrino subsample does not contradict to that from the global fit, while it is not so for the antineutrino data subsample where the discrepancy is essential. This discrepancy can be attributed (at least, partially) to the vagueness of the model for nuclear effects used in the analyses of the GGM data. Since the details of the GGM nuclear Monte Carlo are not available, we do not include this data sample into the global fit. We note, however, that the inclusion of these data (also without the low-energy datapoints) into the fit only leads to a small decrease of the output values of M ν A , M ν A , and M ν,ν A -by, respectively, 0.4, 2.2, and 0.9% for BBBA(07) and 0.3, 2.0, and 0.8% for GKex(05) form factors. The corresponding χ 2 /NDF values remain nearly the same.

Q 2 distributions
An additional fruitful set of available data is the Q 2 distributions dN/dQ 2 of the QES events measured in several experiments with different nuclear targets. Usually just dN/dQ 2 is considered as the observable most appropriate for extracting axial mass value, since it is less dependent of the flux and spectrum uncertainties in comparison with the differential or total cross sections. However, in comparison with the differential cross section, the Q 2 distribution has two drawbacks: it contains an uncertainty due to normalization, and it is generally less responsive to variations of M A at high Q 2 . Figure 13 illustrates the second point. It shows the Q 2 distributions and differential cross sections for ν µ and ν µ quasielastic scattering off free nucleons, evaluated with different values of M A and normalized to the corresponding quantities calculated with M A = 1 GeV. The calculations are done with the fixed values of energy corresponding to the mean (anti)neutrino beam energies in experiments [34,43,56,62]. It is seen from the figure that the region Q 2 0.15 GeV 2 strongly affected by the nuclear effects, is sensitive to M A for dN/dQ 2 and less sensitive for dσ/dQ 2 ; the situation is opposite for the high Q 2 region for which the nuclear corrections are less important.
We use the measured Q 2 distributions for a consistency test of our analysis. For illustration, we show the four sets of data on Q 2 distributions measured in experiments HLBC 1969 (propane) [43] (Fig. 14), IHEP SKAT 1981 (freon) [62] (Fig. 15), CERN GGM 1979 (propane-freon mixture) [56] (Fig. 16), and FNAL E180 (neon-hydrogen mixture) [33,34] (Fig. 17). The curves shown in the figures are calculated with the global-fit M A and normalized to the data after fitting of the normalization factor N . The shaded bands indicate the uncertainty due mainly to indetermination of this factor. The obtained best-fit values of N should be compared with these evaluated directly from the experimental data (all values are shown in the legends of the figures). One can see that the agreement is excellent everywhere. So, we may conclude that this test was quite successful.
Another important confirmation of our result is a reasonably good agreement with the M A value extracted in our earlier analysis of the data on total inelastic ν µ N and ν µ N CC cross sections and relevant observables [117].
Finally, Fig. 18 presents a comparison of the total QES cross sections for ν e , ν µ , ν τ , ν e , ν µ , and ν τ interactions with free nucleons, calculated with the obtained best-fit value of M A = 0.999 ± 0.011 GeV by using the BBBA(07) model of vector form factors. The shaded bands reproduce the uncertainty due to the 1σ error in M A .

Discussion and conclusions
We performed a statistical study of the QES total and differential cross section data in order to extract the best-fit values of the parameters M A . Our main results are summarized in Table 2 are, of course, model dependent and can be recommended for use only within the same (or numerically equivalent) model assumptions as in the present analysis. The best-fit values of the axial mass obtained by different fits do not contradict to each other and agree with the recent re-extraction of M A from ν µ d, ν µ H, and pion electroproduction experiments, reported in Ref. [3]. They are also in agreement with the preliminary result of high-statistical NOMAD experiment at CERN, as well as with the numerous earlier data which were not included into the likelihood analysis. It has been demonstrated that removing the data subsets obtained in experiments with non-active targets, particularly the NuTeV dataset, leads to a further decrease of the extracted values of M A (see Table 3). In other words, there is no way to increase the M A value which follows from essentially all (anti)neutrino data on total and differential QES cross sections.
On the other hand, our best-fit value of M A is in a conflict with the mean values of M A reported by K2K and MiniBooNE Collaborations [4,5], even after accounting for the maximum possible systematic error of our analysis related primarily to its susceptibility to the choice of the data subsets. To expound the problem, let us consider the representative K2K result with more details.
The M A value reported in Ref. [4] has been obtained with a water target by fitting the Q 2 distributions of muon tracks reconstructed from neutrino-oxygen quasielastic interactions by using the combined K2K-I and K2K-IIa data from the Scintillating Fiber detector (SciFi) in the KEK accelerator to Kamioka muon neutrino beam. The experimental data from the continuation of the K2K-II period were not used in the analysis of Ref. [4]. The best-fit values of M A obtained from the K2K-I and K2K-IIa data subsets separately are, respectively, 1.12 ± 0.12 GeV (χ 2 /NDF = 150/127) and 1.25 ± 0.18 GeV (χ 2 /NDF = 109/101). Figure 19 shows the ν µ n → µ − p total cross section per neutron bound in oxygen, recalculated from the fitted values of M A derived in Ref. [4] from the Q 2 distribution shape for each reconstructed neutrino energy. It is necessary to underline here that the authors do not consider their result for each energy bin as a measurement, but rather a consistency test. All calculations represented in Fig. 19 were done with our default inputs that introduces an uncertainty of at most 2%; this uncertainty is added quadratically to the quoted error bars. Also shown are the cross sections evaluated by using our best fit value (12), the K2K value of 1.20±0.12 GeV, and the value of 1.1 GeV used as a default in the recent neutrino oscillation analyses to the data from K2K [118,119] and Super-Kamiokande I [120]. A significant systematic discrepancy is clearly seen at E ν > 1 GeV. Since the energy region covered by the K2K analysis extends to about 4 GeV, it seems problematic to explain this discrepancy by the inapplicability of the RFG model alone.
Considering that the low-energy K2K and MiniBooNE data are in agreement with each other and do not contradict to the high-energy NuTeV results, we may conclude that the new generation experiments for studying the quasielastic neutrino and antineutrino interactions with nucleons and nuclei are of urgent necessity, in order to resolve the inconsistencies between the old and new measurements of the axial-vector mass.

Acknowledgements
This study is currently supported by the Russian Foundation for Basic Research under Grant No. 07-02-00215-a. The authors would like to thank Krzysztof M. Graczyk, Sergey A. Kulagin, Dmitry V. Naumov, Jan T. Sobczyk, and Oleg V. Teryaev for helpful discussions. We thank the NOMAD Collaboration for permission to use their data prior to publication and Antonio Bueno for explaining us some points of LAr TPC experiment. We are especially grateful to Arie Bodek for his constructive comments and suggestions. V. V. L. is very thankful to LPNHE (Paris) for warm hospitality and financial support during a stage of this work.    Fig. 4. The ratio R = σ`Q 2 < Q 2 1´/ σ vs. Q 2 1 , evaluated for νµ and νν quasielastic interactions with carbon target at several (anti)neutrino energies. The MA value is taken to be 1 GeV.  5. Flux-weighted differential cross section for νµn → µ − p measured in the WA25 experiment with the CERN bubble chamber BEBC filled with deuterium and exposed to high-energy νµ beam at the CERN-SPS [57]. The data were converted to a free neutron target by the authors of the experiment. The curves are the calculated cross sections averaged over the experimental νµ energy spectrum borrowed from Ref. [112]. The energy range and estimated mean energy are given in the legend. The dashed curves are for the best fit to the WA25 data, while the solid curves correspond to the global fit to all QES data. Shaded band represents 1σ deviation from the best-fitted value of MA given in the legend.    Fig. 8. Flux-weighted differential cross sections for νµn → µ − p (a) and νµp → µ + n (b) measured with the heavy-liquid bubble chamber Gargamelle filled with heavy freon and exposed to the CERN-PS νµ and νµ beams [50,101]. The error bars contain the statistical fluctuation and the indetermination on the νµ and νµ fluxes. The curves are the calculated cross sections averaged over the experimental νµ and νµ energy spectra given in Ref. [50]. Only the events with Eν, ν > 1.5 GeV were accepted. The dashed curves are for the best fit to the GGM 1977 data, while the solid curves correspond to the global fit to all QES data. The points shown by grey symbols are excluded from the fits (see text). Shaded bands represent 1σ deviations from the best-fitted values of MA given in the legends.  Fig. 9. Flux-weighted differential cross sections for νµn → µ − p (a) and νµp → µ + n (b) measured with the freon filled bubble chamber SKAT exposed to the U70 broad-band νµ and νµ beams of the Serpukhov PS [70,104] (see also Refs. [69] for the earlier analyses of the same data sample). The data were converted to a free nucleon target by the authors of the experiment. The inner and outer bars indicate statistical and total errors, respectively; the systematic error includes the uncertainties due to the cross section normalization and nuclear Monte Carlo. The curves are the calculated cross sections averaged over the experimental νµ and νµ energy spectra borrowed from Ref. [104]. The energy range and estimated mean energies are given in the legends. The dashed curves are for the best fit to the SKAT 1988 data, while the solid curves correspond to the global fit to all QES data (the SKAT 1988 data are excluded from the global fit). Shaded bands represent 1σ deviations from the best-fitted values of MA given in the legends.   Fig. 10. Flux-weighted differential cross sections for νµn → µ − p (a) and νµp → µ + n (b) measured with the freon filled bubble chamber SKAT exposed to the U70 broad-band νµ and νµ beams of the Serpukhov PS [71]. The data were converted to a free nucleon target by the authors of the experiment. The inner and outer bars indicate statistical and total errors, respectively; the systematic error includes the uncertainties due to the cross section normalization and nuclear Monte Carlo. The curves are the calculated cross sections averaged over the experimental νµ and νµ energy spectra borrowed from Ref. [104]. The energy range and estimated mean energies are given in the legends. The dashed curves are for the best fit to the SKAT 1990 data, while the solid curves correspond to the global fit to all QES data. The points shown by grey symbols are excluded from the fits (see text). Shaded bands represent 1σ deviations from the best-fitted values of MA given in the legends.   Fig. 11. Flux-weighted differential cross sections for νµn → µ − p (a) and νµp → µ + n (b) measured with the bubble chamber Gargamelle filled with light propane-freon mixture and exposed to the CERN-PS νµ and νµ beams [55,52]. The inner and outer bars in panel (a) indicate statistical and total errors, respectively; the error bars in panel (b) contain the statistical fluctuation and the indetermination on the νµ flux. The curves are the calculated cross sections averaged over the experimental νµ and νµ energy spectra given in Refs. [ Fig. 14. Flux-weighted Q 2 distribution for νµn → µ − p measured with the CERN heavy-liquid bubble chamber (HLBC) filled with propane and exposed to the CERN PS νµ beam [43]. The curve is the distribution calculated with MA obtained from the global fit, averaged over the experimental νµ energy spectrum from Ref. [111], and normalized to the HLBC 1969 data. The spectrum is estimated to be accurate within ±15% (the error includes an estimate of systematic effects). The energy range and estimated mean energy are given in the legends. Shaded band represents 1σ variation from the average due to uncertainties in MA and normalization factor N .  Fig. 15. Flux-weighted Q 2 distribution for νµn → µ − p measured with the freon filled bubble chamber SKAT exposed to the U70 broad-band νµ beam of the Serpukhov PS [62]. The data were converted to a free nucleon target by the authors of the experiment. The inner and outer bars indicate statistical and total errors, respectively; the systematic error includes the uncertainties due to the flux normalization and nuclear Monte Carlo. The curve is the distribution calculated with MA obtained from the global fit, averaged over the experimental νµ energy spectrum from Ref. [104], and normalized to the SKAT 1981 data. The energy range and estimated mean energy are given in the legends. Shaded band represents 1σ variation from the average due to uncertainties in MA and normalization factor N .   Fig. 17. Flux-weighted Q 2 distribution for νµp → µ + n measured in the FNAL E180 experiment with a 15' bubble chamber filled with heavy neon-hydrogen mixture (64% of neon atoms) and exposed to the FNAL wide-band νµ beam [33,34] (see also Ref. [31] for an earlier version). The curve is the distribution calculated at the mean antineutrino energy of 12.7 ± 0.2 GeV, with MA obtained from the global fit and then normalized to the E180 data. [The spectrum averaging procedure cannot be applied here, since the νµ spectrum has been evaluated just from the quoted Q 2 distribution.] Shaded band represents 1σ variation from the average due to uncertainties in MA and normalization factor N .