22 #include "CLHEP/Random/RandFlat.h"
66 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Init * START *" <<
endmsg;
68 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Init * initializing muon selector *" <<
endmsg;
70 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Init * initializing electron selector *" <<
endmsg;
73 if (
m_workAsFourMuons) {(*m_msgStream) << MSG::INFO <<
" * FourMuonEvent::Init * working mode: 4 muons" <<
endmsg; }
78 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Init * Completed * " <<
endmsg;
99 if (pxMuonContainer !=
nullptr) {
101 <<
" track list has "<< pxMuonContainer->
size()
109 while ( xMuonItr != xMuonItrE ){
115 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Reco ** muon " << theCount <<
" is accepted " <<
endmsg;
123 if (!pxMuonContainer) {
124 std::cout <<
" * FourMuonEvent::Reco * Can't retrieve combined muon collection (container: " <<
m_container <<
") " << std::endl;
139 if (pxElecContainer !=
nullptr) {
140 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Reco * retrieving xAOD::ElectronContainer SUCCESS. "
142 <<
" size: " << pxElecContainer->
size()
144 if (pxElecContainer->
size() > 0 ){
151 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Reco * retrieving xAOD::ElectronContainer -- FAILED -- eventcount: " <<
m_eventCount <<
endmsg;
170 if (
m_doDebug) std::cout <<
" * FourMuonEvent::Reco * === Event " <<
m_eventCount <<
" is a GOOD 4-muon event " << std::endl;
173 if (
m_doDebug) std::cout <<
" * FourMuonEvent::Reco * === Event " <<
m_eventCount <<
" FAILS the 4-muon event selection " << std::endl;
179 if (
m_doDebug) std::cout <<
" * FourMuonEvent::Reco * applying 4 electron selection " << std::endl;
188 if (
m_doDebug) std::cout <<
" * FourMuonEvent::Reco * === Event " <<
m_eventCount <<
" is a GOOD 4-electrom event " << std::endl;
191 if (
m_doDebug) std::cout <<
" * FourMuonEvent::Reco * === Event " <<
m_eventCount <<
" FAILS the 4-electron event selection " << std::endl;
195 if (
m_doDebug) std::cout <<
" * FourMuonEvent::Reco * === Event " <<
m_eventCount <<
" is not a 4-electron event " << std::endl;
201 if (
m_doDebug) std::cout <<
" * FourMuonEvent::Reco * applying 4 lepton selection " << std::endl;
204 bool enoughleptons =
false;
209 if ( enoughleptons) {
211 std::cout <<
" * FourMuonEvent::Reco * Global statistics: Total number of accepted muons so far: " <<
m_acceptedMuonCount
227 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::Reco * === Event " <<
m_eventCount <<
" -- FAILS -- the 4-lepton event selection " <<
endmsg;
266 bool eventisgood =
true;
267 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent::EventSelection( type= " << eType <<
") ** started ** " << std::endl
278 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * Failing number of good muons == 4 :( "
287 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * Failing number of good electrons == 4 :( "
294 bool thisselection =
false;
298 if (!thisselection) {
299 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * Failing number of good muons >= 2 && electrons >= 2 :( "
304 eventisgood = thisselection;
310 unsigned int npassleadingpt = 0;
311 unsigned int npasssecondpt = 0;
317 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * using muon " <<
i <<
" with pt: " <<
m_pxIDTrack[
i]->
pt() << std::endl;}
321 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * #muons with pt > leading pt: " <<
m_LeadingMuonPtCut*
CLHEP::GeV <<
" = " << npassleadingpt << std::endl;}
322 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * #muons with pt > second pt: " <<
m_SecondMuonPtCut*
CLHEP::GeV <<
" = " << npasssecondpt << std::endl;}
323 if (npassleadingpt == 0) {
341 if (
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * using electron " <<
i <<
" with pt: " <<
m_pxELTrack[
i]->
pt() << std::endl;}
345 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * #elecs with pt > leading pt: " <<
m_LeadingMuonPtCut*
CLHEP::GeV <<
" = " << npassleadingpt << std::endl;}
346 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * #elecs with pt > second pt: " <<
m_SecondMuonPtCut*
CLHEP::GeV <<
" = " << npasssecondpt << std::endl;}
347 if (npassleadingpt == 0) {
362 std::cout <<
" * FourMuonEvent::EventSelection * Failing mass window low cut: reco m= " <<
m_fInvariantMass[eType] <<
" > " <<
m_MassWindowLow << std::endl;
368 std::cout <<
" * FourMuonEvent::EventSelection * Failing mass window high cut: reco m= " <<
m_fInvariantMass[eType] <<
" > " <<
m_MassWindowHigh << std::endl;
376 bool vertexstatus =
false;
378 std::vector <float> vtxListX;
379 std::vector <float> vtxListY;
380 std::vector <float> vtxListZ;
381 std::vector <int> vtxNpart;
383 int noVertexCountMuon = 0;
384 int noVertexCountElec = 0;
387 for (
unsigned int imu = 0; imu <
NUM_MUONS; imu++) {
425 for (
unsigned int iel = 0; iel <
NUM_MUONS; iel++) {
472 if (vtxListX.size()>0) vertexstatus =
true;
476 if (noVertexCountMuon > 0) vertexstatus =
false;
478 std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * vertices ID of the muons = " << std::endl
489 if (vtxListX.size()>0) vertexstatus =
true;
493 if (noVertexCountElec > 0) vertexstatus =
false;
495 std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * vertices ID of the electrons = " << std::endl
506 if (vtxListX.size()>0) vertexstatus =
true;
511 std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * vertices in event = " << vtxListX.size() << std::endl;
512 for (
size_t imu=0; imu <
NUM_MUONS; imu++) {
515 for (
size_t iel=0; iel <
NUM_MUONS; iel++) {
524 std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * Number of vertex found = " << vtxListX.size()
525 <<
" and mu without vertex: " << noVertexCountMuon
526 <<
" and elec without vertex: " << noVertexCountElec << std::endl;
527 for (
unsigned int ivtx=0; ivtx < vtxListX.size(); ivtx++) {
528 std::cout <<
" vtx[" << ivtx <<
"]= "
529 <<
"( " << vtxListX.at(ivtx)
530 <<
", " << vtxListY.at(ivtx)
531 <<
", " << vtxListZ.at(ivtx)
532 <<
") nparticles: " << vtxNpart.at(ivtx)
538 std::cout <<
" * FourMuonEvent::EventSelection * FAILED * number of vertex found = " << vtxListX.size()
539 <<
" mu without vertex: " << noVertexCountMuon
540 <<
" elec without vertex: " << noVertexCountElec << std::endl;
543 eventisgood = vertexstatus;
547 if(
m_doDebug){ std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") ** completed ** result= " << eventisgood << std::endl;}
556 if(
m_doDebug){ std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") ** started ** " << std::endl
565 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * Failing number of good muons and electrons == 4 :( "
572 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * Too many muons !! Failing number of good muons == 4 :( "
578 if(
m_doDebug) {std::cout <<
" * FourMuonEvent::EventSelection(" << eType <<
") * Too many electrons !! Failing number of good electrons == 4 :( "
584 double leadingMuonPt = -1., secondMuonPt=-1., thirdMuonPt=-1, fourthMuonPt=-1.;
628 double theLeadingPt = leadingMuonPt;
629 if (secondMuonPt > theLeadingPt) { theLeadingPt = secondMuonPt;}
630 if (thirdMuonPt > theLeadingPt) { theLeadingPt = thirdMuonPt;}
631 if (fourthMuonPt > theLeadingPt) { theLeadingPt = fourthMuonPt;}
633 double theTrailingPt = leadingMuonPt;
634 if (secondMuonPt < theTrailingPt && secondMuonPt > 0) { theTrailingPt = secondMuonPt;}
635 if (thirdMuonPt < theTrailingPt && thirdMuonPt > 0) { theTrailingPt = thirdMuonPt;}
636 if (fourthMuonPt < theTrailingPt && fourthMuonPt > 0) { theTrailingPt = fourthMuonPt;}
641 std::cout <<
" Pt of muons in this event 1: " << leadingMuonPt << std::endl
642 <<
" 2: " << secondMuonPt << std::endl
643 <<
" 3: " << thirdMuonPt << std::endl
644 <<
" 4: " << fourthMuonPt << std::endl
645 <<
" leading Pt: " << theLeadingPt << std::endl
646 <<
" trailing Pt: " << theTrailingPt << std::endl;
652 std::cout <<
" * FourMuonEvent::EventSelection * Failing pt cut * Reco Pt: leading " << theLeadingPt <<
" --> trailing " << theTrailingPt << std::endl;
657 std::cout <<
" * FourMuonEvent::EventSelection * Event passed the pt cuts: leading muon pt: " << leadingMuonPt << std::endl;
658 std::cout <<
" trailing muon pt: " << theTrailingPt << std::endl;
665 std::cout <<
" * FourMuonEvent::EventSelection * Failing mass window low cut: reco m= " <<
m_fInvariantMass[eType] <<
" > " <<
m_MassWindowLow << std::endl;
676 std::cout <<
" * FourMuonEvent::EventSelection * Event passed the mass window: " <<
m_fInvariantMass[eType] << std::endl;
865 std::cout <<
" * FourMuonEvent::EventSelection( type= " << eType <<
")* Good 4-muon set: pt range from " << leadingMuonPt/1000
866 <<
" to " << secondMuonPt/1000
867 <<
" GeV 4-muon invariant mass = " <<
m_fInvariantMass[eType] <<
" GeV " << std::endl;
868 std::cout <<
" * FourMuonEvent::EventSelection( type= " << eType <<
")* completed * " << std::endl;
924 bool thisdebug =
false;
928 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * RecordMuon * bad pxMuon --> EXIT "<< std::endl;}
943 if (
m_doDebug){ std::cout <<
" * FourMuonEvent * RecordMuon * bad pxMSmuon --> EXIT "<< std::endl;}
956 if (
m_doDebug){ std::cout <<
" * FourMuonEvent * RecordMuon * bad pxIDTrack for this muon--> EXIT "<< std::endl;}
977 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * ReconstructKinematics * -- start -- " << std::endl; }
978 bool kinematicscomputed =
false;
983 bool goodidentifiers =
true;
989 if (goodidentifiers) {
991 bool goodtracks =
true;
1009 kinematicscomputed =
true;
1012 std::cout <<
" * FourMuonEvent * ReconstructKinematics4Elec * -- Muon ID Tracks -- new -- " << std::endl
1025 if (!kinematicscomputed) {
1026 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * ReconstructKinematics * -- FAILED -- " << std::endl; }
1029 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * ReconstructKinematics * -- completed -- status: " << kinematicscomputed << std::endl; }
1030 return kinematicscomputed;
1036 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * ReconstructKinematics4Elec * -- start -- " << std::endl; }
1037 bool kinematicscomputed =
false;
1042 bool goodtracks =
true;
1047 if (
m_pxELTrack[0] ==
nullptr) goodtracks =
false;
1048 if (
m_pxELTrack[1] ==
nullptr) goodtracks =
false;
1049 if (
m_pxELTrack[2] ==
nullptr) goodtracks =
false;
1050 if (
m_pxELTrack[3] ==
nullptr) goodtracks =
false;
1059 kinematicscomputed =
true;
1062 std::cout <<
" * FourMuonEvent * ReconstructKinematics4Elec * -- Electron Tracks -- " << std::endl
1073 if (!kinematicscomputed) {
1074 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * ReconstructKinematics4Elec * -- FAILED -- " << std::endl; }
1077 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * ReconstructKinematics4Elec * -- completed -- status: " << kinematicscomputed << std::endl; }
1078 return kinematicscomputed;
1084 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent * ReconstructKinematicsNew * -- START -- " <<
endmsg;
1086 bool kinematicscomputed =
false;
1110 bool goodtracks =
true;
1111 if (
m_pxMUTrack[0] ==
nullptr) goodtracks =
false;
1112 if (
m_pxMUTrack[1] ==
nullptr) goodtracks =
false;
1113 if (
m_pxMUTrack[2] ==
nullptr) goodtracks =
false;
1114 if (
m_pxMUTrack[3] ==
nullptr) goodtracks =
false;
1122 kinematicscomputed =
true;
1127 std::cout <<
" * FourMuonEvent * ReconstructKinematicsNew * -- Muon ID Tracks -- new -- " << std::endl;
1132 if (kinematicscomputed) std::cout <<
" invariant mass (4mu) = " <<
m_fInvariantMass[
ID] << std::endl;
1135 double invmass_test = -1.;
1138 bool goodtracks =
true;
1139 if (
m_pxELTrack[0] ==
nullptr) goodtracks =
false;
1140 if (
m_pxELTrack[1] ==
nullptr) goodtracks =
false;
1141 if (
m_pxELTrack[2] ==
nullptr) goodtracks =
false;
1142 if (
m_pxELTrack[3] ==
nullptr) goodtracks =
false;
1144 if (goodtracks && !kinematicscomputed) {
1152 kinematicscomputed =
true;
1156 std::cout <<
" * FourMuonEvent * ReconstructKinematicsNew * -- Electron Tracks -- " << std::endl;
1161 std::cout <<
" invariant mass (4e) = " << invmass_test << std::endl;
1166 bool goodtracks =
true;
1167 if (
m_pxMUTrack[0] ==
nullptr) goodtracks =
false;
1168 if (
m_pxMUTrack[2] ==
nullptr) goodtracks =
false;
1169 if (
m_pxELTrack[0] ==
nullptr) goodtracks =
false;
1170 if (
m_pxELTrack[2] ==
nullptr) goodtracks =
false;
1172 if (goodtracks && !kinematicscomputed) {
1180 kinematicscomputed =
true;
1184 std::cout <<
" * FourMuonEvent * ReconstructKinematicsNew * -- Muon and Electron Tracks -- " << std::endl;
1194 std::cout <<
" invariant mass test = " << invmass_test << std::endl;
1197 if (!kinematicscomputed) {
1198 if(
m_doDebug){ std::cout <<
" * FourMuonEvent * ReconstructKinematicsNew * -- FAILED -- " << std::endl; }
1201 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonEvent * ReconstructKinematicsNew * -- COMPLETED -- status " << kinematicscomputed <<
endmsg;
1202 return kinematicscomputed;
1290 if ( pxTrackContainer )
1294 while ( xTrkItr != xTrkItrE )
1297 if ( !pxTrack )
continue;
1299 if(!pxTrkTrack)
continue;
1301 if ( !pxPerigee )
continue;
1302 const float fTrkPhi = pxPerigee->parameters()[
Trk::phi];
1303 const float fTrkEta = pxPerigee->eta();
1307 float fDR = sqrt( fDPhi*fDPhi + fDEta*fDEta );
1334 std::cout <<
" * FourMuonEvent * SetLeadingMuonPtCut * new Pt cuts: " <<
m_LeadingMuonPtCut <<
" & "
1356 std::cout <<
" * FourMuonEvent * SetSecondMuonPtCut * new Pt cuts: " <<
m_LeadingMuonPtCut
1368 bool thisdebug =
false;
1370 if (
m_doDebug || thisdebug) {std::cout <<
" * FourMuonEvent::OrderMuonList * -- start -- " << std::endl
1375 int muMinus1Id = -9;
1376 int muMinus2Id = -9;
1377 double muPlus1Pt = 0.;
1378 double muPlus2Pt = 0.;
1379 double muMinus1Pt = 0.;
1380 double muMinus2Pt = 0.;
1388 if (nMuonsAtEntry >= 2) {
1389 for (
int imuon=0; imuon < (
int) nMuonsAtEntry; imuon++) {
1390 if(
m_doDebug &&
false ){ std::cout <<
" * FourMuonEvent::OrderMuonList * testing imuon= " << imuon
1401 muPlus2Pt = muPlus1Pt;
1402 muPlus2Id = muPlus1Id;
1418 muMinus2Pt = muMinus1Pt;
1419 muMinus2Id = muMinus1Id;
1433 if (nMuonsAtEntry >= 2 && (muposcount == 0 || munegcount == 0)) {
1434 if (
m_doDebug) std::cout <<
" -- FourMuonEvent::OrderMuonList -- No opposite charge muons in the " << nMuonsAtEntry <<
" input muons"
1435 <<
" #mu+ " << muposcount
1436 <<
" #mu- " << munegcount
1437 <<
" --> DISCARD ALL MUONS -- " << std::endl;
1455 std::cout <<
" * FourMuonEvent::OrderMuonList * taking " <<
m_numberOfFullPassMuons <<
" muons from the input list of " << nMuonsAtEntry <<
" muons: " << std::endl;
1456 if (muMinus1Id >= 0) std::cout <<
" leading mu-: " << muMinus1Id <<
" Pt = " << muMinus1Pt << std::endl;
1457 if (muMinus2Id >= 0) std::cout <<
" second mu-: " << muMinus2Id <<
" Pt = " << muMinus2Pt << std::endl;
1458 if (muPlus1Id >= 0) std::cout <<
" leading mu+: " << muPlus1Id <<
" Pt = " << muPlus1Pt << std::endl;
1459 if (muPlus2Id >= 0) std::cout <<
" second mu+: " << muPlus2Id <<
" Pt = " << muPlus2Pt << std::endl;
1462 if (
m_doDebug) std::cout <<
" * FourMuonEvent::OrderMuonList * This event has less than 2 muons :(" << std::endl;
1465 if (
m_doDebug || thisdebug) std::cout <<
" * FourMuonEvent::OrderMuonList * completed * m_numberOfFullPassMuons= " <<
m_numberOfFullPassMuons << std::endl;
1472 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonsEvents::CheckMuonVertices * -- START --" <<
endmsg;
1474 bool goodvertices =
false;
1475 goodvertices =
true;
1476 int nverticesfound = 0;
1480 for (
unsigned int imuneg = 0; imuneg <= 1; imuneg++) {
1516 if (nverticesfound >= 1) goodvertices =
true;
1518 if (nverticesfound == 0)
if (
m_doDebug) std::cout <<
" -- FourMuonEvent::CheckMuonVertices -- WARNING -- MUONS DO NOT COME FROM SAME VERTEX " << std::endl;
1520 (*m_msgStream) <<
MSG::DEBUG <<
" * FourMuonsEvents::CheckMuonVertices * -- COMPLETED -- status: " << goodvertices <<
endmsg;
1521 return goodvertices;