#include <CSCSegAlgoDF.h>

Inheritance diagram for CSCSegAlgoDF:

Public Types
typedef std::deque< bool >	BoolContainer

typedef std::vector< const CSCRecHit2D * >	ChamberHitContainer

typedef std::vector< const CSCRecHit2D * >::const_iterator	ChamberHitContainerCIt

typedef std::vector< int >	LayerIndex
	Typedefs. More...

Public Member Functions
std::vector< CSCSegment >	buildSegments (const ChamberHitContainer &rechits)

	CSCSegAlgoDF (const edm::ParameterSet &ps)
	Constructor. More...

std::vector< CSCSegment >	run (const CSCChamber *aChamber, const ChamberHitContainer &rechits) override

	~CSCSegAlgoDF () override
	Destructor. More...

Public Member Functions inherited from CSCSegmentAlgorithm
	CSCSegmentAlgorithm (const edm::ParameterSet &)
	Constructor. More...

virtual std::vector< CSCSegment >	run (const CSCChamber chamber, const std::vector< const CSCRecHit2D > &rechits)=0

virtual	~CSCSegmentAlgorithm ()
	Destructor. More...

Private Member Functions
bool	addHit (const CSCRecHit2D *hit, int layer)

void	compareProtoSegment (const CSCRecHit2D *h, int layer)

void	dumpSegment (const CSCSegment &seg) const

void	flagHitsAsUsed (const ChamberHitContainer &rechitsInChamber)

bool	hasHitOnLayer (int layer) const

bool	isHitNearSegment (const CSCRecHit2D *h) const

void	pruneFromResidual (void)

void	tryAddingHitsToSegment (const ChamberHitContainer &rechitsInChamber, const ChamberHitContainerCIt i1, const ChamberHitContainerCIt i2, const LayerIndex &layerIndex)
	Utility functions. More...

void	updateParameters (void)

Private Attributes
float	chi2Max

ChamberHitContainer	closeHits

bool	debug

double	dPhiFineMax

double	dRPhiFineMax

float	maxRatioResidual

int	minHitsForPreClustering

int	minHitsPerSegment

int	minLayersApart

const std::string	myName

int	nHitsPerClusterIsShower

CSCSegAlgoPreClustering *	preCluster_

bool	preClustering

ChamberHitContainer	protoSegment

bool	Pruning

ChamberHitContainer	secondSeedHits

CSCSegFit *	sfit_

CSCSegAlgoShowering *	showering_

float	tanPhiMax

float	tanThetaMax

const CSCChamber *	theChamber

BoolContainer	usedHits

Detailed Description

This is a modified version of the SK algorithm for building endcap muon track segments out of the rechit's in a CSCChamber.

A CSCSegment is a RecSegment4D, and is built from CSCRecHit2D objects, each of which is a RecHit2DLocalPos.

This builds segments by first creating proto-segments from at least 3 hits. We intend to try all possible pairs of hits to start segment building. 'All possible' means each hit lies on different layers in the chamber. Once a hit has been assigned to a segment, we don't consider it again, THAT IS, FOR THE FIRST PASS ONLY ! In fact, this is one of the possible flaw with the SK algorithms as it sometimes manages to build segments with the wrong starting points. In the DF algorithm, the endpoints are tested as the best starting points in a 2nd and 3rd loop.

Another difference with the from the SK algorithm is that rechits can be added to proto segments if they fall within n sigmas of the projected track within a given layer. Hence, a cylinder isn't used as in the SK algorimthm, which allows for pseudo 2D hits built from wire or strip only hits to be used in segment reconstruction.

Also, only a certain muonsPerChamberMax maximum number of segments can be produced in the chamber. [Seems to be hardwired rather than using this variable?]

Alternative algorithms can be used for the segment building by writing classes like this, and then selecting which one is actually used via the CSCSegmentBuilder.

Author: Dominique Fortin - UCR

Definition at line 48 of file CSCSegAlgoDF.h.

Member Typedef Documentation

typedef std::deque<bool> CSCSegAlgoDF::BoolContainer

Definition at line 55 of file CSCSegAlgoDF.h.

typedef std::vector<const CSCRecHit2D*> CSCSegAlgoDF::ChamberHitContainer

Definition at line 53 of file CSCSegAlgoDF.h.

typedef std::vector<const CSCRecHit2D*>::const_iterator CSCSegAlgoDF::ChamberHitContainerCIt

Definition at line 54 of file CSCSegAlgoDF.h.

typedef std::vector<int> CSCSegAlgoDF::LayerIndex

Typedefs.

Definition at line 52 of file CSCSegAlgoDF.h.

Constructor & Destructor Documentation

CSCSegAlgoDF::CSCSegAlgoDF ( const edm::ParameterSet & ps )

explicit

Constructor.

Definition at line 34 of file CSCSegAlgoDF.cc.

References chi2Max, debug, dPhiFineMax, dRPhiFineMax, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), maxRatioResidual, minHitsForPreClustering, minHitsPerSegment, minLayersApart, nHitsPerClusterIsShower, preCluster_, preClustering, Pruning, showering_, tanPhiMax, and tanThetaMax.

     : CSCSegmentAlgorithm(ps), myName("CSCSegAlgoDF"), sfit_(nullptr) {
   debug = ps.getUntrackedParameter<bool>("CSCSegmentDebug");
   minLayersApart = ps.getParameter<int>("minLayersApart");
   minHitsPerSegment = ps.getParameter<int>("minHitsPerSegment");
   dRPhiFineMax = ps.getParameter<double>("dRPhiFineMax");
   dPhiFineMax = ps.getParameter<double>("dPhiFineMax");
   tanThetaMax = ps.getParameter<double>("tanThetaMax");
   tanPhiMax = ps.getParameter<double>("tanPhiMax");
   chi2Max = ps.getParameter<double>("chi2Max");
   preClustering = ps.getUntrackedParameter<bool>("preClustering");
   minHitsForPreClustering = ps.getParameter<int>("minHitsForPreClustering");
   nHitsPerClusterIsShower = ps.getParameter<int>("nHitsPerClusterIsShower");
   Pruning = ps.getUntrackedParameter<bool>("Pruning");
   maxRatioResidual = ps.getParameter<double>("maxRatioResidualPrune");
 
   preCluster_ = new CSCSegAlgoPreClustering(ps);
   showering_ = new CSCSegAlgoShowering(ps);
 }

CSCSegAlgoDF::~CSCSegAlgoDF ( )

override

Destructor.

Definition at line 57 of file CSCSegAlgoDF.cc.

References preCluster_, and showering_.

                             {
   delete preCluster_;
   delete showering_;
 }

Member Function Documentation

bool CSCSegAlgoDF::addHit	(	const CSCRecHit2D *	hit,
		int	layer
	)

private

Definition at line 423 of file CSCSegAlgoDF.cc.

References protoSegment.

Referenced by compareProtoSegment(), and tryAddingHitsToSegment().

                                                             {
   //  std::cout << "[CSCSegAlgoDF::addHit] on layer " << layer << " to protoSegment.size() = "
   //        << protoSegment.size() << std::endl;
 
   // Return true if hit was added successfully and then parameters are updated.
   // Return false if there is already a hit on the same layer, or insert failed.
 
   if (protoSegment.size() > 5)
     return false;  //@@ can only have 6 hits at most
 
   // Test that we are not trying to add the same hit again
   for (ChamberHitContainer::const_iterator it = protoSegment.begin(); it != protoSegment.end(); ++it)
     if (aHit == (*it))
       return false;
 
   protoSegment.push_back(aHit);
 
   return true;
 }

std::vector< CSCSegment > CSCSegAlgoDF::buildSegments ( const ChamberHitContainer & rechits )

Build track segments in this chamber (this is where the actual segment-building algorithm hides.)

Definition at line 109 of file CSCSegAlgoDF.cc.

References CSCSegFit::chi2(), chi2Max, CSCSegFit::covarianceMatrix(), debug, dumpSegment(), PVValHelper::dz, flagHitsAsUsed(), CSCSegFit::hits(), mps_fire::i, testProducerWithPsetDescEmpty_cfi::i1, testProducerWithPsetDescEmpty_cfi::i2, cuy::ib, CSCSegFit::intercept(), CSCChamber::layer(), CSCSegFit::localdir(), CSCRecHit2D::localPosition(), minHitsPerSegment, minLayersApart, nHitsPerClusterIsShower, ticl::constants::nLayers, CSCSegment::nRecHits(), GeomDet::position(), preClustering, protoSegment, pruneFromResidual(), Pruning, TrackInfoProducer_cfi::rechits, groupFilesInBlocks::reverse, sfit_, showering_, CSCSegAlgoShowering::showerSeg(), tanPhiMax, tanThetaMax, groupFilesInBlocks::temp, theChamber, GeomDet::toGlobal(), GeomDet::toLocal(), tryAddingHitsToSegment(), updateParameters(), usedHits, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by run().

                                                                                      {
   ChamberHitContainer rechits = _rechits;
   // Clear buffer for segment vector
   std::vector<CSCSegment> segmentInChamber;
   segmentInChamber.clear();
 
   unsigned nHitInChamber = rechits.size();
 
   //  std::cout << "[CSCSegAlgoDF::buildSegments] address of chamber = " << theChamber << std::endl;
   //  std::cout << "[CSCSegAlgoDF::buildSegments] starting in " << theChamber->id()
   //         << " with " << nHitInChamber << " rechits" << std::endl;
 
   // Return #1 - OK, there aren't enough rechits to build a segment
   if (nHitInChamber < 3)
     return segmentInChamber;
 
   LayerIndex layerIndex(nHitInChamber);
 
   size_t nLayers = 0;
   size_t old_layer = 0;
   for (size_t i = 0; i < nHitInChamber; ++i) {
     size_t this_layer = rechits[i]->cscDetId().layer();
     //    std::cout << "[CSCSegAlgoDF::buildSegments] this_layer = " << this_layer << std::endl;
     layerIndex[i] = this_layer;
     //    std::cout << "[CSCSegAlgoDF::buildSegments] layerIndex[" << i << "] = " << layerIndex[i] << std::endl;
     if (this_layer != old_layer) {
       old_layer = this_layer;
       ++nLayers;
     }
   }
 
   //  std::cout << "[CSCSegAlgoDF::buildSegments] layers with rechits = " << nLayers << std::endl;
 
   // Return #2 - OK, there aren't enough hit layers to build a segment
   if (nLayers < 3)
     return segmentInChamber;
 
   double z1 = theChamber->layer(1)->position().z();
   double z6 = theChamber->layer(6)->position().z();
 
   if (z1 > 0.) {
     if (z1 > z6) {
       reverse(layerIndex.begin(), layerIndex.end());
       reverse(rechits.begin(), rechits.end());
     }
   } else if (z1 < 0.) {
     if (z1 < z6) {
       reverse(layerIndex.begin(), layerIndex.end());
       reverse(rechits.begin(), rechits.end());
     }
   }
 
   //  std::cout << "[CSCSegAlgoDF::buildSegments] rechits have been ordered" << std::endl;
 
   // Showering muon
   if (preClustering && int(nHitInChamber) > nHitsPerClusterIsShower && nLayers > 2) {
     // std::cout << "[CSCSegAlgoDF::buildSegments] showering block" << std::endl;
 
     CSCSegment segShower = showering_->showerSeg(theChamber, rechits);
 
     // Return #3 - OK, this is now 'effectve' rechits
     // Make sure have at least 3 hits...
     if (segShower.nRecHits() < 3)
       return segmentInChamber;
 
     segmentInChamber.push_back(segShower);
     if (debug)
       dumpSegment(segShower);
 
     // Return #4 - OK, only get one try at building a segment from shower
     return segmentInChamber;
   }
 
   // Initialize flags that a given hit has been allocated to a segment
   BoolContainer used_ini(rechits.size(), false);
   usedHits = used_ini;
 
   ChamberHitContainerCIt ib = rechits.begin();
   ChamberHitContainerCIt ie = rechits.end();
 
   //  std::cout << "[CSCSegAlgoDF::buildSegments] entering rechit loop" << std::endl;
 
   // Now Loop over hits within the chamber to find 1st seed for segment building
   for (ChamberHitContainerCIt i1 = ib; i1 < ie; ++i1) {
     if (usedHits[i1 - ib])
       continue;
 
     const CSCRecHit2D* h1 = *i1;
     int layer1 = layerIndex[i1 - ib];
     const CSCLayer* l1 = theChamber->layer(layer1);
     GlobalPoint gp1 = l1->toGlobal(h1->localPosition());
     LocalPoint lp1 = theChamber->toLocal(gp1);
 
     // Loop over hits backward to find 2nd seed for segment building
     for (ChamberHitContainerCIt i2 = ie - 1; i2 > ib; --i2) {
       if (usedHits[i2 - ib])
         continue;  // Hit has been used already
 
       int layer2 = layerIndex[i2 - ib];
       if ((layer2 - layer1) < minLayersApart)
         continue;
 
       const CSCRecHit2D* h2 = *i2;
       const CSCLayer* l2 = theChamber->layer(layer2);
       GlobalPoint gp2 = l2->toGlobal(h2->localPosition());
       LocalPoint lp2 = theChamber->toLocal(gp2);
 
       // Clear proto segment so it can be (re)-filled
       protoSegment.clear();
 
       // We want hit wrt chamber (and local z will be != 0)
       float dz = gp2.z() - gp1.z();
       float slope_u = (lp2.x() - lp1.x()) / dz;
       float slope_v = (lp2.y() - lp1.y()) / dz;
 
       // Test if entrance angle is roughly pointing towards IP
       if (fabs(slope_v) > tanThetaMax)
         continue;
       if (fabs(slope_u) > tanPhiMax)
         continue;
 
       protoSegment.push_back(h1);
       protoSegment.push_back(h2);
 
       //      std::cout << "[CSCSegAlgoDF::buildSegments] about to fit 2 hits on layers "
       //        << layer1 << " and " << layer2 << std::endl;
 
       // protoSegment has just 2 hits - but need to create a CSCSegFit to hold it in case
       // we fail to add any more hits
       updateParameters();
 
       // Try adding hits to proto segment
       tryAddingHitsToSegment(rechits, i1, i2, layerIndex);
 
       // Check no. of hits on segment to see if segment is large enough
       bool segok = true;
       unsigned iadd = 0;
 
       if (protoSegment.size() < minHitsPerSegment + iadd)
         segok = false;
 
       if (Pruning && segok)
         pruneFromResidual();
 
       // Check if segment satisfies chi2 requirement
       if (sfit_->chi2() > chi2Max)
         segok = false;
 
       if (segok) {
         // Create an actual CSCSegment - retrieve all info from the current fit
         CSCSegment temp(sfit_->hits(), sfit_->intercept(), sfit_->localdir(), sfit_->covarianceMatrix(), sfit_->chi2());
         //  std::cout << "[CSCSegAlgoDF::buildSegments] about to delete sfit= = " << sfit_ << std::endl;
         delete sfit_;
         sfit_ = nullptr;  // avoid possibility of attempting a second delete later
 
         segmentInChamber.push_back(temp);
         if (debug)
           dumpSegment(temp);
 
         // Return #5 - OK, fewer than 3 rechits not on this segment left in chamber
         if (nHitInChamber - protoSegment.size() < 3)
           return segmentInChamber;
         // Return $6 - already have more than 4 segments in this chamber
         if (segmentInChamber.size() > 4)
           return segmentInChamber;
 
         // Flag used hits
         flagHitsAsUsed(rechits);
       }
     }
   }
   // Return #7
   return segmentInChamber;
 }

void CSCSegAlgoDF::compareProtoSegment	(	const CSCRecHit2D *	h,
		int	layer
	)

private

Definition at line 483 of file CSCSegAlgoDF.cc.

References addHit(), CSCSegFit::chi2(), CSCSegFit::fit(), convertSQLiteXML::ok, protoSegment, sfit_, and theChamber.

Referenced by tryAddingHitsToSegment().

                                                                       {
   //  std::cout << "[CSCSegAlgoDF::compareProtoSegment] for hit on layer " << layer
   //        << " with protoSegment.size() = " << protoSegment.size() << std::endl;
 
   // Try adding the hit to existing segment, and remove old one existing in same layer
   ChamberHitContainer::iterator it;
   for (it = protoSegment.begin(); it != protoSegment.end();) {
     if ((*it)->cscDetId().layer() == layer) {
       it = protoSegment.erase(it);
     } else {
       ++it;
     }
   }
 
   //  std::cout << "[CSCSegAlgoDF::compareProtoSegment] about to add hit on layer " << layer
   //        << " with protoSegment.size() = " << protoSegment.size() << std::endl;
 
   bool ok = addHit(h, layer);
 
   CSCSegFit* newfit = nullptr;
   if (ok) {
     newfit = new CSCSegFit(theChamber, protoSegment);
     //    std::cout << "[CSCSegAlgoDF::compareProtoSegment] newfit = " << newfit << std::endl;
     newfit->fit();
   }
   if (!ok || (newfit->chi2() > sfit_->chi2())) {
     //    std::cout << "[CSCSegAlgoDF::compareProtoSegment] about to delete newfit = " << newfit << std::endl;
     delete newfit;  // failed to add a hit or new fit is worse
   } else {
     //    std::cout << "[CSCSegAlgoDF::compareProtoSegment] about to delete sfit_ = " << sfit_ << std::endl;
     delete sfit_;  // new fit is better
     sfit_ = newfit;
     //    std::cout << "[CSCSegAlgoDF::compareProtoSegment] reset sfit_ = " << sfit_ << std::endl;
   }
 }

void CSCSegAlgoDF::dumpSegment ( const CSCSegment & seg ) const

private

Definition at line 608 of file CSCSegAlgoDF.cc.

References CSCSegment::chi2(), CSCSegment::degreesOfFreedom(), CSCChamber::id(), CSCSegment::localDirection(), CSCSegment::localDirectionError(), CSCSegment::localPosition(), CSCSegment::localPositionError(), CSCSegment::parametersError(), CSCSegment::specificRecHits(), theChamber, and CSCSegment::time().

Referenced by buildSegments().

                                                           {
   edm::LogVerbatim("CSCSegment") << "CSCSegment in " << theChamber->id() << "\nlocal position = " << seg.localPosition()
                                  << "\nerror = " << seg.localPositionError()
                                  << "\nlocal direction = " << seg.localDirection()
                                  << "\nerror =" << seg.localDirectionError() << "\ncovariance matrix"
                                  << seg.parametersError() << "chi2/ndf = " << seg.chi2() << "/"
                                  << seg.degreesOfFreedom() << "\n#rechits = " << seg.specificRecHits().size()
                                  << "\ntime = " << seg.time();
 }

void CSCSegAlgoDF::flagHitsAsUsed ( const ChamberHitContainer & rechitsInChamber )

private

Flag hits on segment as used

Definition at line 524 of file CSCSegAlgoDF.cc.

References closeHits, cuy::ib, protoSegment, and usedHits.

Referenced by buildSegments().

                                                                              {
   // Flag hits on segment as used
   ChamberHitContainerCIt ib = rechitsInChamber.begin();
   ChamberHitContainerCIt hi, iu;
 
   for (hi = protoSegment.begin(); hi != protoSegment.end(); ++hi) {
     for (iu = ib; iu != rechitsInChamber.end(); ++iu) {
       if (*hi == *iu)
         usedHits[iu - ib] = true;
     }
   }
   // Don't reject hits marked as "nearby" for now.
   // So this is bypassed at all times for now !!!
   // Perhaps add back to speed up algorithm some more
   if (!closeHits.empty())
     return;
   for (hi = closeHits.begin(); hi != closeHits.end(); ++hi) {
     for (iu = ib; iu != rechitsInChamber.end(); ++iu) {
       if (*hi == *iu)
         usedHits[iu - ib] = true;
     }
   }
 }

bool CSCSegAlgoDF::hasHitOnLayer ( int layer ) const

private

Definition at line 466 of file CSCSegAlgoDF.cc.

References protoSegment.

Referenced by tryAddingHitsToSegment().

                                                 {
   //  std::cout << "[CSCSegAlgoDF::hasHitOnLayer] on layer " << layer << std::endl;
 
   // Is there already a hit on this layer?
   for (ChamberHitContainerCIt it = protoSegment.begin(); it != protoSegment.end(); it++)
     if ((*it)->cscDetId().layer() == layer)
       return true;
 
   return false;
 }

bool CSCSegAlgoDF::isHitNearSegment ( const CSCRecHit2D * h ) const

private

Definition at line 381 of file CSCSegAlgoDF.cc.

References CSCRecHit2D::cscDetId(), SiPixelRawToDigiRegional_cfi::deltaPhi, dPhiFineMax, dRPhiFineMax, CSCDetId::layer(), CSCChamber::layer(), CSCRecHit2D::localPosition(), M_PI, PV3DBase< T, PVType, FrameType >::phi(), dttmaxenums::R, sfit_, mathSSE::sqrt(), theChamber, GeomDet::toGlobal(), GeomDet::toLocal(), PV3DBase< T, PVType, FrameType >::x(), CSCSegFit::xfit(), PV3DBase< T, PVType, FrameType >::y(), CSCSegFit::yfit(), z, and PV3DBase< T, PVType, FrameType >::z().

Referenced by tryAddingHitsToSegment().

                                                                 {
   const CSCLayer* layer = theChamber->layer(hit->cscDetId().layer());
 
   // hit phi position in global coordinates
   GlobalPoint Hgp = layer->toGlobal(hit->localPosition());
   float Hphi = Hgp.phi();
   if (Hphi < 0.)
     Hphi += 2. * M_PI;
   LocalPoint Hlp = theChamber->toLocal(Hgp);
   float z = Hlp.z();
 
   float LocalX = sfit_->xfit(z);
   float LocalY = sfit_->yfit(z);
 
   LocalPoint Slp(LocalX, LocalY, z);
   GlobalPoint Sgp = theChamber->toGlobal(Slp);
   float Sphi = Sgp.phi();
   if (Sphi < 0.)
     Sphi += 2. * M_PI;
   float R = sqrt(Sgp.x() * Sgp.x() + Sgp.y() * Sgp.y());
 
   float deltaPhi = Sphi - Hphi;
   if (deltaPhi > 2. * M_PI)
     deltaPhi -= 2. * M_PI;
   if (deltaPhi < -2. * M_PI)
     deltaPhi += 2. * M_PI;
   if (deltaPhi < 0.)
     deltaPhi = -deltaPhi;
 
   float RdeltaPhi = R * deltaPhi;
 
   if (RdeltaPhi < dRPhiFineMax && deltaPhi < dPhiFineMax)
     return true;
 
   return false;
 }

void CSCSegAlgoDF::pruneFromResidual ( void )

private

Prune bad segment from the worse hit based on residuals

Definition at line 549 of file CSCSegAlgoDF.cc.

References CSCRecHit2D::cscDetId(), runTauDisplay::gp, dqmiolumiharvest::j, CSCChamber::layer(), CSCDetId::layer(), CSCRecHit2D::localPosition(), maxRatioResidual, protoSegment, CSCSegFit::Rdev(), sfit_, theChamber, GeomDet::toGlobal(), GeomDet::toLocal(), updateParameters(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by buildSegments().

                                          {
   // Only prune if have at least 5 hits
   if (protoSegment.size() < 5)
     return;
 
   // Now Study residuals
 
   float maxResidual = 0.;
   float sumResidual = 0.;
   int nHits = 0;
   int badIndex = -1;
   int j = 0;
 
   ChamberHitContainer::const_iterator ih;
 
   for (ih = protoSegment.begin(); ih != protoSegment.end(); ++ih) {
     const CSCRecHit2D& hit = (**ih);
     const CSCLayer* layer = theChamber->layer(hit.cscDetId().layer());
     GlobalPoint gp = layer->toGlobal(hit.localPosition());
     LocalPoint lp = theChamber->toLocal(gp);
 
     float residual = sfit_->Rdev(lp.x(), lp.y(), lp.z());
 
     sumResidual += residual;
     nHits++;
     if (residual > maxResidual) {
       maxResidual = residual;
       badIndex = j;
     }
     j++;
   }
 
   float corrAvgResidual = (sumResidual - maxResidual) / (nHits - 1);
 
   // Keep all hits
   if (maxResidual / corrAvgResidual < maxRatioResidual)
     return;
 
   // Drop worse hit and recompute segment properties + fill
 
   ChamberHitContainer newProtoSegment;
 
   j = 0;
   for (ih = protoSegment.begin(); ih != protoSegment.end(); ++ih) {
     if (j != badIndex)
       newProtoSegment.push_back(*ih);
     j++;
   }
 
   protoSegment.clear();
 
   for (ih = newProtoSegment.begin(); ih != newProtoSegment.end(); ++ih) {
     protoSegment.push_back(*ih);
   }
 
   // Compute segment parameters
   updateParameters();
 }

std::vector< CSCSegment > CSCSegAlgoDF::run	(	const CSCChamber *	aChamber,
		const ChamberHitContainer &	rechits
	)

override

Here we must implement the algorithm

Definition at line 65 of file CSCSegAlgoDF.cc.

References buildSegments(), CSCSegAlgoPreClustering::clusterHits(), minHitsForPreClustering, preCluster_, preClustering, and theChamber.

                                                                                                       {
   // Store chamber info in temp memory
   theChamber = aChamber;
 
   int nHits = rechits.size();
 
   // Segments prior to pruning
   std::vector<CSCSegment> segments_temp;
 
   if (preClustering && nHits > minHitsForPreClustering) {
     // This is where the segment origin is in the chamber on avg.
     std::vector<CSCSegment> testSegments;
     std::vector<ChamberHitContainer> clusteredHits = preCluster_->clusterHits(theChamber, rechits);
     // loop over the found clusters:
     for (std::vector<ChamberHitContainer>::iterator subrechits = clusteredHits.begin();
          subrechits != clusteredHits.end();
          ++subrechits) {
       // build the subset of segments:
       std::vector<CSCSegment> segs = buildSegments((*subrechits));
       // add the found subset of segments to the collection of all segments in this chamber:
       segments_temp.insert(segments_temp.end(), segs.begin(), segs.end());
     }
   } else {
     std::vector<CSCSegment> segs = buildSegments(rechits);
     // add the found subset of segments to the collection of all segments in this chamber:
     segments_temp.insert(segments_temp.end(), segs.begin(), segs.end());
   }
 
   return segments_temp;
 }

void CSCSegAlgoDF::tryAddingHitsToSegment	(	const ChamberHitContainer &	rechitsInChamber,
		const ChamberHitContainerCIt	i1,
		const ChamberHitContainerCIt	i2,
		const LayerIndex &	layerIndex
	)

private

Utility functions.

Try adding non-used hits to segment
Skip the layers containing the segment endpoints on first 2 passes, but then
try hits on layer containing the segment starting points on 2nd and/or 3rd pass
if segment has >2 hits. Test each hit on the other layers to see if it is near
the segment using rechit error matrix.
If it is, see whether there is already a hit on the segment from the same layer

if so, and there are more than 2 hits on the segment, copy the segment,
replace the old hit with the new hit. If the new segment chi2 is better
then replace the original segment with the new one
if not, copy the segment, add the hit if it's within a certain range.

Definition at line 290 of file CSCSegAlgoDF.cc.

References addHit(), closeHits, compareProtoSegment(), h, hasHitOnLayer(), mps_fire::i, cuy::ib, isHitNearSegment(), protoSegment, updateParameters(), and usedHits.

Referenced by buildSegments().

                                                                         {
   /* Iterate over the layers with hits in the chamber
  * Skip the layers containing the segment endpoints on first pass, but then
  * try hits on layer containing the segment starting points on 2nd pass
  * if segment has >2 hits.  Once a hit is added to a layer, don't replace it 
  * until 2nd iteration
  */
 
   //  std::cout << "[CSCSegAlgoDF::tryAddingHitsToSegment] entering"
   //        << " with rechits.size() = " << rechits.size() << std::endl;
 
   ChamberHitContainerCIt ib = rechits.begin();
   ChamberHitContainerCIt ie = rechits.end();
   closeHits.clear();
 
   //  int counter1 = 0;
   //  int counter2 = 0;
 
   for (ChamberHitContainerCIt i = ib; i != ie; ++i) {
     //    std::cout << "counter1 = " << ++counter1 << std::endl;
     if (i == i1 || i == i2)
       continue;
     if (usedHits[i - ib])
       continue;  // Don't use hits already part of a segment.
 
     //    std::cout << "counter2 = " << ++counter2 << std::endl;
     const CSCRecHit2D* h = *i;
     int layer = layerIndex[i - ib];
     int layer1 = layerIndex[i1 - ib];
     int layer2 = layerIndex[i2 - ib];
 
     //    std::cout << "[CSCSegAlgoDF::tryAddingHitsToSegment] layer, layer1, layer2 = "
     //        << layer << ", " << layer1 << ", " << layer2 << std::endl;
 
     // Low multiplicity case
     // only adds hit to protoSegment if no hit on layer already; otherwise adds to closeHits
     if (rechits.size() < 9) {
       //      std::cout << "low mult" << std::endl;
       if (isHitNearSegment(h)) {
         if (!hasHitOnLayer(layer)) {
           addHit(h, layer);
         } else {
           closeHits.push_back(h);
         }
       }
 
       // High multiplicity case
       // only adds hit to protoSegment if no hit on layer already AND then refits; otherwise adds to closeHits
     } else {
       //      std::cout << "high mult" << std::endl;
       if (isHitNearSegment(h)) {
         //  std::cout << "near seg" << std::endl;
         if (!hasHitOnLayer(layer)) {
           //      std::cout << "no hit on layer" << std::endl;
           if (addHit(h, layer)) {
             //      std::cout << "update fit" << std::endl;
             updateParameters();
           }
           // Don't change the starting points at this stage !!!
         } else {
           //      std::cout << "already hit on layer" << std::endl;
           closeHits.push_back(h);
           if (layer != layer1 && layer != layer2)
             compareProtoSegment(h, layer);
         }
       }
     }
   }
 
   if (int(protoSegment.size()) < 3)
     return;
   //  std::cout << "final fit" << std::endl;
   updateParameters();
 
   // 2nd pass to remove biases
   // This time, also consider changing the endpoints
   for (ChamberHitContainerCIt i = closeHits.begin(); i != closeHits.end(); ++i) {
     //    std::cout << "2nd pass" << std::endl;
     const CSCRecHit2D* h = *i;
     int layer = (*i)->cscDetId().layer();
     compareProtoSegment(h, layer);
   }
 }

void CSCSegAlgoDF::updateParameters ( void )

private

Definition at line 448 of file CSCSegAlgoDF.cc.

References CSCSegFit::fit(), protoSegment, sfit_, and theChamber.

Referenced by buildSegments(), pruneFromResidual(), and tryAddingHitsToSegment().

                                     {
   // Delete existing CSCSegFit, create a new one and make the fit
   // Uses internal variables - theChamber & protoSegment
 
   //  std::cout << "[CSCSegAlgoDF::updateParameters] about to delete sfit_ = " << sfit_ << std::endl;
   delete sfit_;
   //  std::cout << "[CSCSegAlgoDF::updateParameters] protoSegment.size() = "
   //                                 << protoSegment.size() << std::endl;
   //  std::cout  << "[CSCSegAlgoDF::updateParameters] theChamber = " << theChamber << std::endl;
   sfit_ = new CSCSegFit(theChamber, protoSegment);
   //  std::cout << "[CSCSegAlgoDF::updateParameters] new sfit_ = " << sfit_ << std::endl;
   sfit_->fit();
 }