TkHitPairs/src/CosmicHitPairGeneratorFromLayerPair.cc

0001 #include "RecoTracker/TkHitPairs/interface/CosmicHitPairGeneratorFromLayerPair.h"
0002 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegion.h"
0003 #include "TrackingTools/DetLayers/interface/DetLayer.h"
0004 #include "RecoTracker/TkHitPairs/interface/OrderedHitPair.h"
0005 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "TrackingTools/DetLayers/interface/BarrelDetLayer.h"
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h"
0010 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
0011
0012 using namespace std;
0013 // typedef TransientTrackingRecHit::ConstRecHitPointer TkHitPairsCachedHit;
0014
0015 CosmicHitPairGeneratorFromLayerPair::CosmicHitPairGeneratorFromLayerPair(const LayerWithHits* inner,
0016                                                                          const LayerWithHits* outer,
0017                                                                          const TrackerGeometry& geom)
0018     : trackerGeometry(&geom), theOuterLayer(outer), theInnerLayer(inner) {}
0019 CosmicHitPairGeneratorFromLayerPair::~CosmicHitPairGeneratorFromLayerPair() {}
0020
0021 void CosmicHitPairGeneratorFromLayerPair::hitPairs(const TrackingRegion& region, OrderedHitPairs& result) {
0022   //  static int NSee = 0; static int Ntry = 0; static int Nacc = 0;
0023
0024   typedef OrderedHitPair::InnerRecHit InnerHit;
0025   typedef OrderedHitPair::OuterRecHit OuterHit;
0026
0027   if (theInnerLayer->recHits().empty())
0028     return;
0029
0030   if (theOuterLayer->recHits().empty())
0031     return;
0032   //  const DetLayer* innerlay=theOuterLayer->layer();
0033   // const BarrelDetLayer *pippo=dynamic_cast<const BarrelDetLayer*>(theOuterLayer->layer());
0034
0035   // ************ Daniele
0036
0037   const DetLayer* blay1;
0038   const DetLayer* blay2;
0039   blay1 = dynamic_cast<const BarrelDetLayer*>(theInnerLayer->layer());
0040   blay2 = dynamic_cast<const BarrelDetLayer*>(theOuterLayer->layer());
0041
0042   bool seedfromoverlaps = false;
0043   bool InTheBarrel = false;
0044   bool InTheForward = false;
0045   if (blay1 && blay2) {
0046     InTheBarrel = true;
0047   } else
0048     InTheForward = true;
0049
0050   if (InTheBarrel) {
0051     float radius1 = dynamic_cast<const BarrelDetLayer*>(theInnerLayer->layer())->specificSurface().radius();
0052     float radius2 = dynamic_cast<const BarrelDetLayer*>(theOuterLayer->layer())->specificSurface().radius();
0053     seedfromoverlaps = (abs(radius1 - radius2) < 0.1) ? true : false;
0054   }
0055
0056   vector<OrderedHitPair> allthepairs;
0057
0058   for (auto ohh = theOuterLayer->recHits().begin(); ohh != theOuterLayer->recHits().end(); ohh++) {
0059     for (auto ihh = theInnerLayer->recHits().begin(); ihh != theInnerLayer->recHits().end(); ihh++) {
0060       auto oh = static_cast<BaseTrackerRecHit const* const>(*ohh);
0061       auto ih = static_cast<BaseTrackerRecHit const* const>(*ihh);
0062
0063       float z_diff = ih->globalPosition().z() - oh->globalPosition().z();
0064       float inny = ih->globalPosition().y();
0065       float outy = oh->globalPosition().y();
0066       float innx = ih->globalPosition().x();
0067       float outx = oh->globalPosition().x();
0068       ;
0069       float dxdy = abs((outx - innx) / (outy - inny));
0070       float DeltaR = oh->globalPosition().perp() - ih->globalPosition().perp();
0071       ;
0072
0073       if (InTheBarrel &&
0074           (abs(z_diff) < 30)  // && (outy > 0.) && (inny > 0.)
0075           //&&((abs(inny-outy))<30)
0076           && (dxdy < 2) && (inny * outy > 0) && (abs(DeltaR) > 0)) {
0077         //  cout << " ******** sono dentro inthebarrel *********** " << endl;
0078         if (seedfromoverlaps) {
0079           //this part of code works for MTCC
0080           // for the other geometries must be verified
0081           //Overlaps in the difference in z is decreased and the difference in phi is
0082           //less than 0.05
0083           if ((DeltaR < 0) && (abs(z_diff) < 18) &&
0084               (abs(ih->globalPosition().phi() - oh->globalPosition().phi()) < 0.05) && (dxdy < 2))
0085             result.push_back(OrderedHitPair(ih, oh));
0086         } else
0087           result.push_back(OrderedHitPair(ih, oh));
0088       }
0089       if (InTheForward && (abs(z_diff) > 1.)) {
0090         //  cout << " ******** sono dentro intheforward *********** " << endl;
0091         result.push_back(OrderedHitPair(ih, oh));
0092       }
0093     }
0094   }
0095
0096   /*
0097   // uncomment if the sort operation below needs to be called
0098   class CompareHitPairsY {
0099   public:
0100     CompareHitPairsY(const TrackerGeometry& t): tracker{&t} {}
0101     bool operator()(const OrderedHitPair& h1, const OrderedHitPair& h2) {
0102       const TrackingRecHit* trh1i = h1.inner()->hit();
0103       const TrackingRecHit* trh2i = h2.inner()->hit();
0104       const TrackingRecHit* trh1o = h1.outer()->hit();
0105       const TrackingRecHit* trh2o = h2.outer()->hit();
0106       GlobalPoint in1p = tracker->idToDet(trh1i->geographicalId())->surface().toGlobal(trh1i->localPosition());
0107       GlobalPoint in2p = tracker->idToDet(trh2i->geographicalId())->surface().toGlobal(trh2i->localPosition());
0108       GlobalPoint ou1p = tracker->idToDet(trh1o->geographicalId())->surface().toGlobal(trh1o->localPosition());
0109       GlobalPoint ou2p = tracker->idToDet(trh2o->geographicalId())->surface().toGlobal(trh2o->localPosition());
0110       if (ou1p.y() * ou2p.y() < 0)
0111         return ou1p.y() > ou2p.y();
0112       else {
0113         float dist1 = 100 * std::abs(ou1p.z() - in1p.z()) - std::abs(ou1p.y()) - 0.1 * std::abs(in1p.y());
0114         float dist2 = 100 * std::abs(ou2p.z() - in2p.z()) - std::abs(ou2p.y()) - 0.1 * std::abs(in2p.y());
0115         return dist1 < dist2;
0116       }
0117     }
0118
0119   private:
0120     const TrackerGeometry* tracker;
0121   };
0122   stable_sort(allthepairs.begin(),allthepairs.end(),CompareHitPairsY(*trackerGeometry));
0123   //Seed from overlaps are saved only if
0124   //no others have been saved
0125
0126   if (allthepairs.size()>0) {
0127     if (seedfromoverlaps) {
0128       if (result.size()==0) result.push_back(allthepairs[0]);
0129     }
0130     else result.push_back(allthepairs[0]);
0131   }
0132 */
0133 }