Mass distributions marginalized over per-event errors

Nuclear Instruments and Methods in Physics Research A 764 (2014) 150–155

Contents lists available at ScienceDirect

Nuclear Instruments and Methods in Physics Research A journal homepage: www.elsevier.com/locate/nima

Mass distributions marginalized over per-event errors D. Martínez Santos a,n, F. Dupertuis b a b

NIKHEF and VU University Amsterdam, Amsterdam, The Netherlands Ecole Polytechnique Fédérale de Lausanne (EPFL), Lausanne, Switzerland

art ic l e i nf o

a b s t r a c t

Article history: Received 16 January 2014 Received in revised form 22 June 2014 Accepted 27 June 2014 Available online 15 July 2014

We present generalizations of the Crystal Ball function to describe mass peaks in which the per-event mass resolution is unknown and marginalized over. The presented probability density functions are tested using a series of toy Monte Carlo samples generated with Pythia and smeared with different amounts of multiple scattering and for different detector resolutions. & 2014 Elsevier B.V. All rights reserved.

Keywords: Statistics Invariant mass peaks Resolution modelling

1. Introduction A very common probability density function (p:d:f :) used to fit the mass peak of a resonance in experimental particle physics is the so-called Crystal Ball (CB) function [1–3]: 8 mμ >  ð1=2Þððm  μÞ=σ Þ2 > if 4 a > otherwise : A B

σ

where m is the free variable (the measured mass), μ is the most probable value (the resonance mass), σ the resolution, a is called the transition point and n the power-law exponent. A and B are calculated by imposing the continuity of the function and its derivative at the transition point a. This function consists of a Gaussian core, that models the detector resolution, with a tail on the left-hand side that parametrizes the effect of photon radiation by the final state particles in the decay. In data analysis, one may deal with events which have different uncertainties on the measured mass, therefore distorting the core of the Crystal Ball function, which will not be a Gaussian any more. This can be very relevant if one wants to estimate signal yields or efficiencies at the percent precision, and it becomes even more important if we study a small signal that sits on top of the resolution tail of another resonance (see for example the analysis of B0s -J=ψ K π [4], where a small B0s signal sits on top of the B0d -J=ψ K π mass tails). The nonGaussian resolution tails are sometimes modelled by the sum of two or three Crystal Ball functions, which is equivalent to assume n Corresponding author at: NIKHEF and VU University Amsterdam, Science Park 105 1098 XG, Amsterdam, The Netherlands. E-mail address: [email protected] (D. Martínez Santos).

http://dx.doi.org/10.1016/j.nima.2014.06.081 0168-9002/& 2014 Elsevier B.V. All rights reserved.

that the per-event uncertainty is a sum of two or three delta functions. However, per-event uncertainties are usually continuous functions very different from a sum of a small number of delta functions. One way of dealing with per-event uncertainties that follow a certain distribution, is to either make a p:d:f : conditional on the per-event uncertainty (if its distribution is known) or perform the analysis in bins of the quantities that affect the perevent uncertainties (for example, particle momenta) and combine them afterwards. However, those procedures can significantly complicate the analysis. For example, if one wants to build a p:d:f : conditional on the per-event uncertainty, one needs to make sure that such uncertainty has been accurately evaluated or, if not, do a very careful calibration. Thus, in some cases one may prefer to simply marginalize over the mass error (for which only the functional form is needed, neither the values of its parameters nor an actual evaluation of the uncertainty event-by-event) and have a p:d:f : that describes the final mass peak, as Z 1 2 1 pffiffiffie  ð1=2vÞðm  μÞ ρðvÞ dv pðmÞ p ð2Þ v 0 where v is the variance and ρðvÞ the prior density of the variance. In this paper, we will define some extensions of the Crystal Ball distribution for different assumptions on ρðvÞ. We will fit the proposed mass models to J=ψ -μ þ μ toy MC samples where we can modify the relative importance of multiple scattering (MS) and detector spatial resolution (hereafter SR). Section 2 describes the generation of the toy MC samples. Section 3 defines an extension of the CB function using a hyperbolic distribution core. Sections 4 and 5 generalize the function defined in Section 3. Section 6 gives a brief discussion of the meaning of the fit parameters. Section 7 discusses other effects on the invariant mass line-shape that are not directly related to resolution. Conclusions are drawn in Section 8.

D. Martínez Santos, F. Dupertuis / Nuclear Instruments and Methods in Physics Research A 764 (2014) 150–155

2. Simulation of J=ψ -μ þ μ decays pffiffi We generate J=ψ events at s ¼ 8 TeV using the main17.cc script of (Pythia8.176) [5]. The J/ψ's are then isotropically decayed into two muons. No photons are added, as the radiative tail of the mass distribution should be well accounted by the Crystal Ball tail. The generated muon momenta are smeared with a Gaussian resolution which has a momentum dependence:

σ ðpÞ p

¼ a þbp

ð3Þ

where a mimics the multiple scattering (MS) and b mimics the effect of the hit resolution. We take as typical values a ¼ 3  10  3 and b ¼ 2  10  5 GeV  1 c inspired by [6], although we will vary them for different tests.

3. Hyperbolic resolution model A very flexible function that describes asymmetric unimodal p:d:f :'s defined above a certain threshold (like per-event error distributions usually look like) is the so-called Amoroso distribution [7,8] (see Fig. 1 for an example). If we consider the Amoroso distribution as a potential implementation for ρðvÞ     v  v0 αβ  1  ððv  v0 Þ=θÞβ ð4Þ ρðvÞ ¼ Amorosoðvv0 ; α; β; θÞ p e  θ then the corresponding core of the invariant mass p:d:f : will be the following:

ΦðmÞ p

Z

1

0

  2 1 v  v0 αβ  1  ððv  v0 Þ=θÞβ pffiffiffie  ð1=2vÞðm  μÞ e dv θ v

ð5Þ

Unfortunately, the above integral cannot be solved analytically. It would require a numerical implementation of the core and its derivative. This would make the matching with the radiative tail difficult. Evaluating Eq. (5) numerically for different values of Amoroso parameters, we find distributions that exhibit a hyperbolic profile when plotted in logarithmic scale. Based on that observation, we define a possible core: pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 cðmÞ p e  b 1 þ ðm  μÞ =δ ð6Þ i.e., c(m) is the symmetric limit of the hyperbolic distribution. It can also be rewritten in such a way that the mass resolution σ appears explicitly, as it will be discussed in Section 4. Adding a CB-like tail to Eq. (6), we obtain the following p:d:f :: Aðm; μ; b; δ; a; nÞ

300 250 200 150 100 50 0 40

50

60

70 2

80

4

v (MeV /c ) Fig. 1. Example of the Amoroso distribution fitted to the distribution of the perevent variance v in one of our toy MC samples.

if

m μ

δ

Z a

otherwise

ð7Þ hereafter referred to as the Apollonios distribution, currently being used for data analysis of the LHCb experiment (see Refs. [9,10]). The core (6) can be obtained analytically for a variance prior density:

ρðv; b; δÞ p e  ðb

2

v=δ þ δ =vÞ 2

2

ð8Þ J=ψ pT

A We fit the mass peak for J=ψ -μ μ decays satisfying 7 μ7 μ7 ½0; 14 GeV=c, θ A ½20; 300 mrad, pμ 46 GeV=c, pT 40:5 GeV=c (which mimics LHCb-like conditions) to the core of the Apollonios distribution,1 and find a very good agreement as can be seen in Fig. 2. Now, the good agreement between this model and the toy MC samples used for testing can be broken without too much effort. For example, we now repeat the exercise releasing all kinematic and acceptance cuts, and switching off the MS term. These changes modify the distribution and fit results are shown in Fig. 3, where we see that Eq. (6) cannot fit the generated data. However, it is also interesting to note that, even in this extreme case, Eq. (6) can do a good job if the fit is performed in a region of about two standard deviations around the peak. þ



4. Generalized hyperbolic resolution model The core of Eq. (6) is a limit case of the generalized hyperbolic distribution [11]: Gðm; μ; λ; α; β; δÞ ¼ ððm  μÞ2 þ δ Þð1=2Þλ  1=4 eβðm  μÞ K λ  1=2 ðα 2

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 ðm μÞ2 þ δ Þ

ð9Þ where K λ are the cylindrical harmonics or special Bessel functions of the second kind [12]. The parameter β is related to a variancedependency of the Gaussian mean, and not to the per event variance distribution. For the purpose of this paper β can be considered zero. In principle, β2 is constrained to be smaller than α2. In practice that condition can be ignored if the fitting range is 2 finite, but one has to be careful that if β 4 α2 one of the tails will start rising at some point. The generalized hyperbolic distribution also has an important limit case, the Student's-t distribution, as indicated in Table 1. The p:d:f : in Eq. (9) can also be obtained by marginalizing over a variance density:

ρðv; λ; α; δÞ p vλ  1 e  ð1=2Þðα

350

Events / ( 0.2 MeV2/c4 )

8 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 > > e  b 1 þ ðm  μÞ =δ > > < !n pffiffiffiffiffiffiffiffiffiffiffiffiffi p pffiffiffiffiffiffiffiffiffi n 1 þ a2 >  b 1 þ a2 > ffiffiffiffiffiffiffiffiffiffiffiffi ffi p e > > : baðn 1 þ a2 =ðbaÞ  a  ðm  μÞ=δÞ

151

2

v þ δ =vÞ 2

ð10Þ

The distribution (10) is the generalized inverse Gaussian distribution and describes very well the density we find for the per-event mass error squared (hereafter σ 2μμ ), for the example in Fig. 3. This is shown in Fig. 4, together with the good agreement between the simulated data and the generalized hyperbolic distribution. We find that Eq. (10) fits well the mass variance distribution for all the generated J=ψ -μ þ μ samples that we have tested, although one needs to add an overall offset to the per-event error, i.e., to change v by v v0 in Eq. (10). The effect of an overall displacement of the per-event error distribution is further discussed in Section 5. The following re-parametrization fα; δg-fσ ; ζ g:

ζ ¼ αδ 1

ð11Þ

In other words, we set a-1, since no final state radiation was simulated.

152

D. Martínez Santos, F. Dupertuis / Nuclear Instruments and Methods in Physics Research A 764 (2014) 150–155

104

7000

Events / ( 0.6 MeV/c2)

Events / ( 0.6 MeV/c2)

8000 6000 5000 4000 3000 2000 1000 0

3050

3100

103 102 10 1 10-1

3150

3050

3100

Mμμ (MeV/c2)

3150

Mμμ (MeV/c2) 7

Fig. 2. Fit of the invariant mass distribution of a J=ψ-μ þ μ generated sample with pT A ½0; 14 GeV/c, θμ A ½20; 300 mrad, pμ 4 6 GeV/c and pμT 40:5 GeV=c. The solid pink line corresponds to the fit to a hyperbolic distribution. The dashed black line corresponds to the fit to a Gaussian. Left: linear scale. Right: Logarithmic scale. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.) 7

J=ψ

3000

Events / ( 0.02 MeV/c2)

Events / ( 0.6 MeV/c2)

105

7

104 103 102 10 1

3050

3100

3150

Mμμ(MeV/c2) þ

2500 2000 1500 1000 500 0 3095

3096

3097

Mμμ

3098

3099

(MeV/c2)



Fig. 3. Fit of the invariant mass distribution of a J=ψ-μ μ generated sample without any phase space restriction and without the multiple scattering term in the momentum resolution. The pink line corresponds to the fit to an hyperbolic distribution. The dashed black line corresponds to the fit to a Gaussian. Left: Fit in the full mass range. Right: Fit in a region of about two standard deviations around the mean. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)

Table 1 Limit cases of the generalized hyperbolic distribution. Distributions

Gðm; μ; λ; α; β; δÞ

Hyperbolic Symmetric hyperbolic Student's t

λ ¼ 1; αδ ¼ b λ ¼ 1; β ¼ 0; αδ ¼ b pffiffiffi ν λ ¼  ; α ¼ 0; β ¼ 0; δ ¼ ν 2 ν λ ¼  ; α ¼ 0; β ¼ 0 2

Non-standardized Student's t

5. Effect of the offset

ζÞ 2 ¼ δ Aλ 2 ðζ Þ σ ¼δ ζ K λ ðζ Þ 2

2K λ þ 1ð

ð12Þ

is more suitable for fitting purposes as it allows us to specify the rms (σ) of the distribution in the symmetric case (β ¼ 0) as an explicit parameter. A2λ ¼ ζ K λ ðζ Þ=K λ þ 1 ðζ Þ is introduced for further convenience. In that parametrization Gðm; μ; σ ; λ; ζ ; βÞ

p ððm  μÞ2

2 ð1=2Þλ  1=4 βðm  μÞ

þ Aλ ðζ Þσ Þ 2

e

hereafter referred to as Hypatia distribution, where G0 is the derivative of the G defined in Eq. (9). The generalized hyperbolic core can describe most of the examples that were generated (see Fig. 4, bottom), but can also be broken with high statistics samples if J=ψ -μ þ μ events are taken all over the phase space without any kinematic or acceptance requirement, as can be seen in Fig. 6.

0 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi1   m  μ 2A K λ  ð1=2Þ @ζ 1 þ A λ ðζ Þσ

ð13Þ

Fig. 5 shows Gðm; μ; σ ; λ; ζ ; βÞ for different values of ζ and λ. Using Eq. (13) as the core of a CB-like function, we define Iðm; μ; σ ; λ; ζ ; β; a; nÞ

8 0 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi1   > m  μ 2A mμ > > if 4 a < ððm  μÞ2 þ A2λ ðζ Þσ 2 Þð1=2Þλ  1=4 eβðm  μÞ K ðλ  ð1=2ÞÞ @ζ 1 þ σ Aλ ðζ Þσ p   n > > Gðμ aσ ; μ; σ ; λ; ζ ; βÞ > : Gðμ  aσ ; μ; σ ; λ; ζ; βÞ= 1 m n 0  aσ Þ otherwise G ðμ aσ ; μ; σ ; λ; ζ ; βÞ

ð14Þ

We have seen that Eq. (10) is a flexible function that can parametrize mass variance distributions if an offset is added to it. Yet, by adding the offset, the marginalization does not yield a generalized hyperbolic distribution for the most general case. We can see that adding an offset to the per-event error distribution is equivalent to performing a convolution: Z 1 2 1 pffiffiffie  ð1=2vÞðm  μÞ ρðv  v0 Þ dv v 0 Z 1 2 1 pffiffiffiffiffiffiffiffiffiffiffiffiffie  ð1=2ðv0 þ ΔÞÞðm  μÞ ρðΔÞ dΔ ¼ v0 þ Δ 0  Z 1 Z þ 1 2 1 2 1 ¼ pffiffiffiffiffie  ð1=2ðv0 ÞÞðm  tÞ pffiffiffiffie  ð1=2ðΔÞÞðt  μÞ dt ρðΔÞ dΔ v0 Δ 0 1 Z 1  Z þ1 2 2 1 1 pffiffiffiffie  ð1=2ðΔÞÞðt  μÞ ρðΔÞ dΔ dt ¼ pffiffiffiffiffie  ð1=2ðv0 ÞÞðm  tÞ v0 Δ 1 0 Z 1  2 2 1 1 pffiffiffiffie  ð1=2ðΔÞÞðt  μÞ ρðΔÞ dΔ ¼ pffiffiffiffiffie  ð1=2ðv0 ÞÞðmÞ n ð15Þ v0 Δ 0

D. Martínez Santos, F. Dupertuis / Nuclear Instruments and Methods in Physics Research A 764 (2014) 150–155

Events / ( 0.05 MeV2/c4)

The convolution of a generalized hyperbolic distribution with a Gaussian is not in general another generalized hyperbolic. However, we can argue that if v0 b Δ, we will have a single Gaussian (that is a limit case of the generalized hyperbolic) and, on the contrary, that if v0 5 Δ in most of the Δ range, we will recover the generalized hyperbolic distribution. One can also argue that as we are looking for corrections to the Gaussian distribution, the

102

2

4

6

8

10

σ2μμ (MeV2/c4) 105

1 v0

Υ ðm; μ; σ SR ; λ; ζ ; β; a; n; v0 Þ ¼ pffiffiffiffiffie  ð1=2ðv0 ÞÞðmÞ nIðm; μ; σ SR ; λ; ζ ; β; a; nÞ 2

ð16Þ

104

Events / ( 1 MeV/c2 )

Events / ( 0.6 MeV/c2)

convolution properties of the Gaussian function still hold approximately. Yet, it will not be exact, and therefore a smeared Hypatia distribution

can provide a better fit than Iðm; μ; σ ; λ; ζ ; β; a; nÞ for some complicated cases with high statistics, without a real increase in the pffiffiffiffiffi number of fit parameters ( v0 can be fixed in a somewhat arbitrary point at the start-up of the mass error distribution), although at the cost of a numerical convolution. The latter can be done in RooFit [13] by calling the RooFFTConvPdfclass on top of the implementation of Iðm; μ; σ ; λ; ζ ; β; a; nÞ. If written this way, pffiffiffiffiffi v0 can be interpreted as an estimate of the mass resolution due to multiple scattering, σ SR the dispersion caused by the spatial resolution of the detector given the kinematics of the final state, qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi and the total resolution would be σ ¼ v0 þ ðσ SR Þ2 . Fig. 6 shows a fit of Iðm; μ; σ ; λ; ζ ; 0; 1; 1Þ and Υ ðm; μ; σ SR ; λ; ζ ; 0; 1; 1; 6:5 MeV=c2 Þ to the simulated J=ψ -μ þ μ data, for the full sample without any kinematic constraint, i.e., where very low momentum (MS dominated) and very high momentum (hit resolution

103

10

153

103 102 10 1

3050

3100

3150

105 104 103 102

Mμμ (MeV/c2)

1 3000

G(m,0,1,λ,0.5,0)

1

λ = -0.5 λ = -1 λ = -10

10-1

10-2

3050

3100

3150

3200

2

Mμμ(MeV/c ) Fig. 6. Fit of Iðm; μ; σ; λ; ζ; 0; 1; 1Þ (red dashed) and Υ ðm; μ; σ SR ; λ; ζ; 0; 1; 1; 6:5 MeV=c2 Þ (solid blue) to the simulated J=ψ-μ þ μ data, for the full sample without any kinematic constraint. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)

G(m,0,1,λ,0.5,0)

Fig. 4. Top: Per-event mass error squared fitted to Eq. (10) in a J=ψ-μ þ μ sample generated without MS. Bottom: Fit to the mass distribution on the same sample. The pink solid line shows the generalized hyperbolic. The dot-dashed blue line the Student's-t case, and the dashed red line the hyperbolic distribution. (For interpretation of the references to color in this figure caption, the reader is referred to the web version of this paper.)

λ = 0.5 λ=1 λ = 10

10-1

10-2

10-3

10-3 -2

0

m

2

4

ζ = 10-2 ζ = 0.5 ζ = 10

1

10-1

10-2

-4

-2

0

m

2

4

-4

-2

0

m

2

1

G(m,0,1,0.5,ζ,0)

-4

G(m,0,1,-0.5,ζ,0)

10

4

ζ = 10-2 ζ = 0.5 ζ = 10

10-1

10-2

10-3

-4

-2

0

2

4

m

Fig. 5. Gðm; μ; σ; λ; ζ; βÞ is plotted for standard values of μ, σ, β and different values of ζ and λ. Top: ζ is fixed to 0.1 and λ is varied. Bottom: λ is fixed to  0.5 for the left plot and 0.5 for the right plot and ζ is varied.

154

D. Martínez Santos, F. Dupertuis / Nuclear Instruments and Methods in Physics Research A 764 (2014) 150–155

dominated) coexist. The fitting range corresponds to about 11 standard deviations. An excellent agreement between Υ ðm; μ; σ SR ; λ; ζ ; 0; 1; 1; 6:5 MeV=c2 Þ and the simulated data is found. We failed to find any subsample of the J=ψ -μ þ μ data that could not be fitted by Υ ðm; μ; σ SR ; λ; ζ ; 0; 1; 1; v0 Þ, using as free parameters μ; σ SR ; λ, and ζ.

6. Properties of

λ

We have seen that using the Hypatia Υ distribution we can factorize the mass resolution modelling into MS and SR. The first pffiffiffiffiffi part is governed by a resolution parameter v0 that can be estimated from the start-up of the per-event variance distribution. The second part is governed by the parameters ζ ; λ; σ of the generalized hyperbolic distribution, where σ corresponds to the resolution introduced by SR and where, empirically, we have found that ζ is in most cases small. In this section we will derive a physical meaning for λ, at least in the small ζ limit. In the α ¼ 0 (-ζ ¼ 0) limit case, the generalized hyperbolic distribution becomes a Student's-t distribution, which can be understood as a marginalization over a per-event variance density:

ρðvÞ p vλ  1 e  b=v

ð17Þ

The mean (M) and mode (μ) of Eq. (17) are M¼

b ; λ1

μ¼

b λþ1

ð18Þ

thus

λ¼

1 þMðvÞ=μðvÞ o0 1 MðvÞ=μðvÞ

ð19Þ

and we can get an estimate of λ by looking at the per-event error (squared) distribution, and making the ratio of its mean and mode after shifting it to start at zero. But, we can further exploit this relation. From Eq. (3) we can suppose that the per-event uncertainty will be strongly correlated with the particle momenta. Indeed, Fig. 7 supports this.

σ SR i  cte  pJ=ψ ;i

ð20Þ

If this is the case, then MðvSR Þ=μðvSR Þ  Mðp2J=ψ Þ=μðp2J=ψ Þ and λ does not depend on detector effects, only on particle kinematics. This is an interesting result, because if we have a MC simulation with a good description of the momentum distribution of the particles in the lab frame, then the values of λ obtained in simulation should be reasonably valid for data, regardless of having an accurate description of detector simulation. A numeric test of this approach can be seen in Table 2.

7. Mass constraints on intermediate resonances Up to now we have described resolution effects. In more complicated cases, it is sometimes very useful to apply constraints on the decay products. For example, one can significantly improve the invariant mass resolution of B0s -J=ψϕ by constraining the two muons to have the PDG J=ψ mass [14]. This kind of approach, although great at improving the overall resolution, can also generate tails on the mass distribution, due to the photon energy radiated in the J=ψ -μ þ μ decay. Let us consider a simple case in which the constraint is just applied by substituting the mass of the dimuon by the mass of the J=ψ . qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiqffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  m2c ¼ m2J=ψ þ m2KK þ 2 m2J=ψ þ p2μμ m2KK þ p2KK  pμμ pKK cos ðθÞ ð21Þ while ideally one would have wanted to implement: qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiqffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  m2true ¼ m2J=ψ þ m2KK þ 2 m2J=ψ þ p2J=ψ m2KK þ p2KK  pJ=ψ pKK cos ðθÞ ð22Þ m2c  m2true

The difference is not zero but rather a function of the energy of the photons generated in the J=ψ -μ þ μ decay. This difference can be greater than zero, generating a tail on the righthand side. Hence, even with a perfect detector resolution, the combination of the mass constraint and the photon radiation will generate non-Gaussian tails. In practice, this effect is expected to be small because the J=ψ -μ þ μ decays are selected with a mass window cut that allows only low energy photons. Otherwise, it can be partially accommodated either by the resolution model (e.g. Eq. (9)) or by using a CB-like tail on the right-hand side (i.e., using a double-sided Hypatia). A further discussion of such effects goes beyond the scope of this paper, which is to provide models marginalized over per-event errors. Table 2 Results of a fit to Υ ðm; μ; σ SR ; λ; 0; 0; 1; 1; σ MS Þ for toy MC J=ψ events smeared with different values of a and b in (3). The parameter σ MS is fixed at the start-up of the per-event variance. The parameter λ is found to be very stable with respect to smearing parameters, which are varied by 100%. However the uncertainty on λ varies significantly, and increases with a=b. a ½10  3 

b ð GeV=c  1 Þ

λ

σ MS ð MeV=c2 Þ

σ SR ð MeV=c2 Þ

3 1.5 6 3 3

2  10  5 2  10  5 2  10  5 4  10  5 1  10  5

 2.40 7 0.06  2.107 0.03  2.677 0.16  2.11 70.03 2.65 70.15

6.81 3.53 13.3 7.07 6.67

4.75 7 0.02 3.71 70.01 6.3 7 0.05 7.53 70.03 3.06 70.03

Fig. 7. Per-event mass uncertainty versus J=ψ momentum. Left: with multiple scattering. Right: only detector resolution.

D. Martínez Santos, F. Dupertuis / Nuclear Instruments and Methods in Physics Research A 764 (2014) 150–155

8. Conclusions We have presented a generalization of the Crystal Ball function that gives an excellent description of mass resolution non-Gaussian tails. This function, that we name the Hypatia distribution, I, corresponds to a CB-like tail with a generalized hyperbolic core. The smeared Hypatia distribution, Υ provides an improved description of mass peaks and its fit parameters have clearer fundamental meaning, although the price to pay is a numeric convolution. A second, right-hand side CB-like tail can be added in cases where one has other non-resolution effects, such as those coming from constraining the mass of intermediate resonances of a decay. Acknowledgments We would like to thank L. Carson and V. Gligorov for their useful comments on our paper draft. We would also like to thank W. Hulsbergen for helpful discussions during the preparation of this work. References [1] M. Oreglia, A study of the reactions ψ 0 -γγψ (Ph.D. thesis), 1980. [2] J. Gaiser, Charmonium spectroscopy from radiative decays of the J=ψ and ψ 0 (Ph.D. thesis), 1982.

155

[3] T. Skwarnicki, A study of the radiative cascade transitions between the Upsilon-prime and Upsilon resonances (Ph.D. thesis), Institute of Nuclear Physics, Krakow, DESY-F31-86-02, 1986. [4] R. Aaij, et al., Physical Review D 86 (2012) 071102, http://dx.doi.org/10.1103/ PhysRevD.86.071102 arXiv:1208.0738. [5] T. Sjöstrand, S. Mrenna, P. Skands, Journal of High Energy Physics 05 (2006) 026, http://dx.doi.org/10.1088/1126-6708/2006/05/026 arXiv:hep-ph/0603175. [6] A.A. Alves Jr., et al., Journal of Instrumentation 3 (2008) S08005, http://dx.doi. org/10.1088/1748-0221/3/08/S08005. [7] L. Amoroso, Annali Di Matematica Pura Ed Applicata 21 (1925) 123. [8] G.E. Crooks, The Amoroso Distribution arXiv:1005.3274. [9] R. Aaij, et al., Measurements of the B þ , B 0, B0s meson and Λ0b baryon lifetimes arXiv:1402.2554. [10] R. Aaij, et al., Study of beauty hadron decays into pairs of charm hadrons arXiv:1403.3606. [11] A.J. McNeil, R. Frey, P. Embrechts, Quantitative Risk Management Concepts, Techniques and Tools, Princeton University Press. ISBN: 9780691122557. [12] Mathworld, 〈http://mathworld.wolfram.com/ModifiedBesselFunctionoftheSe condKind.html〉. [13] Roofit, 〈http://root.cern.ch/drupal/content/roofit〉. [14] J. Beringer, et al., Physical Review D 86 (2012) 010001, http://dx.doi.org/ 10.1103/PhysRevD.86.010001.