CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Private Attributes
Onia2MuMuPAT Class Reference

#include <Onia2MuMuPAT.h>

Inheritance diagram for Onia2MuMuPAT:
edm::EDProducer edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Public Member Functions

 Onia2MuMuPAT (const edm::ParameterSet &)
 
 ~Onia2MuMuPAT ()
 
- Public Member Functions inherited from edm::EDProducer
 EDProducer ()
 
ModuleDescription const & moduleDescription () const
 
virtual ~EDProducer ()
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
std::vector< edm::ProductResolverIndex > const & indiciesForPutProducts (BranchType iBranchType) const
 
 ProducerBase ()
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription const &)> registrationCallback () const
 used by the fwk to register list of products More...
 
void resolvePutIndicies (BranchType iBranchType, std::unordered_multimap< std::string, edm::ProductResolverIndex > const &iIndicies, std::string const &moduleLabel)
 
virtual ~ProducerBase () noexcept(false)
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

virtual void beginJob ()
 
virtual void endJob ()
 
std::pair< int, float > findJpsiMCInfo (reco::GenParticleRef genJpsi)
 
bool isAbHadron (int pdgID)
 
bool isAMixedbHadron (int pdgID, int momPdgID)
 
virtual void produce (edm::Event &, const edm::EventSetup &)
 

Private Attributes

bool addCommonVertex_
 
bool addMCTruth_
 
bool addMuonlessPrimaryVertex_
 
StringCutObjectSelector< reco::Candidate, true > dimuonSelection_
 
StringCutObjectSelector< pat::MuonhigherPuritySelection_
 
StringCutObjectSelector< pat::MuonlowerPuritySelection_
 
InvariantMassFromVertex massCalculator
 
edm::EDGetTokenT< edm::View< pat::Muon > > muons_
 
bool resolveAmbiguity_
 
edm::EDGetTokenT< reco::BeamSpotrevtxbs_
 
edm::EDGetTokenT< reco::TrackCollectionrevtxtrks_
 
edm::EDGetTokenT< reco::BeamSpotthebeamspot_
 
edm::EDGetTokenT< reco::VertexCollectionthePVs_
 
GreaterByVProb< pat::CompositeCandidatevPComparator_
 

Additional Inherited Members

- Public Types inherited from edm::EDProducer
typedef EDProducer ModuleType
 
- Public Types inherited from edm::ProducerBase
typedef ProductRegistryHelper::TypeLabelList TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::EDProducer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

Definition at line 37 of file Onia2MuMuPAT.h.

Constructor & Destructor Documentation

Onia2MuMuPAT::Onia2MuMuPAT ( const edm::ParameterSet iConfig)
explicit

Definition at line 29 of file Onia2MuMuPAT.cc.

References revtxbs_, and revtxtrks_.

29  :
31  thebeamspot_(consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("beamSpotTag"))),
32  thePVs_(consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("primaryVertexTag"))),
33  higherPuritySelection_(iConfig.getParameter<std::string>("higherPuritySelection")),
34  lowerPuritySelection_(iConfig.getParameter<std::string>("lowerPuritySelection")),
35  dimuonSelection_(iConfig.existsAs<std::string>("dimuonSelection") ? iConfig.getParameter<std::string>("dimuonSelection") : ""),
36  addCommonVertex_(iConfig.getParameter<bool>("addCommonVertex")),
37  addMuonlessPrimaryVertex_(iConfig.getParameter<bool>("addMuonlessPrimaryVertex")),
38  resolveAmbiguity_(iConfig.getParameter<bool>("resolvePileUpAmbiguity")),
39  addMCTruth_(iConfig.getParameter<bool>("addMCTruth"))
40 {
41  revtxtrks_ = consumes<reco::TrackCollection>((edm::InputTag)"generalTracks"); //if that is not true, we will raise an exception
42  revtxbs_ = consumes<reco::BeamSpot>((edm::InputTag)"offlineBeamSpot");
43  produces<pat::CompositeCandidateCollection>();
44 }
T getParameter(std::string const &) const
bool existsAs(std::string const &parameterName, bool trackiness=true) const
checks if a parameter exists as a given type
Definition: ParameterSet.h:186
bool addCommonVertex_
Definition: Onia2MuMuPAT.h:61
StringCutObjectSelector< pat::Muon > higherPuritySelection_
Definition: Onia2MuMuPAT.h:58
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
edm::EDGetTokenT< reco::BeamSpot > revtxbs_
Definition: Onia2MuMuPAT.h:57
StringCutObjectSelector< reco::Candidate, true > dimuonSelection_
Definition: Onia2MuMuPAT.h:60
edm::EDGetTokenT< reco::TrackCollection > revtxtrks_
Definition: Onia2MuMuPAT.h:56
bool addMuonlessPrimaryVertex_
Definition: Onia2MuMuPAT.h:61
bool resolveAmbiguity_
Definition: Onia2MuMuPAT.h:62
StringCutObjectSelector< pat::Muon > lowerPuritySelection_
Definition: Onia2MuMuPAT.h:59
bool addMCTruth_
Definition: Onia2MuMuPAT.h:63
edm::EDGetTokenT< reco::BeamSpot > thebeamspot_
Definition: Onia2MuMuPAT.h:54
edm::EDGetTokenT< reco::VertexCollection > thePVs_
Definition: Onia2MuMuPAT.h:55
edm::EDGetTokenT< edm::View< pat::Muon > > muons_
Definition: Onia2MuMuPAT.h:53
Onia2MuMuPAT::~Onia2MuMuPAT ( )

Definition at line 47 of file Onia2MuMuPAT.cc.

48 {
49 
50  // do anything here that needs to be done at desctruction time
51  // (e.g. close files, deallocate resources etc.)
52 
53 }

Member Function Documentation

void Onia2MuMuPAT::beginJob ( void  )
privatevirtual

Reimplemented from edm::EDProducer.

Definition at line 487 of file Onia2MuMuPAT.cc.

488 {
489 }
void Onia2MuMuPAT::endJob ( void  )
privatevirtual

Reimplemented from edm::EDProducer.

Definition at line 493 of file Onia2MuMuPAT.cc.

References DEFINE_FWK_MODULE.

493  {
494 }
std::pair< int, float > Onia2MuMuPAT::findJpsiMCInfo ( reco::GenParticleRef  genJpsi)
private

Definition at line 427 of file Onia2MuMuPAT.cc.

References isAbHadron(), isAMixedbHadron(), edm::Ref< C, T, F >::isNonnull(), edm::Ref< C, T, F >::isNull(), and mps_fire::result.

Referenced by produce().

427  {
428 
429  int momJpsiID = 0;
430  float trueLife = -99.;
431 
432  if (genJpsi->numberOfMothers()>0) {
433 
434  TVector3 trueVtx(0.0,0.0,0.0);
435  TVector3 trueP(0.0,0.0,0.0);
436  TVector3 trueVtxMom(0.0,0.0,0.0);
437 
438  trueVtx.SetXYZ(genJpsi->vertex().x(),genJpsi->vertex().y(),genJpsi->vertex().z());
439  trueP.SetXYZ(genJpsi->momentum().x(),genJpsi->momentum().y(),genJpsi->momentum().z());
440 
441  bool aBhadron = false;
442  reco::GenParticleRef Jpsimom = genJpsi->motherRef(); // find mothers
443  if (Jpsimom.isNull()) {
444  std::pair<int, float> result = std::make_pair(momJpsiID, trueLife);
445  return result;
446  } else {
447  reco::GenParticleRef Jpsigrandmom = Jpsimom->motherRef();
448  if (isAbHadron(Jpsimom->pdgId())) {
449  if (Jpsigrandmom.isNonnull() && isAMixedbHadron(Jpsimom->pdgId(),Jpsigrandmom->pdgId())) {
450  momJpsiID = Jpsigrandmom->pdgId();
451  trueVtxMom.SetXYZ(Jpsigrandmom->vertex().x(),Jpsigrandmom->vertex().y(),Jpsigrandmom->vertex().z());
452  } else {
453  momJpsiID = Jpsimom->pdgId();
454  trueVtxMom.SetXYZ(Jpsimom->vertex().x(),Jpsimom->vertex().y(),Jpsimom->vertex().z());
455  }
456  aBhadron = true;
457  } else {
458  if (Jpsigrandmom.isNonnull() && isAbHadron(Jpsigrandmom->pdgId())) {
459  reco::GenParticleRef JpsiGrandgrandmom = Jpsigrandmom->motherRef();
460  if (JpsiGrandgrandmom.isNonnull() && isAMixedbHadron(Jpsigrandmom->pdgId(),JpsiGrandgrandmom->pdgId())) {
461  momJpsiID = JpsiGrandgrandmom->pdgId();
462  trueVtxMom.SetXYZ(JpsiGrandgrandmom->vertex().x(),JpsiGrandgrandmom->vertex().y(),JpsiGrandgrandmom->vertex().z());
463  } else {
464  momJpsiID = Jpsigrandmom->pdgId();
465  trueVtxMom.SetXYZ(Jpsigrandmom->vertex().x(),Jpsigrandmom->vertex().y(),Jpsigrandmom->vertex().z());
466  }
467  aBhadron = true;
468  }
469  }
470  if (!aBhadron) {
471  momJpsiID = Jpsimom->pdgId();
472  trueVtxMom.SetXYZ(Jpsimom->vertex().x(),Jpsimom->vertex().y(),Jpsimom->vertex().z());
473  }
474  }
475 
476  TVector3 vdiff = trueVtx - trueVtxMom;
477  //trueLife = vdiff.Perp()*3.09688/trueP.Perp();
478  trueLife = vdiff.Perp()*genJpsi->mass()/trueP.Perp();
479  }
480  std::pair<int, float> result = std::make_pair(momJpsiID, trueLife);
481  return result;
482 
483 }
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:252
bool isAbHadron(int pdgID)
bool isNull() const
Checks for null.
Definition: Ref.h:249
bool isAMixedbHadron(int pdgID, int momPdgID)
bool Onia2MuMuPAT::isAbHadron ( int  pdgID)
private

Definition at line 409 of file Onia2MuMuPAT.cc.

References funct::abs().

Referenced by findJpsiMCInfo().

409  {
410 
411  if (abs(pdgID) == 511 || abs(pdgID) == 521 || abs(pdgID) == 531 || abs(pdgID) == 5122) return true;
412  return false;
413 
414 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool Onia2MuMuPAT::isAMixedbHadron ( int  pdgID,
int  momPdgID 
)
private

Definition at line 417 of file Onia2MuMuPAT.cc.

References funct::abs().

Referenced by findJpsiMCInfo().

417  {
418 
419  if ((abs(pdgID) == 511 && abs(momPdgID) == 511 && pdgID*momPdgID < 0) ||
420  (abs(pdgID) == 531 && abs(momPdgID) == 531 && pdgID*momPdgID < 0))
421  return true;
422  return false;
423 
424 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
void Onia2MuMuPAT::produce ( edm::Event iEvent,
const edm::EventSetup iSetup 
)
privatevirtual

DCA

end DCA

Definition at line 62 of file Onia2MuMuPAT.cc.

References IPTools::absoluteImpactParameter3D(), addCommonVertex_, reco::CompositeCandidate::addDaughter(), addMCTruth_, addMuonlessPrimaryVertex_, pat::PATObject< ObjectType >::addUserData(), pat::PATObject< ObjectType >::addUserFloat(), pat::PATObject< ObjectType >::addUserInt(), align::BeamSpot, TransientTrackBuilder::build(), TwoTrackMinimumDistance::calculate(), ClosestApproachInRPhi::calculate(), gather_cfg::cout, TransientVertex::degreesOfFreedom(), dimuonSelection_, VertexDistanceXY::distance(), ClosestApproachInRPhi::distance(), pat::PATObject< ObjectType >::embedGenParticle(), Measurement1D::error(), relativeConstraints::error, reco::Vertex::error(), cppFunctionSkipper::exception, findJpsiMCInfo(), edm::EventSetup::get(), edm::Event::getByToken(), reco::Vertex::hasRefittedTracks(), higherPuritySelection_, reco::TrackBase::highPurity, edm::HandleBase::id(), edm::Ref< C, T, F >::id(), reco::Muon::innerTrack(), InvariantMassFromVertex::invariantMass(), edm::Ref< C, T, F >::isNonnull(), edm::HandleBase::isValid(), TrajectoryStateClosestToPoint::isValid(), TransientVertex::isValid(), jpsi, edm::RefToBase< T >::key(), edm::Ref< C, T, F >::key(), lowerPuritySelection_, seedCreatorFromRegionConsecutiveHitsEDProducer_cff::magneticField, OniaVtxReProducer::makeVertices(), reco::LeafCandidate::mass(), massCalculator, GlobalErrorBase< T, ErrorWeightType >::matrix(), Validation_cff::minDz, eostools::move(), electronCleaner_cfi::muons, muons_, MagneticField::nominalValue(), reco::Vertex::originalTrack(), reco::Candidate::pdgId(), TwoTrackMinimumDistance::points(), reco::Vertex::position(), TransientVertex::position(), edm::ESHandle< T >::product(), reco::TrackBase::pt(), reco::TrackBase::ptError(), edm::Event::put(), reco::LeafCandidate::px(), reco::LeafCandidate::py(), reco::LeafCandidate::pz(), reco::TrackBase::quality(), reco::Vertex::refittedTracks(), resolveAmbiguity_, reco::return(), revtxbs_, revtxtrks_, reco::LeafCandidate::setCharge(), pat::PATObject< ObjectType >::setGenParticleRef(), reco::LeafCandidate::setP4(), mathSSE::sqrt(), ClosestApproachInRPhi::status(), mps_update::status, thebeamspot_, thePVs_, TrajectoryStateClosestToPoint::theState(), TransientVertex::totalChiSquared(), HiIsolationCommonParameters_cff::track, reco::Muon::track(), reco::Vertex::tracks_begin(), reco::Vertex::tracks_end(), reco::Vertex::tracksSize(), reco::Vertex::trackWeight(), groupFilesInBlocks::tt, Measurement1D::value(), KalmanVertexFitter::vertex(), vPComparator_, badGlobalMuonTaggersAOD_cff::vtx, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

63 {
64  using namespace edm;
65  using namespace std;
66  using namespace reco;
68 
69  vector<double> muMasses;
70  muMasses.push_back( 0.1056583715 );
71  muMasses.push_back( 0.1056583715 );
72 
73  std::unique_ptr<pat::CompositeCandidateCollection> oniaOutput(new pat::CompositeCandidateCollection);
74 
75  Vertex thePrimaryV;
76  Vertex theBeamSpotV;
77 
79  iSetup.get<IdealMagneticFieldRecord>().get(magneticField);
80  const MagneticField* field = magneticField.product();
81 
82  Handle<BeamSpot> theBeamSpot;
83  iEvent.getByToken(thebeamspot_,theBeamSpot);
84  BeamSpot bs = *theBeamSpot;
85  theBeamSpotV = Vertex(bs.position(), bs.covariance3D());
86 
88  iEvent.getByToken(thePVs_, priVtxs);
89  if ( priVtxs->begin() != priVtxs->end() ) {
90  thePrimaryV = Vertex(*(priVtxs->begin()));
91  }
92  else {
93  thePrimaryV = Vertex(bs.position(), bs.covariance3D());
94  }
95 
97  iEvent.getByToken(muons_,muons);
98 
100  iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder",theTTBuilder);
101  KalmanVertexFitter vtxFitter(true);
102  TrackCollection muonLess;
103 
104  // JPsi candidates only from muons
105  for(View<pat::Muon>::const_iterator it = muons->begin(), itend = muons->end(); it != itend; ++it){
106  // both must pass low quality
107  if(!lowerPuritySelection_(*it)) continue;
108  for(View<pat::Muon>::const_iterator it2 = it+1; it2 != itend;++it2){
109  // both must pass low quality
110  if(!lowerPuritySelection_(*it2)) continue;
111  // one must pass tight quality
112  if (!(higherPuritySelection_(*it) || higherPuritySelection_(*it2))) continue;
113 
115  vector<TransientVertex> pvs;
116 
117  // ---- no explicit order defined ----
118  myCand.addDaughter(*it, "muon1");
119  myCand.addDaughter(*it2,"muon2");
120 
121  // ---- define and set candidate's 4momentum ----
122  LorentzVector jpsi = it->p4() + it2->p4();
123  myCand.setP4(jpsi);
124  myCand.setCharge(it->charge()+it2->charge());
125 
126  // ---- apply the dimuon cut ----
127  if(!dimuonSelection_(myCand)) continue;
128 
129  // ---- fit vertex using Tracker tracks (if they have tracks) ----
130  if (it->track().isNonnull() && it2->track().isNonnull()) {
131 
132  //build the dimuon secondary vertex
133  vector<TransientTrack> t_tks;
134  t_tks.push_back(theTTBuilder->build(*it->track())); // pass the reco::Track, not the reco::TrackRef (which can be transient)
135  t_tks.push_back(theTTBuilder->build(*it2->track())); // otherwise the vertex will have transient refs inside.
136  TransientVertex myVertex = vtxFitter.vertex(t_tks);
137 
138  CachingVertex<5> VtxForInvMass = vtxFitter.vertex( t_tks );
139 
140  Measurement1D MassWErr(jpsi.M(),-9999.);
141  if ( field->nominalValue() > 0 ) {
142  MassWErr = massCalculator.invariantMass( VtxForInvMass, muMasses );
143  } else {
144  myVertex = TransientVertex(); // with no arguments it is invalid
145  }
146 
147  myCand.addUserFloat("MassErr",MassWErr.error());
148 
149  if (myVertex.isValid()) {
150  float vChi2 = myVertex.totalChiSquared();
151  float vNDF = myVertex.degreesOfFreedom();
152  float vProb(TMath::Prob(vChi2,(int)vNDF));
153 
154  myCand.addUserFloat("vNChi2",vChi2/vNDF);
155  myCand.addUserFloat("vProb",vProb);
156 
157  TVector3 vtx;
158  TVector3 pvtx;
159  VertexDistanceXY vdistXY;
160 
161  vtx.SetXYZ(myVertex.position().x(),myVertex.position().y(),0);
162  TVector3 pperp(jpsi.px(), jpsi.py(), 0);
163  AlgebraicVector3 vpperp(pperp.x(),pperp.y(),0);
164 
165  if (resolveAmbiguity_) {
166 
167  float minDz = 999999.;
170  GlobalPoint(myVertex.position().x(), myVertex.position().y(), myVertex.position().z()),
171  GlobalVector(myCand.px(),myCand.py(),myCand.pz()),TrackCharge(0),&(*magneticField)),
173  GlobalPoint(bs.position().x(), bs.position().y(), bs.position().z()),
174  GlobalVector(bs.dxdz(), bs.dydz(), 1.),TrackCharge(0),&(*magneticField)));
175  float extrapZ=-9E20;
176  if (status) extrapZ=ttmd.points().first.z();
177 
178  for(VertexCollection::const_iterator itv = priVtxs->begin(), itvend = priVtxs->end(); itv != itvend; ++itv){
179  float deltaZ = fabs(extrapZ - itv->position().z()) ;
180  if ( deltaZ < minDz ) {
181  minDz = deltaZ;
182  thePrimaryV = Vertex(*itv);
183  }
184  }
185  }
186 
187  Vertex theOriginalPV = thePrimaryV;
188 
189  muonLess.clear();
190  muonLess.reserve(thePrimaryV.tracksSize());
191  if( addMuonlessPrimaryVertex_ && thePrimaryV.tracksSize()>2) {
192  // Primary vertex matched to the dimuon, now refit it removing the two muons
193  OniaVtxReProducer revertex(priVtxs, iEvent);
195  iEvent.getByToken(revtxtrks_, pvtracks);
196  if( !pvtracks.isValid()) { std::cout << "pvtracks NOT valid " << std::endl; }
197  else {
198  edm::Handle<reco::BeamSpot> pvbeamspot;
199  iEvent.getByToken(revtxbs_, pvbeamspot);
200  if (pvbeamspot.id() != theBeamSpot.id()) edm::LogWarning("Inconsistency") << "The BeamSpot used for PV reco is not the same used in this analyzer.";
201  // I need to go back to the reco::Muon object, as the TrackRef in the pat::Muon can be an embedded ref.
202  const reco::Muon *rmu1 = dynamic_cast<const reco::Muon *>(it->originalObject());
203  const reco::Muon *rmu2 = dynamic_cast<const reco::Muon *>(it2->originalObject());
204  // check that muons are truly from reco::Muons (and not, e.g., from PF Muons)
205  // also check that the tracks really come from the track collection used for the BS
206  if (rmu1 != 0 && rmu2 != 0 && rmu1->track().id() == pvtracks.id() && rmu2->track().id() == pvtracks.id()) {
207  // Save the keys of the tracks in the primary vertex
208  // std::vector<size_t> vertexTracksKeys;
209  // vertexTracksKeys.reserve(thePrimaryV.tracksSize());
210  if( thePrimaryV.hasRefittedTracks() ) {
211  // Need to go back to the original tracks before taking the key
212  std::vector<reco::Track>::const_iterator itRefittedTrack = thePrimaryV.refittedTracks().begin();
213  std::vector<reco::Track>::const_iterator refittedTracksEnd = thePrimaryV.refittedTracks().end();
214  for( ; itRefittedTrack != refittedTracksEnd; ++itRefittedTrack ) {
215  if( thePrimaryV.originalTrack(*itRefittedTrack).key() == rmu1->track().key() ) continue;
216  if( thePrimaryV.originalTrack(*itRefittedTrack).key() == rmu2->track().key() ) continue;
217  // vertexTracksKeys.push_back(thePrimaryV.originalTrack(*itRefittedTrack).key());
218  muonLess.push_back(*(thePrimaryV.originalTrack(*itRefittedTrack)));
219  }
220  }
221  else {
222  std::vector<reco::TrackBaseRef>::const_iterator itPVtrack = thePrimaryV.tracks_begin();
223  for( ; itPVtrack != thePrimaryV.tracks_end(); ++itPVtrack ) if (itPVtrack->isNonnull()) {
224  if( itPVtrack->key() == rmu1->track().key() ) continue;
225  if( itPVtrack->key() == rmu2->track().key() ) continue;
226  // vertexTracksKeys.push_back(itPVtrack->key());
227  muonLess.push_back(**itPVtrack);
228  }
229  }
230  if (muonLess.size()>1 && muonLess.size() < thePrimaryV.tracksSize()){
231  pvs = revertex.makeVertices(muonLess, *pvbeamspot, iSetup) ;
232  if (!pvs.empty()) {
233  Vertex muonLessPV = Vertex(pvs.front());
234  thePrimaryV = muonLessPV;
235  }
236  }
237  }
238  }
239  }
240 
241  // count the number of high Purity tracks with pT > 900 MeV attached to the chosen vertex
242  double vertexWeight = -1., sumPTPV = -1.;
243  int countTksOfPV = -1;
244  const reco::Muon *rmu1 = dynamic_cast<const reco::Muon *>(it->originalObject());
245  const reco::Muon *rmu2 = dynamic_cast<const reco::Muon *>(it2->originalObject());
246  try{
247  for(reco::Vertex::trackRef_iterator itVtx = theOriginalPV.tracks_begin(); itVtx != theOriginalPV.tracks_end(); itVtx++) if(itVtx->isNonnull()){
248  const reco::Track& track = **itVtx;
249  if(!track.quality(reco::TrackBase::highPurity)) continue;
250  if(track.pt() < 0.5) continue; //reject all rejects from counting if less than 900 MeV
251  TransientTrack tt = theTTBuilder->build(track);
252  pair<bool,Measurement1D> tkPVdist = IPTools::absoluteImpactParameter3D(tt,thePrimaryV);
253  if (!tkPVdist.first) continue;
254  if (tkPVdist.second.significance()>3) continue;
255  if (track.ptError()/track.pt()>0.1) continue;
256  // do not count the two muons
257  if (rmu1 != 0 && rmu1->innerTrack().key() == itVtx->key())
258  continue;
259  if (rmu2 != 0 && rmu2->innerTrack().key() == itVtx->key())
260  continue;
261 
262  vertexWeight += theOriginalPV.trackWeight(*itVtx);
263  if(theOriginalPV.trackWeight(*itVtx) > 0.5){
264  countTksOfPV++;
265  sumPTPV += track.pt();
266  }
267  }
268  } catch (std::exception & err) {std::cout << " muon Selection%G�%@failed " << std::endl; return ; }
269 
270  myCand.addUserInt("countTksOfPV", countTksOfPV);
271  myCand.addUserFloat("vertexWeight", (float) vertexWeight);
272  myCand.addUserFloat("sumPTPV", (float) sumPTPV);
273 
275  TrajectoryStateClosestToPoint mu1TS = t_tks[0].impactPointTSCP();
276  TrajectoryStateClosestToPoint mu2TS = t_tks[1].impactPointTSCP();
277  float dca = 1E20;
278  if (mu1TS.isValid() && mu2TS.isValid()) {
280  cApp.calculate(mu1TS.theState(), mu2TS.theState());
281  if (cApp.status() ) dca = cApp.distance();
282  }
283  myCand.addUserFloat("DCA", dca );
285 
287  myCand.addUserData("muonlessPV",Vertex(thePrimaryV));
288  else
289  myCand.addUserData("PVwithmuons",thePrimaryV);
290 
291  // lifetime using PV
292  pvtx.SetXYZ(thePrimaryV.position().x(),thePrimaryV.position().y(),0);
293  TVector3 vdiff = vtx - pvtx;
294  double cosAlpha = vdiff.Dot(pperp)/(vdiff.Perp()*pperp.Perp());
295  Measurement1D distXY = vdistXY.distance(Vertex(myVertex), thePrimaryV);
296  //double ctauPV = distXY.value()*cosAlpha*3.09688/pperp.Perp();
297  double ctauPV = distXY.value()*cosAlpha * myCand.mass()/pperp.Perp();
298  GlobalError v1e = (Vertex(myVertex)).error();
299  GlobalError v2e = thePrimaryV.error();
300  AlgebraicSymMatrix33 vXYe = v1e.matrix()+ v2e.matrix();
301  //double ctauErrPV = sqrt(vXYe.similarity(vpperp))*3.09688/(pperp.Perp2());
302  double ctauErrPV = sqrt(ROOT::Math::Similarity(vpperp,vXYe))*myCand.mass()/(pperp.Perp2());
303 
304  myCand.addUserFloat("ppdlPV",ctauPV);
305  myCand.addUserFloat("ppdlErrPV",ctauErrPV);
306  myCand.addUserFloat("cosAlpha",cosAlpha);
307 
308  // lifetime using BS
309  pvtx.SetXYZ(theBeamSpotV.position().x(),theBeamSpotV.position().y(),0);
310  vdiff = vtx - pvtx;
311  cosAlpha = vdiff.Dot(pperp)/(vdiff.Perp()*pperp.Perp());
312  distXY = vdistXY.distance(Vertex(myVertex), theBeamSpotV);
313  //double ctauBS = distXY.value()*cosAlpha*3.09688/pperp.Perp();
314  double ctauBS = distXY.value()*cosAlpha*myCand.mass()/pperp.Perp();
315  GlobalError v1eB = (Vertex(myVertex)).error();
316  GlobalError v2eB = theBeamSpotV.error();
317  AlgebraicSymMatrix33 vXYeB = v1eB.matrix()+ v2eB.matrix();
318  //double ctauErrBS = sqrt(vXYeB.similarity(vpperp))*3.09688/(pperp.Perp2());
319  double ctauErrBS = sqrt(ROOT::Math::Similarity(vpperp,vXYeB))*myCand.mass()/(pperp.Perp2());
320 
321  myCand.addUserFloat("ppdlBS",ctauBS);
322  myCand.addUserFloat("ppdlErrBS",ctauErrBS);
323 
324  if (addCommonVertex_) {
325  myCand.addUserData("commonVertex",Vertex(myVertex));
326  }
327  } else {
328  myCand.addUserFloat("vNChi2",-1);
329  myCand.addUserFloat("vProb", -1);
330  myCand.addUserFloat("ppdlPV",-100);
331  myCand.addUserFloat("ppdlErrPV",-100);
332  myCand.addUserFloat("cosAlpha",-100);
333  myCand.addUserFloat("ppdlBS",-100);
334  myCand.addUserFloat("ppdlErrBS",-100);
335  myCand.addUserFloat("DCA", -1 );
336  if (addCommonVertex_) {
337  myCand.addUserData("commonVertex",Vertex());
338  }
340  myCand.addUserData("muonlessPV",Vertex());
341  } else {
342  myCand.addUserData("PVwithmuons",Vertex());
343  }
344  }
345  }
346 
347  // ---- MC Truth, if enabled ----
348  if (addMCTruth_) {
349  reco::GenParticleRef genMu1 = it->genParticleRef();
350  reco::GenParticleRef genMu2 = it2->genParticleRef();
351  if (genMu1.isNonnull() && genMu2.isNonnull()) {
352  if (genMu1->numberOfMothers()>0 && genMu2->numberOfMothers()>0){
353  reco::GenParticleRef mom1 = genMu1->motherRef();
354  reco::GenParticleRef mom2 = genMu2->motherRef();
355  if (mom1.isNonnull() && (mom1 == mom2)) {
356  myCand.setGenParticleRef(mom1); // set
357  myCand.embedGenParticle(); // and embed
358  std::pair<int, float> MCinfo = findJpsiMCInfo(mom1);
359  myCand.addUserInt("momPDGId",MCinfo.first);
360  myCand.addUserFloat("ppdlTrue",MCinfo.second);
361  } else {
362  myCand.addUserInt("momPDGId",0);
363  myCand.addUserFloat("ppdlTrue",-99.);
364  }
365  } else {
367  edm::EDGetTokenT<reco::GenParticleCollection> genCands_ = consumes<reco::GenParticleCollection>((edm::InputTag)"genParticles");
368  iEvent.getByToken(genCands_, theGenParticles);
369  if (theGenParticles.isValid()){
370  for(size_t iGenParticle=0; iGenParticle<theGenParticles->size();++iGenParticle) {
371  const Candidate & genCand = (*theGenParticles)[iGenParticle];
372  if (genCand.pdgId()==443 || genCand.pdgId()==100443 ||
373  genCand.pdgId()==553 || genCand.pdgId()==100553 || genCand.pdgId()==200553) {
374  reco::GenParticleRef mom1(theGenParticles,iGenParticle);
375  myCand.setGenParticleRef(mom1);
376  myCand.embedGenParticle();
377  std::pair<int, float> MCinfo = findJpsiMCInfo(mom1);
378  myCand.addUserInt("momPDGId",MCinfo.first);
379  myCand.addUserFloat("ppdlTrue",MCinfo.second);
380  }
381  }
382  } else {
383  myCand.addUserInt("momPDGId",0);
384  myCand.addUserFloat("ppdlTrue",-99.);
385  }
386  }
387  } else {
388  myCand.addUserInt("momPDGId",0);
389  myCand.addUserFloat("ppdlTrue",-99.);
390  }
391  }
392 
393 
394 
395 
396  // ---- Push back output ----
397  oniaOutput->push_back(myCand);
398  }
399  }
400 
401  std::sort(oniaOutput->begin(),oniaOutput->end(),vPComparator_);
402 
403  iEvent.put(std::move(oniaOutput));
404 
405 }
Analysis-level particle class.
virtual bool calculate(const TrajectoryStateOnSurface &sta, const TrajectoryStateOnSurface &stb)
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:122
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:252
virtual double mass() const final
mass
Covariance3DMatrix covariance3D() const
return only 3D position covariance matrix
Definition: BeamSpot.h:118
bool addCommonVertex_
Definition: Onia2MuMuPAT.h:61
trackRef_iterator tracks_end() const
last iterator over tracks
Definition: Vertex.cc:81
StringCutObjectSelector< pat::Muon > higherPuritySelection_
Definition: Onia2MuMuPAT.h:58
ProductID id() const
Definition: HandleBase.cc:15
virtual TrackRef innerTrack() const
Definition: Muon.h:48
int nominalValue() const
The nominal field value for this map in kGauss.
Definition: MagneticField.h:56
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:460
TrackBaseRef originalTrack(const Track &refTrack) const
Definition: Vertex.cc:108
const FreeTrajectoryState & theState() const
const AlgebraicSymMatrix33 matrix() const
bool hasRefittedTracks() const
Checks whether refitted tracks are stored.
Definition: Vertex.h:149
float totalChiSquared() const
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:14
reco::TransientTrack build(const reco::Track *p) const
std::pair< bool, Measurement1D > absoluteImpactParameter3D(const reco::TransientTrack &transientTrack, const reco::Vertex &vertex)
Definition: IPTools.cc:37
T y() const
Definition: PV3DBase.h:63
std::pair< int, float > findJpsiMCInfo(reco::GenParticleRef genJpsi)
virtual TrackRef track() const
reference to a Track
Definition: Muon.h:49
InvariantMassFromVertex massCalculator
Definition: Onia2MuMuPAT.h:66
void addUserFloat(const std::string &label, float data, const bool overwrite=false)
Set user-defined float.
Definition: PATObject.h:813
key_type key() const
Accessor for product key.
Definition: Ref.h:264
const std::vector< Track > & refittedTracks() const
Returns the container of refitted tracks.
Definition: Vertex.h:164
const Point & position() const
position
Definition: Vertex.h:109
void setGenParticleRef(const reco::GenParticleRef &ref, bool embed=false)
Set the generator level particle reference.
Definition: PATObject.h:671
ProductID id() const
Accessor for product ID.
Definition: Ref.h:258
double dydz() const
dydz slope
Definition: BeamSpot.h:84
int TrackCharge
Definition: TrackCharge.h:4
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
edm::EDGetTokenT< reco::BeamSpot > revtxbs_
Definition: Onia2MuMuPAT.h:57
float degreesOfFreedom() const
size_t key() const
Definition: RefToBase.h:250
void embedGenParticle()
Definition: PATObject.h:694
virtual int pdgId() const =0
PDG identifier.
virtual double px() const final
x coordinate of momentum vector
T sqrt(T t)
Definition: SSEVec.h:18
GlobalPoint position() const
double pt() const
track transverse momentum
Definition: TrackBase.h:616
T z() const
Definition: PV3DBase.h:64
double ptError() const
error on Pt (set to 1000 TeV if charge==0 for safety)
Definition: TrackBase.h:758
virtual void setCharge(Charge q) final
set electric charge
Definition: LeafCandidate.h:93
StringCutObjectSelector< reco::Candidate, true > dimuonSelection_
Definition: Onia2MuMuPAT.h:60
ROOT::Math::SVector< double, 3 > AlgebraicVector3
edm::EDGetTokenT< reco::TrackCollection > revtxtrks_
Definition: Onia2MuMuPAT.h:56
float trackWeight(const TREF &r) const
returns the weight with which a Track has contributed to the vertex-fit.
Definition: Vertex.h:81
bool addMuonlessPrimaryVertex_
Definition: Onia2MuMuPAT.h:61
virtual Measurement1D distance(const GlobalPoint &vtx1Position, const GlobalError &vtx1PositionError, const GlobalPoint &vtx2Position, const GlobalError &vtx2PositionError) const
bool isValid() const
Definition: HandleBase.h:74
void addUserInt(const std::string &label, int32_t data, const bool overwrite=false)
Set user-defined int.
Definition: PATObject.h:855
void addDaughter(const Candidate &, const std::string &s="")
add a clone of the passed candidate as daughter
double dxdz() const
dxdz slope
Definition: BeamSpot.h:82
virtual double pz() const final
z coordinate of momentum vector
bool resolveAmbiguity_
Definition: Onia2MuMuPAT.h:62
StringCutObjectSelector< pat::Muon > lowerPuritySelection_
Definition: Onia2MuMuPAT.h:59
bool addMCTruth_
Definition: Onia2MuMuPAT.h:63
virtual bool calculate(const TrajectoryStateOnSurface &sta, const TrajectoryStateOnSurface &stb)
std::vector< CompositeCandidate > CompositeCandidateCollection
const T & get() const
Definition: EventSetup.h:56
virtual void setP4(const LorentzVector &p4) final
set 4-momentum
double value() const
Definition: Measurement1D.h:28
Error error() const
return SMatrix
Definition: Vertex.h:139
edm::EDGetTokenT< reco::BeamSpot > thebeamspot_
Definition: Onia2MuMuPAT.h:54
bool quality(const TrackQuality) const
Track quality.
Definition: TrackBase.h:505
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:37
return(e1-e2)*(e1-e2)+dp *dp
#define jpsi
fixed size matrix
HLT enums.
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Definition: View.h:86
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector<TrackRef>
Definition: Vertex.h:37
edm::EDGetTokenT< reco::VertexCollection > thePVs_
Definition: Onia2MuMuPAT.h:55
const Point & position() const
position
Definition: BeamSpot.h:62
GreaterByVProb< pat::CompositeCandidate > vPComparator_
Definition: Onia2MuMuPAT.h:64
trackRef_iterator tracks_begin() const
first iterator over tracks
Definition: Vertex.cc:76
virtual float distance() const
T x() const
Definition: PV3DBase.h:62
T const * product() const
Definition: ESHandle.h:86
bool isValid() const
virtual std::pair< GlobalPoint, GlobalPoint > points() const
virtual double py() const final
y coordinate of momentum vector
void addUserData(const std::string &label, const T &data, bool transientOnly=false, bool overwrite=false)
Definition: PATObject.h:309
def move(src, dest)
Definition: eostools.py:510
size_t tracksSize() const
number of tracks
Definition: Vertex.cc:71
virtual bool status() const
Measurement1D invariantMass(const CachingVertex< 5 > &vertex, const std::vector< double > &masses) const
edm::EDGetTokenT< edm::View< pat::Muon > > muons_
Definition: Onia2MuMuPAT.h:53
math::PtEtaPhiELorentzVectorF LorentzVector
Global3DVector GlobalVector
Definition: GlobalVector.h:10

Member Data Documentation

bool Onia2MuMuPAT::addCommonVertex_
private

Definition at line 61 of file Onia2MuMuPAT.h.

Referenced by produce().

bool Onia2MuMuPAT::addMCTruth_
private

Definition at line 63 of file Onia2MuMuPAT.h.

Referenced by produce().

bool Onia2MuMuPAT::addMuonlessPrimaryVertex_
private

Definition at line 61 of file Onia2MuMuPAT.h.

Referenced by produce().

StringCutObjectSelector<reco::Candidate, true> Onia2MuMuPAT::dimuonSelection_
private

Definition at line 60 of file Onia2MuMuPAT.h.

Referenced by produce().

StringCutObjectSelector<pat::Muon> Onia2MuMuPAT::higherPuritySelection_
private

Definition at line 58 of file Onia2MuMuPAT.h.

Referenced by produce().

StringCutObjectSelector<pat::Muon> Onia2MuMuPAT::lowerPuritySelection_
private

Definition at line 59 of file Onia2MuMuPAT.h.

Referenced by produce().

InvariantMassFromVertex Onia2MuMuPAT::massCalculator
private

Definition at line 66 of file Onia2MuMuPAT.h.

Referenced by produce().

edm::EDGetTokenT<edm::View<pat::Muon> > Onia2MuMuPAT::muons_
private

Definition at line 53 of file Onia2MuMuPAT.h.

Referenced by produce().

bool Onia2MuMuPAT::resolveAmbiguity_
private

Definition at line 62 of file Onia2MuMuPAT.h.

Referenced by produce().

edm::EDGetTokenT<reco::BeamSpot> Onia2MuMuPAT::revtxbs_
private

Definition at line 57 of file Onia2MuMuPAT.h.

Referenced by Onia2MuMuPAT(), and produce().

edm::EDGetTokenT<reco::TrackCollection> Onia2MuMuPAT::revtxtrks_
private

Definition at line 56 of file Onia2MuMuPAT.h.

Referenced by Onia2MuMuPAT(), and produce().

edm::EDGetTokenT<reco::BeamSpot> Onia2MuMuPAT::thebeamspot_
private

Definition at line 54 of file Onia2MuMuPAT.h.

Referenced by produce().

edm::EDGetTokenT<reco::VertexCollection> Onia2MuMuPAT::thePVs_
private

Definition at line 55 of file Onia2MuMuPAT.h.

Referenced by produce().

GreaterByVProb<pat::CompositeCandidate> Onia2MuMuPAT::vPComparator_
private

Definition at line 64 of file Onia2MuMuPAT.h.

Referenced by produce().