Skip to content

Commit

Permalink
Merge pull request #12680 from cms-btv-pog/BoostedDoubleSVTaggerV2-Th…
Browse files Browse the repository at this point in the history
…readSafe_from-CMSSW_7_6_1

Making the boosted double SV tagger thread-safe
  • Loading branch information
cmsbuild committed Dec 8, 2015
2 parents ec5722e + b7970d9 commit 97777d9
Showing 2 changed files with 7 additions and 7 deletions.
Original file line number Diff line number Diff line change
@@ -15,7 +15,6 @@
#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"

#include "fastjet/PseudoJet.hh"
#include "fastjet/contrib/Njettiness.hh"

class CandidateBoostedDoubleSecondaryVertexComputer : public JetTagComputer {

@@ -34,8 +33,6 @@ class CandidateBoostedDoubleSecondaryVertexComputer : public JetTagComputer {

const double beta_;
const double R0_;
// N-subjettiness calculator
fastjet::contrib::Njettiness njettiness_;

const double maxSVDeltaRToJet_;
const bool useCondDB_;
Original file line number Diff line number Diff line change
@@ -14,11 +14,11 @@
#include "RecoBTau/JetTagComputer/interface/JetTagComputer.h"
#include "TrackingTools/IPTools/interface/IPTools.h"

#include "fastjet/contrib/Njettiness.hh"

CandidateBoostedDoubleSecondaryVertexComputer::CandidateBoostedDoubleSecondaryVertexComputer(const edm::ParameterSet & parameters) :
beta_(parameters.getParameter<double>("beta")),
R0_(parameters.getParameter<double>("R0")),
njettiness_(fastjet::contrib::OnePass_KT_Axes(), fastjet::contrib::NormalizedMeasure(beta_,R0_)),
maxSVDeltaRToJet_(parameters.getParameter<double>("maxSVDeltaRToJet")),
useCondDB_(parameters.getParameter<bool>("useCondDB")),
gbrForestLabel_(parameters.existsAs<std::string>("gbrForestLabel") ? parameters.getParameter<std::string>("gbrForestLabel") : ""),
@@ -583,10 +583,13 @@ void CandidateBoostedDoubleSecondaryVertexComputer::calcNsubjettiness(const reco
edm::LogWarning("MissingJetConstituent") << "Jet constituent required for N-subjettiness computation is missing!";
}

// N-subjettiness calculator
fastjet::contrib::Njettiness njettiness(fastjet::contrib::OnePass_KT_Axes(), fastjet::contrib::NormalizedMeasure(beta_,R0_));

// calculate N-subjettiness
tau1 = njettiness_.getTau(1, fjParticles);
tau2 = njettiness_.getTau(2, fjParticles);
currentAxes = njettiness_.currentAxes();
tau1 = njettiness.getTau(1, fjParticles);
tau2 = njettiness.getTau(2, fjParticles);
currentAxes = njettiness.currentAxes();
}


0 comments on commit 97777d9

Please sign in to comment.