CMS 3D CMS Logo

Public Member Functions | Private Member Functions | Private Attributes

QualityFilter Class Reference

#include <QualityFilter.h>

Inheritance diagram for QualityFilter:
edm::EDProducer edm::ProducerBase edm::ProductRegistryHelper

List of all members.

Public Member Functions

 QualityFilter (const edm::ParameterSet &)
 ~QualityFilter ()

Private Member Functions

virtual void beginRun (edm::Run &run, const edm::EventSetup &)
virtual void endJob ()
virtual void produce (edm::Event &, const edm::EventSetup &)

Private Attributes

bool copyExtras_
edm::InputTag tkTag
reco::TrackBase::TrackQuality trackQuality_

Detailed Description

Definition at line 10 of file QualityFilter.h.


Constructor & Destructor Documentation

QualityFilter::QualityFilter ( const edm::ParameterSet iConfig) [explicit]

Definition at line 17 of file QualityFilter.cc.

References edm::ParameterSet::getParameter(), and edm::ParameterSet::getUntrackedParameter().

{
  copyExtras_ = iConfig.getUntrackedParameter<bool>("copyExtras", false);

  produces<reco::TrackCollection>();
  if (copyExtras_) {
      produces<TrackingRecHitCollection>();
      produces<reco::TrackExtraCollection>();
  }
  produces<std::vector<Trajectory> >();
  produces<TrajTrackAssociationCollection>();

  tkTag  = iConfig.getParameter<edm::InputTag>("recTracks");
  trackQuality_=TrackBase::qualityByName(iConfig.getParameter<std::string>("TrackQuality"));
}
QualityFilter::~QualityFilter ( )

Definition at line 34 of file QualityFilter.cc.

{
 
   // do anything here that needs to be done at desctruction time
   // (e.g. close files, deallocate resources etc.)

}

Member Function Documentation

void QualityFilter::beginRun ( edm::Run run,
const edm::EventSetup  
) [private, virtual]

Reimplemented from edm::EDProducer.

Definition at line 168 of file QualityFilter.cc.

{
}
void QualityFilter::endJob ( void  ) [private, virtual]

Reimplemented from edm::EDProducer.

Definition at line 174 of file QualityFilter.cc.

                      {
}
void QualityFilter::produce ( edm::Event iEvent,
const edm::EventSetup iSetup 
) [private, virtual]

Implements edm::EDProducer.

Definition at line 49 of file QualityFilter.cc.

References edm::Event::getByLabel(), edm::Event::getRefBeforePut(), getHLTprescales::index, reco::Track::innerDetId(), reco::Track::innerMomentum(), reco::Track::innerOk(), reco::Track::innerPosition(), reco::Track::innerStateCovariance(), reco::Track::outerDetId(), reco::Track::outerMomentum(), reco::Track::outerOk(), reco::Track::outerPosition(), reco::Track::outerStateCovariance(), edm::Handle< T >::product(), edm::Event::put(), reco::Track::recHitsSize(), reco::Track::seedDirection(), and reco::Track::seedRef().

{

  
  auto_ptr<TrackCollection> selTracks(new TrackCollection);
  auto_ptr<TrackingRecHitCollection> selHits(copyExtras_ ? new TrackingRecHitCollection() : 0);
  auto_ptr<TrackExtraCollection> selTrackExtras(copyExtras_ ? new TrackExtraCollection() : 0);
  auto_ptr<vector<Trajectory> > outputTJ(new vector<Trajectory> );
  auto_ptr<TrajTrackAssociationCollection> trajTrackMap( new TrajTrackAssociationCollection() );
  
  
  TrackExtraRefProd rTrackExtras; 
  TrackingRecHitRefProd rHits;
  if (copyExtras_) {
    rTrackExtras = iEvent.getRefBeforePut<TrackExtraCollection>();
    rHits = iEvent.getRefBeforePut<TrackingRecHitCollection>();
  }

  Handle<std::vector<Trajectory> > TrajectoryCollection;
  Handle<TrajTrackAssociationCollection> assoMap;
  
  Handle<TrackCollection> tkCollection;  
  iEvent.getByLabel( tkTag, tkCollection);
  const reco::TrackCollection*  tC = tkCollection.product();
  TrackCollection::const_iterator itxc;
  TrackCollection::const_iterator firstTrack = tC->begin();
  TrackCollection::const_iterator lastTrack = tC->end();

    
  iEvent.getByLabel(tkTag,TrajectoryCollection);
  iEvent.getByLabel(tkTag,assoMap);



  TrajTrackAssociationCollection::const_iterator it = assoMap->begin();
  TrajTrackAssociationCollection::const_iterator lastAssoc = assoMap->end();
  for( ; it != lastAssoc; ++it ) {
    const Ref<vector<Trajectory> > traj = it->key;
    const reco::TrackRef itc = it->val;
    bool goodTk = (itc->quality(trackQuality_));
 
    
    if (goodTk){
      Track track =(*itc);
      //tracks and trajectories
      selTracks->push_back( track );
      outputTJ->push_back( *traj );
      if (copyExtras_) {
          //TRACKING HITS
          trackingRecHit_iterator irhit   =(*itc).recHitsBegin();
          trackingRecHit_iterator lasthit =(*itc).recHitsEnd();
          for (; irhit!=lasthit; ++irhit) {
            selHits->push_back((*irhit)->clone() );
          }
      }
          
    }

  }

  unsigned nTracks = selTracks->size();
  if (copyExtras_) {
      //PUT TRACKING HITS IN THE EVENT
      OrphanHandle<TrackingRecHitCollection> theRecoHits = iEvent.put(selHits );
  
      //PUT TRACK EXTRA IN THE EVENT
      selTrackExtras->reserve(nTracks);
      unsigned hits=0;

      for ( unsigned index = 0; index<nTracks; ++index ) { 

        reco::Track& aTrack = selTracks->at(index);
        TrackExtra aTrackExtra(aTrack.outerPosition(),
                               aTrack.outerMomentum(),
                               aTrack.outerOk(),
                               aTrack.innerPosition(),
                               aTrack.innerMomentum(),
                               aTrack.innerOk(),
                               aTrack.outerStateCovariance(),
                               aTrack.outerDetId(),
                               aTrack.innerStateCovariance(),
                               aTrack.innerDetId(),
                               aTrack.seedDirection(),
                               aTrack.seedRef());
            

        //unsigned nHits = aTrack.numberOfValidHits();
        unsigned nHits = aTrack.recHitsSize();
        for ( unsigned int ih=0; ih<nHits; ++ih) {
          aTrackExtra.add(TrackingRecHitRef(theRecoHits,hits++));
        }
        selTrackExtras->push_back(aTrackExtra);
      }

      //CORRECT REF TO TRACK
      OrphanHandle<TrackExtraCollection> theRecoTrackExtras = iEvent.put(selTrackExtras); 
      for ( unsigned index = 0; index<nTracks; ++index ) { 
        const reco::TrackExtraRef theTrackExtraRef(theRecoTrackExtras,index);
        (selTracks->at(index)).setExtra(theTrackExtraRef);
      }
  } // END IF COPY EXTRAS

  //TRACKS AND TRAJECTORIES
  OrphanHandle<TrackCollection> theRecoTracks = iEvent.put(selTracks);
  OrphanHandle<vector<Trajectory> > theRecoTrajectories = iEvent.put(outputTJ);

  //TRACKS<->TRAJECTORIES MAP 
  nTracks = theRecoTracks->size();
  for ( unsigned index = 0; index<nTracks; ++index ) { 
    Ref<vector<Trajectory> > trajRef( theRecoTrajectories, index );
    Ref<TrackCollection>    tkRef( theRecoTracks, index );
    trajTrackMap->insert(trajRef,tkRef);
  }
  //MAP IN THE EVENT
  iEvent.put( trajTrackMap );
}

Member Data Documentation

Definition at line 25 of file QualityFilter.h.

Definition at line 23 of file QualityFilter.h.

Definition at line 24 of file QualityFilter.h.