PatMuonEDMAnalyzer.cc

Go to the documentation of this file.

#include <memory>
#include <string>
#include <vector>
#include <sstream>
#include <fstream>
#include <iostream>

#include <TH1F.h>
#include <TROOT.h>
#include <TFile.h>
#include <TSystem.h>

#include "DataFormats/FWLite/interface/Event.h"
#include "DataFormats/Common/interface/Handle.h"
#include "FWCore/FWLite/interface/AutoLibraryLoader.h"

#include "DataFormats/MuonReco/interface/Muon.h"
#include "DataFormats/PatCandidates/interface/Muon.h"
#include "PhysicsTools/FWLite/interface/TFileService.h"
#include "PhysicsTools/FWLite/interface/CommandLineParser.h"

int main(int argc, char* argv[]) 
{
  // ----------------------------------------------------------------------
  // First Part: 
  //
  //  * enable the AutoLibraryLoader 
  //  * book the histograms of interest 
  //  * open the input file
  // ----------------------------------------------------------------------

  // load framework libraries
  gSystem->Load( "libFWCoreFWLite" );
  AutoLibraryLoader::enable();

  // initialize command line parser
  optutl::CommandLineParser parser ("Analyze FWLite Histograms");

  // set defaults
  parser.integerValue ("maxEvents"  ) = 1000;
  parser.integerValue ("outputEvery") =   10;
  parser.stringValue  ("outputFile" ) = "analyzeEdmTuple.root";

  // parse arguments
  parser.parseArguments (argc, argv);
  int maxEvents_ = parser.integerValue("maxEvents");
  unsigned int outputEvery_ = parser.integerValue("outputEvery");
  std::string outputFile_ = parser.stringValue("outputFile");
  std::vector<std::string> inputFiles_ = parser.stringVector("inputFiles");

  // book a set of histograms
  fwlite::TFileService fs = fwlite::TFileService(outputFile_.c_str());
  TFileDirectory dir = fs.mkdir("analyzePatMuon");
  TH1F* muonPt_  = dir.make<TH1F>("muonPt"  , "pt"  ,   100,   0., 300.);
  TH1F* muonEta_ = dir.make<TH1F>("muonEta" , "eta" ,   100,  -3.,   3.);
  TH1F* muonPhi_ = dir.make<TH1F>("muonPhi" , "phi" ,   100,  -5.,   5.);  

  // loop the events
  int ievt=0;  
  for(unsigned int iFile=0; iFile<inputFiles_.size(); ++iFile){
    // open input file (can be located on castor)
    TFile* inFile = TFile::Open(inputFiles_[iFile].c_str());
    if( inFile ){
      // ----------------------------------------------------------------------
      // Second Part: 
      //
      //  * loop the events in the input file 
      //  * receive the collections of interest via fwlite::Handle
      //  * fill the histograms
      //  * after the loop close the input file
      // ----------------------------------------------------------------------      
      fwlite::Event ev(inFile);
      for(ev.toBegin(); !ev.atEnd(); ++ev, ++ievt){
    edm::EventBase const & event = ev;
    // break loop if maximal number of events is reached 
    if(maxEvents_>0 ? ievt+1>maxEvents_ : false) break;
    // simple event counter
    if(outputEvery_!=0 ? (ievt>0 && ievt%outputEvery_==0) : false) 
      std::cout << "  processing event: " << ievt << std::endl;

    // Handle to the muon pt
    edm::Handle<std::vector<float> > muonPt;
    event.getByLabel(std::string("patMuonAnalyzer:pt"), muonPt);
    // loop muon collection and fill histograms
    for(std::vector<float>::const_iterator mu1=muonPt->begin(); mu1!=muonPt->end(); ++mu1){
      muonPt_ ->Fill( *mu1 );
    }
    // Handle to the muon eta
    edm::Handle<std::vector<float> > muonEta;
    event.getByLabel(std::string("patMuonAnalyzer:eta"), muonEta);
    for(std::vector<float>::const_iterator mu1=muonEta->begin(); mu1!=muonEta->end(); ++mu1){
      muonEta_ ->Fill( *mu1 );
    }
    // Handle to the muon phi
    edm::Handle<std::vector<float> > muonPhi;
    event.getByLabel(std::string("patMuonAnalyzer:phi"), muonPhi);
    for(std::vector<float>::const_iterator mu1=muonPhi->begin(); mu1!=muonPhi->end(); ++mu1){
      muonPhi_ ->Fill( *mu1 );
    }
      }  
      // close input file
      inFile->Close();
    }
    // break loop if maximal number of events is reached:
    // this has to be done twice to stop the file loop as well
    if(maxEvents_>0 ? ievt+1>maxEvents_ : false) break;
  }
  return 0;
}