BHLS$_2$, a New Breaking of the HLS Model and its Phenomenology

Previous studies have shown that the Hidden Local Symmetry (HLS) Model, supplied with appropriate symmetry breaking, provides an Effective Lagrangian (BHLS) able to encompass a large number of processes within a unified framework. A global fit provides a simultaneous description of the $e^+ e^-$ annihilation into $\pi^+\pi^-$, $\pi^0\gamma$, $\eta \gamma$, $\pi^+\pi^-\pi^0$, $K^+K^-$, $K_L K_S$, the dipion spectrum in the $\tau$ decay and other light meson partial widths. Additional breaking schemes extend its scope to spacelike processes within a new framework (BHLS$_2$). The phenomenology previously explored with BHLS is revisited in the BHLS$_2$ context with emphasis on the $\Phi$ mass region. It is shown that BHLS$_2$ addresses the close spacelike region covered by NA7 and Fermilab data; also,the Lattice QCD pion form factor is accurately predicted by BHLS$_2$ using fits to annihilation data only. These channels contribution to $a_\mu$ over the range of validity of BHLS$_2$ is shown to strongly reduce the former BHLS systematics. The uncertainty on $a_\mu^{\rm th}(\sqrt{s}<1.05$ GeV) is much improved compared to standard approaches relying on direct integration methods of measured spectra. Using the BHLS$_2$ results, the leading order HVP contribution to the muon anomalous moment is $a_\mu^{HVP-LO}= 686.65 \pm 3.01 +(+1.16,-0.75)_{syst.}$ in units of $10^{-10}$. Using a conservative estimate for the light--by--light contribution, our evaluation is $a_\mu^{\rm th}=\left [11\,659\,175.96 \pm 4.17 +(+1.16,-0.75)_{syst.}\right]\times 10^{-10}$. The $\rho^0\gamma$ mixing in dispersive and LQCD approaches is discussed and may amount to a shift $\delta a_\mu[\pi\pi]_{\rho\gamma} = +(3.10\pm0.31)\times 10^{-10}$ at LO+NLO, presently treated as additional systematics. Taking this shift into account, $a_\mu^{\rm th}-a_\mu^{\rm BNL}$ exhibits a significance not smaller than $3.8 \sigma$.

µ of these annihilation channels over the range of validity of BHLS 2 (up to ≃ 1.05 GeV) is updated within the new BHLS 2 framework and shown to strongly reduce the former BHLS systematics. The uncertainty on a th µ ( √ s < 1.05 GeV) is much improved compared to standard approaches relying on direct integration methods of measured spectra. Using the BHLS 2 results, the leading order HVP contribution to the muon anomalous moment is a HV P −LO

Introduction
The Standard Model is widely recognized as the (gauge) theory which unifies the whole realm of weak, electromagnetic and strong interactions among quarks, leptons and the various gauge bosons (gluons, photons, W ± , Z 0 ). For the physics processes -and quantities -involving strong interactions, QCD is at work under two different regimes tightly connected with the energy involved and the onset of the perturbative regime is a priori expected to occur at high energies. However, as clearly shown by the data recently collected by KEDR [1,2] on the ratio R(s) = σ(e + e − → hadrons)/σ(e + e − → µ + µ − ), the {u, d, s} sector of this ratio reaches the perturbative regime at energies as low as ≃ 2.0 GeV (see [3,4], for instance); above this energy, the observed departures from perturbative QCD predictions -including, of course, the cc and bb threshold effects therein -are only spikes and narrow bumps associated with the charmonium and bottomium states which require an additional specific treatment in the mass range where they are located.
Therefore, above a relatively low energy threshold, perturbative QCD 1 can be expected to provide accurate predictions. However, in the low energy region where the non-perturbative regime of QCD is involved, getting theoretical predictions able to compete with the accuracy of some important experimental measurements may be challenging. This is, in particular, the issue met with the photon hadronic vacuum polarization (HVP), which contributes importantly to the muon anomalous magnetic moment a µ , one of the best measured particle properties. The hadronic part of a µ is related to the so-called R(s) ratio and is given (at leading order) by : where s 0 = m 2 π 0 is the lowest hadronic threshold andK(s) is a known smoothly varying positive function [5]. The 1/s 2 factor strongly enhances the low energy contribution of the spectrum, i.e. just in the validity domain of non-perturbative QCD. Fortunately, as obvious from the KEDR data [1,2], one can undoubtfully consider that the non-perturbative regime does not extend to energies larger than ≃ 2.0 GeV. So, the real issue with predictions for objects such as a µ is to get precise estimates of the effects covered by the non-perturbative regime of QCD. Chiral Perturbation Theory (ChPT) [6,7], the low energy limit of QCD, is of limited help for the present purpose as its realm does not extend much beyond the 400 ÷ 500 MeV region and misses the quite important meson resonance mass region.
The most promising approach to the non-perturbative regime of QCD is certainly Lattice QCD (LQCD) which already provides valuable information at low energies [8]. The recently derived LQCD evaluations of a Hadr µ [9,10,11,12] have been found in accord with the experimental measurement performed at BNL [13,14]. However, the magnitude of the reported uncertainties is by far too large to fruitfully compare with the already existing BNL datum and, a fortiori, with the measurements expected from the Fermilab experiment [15,16], already running, or from the experiment planned to start later on at J-PARC [17], as both are expected to improve the uncertainties by a factor of 4. So much progress remains to be done before getting satisfactory uncertainties from LQCD. This leaves room for low energy Effective Resonance Lagrangian Approaches, which can contribute to improve the knowledge of a Hadr µ by providing a good description of the rich amount of experimental data collected in the timelike region. Among the richest possible Lagrangians, the elegant Hidden Local Symmetry (HLS) Model [18] is worth to be considered; it has been proven to be equivalent to RχPT [19,20], provided consistency with the QCD asymptotic behavior is incorporated. It thus follows that the HLS model is a motivated and constraining QCD rooted framework, moreover, easy enough to handle in phenomenological applications.
The original HLS Model deals with the lowest vector meson nonet and provides a framework for the e + e − annihilations naturally bounded by the Φ mass region -i.e. up to ≃ 1.05 GeV. It thus represents a tool giving a handle on a mass region contributing for ≃ 83% of the total muon hadronic VP. The region extending from just above the Φ meson mass to 2 GeV only contributes for 7%, slightly less than the [2 GeV, ∞] region (10%).
As such, the non-anomalous HLS Lagrangian [21] sets up a unified framework which encompasses the e + e − → π + π − , e + e − → K + K − and e + e − → K 0 K 0 annihilation channels; the τ ± → π ± π 0 ν spectrum also belongs to the same framework, allowing naturally a simultaneous treatment of all the mentioned annihilation and decay processes. On the other hand, the HLS Model possesses an anomalous sector [22,18] which also brings the e + e − → π 0 γ, e + e − → ηγ and e + e − → π + π − π 0 annihilation channels inside the same framework, together with some radiative decay modes and the dipion spectra in the η/η ′ → π + π − γ decays as shown in former studies [23,24]. The annihilation channels just listed exhaust almost completly the processes contributing to a Hadr µ below 1.05 GeV; indeed the missing channels 2 (4π, 2πη, . . . ) contributes only ≃ 2 of the full HVP. Such a broken HLS (BHLS) Model has already been built up [25,26] and shown to provide a pretty good simultaneous description of almost all data samples covering the six channels listed above up to ≃ 1.05 GeV. This work has also proved that there was no contradiction between the e + e − → π + π − and τ ± → π ± π 0 ν spectra.
The non-anomalous [21] and anomalous [22] sectors of the HLS Model thus open a unified framework able to encompass a large corpus of data and physics processes. However, as such, the HLS framework -with only the universal vector coupling g as free parameter-cannot pretend to provide a satisfactory simultaneous description of the wide ensemble of high statistics data samples collected by several sophisticated experiments in several annihilation channels.
In order to achieve such a program, the HLS model must be supplied with appropriate symmetry breaking mechanisms not present in its original formulation [18]. A first successful attempt has been done in [25,26,27] which has set up a model (BHLS) based on a breaking mechanism originally proposed in [28], hereafter named BKY.
The present work reports on a new breaking scheme for the HLS model which aims at improving the behavior of the form factors in the dipion threshold region and also in the Φ region where the original BHLS [25] meets some difficulty leading to introduce additional systematic uncertainties in the evaluation of the muon g − 2. The model (BHLS 2 ) presented here will be shown to widen the scope of BHLS, in particular to the spacelike domain and improve the BHLS prediction for the muon g − 2. Indeed, using variants of BHLS 2 , together with the original BHLS should allow to evaluate the model dependence effects in the estimates 2 The e + e − → η ′ γ channel has been considered in [25,26] and found to contributes to a Hadr µ at the 10 −12 level only and can thus be safely discarded. of derived physics quantities, noticeably the photon HVP.
The layout of the paper is as follows. Sections 2 gives a brief reminder of the original HLS model while Section 3 briefly reminds the BKY mechanism, still involved in BHLS 2 . Sections 4 and 5 introduce a new breaking mechanism at the level of the covariant derivative which is the basic object leading to the original HLS model. Section 6 describes the O(p 4 ) terms [18] of the HLS framework which are included in the BHLS 2 Lagrangian.
At this step, one has at hand the first variant of the BHLS 2 model, named Basic Solution (BS), and, relying on its vector meson mass term, Section 7 defines its parameter properties. Remarking that there is no fundamental reason why the neutral vector fields involved in physics processes should be their ideal combinations (ρ I 0 , ω I , Φ I ), one allows for combinations of these via a mechanism called Primordial Mixing described in Section 8. This gives rise to a second variant of the BHLS 2 model, named Reference Solution (RS). In Section 9, one first reminds the parametrization of the propagators for the ω and Φ fields [25] in connection with form factor behaviors at s = 0; this part of the work also shows that the breaking of nonet symmetry in the vector meson sector must be accompanied by a breaking of SU(3) of same intensity.
The dynamical breaking of vector meson [29,25] is revisited in Section 10 and the role of anomalous loop corrections is emphasized, especially in order to connect smoothly the timelike and spacelike branches of the pion form factor. The derivation of the pion and kaon form factors is the subject of Section 11. The anomalous sector is fully analyzed in Section 12 and the cross sections for the e + e − → (π 0 /η)γ and e + e − → π + π − π 0 annihilations are derived.
Altogether, the above listed Sections and the Appendix fully report on the BHLS 2 model properties and tools. Section 13 provides a comprehensive discussion of the available data samples falling into the BHLS 2 scope and their peculiarities. Relying on several preliminary studies, one comments here on the 3 discarded data samples; one should note that the number of data samples found to accomodate with each other within the BHLS 2 framework exceeds now 50 and the total number of data points reaches 1237.
Sections 15 and 16 give a full description of all aspects of the fits to the 6 annihilation channels performed within the global BHLS 2 framework. Moreover, Section 16 also provides a comprehensive analysis within our global framework of the existing KK data samples. The aim is to substantiate the issue between the BaBar, CMD-3 and CMD-2 spectra revealed by fits. A summary of the various fit properties is made in Section 17.
The description of the spacelike region for the pion and kaon form factors as coming from BHLS 2 is the purpose of Section 18. The existing model-independent data for the π ± and K ± form factors are found to naturally accomodate our global framework, giving support to the low energy behavior predicted by BHLS 2 for all meson form factors. This fair agreement extends to the π ± form factor data provided by several Lattice QCD (LQCD) groups. One may infer from these data and from other model-dependent data that BHLS 2 should almost certainly apply down to ≃ −1 GeV 2 . The prediction for the K 0 form factor is shown to cope with expectations at s ≃ 0. Information on the charge radii of the π ± and K ± mesons is provided. Section 19 examines the P -wave ππ phase shift prediction and shows its agreement with data and other predictions, including our former BHLS, over an energy range extending up to the Φ mass region.
The values found for some parameters, noticeably the a HLS parameter and some BKY breaking parameters, deserve a specific account given in Section 20.
The main motivation of the present work is the study of the BHLS 2 predictions and their consequences on the muon g − 2. This topic is addressed in Section 21. The part (≃ 83%) of the muon HVP contribution to a µ covered by BHLS 2 is estimated in several contexts to derive a HV P −LO µ ( √ s < 1.05 GeV) and estimate additional systematics possibly due to modeling effects and to observed tension among some data samples. Complementing this piece by the non-HLS part of the HVP derived by more usual means, one gives our best evaluation of the full HVP and our estimate of the muon anomalous magnetic moment. Comparison is done with other published evaluations. Finally, Section 22 summarizes our results and conclusions.

Outline of the Hidden Local Symmetry (HLS) Model
Let us briefly outline the derivation of the HLS model 3 . The usual ChPT Lagrangian [6,7] can be written in two different manners [18,30] : where f π (=92.42 MeV) is the pion decay constant and : when working in the so-called unitary gauge which removes a scalar field term in the definition of ξ R/L (x); P (x) is the pseudoscalar (PS) field matrix. The hidden local symmetry which gives its name to the HLS Lagrangian has transformation properties which affect ξ R/L (x) letting U unchanged. Ignoring for the moment the weak sector to ease the discussion, the HLS Lagrangian is derived by replacing in Equation (2) the usual derivative by the covariant derivative : where A µ is the photon field, Q = Diag[2/3, −1/3, −1/3] the quark charge matrix and V µ is the vector field matrix; the expressions for P and V are the usual ones -fulfilling the U(3) flavor symmetry -and can be found in [18,30,25]. In the expressions for D µ ξ R/L , the universal vector coupling g occurs beside the unit electric charge e. The ω and Φ fields occuring in the diagonal of the V matrix, also denoted V I below, are the so-called ideal combinations generally denoted ω I and Φ I . Substituting the covariant derivative Equation (4) to the usual derivative in Equation (2), L chiral becomes the first HLS Lagrangian piece denoted L A while another piece L V shows up which vanishes in the reversed substitution D µ ⇒ ∂ µ : The full HLS Lagrangian is then defined by : where a is a parameter specific of the HLS approach [18]. The usual VMD framework is obtained by setting a = 2 [30] while the presently reported phenomenology prefers larger values [24,29]. The explicit expression for this unbroken L HLS can be found fully developed in [30]. Let us note that L HLS contains a photon mass term; however this Reference also showed that the loop dressing (in the HLS context) cancels out the physical photon mass and thus, this mass term can be safely ignored in phenomenological studies. So, the unperturbed (i.e. unbroken) HLS Lagrangian depends on only two parameters (g and a) which can be ajusted using data; it is obviously unrealistic to expect this basic L HLS to describe precisely the large amount of data covering its scope with so few parameters. To have any chance to account for real data, additional input should enter the HLS Lagrangian to make it more flexible; this is what motivates the introduction of breaking procedures within the HLS framework.
As already noted, both HLS Lagrangian pieces fulfill a U(N f ) × U(N f ) symmetry and not SU(N f )×SU(N f ). To reduce this symmetry which introduces a ninth PS meson, one includes [31,25] the 'tHooft determinant terms [32] which break the axial U(1) symmetry; this turns out to complement the L A Lagrangian piece with : The singlet mass term manifestly breaks nonet symmetry in the PS sector. Nothing analogous affects the vector sector.

Breaking the HLS Lagrangian I : The BKY Mechanism
Up to now, a single breaking mechanism for the HLS Lagrangian has been proposed [28]; it is named here BKY after its proponent names. However, undesirable properties of the original BKY proposal have led to modifications of this breaking scheme [33,34] originally proposed to break (only) the flavor symmetry of the Lagrangian. The "new scheme" variant from [34] has been analyzed with data and its results shown to compare fairly well to ChPT expectations [31]. It has been extended to include isospin breaking effects in [25] following the lines of [35]. This (modified) BKY breaking scheme, intended to cover isospin and SU(3) breaking effects within the same framework, has been proved fairly successful, allowing for a high quality global fit of (almost) all available data covering the validity range of the HLS model [25,26,27] : e + e − annhiliations, τ decay spectra and radiative decays of light flavor mesons; in particular, the dipion spectra in the η/η ′ → π + π − γ decays are fairly well predicted [23,24], especially the detailed form of the drop-off in the ρ − ω interference region as recently measured at BESIII [36,37]. Nevertheless, as already stated in [25,26,27], the description of the threshold region in the dipion spectrum and of the φ mass region in the three pion annihilation channel deserves improvements; this issue is addressed by the present paper.
As BKY is one of the breaking mechanisms used in the present work, let us briefly remind how it is implemented. To lighten writing, we use the notations L = D µ ξ L ξ † L and R = D µ ξ R ξ † R . The (modified and extended) BKY breaking which also underlies the broken HLS (BHLS) model developed in [25] turns out to modify Equations (5) as follows : where X A/V =Diag(q A/V , y A/V , z A/V ) are constant real matrices. In practice, one prefers setting q A/V = 1 + (Σ A/V + ∆ A/V )/2 and y A/V = 1 + (Σ A/V − ∆ A/V )/2. As z A/V are affecting the ss entries, their departure from 1 can be large compared to q A/V and y A/V -which refer to resp. the uu and dd entries; so, the Σ's and ∆'s are expected small. Within the previous broken HLS framework -hereafter named BHLS - [25,26], one got z A ≃ [f K /f π ] ≃ 1.5, z V ≃ 1.2 while the Σ's and ∆'s were found at the few percent level. Once these breakings are applied, the PS kinetic energy term contained in L A is no longer diagonal and a PS field redefinition has to be performed in order to restore the kinetic energy term to canonical form. The procedure is fully described in Section 4 of [25] and is valid unchanged in the present work; it is not repeated here.
In the BHLS model [25], after the X V breaking, the vector meson mass term is no longer diagonal and writes (m 2 = ag 2 f 2 π ) : discarding the unessential K * sector and neglecting terms of degree greater than 1 in the breaking parameters. This mass term can be diagonalized by a 45 • rotation, a quite unacceptable solution as it does not give a smooth limit in the symmetry limit ∆ V → 0. A smooth solution is obtained by the following transform to renormalized (R) fields : where h V is a parameter submitted to fit, together with z V , Σ V and ∆ V . This solution, adopted in [25,26,27], provides a quite satisfactory description of the data. This transformation introduces a ∆ V ∂ρ R ∂ω R term in the kinetic energy term of the vector mesons; this issue is known to imply the occurence of wave-function renormalization factors absorbed in the effective couplings [38]; in the case of BHLS, they can be considered absorbed in the breaking parameters. The full account of the BHLS model also requires the dynamical breaking [29,25] generated by PS loop effects which also calls for an additional s-dependent renormalization step. It is rediscussed below within the context of the present work.
4 Breaking the HLS Lagrangian II :: The Covariant Derivative (CD) Handle As clear in Equations (8), BKY breaks the symmetry of the HLS Model in the very definition of its Lagrangian. This is quite legitimate and strongly validated by its remarquable description of a considerable amount of data [25,26,27]. Nevertheless, it is of interest to explore other possibilities to improve the description of the data in some specific mass regions and, then, reduce or cancel out the additional systematic uncertainties reported in [27].
Let us focus on Equation (4) which defines the original HLS covariant derivative : The superscript I assigned to the vector field matrix V aims at reminding that the Isospin zero mesons are the so-called ideal combinations ω I and Φ I . These naturally occur in the U(3) symmetric expression of V .
The vector meson term in D µ ξ R/L can be written gV I = g a=0, 8 V a T a where the V a 's are the vector meson fields and T a , (a = 1, · · · 8) are the Gell-Mann matrices normalized such that Tr[T a T b ] = δ ab /2; T 0 is the unit matrix appropriately normalized : T 0 = I/ √ 6. Thus, the (nonet) U(3) V symmetry of HLS is obtained by : • Plugging the appropriately constructed V I into the covariant derivative D µ ξ R/L , • Assuming the universality of the vector coupling g to the external world.
Therefore, we propose a direct breaking of the covariant derivative, a new tool independent of the BKY mechanism. Such a breaking mechanism is expressed by the modified covariant derivative : where δV µ can be chosen to break the U(3) V symmetry in a controlled way. For instance, breaking solely the nonet symmetry of V I turns out to allow the coupling to the singlet V 0 µ to differ from those of the octet fields, preserving in this way the SU(3) V symmetry. It will become clearer below why a more systematic approach must be preferred. Identifying the field combinations associated with each of the canonical T a matrices, one is led to define the following components which can participate to δV µ separately or together : in terms of the usual ideal field combinations. The (free) breaking parameters ξ 0 , ξ 8 and ξ 3 are only requested to be real in order that δV µ is hermitian as V I µ itself. Clearly, δV 0 µ defines a breaking of the nonet symmetry down to SU(3) V × U(1) V , δV 8 µ rather expresses the breaking of the SU(3) V symmetry, while δV 3 µ is related to a direct breaking of Isospin symmetry in the vector sector.
One could also introduce some breaking affecting the ρ ± entries of the V I matrix 4 . The ρ ± entries are given by the sum of the Gell-Mann matrix terms T 1 and T 2 ; forcing a breaking for these entries requires two real parameters which should be equal (ξ 1,2 ) in order to preserve hermiticity. So, ξ 3 and this ξ 1,2 could summarize the whole isospin breaking effects in the vector meson side; however, one can choose ξ 1,2 = 0 as no theoretical value for the (unbroken) universal coupling g is presently available. Indeed, all vector couplings in the HLS Lagrangian could then be reexpressed in terms of g ′ = g(1 + ξ 1,2 ); in this case, all physics quantities will depend on g ′ and ξ ′ i = ξ i /(1 + ξ 1,2 ) (i = 0, 8, 3) without any other occurence of ξ 1,2 dependency. Therefore, phenomenologically, the ξ i 's and g itself are defined up to a normalization factor presently out of reach. This scaling property has obviously no consequence on physics observables like cross-sections or form factors.
So, from now on, one assumes the maximal breaking of U(3) V experimentally accessible : δV µ = δV 0 µ + δV 8 µ + δV 3 µ . Compared with the original V I entries, this turns out to modify only the diagonal entries of V I by the following substitutions 5 : Then, the U(3) V breaking of the covariant derivative generates a breaking of the vector coupling universality. For this purpose, one should note that a vector mixing is generated, except if ξ 0 = ξ 8 . For the sake of conciseness, this mechanism is referred to below as CD breaking.

BHLS 2 : A New Broken HSL Lagrangian
Here, we define a new version of the broken HLS Lagrangian, merging the two breaking schemes just presented -BKY and CD -as these two ways of breaking are conceptually independent. Their contributions are thus complementing each other. In order to fully take into account the electroweak sector, one should modify correspondingly the covariant derivative to : where, as usual [18] : discarding the Z 0 boson terms of no concern for the phenomenology we address. The T ± matrices are constructed out of the matrix elements V ud and V us of the Cabibbo-Kobayashi-Maskawa matrix and can be found in [18]. Finally, the weak coupling g 2 is related to the Fermi constant by g 2 = 2m W G F √ 2. The expression for the non-anomalous HLS Lagrangian pieces given in Section 2 remains formally valid, being understood that the covariant derivatives are modified according to Equations (14)(15). As a result, the L V piece substantially differs from its partner in [25] while the L A piece in the present scheme is strictly identical to those given in [25]. The non-anomalous BHLS 2 Lagrangian is given expanded in Appendix A.

The order O(p 4 ) Terms of the HLS Lagrangian
Beside the L A and L V pieces which are O(p 2 ), the HLS approach also possesses terms of order O(p 4 )(see Section 4.3 in [18]), which modify the V − A/W couplings in a specific 5 Remind that a 1/ √ 2 term is factored out in the definition of V I . way. As the role of such terms has never been really examined in phenomenology 6 , it looks worthwhile examining their relevance when dealing with data of high accuracy. Using L µ and R µ just given, one first defines : and also : where Defining also : the O(p 4 ) Lagrangian writes [18] : where one has generically defined : and where the z i are constants not theoretically constrained. The effect of the z 1 and z 2 terms can be absorbed by a redefinition of the electromagnetic and weak couplings e and g 2 . The z 3 term is more involved in the phenomenology one addresses and writes : at lowest order in the expansion of the ξ L/R fields; here one has kept the unbroken V µ matrix for clarity. By integrating by part and fixing the gauge condition to ∂ µ X µ = 0 for all vector fields, this piece becomes : Thus, the L z 3 piece exhibits s-dependent parts for the V µ − A µ and V µ − W µ transition amplitudes which complement their constant parts provided by the usual O(p 2 ) HLS Lagrangian. The 'broken' version of Equation (22) is given in Appendix A.3.

The Basic Solution (BS) and the Vector Mass Term
In the context where both the BKY and CD breaking mechanisms are at work, the vector mass term in L V becomes : 6 See, however, the discussion in [23] which is revisited here.
Ignoring the K * sector and setting m 2 = ag 2 f 2 π , it writes : 24) Therefore the CD breaking, via its δV 3 component, allows for a ρ 0 − ρ ± (HK) mass difference, provided fits favor non-vanishing ξ 3 values. Moreover, all the components of δV defined above contribute to generate different HK masses for the ρ 0 and ω mesons as ∆m 2 This contrasts with the previous BHLS model [25] where m 2 ρ 0 = m 2 ρ ± = m 2 ω is fulfilled. Because of the non-vanishing HK mass difference between the ρ 0 and ω, one can already expect an improved treatment of the dipion threshold and spacelike regions in the BHLS 2 framework compared to BHLS.
Equation (24) is manifestly diagonalized by setting : This solution -which lets z V unconstrained-defines our Basic Solution (BS). Within the framework of this solution, the breaking parameters to be derived from fits are Here one feels the issue with assuming solely nonet symmetry breaking; indeed, as this turns out to state ξ 8 = 0, the CD breaking implies ξ 0 = 0 and thus breaking HLS via the covariant derivative intrinsically implies that nonet symmetry cannot be solely broken.
At leading order in breaking parameters, the vector meson mass term in L V becomes diagonal and one has : So, BHLS 2 a priori yields different HK masses for all the vector mesons.

The Primordial Mixing (PM) and the Reference Solution (RS)
Another unused mechanism can be invoked; indeed, independently of the BKY and CD mechanisms just defined, one can always consider that the neutral vector fields ρ 0 , ω, Φ involved in physical processes are not directly the ideal ones but combinations of these. For this purpose, let us define an infinitesimal rotation matrix R(U 3 ) = 1 + O(ǫ), which associates to the ideal field vector V I = (ρ I , ω I , Φ I ) a (first step) renormalized field vector This is also a quite legitimate tool to extend the flexibility of BHLS 2 , , up to terms of order O(ǫ 2 ) which are discarded in our O(ǫ) approach. As it is a rotation, this transformation preserves the canonical structure of the vector field kinetic term provided by the Yang-Mills Lagrangian (up to O(ǫ 2 ) terms).
However, the (real) ψ (Euler) angles from tranformation in Equation (27) should be chosen in such a way that L mass (V I ) remains canonical in the change of fields V I → V R , i.e. the crossed terms in L mass (V R ) should be cancelled out. Restarting from the mass term given by Equation (24), three conditions should be fulfilled : The last bracketed terms in these expressions can be discarded as clearly of second order in the breaking parameters {ξ i , ψ j }; this already implies that ∆ V = 0 and that ψ ω becomes unconstrained.
As for the parameter z V -generated by the BKY breaking in its ss entry -the situation deserves further comments : • Either z V is found such that 7 1 − z V ≃ O(ǫ), then none of the ψ angles is constrained at order O(ǫ) and, moreover, the last two Equations (28) then imposes ξ 0 = ξ 8 .
• Or z V is returned by fits such that 1 − z V is not O(ǫ). Then, at first order in breakings, the last two equations imply ψ Φ = ψ 0 and that they are proportional to ξ 0 − ξ 8 .
Preliminary fits, performed using this mechanism -named from now on -Primordial Mixing (PM) -indicate that 1 − z V is in the range (1 ÷ 2)ǫ and thus, in this case also ξ 0 = ξ 8 can be imposed. Therefore, beside the Basic Solution (BS), one gets an additional one, we name it Reference Solution (RS), which also includes the Primordial Mixing (PM). Also considering the former BHLS model [25,26,27], one thus has at disposal three different models. This allows for a better insight into the systematics and model dependence effects.

Vector Meson Propagators And Form Factors
Let us anticipate on form factor calculations using the Lagrangian given in Appendix A and have a look at the form factor values at s = 0 -before introducing mixing effects due to loop corrections. For this purpose, let us start with a preliminary digression on vector meson propagators, especially those for the ω and Φ mesons, as these play an important role at the chiral point and in the physics of the (close) spacelike region which is also addressed in the present work.
For the ρ meson, following the pioneering work [39], the inverse propagator at one loop can be written : [23,24,25] where the self-energy Π(s) is the sum of the loops allowed by the non-anomalous Effective Lagrangian given in Appendix A; in our context, these are essentially pion and kaon loops 8 . Π(s) is a real analytic function (Π(s * ) = [Π(s)] * ) which vanishes at s = 0 because of current conservation, and is real for negative s -in fact, this property holds already below m 2 π 0 , the lowest energy hadronic threshold. Hence, this implies that D ρ (0) = −m 2 ρ , where m 2 ρ is displayed in Equations (26).
In principle, this also applies to the ω and Φ propagators; however, taking into account their narrow character, it looks unmotivated to enter into such complications when fitting objects like e + e − annihilation cross sections. In this case, phenomenology has lengthily experienced a broad success using Breit-Wigner (BW) lineshapes. However, as discussed in [25], some physics quantities like the contribution of the e + e − → π + π − cross-section to the estimate for the muon HVP, may require some care about the behavior of approximations close to the chiral point, as this region gives an enhanced contribution to the muon HVP. For this purpose, [25] proposed the following modified BW lineshape for the ω and Φ mesons : where ( m V , Γ V ) are parameters to be fitted and the m 2 V 's are the relevant HK square masses as displayed in Equations (26). Numerically, this expression for the D V (s)'s is close to usual BW lineshapes 9 . This BW-like parametrization imposes the ω and Φ meson inverse propagators to fulfill D V (0) = −m 2 V as they should -and as already does the ρ meson inverse propagator; this property of the propagator at s = 0 is essential to recover the expected values for the pion and kaon form factors at the origin -within the HLS framework.
On the other hand, if one wishes to examine the analytic continuation of the ω and Φ propagators somewhat inside the spacelike region, it might be desirable to make it real there; this can be achieved by simply replacing Γ V in Equation (29) by θ(s) Γ V -or rather, θ(s − m 2 π 0 ) Γ V . 8 Loops generated by the anomalous HLS Lagrangian pieces [18,22] also contribute but can be discarded [25] as reminded in the next Section; however, some anomalous loops, suppressed by a factor of e 2 , can have to be kept as they may play some role in specific expressions as will be emphasized below. 9 Setting m V = m V in this approximation gives D V (s) = s − m 2 V + is Γ V /m V ; more common choices for the imaginary part here are m V Γ V or √ s Γ V , quite analogous to our own choice.
Using the Basic Solution (BS) or the Reference Solution (RS) to redefine the physical vector meson fields, one can show, using the Lagrangian given in the Appendix A, that : up to terms of order O(ǫ 2 ). One can check that the conditions ∆ V = 0 and ξ 0 = ξ 8 are essential in the derivation of these constraints. So, the breaking parameters to be derived from fits are Σ V , ξ 3 , ξ 0 (= ξ 8 ) and z V when working within the BS framework; one should also include the ψ rotation angles when extending BHLS 2 to the Reference Solution.
It is worthwhile reminding that the effects of δV 0 µ and δV 8 µ are intimately intricated within our new Lagrangian frameworks. Here, indeed, a breaking of solely nonet symmetry (i.e. ξ 0 = 0) cannot occur if not accompanied by a breaking of the SU(3) V symmetry of the same intensity.
In conclusion, the present model (BHLS 2 ) is not a trivial variant of BHLS [25]. Some important properties of BHLS 2 versus BHLS mentioned below, will further substantiate this statement. As previously noted [29,25] and reminded above, all variants of the HLS Model, including BHLS 2 (see Appendix A.1), exhibit ρ 0 /ω/Φ → KK couplings. This implies that, at one loop, the squared mass matrix of the ρ 0 /ω/Φ system receives non-diagonal entries, i.e. the vector fields occuring in the Lagrangian Eq (116) are no longer mass eigenstates. Mass eigenstates are constructed using perturbative methods as performed in [29,25]. As the loops are (real analytic) functions of s, the mass eigenstates become also s-dependent; actually, it is this property which has given the solution [29,25] to the long standing e + e − versus τ puzzle [40,41,42].
Basically, the Dynamical (or Loop) Mixing of Vector Mesons has been first defined 10 and motivated in [29]. In order to ease the reading, we remind it and emphasize the new features provided by the BHLS 2 context. The (squared) mass matrix of the ρ 0 /ω/Φ sector can be written : where : The Higgs-Kibble masses m V are displayed in Equations (26), and Π ππ (s) is the pion loop weighted by the square of the ρ 0 π + π − coupling constant. Because all the loop functions are real analytic function of s, M 2 0 (s), δM 2 (s), and hence M 2 (s), are hermitian analytic matrices (e.g. fulfilling [X(s)] † = X(s * )).
In order to construct explicitly δM 2 (s), it is worth reexpressing some V P P coupling constants in a suitable manner. One can write : having defined G = ag/(4z A ). Using Equations (25), common to both the BS and RS variants, one finds : and : When dealing specifically with the BS variant, these two sets of equations are used by simply dropping out the ψ α parameters. Denoting by resp. Π ± (s) and Π 0 (s) the amputated charged and neutral kaon loops, the V i R → V j R transition amplitudes (i, j = ρ 0 , ω, Φ) are given by : using the notations just defined. Then, the elements of the δM 2 (s) matrix are : As in the former BHLS [25], δM 2 (s) is non-diagonal for non-zero values of s. Therefore, at one-loop order, the field combinations defined in both the BS and RS variants are not mass eigenstates, as in the previous BHLS release. This calls for a mass-dependent diagonalization which is reminded in the next Subsection.
Let us make a few more remarks about the non-diagonal entries of δM 2 . In the no breaking limit, where the neutral and charged kaon loops coincide, the ρ 0 − ω and ρ 0 − Φ entries identically vanish; however, the ω − Φ entry, proportional to the sum of the neutral and charged kaon loops, does not vanish, indicating that the ω − Φ mixing is always at work at one loop order within the HLS framework.
In the former BHLS, however, if no breaking is applied in the L A sector and if there is no mass breaking among the kaons, there was no ρ 0 − ω and ρ 0 − Φ transitions. This remains true within the BS framework defined above. However, within the RS framework, the ψ α 's come in such a way that none of the entries of δM 2 (s) vanishes.

. . . And its Diagonalization
The diagonalization of L mass at one loop order is performed by means of another (sdependent) rotation matrix to a second step of vector field renormalization R ⇒ R ′ : so that the full renormalization procedure is obtained by : The R ′ fields are the physical vector meson fields. When working within the Basic Solution variant, V R ≡ V I and R(U 3 ) can be dropped out from Equation (39); otherwise, within the RS framework, R(U 3 ) is given in Equation (27). The complex "angles" in Equation (38) can be derived from the δM 2 (s) matrix elements : where the numerators manifestly depend on the kaon loops only and the λ's are the eigenvalues of M 2 (s) matrix. At leading order in breaking parameters, these write [25] : the vector meson masses being those displayed in Equations (26). As the three m 2 V 's occuring here are different, the three angles defined above vanish at s = 0; in constrast, within the former BHLS context [25], because the HK masses for ρ 0 and ω were equal, α(0) has a non-zero limit at s = 0, compromising a proper analytic continuation of the form factors downwards into the spacelike region 11 .
In summary the full transform from ideal to fully renormalized ρ 0 , ω, Φ fields is given either by (BS solution) : or by (RS solution) : by discarding breaking terms of order greater than 1. We will refer from now on to Equation (43) for both the BS and RS variants, being understood that in the former case, the ψ i 's are zero. 11 See, in particular, the comments in Section 6.3 of [25].
Moreover, in the following, the fully renormalized fields can be either indexed by R ′ -if useful -or deprived of indexation to lighten expressions. For clarity, in all displayed Lagrangian pieces, one always writes the couplings for the ideal ρ 0 , ω, Φ fields. To go to physical fields, one has to apply the appropriate transformation R(I ⇒ R ′ ) and collects the contributions in order to yield the physical ρ 0 , ω, Φ couplings.

The Effects of Anomalous Loops
When defining the loops contributing to δM 2 (s), we have only considered those generated by the non-anomalous HLS Lagrangian. However, the anomalous (FKTUY) HLS Lagrangian pieces [22,18] allow for (anomalous) couplings also generating loop contributions 12 to δM 2 (s) as, for instance, K * K or 3-pion loops. As all the two particle channels have thresholds above our fitting region, their loops -of order O(g 2 ) -are real in our fitting range and, thus, can be discarded, assuming their effects numerically absorbed by the subtraction polynomials of the other (pion and kaon) loops involved [25].
However, one should also note that the anomalous (FKTUY) sector of the HLS Lagrangian generates couplings of the neutral vector mesons to the π 0 γ, ηγ and η ′ γ final states, even when setting c 3 = c 4 [25,27]. Then the FKTUY sector generates the corresponding loops, which have couplings of order O(g 2 e 2 ) and develop tiny imaginary parts far inside our fitting region for the first two and close to its upper border for the third one. These contributions, and their imaginary parts, are higher order and can generally be discarded; nevertheless, their tiny imaginary parts were accounted for in our previous studies [25,26,27] by simply adding some fixed iε to the eigenvalue differences λ i (s) − λ j (s) which occur in the denominators of the mixing "angles" shown in Equations (40).
So, these (complex) mixing "angles" exhibit a dependence upon differences of the M 2 (s) eigenvalues. However, these eigenvalues -see Equations (41) -being s-dependent, their differences may exhibit zeros for real values of s which, accordingly, generate real poles for the mixing "angles". Poles occuring at real negative s are not a real issue as they can be dropped out by means of customary dispersive technics. For the others, the iε trick emphasized just above, permits to avoid unwanted singularities on the physical region s ≥ m π 0 .
In the previous BHLS version [25], the main purpose of this ad hoc iε was to force 13 α(s = 0) = 0, which is no longer necessary in the present version of BHLS. A more natural way to proceed -giving the same results -is, however, to take into account the π 0 γ loop which introduces a small imaginary part to the eigenvalue differences, for s 0 ≥ m 2 π 0 upwards. For this purpose, we use the loop expression derived in [43] (see its Appendix C) which can be written (s 0 = m 2 π 0 ) : 12 The Yang-Mills part of the full Lagrangian may also contribute with other loops like K * K * . The V P P P Lagrangian generate 2-loop contributions with coefficients of O(ǫ) or O(ǫ 2 ) to the transitions among the vector mesons; they are presently ignored. 13 As clear from the first Equations (40), the numerator of α(s) behaves as s close to the origin and, as m ρ 0 = m ω in the previous BHLS [25], λ ρ (s) − λ ω (s) exhibits a similar s-behavior.
Along the real axis in the complex s-plane, one has : The derivation of the g R V π 0 γ couplings is given in Subsection 12.2. The constant term inside the bracket in Equation (44) ensures that ǫ π 0 (0) = 0. Full consistency would impose to add (using obvious notations) : to δM 2 (s) as expressed above; however, because of the e 2 factor, they are higher order and can be neglected, except in the diagonal where they avoid possible real poles in the physical region. Additionally, the ǫ η (s) and ǫ η ′ (s) contributions displayed for completeness will also be discarded. In order to substantiate the effect of the π 0 γ loops, let us anticipate on the global fit information. The main effect of these loops is to prevent naturally the occurence of poles on the physical region. This effect is prominent for the "angle" α(s) which produces a ρ 0 − ω mixing; Figure 1 displays the behavior of Re[α(s)] for the BHLS 2 BS and RS variants and also, for comparison, its behavior within the former BHLS, where a iǫ trick prevents a real pole close to the 2-pion threshold [25]. Even if the overall lineshapes are similar, one observes significant differences between the various solutions ; for instance, α(s) vanishes at s = 0 within the two variants of BHLS2, as the denominator is non-zero at the chiral point, allowing for a smooth connection between the spacelike and timelike regions.
Over the range [−3.5 ÷ 1.2] GeV 2 , the behavior observed for the two other mixing angles β(s) and γ(s) is much smoother (see also [23]); nevertheless, specific parametrizations may lead to a real pole 14 for β(s) at negative s below ≃ −1.5 GeV 2 . As noted above, such poles can be eliminated by means of usual dipersive methods.

The Model Pion and Kaon Form Factors
Using the non-anomalous sectors of the HLS model, broken as pointed out in the previous Sections (see Appendix A), one can already derive some of the form factors expressions needed for our global fit framework.

The Pion Form factor in the τ Decay
The τ sector of BHLS 2 is almost identical to that of BHLS [25]; indeed the main L τ piece reminded in Appendix A.2 is unchanged, but the newly introduced L z 3 contributes. This yields : with : where z 3 is a new HLS parameter [18] introduced by the O(p 4 ) terms of the HLS Lagrangian (see Appendix A.3). Π τ ρρ (s) is the loop correction to the ρ ± propagator and Π W (s) the loop correction to the W ± − ρ ± transition amplitude; both are defined just below.
Using the following short-hand notations : one first defines the pion and kaon loop contributions to the ρ ± self-energy Π τ ρρ (s) : where ℓ τ π (s) and ℓ τ K (s) denote the amputated loop functions for resp. π ± π 0 and K ± K 0 , each having absorbed only a factor of resp. G 2 π and G 2 K ; in this way, only the breaking parameters affecting the V sector of BHLS 2 are factored out and displayed. P τ π (s) and P τ K (s) are subtraction 14 One should note that a pole at such a location is exhibited by the usual Gounaris-Sakurai formula [44]. polynomials chosen to vanish at s = 0 because of current conservation. Then one has 15 : As the loop functions vanish at s = 0, one clearly has F τ π (0) = 1. One should note that the same subtraction polynomials occur in the W − ρ transition amplitude F τ ρ (s) and in the inverse propagator D ρ (s) but differently weighted.
The expression for the loop functions ℓ τ π (s) and ℓ τ K (s) along the real s-axis can be found in Subsection F.2 of the Appendix to [29]; the subtraction poynomial needed to ensure that they vanish at s = 0 has coefficients given in Equations (122) of this reference. The coefficients for the constant term and for the second degree one given therein are correct, but the coefficient for the term linear in s should be corrected to : using the notations 16 of [29].

Parametrization of Loops in e + e − Annihilations and τ Decays
The ρ 0 self-energy and the γ ⇒ (ρ 0 /ω/Φ) transition amplitudes involve important loop corrections; these already play a central role in the dynamical mixing of vector mesons as shown in Section 10. Let us define a parametrization common to the computation of the "angles" α(s), β(s) and γ(s), to the self-energies and to the A ⇒ V transition amplitudes as these are involved in all form factors addressed in our global fitting code.
One first defines the following loops : where ℓ e π (s), ℓ e K ± (s) and ℓ e K 0 (s) are resp. the π + π − , K + K − and K 0 K 0 loops. These are defined as the amputated loops, including only a factor of G 2 π for the former loop and of G 2 K for the two others (see Equations (49)). The expression for these loops along the real s-axis can be found in Section F.1 of the Appendix to [29]. The loop parts being defined this way, P e π (s), P e K ± (s) and P e K 0 (s) are the subtraction polynomials, chosen of second degree and vanishing at s = 0.
However, as already observed in Section 10, besides the pion loop, it is rather the sum and difference of the kaon loops which are relevant : 15 Here and in the following, one can be led to keep breaking parameter dependencies higher than first order to simplify the expression writing. 16  denoting P e S (s) and P e D (s) their respective subtraction polynomials. As already noted, one can consider Π e K dif f (s) and P e D (s) to be first order in breaking (O(ǫ)). So, it is certainly more appropriate to parametrize P e S (s) and P e D (s) for fits and propagate them into the subtraction polynomials of Π e K ± (s) and Π e K 0 (s) by stating : On the other hand, as one certainly has : it is appropriate to impose : and submit to fit the coefficients of the δP τ (s) polynomial rather than those for P τ π (s) directly. This turns out to attribute the full breaking in e + e − versus τ to the pion loop subtraction term; as our fitting region is mostly below the KK thresholds, this looks a safe assumption. The second relation, additionnally, is an assumption which allows to reduce the number of free parameters in the minimization procedure by two units without any degradation of the fit quality.
As a general statement, all our subtraction polynomials have been chosen of second degree and vanishing at s = 0. An exception is made , however, by choosing the third degree for δP τ (s) which is found to provide a significantly better description of the τ spectra, especially for the Belle dipion spectrum [45].

The Pion Form factor in the e + e − Annihilations
Because of the vector meson mixing, deriving the pion form factor in the e + e − annihilation is slightly more involved than in the τ decay. To start, one needs to propagate the transformation in Equation (43) into the Lagrangian in Appendix A.1 and collect all contributions to the fully renormalized ρ, ω and Φ field couplings. The derived pion form factor in e + e − annihilations F e π (s) writes : Using the ρ I π + π − coupling in the Lagrangian displayed in Appendix A.1 and the I ⇒ R transform given by Equation (43), one readily gets : keeping only the leading O(ǫ) terms in breaking. So, BHLS 2 generates ω and Φ direct couplings to ππ. The inverse ρ 0 propagator writes : where m 2 ρ 0 is defined in Equations (26) and Π e ρρ (s) is the self-energy of the physical (i.e. fully renormalized) ρ 0 . This can be defined in terms of the loop functions constructed in the Subsection just above. Π e ρρ (s) is made up of two pieces : where the couplings g ρ K± and g ρ K0 are given in Appendix A.4. A term g ρ 2 P e D (s)/2 has been omitted in the expression for Π e K (s) as being O(ǫ 2 ). As for the inverse propagators D ω (s) and D Φ (s), one uses the expression given by Equation (29) with the appropriate values for HK masses m 2 V , and having ([ m V , Γ V ], V = ω, Φ) as parameters to be determined by fit.
The loop corrected V − γ transitions amplitudes F e V γ (s) are given by : where the new-comer z 3 s terms occur as in F τ ρ (s) above. Using together the ideal couplings Equations (A.3), the renormalization matrix Equation (43) and the RS conditions (see Section 8), one derives : which become s-dependent thanks to the dynamical vector mixing. Similarly, using the couplings which can be read off L z 3 in Equation (119) and the tranformation matrix given by Equation (43), one also derives : neglecting terms of breaking order greater than 1. Finally, the relevant loop transition terms Π V γ (s) are collected in Appendix A.5.

The Charged and Neutral Kaon Form Factors
The charged and neutral kaon electromagnetic form-factors can also be derived from the BHLS 2 Lagrangian given in Appendix A.1 : The experimental spectra measured in the spacelike region, in particular by the NA7 Collaboration [46,47], being squared form factors, the modulus squared of the form factors given in Equation (57) and in Equations (64) apply directly. In the timelike region, experimentalists rather publish cross-sections : where q P = s − 4m 2 P /2 is the meson momentum in the center-of-mass system (P = K + , K 0 ). In the case of the K + K − final state, the cross-section should take into account the significant interaction between the charged kaons emerging from the e + e − annihilation. Traditionally, this is done by multiplying the e + e − → K + K − cross-section by the leading term describing the Coulomb interaction as formulated by [48,49] : The data having become more and more precise, it looks more convenient to use the exponentiated expression derived in [50] which includes also the Final State Radiation (FSR) correction factor; indeed, the K + K − spectra provided by the various experiments are not amputated from FSR effects.

Form Factors in BHLS versus BHLS 2
It looks worth pointing out the difference between the form factor expressions in BHLS and in BHLS 2 . As already noted, the L A Lagrangian is common to both frameworks. In contrast, the L V Lagrangian is significantly different. Obvious differences for L V have already been noted. For instance, within BHLS 2 , ∆ V = 0 and its sharing parameter h V drop out whereas they are a key ingredient in the diagonalization procedure of BHLS.
Three parameters 17 arise from the Covariant Derivative breaking (ξ 3 , ξ 0 , ξ 8 ). Moreover, a key role is played by the ψ rotation angles within the RS solution of BHLS 2 . None of these parameters occurs in the BHLS framework of [25].
The effect of a significant ξ 3 value is obviously equivalent to shifting apart the Higgs-Kibble masses of both ρ's and their couplings to ππ in a correlated way. As this kind of exercise performed within BHLS concluded to a small effect [25,26], one expects a small value for ξ 3 . Also, BHLS 2 deals with different subtraction polynomials -the δP τ (s) function already referred to above -for the π ± π 0 and π + π − loops, a new (loop) mechanism able alone to contribute to physical mass and width differences between the ρ 0 and ρ ± mesons.
One should also note another difference : In the present framework, as already noted in Subsection 11.1, the pion and kaon loops entering the amplitudes F τ ρ (s) and F e ρ (s) carry fitted subtraction polynomials already involved in the charged and neutral ρ inverse propagators D ρ (s). In BHLS, one allowed the subtraction polynomial in F e/τ ρ (s) to carry a piece independent of both D ρ (s)'s. Because of the O(p 4 ) term still introduced, it means that a second degree term in s has been removed from F e/τ ρ (s) in the BHLS 2 framework. We should see below the effect of this on the ππ P -wave phase-shift above the Φ mass.

The Anomalous Sector of BHLS 2
As in [25], the anomalous FKTUY Lagrangian [22,18] is used to address the e + e − → (π 0 /η)γ and e + e − → π 0 π + π − processes, and also the radiative decays with couplings of the form V P γ and P γγ which are parts of our global fitting framework. One can display its various pieces : where, temporarily, the weight of each piece has been extracted out to exhibit its dependence upon the unconstrained c i constants of the anomalous HLS Lagrangian [22]. If one imposes to yield the Wess-Zumino-Witten terms [51,52] at the chiral point -the so-called triangle and box anomalies -the condition c 3 = c 4 is known to be mandatory [21,18]. Reference [25] having also shown that this condition is indeed favored by fits to the annihilation data, one endorses the c 3 = c 4 constraint and, thus, cancels out the last term in Equation (67). Hence, the V P γ couplings become entirely generated by the combination of V V P couplings with the V − γ transitions provided by the non-anomalous HLS Lagrangian. The surviving Lagrangian pieces displayed in Equation (67) are given by : where N c is the number of colors and the c i factors have been reabsorbed in the normalization of the various Lagrangian pieces. P and V denote the bare pseudoscalar and vector field matrices.
In the present context, one has, however : where V I µ is the (usual) U(3) symmetric vector field matrix and δV µ = δV 3 µ + δV 0 µ + δV 8 µ which has been expressed in term of the ideal field combinations in Equations (12).
For the reader convenience, Appendix B reminds briefly the renormalization procedure which leads from the bare PS fields to their renormalized partners [25].

Renormalization of the V V P Couplings
As displayed in Equation (39), referring to : the relation between the fully renormalized vector fields (denoted R here and from now on) and their ideal combination is given by : where the product R(U 3 ) · R(Loop), given by Equation (43), will be denoted from now on : Let us point out that R(s) is a real analytic 3 × 3 matrix function and fulfills at order O(ǫ) : along the right-hand cut in the physical sheet [29]. The V V P Lagrangian displayed in Appendix C can be symbolically written : in terms of ideal field combinations, P 0 being any of the π 0 , η or η ′ fields; summation over greek and latin indices is understood. The 3-vector G I (π ± ) and the 3×3 matrix G I (P 0 ) can easily be constructed using the relevant pieces of information given in Appendix C. One thus has : with the obvious (ρ, ω, Φ) component indexing. On the other hand, one has defined : where, for each P 0 , g ij is the coupling of the form P 0 V i V j which can be read off the relevant expressions given in Appendix C.
Performing the replacement V I = R(s)V R into Equation (72), one can derive the fully renormalized V V P Lagrangian : where, the coupling vector G R (π ± ) and the coupling matrix G R (P 0 ) inherit from R a sdependence; they are related to their ideal partners by : which is a concise way to express the renormalized couplings.

The Renormalized AVP Effective Lagrangian
As already noted, setting c 3 = c 4 turns out to cancel the direct FKTUY AV P Lagrangian piece; hence, all AV P couplings inside BHLS 2 are generated by V V P couplings followed by one V → A transition. Let us consider the V −A transition term of the non-anomalous BHLS 2 Lagrangian (see Equation (116)). Having defined the transposed vector , one can rewrite it : which defines the 3-vector f R V γ as : The effective L AV P Lagrangian can be derived from Equation (75) by replacing one neutral vector meson, say V R i , using the rule : where m 2 i is the renormalized squared mass of the V i vector meson as given in Equations (26); let us note that each H i γ carries a hidden factor of 1/g. Starting with the charged pion part of Equation (75), one gets : where H γV is the 3-vector constructed from the H i γ defined just above : For the P 0 part of the Lagrangian Equation (75), the algebra is slightly more involved and gives : This effective Lagrangian is the tool to parametrize the e + e − → (π 0 /η)γ cross-sections and the V P γ radiative decays.

12.3
The e + e − → (π 0 /η)γ Cross-Sections Using the AAP and AV P Lagrangians derived above in terms of renormalized vector and pseudoscalar fields, one can construct the amplitudes for γ * → (π 0 /η)γ transitions. These can be written : where the (co)vector H γV is defined in Equation (81), the G R (P 0 ) matrix in Equation (76) (and Equation (74)). One has also defined : the functions F R γV (s) and D V (s) having been defined in Subsection 11.3. Finally, one has : the g P 0 γγ , derived in [25], are reminded in Appendix D.
In the case of the e + e − → (π 0 /η)γ annihilations, the available data are always crosssections which are related to F P 0 γ (s) just defined by : 12.4 The e + e − → π 0 π + π − Cross-section The amplitude for the γ * → π 0 π + π − transition involves most of the FKTUY Lagrangian pieces; it can be written : labeling each term by the particular piece of the FKTUY Lagrangian from which it originates. As already noted, because c 3 = c 4 is assumed, there is no T AV P piece.
The T AP P P contribution to the full T (γ * → π 0 π + π − ) amplitude can be derived from the information given in Appendix D; it writes : where one has defined : Three pieces are coming from the V P P P Lagrangian also displayed in Appendix E; they are collected in : where the renormalized vector couplings g R V π (s) to 3 pions have to be derived from using Equations (148) and Equation (151), R(s) being given by Equation (43). The V −γ amplitudes F e V γ (s) have been constructed in Subsection 11.3 and the inverse propagators in Subsection 11.3 for the ρ 0 meson. The inverse propagators for the ω and Φ mesons have been discussed and defined in Section 7. One has also defined : The V V P Lagrangian piece in Equation (140) has been rewritten in terms of renormalized fields in Equation (75) with (renormalized) couplings derivable by means of Equation (76). The simplest way to write down T (γ * → π 0 π + π − ) in a way easy to code within our global fit procedure is displayed just below.
One first defines the H i (s) functions : where s is the incoming squared energy and the s ij 's indicate the invariant mass squared of the corresponding outgoing (i, j) pairs. T V V P depends on the 5 functions (H i (s), i = 1 · · · 5) with sdependent coefficients F i (s) given below.
Collecting all terms, the full amplitude writes : with : In this way, to write down the full amplitude, the various F i (s) functions only depend on the incoming off-shell photon energy squared s; the dependence upon the various sub-energies s ij is, instead, only carried by the H i (s) functions as clear from Equations (93). One has : One should note, as in BHLS [25] and earlier in [23], that all parameters already fitted with the (five) e + e − → (ππ/KK/π 0 γ/ηγ) processes fully determine all the (F i (s), i = 1 · · · 5). The only new parameter coming with e + e − → π 0 π + π − is c 1 − c 2 , only affecting F 0 (s); for practical purposes, it is handled in a specific manner, as in the previous versions of the (broken) HLS model just quoted.
Indeed, a global fit to all cross-sections but e + e − → π 0 π + π − allows to yield the relevant parameters with a good approximation; thus, having at hand all ingredients defining the F i (s)'s, a first minimization run including the e + e − → π 0 π + π − cross-section can be performed to also derive a first estimate for c 1 − c 2 . The output of this MINUIT minimization run is then used as input for a next minimization step. This initiates an iteration procedure involving all crosssections which is carried on up to convergence -when some criterium is met.
This method has been proven to converge in a very few minimization steps [25]. What makes such a minimization procedure unavoidable is that the e + e − → π 0 π + π − cross-section expression implies to integrate over the s ij 's at each step. This is obviously prohibitively computer time consuming for a negligible gain. Hence, at each step, one starts by tabulating coefficient functions, exhibited in the next expression between curly brackets : and these tables are used all along this step. Equation (97) uses the Kuraev-Silagadze parametrization [53] and its kernel G(x, y) function; these are reminded in Appendix H of [25].

The Data Samples Submitted to Global Fits
Basically, 48 experimental data samples are presently available which enter the scope of the HLS Model. Relying on the global fits performed in [25,26], one has been led to discard a few of them from our fitting procedure; this situation happens again within the present framework. Before going on, let us remind the most important part of the data samples and list the newcomers. Newcomers encompass either newly collected samples or existing form factor spectra covering spacelike momenta not treated within the former BHLS framework [25].
Substantially, the set of available data samples includes : • i/ For the π + π − annihilation channel : The earliest data samples date back to the beginning of the eighties [54,55]. These have been followed by quite precise data samples collected in scan mode by CMD-2 [56,57,58] and SND [59] on the VEPP-2M collider at Novosibirk. The published spectra cover the energy region from 370 MeV to 970 MeV; presently, there is still no published spectrum covering the Φ(1020) energy region.
These measurements, referred to globally as NSK, have been followed by large statistics data collected using the so-called initial state radiation (ISR) mechanism (see [60], for instance). In this way, three data samples [61,62,63] have been collected with the KLOE detector running on the DAPHNE collider up to ≃ 1 GeV, just below the Φ meson mass.
In the same period of time, another π + π − sample has been collected by BaBar [64,65] covering the energy range up to 3 GeV. By 2015, BESIII has also published a π + π − sample [66] and, by end of 2017, a group from the CLEO-c Collaboration has published the latest to date π + π − sample [67] collected with the CESR detector.
As a whole, the most precise data samples represent now 8 spectra. One has been led to discard from our global framework KLOE08 [61] and BaBar [64,65] for different reasons 18 discussed in [26,27]. A reanalysis unifying the three KLOE data samples has been recently performed [69,70] and a 85-data point spectrum has been derived which will be commented at the appropriate place. Moreover, the full error covariance matrix including the correlations between KLOE08, KLOE10 and KLOE12 is now available; this is used here to account for the (weak) correlations between KLOE10 and KLOE12.
• ii/ For the K + K − annihilation channel : Up to very recently, the available data samples covering this channel were the three scans collected by CMD-2 [71,72] and the two scans of SND [73]. They were included within the global fits performed in the (previous) BHLS framework (see [25,27]). Quite recently, the CMD-3 Collaboration has found that the two CMD-2 scans [72] should undergo a rescaling by [74] 1.094 ± 0.04 which represents an important correction of the spectrum absolute normalizations and of its (correlated) systematics (2.2% → 4.45%).
The influence of this correction on the the physics information previously derived within BHLS in [25,26,27] imposes the quite significant update performed in the present study.
Babar has published in 2013 the first measurement [75] of this cross section performed in the ISR mode. This high statistics spectrum covers the full energy range from threshold to 3 GeV. Moreover, very recently, CMD-3 has also produced a new high statistics measurement of this spectrum in scan mode [76]. However, as emphasized in the CMD-3 publication, this measurement is inconsistent with the BaBar data sample [75]. This issue is addressed below.
• iii/ For the K L K S annihilation channel : To our knowledge, the first reported measurement of the e + e − → K L K S cross section has been provided by CMD-2 in 1995 [71] and has been followed in 2004 by the 4 scans reported in [77]. In the mean time SND also produced 4 data samples , 2 in the charged decay mode of the K S meson, 2 in its neutral mode [73]. These were the data encompassed in our former studies. The present work also includes the CMD-3 data sample [78] recently published.
• iv/ For the (π 0 /η)γ annihilation channels : As for the e + e − → π 0 γ cross section, in our energy range, the available data samples are scarce : One sample has been provided by CMD-2 [79] and another one by SND [80] which superseds the former [81] used in our previous studies [25,26,27]; the older SND π 0 γ spectrum given in [82] is also considered.
Regarding the e + e − → ηγ cross section, the situation is more favorable as six data samples corresponding to different η decay modes have been collected by CMD-2 [83,84,79] and SND [82,85]. Nevertheless, no newly collected data sample for this channel has been reported.
• v/ For the π + π − π 0 annihilation channel : As no newly collected sample covering the energy region up to the Φ mass has been reported, we are left with only the samples already examined in our previous analyses [25,26,27]. These have been collected by CMD-2 [71,77,86,87] and SND [88,89] and cover either the ω peak region or the Φ peak region. We also consider the 3-pion spectrum given in [90] as it mostly deals with the region in between the ω and φ peaks, allowing for a valuable constraint on the ρ 0 meson physics background.
However, preliminary fits showing that the sample in [86] exhibits an average χ 2 per data point much above 2 -as also yielded by [91] in a different theoretical context -we have been led to discard it from our fits. Overall, one is left with nine data samples, including an old sample from the former CMD Collaboration.
• vi/ The pion and kaon form factors in the spacelike region : One of the motivations to develop a new breaking procedure for the HLS Model, is to design it in such a way that the resulting broken model (BHLS 2 ) could encompass accurately the close spacelike region of the pion and kaon form factors. Indeed, a fair account of F π (s) in a region which stretches on both sides of s = 0 is the best way to ascertain the low energy behavior of the e + e − → π + π − cross section which provides an enhanced contribution to (g − 2) µ .
The NA7 Collaboration has measured the pion form factor [46] in the spacelike range Q 2 ∈ [−0.253, −0.015] GeV 2 with a good statistics and a small reported normalization error (0.9%); two years later, NA7 has also published a measurement of the charged kaon form factor [47] [47,93] as well -report on results derived in resp. π and K scattering on atomic electrons and are, thus, model independent.
As a matter of fact, one should also note that there is no K 0 form factor spectrum measurement in the spacelike region, which would be a real challenge; nevertheless, measurements of < r 2 K 0 > exist [101,102] which report negative values.

Comments on the Model Parameters
Before going on, it is worthwhile commenting on the fit parameters of the broken HLS Models. Basically, the number of parameters of the original (unbroken) HLS Model [21,18,22] to be fixed from data is small. As for its non-anomalous sector, besides f π , these are the universal vector coupling g and the specific HLS parameter a, generally found to slightly depart from 2, the value expected from standard VMD approaches. Its anomalous sector [22] involves a few more parameters already displayed in the Lagrangian pieces presented in Section 12 : c 1 − c 2 , c 3 and c 4 . There are strong motivations [21,18] to impose c 3 = c 4 . So, when unbroken, the HLS model depends on only four unconstrained parameters.
There are two parameter categories introduced in our approach : the coefficients of the subtraction polynomials which supplement the corresponding loops and, on the other hand, the breaking parameters inherently affecting the present Lagrangian -and those defined in [25].
Traditionally, as illustrated by the Gounaris-Sakurai formula [44], the coefficients of the subtraction polynomials are reexpressed in terms of the resonance mass and width, which are determined by fit to the data. In a global approach like ours, the subtraction polynomials affect loops which come simultaneously in the various amplitudes associated with the various annihilation processes simultaneously under examination. In this case, as there is no reason to favor a specific channel to fix the meaning and the value of parameters which come in all channels, the strategy adopted is to let the data, as a whole, determine the subtraction parameter values via a global fit. Even if the parameter content may look less intuitive, this approachalso adopted in BHLS [25] -looks the most motivated. This does not prevent to extract usual physics quantities, relying on the fit results, as done since [24,25], for instance.
In the second category of parameters, besides the parameters of the unbroken HLS model reminded just above, one finds all those introduced by the various aspects of the breaking procedure at work. As for the parameters generated by the BKY mechanism and reminded in Section 3, six in total, the condition Σ A = 0 is already stated since the original BHLS Model [25]. It has been shown in Section 7 that the diagonalization procedure of BHLS 2 imposes ∆ V = 0. Therefore, four BKY parameters only (Σ V , ∆ A , z A and z V ) remain free.
In Sections 7 and 9, it was also shown that ξ 0 = ξ 8 comes as a natural constraint on the parameters of the CD breaking defined in Section 4. Therefore the CD breaking involves two more parameters 21 ξ 0 (= ξ 8 ) and ξ 3 . Finally, switching on the Primordial Mixing introduces three more parameters (ψ 0 , ψ ω and ψ Φ ) into the fit procedure.

The BHLS 2 Global Fits : General Features
Various kinds of global fits have been performed using the data samples and channels listed in Section 13, with or without the spacelike data, with or without the τ data, fixing some of the model parameters (the ψ α 's angles and z V ) or leaving all of them free. When mixing as many data samples carrying important -sometimes dominant -overall normalization uncertainties, the effect of the so-called Peelle's Pertinent Puzzle (PPP) [103] cannot be ignored to avoid getting biased estimates for physics quantities. The use of χ 2 minimization methods, first questioned, has been finally justified [104,105]; however, when dealing with distributions or form factors, the pertinent method turns out to invoke the underlying (theoretical) model function which is just what fits are supposed to provide. Following [106], [27] has defined an iterative procedure proved to cancel out biases (see the Appendix in [27]).
It is worthwhile briefly sketching the conclusions of 22 [27] concerning the issue just reminded : 1/ The normalization uncertainty (denoted σ, possibly s-dependent) of any spectrum T is absorbed within the total covariance matrix contributing to the general χ 2 , 2/ The correction which affects the normalization of T is derived from the data, the theoretical function provided by the fit, σ and the statistical covariance matrix -which may absorb the uncorrelated systematics.
Let us make a few general statements about fit results and data before focusing on specific topics, mainly the KK and τ sectors and the spacelike region behavior of the the pion and kaon form factors. The detailed properties will be emphasized in the next Sections. Figure 3: BHLS 2 fit to the π + π − π 0 spectra, the RS solution : Top panels show the fit and data in, resp. the ω and φ regions, the bottom one focuses on the intermediate region.
• As for the behavior of the data samples within our global fits, one can state that the issues met with the dipion spectra from BaBar [64,65] and KLOE08 [61] within the BHLS framework [25] are confirmed within the BS and RS variants of the BHLS 2 framework studied here; they are also discarded in the present study. Whatever the specific con- Figure 4: BHLS 2 fit to the (π 0 /η)γ spectra, the RS solution : Top panels show the case for the π 0 γ spectra in the ω and Φ regions, the bottom panels display the corresponding plots for the ηγ spectra.
ditions of the various fits performed, the π + π − data samples from NSK (i.e. CMD-2 and SND), KLOE10/12, BESIII [66] and CLEO-c [67] are as well described within the BHLS 2 framework as they were already within BHLS [25]; more precisely, NSK and KLOE10/12 contribute an average χ 2 per data point of ≃ 1, while the BESIII and CLEO-c samples yield χ 2 /N pts ≃ 0.6. Figure 2 displays the spectra and fits for the pion form factor in the e + e − annihilation and in the τ decay.
• In order to account for the difficulty to fully address the π + π − π 0 channel in the Φ mass region within BHLS, a so-called B model was defined which simply turns out to discard from fit this mass region for (only) this channel [25]. In contrast, BHLS 2 does not meet any issue 23 and achieves a fair description of the π + π − π 0 data up to, and including, the Φ mass region. The spectra are shown in Figure 3.
• Finally, the e + e − → (π 0 /η)γ spectra are also nicely fitted within BHLS 2 -as well as in 23 As for the discarded data sample [86], see the above Section 13 and also [91].
In order to perform global fits, one obviously should make an assumption on the energy calibration of each of the data samples considered. For the dipion spectra, one has gathered into the fit data coming from the VEPP-2M machine at Novosibirsk (CMD-2 and SND), data from DAPHNE (KLOE10 and KLOE12), from BESIII and from CESR. A significant mismatch of the relative energy calibrations of these would have produced a failure of their common fit because of the ρ − ω drop off region; this is clearly not seen neither in the χ 2 contributions of the various samples, nor in Figure 2 (see the inset in its leftmost panel). Moreover, comparing both panels in this Figure excludes a significant energy calibration mismatch between e + e − and τ data.
Moreover, the data collected by CMD-2 and SND around the ω and Φ peaks in the 3-pion channel as well as in the (π 0 /η)γ final states indicate a good consistency between the energy calibrations of the various data samples collected in different runs distributed over several years. This is visible in the fit properties as well as in the Figures where the fit function and the various data are displayed. 16 The BHLS 2 Global Fits : The e + e − → KK Channels The available data covering the K + K − channel 24 come from SND, CMD-2 and BaBaR. Recently CMD-3 has published a new K + K − measurement [76]. As for the neutral channel, most data come from CMD-2 and SND and have been listed in Section 13. A new data sample covering our fitting range has also been produced by CMD-3 [78]; unfortunately, the published K L K S spectrum from BaBar starts at 1.06 GeV, just above our fitting range [107].

Preliminary Remarks
As for the K + K − channel : In its [76], the CMD-3 Collaboration points toward disagreements with both the BaBar [75] and the original CMD-2 data samples [72]. Once the correction factor (1.094 ± 0.040) for the trigger efficiency of CMD-2 is applied, the discrepancy between the CMD-3 and CMD-2 K + K − samples looks much reduced [74] -possibly washed out. Anyway, a full examination of the various e + e − → KK data samples 25 is worth addressing.
On the other hand, preliminary BHLS 2 analyses based on the CMD-3 and BaBaR data, made us aware of a visible degradation of both KK channel data description above 1.025 ÷ 1.030 GeV; therefore, in the present work, one limits the fit region for the KK data to the √ s ∈ [2m K , 1.025 GeV] range. The explanation for this effect is unclear. It could be an intrinsic weakness of the BHLS/BHLS 2 models 26 which, as it stands, does address a possible 24 In the following, as already noted in Section 13, the two CMD-2 scans [72] are corrected as stated in [74]. As for the BHLS predictions in [25,26,27], they are updated correspondingly. 25 As for the contribution of the energy region from s ∈ [4m 2 K , 1.06 GeV 2 ] to a µ (K + K − ), BaBaR obtains 18.64 ± 0.16 stat ± 0.13 syst ± 0.03 V P while CMD-3 gets 19.33 ± 0.040 in units of 10 −10 , a difference of ≃ 0.7 × 10 −10 between the central values, a 3.25σ effect. 26 The standard HLS framework [18], which underlies the BHLS/BHLS 2 models, does not account for vector particles not belonging to the fundamental vector nonet. onset of the high mass vector meson contributions; it could also reflect an experimental issue. For the sake of consistency, we have chosen to impose the same fitting range also to the former CMD-2 and SND data samples in both the charged and neutral modes,

A Specific Issue of the KK data
As already pointed out, one does not detect any indication for a mismatch in the energy calibration of the various data samples covering the 3-pion and (π 0 /η)γ channels, all collected with the CMD-2 and SND detectors on the same VEPP-2M collider. One also does not observe any mismatch between the energy calibration of the 3-pion and (π 0 /η)γ data and those collected by the same collaborations in the KK channels. Nevertheless 27 , CMD-2 and SND carry a (correlated) absolute energy uncertainty in the range of 30 ÷ 40 keV, while for CMD-3, the absolute energy uncertainty is ≃ 60 keV, statistically independent of that of CMD-2 & SND. Moreover, BaBaR also reports a 50 keV correlated energy scale uncertainty [75]. When dealing with an object as narrow as the Φ meson within a global framework, these uncertainties should be accounted for otherwise the fits cannot succeed.
As the bulk of data samples for both KK channels has been provided by CMD-2 and SND, it looks quite natural to choose their common energy calibration as reference; this will be denoted E N SK in the following. For this purpose, one defines E CM D3 = E N SK + δE CM D3 and E BaBar = E N SK + δE BaBar to recalibrate the CMD-3 and BaBaR data point energies in a fully correlated way and derive their E N SK equivalent energy.
Therefore, in order to deal with the KK data, we are faced with a case where two kinds of correlated uncertainties have to be accounted for : Energy calibration and absolute normalization of the cross sections. Hence, one merges in the minimization procedure the method developped in [27] and the two calibration parameters δE CM D3 and δE BaBar to be also determined within the same fits.
For this purpose, one assumes that there is no functional relation between the normalization uncertainty λ i -of an experiment labelled i -and the energy scale uncertainty δE i ; then, as ∂λ i /∂[δE i ] = 0 for each experiment, substituting, as emphasized in [27], the value for λ i derived from ∂χ 2 /∂λ i = 0 into the expression for χ 2 is unchanged compared to [27].

Fitting the KK Channels : NSK, BaBaR & CMD-3 Separately
In the following, the non-kaon sector data and channnels used in our fits include, unless otherwise stated, the π + π − , π + π − π 0 , π 0 γ, ηγ annihilations, the τ decay to π ± π 0 and a few decay partial widths (see [25]). One starts by discussing the fits where the CMD-2 & SND (NSK), BaBaR and CMD-3 kaon data are not influencing each other within the minimization procedure. Table 1 gathers most of our results in this case; these are commented right now.
• With regard to NSK : < χ 2 >, the average (partial) χ 2 per NSK data points, is of order 1 in both the neutral and charged kaon channels; thanks to the fair description of the other channels involved in the fit, the global fit probability exceeds the 90% level.
Therefore the CMD-2 and SND data samples successfully fit the BHLS 2 framework and, moreover, the derived scale corrections meet the reported scale uncertainty expectations.
• With regard to BaBaR : One has chosen to supply the NSK neutral data to the fits in order to populate each dikaon channel as the BaBaR K L K S data [107] do not cover our fitting energy range, As the BaBaR charged spectrum [75] has been collected in the ISR mode, it looks appropriate to compare, within the fit, each spectrum datum in a bin to an average value of the model spectrum over the bin width. The second data column in Table 1 shows that the BaBaR sample yields < χ 2 >= 1.6, somewhat large, and a global fit probability also above the 90% level.
The minimization procedure returns as energy shift vs NSK δE BaBar = −122.3 ± 18.8 keV; this looks fair when compared to the independent 30 ÷ 40 keV (NSK) and 50 keV (BaBaR) energy calibration uncertainties. Moreover, the spectrum scale correction is at the 10 −4 level -so that the absolute scale of the BaBar spectrum has not to be corrected. On the other hand, the scale corrections for the neutral NSK data included in the fit are -2.3% (SND) and +0.5% (CMD-2), almost identical to what has been obtained in the fit of the NSK data in isolation mode.
Interestingly one can also derive the scale correction for the charged data of SND (+5.8%) and CMD-2 (-3.2 %) even if they are not submitted to the fit; these predictions are clearly in the expected ballpark. This indicates, beforehand, that NSK and BaBar data in the kaon sectors are consistent with each other.
• As for CMD-3 : This Collaboration has provided data samples in both dikaon final states, When fitting these data, the correspondence is made -as with NSK -between each datum and the theoretical model computed at the nominal energy. As for NSK and BaBaR, the systematic uncertainties are treated as fully correlated and serve to build up the systematic covariance matrices to be added up to the statistical ones; these are used in the fitting of the neutral and charged CMD-3 data within the global BHLS 2 framework. The fit results are displayed in Table 1.
The < χ 2 > for both CMD-3 data samples are clearly very large and the global fit probability poor. The scale corrections found are +2.6% (charged channel) and −2.6% (neutral channel) to be compared to the reported r.m.s., resp. 2.0% and 1.8%. The possibility of treating the bin information as in ISR experiments 29 has been examinedspecifically, only in the neutral channel -and in view of the outcome, we gave up. 28 One should stress that such a "renormalization" is the natural outcome to the scale uncertainty problem and that the main avantage of any kind of global fit is to allow determining such a correction in consistency with a (large) set of data samples coming from different sources. 29 Comparing the datum information to the fit function averaged over the bin width, as just sketched for BaBaR.
For the sake of completeness, we have nevertheless also performed the fit of the CMD-3 data in isolation by discarding the non-diagonal part of the total covariance matrix 30 , as seemingly done in [75] for the BaBar data; the fit results change dramatically as one reaches < χ 2 >= 14/16 and < χ 2 >= 17/17 for resp. the charged and neutral modes of the CMD-3 data and with a fit probability exceeding 90%. Moreover, in this case (CMD-3 data in isolation) the corrections to the original normalization of the crosssections become negligible (≃ 10 −4 ) for both the charged and neutral channels 31 . As for the energy shift vs NSK it moves from −60.6 ± 13.3 keV (see Table 1) to −41.6 ± 15.1 keV, still in the expected ballpark.  Table 1: BHLS 2 fit properties in the kaon sectors. The kaon data from NSK, BaBaR and CMD-3 are fitted in isolation. Running BaBaR data "alone" is complemented by the NSK data for the K L K S channel. The data and channels other than kaon are identical in all cases; the full χ 2 value and the total number of data points (N points ) refer to the global fit. The last line displays the shift values δE CM D3 or δE Babar relative to the NSK data energy calibration.

Fitting the KK Channels : NSK, BaBaR & CMD-3 Combined
This Subsection is devoted to analyzing the consistency of the three available groups of data. Quite generally, the covariance matrices used within the fit procedure assume full correlation for the systematics, as in the previous paragraph. Nevertheless, one also briefly reports on additional information about a fit performed assuming the systematics uncorrelated for CMD-3 data. It is worthwhile proceeding with pairwise combinations. 30 This turns out to consider that the systematic uncertainties are fully uncorrelated; some configurations intermediate between full correlation and no correlation can lead to similar results. 31 The < χ 2 >'s of the fit prediction to NSK data remain in reasonable correspondence with the NSK fit in isolation (49/49 and 117/92) while the distance of this prediction to BaBaR data degrades more severely (64/27).
• As clear from the first data column in Table 2 -and as could be inferred from the results already reported -one observes a fair consistency between NSK and BaBaR samples within the global BHLS 2 framework. Comparing the present NSK and BaBar contributions to the (total) χ 2 with their analogs in their fits in isolation clearly shows that there is no significant tension between them. The only noticeable difference is the central value for δE BaBar found closer to 100 keV by 1σ.
Other pieces of information concerning the absolute spectrum normalizations 32 are worth to be mentioned. The fit outcome indicates that the BaBar data normalization should be downscaled by ≃ 0.4%, the (updated) CMD-2 charged spectrum by 3% while the SND charged data normalization should be increased by 6.1%. These corrections are always consistent with the respective expectations reminded above.  Table 2: BHLS 2 fit properties in the kaon sectors. The kaon data from NSK, BaBaR and CMD-3 are fitted pairwise combined. The χ 2 values are displayed for individual channels or for the total number of data points N points . In the two rightmost data columns, the top entry is δE CM D3 and the bottom one δE Babar . The rightmost column displays the fit results when treating the CMD-3 systematics as uncorrelated.

KK
• The second data column in Table 2 displays results from fits combining NSK and CMD-3. Here also one observes a striking resemblance between the individual partial χ 2 's in the combined fit and in their fits in isolation, i.e. the NSK data samples go on yielding < χ 2 >'s of order 1 while CMD-3 are much larger; the global fit probability is also poor.
Besides the energy shift δE CM D3 = −47.9 ± 6.7 keV, the spectrum rescalings derived from the fit are interesting pieces of information. In this way, the fit returns an increase of the charged CMD-3 data normalization by 0.6% while the neutral CMD-3 ones should be downscaled by 2.8% (a 1.5σ effect). One also finds that the SND charged data normalization should undergo an increase by 10.4%, while for CMD-2 the increase is only 1.5%. As for the neutral kaon normalizations SND should undergo an increase by only 0.3%, while CMD-2 is downscaled by 2.1%. So, one does not observe major changes in the NSK normalizations.
• The third data column in Table 2 displays results from fits combining BaBar and CMD-3 kaon data, excluding their NSK analogs. The global fit probability clearly indicates that such a combination is not really favored; one should remark, nevertheless, that all properties exhibited by CMD-2 and SND are in fair correspondence with their analogs yielded in their fits in isolation; additionally, even the χ 2 distances of the CMD-3+BaBaR fit function to the NSK data (not under fit) -given within parentheses -are analog to their values derived in NSK alone or NSK+BaBaR fits. Quite interestingly, the tension observed between CMD-3 one the one hand and BaBar and NSK on the other hand is not manifest in plots, as soon as the energy calibration of CMD-3 vs NSK is applied. For instance Figure (5), shows the KK cross sections derived from fit to NSK + BaBaR; the corresponding data are rescaled and appropriately shifted as emphasized in the previous item; one has superimposed the (unfitted) data from CMD-3 shifted by −61.8 keV and rescaled as predicted by the NSK + BaBaR fit : −3.9% and 2.9% for resp. the charged and neutral CMD-3 data samples.
On the other hand, one has found worth to also perform the fit sketched in the rightmost data column in Table 2, i.e. treating the NSK and BaBar uncertainties as fully correlated and both CMD-3 data samples as carrying uncorrelated uncertainties. Using the corrections to the absolute normalization of each sample as derived from this fit 33  This clearly indicates that the observed tension is essentially due to the non diagonalparts of the CMD-3 error covariance matrices constructed from experimental information, as also done for the other data samples, CMD-2 in particular.

Fit Properties of the Various Broken HLS Models
All fits referred to in this Section involve the collection of data samples and channels already described and listed in Section 13; as for the KK channels, one limits oneself to using only the samples from CMD-2 (updated), SND and BaBaR which do not show any obvious tension. The motivation for this choice is to avoid disentangling model dependence effects, on the one hand, and identified tensions among data samples 34 on the other hand. There is indeed no gain in increasing the statistics if this turns out to increase the systematics in an uncontrollable way.
The BHLS 2 fits presented up to now have been derived using the so-called Reference Solution (RS); this solution combines the so-called Basic Solution (BS) defined in Section 7 and the Primordial Mixing defined in Section 8.
Nevertheless, in order to explore model dependence effects, analyzing also the fit results derived using BS is certainly relevant. On the other hand, for the same concern, we have rerun the BHLS model [27] using the SND together with the updated CMD-2 data; it then looked appropriate to also include the BaBaR dikaon sample. Therefore, with the BS and RS variants of BHLS 2 and with BHLS, one actually has at disposal three approaches to the same data; as for BHLS, one concentrates on the so-called model B which circumvents the issue met at the Φ mass in the three pion annihilation channel which is one of the motivations having led to BHLS 2 . The other one was to ascertain the description of the dipion threshold region and to improve the account of the spacelike region. The main fit properties are sketched in Table 3 and deserve some comments.
The leftmost data column -dealing with BHLS -shows the various χ 2 for the displayed channels and the sum of the contributions of all data submitted to the fit. In all channels the average χ 2 is of the order 1 or better 35 and the (global) fit probability is fairly good.
The rightmost data column displays the corresponding information derived with the RS variant of BHLS 2 . Here also, the fit quality is good and one yields a fair probability. One should note the quite favorable < χ 2 >'s for the pion and kaon form factors in the spacelike energy region in conjonction with good χ 2 's in their timelike regions; the spacelike data samplescoming from NA7 and Fermilab -are discussed in detail in the next Section. As far as we know, 33 This amounts to the following rescalings w.r.t. the original normalizations : -1.2% (CMD3 K + K − ), -0.7% (CMD3 K L K S ), +8.5% and -0.3% (K + K − from resp. SND & CMD2), +0.3% and 1.9% (K L K S from resp. SND & CMD2) and ≃ −0.3% for BaBaR K + K − . 34 The issue between BaBaR and CMD-3 charged kaon data has already been noted in [76]. 35 One should note the smaller number of data points considered in the fit procedure for the π + π − π 0 annihilation channel, because of excluding the data points from the Φ region.  the fits presented here are the first ones involving the pion and kaon form factors simultaneously in their timelike and spacelike regions with such a quality. The third data column displays the fit properties derived from the BS variant of BHLS 2 using the same data as those just discussed for the RS variant. The fit probability is not disastrous but unusually low in the HLS context. Scrutinizing the various items in Table 3 one clearly observes a significant tension in the τ decay and for the 3-pion annihilation data where the χ 2 is larger than obtained by the RS variant by resp. 20 and 13 units; the pion form factor in the spacelike region yields also a 8 units excess.
In view of this observation, one has redone the fit by discarding the τ data from the BS fit; the results are reported in the second data column of Table 3 and are in clear correspondence with the results obtained by fitting with the RS variant of the BHLS 2 model. Nevertheless, some tension survives in the 3-pion sector while the spacelike data recover an optimum fit quality. This may indicate that some improvement could come from an improved treatment of the breaking in the L A sector of the HLS Model 36 .
In conclusion, despite their different structures, these models/variants open the possibility to evaluate motivatedly model dependence effects in the evaluation of some reconstructed physical quantities; it could then be interesting to compare these to the corresponding effects generated when dealing with contradictory samples.   Table 4: BS variant global fits excluding τ data: The leftmost data column shows the partial χ 2 's in the fit including both kinds of spacelike data, the rightmost one shows the corresponding results in the fit excluding all spacelike data. Gobal fit information is displayed in the last lines. The numbers of fitted data points in each channel are given within parentheses. 36 Issue under investigation. Figure 6 displays the pion and kaon form factors in the spacelike region coming out from fit. The spacelike data supplied to the global fit procedure are only the (model-independent) data samples collected by the NA7 [46,47] and Fermilab experiments [92,93]. Altogether, these experiments cover photon virtualities down to ≃ −0.25 GeV 2 and ≃ −0.12 GeV 2 for resp. the pion and kaon form factors.

The Pion and Kaon Form Factors in the Close Spacelike Region
The pion and kaon form factors derived from the global fit including the spacelike data are the red curves in Figure 6. Global fits have also been performed using the BS and RS variants of BHLS 2 , excluding the spacelike data. The relevant fit results obtained using, for instance, the BS variant are displayed in Table 4. They correspond to the green curve in Figure 6; as for the kaon form factor the predicted (green) curve superimposes exactly to the fitted (red) one. Table 4 and Figure 6 obviously prove that the spacelike data are not a real constraint in the BHLS 2 global fit approach; stated otherwise, BHLS 2 supplied with only timelike data fairly predicts the spacelike behaviors far inside the spacelike region for both the pion and kaon form factors. It is worthwhile to note that these analytic continuations are performed across energy squared distances which can be as large as 1 GeV 2 . For instance, the kaon form factor displayed here is the analytic continuation of the charged kaon form factor fitted in the Φ mass region, i.e. the curve shown on the rightmost panel of Figure 5. Such a property is rare enough to be underlined. Figure 6 gives the opportunity to comment on the representation of spectra affected by a normalization (scale) uncertainty; the NA7 spacelike data are a quite simple illustration of handling the rescaling effects. When a scale uncertainty affects a spectrum, there are obviously two possible ways to display the result : Either rescaling the data, keeping the fit/prediction untouched or performing the converse.
So, in Figure 6 one has chosen to display the NA7 data rescaled, keeping unmodified the fit (and prediction) function(s), the plotted uncertainties being the statistical ones. As for the Fermilab data samples, taking into account their low accuracy, no rescaling has been performed and the plotted uncertainties are just the reported total errors. One can, for once, go into a few plotting details in simple wording. Following the method given in [27], the scale correction is given by : where (NA7) σ = 0.9% or 1.0%, f i = |F P (s i )| 2 (P = π ± , K ± ), V is the statistical error covariance matrix and m i is the measured datum at s = s i . The plotted NA7 data m ′ are related with the original ones m by 37 : Looking at Figure 6, the precise data (NA7 pion form factor) are indeed observed to tightly follow the fit function, once the normalization correction is applied. 37 Numerically, one finds λ π = −1.2% and λ K = −1.0%.

The Pion Form Factor : The Model-Dependent Data
As one cannot measure the pion form factor in a model-independent way for photon virtualities Q 2 = −s > 0.3 GeV 2 , this physics has been performed using the scattering process pe → eπ + n, referred to in the literature as 1 H(e, e ′ π + )n. In this approach, two questions should be addressed : i/ Separating out the longitudinal cross section σ L , ii/ Extracting F π from σ L .

JLAB 2008 DESY (updt 2008)
BeBek 1978 BeBek 1976 BeBek 1974 CEA 1973 Figure 7: The pion form factor in the spacelike region; the curve displays the function F π (s) coming from our global fit of the timelike data and of the spacelike data from NA7 [46] and Fermilab [92]. The (unfitted) data flagged by JLAB and DESY are extracted from [98], the other can be found in [94]. See also Figure 6.
One can refer to [96,97,98] for comprehensive accounts of the issues. As for item i/, the method has evolved and improved in the course of time from the early experiments (also reported in [94]) to the more recent ones [99,100,95]. Item ii/ addresses more deeply the nuclear physics content of the process in order to account at best for the π + pn vertex and for the off-shell character of the intermediate π + . Presently, the favored approach is the Regge Model developed by Vanderhaeghen, Guidal and Laget [108,109]. The data extracted for F π (s < 0) have thus improved in time but remain unavoidably model-dependent. These data samples have not been included in our fitting approach, but simply compared with the extrapolation of our fit function at s < 0, a pure prediction below s ≃ −0.25 GeV 2 .
In Figure 6, besides the fitted data from NA7 and Fermilab, one has displayed the (unfitted) data from the DESY experiments [99,100] and the lowest photon virtuality data collected at Jefferson Lab by the F π Collaboration. The data points flagged by JLAB 2008 and DESY have been extracted from the latest -to our knowledge -publication of the F π Collaboration [98]. One can observe that the updated [98] DESY data points fall exactly on the F π (s) function derived from a fit using also the NA7 [46] and Fermilab [92] data. The three JLab plotted data points are found systematically below the fit function; however, most of these JLab data points can also be found in a previous JLab publication [95]. The difference between the values in [98] and [95] seems mostly coming from using different extraction methods of σ L (Q 2 ) [96]. The data points from [95] seem in closer agreement with the continuation of our F π (s) than the estimates given in [98]. They are also in better agreement with what can be expected from the DESY data points as updated in [98]; this updating seems to deal with the extraction of F π from the longitudinal cross section σ L . The behavior of (the extrapolated) F π (s) beyond −1 GeV 2 is shown in Figure 7 together the other reported data.
The content of Figure 6 gives a hint that the extrapolation of the BHLS 2 F π (s) could well be valid down to s ≃ −1 GeV 2 . This also reflects a ρ dominance well inside the non-perturbative spacelike region. This dominance is illustrated by fit properties : When fitting the pion spacelike data, beside the ρ contribution, the fit function contains also those from the ω and Φ mesons; the fit returns < χ 2 > /N = 55.4/59; however, computing the χ 2 distance of the fit function amputated from the ω and Φ contributions to the pion data from NA7 and Fermilab, one obtains < χ 2 > /N = 53.0/59. This indicates that the (tiny) contributions of the ω and Φ mesons introduce a (tiny) noise and that nothing else than the ρ contribution is actually required.
This ρ dominance phenomenon does not happen for the kaon form factor as, truncating the prediction for F K (s) from its ω and Φ contributions returns a disastreous χ 2 distance; for instance, the full χ 2 distance of the NA7 kaon data to F K (s) which is ≃ 14, becomes ≃ 2300 when F K (s) is amputated from its ω and Φ contributions. Such a property is quite unexpected as the parametrisation for the ω and Φ [25], reminded in Section 9, is nothing but a Breit-Wigner lineshape modified in such a way that the inverse propagators fulfill D V (0) = −m 2 V , these HK masses being those in Equations (26). That the analytic continuation, over more than ≃ 1 GeV 2 , of such a simple parametrization, with parameter values fixed at the Φ mass, can provide a good description of the existing kaon spacelike data may look quite amazing.

The Pion Form Factor : The Lattice QCD Data
Some Collaborations have published information on the pion form factor derived from Lattice QCD (LQCD) data at pion masses close to the physical pion mass value [3]. We compare these LQCD form factors to our fit function F π (s < 0) which, already, fairly well describes the accurate NA7 experimental data as obvious from Figure 6 and Tables 3 or 4.
The HPQCD Collaboration has thus produced a parametrization 38 of the pion form factor [110] F π (Q 2 ) for virtualites in the range Q 2 ∈ [−0.1, 0.0] GeV 2 . In order to perform a numerical comparison with our fit F π (s < 0), we have plotted their pole parametrization, setting the lattice spacing a at zero and using their pion radius squared < r 2 π >= 0.403(18)(6) fm 2 . In Figure 8, our fit function for F π (s < 0) is shown by the red curve; the blue curve is the HPQCD parametrization at the central value for < r 2 π >. The red and blue curves are clearly very close from each other, if not overlapping. The green curves display the HPQCD parametrization with < r 2 π > at ±1σ from its central value. This indicates that the LQCD parametrization 38 We thanks L. Lellouch for having drawn our attention on this paper.  Figure 8: The red curve is |F π (s)| 2 derived from our global fit, the blue curve is the HPQCD parametrization and the green curves display this parametrization at ±1σ of the central value for their < r 2 π >. The LQCD spectra from ETMC and the experimental data shown are commented in Subsection 18.3; the data point from DESY is not fitted. [110] is in fair agreement with the BHLS 2 fit function F π (s < 0) much beyond the expected HPQCD range of validity and validated by the DESY data points [99,100]; Figure 6 clearly shows that the agreement between HPQCD [110] and BHLS 2 extends to the prediction, shown in this Figure by the green curve, derived when fitting with discarding the spacelike data. Six spectra for F π (Q 2 ) have been published by the ETM Collaboration [111]; two of them correspond to a pion mass close to its physical value, namely the gauge ensembles cA2.09.48 and cA2.09.64. These spectra are given by the form factor values F π (Q 2 ) and uncertainties as functions of Q 2 /M 2 π . In order to compare with phenomenology and experimental data, one has to restore the values for Q 2 from the provided values 39 for Q 2 /M 2 π . Figure 8 shows the |F π (Q 2 )| 2 spectra derived from the cA2.09.48 and cA2.09.64 ensembles in [111]. They are displayed with their reported statistical errors. The χ 2 distances of cA2.09.48 and cA2.09.64 to the BHLS 2 form factor squared |F π (s)| 2 have been computed, assuming the systematics negligible compared to the statistical errors [111]; they are displayed inside Figure  8. One can remark that the average χ 2 distance per datum is 0.91 for cA2.09.48, quite a good value for data outside the fit procedure; it is only 0.067 for the cA2.09.64 form factor. This can be considered as an independent confirmation concerning the smallness of the systematics versus the statistical errors for the considered ETM spectra.
Stated otherwise, the data points of the cA2.09.64 form factor are almost exactly on the fit function. Out of curiosity, we have redone the fit in the RS variant of BHLS 2 discarding the spacelike experimental data in order to get, online, the predicted form factor for any s < 0 using only timelike data; the average χ 2 distances become 0.93 and 0.072 for resp. cA2.09.48 and cA2.09.64. So, the χ 2 distance of the ETM samples to the BHLS 2 prediction indicates that the form factor lineshape is quite consistent with LQCD expectations. In other words, the ETM spectra provide a strong support to the BHLS 2 fits; this BHLS 2 model spectrum is also intimately related to the dealing with the reported correlated systematics as promoted in [27].
For illustration, one has also displayed in Figure 8 the closest among the two data points from DESY and a part of the NA7 spectrum rescaled and not rescaled. As expected from their χ 2 distances to the fit function, the ETM spectra are observed closer to the rescaled NA7 spectrum than to the original one.

The Neutral Kaon Form Factor at s < 0
As seen in the previous Subsections, BHLS 2 predicts the charged pion and kaon form factors in the close spacelike region with fair accuracy. Likewise, BHLS 2 also predicts the neutral kaon form factor F K0 (s); it is displayed in Figure 9. This is the analytic continuation of the form factor given in Subsection 11.4 and numerically determined within the fit procedure by the e + e − → K L K S annihilation data. One clearly observes a behavior quite different of those of F K ± (s), with a negative slope down to almost s = 0, when coming from the spacelike region.
Actually, as already stated in Section 9 one has F K0 (s = 0) = 0 up to contributions of second order in breaking parameters. What is displayed in Figure 9 is F K0 (s) computed at the fitted central value of the fit parameters 40 . To take into account their uncertainties, a gaussian Monte-Carlo has been run which returns the averaged central value F K0 (0) ≃ 10 −6 or less, depending on the fit conditions. This indicates that the negative excursion close to zero observed in Figure  9 could well be a numerical artifact.
Among the relevant physics properties of the pion and kaon form factors, their behavior close to s = 0 deserves special interest. As clear from Section 11, the charged pion and kaon vector form factors are analytic functions of s and, so, can be expanded in Taylor series around the origin. Neglecting terms of degree higher than 2 in s, one can write : where the a P 's are expected equal to 1 (up to second order terms in breakings) and the b P 's are related with the so-called charged radii by b P =< r 2 P > /6 for which several values have been reported.
In order to have a reasonable lever arm, one has chosen a s interval extending on both sides of s = 0 and bounded 41 As it is more customary to express the slope term in units of fm 2 , one displays our values together with external information in Table 5. Our results are given with two uncertainties; 41 Considering the real part of F P (s) permits to avoid (minor) issues with the π 0 γ threshold effects. one is derived from the sampling on the error covariance matrix and is indexed by "fit". One has also varied the content of the data sample set submitted to the global fit and compared the results to those derived with the maximal sample set. It is observed that the central values for < r 2 π ± > and < r 2 K ± > get shifted while the uncertainties are almost unchanged. One has also switched from the RS to the BS variant of BHLS 2 ; in this case, the central values get marginally shifted. Merged together, these shifts contribute an uncertainty denoted "mod" for "model" which, however, reflects as much the (weak) tensions between the various data samples of the largest sample set -defined in Section 17 -than the model properties by themselves.
As for < r 2 π ± >, Table 5 shows an agreement of our result with the experimental data from NA7 and Fermilab at the 1σ level or better. The result of G. Colangelo et al. [113] (CHS) is derived from a global fit of the timelike and spacelike pion form factor data based on a parametrization through the Omnès representation of F π (s) taking into account isospin breaking effects (ω) and prescriptions for the asymptotics of the ππ P -wave phase shift; this estimate is in perfect agreement with our own fit result. Using an elaborate method mixing the low energy ππ phase information and some form factor modulus information, B. Ananthanarayan et al. [114] recently derived the value flagged by ACD in Table 5 which is very close to ours. The two flavor result of J. Bijnens et al. [115] flagged as "ChPT 2 flavors" as well as the three flavor evaluation in [115] and the most recent Lattice QCD derivations of the HPQCD [110] and ETM [112] Collaborations are also in good correspondence with ours. Former LQCD evaluations for < r 2 π ± > can also be found in Table 22 of [117]. As for the curvature term, varying the fit conditions as just indicated for < r 2 π ± >, one can assess the effects observed when switching from the RS to the BS variant of BHLS 2 and when modifying the set of data samples submitted to fit. Altogether, this leads to our final result : The 2 flavor ChPT estimate based on all data samples available at that time (1998), [115] reports c π = 3.85 ± 0.60 GeV 2 and the 3 flavor evaluation [116] (c π = 4.49 ± 0.28 GeV 2 ) look in fair agreement with the BHLS 2 evaluation given in Equation (103).
As the pion and kaon form factors, in the close spacelike region, are frequently parametrized by a pole expression : we performed alike without, however, introducing additional scale factors as commonly done as indeed we believe in the relevance of the absolute scale of the F P (s)'s returned by the BHLS 2 fits. Such a parametrization turns out to impose the constraint c P = b 2 P to the curvature term. This leads to : < r 2 π ± >= 0.424 ± 0.001 and < r 2 K ± >= 0.274 ± 0.001, both exhibiting a shift of 0.005 fm 2 vs the second order Taylor expansion.
For what concerns < r 2 K ± >, the numbers reported in Table 5, show that our evaluation is at 1σ from the NA7 measurement [47] and coincides with those from [93]. As for the ChPT prediction from Bijnens and Talavera [116], it differs from the value returned by BHLS 2 at the 1σ level. There is no reported evaluation of the curvature term of the kaon from factor to which one could compare our own evaluation in Equations (102). Therefore, concerning the low energy behavior of the pion and kaon form factors, the overall picture is a fair agreement of BHLS 2 with the corresponding experimental measurements as well as with the exisiting LQCD and ChPT evaluations. Figure 10: The δ 11 phase shift; In both panels, one displays the Cern-Munich data [118], the phase shift from [119] (JS11) is shown by the green curves and the blue curves are the phase shift from [120] based on the Roy Equations. The red curves show the predictions from BHLS 2 (top panel) and from BHLS [25,27] (bottom panel).

The Isospin 1 P -wave ππ Phase Shift
Another relevant piece of information is the I = 1 dipion P -wave phase shift δ 11 . Both panels of Figure 10 display the experimental data collected by the Cern-Munich Collaboration [118], the phase shift derived using the Roy Equations [120] and the phase shift derived in [119] in the context of a common fit of the e + e − → π + π − cross section and of the dipion spectrum in the τ decay 42 . The phase shift coming from BHLS 2 is displayed in the top panel while the phase shift from the previous BHLS [25] is displayed in the bottom panel. In both cases, one observes tiny blips at the mass locations of the ω and Φ mesons.
In the bottom panel, the BHLS phase shift (red curve) exhibits a small structure close to threshold already encountered and discussed in our [27]; this was identified as an effect due to having a vanishing ρ 0 -ω HK mass difference in the BHLS framework [25]. One should also note that the BHLS phase and the Cern-Munich data are in fair agreement much above the validity range of the BHLS model (the Φ mass region); this can be traced back to having allowed the pion loop in the A − V transition amplitude F e ρ and in the ρ self-mass amplitude to carry different subtraction polynomials (see Subsection 11 above). This clearly accounts for the onset of the high mass vector mesons as the agreement remains valid up to 1.3 GeV and somewhat above.
In the case of BHLS 2 (red curve in top panel), the phase shift superimposes almost exactly on the data [118] and on the JS11 curve up to the Φ mass region. Both JS11 and BHLS 2 then start to depart from the data increasingly with increasing energies. Put together with the remarks already stated about the behavior of the dikaon cross sections slightly above the Φ mass, this gives a hint for an onset of higher mass vector mesons presently not accounted for in the HLS framework. The different behavior of BHLS and BHLS 2 should, anyway, be noted. Table 6 reports on the Lagrangian parameter values 43 derived when running BHLS with updated kaon data (first data column). Running the basic variant BS of BHLS 2 with τ data discarded gives parameters values displayed in the second data column, while the third data column reports on the numerical results from the BHLS 2 Reference Solution (RS) using the largest set of consistent data samples.

Lagrangian Parameter Values in BHLS and BHLS 2
Clearly, the Lagrangian parameters common to the various modelings cover a wide range of values; this does not prevent the fit properties (reminded in the last two data lines) to always exhibit a fairly good description of the data.

The (z V ÷ a) Correlation within BHLS 2
The value obtained for z V from fitting with the RS variant of BHLS 2 may raise questions -see Section 8 -as z V − 1 is numerically large 44 , even larger than z A − 1. For this purpose, it has been worthwile to explore the RS variant behavior in fits performed discarding the 3π final state data; this has revealed a numerical anti-correlation between the HLS parameter a and the BKY breaking parameter z V , undetectable in the parameter error correlation matrix as < δa δz V >= −0.22. Several solutions with comparable total χ 2 were found and the most striking difference among them just concerns the values for a and z V , as shown by :  43 The full list of parameter values and uncertainties, especially the subtraction polynomial coefficients, are not given; they can be provided upon request. 44 Even if only a 2.4ǫ effect, referring to the numerical estimate for generic ǫ's given in Footnote 7.
Channels  exhibiting a clear anti-correlation between the numerical values for a and z V . The origin of this anti-correlation can be traced back to the coupling of kaon pairs to the Φ; indeed, setting ξ 8 = ξ 0 in the Lagrangian Equation (116), one gets : and, using Equations (105), one can check that az V = (2.59 ÷ 2.61).
One has also run the BHLS 2 RS variant fixing z V = 1 and found the results shown in the last data column of Table 6. The 6 units difference between the minimum χ 2 obtained fixing z V = 1 and the best fit leaving z V free is almost entirely concentrated in χ 2 (e + e − → 3π) which becomes 146.3 (compared to 140.4 for 158 data points). Running also the BHLS 2 RS variant imposing via MINUIT, a lower bound a ≥ 2 returns : a = 2.001 ± 0.004 , z V = 1.305 ± 0.002 , χ 2 /N pts = 1131.9/1237 (107) where the 4 χ 2 units difference with the best fit is also concentrated in the e + e − → 3π channel. These fits, as all fits performed within the BHLS 2 context, return az V ≃ 2.6 and show fairly good quality. Moreover, all solutions return consistent information about the physics observables as will be exemplified below. Therefore, the relatively large value for z V − 1 obtained in the best fit with the (unconstrained) RS variant of BHLS 2 is not a real issue. One should also note that the multiplicity of solutions reflected by Equations (105) reduces to a single one when the 3-pion data are accounted for.

A Remark on the HLS Parameter a
Correlated with the remark on z V , one should note the value for the standard HLS model parameter a coming out of the best RS fit : 1.534 ± 0.001. This is much smaller than any of its previous determinations and the BS or BHLS best fit results (see Table 6) which also follow the trend observed previously.
It is commonly stated -see [18], for instance -that the standard VMD approach corresponds to the HLS model with a = 2 imposed. This corresponds to cancelling out any direct coupling of the photon with PS meson pairs in the Lagrangian. However, a look at the Lagrangian displayed in Appendix A.1 clearly indicates that there is no likely way to cancel out simultaneously the γπ + π − , γK + K − and γK L K S couplings when certain breaking mechanisms are at work. Therefore the properties usually attributed to the a = 2 constraint may vanish when the HLS Lagrangian is broken, as exemplified by the BKY mechanism (see Section 3).

Other Physics Results
As for other model parameters, focusing on BHLS and the RS variant of BHLS 2 : • One may wonder that the value for the vector coupling g increases by more than 20% when going from BHLS to BHLS 2 , keeping good fit qualities. However, the relevant piece of information, closer to observables, being the coupling g ρ 0 ππ = ag(1 + Σ V )(1 + ξ 3 )/2, it is interesting to compare BHLS to the RS variant of BHLS 2 : BHLS returns g ρ 0 ππ = 6.156 ± 0.005 while the full BHLS 2 fits gets g ρ 0 ππ = 6.133 ± 0.012 in reasonable correspondence with each other.
• The parameters associated with the anomalous sectors show that 45 (c 3 + c 4 )/2 yields comparable values with BHL and with BHLS 2 . However, for c 1 − c 2 , the difference is large; in this case, however, one should remind that the so-called "model B" variant of BHLS excludes from the fit the 3π data in the Φ mass region while the whole 3π spectra are considered within BHLS 2 .
• The Lagrangian term of order p 4 driven by z 3 is found small (z 3 ≃ −0.35%) but highly significant. This is its first experimental determination.
• The values for z A are also interesting as they are related with f K . Indeed, one knows that [25] : is still valid within the BHLS 2 framework. For this ratio, BHLS returns r = • The additional parameters introduced by the CD breaking and the Primordial Mixing are all found at the few percent level. The difference for ξ 3 between the (full) RS variant and the BS variant should be noted; it may have to be revisited later on, when a complementary mechanism may allow to absorb satisfactorily the τ data within the BS global fit procedure. This may also influence the 3π sector.

The muon LO-HVP Evaluation
BHLS 2 (or BHLS) fits provide the contributions of the various hadronic channels 46 it covers to a µ ( √ s ≤ 1.05 GeV). An overview of the main results, obtained in different running conditions is given in Table 7. These have been derived by submitting to fit most of the data samples presented in Section 13 which define a set of consistent samples covering the six annihilations considered; altogether, these 6 channels and the radiative meson partial width decay represent 1237 data points.
The recent [113], in a quite different theoretical context, provides the estimate G = a ππ µ ( √ s <  [113] and ours are different, but more importantly for this energy range, the dipion BaBar data included in the study [113] and discarded in ours, certainly contribute to the difference.

Estimates for a HV P −LO
In order to get the full contribution to a HV P −LO µ of the energy region up to 1.05 GeV, one should complement the contribution from the channels covered by BHLS/BHLS 2 by those of the channels presently outside this framework (4π, 2πη, . . . ). The most recent evaluation of this quantity, estimated by direct numerical integration of data is : In view of model effect studies, it has been found motivated to update and correct the data samples submitted to the BHLS fit in the way stated in Section 13. One has thus derived the updated BHLS evaluations given in the first data data column of Table 7; this supersedes the BHLS results formerly given in Table 4 of [27]. Differences for some channel contributions are observed and, in fine, a µ (HLS, √ s ≤ 1.05 GeV) increases by ≃ 2 × 10 −10 while its uncertainty slightly decreases (1.08 → 0.92, in units of 10 −10 ). An evaluation of possible additional systematics, specific of the BHLS modeling, has been performed and summarizes by [27] : Each piece shown here was found to act by shifting the central value for a µ rather than enlarging its uncertainty. The first term refers to the uncertainty in the treatment of Φ mass region in the 3-pion spectra, the second term takes into account the uncertainty on the ππ threshold behavior within BHLS (see Subsection 10.3 above and, especially, Figure 1); finally, the third piece reflects differences observed by running the BHLS fit procedure [25,26] with the τ data included and excluded.
As already stated, the BHLS 2 modeling has been motivated by the aim to cancel out the first two sources of uncertainty reported in Equation (108) by a better understanding of the s = 0 region and a full account of the 3-pion data up to the Φ region.
For this purpose, the three data columns in Table 7 referring to BHLS 2 show that : 47 As the (0.630 < √ s < 0.958) GeV interval is relatively far from both the threshold and the Φ mass regions, one could indeed expect similar results for BHLS and BHLS 2 .
• The evaluation for the various contributions to a µ ( √ s ≤ 1.05 GeV) derived using the BS and both RS variants of BHLS 2 carries differences generally at a few 10 −12 level, except for the 3π channel rather found at the 10 −11 level.
• Using the τ data has a significant effect. Indeed, Table 7 shows that running the RS variant of BHLS 2 including the τ data improves the uncertainty in the dipion channel while leaving the central value almost unchanged. This can be expected if the corresponding e + e − and τ data are not conflicting.
The data update having been performed, BHLS and BHLS 2 have been run on the same data and, thus, the differences between their respective evaluations reflect modelling effects. One can also compare these differences to the numbers listed in Equation (108). As shown by the various topics examined in Section 18, one can legitimately consider the π + π − threshold region accurately treated. Therefore, the difference 494.51 − 493.73 = 0.78 between the BHLS 2 and BHLS estimates for a µ (π + π − , √ s ≤ 1.05 GeV) is a motivated evaluation of the unaccounted for nonet symmetry breaking effects in the BHLS framework; it can thus motivatedly replace the former +1.4 estimate in Equation (108). Table 3 shows that all variants of BHLS 2 take well into account the 3π annihilation data over the whole energy domain up to the Φ mass region, the fit quality of the RS variant being noticeable. The difference between BHLS 2 and BHLS for this channel ≃ 0.5 × 10 −10 is in the range expected from the estimate [27] reported by the first term in Equation (108). Nevetheless, the 0.13 × 10 −10 difference for the 3π entry between the RS variant of BHLS 2 runnings including and excluding the τ data might be considered as additional systematics 48 .
It is interesting to compare the outcomes from fits to the values derived by the direct integration of the same data sample set treated in the BLS 2 framework. One may note that the sum of all HLS channel contributions corresponds to the fit expectations; however, this hides the fact that the dipion channel is fitted 3.3 × 10 −10 smaller and the 3-pion channel is fitted 1.74 × 10 −10 larger.
In summary : The issue is to examine whether additional sources of systematics are not at work. For this purpose, the consequence on a µ (HLS) of the various fits referred to in the previous Section carry a relevant piece of information. as also for the fits underlying the results given in Equation (107) or those displayed in the last data column in Table 6. These fits, covering a wide range of parameter values, are relevant for systematics estimates.
The runs leaving aside the 3-pion data briefly reported in Equations (105), show that the dominant ππ contribution can be as low as (494.04 ± 0.89) × 10 −10 , i.e. 0.47 units of 10 −10 smaller than the corresponding reference entry in Table 7. Nevertheless, the middle entry line in Equations (105) is very close to Equations (107), where the 3-pion data have been included and exhibits negligible differences with respect to our reference in Table 7. Finally, the fit corresponding to the last data column in Table 6 has also been analyzed and returns a value for a µ (π + π − , √ s ≤ 1.05 GeV) = (494.02 ± 0.93) × 10 −10 . This may indicate another additional systematic able to shift the central value by, at most, −0.49 × 10 −10 affecting our estimate for a HV P −LO µ ( √ s ≤ 1.05 GeV) given in Equation (109) Table 7: HLS contributions to 10 10 × a HV P −LO integrated up to 1.05 GeV, including FSR. The first data column displays the results using the former BHLS [25,27] and, the second one, those derived from the Basic Solution for BHLS 2 , the τ decay data being discarded. The next two data columns refer to the results obtained using the Reference Solution for BHLS 2 using the largest set of data samples, keeping or discarding the τ data. The last data column refers to the numerical integration for each channel of the same set of data which are used in the BHLS/BHLS 2 fits.
On the other hand, the question of using the CMD-3 data can be raised. We have run our code using the CMD-3 data [76,78] for both KK decay channels, discarding the kaon data from BaBar, CMD-2 and SND and using only the diagonal part of their error covariance matrices (see Section 16 above). In total, the change for a µ ( √ s ≤ 1.05 GeV) is noticeable where the additional systematics come from both model variations and data sample consistency issues. The model uncertainties include the marginal effect attributable to the τ data as this could reflect some (marginal) isospin breaking shortcoming. We have chosen conservatively as reference the RS variant fit excluding the τ data, One should also note that the tension between dikaon data samples introduces non-negligible systematics which contribute a bias. This may have to be revisited when new data will arise. : On top, the result derived by direct integration of the data combined with perturbative QCD; the next five points display some recent evaluations derived by LQCD methods and reported in resp. [124], [125], [126], [127] and [11] with N f = 2 + 1 + 1. The second point from [126] displayed has been derived by supplementing lattice data with some phenomenological information. These are followed by the evaluations from [70], [128] and [129]. The value derived using BHLS [27] -updated with the presently available data -and the evaluation from BHLS 2 are given with their systematics.
Taking into account the data upgrade, the corresponding quantity fo BHLS is : taking into account the present findings to go beyond the former systematics evaluation reminded in Equation (108).  (110) and (111) and the information given in the left part of

BHLS
In the case of BHLS, the systematics due to the τ channel have been added linearly to the others. Indeed, as stated several times, the observed systematics rather come as biases and so, should be treated separately from the fit error which fully absorbs the reported statistical and systematics from the various experiments encompassed inside each of our HLS frameworks. The additional systematics absorb uncertainties which can be attributed to the model workings and to the (marginal) tension observed between some of the selected experimental data samples.
Our evaluations for a HV P −LO µ , given by Equations (112), are shown in Figure 11 together with some other recent results derived in different phenomenological contexts or via LQCD data samples with N f = 2 + 1 + 1. As for their specific HLS modeling parts, the corresponding fit probabilities shown by the downmost entries in Table 7 are always fairly good.
As for modeling effects, the central values for BHLS and BHLS 2 -running on the same set of data samples -are shifted apart by ≃ 1.5 × 10 −10 ; however, if one takes into account the magnitude of the additional systematics, Figure 11 clearly indicates that this shift is not significant and so no noticeable modelling effects show up. Moreover, Figure 11 clearly shows that both HLS based estimates for a HV P −LO µ are consistent with the evaluation derived by numerical integration of the data.

BHLS 2 Evaluation of Muon Anomalous Magnetique Moment
Having at hand our estimate Equation (110), one can derive the BHLS 2 evaluation for a µ = (g − 2)/2, the muon anomalous magnetic moment. The left part of Table 8 collects the pieces to add up with a µ (BHLS 2 ) to get the full leading order HVP, its value is shown as the top entry in the right-hand part of Table 8 where the other contributions are included. The value for ∆a µ = a µ (exp.) − a µ (th.) and its statistical significance are the bottom entries there. Taking into account the additional systematics, one thus gets : and, taking into the systematics : 31.78 ± 4.17 ≤ 10 10 × ∆a µ (= 11 659 176.42 ± 4.17) ≤ 33.26 ± 4.17 (114) which means that the statistical significance of ∆a µ is predicted greater than 4.2σ. The shift just reported is partly due to an estimate of the modeling effects and partly due to somewhat conflicting aspects of some data samples.  Our estimate for a µ given in Equation (113) relies on a variant of the traditional estimate [131] for the Light-by-Light (LBL) contribution ((10.34 ± 2.88) × 10 −10 ). Danilkin, Redmer and Vanderhaegen (DRV) have quite recently published a comprehensive analysis of all ingredients participating to the LBL term and proposed for it the much more precise value [136] (8.7 ± 1.3) × 10 −10 which leads to :

Summary and Conclusions
The present paper has addressed several topics covering modeling, phenomenology and the evaluation of the muon anomalous magnetic moment which puts forward, as well known, a hint for a physics effect beyond the Standard Model expectations. For the sake of clarity, each of these topics is considered separately.

From BHLS to BHLS 2
The BHLS framework [25], although quite successfull [26,27,137], has been shown to exhibit a non-optimum behavior in limited and well identified kinematical regions : The threshold region and some aspects of Φ mass region; the former issue renders delicate the continuation of form factors across the chiral point, the latter has led to discard the Φ mass region for 3-pion final state -the so-called Model B [25].
In order to cure theses diseases, one has initiated a new breaking procedure of the general HLS framework [18]. For this purpose, beside the BKY mechanism [28,34,35], one has introduced a breaking mechanism taking place at the level of the covariant derivative (CD) itself, i.e. inside the basic ingredient of the HLS Model. This (seemingly novel) CD breaking mechanism allows naturally the ρ 0 , ω and ρ ± fields to carry different Higgs-Kibble masses. This takes its importance when going to loop corrections to the vector meson mass matrix [29,25] which generate a dynamical mixing of the vector mesons calling for a s-dependent redefinition of the vector meson fields to recover mass eigenstates.
Thanks to its nonet symmetry breaking piece, the CD breaking scheme generates a nonzero m ρ 0 −m ω difference which naturally makes the three s-dependent mixing angles to vanish at the chiral point. In this way, the timelike to spacelike analytic continuations become smooth without any ad hoc trick.
Relying on BKY and the CD breaking, the Basic Solution (BS) is defined and represents a first variant of BHLS 2 . Another mechanism can also be invoked which turns out to state that the neutral vector fields involved in physical processes are not their ideal combinations but mixtures of these, constructed via a rotation. This Primordial Mixing (PM) mechanism 49 of vector mesons assumes a 1 + O(ǫ) rotation. Working at first order in breaking parameters, the merging the BKY and CD breakings together with the PM mechanism define our Reference Solution (RS), a second variant of BHLS 2 covering the full HLS scope more easily than the BS variant (as it stands presently).
In both variants of BHLS 2 , the L V piece of the non-anomalous Lagrangian is modified while its L A piece remains unchanged 50 compared to the former BHLS [25].
The e + e − Annihilation Phenomenology The BKY and CD breaking schemes, together with the Primordial Mixing, allow to construct the non-anomalous and anomalous pieces of BHLS 2 , a new (broken) HLS Lagrangian. This allows for a simultaneous study of the e + e − annihilation into 6 final states (π + π − , π 0 γ, ηγ, π + π − π 0 , K + K − , K L K S ) within a unified framework where global fits can be worked out.
The fits have been performed by taking into account the fact that most of the data samples carry normalization uncertainties, frequently dominant; these are appropriately considered in order to avoid biased results for physics quantities. The iterative method developped in [27] is thus intensively used in the global fits reported here; the convergence criterium chosen to stop the iteration at order n + 1 was δχ 2 = |χ 2 tot (n + 1) − χ 2 tot (n)| < 0.3. As χ 2 tot ≃ 1000 ÷ 1200 for ≃ 1150 ÷ 1240 data points, this convergence criterium is clearly constraining.
The fitting procedure has allowed to detect (and discard) a very few data samples which hardly accomodate the global framework fed with almost all the existing data samples. Beside the already identified [27] two π + π − samples, one also found that, out of the 7 available π + π − π 0 data samples, one [86] hardly fits the global context, as also reported by others [91].  Table 3 prove the fairly good quality of the global description for the π + π − , (π 0 /η)γ and π + π − π 0 channels. As for the 3-pion channel, one should note that BHLS 2 encompasses easily the samples collected in the Φ mass region and so, solves the difficulty encountered there by the former BHLS.
A thorough study of the available data samples has been performed concerning both the K + K − and K L K S channels for which some tension is expected, as reported in [76]. The analysis in Section 16 has shown that the SND, BaBar [75] and the updated [74] CMD-2 data are in fairly good agreement with each other and accomodate easily the HLS context where all other kinds of data are encompassed. Tension has indeed been observed between the CMD-3 data [76,78] and the other samples. This tension has been traced back to the offdiagonal part of the CMD-3 error covariance matrices for both the K + K − and K L K S spectra; this analysis is summarized by Table 2. In order to include the BaBar (and/or CMD-3 data) within the global framework, energy shits relative to the CMD-2/SND energy scale (chosen as energy scale reference) are mandatory and have been fitted at values consistent with the reported expectations [75,76]. To perform this exercise, when using the CMD-3 data, the error covariance matrices of these (only) have been amputated from their off-diagonal part.
One can finally remark that, once the energy scales of the CMD-3 and BaBar [75] data relative to CMD-2/SND are applied, one does not observe a visible tension among them any longer as illustrated by Figure 5 where all data samples go nicely superimposing onto the same fit function.

Meson Form Factors in the Spacelike Region
The BHLS 2 form factors supplied with the parameter values derived from global fits to the timelike region data only have been shown to nicely predict the behavior of the π ± and K ± form factors in the close spacelike region, as they can be expected from NA7 and Fermilab data. Moreover, including these (model-independent) spacelike data within the BHLS 2 fit procedure does not exhibit any gain in the description of the spacelike region as obvious from Table 4 and from Figure 6.
Our derived π ± and K ± form factors, which perform fairly well simultaneously in the spacelike and timelike regions, strongly supports the validity of the BHLS 2 predictions in the s-region located in between. If the fairly good continuation of F π ± (s) is actually performed over a limited gap in s to reach the spacelike region, one may wonder about the so successfull continuation of F K ± (s) over a s-gap greater than 1 GeV 2 . This indicates that our modified Breit-Wigner formulae are performing (unexpectedly) well far beyond the Φ energy region where their parameters are determined.
The LQCD Collaborations HPQCD and ETM have provided, the former a parametrization, the latter spectra for F π ± (s) covering the same spacelike s-region than the NA7 experiment. Figure 8 shows the superposition of the BHLS 2 prediction (without any free parameter) together with the parametrization provided by HPQCD and the two ETM spectra; the accord is very good for HPQCD and simply perfect for both the ETM spectra. There is, unfortunately, no reported LQCD data for the kaon form factors to which one could compare.
The success just emphasized for the BHLS 2 π ± and K ± form factors compared to spacelike experimental or LQCD data gives confidence in other physics information involving the chiral point region. Let us limit oneself, in this Subsection, to briefly comment on the neutral kaon electromagnetic form factor. As reminded in the main text, the slope for the F K 0 (s < 0) form factor at s → 0 is reported negative. This is qualitatively consistent with the BHLS 2 finding for F K 0 (s) shown in Figure 9, as it is found to carry a negative slope at negative s and a positive slope at positive s. However, the effects of the neglected O(ǫ 2 ) corrections actually prevent to provide confidently accurate slope predictions in the neighborhood of s = 0 which is a delicate region for our O(ǫ) approximations.
The pion form factor measured by the NA7 and Fermilab Collaborations are model independent and, actually, it cannot be measured in a model-independent way above |s| ≃ 0.25 GeV 2 . To go beyond, one should accept to consider data where the extraction of F π ± (s) is model-dependent. Figure 6 already indicates that the two unfitted DESY measurements fall (almost) exactly onto the fitting curve together with the model-independent NA7 and Fermilab measurements. The JLAB data points are more or less closer to the BHLS 2 expectation depending on the method used to extract F π ± (s). Nevetheless, our analysis indicates that the BHLS 2 prediction could well be valid down to ≃ −1 GeV 2 , possibly slightly beyond.
In summary, our analysis clearly indicates that, inherently, BHLS 2 only supplied with consistent e + e − annihilation data, leads precisely to the spacelike expectations supported by the precise NA7 data or the available LQCD input as well. Nevertheless, LQCD information on F K ± (s) (and F K 0 (s)) would be welcome to better check the open strangeness sector.

Physical Quantities and Model Dependence Effects
BHLS 2 has derived evaluations of several physics quantities (see Sections 18,19 and 20) and some deserve special emphasis. As shown in Table 5, one reaches precise values for the charged pion and kaon charge radii (fm 2 ) : < r 2 π ± >= 0.429 ± 0.001 mod. ± 0.001 f it , < r 2 K ± >= 0.279 ± 0.002 mod. ± 0.001 f it ; our estimate for < r 2 π ± > is in fairly good agreement with the most recent evaluations [113,114]. Our evaluation for < r 2 K ± > is in accord with the reported values shown in Table 5 and is more precise. The Reference Solution of BHLS 2 also leads to : where the quoted systematic error reflects different variants of the modelling. Table 6 collects the values of the BHLS 2 model parameters obtained under several fit conditions and discussed in Section 20. An interesting numerical correlation is observed between the fundamental HLS parameter a and the BKY breaking parameter z V within the BHLS 2 framework. This can be approximated by az V ≃ 2.6; good fit quality are thus observed with a varying over a large interval, namely from 1.5 to 2.6. The rightmost two data columns in Table 6, which collect the fit values of a large sample of breaking parameters in two extreme cases (a ≃ 1.5 and a ≃ 2.6), also show correlatedly important changes in some of the other breaking parameter values. Table 6 also shows that the fit qualities reached in all cases are almost identical.
One can conclude herefrom that the various model parameter sets efficiently conspire to provide almost identical cross-sections and form factors which are, actually, the real observables feeding BHLS 2 . This has been substantiated by comparing the running of BHLS 2 fed with our consistent set of data samples under the various conditions underlying the fits reported in Table 6.
A motivated piece of information to detect parametrization effects are the integrals a µ ( √ s ≤ 1.05 GeV) given in Table 7. Running the BS variant (excluding the τ data) provides a µ (BS, √ s ≤ which is, at least, at 4.2σ from the BNL measurement. The systematics are actually upper bounds to a possible bias affecting a µ ; they come from both model effects and data sample tension.

Appendices A The Full HLS Non-Anomalous Lagrangian In The New Scheme
For clarity, the full non-anomalous HLS Lagrangian is written : and the various pieces will be given just below.
A.1 The L A + aL V Lagrangian piece First, one displays the part of L V M D = L A + L V relevant for the physics we address, especially e + e − annihilations. This is : The pseudoscalar fields occuring there are the renormalized ones defined exactly as in [25]. For simplicity, we use ideal vector fields. The V − A couplings are : One has discarded the (irrelevant) photon mass term. The vector meson masses are given in Eq. (26). ∆ V and Σ V have been defined in Section 3.
In practical use, one should perform the change of field given in Eq. (27) and collect the terms corresponding to each (neutral) vector meson coupling.

A.2 The L τ Lagrangian piece
The new L τ is given below (after performing an integration by parts to remove the terms with ∂W ± (= 0) : where all fields are renormalized; one should note the appearance of W ± (π 0 /η/η ′ )π ± couplings. One should also note, about the W ± (η/η ′ )π ± couplings that the contributions from L A and L V are not proportional in contrast with the W ± π 0 π ± case. We also note that the ρ ± (η/η ′ )π ± couplings listed above are needed in the transitions amplitudes for η/η ′ → π + π − γ in order to yield the right benhavior at the chiral point. ***** δ P is related with the PS mixing in the one angle mixing scheme [31] by Eq (YYYYY) below where also the ǫ i will be reminded. **** A.3 The L p 4 Lagrangian piece L z 3 As for the L z 3 Lagrangian, it writes in terms of ideal fields : It clearly allows for different V A and W V couplings, especially for the ρ meson.

A.4 The Renormalized Couplings of Vector Mesons to KK
Let us first remind the coupling of the fully renormalized vector mesons to π + π − ; having defined these as g V π ± = ag g π± /2, one has : Using the Lagrangian in Section A.1 and Eq. (38) one can derive the renormalized coupling of vector mesons to KK pairs. Defining: (121) the couplings involved in form factors and cross-sections are given by G V K + K − = ag/[4z A ] g V K± . One also has : and the full couplings are given by G V One cancels out here and in the following the dependence upon ∆ V . One should also note that all these couplings depend on s, the square energy flowing through the relevant vector line.

A.5 The V − γ Transition Loop Terms
As seen in the main text, the V − γ transitions amplitudes write : for each of V = ρ 0 R , ω R , Φ R . In this Section, one displays the expression for the loop terms Π V γ (s) in terms of the basic loop functions denoted Π e π (s), Π e K± (s) and Π e K 0 (s) and defined in Subsection 11.2. For the sake of conciseness, it is convenient to define : neglecting the O(ǫ 2 ) terms. These coupings do not depend on the "angles" α(s), β(s) and γ(s).
Using these definitions and those displayed in the Subsection just above, on finds : Expanding Equations (124) in order to rather parametrize in terms of Π e K sum and Π e K dif f (see Subsection 11.2) allows to exhibit the terms of order O(ǫ 2 ) in breakings which can be dropped out.

B Pseudoscalar Field Renormalization : A Brief Reminder
As already stated, one goes on using unchanged the breaking procedure in the L A sector of the HLS model as defined in [25]. The broken L A Lagrangian displayed in Eq.(8) implies a first step field redefinition in order to make the pseudoscalar kinetic energy term canonical. Additionally, as BHLS and BHLS 2 also include the 'tHooft determinant terms [32], the pseudoscalar kinetic energy receives an additionnal piece : where η 0 is the singlet PS field and λ is a parameter to be fixed. This term imposes a second step renormalization [31] of the (π 0 , η, η ′ ) sector. The first step renormalization turns out to define the (step one) renormalized pseudoscalar field matrix P R 1 in term of the bare one P by [25] : Some pseudoscalar (PS) fields happen to be fully renormalized at first step; they are related to their bare partners by : These renormalized fileds have already been used to derive the Lagrangians given in Appendix A. As for the other PS fields, the first step renormalization leads to the R1 renormalized PS fields : The second step renormalization required by the 'tHoof term gives the R renormalized fields [31] : where : The parameters affecting the PS fields submitted to fit are thus ∆ A , z A and λ (or alternatively v).
x, x ′ and x ′′ reflect the various ways by which nonet symmetry breaking ocurs in the PS sector of broken HLS Lagrangians.

C The V V P Langrangian in Terms of Ideal Vector Fields
The π 0 , η and η ′ pieces of the L V V P Lagrangian are given separately just below. If the PS fields occuring in this Appendix are renormaized by the procedure remindid in Appendix B, it is simpler to express the various V V P Langrangian pieces in terms of the ideal vector fields. The derivation of the couplings involving the renormalized fields is given in the main text.
We split up below the V V P Langrangian into the pieces involving the renormalized π 0 , π ± , η and η ′ pieces. As for the π 0 , we have : To stay the most concise possible, we have found appropriate to let separate two contributions for the π 0 ∂ µ ρ I ν ∂ α ρ I β vertex (with respective relative weights ξ 3 /2 and g 0 ρπ 0 ). In any physical application, they should obviously be summed up.
The V V η Lagrangian is given by : with : The V V η ′ Lagrangian is given by : Finally, it is worth extracting out from Eq. (134) the part involving a charged ρ meson for illustrative purposes. This write : D The AAP and AP P P Lagrangians in the HLS framework For the reader convenience, we have found appropriate to remind the expressions for the AAP Lagrangian and for the part of the AP P P Lagrangian relevant for the annihilation process e + e − → π 0 π + π − . These have been derived in [25]. One first have : with : using the renormalized PS fields defined in Appendix B and their specific parameters.

E The V P P P Langrangian And Its Renormalization
The V P P P anomalous Lagrangian in terms of ideal vector fields can be written : where one has limited oneself to display the P 0 π + π − sector. The leading terms of these couplings are : One may note, once again the clarity reached by using δ P instead of θ P . The couplings to the renormalized fields can be derived from those to the ideal fields by means of the R(s) matrix defined in the body of the text : the elements of the g 0 P 0 vectors being defined in Equations (148), Equations (149) or Equations (150) for resp. the P 0 = π 0 or η, η ′ mesons.
For the couplings of the Φ meson, at leading order in breaking, one gets : One clearly sees that, at the Φ mass, one obtains a non-vanishing direct coupling of the Φ to 3 pions, additionally s-dependent.