#include <RivetAnalyzer.h>

Inheritance diagram for RivetAnalyzer:

Public Member Functions
template<>
TH1F *	aida2root (const AIDA::IHistogram1D *aidah, const std::string &name)

template<>
TProfile *	aida2root (const AIDA::IProfile1D *aidah, const std::string &name)

virtual void	analyze (const edm::Event &, const edm::EventSetup &)

virtual void	beginJob ()

virtual void	beginRun (const edm::Run &, const edm::EventSetup &)

virtual void	endJob ()

virtual void	endRun (const edm::Run &, const edm::EventSetup &)

	RivetAnalyzer (const edm::ParameterSet &)

virtual	~RivetAnalyzer ()

Public Member Functions inherited from edm::EDAnalyzer
	EDAnalyzer ()

std::string	workerType () const

virtual	~EDAnalyzer ()

Private Member Functions
template<class AIDATYPE , class ROOTTYPE >
ROOTTYPE *	aida2root (const AIDATYPE *, const std::string &)

void	normalizeTree (AIDA::ITree &tree)

template<class AIDATYPE , class ROOTTYPE >
ROOTTYPE *	prebook (const AIDATYPE *, const std::string &)

Private Attributes
Rivet::AnalysisHandler	_analysisHandler

bool	_doFinalize

edm::InputTag	_genEventInfoCollection

edm::InputTag	_hepmcCollection

bool	_isFirstEvent

std::vector< MonitorElement * >	_mes

std::string	_outFileName

bool	_produceDQM

bool	_useExternalWeight

DQMStore *	dbe

Additional Inherited Members
Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer	ModuleType

typedef WorkerT< EDAnalyzer >	WorkerType

Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string &	baseType ()

static void	fillDescriptions (ConfigurationDescriptions &descriptions)

static void	prevalidate (ConfigurationDescriptions &)

Protected Member Functions inherited from edm::EDAnalyzer
CurrentProcessingContext const *	currentContext () const

Detailed Description

Definition at line 28 of file RivetAnalyzer.h.

Constructor & Destructor Documentation

RivetAnalyzer::RivetAnalyzer ( const edm::ParameterSet & pset )

Definition at line 25 of file RivetAnalyzer.cc.

References _analysisHandler, _genEventInfoCollection, _hepmcCollection, _produceDQM, _useExternalWeight, dbe, edm::hlt::Exception, edm::ParameterSet::exists(), edm::ParameterSet::getParameter(), cppFunctionSkipper::operator, and DQMStore::setVerbose().

                                                         : 
 _analysisHandler(),
 _isFirstEvent(true),
 _outFileName(pset.getParameter<std::string>("OutputFile")),
 //decide whether to finlaize tthe plots or not.
 //deciding not to finalize them can be useful for further harvesting of many jobs
 _doFinalize(pset.getParameter<bool>("DoFinalize")),
 _produceDQM(pset.getParameter<bool>("ProduceDQMOutput"))
 {
   //retrive the analysis name from paarmeter set
   std::vector<std::string> analysisNames = pset.getParameter<std::vector<std::string> >("AnalysisNames");
   
   _hepmcCollection = pset.getParameter<edm::InputTag>("HepMCCollection");
 
   _useExternalWeight = pset.getParameter<bool>("UseExternalWeight");
   if (_useExternalWeight) {
     if (!pset.exists("GenEventInfoCollection")){
       throw cms::Exception("RivetAnalyzer") << "when using an external event weight you have to specify the GenEventInfoProduct collection from which the weight has to be taken " ; 
     }
     _genEventInfoCollection = pset.getParameter<edm::InputTag>("GenEventInfoCollection");
   }
 
   //get the analyses
   _analysisHandler.addAnalyses(analysisNames);
 
   //go through the analyses and check those that need the cross section
   const std::set< AnaHandle, AnaHandleLess > & analyses = _analysisHandler.analyses();
 
   std::set< AnaHandle, AnaHandleLess >::const_iterator ibeg = analyses.begin();
   std::set< AnaHandle, AnaHandleLess >::const_iterator iend = analyses.end();
   std::set< AnaHandle, AnaHandleLess >::const_iterator iana; 
   double xsection = -1.;
   xsection = pset.getParameter<double>("CrossSection");
   for (iana = ibeg; iana != iend; ++iana){
     if ((*iana)->needsCrossSection())
       (*iana)->setCrossSection(xsection);
   }
   if (_produceDQM){
     // book stuff needed for DQM
     dbe = 0;
     dbe = edm::Service<DQMStore>().operator->();
     dbe->setVerbose(50);
   }  
 
 }

RivetAnalyzer::~RivetAnalyzer ( )

virtual

Definition at line 71 of file RivetAnalyzer.cc.

71 {

72 }

Member Function Documentation

template<class AIDATYPE , class ROOTTYPE >

ROOTTYPE* RivetAnalyzer::aida2root	(	const AIDATYPE *	,
		const std::string &
	)

private

template<>

TH1F* RivetAnalyzer::aida2root	(	const AIDA::IHistogram1D *	aidah,
		const std::string &	name
	)

Definition at line 91 of file RivetAnalyzer.h.

References h, i, and mergeVDriftHistosByStation::name.

                                                                                   {
   /*TH1F* h = 0;
   if (aidah->axis().isFixedBinning() ) {//equidistant binning (easier case)
     int nbins = aidah->axis().bins();
     h = new TH1F(name.c_str(), name.c_str(), nbins, aidah->axis().lowerEdge(), aidah->axis().upperEdge());
   } else {
     int nbins = aidah->axis().bins();
     //need to dyn cast, IAxis lacks polymorfism
     const LWH::VariAxis* vax = dynamic_cast<const LWH::VariAxis*>(&aidah->axis());
     if (! vax ){
       throw cms::Exception("RivetAnalyzer") << "Cannot dynamix cast an AIDA axis to VariAxis ";
     }
     double* bins = new double[nbins+1];
     for (int i=0; i<nbins; ++i) {
       bins[i] = vax->binEdges(i).first;
     }
     bins[nbins] = vax->binEdges(nbins-1).second; //take last bin right border
     h = new TH1F(name.c_str(), name.c_str(), nbins, bins);
     delete bins;
   }
   */
   TH1F* h = prebook<AIDA::IHistogram1D, TH1F>(aidah, name);
   for (int i = 0; i < aidah->axis().bins(); ++i){
     h->SetBinContent(i+1, aidah->binHeight(i));
     h->SetBinError(i+1, aidah->binError(i));
   }  
   return h;
 }

template<>

TProfile* RivetAnalyzer::aida2root	(	const AIDA::IProfile1D *	aidah,
		const std::string &	name
	)

Definition at line 121 of file RivetAnalyzer.h.

References h, i, and mergeVDriftHistosByStation::name.

                                                                                     {
   TProfile* h = prebook<AIDA::IProfile1D, TProfile>(aidah, name);
   for (int i = 0; i < aidah->axis().bins(); ++i){
     h->SetBinContent(i+1, aidah->binMean(i));
     h->SetBinError(i+1, aidah->binRms(i));
   }
   return h;
 }

void RivetAnalyzer::analyze	(	const edm::Event &	iEvent,
		const edm::EventSetup &	iSetup
	)

virtual

Implements edm::EDAnalyzer.

Definition at line 89 of file RivetAnalyzer.cc.

References _analysisHandler, _genEventInfoCollection, _hepmcCollection, _isFirstEvent, _useExternalWeight, edm::hlt::Exception, and edm::Event::getByLabel().

                                                                              {
   
   //get the hepmc product from the event
   edm::Handle<HepMCProduct> evt;
   iEvent.getByLabel(_hepmcCollection, evt);
 
   // get HepMC GenEvent
   const HepMC::GenEvent *myGenEvent = evt->GetEvent();
   //if you want to use an external weight we have to clene the GenEvent and change the weight  
   if ( _useExternalWeight ){
     HepMC::GenEvent * tmpGenEvtPtr = new HepMC::GenEvent( *(evt->GetEvent()) );
     if (tmpGenEvtPtr->weights().size() == 0) {
       throw cms::Exception("RivetAnalyzer") << "Original weight container has 0 size ";
     }
     if (tmpGenEvtPtr->weights().size() > 1) {
       edm::LogWarning("RivetAnalyzer") << "Original event weight size is " << tmpGenEvtPtr->weights().size() << ". Will change only the first one ";  
     }
     edm::Handle<GenEventInfoProduct> genEventInfoProduct;
     iEvent.getByLabel(_genEventInfoCollection, genEventInfoProduct);
     tmpGenEvtPtr->weights()[0] = genEventInfoProduct->weight();
     myGenEvent = tmpGenEvtPtr; 
   }
     
 
   //aaply the beams initialization on the first event
   if (_isFirstEvent){
     _analysisHandler.init(*myGenEvent);
     _isFirstEvent = false;
   }
 
   //run the analysis
   _analysisHandler.analyze(*myGenEvent);
 
   //if we have cloned the GenEvent, we delete it
   if ( _useExternalWeight ) 
     delete myGenEvent;
 }

void RivetAnalyzer::beginJob ( void )

virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 74 of file RivetAnalyzer.cc.

                             {
   //set the environment, very ugly but rivet is monolithic when it comes to paths
   char * cmsswbase    = getenv("CMSSW_BASE");
   char * cmsswrelease = getenv("CMSSW_RELEASE_BASE");
   std::string rivetref, rivetinfo;
   rivetref = "RIVET_REF_PATH=" + string(cmsswbase) + "/src/GeneratorInterface/RivetInterface/data:" + string(cmsswrelease) + "/src/GeneratorInterface/RivetInterface/data";
   rivetinfo = "RIVET_INFO_PATH=" + string(cmsswbase) + "/src/GeneratorInterface/RivetInterface/data:" + string(cmsswrelease) + "/src/GeneratorInterface/RivetInterface/data";
   putenv(strdup(rivetref.c_str()));
   putenv(strdup(rivetinfo.c_str()));
 }

void RivetAnalyzer::beginRun	(	const edm::Run &	iRun,
		const edm::EventSetup &	iSetup
	)

virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 85 of file RivetAnalyzer.cc.

                                                                           {
   return;
 }

void RivetAnalyzer::endJob ( void )

virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 140 of file RivetAnalyzer.cc.

140 {

141 }

void RivetAnalyzer::endRun	(	const edm::Run &	iRun,
		const edm::EventSetup &	iSetup
	)

virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 128 of file RivetAnalyzer.cc.

References _analysisHandler, _doFinalize, _outFileName, and normalizeTree().

                                                                         {
   if (_doFinalize)
     _analysisHandler.finalize();
   else {
     //if we don't finalize we just want to do the transformation from histograms to DPS
     normalizeTree(_analysisHandler.tree());
   }
   _analysisHandler.writeData(_outFileName);
 
   return;
 }

void RivetAnalyzer::normalizeTree ( AIDA::ITree & tree )

private

Definition at line 143 of file RivetAnalyzer.cc.

References _analysisHandler, _mes, _produceDQM, DQMStore::book1D(), DQMStore::bookProfile(), gather_cfg::cout, dbe, h, interpolateCardsSimple::histo, nevent, AlCaHLTBitMon_ParallelJobs::p, scaleCards::path, EgammaValidation_cff::paths, and DQMStore::setCurrentFolder().

Referenced by endRun().

                                                     {
   using namespace AIDA;
   std::vector<string> analyses = _analysisHandler.analysisNames();
   tree.ls(".", true);
   const string tmpdir = "/RivetNormalizeTmp";
   tree.mkdir(tmpdir);
   foreach (const string& analysis, analyses) {
     if (_produceDQM){
       dbe->setCurrentFolder(("Rivet/"+analysis).c_str());
       //global variables that are always present
       //sumOfWeights
       TH1F nevent("nEvt", "n analyzed Events", 1, 0., 1.);
       nevent.SetBinContent(1,_analysisHandler.sumOfWeights());
       _mes.push_back(dbe->book1D("nEvt",&nevent));
     }  
     //cross section
     //TH1F xsection("xSection", "Cross Section", 1, 0., 1.);
     //xsection.SetBinContent(1,_analysisHandler.crossSection());
     //_mes.push_back(dbe->book1D("xSection",&xsection)); 
     //now loop over the histograms
     const vector<string> paths = tree.listObjectNames("/"+analysis, true); // args set recursive listing
     std::cout << "Number of objects in AIDA tree for analysis " << analysis << " = " << paths.size() << std::endl;
     foreach (const string& path, paths) {
       IManagedObject* hobj = tree.find(path);
       if (hobj) {
         // Weird seg fault on SLC4 when trying to dyn cast an IProfile ptr to a IHistogram
         // Fix by attempting to cast to IProfile first, only try IHistogram if it fails.
         IHistogram1D* histo = 0;
         IProfile1D* prof = dynamic_cast<IProfile1D*>(hobj);
         if (!prof) histo = dynamic_cast<IHistogram1D*>(hobj);
 
         std::cout << "Converting histo " << path << " to DPS" << std::endl;
         tree.mv(path, tmpdir);
         const size_t lastslash = path.find_last_of("/");
         const string basename = path.substr(lastslash+1, path.length() - (lastslash+1));
         const string tmppath = tmpdir + "/" + basename;
 
         // If it's a normal histo:
         if (histo) {
           IHistogram1D* tmphisto = dynamic_cast<IHistogram1D*>(tree.find(tmppath));
           if (tmphisto) {
             _analysisHandler.datapointsetFactory().create(path, *tmphisto);
           }
           //now convert to root and then ME
           TH1F* h = aida2root<IHistogram1D, TH1F>(histo, basename);
           if (_produceDQM)
             _mes.push_back(dbe->book1D(h->GetName(), h));
           delete h;
           tree.rm(tmppath);
         }
         // If it's a profile histo:
         else if (prof) {
           IProfile1D* tmpprof = dynamic_cast<IProfile1D*>(tree.find(tmppath));
           if (tmpprof) {
             _analysisHandler.datapointsetFactory().create(path, *tmpprof);
           }
           //now convert to root and then ME
           TProfile* p = aida2root<IProfile1D, TProfile>(prof, basename);
           if (_produceDQM)
             _mes.push_back(dbe->bookProfile(p->GetName(), p));
           delete p;
           tree.rm(tmppath);
         }
       }
     }
   }
   tree.rmdir(tmpdir);  
 }

template<class AIDATYPE , class ROOTTYPE >

ROOTTYPE * RivetAnalyzer::prebook	(	const AIDATYPE *	aidah,
		const std::string &	name
	)

private

Definition at line 67 of file RivetAnalyzer.h.

References combineCards::bins, edm::hlt::Exception, h, i, and pileupCalc::nbins.

                                                                   {
   ROOTTYPE* h = 0;
   if (aidah->axis().isFixedBinning() ) {//equidistant binning (easier case)
     int nbins = aidah->axis().bins();
     h = new ROOTTYPE(name.c_str(), name.c_str(), nbins, aidah->axis().lowerEdge(), aidah->axis().upperEdge());
   } else {
     int nbins = aidah->axis().bins();
     //need to dyn cast, IAxis lacks polymorfism
     const LWH::VariAxis* vax = dynamic_cast<const LWH::VariAxis*>(&aidah->axis());
     if (! vax ){
       throw cms::Exception("RivetAnalyzer") << "Cannot dynamix cast an AIDA axis to VariAxis ";
     }
     double* bins = new double[nbins+1];
     for (int i=0; i<nbins; ++i) {
       bins[i] = vax->binEdges(i).first;
     }
     bins[nbins] = vax->binEdges(nbins-1).second; //take last bin right border
     h = new ROOTTYPE(name.c_str(), name.c_str(), nbins, bins);
     delete bins;
   }
   return h; 
 }