#include <DTCombinatorialExtendedPatternReco.h>

Inheritance diagram for DTCombinatorialExtendedPatternReco:

Public Member Functions
virtual std::string	algoName () const
	return the algo name
	DTCombinatorialExtendedPatternReco (const edm::ParameterSet &pset)
	Constructor.
virtual edm::OwnVector < DTSLRecSegment2D >	reconstruct (const DTSuperLayer *sl, const std::vector< DTRecHit1DPair > &hits)
	this function is called in the producer
void	setClusters (std::vector< DTSLRecCluster > clusters)
virtual void	setES (const edm::EventSetup &setup)
virtual	~DTCombinatorialExtendedPatternReco ()
	Destructor.
Private Types
typedef std::pair < DTHitPairForFit *, DTEnums::DTCellSide >	AssPoint
Private Member Functions
DTSegmentExtendedCand *	buildBestSegment (std::vector< AssPoint > &assHits, const DTSuperLayer *sl)
void	buildPointsCollection (std::vector< AssPoint > &points, std::deque< DTHitPairForFit * > &pointsNoLR, std::vector< DTSegmentCand * > &candidates, const DTSuperLayer *sl)
std::vector< DTSegmentCand * >	buildSegments (const DTSuperLayer sl, const std::vector< DTHitPairForFit > &hits)
bool	checkDoubleCandidates (std::vector< DTSegmentCand * > &segs, DTSegmentCand *seg)
bool	closeSL (const DTSuperLayerId &id1, const DTSuperLayerId &id2)
std::vector < DTSegmentExtendedCand * >	extendCandidates (std::vector< DTSegmentCand * > &candidates, const DTSuperLayer *sl)
std::vector< AssPoint >	findCompatibleHits (const LocalPoint &pos, const LocalVector &dir, const std::vector< DTHitPairForFit * > &hits)
std::vector< DTHitPairForFit * >	initHits (const DTSuperLayer *sl, const std::vector< DTRecHit1DPair > &hits)
Private Attributes
bool	debug
std::string	theAlgoName
double	theAlphaMaxPhi
double	theAlphaMaxTheta
DTSegmentCleaner *	theCleaner
std::vector< DTSLRecCluster >	theClusters
edm::ESHandle< DTGeometry >	theDTGeometry
unsigned int	theMaxAllowedHits
std::vector< std::vector< int > >	theTriedPattern
DTSegmentUpdator *	theUpdator
bool	usePairs

Detailed Description

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

Date:: 2008/12/03 12:52:22

Revision:: 1.1

Author:: Stefano Lacaprara - INFN Legnaro <stefano.lacaprara@pd.infn.it>; Riccardo Bellan - INFN TO <riccardo.bellan@cern.ch>

Definition at line 43 of file DTCombinatorialExtendedPatternReco.h.

Member Typedef Documentation

typedef std::pair<DTHitPairForFit*, DTEnums::DTCellSide> DTCombinatorialExtendedPatternReco::AssPoint [private]

Definition at line 73 of file DTCombinatorialExtendedPatternReco.h.

Constructor & Destructor Documentation

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

Constructor.

Definition at line 37 of file DTCombinatorialExtendedPatternReco.cc.

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

                                                                                                  : 
DTRecSegment2DBaseAlgo(pset), theAlgoName("DTCombinatorialExtendedPatternReco")
{
  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");
}

DTCombinatorialExtendedPatternReco::~DTCombinatorialExtendedPatternReco ( ) [virtual]

Destructor.

Definition at line 51 of file DTCombinatorialExtendedPatternReco.cc.

{
}

Member Function Documentation

virtual std::string DTCombinatorialExtendedPatternReco::algoName ( void ) const [inline, virtual]

return the algo name

Implements DTRecSegment2DBaseAlgo.

Definition at line 61 of file DTCombinatorialExtendedPatternReco.h.

References theAlgoName.

{ return theAlgoName; }

DTSegmentExtendedCand * DTCombinatorialExtendedPatternReco::buildBestSegment	(	std::vector< AssPoint > &	assHits,
		const DTSuperLayer *	sl
	)		`[private]`

Definition at line 304 of file DTCombinatorialExtendedPatternReco.cc.

References buildPointsCollection(), gather_cfg::cout, debug, extendCandidates(), and DTEnums::undefLR.

Referenced by buildSegments().

                                                                             {
  if (debug) cout << "DTCombinatorialExtendedPatternReco::buildBestSegment " <<
    hits.size()  << endl;
  if (hits.size()<3) {
    //cout << "buildBestSegment: hits " << hits.size()<< endl;
    return 0; // a least 3 point
  }

  // hits with defined LR
  vector<AssPoint> points;

  // without: I store both L and R, a deque since I need front insertion and
  // deletion
  deque<DTHitPairForFit* > pointsNoLR; 

  // first add only the hits with LR assigned
  for (vector<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);

  // here I try to add the external clusters and build a set of "extended
  // segment candidate
  vector<DTSegmentExtendedCand*> extendedCands = extendCandidates(candidates,
                                                                  sl); 
  if (debug) cout << "extended candidates " << extendedCands.size() << endl;

  // so now I have build a given number of segments, I should find the best one,
  // by #hits and chi2.
  vector<DTSegmentExtendedCand*>::const_iterator bestCandIter = extendedCands.end();
  double minChi2=999999.;
  unsigned int maxNumHits=0;
  for (vector<DTSegmentExtendedCand*>::const_iterator iter=extendedCands.begin();
       iter!=extendedCands.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<DTSegmentExtendedCand*>::iterator iter=extendedCands.begin();
       iter!=extendedCands.end(); ++iter) 
    if (iter!=bestCandIter) delete (*iter);

  // return the best candate if any
  if (bestCandIter != extendedCands.end()) {
    return (*bestCandIter);
  }
  return 0;
}

void DTCombinatorialExtendedPatternReco::buildPointsCollection	(	std::vector< AssPoint > &	points,
		std::deque< 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

Referenced by buildBestSegment(), and DTCombinatorialPatternReco::buildPointsCollection().

vector< DTSegmentCand * > DTCombinatorialExtendedPatternReco::buildSegments	(	const DTSuperLayer *	sl,
		const std::vector< DTHitPairForFit * > &	hits
	)		`[private]`

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

Definition at line 111 of file DTCombinatorialExtendedPatternReco.cc.

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

Referenced by reconstruct().

                                                                                  {

  typedef vector<DTHitPairForFit*> hitCont;
  typedef hitCont::const_iterator  hitIter;
  vector<DTSegmentCand*> result;
  
  if(debug) {
    cout << "DTCombinatorialExtendedPatternReco::buildSegments: " << sl->id() << " nHits " << hits.size() << endl;
    for (vector<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) {
      // 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<AssPoint> assHits = findCompatibleHits(posIni, dirIni, hits);
            if(debug) 
              cout << "compatible hits " << assHits.size() << endl;

            // here return an extended candidate (which _has_ the original
            // segment)
            DTSegmentExtendedCand* seg = buildBestSegment(assHits, sl);

            if (seg) {
              if(debug) 
                cout << "segment " << *seg<< endl;

              // check if the chi2 and #hits are ok
              if (!seg->good()) { // good is reimplmented in extended segment
                delete seg;
              } else { 

                // remove duplicated segments 
                if (checkDoubleCandidates(result,seg)) {
                  // add to the vector of hypotesis
                  // still work with extended segments
                  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
  // still with extended candidates
  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;
  }

  // here, finally, I have to return the set of _original_ segments, not the
  // extended ones.
  return result;
}

bool DTCombinatorialExtendedPatternReco::checkDoubleCandidates	(	std::vector< DTSegmentCand * > &	segs,
		DTSegmentCand *	seg
	)		`[private]`

Referenced by buildSegments().

bool DTCombinatorialExtendedPatternReco::closeSL	(	const DTSuperLayerId &	id1,
		const DTSuperLayerId &	id2
	)		`[private]`

Definition at line 513 of file DTCombinatorialExtendedPatternReco.cc.

References abs, DTChamberId::sector(), and DTChamberId::wheel().

Referenced by extendCandidates().

                                                                            {
  if (id1==id2) return false;
  if (abs(id1.wheel()-id2.wheel())>1 ) return false;
  // take into account also sector 13 and 14
  int sec1 = ( id1.sector()==13 ) ? 4: id1.sector();
  sec1=(sec1==14)? 10: sec1;
  int sec2 = ( id2.sector()==13 ) ? 4: id2.sector();
  sec2=(sec2==14)? 10: sec2;
  // take into account also sector 1/12
  if (abs(sec1-sec2)>1 && abs(sec1-sec2)!=11 ) return false;
  //if (abs(id1.station()-id2.station())>1 ) return false;
  return true;
}

vector< DTSegmentExtendedCand * > DTCombinatorialExtendedPatternReco::extendCandidates	(	std::vector< DTSegmentCand * > &	candidates,
		const DTSuperLayer *	sl
	)		`[private]`

extend the candidates with clusters from external SL

Definition at line 442 of file DTCombinatorialExtendedPatternReco.cc.

References DTSegmentExtendedCand::addClus(), closeSL(), gather_cfg::cout, debug, DTSegmentUpdator::fit(), DTSuperLayer::id(), DTSegmentExtendedCand::isCompatible(), DTSegmentExtendedCand::nHits(), pos, query::result, DTSuperLayerId::superLayer(), DTSegmentCand::superLayer(), theClusters, theDTGeometry, theUpdator, GeomDet::toGlobal(), and GeomDet::toLocal().

Referenced by buildBestSegment().

                                                                             {
  if (debug) cout << "extendCandidates " << candidates.size() << endl;
  vector<DTSegmentExtendedCand*> result;

  // in case of phi SL just return
  if (sl->id().superLayer() != 2 ) {
    for (vector<DTSegmentCand*>:: const_iterator cand=candidates.begin();
       cand!=candidates.end(); ++cand) {
      DTSegmentExtendedCand* extendedCand = new DTSegmentExtendedCand(*cand);
      // and delete the original candidate
      delete *cand;
      result.push_back(extendedCand);
    }
    return result;
  }

  // first I have to select the cluster which are compatible with the actual
  // candidate, namely +/-1 sector/station/wheel 
  vector<DTSegmentExtendedCand::DTSLRecClusterForFit> clustersWithPos;
  if (debug) cout << "AllClustersWithPos " << theClusters.size() << endl;
  if(debug) cout << "SL:   " << sl->id() << endl;
  for (vector<DTSLRecCluster>::const_iterator clus=theClusters.begin();
       clus!=theClusters.end(); ++clus) {
    if(debug) cout << "CLUS: " << (*clus).superLayerId() << endl;
    if ((*clus).superLayerId().superLayer()==2 && closeSL(sl->id(),(*clus).superLayerId())) {
      // and then get their pos in the actual SL frame
      const DTSuperLayer* clusSl =
        theDTGeometry->superLayer((*clus).superLayerId());
      LocalPoint pos=sl->toLocal(clusSl->toGlobal((*clus).localPosition()));
      //LocalError err=sl->toLocal(clusSl->toGlobal((*clus).localPositionError()));
      LocalError err=(*clus).localPositionError();
      clustersWithPos.push_back(DTSegmentExtendedCand::DTSLRecClusterForFit(*clus, pos, err) );
    }
  }
  if (debug) cout << "closeClustersWithPos " << clustersWithPos.size() << endl;

  for (vector<DTSegmentCand*>:: const_iterator cand=candidates.begin();
       cand!=candidates.end(); ++cand) {
    // create an extended candidate
    DTSegmentExtendedCand* extendedCand = new DTSegmentExtendedCand(*cand);
    // and delete the original candidate
    delete *cand;
    // do this only for theta SL
    if (extendedCand->superLayer()->id().superLayer() == 2 ) {
      // first check compatibility between cand and clusForFit
      for (vector<DTSegmentExtendedCand::DTSLRecClusterForFit>::const_iterator
           exClus=clustersWithPos.begin(); exClus!=clustersWithPos.end(); ++exClus) {
        if (extendedCand->isCompatible(*exClus)) {
          if (debug) cout << "is compatible " << endl;
          // add compatible cluster
          extendedCand->addClus(*exClus);
        }
      }
      // fit the segment
      if (debug) cout << "extended cands nHits: " << extendedCand->nHits() <<endl;
      if (theUpdator->fit(extendedCand)) {
        // add to result
        result.push_back(extendedCand);
      } else {
        cout << "Bad fit" << endl;
        delete extendedCand;
      }
    } else { // Phi SuperLayer
      result.push_back(extendedCand);
    }
  }

  return result;
}

std::vector<AssPoint> DTCombinatorialExtendedPatternReco::findCompatibleHits	(	const LocalPoint &	pos,
		const LocalVector &	dir,
		const std::vector< DTHitPairForFit * > &	hits
	)		`[private]`

Referenced by buildSegments().

vector< DTHitPairForFit * > DTCombinatorialExtendedPatternReco::initHits	(	const DTSuperLayer *	sl,
		const std::vector< DTRecHit1DPair > &	hits
	)		`[private]`

Definition at line 99 of file DTCombinatorialExtendedPatternReco.cc.

References query::result, and theDTGeometry.

Referenced by reconstruct().

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

edm::OwnVector< DTSLRecSegment2D > DTCombinatorialExtendedPatternReco::reconstruct	(	const DTSuperLayer *	sl,
		const std::vector< DTRecHit1DPair > &	hits
	)		`[virtual]`

this function is called in the producer

Implements DTRecSegment2DBaseAlgo.

Definition at line 56 of file DTCombinatorialExtendedPatternReco.cc.

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

                                                                                       {

  if(debug) cout << "DTCombinatorialExtendedPatternReco::reconstruct" << endl;
  theTriedPattern.clear();
  edm::OwnVector<DTSLRecSegment2D> result;
  vector<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);

    result.push_back(segment);

    if (debug) {
      cout<<"Reconstructed 2D extended segments "<< result.back() <<endl;
    }

    delete *(cand++); // delete the candidate!
  }

  for (vector<DTHitPairForFit*>::iterator it = hitsForFit.begin(), ed = hitsForFit.end(); 
        it != ed; ++it) delete *it;

  return result;
}

void DTCombinatorialExtendedPatternReco::setClusters ( std::vector< DTSLRecCluster > clusters )

Definition at line 94 of file DTCombinatorialExtendedPatternReco.cc.

References theClusters.

                                                                                    {
  theClusters = clusters;
}

void DTCombinatorialExtendedPatternReco::setES ( const edm::EventSetup & setup ) [virtual]

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

Implements DTRecSegment2DBaseAlgo.

Definition at line 88 of file DTCombinatorialExtendedPatternReco.cc.

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

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