Skip to content

Commit

Permalink
improve GeneralPurposeVertexAnalyzer
Browse files Browse the repository at this point in the history
  • Loading branch information
mmusich committed Jan 8, 2025
1 parent b5cf409 commit a3a7371
Showing 1 changed file with 109 additions and 4 deletions.
113 changes: 109 additions & 4 deletions Alignment/OfflineValidation/plugins/GeneralPurposeVertexAnalyzer.cc
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,50 @@

using reco::TrackCollection;

namespace gpVertexAnalyzer {
void setBinLog(TAxis *axis) {
int bins = axis->GetNbins();
float from = axis->GetXmin();
float to = axis->GetXmax();
float width = (to - from) / bins;
std::vector<float> new_bins(bins + 1, 0);
for (int i = 0; i <= bins; i++) {
new_bins[i] = TMath::Power(10, from + i * width);
}
axis->Set(bins, new_bins.data());
}

void setBinLogX(TH1 *h) {
TAxis *axis = h->GetXaxis();
setBinLog(axis);
}

void setBinLogY(TH1 *h) {
TAxis *axis = h->GetYaxis();
setBinLog(axis);
}

template <typename... Args>
TProfile *makeProfileIfLog(const edm::Service<TFileService> &fs, bool logx, bool logy, Args &&...args) {
auto prof = fs->make<TProfile>(std::forward<Args>(args)...);
if (logx)
setBinLogX(prof);
if (logy)
setBinLogY(prof);
return prof;
}

template <typename... Args>
TH1D *makeTH1IfLog(const edm::Service<TFileService> &fs, bool logx, bool logy, Args &&...args) {
auto h1 = fs->make<TH1D>(std::forward<Args>(args)...);
if (logx)
setBinLogX(h1);
if (logy)
setBinLogY(h1);
return h1;
}
} // namespace gpVertexAnalyzer

class GeneralPurposeVertexAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
public:
explicit GeneralPurposeVertexAnalyzer(const edm::ParameterSet &);
Expand All @@ -73,9 +117,9 @@ class GeneralPurposeVertexAnalyzer : public edm::one::EDAnalyzer<edm::one::Share
struct IPMonitoring {
std::string varname_;
float pTcut_;
TH1D *IP_, *IPErr_;
TProfile *IPVsPhi_, *IPVsEta_;
TProfile *IPErrVsPhi_, *IPErrVsEta_;
TH1D *IP_, *IPErr_, *IPPull_;
TProfile *IPVsPhi_, *IPVsEta_, *IPVsPt_;
TProfile *IPErrVsPhi_, *IPErrVsEta_, *IPErrVsPt_;
TProfile2D *IPVsEtaVsPhi_, *IPErrVsEtaVsPhi_;

void bookIPMonitor(const edm::ParameterSet &, const edm::Service<TFileService> fs);
Expand Down Expand Up @@ -128,7 +172,7 @@ class GeneralPurposeVertexAnalyzer : public edm::one::EDAnalyzer<edm::one::Share
TH1D *type[2];
TH1D *bsX, *bsY, *bsZ, *bsSigmaZ, *bsDxdz, *bsDydz, *bsBeamWidthX, *bsBeamWidthY, *bsType;

TH1D *sumpt, *ntracks, *weight, *chi2ndf, *chi2prob;
TH1D *trackpt, *sumpt, *ntracks, *weight, *chi2ndf, *chi2prob;
TH1D *dxy2;
TH1D *phi_pt1, *eta_pt1;
TH1D *phi_pt10, *eta_pt10;
Expand Down Expand Up @@ -157,6 +201,10 @@ void GeneralPurposeVertexAnalyzer::IPMonitoring::bookIPMonitor(const edm::Parame
double EtaMin = config.getParameter<double>("EtaMin");
double EtaMax = config.getParameter<double>("EtaMax");

int PtBin = config.getParameter<int>("PtBin");
double PtMin = config.getParameter<double>("PtMin") * pTcut_;
double PtMax = config.getParameter<double>("PtMax") * pTcut_;

IP_ = fs->make<TH1D>(fmt::format("d{}_pt{}", varname_, pTcut_).c_str(),
fmt::format("PV tracks (p_{{T}} > {} GeV) d_{{{}}} (#mum)", pTcut_, varname_).c_str(),
VarBin,
Expand All @@ -169,6 +217,13 @@ void GeneralPurposeVertexAnalyzer::IPMonitoring::bookIPMonitor(const edm::Parame
0.,
(varname_.find("xy") != std::string::npos) ? 2000. : 10000.);

IPPull_ = fs->make<TH1D>(
fmt::format("d{}Pull_pt{}", varname_, pTcut_).c_str(),
fmt::format("PV tracks (p_{{T}} > {} GeV) d_{{{}}}/#sigma_{{d_{{{}}}}}", pTcut_, varname_, varname_).c_str(),
100,
-5.,
5.);

IPVsPhi_ =
fs->make<TProfile>(fmt::format("d{}VsPhi_pt{}", varname_, pTcut_).c_str(),
fmt::format("PV tracks (p_{{T}} > {}) d_{{{}}} VS track #phi", pTcut_, varname_).c_str(),
Expand All @@ -193,6 +248,21 @@ void GeneralPurposeVertexAnalyzer::IPMonitoring::bookIPMonitor(const edm::Parame
IPVsEta_->SetXTitle("PV track (p_{T} > 1 GeV) #eta");
IPVsEta_->SetYTitle(fmt::format("PV tracks (p_{{T}} > {} GeV) d_{{{}}} (#mum)", pTcut_, varname_).c_str());

IPVsPt_ = gpVertexAnalyzer::makeProfileIfLog(
fs,
true, /* x-axis */
false, /* y-axis */
fmt::format("d{}VsPt_pt{}", varname_, pTcut_).c_str(),
fmt::format("PV tracks (p_{{T}} > {}) d_{{{}}} VS track p_{{T}}", pTcut_, varname_).c_str(),
PtBin,
log10(PtMin),
log10(PtMax),
VarMin,
VarMax,
"");
IPVsPt_->SetXTitle("PV track (p_{T} > 1 GeV) p_{T} [GeV]");
IPVsPt_->SetYTitle(fmt::format("PV tracks (p_{{T}} > {} GeV) d_{{{}}} (#mum)", pTcut_, varname_).c_str());

IPErrVsPhi_ =
fs->make<TProfile>(fmt::format("d{}ErrVsPhi_pt{}", varname_, pTcut_).c_str(),
fmt::format("PV tracks (p_{{T}} > {}) d_{{{}}} error VS track #phi", pTcut_, varname_).c_str(),
Expand All @@ -217,6 +287,21 @@ void GeneralPurposeVertexAnalyzer::IPMonitoring::bookIPMonitor(const edm::Parame
IPErrVsEta_->SetXTitle("PV track (p_{T} > 1 GeV) #eta");
IPErrVsEta_->SetYTitle(fmt::format("PV tracks (p_{{T}} > {} GeV) d_{{{}}} error (#mum)", pTcut_, varname_).c_str());

IPErrVsPt_ = gpVertexAnalyzer::makeProfileIfLog(
fs,
true, /* x-axis */
false, /* y-axis */
fmt::format("d{}ErrVsPt_pt{}", varname_, pTcut_).c_str(),
fmt::format("PV tracks (p_{{T}} > {}) d_{{{}}} error VS track p_{{T}}", pTcut_, varname_).c_str(),
PtBin,
log10(PtMin),
log10(PtMax),
VarMin,
VarMax,
"");
IPErrVsPt_->SetXTitle("PV track (p_{T} > 1 GeV) p_{T} [GeV]");
IPErrVsPt_->SetYTitle(fmt::format("PV tracks (p_{{T}} > {} GeV) d_{{{}}} error (#mum)", pTcut_, varname_).c_str());

IPVsEtaVsPhi_ = fs->make<TProfile2D>(
fmt::format("d{}VsEtaVsPhi_pt{}", varname_, pTcut_).c_str(),
fmt::format("PV tracks (p_{{T}} > {}) d_{{{}}} VS track #eta VS track #phi", pTcut_, varname_).c_str(),
Expand Down Expand Up @@ -294,6 +379,7 @@ GeneralPurposeVertexAnalyzer::GeneralPurposeVertexAnalyzer(const edm::ParameterS
bsBeamWidthX(nullptr),
bsBeamWidthY(nullptr),
bsType(nullptr),
trackpt(nullptr),
sumpt(nullptr),
ntracks(nullptr),
weight(nullptr),
Expand Down Expand Up @@ -411,6 +497,8 @@ void GeneralPurposeVertexAnalyzer::pvTracksPlots(const reco::Vertex &v) {
}

const float pt = t->pt();
trackpt->Fill(pt);

if (pt < 1.f) {
continue;
}
Expand All @@ -435,22 +523,28 @@ void GeneralPurposeVertexAnalyzer::pvTracksPlots(const reco::Vertex &v) {
phi_pt1->Fill(phi);
eta_pt1->Fill(eta);

// dxy pT>1

dxy_pt1.IP_->Fill(Dxy);
dxy_pt1.IPVsPhi_->Fill(phi, Dxy);
dxy_pt1.IPVsEta_->Fill(eta, Dxy);
dxy_pt1.IPVsEtaVsPhi_->Fill(eta, phi, Dxy);

dxy_pt1.IPErr_->Fill(DxyErr);
dxy_pt1.IPPull_->Fill(Dxy / DxyErr);
dxy_pt1.IPErrVsPhi_->Fill(phi, DxyErr);
dxy_pt1.IPErrVsEta_->Fill(eta, DxyErr);
dxy_pt1.IPErrVsEtaVsPhi_->Fill(eta, phi, DxyErr);

// dz pT>1

dz_pt1.IP_->Fill(Dz);
dz_pt1.IPVsPhi_->Fill(phi, Dz);
dz_pt1.IPVsEta_->Fill(eta, Dz);
dz_pt1.IPVsEtaVsPhi_->Fill(eta, phi, Dz);

dz_pt1.IPErr_->Fill(DzErr);
dz_pt1.IPPull_->Fill(Dz / DzErr);
dz_pt1.IPErrVsPhi_->Fill(phi, DzErr);
dz_pt1.IPErrVsEta_->Fill(eta, DzErr);
dz_pt1.IPErrVsEtaVsPhi_->Fill(eta, phi, DzErr);
Expand All @@ -459,22 +553,28 @@ void GeneralPurposeVertexAnalyzer::pvTracksPlots(const reco::Vertex &v) {
phi_pt10->Fill(phi);
eta_pt10->Fill(eta);

// dxy pT>10

dxy_pt10.IP_->Fill(Dxy);
dxy_pt10.IPVsPhi_->Fill(phi, Dxy);
dxy_pt10.IPVsEta_->Fill(eta, Dxy);
dxy_pt10.IPVsEtaVsPhi_->Fill(eta, phi, Dxy);

dxy_pt10.IPErr_->Fill(DxyErr);
dxy_pt10.IPPull_->Fill(Dxy / DxyErr);
dxy_pt10.IPErrVsPhi_->Fill(phi, DxyErr);
dxy_pt10.IPErrVsEta_->Fill(eta, DxyErr);
dxy_pt10.IPErrVsEtaVsPhi_->Fill(eta, phi, DxyErr);

// dz pT>10

dz_pt10.IP_->Fill(Dz);
dz_pt10.IPVsPhi_->Fill(phi, Dz);
dz_pt10.IPVsEta_->Fill(eta, Dz);
dz_pt10.IPVsEtaVsPhi_->Fill(eta, phi, Dz);

dz_pt10.IPErr_->Fill(DzErr);
dz_pt10.IPPull_->Fill(Dz / DzErr);
dz_pt10.IPErrVsPhi_->Fill(phi, DzErr);
dz_pt10.IPErrVsEta_->Fill(eta, DzErr);
dz_pt10.IPErrVsEtaVsPhi_->Fill(eta, phi, DzErr);
Expand Down Expand Up @@ -649,6 +749,8 @@ void GeneralPurposeVertexAnalyzer::beginJob() {

dxy2 = book<TH1D>("dxyzoom", fmt::sprintf("%s d_{xy} (#mum)", s_1).c_str(), dxyBin_, dxyMin_ / 5., dxyMax_ / 5.);

trackpt = gpVertexAnalyzer::makeTH1IfLog(
fs_, true, false, "pt_track", "PV tracks p_{T};PV tracks p_{T} [GeV];#tracks", 49, log10(1.), log10(50.));
phi_pt1 =
book<TH1D>("phi_pt1", fmt::sprintf("%s #phi; PV tracks #phi;#tracks", s_1).c_str(), phiBin_, phiMin_, phiMax_);
eta_pt1 =
Expand Down Expand Up @@ -701,6 +803,9 @@ void GeneralPurposeVertexAnalyzer::fillDescriptions(edm::ConfigurationDescriptio
desc.add<int>("EtaBin2D", 8);
desc.add<double>("EtaMin", -2.7);
desc.add<double>("EtaMax", 2.7);
desc.add<int>("PtBin", 49);
desc.add<double>("PtMin", 1.);
desc.add<double>("PtMax", 50.);
descriptions.addWithDefaultLabel(desc);
}

Expand Down

0 comments on commit a3a7371

Please sign in to comment.