Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-04-21 02:53:31

0001 #include "CLHEP/Units/GlobalPhysicalConstants.h"
0002 #include "SimFastTiming/FastTimingCommon/interface/ETLDeviceSim.h"
0003 #include "DataFormats/Math/interface/GeantUnits.h"
0004 
0005 #include "DataFormats/DetId/interface/DetId.h"
0006 #include "DataFormats/ForwardDetId/interface/MTDDetId.h"
0007 #include "DataFormats/ForwardDetId/interface/ETLDetId.h"
0008 #include "FWCore/Framework/interface/ConsumesCollector.h"
0009 
0010 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0011 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0012 #include "Geometry/MTDGeometryBuilder/interface/ProxyMTDTopology.h"
0013 #include "Geometry/MTDGeometryBuilder/interface/RectangularMTDTopology.h"
0014 
0015 ETLDeviceSim::ETLDeviceSim(const edm::ParameterSet& pset, edm::ConsumesCollector iC)
0016     : geomToken_(iC.esConsumes()),
0017       geom_(nullptr),
0018       MIPPerMeV_(1.0 / pset.getParameter<double>("meVPerMIP")),
0019       bxTime_(pset.getParameter<double>("bxTime")),
0020       tofDelay_(pset.getParameter<double>("tofDelay")) {}
0021 
0022 void ETLDeviceSim::getEventSetup(const edm::EventSetup& evs) { geom_ = &evs.getData(geomToken_); }
0023 
0024 void ETLDeviceSim::getHitsResponse(const std::vector<std::tuple<int, uint32_t, float> >& hitRefs,
0025                                    const edm::Handle<edm::PSimHitContainer>& hits,
0026                                    mtd_digitizer::MTDSimHitDataAccumulator* simHitAccumulator,
0027                                    CLHEP::HepRandomEngine* hre) {
0028   using namespace geant_units::operators;
0029 
0030   //loop over sorted hits
0031   const int nchits = hitRefs.size();
0032   for (int i = 0; i < nchits; ++i) {
0033     const int hitidx = std::get<0>(hitRefs[i]);
0034     const uint32_t id = std::get<1>(hitRefs[i]);
0035     const MTDDetId detId(id);
0036 
0037     // Safety check (this should never happen, it should be an exception
0038     if (detId.det() != DetId::Forward || detId.mtdSubDetector() != 2) {
0039       throw cms::Exception("ETLDeviceSim")
0040           << "got a DetId that was not ETL: Det = " << detId.det() << "  subDet = " << detId.mtdSubDetector();
0041     }
0042 
0043     if (id == 0)
0044       continue;  // to be ignored at RECO level
0045 
0046     ETLDetId etlid(detId);
0047     DetId geoId = ETLDetId(etlid.mtdSide(), etlid.mtdRR(), etlid.module(), etlid.modType());
0048     const MTDGeomDet* thedet = geom_->idToDet(geoId);
0049     if (thedet == nullptr) {
0050       throw cms::Exception("ETLDeviceSim") << "GeographicalID: " << std::hex << geoId.rawId() << " (" << detId.rawId()
0051                                            << ") is invalid!" << std::dec << std::endl;
0052     }
0053     const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0054     const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0055 
0056     const float toa = std::get<2>(hitRefs[i]) + tofDelay_;
0057     const PSimHit& hit = hits->at(hitidx);
0058     const float charge = convertGeVToMeV(hit.energyLoss()) * MIPPerMeV_;
0059 
0060     // calculate the simhit row and column
0061     const auto& pentry = hit.entryPoint();
0062     // ETL is already in module-local coordinates so just scale to cm from mm
0063     Local3DPoint simscaled(convertMmToCm(pentry.x()), convertMmToCm(pentry.y()), convertMmToCm(pentry.z()));
0064     //The following lines check whether the pixel point is actually out of the active area.
0065     //If that is the case it simply ignores the point but in the future some more sophisticated function could be applied.
0066     if (!topo.isInPixel(simscaled)) {
0067       continue;
0068     }
0069     const auto& thepixel = topo.pixel(simscaled);
0070     const uint8_t row(thepixel.first), col(thepixel.second);
0071 
0072     auto simHitIt =
0073         simHitAccumulator->emplace(mtd_digitizer::MTDCellId(id, row, col), mtd_digitizer::MTDCellInfo()).first;
0074 
0075     // Accumulate in 15 buckets of 25ns (9 pre-samples, 1 in-time, 5 post-samples)
0076     const int itime = std::floor(toa / bxTime_) + 9;
0077     if (itime < 0 || itime > 14)
0078       continue;
0079 
0080     // Check if time index is ok and store energy
0081     if (itime >= (int)simHitIt->second.hit_info[0].size())
0082       continue;
0083 
0084     (simHitIt->second).hit_info[0][itime] += charge;
0085 
0086     // Store the time of the first SimHit in the right DataFrame bucket
0087     const float tof = toa - (itime - 9) * bxTime_;
0088 
0089     if ((simHitIt->second).hit_info[1][itime] == 0. || tof < (simHitIt->second).hit_info[1][itime]) {
0090       (simHitIt->second).hit_info[1][itime] = tof;
0091     }
0092   }
0093 }