Rivet analyses referenceALEPH_1996_I398193Measurement of the quark to photon fragmentation functionExperiment: ALEPH (LEP Run 1) Inspire ID: 398193 Status: VALIDATED Authors:
Beam energies: (45.6, 45.6) GeV Run details:
Earlier measurements at LEP of isolated hard photons in hadronic Z decays, attributed to radiation from primary quark pairs, have been extended in the ALEPH experiment to include hard photon production inside hadron jets. Events are selected where all particles combine democratically to form hadron jets, one of which contains a photon with a fractional energy $z > 0.7$. After statistical subtraction of non-prompt photons, the quark-to-photon fragmentation function, $D(z)$, is extracted directly from the measured 2-jet rate. Source code: ALEPH_1996_I398193.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/IdentifiedFinalState.hh"
6#include "Rivet/Projections/FastJets.hh"
7#include "Rivet/Projections/Thrust.hh"
8
9namespace Rivet {
10
11
12 /// ALEPH measurement of quark-to-photon fragmentation function
13 class ALEPH_1996_I398193 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(ALEPH_1996_I398193);
18
19
20 /// @name Analysis methods
21 /// @{
22
23 void init() {
24 // Set up projections
25 FinalState fs;
26 declare(FastJets(fs, JetAlg::DURHAM, 0.7), "DurhamJets");
27 IdentifiedFinalState ifs; //(Cuts::pT > 0);
28 ifs.acceptId(PID::PHOTON);
29 declare(ifs, "Photons");
30 declare(Thrust(fs), "Thrust");
31 declare(ChargedFinalState(), "CFS");
32
33 // Book histograms
34 book(_h_z_2jet_001 ,1, 1, 1);
35 book(_h_z_2jet_006 ,2, 1, 1);
36 book(_h_z_2jet_01 ,3, 1, 1);
37 book(_h_z_2jet_033 ,4, 1, 1);
38 book(_h_z_3jet_001 ,5, 1, 1);
39 book(_h_z_3jet_006 ,6, 1, 1);
40 book(_h_z_3jet_01 ,7, 1, 1);
41 book(_h_z_4jet_001 ,8, 1, 1);
42 }
43
44
45 /// Perform the per-event analysis
46 void analyze(const Event& event) {
47
48 if (apply<FinalState>(event, "CFS").particles().size() < 2) vetoEvent;
49
50 const Particles allphotons = apply<IdentifiedFinalState>(event, "Photons").particles();
51 Particles photons;
52 for (const Particle& photon : allphotons)
53 if (fabs(cos(photon.theta())) < 0.95 && photon.E() > 5.0*GeV)
54 photons.push_back(photon);
55 if (photons.size() < 1) vetoEvent;
56
57 const Thrust& thrust = apply<Thrust>(event, "Thrust");
58 if (fabs(cos(thrust.thrustAxis().theta()))>0.95) vetoEvent;
59
60 const FastJets& durjet = apply<FastJets>(event, "DurhamJets");
61
62 for (const Particle& photon : photons) {
63
64 PseudoJets jets_001 = durjet.clusterSeq()->exclusive_jets_ycut(0.01);
65 for (const fastjet::PseudoJet& jet : jets_001) {
66 if (particleInJet(photon, jet)) {
67 double zgamma = photon.E()/jet.E();
68 if (jets_001.size() == 2) _h_z_2jet_001->fill(zgamma);
69 else if (jets_001.size() == 3) _h_z_3jet_001->fill(zgamma);
70 else if (jets_001.size() > 3) _h_z_4jet_001->fill(zgamma);
71 break;
72 }
73 }
74
75 PseudoJets jets_006 = durjet.clusterSeq()->exclusive_jets_ycut(0.06);
76 for (const fastjet::PseudoJet& jet : jets_006) {
77 if (particleInJet(photon, jet)) {
78 double zgamma = photon.E()/jet.E();
79 if (jets_006.size() == 2) _h_z_2jet_006->fill(zgamma);
80 else if (jets_006.size() == 3) _h_z_3jet_006->fill(zgamma);
81 break;
82 }
83 }
84
85 PseudoJets jets_01 = durjet.clusterSeq()->exclusive_jets_ycut(0.1);
86 for (const fastjet::PseudoJet& jet : jets_01) {
87 if (particleInJet(photon, jet)) {
88 double zgamma = photon.E()/jet.E();
89 if (jets_01.size() == 2) _h_z_2jet_01->fill(zgamma);
90 else if (jets_01.size() == 3) _h_z_3jet_01->fill(zgamma);
91 break;
92 }
93 }
94
95 PseudoJets jets_033 = durjet.clusterSeq()->exclusive_jets_ycut(0.33);
96 for (const fastjet::PseudoJet& jet : jets_033) {
97 if (particleInJet(photon, jet)) {
98 double zgamma = photon.E()/jet.E();
99 if (jets_033.size() == 2) _h_z_2jet_033->fill(zgamma);
100 break;
101 }
102 }
103
104 }
105 }
106
107
108 bool particleInJet(const Particle& p, const fastjet::PseudoJet& jet) {
109 for (const fastjet::PseudoJet& jetpart : jet.constituents()) {
110 if (fuzzyEquals(jetpart.E(), p.E()) &&
111 fuzzyEquals(jetpart.px(), p.px()) &&
112 fuzzyEquals(jetpart.py(), p.py()) &&
113 fuzzyEquals(jetpart.pz(), p.pz())) {
114 return true;
115 }
116 }
117 return false;
118 }
119
120
121 /// Normalise histograms etc., after the run
122 void finalize() {
123 scale(_h_z_2jet_001, 1000.0/sumOfWeights());
124 scale(_h_z_2jet_006, 1000.0/sumOfWeights());
125 scale(_h_z_2jet_01, 1000.0/sumOfWeights());
126 scale(_h_z_2jet_033, 1000.0/sumOfWeights());
127 scale(_h_z_3jet_001, 1000.0/sumOfWeights());
128 scale(_h_z_3jet_006, 1000.0/sumOfWeights());
129 scale(_h_z_3jet_01, 1000.0/sumOfWeights());
130 scale(_h_z_4jet_001, 1000.0/sumOfWeights());
131 }
132
133 /// @}
134
135
136 private:
137
138 /// @name Histograms
139 /// @{
140 Histo1DPtr _h_z_2jet_001, _h_z_2jet_006, _h_z_2jet_01, _h_z_2jet_033;
141 Histo1DPtr _h_z_3jet_001, _h_z_3jet_006, _h_z_3jet_01;
142 Histo1DPtr _h_z_4jet_001;
143 /// @}
144
145 };
146
147
148
149 RIVET_DECLARE_ALIASED_PLUGIN(ALEPH_1996_I398193, ALEPH_1996_S3196992);
150
151}
|