Heavy Quark Masses from Sum Rules in Four-Loop Approximation

New data for the total cross section $\sigma(e^+e^-\to{hadrons})$ in the charm and bottom threshold region are combined with an improved theoretical analysis, which includes recent four-loop calculations, to determine the short distance $\bar{\rm MS}$ charm and bottom quark masses. A detailed discussion of the theoretical and experimental uncertainties is presented. The final result for the $\bar{\rm MS}$-masses, $m_c(3 {GeV})=0.986(13)$ GeV and $m_b(10 {GeV})=3.609(25)$ GeV, can be translated into $m_c(m_c)=1.286(13)$ GeV and $m_b(m_b)=4.164(25)$ GeV. This analysis is consistent with but significantly more precise than a similar previous study.


Introduction
The strong coupling constant and the quark masses are the fundamental input parameters of the theory of strong interaction. Quark masses are an essential input for the evaluation of weak decay rates of heavy mesons and for quarkonium spectroscopy. Decay rates and branching ratios of a light Higgs boson, suggested by electroweak precision measurements, depend critically on the masses of the charm and bottom quarks, m c and m b . Last not least, confronting the predictions for these masses with experiment is an important task for all variants of Grand Unified Theories. To deduce the values in a consistent way from different experimental investigations and with utmost precision is thus a must for current phenomenology.
Let us recapitulate the main ingredients of the sum rule approach which will be used in the following analysis. Originally the idea has been suggested for the analysis of the charm quark mass by the ITEP group [1] long ago. Subsequently the method has been developed further [2] and frequently applied to the bottom quark. Most of these later analyses concentrated on using relatively high moments which are less sensitive to the continuum contribution and exhibit a very strong quark mass dependence [3,4,5,6,7]. However, this approach requires the proper treatment of the threshold, in part the resummation of the higher order terms from the Coulombic binding and a definition of the quark mass adopted to this situation like the potential-or 1S-mass [8,9], which subsequently has to be converted to the MS-mass. The low moments, in contrast, can be directly expressed by the short-distance mass at a scale around 3 GeV and 10 GeV for charm and bottom, respectively. The extrapolation to even higher scales, required for a number of applications is therefore insensitive to the larger corrections which would appear in the low energy region.
A detailed analysis of m c and m b , based on sum rules, has been performed several years ago [10] (see also Ref. [11]) and lead to m c (m c ) = 1.304 (27) GeV and m b (m b ) = 4.191 (51) GeV. This is still one of the most precise results presently available. During the past years new and more precise data for σ(e + e − → hadrons) have become available in the low energy region, in particular for the parameters of the charmonium and bottomonium resonances. Furthermore, the error in the strong coupling constant with its present value 1 α s (M Z ) = 0.1189 ± 0.0020 [12] (as compared to α s (M Z ) = 0.118 ± 0.003 for the last analysis) has been reduced. Last not least, the vacuum polarization induced by massive quarks has recently been computed in four-loop approximation [14,15]; more precisely: its first derivative at q 2 = 0 has been evaluated, which corresponds to the lowest moment of the familiar R-ratio. A fresh look at the determination of the quark masses based on these new developments is thus appropriate and will be presented below.
The outline of the paper is as follows: The basic assumptions of our approach are presented in Section 2. In Section 3 we discuss in detail the evaluation of the background and recall the theory input required for the quark mass determination. The measurements in the charm threshold region are discussed and applied to the determination of the charm quark mass in Section 4. Similar considerations are used in Section 5 to obtain the bottom quark mass. Section 6 contains our conclusions.

Generalities
The extraction of m Q from low moments of the cross section σ(e + e − → QQ) exploits its sharp rise close to the threshold for open charm and bottom production and the important contributions from the narrow quarkonium resonances. At best the properties of the lowest bottomonium state, Υ(1S), can be evaluated in perturbative QCD. In general a differential description of the cross section close to threshold involves necessarily low scales, of O(α s m Q ) or O(α 2 s m Q ) and even non-perturbative contributions, arising, e.g., from condensates G 2 , will enter. In contrast, by evaluating the moments with low values of n, the long distance contributions are averaged out and M n involves short distance physics only, with a characteristic scale of order E threshold = 2m Q . Through dispersion relations the moments are directly related to derivatives of the vacuum polarization function at q 2 = 0, which can be evaluated in perturbative QCD. From dimensional considerations one obtains for the relative error induced by the experimental uncertainties. (The theoretical and parametric errors will be described below.) Given δR/R ≈ 1.5 − 3%, a precision close to 10 MeV for m c and 20 MeV for m b seems within reach. An analysis at this level obviously requires to control a variety of corrections, e.g., contributions from the nonperturbative condensate in the case of m c . Furthermore, a careful definition of the heavy quark production cross section per se is required. We will discuss this now in detail for charm production, with the generalization to bottom production being obvious.
The narrow charmonium resonances J/Ψ, Ψ(2S) and the higher excitations will obviously contribute to the moments. Open charm production exhibits a sharp rise, nearly like a step function. Beyond the Ψ(3770)-resonance a few oscillations are observed which quickly level out into a fairly flat energy dependence. Around and above approximately 5 GeV the cross section is well approximated by perturbative QCD and mass terms can be considered as small corrections. The sensitivity to m Q is, therefore, concentrated on the small region from J/Ψ up to approximately 5 GeV. At first glance one might try to extract σ(e + e − → cc) by measuring the cross section with tagged charmed hadrons, thus eliminating the contribution from light quarks. However, such data are presently only available for selected narrow energy regions. Even more important, this approach is also not tenable from the conceptual side. A sizeable fraction of Ψ(3770) decays (not to speak of J/ψ or ψ ′ ) proceeds to non-charmed final states, nevertheless its total contribution should be attributed to the moments. On the other hand, processes with "secondary" charm production arising from gluon splitting (cf. Fig. 1(a)) should rather be assigned to the light quark cross section. Also cuts through singlet diagrams (cf. Fig. 2) which first appear in order α 3 s and where the fermion lines represent light or heavy quarks, exhibit thresholds at q 2 = 0, 4m 2 and 16m 2 , and these contributions cannot be attributed to a specific quark flavour. For the determination of m Q we therefore rely on the measurement of the total cross section, which is calculated through the absorptive parts of vacuum polarization diagrams. Various contributions to the cross section, arising from different cuts through the same diagrams, will always be treated together. For the moments M n we only consider those contributions to Π(q 2 ) which are analytic for |q 2 | ≤ 4m 2 Q and exhibit a cut starting at q 2 = 4m 2 Q . They result from diagrams with one massive quark loop connecting both photon vertices and will be denoted by Π Q and R Q , respectively. The remaining diagrams with cuts starting at q 2 = 0 will be treated as "background". This includes all diagrams with one light quark loop connecting both photon vertices (including those with internal massive quark loops) and all singlet diagrams. It will be shown below that the background defined in this way exhibits a fairly flat behaviour across the threshold, despite the fact that individual cuts with branching points at q 2 = 4m 2 Q exist. In the remaining part of this Section we fix the notation and give precise definitions of the quantities used in the paper. It is convenient to normalize the radiatively corrected hadronic cross section relative to the point cross section and to define the ratio where σ pt = 4πα 2 /(3s). As an inclusive quantity R(s) is conveniently obtained via the optical theorem from the imaginary part of the polarization function of two vector currents via where Π(q 2 ) is defined through with j µ being the electromagnetic current. The perturbative expansion of R(s) can be written as where the summation is performed over all active quark flavours i. (The small singlet contribution starts to contribute at order α 3 s .) For a comprehensive compilation of the individual pieces we refer to [16,17,18] where also explicit results are given. The full quark mass dependence is available up to order α 2 s [19,20]. For the non-singlet term R (3) i the first three terms of the high-energy expansion are known [21]. The predictions for R(s) which are used to calculate the background and the continuum term far above charm or bottom threshold, respectively, are based on Eq. (7) where the up, down and strange quark masses are taken as massless. For the charm and bottom quark the respective pole masses are chosen as input. If not stated otherwise we will use the following parameters for the evaluation of the background α which cover the full range of all currently accepted results. In Eq. (8) M c and M b refer to the pole masses. At several places of our analysis the renormalization group functions and the matching conditions for the masses and the strong coupling are needed in order to get relations between different energy scales. The corresponding calculations are performed using the package RunDec [22].

The background
As stated above, we distinguish three energy regions: First, the region of the narrow resonances J/ψ and ψ ′ , second, the "charm threshold region" starting from the D-meson threshold at 3.73 GeV up to approximately 5 GeV, where the cross section exhibits rapid variations and, third, the continuum region where pQCD and local duality are expected to give reliable predictions. In this last region the cross section is mainly sensitive to α s and fairly insensitive to m c .
For the present analysis 2 we use the data from the BES collaboration published in 2001 [23] and 2006 [24]. Whereas the older data cover the whole range from 2.0 GeV to 4.8 GeV the newer ones provide a very precise measurement in the region around √ s ≈ 3770 GeV. In Fig. 3 the BES-results are shown together with the measurements from the MD-1 [25] and CLEO [26] experiments between 4.8 GeV and 10.52 GeV. As is evident from Fig. 3, pQCD with α (5) s (M Z ) = 0.1189 provides an excellent description of all recent results. For example, the recent BES-value of R exp = 2.141±0.025(stat.)±0.085(syst.), determined in the interval [3.650GeV, 3.732GeV] is in excellent agreement with the theoretical prediction of R(3.700 GeV) = 2.161±0.017 and the same is true for R exp (10.52 GeV) = 3.56 ± 0.01 ± 0.07 from CLEO [26] to be compared with R(10.52 GeV) = 3.548 ± 0.012.
Up to O(α 2 s ), the prediction for charm quark production incorporates the full quark mass dependence. Starting from order α 2 s also the charm quark mass dependence of "secondary" charm production has to be taken into account. This includes diagrams of the type in Fig. 1(a) as well as those from Fig. 1(b). In addition we include O(α 3 s ) terms from the expansion in (M 2 c /s) n with n = 0, 1 and 2 [16,21]. Last but not least contributions from virtual c quarks ( √ s ≤ 3.73 GeV) and b quarks ( √ s ≤ 10.52 GeV) are included. Since our formulae are expressed in terms of α (4) s , the latter are suppressed ∼ (α s /π) 2 s/(4M 2 b ) and decouple for s ≪ 4M 2 b . These are included in the forth column of Tab. 1. Pure QED final state radiation is tiny and taken into account for completeness. (For a related analysis at √ s = 10.5 GeV see [27].) Let us now discuss in detail the evaluation of the "background". The results will be used to extrapolate R uds from the region below charm threshold up to 4.8 GeV. Here we provide all formulae relevant for the charm threshold region; the modifications to the J/ψ  through the threshold the effective number of flavours will be chosen to be n f = 4 and virtual charm quark effects are taken into account (for a compilation of the relevant formulae see Refs. [16,17,18]). Specifically we define where we take for the purely light degrees of freedom with L sµ = ln(s/µ 2 ) and α s ≡ α (4) s (µ). Contributions with virtual or secondary charmed quarks appear first in order α 2 s and in this order they are known in closed analytical form [28,29,30]: where ρ V and ρ R can be found in Eqs. (E.15) and (E.16) of Ref. [18]. 3 The first three expansion terms for small and large center-of-mass energy are given by δR uds(c) exhibits a logarithmic divergence for small s, which could be removed by matching to a theory with three effective flavours. However, restricting ourselves to √ s > 2 GeV this logarithm is numerically not large and we can restrict ourselves to n f = 4. The corresponding α 3 s -terms with one or two internal heavy quark loops, which in Eq. (11) is parameterized via δR δR with L sMc = ln(s/M 2 c ). In Fig. 4 the exact result for δR (2) uds(c) is compared with the approximations of Eqs. (12) and (13). As can be seen, the approximations work excellent even fairly close to threshold, above and below. The corresponding plot for δR (3) uds(c) is shown in Fig. 5. The dash-dotted curve corresponds to the approximation in Eq. (15) and the dashed one to Eq. (14). One can see a good agreement of the two approximations for √ s > ∼ 3.5 GeV. Since for √ s < ∼ 3.5 GeV the lower dashed curve should be very close to the exact result we adopt the result in Eq. (14) for the whole interval 2 GeV < ∼ √ s < ∼ 5 GeV. It is interesting to note that the numerical effect of δR (3) uds(c) is comparable to the one of δR (2) uds(c) . However, their overall size is small, as can be seen in Tabs. 1 and 2 (see below).
The singlet contributions start at order α 3 s and are fairly small. In the limit s ≫ 4M 2 c R sing can be cast into the form [18] R sing There are no quadratic mass terms in the singlet contribution. Both terms in Eq. (16) are proportional to Q 2 c , since the sum over the light quark charges happens to vanish.  The expansion for small √ s has been considered in Refs. [31,32] in the context of the Z boson decay. It is straightforward to adopt Eq. (20) of [31] to our situation for R sing . It turns out that the lowest four expansion terms for s ≪ 4M 2 c are proportional to Q u + Q d + Q s = 0. Thus the first non-zero contribution is proportional to Q 2 c and arises at order (s/M 2 c ) 4 . This term can be found in Refs. [33,32]. However, its contribution is numerically very small and can safely be neglected in our context.
Combining the information about the small and large energy region for R sing one obtains that its numerical contribution is between −0.55(α s /π) 3 and zero and can thus be safely neglected. Similar arguments also hold for the bottom case.
The magnitude of the most important terms contributing to R background is shown for several characteristic energies in Tabs. 1 and 2 for both the charm and the bottom region. In Tab. 3 they are compared to a few selected experimental results. Note that in contrast to R background a three-(four-)flavour theory has been used for the √ s values below (above) the charm threshold region. All theory predictions are given for the reference values specified in Eq. (8).  The band in Fig. 3 and the uncertainty in Tab. 3 is obtained from the parametric errors and the theoretical uncertainty, the latter obtained by varying the renormalization scale between µ 2 = s/4 and µ 2 = 4s. These were added quadratically. The excellent agreement between prediction and measurement in the continuum justifies to use the background formula, normalized to data below threshold also in the region above to obtain as remainder R c which will be used for the quark mass determination.
In passing let us mention that the formulae in combination with the R measurement between 3 GeV and 10.5 GeV can be used to determine a fairly precise value for α s , as demonstrated in Ref. [10].

Charm quark mass determination from the threshold region
Our determination of the charmed quark mass follows closely Ref. [10]. It is based on the direct comparison of the theoretical and experimental evaluations of the contributions to the derivatives of the polarization function Π c (q 2 ) (as defined in Eq. (6)), the former evaluated in perturbative QCD, the latter through moments of the measured cross section  for charm production in electron-positron annihilation. In its domain of analyticity Π c (q 2 ) can be cast into the form with Q c = 2/3 and z = q 2 /(4m 2 c ) where m c = m c (µ) is the MS charm quark mass at the scale µ. The perturbative series for the coefficientsC n in order α 2 s was evaluated up to n = 8 in Ref. [19,20] (recently [34] even up to n = 30). The four-loop contributions toC 0 andC 1 were calculated in Ref. [14,15]. The coefficients depend on α s and on the charm quark mass through logarithms of the form l mc ≡ ln(m 2 c (µ)/µ 2 ). They can be written as In Tab. 4 the coefficientsC (ij) n are given in numerical form up to n = 4. The terms C (30) n with n ≥ 2 are still unknown. To estimate the theoretical uncertainly from this source we will use the following bounds on C (30) n : The lower limit C (30) n > −6.0 is based on the observation that for fixed order in α s the absolute values of C (20) n , C (10) n and C (0) n in only one case exhibit a slight increase, in general, however, they tend to decrease with increasing n ≥ 1 quite markedly. The upper limit is estimated from the dependence on the order in α s with n fixed, where we assume at maximum a geometric increase and thus C (30) n ≤ (C (20) n ) 2 /C (10) n . In our numerical analysis we also include the two-loop QED corrections to C n which can easily be obtained from C (20) n . It leads to a negative shift in m c by less than 2 MeV. Tab. 4 and Eq. (19) essentially constitutes our theoretical input. We define the moments  Table 4: Coefficients of the photon polarization function in the MS scheme as defined in Eqs. (17) and (18). n f = 4 has been adopted which is appropriate for the charm threshold.
which leads to As demonstrated in Ref. [10], high moments (n > 4) are less suited for the mass extraction. The analysis will therefore be limited to n < 4 and the results for n = 4 will only be presented for illustration.
With the help of the dispersion relation we establish the connection between the polarization function and the experimentally accessible cross section R c (s). In the MS scheme which allows to determine the experimental moments in Eq. (1). Note, that the last term in Eq. (22) which defines the renormalization scheme disappears after taking derivatives with respect to q 2 . Equating Eqs. (21) and (1) leads to an expression from which the charm quark mass can be obtained: Non-perturbative contributions to the moments have been evaluated in [1,35] and, in view of their smallness 5 , were neglected in Ref. [10]. With the further reduction of theoretical and experimental uncertainties this approximation has to be reconsidered.
where we use the parameters listed in Tab. 5.
In the charm threshold region (which includes Ψ(3770)) we have to identify the contribution from the charm quark, i.e. we have to subtract the parts arising from the light u, d and s quark. (As discussed in Section 3, charm production through secondary gluon splitting is treated as part of the background, the same applies to the singlet contribution.) Technically this is done by determining a mean value for R background (cf. Eq. (9)) from the comparison of theoretical predictions and the BES data between 2 GeV and 3.73 GeV and using the theoretically predicted energy dependence to extrapolate into the region between 3.73 GeV and 4.8 GeV [39]. The mean value at the energy just below the DD-threshold is given byR with where the sums run over the s-values where experimental data, R exp , are available. In our case the upper bound is given by √ s − = 3.73 GeV and R th corresponds to the theoretical predictions from R background (s). In a next step we multiply R background (s) with the normalization factor n − and subtract the result for each s-value in the energy region [3.73 GeV, 4.8 GeV] from the data before the integration is performed. The final result for the charm quark mass m c (3 GeV) changes by about −5 MeV if we do not take this renormalization into account. The effect of the energy dependence of R background has practically no effect on the final result (less than 1 MeV). The normalization factors n − = 1.038 for the 2001 data [23] and n − = 0.991 for the 2006 data [24] can be considered as systematic offset of the experimental R measurement, as determined from the measurement below threshold and is determined under the assumption that the energy dependence of the true R-ratio is indeed given by R th . The magnitude of this factor is consistent with the quoted systematic error of 4.3% [23] and of 4.0% (4.9%) outside (within) the Ψ(3770)-region [24] which is kept in the subsequent analysis. The statistical error of n − is negligible. In the continuum region above √ s = 4.8 GeV data are sparse and imprecise. On the other hand, pQCD provides reliable predictions for R(s), which is essentially due to the knowledge of the complete mass dependence up to order α 2 s [19] and the dominant terms of order α 3 s [21]. Thus in this region we will replace data by the theoretical prediction for R(s) as discussed in the beginning of Section 3. As emphasized above, the cross section is only weakly m c -dependent in this region.
In Tab. 6 we present the results for the moments separated according to the three different contributions discussed above. The error of the resonance contribution is due to the uncertainties of the input parameters. The error of the charm threshold contribution is dominated by the correlated normalization error of approximately 4.0% to 4.9% of the BES data.
To estimate the error on M cont The higher moments (in fact already for n above two) are increasingly dominated by the resonance contributions with their 2.5% uncertainty. The large uncertainty in α s /πG 2 starts to matter for n ≥ 3. We use the results of Tab. 6 together with Eqs. (23) and (24) in order to obtain in a first step m c (3 GeV). The results are listed in Tab. 7 together with the total uncertainty and its decomposition into the individual contributions from the experimental moments, the variation of α s (cf. Eq. (8)), the renormalization scale µ = (3±1) GeV and the condensate which are all added quadratically. To determine the µ-dependence, we evaluate Eq. (23) at µ = 2 GeV and 4 GeV. Subsequently, using the renormalization group equation to four-loop accuracy, we evolve the result for m c (µ) to the reference point µ = 3 GeV and compare with the original value.
The first moment is available in four-loop approximation, and the corresponding µdependence of the result is completely negligible. It is interesting to anticipate the corresponding µ-dependence for the four-loop result also for the higher moments. The corresponding results are displayed in Tab. 7. In the eighth column we show the effect of the unknown coefficientsC (30) n (n ≥ 2) where assumptions were made on their respective numerical growth as discussed above. For the extraction of m c (3 GeV) we assumeC (30) n = 0. For δC (30) n we take the larger of the deviations estimated in Eq. (19). Let us compare the results displayed in Tab. 7 with those from [10]. For the moment with n = 1 the main improvement originates from the new experimental results for the electronic widths of the J/Ψ and Ψ ′ . Some additional improvement is due to the new BESdata [24]. The four-loop term leads to a further reduction of the theoretical uncertainty -in view of the dominant experimental uncertainties it only leads to consolidation of the theoretical prediction. The error composition in the analysis based on n = 2 is similar, as far as the experimental input is concerned. However, in this case the error from the R-measurements is somewhat smaller, the error induced by the uncertainty in α s is larger. Furthermore, the error from δC (30) 2 has to be taken into account. Although at first glance this error looks small, a check of the assumptions on C (30) 2 would be extremely useful. The analysis based on n = 3 leads to a result consistent with the one for n = 1 and n = 2. In this case, however, the theoretical uncertainty estimated through the µ dependence starts to become important and the total error is therefore larger.
For the subtraction of the light-quark continuum we have assumed pQCD to be exact. Deviations from this idealized situation have been estimated in Ref. [40]. Based on theoretical and phenomenological considerations it is assumed, that the oscillations of R(s) around the perturbative result, which are observed in the low-energy region and reflect the low-lying resonances, are damped by a power-law or an exponential factor. The following two correction factors were suggested in [40] with ρ = 3 GeV −1 , δ = 1.32 and, alternatively, with σ 2 = 2 GeV 2 , B = 0.5 and N c = 3. Potential effects from this modification are estimated by including the factors f 1 or f 2 into the analysis. 6 In the case of f 1 we also vary the phase δ between 0 and π. The change in the contribution from the region 3.73 GeV ≤ √ s ≤ 4.8 GeV is about 1.2% for f 1 and n = 1 and significantly smaller in the case of f 2 . The 1.2% translate to a shift of −1 MeV for m c (3 GeV) which is completely negligible in the present context. The results of Tab. 7 are also shown in Fig. 6 where the central value and the uncertainties of m c (3 GeV) are plotted for n = 1, 2, 3 and 4. For each n the results are shown for one-to four-loop theory input. For the O(α 3 s ) term the exact result is used for n = 1 whereas for n = 2, 3 and 4 we use C (30) n = 0 and the error estimates discussed above. It is nicely seen that the results for m c (3 GeV) further stabilize when going from three to four loops. At the same time the uncertainty is considerably reduced. Furthermore, the preference for the first three moments is clearly visible. Also the analysis for n = 2 and n = 3 leads to small errors, even if we include to contribution from the yet uncalculated δC (30) n . We emphasize the remarkable consistency between the three results which we consider as additional confirmation of our approach. The result based on the moment with n = 1 is evidently least sensitive to nonperturbative contributions from condensates, to the Coulombic higher order effects, the variation of µ and the parametric α s dependence. The results for n = 2 and n = 3 are practically identical. Since the error in this case is somewhat larger, and furthermore affected by the estimate for δC (30) n , we adopt m c (3 GeV) = 0.986(13) GeV , Our result agrees within the uncertainties with our last determination [10], but is considerably more precise. Using the three-loop relation [41,42,43] Table 8: Predictions for m c (m c ). In the second column some keywords concerning the used method are given (NNLO: next-to-next-to-leading order; NNNLO: next-to-next-tonext-to-leading order).  Table 9: Coefficients of the photon polarization function in the MS scheme as defined in Eqs. (17) and (18) for n f = 5.
M c since it is well known that the perturbative series between the MS and the pole mass is badly behaved. A comparison with a few selected determinations (published 2001 or later) is shown in Tab. 8 and Fig. 7. Within their respective errors all results are consistent. Let us also mention that a first investigation of the numerical effect of the four-loop results for the moments has been performed in Ref. [15] where basically the analysis of Ref. [10] has been adopted with updated parameters. The result of Ref. [15] reads m c (m c ) = 1.295(15) GeV.

The bottom quark mass
The same approach is also applicable to the determination of m b . The coefficientsC n are listed in Tab where the same criteria as in the charm case have been applied. The coefficients in Tab. 9 determine the theoretical moments through Eq. (21) where Q c has to be replaced by   neglected. We included the charm mass terms [49] of order (m c /m b ) 2 to the momentsC n which induce a shift of about −1 MeV in the bottom quark mass. The experimental results for the moments are listed in Tab. 11. The contribution from the resonances include Υ(1S) up to Υ(4S). The values for the masses of the Υ resonances and their electronic width have been taken from Ref. [13] and are listed in Tab. 10. The errors from the three lowest Υ resonances have been combined linearly, since the PDG values for the electronic widths are dominated by the measurement from CLEO [50]. The result is then combined in quadrature with the contribution from the Υ(4S) resonance. The treatment of the bottom threshold region is similar to the one of the charm region. Measurements of R from threshold up to 11.24 GeV have been performed by the CLEO Collaboration more than 20 years ago [51], with a systematic error of 6%. ×10 (2n+1) ×10 (2n+1) ×10 (2n+1) 1 1.394 (23) 0.296 (32)   No radiative corrections were applied. The average value derived from the four data points below threshold amounts toR = 4.559 ± 0.034(stat.) which is 28% larger than the prediction from pQCD. However, a later result of CLEO [26] at practically the same energy, R(10.52 GeV) = 3.56 ± 0.01 ± 0.07, is significantly more precise and in perfect agreement with theory. Applying a rescaling factor of 1/1.28 to the old CLEO data not only enforces agreement between old and new CLEO data and pQCD in the region below the Υ(4S), it leads, in addition, also to excellent agreement between theory and experiment above threshold around 11.2 GeV where pQCD should be applicable also to bottom production. Further support to our approach is provided by the CLEO measurement of the cross section for bottom quark production at √ s = 10.865 GeV which is given by  [52]. In Fig. 8 the original and the rescaled data from [51] is shown and compared to pQCD and data point from [26]. We thus extract the threshold contribution to the moments from the interval 10.62 GeV ≤ √ s ≤ 11.24 GeV by applying the rescaling factor to the data, subtract the "background" from u, d, s and c quarks and attribute a systematic error of 10% to the result. To evaluate the contributions to the moments from above 11.24 GeV, i.e. M cont n , we use the prediction from pQCD for R b (s).
The results for the moments are listed in Tab. 11, those for m b (10 GeV) and m b (m b ) in Tab. 12. Just as in the charm case, a remarkable consistency and stability is observed. For n = 1 the error is dominated by the experimental input. For n = 3 we obtain ±0.010 from the experimental input, ±0.014 from α s and ±0.006 from the variation of µ.
The sensitivity to the inclusion of higher orders is displayed in Fig. 9. For the lowest moment the error is dominated by the experimental uncertainty -nevertheless the theory error is reduced also in this case and the prediction stabilized. In general a significant improvement and stabilization is observed.
The three results based on n = 1, 2 and 3 are of comparable precision. The relative size of the contributions from the threshold and the continuum region decreases for the moments n = 2 and 3. On the other hand, the theory uncertainty estimated from the variation of µ and δC (30) n is still acceptable. We therefore take the result from n = 2 (which is roughly between the n = 1 and n = 3 values) and add the uncertainty from "total"   [51] (circles) and [26] (triangle). The black curves are the predictions from pQCD outside the resonance region. In the lower plot the older data from [51] are rescaled by a factor 1/1.28.
(cf. Tab. 12) and the one induced by δC (30) n linearly which leads to an error estimate of ±25 MeV: This result can also be converted to a pole mass [41,42,43] .800 GeV. Again we do not display the additional error from the MS-pole-mass conversion, which is significantly larger.
A comparison with a few selected determinations is shown in Tab. 13 and Fig. 10. In Ref. [15] where the first error originates from the combined error listed in Eq.   Table 13: Predictions for m b (m b ). In the second column some keywords concerning the used method are given (NNLO: next-to-next-to-leading order; NNNLO: next-to-next-tonext-to-leading order; NNLL: next-to-next-to-leading-logarithmic order).
In connection with Yukawa unification in Grand Unified Theories the value m b (µ = m t ) in the n f = 6 theory may be of interest. Given M t = 171.4 ± 2.1 [57], which we identify within the present uncertainties with the pole mass, we find m t (m t ) = 161.8 ± 2.0 ± 0.2 , to three-loop accuracy [41,42,43]. (Again we do not specify a theory uncertainty. The difference between the two-and three-loop result amounts to 426 MeV.) Evolving m b from µ = 10 GeV to m t (m t ) and matching to the n f = 6 theory we find m b (161.8) = 2.703 ± 0.018 ± 0.019 , where the first error reflects the combined error from Eq. (37) and the second one the uncertainty due to α s . As stated above, the ratio  should be a useful input for Grand Unified Theories. These values are consistent with our previous determination [10] but considerably more precise. Let us stress that, since m c is relatively small, m c (m c ) shows a less stable behaviour against higher order corrections. Thus we propose to use m c (3 GeV) as the better choice for comparisons of the various methods. The evaluation of the running bottom mass at the scale of m t leads to m t (m t )/m b (m t ) = 59.8 ± 1.3, a result of interest for theories with Yukawa unification. We want to mention that for the charm quark mass our analysis constitutes the only one where NNNLO corrections from theory are incorporated. In the case of the bottom quark there is only one further NNNLO analysis which is based on the mass of the Υ(1S) system. However, the accuracy of this approach is strongly limited due to large non-perturbative contributions.