| 16 |
|
#include "UltraFastSim.h" |
| 17 |
|
#include "UFSDataStore.h" |
| 18 |
|
|
| 19 |
+ |
int poisson(double mean, Rndm *rndmPtr) { |
| 20 |
+ |
static double oldMean = -1; |
| 21 |
+ |
static double g; |
| 22 |
+ |
if(mean != oldMean) { |
| 23 |
+ |
oldMean = mean; |
| 24 |
+ |
if(mean == 0) { |
| 25 |
+ |
g = 0; |
| 26 |
+ |
} |
| 27 |
+ |
else { |
| 28 |
+ |
g = exp(-mean); |
| 29 |
+ |
} |
| 30 |
+ |
} |
| 31 |
+ |
double em = -1; |
| 32 |
+ |
double t = 1; |
| 33 |
+ |
do { |
| 34 |
+ |
em++; |
| 35 |
+ |
t *= rndmPtr->flat(); |
| 36 |
+ |
} while(t > g); |
| 37 |
+ |
return em; |
| 38 |
+ |
} |
| 39 |
+ |
|
| 40 |
|
int main(int argc, char **argv) { |
| 41 |
|
// Generator. Process selection. LHC initialization. Histogram. |
| 42 |
|
Pythia pythia; |
| 207 |
|
|
| 208 |
|
UltraFastSim ufs; |
| 209 |
|
char jobName[256]; |
| 210 |
< |
sprintf(jobName, "%s-%4.4d", argv[1], runNumber); |
| 210 |
> |
sprintf(jobName, "%s-PU%3.3d-%4.4d", argv[1], meanPileupEventCount, runNumber); |
| 211 |
|
UFSDataStore dataStore(jobName, &ufs); |
| 212 |
|
|
| 213 |
|
// Begin event loop |
| 221 |
|
// Add pileup |
| 222 |
|
|
| 223 |
|
int pileupEventCount = 0; |
| 224 |
< |
if(meanPileupEventCount > 0) pileupEventCount = meanPileupEventCount; // * rmdmPtr->pois(); (not available yet) |
| 224 |
> |
if(meanPileupEventCount > 0) pileupEventCount = poisson(meanPileupEventCount, rndmPtr); |
| 225 |
|
for (int puEvent = 0; puEvent < pileupEventCount; ) { |
| 226 |
|
if(!pileupPythia.next()) continue; |
| 227 |
|
double vx = rndmPtr->gauss()*2.; // Mean beam spread is 2 mm in x-y and 75 mm in z |
