22 #include "GaudiKernel/IPartPropSvc.h"
23 #include "AthLinks/ElementLink.h"
24 #include "HepPDT/ParticleDataTable.hh"
55 const HepPDT::ParticleData* pd_el = particleDataTable->particle(MC::ELECTRON);
60 if (
m_doTagAndProbe)
ATH_MSG_WARNING(
"You have requested tag and probe mode. Duplicate mu+trk pairs WILL be allowed, charge ordering WILL NOT be done. Tag track will be first in each candidate");
64 bool illogicalOptions(
false);
66 ATH_MSG_WARNING(
"You are requesting incompatible combinations of muons and tracks in the pairs. JpsiEECandidates will be EMPTY!");
67 illogicalOptions=
true;
70 ATH_MSG_WARNING(
"You are requesting Tag and Probe analysis but have not requested mu+trk mode. This is impossible. JpsiEECandidates will be EMPTY!");
71 illogicalOptions=
true;
74 ATH_MSG_WARNING(
"You are requesting incompatible combinations of charges in the pairs. JpsiEECandidates will be EMPTY!");
75 illogicalOptions=
true;
78 ATH_MSG_WARNING(
"You are requesting same-sign or all-sign combinations in a tag and probe analysis. This doesn't make sense. JpsiEECandidates will be EMPTY!");
79 illogicalOptions=
true;
81 if (illogicalOptions)
return StatusCode::FAILURE;
87 return StatusCode::SUCCESS;
96 m_allElectrons(false),
97 m_useTrackMeasurement(true),
104 m_trkThresholdPt(0.0),
105 m_invMassUpper(100000.0),
107 m_collAngleTheta(0.0),
113 m_electronCollectionKey(
"Electrons"),
114 m_TrkParticleCollection(
"InDetTrackParticles"),
115 m_iVertexFitter(
"Trk::TrkVKalVrtFitter"),
116 m_iV0VertexFitter(
"Trk::V0VertexFitter"),
117 m_trkSelector(
"InDet::TrackSelectorTool"),
118 m_vertexEstimator(
"InDet::VertexPointEstimator"),
120 m_elSelection(
"d0_or_nod0"),
121 m_doTagAndProbe(false)
124 declareInterface<JpsiFinder_ee>(
this);
174 return StatusCode::SUCCESS;;
176 importedElectronCollection = ehandle.
cptr();
186 return StatusCode::SUCCESS;;
188 importedTrackCollection = thandle.
cptr();
192 typedef std::vector<const xAOD::TrackParticle*>
TrackBag;
193 typedef std::vector<const xAOD::Electron*> ElectronBag;
200 for (trkCItr=importedTrackCollection->
begin(); trkCItr!=importedTrackCollection->
end(); ++trkCItr) {
204 theIDTracksAfterSelection.push_back(TP);
206 if (theIDTracksAfterSelection.size() == 0)
return StatusCode::SUCCESS;;
207 ATH_MSG_DEBUG(
"Number of tracks after ID track selection: " << theIDTracksAfterSelection.size());
211 ElectronBag theElectronsAfterSelection;
214 for (elItr=importedElectronCollection->
begin(); elItr!=importedElectronCollection->
end(); ++elItr) {
215 if ( *elItr == NULL )
continue;
216 if (!(*elItr)->trackParticleLink().isValid())
continue;
218 elTrk = (*elItr)->trackParticleLink().cachedElement();
220 if ( elTrk==NULL)
continue;
225 theElectronsAfterSelection.push_back(*elItr);
227 if (theElectronsAfterSelection.size() == 0)
return StatusCode::SUCCESS;;
228 ATH_MSG_DEBUG(
"Number of electrons after selection: " << theElectronsAfterSelection.size());
232 std::vector<JpsiEECandidate> jpsiCandidates;
233 if (
m_elel) jpsiCandidates =
getPairs(theElectronsAfterSelection);
237 ATH_MSG_DEBUG(
"Number of pairs with ee from a B decay: " << jpsiCandidates.size() );
249 for (jpsiItr=jpsiCandidates.begin(); jpsiItr!=jpsiCandidates.end(); ++jpsiItr) {
250 (*jpsiItr).collection1 = importedTrackCollection;
251 (*jpsiItr).collection2 = importedTrackCollection;
256 for (jpsiItr=jpsiCandidates.begin(); jpsiItr!=jpsiCandidates.end(); ++jpsiItr) {
258 (*jpsiItr).trackParticle1 = (*jpsiItr).el1->trackParticleLink().cachedElement();
259 (*jpsiItr).trackParticle2 = (*jpsiItr).el2->trackParticleLink().cachedElement();
260 (*jpsiItr).collection1 = importedTrackCollection;
261 (*jpsiItr).collection2 = importedTrackCollection;
263 ATH_MSG_WARNING(
"Not setup for non-track electron measurements yet....");
272 std::vector<int> listToDelete;
273 std::vector<int>::reverse_iterator ii;
274 for(jpsiItr=jpsiCandidates.begin(); jpsiItr!=jpsiCandidates.end();++jpsiItr,++
index) {
275 if( (fabs((*jpsiItr).trackParticle1->pt()) <
m_higherPt) && (fabs((*jpsiItr).trackParticle2->pt()) <
m_higherPt) ) listToDelete.push_back(
index);
277 for (ii=listToDelete.rbegin(); ii!=listToDelete.rend(); ++ii) {
278 jpsiCandidates.erase(jpsiCandidates.begin() + (*ii) );
280 ATH_MSG_DEBUG(
"Number of candidates after higherPt cut: " << jpsiCandidates.size() );
284 std::vector<JpsiEECandidate> sortedJpsiEECandidates;
288 ATH_MSG_DEBUG(
"Number of candidates after charge selection: " << sortedJpsiEECandidates.size() );
293 std::vector<int> listToDelete;
294 std::vector<int>::reverse_iterator ii;
295 for(jpsiItr=sortedJpsiEECandidates.begin(); jpsiItr!=sortedJpsiEECandidates.end();++jpsiItr,++
index) {
296 double deltatheta = fabs( (*jpsiItr).trackParticle1->theta() - (*jpsiItr).trackParticle2->theta() );
298 double deltaphi = std::abs(
xAOD::P4Helpers::deltaPhi((*jpsiItr).trackParticle1->phi0() , (*jpsiItr).trackParticle2->phi0()));
302 for (ii=listToDelete.rbegin(); ii!=listToDelete.rend(); ++ii) {
303 sortedJpsiEECandidates.erase(sortedJpsiEECandidates.begin() + (*ii) );
305 ATH_MSG_DEBUG(
"Number of collimated candidates: " << sortedJpsiEECandidates.size() );
309 std::vector<double> trkMasses;
314 std::vector<int> listToDelete;
315 std::vector<int>::reverse_iterator ii;
316 for(jpsiItr=sortedJpsiEECandidates.begin(); jpsiItr!=sortedJpsiEECandidates.end(); ++jpsiItr,++
index) {
320 listToDelete.push_back(
index);
323 for (ii=listToDelete.rbegin(); ii!=listToDelete.rend(); ++ii) {
324 sortedJpsiEECandidates.erase(sortedJpsiEECandidates.begin() + (*ii) );
326 ATH_MSG_DEBUG(
"Number of candidates passing invariant mass selection: " << sortedJpsiEECandidates.size() );
329 if (sortedJpsiEECandidates.size() == 0)
return StatusCode::SUCCESS;;
333 for(jpsiItr=sortedJpsiEECandidates.begin(); jpsiItr!=sortedJpsiEECandidates.end(); ++jpsiItr) {
335 std::vector<const xAOD::TrackParticle*> theTracks; theTracks.clear();
336 theTracks.push_back((*jpsiItr).trackParticle1);
337 theTracks.push_back((*jpsiItr).trackParticle2);
338 ATH_MSG_DEBUG(
"theTracks size (should be two!) " << theTracks.size() <<
" being vertexed with tracks " << importedTrackCollection);
339 std::unique_ptr<xAOD::Vertex> myVxCandidate{
fit(theTracks,importedTrackCollection)};
340 if (myVxCandidate != 0) {
342 double chi2 = myVxCandidate->chiSquared();
349 std::vector<const xAOD::Electron*> theStoredElectrons;
350 theStoredElectrons.push_back((*jpsiItr).el1);
351 if (
m_elel) theStoredElectrons.push_back((*jpsiItr).el2);
352 jpsiHelper.
setElectrons(theStoredElectrons,importedElectronCollection);
355 vxContainer.
push_back(std::move(myVxCandidate));
366 return StatusCode::SUCCESS;;
382 if(concreteVertexFitter == 0) {
383 ATH_MSG_FATAL(
"The vertex fitter passed is not a V0 Vertex Fitter");
393 if (errorcode != 0) {startingPoint(0) = 0.0; startingPoint(1) = 0.0; startingPoint(2) = 0.0;}
398 if(myVxCandidate != 0){
399 std::vector<ElementLink<DataVector<xAOD::TrackParticle> > > newLinkVector;
403 newLinkVector.push_back( mylink ); }
411 return myVxCandidate;
414 ATH_MSG_DEBUG(
"Initial fit was a success! " << myVxCandidate);
416 if(myVxCandidate != 0){
417 std::vector<ElementLink<DataVector<xAOD::TrackParticle> > > newLinkVector;
421 newLinkVector.push_back( mylink );
429 return myVxCandidate;
447 std::vector<JpsiEECandidate> myPairs;
449 std::vector<const xAOD::TrackParticle*>::const_iterator outerItr;
450 std::vector<const xAOD::TrackParticle*>::const_iterator innerItr;
452 if(TracksIn.size()>=2){
453 for(outerItr=TracksIn.begin();outerItr<TracksIn.end();++outerItr){
454 for(innerItr=(outerItr+1);innerItr!=TracksIn.end();++innerItr){
458 myPairs.push_back(pair);
474 std::vector<JpsiEECandidate> myPairs;
476 std::vector<const xAOD::Electron*>::const_iterator outerItr;
477 std::vector<const xAOD::Electron*>::const_iterator innerItr;
479 if(electronsIn.size()>=2){
480 for(outerItr=electronsIn.begin();outerItr<electronsIn.end();++outerItr){
481 for(innerItr=(outerItr+1);innerItr!=electronsIn.end();++innerItr){
482 pair.
el1 = *innerItr;
483 pair.
el2 = *outerItr;
485 myPairs.push_back(pair);
501 std::vector<JpsiEECandidate> myPairs;
505 std::vector<const xAOD::TrackParticle*> tracksToKeep;
507 if(tracks.size()>=1 &&
electrons.size()>=1){
509 bool trackIsElectron(
false);
511 if ( ele->trackParticleLink().cachedElement() == trk ) {
512 trackIsElectron=
true;
516 if (!trackIsElectron) tracksToKeep.push_back(trk);
519 }
else {tracksToKeep = tracks;}
521 if(tracksToKeep.size()>=1 &&
electrons.size()>=1){
529 myPairs.push_back(pair);
546 double mass1 = massHypotheses.at(0);
547 double mass2 = massHypotheses.at(1);
557 return (mu1+mu2).M();
569 bool opposite(
false),
same(
false),
all(
false);
570 if (
selection==
"OPPOSITE") opposite=
true;
575 std::vector<JpsiEECandidate> jpsis;
578 for(
auto jpsiItr=jpsisIn.cbegin();jpsiItr!=jpsisIn.cend();jpsiItr++){
579 bool oppCh(
false),sameCh(
false);
582 qOverP2=(*jpsiItr).trackParticle2->qOverP();
583 if(qOverP1*qOverP2<0.0) oppCh=
true;
584 if(qOverP1*qOverP2>0.0) sameCh=
true;
590 tmpJpsi.
el1 = (*jpsiItr).el2;
591 tmpJpsi.
el2 = (*jpsiItr).el1;
595 if (oppCh && (opposite ||
all) ) jpsis.push_back(tmpJpsi);
596 if (sameCh && (
same ||
all) ) jpsis.push_back(tmpJpsi);
613 bool passesLHVLoose = isLHVeryLoose(*
electron);
614 bool passesLHVLoosenod0 = isLHVeryLoosenod0(*
electron);
632 return std::find(theCollection->
begin(), theCollection->
end(), theTrack) != theCollection->
end();
641 double px = 0.,
py = 0.,
pz = 0.;
642 if (0 != vxCandidate) {