Fork me on GitHub

Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • modules/JetFakeParticle.cc

    r341014c r410f3a2  
    1717 */
    1818
     19
    1920/** \class JetFakeParticle
    2021 *
     
    3233#include "classes/DelphesFormula.h"
    3334
     35#include "ExRootAnalysis/ExRootResult.h"
     36#include "ExRootAnalysis/ExRootFilter.h"
    3437#include "ExRootAnalysis/ExRootClassifier.h"
    35 #include "ExRootAnalysis/ExRootFilter.h"
    36 #include "ExRootAnalysis/ExRootResult.h"
    37 
     38
     39#include "TMath.h"
     40#include "TString.h"
     41#include "TFormula.h"
     42#include "TRandom3.h"
     43#include "TObjArray.h"
    3844#include "TDatabasePDG.h"
    39 #include "TFormula.h"
    4045#include "TLorentzVector.h"
    41 #include "TMath.h"
    42 #include "TObjArray.h"
    43 #include "TRandom3.h"
    44 #include "TString.h"
    4546
    4647#include <algorithm>
     48#include <stdexcept>
    4749#include <iostream>
    4850#include <sstream>
    49 #include <stdexcept>
     51
    5052
    5153using namespace std;
     
    7981  fEfficiencyMap.clear();
    8082
    81   for(i = 0; i < size / 2; ++i)
     83  for(i = 0; i < size/2; ++i)
    8284  {
    8385    formula = new DelphesFormula;
    84     formula->Compile(param[i * 2 + 1].GetString());
    85     pdgCode = param[i * 2].GetInt();
     86    formula->Compile(param[i*2 + 1].GetString());
     87    pdgCode = param[i*2].GetInt();
    8688
    8789    if(TMath::Abs(pdgCode) != 11 && TMath::Abs(pdgCode) != 13 && TMath::Abs(pdgCode) != 22)
     
    9092    }
    9193
    92     fEfficiencyMap[param[i * 2].GetInt()] = formula;
     94    fEfficiencyMap[param[i*2].GetInt()] = formula;
    9395  }
    9496
     
    114116  fPhotonOutputArray = ExportArray(GetString("PhotonOutputArray", "fakePhotons"));
    115117  fJetOutputArray = ExportArray(GetString("JetOutputArray", "jets"));
     118
    116119}
    117120
     
    144147
    145148  fItInputArray->Reset();
    146   while((candidate = static_cast<Candidate *>(fItInputArray->Next())))
     149  while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
    147150  {
    148151    const TLorentzVector &candidatePosition = candidate->Position;
     
    167170      if(total <= r && r < total + p)
    168171      {
    169         fake = static_cast<Candidate *>(candidate->Clone());
     172        fake = static_cast<Candidate*>(candidate->Clone());
    170173
    171174        // convert jet
     
    175178          if(candidate->Charge != 0)
    176179          {
    177             fake->Charge = candidate->Charge / TMath::Abs(candidate->Charge);
     180            fake->Charge = candidate->Charge/TMath::Abs(candidate->Charge);
    178181          }
    179182          else
     
    181184            rs = gRandom->Uniform();
    182185            fake->Charge = (rs < 0.5) ? -1 : 1;
     186           
    183187          }
    184188        }
Note: See TracChangeset for help on using the changeset viewer.