An iterative method to estimate the combinatorial background
aa r X i v : . [ phy s i c s . d a t a - a n ] A ug EPJ manuscript No. (will be inserted by the editor)
An iterative method to estimate the combinatorial background
Georgy Kornakov ,Tetyana Galatyuk , Technische Universit¨at Darmstadt GSI Helmholtzzentrum f¨ur SchwerionenforschungAugust 2018
Abstract.
The reconstruction of broad resonances is important for understanding the dynamics of heavyion collisions. However, large combinatorial background makes this objective very challenging. In this workan innovative iterative method which identifies signal and background contributions without input modelsfor normalization constants is presented. This technique is successfully validated on a simulated thermalcocktail of resonances. This demonstrates that the iterative procedure is a powerful tool to reconstructmulti-differentially inclusive resonant signals in high multiplicity events as produced in heavy ion collisions.
PACS.
29. Experimental methods and instrumentation for elementary-particle and nuclear physics –29.85.-c Computer data analysis
In heavy ion collisions, one of the observables that is com-monly used to characterize the properties of producedmatter is the two-body invariant mass spectrum. The short-lived (unstable) particles appear as a bump with a widthrelated to the lifetime of the particle, ( Γ ∼ /τ ). In thecase of excited baryon states like ∆ , N ∗ or the ρ meson,their broad width and a poor signal-to-background ratiomake the reconstruction difficult. The main challenge isthe presence of a large combinatorial background due torandom combinations of measured particles.The precise measurement of the differential inclusivespectra of short lived baryons and mesons would allow toimprove the understanding of the mechanisms of parti-cle production and to test available theoretical models byproviding precise inclusive spectra.Several approaches have been proposed to solve thisproblem, among them can be named techniques like Monte-Carlo phase-space generators (integrators) [1], a commontool for inclusive analysis in elementary reactions. Thelike-sign method, which is commonly used for dilepton[2,3,4,5] and dipion analysis at high energies [6], can beproven to be exact for certain cases, but it cannot be ap-plied universally. For many-particle events, methods basedon position swapping have been proposed with good re-sults for the diphoton invariant mass reconstruction [7].However, the most common technique is the Event-Mixingmethod [8], introduced originally to generate the refer-ence distribution to study the identical particle correla-tions. The main idea is to combine particles from different Send offprint requests to : a Present address:
[email protected] collisions, generating an uncorrelated spectrum that afternormalization is subtracted from the total same-event dis-tribution, yielding the signal spectrum as their difference.Even though it has been pointed out that the backgroundobtained with this technique has several problems due tothe determination of normalization constants, the parti-cle and event class selection or acceptance [9,10,11], it isthe most commonly applied method [12,13,14,15,16,17,18,19].In general, decay products are used to generate thebackground as well as uncorrelated particles, having com-monly different kinematics. Consequently, contributionsfrom signal particles might distort the generated back-ground. This is a major drawback, since the uncorrelateddistribution will be usually different from the true com-binatorial background [10]. One possible solution to over-come this is to generate the contribution from signal pairsto the combinatorial background by modelling the yield,spectral shape and decay kinematics of the studied reso-nance. However, it introduces a strong model dependencemaking the study of the modification of the line shapes im-possible. In this paper, we describe a new unbiased methodto estimate the correlated distribution directly from data.It is worthwhile to mention that a similar approach hasbeen developed in gamma-spectroscopy in order to sup-press contributions from Compton and backward scatter-ing, where the sought signal is stripped or unfolded fromthe total measured gamma spectrum in a similar iterativeprocedure [20].This paper is organized as follows: in Sec. 2 the iter-ative method is presented. In Sec. 3 a simulated cocktailof resonances and thermal sources is used to demonstratethe capabilities of the developed technique. In Sec. 4 themethod is discussed regarding its capabilities and limita-
Georgy Kornakov,Tetyana Galatyuk: An iterative method to estimate the combinatorial background tions. Also, it is put into context of other techniques devel-oped for background subtraction in other fields. Finally,the work is summarized.
Let T be the total discrete multi-differential spectrum ob-tained after adding all possible same-event four-momentacombinations of two particles of type 1 and 2 over allevents.The signal, S , is comprised of those combinations ofparticles that are daughters of same mother, R → P + P .The remaining combinations belong to the combinatorialbackground B and they relate to each other as T = B + S . (1)The goal of the iterative method is to transform T such that only signal remains at the end. For that pur-pose, uncorrelating operators acting only on signal canbe used. Then, let ˆ U be a non-linear uncorrelating opera-tor, for example, position swapping, event mixing or trackrotation. The latter technique consists in modifying theorientation of single particle three-momentum preservingits module. It is widely used to calculate the combinato-rial background beneath narrow peaks [22,23,24,25] andis the operator of choice in this work.In order to illustrate these definitions a simulated setof three uncorrelated ρ (770) per event decaying into pionpairs has been produced. The proportion of signal pairsto random combinations is 3:6. The invariant mass projec-tions for total, signal, background and the total rotatedspectra are shown in Fig. 1. The invariant mass projec-tions of the difference between T and ˆ U ( T ) is shown inthe inset of the same figure. As pairs from background areinsensitive to rotations, these differences are only due tosignal. After rotating one track from the pair, the signalspreads over a broader phase space region. Hence, posi-tive values are expected in the peak and negative valueson the sides. Provided a model exists for the signal, it ispossible to extract directly the correct signal yield, whose S − ˆ U ( S ), saturates the T − ˆ U ( T ) distribution. Althoughthe reconstruction would be biased [10].Alternatively, the signal can be approached by an iter-ative solution. The positive difference areas of T − ˆ U ( T )are due to the signal. The first iterative solution, S , corre-sponds to these regions. Then, the uncorrelator operatoris applied upon the first solution and the improved ap-proximate solution S = T − ( ˆ U ( T ) − ˆ U ( S )) is obtained.This process repeats until complete saturation of the totaldifference distribution.Mathematically, the proposed algorithm can be con-sidered to be similar to the Landweber Iterative Method[26], widely used to solve ill-conditioned, noisy and non-linear systems [27,28]. From Eq. 1, B = ˆ U ( B ) and the factthat the positive areas of T − ˆ U ( T ) belong to signal, theiterative solution for the S matrix can be written as S k +1 = max {T − ( ˆ U ( T ) − ˆ U ( S k )) , } , (2) ] [MeV/c - π + π inv M
500 1000 1500 a . u . r o t a t ed s i gna l background ] [MeV/c - π + π inv M500 1000 a . u . total - rotated Fig. 1.
Example of an invariant mass distribution (full circles), T , which is comprised of signal S (dashed line) and combina-torial background B (grey area). The solid line is the uncorre-lated ˆ U ( T ) distribution obtained with the rotation technique(see for details in text). The inset shows the invariant mass ofthe T − ˆ U ( T ) distribution. Notice the positive and negativeareas coinciding with the peak location and its sides. where only positive values are chosen and k is the iterationindex. Then, convergence is reached once S k +1 = S k . (3)This is exact for infinite statistics. Otherwise, the bin-to-bin statistical fluctuations introduce a systematic incre-ment of the signal with every iteration. This incrementdepends on the number of pairs contributing to the totalmatrix T and those to ˆ U ( B ). As a consequence, S k willgrow continuously and eventually will reach T . A simpleexample can be a data set without signal in it. In such acase, after applying the uncorrelating operator and sub-tracting the produced distribution from the total approx-imately half of the bins will have a positive value. Thisphony signal will increment further with each iterationthe entries in S k . The stopping criteria for the iteratorin case of finite statistics is given by the change in thesignal increment rate as it moves from the real signal do-main towards the statistical increment regime. There is nocommon stop iteration for all the bins together. Individualvalues have to be extracted by identifying the transitionfrom increase due to real signal to increase due to statisti-cal fluctuations. In the abovementioned example withoutsignal, the increment rate will be constant and only de-pending on the number of entries in each bin. eorgy Kornakov,Tetyana Galatyuk: An iterative method to estimate the combinatorial background 3 ] [MeV/c p - π M c oun t s t o t a l r o t a t ed s i gna l s i gna l Fig. 2.
Signal generated with Pluto [21] imitating the environ-ment of a heavy ion collision. The solid curve is the total invari-ant mass distribution T , the short-dashed line is the resonantsignal S , and the long-dashed line is the uncorrelated resonantsignal distribution ˆ U ( S ). The signal-to-background ratio is be-low 5 %. The ˆ U ( T ) invariant mass distribution is not shown asit would be indistinguishable in this representation and equalto T . Table 1.
Cocktail of particles generated by Pluto [21] to im-itate the environment of a heavy ion collision at intermediateenergies. A thermal source of 65 MeV was considered.N/Nprot Channels
Γ/Γ i ∆ (1232) π − p 1 N (1440) π − p, p ρ − , ∆ + π − N (1535) π − p, p ρ − , ∆ + π − ∆ (1620) π − p, p ρ − , ∆ + π − ∆ (1920) π − p, p ρ − , ∆ + π − π −
10 %
In order to present the complete picture, a simulationconsisting of pions and protons from a thermal source to-gether with a cocktail of resonances decaying into a protonand a negative pion in the final state was produced withPluto [21], mimicking particle production in a heavy ioncollision at 1.25 GeV per nucleon. The resonances weregenerated assuming a thermal source of temperature of65 MeV, their relative abundance, decay channels andbranching ratios are summarized in Table 1. The invariantmass distribution of T , S and ˆ U ( S ) are shown in Fig. 2,having a signal-to-background ratio below 5%.The T − ˆ U ( T ), S − ˆ U ( S ) and B − ˆ U ( B ) difference dis-tributions of the invariant mass of π − p pairs are shown c oun t s − × total - total rotatedsignal - signal rotated ] [MeV/c invp - π M c oun t s −
010 10 × × background - background rotated × Fig. 3.
Top panel, invariant mass distribution
T − ˆ U ( T ) (blackcircles) and S − ˆ U ( S ) distribution (red dots) obtained from thecocktail simulation. The non-zero difference is only due to thesignal contribution. Bottom panel, invariant mass distributionof B − ˆ U ( B ), which is compatible with zero within statisticaluncertainties. in Fig. 3. It is evident from this figure that the differencein the shape is only due to the signal and that the back-ground does not contribute to the total difference. Fig. 4shows the difference distributions of the opening angle α and pair transverse momentum P T . The non-uniformitypoints to the necessity of a differential analysis. In previ-ous works [11], it has already been observed that differ-ent particle or phase-space selections were modifying theshape of the reconstructed background and it has to beaccounted for during the reconstruction. Therefore, the re-construction of the signal is performed in four dimensions:invariant mass, pair transverse momentum, pair rapidityand opening angle. The isotropy of the azimuthal emissionangle and the decay plane are assumed. Otherwise, theseextra two dimensions have to be included as well in thereconstruction.The starting minimization hypothesis is absence of sig-nal, S = 0. The available phase-space is binned as follows:40 in M π − p , 5 in P T , 5 in α and 5 in Y . Then, following theEq. 2, the first solution S is given by max {T − ˆ U ( T ) , } .The next iteration includes the contribution of the signalobtained in the first step after rotating the pairs from it.For the simulated data sample 300 iterations have beenperformed. The increment in the total yield of the sig-nal is shown in Figure 5. The largest increase happens inthe first iterations as expected. In this particular case, anankle around the 80th iteration indicates that the statis-tical increment regime was reached and no more resonantparticle pairs can be incorporated into the reconstructedsignal. Georgy Kornakov,Tetyana Galatyuk: An iterative method to estimate the combinatorial background [MeV/c] p - π T P c oun t s − × T p pair P - π total - total rotated [rad] o.a.p - π θ c oun t s − × -p - π Opening angle total - total rotated
Fig. 4.
Projections of
T − ˆ U ( T ) from cocktail simulationto the transverse momentum (top) and opening angle (bot-tom). These projections point towards the necessity of multi-differential analysis in order to account for differences in thesignal and background in different phase-space bins. iteration k - S ∫ - k S ∫ Fig. 5.
Integrated increment of the signal yield in every iter-ation step. The change in slope indicates a transition betweenthe regime with real signal growth and the rise of signal dueto bin-to-bin statistical fluctuations. In this example the tran-sition happens between the 40th and the 80th iteration.
However, as it was introduced in Sec. 2, it is not enoughto stop the iterator at this value. Due to different signal tobackground ratios over the available phase-space, the yieldevolution has to be studied independently for each bin.Three regions have to be identified, the first with a fastrising, then an increase that linearly continues until theankle found in the total increment evolution and finallythe continuous increase, as shown in Fig. 6. An algorithm iteration c oun t s − Fig. 6.
Signal yield evolution for each iteration step for a given M and P T bin, integrated over the whole rapidity and pairopening angle. The true signal, in this phase-space bin, growsonly until the 28th iteration whereas is incompatible with alinear trend. ] [MeV/c p - π M E n t r i e s
365 MeV/c ≤ T P ] [MeV/c p - π M E n t r i e s
545 MeV/c ≤ T
365 < P ] [MeV/c p - π M E n t r i e s
720 MeV/c ≤ T
545 < P ] [MeV/c p - π M E n t r i e s
940 MeV/c ≤ T
720 < P ] [MeV/c p - π M E n t r i e s > 940 MeV/c T P reconstructed signalinput signal
Fig. 7.
Differential π − p invariant mass distribution in 5 P T bins after bin-by-bin convergence. The input resonance signalfrom Pluto is shown in red. finds the iteration at which the increase becomes compat-ible with linear in the second region. After this procedureis applied for all bins, one more iteration should be per-formed to ensure that the condition T − ˆ U ( T ) = S − ˆ U ( S )is fulfilled.Finally the differential and integrated solutions are ob-tained. In this example, with a signal-to-background ratiosmaller than 5%, the reconstructed and the input differ-ential invariant mass are shown in Fig. 7. The integratedsignal can be resolved with an accuracy better than 10%,within statistical uncertainties, as shown in Fig. 8. eorgy Kornakov,Tetyana Galatyuk: An iterative method to estimate the combinatorial background 5 - ] E n t r i e s [ M e V / c S i gna l / B a ck g r ound ] [MeV/c p - π M R /I Fig. 8.
Final integrated reconstructed invariant mass distribu-tion (open circles) after applying the iterative procedure com-pared to the input distribution (red line). The signal to back-ground ratio (green open circles). The ratio between input andreconstructed signals (R/I), shows agreement within statisticaluncertainties, which are of the order of 10%.
The inclusive reconstruction of short lived states ( τ ≈ − s) is technically challenging due to the existence of alarge combinatorial background. In this work the develop-ment of a model-independent and normalization-free iter-ative method has been addressed, in order to overcome thedifficulties that other techniques have [10]. This methodwill facilitate a precise multi-differential identification ofsignal and background contributions to the reconstructedspectra.The difference of the total and the generated uncor-related distributions from the same data sample preserveinformation about the resonant signal. It was found thatthe positive areas belong to the signal distribution. Thesolution is obtained in an iterative approach, similar to theLandweber Iterative Method used to solve ill-conditioned,noisy and non-linear systems [26,27].As it happened with the Landweber Iterative Method,the first solutions were slow and inefficient. There aremany possibilities to be explored in order to improve theminimization: improved convergence criteria to stop ear-lier the number of iterations, the pre-conditioning of theinitial matrices, and the accelerated or boosted versionswith larger gradients.Most of the effects of a real measurement have notbeen considered such as finite resolution, wrong parti-cle identification, efficiency, acceptance nor occupancy. Anon-uniform probability distribution of the rotation an-gle can account for them. The needed corrections can bederived from data.The algorithm has been validated successfully on asimulated cocktail of thermal resonances. The multi-di- fferential signal has been reconstructed, with a signal-to-background ratio below 5%, in all mass ranges with anaccuracy better that 10% within statistical uncertainties.The method capability is limited by the statistical fluc-tuations of the binned background. The more the broadsignal approaches in shape the combinatorial backgroundthe larger becomes the number of needed iterations. Thisdemands to perform a feasibility pre-study of the differentdistributions before applying the technique. Even contin-uum signals like c ¯ c or QGP signaltures like q ¯ q → e + e − should be accessible.The application of this technique would allow to re-construct differentially inclusive spectra in both elemen-tary and heavy ion collisions for the study of resonanceproduction and their contribution to the yields and differ-ential spectra of stable particles measured by detectors.The reconstruction of π ± p and π + π − channels measuredin Au+Au collisions at √ s = 2 .
42 GeV in HADES is on-going.The authors wish to acknowledge stimulating discus-sions with M. Gumberidze, S. Harabasz, R. Holzmann,C. M¨untz, V. Pechenov, O. Pechenova, B. Ramstein, P.Salabura, H. Str¨obele and J. Stroth. This work has beensupported by VH-NG-823 and Helmholtz Alliance HA216-EMMI.
References
1. F. James, CERN-68-15.2. S. Constantinescu, S. Dita, and D. Jouan, 1996. Report PN-ODER9601.3. M. Gazdzicki and M. I. Gorenstein, hep-ph/0003319.4. M. C. Abreu et al. [NA38 and NA50 Collaborations], Eur.Phys. J. C (2000) 443.5. A. Adare et al. [PHENIX Collaboration], Phys. Rev. C (2010) 0349116. J. Adams et al. [STAR Collaboration], Phys. Rev. Lett. (2004) 0923017. N. van Eijndhoven and W. Wetzels, Nucl. Instrum. Meth.A (2002) 5138. G. I. Kopylov, Phys. Lett. (1974) 472.9. P. D. Higgins et al. , Phys. Rev. D (1979) 65.10. D. L’Hˆote, Nucl. Instrum. Meth. A (1994) 544.11. M. Kaskulov, E. Hern´andez and E. Oset, Eur. Phys. J. A (2010) 22312. G. Jancso et al. , Nucl. Phys. B (1977) 1.13. A. Breakstone et al. , Z. Phys. C (1984) 321.14. D. Drijard, H. G. Fischer and T. Nakada, Nucl. Instrum.Meth. A (1984) 367.15. M. Trzaska et al. , Z. Phys. A (1991) 325.16. M. Eskef and the FOPI collaboration, Eur. Phys. Journ.A (1998) 33517. P. Crochet and P. Braun-Munzinger, Nucl. Instrum. Meth.A (2002) 56418. R. Arnaldi et al. [NA60 Collaboration], Phys. Rev. Lett. (2006) 16230219. B. I. Abelev et al. [STAR Collaboration], Phys. Rev. C (2009) 06490320. D.C. Radford, I. Ahmad, R. Holzmann, R.V.F. Janssensand T.L. Khoo, Nucl. Instrum. Meth. A (1987) 111. Georgy Kornakov,Tetyana Galatyuk: An iterative method to estimate the combinatorial background21. I. Fr¨ohlich et al. , PoS ACAT (2007) 07622. J. Adams et al. [STAR Collaboration], Phys. Rev. Lett. (2005) 12230123. B. I. Abelev et al. [STAR Collaboration], Science (2010) 5824. D. Adamova et al. [CERES Collaboration], Nucl. Phys. A (2012) 4125. B. Abelev et al. [ALICE Collaboration], Phys. Lett. B (2012) 16526. L. Landweber, American Journal of Mathematics, 73(3),615-624 (1951). doi:10.2307/237231327. O. Scherzer, Journal of Mathematical Analysis and Appli-cations (1995) 911.28. Q. Jin, Journal of Mathematical Analysis and Applications253