#include <LumiReWeighting.h>

Public Member Functions
	LumiReWeighting (std::string generatedFile, std::string dataFile, std::string GenHistName="pileup", std::string DataHistName="pileup", const edm::InputTag &PileupSumInfoInputTag=edm::InputTag("addPileupInfo"))

	LumiReWeighting (const std::vector< float > &MC_distr, const std::vector< float > &Lumi_distr, const edm::InputTag &PileupSumInfoInputTag=edm::InputTag("addPileupInfo"))

	LumiReWeighting ()

double	weight (int npv)

double	weight (float npv)

double	weight (const edm::EventBase &e)

double	weightOOT (const edm::EventBase &e)

Protected Attributes
boost::shared_ptr< TH1 >	Data_distr_

boost::shared_ptr< TFile >	dataFile_

std::string	dataFileName_

std::string	DataHistName_

bool	FirstWarning_

boost::shared_ptr< TFile >	generatedFile_

std::string	generatedFileName_

std::string	GenHistName_

boost::shared_ptr< TH1 >	MC_distr_

int	OldLumiSection_

edm::InputTag	pileupSumInfoTag_

boost::shared_ptr< TH1 >	weights_

Detailed Description

Definition at line 220 of file LumiReWeighting.h.

Constructor & Destructor Documentation

LumiReWeighting::LumiReWeighting	(	std::string	generatedFile,
		std::string	dataFile,
		std::string	GenHistName = `"pileup"`,
		std::string	DataHistName = `"pileup"`,
		const edm::InputTag &	PileupSumInfoInputTag = `edm::InputTag( "addPileupInfo" )`
	)

Definition at line 39 of file LumiReWeighting.cc.

References gather_cfg::cout, Data_distr_, dataFile_, dataFileName_, DataHistName_, FirstWarning_, generatedFile_, generatedFileName_, GenHistName_, MC_distr_, EnsembleCalibrationLA_cfg::NBins, OldLumiSection_, and weights_.

                                                       :
       generatedFileName_( generatedFile), 
       dataFileName_     ( dataFile ), 
       GenHistName_      ( GenHistName ), 
       DataHistName_     ( DataHistName ),
       pileupSumInfoTag_ ( PileupSumInfoInputTag )
       {
     generatedFile_ = boost::shared_ptr<TFile>( new TFile(generatedFileName_.c_str()) ); //MC distribution
     dataFile_      = boost::shared_ptr<TFile>( new TFile(dataFileName_.c_str()) );      //Data distribution
 
     Data_distr_ = boost::shared_ptr<TH1>(  (static_cast<TH1*>(dataFile_->Get( DataHistName_.c_str() )->Clone() )) );
     MC_distr_ = boost::shared_ptr<TH1>(  (static_cast<TH1*>(generatedFile_->Get( GenHistName_.c_str() )->Clone() )) );
 
     // MC * data/MC = data, so the weights are data/MC:
 
     // normalize both histograms first
 
     Data_distr_->Scale( 1.0/ Data_distr_->Integral() );
     MC_distr_->Scale( 1.0/ MC_distr_->Integral() );
 
     weights_ = boost::shared_ptr<TH1>( static_cast<TH1*>(Data_distr_->Clone()) );
 
     weights_->SetName("lumiWeights");
 
     TH1* den = dynamic_cast<TH1*>(MC_distr_->Clone());
 
     //den->Scale(1.0/ den->Integral());
 
     weights_->Divide( den );  // so now the average weight should be 1.0
 
     std::cout << " Lumi/Pileup Reweighting: Computed Weights per In-Time Nint " << std::endl;
 
     int NBins = weights_->GetNbinsX();
 
     for(int ibin = 1; ibin<NBins+1; ++ibin){
       std::cout << "   " << ibin-1 << " " << weights_->GetBinContent(ibin) << std::endl;
     }
 
 
     FirstWarning_ = true;
     OldLumiSection_ = -1;
 }

LumiReWeighting::LumiReWeighting	(	const std::vector< float > &	MC_distr,
		const std::vector< float > &	Lumi_distr,
		const edm::InputTag &	PileupSumInfoInputTag = `edm::InputTag( "addPileupInfo" )`
	)

Definition at line 86 of file LumiReWeighting.cc.

References ecal_dqm_sourceclient-live_cfg::cerr, gather_cfg::cout, Data_distr_, FirstWarning_, MC_distr_, EnsembleCalibrationLA_cfg::NBins, OldLumiSection_, and weights_.

                                                                                                                                                    :
       pileupSumInfoTag_ ( PileupSumInfoInputTag )
 {
   // no histograms for input: use vectors
   
   // now, make histograms out of them:
 
   // first, check they are the same size...
 
   if( MC_distr.size() != Lumi_distr.size() ){   
 
     std::cerr <<"ERROR: LumiReWeighting: input vectors have different sizes. Quitting... \n";
     return;
 
   }
 
   Int_t NBins = MC_distr.size();
 
   MC_distr_ = boost::shared_ptr<TH1> ( new TH1F("MC_distr","MC dist",NBins,-0.5, float(NBins)-0.5) );
   Data_distr_ = boost::shared_ptr<TH1> ( new TH1F("Data_distr","Data dist",NBins,-0.5, float(NBins)-0.5) );
 
   weights_ = boost::shared_ptr<TH1> ( new TH1F("luminumer","luminumer",NBins,-0.5, float(NBins)-0.5) );
   TH1* den = new TH1F("lumidenom","lumidenom",NBins,-0.5, float(NBins)-0.5) ;
 
   for(int ibin = 1; ibin<NBins+1; ++ibin ) {
     weights_->SetBinContent(ibin, Lumi_distr[ibin-1]);
     Data_distr_->SetBinContent(ibin, Lumi_distr[ibin-1]);
     den->SetBinContent(ibin,MC_distr[ibin-1]);
     MC_distr_->SetBinContent(ibin,MC_distr[ibin-1]);
   }
 
   // check integrals, make sure things are normalized
 
   float deltaH = weights_->Integral();
   if(fabs(1.0 - deltaH) > 0.02 ) { //*OOPS*...
     weights_->Scale( 1.0/ weights_->Integral() );
     Data_distr_->Scale( 1.0/ Data_distr_->Integral() );
   }
   float deltaMC = den->Integral();
   if(fabs(1.0 - deltaMC) > 0.02 ) {
     den->Scale(1.0/ den->Integral());
     MC_distr_->Scale(1.0/ MC_distr_->Integral());
   }
 
   weights_->Divide( den );  // so now the average weight should be 1.0    
 
   std::cout << " Lumi/Pileup Reweighting: Computed Weights per In-Time Nint " << std::endl;
 
   for(int ibin = 1; ibin<NBins+1; ++ibin){
     std::cout << "   " << ibin-1 << " " << weights_->GetBinContent(ibin) << std::endl;
   }
 
   FirstWarning_ = true;
   OldLumiSection_ = -1;
 }

edm::LumiReWeighting::LumiReWeighting ( )

inline

Definition at line 230 of file LumiReWeighting.h.

230 { } ;

Member Function Documentation

double LumiReWeighting::weight ( int npv )

Definition at line 142 of file LumiReWeighting.cc.

References newFWLiteAna::bin, and weights_.

Referenced by cuy.ValElement::__init__(), cuy.additionElement::__init__(), cuy.superimposeElement::__init__(), cuy.graphElement::__init__(), and weight().

                                         {
   int bin = weights_->GetXaxis()->FindBin( npv );
   return weights_->GetBinContent( bin );
 }

double LumiReWeighting::weight ( float npv )

Definition at line 147 of file LumiReWeighting.cc.

References newFWLiteAna::bin, and weights_.

Referenced by cuy.ValElement::__init__(), cuy.additionElement::__init__(), cuy.superimposeElement::__init__(), and cuy.graphElement::__init__().

                                           {
   int bin = weights_->GetXaxis()->FindBin( npv );
   return weights_->GetBinContent( bin );
 }

double LumiReWeighting::weight ( const edm::EventBase & e )

Definition at line 158 of file LumiReWeighting.cc.

References edm::ProcessHistory::begin(), rpcdqm::BX, ecal_dqm_sourceclient-live_cfg::cerr, edm::ProcessHistory::end(), FirstWarning_, edm::EventBase::getByLabel(), pileupSumInfoTag_, edm::EventBase::processHistory(), edm::ProcessConfiguration::processName(), edm::ProcessConfiguration::releaseVersion(), AlCaHLTBitMon_QueryRunRegistry::string, and weight().

Referenced by cuy.ValElement::__init__(), cuy.additionElement::__init__(), cuy.superimposeElement::__init__(), and cuy.graphElement::__init__().

                                                       {
 
   // find provenance of event objects, just to check at the job beginning if there might be an issue  
 
   if(FirstWarning_) {
 
     edm::ProcessHistory PHist = e.processHistory();
     edm::ProcessHistory::const_iterator PHist_iter = PHist.begin();
 
     for(; PHist_iter<PHist.end() ;++PHist_iter) {
       edm::ProcessConfiguration PConf = *(PHist_iter);
       edm::ReleaseVersion Release =  PConf.releaseVersion() ;
       const std::string Process =  PConf.processName();
 
     }
     //    SetFirstFalse();
     FirstWarning_ = false;
   }
 
   // get pileup summary information
 
   Handle<std::vector< PileupSummaryInfo > >  PupInfo;
   e.getByLabel(pileupSumInfoTag_, PupInfo);
 
   std::vector<PileupSummaryInfo>::const_iterator PVI;
 
   int npv = -1;
   for(PVI = PupInfo->begin(); PVI != PupInfo->end(); ++PVI) {
 
     int BX = PVI->getBunchCrossing();
 
     if(BX == 0) { 
       npv = PVI->getPU_NumInteractions();
       continue;
     }
 
   }
 
   if(npv < 0) std::cerr << " no in-time beam crossing found\n! " ;
 
   return weight(npv);
 
 }

double LumiReWeighting::weightOOT ( const edm::EventBase & e )

Definition at line 206 of file LumiReWeighting.cc.

References edm::ProcessHistory::begin(), newFWLiteAna::bin, rpcdqm::BX, ecal_dqm_sourceclient-live_cfg::cerr, edm::ProcessHistory::end(), edm::EventBase::getByLabel(), edm::EventBase::luminosityBlock(), OldLumiSection_, pileupSumInfoTag_, edm::EventBase::processHistory(), edm::ProcessConfiguration::processName(), edm::ProcessConfiguration::releaseVersion(), AlCaHLTBitMon_QueryRunRegistry::string, and weights_.

                                                          {
 
 
   //int Run = e.run();
   int LumiSection = e.luminosityBlock();
   
   // do some caching here, attempt to catch file boundaries
 
   if(LumiSection != OldLumiSection_) {
 
     edm::ProcessHistory PHist = e.processHistory();
     edm::ProcessHistory::const_iterator PHist_iter = PHist.begin();
 
     for(; PHist_iter<PHist.end() ;++PHist_iter) {
       edm::ProcessConfiguration PConf = *(PHist_iter);
       edm::ReleaseVersion Release =  PConf.releaseVersion() ;
       const std::string Process =  PConf.processName();
 
     }
     OldLumiSection_ = LumiSection;
   }
 
   // find the pileup summary information
 
   Handle<std::vector< PileupSummaryInfo > >  PupInfo;
   e.getByLabel(pileupSumInfoTag_, PupInfo);
 
   std::vector<PileupSummaryInfo>::const_iterator PVI;
 
   int npv = -1;
   int npv50ns = -1;
 
   for(PVI = PupInfo->begin(); PVI != PupInfo->end(); ++PVI) {
 
     int BX = PVI->getBunchCrossing();
 
     if(BX == 0) { 
       npv = PVI->getPU_NumInteractions();
     }
 
     if(BX == 1) { 
       npv50ns = PVI->getPU_NumInteractions();
     }
 
   }
 
   // Note: for the "uncorrelated" out-of-time pileup, reweighting is only done on the 50ns
   // "late" bunch (BX=+1), since that is basically the only one that matters in terms of 
   // energy deposition.  
 
   if(npv < 0) {
     std::cerr << " no in-time beam crossing found\n! " ;
     std::cerr << " Returning event weight=0\n! ";
     return 0.;
   }
   if(npv50ns < 0) {
     std::cerr << " no out-of-time beam crossing found\n! " ;
     std::cerr << " Returning event weight=0\n! ";
     return 0.;
   }
 
   int bin = weights_->GetXaxis()->FindBin( npv );
 
   double inTimeWeight = weights_->GetBinContent( bin );
 
   double TotalWeight = 1.0;
 
   TotalWeight = inTimeWeight;
 
   return TotalWeight;
  
 }