#include <LumiReWeighting.h>

Public Member Functions
	LumiReWeighting (std::string generatedFile, std::string dataFile, std::string GenHistName, std::string DataHistName)
	LumiReWeighting (std::vector< float > MC_distr, std::vector< float > Lumi_distr)
	LumiReWeighting ()
double	weight (float npv)
double	weight (int 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_
boost::shared_ptr< TH1 >	weights_

Detailed Description

Definition at line 219 of file LumiReWeighting.h.

Constructor & Destructor Documentation

LumiReWeighting::LumiReWeighting	(	std::string	generatedFile,
		std::string	dataFile,
		std::string	GenHistName = `"pileup"`,
		std::string	DataHistName = `"pileup"`
	)

Definition at line 38 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 )
      {
        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	(	std::vector< float >	MC_distr,
		std::vector< float >	Lumi_distr
	)

Definition at line 83 of file LumiReWeighting.cc.

References dtNoiseDBValidation_cfg::cerr, gather_cfg::cout, Data_distr_, FirstWarning_, MC_distr_, EnsembleCalibrationLA_cfg::NBins, OldLumiSection_, and weights_.

                                                                                            {
  // 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 228 of file LumiReWeighting.h.

{ } ;

Member Function Documentation

double LumiReWeighting::weight ( int npv )

Definition at line 137 of file LumiReWeighting.cc.

References newFWLiteAna::bin, and weights_.

Referenced by weight().

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

double LumiReWeighting::weight ( float npv )

Definition at line 142 of file LumiReWeighting.cc.

References newFWLiteAna::bin, and weights_.

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

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

Definition at line 153 of file LumiReWeighting.cc.

References edm::ProcessHistory::begin(), dtNoiseDBValidation_cfg::cerr, edm::ProcessHistory::end(), FirstWarning_, edm::EventBase::getByLabel(), edm::EventBase::processHistory(), edm::ProcessConfiguration::processName(), edm::ProcessConfiguration::releaseVersion(), AlCaHLTBitMon_QueryRunRegistry::string, and weight().

                                                      {

  // 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(edm::InputTag("addPileupInfo"), 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 201 of file LumiReWeighting.cc.

References edm::ProcessHistory::begin(), newFWLiteAna::bin, dtNoiseDBValidation_cfg::cerr, edm::ProcessHistory::end(), edm::EventBase::getByLabel(), edm::EventBase::luminosityBlock(), OldLumiSection_, 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(edm::InputTag("addPileupInfo"), 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;
 
}