diff --git a/DPGAnalysis/MuonTools/BuildFile.xml b/DPGAnalysis/MuonTools/BuildFile.xml new file mode 100644 index 0000000000000..13cf66377b213 --- /dev/null +++ b/DPGAnalysis/MuonTools/BuildFile.xml @@ -0,0 +1,32 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/DPGAnalysis/MuonTools/plugins/BuildFile.xml b/DPGAnalysis/MuonTools/plugins/BuildFile.xml new file mode 100644 index 0000000000000..ccd1c6c40c741 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/BuildFile.xml @@ -0,0 +1,18 @@ + + + + + + + + + + + + + + + + + + diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.cc new file mode 100644 index 0000000000000..15cc09dd7d2b9 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.cc @@ -0,0 +1,175 @@ +/** \class MuNtupleBmtfFiller MuNtupleBmtfFiller.cc DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.cc + * + * Helper class : the BMTF filler + * + * \author L. Borgonovi (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.h" + +#include + +MuNtupleBmtfFiller::MuNtupleBmtfFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller(config), + m_tpgPhiToken{config, consumesCollector(), "dtTpTag"}, + m_bmtfToken{config, consumesCollector(), "bmtfTag"} { + produces(); +} + +void MuNtupleBmtfFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "l1tBmtfOut"); + desc.add("dtTpTag", edm::InputTag{"bmtfDigis"}); + desc.add("bmtfTag", edm::InputTag{"bmtfDigis"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleBmtfFiller::fill(edm::Event& ev) { + unsigned int nBmtfCands{0}; + + std::vector wheel; + std::vector sector; + + std::vector pt; + std::vector phi; + std::vector eta; + std::vector bx; + std::vector dxy; + std::vector qual; + std::vector etaFine; + + std::vector matchedTpIdxMB1; + std::vector matchedTpIdxMB2; + std::vector matchedTpIdxMB3; + std::vector matchedTpIdxMB4; + + auto bmtfColl = m_bmtfToken.conditionalGet(ev); + auto tpColl = m_tpgPhiToken.conditionalGet(ev); + + if (bmtfColl.isValid()) { + auto bmtfCandBX = bmtfColl->getFirstBX(); + auto bmtfCandLastBX = bmtfColl->getLastBX(); + + for (; bmtfCandBX <= bmtfCandLastBX; ++bmtfCandBX) { + auto bmtfCand = bmtfColl->begin(bmtfCandBX); + auto bmtfCandLast = bmtfColl->end(bmtfCandBX); + + for (; bmtfCand != bmtfCandLast; ++bmtfCand) { + std::map mapTA = bmtfCand->trackAddress(); + + int wsign = mapTA[0] == 0 ? 1 : -1; + int wh = wsign * mapTA[1]; + int sec = bmtfCand->processor(); + + wheel.push_back(wh); + sector.push_back(sec); + + std::array ts_mb{{DEFAULT_INT_VAL, DEFAULT_INT_VAL, DEFAULT_INT_VAL, DEFAULT_INT_VAL}}; + std::array w_mb{{DEFAULT_INT_VAL, DEFAULT_INT_VAL, DEFAULT_INT_VAL, DEFAULT_INT_VAL}}; + std::array s_mb{{DEFAULT_INT_VAL, DEFAULT_INT_VAL, DEFAULT_INT_VAL, DEFAULT_INT_VAL}}; + + std::array triggerIndex{{0, 0, 0, 0}}; + + if (mapTA[2] != 3) { + w_mb[0] = wh; + s_mb[0] = sec; + ts_mb[0] = (mapTA[2] & 1) ? 1 : 0; + } + + for (int iSt = 1; iSt <= N_STAT; ++iSt) { + if (mapTA[iSt + 2] != 15) { + ts_mb[iSt] = (mapTA[iSt + 2] & 1) ? 1 : 0; // 0 for ts1 , 1 for ts2 + w_mb[iSt] = + (mapTA[iSt + 2] & 8) + ? wh + : wh + + wsign; // own wheel if true, nex wheel if false (depends on sign of the wheel -> case 0+, +1: +1, case 0-, -1: -1) + int tmpMapValue = mapTA[iSt + 2] >> 1; // temp value to remove less significant bit + if (tmpMapValue & 1) + s_mb[iSt] = + sec != 0 ? sec - 1 : 11; // if last two remained bit == 01 -> sector-1 (N+1 column of ref table) + else if (tmpMapValue & 2) + s_mb[iSt] = + sec != 11 ? sec + 1 : 0; // if last two remained bit == 10 -> sector+1 (N-1 column of ref table) + else + s_mb[iSt] = sec; // if last two remained bit == 00 -> sector (N column of ref table) + } + } + + int tmpIdx = 0; + + for (int iSt = 0; iSt < N_STAT; ++iSt) { + int iTP = 0; + + if (tpColl.isValid()) { + const auto trigs = tpColl->getContainer(); + + for (const auto& trig : (*trigs)) { + if (trig.code() != 7) { + if (bmtfCandBX == trig.bxNum()) { + if ((w_mb[iSt] == trig.whNum()) && (ts_mb[iSt] == trig.Ts2Tag()) && (s_mb[iSt] == trig.scNum()) && + (iSt + 1 == trig.stNum())) { + triggerIndex[tmpIdx] = iTP; + tmpIdx++; + } + } + } + + iTP++; + } + } + } + + int iPhi = bmtfCand->hwPhi() + bmtfCand->processor() * 48 - 15; + if (iPhi < 0) + iPhi += 576; + + bx.push_back(bmtfCandBX); + pt.push_back((bmtfCand->hwPt()) * PT_SCALE); + phi.push_back(iPhi * PHI_SCALE); //no conversion yet + eta.push_back((bmtfCand->hwEta()) * ETA_SCALE); + dxy.push_back(bmtfCand->hwDXY()); + qual.push_back(bmtfCand->hwQual()); + etaFine.push_back(bmtfCand->hwHF()); + + matchedTpIdxMB1.push_back(triggerIndex[0]); + matchedTpIdxMB2.push_back(triggerIndex[1]); + matchedTpIdxMB3.push_back(triggerIndex[2]); + matchedTpIdxMB4.push_back(triggerIndex[3]); + + ++nBmtfCands; + } + } + } + + auto table = std::make_unique(nBmtfCands, m_label, false, false); + + table->setDoc("BMTF information"); + + addColumn(table, "pt", pt, "BMTF cand pt - GeV/c "); + addColumn(table, "phi", phi, "BMTF cand phi - rad"); + addColumn(table, "eta", eta, "BMTF cand eta"); + + addColumn(table, "wheel", pt, "BMTF cand wheel"); + addColumn(table, "sector ", pt, "BMTF cand sector"); + + addColumn(table, "dxy", dxy, "BMTF cand dxy - units?"); + addColumn(table, "qual", qual, "BMTF cand quality"); + addColumn(table, "etaFine", etaFine, "BMTF cand fine eta bit"); + + addColumn(table, "matchedTpIdxMB1", matchedTpIdxMB1, "BMTF link to TP index (MB1)"); + addColumn(table, "matchedTpIdxMB2", matchedTpIdxMB2, "BMTF link to TP index (MB2)"); + addColumn(table, "matchedTpIdxMB3", matchedTpIdxMB3, "BMTF link to TP index (MB3)"); + addColumn(table, "matchedTpIdxMB4", matchedTpIdxMB4, "BMTF link to TP index (MB4)"); + + ev.put(std::move(table)); +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleBmtfFiller); \ No newline at end of file diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.h new file mode 100644 index 0000000000000..d9b5099d5936d --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.h @@ -0,0 +1,51 @@ +#ifndef MuNtuple_MuNtupleBmtfFiller_h +#define MuNtuple_MuNtupleBmtfFiller_h + +/** \class MuNtupleBmtfFiller MuNtupleBmtfFiller.h DPGAnalysis/MuonTools/plugins/MuNtupleBmtfFiller.h + * + * Helper class : the BMTF filler + * + * \author L. Borgonovi (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambPhContainer.h" +#include "DataFormats/L1TMuon/interface/RegionalMuonCand.h" + +class MuNtupleBmtfFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleBmtfFiller(const edm::ParameterSet &); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); + +protected: + /// Fill tree branches for a given events + void fill(edm::Event &ev) final; + +private: + /// Scale to convert HW pt to GeV + static constexpr double PT_SCALE = 0.5; + + /// Scale to convert HW eta to phisical value + static constexpr double ETA_SCALE = 0.010875; + + /// Scale to convert HW eta to rad + static constexpr double PHI_SCALE = 0.010908308; // 2 * pi / 576 + + /// Number of DT stations + static constexpr int N_STAT = 4; + + /// The trigger tokens + nano_mu::EDTokenHandle m_tpgPhiToken; + nano_mu::EDTokenHandle m_bmtfToken; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.cc new file mode 100644 index 0000000000000..9ca13cc2056c2 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.cc @@ -0,0 +1,98 @@ +/** \class MuNtupleDigiFiller MuNtupleDigiFiller.cc DPGAnalysis/MuonTools/src/MuNtupleDTDigiFiller.cc + * + * Helper class : the digi filler for Phase-1 / Phase2 DT digis (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.h" + +#include + +MuNtupleDTDigiFiller::MuNtupleDTDigiFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller{config}, m_token{config, consumesCollector(), "dtDigiTag"} { + produces(); +} + +void MuNtupleDTDigiFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "dtDigis"); + desc.add("dtDigiTag", edm::InputTag{"muonDTDigis"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleDTDigiFiller::fill(edm::Event& ev) { + unsigned int nDigis{0}; + + std::vector wheel; + std::vector sector; + std::vector station; + + std::vector superLayer; + std::vector layer; + std::vector wire; + + std::vector time; + + auto dtDigis = m_token.conditionalGet(ev); + + if (dtDigis.isValid()) { + auto dtLayerIdIt = dtDigis->begin(); + auto dtLayerIdEnd = dtDigis->end(); + + for (; dtLayerIdIt != dtLayerIdEnd; ++dtLayerIdIt) { + const auto& [dtLayerId, range] = (*dtLayerIdIt); + + for (auto digiIt = range.first; digiIt != range.second; ++digiIt) { + wheel.push_back(dtLayerId.wheel()); + sector.push_back(dtLayerId.sector()); + station.push_back(dtLayerId.station()); + + superLayer.push_back(dtLayerId.superLayer()); + layer.push_back(dtLayerId.layer()); + wire.push_back(digiIt->wire()); + + time.push_back(digiIt->time()); + + ++nDigis; + } + } + } + + auto table = std::make_unique(nDigis, m_label, false, false); + + table->setDoc("DT digi information"); + + addColumn(table, "wheel", wheel, "wheel - [-2:2] range"); + addColumn(table, + "sector", + sector, + "sector - [1:14] range" + "
sector 13 used for the second MB4 of sector 4" + "
sector 14 used for the second MB4 of sector 10"); + addColumn(table, "station", station, "station - [1:4] range"); + addColumn(table, + "superLayer", + superLayer, + "superlayer - [1:3] range" + "
SL 1 and 3 are phi SLs" + "
SL 2 is theta SL"); + addColumn(table, "layer", layer, "station - [1:4] range"); + addColumn(table, + "wire", + wire, + "wire - [1:X] range" + "
X varies for different chambers SLs and layers"); + addColumn(table, "time", time, "digi time in ns (no pedestal subtraction)"); + + ev.put(std::move(table)); +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleDTDigiFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.h new file mode 100644 index 0000000000000..5b24cf433c653 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.h @@ -0,0 +1,37 @@ +#ifndef MuNtuple_MuNtupleDTDigiFiller_h +#define MuNtuple_MuNtupleDTDigiFiller_h + +/** \class MuNtupleDTDigiFiller MuNtupleDTDigiFiller.h DPGAnalysis/MuonTools/plugins/MuNtupleDTDigiFiller.h + * + * Helper class : the digi filler for Phase-1 / Phase2 DT digis (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/DTDigi/interface/DTDigiCollection.h" + +class MuNtupleDTDigiFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleDTDigiFiller(const edm::ParameterSet&); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + /// Fill tree branches for a given event + void fill(edm::Event& ev) final; + +private: + /// The digi token + nano_mu::EDTokenHandle m_token; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTSegmentFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleDTSegmentFiller.cc new file mode 100644 index 0000000000000..2331888736445 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTSegmentFiller.cc @@ -0,0 +1,429 @@ +/** \class MuNtupleDTSegmentFiller MuNtupleDTSegmentFiller.cc DTDPGAnalysis/MuonDPGNtuples/src/MuNtupleDTSegmentFiller.cc + * + * Helper class : the segment filler for Phase-1 / Phase2 segments (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleDTSegmentFiller.h" + +#include "TrackingTools/GeomPropagators/interface/StraightLinePlaneCrossing.h" +#include "TrackingTools/GeomPropagators/interface/Propagator.h" + +#include "CalibMuon/DTDigiSync/interface/DTTTrigSyncFactory.h" + +#include + +MuNtupleDTSegmentFiller::MuNtupleDTSegmentFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller{config}, + m_token{config, consumesCollector(), "dtSegmentTag"}, + m_fillHits{config.getParameter("fillHits")}, + m_fillExtr{config.getParameter("fillExtrapolation")}, + m_dtGeometry{consumesCollector()}, + m_trackingGeometry{consumesCollector()} { + produces(); + if (m_fillHits) { + produces("hits"); + } + if (m_fillExtr) { + produces("extr"); + } + + m_dtSync = DTTTrigSyncFactory::get()->create(config.getParameter("tTrigMode"), + config.getParameter("tTrigModeConfig"), + consumesCollector()); +} + +void MuNtupleDTSegmentFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "dtSegments"); + desc.add("dtSegmentTag", edm::InputTag{"dt4DSegments"}); + desc.add("fillExtrapolation", true); + desc.add("fillHits", true); + + desc.add("tTrigMode", "DTTTrigSyncFromDB"); + + edm::ParameterSetDescription tTrigModeParams; + tTrigModeParams.add("doTOFCorrection", true); + tTrigModeParams.add("tofCorrType", 2); + + tTrigModeParams.add("vPropWire", 24.4); + tTrigModeParams.add("doWirePropCorrection", true); + tTrigModeParams.add("wirePropCorrType", 0); + + tTrigModeParams.add("tTrigLabel", ""); + tTrigModeParams.add("doT0Correction", true); + tTrigModeParams.add("t0Label", ""); + tTrigModeParams.addUntracked("debug", false); + + desc.add("tTrigModeConfig", tTrigModeParams); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleDTSegmentFiller::getFromES(const edm::Run& run, const edm::EventSetup& environment) { + m_dtGeometry.getFromES(environment); +} + +void MuNtupleDTSegmentFiller::getFromES(const edm::EventSetup& environment) { + m_trackingGeometry.getFromES(environment); + m_dtSync->setES(environment); +} + +void MuNtupleDTSegmentFiller::fill(edm::Event& ev) { + unsigned int nSegments{0}; + + std::vector seg4D_wheel; + std::vector seg4D_sector; + std::vector seg4D_station; + + std::vector seg4D_hasPhi; + std::vector seg4D_hasZed; + + std::vector seg4D_posLoc_x; + std::vector seg4D_posLoc_y; + std::vector seg4D_posLoc_z; + std::vector seg4D_dirLoc_x; + std::vector seg4D_dirLoc_y; + std::vector seg4D_dirLoc_z; + + std::vector seg4D_posLoc_x_SL1; + std::vector seg4D_posLoc_x_SL3; + std::vector seg4D_posLoc_x_midPlane; + + std::vector seg4D_posGlb_phi; + std::vector seg4D_posGlb_eta; + std::vector seg4D_dirGlb_phi; + std::vector seg4D_dirGlb_eta; + + std::vector seg4D_extr_begin; + std::vector seg4D_extr_end; + + std::vector seg2D_phi_t0; + std::vector seg2D_phi_vDrift; + std::vector seg2D_phi_normChi2; + std::vector seg2D_phi_nHits; + + std::vector seg2D_z_normChi2; + std::vector seg2D_z_nHits; + + std::vector seg2D_hits_begin; + std::vector seg2D_hits_end; + + // segment extrapolation to DT layers filled, if m_fillExtr == true + unsigned int nExtr{0}; + + std::vector seg4D_hitsExpPos; + std::vector seg4D_hitsExpPosCh; + std::vector seg4D_hitsExpWire; + + // rec-hits vectors, filled if m_fillHits == true + unsigned int nHits{0}; + + std::vector seg2D_hits_pos; + std::vector seg2D_hits_posCh; + std::vector seg2D_hits_posErr; + std::vector seg2D_hits_side; + std::vector seg2D_hits_wire; + std::vector seg2D_hits_wirePos; + std::vector seg2D_hits_layer; + std::vector seg2D_hits_superLayer; + std::vector seg2D_hits_time; + std::vector seg2D_hits_timeCali; + + auto fillHits = [&](const DTRecSegment2D* seg, const GeomDet* chamb) { + const auto& recHits = seg->specificRecHits(); + + for (const auto& recHit : recHits) { + ++nHits; + + const auto wireId = recHit.wireId(); + const auto layerId = wireId.layerId(); + const auto* layer = m_dtGeometry->layer(layerId); + + seg2D_hits_pos.push_back(recHit.localPosition().x()); + seg2D_hits_posCh.push_back(chamb->toLocal(layer->toGlobal(recHit.localPosition())).x()); + seg2D_hits_posErr.push_back(recHit.localPositionError().xx()); + + seg2D_hits_side.push_back(recHit.lrSide()); + seg2D_hits_wire.push_back(wireId.wire()); + seg2D_hits_wirePos.push_back(layer->specificTopology().wirePosition(wireId.wire())); + seg2D_hits_layer.push_back(layerId.layer()); + seg2D_hits_superLayer.push_back(layerId.superlayer()); + + auto digiTime = recHit.digiTime(); + + auto tTrig = m_dtSync->offset(wireId); + + seg2D_hits_time.push_back(digiTime); + seg2D_hits_timeCali.push_back(digiTime - tTrig); + } + }; + + auto segments4D = m_token.conditionalGet(ev); + + if (segments4D.isValid()) { + auto chambIt = segments4D->id_begin(); + const auto chambEnd = segments4D->id_end(); + + for (; chambIt != chambEnd; ++chambIt) { + const auto& range = segments4D->get(*chambIt); + + for (auto segment4D = range.first; segment4D != range.second; ++segment4D) { + auto wheel = (*chambIt).wheel(); + auto sector = (*chambIt).sector(); + auto station = (*chambIt).station(); + + seg4D_wheel.push_back(wheel); + seg4D_sector.push_back(sector); + seg4D_station.push_back(station); + + bool hasPhi = segment4D->hasPhi(); + bool hasZed = segment4D->hasZed(); + + seg4D_hasPhi.push_back(hasPhi); + seg4D_hasZed.push_back(hasZed); + + auto pos = segment4D->localPosition(); + auto dir = segment4D->localDirection(); + + seg4D_posLoc_x.push_back(pos.x()); + seg4D_posLoc_y.push_back(pos.y()); + seg4D_posLoc_z.push_back(pos.z()); + + seg4D_dirLoc_x.push_back(dir.x()); + seg4D_dirLoc_y.push_back(dir.y()); + seg4D_dirLoc_z.push_back(dir.z()); + + double xPosLocSL[2] = {DEFAULT_DOUBLE_VAL, DEFAULT_DOUBLE_VAL}; + bool hasPptSL[2] = {false, false}; + auto xPosLocMidPlane = DEFAULT_DOUBLE_VAL; + + const auto* chamb = m_dtGeometry->chamber(*chambIt); + + StraightLinePlaneCrossing segmentPlaneCrossing{ + chamb->toGlobal(pos).basicVector(), chamb->toGlobal(dir).basicVector(), anyDirection}; + + if (hasPhi) { + for (int iSL = 0; iSL < 2; ++iSL) { + const auto* sl = chamb->superLayer(1 + iSL * 2); + + auto pptSL = segmentPlaneCrossing.position(sl->surface()); + hasPptSL[iSL] = pptSL.first; + + if (hasPptSL[iSL]) { + GlobalPoint segExrapolationToSL(pptSL.second); + xPosLocSL[iSL] = chamb->toLocal(segExrapolationToSL).x(); + } + } + } + + seg4D_posLoc_x_SL1.push_back(xPosLocSL[0]); + seg4D_posLoc_x_SL3.push_back(xPosLocSL[1]); + + if (hasPptSL[0] && hasPptSL[1]) + xPosLocMidPlane = (xPosLocSL[0] + xPosLocSL[1]) * 0.5; + + seg4D_posLoc_x_midPlane.push_back(xPosLocMidPlane); + + const auto begin = seg4D_hitsExpPos.size(); + + auto size{station == 4 ? 8 : 12}; + + nExtr += size; + seg4D_extr_begin.push_back(begin); + seg4D_extr_end.push_back(begin + size); + + std::vector iSLs{1, 2, 3}; + + // CB find a better way + if (station == 4) + iSLs.erase(++iSLs.begin()); + + for (int iL = 1; iL < 5; ++iL) { + for (const auto iSL : iSLs) { + auto* layer = m_dtGeometry->layer(DTWireId{wheel, station, sector, iSL, iL, 2}); + auto ppt = segmentPlaneCrossing.position(layer->surface()); + + bool success{ppt.first}; // check for failure + + auto expPos{DEFAULT_DOUBLE_VAL}; + auto expPosCh{DEFAULT_DOUBLE_VAL}; + auto expWire{DEFAULT_INT_VAL_POS}; + + if (success) { + GlobalPoint segExrapolationToLayer(ppt.second); + LocalPoint segPosAtZWireLayer = layer->toLocal(segExrapolationToLayer); + LocalPoint segPosAtZWireChamber = chamb->toLocal(segExrapolationToLayer); + + if (hasPhi && iSL != 2) { + expPos = segPosAtZWireLayer.x(); + expPosCh = segPosAtZWireChamber.x(); + expWire = layer->specificTopology().channel(segPosAtZWireLayer); + } else if (hasZed && iSL == 2) { + expPos = segPosAtZWireLayer.x(); + expPosCh = segPosAtZWireChamber.y(); + expWire = layer->specificTopology().channel(segPosAtZWireLayer); + } + } + + seg4D_hitsExpPos.push_back(expPos); + seg4D_hitsExpPosCh.push_back(expPosCh); + seg4D_hitsExpWire.push_back(expWire); + } + } + + const GeomDet* geomDet = m_trackingGeometry->idToDet(segment4D->geographicalId()); + auto posGlb = geomDet->toGlobal(pos); + auto dirGlb = geomDet->toGlobal(dir); + + seg4D_posGlb_phi.push_back(posGlb.phi()); + seg4D_posGlb_eta.push_back(posGlb.eta()); + + seg4D_dirGlb_phi.push_back(dirGlb.phi()); + seg4D_dirGlb_eta.push_back(dirGlb.eta()); + + seg2D_hits_begin.push_back(seg2D_hits_pos.size()); + + if (hasPhi) { + const auto* phiSeg = segment4D->phiSegment(); + + seg2D_phi_t0.push_back(phiSeg->t0()); + seg2D_phi_vDrift.push_back(phiSeg->vDrift()); + seg2D_phi_normChi2.push_back(phiSeg->chi2() / phiSeg->degreesOfFreedom()); + seg2D_phi_nHits.push_back(phiSeg->specificRecHits().size()); + + fillHits(phiSeg, geomDet); + } else { + seg2D_phi_t0.push_back(DEFAULT_DOUBLE_VAL); + seg2D_phi_vDrift.push_back(DEFAULT_DOUBLE_VAL); + seg2D_phi_normChi2.push_back(DEFAULT_DOUBLE_VAL_POS); + seg2D_phi_nHits.push_back(0); + } + + if (hasZed) { + const auto* zSeg = segment4D->zSegment(); + + seg2D_z_normChi2.push_back(zSeg->chi2() / zSeg->degreesOfFreedom()); + seg2D_z_nHits.push_back(zSeg->specificRecHits().size()); + + fillHits(zSeg, geomDet); + } else { + seg2D_z_normChi2.push_back(DEFAULT_DOUBLE_VAL_POS); + seg2D_z_nHits.push_back(0); + } + + seg2D_hits_end.push_back(seg2D_hits_pos.size()); + ++nSegments; + } + } + } + + auto table = std::make_unique(nSegments, m_label, false, false); + + table->setDoc("DT segment information"); + + addColumn(table, "wheel", seg4D_wheel, "wheel - [-2:2] range"); + addColumn(table, + "sector", + seg4D_sector, + "sector - [1:14] range" + "
sector 13 used for the second MB4 of sector 4" + "
sector 14 used for the second MB4 of sector 10"); + addColumn(table, "station", seg4D_station, "station - [1:4] range"); + + addColumn(table, "hasPhi", seg4D_hasPhi, "has segment phi view - 0/1 = no/yes"); + addColumn(table, "hasZed", seg4D_hasZed, "has segment zed view - 0/1 = no/yes"); + + addColumn(table, "seg4D_posLoc_x", seg4D_posLoc_x, "position x in local coordinates - cm"); + addColumn(table, "seg4D_posLoc_y", seg4D_posLoc_y, "position y in local coordinates - cm"); + addColumn(table, "seg4D_posLoc_z", seg4D_posLoc_z, "position z in local coordinates - cm"); + addColumn(table, "seg4D_dirLoc_x", seg4D_dirLoc_x, "direction x in local coordinates"); + addColumn(table, "seg4D_dirLoc_y", seg4D_dirLoc_y, "direction y in local coordinates"); + addColumn(table, "seg4D_dirLoc_z", seg4D_dirLoc_z, "direction z in local coordinates"); + + addColumn(table, "seg4D_posLoc_x_SL1", seg4D_posLoc_x_SL1, "position x at SL1 in local coordinates - cm"); + addColumn(table, "seg4D_posLoc_x_SL3", seg4D_posLoc_x_SL3, "position x at SL3 in local coordinates - cm"); + addColumn(table, + "seg4D_posLoc_x_midPlane", + seg4D_posLoc_x_midPlane, + "position x at SL1 - SL3 mid plane in local coordinates - cm"); + + addColumn(table, "seg4D_posGlb_phi", seg4D_posGlb_phi, "position phi in global coordinates - radians [-pi:pi]"); + addColumn(table, "seg4D_posGlb_eta", seg4D_posGlb_eta, "position eta in global coordinates"); + addColumn(table, "seg4D_dirGlb_phi", seg4D_dirGlb_phi, "direction phi in global coordinates - radians [-pi:pi]"); + addColumn(table, "seg4D_dirGlb_eta", seg4D_dirGlb_eta, "direction eta in global coordinates"); + + addColumn(table, "seg2D_phi_t0", seg2D_phi_t0, "t0 from segments with phi view - ns"); + + addColumn(table, "seg2D_phi_vDrift", seg2D_phi_vDrift, "v_drift from segments with phi view"); + addColumn(table, "seg2D_phi_normChi2", seg2D_phi_normChi2, "chi2/n.d.o.f. from segments with phi view"); + addColumn(table, "seg2D_phi_nHits", seg2D_phi_nHits, "# hits in phi view - [0:8] range"); + + addColumn(table, "seg2D_z_normChi2", seg2D_z_normChi2, "chi2/n.d.o.f. from segments with z view"); + addColumn(table, "seg2D_z_nHits", seg2D_z_nHits, "# hits in z view - [0:4] range"); + + addColumn(table, "seg2D_hits_begin", seg2D_hits_begin, "begin of range of quantities in the hits_* vectors"); + addColumn(table, "seg2D_hits_end", seg2D_hits_end, "end of range of quantities in the hits_* vectors"); + + addColumn(table, "seg4D_extr_begin", seg4D_extr_begin, "begin of range of quantities in the extr_* vectors"); + addColumn(table, "seg4D_extr_end", seg4D_extr_end, "end of range of quantities in the extr_* vectors"); + + ev.put(std::move(table)); + + if (m_fillHits) { + auto tabHits = std::make_unique(nHits, m_label + "_hits", false, false); + + tabHits->setDoc("Size of DT segment hits_* vectors"); + + addColumn(tabHits, "pos", seg2D_hits_pos, "local position of a hit in layer local coordinates - x coordinate"); + addColumn( + tabHits, "posCh", seg2D_hits_posCh, "local position of a hit in chamber local coordinates - x coordinate"); + addColumn(tabHits, + "posErr", + seg2D_hits_posErr, + "local position error of a hit in layer local coordinates - xx component of error matrix"); + addColumn(tabHits, "side", seg2D_hits_side, "is hit on L/R side of a cell wire - 1/2 is R/L"); + addColumn(tabHits, "wire", seg2D_hits_wire, "hit wire number - range [1:X] X"); + addColumn(tabHits, "wirePos", seg2D_hits_wirePos, "hit wire position in layer local coordinates - x coordinate"); + addColumn(tabHits, "layer", seg2D_hits_layer, "hit layer number - range [1:4]"); + addColumn(tabHits, "superLayer", seg2D_hits_superLayer, "hit SL number - [1 or 3] the two phi SLs"); + addColumn(tabHits, "time", seg2D_hits_time, "digi time - ns, pedestal not subtracted"); + addColumn(tabHits, "timeCali", seg2D_hits_timeCali, "digi time - ns, pedestal subtracted"); + + ev.put(std::move(tabHits), "hits"); + } + + if (m_fillExtr) { + auto tabExtr = std::make_unique(nExtr, m_label + "_extr", false, false); + + tabExtr->setDoc("Size of DT segment extr_* vectors"); + addColumn(tabExtr, + "ExpPos", + seg4D_hitsExpPos, + "expected position of segment extrapolated" + "
to a given layer in layer local coordinates - cm"); + + addColumn(tabExtr, + "ExpPosCh", + seg4D_hitsExpPosCh, + "expected position of segment extrapolated" + "
to a given layer in chhamber local coordinates - cm"); + + addColumn(tabExtr, + "ExpWire", + seg4D_hitsExpWire, + "expected wire crossed by segment extrapolated" + "
to a given layer - range depends on chamber size"); + + ev.put(std::move(tabExtr), "extr"); + } +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleDTSegmentFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTSegmentFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleDTSegmentFiller.h new file mode 100644 index 0000000000000..8e14b2a078842 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTSegmentFiller.h @@ -0,0 +1,65 @@ +#ifndef MuNtuple_MuNtupleDTSegmentFiller_h +#define MuNtuple_MuNtupleDTSegmentFiller_h + +/** \class MuNtupleDTSegmentFiller MuNtupleDTSegmentFiller.h DPGAnalysis/MuonTools/src/MuNtupleDTSegmentFiller.h + * + * Helper class : the segment filler for Phase-1 / Phase2 DT segments (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "CalibMuon/DTDigiSync/interface/DTTTrigBaseSync.h" +#include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h" + +#include "Geometry/DTGeometry/interface/DTGeometry.h" +#include "Geometry/Records/interface/MuonGeometryRecord.h" + +#include "Geometry/CommonDetUnit/interface/GlobalTrackingGeometry.h" +#include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h" + +class MuNtupleDTSegmentFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleDTSegmentFiller(const edm::ParameterSet &); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); + +protected: + /// Fill tree branches for a given event + void fill(edm::Event &ev) final; + + /// Get info from the ES by run + void getFromES(const edm::Run &run, const edm::EventSetup &environment) final; + + /// Get info from the ES for a given event + void getFromES(const edm::EventSetup &environment) final; + +private: + /// The segment token + nano_mu::EDTokenHandle m_token; + + /// Fill rec-hit table + bool m_fillHits; + + /// Fill segment extrapolation table + bool m_fillExtr; + + /// DT Geometry + nano_mu::ESTokenHandle m_dtGeometry; + + /// Tracking Geometry + nano_mu::ESTokenHandle m_trackingGeometry; + + /// Handle DT trigger time pedestals + std::unique_ptr m_dtSync; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.cc new file mode 100644 index 0000000000000..616f267e4e92f --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.cc @@ -0,0 +1,168 @@ +/** \class MuNtupleDTTPGPhiFiller MuNtupleDTTPGPhiFiller.cc DPGAnalysis/MuonTools/src/MuNtupleDTTPGPhiFiller.cc + * + * Helper class : the Phase-1 local trigger filler for TwinMux in/out and BMTF in (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.h" +#include "FWCore/ParameterSet/interface/allowedValues.h" + +#include +#include + +MuNtupleDTTPGPhiFiller::MuNtupleDTTPGPhiFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller{config}, + m_tag{getTag(config)}, + m_token{config, consumesCollector(), "dtTpTag"}, + m_trigGeomUtils{consumesCollector()} { + produces(); +} + +void MuNtupleDTTPGPhiFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "ltBmtfIn"); + desc.ifValue(edm::ParameterDescription("tag", "BMTF_IN", true), + edm::allowedValues("BMTF_IN", "TM_IN", "TM_OUT")); + desc.add("dtTpTag", edm::InputTag{"bmtfDigis"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleDTTPGPhiFiller::getFromES(const edm::Run& run, const edm::EventSetup& environment) { + m_trigGeomUtils.getFromES(run, environment); +} + +void MuNtupleDTTPGPhiFiller::fill(edm::Event& ev) { + unsigned int nTrigs{0}; + + std::vector wheel; + std::vector sector; + std::vector station; + + std::vector quality; + std::vector rpcBit; + + std::vector phi; + std::vector phiB; + + std::vector posLoc_x; + std::vector dirLoc_phi; + + std::vector bx; + std::vector is2nd; + + auto trigColl = m_token.conditionalGet(ev); + + if (trigColl.isValid()) { + const auto trigs = trigColl->getContainer(); + for (const auto& trig : (*trigs)) { + if (trig.code() != 7) { + wheel.push_back(trig.whNum()); + sector.push_back(trig.scNum() + (m_tag != TriggerTag::BMTF_IN ? 1 : 0)); + station.push_back(trig.stNum()); + + quality.push_back(trig.code()); + + if (m_tag == TriggerTag::TM_OUT) + rpcBit.push_back(trig.RpcBit()); + + phi.push_back(trig.phi()); + phiB.push_back(trig.phiB()); + + auto [x, dir] = m_trigGeomUtils.trigToReco(&trig); + + posLoc_x.push_back(x); + dirLoc_phi.push_back(dir); + + bx.push_back(trig.bxNum() - (m_tag == TriggerTag::TM_IN && trig.Ts2Tag() ? 1 : 0)); + is2nd.push_back(trig.Ts2Tag()); + + ++nTrigs; + } + } + } + + auto table = std::make_unique(nTrigs, m_label, false, false); + + table->setDoc("DT trigger primitive information (phi view)"); + + addColumn(table, "wheel", wheel, "wheel - [-2:2] range"); + addColumn(table, + "sector", + sector, + "sector" + "
- [1:12] range for TwinMux" + "
- [0:11] range for BMTF input" + "
- double MB4 stations are part of S4 and S10 in TwinMux" + "
- double MB4 stations are part of S3 and S9 in BMTF input"); + addColumn(table, "station", station, "station - [1:4] range"); + addColumn(table, + "quality", + quality, + "quality - [0:6] range" + "
- [0:1] : uncorrelated L triggers" + "
- [2:3] : uncorrelated H triggers" + "
- 4 : correlated LL triggers" + "
- 5 : correlated HL triggers" + "
- 6 : correlated HH triggers"); + if (m_tag == TriggerTag::TM_OUT) { + addColumn(table, + "rpcBit", + rpcBit, + "use of RPC - [0:2] range" + "
- 0 : RPC not used" + "
- 1 : RPC+DT combined trigger" + "
- 2 : RPC-only trigger"); + } + + addColumn(table, + "phi", + phi, + "phi - scale and range:" + "
- 4096 correstpond to 1 rad" + "
- 0 is @ (DT sector - 1) * 30 deg in global CMS phi"); + addColumn(table, + "phiB", + phiB, + "phiB - scale and range:" + "
- 512 correstpond to 1 rad" + "
- 0 is a muon with infinite pT (straight line)"); + addColumn(table, "posLoc_x", posLoc_x, "position x in chamber local coordinates - cm"); + addColumn(table, "dirLoc_phi", dirLoc_phi, "direction phi angle in chamber local coordinates - deg"); + addColumn(table, + "bx", + bx, + "bx:" + "
- BX = 0 is the one where the event is collected" + "
- TwinMux range [X:Y]" + "
- BMT input range [X:Y]"); + addColumn(table, "is2nd", is2nd, "1st/2nd track flag - [0:1]"); + + ev.put(std::move(table)); +} + +MuNtupleDTTPGPhiFiller::TriggerTag MuNtupleDTTPGPhiFiller::getTag(const edm::ParameterSet& config) { + auto tag{TriggerTag::TM_IN}; + + auto tagName = config.getParameter("tag"); + + if (tagName != "TM_IN" && tagName != "TM_OUT" && tagName != "BMTF_IN") + edm::LogError("") << "[MuNtupleDTTPGPhiFiller]::getTag: " << tagName << " is not a valid tag, defaulting to TM_IN"; + + if (tagName == "TM_OUT") { + tag = TriggerTag::TM_OUT; + } else if (tagName == "BMTF_IN") { + tag = TriggerTag::BMTF_IN; + } + + return tag; +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleDTTPGPhiFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.h new file mode 100644 index 0000000000000..f9dd020e94548 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.h @@ -0,0 +1,52 @@ +#ifndef MuNtuple_MuNtupleDTTPGPhiFiller_h +#define MuNtuple_MuNtupleDTTPGPhiFiller_h + +/** \class MuNtupleDTTPGPhiFiller MuNtupleDTTPGPhiFiller.h DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGPhiFiller.h + * + * Helper class : the Phase-1 local trigger filler for TwinMux in/out and BMTF in (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambPhContainer.h" + +class MuNtupleDTTPGPhiFiller : public MuNtupleBaseFiller { +public: + enum class TriggerTag { TM_IN = 0, TM_OUT, BMTF_IN }; + + /// Constructor + MuNtupleDTTPGPhiFiller(const edm::ParameterSet&); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + /// Fill tree branches for a given event + void fill(edm::Event& ev) final; + + /// Get info from the ES by run + void getFromES(const edm::Run& run, const edm::EventSetup& environment) final; + +private: + /// Enum to activate "flavour-by-flavour" + /// changes in the filling logic + TriggerTag m_tag; + + /// The trigger-primitive token + nano_mu::EDTokenHandle m_token; + + /// The class to perform DT local trigger coordinate conversions + nano_mu::DTTrigGeomUtils m_trigGeomUtils; + + /// Helper function translating config parameter into TriggerTag + TriggerTag getTag(const edm::ParameterSet& config); +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.cc new file mode 100644 index 0000000000000..1af9dc3ab3869 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.cc @@ -0,0 +1,124 @@ +/** \class MuNtupleDTTPGThetaFiller MuNtupleDTTPGThetaFiller.cc DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.cc + * + * Helper class : the Phase-1 local trigger filler for TwinMux in/out and BMTF in (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.h" +#include "FWCore/ParameterSet/interface/allowedValues.h" + +#include +#include + +MuNtupleDTTPGThetaFiller::MuNtupleDTTPGThetaFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller{config}, m_tag{getTag(config)}, m_token{config, consumesCollector(), "dtTpTag"} { + produces(); +} + +void MuNtupleDTTPGThetaFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "ltBmtfInTh"); + desc.ifValue(edm::ParameterDescription("tag", "BMTF_IN", true), + edm::allowedValues("BMTF_IN", "TM_IN")); + desc.add("dtTpTag", edm::InputTag{"bmtfDigis"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleDTTPGThetaFiller::fill(edm::Event& ev) { + unsigned int nTrigs{0}; + + std::vector wheel; + std::vector sector; + std::vector station; + + std::vector bx; + std::vector hitMap; + + auto trigColl = m_token.conditionalGet(ev); + + if (trigColl.isValid()) { + const auto trigs = trigColl->getContainer(); + for (const auto& trig : (*trigs)) { + bool hasData = false; + for (int pos = 0; pos < 7; ++pos) { + if (trig.code(pos)) { + hasData = true; + break; + } + } + + if (!hasData) + continue; + + wheel.push_back(trig.whNum()); + sector.push_back(trig.scNum() + (m_tag != TriggerTag::BMTF_IN ? 1 : 0)); + station.push_back(trig.stNum()); + + bx.push_back(trig.bxNum()); + + uint32_t hitMapCh = 0; + for (int pos = 0; pos < 7; ++pos) + if (trig.code(pos)) + hitMapCh = hitMapCh | (0x1 << pos); + + hitMap.push_back(hitMapCh); + + ++nTrigs; + } + } + + auto table = std::make_unique(nTrigs, m_label, false, false); + + table->setDoc("DT trigger primitive information (theta view)"); + + addColumn(table, "wheel", wheel, "wheel - [-2:2] range"); + addColumn(table, + "sector", + sector, + "sector" + "
- [1:12] range for TwinMux" + "
- [0:11] range for BMTF input" + "
- double MB4 stations are part of S4 and S10 in TwinMux" + "
- double MB4 stations are part of S3 and S9 in BMTF input"); + addColumn(table, "station", station, "station - [1:3] range"); + addColumn(table, + "bx", + bx, + "bx:" + "
- BX = 0 is the one where the event is collected" + "
- TwinMux range [X:Y]" + "
- BMTF input range [X:Y]"); + addColumn(table, + "hitMap", + hitMap, + "Map groups of BTIs that fired (unsigned int):" + "
there are 7 groups of BTI per chamber, the first one" + "
being the less significant bit of the map [CHECK]"); + + ev.put(std::move(table)); +} + +MuNtupleDTTPGThetaFiller::TriggerTag MuNtupleDTTPGThetaFiller::getTag(const edm::ParameterSet& config) { + auto tag{TriggerTag::TM_IN}; + + auto tagName = config.getParameter("tag"); + + if (tagName != "TM_IN" && tagName != "BMTF_IN") + edm::LogError("") << "[MuNtupleDTTPGThetaFiller]::getTag: " << tagName + << " is not a valid tag, defaulting to TM_IN"; + + if (tagName == "BMTF_IN") + tag = TriggerTag::BMTF_IN; + + return tag; +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleDTTPGThetaFiller); \ No newline at end of file diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.h new file mode 100644 index 0000000000000..1df81ceec51ef --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.h @@ -0,0 +1,46 @@ +#ifndef MuNtuple_MuNtupleDTTPGThetaFiller_h +#define MuNtuple_MuNtupleDTTPGThetaFiller_h + +/** \class MuNtupleDTTPGThetaFiller MuNtupleDTTPGThetaFiller.h DPGAnalysis/MuonTools/plugins/MuNtupleDTTPGThetaFiller.h + * + * Helper class : the Phase-1 local trigger filler for TwinMux and BMTF in (the DataFormat is the same) + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambThContainer.h" + +class MuNtupleDTTPGThetaFiller : public MuNtupleBaseFiller { +public: + enum class TriggerTag { TM_IN = 0, BMTF_IN }; + + /// Constructor + MuNtupleDTTPGThetaFiller(const edm::ParameterSet&); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + /// Fill tree branches for a given events + void fill(edm::Event& ev) final; + +private: + /// Enum to activate "flavour-by-flavour" + /// changes in the filling logic + TriggerTag m_tag; + + /// The trigger-primitive token + nano_mu::EDTokenHandle m_token; + + /// Helper function translating config parameter into TriggerTag + TriggerTag getTag(const edm::ParameterSet& config); +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMDigiFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMDigiFiller.cc new file mode 100644 index 0000000000000..6a338198e39b2 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMDigiFiller.cc @@ -0,0 +1,129 @@ +#include "DPGAnalysis/MuonTools/plugins/MuNtupleGEMDigiFiller.h" + +#include + +MuNtupleGEMDigiFiller::MuNtupleGEMDigiFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller(config), m_token{config, consumesCollector(), "gemDigiTag"}, m_gemGeometry{consumesCollector()} { + produces(); +} + +void MuNtupleGEMDigiFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "gemDigis"); + desc.add("gemDigiTag", edm::InputTag{"muonGEMDigis"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleGEMDigiFiller::getFromES(const edm::Run& run, const edm::EventSetup& environment) { + m_gemGeometry.getFromES(environment); +} + +void MuNtupleGEMDigiFiller::fill(edm::Event& ev) { + unsigned int nDigis{0}; + + std::vector station; + std::vector roll; + std::vector strip; + std::vector bx; + std::vector region; + + std::vector g_r; + std::vector g_phi; + std::vector g_eta; + std::vector g_x; + std::vector g_y; + std::vector g_z; + + auto gemDigis = m_token.conditionalGet(ev); + + if (gemDigis.isValid()) { + auto gemDetIdIt = gemDigis->begin(); + auto gemDetIdEnd = gemDigis->end(); + + for (; gemDetIdIt != gemDetIdEnd; ++gemDetIdIt) { + const auto& [gemDetId, range] = (*gemDetIdIt); + + const auto etaPart = m_gemGeometry->etaPartition(gemDetId); + const auto& surface = etaPart->surface(); + + for (auto digi = range.first; digi != range.second; ++digi) { + station.push_back(gemDetId.station()); + roll.push_back(gemDetId.roll()); + strip.push_back(digi->strip()); + bx.push_back(digi->bx()); + + region.push_back(gemDetId.region()); + + const auto& localPos = etaPart->centreOfStrip(digi->strip()); + const auto& globalPos = surface.toGlobal(localPos); + + g_r.push_back(globalPos.perp()); + g_phi.push_back(globalPos.phi()); + g_eta.push_back(globalPos.eta()); + g_x.push_back(globalPos.x()); + g_y.push_back(globalPos.y()); + g_z.push_back(globalPos.z()); + + ++nDigis; + } + } + } + + auto table = std::make_unique(nDigis, m_label, false, false); + + table->setDoc("GEM digi information"); + + addColumn(table, "station", station, "GEM station
always 1 for GE1/1)"); + addColumn(table, + "roll", + roll, + "roll id (also known as eta partition)" + "
(partitions numbered from 1 to 8"); + addColumn(table, "strip", strip, "index of the readout strip associated to the digi"); + addColumn(table, + "region", + region, + "GE11 region where the digi is detected" + "
(int, positive endcap: +1, negative endcap: -1"); + addColumn(table, "bx", bx, "bunch crossing associated to the digi"); + + addColumn(table, + "g_r", + g_r, + "Global digi position" + "
(global radial coordinates, cm)"); + addColumn(table, + "g_phi", + g_phi, + "Global digi position" + "
(global azimuthal coordinates, rad)"); + addColumn(table, + "g_eta", + g_eta, + "Global digi position" + "
(global pseudorapidity coordinates)"); + addColumn(table, + "g_x", + g_x, + "Global digi position" + "
(global x coordinates, cm)"); + addColumn(table, + "g_y", + g_y, + "Global digi position" + "
(global y coordinates, cm)"); + addColumn(table, + "g_z", + g_z, + "Global digi position" + "
(global z coordinates, cm)"); + + ev.put(std::move(table)); +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleGEMDigiFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMDigiFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMDigiFiller.h new file mode 100644 index 0000000000000..29e8bea7d6b62 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMDigiFiller.h @@ -0,0 +1,36 @@ +#ifndef MuNtuple_MuNtupleGEMDigiFiller_h +#define MuNtuple_MuNtupleGEMDigiFiller_h + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/GEMDigi/interface/GEMDigiCollection.h" + +#include "Geometry/GEMGeometry/interface/GEMGeometry.h" +#include "Geometry/Records/interface/MuonGeometryRecord.h" + +class MuNtupleGEMDigiFiller : public MuNtupleBaseFiller { +public: + //Constructor + MuNtupleGEMDigiFiller(const edm::ParameterSet &); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); + +protected: + void fill(edm::Event &ev) final; + + /// Get info from the ES by run + void getFromES(const edm::Run &run, const edm::EventSetup &environment) final; + +private: + /// The digi token + nano_mu::EDTokenHandle m_token; + + /// GEM Geometry + nano_mu::ESTokenHandle m_gemGeometry; +}; + +#endif \ No newline at end of file diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMMuonFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMMuonFiller.cc new file mode 100644 index 0000000000000..45badbf846b9b --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMMuonFiller.cc @@ -0,0 +1,600 @@ +#include "DPGAnalysis/MuonTools/plugins/MuNtupleGEMMuonFiller.h" + +#include "DataFormats/MuonDetId/interface/MuonSubdetId.h" +#include "DataFormats/MuonDetId/interface/GEMDetId.h" +#include "DataFormats/MuonDetId/interface/CSCDetId.h" + +#include "DataFormats/GeometryCommonDetAlgo/interface/ErrorFrameTransformer.h" +#include "TrackPropagation/SteppingHelixPropagator/interface/SteppingHelixStateInfo.h" + +#include + +MuNtupleGEMMuonFiller::MuNtupleGEMMuonFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller{config}, + m_token{config, consumesCollector(), "muonToken"}, + m_fillPropagated{config.getParameter("fillPropagated")}, + m_gemGeometry{consumesCollector()}, + m_transientTrackBuilder{consumesCollector(), "TransientTrackBuilder"}, + m_muonSP{std::make_unique(config.getParameter("ServiceParameters"), + consumesCollector())} { + produces(); + + if (m_fillPropagated) { + produces("propagated"); + } +} + +void MuNtupleGEMMuonFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "gemMuon"); + desc.add("muonToken", edm::InputTag{"muons"}); + + desc.add("fillPropagated", true); + desc.setAllowAnything(); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleGEMMuonFiller::getFromES(const edm::Run& run, const edm::EventSetup& environment) { + m_gemGeometry.getFromES(environment); +} + +void MuNtupleGEMMuonFiller::getFromES(const edm::EventSetup& environment) { + m_transientTrackBuilder.getFromES(environment); + m_muonSP->update(environment); +} + +void MuNtupleGEMMuonFiller::fill(edm::Event& ev) { + unsigned int nMuons{0}; + + std::vector pt; + std::vector phi; + std::vector eta; + std::vector charge; + + std::vector isGlobal; + std::vector isStandalone; + std::vector isTracker; + std::vector isGEM; + std::vector isCSC; + std::vector isME11; + + std::vector isLoose; + std::vector isMedium; + std::vector isTight; + + std::vector numberOfValidPixelHits; + std::vector innerTracker_ValidFraction; + std::vector numberOfValidTrackerHits; + + std::vector trackNormChi2; + + std::vector innermost_x; + std::vector innermost_y; + std::vector innermost_z; + + std::vector outermost_x; + std::vector outermost_y; + std::vector outermost_z; + + unsigned int nProp{0}; + + std::vector propagated_muIdx; + + std::vector propagated_isincoming; + std::vector propagated_isinsideout; + std::vector propagated_region; + std::vector propagated_layer; + std::vector propagated_chamber; + std::vector propagated_etaP; + + std::vector propagatedLoc_x; + std::vector propagatedLoc_y; + std::vector propagatedLoc_z; + std::vector propagatedLoc_r; + std::vector propagatedLoc_phi; + std::vector propagatedLoc_dirX; + std::vector propagatedLoc_dirY; + std::vector propagatedLoc_dirZ; + std::vector propagatedLoc_errX; + std::vector propagatedLoc_errY; + + std::vector propagatedGlb_x; + std::vector propagatedGlb_y; + std::vector propagatedGlb_z; + std::vector propagatedGlb_r; + std::vector propagatedGlb_phi; + std::vector propagatedGlb_errX; + std::vector propagatedGlb_errY; + std::vector propagatedGlb_phierr; + std::vector propagatedGlb_rerr; + + std::vector propagated_EtaPartition_centerX; + std::vector propagated_EtaPartition_centerY; + std::vector propagated_EtaPartition_phiMax; + std::vector propagated_EtaPartition_phiMin; + std::vector propagated_EtaPartition_rMax; + std::vector propagated_EtaPartition_rMin; + + std::vector propagated_nME1hits; + std::vector propagated_nME2hits; + std::vector propagated_nME3hits; + std::vector propagated_nME4hits; + + auto muons = m_token.conditionalGet(ev); + + // edm::ESHandle + auto&& propagator_any = m_muonSP->propagator("SteppingHelixPropagatorAny"); + auto&& propagator_along = m_muonSP->propagator("SteppingHelixPropagatorAlong"); + auto&& propagator_opposite = m_muonSP->propagator("SteppingHelixPropagatorOpposite"); + + if (!propagator_any.isValid() || !propagator_along.isValid() || !propagator_opposite.isValid()) { + return; + } + + if (muons.isValid() && m_transientTrackBuilder.isValid()) { + //loop on recoMuons + for (const auto& muon : (*muons)) { + pt.push_back(muon.pt()); + eta.push_back(muon.eta()); + phi.push_back(muon.phi()); + charge.push_back(muon.charge()); + + isGlobal.push_back(muon.isGlobalMuon()); + isStandalone.push_back(muon.isStandAloneMuon()); + isTracker.push_back(muon.isTrackerMuon()); + isGEM.push_back(muon.isGEMMuon()); + + isLoose.push_back(muon.passed(reco::Muon::CutBasedIdLoose)); + isMedium.push_back(muon.passed(reco::Muon::CutBasedIdMedium)); + isTight.push_back(muon.passed(reco::Muon::CutBasedIdTight)); + + if (!muon.innerTrack().isNull()) { + numberOfValidPixelHits.push_back(muon.innerTrack()->hitPattern().numberOfValidPixelHits()); + innerTracker_ValidFraction.push_back(muon.innerTrack()->validFraction()); + numberOfValidTrackerHits.push_back(muon.innerTrack()->hitPattern().numberOfValidTrackerHits()); + } else { + numberOfValidPixelHits.push_back(DEFAULT_DOUBLE_VAL); + innerTracker_ValidFraction.push_back(DEFAULT_DOUBLE_VAL); + numberOfValidTrackerHits.push_back(DEFAULT_DOUBLE_VAL); + } + ++nMuons; + + bool is_csc = false; + bool is_me11 = false; + + if (!muon.outerTrack().isNull()) { + const auto track = muon.outerTrack().get(); + const auto outerTrackRef = muon.outerTrack(); + + float p2_in = track->innerMomentum().mag2(); + float p2_out = track->outerMomentum().mag2(); + float pos_out = track->outerPosition().mag2(); + float pos_in = track->innerPosition().mag2(); + + bool is_insideout = pos_in > pos_out; + + if (is_insideout) { + std::swap(pos_in, pos_out); + std::swap(p2_in, p2_out); + } + + bool is_incoming = p2_out > p2_in; + + const auto&& transient_track = m_transientTrackBuilder->build(track); + const auto& htp = transient_track.hitPattern(); + + if (transient_track.isValid()) { + trackNormChi2.push_back(transient_track.normalizedChi2()); + innermost_x.push_back(transient_track.innermostMeasurementState().globalPosition().x()); + innermost_y.push_back(transient_track.innermostMeasurementState().globalPosition().y()); + innermost_z.push_back(transient_track.innermostMeasurementState().globalPosition().z()); + outermost_x.push_back(transient_track.outermostMeasurementState().globalPosition().x()); + outermost_y.push_back(transient_track.outermostMeasurementState().globalPosition().y()); + outermost_z.push_back(transient_track.outermostMeasurementState().globalPosition().z()); + + } else { + trackNormChi2.push_back(DEFAULT_DOUBLE_VAL); + innermost_x.push_back(DEFAULT_DOUBLE_VAL); + innermost_y.push_back(DEFAULT_DOUBLE_VAL); + innermost_z.push_back(DEFAULT_DOUBLE_VAL); + outermost_x.push_back(DEFAULT_DOUBLE_VAL); + outermost_y.push_back(DEFAULT_DOUBLE_VAL); + outermost_z.push_back(DEFAULT_DOUBLE_VAL); + continue; + } + + const auto&& start_state = + is_insideout ? transient_track.outermostMeasurementState() : transient_track.innermostMeasurementState(); + auto& propagator = is_incoming ? propagator_along : propagator_opposite; + + auto recHitMu = outerTrackRef->recHitsBegin(); + auto recHitMuEnd = outerTrackRef->recHitsEnd(); + + //loop on recHits which form the outerTrack + for (; recHitMu != recHitMuEnd; ++recHitMu) { + DetId detId{(*recHitMu)->geographicalId()}; + + if (detId.det() == DetId::Muon && detId.subdetId() == MuonSubdetId::CSC) { + is_csc = true; + + const CSCDetId csc_id{detId}; + // ME11 chambers consist of 2 subchambers: ME11a, ME11b. + // In CMSSW they are referred as Stat. 1 Ring 1, Stat. 1 Ring. 4 respectively + if (csc_id.station() == 1 && ((csc_id.ring() == 1) || (csc_id.ring() == 4))) + is_me11 = true; + } + } //loop on recHits + + //if at least one CSC hit is found, perform propagation + if (is_csc) { + // CSC Hits + int8_t nME1_hits = 0; + int8_t nME2_hits = 0; + int8_t nME3_hits = 0; + int8_t nME4_hits = 0; + + int nHits{htp.numberOfAllHits(htp.TRACK_HITS)}; + + for (int i = 0; i != nHits; ++i) { + uint32_t hit = htp.getHitPattern(htp.TRACK_HITS, i); + int substructure = htp.getSubStructure(hit); + int hittype = htp.getHitType(hit); + + if (substructure == 2 && hittype == 0) { // CSC Hits + int CSC_station = htp.getMuonStation(hit); + + switch (CSC_station) { + case 1: + ++nME1_hits; + break; + case 2: + ++nME2_hits; + break; + case 3: + ++nME3_hits; + break; + case 4: + ++nME4_hits; + break; + default: + // do nothing + break; + } + } + } + //loop on GEM etaPartitions + for (const auto& eta_partition : m_gemGeometry->etaPartitions()) { + if (eta_partition->id().station() != 1) { + continue; //Only takes GE1/1 + } + const GEMDetId&& gem_id = eta_partition->id(); + + bool is_opposite_region = muon.eta() * gem_id.region() < 0; + if (is_incoming xor is_opposite_region) { + continue; //Check on muon direction + } + const BoundPlane& bound_plane = eta_partition->surface(); + + const auto& dest_state = propagator->propagate(start_state, bound_plane); + if (!dest_state.isValid()) { + // std::cout << "failed to propagate" << std::endl; + continue; + } + const GlobalPoint&& dest_global_pos = dest_state.globalPosition(); + const LocalPoint&& local_point = eta_partition->toLocal(dest_global_pos); + const LocalPoint local_point_2d{local_point.x(), local_point.y(), 0.0f}; + + if (eta_partition->surface().bounds().inside(local_point_2d)) { + //// PROPAGATED HIT ERROR EVALUATION + // X,Y + double xx = dest_state.curvilinearError().matrix()(3, 3); + double yy = dest_state.curvilinearError().matrix()(4, 4); + double xy = dest_state.curvilinearError().matrix()(4, 3); + double dest_glob_error_x = sqrt(0.5 * (xx + yy - sqrt((xx - yy) * (xx - yy) + 4 * xy * xy))); + double dest_glob_error_y = sqrt(0.5 * (xx + yy + sqrt((xx - yy) * (xx - yy) + 4 * xy * xy))); + + // R,Phi + const LocalPoint&& dest_local_pos = eta_partition->toLocal(dest_global_pos); + const LocalError&& dest_local_err = dest_state.localError().positionError(); + const GlobalError& dest_global_err = + ErrorFrameTransformer().transform(dest_local_err, eta_partition->surface()); + const double dest_global_r_err = std::sqrt(dest_global_err.rerr(dest_global_pos)); + const double dest_global_phi_err = std::sqrt(dest_global_err.phierr(dest_global_pos)); + + ++nProp; + + propagated_muIdx.push_back(nMuons - 1); + + propagated_nME1hits.push_back(nME1_hits); + propagated_nME2hits.push_back(nME2_hits); + propagated_nME3hits.push_back(nME3_hits); + propagated_nME4hits.push_back(nME4_hits); + + propagated_EtaPartition_centerX.push_back(eta_partition->position().x()); + propagated_EtaPartition_centerY.push_back(eta_partition->position().y()); + propagated_EtaPartition_rMin.push_back(eta_partition->surface().rSpan().first); + propagated_EtaPartition_rMax.push_back(eta_partition->surface().rSpan().second); + propagated_EtaPartition_phiMin.push_back(eta_partition->surface().phiSpan().first); + propagated_EtaPartition_phiMax.push_back(eta_partition->surface().phiSpan().second); + + propagatedGlb_x.push_back(dest_global_pos.x()); + propagatedGlb_y.push_back(dest_global_pos.y()); + propagatedGlb_z.push_back(dest_global_pos.z()); + propagatedGlb_r.push_back(dest_global_pos.perp()); + propagatedGlb_phi.push_back(dest_global_pos.phi()); + + propagatedLoc_x.push_back(dest_local_pos.x()); + propagatedLoc_y.push_back(dest_local_pos.y()); + propagatedLoc_z.push_back(dest_local_pos.z()); + propagatedLoc_r.push_back(dest_local_pos.perp()); + propagatedLoc_phi.push_back(dest_local_pos.phi()); + propagatedLoc_dirX.push_back(dest_state.localDirection().x()); + propagatedLoc_dirY.push_back(dest_state.localDirection().y()); + propagatedLoc_dirZ.push_back(dest_state.localDirection().z()); + + propagatedLoc_errX.push_back(dest_local_err.xx()); + propagatedLoc_errY.push_back(dest_local_err.yy()); + + propagatedGlb_errX.push_back(dest_glob_error_x); + propagatedGlb_errY.push_back(dest_glob_error_y); + propagatedGlb_rerr.push_back(dest_global_r_err); + propagatedGlb_phierr.push_back(dest_global_phi_err); + + propagated_region.push_back(gem_id.region()); + propagated_layer.push_back(gem_id.layer()); + propagated_chamber.push_back(gem_id.chamber()); + propagated_etaP.push_back(gem_id.roll()); + + propagated_isinsideout.push_back(is_insideout); + propagated_isincoming.push_back(is_incoming); + + } //propagation is inside boundaries + } //loop on EtaPartitions + } //is_csc therefore perform propagation + } else { //!muon.outerTrack().isNull() + trackNormChi2.push_back(DEFAULT_DOUBLE_VAL); + innermost_x.push_back(DEFAULT_DOUBLE_VAL); + innermost_y.push_back(DEFAULT_DOUBLE_VAL); + innermost_z.push_back(DEFAULT_DOUBLE_VAL); + outermost_x.push_back(DEFAULT_DOUBLE_VAL); + outermost_y.push_back(DEFAULT_DOUBLE_VAL); + outermost_z.push_back(DEFAULT_DOUBLE_VAL); + } + isCSC.push_back(is_csc); + isME11.push_back(is_me11); + + } //loop on reco muons + } + + auto table = std::make_unique(nMuons, m_label, false, false); + + table->setDoc("RECO muon information (GEM specific object)"); + + addColumn(table, "pt", pt, ""); + addColumn(table, "phi", phi, ""); + addColumn(table, "eta", eta, ""); + addColumn(table, "charge", charge, ""); + + addColumn(table, "isGlobal", isGlobal, ""); + addColumn(table, "isStandalone", isStandalone, ""); + addColumn(table, "isTracker", isTracker, ""); + addColumn(table, "isGEM", isGEM, ""); + addColumn(table, "isCSC", isCSC, ""); + addColumn(table, "isME11", isME11, ""); + + addColumn(table, "isLoose", isLoose, ""); + addColumn(table, "isMedium", isMedium, ""); + addColumn(table, "isTight", isTight, ""); + + addColumn(table, "numberOfValidPixelHits", numberOfValidPixelHits, ""); + addColumn(table, "innerTracker_ValidFraction", innerTracker_ValidFraction, ""); + addColumn(table, "numberOfValidTrackerHits", numberOfValidTrackerHits, ""); + + addColumn(table, "trackNormChi2", trackNormChi2, ""); + + addColumn(table, "innermost_x", innermost_x, ""); + addColumn(table, "innermost_y", innermost_y, ""); + addColumn(table, "innermost_z", innermost_z, ""); + + addColumn(table, "outermost_x", outermost_x, ""); + addColumn(table, "outermost_y", outermost_y, ""); + addColumn(table, "outermost_z", outermost_z, ""); + ev.put(std::move(table)); + + if (m_fillPropagated) { + auto tabProp = std::make_unique(nProp, m_label + "_propagated", false, false); + + addColumn(tabProp, "propagated_muIdx", propagated_muIdx, ""); + + addColumn(tabProp, + "propagated_nME1hits", + propagated_nME1hits, + "number of hits in the CSC ME1 station" + "in the STA muon track extrapolated to GE11"); + addColumn(tabProp, + "propagated_nME2hits", + propagated_nME2hits, + "number of hits in the CSC ME2 station" + "in the STA muon track extrapolated to GE11"); + addColumn(tabProp, + "propagated_nME3hits", + propagated_nME3hits, + "number of hits in the CSC ME3 station" + "in the STA muon track extrapolated to GE11"); + addColumn(tabProp, + "propagated_nME4hits", + propagated_nME4hits, + "number of hits in the CSC ME4 station" + "in the STA muon track extrapolated to GE11"); + + addColumn( + tabProp, "propagated_isincoming", propagated_isincoming, "bool, condition on the muon STA track direction"); + addColumn( + tabProp, "propagated_isinsideout", propagated_isinsideout, "bool, condition on the muon STA track direction"); + addColumn(tabProp, + "propagated_region", + propagated_region, + "GE11 region where the extrapolated muon track falls" + "
(int, positive endcap: +1, negative endcap: -1"); + addColumn(tabProp, + "propagated_layer", + propagated_layer, + "GE11 layer where the extrapolated muon track falls" + "
(int, layer1: 1, layer2: 2"); + addColumn(tabProp, + "propagated_chamber", + propagated_chamber, + "GE11 superchamber where the extrapolated muon track falls" + "
(int, chambers numbered from 0 to 35"); + addColumn(tabProp, + "propagated_etaP", + propagated_etaP, + "GE11 eta partition where the extrapolated muon track falls" + "
(int, partitions numbered from 1 to 8"); + + addColumn(tabProp, + "propagatedLoc_x", + propagatedLoc_x, + "expected position of muon track extrapolated to GE11 surface" + "
(float, local layer x coordinates, cm)"); + addColumn(tabProp, + "propagatedLoc_y", + propagatedLoc_y, + "expected position of muon track extrapolated to GE11 surface" + "
(float, local layer y coordinates, cm)"); + addColumn(tabProp, + "propagatedLoc_z", + propagatedLoc_z, + "expected position of muon track extrapolated to GE11 surface" + "
(float, local layer z coordinates, cm)"); + addColumn(tabProp, + "propagatedLoc_r", + propagatedLoc_r, + "expected position of muon track extrapolated to GE11 surface" + "
(float, local layer radial coordinate, cm)"); + addColumn(tabProp, + "propagatedLoc_phi", + propagatedLoc_phi, + "expected position of muon track extrapolated to GE11 surface" + "
(float, local layer phi coordinates, rad)"); + + addColumn(tabProp, + "propagatedLoc_dirX", + propagatedLoc_dirX, + "direction cosine of angle between local x axis and GE11 plane" + "
(float, dir. cosine)"); + addColumn(tabProp, + "propagatedLoc_dirY", + propagatedLoc_dirY, + "direction cosine of angle between local y axis and GE11 plane" + "
(float, dir. cosine)"); + addColumn(tabProp, + "propagatedLoc_dirZ", + propagatedLoc_dirZ, + "direction cosine of angle between local z axis and GE11 plane" + "
(float, dir. cosine)"); + + addColumn(tabProp, + "propagatedLoc_errX", + propagatedLoc_errX, + "uncertainty on expected position of muon track extrapolated to GE11 surface" + "
(float, local layer x coordinates, cm)"); + addColumn(tabProp, + "propagatedLoc_errY", + propagatedLoc_errY, + "uncertainty on expected position of muon track extrapolated to GE11 surface" + "
(float, local layer y coordinates, cm)"); + + addColumn(tabProp, + "propagatedGlb_x", + propagatedGlb_x, + "expected position of muon track extrapolated to GE11 surface" + "
(float, global x coordinates, cm)"); + addColumn(tabProp, + "propagatedGlb_y", + propagatedGlb_y, + "expected position of muon track extrapolated to GE11 surface" + "
(float, global y coordinates, cm)"); + addColumn(tabProp, + "propagatedGlb_z", + propagatedGlb_z, + "expected position of muon track extrapolated to GE11 surface" + "
(float, global z coordinates, cm)"); + addColumn(tabProp, + "propagatedGlb_r", + propagatedGlb_r, + "expected position of muon track extrapolated to GE11 surface" + "
(float, global radial (r) coordinates, cm)"); + addColumn(tabProp, + "propagatedGlb_phi", + propagatedGlb_phi, + "expected position of muon track extrapolated to GE11 surface" + "
(float, global phi coordinates, rad)"); + addColumn(tabProp, + "propagatedGlb_errX", + propagatedGlb_errX, + "uncertainty on position of muon track extrapolated to GE11 surface" + "
(float, global x coordinates, cm)"); + addColumn(tabProp, + "propagatedGlb_errY", + propagatedGlb_errY, + "uncertainty on position of muon track extrapolated to GE11 surface" + "
(float, global y coordinates, cm)"); + addColumn(tabProp, + "propagatedGlb_rerr", + propagatedGlb_rerr, + "uncertainty on position of muon track extrapolated to GE11 surface" + "
(float, global radial (r) coordinates, cm)"); + addColumn(tabProp, + "propagatedGlb_phierr", + propagatedGlb_phierr, + "uncertainty on position of muon track extrapolated to GE11 surface" + "
(float, global phi coordinates, rad)"); + + addColumn(tabProp, + "propagated_EtaPartition_centerX", + propagated_EtaPartition_centerX, + "global X coordinate of the center of the etaPartition" + "
where the extrapolated muon track position falls" + "
(float, global x coordinates, cm)"); + addColumn(tabProp, + "propagated_EtaPartition_centerY", + propagated_EtaPartition_centerY, + "global Y coordinate of the center of the etaPartition" + "
where the extrapolated muon track position falls" + "
(float, global x coordinates, cm)"); + addColumn(tabProp, + "propagated_EtaPartition_phiMax", + propagated_EtaPartition_phiMax, + "upper edge in phi global coordinates of the etaPartition" + "
where the extrapolated muon track position falls" + "
(float, global phi coordinates, rad)"); + addColumn(tabProp, + "propagated_EtaPartition_phiMin", + propagated_EtaPartition_phiMin, + "lower edge in phi global coordinates of the etaPartition" + "
where the extrapolated muon track position falls" + "
(float, global phi coordinates, rad)"); + addColumn(tabProp, + "propagated_EtaPartition_rMax", + propagated_EtaPartition_rMax, + "upper edge in r global coordinates of the etaPartition" + "
where the extrapolated muon track position falls" + "
(float, global radial (r) coordinates, cm)"); + addColumn(tabProp, + "propagated_EtaPartition_rMin", + propagated_EtaPartition_rMin, + "lower edge in r global coordinates of the etaPartition" + "
where the extrapolated muon track position falls" + "
(float, global radial (r) coordinates, cm)"); + + ev.put(std::move(tabProp), "propagated"); + } +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleGEMMuonFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMMuonFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMMuonFiller.h new file mode 100644 index 0000000000000..b69d9656869cb --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMMuonFiller.h @@ -0,0 +1,55 @@ +#ifndef MuNtuple_MuNtupleGEMMuonFiller_h +#define MuNtuple_MuNtupleGEMMuonFiller_h + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/MuonReco/interface/Muon.h" +#include "DataFormats/MuonReco/interface/MuonFwd.h" + +#include "Geometry/GEMGeometry/interface/GEMGeometry.h" +#include "Geometry/Records/interface/MuonGeometryRecord.h" + +#include "RecoMuon/TrackingTools/interface/MuonServiceProxy.h" + +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "TrackingTools/Records/interface/TransientTrackRecord.h" + +class MuNtupleGEMMuonFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleGEMMuonFiller(const edm::ParameterSet &); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); + +protected: + /// Fill tree branches for a given event + void fill(edm::Event &ev) final; + + /// Get info from the ES by run + void getFromES(const edm::Run &run, const edm::EventSetup &environment) final; + + /// Get info from the ES for a given event + void getFromES(const edm::EventSetup &environment) final; + +private: + /// The RECO mu token + nano_mu::EDTokenHandle m_token; + + /// Fill matches table + bool m_fillPropagated; + + /// GEM Geometry + nano_mu::ESTokenHandle m_gemGeometry; + + /// Transient Track Builder + nano_mu::ESTokenHandle m_transientTrackBuilder; + + /// Muon service proxy + std::unique_ptr m_muonSP; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMRecHitFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMRecHitFiller.cc new file mode 100644 index 0000000000000..5222bdc53c2a5 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMRecHitFiller.cc @@ -0,0 +1,175 @@ +#include "DPGAnalysis/MuonTools/plugins/MuNtupleGEMRecHitFiller.h" + +#include + +MuNtupleGEMRecHitFiller::MuNtupleGEMRecHitFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller{config}, + m_token{config, consumesCollector(), "gemRecHitTag"}, + m_gemGeometry{consumesCollector()} { + produces(); +} + +void MuNtupleGEMRecHitFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "gemRecHits"); + desc.add("gemRecHitTag", edm::InputTag{"gemRecHits"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleGEMRecHitFiller::getFromES(const edm::Run& run, const edm::EventSetup& environment) { + m_gemGeometry.getFromES(environment); +} + +void MuNtupleGEMRecHitFiller::fill(edm::Event& ev) { + unsigned int nRecHits{0}; + + std::vector cluster_size; + std::vector firstClusterStrip; + std::vector bx; + + std::vector region; + std::vector chamber; + std::vector layer; + std::vector etaPartition; + + std::vector loc_r; + std::vector loc_phi; + std::vector loc_x; + std::vector loc_y; + std::vector loc_z; + + std::vector g_r; + std::vector g_phi; + std::vector g_x; + std::vector g_y; + std::vector g_z; + + auto gemRecHits = m_token.conditionalGet(ev); + + if (gemRecHits.isValid()) { + auto rechitIt = gemRecHits->begin(); + auto recHitEnd = gemRecHits->end(); + + for (; rechitIt != recHitEnd; ++rechitIt) { + GEMDetId gemDetId{rechitIt->gemId()}; + + etaPartition.push_back(gemDetId.roll()); + region.push_back(gemDetId.region()); + chamber.push_back(gemDetId.chamber()); + layer.push_back(gemDetId.layer()); + + const BoundPlane& surface = m_gemGeometry->idToDet(gemDetId)->surface(); + const auto rechit_local_pos = rechitIt->localPosition(); + const auto rechit_global_pos = surface.toGlobal(rechit_local_pos); + + loc_r.push_back(rechit_local_pos.perp()); + loc_phi.push_back(rechit_local_pos.phi()); + loc_x.push_back(rechit_local_pos.x()); + loc_y.push_back(rechit_local_pos.y()); + loc_z.push_back(rechit_local_pos.z()); + + g_r.push_back(rechit_global_pos.perp()); + g_phi.push_back(rechit_global_pos.phi()); + g_x.push_back(rechit_global_pos.x()); + g_y.push_back(rechit_global_pos.y()); + g_z.push_back(rechit_global_pos.z()); + + bx.push_back(rechitIt->BunchX()); + cluster_size.push_back(rechitIt->clusterSize()); + firstClusterStrip.push_back(rechitIt->firstClusterStrip()); + + ++nRecHits; + } + } + + auto table = std::make_unique(nRecHits, m_label, false, false); + + table->setDoc("GEM rec-hit information"); + + addColumn(table, "cluster_size", cluster_size, "recHit cluster size (strip multiplicity)"); + addColumn( + table, "firstClusterStrip", firstClusterStrip, "index of the first readout strip associated to the cluster"); + addColumn(table, "bx", bx, "bunch crossing associated to the recHit"); + + addColumn(table, + "region", + region, + "GE11 region where the hit is reconstructed" + "
(positive endcap: +1, negative endcap: -1"); + addColumn(table, + "chamber", + chamber, + "GE11 superchamber where the hit is reconstructed" + "
(chambers numbered from 0 to 35"); + addColumn(table, + "layer", + layer, + "GE11 layer where the hit is reconstructed" + "
(layer1: 1, layer2: 2"); + addColumn(table, + "etaPartition", + etaPartition, + "GE11 eta partition where the hit is reconstructed" + "
(partitions numbered from 1 to 8"); + + addColumn(table, + "loc_r", + loc_r, + "Local recHit position" + "
(local layer radial coordinates, cm)"); + addColumn(table, + "loc_phi", + loc_phi, + "Local recHit position" + "
(local layer azimuthal coordinates, rad)"); + addColumn(table, + "loc_x", + loc_x, + "Local recHit position" + "
(local layer x coordinates, cm)"); + addColumn(table, + "loc_y", + loc_y, + "Local recHit position" + "
(local layer y coordinates, cm)"); + addColumn(table, + "loc_z", + loc_z, + "Local recHit position" + "
(local layer z coordinates, cm)"); + + addColumn(table, + "g_r", + g_r, + "Global recHit position" + "
(global radial coordinates, cm)"); + addColumn(table, + "g_phi", + g_phi, + "Global recHit position" + "
(global azimuthal coordinates, rad)"); + addColumn(table, + "g_x", + g_x, + "Global recHit position" + "
(global x coordinates, cm)"); + addColumn(table, + "g_y", + g_y, + "Global recHit position" + "
(global y coordinates, cm)"); + addColumn(table, + "g_z", + g_z, + "Global recHit position" + "
(global z coordinates, cm)"); + + ev.put(std::move(table)); +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleGEMRecHitFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMRecHitFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMRecHitFiller.h new file mode 100644 index 0000000000000..3c57ab361f577 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMRecHitFiller.h @@ -0,0 +1,37 @@ +#ifndef MuNtuple_MuNtupleGEMRecHitFiller_h +#define MuNtuple_MuNtupleGEMRecHitFiller_h + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "Geometry/GEMGeometry/interface/GEMGeometry.h" +#include "Geometry/Records/interface/MuonGeometryRecord.h" + +#include "DataFormats/GEMRecHit/interface/GEMRecHitCollection.h" + +class MuNtupleGEMRecHitFiller : public MuNtupleBaseFiller { +public: + //Constructor + MuNtupleGEMRecHitFiller(const edm::ParameterSet &); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); + +protected: + /// Fill tree branches for a given events + void fill(edm::Event &ev) final; + + /// Get info from the ES by run + void getFromES(const edm::Run &run, const edm::EventSetup &environment) final; + +private: + /// The rec-hit token + nano_mu::EDTokenHandle m_token; + + /// GEM Geometry + nano_mu::ESTokenHandle m_gemGeometry; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMSegmentFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMSegmentFiller.cc new file mode 100644 index 0000000000000..036b902a0695a --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMSegmentFiller.cc @@ -0,0 +1,126 @@ +#include "DPGAnalysis/MuonTools/plugins/MuNtupleGEMSegmentFiller.h" + +#include + +MuNtupleGEMSegmentFiller::MuNtupleGEMSegmentFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller(config), + m_token{config, consumesCollector(), "gemSegmentTag"}, + m_trackingGeometry{consumesCollector()} { + produces(); +} + +void MuNtupleGEMSegmentFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "gemSegments"); + desc.add("gemSegmentTag", edm::InputTag{"gemSegments"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleGEMSegmentFiller::getFromES(const edm::EventSetup& environment) { + m_trackingGeometry.getFromES(environment); +} + +void MuNtupleGEMSegmentFiller::fill(edm::Event& ev) { + unsigned int nSegments{0}; + + std::vector region; + std::vector ring; + std::vector station; + std::vector posLoc_x; + std::vector posLoc_y; + std::vector posLoc_z; + std::vector dirLoc_x; + std::vector dirLoc_y; + std::vector dirLoc_z; + + std::vector posGlb_x; + std::vector posGlb_y; + std::vector posGlb_z; + + std::vector posGlb_phi; + std::vector posGlb_eta; + std::vector dirGlb_phi; + std::vector dirGlb_eta; + + std::vector time; + std::vector time_err; + std::vector chi2; + + auto segments = m_token.conditionalGet(ev); + + if (segments.isValid()) { + for (auto range_iter = segments->id_begin(); range_iter != segments->id_end(); range_iter++) { + const auto range = segments->get(*range_iter); + + for (auto segment = range.first; segment != range.second; ++segment) { + region.push_back((*range_iter).region()); + ring.push_back((*range_iter).ring()); + station.push_back((*range_iter).station()); + + auto pos = segment->localPosition(); + auto dir = segment->localDirection(); + + posLoc_x.push_back(pos.x()); + posLoc_y.push_back(pos.y()); + posLoc_z.push_back(pos.z()); + + dirLoc_x.push_back(dir.x()); + dirLoc_y.push_back(dir.y()); + dirLoc_z.push_back(dir.z()); + + const GeomDet* geomDet = m_trackingGeometry->idToDet(segment->geographicalId()); + auto posGlb = geomDet->toGlobal(pos); + auto dirGlb = geomDet->toGlobal(dir); + + posGlb_x.push_back(posGlb.x()); + posGlb_y.push_back(posGlb.y()); + posGlb_z.push_back(posGlb.z()); + posGlb_phi.push_back(posGlb.phi()); + posGlb_eta.push_back(posGlb.eta()); + + dirGlb_phi.push_back(dirGlb.phi()); + dirGlb_eta.push_back(dirGlb.eta()); + + time.push_back(segment->bunchX()); + chi2.push_back(segment->chi2()); + + ++nSegments; + } + } + } + + auto table = std::make_unique(nSegments, m_label, false, false); + + table->setDoc("GEM segment information"); + + addColumn(table, "region", region, ""); + addColumn(table, "ring", ring, ""); + addColumn(table, "station", station, ""); + + addColumn(table, "posLoc_x", posLoc_x, "position x in local coordinates - cm"); + addColumn(table, "posLoc_y", posLoc_y, "position y in local coordinates - cm"); + addColumn(table, "posLoc_z", posLoc_z, "position z in local coordinates - cm"); + + addColumn(table, "dirLoc_x", dirLoc_x, "direction x in local coordinates"); + addColumn(table, "dirLoc_y", dirLoc_y, "direction y in local coordinates"); + addColumn(table, "dirLoc_z", dirLoc_z, "direction z in local coordinates"); + + addColumn(table, "posGlb_x", posGlb_x, "position x in global coordinates - cm"); + addColumn(table, "posGlb_y", posGlb_y, "position y in global coordinates - cm"); + addColumn(table, "posGlb_z", posGlb_z, "position z in global coordinates - cm"); + + addColumn(table, "dirGlb_phi", dirGlb_phi, "direction phi in global coordinates - radians [-pi:pi]"); + addColumn(table, "dirGlb_eta", dirGlb_eta, "direction eta in global coordinates"); + + addColumn(table, "time", time, ""); + addColumn(table, "chi2", chi2, ""); + + ev.put(std::move(table)); +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleGEMSegmentFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleGEMSegmentFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMSegmentFiller.h new file mode 100644 index 0000000000000..bc24fcc235a50 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleGEMSegmentFiller.h @@ -0,0 +1,37 @@ +#ifndef MuNtuple_MuNtupleGEMSegmentFiller_h +#define MuNtuple_MuNtupleGEMSegmentFiller_h + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/GEMRecHit/interface/GEMSegmentCollection.h" + +#include "Geometry/CommonDetUnit/interface/GlobalTrackingGeometry.h" +#include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h" + +class MuNtupleGEMSegmentFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleGEMSegmentFiller(const edm::ParameterSet&); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + /// Fill tree branches for a given event + void fill(edm::Event& ev) final; + + /// Get info from the ES for a given event + void getFromES(const edm::EventSetup& environment) final; + +private: + /// The segment token + nano_mu::EDTokenHandle m_token; + + /// Tracking Geometry + nano_mu::ESTokenHandle m_trackingGeometry; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.cc new file mode 100644 index 0000000000000..f0d57883aabb8 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.cc @@ -0,0 +1,494 @@ +/** \class MuNtupleMuonFiller MuNtupleMuonFiller.cc DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.cc + * + * Helper class : the muon filler + * + * \author L. Lunerti (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.h" + +#include "DataFormats/DetId/interface/DetId.h" +#include "DataFormats/MuonDetId/interface/DTChamberId.h" +#include "DataFormats/MuonReco/interface/MuonChamberMatch.h" +#include "DataFormats/MuonReco/interface/MuonSegmentMatch.h" +#include "DataFormats/MuonReco/interface/MuonSelectors.h" + +#include "DataFormats/DTRecHit/interface/DTRecSegment4D.h" +#include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h" +#include "DataFormats/TrackingRecHit/interface/TrackingRecHitFwd.h" +#include "DataFormats/MuonDetId/interface/DTChamberId.h" + +#include "DataFormats/Math/interface/deltaR.h" + +#include "TString.h" +#include "TRegexp.h" + +#include +#include + +MuNtupleMuonFiller::MuNtupleMuonFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller(config), + m_muToken{config, consumesCollector(), "muonTag"}, + m_dtSegmentToken{config, consumesCollector(), "dtSegmentTag"}, + m_primaryVerticesToken{config, consumesCollector(), "primaryVerticesTag"}, + m_trigResultsToken{config, consumesCollector(), "trigResultsTag"}, + m_trigEventToken{config, consumesCollector(), "trigEventTag"}, + m_fillMatches{config.getParameter("fillMatches")}, + m_trigName{config.getParameter("trigName")}, + m_isoTrigName{config.getParameter("isoTrigName")}, + m_dtGeometry{consumesCollector()} { + produces(); + if (m_fillMatches) { + produces("matches"); + produces("staMatches"); + } +} + +void MuNtupleMuonFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "muon"); + desc.add("muonTag", edm::InputTag{"muons"}); + desc.add("dtSegmentTag", edm::InputTag{"dt4DSegments"}); + desc.add("primaryVerticesTag", edm::InputTag{"offlinePrimaryVertices"}); + + desc.add("fillMatches", true); + + desc.add("trigEventTag", edm::InputTag{"none"}); + desc.add("trigResultsTag", edm::InputTag{"none"}); + + desc.add("trigName", "none"); + desc.add("isoTrigName", "none"); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleMuonFiller::getFromES(const edm::Run& run, const edm::EventSetup& environment) { + m_dtGeometry.getFromES(environment); + + bool changed{true}; + m_hltConfig.init(run, environment, "HLT", changed); + + const bool enableWildcard{true}; + + TString tName = TString(m_trigName); + TRegexp tNamePattern = TRegexp(tName, enableWildcard); + + for (unsigned iPath = 0; iPath < m_hltConfig.size(); ++iPath) { + TString pathName = TString(m_hltConfig.triggerName(iPath)); + if (pathName.Contains(tNamePattern)) + m_trigIndices.push_back(static_cast(iPath)); + } + + tName = TString(m_isoTrigName); + tNamePattern = TRegexp(tName, enableWildcard); + + for (unsigned iPath = 0; iPath < m_hltConfig.size(); ++iPath) { + TString pathName = TString(m_hltConfig.triggerName(iPath)); + if (pathName.Contains(tNamePattern)) + m_isoTrigIndices.push_back(static_cast(iPath)); + } +} + +void MuNtupleMuonFiller::fill(edm::Event& ev) { + unsigned int nMuons{0}; + + std::vector pt; + std::vector phi; + std::vector eta; + std::vector charge; + + std::vector isGlobal; + std::vector isStandalone; + std::vector isTracker; + std::vector isTrackerArb; + std::vector isRPC; + + std::vector firesIsoTrig; + std::vector firesTrig; + + std::vector isLoose; + std::vector isMedium; + std::vector isTight; + + std::vector trkIso03; + std::vector pfIso04; + + std::vector trk_dxy; + std::vector trk_dz; + + std::vector trk_algo; + std::vector trk_origAlgo; + + std::vector trk_numberOfValidPixelHits; + std::vector trk_numberOfValidTrackerLayers; + + std::vector trkMu_stationMask; + std::vector trkMu_numberOfMatchedStations; + std::vector trkMu_numberOfMatchedRPCLayers; + + std::vector staMu_numberOfValidMuonHits; + std::vector staMu_normChi2; + std::vector glbMu_normChi2; + + std::vector nMatches; + std::vector staMu_nMatchSeg; + + std::vector matches_begin; + std::vector matches_end; + + std::vector staMatches_begin; + std::vector staMatches_end; + + std::vector matches_wheel; + std::vector matches_sector; + std::vector matches_station; + + std::vector matches_x; + std::vector matches_y; + + std::vector matches_phi; + std::vector matches_eta; + std::vector matches_edgeX; + std::vector matches_edgeY; + + std::vector matches_dXdZ; + std::vector matches_dYdZ; + + std::vector staMatches_MuSegIdx; + + auto muons = m_muToken.conditionalGet(ev); + auto segments = m_dtSegmentToken.conditionalGet(ev); + auto vtxs = m_primaryVerticesToken.conditionalGet(ev); + + auto triggerResults = m_trigResultsToken.conditionalGet(ev); + auto triggerEvent = m_trigEventToken.conditionalGet(ev); + + if (muons.isValid() && segments.isValid() && vtxs.isValid()) { + for (const auto& muon : (*muons)) { + pt.push_back(muon.pt()); + eta.push_back(muon.eta()); + phi.push_back(muon.phi()); + charge.push_back(muon.charge()); + + isGlobal.push_back(muon.isGlobalMuon()); + isStandalone.push_back(muon.isStandAloneMuon()); + isTracker.push_back(muon.isTrackerMuon()); + isTrackerArb.push_back(muon::isGoodMuon(muon, muon::TrackerMuonArbitrated)); + isRPC.push_back(muon.isRPCMuon()); + + if (triggerResults.isValid() && triggerEvent.isValid()) { + const auto& triggerObjects = triggerEvent->getObjects(); + + bool hasIsoTrig = hasTrigger(m_isoTrigIndices, triggerObjects, triggerEvent, muon); + bool hasTrig = hasTrigger(m_trigIndices, triggerObjects, triggerEvent, muon); + + firesIsoTrig.push_back(hasIsoTrig); + firesTrig.push_back(hasTrig); + + } else { + firesIsoTrig.push_back(false); + firesTrig.push_back(false); + } + + isLoose.push_back(muon.passed(reco::Muon::CutBasedIdLoose)); + isMedium.push_back(muon.passed(reco::Muon::CutBasedIdMedium)); + isTight.push_back(muon.passed(reco::Muon::CutBasedIdTight)); + + trkIso03.push_back(computeTrkIso(muon.isolationR03(), muon.pt())); + pfIso04.push_back(computePFIso(muon.pfIsolationR04(), muon.pt())); + + //INNER TRACK VARIABLES + if (!muon.innerTrack().isNull()) { + const reco::TrackRef innerTrackRef = muon.innerTrack(); + const reco::Vertex& vertex = vtxs->at(0); + + trk_dxy.push_back(innerTrackRef->dxy(vertex.position())); + trk_dz.push_back(innerTrackRef->dz(vertex.position())); + trk_algo.push_back(innerTrackRef->algo()); + trk_origAlgo.push_back(innerTrackRef->originalAlgo()); + trk_numberOfValidPixelHits.push_back(innerTrackRef->hitPattern().numberOfValidPixelHits()); + trk_numberOfValidTrackerLayers.push_back(innerTrackRef->hitPattern().trackerLayersWithMeasurement()); + + } else { + trk_dxy.push_back(MuNtupleBaseFiller::DEFAULT_DOUBLE_VAL); + trk_dz.push_back(MuNtupleBaseFiller::DEFAULT_DOUBLE_VAL); + trk_algo.push_back(MuNtupleBaseFiller::DEFAULT_INT_VAL_POS); + trk_origAlgo.push_back(MuNtupleBaseFiller::DEFAULT_INT_VAL_POS); + trk_numberOfValidPixelHits.push_back(MuNtupleBaseFiller::DEFAULT_INT_VAL_POS); + trk_numberOfValidTrackerLayers.push_back(MuNtupleBaseFiller::DEFAULT_INT_VAL_POS); + } + + //TRACKER / RPC MUON VARIABLES + trkMu_stationMask.push_back(muon.stationMask()); + trkMu_numberOfMatchedStations.push_back(muon.numberOfMatchedStations()); + trkMu_numberOfMatchedRPCLayers.push_back(muon.numberOfMatchedRPCLayers()); + + //STANDALONE MUON VARIABLES + if (muon.isStandAloneMuon()) { + const reco::TrackRef outerTrackRef = muon.outerTrack(); + staMu_numberOfValidMuonHits.push_back(outerTrackRef->hitPattern().numberOfValidMuonHits()); + staMu_normChi2.push_back(outerTrackRef->chi2() / outerTrackRef->ndof()); + + } else { + staMu_numberOfValidMuonHits.push_back(MuNtupleBaseFiller::DEFAULT_INT_VAL_POS); + staMu_normChi2.push_back(MuNtupleBaseFiller::DEFAULT_DOUBLE_VAL_POS); + } + + //GLOBAL MUON VARIABLES + if (muon.isGlobalMuon()) { + const reco::TrackRef globalTrackRef = muon.globalTrack(); + glbMu_normChi2.push_back(globalTrackRef->chi2() / globalTrackRef->ndof()); + + } else { + glbMu_normChi2.push_back(MuNtupleBaseFiller::DEFAULT_DOUBLE_VAL_POS); + } + + size_t iMatches = 0; + size_t iSegMatches = 0; + + if (m_fillMatches) { + matches_begin.push_back(matches_wheel.size()); + + if (muon.isMatchesValid()) { + for (const auto& match : muon.matches()) { + if (iMatches < 16 && match.id.det() == DetId::Muon && match.id.subdetId() == MuonSubdetId::DT) { + DTChamberId dtId(match.id.rawId()); + const auto chamb = m_dtGeometry->chamber(static_cast(match.id)); + + matches_wheel.push_back(dtId.wheel()); + matches_sector.push_back(dtId.sector()); + matches_station.push_back(dtId.station()); + + matches_x.push_back(match.x); + matches_y.push_back(match.y); + + matches_phi.push_back(chamb->toGlobal(LocalPoint(match.x, match.y, 0.)).phi()); + matches_eta.push_back(chamb->toGlobal(LocalPoint(match.x, match.y, 0.)).eta()); + + matches_edgeX.push_back(match.edgeX); + matches_edgeY.push_back(match.edgeY); + + matches_dXdZ.push_back(match.dXdZ); + matches_dYdZ.push_back(match.dYdZ); + + ++iMatches; + } + } + } + + matches_end.push_back(matches_wheel.size()); + + //SEGMENT MATCHING VARIABLES + + staMatches_begin.push_back(staMatches_MuSegIdx.size()); + + if (!muon.outerTrack().isNull()) { + reco::TrackRef outerTrackRef = muon.outerTrack(); + + auto recHitIt = outerTrackRef->recHitsBegin(); + auto recHitEnd = outerTrackRef->recHitsEnd(); + + for (; recHitIt != recHitEnd; ++recHitIt) { + DetId detId = (*recHitIt)->geographicalId(); + + if (detId.det() == DetId::Muon && detId.subdetId() == MuonSubdetId::DT) { + const auto dtSegmentSta = dynamic_cast((*recHitIt)); + int iSeg = 0; + + for (const auto& segment : (*segments)) { + if (dtSegmentSta && dtSegmentSta->chamberId().station() == segment.chamberId().station() && + dtSegmentSta->chamberId().wheel() == segment.chamberId().wheel() && + dtSegmentSta->chamberId().sector() == segment.chamberId().sector() && + std::abs(dtSegmentSta->localPosition().x() - segment.localPosition().x()) < 0.01 && + std::abs(dtSegmentSta->localPosition().y() - segment.localPosition().y()) < 0.01 && + std::abs(dtSegmentSta->localDirection().x() - segment.localDirection().x()) < 0.01 && + std::abs(dtSegmentSta->localDirection().y() - segment.localDirection().y()) < 0.01) { + staMatches_MuSegIdx.push_back(iSeg); + ++iSegMatches; + } + + ++iSeg; + } //loop over segments + } + + } //loop over recHits + } + + staMatches_end.push_back(staMatches_MuSegIdx.size()); + } + + nMatches.push_back(iMatches); + staMu_nMatchSeg.push_back(iSegMatches); + + ++nMuons; + } + } + + auto table = std::make_unique(nMuons, m_label, false, false); + + table->setDoc("RECO muon information"); + + addColumn(table, "pt", pt, "muon pT - GeV/c"); + addColumn(table, "phi", phi, "muon phi - rad"); + addColumn(table, "eta", eta, "muon eta"); + addColumn(table, "charge", charge, "muon charge"); + + addColumn(table, "isGlobal", isGlobal, ""); + addColumn(table, "isStandalone", isStandalone, ""); + addColumn(table, "isTracker", isTracker, ""); + addColumn(table, "isTrackerArb", isTrackerArb, ""); + addColumn(table, "isRPC", isRPC, ""); + + addColumn(table, + "firesIsoTrig", + firesIsoTrig, + "True if the muon is matched to an isolated trigger" + "
specified in the ntuple config file"); + + addColumn(table, + "firesTrig", + firesTrig, + "True if the muon is matched to a (non isolated)trigger" + "
specified in the ntuple config file"); + + addColumn(table, "isLoose", isLoose, "Loose muon ID"); + addColumn(table, "isMedium", isMedium, "Medium muon ID"); + addColumn(table, "isTight", isTight, "Tight muon ID"); + + addColumn(table, "trkIso03", trkIso03, "Relative tracker isolation (0.3 cone)"); + addColumn(table, "pfIso04", pfIso04, "Relative PF-isolation (delta beta corrected, 0.4 cone)"); + + addColumn(table, "trk_dxy", trk_dxy, "Inner track dxy parameter with respect to the primary vertex - cm"); + addColumn(table, "trk_dz", trk_dz, "Inner track dz parameter with respect to the primary vertex - cm"); + + addColumn(table, "trk_algo", trk_algo, "iterative tracking algorithm used to build the inner track"); + addColumn(table, + "trk_origAlgo", + trk_origAlgo, + "Original (pre muon iterations) iterative tracking algorithm used to build the inner track"); + + addColumn(table, "trk_numberOfValidPixelHits", trk_numberOfValidPixelHits, "Number of valid pixel hits"); + addColumn(table, "trk_numberOfValidTrackerLayers", trk_numberOfValidTrackerLayers, "Number of valid tracker layers"); + + addColumn(table, + "trkMu_stationMask", + trkMu_stationMask, + "Bit map of stations with tracks within given distance (in cm) of chamber edges"); + addColumn(table, "trkMu_numberOfMatchedStations", trkMu_numberOfMatchedStations, "Number of matched DT/CSC stations"); + addColumn(table, "trkMu_numberOfMatchedRPCLayers", trkMu_numberOfMatchedRPCLayers, "Number of matched RPC layers"); + + addColumn(table, "staMu_numberOfValidMuonHits", staMu_numberOfValidMuonHits, "Number of valid muon hits"); + + addColumn(table, "staMu_normChi2", staMu_normChi2, "chi2/ndof (standalone track)"); + + addColumn(table, "glbMu_normChi2", glbMu_normChi2, "chi2/ndof (global track)"); + + addColumn(table, "nMatches", nMatches, "Number of muon chamber matches (DT only)"); + addColumn(table, "staMu_nMatchSeg", staMu_nMatchSeg, "Number of segments used in the standalone track (DT only)"); + + addColumn( + table, "matches_begin", matches_begin, "begin of range of quantities for a given muon in the matches_* vectors"); + addColumn(table, "matches_end", matches_end, "end of range of quantities for a given muon in the matches_* vectors"); + + addColumn(table, + "staMatches_begin", + staMatches_begin, + "begin of range of quantities for a given muon in the matches_staMuSegIdx vector"); + addColumn(table, + "staMatches_end", + staMatches_end, + "end of range of quantities for a given muon in the matches_staMuSegIdx vector"); + + ev.put(std::move(table)); + + if (m_fillMatches) { + auto sum = [](std::vector v) { return std::accumulate(v.begin(), v.end(), 0); }; + + auto tabMatches = std::make_unique(sum(nMatches), m_label + "_matches", false, false); + + tabMatches->setDoc("Size of RECO muon matches_* vectors"); + + addColumn(tabMatches, "x", matches_x, "X position of the extrapolated track on the matched chamber"); + addColumn(tabMatches, "y", matches_y, "Y position of the extrapolated track on the matched chamber"); + + addColumn( + tabMatches, "phi", matches_phi, "Phi of the (x,y) position on the matched chamber (global reference frame)"); + addColumn( + tabMatches, "eta", matches_eta, " Eta of the (x,y) position on the matched chamber (global reference frame)"); + + addColumn(tabMatches, "dXdZ", matches_dXdZ, "dXdZ of the extrapolated track on the matched chamber"); + addColumn(tabMatches, "dYdZ", matches_dYdZ, "dYdZ of the extrapolated track on the matched chamber"); + + ev.put(std::move(tabMatches), "matches"); + + auto tabStaMatches = + std::make_unique(sum(staMu_nMatchSeg), m_label + "_staMatches", false, false); + + tabStaMatches->setDoc("Size of RECO muon staMatches_* vector"); + + addColumn(tabStaMatches, + "MuSegIdx", + staMatches_MuSegIdx, + "Index of DT segments used in the standalone track it corresponds" + "
to the index of a given segment in the ntuple seg_* collection"); + + ev.put(std::move(tabStaMatches), "staMatches"); + } +} + +float MuNtupleMuonFiller::computeTrkIso(const reco::MuonIsolation& isolation, float muonPt) { + return isolation.sumPt / muonPt; +} + +float MuNtupleMuonFiller::computePFIso(const reco::MuonPFIsolation& pfIsolation, float muonPt) { + return (pfIsolation.sumChargedHadronPt + + std::max(0., pfIsolation.sumNeutralHadronEt + pfIsolation.sumPhotonEt - 0.5 * pfIsolation.sumPUPt)) / + muonPt; +} + +bool MuNtupleMuonFiller::hasTrigger(std::vector& trigIndices, + const trigger::TriggerObjectCollection& trigObjs, + edm::Handle& trigEvent, + const reco::Muon& muon) { + double matchDeltaR = 999.; + + for (const auto& trigIndex : trigIndices) { + const std::vector moduleLabels(m_hltConfig.moduleLabels(trigIndex)); + + // find index of the last module: + const unsigned moduleIndex = m_hltConfig.size(trigIndex) - 2; + + // find index of HLT trigger name: + const unsigned hltFilterIndex = trigEvent->filterIndex(edm::InputTag(moduleLabels[moduleIndex], "", "HLT")); + + if (hltFilterIndex < trigEvent->sizeFilters()) { + const trigger::Keys triggerKeys(trigEvent->filterKeys(hltFilterIndex)); + const trigger::Vids triggerVids(trigEvent->filterIds(hltFilterIndex)); + + const size_t nTriggers = triggerVids.size(); + + for (size_t iTrig = 0; iTrig < nTriggers; ++iTrig) { + // loop over all trigger objects: + const trigger::TriggerObject trigObject = trigObjs[triggerKeys[iTrig]]; + + double dR = deltaR(muon, trigObject); + + if (dR < matchDeltaR) + matchDeltaR = dR; + + } // loop over different trigger objects + + } // if trigger is in event (should apply hltFilter with used trigger...) + + } // loop over muon candidates + + return matchDeltaR < 0.1; //CB should get it programmable +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleMuonFiller); \ No newline at end of file diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.h new file mode 100644 index 0000000000000..93726c49f6ea2 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.h @@ -0,0 +1,80 @@ +#ifndef MuNtuple_MuNtupleMuonFiller_h +#define MuNtuple_MuNtupleMuonFiller_h + +/** \class MuNtupleMuonFiller MuNtupleMuonFiller.h DPGAnalysis/MuonTools/plugins/MuNtupleMuonFiller.h + * + * Helper class : the muon filler + * + * \author L. Lunerti (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/MuonReco/interface/Muon.h" +#include "DataFormats/VertexReco/interface/Vertex.h" +#include "DataFormats/MuonReco/interface/MuonFwd.h" +#include "DataFormats/MuonReco/interface/MuonIsolation.h" +#include "DataFormats/MuonReco/interface/MuonPFIsolation.h" +#include "DataFormats/DTRecHit/interface/DTRecSegment4D.h" + +#include "DataFormats/Common/interface/TriggerResults.h" +#include "DataFormats/HLTReco/interface/TriggerEvent.h" +#include "DataFormats/HLTReco/interface/TriggerObject.h" +#include "HLTrigger/HLTcore/interface/HLTConfigProvider.h" + +class MuNtupleMuonFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleMuonFiller(const edm::ParameterSet &l); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); + +protected: + /// Fill tree branches for a given events + void fill(edm::Event &ev) final; + + /// Get info from the ES by run + void getFromES(const edm::Run &run, const edm::EventSetup &environment) final; + +private: + /// Tokens + nano_mu::EDTokenHandle m_muToken; + nano_mu::EDTokenHandle m_dtSegmentToken; + nano_mu::EDTokenHandle> m_primaryVerticesToken; + + nano_mu::EDTokenHandle m_trigResultsToken; + nano_mu::EDTokenHandle m_trigEventToken; + + /// Fill matches table + bool m_fillMatches; + + /// Name of the triggers used by muon filler for trigger matching + std::string m_trigName; + std::string m_isoTrigName; + + /// DT Geometry + nano_mu::ESTokenHandle m_dtGeometry; + + /// HLT config provider + HLTConfigProvider m_hltConfig; + + /// Indices of the triggers used by muon filler for trigger matching + std::vector m_trigIndices; + std::vector m_isoTrigIndices; + + bool hasTrigger(std::vector &trigIndices, + const trigger::TriggerObjectCollection &trigObjs, + edm::Handle &trigEvent, + const reco::Muon &muon); + + float computeTrkIso(const reco::MuonIsolation &isolation, float muonPt); + float computePFIso(const reco::MuonPFIsolation &PFisolation, float muonPt); +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.cc new file mode 100644 index 0000000000000..a42f3f932a029 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.cc @@ -0,0 +1,108 @@ +/** \class MuNtupleDigiFiller MuNtupleDigiFiller.cc DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.cc + * + * Helper class : the digi filler for RPC digis + * + * \author Eliza Melo Da Costa (UERJ) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.h" + +#include "DataFormats/MuonDetId/interface/RPCDetId.h" + +#include + +MuNtupleRPCDigiFiller::MuNtupleRPCDigiFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller(config), m_token{config, consumesCollector(), "rpcDigiTag"} { + produces(); +} + +void MuNtupleRPCDigiFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "rpcDigis"); + desc.add("rpcDigiTag", edm::InputTag{"muonRPCDigis"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleRPCDigiFiller::fill(edm::Event& ev) { + unsigned int nDigis{0}; + + std::vector strip; + std::vector bx; + + std::vector region; + std::vector ring; + std::vector station; + std::vector layer; + std::vector sector; + std::vector subsector; + std::vector roll; + + std::vector rawId; + + auto rpcDigis = m_token.conditionalGet(ev); + + if (rpcDigis.isValid()) { + auto detUnitIt = rpcDigis->begin(); + auto detUnitEnd = rpcDigis->end(); + + for (; detUnitIt != detUnitEnd; ++detUnitIt) { + const auto& [rpcDetId, range] = (*detUnitIt); + + for (auto digiIt = range.first; digiIt != range.second; ++digiIt) { + strip.push_back(digiIt->strip()); + bx.push_back(digiIt->bx()); + + region.push_back(rpcDetId.region()); + ring.push_back(rpcDetId.ring()); + station.push_back(rpcDetId.station()); + layer.push_back(rpcDetId.layer()); + sector.push_back(rpcDetId.sector()); + subsector.push_back(rpcDetId.subsector()); + roll.push_back(rpcDetId.roll()); + + rawId.push_back(rpcDetId.rawId()); + + ++nDigis; + } + } + } + + auto table = std::make_unique(nDigis, m_label, false, false); + + table->setDoc("RPC digi information"); + + addColumn(table, "strip", strip, "Strip number fired in a given BX"); + addColumn(table, "bx", bx, "Bunch crossing number"); + + addColumn(table, "region", region, "0: Barrel, +-1: Endcap"); + addColumn(table, "ring", ring, "Ring id: Wheel number in Barrel (from -2 to +2) Ring Number in Endcap (from 1 to 3)"); + addColumn(table, "station", station, "Chambers at same R in barrel, chamber at same Z ion endcap"); + addColumn(table, + "layer", + layer, + "Layer id: in station 1 and 2 for barrel, we have two layers of chambers: " + "layer 1 is the inner chamber and layer 2 is the outer chamber"); + addColumn(table, "sector", sector, "Group of chambers at same phi"); + addColumn(table, + "subsector", + subsector, + "Some sectors are divided along the phi direction in subsectors " + "(from 1 to 4 in Barrel, from 1 to 6 in Endcap)"); + addColumn(table, + "roll", + roll, + "Roll id (also known as eta partition) each chamber is divided along the strip direction in"); + + addColumn(table, "rawId", rawId, "Unique detector unit ID"); + + ev.put(std::move(table)); +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleRPCDigiFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.h new file mode 100644 index 0000000000000..8eed2df79171c --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.h @@ -0,0 +1,37 @@ +#ifndef MuNtuple_MuNtupleRPCDigiFiller_h +#define MuNtuple_MuNtupleRPCDigiFiller_h + +/** \class MuNtupleRPCDigiFiller MuNtupleRPCDigiFiller.h DPGAnalysis/MuonTools/plugins/MuNtupleRPCDigiFiller.h + * + * Helper class : the digi filler for RPC digis + * + * \author Eliza Melo Da Costa (UERJ) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "DataFormats/RPCDigi/interface/RPCDigiCollection.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +class MuNtupleRPCDigiFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleRPCDigiFiller(const edm::ParameterSet&); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + /// Fill tree branches for a given events + void fill(edm::Event& ev) final; + +private: + /// The digi token + nano_mu::EDTokenHandle m_token; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.cc b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.cc new file mode 100644 index 0000000000000..85a6b385e1fd1 --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.cc @@ -0,0 +1,121 @@ +/** \class MuNtupleDigiFiller MuNtupleDigiFiller.cc DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.cc + * + * Helper class : the digi filler for RPC RecHits + * + * \author Eliza Melo Da Costa (UERJ) + * + * + */ + +#include "DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.h" + +#include "DataFormats/MuonDetId/interface/RPCDetId.h" + +#include + +MuNtupleRPCRecHitFiller::MuNtupleRPCRecHitFiller(const edm::ParameterSet& config) + : MuNtupleBaseFiller(config), m_token{config, consumesCollector(), "rpcRecHitTag"} { + produces(); +} + +void MuNtupleRPCRecHitFiller::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + + desc.add("label", "rpcRecHits"); + desc.add("rpcRecHitTag", edm::InputTag{"rpcRecHits"}); + + descriptions.addWithDefaultLabel(desc); +} + +void MuNtupleRPCRecHitFiller::fill(edm::Event& ev) { + unsigned int nRecHit{0}; + + std::vector firstClusterStrip; + std::vector clusterSize; + std::vector bx; + std::vector time; + std::vector coordinateX; + std::vector coordinateY; + std::vector coordinateZ; + + std::vector region; + std::vector ring; + std::vector station; + std::vector layer; + std::vector sector; + std::vector subsector; + std::vector roll; + + std::vector rawId; + + auto rpcRecHits = m_token.conditionalGet(ev); + + if (rpcRecHits.isValid()) { + auto recHitIt = rpcRecHits->begin(); + auto recHitEnd = rpcRecHits->end(); + + for (; recHitIt != recHitEnd; recHitIt++) { + firstClusterStrip.push_back(recHitIt->firstClusterStrip()); + clusterSize.push_back(recHitIt->clusterSize()); + bx.push_back(recHitIt->BunchX()); + time.push_back(recHitIt->time()); + coordinateX.push_back(recHitIt->localPosition().x()); + coordinateY.push_back(recHitIt->localPosition().y()); + coordinateZ.push_back(recHitIt->localPosition().z()); + + auto rpcDetId = (RPCDetId)(*recHitIt).rpcId(); + region.push_back(rpcDetId.region()); + ring.push_back(rpcDetId.ring()); + station.push_back(rpcDetId.station()); + layer.push_back(rpcDetId.layer()); + sector.push_back(rpcDetId.sector()); + subsector.push_back(rpcDetId.subsector()); + roll.push_back(rpcDetId.roll()); + rawId.push_back(rpcDetId.rawId()); + ++nRecHit; + } + } + auto table = std::make_unique(nRecHit, m_label, false, false); + + table->setDoc("RPC rec-hit information"); + + addColumn(table, "bx", bx, "Bunch crossing number"); + addColumn(table, "time", time, "Float with time information in ns"); + addColumn(table, "firstClusterStrip", firstClusterStrip, "Lowest-numbered strip in the cluster"); + addColumn(table, "clusterSize", clusterSize, "Number of strips in the cluster"); + addColumn(table, "coordinateX", coordinateX, "Local x coordinate in cm"); + addColumn(table, "coordinateY", coordinateY, "Local y coordinate in cm"); + addColumn(table, "coordinateZ", coordinateZ, "Local z coordinate in cm"); + + addColumn(table, "region", region, "0: Barrel, +-1: Endcap"); + addColumn(table, + "ring", + ring, + "Ring id: Wheel number in Barrel (from -2 to +2) " + "Ring Number in Endcap (from 1 to 3)"); + addColumn(table, "station", station, "Chambers at same R in barrel, chamber at same Z ion endcap"); + addColumn(table, + "layer", + layer, + "Layer id: in station 1 and 2 for barrel, we have two layers of chambers: " + "layer 1 is the inner chamber and layer 2 is the outer chamber"); + addColumn(table, "sector", sector, "Group of chambers at same phi"); + addColumn(table, + "subsector", + subsector, + "Some sectors are divided along the phi direction in subsectors " + "(from 1 to 4 in Barrel, from 1 to 6 in Endcap)"); + addColumn(table, + "roll", + roll, + "Roll id (also known as eta partition) each chamber is divided along the strip direction in"); + + addColumn(table, "rawId", rawId, "Unique detector unit ID"); + + ev.put(std::move(table)); +} + +#include "FWCore/PluginManager/interface/ModuleDef.h" +#include "FWCore/Framework/interface/MakerMacros.h" + +DEFINE_FWK_MODULE(MuNtupleRPCRecHitFiller); diff --git a/DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.h b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.h new file mode 100644 index 0000000000000..2996fb734bd1a --- /dev/null +++ b/DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.h @@ -0,0 +1,37 @@ +#ifndef MuNtuple_MuNtupleRPCRecHitFiller_h +#define MuNtuple_MuNtupleRPCRecHitFiller_h + +/** \class MuNtupleRPCRecHitFiller MuNtupleRPCRecHitFiller.h DPGAnalysis/MuonTools/plugins/MuNtupleRPCRecHitFiller.h + * + * Helper class : the digi filler for RPC recHits + * + * \author Eliza Melo Da Costa (UERJ) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" + +#include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h" + +class MuNtupleRPCRecHitFiller : public MuNtupleBaseFiller { +public: + /// Constructor + MuNtupleRPCRecHitFiller(const edm::ParameterSet&); + + /// Fill descriptors + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + +protected: + /// Fill tree branches for a given events + void fill(edm::Event& ev) final; + +private: + /// The rec-hit token + nano_mu::EDTokenHandle m_token; +}; + +#endif diff --git a/DPGAnalysis/MuonTools/python/customiseMuNtuples_cff.py b/DPGAnalysis/MuonTools/python/customiseMuNtuples_cff.py new file mode 100644 index 0000000000000..419c58afd4cc2 --- /dev/null +++ b/DPGAnalysis/MuonTools/python/customiseMuNtuples_cff.py @@ -0,0 +1,12 @@ +import FWCore.ParameterSet.Config as cms + +def customiseForRunningOnMC(process, pathName) : + + if hasattr(process,"muNtupleProducer") : + print("[customiseForRunningOnMC]: updating ntuple input tags") + + process.muNtupleTwinMuxInFiller.dtTpTag = "none" + process.muNtupleTwinMuxOutFiller.dtTpTag = "none" + process.muNtupleTwinMuxInThFiller.dtTpTag = "none" + + return process diff --git a/DPGAnalysis/MuonTools/python/muNtupleProducer_cff.py b/DPGAnalysis/MuonTools/python/muNtupleProducer_cff.py new file mode 100644 index 0000000000000..1574db80464f0 --- /dev/null +++ b/DPGAnalysis/MuonTools/python/muNtupleProducer_cff.py @@ -0,0 +1,43 @@ +import FWCore.ParameterSet.Config as cms +from RecoMuon.TrackingTools.MuonServiceProxy_cff import MuonServiceProxy + +from DPGAnalysis.MuonTools.muNtupleDTDigiFiller_cfi import muNtupleDTDigiFiller +from DPGAnalysis.MuonTools.muNtupleRPCDigiFiller_cfi import muNtupleRPCDigiFiller +from DPGAnalysis.MuonTools.muNtupleGEMDigiFiller_cfi import muNtupleGEMDigiFiller + +from DPGAnalysis.MuonTools.muNtupleRPCRecHitFiller_cfi import muNtupleRPCRecHitFiller +from DPGAnalysis.MuonTools.muNtupleGEMRecHitFiller_cfi import muNtupleGEMRecHitFiller + +from DPGAnalysis.MuonTools.muNtupleDTSegmentFiller_cfi import muNtupleDTSegmentFiller +from DPGAnalysis.MuonTools.muNtupleGEMSegmentFiller_cfi import muNtupleGEMSegmentFiller + +from DPGAnalysis.MuonTools.muNtupleMuonFiller_cfi import muNtupleMuonFiller +from DPGAnalysis.MuonTools.muNtupleGEMMuonFiller_cfi import muNtupleGEMMuonFiller +muNtupleGEMMuonFiller.ServiceParameters = MuonServiceProxy.ServiceParameters + +from DPGAnalysis.MuonTools.muNtupleDTTPGPhiFiller_cfi import muNtupleDTTPGPhiFiller + +muNtupleBmtfInFiller = muNtupleDTTPGPhiFiller.clone() +muNtupleTwinMuxInFiller = muNtupleDTTPGPhiFiller.clone(tag = 'TM_IN', label = 'ltTwinMuxIn', dtTpTag = cms.InputTag('twinMuxStage2Digis','PhIn')) +muNtupleTwinMuxOutFiller = muNtupleDTTPGPhiFiller.clone(tag = 'TM_OUT', label = 'ltTwinMuxOut', dtTpTag = cms.InputTag('twinMuxStage2Digis','PhOut')) + +from DPGAnalysis.MuonTools.muNtupleDTTPGThetaFiller_cfi import muNtupleDTTPGThetaFiller + +muNtupleBmtfInThFiller = muNtupleDTTPGThetaFiller.clone() +muNtupleTwinMuxInThFiller = muNtupleDTTPGPhiFiller.clone(tag = 'TM_IN', label = 'ltTwinMuxInTh', dtTpTag = cms.InputTag('twinMuxStage2Digis','ThIn')) + +muNtupleProducer = cms.Sequence(muNtupleDTDigiFiller + + muNtupleRPCDigiFiller + + muNtupleGEMDigiFiller + + muNtupleRPCRecHitFiller + + muNtupleGEMRecHitFiller + + muNtupleDTSegmentFiller + + muNtupleGEMSegmentFiller + + muNtupleTwinMuxInFiller + + muNtupleTwinMuxOutFiller + + muNtupleBmtfInFiller + + muNtupleTwinMuxInThFiller + + muNtupleBmtfInThFiller + + muNtupleMuonFiller + + muNtupleGEMMuonFiller + ) \ No newline at end of file diff --git a/DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.cc b/DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.cc new file mode 100644 index 0000000000000..bc3104bc87d12 --- /dev/null +++ b/DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.cc @@ -0,0 +1,22 @@ +/** \class MuNtupleBaseFiller MuNtupleBaseFiller.cc DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.cc + * + * Helper class defining the generic interface of a filler + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h" + +MuNtupleBaseFiller::MuNtupleBaseFiller(const edm::ParameterSet &config) + : m_label{config.getParameter("label")} {} + +void MuNtupleBaseFiller::beginRun(const edm::Run &run, const edm::EventSetup &environment) { + getFromES(run, environment); +} + +void MuNtupleBaseFiller::produce(edm::Event &ev, const edm::EventSetup &environment) { + getFromES(environment); + fill(ev); +} diff --git a/DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h b/DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h new file mode 100644 index 0000000000000..8ab1a68e2c8a9 --- /dev/null +++ b/DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h @@ -0,0 +1,72 @@ +#ifndef MuNtuple_MuNtupleBaseFiller_h +#define MuNtuple_MuNtupleBaseFiller_h + +/** \class MuNtupleBaseFiller MuNtupleBaseFiller.h DPGAnalysis/MuonTools/src/MuNtupleBaseFiller.h + * + * Helper class defining the generic interface of a filler + * + * \author C. Battilana (INFN BO) + * + * + */ + +#include "DPGAnalysis/MuonTools/src/MuNtupleUtils.h" + +#include "FWCore/Framework/interface/stream/EDProducer.h" +#include "FWCore/Framework/interface/ConsumesCollector.h" + +#include "DataFormats/NanoAOD/interface/FlatTable.h" + +#include + +class MuNtupleBaseFiller : public edm::stream::EDProducer<> { +public: + /// Constructor + explicit MuNtupleBaseFiller(const edm::ParameterSet &); + + /// Configure event setup for each run + void beginRun(const edm::Run &run, const edm::EventSetup &config) final; + + /// Fill ntuples event by event + void produce(edm::Event &, const edm::EventSetup &) final; + + /// Empty, needed by interface + void endRun(const edm::Run &, const edm::EventSetup &) final {} + +protected: + /// The label name of the filler + std::string m_label; + + /// Get info from the ES by run + virtual void getFromES(const edm::Run &run, const edm::EventSetup &environment) {} + + /// Get info from the ES for a given event + virtual void getFromES(const edm::EventSetup &environment) {} + + /// Fill ntuple + virtual void fill(edm::Event &ev) = 0; + + /// Definition of default values for int variables + static constexpr int DEFAULT_INT_VAL = -999; + + /// Definition of default values for int8 variables + static constexpr int8_t DEFAULT_INT8_VAL = -99; + + /// Definition of default values for positive int variables + static constexpr int DEFAULT_INT_VAL_POS = -1; + + /// Definition of default values for float variables + static constexpr double DEFAULT_DOUBLE_VAL = -999.; + + /// Definition of default values for positive float variables + static constexpr double DEFAULT_DOUBLE_VAL_POS = -1.; + + template + void addColumn(std::unique_ptr &table, + const std::string name, + const std::vector &vec, + const std::string descr) { + table->template addColumn>(name, vec, descr); + } +}; +#endif diff --git a/DPGAnalysis/MuonTools/src/MuNtupleUtils.cc b/DPGAnalysis/MuonTools/src/MuNtupleUtils.cc new file mode 100644 index 0000000000000..ca536949c4e03 --- /dev/null +++ b/DPGAnalysis/MuonTools/src/MuNtupleUtils.cc @@ -0,0 +1,101 @@ +/* + * \file MuNtupleUtils.cc + * + * \author C. Battilana - INFN (BO) + * \author L. Giuducci - INFN (BO) +*/ + +#include "DPGAnalysis/MuonTools/src/MuNtupleUtils.h" + +#include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambPhDigi.h" +#include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTPhDigi.h" + +#include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h" + +// CB +// formulas to be re-checked +// can use template is_same, static_assert, enable_if ... + +nano_mu::DTTrigGeomUtils::DTTrigGeomUtils(edm::ConsumesCollector&& collector, bool dirInDeg) + : m_dtGeom{std::move(collector), "idealForDigi"} {} + +nano_mu::DTTrigGeomUtils::chambCoord nano_mu::DTTrigGeomUtils::trigToReco(const L1MuDTChambPhDigi* trig) { + auto wh{trig->whNum()}; + auto sec{trig->scNum() + 1}; + auto st{trig->stNum()}; + auto phi{trig->phi()}; + auto phib{trig->phiB()}; + + auto recoChamb = [&]() { + if (st != 4) { + return DTChamberId(wh, st, sec); + } + int reco_sec{(sec == 4 && phi > 0) ? 13 : (sec == 10 && phi > 0) ? 14 : sec}; + return DTChamberId(wh, st, reco_sec); + }; + + auto gpos{m_dtGeom->chamber(recoChamb())->position()}; + auto r{gpos.perp()}; + + auto delta_phi{gpos.phi() - (sec - 1) * Geom::pi() / 6}; + + // zcn is in local coordinates -> z invreases approching to vertex + // LG: zcn offset was removed <- CB do we need to fix this? + float x = r * tan((phi - (phi < 0 ? 1 : 0)) / PH1_PHI_R) * cos(delta_phi) - r * sin(delta_phi); + float dir = (phib / PH1_PHIB_R + phi / PH1_PHI_R); + + // change sign in case of positive wheels + if (hasPosRF(wh, sec)) { + x = -x; + } else { + dir = -dir; + } + + return {x, dir}; +} + +nano_mu::DTTrigGeomUtils::chambCoord nano_mu::DTTrigGeomUtils::trigToReco(const L1Phase2MuDTPhDigi* trig) { + auto wh{trig->whNum()}; + auto sec{trig->scNum() + 1}; + auto st{trig->stNum()}; + auto phi{trig->phi()}; + auto phib{trig->phiBend()}; + auto quality{trig->quality()}; + auto sl{trig->slNum()}; + + auto recoChamb = [&]() { + if (st != 4) { + return DTChamberId(wh, st, sec); + } + int reco_sec{(sec == 4 && phi > 0) ? 13 : (sec == 10 && phi > 0) ? 14 : sec}; + return DTChamberId(wh, st, reco_sec); + }; + + auto gpos{m_dtGeom->chamber(recoChamb())->position()}; + auto r{gpos.perp()}; + + auto delta_phi{gpos.phi() - (sec - 1) * Geom::pi() / 6}; + + // CB to be potentially updated based on Silvia's results + double zRF = 0; + if (quality >= 6 && quality != 7) + zRF = m_zcn[st - 1]; + if ((quality < 6 || quality == 7) && sl == 1) + zRF = m_zsl1[st - 1]; + if ((quality < 6 || quality == 7) && sl == 3) + zRF = m_zsl3[st - 1]; + + // zcn is in local coordinates -> z invreases approching to vertex + // LG: zcn offset was removed <- CB Mist confirm it is tryly accurate? + float x = r * tan((phi - (phi < 0 ? 1 : 0)) / PH1_PHI_R) * (cos(delta_phi) - zRF) - r * sin(delta_phi); + float dir = (phib / PH2_PHIB_R + phi / PH2_PHI_R); + + // change sign in case of positive wheels + if (hasPosRF(wh, sec)) { + x = -x; + } else { + dir = -dir; + } + + return {x, dir}; +} diff --git a/DPGAnalysis/MuonTools/src/MuNtupleUtils.h b/DPGAnalysis/MuonTools/src/MuNtupleUtils.h new file mode 100644 index 0000000000000..71116cf55956d --- /dev/null +++ b/DPGAnalysis/MuonTools/src/MuNtupleUtils.h @@ -0,0 +1,135 @@ +#ifndef MuNtuple_MuNtupleUtils_h +#define MuNtuple_MuNtupleUtils_h + +/** \class MuNtupleUtils MuNtupleUtils.h MuDPGAnalysis/MuNtuples/src/MuNtupleUtils.h + * + * A set of helper classes class to handle : + * - Handing of InputTags and tokens + * - DB information from edm::EventSetup + * - Conversion between L1T trigger primitive coordinates and CMS global ones + * + * \author C. Battilana - INFN (BO) + * \author L. Giuducci - INFN (BO) + * + * + */ + +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Utilities/interface/InputTag.h" +#include "FWCore/Framework/interface/ESHandle.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/Framework/interface/ConsumesCollector.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" + +#include "Geometry/DTGeometry/interface/DTGeometry.h" +#include "Geometry/Records/interface/MuonGeometryRecord.h" + +#include +#include + +class L1MuDTChambPhDigi; +class L1Phase2MuDTPhDigi; + +namespace nano_mu { + + template + class EDTokenHandle { + public: + /// Constructor + EDTokenHandle(const edm::ParameterSet& config, edm::ConsumesCollector&& collector, std::string name) + : m_name{name}, m_inputTag{config.getParameter(name)} { + if (m_inputTag.label() != "none") { + m_token = collector.template consumes(m_inputTag); + } + } + + /// Conditional getter + /// checks whether a token is valid and if + /// retireving the data collection succeded + auto conditionalGet(const edm::Event& ev) const { + edm::Handle collection; + + if (!m_token.isUninitialized() && !ev.getByToken(m_token, collection)) + edm::LogError("") << "[EDTokenHandle]::conditionalGet: " << m_inputTag.label() + << " collection does not exist !!!"; + + return collection; + } + + private: + std::string m_name; + edm::InputTag m_inputTag; + edm::EDGetTokenT m_token; + }; + + template + class ESTokenHandle { + public: + /// Constructor + ESTokenHandle(edm::ConsumesCollector&& collector, const std::string& label = "") + : m_token{collector.template esConsumes(edm::ESInputTag{"", label})} {} + + ///Get Handle from ES + void getFromES(const edm::EventSetup& environment) { m_handle = environment.getHandle(m_token); } + + /// Check validity + bool isValid() { return m_handle.isValid(); } + + /// Return handle + T const* operator->() { return m_handle.product(); } + + private: + edm::ESGetToken m_token; + edm::ESHandle m_handle; + }; + + class DTTrigGeomUtils { + public: + struct chambCoord { + double pos{}; + double dir{}; + }; + + /// Constructor + DTTrigGeomUtils(edm::ConsumesCollector&& collector, bool dirInDeg = true); + + /// Return local position and direction in chamber RF - legacy + chambCoord trigToReco(const L1MuDTChambPhDigi* trig); + + /// Return local position and direction in chamber RF - analytical method + chambCoord trigToReco(const L1Phase2MuDTPhDigi* trig); + + /// Checks id the chamber has positive RF; + bool hasPosRF(int wh, int sec) { return wh > 0 || (wh == 0 && sec % 4 > 1); }; + + /// Update EventSetup information + void getFromES(const edm::Run& run, const edm::EventSetup& environment) { + m_dtGeom.getFromES(environment); + for (int i_st = 0; i_st != 4; ++i_st) { + const DTChamberId chId(-2, i_st + 1, 4); + const DTChamber* chamb = m_dtGeom->chamber(chId); + const DTSuperLayer* sl1 = chamb->superLayer(DTSuperLayerId(chId, 1)); + const DTSuperLayer* sl3 = chamb->superLayer(DTSuperLayerId(chId, 3)); + m_zsl1[i_st] = chamb->surface().toLocal(sl1->position()).z(); + m_zsl3[i_st] = chamb->surface().toLocal(sl3->position()).z(); + m_zcn[i_st] = 0.5 * (m_zsl1[i_st] + m_zsl3[i_st]); + } + }; + + private: + ESTokenHandle m_dtGeom; + + std::array m_zcn; + std::array m_zsl1; + std::array m_zsl3; + + static constexpr double PH1_PHI_R = 4096.; + static constexpr double PH1_PHIB_R = 512.; + + static constexpr double PH2_PHI_R = 65536. / 0.8; + static constexpr double PH2_PHIB_R = 2048. / 1.4; + }; + +} // namespace nano_mu + +#endif diff --git a/DPGAnalysis/MuonTools/test/muDpgNtuples_cfg.py b/DPGAnalysis/MuonTools/test/muDpgNtuples_cfg.py new file mode 100644 index 0000000000000..4c98172e96e55 --- /dev/null +++ b/DPGAnalysis/MuonTools/test/muDpgNtuples_cfg.py @@ -0,0 +1,118 @@ +import FWCore.ParameterSet.Config as cms +import FWCore.ParameterSet.VarParsing as VarParsing +from Configuration.StandardSequences.Eras import eras +from Configuration.AlCa.GlobalTag import GlobalTag + +import subprocess +import sys + +options = VarParsing.VarParsing() + +options.register('globalTag', + '124X_mcRun3_2021_realistic_v1', #default value + VarParsing.VarParsing.multiplicity.singleton, + VarParsing.VarParsing.varType.string, + "Global Tag") + +options.register('nEvents', + 1000, #default value + VarParsing.VarParsing.multiplicity.singleton, + VarParsing.VarParsing.varType.int, + "Maximum number of processed events") + +options.register('inputFolder', + '/RelValZMM_14/CMSSW_12_4_0_pre4-PU_124X_mcRun3_2021_realistic_v1-v1/GEN-SIM-RECO', #default value + VarParsing.VarParsing.multiplicity.singleton, + VarParsing.VarParsing.varType.string, + "EOS folder with input files") + +options.register('secondaryInputFolder', + '/RelValZMM_14/CMSSW_12_4_0_pre4-PU_124X_mcRun3_2021_realistic_v1-v1/GEN-SIM-DIGI-RAW', #default value + VarParsing.VarParsing.multiplicity.singleton, + VarParsing.VarParsing.varType.string, + "EOS folder with input files for secondary files") + +options.register('ntupleName', + './MuDPGNtuple_nanoAOD_ZMMPU.root', #default value + VarParsing.VarParsing.multiplicity.singleton, + VarParsing.VarParsing.varType.string, + "Folder and name ame for output ntuple") + +options.register('runOnMC', + True, #default value + VarParsing.VarParsing.multiplicity.singleton, + VarParsing.VarParsing.varType.bool, + "Apply customizations to run on MC") + +options.parseArguments() + +process = cms.Process("MUNTUPLES",eras.Run3) + +process.load('FWCore.MessageService.MessageLogger_cfi') + +process.options = cms.untracked.PSet( wantSummary = cms.untracked.bool(True), + numberOfThreads = cms.untracked.uint32(4)) +process.MessageLogger.cerr.FwkReport.reportEvery = 100 +process.maxEvents = cms.untracked.PSet(input = cms.untracked.int32(options.nEvents)) + +process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff') +process.GlobalTag.globaltag = cms.string(options.globalTag) + +process.source = cms.Source("PoolSource", + fileNames = cms.untracked.vstring(), + secondaryFileNames = cms.untracked.vstring() +) + +files = subprocess.check_output(['dasgoclient', '--query', + f'file dataset={options.inputFolder}']) +process.source.fileNames = [f.decode("utf-8") for f in files.split() if b".root" in f] + +print(process.source.fileNames) + +if options.secondaryInputFolder != "" : + files = subprocess.check_output(['dasgoclient', '--query', + f'file dataset={options.secondaryInputFolder}']) + process.source.secondaryFileNames = [f.decode("utf-8") for f in files.split() if b".root" in f] + + print(process.source.secondaryFileNames) + +process.load('Configuration/StandardSequences/GeometryRecoDB_cff') +process.load("Configuration.StandardSequences.MagneticField_cff") + +process.load("TrackingTools/TransientTrack/TransientTrackBuilder_cfi") +process.load('TrackPropagation.SteppingHelixPropagator.SteppingHelixPropagatorAny_cfi') +process.load('TrackPropagation.SteppingHelixPropagator.SteppingHelixPropagatorAlong_cfi') +process.load('TrackPropagation.SteppingHelixPropagator.SteppingHelixPropagatorOpposite_cfi') + +process.load('Configuration.StandardSequences.RawToDigi_Data_cff') +process.load('DPGAnalysis.MuonTools.muNtupleProducer_cff') + +import EventFilter.RPCRawToDigi.rpcUnpacker_cfi +process.muonRPCDigis = EventFilter.RPCRawToDigi.rpcUnpacker_cfi.rpcunpacker.clone() + +process.nanoMuDPGPath = cms.Path(process.muonDTDigis + + process.muonRPCDigis + + process.muonGEMDigis + + process.twinMuxStage2Digis + + process.bmtfDigis + + process.muNtupleProducer) + +process.load("PhysicsTools.NanoAOD.NanoAODEDMEventContent_cff") + +process.out = cms.OutputModule("NanoAODOutputModule", + fileName = cms.untracked.string(options.ntupleName), + outputCommands = process.NANOAODEventContent.outputCommands \ + + ["keep nanoaodFlatTable_*_*_*", + "drop edmTriggerResults_*_*_*"], + SelectEvents = cms.untracked.PSet( + SelectEvents=cms.vstring("nanoMuDPGPath") + ) +) + +process.end = cms.EndPath(process.out) + +process.schedule = cms.Schedule(process.nanoMuDPGPath, process.end) + +if options.runOnMC : + from DPGAnalysis.MuonTools.customiseMuNtuples_cff import customiseForRunningOnMC + customiseForRunningOnMC(process,"p")