Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-01 02:23:32

0001 #include "CUDADataFormats/Vertex/interface/ZVertexHeterogeneous.h"
0002 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0003 #include "DataFormats/Common/interface/OrphanHandle.h"
0004 #include "DataFormats/TrackReco/interface/Track.h"
0005 #include "DataFormats/TrackReco/interface/TrackExtra.h"
0006 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0007 #include "DataFormats/VertexReco/interface/Vertex.h"
0008 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0009 #include "FWCore/Framework/interface/ConsumesCollector.h"
0010 #include "FWCore/Framework/interface/ESHandle.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/EventSetup.h"
0013 #include "FWCore/Framework/interface/MakerMacros.h"
0014 #include "FWCore/Framework/interface/global/EDProducer.h"
0015 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0018 #include "FWCore/PluginManager/interface/ModuleDef.h"
0019 #include "FWCore/Utilities/interface/EDGetToken.h"
0020 #include "FWCore/Utilities/interface/InputTag.h"
0021 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0022 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0023 
0024 #undef PIXVERTEX_DEBUG_PRODUCE
0025 
0026 class PixelVertexProducerFromSoA : public edm::global::EDProducer<> {
0027 public:
0028   using IndToEdm = std::vector<uint16_t>;
0029 
0030   explicit PixelVertexProducerFromSoA(const edm::ParameterSet &iConfig);
0031   ~PixelVertexProducerFromSoA() override = default;
0032 
0033   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0034 
0035 private:
0036   void produce(edm::StreamID streamID, edm::Event &iEvent, const edm::EventSetup &iSetup) const override;
0037 
0038   edm::EDGetTokenT<ZVertexHeterogeneous> tokenVertex_;
0039   edm::EDGetTokenT<reco::BeamSpot> tokenBeamSpot_;
0040   edm::EDGetTokenT<reco::TrackCollection> tokenTracks_;
0041   edm::EDGetTokenT<IndToEdm> tokenIndToEdm_;
0042 };
0043 
0044 PixelVertexProducerFromSoA::PixelVertexProducerFromSoA(const edm::ParameterSet &conf)
0045     : tokenVertex_(consumes<ZVertexHeterogeneous>(conf.getParameter<edm::InputTag>("src"))),
0046       tokenBeamSpot_(consumes<reco::BeamSpot>(conf.getParameter<edm::InputTag>("beamSpot"))),
0047       tokenTracks_(consumes<reco::TrackCollection>(conf.getParameter<edm::InputTag>("TrackCollection"))),
0048       tokenIndToEdm_(consumes<IndToEdm>(conf.getParameter<edm::InputTag>("TrackCollection"))) {
0049   produces<reco::VertexCollection>();
0050 }
0051 
0052 void PixelVertexProducerFromSoA::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0053   edm::ParameterSetDescription desc;
0054 
0055   desc.add<edm::InputTag>("TrackCollection", edm::InputTag("pixelTracks"));
0056   desc.add<edm::InputTag>("beamSpot", edm::InputTag("offlineBeamSpot"));
0057   desc.add<edm::InputTag>("src", edm::InputTag("pixelVerticesSoA"));
0058 
0059   descriptions.add("pixelVertexFromSoA", desc);
0060 }
0061 
0062 void PixelVertexProducerFromSoA::produce(edm::StreamID streamID, edm::Event &iEvent, const edm::EventSetup &) const {
0063   auto vertexes = std::make_unique<reco::VertexCollection>();
0064 
0065   auto tracksHandle = iEvent.getHandle(tokenTracks_);
0066   auto tracksSize = tracksHandle->size();
0067   auto const &indToEdm = iEvent.get(tokenIndToEdm_);
0068   auto bsHandle = iEvent.getHandle(tokenBeamSpot_);
0069 
0070   float x0 = 0, y0 = 0, z0 = 0, dxdz = 0, dydz = 0;
0071   std::vector<int32_t> itrk;
0072   itrk.reserve(64);  // avoid first relocations
0073   if (!bsHandle.isValid()) {
0074     edm::LogWarning("PixelVertexProducer") << "No beamspot found. returning vertexes with (0,0,Z) ";
0075   } else {
0076     const reco::BeamSpot &bs = *bsHandle;
0077     x0 = bs.x0();
0078     y0 = bs.y0();
0079     z0 = bs.z0();
0080     dxdz = bs.dxdz();
0081     dydz = bs.dydz();
0082   }
0083 
0084   auto const &soa = *(iEvent.get(tokenVertex_).get());
0085 
0086   int nv = soa.nvFinal;
0087 
0088 #ifdef PIXVERTEX_DEBUG_PRODUCE
0089   std::cout << "converting " << nv << " vertices "
0090             << " from " << indToEdm.size() << " tracks" << std::endl;
0091 #endif  // PIXVERTEX_DEBUG_PRODUCE
0092 
0093   std::set<uint16_t> uind;  // for verifing index consistency
0094   for (int j = nv - 1; j >= 0; --j) {
0095     auto i = soa.sortInd[j];  // on gpu sorted in ascending order....
0096     assert(i < nv);
0097     uind.insert(i);
0098     assert(itrk.empty());
0099     auto z = soa.zv[i];
0100     auto x = x0 + dxdz * z;
0101     auto y = y0 + dydz * z;
0102     z += z0;
0103     reco::Vertex::Error err;
0104     err(2, 2) = 1.f / soa.wv[i];
0105     err(2, 2) *= 2.;  // artifically inflate error
0106     //Copy also the tracks (no intention to be efficient....)
0107     for (auto k = 0U; k < indToEdm.size(); ++k) {
0108       if (soa.idv[k] == int16_t(i))
0109         itrk.push_back(k);
0110     }
0111     auto nt = itrk.size();
0112     if (nt == 0) {
0113 #ifdef PIXVERTEX_DEBUG_PRODUCE
0114       std::cout << "vertex " << i << " with no tracks..." << std::endl;
0115 #endif  // PIXVERTEX_DEBUG_PRODUCE
0116       continue;
0117     }
0118     if (nt < 2) {
0119       itrk.clear();
0120       continue;
0121     }  // remove outliers
0122     (*vertexes).emplace_back(reco::Vertex::Point(x, y, z), err, soa.chi2[i], soa.ndof[i], nt);
0123     auto &v = (*vertexes).back();
0124     v.reserve(itrk.size());
0125     for (auto it : itrk) {
0126       assert(it < int(indToEdm.size()));
0127       auto k = indToEdm[it];
0128       if (k > tracksSize) {
0129         edm::LogWarning("PixelVertexProducer") << "oops track " << it << " does not exists on CPU " << k;
0130         continue;
0131       }
0132       auto tk = reco::TrackRef(tracksHandle, k);
0133       v.add(tk);
0134     }
0135     itrk.clear();
0136   }
0137 
0138   LogDebug("PixelVertexProducer") << ": Found " << vertexes->size() << " vertexes\n";
0139   for (unsigned int i = 0; i < vertexes->size(); ++i) {
0140     LogDebug("PixelVertexProducer") << "Vertex number " << i << " has " << (*vertexes)[i].tracksSize()
0141                                     << " tracks with a position of " << (*vertexes)[i].z() << " +- "
0142                                     << std::sqrt((*vertexes)[i].covariance(2, 2));
0143   }
0144 
0145   // legacy logic....
0146   if (vertexes->empty() && bsHandle.isValid()) {
0147     const reco::BeamSpot &bs = *bsHandle;
0148 
0149     GlobalError bse(bs.rotatedCovariance3D());
0150     if ((bse.cxx() <= 0.) || (bse.cyy() <= 0.) || (bse.czz() <= 0.)) {
0151       AlgebraicSymMatrix33 we;
0152       we(0, 0) = 10000;
0153       we(1, 1) = 10000;
0154       we(2, 2) = 10000;
0155       vertexes->push_back(reco::Vertex(bs.position(), we, 0., 0., 0));
0156 
0157       edm::LogInfo("PixelVertexProducer") << "No vertices found. Beamspot with invalid errors " << bse.matrix()
0158                                           << "\nWill put Vertex derived from dummy-fake BeamSpot into Event.\n"
0159                                           << (*vertexes)[0].x() << "\n"
0160                                           << (*vertexes)[0].y() << "\n"
0161                                           << (*vertexes)[0].z() << "\n";
0162     } else {
0163       vertexes->push_back(reco::Vertex(bs.position(), bs.rotatedCovariance3D(), 0., 0., 0));
0164 
0165       edm::LogInfo("PixelVertexProducer") << "No vertices found. Will put Vertex derived from BeamSpot into Event:\n"
0166                                           << (*vertexes)[0].x() << "\n"
0167                                           << (*vertexes)[0].y() << "\n"
0168                                           << (*vertexes)[0].z() << "\n";
0169     }
0170   } else if (vertexes->empty() && !bsHandle.isValid()) {
0171     edm::LogWarning("PixelVertexProducer") << "No beamspot and no vertex found. No vertex returned.";
0172   }
0173 
0174   iEvent.put(std::move(vertexes));
0175 }
0176 
0177 DEFINE_FWK_MODULE(PixelVertexProducerFromSoA);