Fork me on GitHub

Changeset fe0273c in git for modules


Ignore:
Timestamp:
Jun 26, 2015, 12:08:48 PM (10 years ago)
Author:
Pavel Demin <pavel.demin@…>
Branches:
ImprovedOutputFile, Timing, dual_readout, llp, master
Children:
28027d5
Parents:
d38348d
Message:

merge BTagging and BTaggingCMS

Location:
modules
Files:
2 deleted
4 edited

Legend:

Unmodified
Added
Removed
  • modules/BTagging.cc

    rd38348d rfe0273c  
    2222 *  Determines origin of jet,
    2323 *  applies b-tagging efficiency (miss identification rate) formulas
    24  *  and sets b-tagging flags 
     24 *  and sets b-tagging flags
    2525 *
    2626 *  \author P. Demin - UCL, Louvain-la-Neuve
     
    4646#include "TLorentzVector.h"
    4747
    48 #include <algorithm> 
     48#include <algorithm>
    4949#include <stdexcept>
    5050#include <iostream>
     
    6262
    6363  Int_t GetCategory(TObject *object);
    64  
     64
    6565  Double_t fEtaMax, fPTMin;
    6666};
     
    7575
    7676  if(momentum.Pt() <= fPTMin || TMath::Abs(momentum.Eta()) > fEtaMax) return -1;
    77  
     77
    7878  pdgCode = TMath::Abs(parton->PID);
    7979  if(pdgCode != 21 && pdgCode > 5) return -1;
     
    117117  param = GetParam("EfficiencyFormula");
    118118  size = param.GetSize();
    119  
     119
    120120  fEfficiencyMap.clear();
    121121  for(i = 0; i < size/2; ++i)
     
    143143
    144144  fFilter = new ExRootFilter(fPartonInputArray);
    145  
     145
    146146  fJetInputArray = ImportArray(GetString("JetInputArray", "FastJetFinder/jets"));
    147147  fItJetInputArray = fJetInputArray->MakeIterator();
     
    170170void BTagging::Process()
    171171{
    172   Candidate *jet, *parton;
     172  Candidate *jet;
    173173  Double_t pt, eta, phi, e;
    174174  TObjArray *partonArray;
    175175  map< Int_t, DelphesFormula * >::iterator itEfficiencyMap;
    176176  DelphesFormula *formula;
    177   Int_t pdgCode, pdgCodeMax;
    178177
    179178  // select quark and gluons
    180179  fFilter->Reset();
    181180  partonArray = fFilter->GetSubArray(fClassifier, 0);
    182  
     181
    183182  if(partonArray == 0) return;
    184183
    185184  TIter itPartonArray(partonArray);
    186  
     185
    187186  // loop over all input jets
    188187  fItJetInputArray->Reset();
     
    190189  {
    191190    const TLorentzVector &jetMomentum = jet->Momentum;
    192     pdgCodeMax = -1;
    193191    eta = jetMomentum.Eta();
    194192    phi = jetMomentum.Phi();
    195193    pt = jetMomentum.Pt();
    196194    e = jetMomentum.E();
    197    
    198     // loop over all input partons
    199     itPartonArray.Reset();
    200     while((parton = static_cast<Candidate*>(itPartonArray.Next())))
    201     {
    202       pdgCode = TMath::Abs(parton->PID);
    203       if(pdgCode == 21) pdgCode = 0;
    204       if(jetMomentum.DeltaR(parton->Momentum) <= fDeltaR)
    205       {
    206         if(pdgCodeMax < pdgCode) pdgCodeMax = pdgCode;
    207       }
    208     }
    209     if(pdgCodeMax == 0) pdgCodeMax = 21;
    210     if(pdgCodeMax == -1) pdgCodeMax = 0;
    211195
    212196    // find an efficency formula
    213     itEfficiencyMap = fEfficiencyMap.find(pdgCodeMax);
     197    itEfficiencyMap = fEfficiencyMap.find(jet->Flavor);
    214198    if(itEfficiencyMap == fEfficiencyMap.end())
    215199    {
     
    220204    // apply an efficency formula
    221205    jet->BTag |= (gRandom->Uniform() <= formula->Eval(pt, eta, phi, e)) << fBitNumber;
    222   }
    223 }
    224 
    225 //------------------------------------------------------------------------------
     206
     207    // find an efficency formula for algo flavor definition
     208    itEfficiencyMap = fEfficiencyMap.find(jet->FlavorAlgo);
     209    if(itEfficiencyMap == fEfficiencyMap.end())
     210    {
     211      itEfficiencyMap = fEfficiencyMap.find(0);
     212    }
     213    formula = itEfficiencyMap->second;
     214
     215    // apply an efficency formula
     216    jet->BTagAlgo |= (gRandom->Uniform() <= formula->Eval(pt, eta, phi, e)) << fBitNumber;
     217
     218    // find an efficency formula for phys flavor definition
     219    itEfficiencyMap = fEfficiencyMap.find(jet->FlavorPhys);
     220    if(itEfficiencyMap == fEfficiencyMap.end())
     221    {
     222      itEfficiencyMap = fEfficiencyMap.find(0);
     223    }
     224    formula = itEfficiencyMap->second;
     225
     226    // apply an efficency formula
     227    jet->BTagPhys |= (gRandom->Uniform() <= formula->Eval(pt, eta, phi, e)) << fBitNumber;
     228  }
     229}
     230
     231//------------------------------------------------------------------------------
  • modules/JetFlavorAssociation.cc

    rd38348d rfe0273c  
    216216{
    217217  float maxPt = 0;
    218   float minDr = 1000;
    219218  bool isGoodParton = true;
    220219  int daughterCounter = 0;
    221220  Candidate *parton, *partonLHEF;
    222   Candidate *tempParton = 0, *tempPartonNearest = 0, *tempPartonHighestPt = 0;
     221  Candidate *tempParton = 0, *tempPartonHighestPt = 0;
    223222  int pdgCode, pdgCodeMax = -1;
    224223
     
    263262      if(jet->Momentum.DeltaR(parton->Momentum) <= fDeltaR)
    264263      {
    265         if(jet->Momentum.DeltaR(parton->Momentum) < minDr)
    266         {
    267           minDr = jet->Momentum.DeltaR(parton->Momentum);
    268           tempPartonNearest = parton;
    269         }
    270 
    271264        // if not yet found && pdgId is a c, take as c
    272265        if(TMath::Abs(parton->PID) == 4) tempParton = parton;
     
    281274  }
    282275
    283   jet->FlavorHeaviest = tempParton ? TMath::Abs(tempParton->PID) : 0;
    284   jet->FlavorHighestPt = tempPartonHighestPt ? TMath::Abs(tempPartonHighestPt->PID) : 0;
    285   jet->FlavorNearest2 = tempPartonNearest ? TMath::Abs(tempPartonNearest->PID) : 0;
    286276  if(!tempParton) tempParton = tempPartonHighestPt;
    287277  jet->FlavorAlgo = tempParton ? TMath::Abs(tempParton->PID) : 0;
     
    290280  if(pdgCodeMax == -1) pdgCodeMax = 0;
    291281
    292   jet->FlavorDefault = pdgCodeMax;
     282  jet->Flavor = pdgCodeMax;
    293283}
    294284
     
    297287void JetFlavorAssociation::GetPhysicsFlavor(Candidate *jet, TIter &itPartonArray, TIter &itPartonLHEFArray)
    298288{
    299   float minDr = 1000;
    300289  int partonCounter = 0;
    301290  float biggerConeSize = 0.7;
     
    305294  int motherCounter = 0;
    306295  Candidate *parton, *partonLHEF, *mother1, *mother2;
    307   Candidate *tempParton = 0, *tempPartonNearest = 0;
     296  Candidate *tempParton = 0;
    308297  vector<Candidate *> contaminations;
    309298  vector<Candidate *>::iterator itContaminations;
     
    315304  {
    316305    dist = jet->Momentum.DeltaR(partonLHEF->Momentum); // take the DR
    317     if(partonLHEF->Status == 1 && dist < minDr)
    318     {
    319       tempPartonNearest = partonLHEF;
    320       minDr = dist;
    321     }
    322306
    323307    if(partonLHEF->Status == 1 && dist <= fDeltaR)
     
    354338  }
    355339
    356   jet->FlavorNearest3 = tempPartonNearest ? TMath::Abs(tempPartonNearest->PID) : 0;
    357 
    358340  if(partonCounter != 1)
    359341  {
    360     jet->FlavorPhysics = 0;
     342    jet->FlavorPhys = 0;
    361343  }
    362344  else if(contaminations.size() == 0)
    363345  {
    364     jet->FlavorPhysics = TMath::Abs(tempParton->PID);
     346    jet->FlavorPhys = TMath::Abs(tempParton->PID);
    365347  }
    366348  else if(contaminations.size() > 0)
    367349  {
    368     jet->FlavorPhysics = TMath::Abs(tempParton->PID);
     350    jet->FlavorPhys = TMath::Abs(tempParton->PID);
    369351
    370352    for(itContaminations = contaminations.begin(); itContaminations != contaminations.end(); ++itContaminations)
     
    391373        // keep association --> the initialParton is a c --> the contaminated parton is a c
    392374        if(contaminatingFlavor == 4) continue;
    393         jet->FlavorPhysics = 0; // all the other cases reject!
     375        jet->FlavorPhys = 0; // all the other cases reject!
    394376        break;
    395377      }
  • modules/ModulesLinkDef.h

    rd38348d rfe0273c  
    4242#include "modules/UniqueObjectFinder.h"
    4343#include "modules/TrackCountingBTagging.h"
    44 #include "modules/BTaggingCMS.h"
    4544#include "modules/BTagging.h"
    4645#include "modules/TauTagging.h"
     
    8584#pragma link C++ class UniqueObjectFinder+;
    8685#pragma link C++ class TrackCountingBTagging+;
    87 #pragma link C++ class BTaggingCMS+;
    8886#pragma link C++ class BTagging+;
    8987#pragma link C++ class TauTagging+;
  • modules/TreeWriter.cc

    rd38348d rfe0273c  
    490490    const TLorentzVector &position = candidate->Position;
    491491
    492 
    493492    pt = momentum.Pt();
    494493    cosTheta = TMath::Abs(momentum.CosTheta());
     
    567566    entry->DeltaPhi = candidate->DeltaPhi;
    568567
     568    entry->Flavor = candidate->Flavor;
     569    entry->FlavorAlgo = candidate->FlavorAlgo;
     570    entry->FlavorPhys = candidate->FlavorPhys;
     571
    569572    entry->BTag = candidate->BTag;
    570573
    571574    entry->BTagAlgo = candidate->BTagAlgo;
    572     entry->BTagDefault = candidate->BTagDefault;
    573     entry->BTagPhysics = candidate->BTagPhysics;
    574     entry->BTagNearest2 = candidate->BTagNearest2;
    575     entry->BTagNearest3 = candidate->BTagNearest3;
    576     entry->BTagHeaviest = candidate->BTagHeaviest;
    577     entry->BTagHighestPt = candidate->BTagHighestPt;
    578 
    579     entry->FlavorAlgo = candidate->FlavorAlgo;
    580     entry->FlavorDefault = candidate->FlavorDefault;
    581     entry->FlavorPhysics = candidate->FlavorPhysics;
    582     entry->FlavorNearest2 = candidate->FlavorNearest2;
    583     entry->FlavorNearest3 = candidate->FlavorNearest3;
    584     entry->FlavorHeaviest = candidate->FlavorHeaviest;
    585     entry->FlavorHighestPt = candidate->FlavorHighestPt;
     575    entry->BTagPhys = candidate->BTagPhys;
    586576
    587577    entry->TauTag = candidate->TauTag;
Note: See TracChangeset for help on using the changeset viewer.