CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/RecoVertex/KalmanVertexFit/plugins/KVFTest.cc

Go to the documentation of this file.
00001 #include "RecoVertex/KalmanVertexFit/plugins/KVFTest.h"
00002 
00003 #include "DataFormats/VertexReco/interface/Vertex.h"
00004 #include "DataFormats/VertexReco/interface/VertexFwd.h"
00005 #include "DataFormats/TrackReco/interface/Track.h"
00006 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00007 #include "DataFormats/Common/interface/Handle.h"
00008 #include "FWCore/Framework/interface/MakerMacros.h"
00009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00010 #include "FWCore/Framework/interface/ESHandle.h"
00011 
00012 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
00013 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
00014 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
00015 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
00016 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
00017 #include "SimTracker/Records/interface/TrackAssociatorRecord.h"
00018 
00019 #include <iostream>
00020 
00021 using namespace reco;
00022 using namespace edm;
00023 using namespace std;
00024 
00025 KVFTest::KVFTest(const edm::ParameterSet& iConfig)
00026   : theConfig(iConfig), associatorForParamAtPca(0), tree(0)
00027 {
00028   trackLabel_ = iConfig.getParameter<std::string>("TrackLabel");
00029   outputFile_ = iConfig.getUntrackedParameter<std::string>("outputFile");
00030   kvfPSet = iConfig.getParameter<edm::ParameterSet>("KVFParameters");
00031   rootFile_ = TFile::Open(outputFile_.c_str(),"RECREATE"); 
00032   edm::LogInfo("RecoVertex/KVFTest") 
00033     << "Initializing KVF TEST analyser  - Output file: " << outputFile_ <<"\n";
00034 }
00035 
00036 
00037 KVFTest::~KVFTest() {
00038   delete rootFile_;
00039 }
00040 
00041 void KVFTest::beginJob(){
00042 }
00043 
00044 
00045 void KVFTest::endJob() {
00046   delete tree;
00047 }
00048 
00049 //
00050 // member functions
00051 //
00052 
00053 void
00054 KVFTest::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00055 {
00056   if ( associatorForParamAtPca==0 ) {
00057     edm::ESHandle<TrackAssociatorBase> theAssociatorForParamAtPca;
00058     iSetup.get<TrackAssociatorRecord>().get("TrackAssociatorByChi2",theAssociatorForParamAtPca);
00059     associatorForParamAtPca = (TrackAssociatorByChi2 *) theAssociatorForParamAtPca.product();
00060 
00061     tree = new SimpleVertexTree("VertexFitter", associatorForParamAtPca);
00062   }
00063 
00064 
00065 
00066   edm::LogInfo("RecoVertex/KVFTest") 
00067     << "Reconstructing event number: " << iEvent.id() << "\n";
00068     
00069   // get RECO tracks from the event
00070   // `tks` can be used as a ptr to a reco::TrackCollection
00071   edm::Handle<edm::View<reco::Track> > tks;
00072   iEvent.getByLabel(trackLabel_, tks);
00073   if (!tks.isValid()) {
00074     edm::LogInfo("RecoVertex/KVFTest") 
00075       << "Exception during event number: " << iEvent.id()
00076       << "\n";
00077   } else {
00078     edm::LogInfo("RecoVertex/KVFTest") 
00079       << "Found: " << (*tks).size() << " reconstructed tracks" << "\n";
00080     std::cout << "got " << (*tks).size() << " tracks " << std::endl;
00081     
00082     // Transform Track to TransientTrack
00083 
00084     //get the builder:
00085     edm::ESHandle<TransientTrackBuilder> theB;
00086     iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder",theB);
00087     //do the conversion:
00088     std::vector<TransientTrack> t_tks = (*theB).build(tks);
00089 
00090     edm::LogInfo("RecoVertex/KVFTest") 
00091       << "Found: " << t_tks.size() << " reconstructed tracks" << "\n";
00092     
00093     // Call the KalmanVertexFitter if more than 1 track
00094     if (t_tks.size() > 1) {
00095       //      KalmanVertexFitter kvf(kvfPSet);
00096       KalmanVertexFitter kvf(true);
00097       TransientVertex tv = kvf.vertex(t_tks);
00098 
00099       std::cout << "Position: " << Vertex::Point(tv.position()) << "\n";
00100 
00101       // For the analysis: compare to your SimVertex
00102       TrackingVertex sv = getSimVertex(iEvent);
00103       edm::Handle<TrackingParticleCollection>  TPCollectionH ;
00104       iEvent.getByLabel("trackingtruth","TrackTruth",TPCollectionH);
00105       if (!TPCollectionH.isValid()) {
00106         edm::LogInfo("RecoVertex/KVFTest") 
00107           << "Exception during event number: " << iEvent.id() 
00108           << "\n";
00109       } else {
00110         const TrackingParticleCollection tPC = *(TPCollectionH.product());
00111         reco::RecoToSimCollection recSimColl=associatorForParamAtPca->associateRecoToSim(tks,
00112                                                                                          TPCollectionH,
00113                                                                                          &iEvent);    
00114         tree->fill(tv, &sv, &recSimColl);
00115       }
00116     }
00117   }  
00118 }
00119 
00120 //Returns the first vertex in the list.
00121 
00122 TrackingVertex KVFTest::getSimVertex(const edm::Event& iEvent) const
00123 {
00124    // get the simulated vertices
00125   edm::Handle<TrackingVertexCollection>  TVCollectionH ;
00126   iEvent.getByLabel("trackingtruth","VertexTruth",TVCollectionH);
00127   const TrackingVertexCollection tPC = *(TVCollectionH.product());
00128 
00129 //    Handle<edm::SimVertexContainer> simVtcs;
00130 //    iEvent.getByLabel("g4SimHits", simVtcs);
00131 //    std::cout << "SimVertex " << simVtcs->size() << std::endl;
00132 //    for(edm::SimVertexContainer::const_iterator v=simVtcs->begin();
00133 //        v!=simVtcs->end(); ++v){
00134 //      std::cout << "simvtx "
00135 //             << v->position().x() << " "
00136 //             << v->position().y() << " "
00137 //             << v->position().z() << " "
00138 //             << v->parentIndex() << " "
00139 //             << v->noParent() << " "
00140 //               << std::endl;
00141 //    }
00142    return *(tPC.begin());
00143 }
00144 DEFINE_FWK_MODULE(KVFTest);