File indexing completed on 2023-03-17 11:18:08
0001
0002
0003
0004
0005
0006
0007
0008 #include <iostream>
0009
0010 #include "FWCore/Framework/interface/Frameworkfwd.h"
0011 #include "FWCore/Framework/interface/stream/EDProducer.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/Framework/interface/ESHandle.h"
0015 #include "FWCore/Framework/interface/MakerMacros.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018 #include "FWCore/Utilities/interface/ESGetToken.h"
0019
0020 #include "DataFormats/Common/interface/Handle.h"
0021 #include "DataFormats/GeometryVector/interface/LocalPoint.h"
0022 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0023 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0024 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
0025 #include "DataFormats/HeavyIonEvent/interface/ClusterCompatibility.h"
0026
0027 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0028 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0029 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0030 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0031 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0032
0033
0034
0035
0036
0037 class ClusterCompatibilityProducer : public edm::stream::EDProducer<> {
0038 public:
0039 explicit ClusterCompatibilityProducer(const edm::ParameterSet &);
0040 ~ClusterCompatibilityProducer() override;
0041
0042 void produce(edm::Event &, const edm::EventSetup &) override;
0043
0044 private:
0045 edm::EDGetTokenT<SiPixelRecHitCollection> inputToken_;
0046 edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> trackerToken_;
0047 edm::InputTag inputTag_;
0048 double minZ_;
0049 double maxZ_;
0050 double zStep_;
0051
0052 struct VertexHit {
0053 float z;
0054 float r;
0055 float w;
0056 };
0057
0058 struct ContainedHits {
0059 float z0;
0060 int nHit;
0061 float chi;
0062 };
0063
0064 ContainedHits getContainedHits(const std::vector<VertexHit> &hits, double z0) const;
0065 };
0066
0067 ClusterCompatibilityProducer::ClusterCompatibilityProducer(const edm::ParameterSet &config)
0068 : inputTag_(config.getParameter<edm::InputTag>("inputTag")),
0069 minZ_(config.getParameter<double>("minZ")),
0070 maxZ_(config.getParameter<double>("maxZ")),
0071 zStep_(config.getParameter<double>("zStep")) {
0072 inputToken_ = consumes<SiPixelRecHitCollection>(inputTag_);
0073 trackerToken_ = esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>();
0074 LogDebug("") << "Using the " << inputTag_ << " input collection";
0075 produces<reco::ClusterCompatibility>();
0076 }
0077
0078 ClusterCompatibilityProducer::~ClusterCompatibilityProducer() {}
0079
0080 void ClusterCompatibilityProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0081 auto creco = std::make_unique<reco::ClusterCompatibility>();
0082
0083
0084 edm::Handle<SiPixelRecHitCollection> hRecHits;
0085 iEvent.getByToken(inputToken_, hRecHits);
0086
0087
0088 if (hRecHits.isValid()) {
0089 edm::ESHandle<TrackerGeometry> trackerHandle = iSetup.getHandle(trackerToken_);
0090 const TrackerGeometry *tgeo = trackerHandle.product();
0091 const SiPixelRecHitCollection *hits = hRecHits.product();
0092
0093
0094 int nPxlHits = 0;
0095 std::vector<VertexHit> vhits;
0096 for (SiPixelRecHitCollection::DataContainer::const_iterator hit = hits->data().begin(), end = hits->data().end();
0097 hit != end;
0098 ++hit) {
0099 if (!hit->isValid())
0100 continue;
0101 ++nPxlHits;
0102 DetId id(hit->geographicalId());
0103 if (id.subdetId() != int(PixelSubdetector::PixelBarrel))
0104 continue;
0105 const PixelGeomDetUnit *pgdu = static_cast<const PixelGeomDetUnit *>(tgeo->idToDet(id));
0106 const PixelTopology *pixTopo = &(pgdu->specificTopology());
0107 std::vector<SiPixelCluster::Pixel> pixels(hit->cluster()->pixels());
0108 bool pixelOnEdge = false;
0109 for (std::vector<SiPixelCluster::Pixel>::const_iterator pixel = pixels.begin(); pixel != pixels.end(); ++pixel) {
0110 int pixelX = pixel->x;
0111 int pixelY = pixel->y;
0112 if (pixTopo->isItEdgePixelInX(pixelX) || pixTopo->isItEdgePixelInY(pixelY)) {
0113 pixelOnEdge = true;
0114 break;
0115 }
0116 }
0117 if (pixelOnEdge)
0118 continue;
0119
0120 LocalPoint lpos = LocalPoint(hit->localPosition().x(), hit->localPosition().y(), hit->localPosition().z());
0121 GlobalPoint gpos = pgdu->toGlobal(lpos);
0122 VertexHit vh;
0123 vh.z = gpos.z();
0124 vh.r = gpos.perp();
0125 vh.w = hit->cluster()->sizeY();
0126 vhits.push_back(vh);
0127 }
0128
0129 creco->setNValidPixelHits(nPxlHits);
0130
0131
0132 for (double z0 = minZ_; z0 <= maxZ_; z0 += zStep_) {
0133 ContainedHits c = getContainedHits(vhits, z0);
0134 creco->append(c.z0, c.nHit, c.chi);
0135 }
0136 }
0137 iEvent.put(std::move(creco));
0138 }
0139
0140 ClusterCompatibilityProducer::ContainedHits ClusterCompatibilityProducer::getContainedHits(
0141 const std::vector<VertexHit> &hits, double z0) const {
0142
0143 int n = 0;
0144 double chi = 0.;
0145
0146 for (std::vector<VertexHit>::const_iterator hit = hits.begin(); hit != hits.end(); hit++) {
0147
0148
0149
0150
0151 double p = 2 * fabs(hit->z - z0) / hit->r + 0.5;
0152 if (fabs(p - hit->w) <= 1.) {
0153 chi += fabs(p - hit->w);
0154 n++;
0155 }
0156 }
0157 ClusterCompatibilityProducer::ContainedHits output;
0158 output.z0 = z0;
0159 output.nHit = n;
0160 output.chi = chi;
0161 return output;
0162 }
0163
0164
0165 DEFINE_FWK_MODULE(ClusterCompatibilityProducer);