Fork me on GitHub

source: git/modules/PileUpMerger.cc@ b20af7f

ImprovedOutputFile Timing dual_readout llp
Last change on this file since b20af7f was 76d3973, checked in by pavel <pavel@…>, 11 years ago

added PileUpDistribution parameter to the PileUpMerger module

  • Property mode set to 100644
File size: 3.7 KB
Line 
1/** \class PileUpMerger
2 *
3 * Merges particles from pile-up sample into event
4 *
5 *
6 * $Date: 2013-02-12 15:13:59 +0100 (Tue, 12 Feb 2013) $
7 * $Revision: 907 $
8 *
9 *
10 * \author M. Selvaggi - UCL, Louvain-la-Neuve
11 *
12 */
13
14#include "modules/PileUpMerger.h"
15
16#include "classes/DelphesClasses.h"
17#include "classes/DelphesFactory.h"
18#include "classes/DelphesFormula.h"
19#include "classes/DelphesPileUpReader.h"
20
21#include "ExRootAnalysis/ExRootResult.h"
22#include "ExRootAnalysis/ExRootFilter.h"
23#include "ExRootAnalysis/ExRootClassifier.h"
24
25#include "TMath.h"
26#include "TString.h"
27#include "TFormula.h"
28#include "TRandom3.h"
29#include "TObjArray.h"
30#include "TDatabasePDG.h"
31#include "TLorentzVector.h"
32
33#include <algorithm>
34#include <stdexcept>
35#include <iostream>
36#include <sstream>
37
38using namespace std;
39
40//------------------------------------------------------------------------------
41
42PileUpMerger::PileUpMerger() :
43 fReader(0), fItInputArray(0)
44{
45}
46
47//------------------------------------------------------------------------------
48
49PileUpMerger::~PileUpMerger()
50{
51}
52
53//------------------------------------------------------------------------------
54
55void PileUpMerger::Init()
56{
57 const char *fileName;
58
59 fPileUpDistribution = GetInt("PileUpDistribution", 0);
60
61 fMeanPileUp = GetDouble("MeanPileUp", 10);
62 fZVertexSpread = GetDouble("ZVertexSpread", 0.05)*1.0E3;
63
64 fileName = GetString("PileUpFile", "MinBias.pileup");
65 fReader = new DelphesPileUpReader(fileName);
66
67 // import input array
68 fInputArray = ImportArray(GetString("InputArray", "Delphes/stableParticles"));
69 fItInputArray = fInputArray->MakeIterator();
70
71 // create output arrays
72 fOutputArray = ExportArray(GetString("OutputArray", "stableParticles"));
73}
74
75//------------------------------------------------------------------------------
76
77void PileUpMerger::Finish()
78{
79 if(fReader) delete fReader;
80}
81
82//------------------------------------------------------------------------------
83
84void PileUpMerger::Process()
85{
86 TDatabasePDG *pdg = TDatabasePDG::Instance();
87 TParticlePDG *pdgParticle;
88 Int_t pid;
89 Float_t x, y, z, t;
90 Float_t px, py, pz, e;
91 Double_t dz, dphi;
92 Int_t numberOfEvents, event;
93 Long64_t allEntries, entry;
94 Candidate *candidate;
95 DelphesFactory *factory;
96
97 fItInputArray->Reset();
98 while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
99 {
100 fOutputArray->Add(candidate);
101 }
102
103 factory = GetFactory();
104
105 switch(fPileUpDistribution)
106 {
107 case 0:
108 numberOfEvents = gRandom->Poisson(fMeanPileUp);
109 break;
110 case 1:
111 numberOfEvents = gRandom->Integer(2*fMeanPileUp + 1);
112 break;
113 default:
114 numberOfEvents = gRandom->Poisson(fMeanPileUp);
115 break;
116 }
117
118 allEntries = fReader->GetEntries();
119
120 for(event = 0; event < numberOfEvents; ++event)
121 {
122 do
123 {
124 entry = TMath::Nint(gRandom->Rndm()*allEntries);
125 }
126 while(entry >= allEntries);
127
128 fReader->ReadEntry(entry);
129
130 dz = gRandom->Gaus(0.0, fZVertexSpread);
131 dphi = gRandom->Uniform(-TMath::Pi(), TMath::Pi());
132
133 while(fReader->ReadParticle(pid, x, y, z, t, px, py, pz, e))
134 {
135 candidate = factory->NewCandidate();
136
137 candidate->PID = pid;
138
139 candidate->Status = 1;
140
141 pdgParticle = pdg->GetParticle(pid);
142 candidate->Charge = pdgParticle ? Int_t(pdgParticle->Charge()/3.0) : -999;
143 candidate->Mass = pdgParticle ? pdgParticle->Mass() : -999.9;
144
145 candidate->IsPU = 1;
146
147 candidate->Momentum.SetPxPyPzE(px, py, pz, e);
148 candidate->Momentum.RotateZ(dphi);
149
150 candidate->Position.SetXYZT(x, y, z + dz, t);
151 candidate->Position.RotateZ(dphi);
152
153 fOutputArray->Add(candidate);
154 }
155 }
156}
157
158//------------------------------------------------------------------------------
159
Note: See TracBrowser for help on using the repository browser.