Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-05-26 22:39:47

0001 #include "FWCore/Framework/interface/global/EDProducer.h"
0002 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0003 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/EventSetup.h"
0006 
0007 #include "DataFormats/Common/interface/Handle.h"
0008 #include "DataFormats/FTLRecHit/interface/FTLClusterCollections.h"
0009 
0010 #include "Geometry/Records/interface/MTDDigiGeometryRecord.h"
0011 #include "Geometry/CommonTopologies/interface/Topology.h"
0012 #include "Geometry/MTDGeometryBuilder/interface/MTDGeometry.h"
0013 #include "DataFormats/TrackerRecHit2D/interface/MTDTrackingRecHit.h"
0014 #include "DataFormats/GeometryCommonDetAlgo/interface/MeasurementError.h"
0015 #include "DataFormats/GeometryCommonDetAlgo/interface/MeasurementPoint.h"
0016 
0017 #include "RecoLocalFastTime/Records/interface/MTDCPERecord.h"
0018 #include "RecoLocalFastTime/FTLClusterizer/interface/MTDClusterParameterEstimator.h"
0019 
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021 
0022 class MTDTrackingRecHitProducer : public edm::global::EDProducer<> {
0023 public:
0024   explicit MTDTrackingRecHitProducer(const edm::ParameterSet& ps);
0025   ~MTDTrackingRecHitProducer() override = default;
0026   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0027 
0028   void produce(edm::StreamID, edm::Event& evt, const edm::EventSetup& es) const override;
0029 
0030 private:
0031   const edm::EDGetTokenT<FTLClusterCollection> ftlbClusters_;  // collection of barrel digis
0032   const edm::EDGetTokenT<FTLClusterCollection> ftleClusters_;  // collection of endcap digis
0033 
0034   const edm::ESGetToken<MTDGeometry, MTDDigiGeometryRecord> mtdgeoToken_;
0035   const edm::ESGetToken<MTDClusterParameterEstimator, MTDCPERecord> cpeToken_;
0036 };
0037 
0038 MTDTrackingRecHitProducer::MTDTrackingRecHitProducer(const edm::ParameterSet& ps)
0039     : ftlbClusters_(consumes<FTLClusterCollection>(ps.getParameter<edm::InputTag>("barrelClusters"))),
0040       ftleClusters_(consumes<FTLClusterCollection>(ps.getParameter<edm::InputTag>("endcapClusters"))),
0041       mtdgeoToken_(esConsumes<MTDGeometry, MTDDigiGeometryRecord>()),
0042       cpeToken_(esConsumes<MTDClusterParameterEstimator, MTDCPERecord>(edm::ESInputTag("", "MTDCPEBase"))) {
0043   produces<MTDTrackingDetSetVector>();
0044 }
0045 
0046 // Configuration descriptions
0047 void MTDTrackingRecHitProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0048   edm::ParameterSetDescription desc;
0049   desc.add<edm::InputTag>("barrelClusters", edm::InputTag("mtdClusters:FTLBarrel"));
0050   desc.add<edm::InputTag>("endcapClusters", edm::InputTag("mtdClusters:FTLEndcap"));
0051   descriptions.add("mtdTrackingRecHitProducer", desc);
0052 }
0053 
0054 void MTDTrackingRecHitProducer::produce(edm::StreamID, edm::Event& evt, const edm::EventSetup& es) const {
0055   auto const& geom = es.getData(mtdgeoToken_);
0056 
0057   auto const& cpe = es.getData(cpeToken_);
0058 
0059   edm::Handle<FTLClusterCollection> inputBarrel;
0060   evt.getByToken(ftlbClusters_, inputBarrel);
0061 
0062   edm::Handle<FTLClusterCollection> inputEndcap;
0063   evt.getByToken(ftleClusters_, inputEndcap);
0064 
0065   std::array<edm::Handle<FTLClusterCollection>, 2> inputHandle{{inputBarrel, inputEndcap}};
0066 
0067   auto outputhits = std::make_unique<MTDTrackingDetSetVector>();
0068   auto& theoutputhits = *outputhits;
0069 
0070   //---------------------------------------------------------------------------
0071   //!  Iterate over DetUnits, then over Clusters and invoke the CPE on each,
0072   //!  and make a RecHit to store the result.
0073   //---------------------------------------------------------------------------
0074 
0075   for (auto const& theInput : inputHandle) {
0076     if (!theInput.isValid()) {
0077       edm::LogWarning("MTDTrackingRecHitProducer") << "MTDTrackingRecHitProducer: Invalid collection";
0078       continue;
0079     }
0080     const edmNew::DetSetVector<FTLCluster>& input = *theInput;
0081 
0082     LogDebug("MTDTrackingRecHitProducer") << "inputCollection " << input.size();
0083     for (const auto& DSVit : input) {
0084       unsigned int detid = DSVit.detId();
0085       DetId detIdObject(detid);
0086       const auto genericDet = geom.idToDetUnit(detIdObject);
0087       if (genericDet == nullptr) {
0088         throw cms::Exception("MTDTrackingRecHitProducer")
0089             << "GeographicalID: " << std::hex << detid << " is invalid!" << std::dec << std::endl;
0090       }
0091 
0092       MTDTrackingDetSetVector::FastFiller recHitsOnDet(theoutputhits, detid);
0093 
0094       LogDebug("MTDTrackingRecHitProducer") << "MTD cluster DetId " << detid << " # cluster " << DSVit.size();
0095 
0096       for (const auto& clustIt : DSVit) {
0097         LogDebug("MTDTrackingRecHitProducer") << "Cluster: size " << clustIt.size() << " " << clustIt.x() << ","
0098                                               << clustIt.y() << " " << clustIt.energy() << " " << clustIt.time();
0099         MTDClusterParameterEstimator::ReturnType tuple = cpe.getParameters(clustIt, *genericDet);
0100         LocalPoint lp(std::get<0>(tuple));
0101         LocalError le(std::get<1>(tuple));
0102 
0103         // Create a persistent edm::Ref to the cluster
0104         edm::Ref<edmNew::DetSetVector<FTLCluster>, FTLCluster> cluster = edmNew::makeRefTo(theInput, &clustIt);
0105         // Make a RecHit and add it to the DetSet
0106         MTDTrackingRecHit hit(lp, le, *genericDet, cluster);
0107         LogDebug("MTDTrackingRecHitProducer")
0108             << "MTD_TRH: " << hit.localPosition().x() << "," << hit.localPosition().y() << " : "
0109             << hit.localPositionError().xx() << "," << hit.localPositionError().yy() << " : " << hit.time() << " : "
0110             << hit.timeError();
0111         // Now save it =================
0112         recHitsOnDet.push_back(hit);
0113       }  //  <-- End loop on Clusters
0114     }    //    <-- End loop on DetUnits
0115     LogDebug("MTDTrackingRecHitProducer") << "outputCollection " << theoutputhits.size();
0116   }
0117 
0118   evt.put(std::move(outputhits));
0119 }
0120 
0121 #include "FWCore/Framework/interface/MakerMacros.h"
0122 DEFINE_FWK_MODULE(MTDTrackingRecHitProducer);