Run.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 #include "Rivet/Run.hh"
00003 #include "Rivet/AnalysisHandler.hh"
00004 #include "HepMC/IO_GenEvent.h"
00005 #include "Rivet/Projections/Beam.hh"
00006 #include <limits>
00007 
00008 namespace Rivet {
00009 
00010 
00011   Run::Run(AnalysisHandler& ah) 
00012     : _ah(ah), _xs(-1.0)
00013   { }
00014 
00015 
00016   Run::~Run() { }
00017 
00018 
00019   Run& Run::setCrossSection(const double xs) {
00020     _xs = xs;
00021     return *this;
00022   }
00023 
00024 
00025   Run& Run::setListAnalyses(const bool dolist) {
00026     _listAnalyses = dolist;
00027     return *this;
00028   }
00029 
00030 
00031   // Fill event and check for a bad read state
00032   bool Run::readEvent() {
00033     /// @todo Clear rather than new the GenEvent object per-event?
00034     _evt.reset(new GenEvent());
00035     if (_io->rdstate() != 0 || !_io->fill_next_event(_evt.get()) ) {
00036       Log::getLog("Rivet.Run") << Log::DEBUG
00037             << "Read failed. End of file?" << endl;
00038       return false;
00039     }
00040     return true;
00041   }
00042   
00043   
00044   bool Run::openFile(const std::string& evtfile) {
00045     // Set up HepMC input reader objects
00046     if (evtfile == "-") {
00047       _io.reset(new HepMC::IO_GenEvent(std::cin));
00048     } else {
00049       // Ignore the HepMC::IO_GenEvent(filename, ios) constructor, since it's only available from HepMC 2.4
00050       _istr.reset(new std::fstream(evtfile.c_str(), std::ios::in));
00051       _io.reset(new HepMC::IO_GenEvent(*_istr));
00052     }
00053     if (_io->rdstate() != 0) {
00054       Log::getLog("Rivet.Run") << Log::ERROR << "Read error on file " << evtfile << endl;
00055       return false;
00056     }
00057     return true;
00058   }
00059 
00060 
00061   bool Run::init(const std::string& evtfile) {
00062     if (!openFile(evtfile)) return false;
00063 
00064     // Read first event to define run conditions
00065     bool ok = readEvent();
00066     if (!ok) return false;
00067     if (_evt->particles_size() == 0) {
00068       Log::getLog("Rivet.Run") << Log::ERROR << "Empty first event." << endl;
00069       return false;
00070     }
00071 
00072     // Initialise AnalysisHandler with beam information from first event
00073     _ah.init(*_evt);
00074 
00075     // Set cross-section from command line
00076     if (_xs >= 0.0) {
00077       Log::getLog("Rivet.Run") 
00078         << Log::DEBUG << "Setting user cross-section = " << _xs << " pb" << endl;
00079       _ah.setCrossSection(_xs);
00080     }
00081 
00082     // List the chosen & compatible analyses if requested
00083     if (_listAnalyses) {
00084       foreach (const std::string& ana, _ah.analysisNames()) {
00085         cout << ana << endl;
00086       }
00087     }
00088 
00089     return true;
00090   }
00091 
00092 
00093   bool Run::processEvent() {
00094     // Set cross-section if found in event and not from command line
00095     #ifdef HEPMC_HAS_CROSS_SECTION
00096     if (_xs < 0.0 && _evt->cross_section()) {
00097       const double xs = _evt->cross_section()->cross_section(); //< in pb
00098       Log::getLog("Rivet.Run")
00099         << Log::DEBUG << "Setting cross-section = " << xs << " pb" << endl;
00100       _ah.setCrossSection(xs);
00101     }
00102     #endif
00103     // Complain about absence of cross-section if required!
00104     if (_ah.needCrossSection() && !_ah.hasCrossSection()) {
00105       Log::getLog("Rivet.Run") 
00106         << Log::ERROR
00107         << "Total cross-section needed for at least one of the analyses. "
00108         << "Please set it (on the command line with '-x' if using the 'rivet' program)" << endl;
00109       return false;
00110     }
00111      
00112     // Analyze event
00113     _ah.analyze(*_evt);
00114  
00115     return true;
00116   }
00117 
00118 
00119   bool Run::finalize() {
00120     _evt.reset();
00121     _istr.reset();
00122     _io.reset();
00123     return true;
00124   }
00125 
00126 
00127 
00128 
00129 }