#include <SiStripElectronSeedGenerator.h>

Public Types
typedef TransientTrackingRecHit::ConstRecHitPointer	ConstRecHitPointer
typedef edm::OwnVector < TrackingRecHit >	PRecHitContainer
typedef TransientTrackingRecHit::RecHitContainer	RecHitContainer
typedef TransientTrackingRecHit::RecHitPointer	RecHitPointer
Public Member Functions
void	run (edm::Event &, const edm::EventSetup &setup, const edm::Handle< reco::SuperClusterCollection > &, reco::ElectronSeedCollection &)
void	setupES (const edm::EventSetup &setup)
	SiStripElectronSeedGenerator (const edm::ParameterSet &)
	~SiStripElectronSeedGenerator ()
Private Member Functions
bool	altCheckHitsAndTSOS (std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit1, std::vector< const SiStripRecHit2D * >::const_iterator hit2, double scr, double scz, double pT, double scEta)
const SiStripRecHit2D *	backupHitConverter (ConstRecHitPointer crhp)
bool	checkHitsAndTSOS (std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit1, std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit2, double scr, double scz, double pT, double scEta)
void	findSeedsFromCluster (edm::Ref< reco::SuperClusterCollection >, edm::Handle< reco::BeamSpot >, reco::ElectronSeedCollection &)
const SiStripMatchedRecHit2D *	matchedHitConverter (ConstRecHitPointer crhp)
double	normalPhi (double phi) const
double	phiDiff (double phi1, double phi2)
bool	preselection (GlobalPoint position, GlobalPoint superCluster, double phiVsRSlope, int hitLayer)
double	unwrapPhi (double phi) const
std::vector< bool >	useDetLayer (double scEta)
int	whichSubdetector (std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit)
Private Attributes
std::vector< const SiStripRecHit2D * >	backupLayer2Hits_
edm::InputTag	beamSpotTag_
unsigned long long	cacheIDCkfComp_
unsigned long long	cacheIDMagField_
unsigned long long	cacheIDTrkGeom_
std::vector< const SiStripMatchedRecHit2D * >	layer1Hits_
std::vector< const SiStripMatchedRecHit2D * >	layer2Hits_
int	maxSeeds_
edm::ESHandle< MeasurementTracker >	measurementTrackerHandle
double	monoDeltaPsiCut_
int	monoMaxHits_
double	monoOriginZCut_
double	monoPhiMissHit2Cut_
PTrajectoryStateOnDet	pts_
PRecHitContainer	recHits_
double	tecDeltaPsiCut_
int	tecMaxHits_
double	tecOriginZCut_
double	tecPhiMissHit2Cut_
double	tecRMissHit2Cut_
edm::Handle< reco::BeamSpot >	theBeamSpot
Chi2MeasurementEstimator *	theEstimator
edm::ESHandle< MagneticField >	theMagField
const SiStripRecHitMatcher *	theMatcher_
const MeasurementTracker *	theMeasurementTracker
std::string	theMeasurementTrackerName
PropagatorWithMaterial *	thePropagator
const edm::EventSetup *	theSetup
KFUpdator *	theUpdator
double	tibDeltaPsiCut_
double	tibOriginZCut_
double	tibPhiMissHit2Cut_
double	tibZMissHit2Cut_
double	tidDeltaPsiCut_
double	tidEtaUsage_
int	tidMaxHits_
double	tidOriginZCut_
double	tidPhiMissHit2Cut_
double	tidRMissHit2Cut_
edm::ESHandle< TrackerGeometry >	trackerGeometryHandle

Detailed Description

Class to generate the trajectory seed from two Si Strip hits.

Author:: Chris Macklin, Avishek Chatterjee

Version:: March 2009 (Adapt code to simplify call to SetupES)

Description: SiStrip-driven electron seed finding algorithm.

Definition at line 64 of file SiStripElectronSeedGenerator.h.

Member Typedef Documentation

typedef TransientTrackingRecHit::ConstRecHitPointer SiStripElectronSeedGenerator::ConstRecHitPointer

Definition at line 69 of file SiStripElectronSeedGenerator.h.

typedef edm::OwnVector<TrackingRecHit> SiStripElectronSeedGenerator::PRecHitContainer

Definition at line 68 of file SiStripElectronSeedGenerator.h.

typedef TransientTrackingRecHit::RecHitContainer SiStripElectronSeedGenerator::RecHitContainer

Definition at line 71 of file SiStripElectronSeedGenerator.h.

typedef TransientTrackingRecHit::RecHitPointer SiStripElectronSeedGenerator::RecHitPointer

Definition at line 70 of file SiStripElectronSeedGenerator.h.

Constructor & Destructor Documentation

SiStripElectronSeedGenerator::SiStripElectronSeedGenerator ( const edm::ParameterSet & pset )

Definition at line 49 of file SiStripElectronSeedGenerator.cc.

References beamSpotTag_, Chi2MeasurementEstimatorESProducer_cfi::Chi2MeasurementEstimator, edm::ParameterSet::exists(), edm::ParameterSet::getParameter(), theEstimator, theMeasurementTrackerName, and theUpdator.

 : beamSpotTag_("offlineBeamSpot"),
   theUpdator(0),thePropagator(0),theMeasurementTracker(0),
   theSetup(0), theMatcher_(0),
   cacheIDMagField_(0),cacheIDCkfComp_(0),cacheIDTrkGeom_(0),
   tibOriginZCut_(pset.getParameter<double>("tibOriginZCut")),
   tidOriginZCut_(pset.getParameter<double>("tidOriginZCut")),
   tecOriginZCut_(pset.getParameter<double>("tecOriginZCut")),
   monoOriginZCut_(pset.getParameter<double>("monoOriginZCut")),
   tibDeltaPsiCut_(pset.getParameter<double>("tibDeltaPsiCut")),
   tidDeltaPsiCut_(pset.getParameter<double>("tidDeltaPsiCut")),
   tecDeltaPsiCut_(pset.getParameter<double>("tecDeltaPsiCut")),
   monoDeltaPsiCut_(pset.getParameter<double>("monoDeltaPsiCut")),
   tibPhiMissHit2Cut_(pset.getParameter<double>("tibPhiMissHit2Cut")),
   tidPhiMissHit2Cut_(pset.getParameter<double>("tidPhiMissHit2Cut")),
   tecPhiMissHit2Cut_(pset.getParameter<double>("tecPhiMissHit2Cut")),
   monoPhiMissHit2Cut_(pset.getParameter<double>("monoPhiMissHit2Cut")),
   tibZMissHit2Cut_(pset.getParameter<double>("tibZMissHit2Cut")),
   tidRMissHit2Cut_(pset.getParameter<double>("tidRMissHit2Cut")),
   tecRMissHit2Cut_(pset.getParameter<double>("tecRMissHit2Cut")),
   tidEtaUsage_(pset.getParameter<double>("tidEtaUsage")),
   tidMaxHits_(pset.getParameter<int>("tidMaxHits")),
   tecMaxHits_(pset.getParameter<int>("tecMaxHits")),
   monoMaxHits_(pset.getParameter<int>("monoMaxHits")),
   maxSeeds_(pset.getParameter<int>("maxSeeds"))
{
  // use of a theMeasurementTrackerName
  if (pset.exists("measurementTrackerName"))
   { theMeasurementTrackerName = pset.getParameter<std::string>("measurementTrackerName") ; }

  // new beamSpot tag
  if (pset.exists("beamSpot"))
   { beamSpotTag_ = pset.getParameter<edm::InputTag>("beamSpot") ; }

  theUpdator = new KFUpdator();
  theEstimator = new Chi2MeasurementEstimator(30,3);
}

SiStripElectronSeedGenerator::~SiStripElectronSeedGenerator ( )

Definition at line 88 of file SiStripElectronSeedGenerator.cc.

References thePropagator, and theUpdator.

                                                            {
  delete thePropagator;
  delete theUpdator;
}

Member Function Documentation

bool SiStripElectronSeedGenerator::altCheckHitsAndTSOS	(	std::vector< const SiStripMatchedRecHit2D * >::const_iterator	hit1,
		std::vector< const SiStripRecHit2D * >::const_iterator	hit2,
		double	scr,
		double	scz,
		double	pT,
		double	scEta
	)		`[private]`

Definition at line 607 of file SiStripElectronSeedGenerator.cc.

References a, abs, b, newFWLiteAna::build, FastHelix::isValid(), monoPhiMissHit2Cut_, trajectoryStateTransform::persistentState(), phiDiff(), PropagatorWithMaterial::propagate(), pts_, diffTwoXMLs::r1, diffTwoXMLs::r2, mathSSE::sqrt(), FastHelix::stateAtVertex(), theMatcher_, thePropagator, theSetup, theUpdator, trackerGeometryHandle, and KFUpdator::update().

                                                                                                   {

  bool seedCutSatisfied = false;

  using namespace std;

  GlobalPoint hit1Pos = trackerGeometryHandle->idToDet((*hit1)->geographicalId())->surface().toGlobal((*hit1)->localPosition());
  double r1 = sqrt(hit1Pos.x()*hit1Pos.x() + hit1Pos.y()*hit1Pos.y());
  double phi1 = hit1Pos.phi();
  double z1=hit1Pos.z();

  GlobalPoint hit2Pos = trackerGeometryHandle->idToDet((*hit2)->geographicalId())->surface().toGlobal((*hit2)->localPosition());
  double r2 = sqrt(hit2Pos.x()*hit2Pos.x() + hit2Pos.y()*hit2Pos.y());
  double phi2 = hit2Pos.phi();
  double z2 = hit2Pos.z();

  if(r2 > r1 && std::abs(z2) > std::abs(z1)) {

    //Consider the circle made of IP and Hit 1; Calculate it's radius using pT

    double curv = pT*100*.877;

    //Predict phi of hit 2
    double a = (r2-r1)/(2*curv);
    double b = phiDiff(phi2,phi1);
    double phiMissHit2 = 0;
    if(std::abs(b - a)<std::abs(b + a)) phiMissHit2 = b - a;
    if(std::abs(b - a)>std::abs(b + a)) phiMissHit2 = b + a;

    if(std::abs(phiMissHit2) < monoPhiMissHit2Cut_) seedCutSatisfied = true;

  }

  if(!seedCutSatisfied) return false;

  // seed checks borrowed from pixel-based algoritm

 

  /* Some of this code could be better optimized.  The Pixel algorithm natively
     takes Transient rec hits, so to recycle code we have to build them.
  */

  RecHitPointer hit1Trans = TSiStripMatchedRecHit::build(trackerGeometryHandle->idToDet((*hit1)->geographicalId()), *hit1, theMatcher_);
  RecHitPointer hit2Trans = TSiStripMatchedRecHit::build(trackerGeometryHandle->idToDet((*hit2)->geographicalId()), *hit2, theMatcher_);

  typedef TrajectoryStateOnSurface TSOS;

  double vertexZ = z1 - (r1 * (zc - z1) ) / (rc - r1);
  GlobalPoint eleVertex(0.,0.,vertexZ);

  // make a spiral
  FastHelix helix(hit2Pos,hit1Pos,eleVertex,*theSetup);
  if (!helix.isValid()) return false;

  FreeTrajectoryState fts = helix.stateAtVertex();
  TSOS propagatedState = thePropagator->propagate(fts,hit1Trans->det()->surface());

  if (!propagatedState.isValid()) return false;

  TSOS updatedState = theUpdator->update(propagatedState, *hit1Trans);
  TSOS propagatedState_out = thePropagator->propagate(fts,hit2Trans->det()->surface()) ;

  if (!propagatedState_out.isValid()) return false;

  // the seed has now passed all the cuts

  TSOS updatedState_out = theUpdator->update(propagatedState_out, *hit2Trans);

  pts_ =  trajectoryStateTransform::persistentState(updatedState_out, hit2Trans->geographicalId().rawId());

  return true;
}

const SiStripRecHit2D * SiStripElectronSeedGenerator::backupHitConverter ( ConstRecHitPointer crhp ) [private]

Definition at line 737 of file SiStripElectronSeedGenerator.cc.

                                                                                              {
  const TrackingRecHit* trh = crhp->hit();
  const SiStripRecHit2D* backupHit = dynamic_cast<const SiStripRecHit2D*>(trh);
  return backupHit;
}

bool SiStripElectronSeedGenerator::checkHitsAndTSOS	(	std::vector< const SiStripMatchedRecHit2D * >::const_iterator	hit1,
		std::vector< const SiStripMatchedRecHit2D * >::const_iterator	hit2,
		double	scr,
		double	scz,
		double	pT,
		double	scEta
	)		`[private]`

Definition at line 505 of file SiStripElectronSeedGenerator.cc.

References a, abs, b, newFWLiteAna::build, FastHelix::isValid(), trajectoryStateTransform::persistentState(), phiDiff(), PropagatorWithMaterial::propagate(), pts_, diffTwoXMLs::r1, diffTwoXMLs::r2, mathSSE::sqrt(), FastHelix::stateAtVertex(), tecPhiMissHit2Cut_, tecRMissHit2Cut_, theMatcher_, thePropagator, theSetup, theUpdator, tibPhiMissHit2Cut_, tibZMissHit2Cut_, tidPhiMissHit2Cut_, tidRMissHit2Cut_, trackerGeometryHandle, KFUpdator::update(), and whichSubdetector().

                                                                                                {

  bool seedCutsSatisfied = false;

  using namespace std;

  GlobalPoint hit1Pos = trackerGeometryHandle->idToDet((*hit1)->geographicalId())->surface().toGlobal((*hit1)->localPosition());
  double r1 = sqrt(hit1Pos.x()*hit1Pos.x() + hit1Pos.y()*hit1Pos.y());
  double phi1 = hit1Pos.phi();
  double z1=hit1Pos.z();

  GlobalPoint hit2Pos = trackerGeometryHandle->idToDet((*hit2)->geographicalId())->surface().toGlobal((*hit2)->localPosition());
  double r2 = sqrt(hit2Pos.x()*hit2Pos.x() + hit2Pos.y()*hit2Pos.y());
  double phi2 = hit2Pos.phi();
  double z2 = hit2Pos.z();

  if(r2 > r1 && (std::abs(z2) > std::abs(z1) || std::abs(scEta) < 0.25)) {

    //Consider the circle made of IP and Hit 1; Calculate it's radius using pT

    double curv = pT*100*.877;

    //Predict phi of hit 2
    double a = (r2-r1)/(2*curv);
    double b = phiDiff(phi2,phi1);
    //UB added '=0' to avoid compiler warning
    double phiMissHit2=0;
    if(std::abs(b - a)<std::abs(b + a)) phiMissHit2 = b - a;
    if(std::abs(b - a)>std::abs(b + a)) phiMissHit2 = b + a;

    double zMissHit2 = z2 - (r2*(zc-z1)-r1*zc+rc*z1)/(rc-r1);

    double rPredHit2 = r1 + (rc-r1)/(zc-z1)*(z2-z1);
    double rMissHit2 = r2 - rPredHit2;

    int subdetector = whichSubdetector(hit2);

    bool zDiff = true;
    double zVar1 = std::abs(z1);
    double zVar2 = std::abs(z2 - z1);
    if(zVar1 > 75 && zVar1 < 95 && (zVar2 > 18 || zVar2 < 5)) zDiff = false;
    if(zVar1 > 100 && zVar1 < 110 && (zVar2 > 35 || zVar2 < 5)) zDiff = false;
    if(zVar1 > 125 && zVar1 < 150 && (zVar2 > 18 || zVar2 < 5)) zDiff = false;

    if(subdetector == 1){
      int tibExtraCut = 0;
      if(r1 > 23 && r1 < 28 && r2 > 31 && r2 < 37) tibExtraCut = 1;
      if(std::abs(phiMissHit2) < tibPhiMissHit2Cut_ && std::abs(zMissHit2) < tibZMissHit2Cut_ && tibExtraCut == 1) seedCutsSatisfied = true;
    }else if(subdetector == 2){
      int tidExtraCut = 0;
      if(r1 > 23 && r1 < 34 && r2 > 26 && r2 < 42) tidExtraCut = 1;
      if(std::abs(phiMissHit2) < tidPhiMissHit2Cut_ && std::abs(rMissHit2) < tidRMissHit2Cut_ && tidExtraCut == 1 && zDiff) seedCutsSatisfied = true;
    }else if(subdetector == 3){
      int tecExtraCut = 0;
      if(r1 > 23 && r1 < 32 && r2 > 26 && r2 < 42) tecExtraCut = 1;
      if(std::abs(phiMissHit2) < tecPhiMissHit2Cut_ && std::abs(rMissHit2) < tecRMissHit2Cut_ && tecExtraCut == 1 && zDiff) seedCutsSatisfied = true;
    }

  }

  if(!seedCutsSatisfied) return false;

  // seed checks borrowed from pixel-based algoritm


  /* Some of this code could be better optimized.  The Pixel algorithm natively
     takes Transient rec hits, so to recycle code we have to build them.
  */

  RecHitPointer hit1Trans = TSiStripMatchedRecHit::build(trackerGeometryHandle->idToDet((*hit1)->geographicalId()), *hit1, theMatcher_);
  RecHitPointer hit2Trans = TSiStripMatchedRecHit::build(trackerGeometryHandle->idToDet((*hit2)->geographicalId()), *hit2, theMatcher_);

  typedef TrajectoryStateOnSurface TSOS;

  double vertexZ = z1 - (r1 * (zc - z1) ) / (rc - r1);
  GlobalPoint eleVertex(0.,0.,vertexZ);

  // make a spiral
  FastHelix helix(hit2Pos,hit1Pos,eleVertex,*theSetup);
  if (!helix.isValid()) return false;

  FreeTrajectoryState fts = helix.stateAtVertex();
  TSOS propagatedState = thePropagator->propagate(fts,hit1Trans->det()->surface());

  if (!propagatedState.isValid()) return false;

  TSOS updatedState = theUpdator->update(propagatedState, *hit1Trans);
  TSOS propagatedState_out = thePropagator->propagate(fts,hit2Trans->det()->surface()) ;

  if (!propagatedState_out.isValid()) return false;

  // the seed has now passed all the cuts

  TSOS updatedState_out = theUpdator->update(propagatedState_out, *hit2Trans);

  pts_ =  trajectoryStateTransform::persistentState(updatedState_out, hit2Trans->geographicalId().rawId());

  return true;
}

void SiStripElectronSeedGenerator::findSeedsFromCluster	(	edm::Ref< reco::SuperClusterCollection >	seedCluster,
		edm::Handle< reco::BeamSpot >	bs,
		reco::ElectronSeedCollection &	result
	)		`[private]`

Definition at line 138 of file SiStripElectronSeedGenerator.cc.

References abs, alongMomentum, SiPixelRawToDigiRegional_cfi::beamSpot, funct::cos(), dir, LayerMeasurements::measurements(), PV3DBase< T, PVType, FrameType >::perp(), position, funct::pow(), reco::ElectronSeed::setCaloCluster(), funct::sin(), mathSSE::sqrt(), FastHelix::stateAtVertex(), and GeometricSearchTracker::tibLayers().

Referenced by run().

 {
  // clear the member vectors of good hits
  layer1Hits_.clear() ;
  layer2Hits_.clear() ;
  backupLayer2Hits_.clear() ;

  using namespace std;

  double sCenergy = seedCluster->energy();
  math::XYZPoint sCposition = seedCluster->position();
  double scEta = seedCluster->eta();

  double scz = sCposition.z();
  double scr = sqrt(pow(sCposition.x(),2)+pow(sCposition.y(),2));

  double pT = sCenergy * seedCluster->position().rho()/sqrt(seedCluster->x()*seedCluster->x()+seedCluster->y()*seedCluster->y()+seedCluster->z()*seedCluster->z());

  double magneticField = 3.8;

  // cf Jackson p. 581-2, a little geometry
  double phiVsRSlope = -3.00e-3 * magneticField / pT / 2.;


  //Need to create TSOS to feed MeasurementTracker
  GlobalPoint beamSpot(bs->x0(),bs->y0(),bs->z0());
  GlobalPoint superCluster(sCposition.x(),sCposition.y(),sCposition.z());
  double r0 = beamSpot.perp();
  double z0 = beamSpot.z();

  //We need to pick a charge for the particle we want to reconstruct before hits can be retrieved
  //Choosing both charges improves seeding efficiency by less than 0.5% for signal events
  //If we pick a single charge, this reduces fake rate and CPU time
  //So we pick a charge that is equally likely to be positive or negative

  int chargeHypothesis;
  double chargeSelector = sCenergy - (int)sCenergy;
  if(chargeSelector >= 0.5) chargeHypothesis = -1;
  if(chargeSelector < 0.5) chargeHypothesis = 1;

  //Use BeamSpot and SC position to estimate 3rd point
  double rFake = 25.;
  double phiFake = phiDiff(superCluster.phi(),chargeHypothesis * phiVsRSlope * (scr - rFake));
  double zFake = (rFake*(scz-z0)-r0*scz+scr*z0)/(scr-r0);
  double xFake = rFake * cos(phiFake);
  double yFake = rFake * sin(phiFake);
  GlobalPoint fakePoint(xFake,yFake,zFake);

  //Use 3 points to make helix
  FastHelix initialHelix(superCluster,fakePoint,beamSpot,*theSetup);

  //Use helix to get FTS
  FreeTrajectoryState initialFTS = initialHelix.stateAtVertex();

  //Use FTS and BeamSpot to create TSOS
  TransverseImpactPointExtrapolator* tipe = new TransverseImpactPointExtrapolator(*thePropagator);
  TrajectoryStateOnSurface initialTSOS = tipe->extrapolate(initialFTS,beamSpot);

  //Use GST to retrieve hits from various DetLayers using layerMeasurements class
  const GeometricSearchTracker* gst = theMeasurementTracker->geometricSearchTracker();

  std::vector<BarrelDetLayer*> tibLayers = gst->tibLayers();
  DetLayer* tib1 = tibLayers.at(0);
  DetLayer* tib2 = tibLayers.at(1);

  std::vector<ForwardDetLayer*> tecLayers;
  std::vector<ForwardDetLayer*> tidLayers;
  if(scEta < 0){
    tecLayers = gst->negTecLayers();
    tidLayers = gst->negTidLayers();
  }
  if(scEta > 0){
    tecLayers = gst->posTecLayers();
    tidLayers = gst->posTidLayers();
  }

  DetLayer* tid1 = tidLayers.at(0);
  DetLayer* tid2 = tidLayers.at(1);
  DetLayer* tid3 = tidLayers.at(2);
  DetLayer* tec1 = tecLayers.at(0);
  DetLayer* tec2 = tecLayers.at(1);
  DetLayer* tec3 = tecLayers.at(2);

  //Figure out which DetLayers to use based on SC Eta
  std::vector<bool> useDL = useDetLayer(scEta);
  bool useTID = false;

  //Use counters to restrict the number of hits in TID and TEC layers
  //This reduces seed multiplicity
  int tid1MHC = 0;
  int tid2MHC = 0;
  int tid3MHC = 0;
  int tid1BHC = 0;
  int tid2BHC = 0;
  int tid3BHC = 0;
  int tec1MHC = 0;
  int tec2MHC = 0;
  int tec3MHC = 0;

  //Use counter to limit the allowed number of seeds
  int seedCounter = 0;

  bool hasLay1Hit = false;
  bool hasLay2Hit = false;
  bool hasBackupHit = false;

  LayerMeasurements layerMeasurements(theMeasurementTracker);

  std::vector<TrajectoryMeasurement> tib1measurements;
  if(useDL.at(0)) tib1measurements = layerMeasurements.measurements(*tib1,initialTSOS,*thePropagator,*theEstimator);
  std::vector<TrajectoryMeasurement> tib2measurements;
  if(useDL.at(1)) tib2measurements = layerMeasurements.measurements(*tib2,initialTSOS,*thePropagator,*theEstimator);

  //Basic idea: Retrieve hits from a given DetLayer
  //Check if it is a Matched Hit and satisfies some cuts
  //If yes, accept hit for seed making

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tib1measurements.begin(); tmIter != tib1measurements.end(); ++ tmIter){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 1)){
        hasLay1Hit = true;
        layer1Hits_.push_back(matchedHit);
      }
    }
  }

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tib2measurements.begin(); tmIter != tib2measurements.end(); ++ tmIter){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 1)){
        hasLay2Hit = true;
        layer2Hits_.push_back(matchedHit);
      }
    }
  }

  if(!(hasLay1Hit && hasLay2Hit)) useTID = true;
  if(std::abs(scEta) > tidEtaUsage_) useTID = true;
  std::vector<TrajectoryMeasurement> tid1measurements;
  if(useDL.at(2) && useTID) tid1measurements = layerMeasurements.measurements(*tid1,initialTSOS,*thePropagator,*theEstimator);
  std::vector<TrajectoryMeasurement> tid2measurements;
  if(useDL.at(3) && useTID) tid2measurements = layerMeasurements.measurements(*tid2,initialTSOS,*thePropagator,*theEstimator);
  std::vector<TrajectoryMeasurement> tid3measurements;
  if(useDL.at(4) && useTID) tid3measurements = layerMeasurements.measurements(*tid3,initialTSOS,*thePropagator,*theEstimator);

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tid1measurements.begin(); tmIter != tid1measurements.end(); ++ tmIter){
    if(tid1MHC < tidMaxHits_){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 2)){
        tid1MHC++;
        hasLay1Hit = true;
        layer1Hits_.push_back(matchedHit);
        hasLay2Hit = true;
        layer2Hits_.push_back(matchedHit);
      }
    }else if(useDL.at(8) && tid1BHC < monoMaxHits_){
      const SiStripRecHit2D* backupHit = backupHitConverter(hit);
      if(backupHit){
        GlobalPoint position = trackerGeometryHandle->idToDet(backupHit->geographicalId())->surface().toGlobal(backupHit->localPosition());
        if(preselection(position, superCluster, phiVsRSlope, 4) && position.perp() > 37.){
          tid1BHC++;
          hasBackupHit = true;
          backupLayer2Hits_.push_back(backupHit);
        }
      }
    }
    }
  }

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tid2measurements.begin(); tmIter != tid2measurements.end(); ++ tmIter){
    if(tid2MHC < tidMaxHits_){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 2)){
        tid2MHC++;
        hasLay1Hit = true;
        layer1Hits_.push_back(matchedHit);
        hasLay2Hit = true;
        layer2Hits_.push_back(matchedHit);
      }
    }else if(useDL.at(8) && tid2BHC < monoMaxHits_){
      const SiStripRecHit2D* backupHit = backupHitConverter(hit);
      if(backupHit){
        GlobalPoint position = trackerGeometryHandle->idToDet(backupHit->geographicalId())->surface().toGlobal(backupHit->localPosition());
        if(preselection(position, superCluster, phiVsRSlope, 4) && position.perp() > 37.){
          tid2BHC++;
          hasBackupHit = true;
          backupLayer2Hits_.push_back(backupHit);
        }
      }
    }
    }
  }

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tid3measurements.begin(); tmIter != tid3measurements.end(); ++ tmIter){
    if(tid3MHC < tidMaxHits_){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 2)){
        tid3MHC++;
        hasLay1Hit = true;
        layer1Hits_.push_back(matchedHit);
        hasLay2Hit = true;
        layer2Hits_.push_back(matchedHit);
      }
    }else if(useDL.at(8) && tid3BHC < monoMaxHits_){
      const SiStripRecHit2D* backupHit = backupHitConverter(hit);
      if(backupHit){
        GlobalPoint position = trackerGeometryHandle->idToDet(backupHit->geographicalId())->surface().toGlobal(backupHit->localPosition());
        if(preselection(position, superCluster, phiVsRSlope, 4) && position.perp() > 37.){
          tid3BHC++;
          hasBackupHit = true;
          backupLayer2Hits_.push_back(backupHit);
        }
      }
    }
    }
  }

  std::vector<TrajectoryMeasurement> tec1measurements;
  if(useDL.at(5)) tec1measurements = layerMeasurements.measurements(*tec1,initialTSOS,*thePropagator,*theEstimator);
  std::vector<TrajectoryMeasurement> tec2measurements;
  if(useDL.at(6)) tec2measurements = layerMeasurements.measurements(*tec2,initialTSOS,*thePropagator,*theEstimator);
  std::vector<TrajectoryMeasurement> tec3measurements;
  if(useDL.at(7)) tec3measurements = layerMeasurements.measurements(*tec3,initialTSOS,*thePropagator,*theEstimator);

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tec1measurements.begin(); tmIter != tec1measurements.end(); ++ tmIter){
    if(tec1MHC < tecMaxHits_){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 3)){
        tec1MHC++;
        hasLay1Hit = true;
        layer1Hits_.push_back(matchedHit);
        hasLay2Hit = true;
        layer2Hits_.push_back(matchedHit);
      }
    }
    }
  }

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tec2measurements.begin(); tmIter != tec2measurements.end(); ++ tmIter){
    if(tec2MHC < tecMaxHits_){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 3)){
        tec2MHC++;
        hasLay1Hit = true;
        layer1Hits_.push_back(matchedHit);
        hasLay2Hit = true;
        layer2Hits_.push_back(matchedHit);
      }
    }
    }
  }

  for(std::vector<TrajectoryMeasurement>::const_iterator tmIter = tec3measurements.begin(); tmIter != tec3measurements.end(); ++ tmIter){
    if(tec3MHC < tecMaxHits_){
    ConstRecHitPointer hit = tmIter->recHit();
    const SiStripMatchedRecHit2D* matchedHit = matchedHitConverter(hit);
    if(matchedHit){
      GlobalPoint position = trackerGeometryHandle->idToDet(matchedHit->geographicalId())->surface().toGlobal(matchedHit->localPosition());
      if(preselection(position, superCluster, phiVsRSlope, 3)){
        tec3MHC++;
        hasLay2Hit = true;
        layer2Hits_.push_back(matchedHit);
      }
    }
    }
  }

  // We have 2 arrays of hits, combine them to form seeds
  if( hasLay1Hit && hasLay2Hit ){

    for (std::vector<const SiStripMatchedRecHit2D*>::const_iterator hit1 = layer1Hits_.begin() ; hit1!= layer1Hits_.end(); ++hit1) {
      for (std::vector<const SiStripMatchedRecHit2D*>::const_iterator hit2 = layer2Hits_.begin() ; hit2!= layer2Hits_.end(); ++hit2) {

        if(seedCounter < maxSeeds_){

          if(checkHitsAndTSOS(hit1,hit2,scr,scz,pT,scEta)) {

            seedCounter++;

            recHits_.clear();

            SiStripMatchedRecHit2D *hit;
            hit=new SiStripMatchedRecHit2D(*(dynamic_cast <const SiStripMatchedRecHit2D *> (*hit1) ) );
            recHits_.push_back(hit);
            hit=new SiStripMatchedRecHit2D(*(dynamic_cast <const SiStripMatchedRecHit2D *> (*hit2) ) );
            recHits_.push_back(hit);

            PropagationDirection dir = alongMomentum;
            reco::ElectronSeed seed(pts_,recHits_,dir) ;
            reco::ElectronSeed::CaloClusterRef caloCluster(seedCluster) ;
            seed.setCaloCluster(caloCluster) ;
            result.push_back(seed);

          }

        }

      }// end of hit 2 loop

    }// end of hit 1 loop

  }//end of seed making

  //Make seeds using TID Ring 3 if necessary

  if(hasLay1Hit && hasBackupHit && seedCounter == 0){

    for (std::vector<const SiStripMatchedRecHit2D*>::const_iterator hit1 = layer1Hits_.begin() ; hit1!= layer1Hits_.end(); ++hit1) {
      for (std::vector<const SiStripRecHit2D*>::const_iterator hit2 = backupLayer2Hits_.begin() ; hit2!= backupLayer2Hits_.end(); ++hit2) {

        if(seedCounter < maxSeeds_){

          if(altCheckHitsAndTSOS(hit1,hit2,scr,scz,pT,scEta)) {

            seedCounter++;

            recHits_.clear();

            SiStripMatchedRecHit2D *innerHit;
            innerHit=new SiStripMatchedRecHit2D(*(dynamic_cast <const SiStripMatchedRecHit2D *> (*hit1) ) );
            recHits_.push_back(innerHit);
            SiStripRecHit2D *outerHit;
            outerHit=new SiStripRecHit2D(*(dynamic_cast <const SiStripRecHit2D *> (*hit2) ) );
            recHits_.push_back(outerHit);

            PropagationDirection dir = alongMomentum;
            reco::ElectronSeed seed(pts_,recHits_,dir) ;
            reco::ElectronSeed::CaloClusterRef caloCluster(seedCluster) ;
            seed.setCaloCluster(caloCluster) ;
            result.push_back(seed);

          }

        }

      }// end of hit 2 loop

    }// end of hit 1 loop

  }// end of backup seed making

} // end of findSeedsFromCluster

const SiStripMatchedRecHit2D * SiStripElectronSeedGenerator::matchedHitConverter ( ConstRecHitPointer crhp ) [private]

Definition at line 731 of file SiStripElectronSeedGenerator.cc.

                                                                                                      {
  const TrackingRecHit* trh = crhp->hit();
  const SiStripMatchedRecHit2D* matchedHit = dynamic_cast<const SiStripMatchedRecHit2D*>(trh);
  return matchedHit;
}

double SiStripElectronSeedGenerator::normalPhi ( double phi ) const [inline, private]

Definition at line 83 of file SiStripElectronSeedGenerator.h.

References M_PI, and phi.

Referenced by phiDiff().

                                     {
    while (phi > 2.* M_PI) { phi -= 2.*M_PI; }
    while (phi < 0) { phi += 2.*M_PI; }
    return phi;
  }

double SiStripElectronSeedGenerator::phiDiff	(	double	phi1,
		double	phi2
	)		`[inline, private]`

Definition at line 89 of file SiStripElectronSeedGenerator.h.

References M_PI, normalPhi(), and query::result.

Referenced by altCheckHitsAndTSOS(), checkHitsAndTSOS(), and preselection().

                                          {
    double result = normalPhi(phi1) - normalPhi(phi2);
    if(result > M_PI) result -= 2*M_PI;
    if(result < -M_PI) result += 2*M_PI;
    return result;
  }

bool SiStripElectronSeedGenerator::preselection	(	GlobalPoint	position,
		GlobalPoint	superCluster,
		double	phiVsRSlope,
		int	hitLayer
	)		`[private]`

Definition at line 684 of file SiStripElectronSeedGenerator.cc.

References abs, monoDeltaPsiCut_, monoOriginZCut_, PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), phi, phiDiff(), alignCSCRings::r, query::result, tecDeltaPsiCut_, tecOriginZCut_, tibDeltaPsiCut_, tibOriginZCut_, tidDeltaPsiCut_, tidOriginZCut_, PV3DBase< T, PVType, FrameType >::z(), and z.

                                                                                                                            {
  double r = position.perp();
  double phi = position.phi();
  double z = position.z();
  double scr = superCluster.perp();
  double scphi = superCluster.phi();
  double scz = superCluster.z();
  double psi = phiDiff(phi,scphi);
  double deltaPsi = psi - (scr-r)*phiVsRSlope;
  double antiDeltaPsi = psi - (r-scr)*phiVsRSlope;
  double dP;
  if (std::abs(deltaPsi)<std::abs(antiDeltaPsi)){
    dP = deltaPsi;
  }else{
    dP = antiDeltaPsi;
  }
  double originZ = (scr*z - r*scz)/(scr-r);

  bool result = false;

  if(hitLayer == 1){
    if(std::abs(originZ) < tibOriginZCut_ && std::abs(dP) < tibDeltaPsiCut_) result = true;
  }else if(hitLayer == 2){
    if(std::abs(originZ) < tidOriginZCut_ && std::abs(dP) < tidDeltaPsiCut_) result = true;
  }else if(hitLayer == 3){
    if(std::abs(originZ) < tecOriginZCut_ && std::abs(dP) < tecDeltaPsiCut_) result = true;
  }else if(hitLayer == 4){
    if(std::abs(originZ) < monoOriginZCut_ && std::abs(dP) < monoDeltaPsiCut_) result = true;
  }

  return result;
}

void SiStripElectronSeedGenerator::run	(	edm::Event &	e,
		const edm::EventSetup &	setup,
		const edm::Handle< reco::SuperClusterCollection > &	clusters,
		reco::ElectronSeedCollection &	out
	)

Definition at line 116 of file SiStripElectronSeedGenerator.cc.

References beamSpotTag_, findSeedsFromCluster(), edm::Event::getByLabel(), i, edm::EventBase::id(), LogDebug, HcalObjRepresent::setup(), theBeamSpot, theMeasurementTracker, theSetup, and MeasurementTracker::update().

                                                                          {
  theSetup= &setup;
  e.getByLabel(beamSpotTag_,theBeamSpot);
  theMeasurementTracker->update(e);

  for  (unsigned int i=0;i<clusters->size();++i) {
    edm::Ref<reco::SuperClusterCollection> theClusB(clusters,i);
    // Find the seeds
    LogDebug ("run") << "new cluster, calling findSeedsFromCluster";
    findSeedsFromCluster(theClusB,theBeamSpot,out);
  }

  LogDebug ("run") << ": For event "<<e.id();
  LogDebug ("run") <<"Nr of superclusters: "<<clusters->size()
                   <<", no. of ElectronSeeds found  = " << out.size();
}

void SiStripElectronSeedGenerator::setupES ( const edm::EventSetup & setup )

Definition at line 94 of file SiStripElectronSeedGenerator.cc.

References alongMomentum, cacheIDCkfComp_, edm::eventsetup::EventSetupRecord::cacheIdentifier(), cacheIDMagField_, cacheIDTrkGeom_, edm::EventSetup::get(), measurementTrackerHandle, edm::ESHandle< T >::product(), theMagField, theMeasurementTracker, theMeasurementTrackerName, thePropagator, and trackerGeometryHandle.

                                                                     {

  if (cacheIDMagField_!=setup.get<IdealMagneticFieldRecord>().cacheIdentifier()) {
    setup.get<IdealMagneticFieldRecord>().get(theMagField);
    cacheIDMagField_=setup.get<IdealMagneticFieldRecord>().cacheIdentifier();
    if (thePropagator) delete thePropagator;
    thePropagator = new PropagatorWithMaterial(alongMomentum,.000511,&(*theMagField));
  }

  if (cacheIDCkfComp_!=setup.get<CkfComponentsRecord>().cacheIdentifier()) {
    setup.get<CkfComponentsRecord>().get(theMeasurementTrackerName,measurementTrackerHandle);
    cacheIDCkfComp_=setup.get<CkfComponentsRecord>().cacheIdentifier();
    theMeasurementTracker = measurementTrackerHandle.product();
  }

  if (cacheIDTrkGeom_!=setup.get<TrackerDigiGeometryRecord>().cacheIdentifier()) {
    setup.get<TrackerDigiGeometryRecord>().get(trackerGeometryHandle);
    cacheIDTrkGeom_=setup.get<TrackerDigiGeometryRecord>().cacheIdentifier();
  }

}

double SiStripElectronSeedGenerator::unwrapPhi ( double phi ) const [inline, private]

Definition at line 96 of file SiStripElectronSeedGenerator.h.

References M_PI, and phi.

                                     {
    while (phi > M_PI) { phi -= 2.*M_PI; }
    while (phi < -M_PI) { phi += 2.*M_PI; }
    return phi;
  }

std::vector< bool > SiStripElectronSeedGenerator::useDetLayer ( double scEta ) [private]

Definition at line 743 of file SiStripElectronSeedGenerator.cc.

References abs.

                                                                     {
  std::vector<bool> useDetLayer;
  double variable = std::abs(scEta);
  if(variable > 0 && variable < 1.8){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 0 && variable < 1.5){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 1 && variable < 2.1){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 1 && variable < 2.2){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 1 && variable < 2.3){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 1.8 && variable < 2.5){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 1.8 && variable < 2.5){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 1.8 && variable < 2.5){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  if(variable > 1.2 && variable < 1.6){
    useDetLayer.push_back(true);
  }else{
    useDetLayer.push_back(false);
  }
  return useDetLayer;
}

int SiStripElectronSeedGenerator::whichSubdetector ( std::vector< const SiStripMatchedRecHit2D * >::const_iterator hit ) [private]

Definition at line 719 of file SiStripElectronSeedGenerator.cc.

References query::result, sistripsummary::TEC, sistripsummary::TIB, and sistripsummary::TID.

Referenced by checkHitsAndTSOS().

                                                                                                            {
  int result = 0;
  if(((*hit)->geographicalId()).subdetId() == StripSubdetector::TIB){
    result = 1;
  }else if(((*hit)->geographicalId()).subdetId() == StripSubdetector::TID){
    result = 2;
  }else if(((*hit)->geographicalId()).subdetId() == StripSubdetector::TEC){
    result = 3;
  }
  return result;
}