Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-03-17 11:22:07 

0001 // -*- C++ -*-
0002 //
0003 // Package:    ConversionSeedFilter
0004 // Class:      ConversionSeedFilter
0005 //
0006 /**\class ConversionSeedFilter ConversionSeedFilter.cc RecoTracker/TkSeedGenerator/src/ConversionSeedFilter.cc
0007 
0008  Description: <one line class summary>
0009 
0010  Implementation:
0011      <Notes on implementation>
0012 */
0013 //
0014 // Original Author:  Giuseppe Cerati & Domenico Giordano
0015 //         Created:  Thu Mar 11 10:48:48 CET 2010
0016 //
0017 //
0018 
0019 #include <memory>
0020 #include "DataFormats/Math/interface/deltaPhi.h"
0021 #include "DataFormats/TrackReco/interface/Track.h"
0022 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0023 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0024 #include "FWCore/Framework/interface/Event.h"
0025 #include "FWCore/Framework/interface/Frameworkfwd.h"
0026 #include "FWCore/Framework/interface/MakerMacros.h"
0027 #include "FWCore/Framework/interface/stream/EDProducer.h"
0028 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0029 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0030 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0031 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0032 #include "MagneticField/Engine/interface/MagneticField.h"
0033 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0034 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0035 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0036 #include "TrackingTools/Records/interface/TrackingComponentsRecord.h"
0037 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0038 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0039 
0040 class ConversionSeedFilter : public edm::stream::EDProducer<> {
0041 public:
0042   explicit ConversionSeedFilter(const edm::ParameterSet&);
0043   ~ConversionSeedFilter() override = default;
0044 
0045 private:
0046   void produce(edm::Event&, const edm::EventSetup&) override;
0047   bool isCompatible(double* vars1, double* vars2);
0048   void getKine(const TrajectoryStateOnSurface& tsos, double* vars);
0049   void SearchAmongSeeds(const TrajectorySeedCollection* pInPos,
0050                         const TrajectorySeedCollection* pInNeg,
0051                         TrajectorySeedCollection& selectedColl,
0052                         std::vector<bool>& idxPosColl1,
0053                         std::vector<bool>& idxPosColl2);
0054   void SearchAmongTracks(const TrajectorySeedCollection* pInSeed,
0055                          const reco::TrackCollection* pInTk,
0056                          TrajectorySeedCollection& selectedColl,
0057                          std::vector<bool>& idxPosColl);
0058   void SearchAmongTrajectories(const TrajectorySeedCollection* pInSeed,
0059                                const Trajectory* InTj,
0060                                TrajectorySeedCollection& selectedColl,
0061                                std::vector<bool>& idxPosColl);
0062 
0063   TrajectoryStateOnSurface getTSOS(const TrajectorySeed& ts);
0064   TrajectoryStateOnSurface getTSOS(const reco::Track& tk);
0065   TrajectoryStateOnSurface getTSOS(const Trajectory& tj, const TrajectorySeed& ts);
0066 
0067   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken;
0068   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> mfToken;
0069   const edm::EDGetTokenT<TrajectorySeedCollection> inputCollSeedPos;
0070   const edm::EDGetTokenT<TrajectorySeedCollection> inputCollSeedNeg;
0071   const edm::EDGetTokenT<TrajTrackAssociationCollection> inputTrajectory;
0072   const double deltaPhiCut, deltaCotThetaCut, deltaRCut, deltaZCut;
0073 
0074   const TrackerGeometry* theG;
0075   const MagneticField* theMF;
0076   uint32_t maxInputSeeds;
0077   bool takeAll;
0078 };
0079 
0080 ConversionSeedFilter::ConversionSeedFilter(const edm::ParameterSet& cfg)
0081     : geomToken(esConsumes()),
0082       mfToken(esConsumes()),
0083       //inputCollTkPos(cfg.getParameter<edm::InputTag>("tkCollectionPos")),
0084       inputCollSeedPos(consumes<TrajectorySeedCollection>(cfg.getParameter<edm::InputTag>("seedCollectionPos"))),
0085       //inputCollTkNeg(cfg.getParameter<edm::InputTag>("tkCollectionNeg")),
0086       inputCollSeedNeg(consumes<TrajectorySeedCollection>(cfg.getParameter<edm::InputTag>("seedCollectionNeg"))),
0087       inputTrajectory(consumes<TrajTrackAssociationCollection>(cfg.getParameter<edm::InputTag>("inputTrajectory"))),
0088       deltaPhiCut(cfg.getParameter<double>("deltaPhiCut")),
0089       deltaCotThetaCut(cfg.getParameter<double>("deltaCotThetaCut")),
0090       deltaRCut(cfg.getParameter<double>("deltaRCut")),
0091       deltaZCut(cfg.getParameter<double>("deltaZCut")),
0092       maxInputSeeds(cfg.getParameter<uint32_t>("maxInputSeeds")),
0093       takeAll(cfg.getParameter<bool>("takeAll")) {
0094   produces<TrajectorySeedCollection>();
0095 }
0096 
0097 void ConversionSeedFilter::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0098   using namespace edm;
0099   using namespace std;
0100   Handle<TrajectorySeedCollection> pInPos;
0101   iEvent.getByToken(inputCollSeedPos, pInPos);
0102   Handle<TrajectorySeedCollection> pInNeg;
0103   iEvent.getByToken(inputCollSeedNeg, pInNeg);
0104 
0105   // Handle<reco::TrackCollection> pInTkPos;  iEvent.getByLabel(inputCollTkPos,pInTkPos);
0106   // Handle<reco::TrackCollection> pInTkNeg;  iEvent.getByLabel(inputCollTkNeg,pInTkNeg);
0107 
0108   edm::Handle<TrajTrackAssociationCollection> trajTrackAssociations;
0109   iEvent.getByToken(inputTrajectory, trajTrackAssociations);
0110 
0111   theG = &iSetup.getData(geomToken);
0112   theMF = &iSetup.getData(mfToken);
0113 
0114   auto result = std::make_unique<TrajectorySeedCollection>();
0115 
0116   TrajectorySeedCollection selectedColl;
0117 
0118   if (takeAll) {
0119     result->insert(result->end(), pInPos->begin(), pInPos->end());
0120     result->insert(result->end(), pInNeg->begin(), pInNeg->end());
0121   } else {
0122     edm::LogInfo("ConversionSeedFilter") << "takeAll " << takeAll;
0123     if (pInPos->size() < maxInputSeeds && pInNeg->size() < maxInputSeeds) {
0124       std::vector<bool> idxPosColl1(pInPos->size(), false);
0125       std::vector<bool> idxPosColl2(pInNeg->size(), false);
0126       selectedColl.reserve(pInPos->size());
0127 
0128       edm::LogInfo("ConversionSeedFilter")
0129           << "New Event \t Pos " << pInPos->size() << " \t Neg " << pInNeg->size() << std::endl;
0130 
0131       SearchAmongSeeds(pInPos.product(), pInNeg.product(), selectedColl, idxPosColl1, idxPosColl2);
0132       //SearchAmongTracks(pInPos.product(),pInTkNeg.product(),selectedColl,idxPosColl1);
0133       //SearchAmongTracks(pInNeg.product(),pInTkPos.product(),selectedColl,idxPosColl2);
0134 
0135       if (trajTrackAssociations.isValid()) {
0136         edm::LogInfo("ConversionSeedFilter") << "Reconstructed tracks " << trajTrackAssociations->size() << std::endl;
0137         for (TrajTrackAssociationCollection::const_iterator association = trajTrackAssociations->begin();
0138              association != trajTrackAssociations->end();
0139              association++) {
0140           const Trajectory* traj = association->key.get();
0141           const reco::Track* track = association->val.get();
0142 
0143           //edm::LogInfo("ConversionSeedFilter") << "Traj charge " << track->charge() << std::endl;
0144 
0145           if (track->charge() < 0) {
0146             SearchAmongTrajectories(pInPos.product(), traj, selectedColl, idxPosColl1);
0147           } else {
0148             SearchAmongTrajectories(pInNeg.product(), traj, selectedColl, idxPosColl2);
0149           }
0150         }
0151       }
0152     }
0153     result->insert(result->end(), selectedColl.begin(), selectedColl.end());
0154   }
0155 
0156   edm::LogInfo("ConversionSeedFilter") << "\nNew Event : result size " << result->size() << std::endl;
0157 
0158   iEvent.put(std::move(result));
0159 }
0160 
0161 void ConversionSeedFilter::SearchAmongSeeds(const TrajectorySeedCollection* pInPos,
0162                                             const TrajectorySeedCollection* pInNeg,
0163                                             TrajectorySeedCollection& selectedColl,
0164                                             std::vector<bool>& idxPosColl1,
0165                                             std::vector<bool>& idxPosColl2) {
0166   for (TrajectorySeedCollection::const_iterator iS1 = pInPos->begin(); iS1 != pInPos->end(); ++iS1) {
0167     bool pushed1 = false;
0168 
0169     double vars1[4];
0170     getKine(getTSOS(*iS1), vars1);
0171 
0172     for (TrajectorySeedCollection::const_iterator iS2 = pInNeg->begin(); iS2 != pInNeg->end(); ++iS2) {
0173       double vars2[4];
0174       getKine(getTSOS(*iS2), vars2);
0175 
0176       if (isCompatible(vars1, vars2)) {
0177         edm::LogInfo("ConversionSeedFilter")
0178             << "[SearchAmongSeeds] match in pos " << iS1 - pInPos->begin() << " " << iS2 - pInNeg->begin() << std::endl;
0179         if (!pushed1) {
0180           idxPosColl1[iS1 - pInPos->begin()] = true;
0181           selectedColl.push_back(*iS1);
0182           pushed1 = true;
0183         }
0184         if (!idxPosColl2[iS2 - pInNeg->begin()]) {
0185           selectedColl.push_back(*iS2);
0186           idxPosColl2[iS2 - pInNeg->begin()] = true;
0187         }
0188       }
0189     }
0190   }
0191 }
0192 
0193 void ConversionSeedFilter::SearchAmongTracks(const TrajectorySeedCollection* pInSeed,
0194                                              const reco::TrackCollection* pInTk,
0195                                              TrajectorySeedCollection& selectedColl,
0196                                              std::vector<bool>& idxPosColl) {
0197   for (TrajectorySeedCollection::const_iterator iS1 = pInSeed->begin(); iS1 != pInSeed->end(); ++iS1) {
0198     if (idxPosColl[iS1 - pInSeed->begin()])
0199       continue;
0200 
0201     double vars1[4];
0202     getKine(getTSOS(*iS1), vars1);
0203 
0204     for (reco::TrackCollection::const_iterator iS2 = pInTk->begin(); iS2 != pInTk->end(); ++iS2) {
0205       double vars2[4];
0206       getKine(getTSOS(*iS2), vars2);
0207 
0208       if (isCompatible(vars1, vars2)) {
0209         edm::LogInfo("ConversionSeedFilter")
0210             << "[SearchAmongTracks] match in pos " << iS1 - pInSeed->begin() << std::endl;
0211         idxPosColl[iS1 - pInSeed->begin()] = true;
0212         selectedColl.push_back(*iS1);
0213         continue;
0214       }
0215     }
0216   }
0217 }
0218 
0219 void ConversionSeedFilter::SearchAmongTrajectories(const TrajectorySeedCollection* pInSeed,
0220                                                    const Trajectory* InTj,
0221                                                    TrajectorySeedCollection& selectedColl,
0222                                                    std::vector<bool>& idxPosColl) {
0223   for (TrajectorySeedCollection::const_iterator iS1 = pInSeed->begin(); iS1 != pInSeed->end(); ++iS1) {
0224     if (idxPosColl[iS1 - pInSeed->begin()])
0225       continue;
0226 
0227     double vars1[4];
0228     getKine(getTSOS(*iS1), vars1);
0229 
0230     double vars2[4];
0231     getKine(getTSOS(*InTj, *iS1), vars2);
0232 
0233     if (isCompatible(vars1, vars2)) {
0234       edm::LogInfo("ConversionSeedFilter")
0235           << "[SearchAmongTrajectories] match seed in pos " << iS1 - pInSeed->begin() << " of " << pInSeed->size()
0236           << " seed charge " << iS1->startingState().parameters().charge() << std::endl;
0237       idxPosColl[iS1 - pInSeed->begin()] = true;
0238       selectedColl.push_back(*iS1);
0239     }
0240   }
0241 }
0242 
0243 TrajectoryStateOnSurface ConversionSeedFilter::getTSOS(const TrajectorySeed& ts) {
0244   PTrajectoryStateOnDet state = ts.startingState();
0245   DetId detId(state.detId());
0246   return trajectoryStateTransform::transientState(state, &(theG->idToDet(detId)->surface()), theMF);
0247 }
0248 
0249 TrajectoryStateOnSurface ConversionSeedFilter::getTSOS(const reco::Track& tk) {
0250   return trajectoryStateTransform::innerStateOnSurface(tk, *theG, theMF);
0251 }
0252 
0253 TrajectoryStateOnSurface ConversionSeedFilter::getTSOS(const Trajectory& tj, const TrajectorySeed& ts) {
0254   PTrajectoryStateOnDet state = ts.startingState();
0255   DetId detId(state.detId());
0256   GlobalPoint p = theG->idToDet(detId)->surface().toGlobal(state.parameters().position());
0257   return tj.closestMeasurement(p).updatedState();
0258 }
0259 
0260 bool ConversionSeedFilter::isCompatible(double* vars1, double* vars2) {
0261   double deltaPhi = fabs(reco::deltaPhi(vars1[0], vars2[0]));
0262   double deltaCotTheta = fabs(vars1[1] - vars2[1]);
0263   double deltaR = fabs(vars1[2] - vars2[2]);
0264   double deltaZ = fabs(vars1[3] - vars2[3]);
0265   /*
0266   std::cout << "\nnewCheck" << std::endl;
0267   if(deltaPhi>deltaPhiCut)             std::cout << "\nko deltaphi" << deltaPhi; 
0268   if(deltaCotTheta>deltaCotThetaCut )  std::cout << "\nko deltaCt " << deltaCotTheta; 
0269   if(deltaR>deltaRCut           )  std::cout << "\nko deltaR  " << deltaR; 
0270   if(deltaZ>deltaZCut               )  std::cout << "\nko deltaZ  " << deltaZ; 
0271   if(deltaPhi<deltaPhiCut && deltaCotTheta<deltaCotThetaCut && deltaR<deltaRCut && deltaZ<deltaZCut) std::cout << "\nok :)\n";
0272   */
0273   return deltaPhi < deltaPhiCut && deltaCotTheta < deltaCotThetaCut && deltaR < deltaRCut && deltaZ < deltaZCut;
0274 }
0275 
0276 void ConversionSeedFilter::getKine(const TrajectoryStateOnSurface& tsos, double* vars) {
0277   vars[0] = tsos.globalMomentum().phi();             //phi
0278   vars[1] = 2 / tan(tsos.globalMomentum().theta());  //cotTheta
0279   vars[2] = tsos.globalPosition().perp();            //R
0280   vars[3] = tsos.globalPosition().z();               //Z
0281 }
0282 
0283 DEFINE_FWK_MODULE(ConversionSeedFilter);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!