30 #include "CLHEP/Random/RandPoissonQ.h" 31 #include "CLHEP/Random/RandPoisson.h" 49 static Double_t
GetRandom(TH1* th1, CLHEP::HepRandomEngine* rng) {
50 Int_t
nbinsx = th1->GetNbinsX();
51 Double_t* fIntegral = th1->GetIntegral();
57 Double_t
r1 = rng->flat();
58 Int_t ibin = TMath::BinarySearch(
nbinsx, fIntegral,
r1);
59 Double_t
x = th1->GetBinLowEdge(ibin + 1);
60 if (
r1 > fIntegral[ibin])
61 x += th1->GetBinWidth(ibin + 1) * (
r1 - fIntegral[ibin]) / (fIntegral[ibin + 1] - fIntegral[ibin]);
68 const std::shared_ptr<PileUpConfig>&
config,
70 const bool mixingConfigFromDB)
72 Source_type_(
config->sourcename_),
73 averageNumber_(
config->averageNumber_),
74 intAverage_(static_cast<
int>(averageNumber_)),
75 histo_(
std::make_shared<TH1F>(*
config->histo_)),
76 histoDistribution_(type_ ==
"histo"),
77 probFunctionDistribution_(type_ ==
"probFunction"),
78 poisson_(type_ ==
"poisson"),
79 fixed_(type_ ==
"fixed"),
80 none_(type_ ==
"none"),
84 ->makeVectorInputSource(
93 PoissonDistribution_(),
96 playback_(
config->playback_),
97 sequential_(
pset.getUntrackedParameter<
bool>(
"sequential",
false)) {
98 if (mixingConfigFromDB) {
106 if (
pset.existsAs<std::vector<ParameterSet> >(
"producers",
true)) {
107 std::vector<ParameterSet>
producers =
pset.getParameter<std::vector<ParameterSet> >(
"producers");
115 std::make_shared<BranchIDListHelper>(),
116 std::make_shared<ThinnedAssociationsHelper>(),
122 if (
pset.exists(
"nbPileupEvents")) {
124 edm::LogWarning(
"MixingModule") <<
"Parameter nbPileupEvents.seed is not supported";
131 <<
"PileUp requires the RandomNumberGeneratorService\n" 132 "which is not present in the configuration file. You must add the service\n" 133 "in the configuration file or remove the modules that require it.";
137 throw cms::Exception(
"Illegal parameter value",
"PileUp::PileUp(ParameterSet const& pset)")
138 <<
"'type' parameter (a string) has a value of '" <<
type_ <<
"'.\n" 139 <<
"Legal values are 'poisson', 'fixed', or 'none'\n";
168 throw cms::Exception(
"Illegal parameter value",
"PileUp::PileUp(ParameterSet const& pset)")
169 <<
" Fixed out-of-time pileup requested, but no fixed value given ";
172 throw cms::Exception(
"Illegal parameter value",
"PileUp::PileUp(ParameterSet const& pset)")
173 <<
"'OOT_type' parameter (a string) has a value of '" <<
OOT_type <<
"'.\n" 174 <<
"Legal values are 'poisson' or 'fixed'\n";
177 <<
" distribution. ";
256 edm::LogError(
"MisConfiguration") <<
"type histo cannot be reloaded from DB, yet";
267 const std::vector<int>& dataProbFunctionVar =
config.probFunctionVariable();
268 std::vector<double> dataProb =
config.probValue();
270 int varSize = (
int)dataProbFunctionVar.size();
271 int probSize = (
int)dataProb.size();
273 if ((dataProbFunctionVar[0] != 0) || (dataProbFunctionVar[varSize - 1] != (varSize - 1)))
275 <<
"Please, check the variables of the probability function! The first variable should be 0 and the " 276 "difference between two variables should be 1." 281 if (probSize < varSize) {
282 edm::LogWarning(
"MixingModule") <<
" The probability function data will be completed with " 283 << (varSize - probSize) <<
" values 0.";
285 for (
int i = 0;
i < (varSize - probSize);
i++)
286 dataProb.push_back(0);
288 probSize = dataProb.size();
289 edm::LogInfo(
"MixingModule") <<
" The number of the P(x) data set after adding the values 0 is " << probSize;
293 int xmin = (
int)dataProbFunctionVar[0];
294 int xmax = (
int)dataProbFunctionVar[varSize - 1] + 1;
295 int numBins = varSize;
297 edm::LogInfo(
"MixingModule") <<
"An histogram will be created with " << numBins <<
" bins in the range (" <<
xmin 298 <<
"," <<
xmax <<
")." << std::endl;
300 histo_.reset(
new TH1F(
"h",
"Histo from the user's probability function", numBins,
xmin,
xmax));
302 LogDebug(
"MixingModule") <<
"Filling histogram with the following data:" << std::endl;
304 for (
int j = 0;
j < numBins;
j++) {
305 LogDebug(
"MixingModule") <<
" x = " << dataProbFunctionVar[
j] <<
" P(x) = " << dataProb[
j];
306 histo_->Fill(dataProbFunctionVar[
j] + 0.5,
312 throw cms::Exception(
"BadProbFunction") <<
"The probability function should be normalized!!! " << std::endl;
317 int oot =
config.outOfTime();
324 }
else if (oot == 2) {
334 CLHEP::HepRandomEngine& engine = *
randomEngine(streamID);
342 CLHEP::HepRandomEngine& engine = *
randomEngine(streamID);
358 std::vector<int>& PileupSelection,
359 std::vector<float>& TrueNumInteractions,
364 int nzero_crossing = -1;
365 double Fnzero_crossing = -1;
385 nzero_crossing =
int(
d);
389 for (
int bx = MinBunch;
bx < MaxBunch + 1; ++
bx) {
392 PileupSelection.push_back(nzero_crossing);
393 TrueNumInteractions.push_back(nzero_crossing);
399 PileupSelection.push_back(
poissonDistr_OOT(streamID)->fire(Fnzero_crossing));
401 TrueNumInteractions.push_back(Fnzero_crossing);
409 PileupSelection.push_back(0);
410 TrueNumInteractions.push_back(0.);
426 PileupSelection.push_back(
int(
d));
427 TrueNumInteractions.push_back(
d);
std::unique_ptr< SecondaryEventProvider > provider_
OOT_type
manage out-of-time pileup setting this to True means that the out-of-time pileup will have a differen...
std::shared_ptr< ProcessConfiguration > processConfiguration_
double averageNumber() const
static Double_t GetRandom(TH1 *th1, CLHEP::HepRandomEngine *rng)
void beginStream(edm::StreamID)
Log< level::Error, false > LogError
PileUp(ParameterSet const &pset, const std::shared_ptr< PileUpConfig > &config, edm::ConsumesCollector iC, const bool mixingConfigFromDB)
std::unique_ptr< CLHEP::RandPoissonQ > const & poissonDistribution(StreamID const &streamID)
std::shared_ptr< TH1F > histo_
bool probFunctionDistribution_
T getUntrackedParameter(std::string const &, T const &) const
void endLuminosityBlock(const edm::LuminosityBlock &lumi, const edm::EventSetup &setup)
std::unique_ptr< VectorInputSource > const input_
Abs< T >::type abs(const T &t)
void beginRun(const edm::Run &run, const edm::EventSetup &setup)
void CalculatePileup(int MinBunch, int MaxBunch, std::vector< int > &PileupSelection, std::vector< float > &TrueNumInteractions, StreamID const &)
void beginLuminosityBlock(const edm::LuminosityBlock &lumi, const edm::EventSetup &setup)
Integral< F, X >::type integral(const F &f)
std::unique_ptr< CLHEP::RandPoisson > const & poissonDistr_OOT(StreamID const &streamID)
std::unique_ptr< CLHEP::RandPoisson > PoissonDistr_OOT_
Log< level::Info, false > LogInfo
std::string getReleaseVersion()
CLHEP::HepRandomEngine * randomEngine_
void setupPileUpEvent(const edm::EventSetup &setup)
std::shared_ptr< LuminosityBlockPrincipal > lumiPrincipal_
void endRun(const edm::Run &run, const edm::EventSetup &setup)
std::shared_ptr< ProductRegistry > productRegistry_
void reload(const edm::EventSetup &setup)
T const & get(Event const &event, InputTag const &tag) noexcept(false)
CLHEP::HepRandomEngine * randomEngine(StreamID const &streamID)
std::unique_ptr< CLHEP::RandPoissonQ > PoissonDistribution_
std::shared_ptr< ProcessContext > processContext_
Log< level::Warning, false > LogWarning
std::unique_ptr< EventPrincipal > eventPrincipal_
std::shared_ptr< RunPrincipal > runPrincipal_
static ParameterSetID emptyParameterSetID()
std::shared_ptr< StreamContext > streamContext_
edm::ESGetToken< MixingModuleConfig, MixingRcd > configToken_