#include <PileUp.h>

Public Member Functions
double	averageNumber () const
void	CalculatePileup (int MinBunch, int MaxBunch, std::vector< int > &PileupSelection, std::vector< float > &TrueNumInteractions)
bool	doPileUp ()
void	dropUnwantedBranches (std::vector< std::string > const &wantedBranches)
void	endJob ()
void	input (unsigned int s)
const unsigned int &	input () const
	PileUp (ParameterSet const &pset, double averageNumber, TH1F *const histo, const bool playback)
template<typename T >
void	playPileUp (const std::vector< edm::EventID > &ids, T eventOperator)
bool	poisson () const
template<typename T >
void	readPileUp (edm::EventID const &signal, std::vector< edm::EventID > &ids, T eventOperator, const int NumPU)
void	reload (const edm::EventSetup &setup)
	~PileUp ()
Private Attributes
double	averageNumber_
bool	fixed_
bool	fixed_OOT_
TH1F *	h1f
TH1F *	histo_
bool	histoDistribution_
TH1F *	hprobFunction
VectorInputSource *const	input_
unsigned int	inputType_
int const	intAverage_
int	intFixed_OOT_
bool	manage_OOT_
bool	none_
bool	playback_
bool	poisson_
bool	poisson_OOT_
CLHEP::RandPoisson *	poissonDistr_OOT_
CLHEP::RandPoissonQ *	poissonDistribution_
TFile *	probFileHisto
bool	probFunctionDistribution_
bool	samelumi_
int	seed_
bool	sequential_
std::string	type_

Detailed Description

Definition at line 31 of file PileUp.h.

Constructor & Destructor Documentation

edm::PileUp::PileUp	(	ParameterSet const &	pset,
		double	averageNumber,
		TH1F *const	histo,
		const bool	playback
	)		`[explicit]`

Definition at line 20 of file PileUp.cc.

References averageNumber_, Exception, edm::ParameterSet::exists(), fixed_, fixed_OOT_, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), histoDistribution_, intFixed_OOT_, edm::Service< T >::isAvailable(), manage_OOT_, none_, poisson_, poisson_OOT_, poissonDistr_OOT_, poissonDistribution_, probFunctionDistribution_, seed_, and type_.

                                                                                                        :
    type_(pset.getParameter<std::string>("type")),
    averageNumber_(averageNumber),
    intAverage_(static_cast<int>(averageNumber)),
    histo_(histo),
    histoDistribution_(type_ == "histo"),
    probFunctionDistribution_(type_ == "probFunction"),
    poisson_(type_ == "poisson"),
    fixed_(type_ == "fixed"),
    none_(type_ == "none"),
    input_(VectorInputSourceFactory::get()->makeVectorInputSource(pset, InputSourceDescription()).release()),
    poissonDistribution_(0),
    poissonDistr_OOT_(0),
    playback_(playback),
    sequential_(pset.getUntrackedParameter<bool>("sequential", false)),
    samelumi_(pset.getUntrackedParameter<bool>("sameLumiBlock", false)),
    seed_(0)
   {
     if (pset.exists("nbPileupEvents"))
       seed_=pset.getParameter<edm::ParameterSet>("nbPileupEvents").getUntrackedParameter<int>("seed",0);

    bool DB=type_=="readDB";

    edm::Service<edm::RandomNumberGenerator> rng;
    if (!rng.isAvailable()) {
      throw cms::Exception("Configuration")
        << "PileUp requires the RandomNumberGeneratorService\n"
        "which is not present in the configuration file.  You must add the service\n"
        "in the configuration file or remove the modules that require it.";
    }
    
    CLHEP::HepRandomEngine& engine = rng->getEngine();
    poissonDistribution_ = new CLHEP::RandPoissonQ(engine, averageNumber_);
    
    // Get seed for the case when using user histogram or probability function
    if (histoDistribution_ || probFunctionDistribution_ || DB){ 
      if(seed_ !=0) {
        gRandom->SetSeed(seed_);
        LogInfo("MixingModule") << " Change seed for " << type_ << " mode. The seed is set to " << seed_;
      }
      else {
        gRandom->SetSeed(engine.getSeed());
      }
    } 
     
        
    if (!(histoDistribution_ || probFunctionDistribution_ || poisson_ || fixed_ || none_) && !DB) {
      throw cms::Exception("Illegal parameter value","PileUp::PileUp(ParameterSet const& pset)")
        << "'type' parameter (a string) has a value of '" << type_ << "'.\n"
        << "Legal values are 'poisson', 'fixed', or 'none'\n";
    }

    if (!DB){
    manage_OOT_ = pset.getUntrackedParameter<bool>("manage_OOT", false);

    if(manage_OOT_) { // figure out what the parameters are

      //      if (playback_) throw cms::Exception("Illegal parameter clash","PileUp::PileUp(ParameterSet const& pset)")
      // << " manage_OOT option not allowed with playback ";

      std::string OOT_type = pset.getUntrackedParameter<std::string>("OOT_type");

      if(OOT_type == "Poisson" || OOT_type == "poisson") {
        poisson_OOT_ = true;
        poissonDistr_OOT_ = new CLHEP::RandPoisson(engine);
      }
      else if(OOT_type == "Fixed" || OOT_type == "fixed") {
        fixed_OOT_ = true;
        // read back the fixed number requested out-of-time
        intFixed_OOT_ = pset.getUntrackedParameter<int>("intFixed_OOT", -1);
        if(intFixed_OOT_ < 0) {
          throw cms::Exception("Illegal parameter value","PileUp::PileUp(ParameterSet const& pset)") 
            << " Fixed out-of-time pileup requested, but no fixed value given ";
        }
      }
      else {
        throw cms::Exception("Illegal parameter value","PileUp::PileUp(ParameterSet const& pset)")
          << "'OOT_type' parameter (a string) has a value of '" << OOT_type << "'.\n"
          << "Legal values are 'poisson' or 'fixed'\n";
      }
      edm::LogInfo("MixingModule") <<" Out-of-time pileup will be generated with a " << OOT_type << " distribution. " ;
    }
    }
    
  }

edm::PileUp::~PileUp ( )

Definition at line 207 of file PileUp.cc.

References poissonDistr_OOT_, and poissonDistribution_.

                  {
    delete poissonDistribution_;
    delete poissonDistr_OOT_ ;
  }

Member Function Documentation

double edm::PileUp::averageNumber ( ) const [inline]

Definition at line 42 of file PileUp.h.

References averageNumber_.

Referenced by reload().

{return averageNumber_;}

void edm::PileUp::CalculatePileup	(	int	MinBunch,
		int	MaxBunch,
		std::vector< int > &	PileupSelection,
		std::vector< float > &	TrueNumInteractions
	)

Definition at line 212 of file PileUp.cc.

References averageNumber_, fixed_, histo_, histoDistribution_, intAverage_, intFixed_OOT_, manage_OOT_, none_, poisson_, poisson_OOT_, poissonDistr_OOT_, poissonDistribution_, and probFunctionDistribution_.

                                                                                                                                 {

    // if we are managing the distribution of out-of-time pileup separately, select the distribution for bunch
    // crossing zero first, save it for later.

    int nzero_crossing = -1;
    double Fnzero_crossing = -1;

    if(manage_OOT_) {
      if (none_){
        nzero_crossing = 0;
      }else if (poisson_){
        nzero_crossing =  poissonDistribution_->fire() ;
      }else if (fixed_){
        nzero_crossing =  intAverage_ ;
      }else if (histoDistribution_ || probFunctionDistribution_){
        double d = histo_->GetRandom();
        //n = (int) floor(d + 0.5);  // incorrect for bins with integer edges
        Fnzero_crossing =  d;
      }

    }

    for(int bx = MinBunch; bx < MaxBunch+1; ++bx) {

      if(manage_OOT_) {
        if(bx==0 && !poisson_OOT_) { 
          PileupSelection.push_back(nzero_crossing) ;
          TrueNumInteractions.push_back( nzero_crossing );
        }
        else{
          if(poisson_OOT_) {
            PileupSelection.push_back(poissonDistr_OOT_->fire(Fnzero_crossing)) ;
            TrueNumInteractions.push_back( Fnzero_crossing );
          }
          else {
            PileupSelection.push_back(intFixed_OOT_) ;
            TrueNumInteractions.push_back( intFixed_OOT_ );
          }  
        }
      }
      else {
        if (none_){
          PileupSelection.push_back(0);
          TrueNumInteractions.push_back( 0. );
        }else if (poisson_){
          PileupSelection.push_back(poissonDistribution_->fire());
          TrueNumInteractions.push_back( averageNumber_ );
        }else if (fixed_){
          PileupSelection.push_back(intAverage_);
          TrueNumInteractions.push_back( intAverage_ );
        }else if (histoDistribution_ || probFunctionDistribution_){
          double d = histo_->GetRandom();
          PileupSelection.push_back(int(d));
          TrueNumInteractions.push_back( d );
        }
      }
    
    }
  }

bool edm::PileUp::doPileUp ( ) [inline]

Definition at line 44 of file PileUp.h.

References averageNumber_, and none_.

{return none_ ? false :  averageNumber_>0.;}

void edm::PileUp::dropUnwantedBranches ( std::vector< std::string > const & wantedBranches ) [inline]

Definition at line 45 of file PileUp.h.

References edm::VectorInputSource::dropUnwantedBranches(), and input_.

                                                                          {
      input_->dropUnwantedBranches(wantedBranches);
    }

void edm::PileUp::endJob ( void ) [inline]

Definition at line 48 of file PileUp.h.

References edm::InputSource::doEndJob(), and input_.

                   {
      input_->doEndJob();
    }

void edm::PileUp::input ( unsigned int s ) [inline]

Definition at line 61 of file PileUp.h.

References inputType_, and alignCSCRings::s.

{inputType_=s;}

const unsigned int& edm::PileUp::input ( ) const [inline]

Definition at line 60 of file PileUp.h.

References inputType_.

{return inputType_;}

template<typename T >

void edm::PileUp::playPileUp	(	const std::vector< edm::EventID > &	ids,
		T	eventOperator
	)

Definition at line 173 of file PileUp.h.

References input_, and edm::VectorInputSource::loopSpecified().

                                                                        {
    //TrueNumInteractions.push_back( ids.size() ) ;
    input_->loopSpecified(ids,eventOperator);
  }

bool edm::PileUp::poisson ( ) const [inline]

Definition at line 43 of file PileUp.h.

References poisson_.

{return poisson_;}

template<typename T >

void edm::PileUp::readPileUp	(	edm::EventID const &	signal,
		std::vector< edm::EventID > &	ids,
		T	eventOperator,
		const int	pileEventCnt
	)

Generates events from a VectorInputSource. This function decides which method of VectorInputSource to call: sequential, random, or pre-specified. The ids are either ids to read or ids to store while reading. eventOperator has a type that matches the eventOperator in VectorInputSource::loopRandom.

The "signal" event is optionally used to restrict the secondary events used for pileup and mixing.

Definition at line 132 of file PileUp.h.

References input_, edm::VectorInputSource::loopRandom(), edm::VectorInputSource::loopRandomWithID(), edm::VectorInputSource::loopSequential(), edm::VectorInputSource::loopSequentialWithID(), fjr2json::lumi, edm::EventID::luminosityBlock(), SiPixelLorentzAngle_cfi::read, edm::EventID::run(), samelumi_, and sequential_.

                                                                                                                     {

    // One reason PileUp is responsible for recording event IDs is
    // that it is the one that knows how many events will be read.
    ids.reserve(pileEventCnt);
    RecordEventID<T> recorder(ids,eventOperator);
    int read;
    if (samelumi_) {
      const edm::LuminosityBlockID lumi(signal.run(), signal.luminosityBlock());
      if (sequential_)
        read = input_->loopSequentialWithID(lumi, pileEventCnt, recorder);
      else
        read = input_->loopRandomWithID(lumi, pileEventCnt, recorder);
    } else {
      if (sequential_) {
        // boost::bind creates a functor from recordEventForPlayback
        // so that recordEventForPlayback can insert itself before
        // the original eventOperator.

        read = input_->loopSequential(pileEventCnt, recorder);
        //boost::bind(&PileUp::recordEventForPlayback<T>,
        //                    boost::ref(*this), _1, boost::ref(ids),
        //                             boost::ref(eventOperator))
        //  );
          
      } else  {
        read = input_->loopRandom(pileEventCnt, recorder);
        //               boost::bind(&PileUp::recordEventForPlayback<T>,
        //                             boost::ref(*this), _1, boost::ref(ids),
        //                             boost::ref(eventOperator))
        //                 );
      }
    }
    if (read != pileEventCnt)
      edm::LogWarning("PileUp") << "Could not read enough pileup events: only " << read << " out of " << pileEventCnt << " requested.";
  }

void edm::PileUp::reload ( const edm::EventSetup & setup )

Definition at line 106 of file PileUp.cc.

References MixingInputConfig::averageNumber(), averageNumber(), averageNumber_, HDQMDatabaseProducer::config, alignCSCRings::e, Exception, fixed_, fixed_OOT_, MixingInputConfig::fixedOutOfTime(), edm::EventSetup::get(), histo_, histoDistribution_, i, inputType_, intFixed_OOT_, j, LogDebug, manage_OOT_, none_, MixingInputConfig::outOfTime(), poisson_, poisson_OOT_, poissonDistr_OOT_, poissonDistribution_, probFunctionDistribution_, MixingInputConfig::probFunctionVariable(), MixingInputConfig::probValue(), MixingInputConfig::type(), and type_.

                                                {
    //get the required parameters from DB.
    edm::ESHandle<MixingModuleConfig> configM;
    setup.get<MixingRcd>().get(configM);

    const MixingInputConfig & config=configM->config(inputType_);

    //get the type
    type_=config.type();
    //set booleans
    histoDistribution_=type_ == "histo";
    probFunctionDistribution_=type_ == "probFunction";
    poisson_=type_ == "poisson";
    fixed_=type_ == "fixed";
    none_=type_ == "none";
    
    if (histoDistribution_) edm::LogError("MisConfiguration")<<"type histo cannot be reloaded from DB, yet";
    
    if (fixed_){
      averageNumber_=averageNumber();
    }
    else if (poisson_)
      {
        averageNumber_=config.averageNumber();
        edm::Service<edm::RandomNumberGenerator> rng; 
        CLHEP::HepRandomEngine& engine = rng->getEngine();            
        delete poissonDistribution_;
        poissonDistribution_ = new CLHEP::RandPoissonQ(engine, averageNumber_);  
      }
    else if (probFunctionDistribution_)
      {
        //need to reload the histogram from DB
        const std::vector<int> & dataProbFunctionVar = config.probFunctionVariable();
        std::vector<double> dataProb = config.probValue();
        
        int varSize = (int) dataProbFunctionVar.size();
        int probSize = (int) dataProb.size();
                
        if ((dataProbFunctionVar[0] != 0) || (dataProbFunctionVar[varSize - 1] != (varSize - 1))) 
          throw cms::Exception("BadProbFunction") << "Please, check the variables of the probability function! The first variable should be 0 and the difference between two variables should be 1." << std::endl;
                
        // Complete the vector containing the probability  function data
        // with the values "0"
        if (probSize < varSize){
          edm::LogWarning("MixingModule") << " The probability function data will be completed with " <<(varSize - probSize)  <<" values 0.";
          
          for (int i=0; i<(varSize - probSize); i++) dataProb.push_back(0);
          
          probSize = dataProb.size();
          edm::LogInfo("MixingModule") << " The number of the P(x) data set after adding the values 0 is " << probSize;
        }
        
        // Create an histogram with the data from the probability function provided by the user           
        int xmin = (int) dataProbFunctionVar[0];
        int xmax = (int) dataProbFunctionVar[varSize-1]+1;  // need upper edge to be one beyond last value
        int numBins = varSize;
        
        edm::LogInfo("MixingModule") << "An histogram will be created with " << numBins << " bins in the range ("<< xmin << "," << xmax << ")." << std::endl;

        if (histo_) delete histo_;
        histo_ = new TH1F("h","Histo from the user's probability function",numBins,xmin,xmax); 
        
        LogDebug("MixingModule") << "Filling histogram with the following data:" << std::endl;
        
        for (int j=0; j < numBins ; j++){
          LogDebug("MixingModule") << " x = " << dataProbFunctionVar[j ]<< " P(x) = " << dataProb[j];
          histo_->Fill(dataProbFunctionVar[j]+0.5,dataProb[j]); // assuming integer values for the bins, fill bin centers, not edges 
        }
        
        // Check if the histogram is normalized
        if ( ((histo_->Integral() - 1) > 1.0e-02) && ((histo_->Integral() - 1) < -1.0e-02)){ 
          throw cms::Exception("BadProbFunction") << "The probability function should be normalized!!! " << std::endl;
        }
        averageNumber_=histo_->GetMean();
      }

    int oot=config.outOfTime();
    manage_OOT_=false;
    if (oot==1)
      {
        manage_OOT_=true;
        poisson_OOT_ = false;
        if (poissonDistr_OOT_){delete poissonDistr_OOT_; poissonDistr_OOT_=0; }
        fixed_OOT_ = true;
        intFixed_OOT_=config.fixedOutOfTime();
      }
    else if (oot==2)
      {
        manage_OOT_=true;
        poisson_OOT_ = true;
        fixed_OOT_ = false;
        if (!poissonDistr_OOT_) {
          //no need to trash the previous one if already there
          edm::Service<edm::RandomNumberGenerator> rng; 
          CLHEP::HepRandomEngine& engine = rng->getEngine();                      
          poissonDistr_OOT_ = new CLHEP::RandPoisson(engine);
        }
      }

    
  }