#include <RecoTracker/NuclearSeedGenerator/plugin/NuclearTrackCorrector.cc>

Inheritance diagram for NuclearTrackCorrector:

Public Types
typedef TrackProducerAlgorithm < reco::Track > ::AlgoProductCollection	AlgoProductCollection
typedef TransientTrackingRecHit::ConstRecHitContainer	ConstRecHitContainer
typedef edm::Ref < TrackCandidateCollection >	TrackCandidateRef
typedef edm::Ref < TrajectoryCollection >	TrajectoryRef
typedef edm::RefVector < TrajectorySeedCollection >	TrajectorySeedRefVector
Public Member Functions
	NuclearTrackCorrector (const edm::ParameterSet &)
	~NuclearTrackCorrector ()
Private Member Functions
virtual void	beginRun (edm::Run &run, const edm::EventSetup &)
virtual void	endJob ()
TrajectoryStateOnSurface	getInitialState (const reco::Track theT, TransientTrackingRecHit::RecHitContainer &hits, const TrackingGeometry theG, const MagneticField *theMF)
	Calculate the inital state to be used to buil the track.
reco::TrackExtra	getNewTrackExtra (const AlgoProductCollection &algoresults)
	get a new TrackExtra from an AlgoProductCollection
bool	getTrackFromTrajectory (const Trajectory &newTraj, const TrajectoryRef &initialTrajRef, AlgoProductCollection &algoResults)
	Get the refitted track from the Trajectory.
bool	newTrajNeeded (Trajectory &newtrajectory, const TrajectoryRef &trajRef, const reco::NuclearInteraction &ni)
	check if the trajectory has to be refitted and get the new trajectory
virtual void	produce (edm::Event &, const edm::EventSetup &)
void	swap_map (const edm::Handle< TrajectoryCollection > &trajColl, std::map< reco::TrackRef, edm::Ref< TrajectoryCollection > > &result)
Private Attributes
edm::ParameterSet	conf_
std::vector< std::pair < unsigned int, unsigned int > >	Indice_Map
int	int_Input_Hit_Distance
int	KeepOnlyCorrectedTracks
const TrajTrackAssociationCollection *	m_TrajToTrackCollection
std::string	str_Input_NuclearInteraction
std::string	str_Input_Trajectory
TrackProducerAlgorithm < reco::Track > *	theAlgo
edm::ESHandle< TrajectoryFitter >	theFitter
edm::ESHandle< TrackerGeometry >	theG
TransientInitialStateEstimator *	theInitialState
edm::ESHandle< MagneticField >	theMF
edm::ESHandle< Propagator >	thePropagator
int	verbosity

Detailed Description

Description: <one line="" class="" summary>="">

Implementation: <Notes on="" implementation>="">

Definition at line 75 of file NuclearTrackCorrector.h.

Member Typedef Documentation

typedef TrackProducerAlgorithm<reco::Track>::AlgoProductCollection NuclearTrackCorrector::AlgoProductCollection

Definition at line 82 of file NuclearTrackCorrector.h.

typedef TransientTrackingRecHit::ConstRecHitContainer NuclearTrackCorrector::ConstRecHitContainer

Definition at line 81 of file NuclearTrackCorrector.h.

typedef edm::Ref<TrackCandidateCollection> NuclearTrackCorrector::TrackCandidateRef

Definition at line 80 of file NuclearTrackCorrector.h.

typedef edm::Ref<TrajectoryCollection> NuclearTrackCorrector::TrajectoryRef

Definition at line 79 of file NuclearTrackCorrector.h.

typedef edm::RefVector<TrajectorySeedCollection> NuclearTrackCorrector::TrajectorySeedRefVector

Definition at line 78 of file NuclearTrackCorrector.h.

Constructor & Destructor Documentation

NuclearTrackCorrector::NuclearTrackCorrector ( const edm::ParameterSet & iConfig ) [explicit]

Definition at line 36 of file NuclearTrackCorrector.cc.

References edm::ParameterSet::getParameter(), KeepOnlyCorrectedTracks, str_Input_NuclearInteraction, str_Input_Trajectory, theAlgo, and verbosity.

                                                                           :
conf_(iConfig),
theInitialState(0)
{
     str_Input_Trajectory           = iConfig.getParameter<std::string>     ("InputTrajectory");
     str_Input_NuclearInteraction   = iConfig.getParameter<std::string>     ("InputNuclearInteraction");
     verbosity                      = iConfig.getParameter<int>             ("Verbosity");
     KeepOnlyCorrectedTracks        = iConfig.getParameter<bool>             ("KeepOnlyCorrectedTracks");


     theAlgo = new TrackProducerAlgorithm<reco::Track>(iConfig);

     produces< TrajectoryCollection >();
     produces< TrajectoryToTrajectoryMap >();

     produces< reco::TrackExtraCollection >();
     produces< reco::TrackCollection >();       
     produces< TrackToTrajectoryMap >();

     produces< TrackToTrackMap >();
}

NuclearTrackCorrector::~NuclearTrackCorrector ( )

Definition at line 59 of file NuclearTrackCorrector.cc.

{
}

Member Function Documentation

void NuclearTrackCorrector::beginRun	(	edm::Run &	run,
		const edm::EventSetup &	iSetup
	)		`[private, virtual]`

Reimplemented from edm::EDProducer.

Definition at line 201 of file NuclearTrackCorrector.cc.

{
}

void NuclearTrackCorrector::endJob ( void ) [private, virtual]

Reimplemented from edm::EDProducer.

Definition at line 207 of file NuclearTrackCorrector.cc.

{
}

TrajectoryStateOnSurface NuclearTrackCorrector::getInitialState	(	const reco::Track *	theT,
		TransientTrackingRecHit::RecHitContainer &	hits,
		const TrackingGeometry *	theG,
		const MagneticField *	theMF
	)		`[private]`

Calculate the inital state to be used to buil the track.

Definition at line 318 of file NuclearTrackCorrector.cc.

References TrajectoryStateOnSurface::globalPosition(), trajectoryStateTransform::innerStateOnSurface(), TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localParameters(), mag2(), trajectoryStateTransform::outerStateOnSurface(), TrajectoryStateOnSurface::rescaleError(), TrajectoryStateOnSurface::surface(), and theMF.

Referenced by getTrackFromTrajectory().

                                                                                             {

  TrajectoryStateOnSurface theInitialStateForRefitting;
  //the starting state is the state closest to the first hit along seedDirection.
  
  //avoiding to use transientTrack, it should be faster;
  TrajectoryStateOnSurface innerStateFromTrack=trajectoryStateTransform::innerStateOnSurface(*theT,*theG,theMF);
  TrajectoryStateOnSurface outerStateFromTrack=trajectoryStateTransform::outerStateOnSurface(*theT,*theG,theMF);
  TrajectoryStateOnSurface initialStateFromTrack = 
    ( (innerStateFromTrack.globalPosition()-hits.front()->globalPosition()).mag2() <
      (outerStateFromTrack.globalPosition()-hits.front()->globalPosition()).mag2() ) ? 
    innerStateFromTrack: outerStateFromTrack;       
  
  // error is rescaled, but correlation are kept.
  initialStateFromTrack.rescaleError(100);
  theInitialStateForRefitting = TrajectoryStateOnSurface(initialStateFromTrack.localParameters(),
                                                         initialStateFromTrack.localError(),                  
                                                         initialStateFromTrack.surface(),
                                                         theMF); 
  return theInitialStateForRefitting;
}

reco::TrackExtra NuclearTrackCorrector::getNewTrackExtra ( const AlgoProductCollection & algoresults ) [private]

get a new TrackExtra from an AlgoProductCollection

Definition at line 284 of file NuclearTrackCorrector.cc.

References alongMomentum, TrajectoryStateOnSurface::curvilinearError(), Trajectory::direction(), Trajectory::firstMeasurement(), TrajectoryStateOnSurface::globalParameters(), Trajectory::lastMeasurement(), GlobalTrajectoryParameters::momentum(), AlCaHLTBitMon_ParallelJobs::p, GlobalTrajectoryParameters::position(), TrajectoryMeasurement::recHit(), TrajectoryMeasurement::updatedState(), v, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by produce().

                                                                                               {
                Trajectory* theTraj          = algoResults[0].first;
                PropagationDirection seedDir = algoResults[0].second.second;

                TrajectoryStateOnSurface outertsos;
                TrajectoryStateOnSurface innertsos;
                unsigned int innerId, outerId;
                if (theTraj->direction() == alongMomentum) {
                  outertsos = theTraj->lastMeasurement().updatedState();
                  innertsos = theTraj->firstMeasurement().updatedState();
                  outerId   = theTraj->lastMeasurement().recHit()->geographicalId().rawId();
                  innerId   = theTraj->firstMeasurement().recHit()->geographicalId().rawId();
                } else {
                  outertsos = theTraj->firstMeasurement().updatedState();
                  innertsos = theTraj->lastMeasurement().updatedState();
                  outerId   = theTraj->firstMeasurement().recHit()->geographicalId().rawId();
                  innerId   = theTraj->lastMeasurement().recHit()->geographicalId().rawId();
                }

                GlobalPoint v = outertsos.globalParameters().position();
                GlobalVector p = outertsos.globalParameters().momentum();
                math::XYZVector outmom( p.x(), p.y(), p.z() );
                math::XYZPoint  outpos( v.x(), v.y(), v.z() );
                v = innertsos.globalParameters().position();
                p = innertsos.globalParameters().momentum();
                math::XYZVector inmom( p.x(), p.y(), p.z() );
                math::XYZPoint  inpos( v.x(), v.y(), v.z() );

                return reco::TrackExtra (outpos, outmom, true, inpos, inmom, true,
                                        outertsos.curvilinearError(), outerId,
                                        innertsos.curvilinearError(), innerId, seedDir);

}

bool NuclearTrackCorrector::getTrackFromTrajectory	(	const Trajectory &	newTraj,
		const TrajectoryRef &	initialTrajRef,
		AlgoProductCollection &	algoResults
	)		`[private]`

Get the refitted track from the Trajectory.

Definition at line 244 of file NuclearTrackCorrector.cc.

References TrackProducerAlgorithm< T >::buildTrack(), getInitialState(), h, LogDebug, m_TrajToTrackCollection, edm::ESHandle< T >::product(), Trajectory::seed(), theAlgo, theFitter, theG, theMF, thePropagator, and Trajectory::validRecHits().

Referenced by produce().

                                                                                                                                                      {

        const Trajectory*  it = &newTraj;

        TransientTrackingRecHit::RecHitContainer hits;
        it->validRecHits( hits  );
        

        float ndof=0;
        for(unsigned int h=0 ; h<hits.size() ; h++)
        {
            if( hits[h]->isValid() )
            {
                ndof = ndof + hits[h]->dimension() * hits[h]->weight();
            }
            else {
                 LogDebug("NuclearSeedGenerator") << " HIT IS INVALID ???";
            }
        }


        ndof = ndof - 5;
        reco::TrackRef  theT     = m_TrajToTrackCollection->operator[]( initialTrajRef );
        LogDebug("NuclearSeedGenerator") << " TrackCorrector - number of valid hits" << hits.size() << "\n"
                                         << "                - number of hits from Track " << theT->recHitsSize() << "\n"
                                         << "                - number of valid hits from initial track " << theT->numberOfValidHits();


        if(  hits.size() > 1){

                TrajectoryStateOnSurface theInitialStateForRefitting = getInitialState(&(*theT),hits,theG.product(),theMF.product()
);

           reco::BeamSpot bs;
           return theAlgo->buildTrack(theFitter.product(), thePropagator.product(), algoResults, hits, theInitialStateForRefitting ,it->seed(), ndof, bs, theT->seedRef());
         }

        return false;
}

bool NuclearTrackCorrector::newTrajNeeded	(	Trajectory &	newtrajectory,
		const TrajectoryRef &	trajRef,
		const reco::NuclearInteraction &	ni
	)		`[private]`

check if the trajectory has to be refitted and get the new trajectory

Definition at line 210 of file NuclearTrackCorrector.cc.

References alongMomentum, LogDebug, m, PV3DBase< T, PVType, FrameType >::mag(), reco::Vertex::position(), Trajectory::push(), mathSSE::sqrt(), theG, verbosity, and reco::NuclearInteraction::vertex().

Referenced by produce().

                                                                                                                               {

        bool needNewTraj=false;
        reco::Vertex::Point vtx_pos = ni.vertex().position();
        double vtx_pos_mag = sqrt (vtx_pos.X()*vtx_pos.X()+vtx_pos.Y()*vtx_pos.Y()+vtx_pos.Z()*vtx_pos.Z());
        if(verbosity>=2) printf("Nuclear Interaction pos = %f\n",vtx_pos_mag );


        newtrajectory = Trajectory(trajRef->seed(), alongMomentum);

        // Look all the Hits of the trajectory and keep only Hits before seeds
        Trajectory::DataContainer Measurements = trajRef->measurements();
        if(verbosity>=2)LogDebug("NuclearTrackCorrector")<<"Size of Measurements  = "<<Measurements.size();

        for(unsigned int m=Measurements.size()-1 ;m!=(unsigned int)-1 ; m--){

                if(!Measurements[m].recHit()->isValid() )continue;
                GlobalPoint hit_pos = theG->idToDet(Measurements[m].recHit()->geographicalId())->surface().toGlobal(Measurements[m].recHit()->localPosition());

                if(verbosity>=2)printf("Hit pos = %f",hit_pos.mag() );

                if(hit_pos.mag()>vtx_pos_mag){
                         if(verbosity>=2)printf(" X ");
                         needNewTraj=true;
                }else{
                        newtrajectory.push(Measurements[m]);
                }
                if(verbosity>=2)printf("\n");
        }

        return needNewTraj;
}

void NuclearTrackCorrector::produce	(	edm::Event &	iEvent,
		const edm::EventSetup &	iSetup
	)		`[private, virtual]`

Implements edm::EDProducer.

Definition at line 64 of file NuclearTrackCorrector.cc.

References edm::RefToBase< T >::castTo(), conf_, event(), first, edm::EventSetup::get(), edm::Event::getByLabel(), getNewTrackExtra(), edm::ParameterSet::getParameter(), edm::Event::getRefBeforePut(), getTrackFromTrajectory(), i, Indice_Map, edm::AssociationMap< Tag >::insert(), KeepOnlyCorrectedTracks, reco::NuclearInteraction::likelihood(), LogDebug, m_TrajToTrackCollection, newTrajNeeded(), reco::NuclearInteraction::primaryTrack(), edm::Handle< T >::product(), edm::Event::put(), edm::second(), str_Input_NuclearInteraction, str_Input_Trajectory, swap_map(), theFitter, theG, theMF, thePropagator, and verbosity.

{
  // Create Output Collections
  // --------------------------------------------------------------------------------------------------
  std::auto_ptr<TrajectoryCollection>           Output_traj          ( new TrajectoryCollection );
  std::auto_ptr<TrajectoryToTrajectoryMap>      Output_trajmap       ( new TrajectoryToTrajectoryMap );

  std::auto_ptr<reco::TrackExtraCollection>     Output_trackextra    ( new reco::TrackExtraCollection );
  std::auto_ptr<reco::TrackCollection>          Output_track         ( new reco::TrackCollection );
  std::auto_ptr<TrackToTrajectoryMap>           Output_trackmap      ( new TrackToTrajectoryMap );

  std::auto_ptr<TrackToTrackMap>                Output_tracktrackmap ( new TrackToTrackMap );





  // Load Reccord
  // --------------------------------------------------------------------------------------------------
  std::string fitterName = conf_.getParameter<std::string>("Fitter");   
  iSetup.get<TrajectoryFitter::Record>().get(fitterName,theFitter);

  std::string propagatorName = conf_.getParameter<std::string>("Propagator");   
  iSetup.get<TrackingComponentsRecord>().get(propagatorName,thePropagator);

  iSetup.get<TrackerDigiGeometryRecord>().get(theG);

  reco::TrackExtraRefProd rTrackExtras = iEvent.getRefBeforePut<reco::TrackExtraCollection>();

   iSetup.get<IdealMagneticFieldRecord>().get(theMF); 

  // Load Inputs
  // --------------------------------------------------------------------------------------------------
  edm::Handle< TrajectoryCollection > temp_m_TrajectoryCollection;
  iEvent.getByLabel( str_Input_Trajectory.c_str(), temp_m_TrajectoryCollection );
  const TrajectoryCollection m_TrajectoryCollection = *(temp_m_TrajectoryCollection.product());

  edm::Handle< NuclearInteractionCollection > temp_m_NuclearInteractionCollection;
  iEvent.getByLabel( str_Input_NuclearInteraction.c_str(), temp_m_NuclearInteractionCollection );
  const NuclearInteractionCollection m_NuclearInteractionCollection = *(temp_m_NuclearInteractionCollection.product());

  edm::Handle< TrajTrackAssociationCollection > h_TrajToTrackCollection;
  iEvent.getByLabel( str_Input_Trajectory.c_str(), h_TrajToTrackCollection );
  m_TrajToTrackCollection = h_TrajToTrackCollection.product();


  // Correct the trajectories (Remove trajectory's hits that are located after the nuclear interacion)
  // --------------------------------------------------------------------------------------------------
  if(verbosity>=1){
    LogDebug("NuclearTrackCorrector")
      <<"Number of trajectories                    = "<<m_TrajectoryCollection.size() <<std::endl
      <<"Number of nuclear interactions            = "<<m_NuclearInteractionCollection.size();
  }

  std::map<reco::TrackRef,TrajectoryRef> m_TrackToTrajMap; 
  swap_map(temp_m_TrajectoryCollection, m_TrackToTrajMap);

  for(unsigned int i = 0 ; i < m_NuclearInteractionCollection.size() ; i++)
  {
        reco::NuclearInteraction ni =  m_NuclearInteractionCollection[i];
        if( ni.likelihood()<0.4) continue;

        reco::TrackRef primTrackRef = ni.primaryTrack().castTo<reco::TrackRef>();

        TrajectoryRef  trajRef = m_TrackToTrajMap[primTrackRef];

        Trajectory newTraj;
        if( newTrajNeeded(newTraj, trajRef, ni) ) {

          AlgoProductCollection   algoResults; 
          bool isOK = getTrackFromTrajectory( newTraj , trajRef, algoResults);

          if( isOK ) {

                pair<unsigned int, unsigned int> tempory_pair;
                tempory_pair.first  = Output_track->size();
                tempory_pair.second = i;
                Indice_Map.push_back(tempory_pair);

                reco::TrackExtraRef teref= reco::TrackExtraRef ( rTrackExtras, i );
                reco::TrackExtra newTrackExtra = getNewTrackExtra(algoResults);
                (algoResults[0].second.first)->setExtra( teref ); 

                Output_track->push_back(*algoResults[0].second.first);        
                Output_trackextra->push_back( newTrackExtra );
                Output_traj->push_back(newTraj);

          }
        }
        else {
           if(!KeepOnlyCorrectedTracks) {
                Output_track->push_back(*primTrackRef);
                Output_trackextra->push_back( *primTrackRef->extra() );
                Output_traj->push_back(*trajRef);
           }
        }

  }
  const edm::OrphanHandle<TrajectoryCollection>     Handle_traj = iEvent.put(Output_traj);
  const edm::OrphanHandle<reco::TrackCollection> Handle_tracks = iEvent.put(Output_track);
  iEvent.put(Output_trackextra);

  // Make Maps between elements
  // --------------------------------------------------------------------------------------------------
  if(Handle_tracks->size() != Handle_traj->size() )
  {
     printf("ERROR Handle_tracks->size() != Handle_traj->size() \n");
     return;
  }



  for(unsigned int i = 0 ; i < Indice_Map.size() ; i++)
  {
        TrajectoryRef      InTrajRef    ( temp_m_TrajectoryCollection, Indice_Map[i].second );
        TrajectoryRef      OutTrajRef   ( Handle_traj, Indice_Map[i].first );
        reco::TrackRef     TrackRef     ( Handle_tracks, Indice_Map[i].first );

        Output_trajmap ->insert(OutTrajRef,InTrajRef);
        Output_trackmap->insert(TrackRef,InTrajRef);

        try{
                reco::TrackRef  PrimaryTrackRef     = m_TrajToTrackCollection->operator[]( InTrajRef );
                Output_tracktrackmap->insert(TrackRef,PrimaryTrackRef);
        }catch(edm::Exception event){}
        
  }
  iEvent.put(Output_trajmap);
  iEvent.put(Output_trackmap);
  iEvent.put(Output_tracktrackmap);


  if(verbosity>=3)printf("-----------------------\n");
}

void NuclearTrackCorrector::swap_map	(	const edm::Handle< TrajectoryCollection > &	trajColl,
		std::map< reco::TrackRef, edm::Ref< TrajectoryCollection > > &	result
	)		`[private]`

Definition at line 343 of file NuclearTrackCorrector.cc.

References i, m_TrajToTrackCollection, and query::result.

Referenced by produce().

                                                                                                                                                       {
  for(unsigned int i = 0 ; i < trajColl->size() ; i++)
  {
     TrajectoryRef      InTrajRef    ( trajColl, i);
     reco::TrackRef  PrimaryTrackRef     = m_TrajToTrackCollection->operator[]( InTrajRef );
     result[ PrimaryTrackRef ] = InTrajRef;
  }
}