File indexing completed on 2024-04-06 12:29:59
0001
0002
0003
0004
0005
0006
0007 #include "DataFormats/Common/interface/Handle.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/EventSetup.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/Framework/interface/stream/EDProducer.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 #include "FWCore/PluginManager/interface/ModuleDef.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015
0016 #include "Geometry/EcalTestBeam/interface/EcalTBHodoscopeGeometry.h"
0017 #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
0018 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
0019 #include "SimDataFormats/EcalTestBeam/interface/HodoscopeDetId.h"
0020 #include "SimDataFormats/EcalTestBeam/interface/PEcalTBInfo.h"
0021 #include "TBDataFormats/EcalTBObjects/interface/EcalTBHodoscopePlaneRawHits.h"
0022 #include "TBDataFormats/EcalTBObjects/interface/EcalTBHodoscopeRawInfo.h"
0023
0024 #include <iostream>
0025 #include <memory>
0026
0027 class TBHodoActiveVolumeRawInfoProducer : public edm::stream::EDProducer<> {
0028 public:
0029
0030 explicit TBHodoActiveVolumeRawInfoProducer(const edm::ParameterSet &ps);
0031
0032
0033 ~TBHodoActiveVolumeRawInfoProducer() override = default;
0034
0035
0036 void produce(edm::Event &event, const edm::EventSetup &eventSetup) override;
0037
0038 private:
0039 double myThreshold;
0040 edm::EDGetTokenT<edm::PCaloHitContainer> m_EcalToken;
0041
0042 std::unique_ptr<EcalTBHodoscopeGeometry> theTBHodoGeom_;
0043 };
0044
0045 using namespace cms;
0046 using namespace std;
0047
0048 TBHodoActiveVolumeRawInfoProducer::TBHodoActiveVolumeRawInfoProducer(const edm::ParameterSet &ps)
0049 : myThreshold(0.05E-3),
0050 m_EcalToken(consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalTBH4BeamHits"))) {
0051 produces<EcalTBHodoscopeRawInfo>();
0052
0053 theTBHodoGeom_ = std::make_unique<EcalTBHodoscopeGeometry>();
0054 }
0055
0056 void TBHodoActiveVolumeRawInfoProducer::produce(edm::Event &event, const edm::EventSetup &eventSetup) {
0057 std::unique_ptr<EcalTBHodoscopeRawInfo> product(new EcalTBHodoscopeRawInfo());
0058
0059
0060 const edm::Handle<edm::PCaloHitContainer> &pCaloHit = event.getHandle(m_EcalToken);
0061 const edm::PCaloHitContainer *caloHits = nullptr;
0062 if (pCaloHit.isValid()) {
0063 caloHits = pCaloHit.product();
0064 LogDebug("EcalTBHodo") << "total # caloHits: " << caloHits->size();
0065 } else {
0066 edm::LogError("EcalTBHodo") << "Error! can't get the caloHitContainer ";
0067 }
0068 if (!caloHits) {
0069 return;
0070 }
0071
0072
0073 std::map<unsigned int, double> energyMap;
0074
0075 for (auto &&aHit : *caloHits) {
0076 double thisHitEne = aHit.energy();
0077
0078 std::map<unsigned int, double>::iterator itmap = energyMap.find(aHit.id());
0079 if (itmap == energyMap.end())
0080 energyMap.insert(pair<unsigned int, double>(aHit.id(), thisHitEne));
0081 else {
0082 (*itmap).second += thisHitEne;
0083 }
0084 }
0085
0086
0087 int nPlanes = theTBHodoGeom_->getNPlanes();
0088 int nFibers = theTBHodoGeom_->getNFibres();
0089 product->setPlanes(nPlanes);
0090
0091 bool firedChannels[4][64];
0092 for (int iPlane = 0; iPlane < nPlanes; ++iPlane) {
0093 for (int iFiber = 0; iFiber < nFibers; ++iFiber) {
0094 firedChannels[iPlane][iFiber] = 0.;
0095 }
0096 }
0097 for (std::map<unsigned int, double>::const_iterator itmap = energyMap.begin(); itmap != energyMap.end(); ++itmap) {
0098 if ((*itmap).second > myThreshold) {
0099 HodoscopeDetId myHodoDetId = HodoscopeDetId((*itmap).first);
0100 firedChannels[myHodoDetId.planeId()][myHodoDetId.fibrId()] = true;
0101 }
0102 }
0103 for (int iPlane = 0; iPlane < nPlanes; ++iPlane) {
0104 EcalTBHodoscopePlaneRawHits planeHit(nFibers);
0105
0106 for (int iFiber = 0; iFiber < nFibers; ++iFiber) {
0107 planeHit.setHit(iFiber, firedChannels[iPlane][iFiber]);
0108 }
0109 product->setPlane((unsigned int)iPlane, planeHit);
0110 }
0111
0112 LogDebug("EcalTBHodo") << (*product);
0113
0114 event.put(std::move(product));
0115 }
0116
0117 DEFINE_FWK_MODULE(TBHodoActiveVolumeRawInfoProducer);