|
[Rivet-svn] r4155 - in trunk: . data/anainfo data/plotinfo src/Analysesblackhole at projects.hepforge.org blackhole at projects.hepforge.orgThu Feb 14 13:39:12 GMT 2013
Author: hoeth Date: Thu Feb 14 13:39:12 2013 New Revision: 4155 Log: merge c4069 and c4070 from branches/2012-06-aidarivet Modified: trunk/ChangeLog trunk/data/anainfo/ALICE_2012_I1181770.info trunk/data/plotinfo/ALICE_2012_I1181770.plot trunk/src/Analyses/ALICE_2012_I1181770.cc Modified: trunk/ChangeLog ============================================================================== --- trunk/ChangeLog Thu Feb 14 13:36:46 2013 (r4154) +++ trunk/ChangeLog Thu Feb 14 13:39:12 2013 (r4155) @@ -24,6 +24,8 @@ * Version number bump to 1.8.2 -- release approaching. + * Rewrite of ALICE_2012_I1181770 analysis to make it a bit more sane and acceptable. + * Adding a note on FourVector and FourMomentum that operator- and operator-= invert both the space and time components: use of -= can result in a vector with negative energy. Modified: trunk/data/anainfo/ALICE_2012_I1181770.info ============================================================================== --- trunk/data/anainfo/ALICE_2012_I1181770.info Thu Feb 14 13:36:46 2013 (r4154) +++ trunk/data/anainfo/ALICE_2012_I1181770.info Thu Feb 14 13:39:12 2013 (r4155) @@ -10,18 +10,21 @@ - Martin Poghosyan <Martin.Poghosyan at cern.ch> - Sercan Sen <Sercan.Sen at cern.ch> - Burak Bilki <bbilki at gmail.com> + - Andy Buckley <andy.buckley at cern.ch> References: - arXiv:1208.4968 [hep-ex] RunInfo: - Inelastic events (non-diffractive and inelastic diffractive). + Inelastic events (non-diffractive and inelastic diffractive). NumEvents: 100K Beams: [p+, p+] Energies: [900, 2760, 7000] Description: - 'Measurements of cross sections of inelastic and diffractive processes in proton-proton collisions at $\sqrt{s}$=900, 2760 and 7000 GeV. - The fractions of diffractive processes in inelastic collisions were determined from a study of gaps in charged particle pseudorapidity distributions. - Single-diffractive events are selected with $M_{X} < 200$ GeV $/c^{2}$ and double-diffractive events defined as NSD events with $\Delta\eta >3$. - To measure the inelastic cross section, beam properties were determined with van der Meer scans using a simulation of diffraction adjusted to data.' + 'Measurements of cross-sections of inelastic and diffractive processes in proton-proton collisions at $\sqrt{s} = 900$, 2760 and 7000 GeV. + The fractions of diffractive processes in inelastic collisions were determined from a study of gaps in charged particle pseudorapidity distributions. + Single-diffractive events are selected with $M_{X} < 200\;\GeV/c^2$ and double-diffractive events defined as NSD events with $\Delta\eta > 3$. + To measure the inelastic cross-section, beam properties were determined with van der Meer scans using a simulation of diffraction adjusted to data. + + Note that these are experimental approximations to theoretical concepts -- it is not totally clear whether the data point values are model-independent.' BibKey: :2012sja BibTeX: '@article{:2012sja, author = "Abelev, Betty and others", Modified: trunk/data/plotinfo/ALICE_2012_I1181770.plot ============================================================================== --- trunk/data/plotinfo/ALICE_2012_I1181770.plot Thu Feb 14 13:36:46 2013 (r4154) +++ trunk/data/plotinfo/ALICE_2012_I1181770.plot Thu Feb 14 13:39:12 2013 (r4155) @@ -3,7 +3,7 @@ # END PLOT # BEGIN PLOT /ALICE_2012_I1181770/d01-x01-y0 -Title=Production ratios of SD with $M_{X} < 200$ GeV$/c^{2}$ to INEL +Title=Production ratios of SD with $M_{X} < 200\;\GeV/c^2$ to INEL YLabel=$\sigma_\text{SD} / \sigma_\text{inel}$ # END PLOT @@ -13,16 +13,16 @@ # END PLOT # BEGIN PLOT /ALICE_2012_I1181770/d03-x01-y0 -Title=Single diffraction cross-section for $M_{X} < 200$ GeV$/c^{2}$ +Title=Single diffraction cross-section for $M_{X} < 200\;\GeV/c^2$ YLabel=$\sigma_\text{SD}$ [mb] # END PLOT # BEGIN PLOT /ALICE_2012_I1181770/d04-x01-y0 -Title=Double diffraction cross-section for $\Delta\eta >3$ +Title=Double diffraction cross-section for $\Delta\eta > 3$ YLabel=$\sigma_\text{DD}$ [mb] # END PLOT # BEGIN PLOT /ALICE_2012_I1181770/d05-x01-y0 -Title=Inleastic cross-section +Title=Inelastic cross-section YLabel=$\sigma_\text{inel}$ [mb] # END PLOT Modified: trunk/src/Analyses/ALICE_2012_I1181770.cc ============================================================================== --- trunk/src/Analyses/ALICE_2012_I1181770.cc Thu Feb 14 13:36:46 2013 (r4154) +++ trunk/src/Analyses/ALICE_2012_I1181770.cc Thu Feb 14 13:39:12 2013 (r4155) @@ -13,145 +13,82 @@ : Analysis("ALICE_2012_I1181770") { } - public: void init() { - addProjection(ChargedFinalState(),"CFS"); + // Projection setup + addProjection(ChargedFinalState(), "CFS"); - if (fuzzyEquals(sqrtS()/GeV, 900, 1E-3)) { - _h_frac_sd_inel = bookScatter2D(1, 1, 1); - _h_frac_dd_inel = bookScatter2D(2, 1, 1); - _h_xsec_sd = bookHisto1D (3, 1, 1); - _h_xsec_dd = bookHisto1D (4, 1, 1); - _h_xsec_inel = bookHisto1D (5, 1, 1); - } else if (fuzzyEquals(sqrtS()/GeV, 2760, 1E-3)) { - _h_frac_sd_inel = bookScatter2D(1, 1, 2); - _h_frac_dd_inel = bookScatter2D(2, 1, 2); - _h_xsec_sd = bookHisto1D (3, 1, 2); - _h_xsec_dd = bookHisto1D (4, 1, 2); - _h_xsec_inel = bookHisto1D (5, 1, 2); - } else if (fuzzyEquals(sqrtS()/GeV, 7000, 1E-3)) { - _h_frac_sd_inel = bookScatter2D(1, 1, 3); - _h_frac_dd_inel = bookScatter2D(2, 1, 3); - _h_xsec_sd = bookHisto1D (3, 1, 3); - _h_xsec_dd = bookHisto1D (4, 1, 3); - _h_xsec_inel = bookHisto1D (5, 1, 3); - } + // Book (energy-specific) histograms + int isqrts = -1; + if (fuzzyEquals(sqrtS()/GeV, 900, 1E-3)) isqrts = 1; + else if (fuzzyEquals(sqrtS()/GeV, 2760, 1E-3)) isqrts = 2; + else if (fuzzyEquals(sqrtS()/GeV, 7000, 1E-3)) isqrts = 3; + assert(isqrts > 0); + + _h_frac_sd_inel = bookScatter2D(1, 1, isqrts); + _h_frac_dd_inel = bookScatter2D(2, 1, isqrts); + _h_xsec_sd = bookHisto1D (3, 1, isqrts); + _h_xsec_dd = bookHisto1D (4, 1, isqrts); + _h_xsec_inel = bookHisto1D (5, 1, isqrts); } - void analyze(const Event& event) { - const double weight = event.weight(); + void analyze(const Event& event) { const ChargedFinalState& cfs = applyProjection<ChargedFinalState>(event, "CFS"); + if (cfs.size() < 2) vetoEvent; // need at least two particles to calculate gaps - // fill INEL plots for each event - if (fuzzyEquals(sqrtS()/GeV, 900, 1E-3)) { - _h_xsec_inel->fill( 900/GeV, weight ); - } else if (fuzzyEquals(sqrtS()/GeV, 2760, 1E-3)) { - _h_xsec_inel->fill( 2760/GeV, weight ); - } else if (fuzzyEquals(sqrtS()/GeV, 7000, 1E-3)) { - _h_xsec_inel->fill( 7000/GeV, weight ); - } - - double Eslowest(0.0), pslowest(0.0), yslowest(999.); - double Efastest(0.0), pfastest(0.0), yfastest(-999.); - int pidslowest(0), pidfastest(0); - - double LRG(0.0), etapre(0.0), gapbwd(0.0), gapfwd(0.0); - unsigned int num_p(0); - - FourMomentum leadP(0.,0.,0.,0.); - double Elead(0.0), plead(0.0); - - foreach(const Particle& p, cfs.particlesByEta()) { //sorted from minus to plus - const PdgId pid = p.pdgId(); - double y = p.momentum().rapidity(); - double eta = p.momentum().eta(); - //SD case - if (y < yslowest) { - Eslowest = p.momentum().E(); - pslowest = p.momentum().vector3().mod(); - yslowest = p.momentum().rapidity(); - pidslowest = pid; - } - if (y > yfastest) { - Efastest = p.momentum().E(); - pfastest = p.momentum().vector3().mod(); - yfastest = p.momentum().rapidity(); - pidfastest = pid; - } - - num_p += 1; - // DD case - if (num_p==1) { - etapre = p.momentum().eta(); - } else if (num_p > 1) { - if (num_p==2) gapbwd = fabs(eta-etapre); - - double gap = fabs(eta-etapre); - LRG = (gap > LRG ? gap : LRG); // largest gap - - if (num_p==cfs.size()) gapfwd = fabs(eta-etapre); - etapre = eta; - } - } + const double weight = event.weight(); - // Mx calculation - if (pidslowest==2212 && pidfastest==2212) { - if (fabs(yslowest) > fabs(yfastest)) { - Elead = Eslowest; - plead = pslowest; - } else if (fabs(yslowest) < fabs(yfastest)) { - Elead = Efastest; - plead = pfastest; - } else { - Elead = Eslowest; - plead = pslowest; - if ( (double)rand() / (double)RAND_MAX > 0.5) { // generate random number in [0.,1.] range and make decision randomly - Elead = Efastest; - plead = pfastest; - } - } - } else if (pidslowest==2212) { - Elead = Eslowest; - plead = pslowest; - } else if (pidfastest==2212) { - Elead = Efastest; - plead = pfastest; - } + // Fill INEL plots for each event + _h_xsec_inel->fill(sqrtS()/GeV, weight); - double Mx = sqrt((sqrtS()/GeV-Elead-plead)*(sqrtS()/GeV-Elead+plead)); - bool singleDiff = false; + // Identify particles with most positive/most negative rapidities + const ParticleVector particlesByRap = cfs.particlesByRapidity(); + const Particle pslowest = particlesByRap.front(); + const Particle pfastest = particlesByRap.back(); + + // Find gap sizes + const ParticleVector particlesByEta = cfs.particlesByEta(); // sorted from minus to plus + const size_t num_particles = particlesByEta.size(); + vector<double> gaps; + for (size_t ip = 1; ip < num_particles; ++ip) { + const Particle& p1 = particlesByEta[ip-1]; + const Particle& p2 = particlesByEta[ip]; + const double gap = p2.momentum().eta() - p1.momentum().eta(); + assert(gap >= 0); + gaps.push_back(gap); + } + + // First, last, and largest gaps + const double gapmax = *max_element(gaps.begin(), gaps.end()); + const double gapbwd = gaps.front(); + const double gapfwd = gaps.back(); - // Fill SD - if (Mx < 200.) { - singleDiff = true; - if (fuzzyEquals(sqrtS()/GeV, 900, 1E-3)) { - _h_xsec_sd->fill( 900/GeV, weight); - } else if (fuzzyEquals(sqrtS()/GeV, 2760, 1E-3)) { - _h_xsec_sd->fill( 2760/GeV, weight); - } else if (fuzzyEquals(sqrtS()/GeV, 7000, 1E-3)) { - _h_xsec_sd->fill( 7000/GeV, weight); - } + // Mx calculation + FourMomentum p4lead; + if (pslowest.pdgId() == PROTON && pfastest.pdgId() == PROTON) { + p4lead = (fabs(pslowest.momentum().rapidity()) > fabs(pfastest.momentum().rapidity())) ? pslowest.momentum() : pfastest.momentum(); + } else if (pslowest.pdgId() == PROTON) { + p4lead = pslowest.momentum(); + } else if (pfastest.pdgId() == PROTON) { + p4lead = pfastest.momentum(); + } + const double Mx = sqrt( (sqrtS()-p4lead.E()-p4lead.vector3().mod()) * (sqrtS()-p4lead.E()+p4lead.vector3().mod()) ); + + // Fill SD (and escape) if Mx is sufficiently low + if (Mx < 200*GeV) { + _h_xsec_sd->fill(sqrtS()/GeV, weight); + return; } - if ( singleDiff ) vetoEvent; // DD events are defined as NSD with large gap. - - // also remove SD-like events in NSD events - if ( std::abs(gapbwd-LRG) < std::numeric_limits<double>::epsilon() || std::abs(gapfwd-LRG) < std::numeric_limits<double>::epsilon() ) vetoEvent; + // Also remove SD-like events in NSD events + if (fuzzyEquals(gapbwd, gapmax) || fuzzyEquals(gapfwd, gapmax)) vetoEvent; // Fill DD plots - if (LRG > 3.) { - if (fuzzyEquals(sqrtS()/GeV, 900, 1E-3)) { - _h_xsec_dd->fill( 900/GeV, weight); - } else if (fuzzyEquals(sqrtS()/GeV, 2760, 1E-3)) { - _h_xsec_dd->fill( 2760/GeV, weight); - } else if (fuzzyEquals(sqrtS()/GeV, 7000, 1E-3)) { - _h_xsec_dd->fill( 7000/GeV, weight); - } - } + if (gapmax > 3) _h_xsec_dd->fill(sqrtS()/GeV, weight); } + void finalize() { // get the ratio plots: SD/inel, DD/inel
More information about the Rivet-svn mailing list |