Skip to content

Commit

Permalink
Merge pull request #42973 from mmusich/addUnbalanceForZPointingAngle_…
Browse files Browse the repository at this point in the history
…13_2_X

[TkAl] miscellaneous updates related to analysis of μμ + vertex events (13.2.X)
  • Loading branch information
cmsbuild authored Oct 16, 2023
2 parents a8046d7 + 6d05bd7 commit bbbd24b
Show file tree
Hide file tree
Showing 6 changed files with 325 additions and 15 deletions.
91 changes: 82 additions & 9 deletions Alignment/OfflineValidation/plugins/DiMuonVertexValidation.cc
Original file line number Diff line number Diff line change
Expand Up @@ -70,12 +70,14 @@ class DiMuonVertexValidation : public edm::one::EDAnalyzer<edm::one::SharedResou
private:
void beginJob() override;
void analyze(const edm::Event&, const edm::EventSetup&) override;
const reco::Vertex* findClosestVertex(const TransientVertex aTransVtx, const reco::VertexCollection* vertices) const;
void endJob() override;

// ----------member data ---------------------------
DiLeptonHelp::Counts myCounts;

bool useReco_;
const bool useReco_;
const bool useClosestVertex_;
std::vector<double> pTthresholds_;
float maxSVdist_;

Expand Down Expand Up @@ -146,6 +148,7 @@ static constexpr float mumass2 = 0.105658367 * 0.105658367; //mu mass squared (
//
DiMuonVertexValidation::DiMuonVertexValidation(const edm::ParameterSet& iConfig)
: useReco_(iConfig.getParameter<bool>("useReco")),
useClosestVertex_(iConfig.getParameter<bool>("useClosestVertex")),
pTthresholds_(iConfig.getParameter<std::vector<double>>("pTThresholds")),
maxSVdist_(iConfig.getParameter<double>("maxSVdist")),
CosPhiConfiguration_(iConfig.getParameter<edm::ParameterSet>("CosPhiConfig")),
Expand Down Expand Up @@ -269,6 +272,14 @@ void DiMuonVertexValidation::analyze(const edm::Event& iEvent, const edm::EventS
}
}

#ifdef EDM_ML_DEBUG
for (const auto& track : myTracks) {
edm::LogVerbatim("DiMuonVertexValidation") << __PRETTY_FUNCTION__ << " pT: " << track->pt() << " GeV"
<< " , pT error: " << track->ptError() << " GeV"
<< " , eta: " << track->eta() << " , phi: " << track->phi() << std::endl;
}
#endif

LogDebug("DiMuonVertexValidation") << "selected tracks: " << myTracks.size() << std::endl;

const TransientTrackBuilder* theB = &iSetup.getData(ttbESToken_);
Expand All @@ -288,6 +299,17 @@ void DiMuonVertexValidation::analyze(const edm::Event& iEvent, const edm::EventS
TLorentzVector p4_tplus(tplus->px(), tplus->py(), tplus->pz(), sqrt((tplus->p() * tplus->p()) + mumass2));
TLorentzVector p4_tminus(tminus->px(), tminus->py(), tminus->pz(), sqrt((tminus->p() * tminus->p()) + mumass2));

#ifdef EDM_ML_DEBUG
// Define a lambda function to convert TLorentzVector to a string
auto tLorentzVectorToString = [](const TLorentzVector& vector) {
return std::to_string(vector.Px()) + " " + std::to_string(vector.Py()) + " " + std::to_string(vector.Pz()) + " " +
std::to_string(vector.E());
};

edm::LogVerbatim("DiMuonVertexValidation") << "mu+" << tLorentzVectorToString(p4_tplus) << std::endl;
edm::LogVerbatim("DiMuonVertexValidation") << "mu-" << tLorentzVectorToString(p4_tminus) << std::endl;
#endif

const auto& Zp4 = p4_tplus + p4_tminus;
float track_invMass = Zp4.M();
hTrackInvMass_->Fill(track_invMass);
Expand Down Expand Up @@ -323,25 +345,41 @@ void DiMuonVertexValidation::analyze(const edm::Event& iEvent, const edm::EventS
// fill the VtxProb plots
VtxProbPlots.fillPlots(SVProb, tktk_p4);

math::XYZPoint MainVertex(0, 0, 0);
const reco::Vertex* theClosestVertex = nullptr;
// get collection of reconstructed vertices from event
edm::Handle<reco::VertexCollection> vertexHandle = iEvent.getHandle(vertexToken_);

math::XYZPoint MainVertex(0, 0, 0);
reco::Vertex TheMainVertex = vertexHandle.product()->front();

if (vertexHandle.isValid()) {
const reco::VertexCollection* vertices = vertexHandle.product();
if ((*vertices).at(0).isValid()) {
auto theMainVtx = (*vertices).at(0);
MainVertex.SetXYZ(theMainVtx.position().x(), theMainVtx.position().y(), theMainVtx.position().z());
}
theClosestVertex = this->findClosestVertex(aTransientVertex, vertices);
} else {
edm::LogWarning("DiMuonVertexMonitor") << "invalid vertex collection encountered Skipping event!";
return;
}

reco::Vertex TheMainVertex;
if (!useClosestVertex_ || theClosestVertex == nullptr) {
// if the closest vertex is not available, or explicitly not chosen
TheMainVertex = vertexHandle.product()->front();
} else {
TheMainVertex = *theClosestVertex;
}

MainVertex.SetXYZ(TheMainVertex.position().x(), TheMainVertex.position().y(), TheMainVertex.position().z());
const math::XYZPoint myVertex(
aTransientVertex.position().x(), aTransientVertex.position().y(), aTransientVertex.position().z());
const math::XYZPoint deltaVtx(
MainVertex.x() - myVertex.x(), MainVertex.y() - myVertex.y(), MainVertex.z() - myVertex.z());

#ifdef EDM_ML_DEBUG
edm::LogVerbatim("DiMuonVertexValidation")
<< "mm vertex position:" << aTransientVertex.position().x() << "," << aTransientVertex.position().y() << ","
<< aTransientVertex.position().z();

edm::LogVerbatim("DiMuonVertexValidation") << "main vertex position:" << TheMainVertex.position().x() << ","
<< TheMainVertex.position().y() << "," << TheMainVertex.position().z();
#endif

if (TheMainVertex.isValid()) {
// Z Vertex distance in the xy plane

Expand Down Expand Up @@ -391,6 +429,11 @@ void DiMuonVertexValidation::analyze(const edm::Event& iEvent, const edm::EventS
hCosPhi_->Fill(cosphi);
hCosPhi3D_->Fill(cosphi3D);

#ifdef EDM_ML_DEBUG
edm::LogVerbatim("DiMuonVertexValidation")
<< "distance " << distance * cmToum << " cosphi3D:" << cosphi3D << std::endl;
#endif

// fill the cosphi plots
CosPhiPlots.fillPlots(cosphi, tktk_p4);

Expand Down Expand Up @@ -489,6 +532,35 @@ void DiMuonVertexValidation::endJob() {
}
}

// compute the closest vertex to di-lepton ------------------------------------
const reco::Vertex* DiMuonVertexValidation::findClosestVertex(const TransientVertex aTransVtx,
const reco::VertexCollection* vertices) const {
reco::Vertex* defaultVtx = nullptr;

if (!aTransVtx.isValid())
return defaultVtx;

// find the closest vertex to the secondary vertex in 3D
VertexDistance3D vertTool3D;
float minD = 9999.;
int closestVtxIndex = 0;
int counter = 0;
for (const auto& vtx : *vertices) {
double dist3D = vertTool3D.distance(aTransVtx, vtx).value();
if (dist3D < minD) {
minD = dist3D;
closestVtxIndex = counter;
}
counter++;
}

if ((*vertices).at(closestVtxIndex).isValid()) {
return &(vertices->at(closestVtxIndex));
} else {
return defaultVtx;
}
}

// ------------ method fills 'descriptions' with the allowed parameters for the module ------------
void DiMuonVertexValidation::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
edm::ParameterSetDescription desc;
Expand All @@ -503,6 +575,7 @@ void DiMuonVertexValidation::fillDescriptions(edm::ConfigurationDescriptions& de
desc.add<edm::InputTag>("tracks", edm::InputTag("generalTracks"));
desc.add<edm::InputTag>("vertices", edm::InputTag("offlinePrimaryVertices"));
desc.add<std::vector<double>>("pTThresholds", {30., 10.});
desc.add<bool>("useClosestVertex", true);
desc.add<double>("maxSVdist", 50.);

{
Expand Down
11 changes: 11 additions & 0 deletions Alignment/OfflineValidation/test/DiMuonVertexValidation_cfg.py
Original file line number Diff line number Diff line change
Expand Up @@ -96,6 +96,17 @@ def best_match(rcd):
elif options.era=='2018':
print("===> running era 2018")
process = cms.Process('Analysis',eras.Run2_2018)
elif options.era=='2022':
print("===> running era 2022")
process = cms.Process('Analysis',eras.Run3)
elif options.era=='2023':
print("===> running era 2023")
process = cms.Process('Analysis',eras.Run3_2023)

###################################################################
# Set the process to run multi-threaded
###################################################################
process.options.numberOfThreads = 8

# import of standard configurations
process.load('Configuration.StandardSequences.Services_cff')
Expand Down
23 changes: 23 additions & 0 deletions DQMOffline/Alignment/interface/DiMuonVertexMonitor.h
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
#include <string>

// user includes
#include "DQMOffline/Alignment/interface/DiLeptonPlotHelpers.h"
#include "DQMServices/Core/interface/DQMEDAnalyzer.h"
#include "DQMServices/Core/interface/DQMStore.h"
#include "FWCore/Framework/interface/Event.h"
Expand Down Expand Up @@ -76,9 +77,31 @@ class DiMuonVertexMonitor : public DQMEDAnalyzer {
MonitorElement *hCosPhi3D_;
MonitorElement *hCosPhiInv_;
MonitorElement *hCosPhiInv3D_;
MonitorElement *hCosPhiUnbalance_;
MonitorElement *hCosPhi3DUnbalance_;
MonitorElement *hInvMass_;
MonitorElement *hCutFlow_;

// 2D histograms of pointing angle vs variable
edm::ParameterSet CosPhi3DConfiguration_;
DiLepPlotHelp::PlotsVsKinematics CosPhi3DPlots_ = DiLepPlotHelp::PlotsVsKinematics(DiLepPlotHelp::MM);

// 2D histograms of 3D PV-SV distance vs variable
edm::ParameterSet SVDistConfiguration_;
DiLepPlotHelp::PlotsVsKinematics SVDistPlots_ = DiLepPlotHelp::PlotsVsKinematics(DiLepPlotHelp::MM);

// 2D histograms of 3D PV-SV distance significance vs variable
edm::ParameterSet SVDistSigConfiguration_;
DiLepPlotHelp::PlotsVsKinematics SVDistSigPlots_ = DiLepPlotHelp::PlotsVsKinematics(DiLepPlotHelp::MM);

// 2D histograms of PV-SV transverse distance vs variable
edm::ParameterSet SVDist3DConfiguration_;
DiLepPlotHelp::PlotsVsKinematics SVDist3DPlots_ = DiLepPlotHelp::PlotsVsKinematics(DiLepPlotHelp::MM);

// 2D histograms of PV-SV transverse distance significance vs variable
edm::ParameterSet SVDist3DSigConfiguration_;
DiLepPlotHelp::PlotsVsKinematics SVDist3DSigPlots_ = DiLepPlotHelp::PlotsVsKinematics(DiLepPlotHelp::MM);

// impact parameters information
MonitorElement *hdxy_;
MonitorElement *hdz_;
Expand Down
14 changes: 12 additions & 2 deletions DQMOffline/Alignment/python/ALCARECOTkAlDQM_cff.py
Original file line number Diff line number Diff line change
Expand Up @@ -164,7 +164,12 @@
decayMotherName = "J/#psi",
vertices = 'offlinePrimaryVertices',
FolderName = "AlCaReco/"+__selectionName,
maxSVdist = 50
maxSVdist = 50,
CosPhi3DConfig = dict(maxDeltaEta = 1.3),
SVDistConfig = dict(maxDeltaEta = 1.3),
SVDistSigConfig = dict(maxDeltaEta = 1.3),
SVDist3DConfig = dict(maxDeltaEta = 1.3),
SVDist3DSigConfig = dict(maxDeltaEta = 1.3)
)

ALCARECOTkAlJpsiMassBiasDQM = DQMOffline.Alignment.DiMuonMassBiasMonitor_cfi.DiMuonMassBiasMonitor.clone(
Expand Down Expand Up @@ -234,7 +239,12 @@
decayMotherName = "#Upsilon",
vertices = 'offlinePrimaryVertices',
FolderName = "AlCaReco/"+__selectionName,
maxSVdist = 50
maxSVdist = 50,
CosPhi3DConfig = dict(maxDeltaEta = 1.6),
SVDistConfig = dict(maxDeltaEta = 1.6),
SVDistSigConfig = dict(maxDeltaEta = 1.6),
SVDist3DConfig = dict(maxDeltaEta = 1.6),
SVDist3DSigConfig = dict(maxDeltaEta = 1.6)
)

ALCARECOTkAlUpsilonMassBiasDQM = DQMOffline.Alignment.DiMuonMassBiasMonitor_cfi.DiMuonMassBiasMonitor.clone(
Expand Down
52 changes: 51 additions & 1 deletion DQMOffline/Alignment/python/DiMuonVertexMonitor_cfi.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,4 +6,54 @@
decayMotherName = cms.string('Z'),
vertices = cms.InputTag('offlinePrimaryVertices'),
FolderName = cms.string('DiMuonVertexMonitor'),
maxSVdist = cms.double(50))
maxSVdist = cms.double(50),
CosPhi3DConfig = cms.PSet(
name = cms.string('CosPhi3D'),
title = cms.string('cos(#phi_{3D})'),
yUnits = cms.string(''),
NxBins = cms.int32(24),
NyBins = cms.int32(50),
ymin = cms.double(-1),
ymax = cms.double(1),
maxDeltaEta = cms.double(3.7)
),
SVDistConfig = cms.PSet(
name = cms.string('SVDist'),
title = cms.string('PV-SV distance'),
yUnits = cms.string('[#mum]'),
NxBins = cms.int32(24),
NyBins = cms.int32(100),
ymin = cms.double(0),
ymax = cms.double(300),
maxDeltaEta = cms.double(3.7)
),
SVDistSigConfig = cms.PSet(
name = cms.string('SVDistSig'),
title = cms.string('PV-SV distance significance'),
yUnits = cms.string('[#mum]'),
NxBins = cms.int32(24),
NyBins = cms.int32(100),
ymin = cms.double(0),
ymax = cms.double(5),
maxDeltaEta = cms.double(3.7)
),
SVDist3DConfig = cms.PSet(
name = cms.string('SVDist3D'),
title = cms.string('PV-SV 3D distance'),
yUnits = cms.string('[#mum]'),
NxBins = cms.int32(24),
NyBins = cms.int32(100),
ymin = cms.double(0),
ymax = cms.double(300),
maxDeltaEta = cms.double(3.7)
),
SVDist3DSigConfig = cms.PSet(
name = cms.string('SVDist3DSig'),
title = cms.string('PV-SV 3D distance significance'),
yUnits = cms.string('[#mum]'),
NxBins = cms.int32(24),
NyBins = cms.int32(100),
ymin = cms.double(0),
ymax = cms.double(5),
maxDeltaEta = cms.double(3.7)
))
Loading

0 comments on commit bbbd24b

Please sign in to comment.