An Update of the HLS Estimate of the Muon g-2

A global fit of parameters allows us to pin down the Hidden Local Symmetry (HLS) effective Lagrangian, which we apply for the prediction of the leading hadronic vacuum polarization contribution to the muon g-2. The latter is dominated by the annihilation channel e+e- ->pi+pi-, for which data are available by scan (CMD-2 and SND) and ISR (KLOE-2008, KLOE-2010 and BaBar) experiments. It is well known that the different data sets are not in satisfactory agreement. In fact it is possible to fix the model parameters without using the pi+pi- data, by using instead the dipion spectra measured in the tau decays together with experimental spectra for the pi0 gamma, eta gamma, pi+pi-pi0, K+K-, K0bar K0 final states supplemented by specific meson decay properties. Among these, the accepted decay width for rho0 ->e+e- and the partial widths and phase information for the omega/phi ->pi+pi- transitions, are considered. It is then shown that, relying on this global data set, the HLS model, appropriately broken, allows to predict accurately the pion form factor up to 1.05 GeV. It is shown that the data samples provided by CMD-2, SND and KLOE-2010 behave consistently with each other and with the other considered data. Consistency problems with the KLOE-2008 and BaBar data samples are substantiated."All data"global fits are investigated by applying reweighting the conflicting data sets. Constraining to our best fit, the broken HLS model yields a_mu(th) = (11659169.55 [+1.26 -0.59]_phi +[+0.00 -2.00]_tau +/- 5.21_(th))~10**-10 associated with a very good global fit probability. Correspondingly, we find that Delta a_mu=a_mu (exp)- a_mu (th) exhibits a significance ranging between 4.7 and 4.9 sigma.


Introduction
The theoretical value for muon anomalous magnetic moment a µ is an important window in the quest for new phenomena in particle physics. The predicted value is the sum of several contributions and the most prominent ones are already derived from the Standard Model with very high accuracies. The QED contribution is thus estimated with an accuracy of a few 10 −12 [1,2,3] and the precision of the electroweak contribution is now of order 10 −11 [4]. The light-by-light contribution to a µ is currently known with an accepted accuracy of 2.6 × 10 −10 [5].
Presently, the uncertainty of the Standard Model prediction for a µ is driven by the uncertainty on the leading order (LO) hadronic vacuum polarization (HVP) up to ≃ 2 GeV [6,7]. This region is covered by the non-perturbative regime of QCD and the leading order HVP (LO − HV P ) is evaluated by means of : which relates the hadronic intermediate state contributions {H i , i = 1 · · · n} to the annihilation cross sections σ(e + e − → H i ) ≡ σ H i (s). K(s) is a known kernel [4] enhancing the weight of the threshold region s H i and s cut is some energy squared where perturbative QCD starts to be applicable. In the region where perturbative QCD holds 1 , its contribution to a µ carries an uncertainty of the order of a few 10 −11 . Up to very recently, the single method used to get the a µ (H i )'s was to plug the experimental cross sections into Eq. (1). Among the most recent studies based on this method, let us quote [8,6,7,9]. When several data sets cover the same cross section σ H i (s), Eq. (1) is used with some appropriate weighting of the various spectra, allowing to improve the corresponding a µ (H i ).
On the other hand, it is now widely accepted that the Vector Meson Dominance (VMD) concept applies to low energy physics [10,11]. VMD based Effective Lagrangians have been proposed like the Resonance Chiral Perturbation Theory or the Hidden Local Symmetry (HLS) Model; it has been proven [12] that these are essentially equivalent. Intrinsically, this means that there exist physics correlations between the various e + e − → H j annihilation channels. Therefore, it becomes conceptually founded to expect improving each a µ (H i ) by means of the data covering the other channels e + e − → H j (j = i).
This is basically the idea proposed in [13] relying on the HLS model [14,15]. Using a symmetry breaking mechanism based on the simple BKY idea [16] and a vector meson mixing scheme, the model has been developed stepwise [17,18,19,20] and its most recent form [13] has been shown to provide a successful simultaneous description of the e + e − annihilation into the π + π − , π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 final states as well as the τ ± → ν τ π ± π 0 decay spectrum. Some more decays of the form 2 V → P γ or P → γγ are considered. 1 The charmonium and bottomium regions carry uncertainties also in the range of a few 10 −11 . 2 We denote by V or P resp. any meson belonging to the (basic) vector or pseudoscalar lowest mass nonets.
As higher mass meson nonets are absent from the standard HLS model, its energy scope is a priori limited upwards by the φ meson mass region (≃ 1.05 GeV). However, as this region contributes more than 80 % to the total HVP, improvements which can follow from the broken HLS model are certainly valuable 3 .
The global simultaneous fit of the data corresponding to the channels quoted above allows to reconstruct the various cross sections σ H i (s) taking automatically into account the physics correlations inside the set H ≡ {H i } of possible final states and decay processes. The fit parameter values and the parameter error covariance matrix summarize optimally the full knowledge of H. This has two important consequences : • One should get the {a µ (H i ), i = 1, · · · n} with improved uncertainties by integrating the model cross sections instead of the measured ones. Indeed, the functional correlations among the various cross sections turn out to provide (much) larger statistics in each channel and thus yield improved uncertainties for each a µ (H i ).
• When several data samples cover the same process H i , one has a handle to motivatedly examine the behavior of each within the global fit context. Stated otherwise, the issue of the consistency of each data set with all the others can be addressed with the (global) fit probability as a tool to detect data samples carrying problematic properties.
Up to now, the broken HLS model (BHLS) [13] -basically an empty shell -has been fed with all existing data sets 4 for what concerns the annihilation channels π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 , with the spectra from ALEPH [21], CLEO [22] and BELLE [23] for the τ dipion decay 5 and with the V P γ/P γγ partial width information extracted from the Review of Particle Properties (RPP) [24]. This already represents more than 40 data sets collected by different groups with different detectors; one may thus consider that the systematics affecting these data sets wash out to a large extent within a global fit framework.
For what concerns the crucial process e + e − → π + π − , the analysis in [13] only deals with the data sets collected in the scan experiments performed at Novosibirsk and referred to globally hereafter as NSK [25,26,27,28,29]. The main reason was, at this step, to avoid discussing the reported tension [30,8] between the various existing π + π − data sets : the scan data sets just quoted, and the data sets collected using the Initial State Radiation (ISR) method by KLOE [31,32] and BaBar [33,34], not to mention the pion form factor data collected in the spacelike region [35,36].
It has thus been shown that the global fit excluding the ISR data sets, allows to yield a splendid fit quality; this proves that the whole collection of data sets considered in [13] is selfconsistent and may provide a safe reference, i.e. a benchmark, to examine the behavior of other data samples.
Using the fit results, the uncertainty on the contribution to a µ of each of the annihilation channels considered was improved by -at least -a factor of 2, compared to the standard 3 The broken HLS model does not include the 4π, 5π, 6π, ηππ and ωπ annihilation channels. Therefore, the (small) contribution of these missing channels [13] to a µ should be still evaluated by direct integration of the experimental cross sections; up to the φ mass, this amounts to [13] (1.55 ± 0.57 tot ) 10 −10 . 4 The full list of data sets can be found in [19] or [13] together with a critical analysis of their individual behavior. 5 The energy region used in the fits has been limited to the [2m π , 1 GeV] interval where the three data sets are in accord with each other. This should lessen the effect of some systematic effects. estimation method based on the numerical integration of the measured cross sections. For the case of the π + π − channel, the final uncertainty was even found slightly better than those obtained with the standard method by merging scan and ISR data, i.e. a statistics about 4 times larger in the π + π − annihilation channel.
The main purpose of the present study is an update of the work in [13] aiming at confronting all scan (NSK) and ISR (from BaBar and KLOE) -and even spacelike [35, 36] -π + π − data and reexamine the reported issues [30,8]. The framework in which our analysis is performed is the same as the one motivated and developed in [13]. The broken HLS model described in [13] happens to provide a tool allowing to compare the behavior of any of these π + π − data sets when confronted with the π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 annihilation data and with the τ dipion spectra. Indeed, the latter data alone, supplemented with some limited information extracted from the Review of Particle Properties 6 (RPP) [24], allow to predict the pion form factor with a surprisingly good precision. The additional RPP information is supposed to carry the Isospin Breaking (IB) information requested in order to derive reliably the π + π − information from the knowledge of the π ± π 0 spectrum. We also take profit of the present work to update the numerical values for some contributions to the muon anomalous moment a µ , all gathered in Table 10 of [13]. Thus, we update the QED entry by using the recent spectacular progress by Aoyama, Hayakawa, Kinoshita and Nio [1,2]. They have been able to perform a complete numerical calculation of the 5-loop QED corrections to a e and a µ . On the other hand, the electroweak contribution, which depends on the Higgs mass at 2-loops is now better known if we accept that ATLAS [37] and CMS [38] have observed at the LHC the Higgs boson at a mass of about 125 GeV in a narrow window. Using this information slightly changes the central value as well as the uncertainty of the EW entry. We also have reevaluated the higher order HVP contribution (HO) within the standard approach based on all π + π − channels (i.e. all scan and ISR data).
The paper is organized as follows. Section 2 reminds the motivations of the BHLS model and a few basic topics concerning the ππ channel description (from [13]); we also reexamines how the isospin breaking corrections apply. In Section 3, the detailed framework -named "τ +PDG" -used to study the differential behavior of the scan and ISR data is presented. Thanks to the (wider than usual) energy range covered by the BaBar spectrum [33, 34], a detailed study of the π + π − spectrum in the φ region can be performed for the first time. This leads to update the φ → π + π − treatment within our computer code; this is emphasized in Subsection 3.3. In Section 4, one confronts the "τ +PDG" predictions with the available scan (NSK) and ISR data samples; it is shown that the NSK data and both KLOE data samples (referred to hereafter as KLOE08 [31] and KLOE10 [32]) have similar properties while BaBar behaves differently, especially in the ρ − ω interference region. Section 5, especially Subsection 5.1, reports on the global fits performed using the various π + π − data samples each in isolation or combined. Subsection 5.2 collects some topics on various aspects of the physics covered by the HLS model. More precisely, Subsection 5.2.1 is devoted to studying the φ region of the pion form factor and Subsection 5.2.2 gives numerical fit information which may allow to compare with corresponding results available from other studies performed using different methods. In Section 6, we focus on the consequences for the muon anomalous moment a µ of the various scan and ISR π + π − spectra and compare results with the BNL [39,40] measurement. The a main ingredient. Phenomenologically, scalars only play a kind of "next-to-leading" role.
The situation is quite different for the spin 1 bosons, which naturally acquire a Yang-Mills effective structure in a low-energy expansion i..e., they naturally exhibit a leading local gauge symmetry structure with masses as generated by a Higgs-Kibble mechanism. There is one important proviso, however : such a low-energy effective structure is pronounced only to the extent that the effective expansion scale Λ eff is high enough, which is not clear at all for QCD unless we understand why Λ eff >> Λ OCD ∼ 400 MeV. However, there is a different approach, namely, a "derivation" of the Extended Nambu-Jona-Lasinio (ENJL) model [43,44], which has also been proved to be largely equivalent to the Resonance Lagrangian Approach (RLA) [45,46]. Last but not least, large-N c QCD [47,48,49,50] in fact predicts the low energy hadron spectrum to be dominated by spin 1 resonances. These arguments are also the guidelines for the construction of the HLS model [14,15]. It provides a specific way to incorporate the phenomenologically known low energy hadron spectrum into an effective field theory.
Most frequently the RLA is applied to study individual processes. In this paper as in a few previous ones, we attempt to fit the whole HLS Lagrangian by a global fit strategy. This is, in our opinion, the only way to single out a phenomenologically acceptable low-energy effective theory, which allows to make predictions which can be confronted with experiments.

The Broken HLS Lagrangian
The (unbroken) HLS Lagrangian is then given by L HLS = L A + aL V , where a is a basic HLS parameter not fixed by the theory, which should be constrained by confrontation with the data. From standard VMD models, one expects a ≃ 2.
It is well known that the global chiral symmetry G global is not realized as an exact symmetry in nature, which implies that the ideal HLS symmetry is evidently not a symmetry of nature either. Therefore, it has obviously to be broken appropriately in order to provide a realistic low energy effective theory mimicking low energy effective QCD.
Unlike in ChPT where one is performing a systematic low energy expansion in low momenta and the quark masses, here one introduces symmetry breaking as phenomenological parameters to be fixed from appropriate data. Since a systematic low energy expansionsà la ChPT does not converge above about ≃ 400 MeV, this is the only way to model phenomenology up to, and including, the φ resonance region.
In our approach, the Lagrangian pieces in Eqs. (2) are broken in a two step procedure. A first breaking mechanism named BKY is used, originating from [16,15]. In order to avoid some undesirable properties [51,52] of the original BKY mechanism, we have adopted the modified BKY scheme proposed in [17]. In its original form, this modified BKY breaking scheme only covers the breaking of the SU(3) symmetry; following [53], it has been extended in order to include isospin symmetry breaking effects. This turns out to modify Eqs. (2) by introducing two constant diagonal matrices X A/V : and the (non-zero) entries in X A/V are fixed from fit to the data. The final broken HLS Lagrangian can be written : One has, here, included L′ tHoof t which provides determinant terms [54] breaking the nonet symmetry in the pseudoscalar sector and thus allowing an improved account of the π 0 , η, η ′ sector. L ′ HLS can be found expanded in the various Appendices of [13]. However, in order to account successfully for the largest possible set of data, isospin symmetry breakingà la BKY should be completed by a second step involving the kaon loop mixing of the neutral vector mesons (ρ 0 I , ω I and φ I ) outlined just below. This implies a change of fields to be performed in the L ′ HLS Lagrangian.

Mixing of Neutral Vector Mesons Through Kaon Loops
It has been shown [19,20,13] that L ′ HLS is insufficient in order to get a good simultaneous account of the e + e − → π + π − annihilation data and of the dipion spectrum measured in the τ ± → ν τ π ± π 0 decay. A consistent solution to this problem is provided by the vector field mixing mechanism first introduced in [18].
Basically, the vector field mixing is motivated by the one-loop corrections to the vector field squared mass matrix. These are generated by the following term of the broken HLS where g is the universal vector coupling and the subscript I indicates the ideal vector fields originally occurring in the Lagrangian. Therefore, the vector meson squared mass matrix M 2 0 , which is diagonal at tree level, undergoes corrections at one-loop. The perturbation matrix δM 2 (s) [18,19,20] depends on the square of the momentum flowing through the vector meson lines. The diagonal entries acquire self-mass corrections -noticeably the ρ 0 entry absorbs the pion loop -but non-diagonal entries are also generated which correspond to transitions among the ideal ρ 0 , ω and φ meson fields which originally enter the HLS Lagrangian 8 : Π ωφ (s), Π ρω (s) and Π ρφ (s). These are linear combinations of the kaon loops 9 . Denoting by resp. Π c (s) and Π n (s) the charged and the neutral kaon loops (including resp. the ρ 0 K + K − and ρ 0 K 0 K 0 coupling constants squared), one defines two combinations of these : 7 For clarity, we have dropped out the isospin breaking corrections generated by the BKY mechanism; the exact formula can be found in the Appendix A of [13]. The parameter z V corresponds to the breaking of the SU(3) symmetry in the Lagrangian piece L V , while z A is associated with the SU(3) breaking of L A . z V has no really intuitive value, while z A can be expressed in terms of the kaon and pion decay constants as z A = [f K /f π ] 2 . 8 These mixing functions were denoted resp. ε ωφ (s), ε ρω (s) and ε ρφ (s) in Section 6 of [13]. 9 Other contributions than kaon loops, like K * K loops, take place [18,19] which are essentially real in the energy region up to the φ meson mass. These can be considered as numerically absorbed by the subtraction polynomials of the kaon loops.
In term of ε 1 (s) and ε 2 (s), the transition amplitudes write : Therefore, at one-loop order, the ideal vector field V I = [ρ 0 I , ω I , φ I ] originally occurring in L ′ HLS are no longer mass eigenstates; the physical vector fields are then (re)defined as the eigenvectors of M 2 = M 2 0 + δM 2 (s). This change of fields should be propagated into the whole broken HLS Lagrangian L ′ HLS , extended in order to include the anomalous couplings [55] as done in [13]. In terms of the combinations (V R1 ) of the original vector fields V I which diagonalize L ′ HLS (see Section 5 in [13]), the physical vector fields -denoted V R -can be derived by inverting : where α(s), β(s) and γ(s) are the (s-dependent) vector mixing angles and s is the 4-momentum squared flowing through the corresponding vector meson line. These functions are proportional to the transition amplitudes reminded above. In contrast to ε 1 (s) which identically vanishes in the Isospin Symmetry limit, ε 2 (s) is always a (small) non-identically vanishing function. Therefore, within our breaking scheme, the ω − φ mixing is a natural feature following from loop corrections and not from IB effects. In contrast, the ρ − ω and ρ − φ mixings are pure effects of Isospin breaking in the pseudoscalar sector. For brevity, the Lagrangian L ′ HLS expressed in terms of the physical fields is referred to as BHLS.

The V ππ and V − γ/W ± Couplings
As the present study focuses on e + e − → π + π − data, it is worth to briefly remind a few relevant pieces of the L ′ HLS Lagrangian. In terms of physical vector fields, i.e. the eigenstates of M 2 = M 2 0 + δM 2 (s), the V ππ Lagrangian piece writes : where Σ V and (1 − h V )∆ V are isospin breaking parameters generated by the BKY mechanism [13], whereas α(s) and β(s) are (complex) "angles" already defined. Their expressions can be found in [13]. Eq. (9) shows how the IB decays ω/φ → π + π − appear in the BHLS Lagrangian.
Another Lagrangian piece relevant for the present update is : where g 2 is the weak SU(2) L gauge coupling and V ud is the element of the (u, d) entry in the CKM matrix. The f V γ (s) functions and f ρW are the transition amplitudes of the physical vector mesons to the photon and the W boson, respectively. At leading order in the breaking parameters, they are given by [13] : Eqs. (9) and (10) exhibit an important property which should be noted. The functions α(s) and β(s) providing the coupling of the physical ω and φ mesons to a pion pair also enter each of the f V γ (s) transition amplitudes, especially into f ργ (s). Therefore, any change in the conditions used in order to account for the decays ω/φ → π − π + correspondingly affects the whole description of the e + e − → π + π − cross section. Using the ω/φ → π + π − branching fractions in place of the π + π − spectrum in the corresponding regions has, of course, local consequences by affecting the corresponding invariant mass regions; it has also quite global consequences : indeed, it also affects the description of the annihilation cross-sections to π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 final states which all carry the f V γ (s) transition amplitudes.
Another effect, already noted in [13], is exhibited by Eqs. (11) : the ratio f ργ (s)/f ρW becomes s-dependent, which is an important difference between τ decays and e + e − annihilations absent from all previous studies, except for [7]. Figure 11 in [13] shows that the difference between f ργ (s) and f ρW is at the few percent level.

The Pion Form Factor
Here, we only remind the BHLS form of the pion form factor in τ decay and in e + e − annihilation and refer the interested reader to [13] for detailed information on the other channels. The pion form factor in the τ ± decay to π ± π 0 ν τ can be written : where a and g are the basic HLS parameters [15] already encountered; Σ V is one of the isospin breaking parameters introduced by the (extended) BKY breaking scheme. The other quantities are : where Π W (s) and Π ′ ρρ (s) are, respectively, the loop correction to the ρ ± − W ± transition amplitude and the charged ρ self-mass (see [13]).
The pion form factor in e + e − annihilation is more complicated and writes : where g ρππ , g ωππ and g φππ can be read off Eq. (9) and the F e V γ are given by : with the f V γ (s) given by Eqs. (11) above and the Π V γ (s) being loop corrections [13]. D ρ (s) = s−m 2 ρ −Π ρρ (s) is the inverse ρ 0 propagator while D ω (s) and D φ (s) are the modified fixed width Breit-Wigner functions defined in [13]; these have been chosen in order to cure the violation of F π (0) = 1 produced by the usual fixed width Breit-Wigner approximation formulae.

IB Distortions of the Dipion τ and e + e − Spectra
IB effects in the ππ channel are of various kinds. In the breaking model developed in [13] and outlined just above, the IB effects following from the neutral vector meson mixing, together or with the photon (see also [7]), are dynamically generated from the HLS Lagrangian. The most relevant effects have been reminded in Subsections 2.3 and 2.4 for the ππ channel. Indeed, Eq. (9) exhibits the generated coupling of the ω and φ mesons to a pion pair and Eqs. (10) and (11) show how the V − γ couplings are modified by the extended BKY breaking and vector mixing mechanisms. Therefore, in principle, all breaking effects 10 of concern for e + e − annihilations are exhausted.
Some IB effects affecting the dipion τ spectrum are also generated by the breaking mechanism, which modifies the W − ρ ± transition amplitude and the ρ ± π ∓ π 0 coupling. In some sense, the breaking mechanism decorrelates the universal coupling g as it occurs in the anomalous sector from those in the non-anomalous sector, where g appears in combinations reflecting IB effects, like g(1 + Σ V ) for the simplest form [13].
On the other hand, and as a general statement, the effects generated by the pion mass difference m π ± − m π 0 do not call for any specific IB treatment, as the appropriate pion masses are utilized at the corresponding places inside the model formulae derived from BHLS; this concerns, in particular, the pion 3-momentum which appears, for instance, in the phase space terms of the charged and neutral ρ widths.
However, there are IB breaking effects in τ decay, which have not yet been taken into account. Indeed, known distortions of the dipion τ spectrum relative to e + e − are produced by the radiative corrections due to photon emission. The long distance effects have been calculated in [56,57,58,59,60,61] and the short distance contributions in [62,63,64,65].
We have adopted the corresponding corrections, G EM (s) and S EW (= 1.0235±0.0003), respectively, as specified in [13]. In Ref. [59] the contribution of the sub-process τ → ωπ − ν τ (ω → π 0 γ) has been evaluated to substantially shift the correction G EM (s) (see Fig. 2 in [60]). This sub-process has been subtracted in the Belle data [23] which supposes that the corresponding correction has not to apply 11 . Hence, we applied to the three dipion τ spectra the correction as given in [56], as in our previous analysis [13].
These IB corrections distort the dipion spectra from the τ decay. They are accounted for by submitting to the global fit the experimental dipion τ distributions [21, 22, 23] using the HLS expression for dΓ ππ (s)/ds (see Eqs. (73) and (74) in [13] and the present Eq. (14)) correctedas usual -in the following way : where Γ τ is the full τ width and B ππ its branching fraction to π ± π 0 ν τ , both extracted from the RPP [24]. Indeed, as our fitting range is bounded by 1.05 GeV, both pieces of information are beyond the scope of our model. These corrections represent, by far, the most important corrections specific of the τ decay not accounted for within the HLS framework.
Another source of isospin breaking which may distort the τ spectrum compared to e + e − is due to the ρ mass difference δM ρ = m ρ ± − m ρ 0 . We note that the Cottingham formula, which provides a rather precise prediction of the m π ± − m π 0 electromagnetic mass difference, predicts for the ρ an electromagnetic mass difference : In principle, within the HLS model a ρ 0 − ρ ± mass shift is also generated by the Higgs-Kibble mechanism (corresponding to the well known shift in the Electroweak Standard Model). In the HLS model m 2 ρ ± = a f 2 π g 2 while m 2 ρ 0 = a f 2 π (g 2 + e 2 ) due to ρ 0 − γ mixing. This leads to a Higgs-Kibble shift of about m ρ 0 − m ρ ± ≃ e 2 2g √ a f π ∼ 1 MeV (see [15]), which essentially compensates the electromagnetic shift obtained from the Cottingham formula. In addition, the masses are subject to modifications by further ρ − ω − φ mixing effects, obtained from diagonalizing the mass matrix after including self-energy effects. The mentioned effects have been estimated in [57] and lead to : When evaluating the anomalous magnetic moment from τ data, several choices have been made; for instance, the analysis in [8] assumes δM ρ ≃ 1.0 ± 0.9 MeV, while Belle [23] preferred δM ρ ≃ 0.0 ± 1.0 MeV.
In our study, we followed Belle and have adopted δM ρ ≃ 0 ± 1 MeV, consistent with the estimate by Bijnens and Gosdzinsky just reminded and with most experimental values reported in the RPP 12 [24]. As noted elsewhere [20,13], based on the available data, auxiliary HLS fits do not improve by letting δM ρ floating. This justified reducing the model freedom by fixing this additional parameter to 0.
Yet another source of isospin breaking which may somewhat distort the τ dipion spectrum compared to its e + e − partner is the width difference δΓ ρ = Γ ρ 0 − Γ ρ ± between the charged and the neutral ρ; however, this can be expected to be small as the accepted average [24], δΓ ρ = (0.3 ± 1.3) MeV, is consistent with 0.
Usually, the evaluation of the (δM ρ , δΓ ρ ) effects is performed using the Gounaris-Sakurai (GS) parametrization [66] of the pion form factor. However, the GS formula does not parametrize the radiative corrections expected to affect the measurement of the pion form factor. Therefore, the correction for the radiative width may not be well taken into account by just shifting the width in the GS formula. In Refs. [61,8] an effective shift of δΓ γ ρ ≃ 1.82 MeV has been estimated by subtracting ρ + → π + π 0 γ in the τ channel and adding ρ 0 → π + π − γ in the e + e − channel. The question is how this affects |F π (s)| 2 . Usually one adopts the GS formula to parametrize the undressed data, which is not precisely what is measured. If one assumes the GS formula to represent the dressed data as well, one may just modify the width for undressing the τ spectrum and redressing the radiative effects in the e + e − channel, as an IB correction.
An increase of the width in the GS formula has two effects. One is to broaden the ρ shape, which results in an increase of the cross section. The second, working in the opposite direction, is to lower the peak cross section. In the standard form of the GS formula (see e.g. CMD-2 [67] or Belle [23]) the second effect wins and one gets a substantial reduction of the muon g − 2 integral by δa had,LO µ [ππ, τ ] = (−5.91 ± 0.59) × 10 −10 [8], a large reduction of the τ -based evaluation. Looking at the Breit-Wigner peak cross section given by it is not a priori clear, which of the different widths are affected. If one keeps fixed the branching fractions for ρ → ee and ρ → ππ, the peak cross section would not change at all. Therefore, the correction for radiative events via the GS parametrization is not unambiguous. In the standard GS parametrization Γ GS ee = ) is a derived quantity and depends on Γ ρ in an unexpected way. We therefore consider this standard procedure of correcting for radiative decays as not well established.
Auxiliary fits allowing a difference δg between the ρ ± and ρ 0 couplings, in order to generate a floating δΓ ρ , have been performed. One observed slightly more sensitivity to a free δg than to a free δM ρ , but nothing conclusive enough to depart from δg = 0 while increasing the number of fit parameters and their correlations. Indeed, within the BHLS framework, the τ data only play the role of an additional constraint and their use is certainly not mandatory, except for testing the "e + e − vs τ discrepancy" which has been shown to disappear [13].
Within the set of data samples which are studied by means of the global fit framework provided by BHLS, the single place where the charged ρ meson plays a noticeable role is the τ dipion spectrum. Taking into account its relatively small statistical weight within this set of data samples, one does not expect to exhibit from global fits a noticeable sensitivity to mass and width differences with its neutral partner.

The Issue
Although the BHLS Lagrangian should be able to describe more complicated hadron production processes, in a first step one obviously has to focus on low multiplicity states, primarily two particle production but also the simplest three particle production channel e + e − → π + π − π 0 . Four pion production, annihilation to KKπ . . . are beyond the scope of the basic setup of the BHLS model. We expect that available data on the lowest multiplicity channels provide a consistent database which allows us to pin down all relevant parameters, such that our effective resonance Lagrangian is able to simultaneously fit all possible low multiplicity channels. In fact, what is considered are essentially all relevant annihilation channels up to the φ; in this energy range the missing channels (4π, 5π, 6π, ηππ, ωπ) contribute less than 0.3% to a had µ [13]. Our previous study [13] has actually shown that the following groups of complementary data samples and/or RPP [24] accepted particle properties (mainly complementary branching fractions) support our global fit strategy : • (i) All e + e − annihilation data into the π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 final states admit a consistent simultaneous fit 13 , • (iii) Some additional partial width from the P γγ and V P γ decays, which are independent of the annihilation channels listed just above, • (iv) Some information concerning the φ → π − π + decay, especially its accepted partial width Γ(φ → π − π + ) [24]. This piece of information is supposed to partly counterbalance the lack of spectrum for the e + e − → π + π − annihilation in the φ mass region 15 .
• (v) All the e + e − → π + π − data sets (NSK) collected 16  They represent a complete reference collection of data samples and lead to fits which do not exhibit any visible tension between the BHLS model parametrization and the data (see for instance Table 3 in [13]). This is worth being noted, as we are dealing with a large number of different data sets collected by different groups using different detectors and different accelerators. The (statistical & systematic) error covariance matrices used within our fit procedure are cautiously constructed following closely the group claims and recommendations 17 . Therefore, the study in [13] leads to think that the model correlations exhibited by BHLS reflect reasonably well the physics correlations expected to exist between the various channels.
However, beside the (NSK) e + e − → π + π − data sets collected in scan mode, there exists now data sets collected using the Initial State Radiation (ISR) method by the KLOE and BaBar Collaborations. All recent studies (see [30,8], for instance) report upon some "tension" between them. As this issue has important consequences concerning the estimate of the muon anomalous magnetic moment, it is worth examining if the origin of this tension can be identified and, possibly, substantiated. Besides scan and ISR data, it is also interesting to reexamine [18] the pion form factor data collected in the spacelike region [35,36] within the BHLS framework; indeed, if valid, these data provide strong constraints on the threshold behavior of the pion form factor and, therefore, an improved information on the muon g − 2.

The Analysis Method
The BHLS model has many parameters and a global fit has to be guided by fitting those parameters to those channels to which they are the most sensitive. Obviously resonance parameters of a given resonance have to be derived from a fit of the corresponding invariant mass region. Similarly, the anomalous type interaction responsible for π 0 → γγ or the πππ final state . . . are sensitive to very specific channels only.
We also have to distinguish the gross features of the HLS model and the chiral symmetry breaking imposed to it. With this in mind, in our approach to comparing the various e + e − → π + π − data samples, the τ decay spectra play a key role since the charged channel is much simpler than the neutral one where γ, ρ 0 , ω and φ are entangled by substantial mixing of the amplitudes, which are not directly observable. In the low energy region, below the kaon pair thresholds and the φ region, what comes into play is the ρ ± form factor obtained from the τ spectra. Together with the isospin breaking due to ρ 0 − ω mixing -characterized by the branching fractions Br(ω → π + π − ) and Br(ω → e + e − ), which in a first step can be taken from the RPP -the ρ ± form factor should provide a good prediction for the e + e − → π + π − channel. Data from the latter can then be used to refine the global fit. This will be our strategy in the following.
The annihilation channels referred to as (i) in the above Subsection as well as the decay information listed in (iii) have little to do with the e + e − → π + π − annihilation channels, except for the physics correlations implied by the BHLS model. On the other hand, as long as one limits oneself to the region (2m π , 1 GeV), there is no noticeable contradiction between the various dipion spectra extracted from the τ ± → ν τ π ± π 0 decay by the various groups [21, 22, 23]. Therefore, it is motivated to examine the behavior of each of the collected e + e − → π + π − data sets, independently of each other, while keeping as common reference the data corresponding to the channels listed in (i) -(iii). Stated otherwise, the data for the channels listed in (i) -(iii), together with the BHLS model, represent a benchmark, able to examine critically any given e + e − → π + π − data sample.
It then only remains to account for isospin breaking effects specific of the e + e − → π + π − channel, in a clearly identified way. A priori, IB effects specific of the e + e − → π + π − annihilation are threefold and cover : • (jj) Information on the decay ω → π + π − , The importance of decay information on ρ 0 → e + e − to determine IB effects has been emphasized in only a few previous works [7,20,13]. Within the BHLS model, the ratio f ργ (s)/f ρW exhibits non-negligible IB effects for this particular coupling (see Fig. 11 in [13]). They amount to several percents in the threshold region quite important for evaluating g − 2.
There is certainly no piece of information in the data covered by the channels listed in (i) -(iii) above concerning the decay information (jj) or (jjj). In contrast, the vertex ρ 0 e + e − is certainly involved in all the annihilation channels considered. Imposing the RPP [24] information Γ(ρ 0 → e + e − ) = 7.04 ± 0.06 keV is, nevertheless, legitimate because the channels (i) -(iii) do not significantly constrain the decay width ρ 0 → e + e − .
For the following discussion we define the branching ratio products F ω .
Two alternative analysis strategies can be followed : • (k) Use the accepted values [24] for the ρ 0 → e + e − , ω → π + π − and φ → π + π − . These are the least experiment dependent pieces of information 18 . We will be even more constraining by supplementing these ω and φ branching ratios by phase information : The so-called Orsay phase concerning the ω decay 19 and the reported phase 20 of the φ → π + π − amplitude relative to ρ 0 → π + π − .
• (kk) Use directly data when possible. Indeed, all relevant IB information carried by ρ 0 → e + e − and ω → π + π − can be numerically derived within the BHLS model by the difference between the π + π − spectra and the dipion spectrum π ± π 0 in the τ decay; more precisely, using the π + π − spectrum within the tiny energy region 0.76 ÷ 0.82 GeV should be enough to derive the relevant IB pieces of information in full consistency with our model. 18 Nevertheless, one should keep in mind that these accepted values are highly influenced by the e + e − → π + π − scan data samples compared to others. Therefore, this choice could favor the CMD-2 and SND data samples when fitting; however, as these accepted values are certainly not influenced by none of the BaBar or KLOE data samples, the behavior of each of the various ISR data samples becomes a crucial piece of information. 19 We will use as input the value 104.7 • ± 4.1 • found by [70], which is consistent with the results recently derived [71] while using an analogous (HLS) model (see Tables VI-IX therein). 20 The single existing measurement −34 • ± 4 • is reported by the SND Collaboration [72].
As all scan (NSK) e + e − → π + π − data samples [25, 26, 27, 28, 29] and both KLOE data sets (KLOE08 and KLOE10) stop below 1 GeV, the φ information should be taken from somewhere else, namely from the RPP. Fortunately, the φ region is now covered by the BaBar data set [33,34]. Therefore, as soon as the consistency of the φ → π + π − information carried by the BaBar data and by [24,72] is established, this part of the spectrum could supplement the scan and KLOE data sets 21 .
Concerning the φ → π + π − information, this second strategy will be used by either taking the [24, 72] information or the BaBar data points between 1.0 and 1.05 GeV.

How to Implement
The vector meson couplings to π + π − or e + e − depend on the s-dependent "mixing angles" α(s), β(s) and γ(s). This does not give rise to any ambiguity as long as one deals with spectra; however, when using the PDG information for vector meson decays, especially to π + π − or e + e − , one has to specify at which value for s each of the vector meson (model) coupling should be evaluated.
Within the HLS model, there are a priori two legitimate choices for the mass of vector mesons; this can either be the Higgs-Kibble (HK) mass which occurs in the Lagrangian after symmetry breaking or, especially for the ω and φ mesons, the experimental (accepted) mass as given in the RPP. Prior to the availability of the BaBar data [33], the published e + e − → π + π − cross section data did not include the φ mass region and, therefore, there was no criterion to check the quality of each possible choice in the φ mass region 22 . The choice made in the previous studies using the broken HLS model [18,19,20,13] was the φ HK mass.
As already noted, the broken HLS model, fed with the data listed in i-iv (see Subsection 3.1, above), provides predictions for the pion form factor independently of the measured e + e − → π + π − data. This procedure is discussed in detail in the next section. Here we anticipate some results specific to the φ mass issue. Fig. (1a) displays the prediction for the pion form factor in the φ region using the HK mass to estimate the φπ − π + coupling constant with the BaBar data superimposed (not fitted); it is clear that the prediction is quite reasonable up to ≃ 0.98 GeV as well as above ≃ 1.05 GeV. However, it is clearly unacceptable for the mass region in-between. In contrast, using the φ mass as given in the RPP to extract the φπ − π + coupling constant from its accepted partial width [24] provides the spectrum shown in Fig. (1b); this alternate choice is certainly reasonable all along the mass region displayed. Therefore, it is motivated to update our former results [13] by performing the change just emphasized 23 . In order to be complete, it is worth mentioning here a fit result obtained by exchanging the PDG/SND φ decay information with the BaBar pion form factor data 21 in the range (1.0 ÷ 1.05) GeV. The result, given in Fig. (1c), shows that the lineshape of the BaBar pion form factor at the φ mass can be satisfactorily accommodated. Figure 1: The e + e − → π + π − cross section around the φ mass together with BaBar data superimposed. The curve in (a) displays the prediction using the RPP φ decay information computed at the φ Higgs-Kibble mass; the curve in (b) displays the prediction using the PDG φ decay information computed at the experimental φ mass. In (c) the PDG φ decay information is replaced by the five BaBar data points located between 1. and 1.05 GeV.
As the exact pole position of the φ meson is determined by a benchmark independent of the e + e − → π − π + process (see Subsection 3), the drop exhibited by Fig. (1c) in the BaBar is perfectly consistent with an expected φ signal.

τ +PDG Predictions
As mentioned before, the charged isovector τ ± → ν τ π ± π 0 dipion spectra are not affected by γ − ρ − ω − φ mixing and hence are of much simpler structure. Supplemented by the basic ρ − ω − φ mixing effects which derive from SU(2) and SU(3) flavor breaking, one has a good starting point to fix the parameters of the BHLS model to predict the process e + e − → π + π − . Specifically, we are using the data including the channels listed in i-iv of Section 3 together with RPP information relevant to fix the IB effects affecting the pion form factor. This method is named, somewhat abusively 24 τ +PDG.
Specifically, the IB effects encoded in are taken from the RPP. For the missing phase information we adopt the result from the fit [70] for the Orsay phase of the ω → π + π − amplitude and the result from SND [72] for the phase of the φ → π + π − amplitude 25 . Following the discussion in the preceding Subsection, the model branching ratios and phases are computed at the vector boson masses accepted by the RPP.
The fit returns a probability of 89.4% with χ 2 /n dof = 553.4/596. The fit quality (χ 2 /n points ) for each of the fitted channels is almost identical to our results in [13] (see the last column in Table 3 therein). Each of the decay partial width extracted from [24] contributes by ≃ 1 to the total χ 2 . It is also worth mentioning that the dipion spectra from [21, 22, 23] are nicely described up to √ s = 1 GeV and provide residual distributions indistinguishable from those shown in Figure 10 of [13]. From this fit, one derives the (τ +PDG) predictions for the pion form factor which can be compared with the various existing e + e − → π + π − data samples.
The overall view of the comparison is shown in Fig. 2. This clearly indicates that the data associated with the channels listed in i-iv, supplemented with a limited PDG information is indeed able to provide already a satisfactory picture of the pion form factor as reported by all experiments having published e + e − → π + π − spectra.
Let us stress that the predicted pion form factor relies on the π ± π 0 spectra provided by ALEPH [21], Belle [23] and CLEO [22] only up to 1.0 GeV. Therefore, the inset in Fig. 2 actually shows the extrapolation of the prediction into the spacelike region with the NA7 data [35] superimposed; this clearly indicates that there is no a priori reason to discard the spacelike data from our data handling. One should also note that the extrapolation of the prediction above the φ mass is quite reasonable up to ≃ 1.2 ÷ 1.3 GeV. This may indicate that the influence of high mass vector mesons is negligible up to this energy region.
In order to make more precise statements, let us magnify piece wise the information carried by Fig. 2. Thus, Fig. 3 displays the behavior of the various e + e − → π + π − data samples in the (0.3 ÷ 0.7) GeV energy region. As a general statement the behavior expected from the existing data samples looks well predicted by the τ +PDG method. A closer inspection allows to infer that the CMD-2 and SND data points (i.e. NSK when used together) are well spread onto 24 By abusively, we mean, first that the "Orsay" phases for both the ω and φ mesons have no entry in the RPP and, second, that the benchmark represented by the processes listed in Subsection 3 within items (i) to (iv) have little to do with τ or the RPP. 25 A preliminary version of the present work was presented [73] at the Workshop on Meson Transition Form Factors held on May 29-30, 2012 in Krakow, Poland. Some minor differences may occur with the present results due to the fact that the SND phase for φ → π + π − was not imposed in the preliminary work. both sides of the predicted curve; this property is also shared by the KLOE10 sample. Even if reasonably well described, the KLOE08 and BaBar data samples are lying slightly above the τ +PDG expectations; this difference should vanish when including the π + π − spectra inside the fit procedure. Fig. 4 displays the behavior of the various e + e − → π + π − data samples in the (0.85 ÷ 1.2) GeV energy region. Here also the predicted curve accounts well for the data behavior. A closer inspection tells that the sparse NSK data are well described. The BaBar data are also well accounted for all along this energy interval except for the φ region. As shown by Fig. (1c) above, this can be cured and one can show that the difference is mostly due to the phase for the φ → π + π − amplitude which departs significantly 26 from those provided by SND [72]. One could also note that both KLOE data samples look slightly below the τ +PDG expectations in this region.
One may conclude from Figs. 3 and 4 that our "τ +PDG" predictions are in good agreement with the data and that a fit using fully these data samples should provide marginal differences between all π + π − data sets 27 .
However, the picture becomes quite different in the medium energy region (0.70 ÷ 0.85) GeV as illustrated by Figure 5. In this region, our τ +PDG prediction follows almost perfectly expectations 27 from both the KLOE08 and KLOE10 data and the detailed lineshape of the ω − ρ interference region is strikingly reproduced. Paradoxically, the NSK data are slightly less favored -especially around 0.8 GeV -despite their influence on the PDG information used in order to account for IB effects in the ω − ρ interference region; however, taking into account experimental uncertainties, we already know that a global fit using NSK data is highly successful [13].
In contrast, the behavior of the BaBar data looks inconsistent with the τ +PDG prediction, especially on the low mass side of the interference region. Actually, the observed overestimate of the BaBar spectrum affects the whole region from threshold to the ω mass but is more important in the range (0.74 ÷0.78) GeV. At higher energies one observes a reasonable agreement with expectations as well as with both KLOE data sets.
One should also note that the ω mass and (total) width induced by the data for the π 0 γ, ηγ, π + π − π 0 final states are in perfect agreement with all the examined data samples; this indicates that the energy calibration around the ω mass is good for all ISR data samples. Fig. (1c) has already shown that the BaBar energy calibration is also good in the φ region.

Predictions Using the ρ 0 − ω Interference Region From Data
As stated in Subsection 3.1, one can replace the PDG information for ω → π + π − and ρ 0 → e + e − by any pion form factor spectrum limited to the region (0.75 ÷ 0.82) GeV. In particular, this turns out to fit F ω = Br(ω → e + e − ) × Br(ω → π + π − ) and the Orsay phase as it comes out from each of the specified data sets.
Concerning the φ → π + π − mode, one can also check the BaBar data [33, 34] versus the SND datum [72]. This will be discussed below when reporting on fitting the whole spectra (see Subsection 5.2.1).
Using the RPP recommended value for Br(ω → π + π − ) and the Orsay phase information from [70] yields a value for F ω in good correspondence with expectations, as clear from the entry τ + PDG in Table 1. Using NSK data or any of both KLOE samples, instead of the PDG information, does not lead to predicted curves substantially different from their analogue already shown and commented upon in the previous Subsection. Interesting parameter values have been extracted from global fits using only the (0.75 ÷0.82) GeV region from the NSK and KLOE spectra for e + e − → π + π − and are reported in Table 1; they are in reasonable agreement with the reported branching ratio product F ω [29] and the Orsay phase [70,71] as well. Indeed, taking the RPP F ω value as reference, our estimates using the 70 MeV interval surrounding the interference region are at 0.1σ, 1.8σ and 3.0σ for respectively the NSK [27, 29], KLOE08 [31] and KLOE10 [32] data samples. The difference between the RPP recommended value for F ω and our entry for NSK also tells that the BHLS parametrization and the more standard form factor lineshape used by SND [29] provide almost identical values for F ω .
As far as BaBar data are concerned, the situation looks different and the most relevant piece of information is provided in Figure 6. This proves that the largest difference between BaBar data and the other analogous data samples [27, 29, 31, 32] is the F ω information inherent to the BaBar data. Clearly, Table 1 shows that the BaBar value for F ω is off from its recommended value by (7 ÷ 8)σ. This strong disagreement is substantiated by comparing Figures 5 and 6. Table 1 also reports the fit probability for each of the examined configurations. With about 90% probabilities, the "τ +PDG" prediction and the NSK, KLOE08 and KLOE10 (global) fits exhibit a full consistency with the rest of our benchmark (i.e. all other annihilation channel physics). The agreement of this with BaBar data, even limited to such a tiny interval, is found much poorer and exhibits a clear tension between F BaBar ω and the rest of the (non-π + π − ) physics accessible to the HLS model.
Additional pieces of information are provided in the last data column of Table 1 which complements the global fit probabilities. These are the values for χ 2 /n π + π − for each of the various π + π − data samples, n π + π − being the number of data points included in the fitted energy range (i.e. (0.75 ÷ 0.82) GeV).

Isospin Breaking Effects in the BHLS Model : Comments
It follows from the developments just above that the BHLS model fed with a limited number of accepted values for some IB pieces of information is indeed able to provide a quite satisfactory prediction for the e + e − → π + π − cross section once the τ spectra are considered. This gives support to our breaking model, especially to the s-dependent vector meson mixing mechanism.
The prediction is found in accord with the scan (NSK) data samples and with both KLOE data sets 28 . Indeed, the predicted lineshape strikingly follows the central values from both KLOE data samples; for the scan data, the prediction based on PDG information is good but not as good as for the KLOE data. However, changing the PDG requested IB information by less than 1 σ -as following from a mere comparison of the first and second lines in Table 1 leads to a perfect description of the NSK spectra over the whole available energy range [13]. In contrast, the RPP branching fraction product F ω has to be changed by about 7σ in order to yield a comparable description of the BaBar [33] data.
Basically, our approach is a τ based prediction of π + π − spectra; it relies on the consistency of several different physics channels, the τ spectra and on a model of isospin symmetry breaking (IB). In fine, our breaking model does not carry IB parameter values plugged in from start, but yields the numerical IB effects in a data driven mode. It is thus interesting to examine the consequences of this τ (and global) based approach on the muon g − 2 estimated value. For this purpose, it is worth stressing that the τ based estimates given just below -and later -are not computed by integrating the τ spectra (and adding corrections like in [23,8], for instance), but by integrating the e + e − cross sections they allow to reconstruct through the global BHLS framework.
This is what is shown in Figure 7. The first line displays the result derived from the "τ +PDG" fit. The four following lines are obtained by replacing the ρ and ω IB information by the 0.76 ÷ 0.82 GeV region of the quoted data sets. The line named BNL displays the experimental result [39,40] and the last line shows the τ based estimate from [8]. The last data column in Figure 7 gives the probability of the corresponding fits.
It is clear that all methods used to include IB effects within our τ based approach give consistent results, all distant from the BNL measurement at the 4σ level. The associated probabilities indicate the quality of the fits from where they are derived.

Global Fits
Using the e + e − → π + π − Spectra As in our previous analysis [13], we have performed global fits using simultaneously all e + e − annihilation data into the π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 final states, the dipion spectra collected in the decay of the τ lepton [21, 23, 22] and the decay information listed in Subsection 3.1. We also use the RPP φ decay properties in the (updated) way emphasized in Subsection 3.

3.
For what concerns the e + e − → π + π − data included into the global fit procedure, we have performed fits using separately the NSK, KLOE08 and KLOE10 data samples. Global fits have also been performed for the BaBar dipion spectrum restricted to the range of validity of our BHLS model. Two options have been considered, using BaBar data up to 1 GeV supplemented by φ → π + π − decay properties from the RPP or using the BaBar data up to 1.05 GeV, thus including its φ region and avoiding the need of using the RPP information about the φ.
In all cases, the errors (and the χ 2 ) were constructed following the information published/ recommended by the experimental groups who collected these data. For instance, concerning the scan data samples, the full covariance matrix is generally constructed by adding the systematic error covariance matrix V syst to the (diagonal) statistical covariance V stat constructed from the tabulated uncertainties as described in [19] (see Section 6 therein); V syst is constructed assuming the reported systematic errors bin-to-bin correlated 29 . Nevertheless, for the rather imprecise e + e − → (π 0 /η)γ data, taking into account the large magnitude of the reported systematics, the systematic and statistical errors were simply added in quadrature. For the e + e − → π + π − π 0 data, we dealt with depending on the level of precision of the relevant data sample (see Subsection 2.2.3 in [20] for more information). 29 When they are reported, uncorrelated systematics are simply added in quadrature to the statistical errors.
For the τ data, as both V syst and V stat are publicly available [21,22,23], one only has to add them up, as already performed for the study in [13]. For the BaBar sample, the systematic uncertainties on the cross section are given (as a function of √ s) in Table V from [34] and are imposed to be 100 % bin-to-bin correlated (following Section F in [34]). Writing, for definiteness, each of these uncertainty functions as f (s i ), for the data point m i located at √ s i , For both the KLOE08 and KLOE10 data samples [31,32], the KLOE Collaboration provides basically the same information as BaBar and therefore, we simply have to proceed likewise.
Moreover, the various π + π − data sets collected by CMD-2 and SND [27, 28, 29] on the one hand and the older data samples from [25] on the other hand, both carry common bin-tobin and sample-to-sample correlated uncertainties estimated resp. to 0.4% and ≃ 1%. As in our previous studies [20,13], this effect is accounted for in the minimization code. These are the most important reported correlations of this type within the scan data [74].
We have also performed global fits using combinations of these individual π + π − data samples. In this case, the contributions of the NSK and KLOE10 data to the total χ 2 were left unweighted as their own χ 2 /n contribution is always of the order 1 in fits using each of them in "isolation" 30 . In contrast, in such combinations involving the KLOE08 and/or BaBar data, the contribution of each of these to the total χ 2 was weighted by the ratio M is the χ 2 of the M data set obtained in the best fit using only M as e + e − → π + π − data set; n M is the corresponding number of data points. In fits involving the spacelike data [35,36], the corresponding weight was also used 31 .
For definiteness, when relevant, we have used f N SK = f KLOE10 = 1, f KLOE08 = 60/90 ≃ 0.67, f BaBar = 270/346 ≃ 0.78 and f space = 59/85 ≃ 0.69. These weights have been varied and it has been found that the sensitivity of the physics results to their precise value is marginal; the main virtue of these weights is to provide probabilities not too much ridiculous. On the other hand, as a matter of principle, when results are displayed which have been obtained using weights, it is quite generally for the reader's information. We have preferred being conclusive by only relying on the largest data set combinations which do not call for any reweighting. Indeed, this simply reflects that global fit probabilities do not raise any objection to trusting the uncertainties as they are reported together with each of the used data sets.
This method 32 , turns out to consider each data set as a global object, rather than defining local (s-dependent) averages as done by others [6]. This method looks better adapted to the global fit method which provides a quality check reflecting the behavior of each e + e − → π + π − data set within the global context of a large number of physics channels. Indeed, doing local averages would prevent to detect discrepancies originating from some given data set only. On 30 We remind that, in this work, the wording "isolation" or "standalone" is always used to qualify the fits performed using a specific e + e − → π + π − data sample. It is understood that each of the π + π − data sets "fitted in isolation" is always fitted in conjunction with all the other data listed in Subsection 3.1. When several of the e + e − → π + π − data samples are submitted to fit (together with the other channels), we then refer to "combined" fits. In order to warn the reader, we prefer keeping the word isolation between quote marks. 31 In the fits referred to in [73], the spacelike data contributions to the total χ 2 were left unweighted. 32 This method is quite parent from the S-factor technique commonly used in the Review of Particle Properties to account, while averaging, for marginal inconsistencies between the various reported measurement/uncertainty of some physics quantity. the other hand, within a global framework as BHLS, such a method would lead to fit parameter values modified in a completely uncontrolled way. It is the reason why our final results will only rely on sample combinations which do not call for any reweighting (i.e. necessarily going beyond experimentally reported uncertainty information).
For completeness, it is also worth noting that the χ 2 contributions of the -more than 40data sets associated with all the other channels were always left unweighted.
A feature common to all fits using the e + e − → π + π − data sets in "isolation" or combined is that the individual χ 2 contributions associated with the other channels (π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 , . . . ) were only marginally affected by the specific choice of e + e − → π + π − data submitted to fit. Their typical values are almost identical to what can be found in the last data column of Table 3 in [13]; more precisely, the χ 2 value provided by each of these channels never varies by more than a few percents. Let us remind that the number of data points submitted to fit -beside the e + e − → π + π − data -is ≃ 600 when working within 33 the configuration B defined in [13] (≃ 675 within configuration A).

Salient Features of the Various
We have performed several tens of fits of the various e + e − → π + π − spectra in "standalone" mode 30 and/or combined. It does not look useful to report on each fit in detail. Instead of overwhelming the reader with unnecessary information and plots, we have preferred focusing on the salient features of their behavior within the global fit context. Beside the fit properties of the full spectra (up to 1 GeV, generally), this covers the muon anomalous magnetic moment value and the behavior at the ω mass, more precisely the value for F ω = Br(ω → e + e − ) × Br(ω → π + π − ). The first of these topics will be addressed in a separate Section below. Concerning the second topic, we remind that the value for F ω expected from the RPP is (1.225± 0.071) 10 −6 .
On the other hand, we do not enter into much detail concerning the effects of the spacelike data, always used weighted in this paper; we limit ourselves to mentioning that they never modify the fit qualities in a significant way.

"Standalone"
Fits of the e + e − → π + π − Spectra When using only 30 the (unweighted) NSK data, the fit returns 34 χ 2 N SK /N π + π − = 128.30/127 and a global fit probability of 96.3%. The fit residuals are displayed in the top panel of Fig.  8 and the fractional deviations from the fitting function in the bottom panel. Taking into account the uncertainties, both plots do not exhibit any significant departure from flatness. The fit residual distributions of the ALEPH, CLEO and Belle dipion spectra also submitted to fit in the present case are rigorously given by Fig. 10 in [13], where the fit also extends over the energy region from threshold to 1 GeV. Therefore, one indeed observes flat residual distributions and flat fractional deviations from the fitting functions simultaneously in the e + e − and τ channels as already claimed in [13].
Similarly, the KLOE10 data set [32] returns χ 2 KLOE10 /N π + π − = 73.68/75 and a global fit probability of 87.7%. The fit residual distribution is shown in the top panel of Figure 9; the bottom panel in this Figure shows the fractional deviations from the fitting function. Both distributions can be considered as reasonably flat. As in the above Fig. 8, the plotted errors are the square roots of the diagonal elements of the (full) error covariance matrix. Therefore, the value for χ 2 KLOE10 and the flatness of the residuals shown in Figure 9 illustrate that the (reported in [32]) full error covariance matrix looks correctly understood.
Therefore, within the global fit context, each of the NSK and KLOE10 data samples exhibits the same outstanding behavior 35 . On the other hand, the fit parameters and error covariance matrix allows to derive, using obvious notations, F N SK ω = (1.205 ± 0.042) 10 −6 close to the PDG value reported at top of Table 1. One also gets F KLOE10 ω = (1.074 ± 0.051) 10 −6 , a ≃ 2.5σ difference with the central value for F N SK ω .
Under the same conditions, the (unweighted) KLOE08 data set returns χ 2 KLOE08 /N π + π − = 96.55/60 and a global fit probability of 56%, which is rather low. It is worth noting the remarkable flatness of the residual distributions displayed by Figure 10 which does not prevent to yield a relatively large value for χ 2 KLOE08 . The large χ 2 KLOE08 and the flatness of the residual distribution shown in Figure 10, considered together, might indicate an issue with the non-diagonal part of the full error covariance matrix. Anyway, one can conclude that the poor KLOE08 fit probability reflects an issue with the KLOE08 error estimate rather than a distorted lineshape. Figure 10 should be compared with the similar distribution derived formerly using a primitive version of the BHLS model (see Figure 3 in [19]). The clear improvement substantiates the gain provided by the BHLS model in its present form. On the other hand, one should also note that the value for F KLOE08 ω = (1.117 ± 0.042) 10 −6 is consistent with F KLOE10 ω . The fit of the (unweighted) BaBar data leads to χ 2 BaBar /N π + π − = 343.08/270 = 1.27 when limited to 1 GeV (17% probability) and to χ 2 BaBar /N π + π − = 340.77/275 = 1.24 when going up to 1.05 GeV (corresponding to a 22% probability). The residual distribution yielded when fitting from threshold to 1.05 GeV is given in Figure 11 together with the fractional distribution. The residual distribution derived when fitting from threshold to 1 GeV, (not shown) is slightly flatter in the ρ − ω region, which indicates that accomodating simultaneously the BaBar ω and φ regions has some price. This has some consequence on the BaBar estimate of the muon anomalous moment, as will be emphasized in Subsection 6.1 below.
Nevertheless, these distributions look reasonable and are associated with quite reasonable values for χ 2 BaBar /N π + π − ; the bottom panel distribution in Figure 11 is even quite similar to those derived by the BaBar Collaboration fit shown in Fig. 47 of [34] (with no quoted fit quality).
Therefore, the real issue is not the description of the BaBar spectrum stricto sensu, but its consistency with the rest of the data and physics involved in the global fit, especially the τ data. This is well reflected by the poor global fit probabilities of 17% or 22%, poorer than for the KLOE08 data sample. Qualitatively, this result could have been expected from comparing the τ +PDG predictions with "τ +(ρ − ω)" (using this region of the BaBar spectrum instead of the PDG information) already analyzed in Subsection 4.2.

Fits
Combining the e + e − → π + π − Spectra In view of the "standalone" fits just reported, we have done fits of different combinations of the scan (NSK) and ISR data samples. In this case, the contributions of the KLOE08 and BaBar samples to the (minimized) χ 2 function are always weighted as already stated.
• Combining the KLOE08 and KLOE10 data : This returns a consistent picture where χ 2 KLOE10 and χ 2 KLOE08 are almost unchanged compared to their "standalone" values and the fit probability reaches 81.6 %. This, indeed, confirms that they share the same physics content. This is confirmed by the fit result F KLOE08/10 ω = (1.121 ± 0.038) 10 −6 , consistent with both of F KLOE08 ω and F KLOE10 ω .
• Combining all the ISR data sets (KLOE08, KLOE10 and BaBar) : This returns χ 2 ISR /N π + π − = 1.34 when including the weights for KLOE08 and BaBar and χ 2 ISR /N π + π − = 1.64 when the weights are not included. The weighting used does not prevent the probability to remain poor 35 (1.8%) reflecting the level of inconsistency of the BaBar and KLOE(08/10) data samples already noted. In this case, one gets F ISR ω = (1.608 ± 0.010) 10 −6 , exhibiting a large distortion towards the BaBar lineshape, despite the weighting.
• Combining all e + e − → π + π − spectra : Taking the weights into account -which is more favorable -one gets χ 2 Global /N π + π − = 1.35 and a fit probability 35 of 1.3%. Once again, the lineshape of the fit function is highly influenced by the BaBar sample in the ρ − ω interference region as shown by F Global ω = (1.582 ± 0.089) 10 −6 .
• Combining the NSK and KLOE10 data : In this case, there is no weight and one gets χ 2 N SK /N π + π − = 131.37/127 and χ 2 KLOE10 /N π + π − = 72.90/75 close to the "standalone" 30 fit results reported in Subsection 5.1.1 above, and thus, χ 2 KLOE10+N SK /N π + π − = 1.01, and the remarkable global fit probability of 96.9 %, showing that the NSK and KLOE10 data samples are quite consistent with each other and with the rest of the BHLS physics as well. Figures 12 and 13 display the fit residuals of this common fit; in both Figures, the leftmost panels show the usual residual distributions of the NSK and KLOE10 samples (i.e. the differences between each measurement and the corresponding fitting function value), while the rightmost panels display the same information corrected for the bin-to-bin correlated uncertainties (see the Appendix).
In both the low and medium energy regions, the data points look reasonably well distributed on both sides of the fitting function (the zero axis) and, also, the corrected residuals look closer to zero than the usual ones. If the effect of residual corrections (of pure graphical concern) looks marginal at low to medium energies, the rightmost panel of Figure 13 clearly indicates the more appropriate character of the corrected residuals to translate the fit quality. = Br(ω → e + e − )Br(ω → π + π − ) and for the Orsay phase (in degrees) using only the (0.75 ÷ 0.82) GeV energy region of each π + π − data sample. The fit is done following the (local) procedure sketched in Subsection 4.2; we use Br(ω → π + π − ) = (1.53 ± 0.13)% [24] and the Orsay phase from [70] as input values for the "τ + PDG" fit. The probabilities are those of the global fit; the last data column shows the contribution of each π + π − data set to the total χ 2 and gives the corresponding number of data points.  Table 2: Probabilities associated with the estimated χ 2 /N p contributions of the π + π − data samples from CMD-2, SND and KLOE10 for the various energy regions and -last columnfor each π + π − data set. Last line is referring to what has been named NSK+KLOE10 in the text and yields a global fit probability is 96.9%; see the last paragraph in Subsection 5.1.2) for more details.
to give more information on their fit quality beyond the global properties just emphasized. This is the purpose of Table 2. We choose here to present the results in terms of probabilities of the χ 2 π + π − contributions to the total χ 2 associated with the corresponding number of data points 36 N π + π − . The last data column thus gives this probability for the data sample(s) provided by the three experiments. From this exercise one can conclude that the three data samples, each as a whole, behave normally with a remarquable fit quality of the SND data sample (χ 2 /N points ≃ 0.7), while KLOE10 (χ 2 /N points ≃ 1.0) and the (collection) of CMD-2 ((χ 2 /N points ≃ 1.1) data samples are quite reasonable.
Even if more approximate, one can perform the same exercise for various energy subregions 37 . This is displayed in the first four data columns of Table 2. The probabilities of the various SND data subsamples (N π + π − = 45) look reasonably flat, pointing to a very good account by the model function and one yields χ 2 /N points 's ranging between 0.5 and 1.1. The probabilities of the CMD2 collection of data subsamples (N π + π − = 82) are less favorable (the χ 2 /N points 's range between 0.7 and 1.3) but remain quite acceptable. The results for KLOE10 (N π + π − = 75) is more appealing as the probabilities are good for three regions (corresponding to χ 2 /N points values of resp. 0.3, 1.2 and 0.5), while the fourth region is much worse (1.7% probability corresponding to a χ 2 /N points of 1.7). This could have been already inferred by inspecting the residuals of the fit in "isolation" shown 38 in Figure (9). Nevertheless, Figure (15)) shows a reasonable account of the fit pion form factor even in this region.
The last line in Table 2, which corresponds to the combination NSK+KLOE10 exhibits good probabilities up to the maximum energy and, in this case, the number of data points per energy sub-region (resp. 42, 53, 71, 36) is large enough to limit the amplitude of possible fluctuations in the π + π − data samples.
The global fit combining the NSK and KLOE10 data also yields F N SK+KLOE10 ω = (1.166 ± 0.036) 10 −6 , which should supersede the present world average value [24] because of the full consistency it exhibits with the largest set of data ever fitted simultaneously and with a splendid probability.

Fits of the e + e − → π + π − Spectra : Concluding Remarks
The results reported just above have allowed us to show that the NSK data are in fair agreement with the physics represented by the π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 annihilation channels, the dipion spectra collected in the decay of the τ lepton and some more decay information listed in Subsection 3.1. This is not a really new result as this conclusion was already reached in our previous analysis [13]. The single difference with [13] is the new way to include the φ → π + π − information (see Subsection 3.3), possible with the BaBar data.
The new information is that the KLOE10 data sample behaves likewise and, more importantly, that the NSK and KLOE10 data sets are consistent with each other as well as with the rest of the physics considered within the BHLS model and the global fit context.
We have performed global fits which have shown that the data from KLOE08 and BaBar have some difficulty to accommodate the global fit context. For what concerns KLOE08, the problem appears to be related with (underestimated?) systematic errors or, possibly, with their correlations. In the case of BaBar data, the problem looks more serious, as it deals with the form factor lineshape itself in the ρ − ω interference region; this issue manifests itself in the value for F ω , much larger than expected.
For this purpose, it is interesting to see the behavior of the most relevant fits of the pion form factor in the ρ − ω resonance region. These are displayed in Figure 14. The best fit obtained using only the NSK data is shown in the top left panel and is clearly quite satisfactory. The top right panel exhibits the case when the NSK data are complemented with the KLOE10 sample; this fit is also quite successful.
Bottom left panel in Figure 14 shows the behavior of the best fit function when only considering the ISR data (BaBar, KLOE08, KLOE10) and bottom right panel when taking into account all existing scan and ISR data. Both are clearly less satisfactory reflecting mostly the tension between BaBar and the KLOE data sets.
As already stated, the largest favored configuration is to include simultaneously the NSK and KLOE10 samples within the fit procedure. Fig. 15 shows the fractional deviations from the fitting functions for the e + e − and τ data derived simultaneously from the global fit. The top panel is clearly consistent with a flat distribution. The bottom panel shows that CLEO and Belle distributions are flat, but also that ALEPH data may show a slight s-dependence starting around 0.8 GeV (quite similar to the top panel in Fig. 55 from [34] or to Fig.12 in [23]). This is a pure consequence of having introduced the KLOE10 data sample 39 .
One may thus consider that the statistically favored configurations (the so-called NSK and NSK+KLOE10 global fit configurations) exhibit flat residual distributions simultaneously in e + e − and τ ππ channels.
The question is now whether the KLOE08 and BaBar data samples can nevertheless help to improve some physics information of important concern as the muon g − 2. This will be discussed in the forthcoming Section. Anyway, the KLOE10 data sample already allows to confirm the results already derived using the NSK data and, even, helps in getting improved results.

The φ Region in the
Up to now, the φ pieces of information used in our fits/predictions are the RPP value for the product F φ = Br(φ → e + e − )Br(φ → π + π − ) = (2.2 ±0.4) 10 −8 and the "Orsay" phase for the φ → π + π − amplitude provided by SND [72], namely −34.0 ± 4.0 degrees. In order to avoid over interpreting the SND phase as an Orsay phase (i.e. identified with the phase of the product F e ωφ (s) g φππ (s) in Eq. (14) at s = m 2 φ ), we have found worth revisiting this assumption. Using 39 Compare bottom panels in the present Fig. 15 and in Fig. 10 in [13] the NSK and KLOE10 data sets as reference π + π − data sample, we have performed several global fits within the BHLS framework. For definiteness, and for an easy comparison, the first line in Table 3 displays the values extracted from [24,72] and used in all fits referred to up to now. The second line in this Table displays the corresponding values reconstructed using the final parameter values from the appropriate global fit. The goodness of fit (94.5% probability) is reflected by the χ 2 values shown in the last two data columns. The reconstructed value for F φ is in very good accord with the (input) RPP value, while the phase is shifted toward slightly more negative values than the SND datum.
The fit results summarized in the third line are obtained by withdrawing the SND phase [72] from the minimization procedure. The fit quality is almost unchanged (95.3% probability) and the expected F φ is still well reproduced; on the other hand, the φ phase moves by about 3σ but remains significantly negative.
In both cases, the pion form factor lineshape exhibits the behavior shown in Figure (1b), e.g. a tiny peak at the φ mass, resembling what is shown in Figure 4 of [72].
The next step has been to remove the reference values for F φ and for the φ phase and supplement the NSK and KLOE10 data with the φ region spectrum from BaBar data ( √ s ∈ [1.00, 1.05] GeV). The important results are shown in the fourth line of Table 3. The fit probability remains good (91.3%), but the χ 2 of the NSK+KLOE10 data is degraded by ≃ 10 units, which may indicate some tension between these data sets and the coherent background beneath the φ peak. The fit value for F φ increases by 50% but remains perfectly consistent with expectations [24]. In contrast, the Orsay phase changes dramatically (from −[30 ÷ 50] • → 150 • ). The lineshape exhibited by the BaBar pion form factor is almost identical to those shown in Figure (1c).  Table 3: Fit results for F φ = Br(φ → e + e − )Br(φ → π + π − ) and for the "Orsay" phase (in degrees) for the φ amplitude. See Subsection 5.2.1 for comments.
As stated above 21 , the [1.00, 1.05] GeV region of the BaBar spectrum carries information on the φ signal but also on the underlying (ρ coherent) background. In order to substantiate its effect, we have finally performed a fit using as input the F φ and φ phase values extracted from fitting the BaBar data up to 1.05 GeV (see Subsection 5.1.1); the corresponding values are given in the fifth data line of Table 3 and, when imposed to the fit procedure, these numbers lead to the results shown in the last line. The fit probability improves to 94%.
As a summary, Table 3 shows that some ambiguity occurs about the pion form factor behavior in the φ region. It becomes interesting to have new data covering this region, in order to decide which among Figures (1b) and (1c) reflects the right behavior. It is also interesting to examine the consequences of this ambiguity on g − 2 estimates.
Without going into details (addressed in the next Section), we can give the values for ∆a µ = a exp µ − a th µ corresponding to the various cases shown in Table 3. In units of 10 −10 , the values for ∆a µ are respectively 39.91 ± 5.21 (line # 2), 39.35 ± 5.23 (line # 3), 37.18 ± 5.16 (line # 4) and 39.32 ± 5.20 (line # 6). The most likely effect this implies on ∆a µ is a possible shif t by 0.59 × 10 −10 (lines # 2, 3, 6). A maximum possible shift is 2.73 × 10 −10 (lines # 2 and 4). Instead, the uncertainty is practically unchanged. Taking into account the uncertainty on a exp µ [39,40], the effect of this shift on the significance of ∆a µ is quite marginal (4.88σ → 4.57σ). Nevertheless, this shift plays as a systematic effect to be accounted for in the final value which will be proposed for a th µ .

The HLS Favored Solution : Fit Parameters and Physics Quantities
In our previous analysis [13], we encountered some tension in the e + e − annihilation data in the φ region between the e + e − → π + π − π 0 channel and the e + e − → K + K − , K 0 K 0 channels. Including all data in the global fit defined what we called "configuration A". Alternatively, we excluded the 3-pion data in the vicinity of the φ mass from the fit. In this case, the corresponding collection of data included in the fit was called "configuration B", which yields a better fit probability. As [13], the present analysis privileges configuration B as a basis. However, as the "observed" tension could not be abnormal, we also consider fits based on configuration A for completeness, essentially in relation with the muon anomalous moment.
Tables 4 and 5 display the global fit parameter values within configuration B, using the e + e − data from only CMD-2, SND and KLOE10. This corresponds to the so-called NSK+KLOE10 case which is favored by the BHLS model and breaking scheme.
The top part of Table 4 displays information directly related with the standard parameters of the original HLS model. The value for a is found, as usual [18] significantly larger than the standard VMD expectation (a = 2) and the value for the universal vector coupling g stays at accepted values [7]. The HLS parameters 40 giving the weight of the anomalous terms [55,15] in the global HLS model (c 3 + c 4 )/2, c 1 − c 2 , get values close to our previous results [13]. The parameter v, related with the 't Hooft determinant terms [54], will be discussed below.
The middle part of Table 4 essentially displays the SU(3) breaking parameters generated by the original BKY mechanism [16]; their role is crucial in the e + e − → KK cross sections and in the V P γ couplings (see Appendix E in [13]). The values for z A , z V and z T show that the breaking of the SU(3) symmetry generates departures from 1 (unbroken case) by 30% to 60%.
The bottom part of Table 4 shows specifically the values for the BKY parameters generating the (direct) isospin breaking mechanism [53,13] defined in Subsection 2.2 : ∆ V , ∆ A and Σ V which vanish when Isospin Symmetry is not broken. Numerically, they are at the level of a few percents, i.e. much smaller than those corresponding to the breaking of the SU(3) symmetry, as expected. The value for g direct ω→ππ which only depends only on the direct IB breaking mechanism is discussed just below.
The fit results shown in Table 5 deal with the subtraction polynomials coming with the various loops (see Section 6 in [13] or, more deeply, Section 12 in [18] and especially Table 3 therein). The bottom part of Table 5 gives the subtraction polynomials of the functions ε 1 (s) and ε 2 (s) which define the vector meson mixing as reminded in Subsection 2.3. They are given for sake of completeness as, generally, they do not carry an obvious intuitive meaning; they mostly play a role in the definition of the vector meson mixing [13]. Nevertheless, as clear from Eqs. 9, via the mixing angles α(s) and β(s), the functions ε 1 (s) and ε 2 (s) play a fundamental role in reconstructing F ω and F φ discussed at several places above and they have been shown to produce the expected effect either isolatedly (F φ only depends on vector meson mixing) or combined with the direct IB [53,13] discussed above (F ω depends on both). Table 6 display physics parameters which can be compared with corresponding information derived by others : different theoretical frameworks or measured quantities not included in our fit procedure; this may give confidence in the validity of the BHLS model and its structure.
If the universal coupling g is indeed returned by fits at usual values, the reconstructed value for g ρππ compares reasonably well with its value derived from fits using the Gounaris-Sakurai parametrization [7] g ρππ = 6.1559.
In their mini-review [24], Rosner and Stone conclude that the most likely estimate for f K /f π can be summarized by f K /f π = 1.197 ± 0.002 exp ± 0.002 CKM ± 0.001 rad. cor. which compares reasonably well with our own estimate reported in Table 6. This also indicates that this proposed value could be included in our fit procedure as an additional constraint.
The following datum is our estimate for g direct ω→ππ . A piece of information recently studied in [71] is related with the ratio G = g direct ω→ππ /g ρππ . Our previous estimate [13], using only the CMD-2 and SND data within our global framework, was (3.47 ± 0.64)% smaller than the findings of [71]; from Table 6, this ratio is now estimated at (6.27 ± 0.94)%, just in between the two estimates proposed by Maltman and Wolfe 41 [71], showing the influence of the KLOE10 data.
Our estimate for the 't Hooft parameter λ [54], which governs the magnitude of nonet symmetry breaking in the pseudoscalar sector (cf. Eq. (7) in [13]), can be derived from the values for v and z A (cf. Eq. (22) in [13]). We confirm that λ is in the range of 10% with a substantial uncertainty (30%); we are not aware of published external estimate of this quantity.
The following part of Table 6 deals with information concerning mixing properties in the pseudoscalar sector. We have adopted as parametrization of the π 0 − η − η ′ mixing those proposed in [75]. Within this parametrization, the π 0 − η mixing is governed by ǫ, the π 0 − η ′ mixing is governed by ǫ ′ . Both parameters are expected to carry values at the few percent level, and one may expect ǫ ′ < ǫ. These expectations are satisfied by our estimates. Moreover, ǫ and ǫ ′ can be expressed [13] in terms of the pseudoscalar singlet-octet mixing angle θ P S and of some ǫ 0 , function of the quark masses [75]. Recent estimates yield ǫ 0 = (1.31 ± 0.24)% [76] and ǫ 0 = (1.16 ± 0.13)% [77], significantly smaller than our former result (ǫ 0 = (3.16 ± 0.23)%). Our present evaluation (see Table 4), is much closer to these expectations.
Our evaluations of θ P S = −13.54 • ± 0.15 • compares favorably with its KLOE measurement 42 [80] . To end up with this topic, we should mention that the fit values for the partial widths π 0 γγ, ηγγ and η ′ γγ are at resp. 0.3σ, 0.5σ and 0.2σ from RPP expectations [24].
Our evaluation for the ω and φ mass and width are also shown in Table 6. These correspond to the usual definition of the inverse propagator with fixed width : They are found in fair consistency with RPP expectations and our global fit estimates compare favorably with the corresponding data extracted from fits to the individual measured annihilation spectra reported in [24].
The last line in Table 6 displays information on the ρ meson mass and width. Concerning objects as broad as the ρ meson, mass and width values are highly definition-dependent [18] and, depending upon the choice, their spread can easily amount to about 20%. In Table 6, we give the location s ρ of the ρ 0 propagator pole located in the unphysical sheet; it has been derived using the parameter values and their error covariance matrix within a Mathematica code. Defining s ρ = M 2 ρ − iM ρ Γ ρ , one can derive mass and width values, also given in Table 6. As far as we know, one can only compare with [81] (see also [82] for a hint on the parametrization dependence of the pole location), who performed fits by approximating the ρ propagator by a Laurent series; this Reference reports M ρ = 757.5 ± 1.5 MeV and Γ = 142.5 ± 3.5 MeV relying on fitting only the data reported in [25], close to our results. The previous estimate of these within the HLS model [18] where larger (resp. 760.4 and 144.6) with comparable uncertainties; however, the data from Belle and KLOE were not available at that time. Because of our approximations (see Subsection 2.6), the present study is not sensitive to differences between the charged and the neutral ρ mesons.

Hadronic Contribution to the Muon g − 2 Estimates
In order to estimate the hadronic contribution to the muon anomalous magnetic moment a µ , the method followed in the present work is identical to the one used in [13]. In this Section, we examine the global fit solutions provided by the data samples listed in Subsection 3.1; the 42 The value published is actually expressed in the flavor basis; for an easy comparison, we convert to the singlet-octet mixing angle using θ P S = φ P − arctan √ 2.
various cases correspond to varying the e + e − → π + π − data sample combination submitted to the global fit procedure. We mostly work under the configuration B which has been reminded in the header of Subsection 5.2.2 just above.

a µ (ππ) Contribution from the Reference m ππ Region
In order to emphasize what is going on, it is worth examining the contribution to a µ (ππ) provided by the invariant mass interval A ≡ [0.630, 0.958] GeV. As for this reference region, all experimental groups have published the numerical integration leading to their estimate of a µ (ππ, A) (using Eq. (1) with the experimental spectra), one can compare our fit outcome with 6.505 ± 0.003 1.279 ± 0.010 0.408 ± 0.061 782.52 ± 0.03 8.66 ± 0.04 1019.25 ± 0.26 4.18 ± 0.02 753.8 ± 0.5 138.10 ± 0.5 0.5682 ± 0.0007 0.1041 ± 0.0004 Table 6: Physics parameters extracted from the BHLS favored fit.
these. When combining data sets within the fit procedure, we compare to the usual weighted average of the a µ (ππ, A) values provided by each experimental group.
The results for a µ (ππ, A) are collected in Fig. 16 and will be discussed from top to bottom. The first point has no experimental partner, as it comes out from the τ +PDG fit described in Section 4. The dashed-dotted line is drawn through the central value of this prediction, the dotted lines show the corresponding ±1σ band.
It should then be noted that most experimental estimates happen to be inside the 1σ band of the τ +PDG prediction. KLOE10 is at the border, while BaBar is rather larger by about 2σ. One should also note that KLOE08 as well as the combination CMD2+SND+KLOE10 (our preferred combination because of its global fit properties) provide central values which coincide almost exactly with the τ +PDG expectation; even if shifting the τ +PDG prediction by (1 ÷ 2) 10 −10 should be performed (see Subsection 6.4.1 below), the overall agreement remains noticeable. This looks to us an important property taking into account the long-standing discussion about the τ − e + e − issue. Our results tend to prove that this issue vanishes once the breaking mechanism is appropriate.
The next point shows the global fit estimate (square symbol) while using 30 "solely" the NSK data (especially those from [27, 28, 29], but also those from [25,26]). The present fit value is found smaller than those published in [13] by about 1 × 10 −10 . This shows that the effect of our updating the φ data decay information extends over the whole spectrum. The uncertainty is unchanged and amounts to a 40% gain compared with the experimental estimate. This gain is a pure effect of the global fit procedure where the channels others than π + π − allow to improve also the π + π − contribution importantly because of the underlying VMD physics correlations.
The following point is derived using 30 "solely" the KLOE08 data sample. As while fitting "solely" the NSK data, the fit outcome also differs by 1 × 10 −10 from the experimental central value. The fit result associated with the KLOE10 data sample is also found in good agreement with the experimental expectation and one should note that both differ by 1.5×10 −10 in opposite directions. In contrast, the fit result for the mixed KLOE10+(weighted) KLOE08 data sample is found to agree very well with expectation (a 0.6 × 10 −10 difference only).
The data point associated with the CMD2+SND+KLOE10 combined data sample is found in much better accord with the average of the experimental data; in this case, the error is also improved by about a factor 2 compared with the experimental value, but, interestingly, first, the central reconstructed value does not exhibit any significant systematic shift (≃ 0.1 × 10 −10 ) compared to the experimental expectation, second, the central reconstructed value is in perfect accord with the τ +PDG prediction, as already noted.
The next pair of points shows the case for the BaBar data. While using the whole spectrum up to 1 GeV, the fitting function provides an estimate for a µ (ππ, A) in good agreement with the expected experimentally integrated value. However, extending the fitting region up to 1.05 GeV, increases the difference between the fit outcome and the experimental expectation from 1.5 × 10 −10 to 4.4 × 10 −10 , which illustrates some issue in the ω region when also fitting the BaBar φ region.
The last two points (ISR and scan + ISR data, resp.) correspond to fitting the whole spectra up to 1 GeV while weighting the KLOE08 and BaBar contributions to the χ 2 . The results look in good correspondence with the weighted averages of the direct integration results, and one remains within the 1σ band of the τ +PDG expectation.
The contribution to a µ (ππ) from the invariant mass interval A ≡ [0.630, 0.958] GeV and the global fit properties (see Section 5) allow to draw a few important conclusions : • (l) There is no mismatch between the τ +PDG expectations and the estimates derived from the fitting functions or those derived from the numerical integration of the measured e + e − → π + π − cross sections. To be more precise, when some departure is observed, it is always closely associated with poor global fit qualities of the corresponding data sample.
• (ll) One can take as reference e + e − → π + π − data the NSK (CMD-2 & SND) and KLOE10 data samples. Separately and together, they are found in perfect accord with all annihilation data considered (π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 final states) as well as with the published dipion spectra in the decay of the τ lepton and the few additional decay information introduced within the global fit BHLS framework. The uncertainty is improved by a factor of about 2.
• (lll) Comparing the central value for a µ (ππ, A) derived from the fitting function -using the NSK and KLOE10 data samples within the global framework -with its direct estimate indicates that the estimate derived from the fit is almost unbiased (≃ 0.1 × 10 −10 ). In contrast, the fitting functions derived from the fits to each of the NSK and KLOE10 data separately may exhibit a small bias (a shift of ≃ 1 × 10 −10 ).
• (lv) One may remark that the global fit using all existing e + e − → π + π − scan & ISR data also provides reasonable values for a µ (ππ, A) without a significant bias and with a much improved uncertainty. This reduced uncertainty may look too optimistic, taking into account the global fits properties of the corresponding sample combinations. Nevertheless, one may note that the last two entries in Fig. 16, derived by weighting the KLOE08 and BaBar samples with respect to all others, do not depart from the general trend.

The Full Muon Anomalous Magnetic Moment a µ
Here we are interested in the theoretical prediction of the muon anomalous magnetic moment and its comparison with the experimental result [39,40]. Our BHLS based global fit results provides the dominant contribution to the hadronic vacuum polarization contribution.
Several global fits have been performed using different groups of e + e − → π + π − data samples. The resulting fitting functions and the parameter error covariance matrices have allowed to derive the leading order hadronic contribution (LO-HVP) to a µ up to 1.05 GeV in various cases. The full HVP is obtained by combining our present results with the other hadronic contributions listed in Table 9 of our previous work [13]. The various evaluations of the muon g −2 are calculated by summing up a LO−HV P µ with the Higher order HVP, the QED and electroweak (EW) contributions and the so-called light-by-light (LBL) contribution; for this purpose, the values given in Table 10 of [13] have been updated (see Table 7 below).
In Figure 17, we display results for the deviation ∆a µ = a exp µ − a th µ between experiment and theoretical predictions in units of 10 −10 , for various fit results which only differ by the indicated e + e − → π + π − data sample considered within the fit procedure.
Specifically, the fit results included in Fig. 17 are all based on the "configuration B" data set introduced in [13]. Furthermore, always included in the fits are all channels listed in Subsection 3.1, especially the τ data.
At top of the Figure, we provide the τ +PDG entry (i.e., without using any e + e − → π + π − data set in the fit procedure). The twin points drawn correspond to ∆a µ = a exp µ − a th µ in units of 10 −10 ; the lower point is obtained using the spacelike data [35,36], the upper one by excluding them from the fit. The number following the point is the estimate for ∆a µ from fits excluding the spacelike data. The significance of this difference is written [x/y σ], xσ being the significance obtained when excluding the spacelike data, yσ when they are included into the fit procedure. The last column in this Figure provides the χ 2 /n points for the e + e − → π + π − data; this number is given only for (combination of) samples used without any weight.
The first remark which can be drawn is that including the spacelike data does not change noticeably the value for ∆a µ nor its statistical significance. The second important point is that all the features noted while discussing a µ (ππ, [0.630, 0.958] GeV) in the previous Subsection survive when dealing with the whole provided spectra (see (l)-(lv) just above).
It was noted in [13] that the use of the ISR data should not dramatically improve the uncertainty of the theoretical estimate : going from scan data only (the NSK entry) to scan+ISR data (the "Global Fit" entry) improves the uncertainty by only 0.3 × 10 −10 ; this is the major gain of working in a global fit context [13]. In contrast, the central value exhibits slightly more sensitivity, as can be seen by comparing the various entries.
Anyway, the present analysis allows to conclude that the discrepancy between the measurement for a µ [39,40] and the Standard Model prediction is certainly in the range of 4.5σ or slightly higher. Taking into account the already reported fit qualities, the most favored value for a LO−HV P µ is : a LO−HV P µ = (681.23 ± 4.51) 10 −10 .
This value is derived by introducing the e + e − NSK and KLOE10 data samples within our global fit procedure; systematics on this number will be discussed with the final result (Subsection 6.4). Finally, we have redone several of our fits within the configuration A. These confirm the results derived within the configuration B which are the body of the present paper. One also confirms the 'tension' reported in [13] between the 3-pion spectra in the φ region and the data collected for both e + e − → KK cross sections. Whether this is an experimental issue, or if this indicates that our breaking model should be extended 43 is an open issue. Nevertheless, some physics results derived within the configuration A will be presented just below for completeness. Fig. 18 reports on some recent estimates of the muon anomalous magnetic moment a µ together with the BNL average value [39,40]. In this Figure, our (favored) results -derived using only the CMD-2, SND and KLOE10 data samples -are given under the entry tags BHLS::A (for configuration A) and BHLS::B (for configuration B). Some evaluations proposed by other authors are also shown; they are extracted from [6] (DHMZ10), [7] (JS11), [8] (DHea09) and [9] (HLMNT11).

Comparison with Other Estimates of a µ
The analyses reported in [6], [9] and [7] provide g −2 evaluations based on all existing data samples. In particular, all π + π − data samples are used (ISR and scan) and all other existing annihilation data samples, especially the 3-pion data samples. These evaluations of g − 2 may also include or not the τ data. Therefore, it is also instructive to give the results derived from BHLS under the same conditions, i.e. working in the so-called configuration A [13] and using all available ISR data samples beside the scan (NSK) data, whatever are the global fit probabilities. The corresponding results are reported at the "Full BHLS" entries in Fig. 18 More precise information is also provided in Table 7. However, one should note that, while the fit probabilities are good when using only the scan (NSK) data (90.6 % or 72.7 % respectively), they are below the 10 −4 level 35 when all the ISR data samples are included. This reflects the statistical inconsistency between both KLOE data samples on the one hand and the BaBar data sample on the other hand.
We have explored this issue in cases including the τ data samples. Without going into much details, it looks worth to comment about fit qualities. The just reported probabilities bounded by 10 −4 are worrying, taking into account that the problematic data samples have been weighted. For this purpose, we have first performed a fit discarding the BaBar data sample, to check the level of statistical inconsistency of NSK, KLOE10 and (weighted) KLOE08; the obtained fit probability is 61.6 % and the numerical results are almost identical to those already derived using only NSK and KLOE10. In particular, the scan data yield χ 2 (NSK)/n = 136/127 = 1.07 and F ω = (1.11 ± 0.03) 10 −6 , i.e. a limited distortion of the scan data description. The same exercise performed using NSK and (weighted) BaBar yield 12.4% probability but χ 2 (NSK)/n = 171/127 = 1.35 and F ω = (1.52 ± 0.04) 10 −6 ; i.e. a severely degraded description of the scan data.  Table 7: The various contributions to 10 10 a µ . ∆a µ = a exp µ − a th µ is given in units of 10 −10 . These results have been derived using all existing annihilation data samples (configuration A).The naming "scan" and "ISR" in the subtitles refers only to the π + π − channel. For the measured value a exp µ , we have adopted the value reported in the RPP which uses the updated value for λ = µ µ /µ p recommended by the CODATA group [83].
It should be noted that, when using all ππ data samples, the value returned by the BHLS model is close to JS11 [7]. An interesting point here is that [7] uses a mixing procedure (γ − ρ − ω) parent to ours (ρ − ω − φ), while the IB effects are introduced differently in order to get the DHMZ10 result [6]. Comparing the Full BHLS values derived when using or not the τ spectra indicates that BHLS is much less sensitive to the use of the τ spectra that the method from [6] (compare the DHMZ10 entries). Finally, comparing DHea09 [8] and BHLS when no ISR data are considered, one observes a good consistency; this seems to indicate that the problem lays in the amount of isospin breaking compared with τ .
Comparing now all the BHLS results (BHLS::A,BHLS::B and Full BHLS) indicates a remarkable stability for ∆a µ = a exp µ − a th µ under various fit conditions. One should also note that the τ +PDG prediction is always in quite nice agreement with using any kind of configuration for ππ data samples (BHLS::A,BHLS::B and Full BHLS), including those discarding the τ spectra and the ISR data (see also Table 7).
Nevertheless, one should stress that, because of the poor fit quality yielded when using the KLOE08 and BaBar data samples beside the NSK and KLOE10 ones, we do not consider the results shown in Table 7 and in Fig. 18 under the entry tag "Full BHLS" as reliable. As a general statement, we prefer relying on estimates using only the NSK and KLOE10 data samples for which the probabilities are good (≃ 90% ÷ 98% for BHLS::B, ≃ 56% ÷ 73% for BHLS::A). Numerically, the BHLS::B and BHLS::A estimates happen to be always close to each other. Nevertheless, we prefer favoring the BHLS::B result over BHLS::A for reasons already explained.
The statistical significance for ∆a µ = a exp µ − a th µ is displayed on the right-hand side of the Figure for each of the reported analyses. The updating of the φ → ππ treatment within our minimization code turns out to increase the significance for ∆a µ by about 0.5σ without any significant change in the fit qualities. The present study, however, has shown that this update is highly favored by the φ information carried by the BaBar spectrum.
On the other hand, one should note that using all ISR data samples only leads to a ≃ 4% improvement for the uncertainty on a th µ compared with discarding all of them. This is actually, as reminded in the Introduction, a specific consequence of the global fit procedure.
Another remark is worth being expressed : Whatever is the configuration (A or B), whatever is the set (or subset) of π + π − data samples considered, Table 7 and Fig. 18 show that the role of the τ data is not limited to returning improved uncertainties on a th µ . Indeed, one should note that the global fits using the τ data always returns a th µ from fits larger than when discarding them; the difference amounts to 1 ÷ 2 units.
Before closing this Section, it is worth remarking that the central values for a LO−HV P µ derived using scan data only or scan ⊕ ISR are very close to each other. This is due to a numerical conspiracy of the BaBar and KLOE data samples : The former data tend to increase a LO−HV P µ and the latter to decrease it, by a quite comparable amount. Taking into account the above reported fit probabilities, the increase of a LO−HV P µ is not justified on statistical grounds; in contrast the decrease of a LO−HV P µ produced by KLOE10 (and also KLOE08) is supported by quite good probabilities of simultaneous fits of scan and KLOE data.

The HLS Favored Solution : Properties Of Its Estimate For a µ
It follows from the present study that a consistent data sample can be defined which allows to yield a reliable estimation of a th µ . This sample contains all the collected data samples covering the π 0 γ, ηγ, K + K − and K L K S annihilation channels and the dipion τ spectra up to 1 GeV. It also contains the π 0 π + π − annihilation channel. For this last process, we prefer working within the so-called configuration B by excluding the data sets collected around the φ mass. However, as clear from Figure 18, taking them into account (configuration A) does not lead to substantial differences. Concerning the fundamental π + π − annihilation channel, the statistical properties exhibited by the fits lead us to rely on using only the scan (NSK) and the KLOE10 data samples; indeed, for this set of data samples one does not observe any mismatch between the claimed uncertainties and the statistical properties of global fits.

Systematics Due To τ Data
Here we come back to the issue of IB, discussed in Subsection 2.6 before, and how the value for a th µ can be affected by the approximations performed on the τ spectrum versus e + e − , namely δM ρ and δΓ ρ ?
So, in total, if one assumes that the traditional evaluations of the corrections to τ estimates apply without any thinking to the BHLS estimate a (τ ),HLS µ , one may have to consider a shift δa had,LO µ [ππ, τ ] ≡ δa µ (tot) ≃ −1.4 10 −10 , if one relies on Belle [23] estimates. In the "τ +PDG" fit configuration, there is no data on the e + e − → π + π − annihilation process. Therefore, the sensitivity to the ρ 0 meson lineshape is certainly quite marginal, while the use of τ data introduces a strong influence of the ρ ± meson. Hence, the value derived for a th µ is likely to be corrected for the (δM ρ , δΓ ρ ) effects discussed in Subsection 2.6 and just reminded.
In the QFT-driven BHLS parametrization, which is similar to what we know from the electroweak Standard Model, the corresponding checks show that the effects of shifting the width by radiative effects leads to much more stable results (i.e. smaller effects). More precisely, using only the Belle data corrected for IB effects as sketched above and the e + e − data from CMD-2, SND and KLOE10 (the NSK + KLOE10 + τ case), the shift for a µ does not exceed ≃ 0.5 10 −10 .
If one accepts the evaluation method described above, a first estimate of the (δM ρ , δΓ ρ ) effects on a th µ is δa had,LO The difference ∆a µ (τ + PDG) − ∆a µ (NSK + KLOE10 + ✚ τ ), obviously carries the same information as τ + PDG does not depend on e + e − data and NSK + KLOE10 + ✚ τ does not depend on τ data and, thus, does not call for any correction. This evaluation yields −3.81, about twice the preceding estimate and in the same direction.
One may prefer favoring the most constrained result, i.e. ∆a µ (NSK + KLOE10 + τ ). Indeed, in this case, the estimation of IB effects is certainly much more constrained by having simultaneously data strongly affected by IB effects (e + e − ) and data marginally depending on them (τ ). The single motivated uncertainty this provides on a th. µ (NSK + KLOE10 + τ ) is a (possible) negative shift (−2 10 −10 ).

Other Sources Of Systematics
The discussion above leads us to conclude that the most favored value for ∆a µ = a exp µ − a th µ is 39.91 ± 5.21 (in units of 10 −10 ) derived from using NSK + KLOE10 + τ . The comparison with data performed in Subsection 6.1 tends to indicate that this value should be almost unbiased.
However, the ambiguity about the choice of the Orsay phase discussed in Subsection 5.2.1 may also indicate a possible systematic shift, likely of 0.59 × 10 −10 , but certainly limited to 2.73 × 10 −10 . This upper bound supposes that the BaBar spectrum behavior in the φ region should be trusted 44 . Only a devoted measurement of the pion form factor in the φ region may allow to conclude about this possible systematic shift.
Looking at the differences between BHLS::A and BHLS::B, one also observes that the most significant difference between both solutions is concentrated -as could be expected -in the contribution of the 3-pion data around the φ peak. This may lessen ∆a µ by 1.26 × 10 −10 .
On the other hand, it is worth noting the behavior of the predicted pion form factor in the region m ππ < 0.5 GeV down to the negative s region. Figure 19 shows that the pion form factor derived from fits describes quite well the (highly constraining) spacelike data down to about −0.15 GeV 2 . Whether spacelike data are included within the data sample submitted to fit (full curve) or not (dashed curve) does not make any difference (both curves coincide within the thickness of the lines). This is a noticeable property as, following Section 6 in [19] and reminded in the Appendix, the NA7 data should be rescaled. The rescaling factor 1 − λ is such that λ depends on the spacelike data, on their scale uncertainties (0.8 % for NA7 [35]) and on the fitting function. Therefore, by ascertaining the threshold behavior of the pion form factor, Figure 19 proves that the contribution of the ππ threshold region to g −2 predicted using the NSK and KLOE10 data is reliable and confirmed by the existing spacelike data. Indeed, these are well accounted for without any regard to their being included within the data sample submitted to fit (see Figure (19c)).

The HLS Favored Estimate for a th µ
From the discussion above, it follows that the preferred solution is those derived using NSK+KLOE10+τ . Its estimate for ∆a µ carries a statistical significance of 4.88σ. On the other hand, Figure 18 has shown that using or not the τ spectra does not significantly modify the significance for ∆a µ (5.0σ → 4.9σ). In contrast, Figure 18 also shows that using the KLOE10 data sample increases the significance by about 0.5σ.
In view of all these considerations, especially the issues encountered with the φ region physics information, the most motivated estimate we can propose is : The quoted [δa th µ ] φ is an estimate of the possible uncertainties affecting the φ region and emphasized in Subsection 5.2.1 and also just above; the central value is derived from a fit within configuration B and avoiding the use the SND phase constraint. [δa th µ ] φ should not be added in quadrature to the theoretical error but linearly to the central value. Accounting for this possible shift slightly lessens the significance for ∆a µ from 4.81σ → 4.67σ. The uncertainty attributed to τ effects, is also as a possible shift and plays in the opposite direction. Therefore, Figure 18 illustrates that the global fit approach allows for a gain in the significance for ∆a µ of about 1σ or more compared to the traditional methods based on numerical integration of the experimental cross sections [8,6,7,9]. Moreover, our approach has allowed us to define the largest set of data samples which exhibit full consistency with each other.
The results involving the fit of the scan ⊕ (all/selected) ISR data, do not lead to a noticeably smaller uncertainty for ∆a µ compared to using only scan data, as already inferred in [13]. Indeed, Figure 18 indicates that in changing the fit conditions [scan ⊕τ ] → [scan ⊕τ ⊕ KLOE10], the uncertainty marginally improves (5.28 → 5.21), while the central value moves by 3 units (36.88 → 39.91), increasing the significance (4.5σ → 4.9σ).

Conclusion and Perspectives
The e + e − → π + π − annihilation channel has been widely studied and one is faced with several data sets collected by different groups under various conditions. Beside the scan experiments performed by the CMD-2 and SND collaborations which have produced valuable data samples [26,27,28,29], experiments using the Initial State Radiation (ISR) method have also been performed by the KLOE and BaBar collaborations and have produced (much) higher statistics data samples [31, 32, 33]. The noticeable gain in statistics -for what concerns the π + π − intermediate state -is partly balanced by the issues raised by the dominance of systematic uncertainties. Indeed, comparing the properties of these data samples reveals inconsistencies which leads to somewhat contradictory conclusions concerning the predicted value for the muon anomalous moment.
It thus becomes a relevant challenge to find a tool able to examine critically the properties of the various available samples and substantiate their differences. The present study proves that the existing data for : • The e + e − annihilation to the π 0 γ, ηγ, π + π − π 0 , K + K − , K 0 K 0 final states, • The dipion spectra in the τ ± → π ± π 0 ν decay, • Some radiative decays (namely, π 0 /η/η ′ → γγ, ρ ± → π ± γ, K * → Kγ, η ′ → ωγ, φ → η ′ γ) as given in the RPP, supplemented with some limited isospin breaking (IB) information, provide a benchmark able to reconstruct with a noticeable precision the pion form factor measured in the e + e − → π + π − annihilation. As could have been expected, the requested IB pieces should cover the ω/φ → π + π − decays, but also the ρ 0 → e + e − decay, rarely stressed explicitly [7]. The tool for this prediction is the broken HLS model (BHLS) defined and studied in [13]. In this last reference, the scan data for the e + e − → π + π − annihilation channel were used. However, the present study shows that replacing these data by only the ρ 0 /ω/φ decay data just referred to, allows BHLS to pin down all the parameters of the model and provide strikingly precise predictions for e + e − → π + π − . These have been named "τ +PDG" predictions of the pion form factor. As these predictions are clearly a new way to formulate τ based predictions for the HVP, BHLS provides also a new tool to explore the long reported τ versus e + e − discrepancy between the HVP evaluations.
Then, the present study leads to the following conclusions : • 1/ There is no mismatch between τ based predictions and direct e + e − evaluation of the HVP if one relies on BHLS and, especially, on its isospin breaking scheme. This is essentially implemented in two steps : direct breaking at the HLS Lagrangian level followed by vector meson mixing. The vector meson mixing, which generates physical vector fields, is unavoidable because, at first order in breaking parameters, the ideal vector fields are no longer mass eigenstates. This simply confirms the study in [13].
• 2/ The τ +PDG predictions compare astonishingly to the KLOE data samples all along the spectra. Because of their influence on the PDG information, it is somewhat paradoxical to observe a difference between the NSK data and the τ +PDG predictions, more marked than with both KLOE data samples. However, a 1σ modification of the width for ω → π + π − allows to recover a good agreement between the NSK data and the τ +PDG predictions. In contrast, the predicted pion form factor reveals an important disagreement with BaBar data, essentially concentrated in the [0.74, 0.78] GeV region.
• 3/ By varying the fit conditions, it has been found that the KLOE10 and the scan (NSK) data samples together lead to global fits with outstanding statistical properties. Even if the KLOE08 data sample has properties similar to KLOE10 (see Figures 10 and 9), it looks premature to include it within the global fit procedure 45 .
From our analysis, one also gets : which exhibits a significance for ∆a µ = a exp µ − a th µ at a (4.7 ÷ 4.9)σ level, significantly larger than the results fully derived by direct numerical integration of the experimental cross sections. This estimate is free of any reweighting going beyond the reported uncertainties affecting the data samples involved in the derivation.
Some additional remarks are worth being made concerning this result and our approach : • 4/ The cross sections which should be integrated within the BHLS model to evaluate a th µ can be considered as an "effective field theory induced interpolation" between data points. By using a relatively fine energy binning, the interpolation uncertainties are certainly minimized in all regions, and for all cross sections, which exhibit sharp energy variations. This is certainly an advantage over the standard method which should rely on trapezoidal estimation between relatively distant measured data points, possibly improved by taking somewhat into account the local curvatures.
• 5/ One has certainly noted that the value for a th µ we privilege is smaller than all estimates involving all ISR data. As explained above, this is because the effects of the KLOE sample is not balanced by the BaBar data. However, as noted several times above, the value for F ω is a criterion which leads us to conclude that NSK and KLOE10 only have homogeneous properties which justifies to privilege a global simultaneous fit of these.
• 6/ The description of the BaBar data within the global fit framework does not look worse than in the really standalone fit published by this Collaboration. Actually, the main issue is not that much the description of the BaBar data sample stricto sensu than its consistency with the underlying physics correlations with other channels implied by the global approach within the same energy range. As, it is unlikely that BHLS may have some inherent reason to exhibit some tropism towards the scan or KLOE data, we indeed consider the a th µ value we propose as the most motivated one on statistical grounds. Of course, new ISR data one can expect from Belle are clearly valuable.
Nevertheless, the results derived using all π + π − ISR data samples have been examined (see Figure 17) and shown to agree with the favored estimate just quoted. Because of the properties expected from the (BHLS) global framework recalled in the Introduction, using only the KLOE10 and NSK data samples provides already as precise results as those derived using all available ISR data samples and the traditional evaluation method of a LO−HV P µ . This allows to avoid arguing on the relevance of (re)weighting procedures, which are always delicate matters.
Our work tells that the significance for ∆a µ starts to be close to the 5σ level. How close to this value it could be, should be confirmed by more precise annihilation data in the [0.95 , 1.05] GeV mass region, especially in the π + π − , π + π − π 0 final states. Finally, new measurements for the experimental value of the muon g − 2 are planned [85,86]; they are important so as to confirm the central value for a exp µ , and also to lessen the experimental uncertainty which starts now to be dominant when estimating the significance for ∆a µ .

A Graphical Account of Correlated Uncertainties
Most of the ππ data we are dealing with are subject to systematic bin-to-bin correlated uncertainties. In the most simple case, this is only a constant (not s-dependent) scale uncertainty with a supposed standard deviation σ. In this case, the χ 2 to be minimized can be written : where m is the vector of measurements, M is the corresponding vector of model values and V is the statistical error covariance matrix, which generally absorbs the parts of systematics which are not bin-to-bin correlated. An appropriate choice is A = m. λ is the parameter which reflects the scale uncertainty; it can be considered as a random variable and each particular experiment can be viewed as a particular sampling of this. Assuming λ independent of the measurements and of the model parameters, its particular value is fixed by solving ∂χ 2 /∂λ = 0 which yields : Substituting this into Eq. (20) leads to the standard χ 2 to be minimized when dealing with a constant bin-to-bin correlated uncertainty (global scale error) : When plotting superimposed the data and the model function, or their difference (standard residual distribution), the drawing may be in visual contradiction with the final χ 2 obtained when minimizing Eq. (21). This generally happens when the scale uncertainty is not very small. However, Eq. (20) clearly tells 46 that it is not m which should compared to M but m − λA; this has been accounted for while plotting the spacelike data [35,36] in Figure (19). A similar issue is encountered with the (standard) residuals m i − M i and is solved by plotting the "corrected residuals" m i − M i − λA i rather than m i − M i .
With the advent of ISR data, the structure of systematic uncertainties has become much more complicated : there is generally a large number of independent bin-to-bin correlated uncertainties 47 ( So, let us assume for definiteness that some data sample -represented by its measurement vector m of length n i and its n i × n i statistical error covariance matrix V -is subject to n α independent bin-to-bin and s-dependent correlated uncertainties. These are generally considered as independent random variables of zero means and are supposed to carry s-dependent standard deviations (σ α (s), α = 1, · · · n α ). In this case, Eq. (20) should be generalized to 48 : (summation over repeated greek indices is understood), defining the n α vectors B α by their components [B α ] i = σ α (s i )m i and using otherwise obvious notations. In this case, the λ α parameters can be considered as independent random variables each of zero mean and of unit standard deviation. Eq. (23) shows that n α scale parameters are to be fixed, one for each of the bin-to-bin correlated uncertainty functions. A trivial mathematical recurrence, using in sequence the solution for ∂χ 2 /∂λ α = 0 (cf. Eq. (21)), allows to prove that the χ 2 to be minimized becomes : Then, one simply gets W ij = [V ij + α σ α (s i )σ α (s j )m i m j ] as expected for independent bin-to-bin correlated uncertainties. No explicit knowledge of the values for the various λ α is needed for this derivation.
For graphical purpose however, one may be interested in getting numerically the various λ α values. One can easily prove that they are solutions of the linear system : 25) which is certainly always regular for systematics small compared with the data. In this case, plotting together m− α λ α B α and M gives a more appropriate representation than plotting m and M; on the other hand, m− α λ α B α −M better reflects the χ 2 value derived by minimizing Eq. (24) than m − M.
This does not exhaust all issues due to bin-to-bin correlations in the graphical representation of data and fits; indeed, at least the spread of the model function (when derived from fits) adds a usefull piece of information. Anyway, as all other correlations (especially those due the statistical errors) are not considered, plotting m − α λ α B α cannot exhaust the full issue but it certainly goes closer to facts.     of the e + e − → π + π − cross section measured by SND detector in the energy region 400-MeV < √ s < 1000-MeV.
[30] M. Davier, A. Hoecker, B. Malaescu, C. Z. Yuan, and Z. Zhang, (2009), 0908.4300, Reevaluation of the hadronic contribution to the muon magnetic anomaly using new e + e − → π + π − cross section data from BABAR.  Top left panel displays the best fit using only the CMD-2 and SND data as e + e − → π + π − spectra. Top right panel shows the case when the CMD-2, SND and KLOE10 data samples are fitted simultaneously. Bottom left panel shows the best fits when the data samples considered are those from BaBar, KLOE08 and KLOE10. Bottom right panel shows the best fit using simultaneously the CMD-2, SND, KLOE08, KLOE10 and BaBar data samples. Measurement of the e + e − → π + π − (γ) Cross Section with the Initial-State Radiation Method at BABAR. [48] A. Manohar, (2001), Hadrons in the 1/N expansion.