EVGEN/AliGenPileup.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 //-------------------------------------------------------------------------
  17 //                          Class AliGenPileup
  18 //   This is a generator of beam-beam pileup.
  19 //   It generates interactions within 3 orbits (+-1) around
  20 //   the trigger event. The trigger event itself is chosen
  21 //   randomly among the bunch crossings within the central orbit.
  22 //   The user can decide whenever to include in the simulation the
  23 //   "trigger" interaction or not. This is handled by the
  24 //   GenerateTrigInteraction(Bool_t flag) method.
  25 //   In the case the trigger interaction is included, it is
  26 //   generated using the same settings (vertex smear for example) as
  27 //   the pileup events.
  28 //   In case the trigger simulation is not included, the user can make
  29 //   a cocktail of generator used to produce the trigger interaction and
  30 //   AliGenPileup. In this case in order to avoid a fake increase of the rate around the
  31 //   trigger, the number of background events within the bunch
  32 //   crossing of the trigger is readuced by one.
  33 //   The beam profile (the list of the active bunch crossings) can be
  34 //   controlled via the SetBCMask(const char *mask) method. The syntax
  35 //   follows the one in AliTriggerBCMask class. For example:
  36 //   "3564H" would mean that all the bunch corssings within the orbit
  37 //   are aloowed (which is of course unphysical). In case one wants to simulate
  38 //   one-bunch-crossing-per-orbit scenario, the way to do it is to put something like:
  39 //   "1H3563L" or similar.
  40 //   The SetGenerator(AliGenerator *generator, Float_t rate) method is
  41 //   used in order to define the generator to be used. The second argument is the pileup
  42 //   rate in terms of #_of_interactions/bunch-crossing = sigma_tot * luminosity.
  43 //   The pileup generation time window can be set via
  44 //   AliGenerator::SetPileUpTimeWindow(Float_t pileUpTimeW) method. By the default the
  45 //   window is set to 88micros (= TPC readout window).
  46 //
  47 // cvetan.cheshkov@cern.ch  9/12/2008
  48 //-------------------------------------------------------------------------
  49
  50 #include <TParticle.h>
  51 #include <TFormula.h>
  52
  53 #include "AliGenPileup.h"
  54 #include "AliLog.h"
  55 #include "AliGenCocktailEventHeader.h"
  56 #include "AliGenCocktailEntry.h"
  57 #include "AliRun.h"
  58 #include "AliStack.h"
  59
  60 ClassImp(AliGenPileup)
  61
  62 AliGenPileup::AliGenPileup():
  63   AliGenCocktail(),
  64   fBCMask("bcm","3564H"),
  65   fGenTrig(kFALSE),
  66   fFlag(kFALSE)
  67 {
  68 // Constructor
  69 // The pileup time window is by default
  70 // set to the TPC readout one
  71     fName = "Pileup";
  72     fTitle= "Beam-beam pileup";
  73
  74     fPileUpTimeWindow = 88e-6;
  75 }
  76
  77 AliGenPileup::~AliGenPileup()
  78 {
  79 // Destructor
  80 }
  81
  82 void AliGenPileup::SetGenerator(AliGenerator *generator, Float_t rate, Bool_t flag)
  83 {
  84   // The method sets the geenrator to be used
  85   // for pileup simulation.
  86   // The second argument is the pileup rate in terms of
  87   // #_of_interactions/bunch-crossing = sigma_tot * luminosity.
  88   // There is a protection in case the generator was already set.
  89   if (fEntries) {
  90     if (FirstGenerator()) {
  91       AliError("Pileup generator has been already set! Nothing done");
  92       return;
  93     }
  94   }
  95   AddGenerator(generator,"pileup generator",rate);
  96   fFlag = flag;
  97 }
  98
  99 void AliGenPileup::AddGenerator(AliGenerator *Generator,
 100                                 const char* Name,
 101                                 Float_t RateExp , TFormula* /*form*/)
 102 {
 103   // The method used to add the pileup generator
 104   // in the cocktail list.
 105   // The method is protected in order to avoid
 106   // its misusage
 107   AliGenCocktail::AddGenerator(Generator,Name,RateExp);
 108 }
 109
 110 Bool_t AliGenPileup::SetBCMask(const char *mask)
 111 {
 112   // Set the active bunch-crossings that
 113   // will be included in the pileup
 114   // simulation. For more details on the
 115   // syntax of the mask - see
 116   // STEER/AliTriggerBCMask.* and the comments
 117   // in the header of this file
 118   return fBCMask.SetMask(mask);
 119 }
 120
 121 void AliGenPileup::Generate()
 122 {
 123   //
 124   // Generate pileup event
 125   // For details see the coments inline
 126
 127   // Check that the pileup generator is correctly set
 128   AliGenCocktailEntry *entry = FirstGenerator();
 129   if (!entry) {
 130     AliFatal("No pileup generator entry is found!");
 131   }
 132
 133   AliGenerator *gen = entry->Generator();
 134   if (!gen) {
 135     AliFatal("No pileup generator specified!");
 136   }
 137   else if (gen->NeedsCollisionGeometry()) {
 138     AliFatal("No Collision Geometry Provided");
 139   }
 140
 141   // Check that the pileup rate is correctly set
 142   Float_t rate = entry->Rate();
 143   if (rate <= 0) {
 144     AliFatal(Form("Invalid rate value: %f",rate));
 145   }
 146
 147   // Create cocktail header
 148   if (fHeader) delete fHeader;
 149   fHeader = new AliGenCocktailEventHeader("Pileup Cocktail Header");
 150
 151   // Generate time of all
 152   // the collisions within one orbit
 153   Int_t *nIntBC = new Int_t[3*AliTriggerBCMask::kNBits];
 154   Int_t *indexBC = new Int_t[3*AliTriggerBCMask::kNBits];
 155   Int_t nTotBC = 0;
 156   while (nTotBC == 0) {
 157     for(Int_t iBC = 0; iBC <  AliTriggerBCMask::kNBits; iBC++) {
 158
 159       if (!fBCMask.GetMask(iBC)) continue;
 160
 161       //      Int_t nInteractions = gRandom->Poisson(rate);
 162       Int_t nInteractions;
 163       if (!fFlag)
 164         nInteractions = gRandom->Poisson(rate);
 165       else
 166         nInteractions = TMath::Nint(rate) + 1;
 167
 168       if (nInteractions == 0) continue;
 169
 170       nIntBC[nTotBC] = nInteractions;
 171       indexBC[nTotBC] = iBC;
 172       nTotBC++;
 173     }
 174   }
 175
 176   // Select the bunch crossing for triggered event
 177   Int_t iTrgBC = gRandom->Integer(nTotBC);
 178   // Subtract one from the number of events
 179   // generated within this bc (only in case
 180   // the user disabled the generation of the trigger
 181   // interaction)
 182   if (!fGenTrig) nIntBC[iTrgBC]--;
 183
 184   // Remove bunch crossings outside pileup
 185   // time window
 186   for(Int_t iBC = 0; iBC <  nTotBC; iBC++) {
 187    if (TMath::Abs(25e-9*(indexBC[iBC]-indexBC[iTrgBC])) > fPileUpTimeWindow)
 188      nIntBC[iBC] = 0;
 189   }
 190
 191   // Generate the two orbits around the central one
 192   // taking into account the pileup time window
 193   for(Int_t iBC = 0; iBC <  AliTriggerBCMask::kNBits; iBC++) {
 194
 195     if (!fBCMask.GetMask(iBC)) continue;
 196
 197     if (TMath::Abs(25e-9*(iBC-AliTriggerBCMask::kNBits-indexBC[iTrgBC])) > fPileUpTimeWindow) continue;
 198
 199     Int_t nInteractions = gRandom->Poisson(rate);
 200     if (nInteractions == 0) continue;
 201
 202     nIntBC[nTotBC] = nInteractions;
 203     indexBC[nTotBC] = iBC-AliTriggerBCMask::kNBits;
 204     nTotBC++;
 205   }
 206   for(Int_t iBC = 0; iBC <  AliTriggerBCMask::kNBits; iBC++) {
 207
 208     if (!fBCMask.GetMask(iBC)) continue;
 209
 210     if (TMath::Abs(25e-9*(iBC+AliTriggerBCMask::kNBits-indexBC[iTrgBC])) > fPileUpTimeWindow) continue;
 211
 212     Int_t nInteractions = gRandom->Poisson(rate);
 213     if (nInteractions == 0) continue;
 214
 215     nIntBC[nTotBC] = nInteractions;
 216     indexBC[nTotBC] = iBC+AliTriggerBCMask::kNBits;
 217     nTotBC++;
 218   }
 219
 220   // Loop over the generated collision times, call the generator
 221   // and correct the partcile times in the stack
 222   AliStack *stack = AliRunLoader::Instance()->Stack();
 223   Int_t lastpart=0;
 224   entry->SetFirst(lastpart);
 225
 226   for(Int_t iBC = 0; iBC <  nTotBC; iBC++) {
 227     Float_t deltat = 25e-9*(indexBC[iBC] - indexBC[iTrgBC]);
 228     for (Int_t i = 0; i < nIntBC[iBC]; i++) {
 229       //  Generate the vertex position and time
 230       Vertex();
 231       TArrayF eventVertex(3);
 232       for (Int_t j=0; j < 3; j++) eventVertex[j] = fVertex[j];
 233       Double_t vTime = deltat + gRandom->Gaus(0,fOsigma[2]/TMath::Ccgs());
 234
 235       gen->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2));
 236       gen->Generate();
 237
 238       for (Int_t k = lastpart; k < stack->GetNprimary(); k++) {
 239         TLorentzVector v;
 240         stack->Particle(k)->ProductionVertex(v);
 241         v[3] = vTime;
 242         stack->Particle(k)->SetProductionVertex(v);
 243       }
 244       lastpart = stack->GetNprimary();
 245
 246       // Store the interaction header in the container of the headers
 247       ((AliGenEventHeader*) fHeader->GetHeaders()->Last())->SetPrimaryVertex(eventVertex);
 248       ((AliGenEventHeader*) fHeader->GetHeaders()->Last())->SetInteractionTime(vTime);
 249     }
 250   }
 251   delete [] nIntBC;
 252   delete [] indexBC;
 253
 254   entry->SetLast(stack->GetNprimary());
 255
 256   fHeader->CalcNProduced();
 257
 258   if (fContainer) {
 259     fContainer->AddHeader(fHeader);
 260   } else {
 261     gAlice->SetGenEventHeader(fHeader);
 262   }
 263
 264 }
 265
 266 void AliGenPileup::SetRandomise(Bool_t /*flag*/)
 267 {
 268   // This setting is not implemented in
 269   // case of pileup generation
 270   // So the method gives an warning and exits
 271   AliWarning("This setting has no effect on the generator!");
 272 }
 273
 274 void AliGenPileup::UsePerEventRates()
 275 {
 276   // This setting is not implemented in
 277   // case of pileup generation
 278   // So the method gives an warning and exits
 279   AliWarning("This setting has no effect on the generator!");
 280 }