Fork me on GitHub

source: git/modules/Weighter.cc@ 6c6efd1

3.4.3pre02
Last change on this file since 6c6efd1 was 341014c, checked in by Pavel Demin <pavel-demin@…>, 6 years ago

apply .clang-format to all .h, .cc and .cpp files

  • Property mode set to 100644
File size: 4.7 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/** \class Weighter
20 *
21 * Apply a weight depending on PDG code.
22 *
23 * \author P. Demin - UCL, Louvain-la-Neuve
24 *
25 */
26
27#include "modules/Weighter.h"
28
29#include "classes/DelphesClasses.h"
30#include "classes/DelphesFactory.h"
31#include "classes/DelphesFormula.h"
32
33#include "ExRootAnalysis/ExRootClassifier.h"
34#include "ExRootAnalysis/ExRootFilter.h"
35#include "ExRootAnalysis/ExRootResult.h"
36
37#include "TDatabasePDG.h"
38#include "TFormula.h"
39#include "TLorentzVector.h"
40#include "TMath.h"
41#include "TObjArray.h"
42#include "TRandom3.h"
43#include "TString.h"
44
45#include <algorithm>
46#include <iostream>
47#include <sstream>
48#include <stdexcept>
49
50using namespace std;
51
52//------------------------------------------------------------------------------
53
54bool Weighter::TIndexStruct::operator<(const Weighter::TIndexStruct &value) const
55{
56 Int_t i;
57
58 for(i = 0; i < 4; ++i)
59 {
60 if(codes[i] != value.codes[i]) return codes[i] < value.codes[i];
61 }
62
63 return false;
64}
65
66//------------------------------------------------------------------------------
67
68Weighter::Weighter() :
69 fItInputArray(0)
70{
71}
72
73//------------------------------------------------------------------------------
74
75Weighter::~Weighter()
76{
77}
78
79//------------------------------------------------------------------------------
80
81void Weighter::Init()
82{
83 ExRootConfParam param, paramCodes;
84 Int_t i, j, size, sizeCodes;
85 Int_t code;
86 TIndexStruct index;
87 Double_t weight;
88
89 fWeightSet.clear();
90
91 // set default weight value
92 fWeightMap.clear();
93 memset(index.codes, 0, 4 * sizeof(Int_t));
94 fWeightMap[index] = 1.0;
95
96 // read weights
97 param = GetParam("Weight");
98 size = param.GetSize();
99 for(i = 0; i < size / 2; ++i)
100 {
101 paramCodes = param[i * 2];
102 sizeCodes = paramCodes.GetSize();
103 weight = param[i * 2 + 1].GetDouble();
104
105 if(sizeCodes < 1 || sizeCodes > 4)
106 {
107 throw runtime_error("only 1, 2, 3 or 4 PDG codes can be specified per weight");
108 }
109
110 memset(index.codes, 0, 4 * sizeof(Int_t));
111
112 for(j = 0; j < sizeCodes; ++j)
113 {
114 code = paramCodes[j].GetInt();
115 index.codes[j] = code;
116 fWeightSet.insert(code);
117 }
118
119 sort(index.codes, index.codes + 4);
120
121 fWeightMap[index] = weight;
122 }
123
124 // import input array(s)
125
126 fInputArray = ImportArray(GetString("InputArray", "Delphes/allParticles"));
127 fItInputArray = fInputArray->MakeIterator();
128
129 // create output array(s)
130
131 fOutputArray = ExportArray(GetString("OutputArray", "weight"));
132}
133
134//------------------------------------------------------------------------------
135
136void Weighter::Finish()
137{
138 if(fItInputArray) delete fItInputArray;
139}
140
141//------------------------------------------------------------------------------
142
143void Weighter::Process()
144{
145 Candidate *candidate;
146 Int_t i;
147 TIndexStruct index;
148 Double_t weight;
149 set<Int_t>::iterator itCodeSet;
150 map<TIndexStruct, Double_t>::iterator itWeightMap;
151
152 DelphesFactory *factory = GetFactory();
153
154 // loop over all particles
155 fCodeSet.clear();
156 fItInputArray->Reset();
157 while((candidate = static_cast<Candidate *>(fItInputArray->Next())))
158 {
159 if(candidate->Status != 3) continue;
160
161 if(fWeightSet.find(candidate->PID) == fWeightSet.end()) continue;
162
163 fCodeSet.insert(candidate->PID);
164 }
165
166 // find default weight value
167 memset(index.codes, 0, 4 * sizeof(Int_t));
168 itWeightMap = fWeightMap.find(index);
169 weight = itWeightMap->second;
170
171 if(fCodeSet.size() <= 4)
172 {
173 i = 0;
174 for(itCodeSet = fCodeSet.begin(); itCodeSet != fCodeSet.end(); ++itCodeSet)
175 {
176 index.codes[i] = *itCodeSet;
177 ++i;
178 }
179
180 sort(index.codes, index.codes + 4);
181
182 itWeightMap = fWeightMap.find(index);
183 if(itWeightMap != fWeightMap.end())
184 {
185 weight = itWeightMap->second;
186 }
187 }
188
189 candidate = factory->NewCandidate();
190 candidate->Momentum.SetPtEtaPhiE(weight, 0.0, 0.0, weight);
191 fOutputArray->Add(candidate);
192}
193
194//------------------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.