Fork me on GitHub

source: git/modules/TrackPileUpSubtractor.cc@ 9c491e1

ImprovedOutputFile Timing dual_readout llp
Last change on this file since 9c491e1 was cab38f6, checked in by Pavel Demin <pavel.demin@…>, 10 years ago

remove svn tags and fix formatting

  • Property mode set to 100644
File size: 4.1 KB
Line 
1/*
2 * Delphes: a framework for fast simulation of a generic collider experiment
3 * Copyright (C) 2012-2014 Universite catholique de Louvain (UCL), Belgium
4 *
5 * This program is free software: you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation, either version 3 of the License, or
8 * (at your option) any later version.
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program. If not, see <http://www.gnu.org/licenses/>.
17 */
18
19
20/** \class TrackPileUpSubtractor
21 *
22 * Subtract pile-up contribution from tracks.
23 *
24 * \author P. Demin - UCL, Louvain-la-Neuve
25 *
26 */
27
28#include "modules/TrackPileUpSubtractor.h"
29
30#include "classes/DelphesClasses.h"
31#include "classes/DelphesFactory.h"
32#include "classes/DelphesFormula.h"
33
34#include "ExRootAnalysis/ExRootResult.h"
35#include "ExRootAnalysis/ExRootFilter.h"
36#include "ExRootAnalysis/ExRootClassifier.h"
37
38#include "TMath.h"
39#include "TString.h"
40#include "TFormula.h"
41#include "TRandom3.h"
42#include "TObjArray.h"
43#include "TDatabasePDG.h"
44#include "TLorentzVector.h"
45
46#include <algorithm>
47#include <stdexcept>
48#include <iostream>
49#include <sstream>
50
51using namespace std;
52
53//------------------------------------------------------------------------------
54
55TrackPileUpSubtractor::TrackPileUpSubtractor()
56{
57}
58
59//------------------------------------------------------------------------------
60
61TrackPileUpSubtractor::~TrackPileUpSubtractor()
62{
63}
64
65//------------------------------------------------------------------------------
66
67void TrackPileUpSubtractor::Init()
68{
69// import input array
70
71 fVertexInputArray = ImportArray(GetString("VertexInputArray", "PileUpMerger/vertices"));
72 fItVertexInputArray = fVertexInputArray->MakeIterator();
73
74 fZVertexResolution = GetDouble("ZVertexResolution", 0.005)*1.0E3;
75
76 // import arrays with output from other modules
77
78 ExRootConfParam param = GetParam("InputArray");
79 Long_t i, size;
80 const TObjArray *array;
81 TIterator *iterator;
82
83 size = param.GetSize();
84 for(i = 0; i < size/2; ++i)
85 {
86 array = ImportArray(param[i*2].GetString());
87 iterator = array->MakeIterator();
88
89 fInputMap[iterator] = ExportArray(param[i*2 + 1].GetString());
90 }
91}
92
93//------------------------------------------------------------------------------
94
95void TrackPileUpSubtractor::Finish()
96{
97 map< TIterator *, TObjArray * >::iterator itInputMap;
98 TIterator *iterator;
99
100 for(itInputMap = fInputMap.begin(); itInputMap != fInputMap.end(); ++itInputMap)
101 {
102 iterator = itInputMap->first;
103
104 if(iterator) delete iterator;
105 }
106
107 if(fItVertexInputArray) delete fItVertexInputArray;
108}
109
110//------------------------------------------------------------------------------
111
112void TrackPileUpSubtractor::Process()
113{
114 Candidate *candidate, *particle;
115 map< TIterator *, TObjArray * >::iterator itInputMap;
116 TIterator *iterator;
117 TObjArray *array;
118 Double_t z, zvtx=0;
119
120
121 // find z position of primary vertex
122
123 fItVertexInputArray->Reset();
124 while((candidate = static_cast<Candidate*>(fItVertexInputArray->Next())))
125 {
126 if(!candidate->IsPU)
127 {
128 zvtx = candidate->Position.Z();
129 // break;
130 }
131 }
132
133 // loop over all input arrays
134 for(itInputMap = fInputMap.begin(); itInputMap != fInputMap.end(); ++itInputMap)
135 {
136 iterator = itInputMap->first;
137 array = itInputMap->second;
138
139 // loop over all candidates
140 iterator->Reset();
141 while((candidate = static_cast<Candidate*>(iterator->Next())))
142 {
143 particle = static_cast<Candidate*>(candidate->GetCandidates()->At(0));
144 z = particle->Position.Z();
145
146 // apply pile-up subtraction
147 // assume perfect pile-up subtraction for tracks outside fZVertexResolution
148
149 if(candidate->IsPU && TMath::Abs(z-zvtx) > fZVertexResolution) continue;
150
151 array->Add(candidate);
152 }
153 }
154}
155
156//------------------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.