[u/mrichter/AliRoot.git] / PWGLF / FORWARD / analysis2 / sim / EGConfig.C

/**
 * @file   EGConfig.C
 * @author Christian Holm Christensen <cholm@nbi.dk>
 * @date   Thu Oct 16 11:01:38 2014
 * 
 * @brief  Specific configuration for event generator. 
 * 
 * 
 */
struct EGCfg : public VirtualEGCfg
{
  Int_t   hftype;     // Heavy flavour type (random)
  Bool_t  fIsLego;    //
  
  EGCfg() 
    : hftype(-1), fIsLego(false);
  {
    hftype = HFType();
  }
  virtual Bool_t IsLego() const { return fIsLego; }
protected:
  /** 
   * Make a random Heavy Flavour type 
   * 
   * @return Heavy flavour type 
   */
  Int_t HFType() const 
  {
    Int_t type = -1;
    Int_t r = gRandom->Rndm();
    if      (r < 0.16) type = 0; 
    else if (r < 0.32) type = 1; 
    else if (r < 0.48) type = 2; 
    else if (r < 0.64) type = 3; 
    else if (r < 0.72) type = 4; 
    else if (r < 0.80) type = 5; 
    else               type = 6; 
    return r;
  }
  /** 
   * Make the generator 
   * 
   * @param rt   Generator ID
   * @param b1   Least impact parameter 
   * @param b2   Largest impact parameter 
   *
   * @return Point to newly allocated generator or null 
   */
  AliGenerator* CreateGenerator(const TString& rt,
				Float_t b1, Float_t b2)
  {
    Bool_t        asym = grp->IsPA()||grp->IsAP();
    AliGenerator* g    = 0;
    TString t(rt);
    if      (t.EqualTo("default"))            t = DeduceRunType();
    if      (t.EndsWith("perugia0chadr"))     g=PythiaHF(0);
    else if (t.EndsWith("perugia0bchadr"))    g=PythiaHF(1);
    else if (t.EndsWith("perugia0cele"))      g=PythiaHF(2);
    else if (t.EndsWith("perugia0bele"))      g=PythiaHF(3);
    else if (t.EndsWith("perugia0jspi2e"))    g=PythiaHF(4);
    else if (t.EndsWith("perugia0btojspi2e")) g=PythiaHF(5);
    else if (t.BeginsWith("pythia"))          g=Pythia(rt);
    else if (t.BeginsWith("hijing2000hf"))    g=HFCocktail(rt, b1, b2);
    else if (t.BeginsWith("hijing2000"))      g=Hijing(b1, b2, asym, 
							false, 2.3);
    else if (t.BeginsWith("hijing"))          g=Hijing(b1, b2, asym, 
							grp->IsAA(), 0);
    else if (t.BeginsWith("ampthf"))          g=HFCocktail(rt,b1,b2);
    else if (t.BeginsWith("ampt"))            g=Ampt(b1, b2);
    else if (t.BeginsWith("dpmjet"))          g=Dpmjet(b1, b2);
    else if (t.BeginsWith("phojet"))          g=Dpmjet(b1, b2);
    else if (t.BeginsWith("hydjet"))          g=Hydjet(b1, b2);
    else if (t.BeginsWith("lego"))            g=Lego(rt);
    if (!g && !fIsLego)
      Fatal("", "Invalid run type \"%s\" specified", t.Data());
    return g;
  }
  /** 
   * Make our decayer 
   * 
   * @return Newly allocated decayer or null
   */
  TVirtualMCDecayer* CreateDecayer(const TString& runType)
  {
    if (runType.BeginsWith("hydjet")) return 0;

    LoadPythia();
    TVirtualMCDecayer* decayer = new AliDecayerPythia();
    if (runType.EqualTo("hijing2000hf") && hftype < 2) 
      decayer->SetForceDecay(kHadronicD);
    else 
      decayer->SetForceDecay(kAll);
    decayer->Init();
    return decayer;
  }

  // === PYTHIA ========================================================
  // Normal 
  /** 
   * Greate a pythia6 event generator 
   * 
   * @param tune Possible tune 
   * 
   * @return newly allocated generator or null
   */
  AliGenerator* Pythia(const TString & tune)
  {
    // Int_t kCTEQ6l = 8;
    if (!grp->IsPP()) Fatal("Setup", "Pythia6 only works for pp");

    TString t(tune);
    t.ToUpper();
    t.ReplaceAll("PYTHIA6", "");
    t.ReplaceAll("PYTHIA", "");
    Info("Setup", "Making Pythia6 event generator (tune: %s)", t.Data());

    LoadPythia();
    AliGenPythia* pythia = new AliGenPythia(-1); 
    pythia->SetMomentumRange(0, 999999.);
    pythia->SetThetaRange(0., 180.);
    pythia->SetYRange(-12.,12.);
    pythia->SetPtRange(0,1000.);
    pythia->SetProcess(kPyMb);
    pythia->SetEnergyCMS(grp->energy);

    if (t == "D6T") {
      //    Tune
      //    109     D6T : Rick Field's CDF Tune D6T 
      //                  (NB: needs CTEQ6L pdfs externally)
      pythia->SetTune(109); // F I X 
      pythia->SetStrucFunc(kCTEQ6l);
    }
    else if (t == "PERUGIA0") { 
      //    Tune
      //    320     Perugia 0
      pythia->SetTune(320); 
      pythia->UseNewMultipleInteractionsScenario();
    }
    else if (t == "ATLAS") {
      //    Tune
      //    C   306 ATLAS-CSC: Arthur Moraes' (new) ATLAS tune 
      //                        (needs CTEQ6L externally)
      pythia->SetTune(306);
      pythia->SetStrucFunc(kCTEQ6l);
    }
    else if (t == "JETS") { 
      pythia->SetProcess(kPyJets);
      pythia->SetStrucFunc(kCTEQ6l);
      pythia->SetJetEtaRange(-1.5, 1.5); 
      pythia->SetJetEtRange(50., 800.);
      pythia->SetPtHard(45., 1000.);
      pythia->SetPycellParameters(2.2, 300, 432, 0., 4., 5., 0.7);
    }
    else if (t == "ATLAS_FLAT") {
      // set high multiplicity trigger
      // this weight achieves a flat multiplicity distribution
      TH1 *weight = new TH1D("weight","weight",201,-0.5,200.5);
      weight->SetBinContent(1,5.49443);
      weight->SetBinContent(2,8.770816);
      weight->SetBinContent(6,0.4568624);
      weight->SetBinContent(7,0.2919915);
      weight->SetBinContent(8,0.6674189);
      weight->SetBinContent(9,0.364737);
      weight->SetBinContent(10,0.8818444);
      weight->SetBinContent(11,0.531885);
      weight->SetBinContent(12,1.035197);
      weight->SetBinContent(13,0.9394057);
      weight->SetBinContent(14,0.9643193);
      weight->SetBinContent(15,0.94543);
      weight->SetBinContent(16,0.9426507);
      weight->SetBinContent(17,0.9423649);
      weight->SetBinContent(18,0.789456);
      weight->SetBinContent(19,1.149026);
      weight->SetBinContent(20,1.100491);
      weight->SetBinContent(21,0.6350525);
      weight->SetBinContent(22,1.351941);
      weight->SetBinContent(23,0.03233504);
      weight->SetBinContent(24,0.9574557);
      weight->SetBinContent(25,0.868133);
      weight->SetBinContent(26,1.030998);
      weight->SetBinContent(27,1.08897);
      weight->SetBinContent(28,1.251382);
      weight->SetBinContent(29,0.1391099);
      weight->SetBinContent(30,1.192876);
      weight->SetBinContent(31,0.448944);
      for (Int_t i = 32; i <= 201; i++) weight->SetBinContent(i,1);
      weight->SetEntries(526);
      
      Int_t limit = weight->GetRandom();
      pythia->SetTriggerChargedMultiplicity(limit, 1.4);
    }
    return pythia;
  }
  /** 
   * Create a Pythia6 generator for high-flavor physics 
   * 
   * @param type    Which kind 
   * @param harder  If true, make harder processes 
   * 
   * @return Newly allocated generator or null
   */
  AliGenerator* PythiaHF(Int_t type, Bool_t harder=0) 
  { 
    LoadPythia();
    if (type == 6) return Pythia("jets");
    if (type == 4) { 
      AliGenParam *jpsi =  AliGenParam(1, AliGenMUONlib::kJpsi,
				       (harder?"CDF pp 8.8":"CDF pp 7"),"Jpsi");
      jpsi->SetPtRange(0.,999.);
      jpsi->SetYRange(-1.0, 1.0);
      jpsi->SetPhiRange(0.,360.);
      jpsi->SetForceDecay(kDiElectron);
      return jpsi;
    }
    AliGenPythia* pythia = static_cast<AliGenPythia*>(Pythia("PERUGIA0"));
    switch (type) { 
    case 0: // chadr
      pythia->SetProcess(kPyCharmppMNRwmi);
      pythia->SetForceDecay(kHadronicD);
      break;
    case 1: // bchadr
      pythia->SetProcess(kPyBeautyppMNRwmi);
      pythia->SetForceDecay(kHadronicD);
      break;
    case 2: // cele
      pythia->SetProcess(kPyCharmppMNRwmi);
      pythia->SetCutOnChild(1);
      pythia->SetPdgCodeParticleforAcceptanceCut(11);
      pythia->SetChildYRange(-1.2,1.2);
      pythia->SetChildPtRange(0,10000.);
      break;
    case 3: // bele
      pythia->SetProcess(kPyBeautyppMNRwmi);
      pythia->SetCutOnChild(1);
      pythia->SetPdgCodeParticleforAcceptanceCut(11);
      pythia->SetChildYRange(-1.2,1.2);
      pythia->SetChildPtRange(0,10000.);
      break;
    case 5:
      pythia->SetProcess(kPyBeautyppMNRwmi);
      pythia->SetCutOnChild(1);
      pythia->SetPdgCodeParticleforAcceptanceCut(443);
      pythia->SetChildYRange(-2,2);
      pythia->SetChildPtRange(0,10000.);
    }
    return pythia;
  }
  /** 
   * Make a Min-Bias AA, pA, or Ap Hijing generator 
   * 
   * @param minB    Least impact parameter 
   * @param maxB    Largest impact parameter 
   * @param quench  If true, enable quenching 
   * @param slowN   If true, make a cocktail with slow neutrons 
   * @param ptHard  Hard pT cut-off 
   * 
   * @return Generator 
   */
  AliGenerator* Hijing(Float_t minB, 
		       Float_t maxB, 
		       Bool_t  slowN=false, 
		       Bool_t  quench=1, 
		       Float_t ptHard=0) 
  {
    LoadHijing();
    AliGenHijing *gener = new AliGenHijing(-1);
    // centre of mass energy 
    gener->SetEnergyCMS(grp->energy);
    gener->SetImpactParameterRange(minB, maxB);	
    // reference frame
    gener->SetReferenceFrame("CMS");
    // projectil
    gener->SetTarget    (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
    gener->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
    // tell hijing to keep the full parent child chain
    gener->KeepFullEvent();
    // enable jet quenching
    gener->SetJetQuenching(quench);
    // enable shadowing
    gener->SetShadowing(slowN);
    // Don't track spectators
    gener->SetSpectators(!slowN);
    // 
    if (ptHard > 0) hi->SetPtHardMin(ptHard);

    // kinematic selection
    gener->SetSelectAll(0);
    // Boosted CMS 
    gener->SetBoostLHC(grp->IsPA() || grp->IsAP());
    // No need for cocktail 
    if (!slowN || !grp->IsPA() || !grp->IsAP()) return gener;

    AliGenCocktail* cocktail = new AliGenCocktail();
    cocktail->SetTarget    (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
    cocktail->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
    cocktail->SetEnergyCMS(grp->energy);

    AliGenSlowNucleons*     gray  = new AliGenSlowNucleons(1);
    AliCollisionGeometry*   coll  = gener->CollisionGeometry();
    AliSlowNucleonModelExp* model = new AliSlowNucleonModelExp();
    //  Not yet in the release...
    //      model->SetSaturation(kTRUE);
    gray->SetSlowNucleonModel(model);
    gray->SetTarget(grp->beam1.a, grp->beam1.z);
    gray->SetThetaDist(1);
    gray->SetProtonDirection(grp->beam1.IsP() ? 1 : 2);
    //      gray->SetDebug(1);
    gray->SetNominalCmsEnergy(2*grp->beamEnergy);
    gray->NeedsCollisionGeometry();
    gray->SetCollisionGeometry(coll);

    cocktail->AddGenerator(gener, "Hijing pPb", 1);
    cocktail->AddGenerator(gray, "Gray Particles", 1);
    
    return cocktail;
  }
  /** 
   * Make a DPMJet generator for pp, AA, pA, or Ap. 
   * 
   * @param fragments If true, make fragments 
   * 
   * @return Generator 
   */
  AliGenerator* Dpmjet(Float_t minB, Float_t maxB, 
		       Bool_t fragments=0)
  {
    LoadDpmjet();
    AliGenDPMjet* dpmjet = new AliGenDPMjet(-1); 
    dpmjet->SetEnergyCMS(grp->energy);
    dpmjet->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
    dpmjet->SetTarget    (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
    dpmjet->SetImpactParameterRange(minB, maxB);
    dpmjet->SetProjectileBeamEnergy(grp->beam2.z*grp->beamEnergy/grp->beam2.a);
    if (grp->IsAA()) { 
      dpmjet->SetPi0Decay(0);
    }
    else if (grp->IsPA() || grp->IsAP()) { 
      // dpmjet->SetTriggerParticle(3312, 1.2, 2.0);
      dpmjet->SetFragmentProd(false/*fragments*/); // Alwas disabled 
    }
    else if (grp->IsPP()) { // PhoJet
      dpmjet->SetMomentumRange(0, 999999.);
      dpmjet->SetThetaRange(0., 180.);
      dpmjet->SetYRange(-12.,12.);
      dpmjet->SetPtRange(0,1000.);
    }
    return dpmjet;
  }
  /** 
   * Make an AMPT generator for AA collisions 
   * 
   * @return Generator 
   */
  AliGenerator* Ampt(Float_t minB, Float_t maxB)
  {
    LoadAmpt();
    AliGenAmpt *genHi = new AliGenAmpt(-1);
    genHi->SetEnergyCMS(grp->energy);
    genHi->SetReferenceFrame("CMS");
    genHi->SetTarget    (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
    genHi->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
    genHi->SetPtHardMin (2);
    genHi->SetImpactParameterRange(minB,maxB);
    // disable jet quenching
    genHi->SetJetQuenching(0); 
    // enable shadowing
    genHi->SetShadowing(1);    
    // neutral pion and heavy particle decays switched off
    genHi->SetDecaysOff(1);
    genHi->SetSpectators(0);   // track spectators 
    genHi->KeepFullEvent();
    genHi->SetSelectAll(0);
    return genHi;
  }
  /** 
   * Make an HydJet generator for A-A
   * 
   * @return Generator 
   */
  AliGenerator* Hydjet(Float_t minB, Float_t maxB)
  {
    LoadHydjet();
    AliGenUHKM *genHi = new AliGenUHKM(-1);
    genHi->SetAllParametersLHC();
    genHi->SetTarget    (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
    genHi->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
    genHi->SetEcms(grp->energy);
    genHi->SetEnergyCMS(grp->energy);
    genHi->SetBmin(minB);
    genHi->SetBmax(maxB);
    genHi->SetPyquenPtmin(9);
    return genHi;
  }
  // === Lego ========================================================
  /** 
   * Greate a lego event generator 
   * 
   * @param tune Possible tune 
   * 
   * @return newly allocated generator or null
   */
  AliGenerator* Lego(const TString & variant)
  {
    fIsLego = true;
    return 0;
#if 0
    TString v(variant);
    v.ToUpper();
    v.ReplaceAll("LEGO", "");
    Info("Setup", "Making Lego event generator (variant: %s)", v.Data());

    AliLegoGenerator* ret = 0;
    // XYZ varies origin of the particles in two dimensions:
    //  X:  o=(0,t1,t2), p=(1,0,0)
    //  Y:  o=(t1,0,t2), p=(0,1,0)
    //  Z:  o=(t1,t2,0), p=(0,0,1)
    // PhiZ varies the momentum in two dimensions
    //  o=(0,0,t1) p=(cos(t2),sin(t2),0)
    // Eta varies momentum in two dimensions
    //  phi=t1
    //  theta=2*atan(exp(-t2))
    //  o=(0,0,0) p=(cos(phi)*sin(theta),sin(phi)*cos(theta),cos(theta))
    // Base varies in two dimensions
    //  phi=t1
    //  theta=t2
    //  o=(0,0,0) p=(cos(phi)*sin(theta),sin(phi)*cos(theta),cos(theta))    
    if (v.BeginsWith("X") || v.BeginsWith) {
      const char* c[] = { v(0), '\0' };
      ret = new AliLegoGeneratorXYZ(c);
      ret->SetCoor1Range(10,-2,2);   // Y, X
      ret->SetCoor2Range(10,-10,10); // Z
    }
    else if (v.BeginsWith("Z")) {
      ret = new AliLegoGeneratorXYZ("Z");
      ret->SetCoor1Range(10,-2,2);   // X
      ret->SetCoor2Range(10,-2,2);   // Y
    }
    else if (v.BeginsWith("PHIZ")) {
      ret = new AliLegoGeneratorPhiZ();
      ret->SetCoor1Range(10,-10,10); // Z
      ret->SetCoor2Range(360,0,360); // phi
    }
    else if (v.BeginsWith("ETA")) {
      ret = new AliLegoGeneratorEta();
      ret->SetCoor1Range(360,0,360); // phi
      Double_t aEta = 7;
      Double_t dEta = (6--4)/200;
      ret->SetCoor2Range(2*aEta/dEta,-aEta,+aEta); // Eta
    }
    else {
      ret = new AliLegoGenerator();
      ret->SetCoor1Range(180,0,180); // theta
      ret->SetCoor1Range(360,0,360); // phi
    }
    return ret;
#endif
  }

  /** 
   * Make a heavy flavour cocktail 
   * 
   * @param base Underlying event. 
   * 
   * @return Generator 
   */
  AliGeneator* HFCocktail(const TString& base, Float_t minB, Float_t maxB) 
  {
    
    AliGenCocktail *cocktail = new AliGenCocktail();
    cocktail->SetTarget    (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
    cocktail->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
    cocktail->SetEnergyCMS(grp->energy);
    
    // Add underlying event
    if (base.BeginsWith("ampt", TString::kIgnoreCase)) { 
      hi = Ampt(minB, maxB);
      cocktail->AddGenerator(hi,"ampt",1);
    }
    else { 
      hi = Hijing(minB, maxB, grp->IsPA() || grp->IsAP(), false, 2.3);
      cocktail->AddGenerator(hi,"hijing",1);

    }

    // --- Default formula -------------------------------------------
    TForumla* one = new TFormula("one", "1.");

    // --- Pythia ----------------------------------------------------
    AliGenerator* pythia = PythiaHF(hftype);
    switch (hftype) { 
    case 6: 
      cocktail->AddGenerator(pythia, "pythiaJets", 1, one);
      break;
    defualt:
      cocktail
	->AddGenerator(pythia, "pythiaHF", 1, 
		       new TFormula("Signals", 
				    "20.*(x<5.)+80./3.*(1.-x/20.)*(x>5.)"));
      break;
    }
    // --- Dummy -----------------------------------------------------
    AliGenParam* param = 0;

    // --- Phos stuff ------------------------------------------------
    AliGenPHOSlib *plib = new AliGenPHOSlib();
    Double_t lower[] = { 0, 3, 6, 9, 12, -1 };
    Double_t *pLow   = lower;
    for (Int_t i = 0; i < 5; i++) {
      param = new AliGenParam(5, plib, AliGenPHOSlib::kPi0);
      param->SetPhiRange(0., 360.) ;
      param->SetYRange(-1.2, 1.2) ;
      param->SetPtRange(lower[i], 30.) ;
      cocktail->AddGenerator(param,Form("Pi0HagPt%d", i), 1., one);
    }

    // --- Jpsi->mu+ mu-----------------------------------------------
    param = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF pp 3.94", "Jpsi");
    param->SetPtRange(0.,999.);
    param->SetYRange(-4.2, -2.3);
    param->SetPhiRange(0.,360.);
    param->SetForceDecay(kDiMuon);
    param->SetCutOnChild(1);
    param->SetChildPhiRange(0.,360.);
    param->SetChildThetaRange(168.5,178.5);
    cocktail->AddGenerator(param, "Jpsi2M", 1, one); 

    // --- Chi_c -> J/Psi + gamma, J/Psi -> e+e- ---------------------
    Float_t thmin  = (180./TMath::Pi())*2.*atan(exp(-1.2));  
    Float_t thmax  = (180./TMath::Pi())*2.*atan(exp( 1.2));  
    param  = new AliGenParam(1, AliGenMUONlib::kChic,"default","Chic");
    param->SetMomentumRange(0, 999.);        // Wide cut on the momentum
    param->SetPtRange(0, 100.);              // Wide cut on Pt
    param->SetYRange(-2.5, 2.5);
    param->SetCutOnChild(1);                 // Enable cuts on decay products
    param->SetChildPhiRange(0., 360.);
    // In the acceptance of the Central Barrel
    param->SetChildThetaRange(thmin, thmax); 
    // Chi_c -> J/Psi + gamma, J/Psi -> e+e-
    param->SetForceDecay(kChiToJpsiGammaToElectronElectron); 
    cocktail->AddGenerator(param, "Chi_c", 1, one); 

    // --- Dummy -----------------------------------------------------
    AliGenBox* box = 0;

    // --- Some particles --------------------------------------------
    Double_t    boxR   = gRandom->Integer(3)+1;
    Int_t       boxP[] = { 310, 3122,   3312, 3322, 
			   (boxR==1 ?3334: boxR==2 ?-3334: -3312) };
    Int_t       boxN[] = {   1,    1,      3,    3, 1 }
    const char* boxT[] = { "K0s", "Lambda", "Xi-", "Xi0", 
			   (boxR==1 ? "Omega-": boxR==2 ? "Omega+": "Xi+") };
    for (Int_t i = 0; i < 5; i++) {
      box = new AliGenBox(boxN[i]);
      box->SetPart(boxP[i]);
      box->SetPtRange(0,13);
      box->SetThetaRange(45, 135);
      cocktail->AddGenerator(box, boxT[i], 1, one);
    }
    
    // --- High pT charged particle ----------------------------------
    TFormula* hptF = new TFormula("High Pt", 
				  "5.*(x<5.)+20./3.*(1.-x/20.)*(x > 5.)");
    Int_t       hptP[] = { 211, 321, 2212 };
    const char* hptT[] = { "pi", "K", "p" };
    for (Int_t i = 0; i < 3; i++) { 
      for (Int_t j = -1; j <= 1; j++) {
	box->SetPart(j*hptP[i]);
	box->SetPtRange(2., 50.);
	box->SetThetaRange(thmin, thmax);
	cocktail->AddGenerator(box, Form("%s%c",hptT[i],(j<0,'-','+')),1,hptF);
      }
    }
    return cocktail;
  }

};

void 
EGConfig()
{
  ::Info("EGConfig", "Creating EG factory");
  egCfg = new EGCfg;
}

// 
// EOF
//
Commit	Line	Data
66ff3038	1	/**
	2	* @file EGConfig.C
	3	* @author Christian Holm Christensen <cholm@nbi.dk>
	4	* @date Thu Oct 16 11:01:38 2014
	5	*
	6	* @brief Specific configuration for event generator.
	7	*
	8	*
	9	*/
	10	struct EGCfg : public VirtualEGCfg
	11	{
	12	Int_t hftype; // Heavy flavour type (random)
15cb4a2b	13	Bool_t fIsLego; //
15cb4a2b	14
66ff3038	15	EGCfg()
15cb4a2b	16	: hftype(-1), fIsLego(false);
66ff3038	17	{
	18	hftype = HFType();
	19	}
15cb4a2b	20	virtual Bool_t IsLego() const { return fIsLego; }
66ff3038	21	protected:
	22	/**
	23	* Make a random Heavy Flavour type
	24	*
	25	* @return Heavy flavour type
	26	*/
	27	Int_t HFType() const
	28	{
	29	Int_t type = -1;
	30	Int_t r = gRandom->Rndm();
	31	if (r < 0.16) type = 0;
	32	else if (r < 0.32) type = 1;
	33	else if (r < 0.48) type = 2;
	34	else if (r < 0.64) type = 3;
	35	else if (r < 0.72) type = 4;
	36	else if (r < 0.80) type = 5;
	37	else type = 6;
	38	return r;
	39	}
	40	/**
	41	* Make the generator
	42	*
	43	* @param rt Generator ID
	44	* @param b1 Least impact parameter
	45	* @param b2 Largest impact parameter
	46	*
	47	* @return Point to newly allocated generator or null
	48	*/
	49	AliGenerator* CreateGenerator(const TString& rt,
	50	Float_t b1, Float_t b2)
	51	{
	52	Bool_t asym = grp->IsPA()\|\|grp->IsAP();
	53	AliGenerator* g = 0;
15cb4a2b	54	TString t(rt);
	55	if (t.EqualTo("default")) t = DeduceRunType();
	56	if (t.EndsWith("perugia0chadr")) g=PythiaHF(0);
	57	else if (t.EndsWith("perugia0bchadr")) g=PythiaHF(1);
	58	else if (t.EndsWith("perugia0cele")) g=PythiaHF(2);
	59	else if (t.EndsWith("perugia0bele")) g=PythiaHF(3);
	60	else if (t.EndsWith("perugia0jspi2e")) g=PythiaHF(4);
	61	else if (t.EndsWith("perugia0btojspi2e")) g=PythiaHF(5);
	62	else if (t.BeginsWith("pythia")) g=Pythia(rt);
	63	else if (t.BeginsWith("hijing2000hf")) g=HFCocktail(rt, b1, b2);
	64	else if (t.BeginsWith("hijing2000")) g=Hijing(b1, b2, asym,
66ff3038	65	false, 2.3);
15cb4a2b	66	else if (t.BeginsWith("hijing")) g=Hijing(b1, b2, asym,
66ff3038	67	grp->IsAA(), 0);
15cb4a2b	68	else if (t.BeginsWith("ampthf")) g=HFCocktail(rt,b1,b2);
	69	else if (t.BeginsWith("ampt")) g=Ampt(b1, b2);
	70	else if (t.BeginsWith("dpmjet")) g=Dpmjet(b1, b2);
	71	else if (t.BeginsWith("phojet")) g=Dpmjet(b1, b2);
	72	else if (t.BeginsWith("hydjet")) g=Hydjet(b1, b2);
	73	else if (t.BeginsWith("lego")) g=Lego(rt);
f3c8e1f4	74	if (!g && !fIsLego)
15cb4a2b	75	Fatal("", "Invalid run type \"%s\" specified", t.Data());
66ff3038	76	return g;
	77	}
	78	/**
	79	* Make our decayer
	80	*
	81	* @return Newly allocated decayer or null
	82	*/
	83	TVirtualMCDecayer* CreateDecayer(const TString& runType)
	84	{
	85	if (runType.BeginsWith("hydjet")) return 0;
	86
	87	LoadPythia();
	88	TVirtualMCDecayer* decayer = new AliDecayerPythia();
	89	if (runType.EqualTo("hijing2000hf") && hftype < 2)
	90	decayer->SetForceDecay(kHadronicD);
	91	else
	92	decayer->SetForceDecay(kAll);
	93	decayer->Init();
	94	return decayer;
	95	}
	96
	97	// === PYTHIA ========================================================
	98	// Normal
	99	/**
	100	* Greate a pythia6 event generator
	101	*
	102	* @param tune Possible tune
	103	*
	104	* @return newly allocated generator or null
	105	*/
	106	AliGenerator* Pythia(const TString & tune)
	107	{
	108	// Int_t kCTEQ6l = 8;
	109	if (!grp->IsPP()) Fatal("Setup", "Pythia6 only works for pp");
	110
	111	TString t(tune);
	112	t.ToUpper();
	113	t.ReplaceAll("PYTHIA6", "");
	114	t.ReplaceAll("PYTHIA", "");
	115	Info("Setup", "Making Pythia6 event generator (tune: %s)", t.Data());
	116
	117	LoadPythia();
	118	AliGenPythia* pythia = new AliGenPythia(-1);
	119	pythia->SetMomentumRange(0, 999999.);
	120	pythia->SetThetaRange(0., 180.);
	121	pythia->SetYRange(-12.,12.);
	122	pythia->SetPtRange(0,1000.);
	123	pythia->SetProcess(kPyMb);
	124	pythia->SetEnergyCMS(grp->energy);
	125
	126	if (t == "D6T") {
	127	// Tune
	128	// 109 D6T : Rick Field's CDF Tune D6T
	129	// (NB: needs CTEQ6L pdfs externally)
	130	pythia->SetTune(109); // F I X
	131	pythia->SetStrucFunc(kCTEQ6l);
	132	}
	133	else if (t == "PERUGIA0") {
	134	// Tune
	135	// 320 Perugia 0
	136	pythia->SetTune(320);
	137	pythia->UseNewMultipleInteractionsScenario();
	138	}
	139	else if (t == "ATLAS") {
140	// Tune
141	// C 306 ATLAS-CSC: Arthur Moraes' (new) ATLAS tune
142	// (needs CTEQ6L externally)
143	pythia->SetTune(306);
144	pythia->SetStrucFunc(kCTEQ6l);
145	}
146	else if (t == "JETS") {
147	pythia->SetProcess(kPyJets);
148	pythia->SetStrucFunc(kCTEQ6l);
149	pythia->SetJetEtaRange(-1.5, 1.5);
150	pythia->SetJetEtRange(50., 800.);
151	pythia->SetPtHard(45., 1000.);
152	pythia->SetPycellParameters(2.2, 300, 432, 0., 4., 5., 0.7);
153	}
154	else if (t == "ATLAS_FLAT") {
155	// set high multiplicity trigger
156	// this weight achieves a flat multiplicity distribution
157	TH1 *weight = new TH1D("weight","weight",201,-0.5,200.5);
158	weight->SetBinContent(1,5.49443);
159	weight->SetBinContent(2,8.770816);
160	weight->SetBinContent(6,0.4568624);
161	weight->SetBinContent(7,0.2919915);
162	weight->SetBinContent(8,0.6674189);
163	weight->SetBinContent(9,0.364737);
164	weight->SetBinContent(10,0.8818444);
165	weight->SetBinContent(11,0.531885);
166	weight->SetBinContent(12,1.035197);
167	weight->SetBinContent(13,0.9394057);
168	weight->SetBinContent(14,0.9643193);
169	weight->SetBinContent(15,0.94543);
170	weight->SetBinContent(16,0.9426507);
171	weight->SetBinContent(17,0.9423649);
172	weight->SetBinContent(18,0.789456);
173	weight->SetBinContent(19,1.149026);
174	weight->SetBinContent(20,1.100491);
175	weight->SetBinContent(21,0.6350525);
176	weight->SetBinContent(22,1.351941);
177	weight->SetBinContent(23,0.03233504);
178	weight->SetBinContent(24,0.9574557);
179	weight->SetBinContent(25,0.868133);
180	weight->SetBinContent(26,1.030998);
181	weight->SetBinContent(27,1.08897);
182	weight->SetBinContent(28,1.251382);
183	weight->SetBinContent(29,0.1391099);
184	weight->SetBinContent(30,1.192876);
185	weight->SetBinContent(31,0.448944);
186	for (Int_t i = 32; i <= 201; i++) weight->SetBinContent(i,1);
187	weight->SetEntries(526);
188
189	Int_t limit = weight->GetRandom();
190	pythia->SetTriggerChargedMultiplicity(limit, 1.4);
191	}
192	return pythia;
193	}
194	/**
195	* Create a Pythia6 generator for high-flavor physics
196	*
197	* @param type Which kind
198	* @param harder If true, make harder processes
199	*
200	* @return Newly allocated generator or null
201	*/
202	AliGenerator* PythiaHF(Int_t type, Bool_t harder=0)
203	{
204	LoadPythia();
205	if (type == 6) return Pythia("jets");
206	if (type == 4) {
207	AliGenParam *jpsi = AliGenParam(1, AliGenMUONlib::kJpsi,
208	(harder?"CDF pp 8.8":"CDF pp 7"),"Jpsi");
209	jpsi->SetPtRange(0.,999.);
210	jpsi->SetYRange(-1.0, 1.0);
211	jpsi->SetPhiRange(0.,360.);
212	jpsi->SetForceDecay(kDiElectron);
213	return jpsi;
214	}
215	AliGenPythia* pythia = static_cast<AliGenPythia*>(Pythia("PERUGIA0"));
216	switch (type) {
217	case 0: // chadr
218	pythia->SetProcess(kPyCharmppMNRwmi);
219	pythia->SetForceDecay(kHadronicD);
220	break;
221	case 1: // bchadr
222	pythia->SetProcess(kPyBeautyppMNRwmi);
223	pythia->SetForceDecay(kHadronicD);
224	break;
225	case 2: // cele
226	pythia->SetProcess(kPyCharmppMNRwmi);
227	pythia->SetCutOnChild(1);
228	pythia->SetPdgCodeParticleforAcceptanceCut(11);
229	pythia->SetChildYRange(-1.2,1.2);
230	pythia->SetChildPtRange(0,10000.);
231	break;
232	case 3: // bele
233	pythia->SetProcess(kPyBeautyppMNRwmi);
234	pythia->SetCutOnChild(1);
235	pythia->SetPdgCodeParticleforAcceptanceCut(11);
236	pythia->SetChildYRange(-1.2,1.2);
237	pythia->SetChildPtRange(0,10000.);
238	break;
239	case 5:
240	pythia->SetProcess(kPyBeautyppMNRwmi);
241	pythia->SetCutOnChild(1);
242	pythia->SetPdgCodeParticleforAcceptanceCut(443);
243	pythia->SetChildYRange(-2,2);
244	pythia->SetChildPtRange(0,10000.);
245	}
246	return pythia;
247	}
248	/**
249	* Make a Min-Bias AA, pA, or Ap Hijing generator
250	*
251	* @param minB Least impact parameter
252	* @param maxB Largest impact parameter
253	* @param quench If true, enable quenching
254	* @param slowN If true, make a cocktail with slow neutrons
255	* @param ptHard Hard pT cut-off
256	*
257	* @return Generator
258	*/
259	AliGenerator* Hijing(Float_t minB,
260	Float_t maxB,
261	Bool_t slowN=false,
262	Bool_t quench=1,
263	Float_t ptHard=0)
264	{
265	LoadHijing();
266	AliGenHijing *gener = new AliGenHijing(-1);
267	// centre of mass energy
268	gener->SetEnergyCMS(grp->energy);
269	gener->SetImpactParameterRange(minB, maxB);
270	// reference frame
271	gener->SetReferenceFrame("CMS");
272	// projectil
f3c8e1f4	273	gener->SetTarget (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
f3c8e1f4	274	gener->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
66ff3038	275	// tell hijing to keep the full parent child chain
	276	gener->KeepFullEvent();
	277	// enable jet quenching
	278	gener->SetJetQuenching(quench);
	279	// enable shadowing
	280	gener->SetShadowing(slowN);
	281	// Don't track spectators
	282	gener->SetSpectators(!slowN);
	283	//
	284	if (ptHard > 0) hi->SetPtHardMin(ptHard);
	285
	286	// kinematic selection
	287	gener->SetSelectAll(0);
	288	// Boosted CMS
	289	gener->SetBoostLHC(grp->IsPA() \|\| grp->IsAP());
	290	// No need for cocktail
	291	if (!slowN \|\| !grp->IsPA() \|\| !grp->IsAP()) return gener;
	292
	293	AliGenCocktail* cocktail = new AliGenCocktail();
f3c8e1f4	294	cocktail->SetTarget (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
f3c8e1f4	295	cocktail->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
66ff3038	296	cocktail->SetEnergyCMS(grp->energy);
	297
	298	AliGenSlowNucleons* gray = new AliGenSlowNucleons(1);
	299	AliCollisionGeometry* coll = gener->CollisionGeometry();
	300	AliSlowNucleonModelExp* model = new AliSlowNucleonModelExp();
	301	// Not yet in the release...
	302	// model->SetSaturation(kTRUE);
	303	gray->SetSlowNucleonModel(model);
f3c8e1f4	304	gray->SetTarget(grp->beam1.a, grp->beam1.z);
66ff3038	305	gray->SetThetaDist(1);
	306	gray->SetProtonDirection(grp->beam1.IsP() ? 1 : 2);
	307	// gray->SetDebug(1);
	308	gray->SetNominalCmsEnergy(2*grp->beamEnergy);
	309	gray->NeedsCollisionGeometry();
	310	gray->SetCollisionGeometry(coll);
	311
	312	cocktail->AddGenerator(gener, "Hijing pPb", 1);
	313	cocktail->AddGenerator(gray, "Gray Particles", 1);
	314
	315	return cocktail;
	316	}
	317	/**
	318	* Make a DPMJet generator for pp, AA, pA, or Ap.
	319	*
	320	* @param fragments If true, make fragments
	321	*
	322	* @return Generator
	323	*/
	324	AliGenerator* Dpmjet(Float_t minB, Float_t maxB,
	325	Bool_t fragments=0)
	326	{
	327	LoadDpmjet();
	328	AliGenDPMjet* dpmjet = new AliGenDPMjet(-1);
	329	dpmjet->SetEnergyCMS(grp->energy);
15cb4a2b	330	dpmjet->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
15cb4a2b	331	dpmjet->SetTarget (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
66ff3038	332	dpmjet->SetImpactParameterRange(minB, maxB);
15cb4a2b	333	dpmjet->SetProjectileBeamEnergy(grp->beam2.z*grp->beamEnergy/grp->beam2.a);
66ff3038	334	if (grp->IsAA()) {
	335	dpmjet->SetPi0Decay(0);
	336	}
	337	else if (grp->IsPA() \|\| grp->IsAP()) {
f3c8e1f4	338	// dpmjet->SetTriggerParticle(3312, 1.2, 2.0);
701a83f4	339	dpmjet->SetFragmentProd(false/fragments/); // Alwas disabled
66ff3038	340	}
	341	else if (grp->IsPP()) { // PhoJet
	342	dpmjet->SetMomentumRange(0, 999999.);
	343	dpmjet->SetThetaRange(0., 180.);
	344	dpmjet->SetYRange(-12.,12.);
	345	dpmjet->SetPtRange(0,1000.);
	346	}
	347	return dpmjet;
	348	}
	349	/**
	350	* Make an AMPT generator for AA collisions
	351	*
	352	* @return Generator
	353	*/
	354	AliGenerator* Ampt(Float_t minB, Float_t maxB)
	355	{
	356	LoadAmpt();
	357	AliGenAmpt *genHi = new AliGenAmpt(-1);
	358	genHi->SetEnergyCMS(grp->energy);
	359	genHi->SetReferenceFrame("CMS");
f3c8e1f4	360	genHi->SetTarget (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
f3c8e1f4	361	genHi->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
66ff3038	362	genHi->SetPtHardMin (2);
	363	genHi->SetImpactParameterRange(minB,maxB);
	364	// disable jet quenching
	365	genHi->SetJetQuenching(0);
	366	// enable shadowing
	367	genHi->SetShadowing(1);
	368	// neutral pion and heavy particle decays switched off
	369	genHi->SetDecaysOff(1);
	370	genHi->SetSpectators(0); // track spectators
	371	genHi->KeepFullEvent();
	372	genHi->SetSelectAll(0);
	373	return genHi;
	374	}
	375	/**
	376	* Make an HydJet generator for A-A
	377	*
	378	* @return Generator
	379	*/
	380	AliGenerator* Hydjet(Float_t minB, Float_t maxB)
	381	{
	382	LoadHydjet();
	383	AliGenUHKM *genHi = new AliGenUHKM(-1);
	384	genHi->SetAllParametersLHC();
f3c8e1f4	385	genHi->SetTarget (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
f3c8e1f4	386	genHi->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
66ff3038	387	genHi->SetEcms(grp->energy);
	388	genHi->SetEnergyCMS(grp->energy);
	389	genHi->SetBmin(minB);
	390	genHi->SetBmax(maxB);
	391	genHi->SetPyquenPtmin(9);
	392	return genHi;
	393	}
15cb4a2b	394	// === Lego ========================================================
	395	/**
	396	* Greate a lego event generator
	397	*
	398	* @param tune Possible tune
	399	*
	400	* @return newly allocated generator or null
	401	*/
	402	AliGenerator* Lego(const TString & variant)
	403	{
	404	fIsLego = true;
	405	return 0;
	406	#if 0
	407	TString v(variant);
	408	v.ToUpper();
	409	v.ReplaceAll("LEGO", "");
	410	Info("Setup", "Making Lego event generator (variant: %s)", v.Data());
	411
	412	AliLegoGenerator* ret = 0;
	413	// XYZ varies origin of the particles in two dimensions:
	414	// X: o=(0,t1,t2), p=(1,0,0)
	415	// Y: o=(t1,0,t2), p=(0,1,0)
	416	// Z: o=(t1,t2,0), p=(0,0,1)
	417	// PhiZ varies the momentum in two dimensions
	418	// o=(0,0,t1) p=(cos(t2),sin(t2),0)
	419	// Eta varies momentum in two dimensions
	420	// phi=t1
	421	// theta=2*atan(exp(-t2))
	422	// o=(0,0,0) p=(cos(phi)sin(theta),sin(phi)cos(theta),cos(theta))
	423	// Base varies in two dimensions
	424	// phi=t1
	425	// theta=t2
	426	// o=(0,0,0) p=(cos(phi)sin(theta),sin(phi)cos(theta),cos(theta))
	427	if (v.BeginsWith("X") \|\| v.BeginsWith) {
	428	const char* c[] = { v(0), '\0' };
	429	ret = new AliLegoGeneratorXYZ(c);
	430	ret->SetCoor1Range(10,-2,2); // Y, X
	431	ret->SetCoor2Range(10,-10,10); // Z
	432	}
	433	else if (v.BeginsWith("Z")) {
	434	ret = new AliLegoGeneratorXYZ("Z");
	435	ret->SetCoor1Range(10,-2,2); // X
	436	ret->SetCoor2Range(10,-2,2); // Y
	437	}
	438	else if (v.BeginsWith("PHIZ")) {
	439	ret = new AliLegoGeneratorPhiZ();
	440	ret->SetCoor1Range(10,-10,10); // Z
	441	ret->SetCoor2Range(360,0,360); // phi
	442	}
	443	else if (v.BeginsWith("ETA")) {
	444	ret = new AliLegoGeneratorEta();
	445	ret->SetCoor1Range(360,0,360); // phi
	446	Double_t aEta = 7;
	447	Double_t dEta = (6--4)/200;
	448	ret->SetCoor2Range(2*aEta/dEta,-aEta,+aEta); // Eta
	449	}
	450	else {
	451	ret = new AliLegoGenerator();
	452	ret->SetCoor1Range(180,0,180); // theta
	453	ret->SetCoor1Range(360,0,360); // phi
	454	}
	455	return ret;
	456	#endif
	457	}
458
66ff3038	459	/**
	460	* Make a heavy flavour cocktail
	461	*
	462	* @param base Underlying event.
	463	*
	464	* @return Generator
	465	*/
	466	AliGeneator* HFCocktail(const TString& base, Float_t minB, Float_t maxB)
	467	{
	468
	469	AliGenCocktail *cocktail = new AliGenCocktail();
f3c8e1f4	470	cocktail->SetTarget (grp->beam1.Name(), grp->beam1.a, grp->beam1.z);
f3c8e1f4	471	cocktail->SetProjectile(grp->beam2.Name(), grp->beam2.a, grp->beam2.z);
66ff3038	472	cocktail->SetEnergyCMS(grp->energy);
	473
	474	// Add underlying event
	475	if (base.BeginsWith("ampt", TString::kIgnoreCase)) {
	476	hi = Ampt(minB, maxB);
	477	cocktail->AddGenerator(hi,"ampt",1);
	478	}
	479	else {
	480	hi = Hijing(minB, maxB, grp->IsPA() \|\| grp->IsAP(), false, 2.3);
	481	cocktail->AddGenerator(hi,"hijing",1);
	482
	483	}
	484
	485	// --- Default formula -------------------------------------------
	486	TForumla* one = new TFormula("one", "1.");
	487
	488	// --- Pythia ----------------------------------------------------
	489	AliGenerator* pythia = PythiaHF(hftype);
	490	switch (hftype) {
	491	case 6:
	492	cocktail->AddGenerator(pythia, "pythiaJets", 1, one);
	493	break;
	494	defualt:
	495	cocktail
	496	->AddGenerator(pythia, "pythiaHF", 1,
	497	new TFormula("Signals",
	498	"20.(x<5.)+80./3.(1.-x/20.)*(x>5.)"));
	499	break;
	500	}
	501	// --- Dummy -----------------------------------------------------
	502	AliGenParam* param = 0;
	503
	504	// --- Phos stuff ------------------------------------------------
	505	AliGenPHOSlib *plib = new AliGenPHOSlib();
	506	Double_t lower[] = { 0, 3, 6, 9, 12, -1 };
	507	Double_t *pLow = lower;
	508	for (Int_t i = 0; i < 5; i++) {
	509	param = new AliGenParam(5, plib, AliGenPHOSlib::kPi0);
	510	param->SetPhiRange(0., 360.) ;
	511	param->SetYRange(-1.2, 1.2) ;
	512	param->SetPtRange(lower[i], 30.) ;
	513	cocktail->AddGenerator(param,Form("Pi0HagPt%d", i), 1., one);
	514	}
	515
	516	// --- Jpsi->mu+ mu-----------------------------------------------
	517	param = new AliGenParam(1, AliGenMUONlib::kJpsi, "CDF pp 3.94", "Jpsi");
	518	param->SetPtRange(0.,999.);
	519	param->SetYRange(-4.2, -2.3);
	520	param->SetPhiRange(0.,360.);
	521	param->SetForceDecay(kDiMuon);
	522	param->SetCutOnChild(1);
	523	param->SetChildPhiRange(0.,360.);
	524	param->SetChildThetaRange(168.5,178.5);
	525	cocktail->AddGenerator(param, "Jpsi2M", 1, one);
	526
	527	// --- Chi_c -> J/Psi + gamma, J/Psi -> e+e- ---------------------
	528	Float_t thmin = (180./TMath::Pi())2.atan(exp(-1.2));
	529	Float_t thmax = (180./TMath::Pi())2.atan(exp( 1.2));
	530	param = new AliGenParam(1, AliGenMUONlib::kChic,"default","Chic");
	531	param->SetMomentumRange(0, 999.); // Wide cut on the momentum
	532	param->SetPtRange(0, 100.); // Wide cut on Pt
	533	param->SetYRange(-2.5, 2.5);
	534	param->SetCutOnChild(1); // Enable cuts on decay products
	535	param->SetChildPhiRange(0., 360.);
536	// In the acceptance of the Central Barrel
537	param->SetChildThetaRange(thmin, thmax);
538	// Chi_c -> J/Psi + gamma, J/Psi -> e+e-
539	param->SetForceDecay(kChiToJpsiGammaToElectronElectron);
540	cocktail->AddGenerator(param, "Chi_c", 1, one);
541
542	// --- Dummy -----------------------------------------------------
543	AliGenBox* box = 0;
544
545	// --- Some particles --------------------------------------------
546	Double_t boxR = gRandom->Integer(3)+1;
547	Int_t boxP[] = { 310, 3122, 3312, 3322,
548	(boxR==1 ?3334: boxR==2 ?-3334: -3312) };
549	Int_t boxN[] = { 1, 1, 3, 3, 1 }
550	const char* boxT[] = { "K0s", "Lambda", "Xi-", "Xi0",
551	(boxR==1 ? "Omega-": boxR==2 ? "Omega+": "Xi+") };
552	for (Int_t i = 0; i < 5; i++) {
553	box = new AliGenBox(boxN[i]);
554	box->SetPart(boxP[i]);
555	box->SetPtRange(0,13);
556	box->SetThetaRange(45, 135);
557	cocktail->AddGenerator(box, boxT[i], 1, one);
558	}
559
560	// --- High pT charged particle ----------------------------------
561	TFormula* hptF = new TFormula("High Pt",
562	"5.(x<5.)+20./3.(1.-x/20.)*(x > 5.)");
563	Int_t hptP[] = { 211, 321, 2212 };
564	const char* hptT[] = { "pi", "K", "p" };
565	for (Int_t i = 0; i < 3; i++) {
566	for (Int_t j = -1; j <= 1; j++) {
567	box->SetPart(j*hptP[i]);
568	box->SetPtRange(2., 50.);
569	box->SetThetaRange(thmin, thmax);
570	cocktail->AddGenerator(box, Form("%s%c",hptT[i],(j<0,'-','+')),1,hptF);
571	}
572	}
573	return cocktail;
574	}
575
576	};
577
578	void
579	EGConfig()
580	{
581	::Info("EGConfig", "Creating EG factory");
582	egCfg = new EGCfg;
583	}
584
585	//
586	// EOF
587	//
588
589