CMS 3D CMS Logo

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