#include <CSCSegAlgoTC.h>

Inheritance diagram for CSCSegAlgoTC:

Public Types
typedef std::deque< bool >	BoolContainer

typedef std::vector< const CSCRecHit2D * >	ChamberHitContainer

typedef ChamberHitContainer::const_iterator	ChamberHitContainerCIt

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

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

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

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

virtual	~CSCSegAlgoTC ()
	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 *aHit, int layer)
	Utility functions. More...

bool	areHitsCloseInGlobalPhi (const CSCRecHit2D h1, const CSCRecHit2D h2) const

bool	areHitsCloseInLocalX (const CSCRecHit2D h1, const CSCRecHit2D h2) const

void	compareProtoSegment (const CSCRecHit2D *h, int layer)

void	dumpHits (const ChamberHitContainer &rechits) const

void	dumpSegment (const CSCSegment &seg) const

void	flagHitsAsUsed (std::vector< CSCSegFit * >::iterator is, const ChamberHitContainer &rechitsInChamber, BoolContainer &used) const

bool	hasHitOnLayer (int layer) const

void	increaseProtoSegment (const CSCRecHit2D *h, int layer)

bool	isHitNearSegment (const CSCRecHit2D *h) const

bool	isSegmentGood (std::vector< CSCSegFit * >::iterator is, const ChamberHitContainer &rechitsInChamber, BoolContainer &used) const

float	phiAtZ (float z) const

void	pruneTheSegments (const ChamberHitContainer &rechitsInChamber)

bool	replaceHit (const CSCRecHit2D *h, int layer)

void	segmentSort (void)

void	tryAddingHitsToSegment (const ChamberHitContainer &rechits, const ChamberHitContainerCIt i1, const ChamberHitContainerCIt i2)

void	updateParameters (void)

Private Attributes
std::vector< CSCSegFit * >	candidates

float	chi2Max

float	chi2ndfProbMin

bool	debugInfo

float	dPhiFineMax

float	dPhiMax

float	dRPhiFineMax

float	dRPhiMax

int	minLayersApart

const std::string	myName

ChamberHitContainer	proto_segment

int	SegmentSorting

CSCSegFit *	sfit_

const CSCChamber *	theChamber
	Member variables. More...

Detailed Description

This is an alternative algorithm for building endcap muon track segments out of the rechit's in a CSCChamber. cf. CSCSegmentizerSK.
'TC' = 'Tim Cox' = Try (all) Combinations

A CSCSegment isa BasicRecHit4D, and is built from CSCRhit objects, each of which isa BasicRecHit2D.

This class is used by the CSCSegmentRecDet.
Alternative algorithms can be used for the segment building by writing classes like this, and then selecting which one is actually used via the CSCSegmentizerBuilder in CSCDetector.

Ported to CMSSW 2006-04-03: Matte.nosp@m.o.Sa.nosp@m.ni@ce.nosp@m.rn.c.nosp@m.h

Replaced least-squares fit by CSCSegFit - Feb 2015

Definition at line 34 of file CSCSegAlgoTC.h.

Member Typedef Documentation

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

Definition at line 47 of file CSCSegAlgoTC.h.

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

Definition at line 39 of file CSCSegAlgoTC.h.

typedef ChamberHitContainer::const_iterator CSCSegAlgoTC::ChamberHitContainerCIt

Definition at line 40 of file CSCSegAlgoTC.h.

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

Typedefs.

Definition at line 38 of file CSCSegAlgoTC.h.

Constructor & Destructor Documentation

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

explicit

Constructor.

Definition at line 26 of file CSCSegAlgoTC.cc.

References chi2Max, chi2ndfProbMin, debugInfo, dPhiFineMax, dPhiMax, dRPhiFineMax, dRPhiMax, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), LogDebug, minLayersApart, myName, and SegmentSorting.

   : CSCSegmentAlgorithm(ps), sfit_(0), myName("CSCSegAlgoTC") {
   
   debugInfo = ps.getUntrackedParameter<bool>("verboseInfo");
   
   dRPhiMax       = ps.getParameter<double>("dRPhiMax");
   dPhiMax        = ps.getParameter<double>("dPhiMax");
   dRPhiFineMax   = ps.getParameter<double>("dRPhiFineMax");
   dPhiFineMax    = ps.getParameter<double>("dPhiFineMax");
   chi2Max        = ps.getParameter<double>("chi2Max");
   chi2ndfProbMin = ps.getParameter<double>("chi2ndfProbMin");
   minLayersApart = ps.getParameter<int>("minLayersApart");
   SegmentSorting = ps.getParameter<int>("SegmentSorting");
   
   LogDebug("CSC") << myName << " has algorithm cuts set to: \n"
           << "--------------------------------------------------------------------\n"
           << "dRPhiMax     = " << dRPhiMax << '\n'
           << "dPhiMax      = " << dPhiMax << '\n'
           << "dRPhiFineMax = " << dRPhiFineMax << '\n'
           << "dPhiFineMax  = " << dPhiFineMax << '\n'
           << "chi2Max      = " << chi2Max << '\n'
           << "chi2ndfProbMin = " << chi2ndfProbMin << '\n'
           << "minLayersApart = " << minLayersApart << '\n'
           << "SegmentSorting = " << SegmentSorting << std::endl;
 }

virtual CSCSegAlgoTC::~CSCSegAlgoTC ( )

inlinevirtual

Destructor.

Definition at line 52 of file CSCSegAlgoTC.h.

52 {};

Member Function Documentation

bool CSCSegAlgoTC::addHit	(	const CSCRecHit2D *	aHit,
		int	layer
	)

private

Utility functions.

Definition at line 259 of file CSCSegAlgoTC.cc.

References convertSQLiteXML::ok, proto_segment, and updateParameters().

Referenced by buildSegments(), increaseProtoSegment(), and replaceHit().

                                                             {
   
   // 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.
 
   bool ok = true;
   ChamberHitContainer::const_iterator it;
   
   for(it = proto_segment.begin(); it != proto_segment.end(); it++)
     if (((*it)->cscDetId().layer() == layer) && (aHit != *it))
       return false; 
   
   if (ok) {
     proto_segment.push_back(aHit);
     updateParameters();
   } 
   return ok;
 }

bool CSCSegAlgoTC::areHitsCloseInGlobalPhi	(	const CSCRecHit2D *	h1,
		const CSCRecHit2D *	h2
	)		const

private

Return true if the difference in (global) phi of two hits is < dPhiMax

Definition at line 372 of file CSCSegAlgoTC.cc.

References CSCRecHit2D::cscDetId(), dPhiMax, CSCChamber::layer(), CSCDetId::layer(), CSCRecHit2D::localPosition(), LogDebug, M_PI, PV3DBase< T, PVType, FrameType >::phi(), theChamber, and GeomDet::toGlobal().

Referenced by buildSegments().

                                                                                              {
   
   const CSCLayer* l1 = theChamber->layer(h1->cscDetId().layer());
   GlobalPoint gp1 = l1->toGlobal(h1->localPosition());  
   const CSCLayer* l2 = theChamber->layer(h2->cscDetId().layer());
   GlobalPoint gp2 = l2->toGlobal(h2->localPosition());  
   
   float h1p = gp1.phi();
   float h2p = gp2.phi();
   float dphi12 = h1p - h2p;
   
   // Into range [-pi, pi) (phi() returns values in this range)
   if (dphi12 < -M_PI) 
     dphi12 += 2.*M_PI;  
   if (dphi12 >  M_PI) 
     dphi12 -= 2.*M_PI;
   LogDebug("CSC") << "    Hits at global phi= " << h1p << ", " 
           << h2p << " have separation= " << dphi12;
   return (fabs(dphi12) < dPhiMax)? true:false;  // +v
 }

bool CSCSegAlgoTC::areHitsCloseInLocalX	(	const CSCRecHit2D *	h1,
		const CSCRecHit2D *	h2
	)		const

private

Return true if the difference in (local) x of two hits is < dRPhiMax

Definition at line 365 of file CSCSegAlgoTC.cc.

References dRPhiMax, CSCRecHit2D::localPosition(), LogDebug, and x.

Referenced by buildSegments().

                                                                                           {
   float deltaX = (h1->localPosition()-h2->localPosition()).x();
   LogDebug("CSC") << "    Hits at local x= " << h1->localPosition().x() << ", " 
           << h2->localPosition().x() << " have separation= " << deltaX;
   return (fabs(deltaX) < (dRPhiMax))? true:false;   // +v
 }

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

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

Definition at line 58 of file CSCSegAlgoTC.cc.

References funct::abs(), addHit(), areHitsCloseInGlobalPhi(), areHitsCloseInLocalX(), candidates, CSCSegFit::chi2(), CSCSegFit::covarianceMatrix(), CSCRecHit2D::cscDetId(), debugInfo, dumpSegment(), CSCSegFit::hits(), i, cuy::ib, CSCSegFit::intercept(), CSCChamber::layer(), CSCDetId::layer(), CSCSegFit::localdir(), CSCRecHit2D::localPosition(), LogDebug, minLayersApart, myName, GeomDet::position(), proto_segment, pruneTheSegments(), HI_PhotonSkim_cff::rechits, sfit_, groupFilesInBlocks::temp, theChamber, GeomDet::toGlobal(), tryAddingHitsToSegment(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by run().

                                                                                      {
   
   //  ChamberHitContainer rechits = urechits;
   ChamberHitContainer rechits(urechits);
 
   //  edm::LogVerbatim("CSCSegment") << "[CSCSegAlgoTC::buildSegments] start building segments in " << theChamber->id();
   //  edm::LogVerbatim("CSCSegment") << "[CSCSegAlgoTC::buildSegments] size of rechit container = " << rechits.size();
 
   if (rechits.size() < 2) {
     LogDebug("CSC") << myName << ": " << rechits.size() << 
       "  hit(s) in chamber is not enough to build a segment.\n";
     return std::vector<CSCSegment>(); 
   }
 
   LayerIndex layerIndex(rechits.size());
 
   for ( size_t i = 0; i < rechits.size(); ++i ) {
     short ilay = rechits[i]->cscDetId().layer();
     layerIndex[i] = ilay;
     //    edm::LogVerbatim("CSCSegment") << "layerIndex[" << i << "] should   be " << rechits[i]->cscDetId().layer();
     //    edm::LogVerbatim("CSCSegment") << "layerIndex[" << i << "] actually is " << layerIndex[i];
   }
   
   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());
     }    
   }
 
   // Dump rechits after sorting?
   //  if (debugInfo) dumpHits(rechits);
 
   //  if (rechits.size() < 2) {
   //    LogDebug("CSC") << myName << ": " << rechits.size() << 
   //      "  hit(s) in chamber is not enough to build a segment.\n";
   //    return std::vector<CSCSegment>(); 
   //  }
   
   // We have at least 2 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.
   // For now we don't care whether a hit has already been allocated to another segment;
   // we'll sort that out after building all possible segments.
   
   // Choose first hit (as close to IP as possible) h1 and 
   // second hit (as far from IP as possible) h2
   // To do this we iterate over hits in the chamber by layer - pick two layers.
   // @@ Require the two layers are at least 3 layers apart. May need tuning?
   // Then we iterate over hits within each of these layers and pick h1 and h2 from these.
   // If they are 'close enough' we build an empty segment.
   // Then try adding hits to this segment.
   
   // Define buffer for segments we build (at the end we'll sort them somehow, and remove
   // those which share hits with better-quality segments.
   
   
   std::vector<CSCSegment> segments;
 
   sfit_ = 0;
   
   ChamberHitContainerCIt ib = rechits.begin();
   ChamberHitContainerCIt ie = rechits.end();
   
   for (ChamberHitContainerCIt i1 = ib; i1 != ie; ++i1) {
     
     int layer1 = layerIndex[i1-ib];
   
     const CSCRecHit2D* h1 = *i1;
     
     for (ChamberHitContainerCIt i2 = ie-1; i2 != i1; --i2) {
       
       int layer2 = layerIndex[i2-ib];
       
       if (abs(layer2 - layer1) < minLayersApart)
         break;
       
       const CSCRecHit2D* h2 = *i2;
       
       if (areHitsCloseInLocalX(h1, h2) && areHitsCloseInGlobalPhi(h1, h2)) {
         
         proto_segment.clear();
         
         const CSCLayer* l1 = theChamber->layer(h1->cscDetId().layer());
         GlobalPoint gp1 = l1->toGlobal(h1->localPosition());                    
         const CSCLayer* l2 = theChamber->layer(h2->cscDetId().layer());
         GlobalPoint gp2 = l2->toGlobal(h2->localPosition());                    
         LogDebug("CSCSegment") << "start new segment from hits " << "h1: " << gp1 << " - h2: " << gp2 << "\n";
     //  edm::LogVerbatim("CSCSegment") << "start new segment from hits " << "h1: " << gp1 << " - h2: " << gp2;
     //  edm::LogVerbatim("CSCSegment") << "on layers " << layer1 << " and " << layer2;
         
         if ( !addHit(h1, layer1) ) { 
           LogDebug("CSCSegment") << "  failed to add hit h1\n";
           continue;
         }
         
         if ( !addHit(h2, layer2) ) { 
           LogDebug("CSCSegment") << "  failed to add hit h2\n";
           continue;
         }
         
         if ( sfit_ ) tryAddingHitsToSegment(rechits, i1, i2); // can only add hits if there's a segment
         
         // if a segment has been found push back it into the segment collection
         if (proto_segment.empty()) {
           LogDebug("CSCSegment") << "No segment found.\n";
       //      edm::LogVerbatim("CSCSegment") << "No segment found.";
         }   
         else {
       //@@ THIS IS GOING TO BE TRICKY - CREATING MANY FITS ON HEAP
       //@@ BUT MEMBER VARIABLE JUST POINTS TO MOST RECENT ONE
           candidates.push_back( sfit_ ); // store the current fit
           LogDebug("CSCSegment") << "Found a segment.\n";
       //      edm::LogVerbatim("CSCSegment") << "Found a segment.";
         }
       }
     }
   }
 
   //  edm::LogVerbatim("CSCSegment") << "[CSCSegAlgoTC::buildSegments] no. of candidates before pruning = " << candidates.size();
   
   // We've built all possible segments. Now pick the best, non-overlapping subset.
   pruneTheSegments(rechits);
   
   // Create CSCSegments for the surviving candidates
   for(unsigned int i = 0; i < candidates.size(); ++i ) {
     CSCSegFit*sfit = candidates[i];
     //    edm::LogVerbatim("CSCSegment") << "candidate fit " << i+1 << " of " << candidates.size() << " is at " << sfit;
     //    if ( !sfit ) {
     //      edm::LogVerbatim("CSCSegment") << "stored a null pointer for element " << i+1 << " of " << candidates.size();
     //      continue;
     //    }
     CSCSegment temp(sfit->hits(), sfit->intercept(), sfit->localdir(), 
             sfit->covarianceMatrix(), sfit->chi2() );    
     delete sfit;
     segments.push_back(temp);   
     if (debugInfo) dumpSegment( temp );
   }
 
   // reset member variables  
   candidates.clear();
   sfit_ = 0;
 
   //  edm::LogVerbatim("CSCSegment") << "[CSCSegAlgoTC::buildSegments] no. of segments returned = " << segments.size();
 
   return segments;
 }

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

private

Definition at line 319 of file CSCSegAlgoTC.cc.

References CSCSegFit::chi2(), LogDebug, convertSQLiteXML::ok, replaceHit(), and sfit_.

Referenced by tryAddingHitsToSegment().

                                                                       {
 
   // Save copy of current fit
   CSCSegFit* oldfit = new CSCSegFit( *sfit_ );
 
   bool ok = replaceHit(h, layer); // possible new fit
   
   if (ok) {
     LogDebug("CSC") << "    hit in same layer as a hit on segment; try replacing old one..." 
             << " chi2 new: " << sfit_->chi2() << " old: " << oldfit->chi2() << "\n";
   }
   
   if ( ok && (sfit_->chi2() < oldfit->chi2()) ) {
     LogDebug("CSC")  << "    segment with replaced hit is better.\n";
     delete oldfit;  // new fit is better
   }
   else {
     delete sfit_;   // new fit is worse
     sfit_ = oldfit; // restore original fit
   }
 }

void CSCSegAlgoTC::dumpHits ( const ChamberHitContainer & rechits ) const

private

Dump global and local coordinate of each rechit in chamber after sort in z

Definition at line 435 of file CSCSegAlgoTC.cc.

References CSCChamber::layer(), LogDebug, Geom::Phi< T >::phi(), PV3DBase< T, PVType, FrameType >::phi(), theChamber, and GeomDet::toGlobal().

                                                                     {
   
   // Dump positions of RecHit's in each CSCChamber
   ChamberHitContainerCIt it;
   
   for(it=rechits.begin(); it!=rechits.end(); it++) {
     
     const CSCLayer* l1 = theChamber->layer((*it)->cscDetId().layer());
     GlobalPoint gp1 = l1->toGlobal((*it)->localPosition()); 
     
     LogDebug("CSC") << "Global pos.: " << gp1 << ", phi: " << gp1.phi() << ". Local position: "
             << (*it)->localPosition() << ", phi: "
             << (*it)->localPosition().phi() << ". Layer: "
             << (*it)->cscDetId().layer() << "\n";
   } 
 }

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

private

Definition at line 602 of file CSCSegAlgoTC.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 CSCSegAlgoTC::flagHitsAsUsed	(	std::vector< CSCSegFit * >::iterator	is,
		const ChamberHitContainer &	rechitsInChamber,
		BoolContainer &	used
	)		const

private

Flag hits on segment as used

Definition at line 506 of file CSCSegAlgoTC.cc.

References cuy::ib.

Referenced by pruneTheSegments().

                                                                               {
 
   // Flag hits on segment as used
 
   ChamberHitContainerCIt ib = rechitsInChamber.begin();
   ChamberHitContainer hits = (*seg)->hits();
    
   for(size_t ish = 0; ish < hits.size(); ++ish) {
     for(ChamberHitContainerCIt iu = ib; iu != rechitsInChamber.end(); ++iu)
       if( hits[ish] == (*iu)) 
     used[iu-ib] = true;
   }
 }

bool CSCSegAlgoTC::hasHitOnLayer ( int layer ) const

private

Definition at line 423 of file CSCSegAlgoTC.cc.

References proto_segment.

Referenced by tryAddingHitsToSegment().

                                                 {
   
   // Is there is already a hit on this layer?
   ChamberHitContainerCIt it;
   
   for(it = proto_segment.begin(); it != proto_segment.end(); it++)
     if ((*it)->cscDetId().layer() == layer)
       return true; 
   
   return false;
 }

void CSCSegAlgoTC::increaseProtoSegment	(	const CSCRecHit2D *	h,
		int	layer
	)

private

Definition at line 341 of file CSCSegAlgoTC.cc.

References addHit(), CSCSegFit::chi2(), chi2Max, LogDebug, CSCSegFit::ndof(), convertSQLiteXML::ok, and sfit_.

Referenced by tryAddingHitsToSegment().

                                                                        {
   
   // save copy of input fit
   CSCSegFit* oldfit = new CSCSegFit( *sfit_ );
 
   bool ok = addHit(h, layer); // possible new fit
   
   if (ok) {
     LogDebug("CSC") << "    hit in new layer: added to segment, new chi2: " 
             << sfit_->chi2() << "\n";
   }
   
   //  int ndf = 2*proto_segment.size() - 4;
   
   if (ok && ((sfit_->chi2() <= 0) || (sfit_->chi2()/sfit_->ndof() < chi2Max))) {
     LogDebug("CSC") << "    segment with added hit is good.\n" ;        
     delete oldfit; // new fit is better
   } 
   else {
     delete sfit_;   // new fit is worse
     sfit_ = oldfit; // restore original fit
   }         
 }       

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

private

Return true if hit is near segment. 'Near' means deltaphi and rxy*deltaphi are within ranges specified by config parameters dPhiFineMax and dRPhiFineMax, where rxy = sqrt(x**2+y**2) of the hit in global coordinates.

Definition at line 393 of file CSCSegAlgoTC.cc.

References CSCRecHit2D::cscDetId(), dPhiFineMax, dRPhiFineMax, AnalysisDataFormats_SUSYBSMObjects::hp, CSCChamber::layer(), CSCDetId::layer(), CSCRecHit2D::localPosition(), LogDebug, M_PI, PV3DBase< T, PVType, FrameType >::perp(), PV3DBase< T, PVType, FrameType >::phi(), phiAtZ(), theChamber, GeomDet::toGlobal(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by tryAddingHitsToSegment().

                                                               {
   
   // Is hit near segment? 
   // Requires deltaphi and rxy*deltaphi within ranges specified
   // in orcarc, or by default, where rxy=sqrt(x**2+y**2) of hit itself.
   // Note that to make intuitive cuts on delta(phi) one must work in
   // phi range (-pi, +pi] not [0, 2pi
 
   const CSCLayer* l1 = theChamber->layer(h->cscDetId().layer());
   GlobalPoint hp = l1->toGlobal(h->localPosition());    
   
   float hphi = hp.phi();               // in (-pi, +pi]
   if (hphi < 0.) 
     hphi += 2.*M_PI;                   // into range [0, 2pi)
   float sphi = phiAtZ(hp.z());         // in [0, 2*pi)
   float phidif = sphi-hphi;
   if (phidif < 0.) 
     phidif += 2.*M_PI;                 // into range [0, 2pi)
   if (phidif > M_PI) 
     phidif -= 2.*M_PI;                 // into range (-pi, pi]
   
   float dRPhi = fabs(phidif)*hp.perp();
   LogDebug("CSC") << "    is hit at phi_h= " << hphi << " near segment phi_seg= " << sphi 
           << "? is " << dRPhi << "<" << dRPhiFineMax << " ? " 
           << " and is |" << phidif << "|<" << dPhiFineMax << " ?";
   
   return ((dRPhi < dRPhiFineMax) && 
       (fabs(phidif) < dPhiFineMax))? true:false;  // +v
 }

bool CSCSegAlgoTC::isSegmentGood	(	std::vector< CSCSegFit * >::iterator	is,
		const ChamberHitContainer &	rechitsInChamber,
		BoolContainer &	used
	)		const

private

Return true if segment is good. In this algorithm, this means it shares no hits with any other segment. If "SegmentSort=2" also require a minimal chi2 probability of "chi2ndfProbMin".

Definition at line 454 of file CSCSegAlgoTC.cc.

References chi2ndfProbMin, ChiSquaredProbability(), cuy::ib, nhits, and SegmentSorting.

Referenced by pruneTheSegments().

                                                                                          {
   
   // Apply any selection cuts to segment
   
   // 1) Require a minimum no. of hits
   //   (@@ THESE VALUES SHOULD BECOME PARAMETERS?)
   
   // 2) Ensure no hits on segment are already assigned to another segment
   //    (typically of higher quality)
 
   size_t iadd = (rechitsInChamber.size() > 20 )?  1 : 0;  
   
   size_t nhits = (*seg)->nhits();
 
   if (nhits < 3 + iadd)
     return false;
   
   // Additional part of alternative segment selection scheme: reject
   // segments with a chi2 probability of less than chi2ndfProbMin. Relies on list 
   // being sorted with "SegmentSorting == 2", that is first by nrechits and then 
   // by chi2prob in subgroups of same no. of rechits.
 
   if( SegmentSorting == 2 ){
     double chi2t = (*seg)->chi2();
     double ndoft = 2*nhits - 4 ;
     if( chi2t > 0 && ndoft > 0 )  {
       if ( ChiSquaredProbability(chi2t,ndoft) < chi2ndfProbMin ) { 
     return false;
       }
     }
     else {
       return false; 
     }
   }
   
   ChamberHitContainer hits_ = (*seg)->hits();
 
   for(size_t ish = 0; ish < nhits; ++ish) {
     
     ChamberHitContainerCIt ib = rechitsInChamber.begin();
     
     for(ChamberHitContainerCIt ic = ib; ic != rechitsInChamber.end(); ++ic) {
 
       if((hits_[ish] == (*ic)) && used[ic-ib])
     return false;
     }
   }
   
   return true;
 }

float CSCSegAlgoTC::phiAtZ ( float z ) const

private

Definition at line 303 of file CSCSegAlgoTC.cc.

References f, CSCSegFit::intercept(), CSCChamber::layer(), CSCSegFit::localdir(), M_PI, phi, sfit_, theChamber, GeomDet::toGlobal(), x, PV3DBase< T, PVType, FrameType >::x(), y, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by isHitNearSegment().

                                         {
   
   // Returns a phi in [ 0, 2*pi )
   const CSCLayer* l1 = theChamber->layer(1);
   GlobalPoint gp = l1->toGlobal(sfit_->intercept());    
   GlobalVector gv = l1->toGlobal(sfit_->localdir());    
   
   float x = gp.x() + (gv.x()/gv.z())*(z - gp.z());
   float y = gp.y() + (gv.y()/gv.z())*(z - gp.z());
   float phi = atan2(y, x);
   if (phi < 0.f) 
     phi += 2. * M_PI;
   
   return phi ;
 }

void CSCSegAlgoTC::pruneTheSegments ( const ChamberHitContainer & rechitsInChamber )

private

Order segments by quality (chi2/#hits) and select the best, requiring that they have unique hits.

Definition at line 521 of file CSCSegAlgoTC.cc.

References candidates, flagHitsAsUsed(), isSegmentGood(), LogDebug, and segmentSort().

Referenced by buildSegments().

                                                                                {
   
   // Sort the segment store according to segment 'quality' (chi2/#hits ?) and
   // remove any segments which contain hits assigned to higher-quality segments.
 
   if (candidates.empty()) 
     return;
   
   // Initialize flags that a given hit has been allocated to a segment
   BoolContainer used(rechitsInChamber.size(), false);
   
   // Sort by chi2/#hits
   segmentSort();
   
   // Select best quality segments, requiring hits are assigned to just one segment
   // Want to erase the bad segments, so the iterator must be incremented
   // inside the loop, and only when the erase is not called
 
   std::vector<CSCSegFit*>::iterator is;
 
   for (is = candidates.begin();  is !=  candidates.end(); ) {
     
     bool goodSegment = isSegmentGood(is, rechitsInChamber, used);
     
     if (goodSegment) {
       LogDebug("CSC") << "Accepting segment: ";
       flagHitsAsUsed(is, rechitsInChamber, used);
       ++is;
     }
     else {
       LogDebug("CSC") << "Rejecting segment: ";
       delete *is; // delete the CSCSegFit* 
       is = candidates.erase(is); // erase the element in container
     }
   }
 }

bool CSCSegAlgoTC::replaceHit	(	const CSCRecHit2D *	h,
		int	layer
	)

private

Definition at line 279 of file CSCSegAlgoTC.cc.

References addHit(), and proto_segment.

Referenced by compareProtoSegment().

                                                              {
 
   // replace a hit from a layer 
   ChamberHitContainer::iterator it;
   for (it = proto_segment.begin(); it != proto_segment.end();) {
     if ((*it)->cscDetId().layer() == layer)
       it = proto_segment.erase(it);
     else
       ++it;   
   }
   
   return addHit(h, layer);
 
 }

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

Here we must implement the algorithm

Definition at line 52 of file CSCSegAlgoTC.cc.

References buildSegments(), and theChamber.

                                                                                                       {
   theChamber = aChamber; 
 
   return buildSegments(rechits); 
 }

void CSCSegAlgoTC::segmentSort ( void )

private

Sort criterion for segment quality, for use in pruneTheSegments.

Definition at line 558 of file CSCSegAlgoTC.cc.

References candidates, ChiSquaredProbability(), i, j, LogDebug, ndof, SegmentSorting, findQualityFiles::size, and groupFilesInBlocks::temp.

Referenced by pruneTheSegments().

                                {
   
   // The segment collection is sorted according to e.g. chi2/#hits 
   
   for(size_t i=0; i<candidates.size()-1; ++i) {
     for(size_t j=i+1; j<candidates.size(); ++j) {
       
       size_t ni = (candidates[i]->hits()).size();
       size_t nj = (candidates[j]->hits()).size();
       
       // Sort criterion: first sort by no. of rechits, then in groups of rechits by chi2prob
       if( SegmentSorting == 2 ){ 
     if ( nj > ni ) { // sort by no. of rechits
       CSCSegFit* temp = candidates[j];
       candidates[j] = candidates[i];
       candidates[i] = temp;
     }
     // sort by chi2 probability in subgroups with equal nr of rechits
     //  if(chi2s[i] != 0. && 2*n2-4 > 0 ) {
         if( candidates[i]->chi2() > 0 && candidates[i]->ndof() > 0 ) {
       if( nj == ni && 
           ( ChiSquaredProbability( candidates[i]->chi2(),(double)(candidates[i]->ndof()) ) < 
             ChiSquaredProbability( candidates[j]->chi2(),(double)(candidates[j]->ndof()) ) ) 
             ){
       CSCSegFit* temp = candidates[j];
       candidates[j] = candidates[i];
       candidates[i] = temp;
       }
     }
       }
       else if( SegmentSorting == 1 ){
     if ((candidates[i]->chi2()/ni) > (candidates[j]->chi2()/nj)) {
       CSCSegFit* temp = candidates[j];
       candidates[j] = candidates[i];
       candidates[i] = temp;
     }
       }
       else{
           LogDebug("CSC") << "No valid segment sorting specified. Algorithm misconfigured! \n";
       }
     }
   }     
 }

void CSCSegAlgoTC::tryAddingHitsToSegment	(	const ChamberHitContainer &	rechits,
		const ChamberHitContainerCIt	i1,
		const ChamberHitContainerCIt	i2
	)

private

Try adding non-used hits to segment

Definition at line 214 of file CSCSegAlgoTC.cc.

References compareProtoSegment(), CSCRecHit2D::cscDetId(), h, hasHitOnLayer(), i, cuy::ib, increaseProtoSegment(), isHitNearSegment(), CSCChamber::layer(), CSCDetId::layer(), CSCRecHit2D::localPosition(), LogDebug, proto_segment, theChamber, and GeomDet::toGlobal().

Referenced by buildSegments().

                                                                            {
   
   // Iterate over the layers with hits in the chamber
   // Skip the layers containing the segment endpoints
   // Test each hit on the other layers to see if it is near the segment
   // 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 (by swap)
   //    - if not, copy the segment, add the hit. If the new chi2/dof is still satisfactory
   //      then replace the original segment with the new one (by swap)
 
   ChamberHitContainerCIt ib = rechits.begin();
   ChamberHitContainerCIt ie = rechits.end();
   
   for (ChamberHitContainerCIt i = ib; i != ie; ++i) {
     
     if ( i == i1 || i == i2 ) 
       continue;
     
     int layer = (*i)->cscDetId().layer();
     const CSCRecHit2D* h = *i;
     
     if (isHitNearSegment(h)) {
       
       const CSCLayer* l1 = theChamber->layer(h->cscDetId().layer());
       GlobalPoint gp1 = l1->toGlobal(h->localPosition());       
       LogDebug("CSC") << "    hit at global " << gp1 << " is near segment\n.";
       
       if (hasHitOnLayer(layer)) {
     if (proto_segment.size() <= 2) {
       LogDebug("CSC") << "    " << proto_segment.size() 
               << " hits on segment...skip hit on same layer.\n";
       continue;
     }
     
     compareProtoSegment(h, layer);
       } 
       else
     increaseProtoSegment(h, layer);
     }   // h & seg close
   }   // i
 }

void CSCSegAlgoTC::updateParameters ( void )

private

Definition at line 294 of file CSCSegAlgoTC.cc.

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

Referenced by addHit().

                                         {
   
   //@@ DO NOT DELETE EXISTING FIT SINCE WE SAVE IT!!
   //  delete sfit_;
   sfit_ = new CSCSegFit( theChamber, proto_segment );
   sfit_->fit();
 
 }

Member Data Documentation

std::vector<CSCSegFit*> CSCSegAlgoTC::candidates

private

Definition at line 148 of file CSCSegAlgoTC.h.

Referenced by buildSegments(), pruneTheSegments(), and segmentSort().

float CSCSegAlgoTC::chi2Max

private

max segment chi squared

Definition at line 152 of file CSCSegAlgoTC.h.

Referenced by CSCSegAlgoTC(), and increaseProtoSegment().

float CSCSegAlgoTC::chi2ndfProbMin

private

min segment chi squared probability. Used ONLY if SegmentSorting is chosen to be 2

Definition at line 157 of file CSCSegAlgoTC.h.

Referenced by CSCSegAlgoTC(), and isSegmentGood().

bool CSCSegAlgoTC::debugInfo

private

Definition at line 193 of file CSCSegAlgoTC.h.

Referenced by buildSegments(), and CSCSegAlgoTC().

float CSCSegAlgoTC::dPhiFineMax

private

max hit deviation in global phi from the segment axis. Function hitNearSegment requires abs(deltaphi) < dPhiFineMax.

Definition at line 167 of file CSCSegAlgoTC.h.

Referenced by CSCSegAlgoTC(), and isHitNearSegment().

float CSCSegAlgoTC::dPhiMax

private

max distance in global phi between hits in one segment

Definition at line 176 of file CSCSegAlgoTC.h.

Referenced by areHitsCloseInGlobalPhi(), and CSCSegAlgoTC().

float CSCSegAlgoTC::dRPhiFineMax

private

max hit deviation in r-phi from the segment axis. Function hitNearSegment requires rxy*abs(deltaphi) < dRPhiFineMax.

Definition at line 162 of file CSCSegAlgoTC.h.

Referenced by CSCSegAlgoTC(), and isHitNearSegment().

float CSCSegAlgoTC::dRPhiMax

private

max distance in local x between hits in one segment @ The name is historical!

Definition at line 172 of file CSCSegAlgoTC.h.

Referenced by areHitsCloseInLocalX(), and CSCSegAlgoTC().

int CSCSegAlgoTC::minLayersApart

private

Require end-points of segment are at least minLayersApart

Definition at line 180 of file CSCSegAlgoTC.h.

Referenced by buildSegments(), and CSCSegAlgoTC().

const std::string CSCSegAlgoTC::myName

private

Name of this class

Definition at line 192 of file CSCSegAlgoTC.h.

Referenced by buildSegments(), and CSCSegAlgoTC().

ChamberHitContainer CSCSegAlgoTC::proto_segment

private

Definition at line 142 of file CSCSegAlgoTC.h.

Referenced by addHit(), buildSegments(), hasHitOnLayer(), replaceHit(), tryAddingHitsToSegment(), and updateParameters().

int CSCSegAlgoTC::SegmentSorting

private

Select which segment sorting to use (the higher the segment is in the list, the better the segment is supposed to be): if value is ==1: Sort segments by Chi2/(#hits on segment) if value is ==2: Sort segments first by #hits on segment, then by Chi2Probability(Chi2/ndf)

Definition at line 188 of file CSCSegAlgoTC.h.

Referenced by CSCSegAlgoTC(), isSegmentGood(), and segmentSort().

CSCSegFit* CSCSegAlgoTC::sfit_

private

Definition at line 145 of file CSCSegAlgoTC.h.

Referenced by buildSegments(), compareProtoSegment(), increaseProtoSegment(), phiAtZ(), and updateParameters().