#include <HitPairGeneratorFromLayerPairForPhotonConversion.h>

Public Types
typedef SeedingLayerSetsHits::SeedingLayer	Layer

typedef LayerHitMapCache	LayerCacheType

typedef SeedingLayerSetsHits::SeedingLayerSet	Layers

Public Member Functions
bool	checkBoundaries (const DetLayer &layer, const ConversionRegion &convRegion, float maxSearchR, float maxSearchZ)

bool	checkRZCompatibilityWithSeedTrack (const RecHitsSortedInPhi::Hit &hit, const DetLayer &layer, const ConversionRegion &convRegion)

float	getLayerRadius (const DetLayer &layer)

float	getLayerZ (const DetLayer &layer)

bool	getPhiRange (float &Phimin, float &Phimax, const DetLayer &layer, const ConversionRegion &convRegion, const edm::EventSetup &es)

bool	getPhiRange (float &Phimin, float &Phimax, const float &layerR, const ConversionRegion &convRegion, const edm::EventSetup &es)

	HitPairGeneratorFromLayerPairForPhotonConversion (unsigned int inner, unsigned int outer, LayerCacheType *layerCache, unsigned int nSize=30000, unsigned int max=0)

void	hitPairs (const ConversionRegion &convRegion, const TrackingRegion &reg, OrderedHitPairs &prs, const edm::Event &ev, const edm::EventSetup &es)

Layer	innerLayer () const

Layer	outerLayer () const

void	setSeedingLayers (Layers layers)

Private Member Functions
double	getCot (double dz, double dr)

Private Attributes
std::stringstream	ss

const unsigned int	theInnerLayer

LayerCacheType &	theLayerCache

const unsigned int	theMaxElement

const unsigned int	theOuterLayer

Layers	theSeedingLayers

Detailed Description

Definition at line 13 of file HitPairGeneratorFromLayerPairForPhotonConversion.h.

Member Typedef Documentation

typedef SeedingLayerSetsHits::SeedingLayer HitPairGeneratorFromLayerPairForPhotonConversion::Layer

Definition at line 19 of file HitPairGeneratorFromLayerPairForPhotonConversion.h.

typedef LayerHitMapCache HitPairGeneratorFromLayerPairForPhotonConversion::LayerCacheType

Definition at line 17 of file HitPairGeneratorFromLayerPairForPhotonConversion.h.

typedef SeedingLayerSetsHits::SeedingLayerSet HitPairGeneratorFromLayerPairForPhotonConversion::Layers

Definition at line 18 of file HitPairGeneratorFromLayerPairForPhotonConversion.h.

Constructor & Destructor Documentation

HitPairGeneratorFromLayerPairForPhotonConversion::HitPairGeneratorFromLayerPairForPhotonConversion	(	unsigned int	inner,
		unsigned int	outer,
		LayerCacheType *	layerCache,
		unsigned int	nSize = `30000`,
		unsigned int	max = `0`
	)

Definition at line 37 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

   :
     theLayerCache(*layerCache), theOuterLayer(outer), theInnerLayer(inner), theMaxElement(max) {
 
 }

Member Function Documentation

bool HitPairGeneratorFromLayerPairForPhotonConversion::checkBoundaries	(	const DetLayer &	layer,
		const ConversionRegion &	convRegion,
		float	maxSearchR,
		float	maxSearchZ
	)

Definition at line 224 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References GeomDetEnumerators::barrel, ConversionRegion::convPoint(), gather_cfg::cout, GeomDetEnumerators::endcap, getLayerRadius(), getLayerZ(), DetLayer::location(), PV3DBase< T, PVType, FrameType >::perp(), dttmaxenums::R, ss, Gflash::Z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by hitPairs().

                                                                                                                                                                  {
 
   if(layer.location() == GeomDetEnumerators::barrel){
 
     float minZEndCap=130;
     if(fabs(convRegion.convPoint().z()) > minZEndCap){
 #ifdef mydebug_Seed
       ss << "\tthe conversion seems to be in the endcap. Zconv " << convRegion.convPoint().z() << std::endl;
       std::cout << ss.str();
 #endif
       return false;
     }
       
     float R=getLayerRadius(layer);
     
     if(convRegion.convPoint().perp()>R){
 #ifdef mydebug_Seed
       ss << "\tthis layer is before the conversion : R layer " << R << " [ Rconv " << convRegion.convPoint().perp() << " Zconv " << convRegion.convPoint().z()<< std::endl;
       std::cout << ss.str();
 #endif
       return false;
     }
     
     if(R - convRegion.convPoint().perp() > maxSearchR ){
 #ifdef mydebug_Seed
       ss << "\tthis layer is far from the conversion more than cut " << maxSearchR << " cm. R layer " << R << " [ Rconv " << convRegion.convPoint().perp() << " Zconv " << convRegion.convPoint().z()<< std::endl;
       std::cout << ss.str();
 #endif
       return false;
     }
     
   }else if (layer.location() == GeomDetEnumerators::endcap){
     
     float Z=getLayerZ(layer);
     if(
        (convRegion.convPoint().z()>0 && convRegion.convPoint().z()>Z)
        ||
        (convRegion.convPoint().z()<0 && convRegion.convPoint().z()<Z)
        ) {
 #ifdef mydebug_Seed
       ss << "\tthis layer is before the conversion : Z layer " << Z << " [ Rconv " << convRegion.convPoint().perp()<< " Zconv " << convRegion.convPoint().z() << std::endl;
       std::cout << ss.str();
 #endif
       return false;
     }    
 
      
     if(fabs(Z - convRegion.convPoint().z()) > maxSearchZ ){
 #ifdef mydebug_Seed
       ss << "\tthis layer is far from the conversion more than cut " << maxSearchZ << " cm. Z layer " << Z << " [ Rconv " << convRegion.convPoint().perp()<< " Zconv " << convRegion.convPoint().z() << std::endl;
       std::cout << ss.str();
 #endif
       return false;
     }
     
   }
   return true;
 }

bool HitPairGeneratorFromLayerPairForPhotonConversion::checkRZCompatibilityWithSeedTrack	(	const RecHitsSortedInPhi::Hit &	hit,
		const DetLayer &	layer,
		const ConversionRegion &	convRegion
	)

Definition at line 326 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References GeomDetEnumerators::barrel, ConversionRegion::convPoint(), ConversionRegion::cotTheta(), ConversionRegion::errTheta(), f, getCot(), PixelRecoRange< T >::intersection(), DetLayer::location(), nSigmaRZ, ConversionRegion::pvtxPoint(), funct::sqr(), mathSSE::sqrt(), ss, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by hitPairs().

                                                                                                                             {
   static const float nSigmaRZ = std::sqrt(12.f);
   Range hitCotTheta;
   
   double sigmaCotTheta = convRegion.errTheta() * (1+convRegion.cotTheta()*convRegion.cotTheta()); //Error Propagation from sigma theta.
   Range allowedCotTheta(convRegion.cotTheta()-nSigmaRZ*sigmaCotTheta,convRegion.cotTheta()+nSigmaRZ*sigmaCotTheta);
   
   double dz        = hit->globalPosition().z()-convRegion.pvtxPoint().z();
   double r_reduced = std::sqrt( sqr(hit->globalPosition().x()-convRegion.pvtxPoint().x())+sqr(hit->globalPosition().y()-convRegion.pvtxPoint().y()));
 
   if (layer.location() == GeomDetEnumerators::barrel){ 
     float zErr = nSigmaRZ * hit->errorGlobalZ();
     hitCotTheta = Range(getCot(dz-zErr,r_reduced),getCot(dz+zErr,r_reduced));
   }else{
     float rErr = nSigmaRZ * hit->errorGlobalR();
     if(dz>0)  
       hitCotTheta = Range(getCot(dz,r_reduced+rErr),getCot(dz,r_reduced-rErr));
     else
       hitCotTheta = Range(getCot(dz,r_reduced-rErr), getCot(dz,r_reduced+rErr));
   }
 
   Range crossRange = allowedCotTheta.intersection(hitCotTheta);
 
 #ifdef mydebug_Seed      
       ss 
     << "\n\t\t cotTheta allowed Range " << allowedCotTheta.min() << " \t, " << allowedCotTheta.max() 
     << "\n\t\t hitCotTheta   Range " << hitCotTheta.min()   << " \t, " << hitCotTheta.max() 
     << "\n\t\t Cross   Range " << crossRange.min()   << " \t, " << crossRange.max() 
     << "\n\t\t the seed track has origin " << convRegion.convPoint() << " \t cotTheta " << convRegion.cotTheta()  
     << std::endl;
 #endif
 
   return crossRange.empty();
 }

double HitPairGeneratorFromLayerPairForPhotonConversion::getCot	(	double	dz,
		double	dr
	)

private

Definition at line 363 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References funct::abs(), alignCSCRings::e, and f.

Referenced by checkRZCompatibilityWithSeedTrack().

                             {
   if ( std::abs(dr) > 1.e-4f )    return dz/dr;
   else
     if(dz>0) return  99999.f; 
     else     return -99999.f; 
 }

float HitPairGeneratorFromLayerPairForPhotonConversion::getLayerRadius ( const DetLayer & layer )

Definition at line 199 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References GeomDetEnumerators::barrel, Surface::bounds(), DetLayer::location(), BarrelDetLayer::specificSurface(), GeometricSearchDet::surface(), and Bounds::thickness().

Referenced by checkBoundaries(), and getPhiRange().

                                                                                            {
   if (layer.location() == GeomDetEnumerators::barrel){ 
     const BarrelDetLayer& bl = static_cast<const BarrelDetLayer&>(layer); 
     float rLayer = bl.specificSurface().radius(); 
     
     // the maximal delta phi will be for the innermost hits
     float theThickness = layer.surface().bounds().thickness();
     return rLayer + 0.5f*theThickness;
   }
 
   //Fixme
   return 0;
 }

float HitPairGeneratorFromLayerPairForPhotonConversion::getLayerZ ( const DetLayer & layer )

Definition at line 213 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References Surface::bounds(), GeomDetEnumerators::endcap, DetLayer::location(), GeometricSearchDet::position(), GeometricSearchDet::surface(), Bounds::thickness(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by checkBoundaries(), and getPhiRange().

                                                                                       {
   if (layer.location() == GeomDetEnumerators::endcap){ 
     float layerZ = layer.position().z();
     float theThickness = layer.surface().bounds().thickness();
     float layerZmax = layerZ > 0 ? layerZ+0.5f*theThickness: layerZ-0.5f*theThickness;
     return layerZmax;
   }else{
     //Fixme
     return 0;
   }
 }

bool HitPairGeneratorFromLayerPairForPhotonConversion::getPhiRange	(	float &	Phimin,
		float &	Phimax,
		const DetLayer &	layer,
		const ConversionRegion &	convRegion,
		const edm::EventSetup &	es
	)

Definition at line 283 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References GeomDetEnumerators::barrel, ConversionRegion::cotTheta(), GeomDetEnumerators::endcap, getLayerRadius(), getLayerZ(), DetLayer::location(), dttmaxenums::R, and Gflash::Z.

Referenced by hitPairs().

                                                                                                                                                                                   {
    if(layer.location() == GeomDetEnumerators::barrel){
      return getPhiRange(Phimin,Phimax,getLayerRadius(layer),convRegion,es);
    }else if (layer.location() == GeomDetEnumerators::endcap){
      float Z=getLayerZ(layer);
      float R=Z/convRegion.cotTheta();
      return getPhiRange(Phimin,Phimax,R,convRegion,es); //FIXME
    }
    return false;
 }

bool HitPairGeneratorFromLayerPairForPhotonConversion::getPhiRange	(	float &	Phimin,
		float &	Phimax,
		const float &	layerR,
		const ConversionRegion &	convRegion,
		const edm::EventSetup &	es
	)

Definition at line 294 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References PixelRecoUtilities::bendingRadius(), ConversionRegion::charge(), ConversionRegion::convPoint(), reco::deltaPhi(), PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), RecoTau_DiTaus_pt_20-420_cfg::Phimin, ptmin, funct::sqr(), and mathSSE::sqrt().

                                                                                                                                                                                 {
   Phimin = reco::deltaPhi(convRegion.convPoint().phi(),0.);
   
   float dphi;
   float ptmin=0.1;
   float DeltaL=layerR-convRegion.convPoint().perp();
 
   if(DeltaL<0){
     Phimin=0;
     Phimax=0;
     return false;
   }
 
   float theRCurvatureMin = PixelRecoUtilities::bendingRadius(ptmin,es);
 
   if(theRCurvatureMin<DeltaL)
     dphi = atan(DeltaL/layerR);
   else
     dphi = atan(theRCurvatureMin/layerR * ( 1 - sqrt(1-sqr(DeltaL/theRCurvatureMin)) ) );
   
   if(convRegion.charge()>0){
     Phimax=Phimin;
     Phimin=Phimax-dphi;
   }else{
     Phimax=Phimin+dphi;
   }
 
   //std::cout << dphi << " " << Phimin << " " << Phimax << " " << layerR << " " << DeltaL  << " " << convRegion.convPoint().phi() << " " << convRegion.convPoint().perp()<< std::endl;
   return true;
 }  

void HitPairGeneratorFromLayerPairForPhotonConversion::hitPairs	(	const ConversionRegion &	convRegion,
		const TrackingRegion &	reg,
		OrderedHitPairs &	prs,
		const edm::Event &	ev,
		const edm::EventSetup &	es
	)

Definition at line 48 of file HitPairGeneratorFromLayerPairForPhotonConversion.cc.

References Reference_intrackfit_cff::barrel, checkBoundaries(), TrackingRegion::checkRZ(), checkRZCompatibilityWithSeedTrack(), ConversionRegion::convPoint(), ConversionRegion::cotTheta(), gather_cfg::cout, SeedingLayerSetsHits::SeedingLayer::detLayer(), RecHitsSortedInPhi::empty(), f, getPhiRange(), RecHitsSortedInPhi::hits(), innerLayer(), PixelRecoRange< T >::intersection(), DetLayer::location(), SeedingLayerSetsHits::SeedingLayer::name(), nSigmaRZ, TrackingRegion::origin(), outerLayer(), PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), HitRZCompatibility::range(), OrderedHitPairs::size(), funct::sqr(), mathSSE::sqrt(), ss, theLayerCache, theMaxElement, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

 {
  auto oldSize = result.size();
  auto maxNum = result.size()+theMaxElement;
 
 #ifdef mydebug_Seed
   ss.str("");
 #endif
 
   typedef OrderedHitPair::InnerRecHit InnerHit;
   typedef OrderedHitPair::OuterRecHit OuterHit;
   typedef RecHitsSortedInPhi::Hit Hit;
 
   Layer innerLayerObj = innerLayer();
   Layer outerLayerObj = outerLayer();
 
 #ifdef mydebug_Seed
   ss << "In " << innerLayerObj.name() << " Out " << outerLayerObj.name() << std::endl;
 #endif
 
   if(!checkBoundaries(*innerLayerObj.detLayer(),convRegion,40.,60.)) return; //FIXME, the maxSearchR(Z) are not optimized
   if(!checkBoundaries(*outerLayerObj.detLayer(),convRegion,50.,60.)) return; //FIXME, the maxSearchR(Z) are not optimized
 
   /*get hit sorted in phi for each layer: NB: doesn't apply any region cut*/
   const RecHitsSortedInPhi & innerHitsMap = theLayerCache(innerLayerObj, region, event, es);
   if (innerHitsMap.empty()) return;
  
   const RecHitsSortedInPhi& outerHitsMap = theLayerCache(outerLayerObj, region, event, es);
   if (outerHitsMap.empty()) return;
   /*----------------*/
 
   /*This object will check the compatibility of the his in phi among the two layers. */
   //InnerDeltaPhi deltaPhi(innerLayerObj.detLayer(), region, es);
 
   static const float nSigmaRZ = std::sqrt(12.f);
   //  static const float nSigmaPhi = 3.f;
   vector<RecHitsSortedInPhi::Hit> innerHits, outerHits;
   float outerPhimin, outerPhimax;
   float innerPhimin, innerPhimax;
 
   /*Getting only the Hits in the outer layer that are compatible with the conversion region*/
   if(!getPhiRange(outerPhimin,outerPhimax,*outerLayerObj.detLayer(),convRegion,es)) return;
   outerHitsMap.hits( outerPhimin, outerPhimax, outerHits);
 
 #ifdef mydebug_Seed
   ss << "\tophimin, ophimax " << outerPhimin << " " << outerPhimax << std::endl;
 #endif
 
   /* loop on outer hits*/
   for (vector<RecHitsSortedInPhi::Hit>::const_iterator oh = outerHits.begin(); oh!= outerHits.end(); ++oh) { 
     RecHitsSortedInPhi::Hit ohit = (*oh);
 #ifdef mydebug_Seed
     GlobalPoint oPos = ohit->globalPosition();  
     
     ss << "\toPos " << oPos << " r " << oPos.perp() << " phi " << oPos.phi() <<  " cotTheta " << oPos.z()/oPos.perp() << std::endl;
 #endif
 
     /*Check the compatibility of the ohit with the eta of the seeding track*/
     if(checkRZCompatibilityWithSeedTrack(ohit,*outerLayerObj.detLayer(),convRegion)) continue;
 
     /*  
     //Do I need this? it uses a compatibility that probably I wouldn't 
     //Removing for the time being
 
     PixelRecoRange<float> phiRange = deltaPhi( oPos.perp(), oPos.phi(), oPos.z(), nSigmaPhi*(ohit->errorGlobalRPhi()));    
     if (phiRange.empty()) continue;
     */
 
     const HitRZCompatibility *checkRZ = region.checkRZ(innerLayerObj.detLayer(), ohit, es);
     if(!checkRZ) {
 #ifdef mydebug_Seed
       ss << "*******\nNo valid checkRZ\n*******" << std::endl;
 #endif
       continue;
     }
     
     /*Get only the inner hits compatible with the conversion region*/
     innerHits.clear();
     if(!getPhiRange(innerPhimin,innerPhimax,*innerLayerObj.detLayer(),convRegion,es)) continue;
     innerHitsMap.hits(innerPhimin, innerPhimax, innerHits);
 
 #ifdef mydebug_Seed
     ss << "\tiphimin, iphimax " << innerPhimin << " " << innerPhimax << std::endl;
 #endif    
 
     /*Loop on inner hits*/
     for ( vector<RecHitsSortedInPhi::Hit>::const_iterator ih=innerHits.begin(), ieh = innerHits.end(); ih < ieh; ++ih) {  
       GlobalPoint innPos = (*ih)->globalPosition();
       
 
 #ifdef mydebug_Seed
       ss << "\tinnPos " << innPos <<  " r " << innPos.perp() << " phi " << innPos.phi() << " cotTheta " << innPos.z()/innPos.perp() <<  std::endl;
 #endif
 
       /*Check the compatibility of the ohit with the eta of the seeding track*/
       if(checkRZCompatibilityWithSeedTrack(*ih,*innerLayerObj.detLayer(),convRegion)) continue;
       
       float r_reduced = std::sqrt( sqr(innPos.x()-region.origin().x())+sqr(innPos.y()-region.origin().y()));
       Range allowed;
       Range hitRZ;
       if (innerLayerObj.detLayer()->location() == barrel) {
         allowed = checkRZ->range(r_reduced);
         float zErr = nSigmaRZ * (*ih)->errorGlobalZ();
         hitRZ = Range(innPos.z()-zErr, innPos.z()+zErr);
       } else {
         allowed = checkRZ->range(innPos.z());
         float rErr = nSigmaRZ * (*ih)->errorGlobalR();
         hitRZ = Range(r_reduced-rErr, r_reduced+rErr);
       }
       Range crossRange = allowed.intersection(hitRZ);
 
 #ifdef mydebug_Seed      
       ss 
     << "\n\t\t allowed Range " << allowed.min() << " \t, " << allowed.max() 
     << "\n\t\t hitRz   Range " << hitRZ.min()   << " \t, " << hitRZ.max() 
     << "\n\t\t Cross   Range " << crossRange.min()   << " \t, " << crossRange.max() 
     << "\n\t\t the seed track has origin " << convRegion.convPoint() << " \t cotTheta " << convRegion.cotTheta() 
     << std::endl;
 #endif
       
       if (! crossRange.empty() ) {
 #ifdef mydebug_Seed      
       ss 
     << "\n\t\t !!!!ACCEPTED!!! \n\n";
 #endif
     if (theMaxElement!=0 && result.size() >= maxNum ){
       result.resize(oldSize);
 #ifdef mydebug_Seed      
       edm::LogError("TooManySeeds")<<"number of pairs exceed maximum, no pairs produced";
 #endif
       delete checkRZ;
 
 #ifdef mydebug_Seed
       std::cout << ss.str();
 #endif
       return;
     }
         result.push_back( OrderedHitPair( *ih, ohit) );
       }
     }
     delete checkRZ;
   }
 
 #ifdef mydebug_Seed
   std::cout << ss.str();
 #endif
 }