Charming penguins and lepton universality violation in b → s (cid:2) + (cid:2) − decays

The LHCb experiment has recently presented new results on Lepton Universality Violation (LUV) in B → K ( ∗ ) (cid:2) + (cid:2) − decays involving K S in the ﬁnal state, which strengthen the recent evidence of LUV obtained in B + → K + (cid:2) + (cid:2) − decays and the previous measurements of B → K ∗ 0 (cid:2) + (cid:2) − . While LUV observables in the Standard Model are theoretically clean, their predictions in New Physics scenarios are sensitive to the details of the hadronic dynamics, and in particular of the charming penguin contribution. In this work, we show how a conservative treatment of hadronic uncertainties is crucial not only to assess the signiﬁcance of deviations from the Standard Model but also to obtain a conservative picture of the New Physics responsible for


Introduction
The LHCb experiment at the Large Hadron Collider has announced evidence of Lepton Universality Violation (LUV) crowning with success a huge experimental effort aimed at detecting deviations from the Standard Model (SM) in rare B q decays. Recently, another piece was added to the already very rich set of data on (semi)leptonic and radiative B q decays: the measurements of [ complementing the analogous search for LUV in B → K * + − decays, R K * [3,4], 1 the measurement of BR(B s → μ + μ − ) [9][10][11][12][13], the angular analyses and BR measurements of B → K ( * ) μ + μ − [14][15][16], B s → φμ + μ − [17][18][19][20] and B → K * e + e − [21,22]. While hadronic uncertainties make the detection of possible New Physics (NP) contributions to B q → K ( * ) + − differential rates very difficult, at least with current data, any observation of LUV beyond the percent level due to QED corrections [23,24] would be a clean signal of NP.
Eur. Phys. J. C (2023) 83:64 In the SM, b → s + − transitions can only arise at the loop level, as all other Flavour Changing Neutral Current (FCNC) processes, and thus they are particularly sensitive to NP. The leading diagrams giving rise to these transitions are illustrated in Fig. 1. Given the hierarchy in the CKM angles, one has V ub V * us V cb V * cs ∼ V tb V * ts , making the contribution of virtual up-type quarks in the loops negligible. The GIM mechanism works differently for different diagrams: Z -penguins and boxes vanish as m 2 q /m 2 W and are therefore dominated by the top quark, while photonic penguins have a logarithmic dependence on the quark mass, allowing for a large contribution by the charm quark.
Another fundamental difference between the two classes of diagrams is due to the chirality of the weak couplings: Zpenguins and boxes involve both vector and axial couplings to leptons, while photon penguins couple vectorially to leptons. This implies that the top-dominated Z -penguins and boxes give rise to the local operators Q 9V ∼bγ μ P L s¯ γ μ and Q 10 A ∼bγ μ P L s¯ γ μ γ 5 at the electroweak scale. Photonic penguins instead are more complicated: the top quark contributes to Q 9V at the electroweak scale, but the charm quark remains dynamical at the scale m b and therefore contributes to b → s + − transitions both via the local operator Q 9V and via the (potentially nonlocal and nonperturbative) matrix elements of current-current operators involving the charm quark, Qb ccs 1,2 ∼bγ μ P L ccγ μ P L s, denoted by charming penguins [25][26][27]. This complication, however, does not affect axial lepton couplings, which remain purely short-distance.
Computing the matrix element of Qb ccs 1,2 is a formidable task. While the calculation of decay amplitudes for exclusive b → s + − transitions is well-defined in the infinite b and c mass limit [28][29][30], and while in the same limit the uncertainty from decay form factors can be eliminated by taking suitable ratios of observables [31,32], in the real world amplitude calculations must cope with power corrections [33,34], which can be sizable or even dominant in several kinematic regions [35][36][37][38][39]. For example, the Operator Product Expansion is known to fail altogether for resonant B → K ( * ) J/ψ → K ( * ) μ + μ − transitions [40], and its accuracy is questionable close to the cc threshold. For this reason, estimating corrections to QCD factorization (QCDF) in the low dilepton invariant mass (low-q 2 ) region of B → K ( * ) + − and B s → φ + − decay amplitudes is a crucial step towards a reliable assessment of possible deviations from SM predictions in these decay channels. Unfortunately, first-principle calculations of these power corrections are not currently available, and a theoretical breakthrough would be needed to perform such calculations, see, e.g., the discussion in [34,41,42].
Therefore, as of now, a conservative analysis of semileptonic B decays can only rely on the use of data-driven methods to account for the theoretical uncertainties and to quantify possible deviations from the SM. In this regard, it is important to stress that while the contribution of Q 9V to the decay amplitude should depend on helicity and on q 2 according to the form factors, long-distance contributions from the charm loop matrix element should show some additional helicity and q 2 dependence. It is then very interesting to use a q 2 and helicity-dependent parameterization of power corrections [33,34] when analyzing experimental data: A sizable deviation from what expected from purely local matrix elements would be a clear confirmation of the presence of power corrections.
Obviously, the charm loop matrix element cannot generate any Lepton Universality Violation (LUV), so that ratios of decay Branching Ratios (BRs) for different leptons in the final state can be very reliably predicted in the SM [23,24,43,44]. However, once lepton non-universal NP is introduced, the hadronic uncertainty related to the charm loop creeps back in, due to the interference between SM and NP contributions in decay amplitudes, so that the inference of NP parameters from LUV observables is not free from hadronic uncertainties. Thus, the relevance of a careful and conservative treatment of hadronic uncertainties in assessing the compatibility of b → s + − data with the SM, and in inferring what kind of NP could lie behind the evidence of LUV, cannot be overstated.
In this work, we focus on b → s + − transitions in a bottom-up perspective. Building on our previous analyses [35,36,38,[62][63][64] and on the data presented in [1][2][3][4][9][10][11][12][13][14][15][16][17][18][19][20][21][22], we aim at answering two fundamental, and deeply related, questions: Fig. 2 The result for the optimized observable P 5 (see Ref. [73]) from a global fit within the SM in three different approaches for the charming penguin contribution, compared to experimental data. Only the approach fully relying on LCSR estimates (Refs. [71,72]) and their extension by means of dispersion relations in the entire dilepton invariant-mass range may lead to the so-called "P 5 anomaly" 2. What is the overall significance for NP in light of the new data, and how do hadronic uncertainties affect the interpretation of the present evidence of LUV?
To this end, following the strategy we originally proposed in Ref. [35], in the following we consider a generic parameterization for non-factorizable QCD power corrections, and let data determine the q 2 and helicity dependence of charming penguins. In particular, in our analysis we determine short-distance and long-distance contributions in a simultaneous fashion exploiting the constraining power of the up-todate experimental information on b → s + − together with the current theoretical knowledge of hadronic correlators and matrix elements for this class of rare B meson decays.
Notice that the approach recently followed in Ref. [65] to allow for an arbitrary lepton-universal correction ΔC U 9 to C 9 is less general, and therefore less conservative, than our approach, unless ΔC U 9 is promoted from a parameter to a q 2and helicity-dependent function.
The paper is organized as follows: in Sect. 2 we present our parameterization for charming penguins and discuss the implications of current data on QCD long-distance effects; in Sect. 3 we present a global analysis of NP effects using both the Standard Model Effective Field Theory (SMEFT) and the Weak Effective Hamiltonian; in Sect. 4 we wrap up the present study with our conclusions.

Charming penguins from current data
For the convenience of the reader, let us briefly summarize our approach to hadronic uncertainties. We write down the helicity-dependent SM decay amplitudes for B → K * + − in the following way [33,66]: with λ = 0, ± and C SM 7,9,10 the SM Wilson coefficients of the operators Q 9V , Q 10 A and Q 7γ ∼ m bbR σ μν F μν s L normalized as in Ref. [63].
The factorizable part of the amplitudes corresponds to seven independent form factors, V 0,± , T 0,± and S, smooth functions of q 2 [67,68]. Instead, h λ (q 2 ) represents the nonfactorizable part of the amplitude [34,35,69], dominated by the charming penguin contribution. Nonperturbative methods working in Euclidean spacetime such as lattice QCD cannot directly evaluate h λ (q 2 ) as at present there is no way to evade the Maiani-Testa no-go theorem [70], which prevents the computation of rescattering and final state interactions away from the threshold. In particular, rescattering from an intermediate D s −D ( * ) state, or from any other on-shell state with flavour content ccsd, into a K ( * ) and a (virtual) photon is currently not computable. Notice that this intermediate state is always kinematically accessible since D s , irrespective of the dilepton invariant mass. Such contribution is intrinsically nonperturbative (it could become computable, and negligible, in the m b m c limit) and could give a large imaginary part to the charming penguin already at low q 2 . Another kind of charming penguin corresponds instead to the decay of the B meson into a K ( * ) and a virtual cc pair, which decays into the lepton pair through the exchange of a photon. This class of charming penguins, giving rise to singularities in q 2 , can be estimated using light-cone sum rules for q 2 m 2 c or for Euclidean q 2 Eur. Phys. J. C (2023) 83:64 (7) in the three approaches for charming penguins. For comparison, the result obtained in the data-driven approach with 2015 data is also reported, along with the short-distance contribution and the factorizable QCD corrections [41,71,72]. To summarize, we can only count on a modeldependent estimate of a subclass of charming penguins in a limited physical region of small q 2 . Given our ignorance of the charming penguin amplitude, we parameterize the hadronic contribution as follows [38]:

Fig. 3 95% probability contours of the posteriors for the functions
which allows us to write the helicity amplitudes as: The parameterization in Eq. (5) is a variation of a simple Taylor expansion formally in q 2 /(4m 2 c ), see [34], whose radius of convergence (for the leading contribution expected from the cc loop) is up to the J/Ψ resonance. While this is clearly not the only possible way to parameterize the hadronic contribution, such a choice has been originally introduced in Ref. [38] with two specific goals in mind. First, from Eq. (6) above it is evident that h (0) − is equivalent to a shift in C 7 , i.e. ΔC 7 , while h (1) − corresponds to a lepton universal correction ΔC 9 . Second, the remaining h parameters appearing in (6) are not equivalent to a shift in the Wilson coefficients of Q 7γ,9V and thus they represent genuine hadronic effects. One could take into account the presence of cc resonances by adding the corresponding poles to the parameterization. However, when Taylor expanding the full parameterization around q 2 = 0, this would just amount to a redefinition of the h (i) λ coefficients. 2 As discussed in detail in [35,71], the hadronic contributions introduced above correspond to the following q 2 -and helicity-dependent shifts in C 9 : An analogous parameterization can be introduced for charming penguins in B → K + − , as well as for B s → φ + − . When considering hadronic contributions in B → K * decays versus those in B s → φ decays the consideration of flavour SU(3) breaking comes into play. Given that the degree of SU(3) breaking originating from m s m d , m u is not quantifiable ab initio, we performed some tests by adding ad hoc SU(3) breaking parameters that multiplicatively modifies the hadronic terms in B s → φ decays vis-à-vis those in B → K * decays, i.e., h B s →φ λ Eur. Phys. J. C (2023) 83:64 The two parameters, δ R and δ I , are taken for simplicity to be independent of helicity and simply modify the real and imaginary parts of the hadronic terms independently. We set a Gaussian prior on δ R and δ I with μ = 0 and σ = 0.3 representative of 30% SU(3) breaking centered at no SU (3) breaking. The posterior distribution of δ R and δ I are not significantly different from the prior distributions leading us to conclude that the experimental results are not precise enough to draw conclusion about SU (3) breaking within this hypothesis. For the rest of the discussion we assume exact flavour SU (3) symmetry for power corrections, which is justified given the current experimental uncertainties. In the future, realistic departures from the SU (3) limit could be probed by data and the present investigation of QCD effects could be further generalized along these lines in a straightforward manner.
Using the HEPfit code [74,75] and the form factors and input parameters used in Refs. [35,36,38,[62][63][64], we perform a Bayesian fit to the data in Refs. [9][10][11][12][13][14][15][16][17][18][19][20][21][22] within the SM, adopting three distinct approaches to account for power corrections to QCDF: 1. A model-dependent approach in which LCSR results are extended over the full range of q 2 taking analiticity constraints into account (in particular, the singularities of the amplitudes in q 2 ) [41,71,72,76,77], where the results either stem directly from a LCSR computation performed at q 2 ≤ 1 GeV 2 , or from an extension of these results in the whole phenomenological region by means of dispersion relations. Within this approach, the parameterization in Eq. (6) is replaced by the ones in Eqs. (7.14) and (7.7) of Ref. [71] for K * and K respectively. In this case, the size of charming penguins in B → K * + − is comparable to the ones provided by QCDF O(α s ) corrections C QCDF 9 , while the size of charming penguins in B → K + − is too small to be phenomenologically relevant. 2. A model-dependent approach in which the LCSR results of Ref. [71] are used to constrain the h λ -parameters in Eq. (5) only for q 2 ≤ 1 GeV 2 (with no specific structure assumed in relation to the singularities in q 2 ). Using the results of Ref. [72] instead of Ref. [71] would imply an even smaller hadronic contribution at q 2 ≤ 1 GeV 2 . Within this approach, the size of charming penguins in B → K * + − can depart from the LCSR estimate once away from the light-cone region, while charmloop effects in B → K + − are still regarded as negligible in light of the estimates in Refs. [71,72] [71], to the one from Ref. [38], and to the data-driven approach pairs, for which no theoretical estimate is currently available.
In Fig. 2, as a representative example, we show the outcome of the fit in the SM for the optimized observable P 5 [73] in the three different approaches to power corrections. As shown in the plot, no tension from data emerges from this observable in the data-driven approach or in the one from Ref. [38]. In general, as already noted in [35], an excellent fit to the data (except, of course, for LUV ratios and with the notable other exception of the time-integrated BR(B s → μ + μ − )) can be obtained within these two approaches. On the other hand, LCSR results extended to larger values of q 2 yield a poor fit of several BRs as well as tension in some of the angular observables, giving rise to the so-called P 5 anomaly. We note that for q 2 → 0, P 5 is only sensitive to the combinations ΔC 9,1±2 (q 2 ). Also, ΔC 9,1+2 (q 2 ) rapidly vanishes towards q 2 = 0 in our model-dependent approach based on Ref. [71]. This is the result of helicity suppression [33], expected to hold on the light cone, i.e. only at q 2 = 0 [71]. Therefore, in the model-dependent approach based on Ref. [38] we enforce such cancellation only at q 2 = 0. Moreover, in the data-driven approach we do not require it to happen. Consequently, in the first bin of P 5 in Fig. 2 the red and orange bands lie above the green one.
The result shown for P 5 highlights the major role played by long-distance effects. In Fig. 3, we further investigate this aspect showing the q 2 and helicity dependence of the charming penguin contributions. In the plot we show the 95% probability regions of the posteriors for the functions ΔC 9,i (q 2 ) obtained in the global fit in the SM under the three differ-ent approaches. In the same figure, as a guideline, we also show the size of the SM short-distance contribution to Q 9V , labeled by C SD 9 , as well as the size of the factorizable QCD corrections.
The posteriors of ΔC 9,i (q 2 ) in Fig. 3 display nonnegligible hadronic contributions -comparable in size to C SD 9 rather than C QCDF 9 -in the whole region of low dilepton mass probed by current data. This is not surprising since power corrections are naively expected to be larger than perturbative QCD corrections of O(α s /(4π)) [25][26][27]. A departure from LCSR expectations even at very low q 2 is hinted at in the data-driven approach, which matches the outcome of the approach based on Ref. [38] only for q 2 4m 2 c . As can also be seen in Fig. 3, by comparing the Data Driven determination of ΔC 9,i (q 2 ) from current data with the one obtained in our 2015 analysis [35] 3 it is evident how improved data on differential BR's allow for a much better knowledge of the charm contribution. In this respect, it will be interesting to see whether more precise data will bring stronger evidence of such hadronic effects. At present, hints of large hadronic contributions from data are still statistically mild, as can be read from Table 1. There, we report the highest Eur. Phys. J. C (2023) 83:64 probability density intervals (HPDI) for the posteriors of the h parameters adopted in the two more conservative approaches. Some h's corresponding to genuine hadronic contributions deviate from 0, but still only at the 2σ level.
To summarize, the data driven scenario stands out as our most conservative choice for a statistically robust inference on NP contributions in current b → s + − data, accounting also for possible large contributions from QCD rescattering. 4 In the following we take it as a reference, but for completeness we present results on NP also in the other two approaches.

New physics in B decays
While experimental data on BR's and angular distributions can be reproduced within the SM in both the data-driven approach and in the model-dependent one based on Ref. [38], reproducing the central values of the LUV ratios for B → K ( * ) + − , as well as the current measurement of BR(B s → μ + μ − ), undoubtedly requires physics beyond the SM.
Given the bounds from direct searches of NP at the LHC, it is reasonable to assume in this context that NP contributions would arise at energies much larger than the weak scale. Then, a suitable framework to describe such contributions is given by the SMEFT, in particular by adding to the SM the following additional dimension-six operators: 5 O L Q (1) where τ A=1,2,3 are Pauli matrices (a sum over A in the equations above is understood), weak doublets are in upper case and SU (2) L singlets are in lower case, and flavour indices are defined in the basis of diagonal down-type quark Yukawa couplings. Since in our analysis operators O L Q (1,3) 2223 always enter as a sum, we collectively denote their Wilson coefficient as C L Q 2223 . For concreteness, we normalize SMEFT Wilson coefficients to a NP scale Λ NP = 30 TeV and we only consider NP contributions to muons. 6 Matching the SMEFT 4 While large hadronic effects will lead to a conservative estimate of NP in C 9 , in our global analysis -where R K ( * ) ratios signal NP irrespective of the hadronic treatment -such conservative approach will also impact the inference of a non-vanishing NP contribution in C 10,μ,e . 5 Notice that these operators may be further generated at one loop via SM RGE effects, see, e.g., Refs. [78,79]. In addition, here we do not consider the possibility that, integrating out NP, one would generate sizable Qb ccs 1,2 as studied e.g. in [80,81]. 6 The focus on LUV effects in muons is mainly motivated by the ∼ 2.3σ tension of the SM with the current experimental average for the timeintegrated BR(B s → μ + μ − ).
with α e the fine-structure constant, v the vacuum expectation value of the SM Higgs field, λ t = V ts V * tb , and alignment in  (22) the down-quark sector assumed, i.e.
We perform a Bayesian fit to the data in Refs. [1][2][3][4][9][10][11][12][13][14][15][16][17][18][19][20][21][22] in several NP scenarios characterized by different combinations of nonvanishing Wilson coefficients. To perform model comparison of different scenarios, we compute the Information Criterion (IC) [83]: where the first and second terms represent mean and variance of the log likelihood posterior distribution, respectively. The first term measures the quality of the fit, while the second one counts effective degrees of freedom and thus penalizes more complicated models. Models with smaller IC should be preferred according to the canonical scale of evidence of Ref. [84], related in this context to (positive) IC differences. For convenience, we always report ΔI C ≡ I C SM − I C NP . As is evident from the discussion in Sect. 2, different assumptions on charming penguins yield different results on NP Wilson coefficients, since LUV ratios depend on the charm loop through the interference between NP and SM contributions. It goes without saying that a conservative inference on NP requires a conservative, i.e. data driven, estimate of charming penguins. Of course, since the SM reproduces much better experimental data in the data-driven approach (and also in the one based on Ref. [38]) than in the approach from Ref. [71], with the first scenario performing better than the second one, the variation ΔI C due to NP (or equivalently the significance of NP) will be the smallest in the data-driven approach and the largest in the one based on Ref. [71].
Let us first consider three very simple NP scenarios: scenario A, in which deviations can only arise in C 9 , scenario B, in which only C L Q 2223 can be nonvanishing, corresponding to C NP 9 = −C NP 10 , and scenario C, in which only C NP 10 is allowed to float. Already in these simple NP scenarios there are dramatic differences in the fit depending on the assumption on charming penguins. Under the data-driven approach (and the one based on Ref. [38]), scenarios B and C perform much better than scenario A, while the opposite is true in the approach from Ref. [71]. As reported in the top panel of Fig. 4, in the data-driven approach charming penguins can even interfere destructively with C NP 9 , allowing for a second solution for LUV observables with positive C NP 9 , albeit with a smaller ΔI C with respect to the solution with negative C NP 9 (see Table 2). 7 The correlation of C NP 9 with the hadronic parameters Re h (1) − for final states involving K * and K mesons is shown in Fig. 5. In the approach based on Ref. [71], scenario A is ideal since it allows to strongly improve the agreement of both LUV and angular observables, while in scenario B the constraint from B s → μ + μ − limits the improvement  Table 3), and scenario C cannot improve the agreement with angular observables at all. Going back to Fig. 2, note that the global fit of scenario A would be perfectly in line with the hadronic contributions estimated in Refs. [71,72], related to what reported in green in the figure; while for scenario B and, most importantly, for C, the prediction for the observable P 5 would confront well with data only in the instance of sizable charming-penguin amplitudes. Conversely, under the data driven hypothesis, scenarios B and C allow to reproduce all observables, including LUV and B s → μ + μ − , and therefore stand out as the preferred NP scenarios. The case where hadronic uncertainties are treated as in Ref. [38] is in a somewhat intermediate position, with scenario C somewhat disfavoured with respect to scenario B due to the constraints on the charming penguin at low q 2 . Obviously, as can be seen in Fig. 4, the p.d.f. for C NP 10 in scenario C is almost independent of the hadronic uncertainties, while the overall quality of the fit strongly depends on the charming penguins, since in this scenario one needs hadronic contributions to reproduce the angular distributions and BRs.
More general scenarios with two or more nonvanishing NP Wilson coefficients, such as scenario D, where C L Q 2223 and C Qe 2322 are allowed to float, or scenario E, where all the coefficients of the operators in Eq. (9) are turned on, are slightly penalized by the number of degrees of freedom unless the approach from Ref. [71] is considered, as can be seen from Table 2. For the reader's convenience, in Table 2 and in Figs. 6, 7 and 8 we present results for scenarios D and E also in the weak effective Hamiltonian basis through Eq. (10).
It is interesting to look at the shape of the probability density contours for the NP parameters in scenario D reported in Fig. 6. In the data-driven approach, C NP 9 is well compatible with 0, while a nonvanishing NP axial lepton coupling emerges from the experimental information coming from "clean observables". A slight preference for a nonvanishing C NP 9 is present in the model-dependent approach of Ref. [38], while a strong evidence for C NP 9 is obtained in the one based on Ref. [71], together with a slight hint of a nonvanishing C NP 10 . Therefore, Fig. 6 represents a clear example of how the choice of the parameterization for the charming penguins can strongly impact the inference of the underlying NP picture, allowing to go from a purely axial NP coupling to a purely vectorial one in the two extreme cases, with dramatic consequences for the model building related to B anomalies.
Concerning scenario E, it is worth noticing that in the case of the approach from Ref. [71], right-handed operators allow to improve the agreement with R K , given the current experimental hint for R K = R K * at the 1σ level, see the discussion in [63]. In the data-driven approach (and in the one based on Ref. [38]) this can be achieved also through the interplay of hadronic corrections with LUV NP (see Table  3). See Figs. 7 and 8 for a comparison of the posteriors for NP coefficients in scenario E.

Conclusions
We have presented a global analysis of the experimental data on b → s + − transitions from Refs. [1][2][3][4][9][10][11][12][13][14][15][16][17][18][19][20][21][22] under three different assumptions about the size and shape of the charming penguin contribution: a data-driven approach, a modeldependent based on Ref. [38] and a model-dependent one based on Ref. [71]. We have shown how current data point to helicity and q 2 dependence of the charm loop, as evinced in red in Table 1 from the HPDI of some of the key hadronic parameters investigated. We have discussed the interplay of NP and hadronic contributions and the dependence of the inferred NP from the assumptions on the charm loop. More conservative hypotheses point to two simple NP scenarios, either a nonvanishing C L Q 2223 , with a ΔI C with respect to the SM of 33 (53) in the data-driven approach (in the one of Ref. [38]), or a nonvanishing C NP 10 , with a ΔI C with respect to the SM of 34 (48) in the data driven case (in the one based on Ref. [38]). The approach based on Ref. [71], instead, favours a more complex scenario with four nonvanishing NP coefficients with a ΔI C of 98, although a ΔI C of 88 can be achieved in the simple scenario of a nonvanishing C N P 9 . Clearly, more data on both LUV observables and differential decay rates is needed to improve our understanding of the charm loop and to single out the correct interpretation of LUV in terms of NP contributions. Hopefully, the LHC [85,86] and Belle II [87] will provide us with the needed precision to identify the NP at the origin of the current evidence of LUV. and C ,NP 10 . Contours and colours as in Fig. 6 are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/. Funded by SCOAP 3 . SCOAP 3 supports the goals of the International Year of Basic Sciences for Sustainable Development.