#include <MuonErrorMatrixAdjuster.h>

Inheritance diagram for MuonErrorMatrixAdjuster:

Public Member Functions
	MuonErrorMatrixAdjuster (const edm::ParameterSet &)
	constructor
	~MuonErrorMatrixAdjuster ()
	destructor
Private Member Functions
bool	attachRecHits (const reco::Track &recotrack_orig, reco::Track &recotrack, reco::TrackExtra &trackextra, TrackingRecHitCollection &RHcol)
	attached rechits to the newly created reco::Track and reco::TrackExtra
virtual void	beginJob ()
	framework method
bool	divide (reco::TrackBase::CovarianceMatrix &num_matrix, const reco::TrackBase::CovarianceMatrix &denom_matrix)
	divide the num_matrix (first argument) by second matrix, TERM by TERM
virtual void	endJob ()
reco::TrackBase::CovarianceMatrix	fix_cov_matrix (const reco::TrackBase::CovarianceMatrix &error_matrix, const GlobalVector &momentum)
	return a corrected error matrix
reco::Track	makeTrack (const reco::Track &recotrack_orig, const FreeTrajectoryState &PCAstate)
	create a corrected reco::Track from itself and trajectory state (redundant information)
reco::TrackExtra *	makeTrackExtra (const reco::Track &recotrack_orig, reco::Track &recotrack, reco::TrackExtraCollection &TEcol)
	create a track extra for the newly created recotrack, scaling the outer/inner measurment error matrix by the scale matrix recotrack/recotrack_orig
void	multiply (reco::TrackBase::CovarianceMatrix &revised_matrix, const reco::TrackBase::CovarianceMatrix &scale_matrix)
	mutliply revised_matrix (first argument) by second matrix TERM by TERM
virtual void	produce (edm::Event &, const edm::EventSetup &)
bool	selectTrack (const reco::Track &recotrack_orig)
	make a selection on the reco:Track. (dummy for the moment)
Private Attributes
std::string	theCategory
	log category: MuonErrorMatrixAdjuster
edm::ESHandle< MagneticField >	theField
	hold on to the magnetic field
std::string	theInstanceName
	instrance name of the created track collecion. rechit and trackextra have no instance name
MuonErrorMatrix *	theMatrixProvider
edm::ParameterSet	theMatrixProvider_pset
	holds the error matrix parametrization
edm::RefProd < TrackingRecHitCollection >	theRefprodRH
	get reference before put rechit to the event, in order to create edm::Ref
edm::RefProd < reco::TrackExtraCollection >	theRefprodTE
	get reference before put track extra to the event, in order to create edm::Ref
bool	theRescale
	select the rescaling or replacing method to correct the error matrix
edm::Ref < TrackingRecHitCollection > ::key_type	theRHi
edm::Ref < reco::TrackExtraCollection > ::key_type	theTEi
edm::InputTag	theTrackLabel
	input tag of the reco::Track collection to be corrected

Detailed Description

EDProducer which duplicatesa collection of track, adjusting their error matrix

track collection is retrieve from the event, duplicated, while the error matrix is corrected rechit are copied into a new collection track extra is also copied and error matrix are corrected by the same scale factors

Dates:: 2007/09/04 13:28

Revision:: 1.3

Author:: Jean-Roch Vlimant UCSB; Finn Rebassoo UCSB

Definition at line 47 of file MuonErrorMatrixAdjuster.h.

Constructor & Destructor Documentation

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

constructor

Definition at line 16 of file MuonErrorMatrixAdjuster.cc.

References edm::ParameterSet::getParameter(), theCategory, theInstanceName, theMatrixProvider_pset, theRescale, and theTrackLabel.

{
  theCategory="MuonErrorMatrixAdjuster";
  theInstanceName = iConfig.getParameter<std::string>("instanceName");
  //register your products
  produces<reco::TrackCollection>( theInstanceName);
  produces<TrackingRecHitCollection>(); 
  produces<reco::TrackExtraCollection>();
          
  theTrackLabel = iConfig.getParameter<edm::InputTag>("trackLabel");
  theRescale = iConfig.getParameter<bool>("rescale");

  theMatrixProvider_pset = iConfig.getParameter<edm::ParameterSet>("errorMatrix_pset");
}

MuonErrorMatrixAdjuster::~MuonErrorMatrixAdjuster ( )

destructor

Definition at line 32 of file MuonErrorMatrixAdjuster.cc.

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

}

Member Function Documentation

bool MuonErrorMatrixAdjuster::attachRecHits	(	const reco::Track &	recotrack_orig,
		reco::Track &	recotrack,
		reco::TrackExtra &	trackextra,
		TrackingRecHitCollection &	RHcol
	)		`[private]`

attached rechits to the newly created reco::Track and reco::TrackExtra

Definition at line 136 of file MuonErrorMatrixAdjuster.cc.

References reco::TrackExtraBase::add(), edm::OwnVector< T, P >::push_back(), reco::Track::recHitsBegin(), reco::Track::recHitsEnd(), reco::TrackBase::setHitPattern(), theRefprodRH, and theRHi.

Referenced by produce().

                                                                       {
  //loop over the hits of the original track
  trackingRecHit_iterator recHit = recotrack_orig.recHitsBegin();
  unsigned int irh=0;
  for (; recHit!=recotrack_orig.recHitsEnd();++recHit){
    //clone it. this is meandatory
    TrackingRecHit * hit = (*recHit)->clone();
    
    //put it on the new track
    recotrack.setHitPattern(*hit,irh++);
    //copy them in the new collection
    RHcol.push_back(hit);
    //do something with the trackextra 
    trackextra.add(TrackingRecHitRef( theRefprodRH, theRHi++));

  }//loop over original rechits
  
  return true; //if nothing fails
}

void MuonErrorMatrixAdjuster::beginJob ( void ) [private, virtual]

framework method

Reimplemented from edm::EDProducer.

Definition at line 232 of file MuonErrorMatrixAdjuster.cc.

References theMatrixProvider, and theMatrixProvider_pset.

{
  theMatrixProvider = new MuonErrorMatrix(theMatrixProvider_pset);
}

bool MuonErrorMatrixAdjuster::divide	(	reco::TrackBase::CovarianceMatrix &	num_matrix,
		const reco::TrackBase::CovarianceMatrix &	denom_matrix
	)		`[private]`

divide the num_matrix (first argument) by second matrix, TERM by TERM

Definition at line 66 of file MuonErrorMatrixAdjuster.cc.

References i, and j.

Referenced by makeTrackExtra().

                                                                                                                                  {
  //divide term by term the matrix
  // the true type of the matrix is such that [i][j] is the same memory object as [j][i]: looping i:0-5, j:0-5 double multiply the terms
  // need to loop only on i:0-5, j:i-5
  for(int i = 0;i!=5;i++){for(int j = i;j!=5;j++){
      if (denom_matrix(i,j)==0) return false;
      num_matrix(i,j)/=denom_matrix(i,j);
   }}
  return true;
}

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

Reimplemented from edm::EDProducer.

Definition at line 239 of file MuonErrorMatrixAdjuster.cc.

References theMatrixProvider.

                                {
  if ( theMatrixProvider) delete theMatrixProvider;
}

reco::TrackBase::CovarianceMatrix MuonErrorMatrixAdjuster::fix_cov_matrix	(	const reco::TrackBase::CovarianceMatrix &	error_matrix,
		const GlobalVector &	momentum
	)		`[private]`

return a corrected error matrix

Definition at line 46 of file MuonErrorMatrixAdjuster.cc.

References MuonErrorMatrix::get(), multiply(), theMatrixProvider, and theRescale.

Referenced by makeTrack().

{   
  //CovarianceMatrix is template for SMatrix
  reco::TrackBase::CovarianceMatrix revised_matrix( theMatrixProvider->get(momentum));
   
  if( theRescale){
    //rescale old error matrix up by a factor
    multiply(revised_matrix,error_matrix); 
  }
  return revised_matrix;
}

reco::Track MuonErrorMatrixAdjuster::makeTrack	(	const reco::Track &	recotrack_orig,
		const FreeTrajectoryState &	PCAstate
	)		`[private]`

create a corrected reco::Track from itself and trajectory state (redundant information)

Definition at line 77 of file MuonErrorMatrixAdjuster.cc.

References reco::TrackBase::charge(), DeDxDiscriminatorTools::charge(), reco::TrackBase::chi2(), reco::TrackBase::covariance(), fix_cov_matrix(), LogDebug, reco::TrackBase::momentum(), FreeTrajectoryState::momentum(), reco::TrackBase::ndof(), reco::TrackBase::referencePoint(), and theCategory.

Referenced by produce().

                                                                               {
  //get the parameters of the track so I can reconstruct it
  double chi2 = recotrack_orig.chi2();
  double ndof = recotrack_orig.ndof();
  math::XYZPoint refpos = recotrack_orig.referencePoint();
  math::XYZVector mom = recotrack_orig.momentum();
  int charge = recotrack_orig.charge();
  
  reco::TrackBase::CovarianceMatrix covariance_matrix = fix_cov_matrix(recotrack_orig.covariance(),PCAstate.momentum());

  LogDebug(theCategory)<<"chi2: "<<chi2
                       <<"\n ndof: "<<ndof
                       <<"\n refpos: "<<refpos
                       <<"\n mom: "<<mom
                       <<"\n charge: "<<charge
                       <<"\n covariance:\n"<<recotrack_orig.covariance()
                       <<"\n replaced by:\n"<<covariance_matrix;

  return reco::Track(chi2,ndof,refpos,mom,charge,covariance_matrix);
}

reco::TrackExtra * MuonErrorMatrixAdjuster::makeTrackExtra	(	const reco::Track &	recotrack_orig,
		reco::Track &	recotrack,
		reco::TrackExtraCollection &	TEcol
	)		`[private]`

create a track extra for the newly created recotrack, scaling the outer/inner measurment error matrix by the scale matrix recotrack/recotrack_orig

Definition at line 99 of file MuonErrorMatrixAdjuster.cc.

References reco::TrackBase::covariance(), divide(), reco::Track::extra(), edm::Ref< C, T, F >::isNull(), multiply(), reco::Track::setExtra(), theCategory, theRefprodTE, and theTEi.

Referenced by produce().

                                                                                      { 
  //get the 5x5 matrix of recotrack/recotrack_orig
  reco::TrackBase::CovarianceMatrix scale_matrix(recotrack.covariance());
  if (!divide(scale_matrix,recotrack_orig.covariance())){
    edm::LogError( theCategory)<<"original track error matrix has term ==0... skipping.";
    return 0; }
  
  const reco::TrackExtraRef & trackExtra_orig = recotrack_orig.extra();
  if (trackExtra_orig.isNull()) {
    edm::LogError( theCategory)<<"original track has no track extra... skipping.";
    return 0;}

  //copy the outer state. rescaling the error matrix
  reco::TrackBase::CovarianceMatrix outerCov(trackExtra_orig->outerStateCovariance());
  multiply(outerCov,scale_matrix);

  //copy the inner state, rescaling the error matrix
  reco::TrackBase::CovarianceMatrix innerCov(trackExtra_orig->innerStateCovariance());
  multiply(innerCov,scale_matrix);

  //put the trackExtra
  TEcol.push_back(reco::TrackExtra (trackExtra_orig->outerPosition(), trackExtra_orig->outerMomentum(), true,
                                    trackExtra_orig->innerPosition(), trackExtra_orig->innerMomentum(), true,
                                    outerCov, trackExtra_orig->outerDetId(),
                                    innerCov, trackExtra_orig->innerDetId(),
                                    trackExtra_orig->seedDirection()));

  //add a reference to the trackextra on the track
  recotrack.setExtra(edm::Ref<reco::TrackExtraCollection>( theRefprodTE, theTEi++));
  
  //return the reference to the last inserted then
  return &(TEcol.back());
}

void MuonErrorMatrixAdjuster::multiply	(	reco::TrackBase::CovarianceMatrix &	revised_matrix,
		const reco::TrackBase::CovarianceMatrix &	scale_matrix
	)		`[private]`

mutliply revised_matrix (first argument) by second matrix TERM by TERM

Definition at line 58 of file MuonErrorMatrixAdjuster.cc.

References i, and j.

Referenced by fix_cov_matrix(), and makeTrackExtra().

                                                                                                                                        {
  //scale term by term the matrix
  // the true type of the matrix is such that [i][j] is the same memory object as [j][i]: looping i:0-5, j:0-5 double multiply the terms
  // need to loop only on i:0-5, j:i-5
  for(int i = 0;i!=5;i++){for(int j = i;j!=5;j++){
      revised_matrix(i,j)*=scale_matrix(i,j);
    }}
}

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

Implements edm::EDProducer.

Definition at line 163 of file MuonErrorMatrixAdjuster.cc.

References attachRecHits(), edm::EventSetup::get(), edm::Event::getByLabel(), edm::Event::getRefBeforePut(), trajectoryStateTransform::initialFreeState(), LogDebug, PV3DBase< T, PVType, FrameType >::mag(), makeTrack(), makeTrackExtra(), FreeTrajectoryState::position(), edm::ESHandle< T >::product(), edm::Event::put(), Reconstruction_cff::recotrack, selectTrack(), theCategory, theField, theInstanceName, theRefprodRH, theRefprodTE, theRHi, theTEi, theTrackLabel, and testEve_cfg::tracks.

{
  using namespace edm;

  //open a collection of track
  edm::Handle<reco::TrackCollection> tracks;
  iEvent.getByLabel( theTrackLabel,tracks);
  LogDebug( theCategory)<<"considering: "<<tracks->size()<<" uncorrected reco::Track from the event.("<< theTrackLabel<<")";
  
  //get the mag field
  iSetup.get<IdealMagneticFieldRecord>().get( theField);
  
  //prepare the output collection
  std::auto_ptr<reco::TrackCollection> Toutput(new reco::TrackCollection());
  std::auto_ptr<TrackingRecHitCollection> TRHoutput(new TrackingRecHitCollection());
  std::auto_ptr<reco::TrackExtraCollection> TEoutput(new reco::TrackExtraCollection());
  theRefprodTE = iEvent.getRefBeforePut<reco::TrackExtraCollection>();
  theTEi=0;
  theRefprodRH =iEvent.getRefBeforePut<TrackingRecHitCollection>();
  theRHi=0;
  
  
  
  for(unsigned int it=0;it!=tracks->size();it++)
    {
      const reco::Track & recotrack_orig = (*tracks)[it];
      FreeTrajectoryState PCAstate = trajectoryStateTransform::initialFreeState(recotrack_orig, theField.product());
      if (PCAstate.position().mag()==0)
        {edm::LogError( theCategory)<<"invalid state from track initial state in "<< theTrackLabel <<". skipping.";continue;}

      //create a reco::Track
      reco::Track recotrack = makeTrack(recotrack_orig,PCAstate);
      
      //make a selection on the create reco::Track
      if (!selectTrack(recotrack)) continue;
      
      Toutput->push_back(recotrack);
      reco::Track & recotrackref = Toutput->back();
      
      //build the track extra
      reco::TrackExtra * extra = makeTrackExtra(recotrack_orig,recotrackref,*TEoutput);
      if (!extra){
        edm::LogError( theCategory)<<"cannot create the track extra for this track.";
        //pop the inserted track
        Toutput->pop_back();
        continue;}
      
      //attach the collection of rechits
      if (!attachRecHits(recotrack_orig,recotrackref,*extra,*TRHoutput)){
        edm::LogError( theCategory)<<"cannot attach any rechits on this track";
        //pop the inserted track
        Toutput->pop_back();
        //pop the track extra
        TEoutput->pop_back();
        theTEi--;
        continue;}
      
    }//loop over the original reco tracks
  
  LogDebug( theCategory)<<"writing: "<<Toutput->size()<<" corrected reco::Track to the event.";
  
  iEvent.put(Toutput, theInstanceName);
  iEvent.put(TEoutput);
  iEvent.put(TRHoutput);

}