#include <DTCombinatorialPatternReco.h>

Inheritance diagram for DTCombinatorialPatternReco:

Classes
class	TriedPattern

Public Types
typedef boost::unordered_set< TriedPattern >	TriedPatterns

Public Member Functions
virtual std::string	algoName () const
	return the algo name More...

	DTCombinatorialPatternReco (const edm::ParameterSet &pset)
	Constructor. More...

virtual edm::OwnVector< DTSLRecSegment2D >	reconstruct (const DTSuperLayer *sl, const std::vector< DTRecHit1DPair > &hits)
	this function is called in the producer More...

virtual void	setES (const edm::EventSetup &setup)

virtual	~DTCombinatorialPatternReco ()
	Destructor. More...

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

virtual	~DTRecSegment2DBaseAlgo ()
	Destructor. More...

Private Member Functions
DTSegmentCand *	buildBestSegment (std::vector< DTSegmentCand::AssPoint > &assHits, const DTSuperLayer *sl)

void	buildPointsCollection (std::vector< DTSegmentCand::AssPoint > &points, std::deque< std::shared_ptr< DTHitPairForFit > > &pointsNoLR, std::vector< DTSegmentCand * > &candidates, const DTSuperLayer *sl)

std::vector< DTSegmentCand * >	buildSegments (const DTSuperLayer *sl, const std::vector< std::shared_ptr< DTHitPairForFit >> &hits)

bool	checkDoubleCandidates (std::vector< DTSegmentCand * > &segs, DTSegmentCand *seg)

std::vector< DTSegmentCand::AssPoint >	findCompatibleHits (const LocalPoint &pos, const LocalVector &dir, const std::vector< std::shared_ptr< DTHitPairForFit >> &hits)

std::vector< std::shared_ptr< DTHitPairForFit > >	initHits (const DTSuperLayer *sl, const std::vector< DTRecHit1DPair > &hits)

Private Attributes
bool	debug

std::string	theAlgoName

double	theAlphaMaxPhi

double	theAlphaMaxTheta

DTSegmentCleaner *	theCleaner

edm::ESHandle< DTGeometry >	theDTGeometry

unsigned int	theMaxAllowedHits

TriedPatterns	theTriedPattern

DTSegmentUpdator *	theUpdator

bool	usePairs

Friends
class	DTCombinatorialPatternReco4D

Detailed Description

Algo for reconstructing 2d segment in DT using a combinatorial approach

Author: Stefano Lacaprara - INFN Legnaro stefa.nosp@m.no.l.nosp@m.acapr.nosp@m.ara@.nosp@m.pd.in.nosp@m.fn.i.nosp@m.t; Riccardo Bellan - INFN TO ricca.nosp@m.rdo..nosp@m.bella.nosp@m.n@ce.nosp@m.rn.ch

Definition at line 41 of file DTCombinatorialPatternReco.h.

Member Typedef Documentation

typedef boost::unordered_set<TriedPattern> DTCombinatorialPatternReco::TriedPatterns

Definition at line 144 of file DTCombinatorialPatternReco.h.

Constructor & Destructor Documentation

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

Constructor.

Definition at line 34 of file DTCombinatorialPatternReco.cc.

References debug, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), theAlphaMaxPhi, theAlphaMaxTheta, theCleaner, theMaxAllowedHits, theUpdator, and usePairs.

                                                                                   : 
 DTRecSegment2DBaseAlgo(pset), theAlgoName("DTCombinatorialPatternReco")
 {
   theMaxAllowedHits = pset.getParameter<unsigned int>("MaxAllowedHits"); // 100
   theAlphaMaxTheta = pset.getParameter<double>("AlphaMaxTheta");// 0.1 ;
   theAlphaMaxPhi = pset.getParameter<double>("AlphaMaxPhi");// 1.0 ;
   debug = pset.getUntrackedParameter<bool>("debug"); //true;
   theUpdator = new DTSegmentUpdator(pset);
   theCleaner = new DTSegmentCleaner(pset);
   string theHitAlgoName = pset.getParameter<string>("recAlgo");
   usePairs = !(theHitAlgoName=="DTNoDriftAlgo");
 }

DTCombinatorialPatternReco::~DTCombinatorialPatternReco ( )

virtual

Destructor.

Definition at line 48 of file DTCombinatorialPatternReco.cc.

References theCleaner, and theUpdator.

                                                         {
   delete theUpdator;
   delete theCleaner;
 }

Member Function Documentation

virtual std::string DTCombinatorialPatternReco::algoName ( void ) const

inlinevirtual

return the algo name

Implements DTRecSegment2DBaseAlgo.

Definition at line 59 of file DTCombinatorialPatternReco.h.

References GeneralSetup::setup().

59 { return theAlgoName; }

DTCombinatorialPatternReco::theAlgoName

std::string theAlgoName

Definition: DTCombinatorialPatternReco.h:100

DTSegmentCand * DTCombinatorialPatternReco::buildBestSegment	(	std::vector< DTSegmentCand::AssPoint > &	assHits,
		const DTSuperLayer *	sl
	)

private

Definition at line 313 of file DTCombinatorialPatternReco.cc.

References buildPointsCollection(), objects.IsoTrackAnalyzer::candidates, gather_cfg::cout, debug, MTVHistoProducerAlgoForTrackerBlock_cfi::minChi2, and DTEnums::undefLR.

Referenced by buildSegments().

                                                                      {
   if (hits.size()<3) {
     //cout << "buildBestSegment: hits " << hits.size()<< endl;
     return 0; // a least 3 point
   }
 
   // hits with defined LR
   vector<DTSegmentCand::AssPoint> points;
 
   // without: I store both L and R, a deque since I need front insertion and
   // deletion
   deque<std::shared_ptr<DTHitPairForFit> > pointsNoLR; 
 
   // first add only the hits with LR assigned
   for (vector<DTSegmentCand::AssPoint>::const_iterator hit=hits.begin();
        hit!=hits.end(); ++hit) {
     if ((*hit).second != DTEnums::undefLR) {
       points.push_back(*hit);
     } else { // then also for the undef'd one
       pointsNoLR.push_back((*hit).first);
     }
   }
 
   if(debug) {
     cout << "points " << points.size() << endl;
     cout << "pointsNoLR " << pointsNoLR.size() << endl;
   }
 
   // build all possible candidates using L/R ambiguity
   vector<DTSegmentCand*> candidates ;
 
   buildPointsCollection(points, pointsNoLR, candidates, sl);
 
   if(debug)
     cout << "candidates " << candidates.size() << endl;
 
   // so now I have build a given number of segments, I should find the best one,
   // by #hits and chi2.
   vector<DTSegmentCand*>::const_iterator bestCandIter = candidates.end();
   double minChi2=999999.;
   unsigned int maxNumHits=0;
   for (vector<DTSegmentCand*>::const_iterator iter=candidates.begin();
        iter!=candidates.end(); ++iter) {
     if ((*iter)->nHits()==maxNumHits && (*iter)->chi2()<minChi2) {
       minChi2=(*iter)->chi2();
       bestCandIter=iter;
     } else if ((*iter)->nHits()>maxNumHits) {
       maxNumHits=(*iter)->nHits();
       minChi2=(*iter)->chi2();
       bestCandIter=iter;
     }
   }
 
   // delete all candidates but the best one!
   for (vector<DTSegmentCand*>::iterator iter=candidates.begin();
        iter!=candidates.end(); ++iter) if (iter!=bestCandIter) delete *iter;
 
        // return the best candate if any
        if (bestCandIter != candidates.end()) {
          return (*bestCandIter);
        }
        return 0;
 }

void DTCombinatorialPatternReco::buildPointsCollection	(	std::vector< DTSegmentCand::AssPoint > &	points,
		std::deque< std::shared_ptr< DTHitPairForFit > > &	pointsNoLR,
		std::vector< DTSegmentCand * > &	candidates,
		const DTSuperLayer *	sl
	)

private

build collection of compatible hits for L/R hits: the candidates is updated with the segment candidates found

Definition at line 379 of file DTCombinatorialPatternReco.cc.

References popcon2dropbox::copy(), gather_cfg::cout, debug, DTSegmentUpdator::fit(), DTEnums::Left, DTEnums::Right, and theUpdator.

Referenced by buildBestSegment().

                                                                           {
 
   if(debug) {
     cout << "buildPointsCollection " << endl;
     cout << "points: " << points.size() << " NOLR: " << pointsNoLR.size()<< endl;
   }
   if (pointsNoLR.size()>0) { // still unassociated points!
     std::shared_ptr<DTHitPairForFit> unassHit = pointsNoLR.front();
     // try with the right
     if(debug)
       cout << "Right hit" << endl;
     points.push_back(DTSegmentCand::AssPoint(unassHit, DTEnums::Right));
     pointsNoLR.pop_front();
     buildPointsCollection(points, pointsNoLR, candidates, sl);
     pointsNoLR.push_front((unassHit));
     points.pop_back();
 
     // try with the left
     if(debug)
       cout << "Left hit" << endl;
     points.push_back(DTSegmentCand::AssPoint(unassHit, DTEnums::Left));
     pointsNoLR.pop_front();
     buildPointsCollection(points, pointsNoLR, candidates, sl);
     pointsNoLR.push_front((unassHit));
     points.pop_back();
   } else { // all associated
 
     if(debug) {
       cout << "The Hits were" << endl;
       copy(points.begin(), points.end(),
            ostream_iterator<DTSegmentCand::AssPoint>(std::cout));
       cout << "----" << endl;
       cout << "All associated " << endl;
     }
     DTSegmentCand::AssPointCont pointsSet;
 
     // for (vector<DTSegmentCand::AssPoint>::const_iterator point=points.begin();
     //      point!=points.end(); ++point) 
     pointsSet.insert(points.begin(),points.end());
 
     if(debug) {
       cout << "The Hits are" << endl;
       copy(pointsSet.begin(), pointsSet.end(),
            ostream_iterator<DTSegmentCand::AssPoint>(std::cout));
       cout << "----" << endl;
     }
 
     DTSegmentCand* newCand = new DTSegmentCand(pointsSet,sl);
     if (theUpdator->fit(newCand,0,0)) candidates.push_back(newCand);
     else delete newCand; // bad seg, too few hits
   }
 }

vector< DTSegmentCand * > DTCombinatorialPatternReco::buildSegments	(	const DTSuperLayer *	sl,
		const std::vector< std::shared_ptr< DTHitPairForFit >> &	hits
	)

private

get two hits in different layers and see if there are other / hits

Definition at line 101 of file DTCombinatorialPatternReco.cc.

References buildBestSegment(), checkDoubleCandidates(), DTSegmentCleaner::clean(), gather_cfg::cout, debug, findCompatibleHits(), DTSegmentCand::good(), hfClusterShapes_cfi::hits, DTSuperLayer::id(), listHistos::IP, DTEnums::Left, mps_fire::result, DTEnums::Right, theAlphaMaxPhi, theAlphaMaxTheta, theCleaner, theMaxAllowedHits, PV3DBase< T, PVType, FrameType >::theta(), theTriedPattern, GeomDet::toGlobal(), and csvLumiCalc::unit.

Referenced by DTCombinatorialPatternReco4D::buildPhiSuperSegmentsCandidates(), and reconstruct().

                                                                                                 {
   // clear the patterns tried
   if (debug) {
       cout << "theTriedPattern.size is " << theTriedPattern.size() << "\n";
   }
   theTriedPattern.clear();
 
   typedef vector<std::shared_ptr<DTHitPairForFit>> hitCont;
   typedef hitCont::const_iterator  hitIter;
   vector<DTSegmentCand*> result;
   
   if(debug) {
     cout << "buildSegments: " << sl->id() << " nHits " << hits.size() << endl;
     for (vector<std::shared_ptr<DTHitPairForFit>>::const_iterator hit=hits.begin();
          hit!=hits.end(); ++hit) cout << **hit<< endl;
   }
 
   // 10-Mar-2004 SL
   // put a protection against heavily populated chambers, for which the segment
   // building could lead to infinite memory usage...
   if (hits.size() > theMaxAllowedHits ) {
     if(debug) {
       cout << "Warning: this SuperLayer " << sl->id() << " has too many hits : "
         << hits.size() << " max allowed is " << theMaxAllowedHits << endl;
       cout << "Skipping segment reconstruction... " << endl;
     }
     return result;
   }
 
   //  compatible with them
   for (hitCont::const_iterator firstHit=hits.begin(); firstHit!=hits.end();
        ++firstHit) {
     for (hitCont::const_reverse_iterator lastHit=hits.rbegin(); 
          (*lastHit)!=(*firstHit); ++lastHit) {
       //if ( (*lastHit)->id().layerId() == (*firstHit)->id().layerId() ) continue; // hits must be in different layers!
       // hits must nor in the same nor in adiacent layers
       if ( fabs((*lastHit)->id().layerId()-(*firstHit)->id().layerId())<=1 ) continue;
       if(debug) {
         cout << "Selected these two hits pair " << endl;
         cout << "First " << *(*firstHit) << " Layer Id: " << (*firstHit)->id().layerId() << endl;
         cout << "Last "  << *(*lastHit)  << " Layer Id: " << (*lastHit)->id().layerId()  << endl;
       }
 
       GlobalPoint IP;
       float DAlphaMax;
       if ((sl->id()).superlayer()==2)  // Theta SL
         DAlphaMax=theAlphaMaxTheta;
       else // Phi SL
         DAlphaMax=theAlphaMaxPhi;
 
       DTEnums::DTCellSide codes[2]={DTEnums::Right, DTEnums::Left};
       for (int firstLR=0; firstLR<2; ++firstLR) {
         for (int lastLR=0; lastLR<2; ++lastLR) {
           // TODO move the global transformation in the DTHitPairForFit class
           // when it will be moved I will able to remove the sl from the input parameter
           GlobalPoint gposFirst=sl->toGlobal( (*firstHit)->localPosition(codes[firstLR]) );
           GlobalPoint gposLast= sl->toGlobal( (*lastHit)->localPosition(codes[lastLR]) );
 
           GlobalVector gvec=gposLast-gposFirst;
           GlobalVector gvecIP=gposLast-IP;
 
           // difference in angle measured
           float DAlpha=fabs(gvec.theta()-gvecIP.theta());
 
           // cout << "DAlpha " << DAlpha << endl;
           if (DAlpha<DAlphaMax) {
 
             // create a segment hypotesis
             // I don't need a true segment, just direction and position
             LocalPoint posIni = (*firstHit)->localPosition(codes[firstLR]);
             LocalVector dirIni = 
               ((*lastHit)->localPosition(codes[lastLR])-posIni).unit();
 
             // search for other compatible hits, with or without the L/R solved
             vector<DTSegmentCand::AssPoint> assHits = findCompatibleHits(posIni, dirIni, hits);
             if(debug) 
               cout << "compatible hits " << assHits.size() << endl;
 
             // get the best segment with these hits: it's just one! 
             // (is it correct?)
             DTSegmentCand* seg = buildBestSegment(assHits, sl);
 
             if (seg) {
               if(debug) 
                 cout << "segment " << *seg<< endl;
 
               // check if the chi2 and #hits are ok
               if (!seg->good()) {
                 delete seg;
               } else { 
 
                 // remove duplicated segments (I know, would be better to do it before the
                 // fit...)
                 if (checkDoubleCandidates(result,seg)) {
                   // add to the vector of hypotesis
                   result.push_back(seg);
                   if(debug) 
                     cout << "result is now " << result.size() << endl;
                 } else { // delete it!
                   delete seg;
                   if(debug) 
                     cout << "already existing" << endl;
                 }
               }
             }
           }
         }
       }
     }
   }
   if (debug) {
     for (vector<DTSegmentCand*>::const_iterator seg=result.begin();
          seg!=result.end(); ++seg) 
       cout << *(*seg) << endl;
   }
 
   // now I have a couple of segment hypotesis, should check for ghost
   result = theCleaner->clean(result);
   if (debug) {
     cout << "result no ghost  " << result.size() << endl;
     for (vector<DTSegmentCand*>::const_iterator seg=result.begin();
          seg!=result.end(); ++seg) 
       cout << *(*seg) << endl;
   }
 
   return result;
 }

bool DTCombinatorialPatternReco::checkDoubleCandidates	(	std::vector< DTSegmentCand * > &	segs,
		DTSegmentCand *	seg
	)

private

Definition at line 436 of file DTCombinatorialPatternReco.cc.

Referenced by buildSegments().

                                                                       {
   for (vector<DTSegmentCand*>::iterator cand=cands.begin();
        cand!=cands.end(); ++cand) 
     if (*(*cand)==*seg) return false;
   return true;
 }

vector< DTSegmentCand::AssPoint > DTCombinatorialPatternReco::findCompatibleHits	(	const LocalPoint &	pos,
		const LocalVector &	dir,
		const std::vector< std::shared_ptr< DTHitPairForFit >> &	hits
	)

private

Definition at line 233 of file DTCombinatorialPatternReco.cc.

References DTCombinatorialPatternReco::TriedPattern::begin(), popcon2dropbox::copy(), gather_cfg::cout, debug, DTCombinatorialPatternReco::TriedPattern::end(), hfClusterShapes_cfi::hits, DTEnums::Left, DTCombinatorialPatternReco::TriedPattern::push_back(), mps_fire::result, DTEnums::Right, theTriedPattern, printsummarytable::tried, DTEnums::undefLR, and usePairs.

Referenced by buildSegments().

                                                                                                    {
   if (debug) cout << "Pos: " << posIni << " Dir: "<< dirIni << endl;
   vector<DTSegmentCand::AssPoint> result;
 
   // counter to early-avoid double counting in hits pattern
   TriedPattern tried;
   int nCompatibleHits=0;
 
   typedef vector<std::shared_ptr<DTHitPairForFit>> hitCont;
   typedef hitCont::const_iterator  hitIter;
   for (hitIter hit=hits.begin(); hit!=hits.end(); ++hit) {
     pair<bool,bool> isCompatible = (*hit)->isCompatible(posIni, dirIni);
     if (debug) 
       cout << "isCompatible " << isCompatible.first << " " <<
         isCompatible.second << endl;
 
     // if only one of the two is compatible, then the LR is assigned,
     // otherwise is undefined
 
     DTEnums::DTCellSide lrcode;
     if (isCompatible.first && isCompatible.second) {
       usePairs ? lrcode=DTEnums::undefLR : lrcode=DTEnums::Left ; // if not usePairs then only use single side 
       tried.push_back(3);
       nCompatibleHits++;
     }
     else if (isCompatible.first) {
       lrcode=DTEnums::Left;
       tried.push_back(2);
       nCompatibleHits++;
     }
     else if (isCompatible.second) {
       lrcode=DTEnums::Right;
       tried.push_back(1);
       nCompatibleHits++;
     }
     else {
       tried.push_back(0);
       continue; // neither is compatible
     }
     result.push_back(DTSegmentCand::AssPoint(*hit, lrcode));
   }
 
 
   // check if too few associated hits or pattern already tried
   // TODO: two different if for nCompatibleHits < 3 =>printout and find ->
   // printour
   if (nCompatibleHits < 3) {
       if (debug) {
           cout << "Less than 3 hits " ;
           copy(tried.begin(), tried.end(), ostream_iterator<int>(std::cout));
           cout << endl;
       }
       // No need to insert into the list of patterns, as you will never search for it.
       //theTriedPattern.insert(tried);
   } else {
       // try to insert, return bool if already there
       bool isnew = theTriedPattern.insert(tried).second;
       if (isnew) {
           if (debug) {
               cout << "NOT Already tried " ;
               copy(tried.begin(), tried.end(), ostream_iterator<int>(std::cout));
               cout << endl;
           }
       } else {
           if (debug) {
               cout << "Already tried " ;
               copy(tried.begin(), tried.end(), ostream_iterator<int>(std::cout));
               cout << endl;
 
           }
           // empty the result vector
           result.clear();
       }
   }
   return result;
 }

vector< std::shared_ptr< DTHitPairForFit > > DTCombinatorialPatternReco::initHits	(	const DTSuperLayer *	sl,
		const std::vector< DTRecHit1DPair > &	hits
	)

private

Definition at line 89 of file DTCombinatorialPatternReco.cc.

References mps_fire::result, and theDTGeometry.

Referenced by DTCombinatorialPatternReco4D::buildPhiSuperSegmentsCandidates(), and reconstruct().

                                                                            {  
   
   vector<std::shared_ptr<DTHitPairForFit>> result;
   for (vector<DTRecHit1DPair>::const_iterator hit=hits.begin();
        hit!=hits.end(); ++hit) {
     result.push_back(std::make_shared<DTHitPairForFit>(*hit, *sl, theDTGeometry));
   }
   return result;
 }

edm::OwnVector< DTSLRecSegment2D > DTCombinatorialPatternReco::reconstruct	(	const DTSuperLayer *	sl,
		const std::vector< DTRecHit1DPair > &	hits
	)

virtual

this function is called in the producer

Implements DTRecSegment2DBaseAlgo.

Definition at line 55 of file DTCombinatorialPatternReco.cc.

References edm::OwnVector< T, P >::back(), buildSegments(), objects.IsoTrackAnalyzer::candidates, gather_cfg::cout, debug, initHits(), edm::OwnVector< T, P >::push_back(), mps_fire::result, theUpdator, and DTSegmentUpdator::update().

Referenced by DTCombinatorialPatternReco4D::reconstruct().

                                                                                {
 
   edm::OwnVector<DTSLRecSegment2D> result;
   vector<std::shared_ptr<DTHitPairForFit>> hitsForFit = initHits(sl, pairs);
 
   vector<DTSegmentCand*> candidates = buildSegments(sl, hitsForFit);
 
   vector<DTSegmentCand*>::const_iterator cand=candidates.begin();
   while (cand<candidates.end()) {
     
     DTSLRecSegment2D *segment = (**cand);
 
     theUpdator->update(segment,0);
 
     result.push_back(segment);
 
     if (debug) {
       cout<<"Reconstructed 2D segments "<< result.back() <<endl;
     }
 
     delete *(cand++); // delete the candidate!
   }
 
   return result;
 }

void DTCombinatorialPatternReco::setES ( const edm::EventSetup & setup )

virtual

Through this function the EventSetup is percolated to the objs which request it

Implements DTRecSegment2DBaseAlgo.

Definition at line 82 of file DTCombinatorialPatternReco.cc.

References edm::EventSetup::get(), DTSegmentUpdator::setES(), theDTGeometry, and theUpdator.

Referenced by DTCombinatorialPatternReco4D::setES().

                                                                 {
   // Get the DT Geometry
   setup.get<MuonGeometryRecord>().get(theDTGeometry);
   theUpdator->setES(setup);
 }