diff --git a/Validation/Configuration/python/mtdSimValid_cff.py b/Validation/Configuration/python/mtdSimValid_cff.py index 0279b8c25a91a..93b57e08af013 100644 --- a/Validation/Configuration/python/mtdSimValid_cff.py +++ b/Validation/Configuration/python/mtdSimValid_cff.py @@ -7,7 +7,9 @@ from Validation.MtdValidation.etlSimHits_cfi import etlSimHits from Validation.MtdValidation.etlDigiHits_cfi import etlDigiHits from Validation.MtdValidation.etlRecHits_cfi import etlRecHits +from Validation.MtdValidation.btlReco_cfi import btlReco +from Validation.MtdValidation.etlReco_cfi import etlReco mtdSimValid = cms.Sequence(btlSimHits + etlSimHits ) mtdDigiValid = cms.Sequence(btlDigiHits + etlDigiHits) -mtdRecoValid = cms.Sequence(btlRecHits + etlRecHits ) +mtdRecoValid = cms.Sequence(btlRecHits + etlRecHits + btlReco + etlReco) diff --git a/Validation/MtdValidation/plugins/BtlRecoHarvester.cc b/Validation/MtdValidation/plugins/BtlRecoHarvester.cc new file mode 100644 index 0000000000000..39b6985b865b0 --- /dev/null +++ b/Validation/MtdValidation/plugins/BtlRecoHarvester.cc @@ -0,0 +1,123 @@ +#include + +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ServiceRegistry/interface/Service.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" + +#include "DQMServices/Core/interface/DQMEDHarvester.h" +#include "DQMServices/Core/interface/DQMStore.h" + +#include "DataFormats/ForwardDetId/interface/BTLDetId.h" + +class BtlRecoHarvester : public DQMEDHarvester { +public: + explicit BtlRecoHarvester(const edm::ParameterSet& iConfig); + ~BtlRecoHarvester() override; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override; + +private: + const std::string folder_; + + // --- Histograms + MonitorElement* meBtlEtaEff_; + MonitorElement* meBtlPhiEff_; + MonitorElement* meBtlPtEff_; +}; + +// ------------ constructor and destructor -------------- +BtlRecoHarvester::BtlRecoHarvester(const edm::ParameterSet& iConfig) + : folder_(iConfig.getParameter("folder")) {} + +BtlRecoHarvester::~BtlRecoHarvester() {} + +// ------------ endjob tasks ---------------------------- +void BtlRecoHarvester::dqmEndJob(DQMStore::IBooker& ibook, DQMStore::IGetter& igetter) { + // --- Get the monitoring histograms + MonitorElement* meTrackEffEtaTot = igetter.get(folder_ + "TrackEffEtaTot"); + MonitorElement* meTrackEffPhiTot = igetter.get(folder_ + "TrackEffPhiTot"); + MonitorElement* meTrackEffPtTot = igetter.get(folder_ + "TrackEffPtTot"); + MonitorElement* meTrackEffEtaMtd = igetter.get(folder_ + "TrackEffEtaMtd"); + MonitorElement* meTrackEffPhiMtd = igetter.get(folder_ + "TrackEffPhiMtd"); + MonitorElement* meTrackEffPtMtd = igetter.get(folder_ + "TrackEffPtMtd"); + + if (!meTrackEffEtaTot || !meTrackEffPhiTot || !meTrackEffPtTot || !meTrackEffEtaMtd || !meTrackEffPhiMtd || + !meTrackEffPtMtd) { + edm::LogError("BtlRecoHarvester") << "Monitoring histograms not found!" << std::endl; + return; + } + + // --- Book histograms + ibook.cd(folder_); + meBtlEtaEff_ = ibook.book1D("BtlEtaEff", + " Track Efficiency VS Eta;#eta;Efficiency", + meTrackEffEtaTot->getNbinsX(), + meTrackEffEtaTot->getTH1()->GetXaxis()->GetXmin(), + meTrackEffEtaTot->getTH1()->GetXaxis()->GetXmax()); + meBtlPhiEff_ = ibook.book1D("BtlPhiEff", + "Track Efficiency VS Phi;#phi;Efficiency [rad]", + meTrackEffPhiTot->getNbinsX(), + meTrackEffPhiTot->getTH1()->GetXaxis()->GetXmin(), + meTrackEffPhiTot->getTH1()->GetXaxis()->GetXmax()); + meBtlPtEff_ = ibook.book1D("BtlPtEff", + "Track Efficiency VS Pt;Pt;Efficiency [GeV]", + meTrackEffPtTot->getNbinsX(), + meTrackEffPtTot->getTH1()->GetXaxis()->GetXmin(), + meTrackEffPtTot->getTH1()->GetXaxis()->GetXmax()); + meBtlEtaEff_->getTH1()->SetMinimum(0.); + meBtlPhiEff_->getTH1()->SetMinimum(0.); + meBtlPtEff_->getTH1()->SetMinimum(0.); + + // --- Calculate efficiency + for (int ibin = 1; ibin <= meTrackEffEtaTot->getNbinsX(); ibin++) { + double eff = meTrackEffEtaMtd->getBinContent(ibin) / meTrackEffEtaTot->getBinContent(ibin); + double bin_err = sqrt((meTrackEffEtaMtd->getBinContent(ibin) * + (meTrackEffEtaTot->getBinContent(ibin) - meTrackEffEtaMtd->getBinContent(ibin))) / + pow(meTrackEffEtaTot->getBinContent(ibin), 3)); + if (meTrackEffEtaTot->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meBtlEtaEff_->setBinContent(ibin, eff); + meBtlEtaEff_->setBinError(ibin, bin_err); + } + for (int ibin = 1; ibin <= meTrackEffPhiTot->getNbinsX(); ibin++) { + double eff = meTrackEffPhiMtd->getBinContent(ibin) / meTrackEffPhiTot->getBinContent(ibin); + double bin_err = sqrt((meTrackEffPhiMtd->getBinContent(ibin) * + (meTrackEffPhiTot->getBinContent(ibin) - meTrackEffPhiMtd->getBinContent(ibin))) / + pow(meTrackEffPhiTot->getBinContent(ibin), 3)); + if (meTrackEffPhiTot->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meBtlPhiEff_->setBinContent(ibin, eff); + meBtlPhiEff_->setBinError(ibin, bin_err); + } + for (int ibin = 1; ibin <= meTrackEffPtTot->getNbinsX(); ibin++) { + double eff = meTrackEffPtMtd->getBinContent(ibin) / meTrackEffPtTot->getBinContent(ibin); + double bin_err = sqrt((meTrackEffPtMtd->getBinContent(ibin) * + (meTrackEffPtTot->getBinContent(ibin) - meTrackEffPtMtd->getBinContent(ibin))) / + pow(meTrackEffPtTot->getBinContent(ibin), 3)); + if (meTrackEffPtTot->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meBtlPtEff_->setBinContent(ibin, eff); + meBtlPtEff_->setBinError(ibin, bin_err); + } +} + +// ------------ method fills 'descriptions' with the allowed parameters for the module ---------- +void BtlRecoHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("folder", "MTD/BTL/Reco/"); + + descriptions.add("btlRecoPostProcessor", desc); +} + +DEFINE_FWK_MODULE(BtlRecoHarvester); diff --git a/Validation/MtdValidation/plugins/BtlRecoValidation.cc b/Validation/MtdValidation/plugins/BtlRecoValidation.cc new file mode 100644 index 0000000000000..ae3e9bb85e9b6 --- /dev/null +++ b/Validation/MtdValidation/plugins/BtlRecoValidation.cc @@ -0,0 +1,219 @@ +#include + +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" + +#include "DQMServices/Core/interface/DQMEDAnalyzer.h" +#include "DQMServices/Core/interface/DQMStore.h" + +#include "DataFormats/Common/interface/ValidHandle.h" +#include "DataFormats/Math/interface/GeantUnits.h" +#include "DataFormats/ForwardDetId/interface/BTLDetId.h" +#include "DataFormats/FTLRecHit/interface/FTLRecHitCollections.h" +#include "DataFormats/FTLRecHit/interface/FTLClusterCollections.h" + +#include "DataFormats/Common/interface/Ptr.h" +#include "DataFormats/Common/interface/PtrVector.h" +#include "DataFormats/Common/interface/RefProd.h" +#include "DataFormats/Common/interface/Ref.h" +#include "DataFormats/Common/interface/RefVector.h" + +#include "DataFormats/TrackReco/interface/Track.h" +#include "DataFormats/TrackReco/interface/TrackFwd.h" +#include "DataFormats/GsfTrackReco/interface/GsfTrack.h" +#include "DataFormats/GsfTrackReco/interface/GsfTrackFwd.h" +#include "DataFormats/VertexReco/interface/Vertex.h" +#include "DataFormats/VertexReco/interface/VertexFwd.h" +#include "DataFormats/TrackerRecHit2D/interface/MTDTrackingRecHit.h" + +#include "Geometry/Records/interface/MTDDigiGeometryRecord.h" +#include "Geometry/Records/interface/MTDTopologyRcd.h" +#include "Geometry/MTDNumberingBuilder/interface/MTDTopology.h" +#include "Geometry/MTDGeometryBuilder/interface/MTDGeometry.h" +#include "Geometry/MTDGeometryBuilder/interface/ProxyMTDTopology.h" +#include "Geometry/MTDGeometryBuilder/interface/RectangularMTDTopology.h" + +class BtlRecoValidation : public DQMEDAnalyzer { +public: + explicit BtlRecoValidation(const edm::ParameterSet&); + ~BtlRecoValidation() override; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +private: + void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override; + + void analyze(const edm::Event&, const edm::EventSetup&) override; + + // ------------ member data ------------ + + const std::string folder_; + const float hitMinEnergy_; + const float trackMinEnergy_; + const float trackMinEta_; + + edm::EDGetTokenT btlRecCluToken_; + edm::EDGetTokenT btlRecTrackToken_; + edm::EDGetTokenT> RecVertexToken_; + + MonitorElement* meCluTime_; + MonitorElement* meCluEnergy_; + MonitorElement* meCluPhi_; + MonitorElement* meCluEta_; + MonitorElement* meCluHits_; + MonitorElement* meCluZvsPhi_; + MonitorElement* meTrackRPTime_; + MonitorElement* meTrackEffEtaTot_; + MonitorElement* meTrackEffPhiTot_; + MonitorElement* meTrackEffPtTot_; + MonitorElement* meTrackEffEtaMtd_; + MonitorElement* meTrackEffPhiMtd_; + MonitorElement* meTrackEffPtMtd_; + MonitorElement* meVerNumber_; + MonitorElement* meVerZ_; + MonitorElement* meVerTime_; +}; + +// ------------ constructor and destructor -------------- +BtlRecoValidation::BtlRecoValidation(const edm::ParameterSet& iConfig) + : folder_(iConfig.getParameter("folder")), + hitMinEnergy_(iConfig.getParameter("hitMinimumEnergy")), + trackMinEnergy_(iConfig.getParameter("trackMinimumEnergy")), + trackMinEta_(iConfig.getParameter("trackMinimumEta")) { + btlRecCluToken_ = consumes(iConfig.getParameter("inputTagC")); + btlRecTrackToken_ = consumes(iConfig.getParameter("inputTagT")); + RecVertexToken_ = consumes>(iConfig.getParameter("inputTagV")); +} + +BtlRecoValidation::~BtlRecoValidation() {} + +// ------------ method called for each event ------------ +void BtlRecoValidation::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) { + using namespace edm; + using namespace geant_units::operators; + using namespace std; + + edm::ESHandle geometryHandle; + iSetup.get().get(geometryHandle); + const MTDGeometry* geom = geometryHandle.product(); + + edm::ESHandle topologyHandle; + iSetup.get().get(topologyHandle); + + auto btlRecCluHandle = makeValid(iEvent.getHandle(btlRecCluToken_)); + auto btlRecTrackHandle = makeValid(iEvent.getHandle(btlRecTrackToken_)); + auto RecVertexHandle = makeValid(iEvent.getHandle(RecVertexToken_)); + + // --- Loop over the BTL RECO tracks --- + for (const auto& track : *btlRecTrackHandle) { + if (fabs(track.eta()) > trackMinEta_) + continue; + if (track.pt() < trackMinEnergy_) + continue; + + // --- all BTL tracks (with and without hit in MTD) --- + meTrackEffEtaTot_->Fill(track.eta()); + meTrackEffPhiTot_->Fill(track.phi()); + meTrackEffPtTot_->Fill(track.pt()); + + bool MTDBtl = false; + for (const auto hit : track.recHits()) { + MTDDetId Hit = hit->geographicalId(); + if ((Hit.det() == 6) && (Hit.subdetId() == 1) && (Hit.mtdSubDetector() == 1)) + MTDBtl = true; + } + + // --- keeping only tracks with last hit in MTD --- + if (MTDBtl == true) { + meTrackEffEtaMtd_->Fill(track.eta()); + meTrackEffPhiMtd_->Fill(track.phi()); + meTrackEffPtMtd_->Fill(track.pt()); + meTrackRPTime_->Fill(track.t0()); + } + } + + // --- Loop over the RECO vertices --- + int nv = 0; + for (const auto& v : *RecVertexHandle) { + if (v.isValid()) { + meVerZ_->Fill(v.z()); + meVerTime_->Fill(v.t()); + nv++; + } else + cout << "The vertex is not valid" << endl; + } + meVerNumber_->Fill(nv); + + // --- Loop over the BTL RECO clusters --- + for (const auto& DetSetClu : *btlRecCluHandle) { + for (const auto& cluster : DetSetClu) { + if (cluster.energy() < hitMinEnergy_) + continue; + BTLDetId cluId = cluster.id(); + DetId detIdObject(cluId); + const auto& genericDet = geom->idToDetUnit(detIdObject); + if (genericDet == nullptr) { + throw cms::Exception("BtlRecoValidation") + << "GeographicalID: " << std::hex << cluId << " is invalid!" << std::dec << std::endl; + } + + const ProxyMTDTopology& topoproxy = static_cast(genericDet->topology()); + const RectangularMTDTopology& topo = static_cast(topoproxy.specificTopology()); + + Local3DPoint local_point(cluster.x() * 5.7, cluster.y() * 0.3, 0.); + local_point = topo.pixelToModuleLocalPoint(local_point, cluId.row(topo.nrows()), cluId.column(topo.ncolumns())); + const auto& global_point = genericDet->toGlobal(local_point); + + meCluEnergy_->Fill(cluster.energy()); + meCluTime_->Fill(cluster.time()); + meCluPhi_->Fill(global_point.phi()); + meCluEta_->Fill(global_point.eta()); + meCluZvsPhi_->Fill(global_point.z(), global_point.phi()); + meCluHits_->Fill(cluster.size()); + } + } +} + +// ------------ method for histogram booking ------------ +void BtlRecoValidation::bookHistograms(DQMStore::IBooker& ibook, edm::Run const& run, edm::EventSetup const& iSetup) { + ibook.setCurrentFolder(folder_); + + // histogram booking + meCluTime_ = ibook.book1D("BtlCluTime", "BTL cluster time ToA;ToA [ns]", 250, 0, 25); + meCluEnergy_ = ibook.book1D("BtlCluEnergy", "BTL cluster energy;E_{RECO} [MeV]", 100, 0, 20); + meCluPhi_ = ibook.book1D("BtlCluPhi", "BTL cluster #phi;#phi_{RECO} [rad]", 100, -3.2, 3.2); + meCluEta_ = ibook.book1D("BtlCluEta", "BTL cluster #eta;#eta_{RECO}", 100, -1.6, 1.6); + meCluHits_ = ibook.book1D("BtlCluHitNumber", "BTL hits per cluster", 10, 0, 10); + meCluZvsPhi_ = ibook.book2D( + "BtlOccupancy", "BTL cluster Z vs #phi;Z_{RECO} [cm]; #phi_{RECO} [rad]", 100, -260., 260., 100, -3.2, 3.2); + meTrackEffEtaTot_ = ibook.book1D("TrackEffEtaTot", "Track efficiency vs eta D;#eta_{RECO}", 100, -1.6, 1.6); + meTrackEffPhiTot_ = ibook.book1D("TrackEffPhiTot", "Track efficiency vs phi D;#phi_{RECO} [rad]", 100, -3.2, 3.2); + meTrackEffPtTot_ = ibook.book1D("TrackEffPtTot", "Track efficiency vs pt D;pt_{RECO} [GeV]", 50, 0, 10); + meTrackEffEtaMtd_ = ibook.book1D("TrackEffEtaMtd", "Track efficiency vs eta N;#eta_{RECO}", 100, -1.6, 1.6); + meTrackEffPhiMtd_ = ibook.book1D("TrackEffPhiMtd", "Track efficiency vs phi N;#phi_{RECO} [rad]", 100, -3.2, 3.2); + meTrackEffPtMtd_ = ibook.book1D("TrackEffPtMtd", "Track efficiency vs pt N;pt_{RECO} [GeV]", 50, 0, 10); + meTrackRPTime_ = ibook.book1D("TrackRPTime", "Track t0 with respect to R.P.;t0 [ns]", 100, -10, 10); + meVerZ_ = ibook.book1D("VerZ", "RECO Vertex Z;Z_{RECO} [cm]", 180, -18, 18); + meVerTime_ = ibook.book1D("VerTime", "RECO Vertex Time;t0 [ns]", 100, -1, 1); + meVerNumber_ = ibook.book1D("VerNumber", "RECO Vertex Number", 100, 0, 500); +} + +// ------------ method fills 'descriptions' with the allowed parameters for the module ------------ + +void BtlRecoValidation::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("folder", "MTD/BTL/Reco"); + desc.add("inputTagC", edm::InputTag("mtdClusters", "FTLBarrel")); + desc.add("inputTagT", edm::InputTag("trackExtenderWithMTD", "")); + desc.add("inputTagV", edm::InputTag("offlinePrimaryVertices4D", "")); + desc.add("hitMinimumEnergy", 1.); // [MeV] + desc.add("trackMinimumEnergy", 1.); // [GeV] + desc.add("trackMinimumEta", 1.5); + + descriptions.add("btlReco", desc); +} + +DEFINE_FWK_MODULE(BtlRecoValidation); diff --git a/Validation/MtdValidation/plugins/EtlRecoHarvester.cc b/Validation/MtdValidation/plugins/EtlRecoHarvester.cc new file mode 100644 index 0000000000000..a5763815f2329 --- /dev/null +++ b/Validation/MtdValidation/plugins/EtlRecoHarvester.cc @@ -0,0 +1,189 @@ +#include + +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ServiceRegistry/interface/Service.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" + +#include "DQMServices/Core/interface/DQMEDHarvester.h" +#include "DQMServices/Core/interface/DQMStore.h" + +#include "DataFormats/ForwardDetId/interface/ETLDetId.h" + +class EtlRecoHarvester : public DQMEDHarvester { +public: + explicit EtlRecoHarvester(const edm::ParameterSet& iConfig); + ~EtlRecoHarvester() override; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override; + +private: + const std::string folder_; + + // --- Histograms + MonitorElement* meEtlEtaEff_[2]; + MonitorElement* meEtlPhiEff_[2]; + MonitorElement* meEtlPtEff_[2]; +}; + +// ------------ constructor and destructor -------------- +EtlRecoHarvester::EtlRecoHarvester(const edm::ParameterSet& iConfig) + : folder_(iConfig.getParameter("folder")) {} + +EtlRecoHarvester::~EtlRecoHarvester() {} + +// ------------ endjob tasks ---------------------------- +void EtlRecoHarvester::dqmEndJob(DQMStore::IBooker& ibook, DQMStore::IGetter& igetter) { + // --- Get the monitoring histograms + MonitorElement* meTrackEffEtaTotZneg = igetter.get(folder_ + "TrackEffEtaTotZneg"); + MonitorElement* meTrackEffPhiTotZneg = igetter.get(folder_ + "TrackEffPhiTotZneg"); + MonitorElement* meTrackEffPtTotZneg = igetter.get(folder_ + "TrackEffPtTotZneg"); + MonitorElement* meTrackEffEtaMtdZneg = igetter.get(folder_ + "TrackEffEtaMtdZneg"); + MonitorElement* meTrackEffPhiMtdZneg = igetter.get(folder_ + "TrackEffPhiMtdZneg"); + MonitorElement* meTrackEffPtMtdZneg = igetter.get(folder_ + "TrackEffPtMtdZneg"); + MonitorElement* meTrackEffEtaTotZpos = igetter.get(folder_ + "TrackEffEtaTotZpos"); + MonitorElement* meTrackEffPhiTotZpos = igetter.get(folder_ + "TrackEffPhiTotZpos"); + MonitorElement* meTrackEffPtTotZpos = igetter.get(folder_ + "TrackEffPtTotZpos"); + MonitorElement* meTrackEffEtaMtdZpos = igetter.get(folder_ + "TrackEffEtaMtdZpos"); + MonitorElement* meTrackEffPhiMtdZpos = igetter.get(folder_ + "TrackEffPhiMtdZpos"); + MonitorElement* meTrackEffPtMtdZpos = igetter.get(folder_ + "TrackEffPtMtdZpos"); + + if (!meTrackEffEtaTotZneg || !meTrackEffPhiTotZneg || !meTrackEffPtTotZneg || !meTrackEffEtaMtdZneg || + !meTrackEffPhiMtdZneg || !meTrackEffPtMtdZneg || !meTrackEffEtaTotZpos || !meTrackEffPhiTotZpos || + !meTrackEffPtTotZpos || !meTrackEffEtaMtdZpos || !meTrackEffPhiMtdZpos || !meTrackEffPtMtdZpos) { + edm::LogError("EtlRecoHarvester") << "Monitoring histograms not found!" << std::endl; + return; + } + + // --- Book histograms + ibook.cd(folder_); + meEtlEtaEff_[0] = ibook.book1D("EtlEtaEffZneg", + " Track Efficiency VS Eta (-Z);#eta;Efficiency", + meTrackEffEtaTotZneg->getNbinsX(), + meTrackEffEtaTotZneg->getTH1()->GetXaxis()->GetXmin(), + meTrackEffEtaTotZneg->getTH1()->GetXaxis()->GetXmax()); + meEtlPhiEff_[0] = ibook.book1D("EtlPhiEffZneg", + "Track Efficiency VS Phi (-Z);#phi [rad];Efficiency", + meTrackEffPhiTotZneg->getNbinsX(), + meTrackEffPhiTotZneg->getTH1()->GetXaxis()->GetXmin(), + meTrackEffPhiTotZneg->getTH1()->GetXaxis()->GetXmax()); + meEtlPtEff_[0] = ibook.book1D("EtlPtEffZneg", + "Track Efficiency VS Pt (-Z);Pt [GeV];Efficiency", + meTrackEffPtTotZneg->getNbinsX(), + meTrackEffPtTotZneg->getTH1()->GetXaxis()->GetXmin(), + meTrackEffPtTotZneg->getTH1()->GetXaxis()->GetXmax()); + meEtlEtaEff_[1] = ibook.book1D("EtlEtaEffZpos", + " Track Efficiency VS Eta (+Z);#eta;Efficiency", + meTrackEffEtaTotZpos->getNbinsX(), + meTrackEffEtaTotZpos->getTH1()->GetXaxis()->GetXmin(), + meTrackEffEtaTotZpos->getTH1()->GetXaxis()->GetXmax()); + meEtlPhiEff_[1] = ibook.book1D("EtlPhiEffZpos", + "Track Efficiency VS Phi (+Z);#phi [rad];Efficiency", + meTrackEffPhiTotZpos->getNbinsX(), + meTrackEffPhiTotZpos->getTH1()->GetXaxis()->GetXmin(), + meTrackEffPhiTotZpos->getTH1()->GetXaxis()->GetXmax()); + meEtlPtEff_[1] = ibook.book1D("EtlPtEffZpos", + "Track Efficiency VS Pt (+Z);Pt [GeV];Efficiency", + meTrackEffPtTotZpos->getNbinsX(), + meTrackEffPtTotZpos->getTH1()->GetXaxis()->GetXmin(), + meTrackEffPtTotZpos->getTH1()->GetXaxis()->GetXmax()); + meEtlEtaEff_[0]->getTH1()->SetMinimum(0.); + meEtlPhiEff_[0]->getTH1()->SetMinimum(0.); + meEtlPtEff_[0]->getTH1()->SetMinimum(0.); + meEtlEtaEff_[1]->getTH1()->SetMinimum(0.); + meEtlPhiEff_[1]->getTH1()->SetMinimum(0.); + meEtlPtEff_[1]->getTH1()->SetMinimum(0.); + + // --- Calculate efficiency + for (int ibin = 1; ibin <= meTrackEffEtaTotZneg->getNbinsX(); ibin++) { + double eff = meTrackEffEtaMtdZneg->getBinContent(ibin) / meTrackEffEtaTotZneg->getBinContent(ibin); + double bin_err = sqrt((meTrackEffEtaMtdZneg->getBinContent(ibin) * + (meTrackEffEtaTotZneg->getBinContent(ibin) - meTrackEffEtaMtdZneg->getBinContent(ibin))) / + pow(meTrackEffEtaTotZneg->getBinContent(ibin), 3)); + if (meTrackEffEtaTotZneg->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meEtlEtaEff_[0]->setBinContent(ibin, eff); + meEtlEtaEff_[0]->setBinError(ibin, bin_err); + } + + for (int ibin = 1; ibin <= meTrackEffEtaTotZpos->getNbinsX(); ibin++) { + double eff = meTrackEffEtaMtdZpos->getBinContent(ibin) / meTrackEffEtaTotZpos->getBinContent(ibin); + double bin_err = sqrt((meTrackEffEtaMtdZpos->getBinContent(ibin) * + (meTrackEffEtaTotZpos->getBinContent(ibin) - meTrackEffEtaMtdZpos->getBinContent(ibin))) / + pow(meTrackEffEtaTotZpos->getBinContent(ibin), 3)); + if (meTrackEffEtaTotZpos->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meEtlEtaEff_[1]->setBinContent(ibin, eff); + meEtlEtaEff_[1]->setBinError(ibin, bin_err); + } + + for (int ibin = 1; ibin <= meTrackEffPhiTotZneg->getNbinsX(); ibin++) { + double eff = meTrackEffPhiMtdZneg->getBinContent(ibin) / meTrackEffPhiTotZneg->getBinContent(ibin); + double bin_err = sqrt((meTrackEffPhiMtdZneg->getBinContent(ibin) * + (meTrackEffPhiTotZneg->getBinContent(ibin) - meTrackEffPhiMtdZneg->getBinContent(ibin))) / + pow(meTrackEffPhiTotZneg->getBinContent(ibin), 3)); + if (meTrackEffPhiTotZneg->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meEtlPhiEff_[0]->setBinContent(ibin, eff); + meEtlPhiEff_[0]->setBinError(ibin, bin_err); + } + + for (int ibin = 1; ibin <= meTrackEffPhiTotZpos->getNbinsX(); ibin++) { + double eff = meTrackEffPhiMtdZpos->getBinContent(ibin) / meTrackEffPhiTotZpos->getBinContent(ibin); + double bin_err = sqrt((meTrackEffPhiMtdZpos->getBinContent(ibin) * + (meTrackEffPhiTotZpos->getBinContent(ibin) - meTrackEffPhiMtdZpos->getBinContent(ibin))) / + pow(meTrackEffPhiTotZpos->getBinContent(ibin), 3)); + if (meTrackEffPhiTotZpos->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meEtlPhiEff_[1]->setBinContent(ibin, eff); + meEtlPhiEff_[1]->setBinError(ibin, bin_err); + } + + for (int ibin = 1; ibin <= meTrackEffPtTotZneg->getNbinsX(); ibin++) { + double eff = meTrackEffPtMtdZneg->getBinContent(ibin) / meTrackEffPtTotZneg->getBinContent(ibin); + double bin_err = sqrt((meTrackEffPtMtdZneg->getBinContent(ibin) * + (meTrackEffPtTotZneg->getBinContent(ibin) - meTrackEffPtMtdZneg->getBinContent(ibin))) / + pow(meTrackEffPtTotZneg->getBinContent(ibin), 3)); + if (meTrackEffPtTotZneg->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meEtlPtEff_[0]->setBinContent(ibin, eff); + meEtlPtEff_[0]->setBinError(ibin, bin_err); + } + + for (int ibin = 1; ibin <= meTrackEffPtTotZpos->getNbinsX(); ibin++) { + double eff = meTrackEffPtMtdZpos->getBinContent(ibin) / meTrackEffPtTotZpos->getBinContent(ibin); + double bin_err = sqrt((meTrackEffPtMtdZpos->getBinContent(ibin) * + (meTrackEffPtTotZpos->getBinContent(ibin) - meTrackEffPtMtdZpos->getBinContent(ibin))) / + pow(meTrackEffPtTotZpos->getBinContent(ibin), 3)); + if (meTrackEffPtTotZpos->getBinContent(ibin) == 0) { + eff = 0; + bin_err = 0; + } + meEtlPtEff_[1]->setBinContent(ibin, eff); + meEtlPtEff_[1]->setBinError(ibin, bin_err); + } +} + +// ------------ method fills 'descriptions' with the allowed parameters for the module ---------- +void EtlRecoHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("folder", "MTD/ETL/Reco/"); + + descriptions.add("etlRecoPostProcessor", desc); +} + +DEFINE_FWK_MODULE(EtlRecoHarvester); diff --git a/Validation/MtdValidation/plugins/EtlRecoValidation.cc b/Validation/MtdValidation/plugins/EtlRecoValidation.cc new file mode 100644 index 0000000000000..31d0086612be5 --- /dev/null +++ b/Validation/MtdValidation/plugins/EtlRecoValidation.cc @@ -0,0 +1,244 @@ +#include + +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" + +#include "DQMServices/Core/interface/DQMEDAnalyzer.h" +#include "DQMServices/Core/interface/DQMStore.h" + +#include "DataFormats/Common/interface/ValidHandle.h" +#include "DataFormats/Math/interface/GeantUnits.h" +#include "DataFormats/ForwardDetId/interface/ETLDetId.h" +#include "DataFormats/FTLRecHit/interface/FTLRecHitCollections.h" +#include "DataFormats/FTLRecHit/interface/FTLClusterCollections.h" + +#include "DataFormats/Common/interface/Ptr.h" +#include "DataFormats/Common/interface/PtrVector.h" +#include "DataFormats/Common/interface/RefProd.h" +#include "DataFormats/Common/interface/Ref.h" +#include "DataFormats/Common/interface/RefVector.h" + +#include "DataFormats/TrackReco/interface/Track.h" +#include "DataFormats/TrackReco/interface/TrackFwd.h" +#include "DataFormats/GsfTrackReco/interface/GsfTrack.h" +#include "DataFormats/GsfTrackReco/interface/GsfTrackFwd.h" +#include "DataFormats/TrackerRecHit2D/interface/MTDTrackingRecHit.h" + +#include "Geometry/Records/interface/MTDDigiGeometryRecord.h" +#include "Geometry/Records/interface/MTDTopologyRcd.h" +#include "Geometry/MTDNumberingBuilder/interface/MTDTopology.h" +#include "Geometry/MTDGeometryBuilder/interface/MTDGeometry.h" +#include "Geometry/MTDGeometryBuilder/interface/ProxyMTDTopology.h" +#include "Geometry/MTDGeometryBuilder/interface/RectangularMTDTopology.h" + +class EtlRecoValidation : public DQMEDAnalyzer { +public: + explicit EtlRecoValidation(const edm::ParameterSet&); + ~EtlRecoValidation() override; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +private: + void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override; + + void analyze(const edm::Event&, const edm::EventSetup&) override; + + // ------------ member data ------------ + + const std::string folder_; + const float hitMinEnergy_; + const float trackMinEnergy_; + const float trackMinEta_; + const float trackMaxEta_; + + edm::EDGetTokenT etlRecCluToken_; + edm::EDGetTokenT etlRecTrackToken_; + + MonitorElement* meCluTime_[2]; + MonitorElement* meCluEnergy_[2]; + MonitorElement* meCluPhi_[2]; + MonitorElement* meCluEta_[2]; + MonitorElement* meCluHits_[2]; + MonitorElement* meCluZvsPhi_[2]; + MonitorElement* meTrackRPTime_[2]; + MonitorElement* meTrackEffEtaTot_[2]; + MonitorElement* meTrackEffPhiTot_[2]; + MonitorElement* meTrackEffPtTot_[2]; + MonitorElement* meTrackEffEtaMtd_[2]; + MonitorElement* meTrackEffPhiMtd_[2]; + MonitorElement* meTrackEffPtMtd_[2]; +}; + +// ------------ constructor and destructor -------------- +EtlRecoValidation::EtlRecoValidation(const edm::ParameterSet& iConfig) + : folder_(iConfig.getParameter("folder")), + hitMinEnergy_(iConfig.getParameter("hitMinimumEnergy")), + trackMinEnergy_(iConfig.getParameter("trackMinimumEnergy")), + trackMinEta_(iConfig.getParameter("trackMinimumEta")), + trackMaxEta_(iConfig.getParameter("trackMaximumEta")) { + etlRecCluToken_ = consumes(iConfig.getParameter("inputTagC")); + etlRecTrackToken_ = consumes(iConfig.getParameter("inputTagT")); +} + +EtlRecoValidation::~EtlRecoValidation() {} + +// ------------ method called for each event ------------ +void EtlRecoValidation::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) { + using namespace edm; + using namespace geant_units::operators; + using namespace std; + + edm::ESHandle geometryHandle; + iSetup.get().get(geometryHandle); + const MTDGeometry* geom = geometryHandle.product(); + + edm::ESHandle topologyHandle; + iSetup.get().get(topologyHandle); + + auto etlRecCluHandle = makeValid(iEvent.getHandle(etlRecCluToken_)); + auto etlRecTrackHandle = makeValid(iEvent.getHandle(etlRecTrackToken_)); + + // --- Loop over the ETL RECO tracks --- + for (const auto& track : *etlRecTrackHandle) { + if (track.pt() < trackMinEnergy_) + continue; + + // --- all ETL tracks (with and without hit in MTD) --- + if ((track.eta() < -trackMinEta_) && (track.eta() > -trackMaxEta_)) { + meTrackEffEtaTot_[0]->Fill(track.eta()); + meTrackEffPhiTot_[0]->Fill(track.phi()); + meTrackEffPtTot_[0]->Fill(track.pt()); + } + + if ((track.eta() > trackMinEta_) && (track.eta() < trackMaxEta_)) { + meTrackEffEtaTot_[1]->Fill(track.eta()); + meTrackEffPhiTot_[1]->Fill(track.phi()); + meTrackEffPtTot_[1]->Fill(track.pt()); + } + + bool MTDEtlZneg = false; + bool MTDEtlZpos = false; + for (const auto hit : track.recHits()) { + MTDDetId Hit = hit->geographicalId(); + if ((Hit.det() == 6) && (Hit.subdetId() == 1) && (Hit.mtdSubDetector() == 2) && (Hit.zside() == -1)) + MTDEtlZneg = true; + if ((Hit.det() == 6) && (Hit.subdetId() == 1) && (Hit.mtdSubDetector() == 2) && (Hit.zside() == 1)) + MTDEtlZpos = true; + } + + // --- keeping only tracks with last hit in MTD --- + if ((track.eta() < -trackMinEta_) && (track.eta() > -trackMaxEta_)) { + if (MTDEtlZneg == true) { + meTrackEffEtaMtd_[0]->Fill(track.eta()); + meTrackEffPhiMtd_[0]->Fill(track.phi()); + meTrackEffPtMtd_[0]->Fill(track.pt()); + meTrackRPTime_[0]->Fill(track.t0()); + } + } + if ((track.eta() > trackMinEta_) && (track.eta() < trackMaxEta_)) { + if (MTDEtlZpos == true) { + meTrackEffEtaMtd_[1]->Fill(track.eta()); + meTrackEffPhiMtd_[1]->Fill(track.phi()); + meTrackEffPtMtd_[1]->Fill(track.pt()); + meTrackRPTime_[1]->Fill(track.t0()); + } + } + } + + // --- Loop over the ETL RECO clusters --- + for (const auto& DetSetClu : *etlRecCluHandle) { + for (const auto& cluster : DetSetClu) { + if (cluster.energy() < hitMinEnergy_) + continue; + ETLDetId cluId = cluster.id(); + DetId detIdObject(cluId); + const auto& genericDet = geom->idToDetUnit(detIdObject); + if (genericDet == nullptr) { + throw cms::Exception("EtlRecoValidation") + << "GeographicalID: " << std::hex << cluId << " is invalid!" << std::dec << std::endl; + } + + const ProxyMTDTopology& topoproxy = static_cast(genericDet->topology()); + const RectangularMTDTopology& topo = static_cast(topoproxy.specificTopology()); + + Local3DPoint local_point(topo.localX(cluster.x()), topo.localY(cluster.y()), 0.); + const auto& global_point = genericDet->toGlobal(local_point); + + int idet = (cluId.zside() + 1) / 2; + + meCluEnergy_[idet]->Fill(cluster.energy()); + meCluTime_[idet]->Fill(cluster.time()); + meCluPhi_[idet]->Fill(global_point.phi()); + meCluEta_[idet]->Fill(global_point.eta()); + meCluZvsPhi_[idet]->Fill(global_point.x(), global_point.y()); + meCluHits_[idet]->Fill(cluster.size()); + } + } +} + +// ------------ method for histogram booking ------------ +void EtlRecoValidation::bookHistograms(DQMStore::IBooker& ibook, edm::Run const& run, edm::EventSetup const& iSetup) { + ibook.setCurrentFolder(folder_); + + // histogram booking + meCluTime_[0] = ibook.book1D("EtlCluTimeZneg", "ETL cluster ToA (-Z);ToA [ns]", 250, 0, 25); + meCluTime_[1] = ibook.book1D("EtlCluTimeZpos", "ETL cluster ToA (+Z);ToA [ns]", 250, 0, 25); + meCluEnergy_[0] = ibook.book1D("EtlCluEnergyZneg", "ETL cluster energy (-Z);E_{RECO} [MeV]", 100, 0, 10); + meCluEnergy_[1] = ibook.book1D("EtlCluEnergyZpos", "ETL cluster energy (+Z);E_{RECO} [MeV]", 100, 0, 10); + meCluPhi_[0] = ibook.book1D("EtlCluPhiZneg", "ETL cluster #phi (-Z);#phi_{RECO} [rad]", 100, -3.2, 3.2); + meCluPhi_[1] = ibook.book1D("EtlCluPhiZpos", "ETL cluster #phi (+Z);#phi_{RECO} [rad]", 100, -3.2, 3.2); + meCluEta_[0] = ibook.book1D("EtlCluEtaZneg", "ETL cluster #eta (-Z);#eta_{RECO}", 100, -3.2, -1.4); + meCluEta_[1] = ibook.book1D("EtlCluEtaZpos", "ETL cluster #eta (+Z);#eta_{RECO}", 100, 1.4, 3.2); + meCluHits_[0] = ibook.book1D("EtlCluHitNumberZneg", "ETL hits per cluster (-Z)", 10, 0, 10); + meCluHits_[1] = ibook.book1D("EtlCluHitNumberZpos", "ETL hits per cluster (+Z)", 10, 0, 10); + meCluZvsPhi_[0] = ibook.book2D( + "EtlOccupancyZneg", "ETL cluster X vs Y (-Z);X_{RECO} [cm]; Y_{RECO} [cm]", 100, -150., 150., 100, -150, 150); + meCluZvsPhi_[1] = ibook.book2D( + "EtlOccupancyZpos", "ETL cluster X vs Y (+Z);X_{RECO} [cm]; Y_{RECO} [cm]", 100, -150., 150., 100, -150, 150); + meTrackEffEtaTot_[0] = + ibook.book1D("TrackEffEtaTotZneg", "Track efficiency vs eta (Tot) (-Z);#eta_{RECO}", 100, -3.2, -1.4); + meTrackEffEtaTot_[1] = + ibook.book1D("TrackEffEtaTotZpos", "Track efficiency vs eta (Tot) (+Z);#eta_{RECO}", 100, 1.4, 3.2); + meTrackEffPhiTot_[0] = + ibook.book1D("TrackEffPhiTotZneg", "Track efficiency vs phi (Tot) (-Z);#phi_{RECO} [rad]", 100, -3.2, 3.2); + meTrackEffPhiTot_[1] = + ibook.book1D("TrackEffPhiTotZpos", "Track efficiency vs phi (Tot) (+Z);#phi_{RECO} [rad]", 100, -3.2, 3.2); + meTrackEffPtTot_[0] = + ibook.book1D("TrackEffPtTotZneg", "Track efficiency vs pt (Tot) (-Z);pt_{RECO} [GeV]", 50, 0, 10); + meTrackEffPtTot_[1] = + ibook.book1D("TrackEffPtTotZpos", "Track efficiency vs pt (Tot) (+Z);pt_{RECO} [GeV]", 50, 0, 10); + meTrackEffEtaMtd_[0] = + ibook.book1D("TrackEffEtaMtdZneg", "Track efficiency vs eta (Mtd) (-Z);#eta_{RECO}", 100, -3.2, -1.4); + meTrackEffEtaMtd_[1] = + ibook.book1D("TrackEffEtaMtdZpos", "Track efficiency vs eta (Mtd) (+Z);#eta_{RECO}", 100, 1.4, 3.2); + meTrackEffPhiMtd_[0] = + ibook.book1D("TrackEffPhiMtdZneg", "Track efficiency vs phi (Mtd) (-Z);#phi_{RECO} [rad]", 100, -3.2, 3.2); + meTrackEffPhiMtd_[1] = + ibook.book1D("TrackEffPhiMtdZpos", "Track efficiency vs phi (Mtd) (+Z);#phi_{RECO} [rad]", 100, -3.2, 3.2); + meTrackEffPtMtd_[0] = + ibook.book1D("TrackEffPtMtdZneg", "Track efficiency vs pt (Mtd) (-Z);pt_{RECO} [GeV]", 50, 0, 10); + meTrackEffPtMtd_[1] = + ibook.book1D("TrackEffPtMtdZpos", "Track efficiency vs pt (Mtd) (+Z);pt_{RECO} [GeV]", 50, 0, 10); + meTrackRPTime_[0] = ibook.book1D("TrackRPTimeZneg", "Track t0 with respect to R.P. (-Z);t0 [ns]", 100, -10, 10); + meTrackRPTime_[1] = ibook.book1D("TrackRPTimeZpos", "Track t0 with respect to R.P. (+Z);t0 [ns]", 100, -10, 10); +} + +// ------------ method fills 'descriptions' with the allowed parameters for the module ------------ + +void EtlRecoValidation::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("folder", "MTD/ETL/Reco"); + desc.add("inputTagC", edm::InputTag("mtdClusters", "FTLEndcap")); + desc.add("inputTagT", edm::InputTag("trackExtenderWithMTD", "")); + desc.add("hitMinimumEnergy", 1.); // [MeV] + desc.add("trackMinimumEnergy", 0.5); // [GeV] + desc.add("trackMinimumEta", 1.4); + desc.add("trackMaximumEta", 3.2); + + descriptions.add("etlReco", desc); +} + +DEFINE_FWK_MODULE(EtlRecoValidation); diff --git a/Validation/MtdValidation/python/MtdPostProcessor_cff.py b/Validation/MtdValidation/python/MtdPostProcessor_cff.py index 10df413cb3c72..f4caa1acc9b40 100644 --- a/Validation/MtdValidation/python/MtdPostProcessor_cff.py +++ b/Validation/MtdValidation/python/MtdPostProcessor_cff.py @@ -1,5 +1,7 @@ import FWCore.ParameterSet.Config as cms from Validation.MtdValidation.btlSimHitsPostProcessor_cfi import btlSimHitsPostProcessor +from Validation.MtdValidation.btlRecoPostProcessor_cfi import btlRecoPostProcessor +from Validation.MtdValidation.etlRecoPostProcessor_cfi import etlRecoPostProcessor -mtdValidationPostProcessor = cms.Sequence(btlSimHitsPostProcessor) +mtdValidationPostProcessor = cms.Sequence(btlSimHitsPostProcessor + btlRecoPostProcessor + etlRecoPostProcessor)