-- SudhirMalik - 23 July 2009

How to Make and Run FWLite Examples

Overview

TO BE WRITTEN Here is a link to the FWLite code in the CVS

Setting up CMSSW Environment

Setup a Release in CMSSW_3_1_1 ( for the first time)

Make a new directory where you will work

 
source /uscmst1/prod/sw/cms/cshrc uaf
cmscvsroot CMSSW
cd WORKINGAREA
scram p CMSSW CMSSW_3_1_1
cd CMSSW_3_1_1/src
cmsenv
addpkg PhysicsTools/FWLite V02-00-02
addpkg DataFormats/FWLite V00-13-00
scram b

If already have CMSSW_3_1_1 release and do not have TWO packages above do:

source /uscmst1/prod/sw/cms/cshrc uaf
cd WORKINGAREA/CMSSW_3_1_1/src
addpkg PhysicsTools/FWLite V02-00-02
addpkg DataFormats/FWLite V00-13-00
scram b
cmsenv

If you need to update a package

source /uscmst1/prod/sw/cms/cshrc uaf
cd WORKINGAREA/CMSSW_3_1_1/src
cvs update -r  PhysicsTools/FWLite
cvs update -r  DataFormats/FWLite
scram b
cmsenv

If you have already have the updated packages, do:

source /uscmst1/prod/sw/cms/cshrc uaf
cd WORKINGAREA/CMSSW_3_1_1/src
cmsenv

Build FWLite Executable - JetPt.exe

There is an script called newFWLiteAna.py to make the necessary directory structure and BuildFile for a new FWLite executable. Build the executable using the above script as follows

newFWLiteAna.py Analysis/SimpleExamples/myJetPt --copy=jetPt
scram b
This creates an executable myJetPt.exe

If this is the first time compiling this executable you need to run the rehash command. You do not need to run this command if you are recompiling.

Now you are ready to use myJetPt.exe

Run JetPt.exe

myJetPt.exe inputFiles=/uscms_data/d2/malik/FWLITE_PAT_TUPLE_311/RelValTTbar_311.root \
 outputEvery=1000 outputFile=myJetPt

Let us look at the command we just executed. In the above command we have the following:

  • The executable myJetPt.exe
  • The options inputFiles=... , outputEvery=... , outputFile=...

The option:

  • inputFiles=/uscms_data/d2/malik/FWLITE_PAT_TUPLE_311/RelValTTbar_311.root
    tells the executable what input file to use.
  • outputEvery=1000
    tells the executable to print out lines like Processing Event: 1000 and gives you a measure of progress of running the code.
  • outputFile=myJetPt
    tells the executable the name of the output root file where your histograms are written. In this case it is myJetPt.root

For further information on command line options see Command_line_option_parsing

Here is the output you see on the screen as a result of executing the command above:

unix > myJetPt.exe inputFiles=/uscms_data/d2/malik/FWLITE_PAT_TUPLE_311/RelValTTbar_311.root \
 outputEvery=1000 outputFile=myJetPt
------------------------------------------------------------------

Integer options:
    jobid          = -1             - jobID given by CRAB,etc. (-1 means append nothing)
    maxevents      = 0              - Maximum number of events to run over (0 for whole file)
    outputevery    = 1000           - Output something once every N events (0 for never)
    section        = 0              - This section (from 1..totalSections inclusive)
    totalsections  = 0              - Total number of sections

Bool options:
    logname        = false          - Print log name and exit

String options:
    outputfile     =                - Output filename
        'myJetPt.root'
    storeprepend   =                - Prepend location on files starting with '/store/'
        ''
    tag            =                - A 'tag' to append to output file (e.g., 'v2', etc.)
        ''

String Vector options:
    inputfiles     =                - List of input files
        /uscms_data/d2/malik/FWLITE_PAT_TUPLE_311/RelValTTbar_311.root
        
------------------------------------------------------------------
Processing Event: 1000
Processing Event: 2000
Processing Event: 3000
Processing Event: 4000
Processing Event: 5000
Processing Event: 6000
Processing Event: 7000
Processing Event: 8000
Processing Event: 9000
EventContainer Summary: Processed 9000 events.
TH1Store::write(): Successfully written to 'myJetPt.root'.

The first part of the printout is a listing of all the variables you can set from the command line and their values when the program is run. You can use --noPrint if you wish to turn this output off. See Command_line_option_parsing for more details

You should now see a root file called, in this case, myJetPt.root.

Before you open the open the root file, you may want to have rootlogon.C file in your current directory.

In order to open the root file do:

root -l myJetPt.root 
Then on the root prompt do
root [0] new TBrowser; 

A browser window opens that looks like this

sudhir_maxevt100rootBrowser.png

Click on the ROOT Files folder on the left side of the browser and then you see this:

sudhir_maxevt100rootBrowser2.png

Click on the myJetPt.root

sudhir_maxevt100root.png

If you click on the jetPt and choose log Y-scale you should see a plot like this:

jetPt

Details of the Code - JetPt.cc

You can browse the full code JetPt.cc here

The code has several parts:

  • CMS includes
    // -*- C++ -*-
    
    // CMS includes
    #include "DataFormats/FWLite/interface/Handle.h"
    #include "DataFormats/PatCandidates/interface/Jet.h"
    
    #include "PhysicsTools/FWLite/interface/EventContainer.h"
    #include "PhysicsTools/FWLite/interface/CommandLineParser.h" 
    
    // Root includes
    #include "TROOT.h"
    
    using namespace std;
    
    

If you want to use other objects such as electron, photons etc, you would add some these includes files.


       #include "DataFormats/PatCandidates/interface/MET.h"
       #include "DataFormats/PatCandidates/interface/Photon.h"
       #include "DataFormats/PatCandidates/interface/Electron.h"
       #include "DataFormats/PatCandidates/interface/Tau.h"
       #include "DataFormats/TrackReco/interface/Track.h"

  • Main Subroutine

All executables in C++ need an int main() subroutine. This is what is run when the executable starts.

 ///////////////////////////
// ///////////////////// //
// // Main Subroutine // //
// ///////////////////// //
///////////////////////////

int main (int argc, char* argv[]) 
{

* Command Line Options

Declare command line option parser

        ////////////////////////////////
        // ////////////////////////// //
        // // Command Line Options // //
        // ////////////////////////// //
        ////////////////////////////////


       // Tell people what this analysis code does and setup default options.
       optutl::CommandLineParser parser ("Plots Jet Pt");

  • Change any defaults or add any new command line options

For example if you do not specify the option outputFile=myJetPt when running the myJetPt.exe, the default name jetPt.root specified in the code line below is used.

        ////////////////////////////////////////////////
        // Change any defaults or add any new command //
        //      line options you would like here.     //
        ////////////////////////////////////////////////
        parser.stringValue ("outputFile") = "jetPt.root";

After adding any new option or changing any defaults, tells the parser to parse command line option.

        // Parse the command line arguments
        parser.parseArguments (argc, argv);

For more details on command line parsing see Command_line_option_parsing

  • Create Event Container
        //////////////////////////////////
        // //////////////////////////// //
        // // Create Event Container // //
        // //////////////////////////// //
        //////////////////////////////////

       // This object 'event' is used both to get all information from the
       // event as well as to store histograms, etc.
       fwlite::EventContainer eventCont (parser);

This object 'event' is used both to get all information from the event as well as to store histograms, etc.

  • Begin Run (e.g., book histograms, etc) Here is where you create ( "book") any histogram you want to fill later in the code.

        ////////////////////////////////////////
        // ////////////////////////////////// //
        // //         Begin Run            // //
        // // (e.g., book histograms, etc) // //
        // ////////////////////////////////// //
        ////////////////////////////////////////

        // Setup a style
        gROOT->SetStyle ("Plain");

        // Book those histograms!
        eventCont.add( new TH1F( "jetPt", "jetPt", 1000, 0, 1000) );

  • Event Loop

Loop over events in input files you specify at the command line.

        //////////////////////
        // //////////////// //
        // // Event Loop // //
        // //////////////// //
        //////////////////////

        for (eventCont.toBegin(); ! eventCont.atEnd(); ++eventCont) 
        {

Here the code

    • Creates a "handle"
    • Tries to hook the handle upto a branch containing the selectedLayer1Jets
    • Makes sure that handle hook up is successful

Noet: This is very similar to how this works inside cmsRun

           //////////////////////////////////
           // Take What We Need From Event //
           //////////////////////////////////
           fwlite::Handle< vector< pat::Jet > > jetHandle;
           jetHandle.getByLabel (eventCont, "selectedLayer1Jets");
           assert ( jetHandle.isValid() );
               

    • A handle can be treated as pointer to whatever data format it is holding. For example, we can put const vector < pat::Jet > &jetvec = *handle and then use jetvec or we can use the handle directly as a pointer as the code does below:

   
         // Loop over the jets
         const vector< pat::Jet >::const_iterator kJetEnd = jetHandle->end();
         for (vector< pat::Jet >::const_iterator jetIter = jetHandle->begin();
              kJetEnd != jetIter; 
              ++jetIter) 
         {         

The following line od code does two distinct things:

    • get pt() for this jet
    • Fills histogram stored in the event container

         eventCont.hist("jetPt")->Fill (jetIter->pt());
         } // for jetIter
      } // for eventCont

  • Clean Up Job

Nothing much to do since histograms are automatically saved

      ////////////////////////
      // ////////////////// //
      // // Clean Up Job // //
      // ////////////////// //
      ////////////////////////

     // Histograms will be automatically written to the root file
     // specificed by command line options.

     // All done!  Bye bye.
     return 0;
}

FWLite Executable Z peak Example

Build the executable. Here is a link to the zPeak.cc

newFWLiteAna.py Analysis/SimpleExamples/myZPeak --copy=zPeak
scram b
rehash
Run the executable
myZPeak.exe inputFiles=/uscms_data/d2/malik/FWLITE_PAT_TUPLE_311/RelValZMM_311.root \ 
outputEvery=1000 outputFile=myZPeak 

You should now see root file called myZPeak.root If you open this root file you browse to the following plot:

zPeak

Let us look at the the code zPeak.cc to see how we looped over the muons. The snippet to loop over muons is here:

This loop is effectively doing the

    • outer = 0.....N-1 (inclusive )
      • inner = outer..... N (indices).

This ensures that while adding muon masses we do not add mass of the muon to itself nor double count the mass additions.

 // O.k.  Let's loop through our muons and see what we can find.
      const vector< pat::Muon >::const_iterator kEndIter       = muonVec.end();
      const vector< pat::Muon >::const_iterator kAlmostEndIter = kEndIter - 1;
      for (vector< pat::Muon >::const_iterator outerIter = muonVec.begin();
           kAlmostEndIter != outerIter;
           ++outerIter)
      {
         for (vector< pat::Muon >::const_iterator innerIter = outerIter + 1;
              kEndIter != innerIter;
              ++innerIter)
         {

This ensures that only masses of opposite charged muons make a Z.

            // make sure that we have muons of opposite charge
            if (outerIter->charge() * innerIter->charge() >= 0) continue;

            // if we're here then we have one positively charged muon
            // and one negatively charged muon.

Here we get the 4-momentum of two muons, add them, get the invariant mass and fill the histogram.

            eventCont.hist("Zmass")->Fill( (outerIter->p4() + innerIter->p4()).M() );
         } // for innerIter
      } // for outerIter

Modify Z peak

  • First, go to the directory containing myZpeak.cc. We are going to use this version to make a new code

newFWLiteAna.py Analysis/SimpleExamples/myZPeakWithCuts  --copy=myZPeak.cc 

  • Edit zPeakWithCuts.cc Change

parser.stringValue ("outputFile") = "zpeak1.root";

to

parser.stringValue ("outputFile") = "zpeak2.root";

and change

vector< pat::Muon > const & muonVec = *muonHandle;

to

      // Create a new vector only with muons that pass our cuts
      vector< pat::Muon > muonVec;
      for (vector< pat::Muon >::const_iterator iter = muonHandle->begin();
           muonHandle->end() != iter;
           ++iter)
      {
         if (iter->pt() > 20 && std::abs(iter->eta()) < 2.5)
         {
            muonVec.push_back( *iter );
         }
      }
Then do:
scram b
Remember to run rehash the first time you successfully compile.

Run the code:

myZPeakWithCuts.exe inputFiles=/uscms_data/d2/malik/FWLITE_PAT_TUPLE_311/RelValZMM_311.root \
outputEvery=1000 outputFile=myZPeakModified

You should now see root file called myZPeakModified.root

Using the ROOT commands below, you can super impose the modified ZPeak and the ZPeak we had made before.

root -l 
root [0] TFile::Open("myZPeak.root")
(class TFile*)0x8366578
root [1] Zmass->Draw()
<TCanvas::MakeDefCanvas>: created default TCanvas with name c1
root [2] TFile::Open("myZPeakModified.root")
(class TFile*)0x85776a0
root [3] Zmass->SetLineColor(kRed)
root [4] Zmass->Draw("Same")
root [5] 

The superimposed plots are:

zPeak_zPeakModified

References

Command line option parsing

Command line option parsing is a method of letting you set the values of different variables when running your FWLite executable from the command line. By default several options are hooked up for you (e.g., inputFiles is the list (std::vector) of files to run over, outputFile is the name of the root file where your histograms will be stored).

You will be able to add new command line options, change the default values of the default command line options, as well as learn how to easily set these options from the command line.

Defining Variables

To define options, one gives:

  • A name,
  • A type,
    • kInteger,
    • kDouble,
    • kString,
    • kBool,
    • kIntegerVector,
    • kDoubleVector, or
    • kStringVector.
  • A description, and
  • (If desired for a non-vector type,) a default value. If this is not given for a non-vector type, 0/""/false is chosen.

For example, here are two variables I hooked up in btagTemplates.cc above:

   parser.addOption ("mode",         optutl::CommandLineParser::kInteger, 
                      "Normal(0), VQQ(1), LF(2), Wc(3)", 
                      0);   
   parser.addOption ("sampleName",   optutl::CommandLineParser::kString, 
                      "Sample name (e.g., top, Wqq, etc.)");   

Here are the default six options that are automatically hocked up:

   parser.addOption ("inputFiles",    kStringVector,
                      "List of input files");
   parser.addOption ("totalSections", kInteger,
                      "Total number of sections",
                       0);
   parser.addOption ("section",       kInteger,
                      "This section (from 1..totalSections inclusive)",
                      0);
   parser.addOption ("maxEvents",     kInteger,
                      "Maximum number of events to run over (0 for whole file)",
                      0);
   parser.addOption ("outputFile",    kString,
                      "Output filename",
                      "output.root");
   parser.addOption ("outputEvery",   kInteger,
                      "Output something once every N events (0 for never)",
                     100);

Accessing Options

To access these variables, you use one of the following access functions:
   int         &integerValue  (std::string key);
   double      &doubleValue   (std::string key);
   std::string &stringValue   (std::string key);
   bool        &boolValue     (std::string key);
   IVec        &integerVector (std::string key);
   DVec        &doubleVector  (std::string key);
   SVec        &stringVector  (std::string key);

Note that these are references and not const references, so you can use these functions to set variables as well. For example, if you wanted to change the default of maxEvents to 10, you can put:
   optutl::integerValue ("maxEvents") = 10;
before the optutl::parseArguments() function call.

Setting Command Line Options

The general way to set a command line option is -varName=value. Here are the rules:

  • General Idea : -varName=value
  • kBool values should be set with 1 for true, 0.
  • Vector values can be set in any mixture of multiple times (e.g., inputFiles=one.root inputFiles=two.root) or with a comma separated list (e.g., inputFiles=three.root,four.root).
  • Vector values can be loaded from a text file (e.g.,=inputFiles_load=fileListingRootFiles.txt=).
  • If you have defined default values for a vector in the code, you can clear the default from the command line (e.g., inputFiles_clear=true).
  • --help will print out all values, their descriptions, and their default values as well as the usage string (set by optutl::setUsageString()) and then exit.
  • By default, all variables, their descriptions, and their current values will be printed and then continue running. However, using the --noPrint option will suppress these print messages.

Try the options, example

Type
newJetPt.exe --help 
in Analysis/SimpleExamples/newJetPt/bin= directory. You see all the options available as below:

jetPt.exe - Plots Jet Pt
------------------------------------------------------------------

Integer options:
    jobid          = -1             - jobID given by CRAB,etc. (-1 means append nothing)
    maxevents      = 0              - Maximum number of events to run over (0 for whole file)
    outputevery    = 0              - Output something once every N events (0 for never)
    section        = 0              - This section (from 1..totalSections inclusive)
    totalsections  = 0              - Total number of sections

Bool options:
    logname        = false          - Print log name and exit

String options:
    outputfile     =                - Output filename
        'output.root'
    storeprepend   =                - Prepend location on files starting with '/store/'
        ''
    tag            =                - A 'tag' to append to output file (e.g., 'v2', etc.)
        ''

String Vector options:
    inputfiles     =                - List of input files
        
------------------------------------------------------------------

In the options, the most important ones are inputfiles and the outputfiles. Option logname can be used to have in Crab jobs to have the sane name for the log file and output file. Options section and totalsections are used to split the number of input file data files. For example your relval_ttbar_300pre8.filelist [ Get the data files first] below has 50 files. You can split input of these 50 files,say, by putting section equal to 5 and totalsection equal to 10.

Histogram storage:

fwlite::EventContainer not only contains what is necessary to access data in the event, but also what is needed for users to store histograms. These histograms will be automatically saved at the end of the job in the file specified by outputName command line option.

1. eventCont.add() is used to store a newly created histogram:

 eventCont.add( new TH1F("ZMass", "Mass of Z candidate", 100, 50, 150) );

 

2. eventCont.hist() is used to access the histograms:

 eventCont.hist ("myHistogram")->Fill (myVariable);

Sudhir Malik

-- Main.Sudhir Malik- 29 Jul 2009

Topic attachments
I Attachment History Action Size Date Who Comment
PNGpng jetPt.png r3 r2 r1 manage 18.1 K 2009-07-31 - 23:15 SudhirMalik  
C source code filec rootlogon.C r1 manage 0.1 K 2009-07-31 - 14:45 SudhirMalik  
PNGpng sudhir_maxevt100root.png r1 manage 18.7 K 2009-07-24 - 00:25 SudhirMalik  
PNGpng sudhir_maxevt100rootBrowser.png r1 manage 20.9 K 2009-07-30 - 21:36 SudhirMalik  
PNGpng sudhir_maxevt100rootBrowser2.png r1 manage 18.1 K 2009-07-30 - 21:52 SudhirMalik  
PNGpng zPeak.png r3 r2 r1 manage 18.2 K 2009-07-31 - 23:16 SudhirMalik  
PNGpng zPeakModified.png r3 r2 r1 manage 18.0 K 2009-07-31 - 23:17 SudhirMalik  
PNGpng zPeak_zPeakModified.png r3 r2 r1 manage 18.4 K 2009-07-31 - 23:14 SudhirMalik  
Edit | Attach | Watch | Print version | History: r21 < r20 < r19 < r18 < r17 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r21 - 2009-08-02 - SudhirMalik
 
    • Cern Search Icon Cern Search
    • TWiki Search Icon TWiki Search
    • Google Search Icon Google Search

    Sandbox All webs login

This site is powered by the TWiki collaboration platform Powered by PerlCopyright & 2008-2020 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
or Ideas, requests, problems regarding TWiki? use Discourse or Send feedback