/** \class MomentumSmearing * * Performs transverse momentum resolution smearing. * * $Date$ * $Revision$ * * * \author P. Demin - UCL, Louvain-la-Neuve * */ #include "modules/MomentumSmearing.h" #include "classes/DelphesClasses.h" #include "classes/DelphesFactory.h" #include "classes/DelphesFormula.h" #include "ExRootAnalysis/ExRootResult.h" #include "ExRootAnalysis/ExRootFilter.h" #include "ExRootAnalysis/ExRootClassifier.h" #include "TMath.h" #include "TString.h" #include "TFormula.h" #include "TRandom3.h" #include "TObjArray.h" #include "TDatabasePDG.h" #include "TLorentzVector.h" #include #include #include #include using namespace std; //------------------------------------------------------------------------------ MomentumSmearing::MomentumSmearing() : fFormula(0), fItInputArray(0) { fFormula = new DelphesFormula; } //------------------------------------------------------------------------------ MomentumSmearing::~MomentumSmearing() { if(fFormula) delete fFormula; } //------------------------------------------------------------------------------ void MomentumSmearing::Init() { // read resolution formula fFormula->Compile(GetString("ResolutionFormula", "0.0")); // import input array fInputArray = ImportArray(GetString("InputArray", "ParticlePropagator/stableParticles")); fItInputArray = fInputArray->MakeIterator(); // create output array fOutputArray = ExportArray(GetString("OutputArray", "stableParticles")); } //------------------------------------------------------------------------------ void MomentumSmearing::Finish() { if(fItInputArray) delete fItInputArray; } //------------------------------------------------------------------------------ void MomentumSmearing::Process() { Candidate *candidate, *mother; Double_t pt, eta, phi; fItInputArray->Reset(); while((candidate = static_cast(fItInputArray->Next()))) { const TLorentzVector &candidatePosition = candidate->Position; const TLorentzVector &candidateMomentum = candidate->Momentum; eta = candidatePosition.Eta(); phi = candidatePosition.Phi(); pt = candidateMomentum.Pt(); // apply smearing formula pt = gRandom->Gaus(pt, fFormula->Eval(pt, eta) * pt); if(pt <= 0.0) continue; mother = candidate; candidate = static_cast(candidate->Clone()); eta = candidateMomentum.Eta(); phi = candidateMomentum.Phi(); candidate->Momentum.SetPtEtaPhiE(pt, eta, phi, pt*TMath::CosH(eta)); candidate->AddCandidate(mother); fOutputArray->Add(candidate); } } //------------------------------------------------------------------------------