Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-07-22 22:47:28

0001 #include "FWCore/Framework/interface/MakerMacros.h"
0002 
0003 #include "DataFormats/Common/interface/DetSetVector.h"
0004 #include "DataFormats/Common/interface/DetSet.h"
0005 #include "FWCore/Framework/interface/stream/EDProducer.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/EventSetup.h"
0008 #include "DataFormats/Common/interface/Handle.h"
0009 #include "FWCore/Framework/interface/ESHandle.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 
0013 #include "FWCore/Framework/interface/ESWatcher.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016 #include "FWCore/Utilities/interface/ESGetToken.h"
0017 
0018 #include "Geometry/Records/interface/VeryForwardRealGeometryRecord.h"
0019 
0020 #include "DataFormats/CTPPSReco/interface/CTPPSPixelRecHit.h"
0021 #include "DataFormats/CTPPSReco/interface/CTPPSPixelLocalTrack.h"
0022 #include "DataFormats/DetId/interface/DetId.h"
0023 
0024 #include "RecoPPS/Local/interface/RPixDetPatternFinder.h"
0025 #include "RecoPPS/Local/interface/RPixDetTrackFinder.h"
0026 
0027 #include <memory>
0028 
0029 #include <string>
0030 #include <vector>
0031 
0032 #include "FWCore/Framework/interface/Frameworkfwd.h"
0033 
0034 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0035 
0036 #include "DataFormats/CTPPSDetId/interface/CTPPSPixelDetId.h"
0037 
0038 #include "RecoPPS/Local/interface/RPixRoadFinder.h"
0039 #include "RecoPPS/Local/interface/RPixPlaneCombinatoryTracking.h"
0040 
0041 #include "CondFormats/PPSObjects/interface/CTPPSPixelAnalysisMask.h"
0042 #include "CondFormats/DataRecord/interface/CTPPSPixelAnalysisMaskRcd.h"
0043 
0044 namespace {
0045   constexpr int rocMask = 0xE000;
0046   constexpr int rocOffset = 13;
0047   constexpr int rocSizeInPixels = 4160;
0048 }  // namespace
0049 
0050 class CTPPSPixelLocalTrackProducer : public edm::stream::EDProducer<> {
0051 public:
0052   explicit CTPPSPixelLocalTrackProducer(const edm::ParameterSet &parameterSet);
0053 
0054   ~CTPPSPixelLocalTrackProducer() override;
0055 
0056   void produce(edm::Event &, const edm::EventSetup &) override;
0057   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0058 
0059 private:
0060   const int verbosity_;
0061   const int maxHitPerPlane_;
0062   const int maxHitPerRomanPot_;
0063   const int maxTrackPerRomanPot_;
0064   const int maxTrackPerPattern_;
0065 
0066   const edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelRecHit>> tokenCTPPSPixelRecHit_;
0067   const edm::ESGetToken<CTPPSGeometry, VeryForwardRealGeometryRecord> tokenCTPPSGeometry_;
0068   edm::ESWatcher<VeryForwardRealGeometryRecord> geometryWatcher_;
0069 
0070   const edm::ESGetToken<CTPPSPixelAnalysisMask, CTPPSPixelAnalysisMaskRcd> tokenCTPPSPixelAnalysisMask_;
0071 
0072   const uint32_t numberOfPlanesPerPot_;
0073 
0074   std::unique_ptr<RPixDetPatternFinder> patternFinder_;
0075   std::unique_ptr<RPixDetTrackFinder> trackFinder_;
0076 };
0077 
0078 //------------------------------------------------------------------------------------------------//
0079 
0080 CTPPSPixelLocalTrackProducer::CTPPSPixelLocalTrackProducer(const edm::ParameterSet &parameterSet)
0081     : verbosity_(parameterSet.getUntrackedParameter<int>("verbosity")),
0082       maxHitPerPlane_(parameterSet.getParameter<int>("maxHitPerPlane")),
0083       maxHitPerRomanPot_(parameterSet.getParameter<int>("maxHitPerRomanPot")),
0084       maxTrackPerRomanPot_(parameterSet.getParameter<int>("maxTrackPerRomanPot")),
0085       maxTrackPerPattern_(parameterSet.getParameter<int>("maxTrackPerPattern")),
0086       tokenCTPPSPixelRecHit_(
0087           consumes<edm::DetSetVector<CTPPSPixelRecHit>>(parameterSet.getParameter<edm::InputTag>("tag"))),
0088       tokenCTPPSGeometry_(esConsumes<CTPPSGeometry, VeryForwardRealGeometryRecord>()),
0089       tokenCTPPSPixelAnalysisMask_(esConsumes<CTPPSPixelAnalysisMask, CTPPSPixelAnalysisMaskRcd>()),
0090       numberOfPlanesPerPot_(parameterSet.getParameter<int>("numberOfPlanesPerPot")) {
0091   std::string patternFinderAlgorithm = parameterSet.getParameter<std::string>("patternFinderAlgorithm");
0092   std::string trackFitterAlgorithm = parameterSet.getParameter<std::string>("trackFinderAlgorithm");
0093 
0094   // pattern algorithm selector
0095   if (patternFinderAlgorithm == "RPixRoadFinder") {
0096     patternFinder_ = std::make_unique<RPixRoadFinder>(parameterSet);
0097   } else {
0098     throw cms::Exception("CTPPSPixelLocalTrackProducer")
0099         << "Pattern finder algorithm" << patternFinderAlgorithm << " does not exist";
0100   }
0101 
0102   std::vector<uint32_t> listOfAllPlanes;
0103   for (uint32_t i = 0; i < numberOfPlanesPerPot_; ++i) {
0104     listOfAllPlanes.push_back(i);
0105   }
0106 
0107   //tracking algorithm selector
0108   if (trackFitterAlgorithm == "RPixPlaneCombinatoryTracking") {
0109     trackFinder_ = std::make_unique<RPixPlaneCombinatoryTracking>(parameterSet);
0110   } else {
0111     throw cms::Exception("CTPPSPixelLocalTrackProducer")
0112         << "Tracking fitter algorithm" << trackFitterAlgorithm << " does not exist";
0113   }
0114   trackFinder_->setListOfPlanes(listOfAllPlanes);
0115   trackFinder_->initialize();
0116   produces<edm::DetSetVector<CTPPSPixelLocalTrack>>();
0117 }
0118 
0119 //------------------------------------------------------------------------------------------------//
0120 
0121 CTPPSPixelLocalTrackProducer::~CTPPSPixelLocalTrackProducer() {}
0122 
0123 //------------------------------------------------------------------------------------------------//
0124 
0125 void CTPPSPixelLocalTrackProducer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0126   edm::ParameterSetDescription desc;
0127 
0128   desc.add<edm::InputTag>("tag", edm::InputTag("ctppsPixelRecHits"))
0129       ->setComment("inputTag of the RecHits input for the tracking algorithm");
0130   desc.add<std::string>("patternFinderAlgorithm", "RPixRoadFinder")->setComment("algorithm type for pattern finder");
0131   desc.add<std::string>("trackFinderAlgorithm", "RPixPlaneCombinatoryTracking")
0132       ->setComment("algorithm type for track finder");
0133   desc.add<uint>("trackMinNumberOfPoints", 3)->setComment("minimum number of planes to produce a track");
0134   desc.addUntracked<int>("verbosity", 0)->setComment("verbosity for track producer");
0135   desc.add<double>("maximumChi2OverNDF", 5.)->setComment("maximum Chi2OverNDF for accepting the track");
0136   desc.add<double>("maximumXLocalDistanceFromTrack", 0.2)
0137       ->setComment("maximum x distance in mm to associate a point not used for fit to the track");
0138   desc.add<double>("maximumYLocalDistanceFromTrack", 0.3)
0139       ->setComment("maximum y distance in mm to associate a point not used for fit to the track");
0140   desc.add<int>("maxHitPerPlane", 20)
0141       ->setComment("maximum hits per plane, events with higher number will not be fitted");
0142   desc.add<int>("maxHitPerRomanPot", 60)
0143       ->setComment("maximum hits per roman pot, events with higher number will not be fitted");
0144   desc.add<int>("maxTrackPerRomanPot", 10)
0145       ->setComment("maximum tracks per roman pot, events with higher track will not be saved");
0146   desc.add<int>("maxTrackPerPattern", 5)
0147       ->setComment("maximum tracks per pattern, events with higher track will not be saved");
0148   desc.add<int>("numberOfPlanesPerPot", 6)->setComment("number of planes per pot");
0149   desc.add<double>("roadRadius", 1.0)->setComment("radius of pattern search window");
0150   desc.add<int>("minRoadSize", 3)->setComment("minimum number of points in a pattern");
0151   desc.add<int>("maxRoadSize", 20)->setComment("maximum number of points in a pattern");
0152   //parameter for bad pot reconstruction patch 45-220-fr 2022
0153   desc.add<double>("roadRadiusBadPot", 0.5)->setComment("radius of pattern search window for bad Pot");
0154 
0155   descriptions.add("ctppsPixelLocalTracks", desc);
0156 }
0157 
0158 //------------------------------------------------------------------------------------------------//
0159 
0160 void CTPPSPixelLocalTrackProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0161   // Step A: get inputs
0162 
0163   edm::Handle<edm::DetSetVector<CTPPSPixelRecHit>> recHits;
0164   iEvent.getByToken(tokenCTPPSPixelRecHit_, recHits);
0165   edm::DetSetVector<CTPPSPixelRecHit> recHitVector(*recHits);
0166 
0167   // get geometry
0168   edm::ESHandle<CTPPSGeometry> geometryHandler = iSetup.getHandle(tokenCTPPSGeometry_);
0169   const CTPPSGeometry &geometry = *geometryHandler;
0170   geometryWatcher_.check(iSetup);
0171 
0172   // get mask
0173   bool isBadPot_45_220 = false;
0174   if (!recHits->empty()) {
0175     const auto &mask = iSetup.getData(tokenCTPPSPixelAnalysisMask_);
0176 
0177     // Read Mask checking if 45-220-far is masked as bad and needs special treatment
0178     std::map<uint32_t, CTPPSPixelROCAnalysisMask> const &maschera = mask.analysisMask;
0179 
0180     bool mask_45_220[6][6] = {{false}};
0181     for (auto const &det : maschera) {
0182       CTPPSPixelDetId detId(det.first);
0183       unsigned int rocNum = (det.first & rocMask) >> rocOffset;
0184       if (rocNum > 5 || detId.plane() > 5)
0185         throw cms::Exception("InvalidRocOrPlaneNumber") << "roc number from mask: " << rocNum;
0186 
0187       if (detId.arm() == 0 && detId.station() == 2 && detId.rp() == 3) {  // pot 45-220-far
0188         if (det.second.maskedPixels.size() == rocSizeInPixels) {          // roc fully masked
0189           mask_45_220[detId.plane()][rocNum] = true;
0190         }
0191       }
0192     }
0193 
0194     // search for specific pattern that requires special reconstruction (isBadPot)
0195     isBadPot_45_220 = mask_45_220[1][4] && mask_45_220[1][5] && mask_45_220[2][4] && mask_45_220[2][5] &&
0196                       mask_45_220[3][4] && mask_45_220[3][5] && mask_45_220[4][4] && mask_45_220[4][5];
0197   }
0198   std::vector<CTPPSPixelDetId> listOfPotWithHighOccupancyPlanes;
0199   std::map<CTPPSPixelDetId, uint32_t> mapHitPerPot;
0200 
0201   for (const auto &recHitSet : recHitVector) {
0202     if (verbosity_ > 2)
0203       edm::LogInfo("CTPPSPixelLocalTrackProducer")
0204           << "Hits found in plane = " << recHitSet.detId() << " number = " << recHitSet.size();
0205     CTPPSPixelDetId tmpRomanPotId = CTPPSPixelDetId(recHitSet.detId()).rpId();
0206     uint32_t hitOnPlane = recHitSet.size();
0207 
0208     //Get the number of hits per pot
0209     if (mapHitPerPot.find(tmpRomanPotId) == mapHitPerPot.end()) {
0210       mapHitPerPot[tmpRomanPotId] = hitOnPlane;
0211     } else
0212       mapHitPerPot[tmpRomanPotId] += hitOnPlane;
0213 
0214     //check is the plane occupancy is too high and save the corresponding pot
0215     if (maxHitPerPlane_ >= 0 && hitOnPlane > (uint32_t)maxHitPerPlane_) {
0216       if (verbosity_ > 2)
0217         edm::LogInfo("CTPPSPixelLocalTrackProducer")
0218             << " ---> To many hits in the plane, pot will be excluded from tracking cleared";
0219       listOfPotWithHighOccupancyPlanes.push_back(tmpRomanPotId);
0220     }
0221   }
0222 
0223   //remove rechit for pot with too many hits or containing planes with too many hits
0224   for (const auto &recHitSet : recHitVector) {
0225     const auto tmpDetectorId = CTPPSPixelDetId(recHitSet.detId());
0226     const auto tmpRomanPotId = tmpDetectorId.rpId();
0227 
0228     if ((maxHitPerRomanPot_ >= 0 && mapHitPerPot[tmpRomanPotId] > (uint32_t)maxHitPerRomanPot_) ||
0229         find(listOfPotWithHighOccupancyPlanes.begin(), listOfPotWithHighOccupancyPlanes.end(), tmpRomanPotId) !=
0230             listOfPotWithHighOccupancyPlanes.end()) {
0231       edm::DetSet<CTPPSPixelRecHit> &tmpDetSet = recHitVector[tmpDetectorId];
0232       tmpDetSet.clear();
0233     }
0234   }
0235 
0236   edm::DetSetVector<CTPPSPixelLocalTrack> foundTracks;
0237 
0238   // Pattern finder
0239 
0240   patternFinder_->clear();
0241   patternFinder_->setHits(&recHitVector);
0242   patternFinder_->setGeometry(&geometry);
0243   patternFinder_->findPattern(isBadPot_45_220);
0244   std::vector<RPixDetPatternFinder::Road> patternVector = patternFinder_->getPatterns();
0245 
0246   //loop on all the patterns
0247   int numberOfTracks = 0;
0248 
0249   for (const auto &pattern : patternVector) {
0250     CTPPSPixelDetId firstHitDetId = CTPPSPixelDetId(pattern.at(0).detId);
0251     CTPPSPixelDetId romanPotId = firstHitDetId.rpId();
0252 
0253     std::map<CTPPSPixelDetId, std::vector<RPixDetPatternFinder::PointInPlane>>
0254         hitOnPlaneMap;  //hit of the pattern organized by plane
0255 
0256     //loop on all the hits of the pattern
0257     for (const auto &hit : pattern) {
0258       CTPPSPixelDetId hitDetId = CTPPSPixelDetId(hit.detId);
0259       CTPPSPixelDetId tmpRomanPotId = hitDetId.rpId();
0260 
0261       if (tmpRomanPotId != romanPotId) {  //check that the hits belong to the same tracking station
0262         throw cms::Exception("CTPPSPixelLocalTrackProducer")
0263             << "Hits in the pattern must belong to the same tracking station";
0264       }
0265 
0266       if (hitOnPlaneMap.find(hitDetId) ==
0267           hitOnPlaneMap.end()) {  //add the plane key in case it is the first hit of that plane
0268         std::vector<RPixDetPatternFinder::PointInPlane> hitOnPlane;
0269         hitOnPlane.push_back(hit);
0270         hitOnPlaneMap[hitDetId] = hitOnPlane;
0271       } else
0272         hitOnPlaneMap[hitDetId].push_back(hit);  //add the hit to an existing plane key
0273     }
0274 
0275     trackFinder_->clear();
0276     trackFinder_->setRomanPotId(romanPotId);
0277     trackFinder_->setHits(&hitOnPlaneMap);
0278     trackFinder_->setGeometry(&geometry);
0279     trackFinder_->setZ0(geometry.rpTranslation(romanPotId).z());
0280     trackFinder_->findTracks(iEvent.getRun().id().run());
0281     std::vector<CTPPSPixelLocalTrack> tmpTracksVector = trackFinder_->getLocalTracks();
0282 
0283     if (verbosity_ > 2)
0284       edm::LogInfo("CTPPSPixelLocalTrackProducer") << "tmpTracksVector = " << tmpTracksVector.size();
0285     if (maxTrackPerPattern_ >= 0 && tmpTracksVector.size() > (uint32_t)maxTrackPerPattern_) {
0286       if (verbosity_ > 2)
0287         edm::LogInfo("CTPPSPixelLocalTrackProducer") << " ---> To many tracks in the pattern, cleared";
0288       continue;
0289     }
0290 
0291     for (const auto &track : tmpTracksVector) {
0292       ++numberOfTracks;
0293       edm::DetSet<CTPPSPixelLocalTrack> &tmpDetSet = foundTracks.find_or_insert(romanPotId);
0294       tmpDetSet.push_back(track);
0295     }
0296   }
0297 
0298   if (verbosity_ > 1)
0299     edm::LogInfo("CTPPSPixelLocalTrackProducer") << "Number of tracks will be saved = " << numberOfTracks;
0300 
0301   for (const auto &track : foundTracks) {
0302     if (verbosity_ > 1)
0303       edm::LogInfo("CTPPSPixelLocalTrackProducer")
0304           << "Track found in detId = " << track.detId() << " number = " << track.size();
0305     if (maxTrackPerRomanPot_ >= 0 && track.size() > (uint32_t)maxTrackPerRomanPot_) {
0306       if (verbosity_ > 1)
0307         edm::LogInfo("CTPPSPixelLocalTrackProducer") << " ---> Too many tracks in the pot, cleared";
0308       CTPPSPixelDetId tmpRomanPotId = CTPPSPixelDetId(track.detId());
0309       edm::DetSet<CTPPSPixelLocalTrack> &tmpDetSet = foundTracks[tmpRomanPotId];
0310       tmpDetSet.clear();
0311     }
0312   }
0313 
0314   iEvent.put(std::make_unique<edm::DetSetVector<CTPPSPixelLocalTrack>>(foundTracks));
0315 
0316   return;
0317 }
0318 
0319 DEFINE_FWK_MODULE(CTPPSPixelLocalTrackProducer);