source: HiSusy/trunk/Delphes/Delphes-3.0.9/modules/EnergySmearing.cc @ 5

Last change on this file since 5 was 5, checked in by zerwas, 11 years ago

update to Delphes-3.0.9

File size: 2.7 KB
Line 
1
2/** \class EnergySmearing
3 *
4 *  Performs energy resolution smearing.
5 *
6 *  $Date: 2013-02-22 12:03:12 +0100 (Fri, 22 Feb 2013) $
7 *  $Revision: 927 $
8 *
9 *
10 *  \author P. Demin - UCL, Louvain-la-Neuve
11 *
12 */
13
14#include "modules/EnergySmearing.h"
15
16#include "classes/DelphesClasses.h"
17#include "classes/DelphesFactory.h"
18#include "classes/DelphesFormula.h"
19
20#include "ExRootAnalysis/ExRootResult.h"
21#include "ExRootAnalysis/ExRootFilter.h"
22#include "ExRootAnalysis/ExRootClassifier.h"
23
24#include "TMath.h"
25#include "TString.h"
26#include "TFormula.h"
27#include "TRandom3.h"
28#include "TObjArray.h"
29#include "TDatabasePDG.h"
30#include "TLorentzVector.h"
31
32#include <algorithm> 
33#include <stdexcept>
34#include <iostream>
35#include <sstream>
36
37using namespace std;
38
39//------------------------------------------------------------------------------
40
41EnergySmearing::EnergySmearing() :
42  fFormula(0), fItInputArray(0)
43{
44  fFormula = new DelphesFormula;
45}
46
47//------------------------------------------------------------------------------
48
49EnergySmearing::~EnergySmearing()
50{
51  if(fFormula) delete fFormula;
52}
53
54//------------------------------------------------------------------------------
55
56void EnergySmearing::Init()
57{
58  // read resolution formula
59
60  fFormula->Compile(GetString("ResolutionFormula", "0.0"));
61
62  // import input array
63
64  fInputArray = ImportArray(GetString("InputArray", "ParticlePropagator/stableParticles"));
65  fItInputArray = fInputArray->MakeIterator();
66
67  // create output array
68
69  fOutputArray = ExportArray(GetString("OutputArray", "stableParticles"));
70}
71
72//------------------------------------------------------------------------------
73
74void EnergySmearing::Finish()
75{ 
76  if(fItInputArray) delete fItInputArray;
77}
78
79//------------------------------------------------------------------------------
80
81void EnergySmearing::Process()
82{
83  Candidate *candidate, *mother;
84  Double_t energy, eta, phi;
85
86  fItInputArray->Reset();
87  while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
88  {
89    const TLorentzVector &candidatePosition = candidate->Position;
90    const TLorentzVector &candidateMomentum = candidate->Momentum;
91    eta = candidatePosition.Eta();
92    phi = candidatePosition.Phi();
93    energy = candidateMomentum.E();
94 
95    // apply smearing formula
96    energy = gRandom->Gaus(energy, fFormula->Eval(0.0, eta, 0.0, energy));
97     
98    if(energy <= 0.0) continue;
99 
100    mother = candidate;
101    candidate = static_cast<Candidate*>(candidate->Clone());
102    eta = candidateMomentum.Eta();
103    phi = candidateMomentum.Phi();
104    candidate->Momentum.SetPtEtaPhiE(energy/TMath::CosH(eta), eta, phi, energy);
105    candidate->AddCandidate(mother);
106 
107    fOutputArray->Add(candidate);
108  }
109}
110
111//------------------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.