Classes imported from EVGEN
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 15 Mar 2003 14:49:39 +0000 (14:49 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 15 Mar 2003 14:49:39 +0000 (14:49 +0000)
22 files changed:
TEPEMGEN/AliGenEpEmv1.cxx [new file with mode: 0644]
TEPEMGEN/AliGenEpEmv1.h [new file with mode: 0644]
TEPEMGEN/TEPEMGENLinkDef.h
TEPEMGEN/libTEPEMGEN.pkg
THbtp/AliGenHBTprocessor.cxx [new file with mode: 0644]
THbtp/AliGenHBTprocessor.h [new file with mode: 0644]
THbtp/THBTprocessorLinkDef.h
THbtp/libTHbtp.pkg
THijing/AliGenHijing.cxx [new file with mode: 0644]
THijing/AliGenHijing.h [new file with mode: 0644]
THijing/AliGenHijingEventHeader.cxx [new file with mode: 0644]
THijing/AliGenHijingEventHeader.h [new file with mode: 0644]
THijing/THijingLinkDef.h
THijing/libTHijing.pkg
TMEVSIM/AliGenMevSim.cxx [new file with mode: 0644]
TMEVSIM/AliGenMevSim.h [new file with mode: 0644]
TMEVSIM/AliMevSimConfig.cxx [new file with mode: 0644]
TMEVSIM/AliMevSimConfig.h [new file with mode: 0644]
TMEVSIM/AliMevSimParticle.cxx [new file with mode: 0644]
TMEVSIM/AliMevSimParticle.h [new file with mode: 0644]
TMEVSIM/TMevSimLinkDef.h
TMEVSIM/libTMEVSIM.pkg

diff --git a/TEPEMGEN/AliGenEpEmv1.cxx b/TEPEMGEN/AliGenEpEmv1.cxx
new file mode 100644 (file)
index 0000000..f4a7f6a
--- /dev/null
@@ -0,0 +1,196 @@
+/**************************************************************************
+ * Copyright(c) 1998-2002, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ *                                                                        *
+ *                                                                        *
+ * Copyright(c) 1997, 1998, 2002, Adrian Alscher and Kai Hencken          *
+ * See $ALICE_ROOT/EpEmGen/diffcross.f for full Copyright notice          *
+ *                                                                        *
+ *                                                                        *
+ * Copyright(c) 2002 Kai Hencken, Yuri Kharlov, Serguei Sadovsky          *
+ * See $ALICE_ROOT/EpEmGen/epemgen.f for full Copyright notice            *
+ *                                                                        *
+ **************************************************************************/
+
+/*
+$Log$
+Revision 1.1  2002/11/06 10:26:45  hristov
+Event generator for e+e- pair production (Yu.Kharlov)
+
+*/
+
+// Event generator of single e+e- pair production in ultraperipheral PbPb collisions
+// at 5.5 TeV/nucleon.
+// The generator is based on 5-dimentional differential cross section of the process.
+//%
+// References:
+// [1] "Multiple electromagnetic electron positron pair production in
+//      relativistic heavy ion collisions".
+//      Adrian Alscher, Kai Hencken, Dirk Trautmann, and Gerhard Baur,
+//      Phys. Rev. A55 (1997) 396.
+// [2] K.Hencken, Yu.Kharlov, S.Sadovsky, Internal ALICE Note 2002-27.
+//%
+// Usage:
+// Initialization:
+//    AliGenEpEmv1 *gener = new AliGenEpEmv1();
+//    gener->SetXXXRange(); // Set kinematics range
+//    gener->Init();
+// Event generation:
+//    gener->Generate(); // Produce one e+e- pair with the event weight assigned 
+//                       // to each track. The sum of event weights, divided by 
+//                       // the total number of generated events, gives the 
+//                       // integral cross section of the process of e+e- pair 
+//                       // production in the above mentioned kinematics range.
+//                       // Sum of the selected event weights, divided by the total 
+//                       // number of generated events, gives the integral cross 
+//                       // section corresponded to the set of selected events
+//%
+// The generator consists of several modules:
+// 1) $ALICE_ROOT/EpEmGen/diffcross.f:
+//    Exact calculation of the total differential e+ e- -pair production
+//    in Relativistic Heavy Ion Collisions for a point particle in an
+//    external field approach. See full comments in the mentioned file.
+// 2) $ALICE_ROOT/EpEmGen/epemgen.f:
+//    Generator of e+e- pairs produced in PbPb collisions at LHC
+//    it generates events according to the parametrization of the
+//    differential cross section. Produces events have weights calculated
+//    by the exact differential cross section calculation (diffcross.f).
+//    See full comments in the mentioned file.
+// 3) Class TEpEmGen:
+//    Interface from the fortran event generator to ALIROOT
+// 4) Class AliGenEpEmv1:
+//    The event generator to call within ALIROOT
+//%
+// Author of this module: Yuri.Kharlov@cern.ch
+// 9 October 2002
+
+#include "AliGenEpEmv1.h"
+#include <TParticle.h>
+#include <TParticlePDG.h>
+#include <TDatabasePDG.h>
+#include <TEpEmGen.h>
+
+ClassImp(AliGenEpEmv1)
+
+//------------------------------------------------------------
+
+AliGenEpEmv1::AliGenEpEmv1()
+{
+  // Default constructor
+  // Avoid zero pt
+  if (fPtMin == 0) fPtMin = 1.E-04;
+}
+
+//____________________________________________________________
+AliGenEpEmv1::AliGenEpEmv1(const AliGenEpEmv1 & gen)
+{
+  // copy constructor
+  gen.Copy(*this);
+}
+
+//____________________________________________________________
+AliGenEpEmv1::~AliGenEpEmv1()
+{
+  // Destructor
+}
+
+//____________________________________________________________
+void AliGenEpEmv1::Init()
+{
+  // Initialisation:
+  // 1) define a generator
+  // 2) initialize the generator of e+e- pair production
+
+  fMass = TDatabasePDG::Instance()->GetParticle(11)->Mass();
+
+  SetMC(new TEpEmGen());
+  fEpEmGen = (TEpEmGen*) fgMCEvGen;
+  fEpEmGen ->Initialize(fYMin,fYMax,fPtMin,fPtMax);
+  fEvent = 0;
+}
+
+//____________________________________________________________
+void AliGenEpEmv1::Generate()
+{
+  //
+  // Generate one e+e- pair
+  // Gaussian smearing on the vertex is done if selected. 
+  //%
+  // Each produced e+e- pair is defined by the following variables:
+  // rapidities of e-, e+ (yElectron,yPositron)
+  // log10(pt in MeV/c) of e-, e+ (xElectron,xPositron)
+  // azymuth angles between e- and e+ (phi12)
+  //%
+  // On output an event weight is given (weight) which is assigned to each track.
+  // The sum of event weights, divided by the total number of generated events, 
+  // gives the integral cross section of the e+e- pair production in the   
+  // selected kinematics range.          
+  //
+
+  Float_t polar[3]= {0,0,0};
+  Float_t origin[3];
+  Float_t p[3];
+
+  Double_t ptElectron,ptPositron, phiElectron,phiPositron, mt;
+  Double_t phi12=0,xElectron=0,xPositron=0,yElectron=0,yPositron=0,weight=0;
+  Int_t   j, nt, id;
+  Float_t random[6];
+
+  fEpEmGen->GenerateEvent(fYMin,fYMax,fPtMin,fPtMax,
+          yElectron,yPositron,xElectron,xPositron,phi12,weight);
+  if (fDebug == 1)
+    printf("AliGenEpEmv1::Generate(): y=(%f,%f), x=(%f,%f), phi=%f\n",
+          yElectron,yPositron,xElectron,xPositron,phi12);
+
+  for (j=0;j<3;j++) origin[j]=fOrigin[j];
+  if(fVertexSmear==kPerEvent) {
+    Rndm(random,6);
+    for (j=0;j<3;j++) {
+      origin[j]+=fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
+       TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
+    }
+  }
+
+  Rndm(random,1);
+  ptElectron  = TMath::Power(10,xElectron) * 1.e-03;;
+  ptPositron  = TMath::Power(10,xPositron) * 1.e-03;;
+  phiElectron = fPhiMin + random[0] * (fPhiMax-fPhiMin);
+  phiPositron = phiElectron + phi12;
+
+  // Produce electron
+  mt = TMath::Sqrt(ptElectron*ptElectron + fMass*fMass);
+  p[0] = ptElectron*TMath::Cos(phiElectron);
+  p[1] = ptElectron*TMath::Sin(phiElectron);
+  p[2] = mt*TMath::SinH(yElectron);
+  id =  11;
+  if (fDebug == 2)
+    printf("id=%+3d, p = (%+11.4e,%+11.4e,%+11.4e) GeV\n",id,p[0],p[1],p[2]);
+  SetTrack(fTrackIt,-1, id,p,origin,polar,0,kPPrimary,nt,weight);
+
+  // Produce positron
+  mt = TMath::Sqrt(ptPositron*ptPositron + fMass*fMass);
+  p[0] = ptPositron*TMath::Cos(phiPositron);
+  p[1] = ptPositron*TMath::Sin(phiPositron);
+  p[2] = mt*TMath::SinH(yPositron);
+  id = -11;
+  if (fDebug == 2)
+    printf("id=%+3d, p = (%+11.4e,%+11.4e,%+11.4e) GeV\n",id,p[0],p[1],p[2]);
+  SetTrack(fTrackIt,-1, id,p,origin,polar,0,kPPrimary,nt,weight);
+  
+  fEvent++;
+  if (fEvent%1000 == 0) {
+    printf("=====> AliGenEpEmv1::Generate(): \n   Event %d, sigma=%f +- %f kb\n",
+          fEvent,fEpEmGen->GetXsection(),fEpEmGen->GetDsection());
+  }
+}
+
diff --git a/TEPEMGEN/AliGenEpEmv1.h b/TEPEMGEN/AliGenEpEmv1.h
new file mode 100644 (file)
index 0000000..59e9c28
--- /dev/null
@@ -0,0 +1,36 @@
+#ifndef ALIGENEPEMV1_H
+#define ALIGENEPEMV1_H
+/* Copyright(c) 1998-2002, ALICE Experiment at CERN, All rights reserved. *
+ * Copyright(c) 1997, 1998, 2002, Adrian Alscher and Kai Hencken          *
+ * Copyright(c) 2002 Kai Hencken, Yuri Kharlov, Serguei Sadovsky          *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+// Event generator of single e+e- pair production in ultraperipheral PbPb collisions
+// at 5.5 TeV/nucleon.
+// Author: Yuri.Kharlov@cern.ch
+// 9 October 2002
+
+#include "AliGenMC.h"
+class TEpEmGen;
+
+//-------------------------------------------------------------
+class AliGenEpEmv1 : public AliGenMC
+{
+public:
+  AliGenEpEmv1();
+  
+  virtual ~AliGenEpEmv1();
+  void Generate();
+  void Init();
+  void SetDebug(Int_t debug) {fDebug=debug;}
+  
+  AliGenEpEmv1(const AliGenEpEmv1 & gen);
+protected:
+  Float_t    fMass;    // electron mass
+  TEpEmGen * fEpEmGen; // e+e- generator
+  Int_t      fDebug;   // debug level
+  Int_t      fEvent;   // internal event number
+ClassDef(AliGenEpEmv1,1) // Generator of single e+e- pair production in PbPb ultra-peripheral collisions
+};
+#endif
index 1498d48..0de7193 100644 (file)
@@ -5,5 +5,6 @@
 #pragma link off all functions;
  
 #pragma link C++ class TEpEmGen+;
+#pragma link C++ class AliGenEpEmv1+;
 
 #endif
index 7dc7089..f0f0f3b 100644 (file)
@@ -1,6 +1,6 @@
-SRCS= TEpEmGen.cxx
+SRCS= TEpEmGen.cxx AliGenEpEmv1.cxx
 
-HDRS= TEpEmGen.h
+HDRS= TEpEmGen.h AliGenEpEmv1.h
 
 DHDR:=TEPEMGENLinkDef.h
 
diff --git a/THbtp/AliGenHBTprocessor.cxx b/THbtp/AliGenHBTprocessor.cxx
new file mode 100644 (file)
index 0000000..b192097
--- /dev/null
@@ -0,0 +1,1151 @@
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+// Implementation of the interface for THBTprocessor
+// Author: Piotr Krzysztof Skowronski <Piotr.Skowronski@cern.ch>
+//////////////////////////////////////////////////////////////////////////////////
+//Wrapper class for "hbt processor" after burner
+//The origibal code is written in fortran by Lanny Ray
+//and is put in the directory $ALICE_ROOT/HBTprocessor
+//Detailed description is on the top of the file hbt_event_processor.f
+//
+//This class can be used ONLY with AliGenCocktailAfterBurner wrapper generator
+//i.e. (in Config.C)
+// ....
+// AliGenCocktailAfterBurner = gener = new AliGenCocktailAfterBurner();
+// gener->SetPhiRange(0, 360); //Set global parameters
+// gener->SetThetaRange(35., 145.); 
+// AliGenHIJINGpara *hijing = new AliGenHIJINGpara(10000); //Add some generator
+// hijing->SetMomentumRange(0, 999);   
+// gener->AddGenerator(hijing,"HIJING PARAMETRIZATION",1); 
+//
+// AliGenHBTprocessor *hbtp = new AliGenHBTprocessor(); //create object
+// hbtp->SetRefControl(2); //set parameters
+// hbtp->SetSwitch1D(1);
+// hbtp->SetR0(6);//fm - radius
+// hbtp->SetLambda(0.7);//chaocity parameter
+// gener->AddAfterBurner(hbtp,"HBT PROCESSOR",1); //add to the main generator
+// 
+// gener->Init();
+//
+//CAUTIONS: 
+//         A)  HBT PROCESSOR NEEDS MORE THAN ONE EVENT TO WORK
+//             AS MORE AS IT BETTER WORKS
+//         B)  IT IS ABLE TO "ADD" CORRELATIONS ONLY UP TO TWO PARTICLE TYPES AT ONES
+//////////////////////////////////////////////////////////////////////////////////
+
+// 11.11.2001 Piotr Skowronski
+// Setting seed (date) in RNG in the constructor
+
+// 09.10.2001 Piotr Skowronski
+// 
+// Theta - Eta cohernecy facilities added:
+//    AliGenerator::SetThetaRange method overriden
+//    Static methods added
+//    EtaToTheta
+//    ThetaToEta 
+//    DegreesToRadians
+//    RadiansToDegrees
+//
+// Class description comments put on proper place
+
+// 27.09.2001 Piotr Skowronski
+// removing of redefinition of defaults velues 
+// in method's implementation. 
+//  
+// 
+
+#include "AliGenHBTprocessor.h"
+#include "TROOT.h"
+#include <Riostream.h>
+#include <TFile.h>
+#include <TTree.h>
+#include "AliRun.h"
+#include "AliStack.h"
+#include "TParticle.h"
+#include "THBTprocessor.h"
+#include "AliGenCocktailAfterBurner.h"
+
+
+
+ClassImp(AliGenHBTprocessor)
+
+
+
+AliGenCocktailAfterBurner*  GetGenerator();
+/*******************************************************************/
+AliGenHBTprocessor::AliGenHBTprocessor(const AliGenHBTprocessor& in)
+{
+//copy contructor
+  // AliGenHBTprocessor::AliGenHBTprocessor();
+}
+
+AliGenHBTprocessor::AliGenHBTprocessor() : AliGenerator() 
+{
+  //
+  // Standard constructor
+  // Sets default veues of all parameters
+  fHbtPStatCodes = 0x0;
+  fHBTprocessor = 0x0;
+
+  SetName("AliGenHBTprocessor");
+  SetTitle("AliGenHBTprocessor");
+  
+  sRandom = fRandom;
+  fHBTprocessor = new THBTprocessor();
+
+  fNPDGCodes = 0;
+  DefineParticles();
+
+  SetTrackRejectionFactor();
+  SetRefControl();
+  SetPIDs();
+  SetNPIDtypes();
+  SetDeltap();
+  SetMaxIterations();
+  SetDelChi();
+  SetIRand();
+  SetLambda();
+  SetR1d() ;
+  SetRSide();
+  SetROut() ;
+  SetRLong() ;
+  SetRPerp();
+  SetRParallel();
+  SetR0();
+  SetQ0();
+  SetSwitch1D();
+  SetSwitch3D();
+  SetSwitchType();
+  SetSwitchCoherence();
+  SetSwitchCoulomb();
+  SetSwitchFermiBose();
+  //SetMomentumRange();
+  SetPtRange();
+  SetPxRange();
+  SetPyRange(); 
+  SetPzRange(); 
+  SetPhiRange(); 
+  SetEtaRange();  
+  SetNPtBins();  
+  SetNPhiBins();  
+  SetNEtaBins();
+  SetNPxBins(); 
+  SetNPyBins();
+  SetNPzBins(); 
+  SetNBins1DFineMesh();
+  SetBinSize1DFineMesh();
+  SetNBins1DCoarseMesh();
+  SetBinSize1DCoarseMesh();
+  SetNBins3DFineMesh();
+  SetBinSize3DFineMesh();
+  SetNBins3DCoarseMesh();
+  SetBinSize3DCoarseMesh();
+  SetNBins3DFineProjectMesh();
+}
+
+/*******************************************************************/
+
+
+/*******************************************************************/
+
+AliGenHBTprocessor::~AliGenHBTprocessor()
+{
+//destructor
+  CleanStatusCodes();
+  if (fHBTprocessor) delete fHBTprocessor; //delete generator
+  
+}
+
+/*******************************************************************/
+
+void AliGenHBTprocessor::InitStatusCodes()
+{
+ //creates and inits status codes array to zero
+  AliGenCocktailAfterBurner *cab = GetGenerator();
+
+  if(!cab) Fatal("InitStatusCodes()","Can not find AliGenCocktailAfterBurner generator");
+
+  Int_t nev = cab->GetNumberOfEvents();
+
+  fHbtPStatCodes = new Int_t* [nev];
+  for( Int_t i =0; i<nev;i++)
+  { 
+    Int_t nprim = cab->GetStack(i)->GetNprimary();
+    fHbtPStatCodes[i] = new Int_t[nprim];
+    for (Int_t k =0 ;k<nprim;k++)
+      fHbtPStatCodes[i][k] =0;
+    
+  }
+  
+}
+/*******************************************************************/
+
+void AliGenHBTprocessor::CleanStatusCodes()
+{
+ //Cleans up status codes
+  if (fHbtPStatCodes)
+  {
+    for (Int_t i =0; i<GetGenerator()->GetNumberOfEvents(); i++)
+      delete [] fHbtPStatCodes[i];  
+    delete fHbtPStatCodes;
+    fHbtPStatCodes = 0;
+  }
+
+}
+/*******************************************************************/
+
+void AliGenHBTprocessor::Init()
+  {  
+  //sets up parameters in generator
+   
+   THBTprocessor *thbtp = fHBTprocessor;
+   
+
+   thbtp->SetTrackRejectionFactor(fTrackRejectionFactor);
+   thbtp->SetRefControl(fReferenceControl);
+   
+   if ((fPid[0] == fPid[1]) || (fPid[0] == 0) || (fPid[1] == 0))
+    {
+       if (fPid[0] == 0)
+         thbtp->SetPIDs(IdFromPDG(fPid[1]) ,0);
+       else
+         thbtp->SetPIDs(IdFromPDG(fPid[0]) ,0);
+       thbtp->SetNPIDtypes(1);
+       
+       if (fSwitchType !=1)
+          Warning("AliGenHBTprocessor::Init","\nThere is only one particle type set,\n\
+                   and Switch_Type differnt then 1,\n which does not make sense.\n\
+                   Setting it to 1.\n");
+                   
+       SetSwitchType(1);
+    }
+   else
+    {
+       thbtp->SetPIDs(IdFromPDG(fPid[0]) ,IdFromPDG(fPid[1]));
+       SetNPIDtypes(2);
+       thbtp->SetSwitchType(fSwitchType); 
+    }
+   
+   thbtp->SetMaxIterations(fMaxit);
+   thbtp->SetDelChi(fDelchi);
+   thbtp->SetIRand(fIrand);
+   thbtp->SetLambda(fLambda);
+   thbtp->SetR1d(fR1d);
+   thbtp->SetRSide(fRside);
+   thbtp->SetROut(fRout);
+   thbtp->SetRLong(fRlong);
+   thbtp->SetRPerp(fRperp);
+   thbtp->SetRParallel(fRparallel);
+   thbtp->SetR0(fR0);
+   thbtp->SetQ0(fQ0);
+   thbtp->SetSwitch1D(fSwitch1d);
+   thbtp->SetSwitch3D(fSwitch3d);
+   thbtp->SetSwitchType(fSwitchType);
+   thbtp->SetSwitchCoherence(fSwitchCoherence);
+   thbtp->SetSwitchCoulomb(fSwitchCoulomb);
+   thbtp->SetSwitchFermiBose(fSwitchFermiBose);
+   thbtp->SetPtRange(fPtMin,fPtMax);
+   thbtp->SetPxRange(fPxMin,fPxMax);
+   thbtp->SetPyRange(fPyMin,fPyMax);
+   thbtp->SetPzRange(fPzMin,fPzMax);
+   thbtp->SetPhiRange(fPhiMin*180./TMath::Pi(),fPhiMax*180./TMath::Pi());
+   thbtp->SetEtaRange(fEtaMin,fEtaMax);
+   thbtp->SetNPtBins(fNPtBins);
+   thbtp->SetNPhiBins(fNPhiBins);
+   thbtp->SetNEtaBins(fNEtaBins);
+   thbtp->SetNPxBins(fNPxBins);
+   thbtp->SetNPyBins(fNPyBins);
+   thbtp->SetNPzBins(fNPzBins);
+   thbtp->SetNBins1DFineMesh(fN1dFine);
+   thbtp->SetBinSize1DFineMesh(fBinsize1dFine);
+   thbtp->SetNBins1DCoarseMesh(fN1dCoarse);
+   thbtp->SetBinSize1DCoarseMesh(fBinsize1dCoarse);
+   thbtp->SetNBins3DFineMesh(fN3dFine);
+   thbtp->SetBinSize3DFineMesh(fBinsize3dFine);
+   thbtp->SetNBins3DCoarseMesh(fN3dCoarse);
+   thbtp->SetBinSize3DCoarseMesh(fBinsize3dCoarse);
+   thbtp->SetNBins3DFineProjectMesh(fN3dFineProject);
+       
+ }
+/*******************************************************************/
+  
+void AliGenHBTprocessor::Generate()
+ {
+ //starts processig 
+   AliGenCocktailAfterBurner* cab = GetGenerator();
+   if (cab == 0x0)
+    {
+      Fatal("Generate()","AliGenHBTprocessor needs AliGenCocktailAfterBurner to be main generator");
+    }
+   if (cab->GetNumberOfEvents() <2)
+    {
+      Fatal("Generate()",
+            "HBT Processor needs more than 1 event to work properly,\
+             but there is only %d set", cab->GetNumberOfEvents());
+    }
+  
+  
+   fHBTprocessor->Initialize(); //reset all fields of common blocks 
+                                   //in case there are many HBT processors
+                                      //run one after one (i.e. in AliCocktailAfterBurner)
+                                   //between Init() called and Generate there might 
+   Init();                         //be different instance running - be sure that we have our settings
+   
+   InitStatusCodes(); //Init status codes
+   
+   fHBTprocessor->GenerateEvent(); //Generates event
+   
+   CleanStatusCodes(); //Clean Status codes - thet are not needed anymore
+ }
+/*******************************************************************/
+
+
+/*******************************************************************/
+void AliGenHBTprocessor::GetParticles(TClonesArray * particles) const
+ {
+ //practically dumm
+   if (!particles)
+    {
+      cout<<"Particles has to be initialized"<<endl;
+      return;
+    } 
+   fHBTprocessor->ImportParticles(particles);
+ }
+
+/*******************************************************************/
+
+Int_t AliGenHBTprocessor::GetHbtPStatusCode(Int_t part) const
+{
+//returns the status code of the given particle in the active event
+//see SetActiveEvent in the bottom of AliGenHBTprocessor.cxx
+//and in AliCocktailAfterBurner
+  Int_t activeEvent = GetGenerator()->GetActiveEventNumber();
+  return fHbtPStatCodes[activeEvent][part];
+}
+
+/*******************************************************************/
+void  AliGenHBTprocessor::SetHbtPStatusCode(Int_t hbtstatcode, Int_t part)
+{
+ //Sets the given status code to given particle number (part) in the active event
+  Int_t activeEvent = GetGenerator()->GetActiveEventNumber();
+  fHbtPStatCodes[activeEvent][part] = hbtstatcode;
+}
+
+/*******************************************************************/
+
+void AliGenHBTprocessor::SetTrackRejectionFactor(Float_t trf) //def 1.0
+  {
+   //Sets the Track Rejection Factor
+   //variates in range 0.0 <-> 1.0
+   //Describes the factor of particles rejected from the output.
+   //Used only in case of low muliplicity particles e.g. lambdas.
+   //Processor generates addisional particles and builds the 
+   //correletions on such a statistics.
+   //At the end these particels are left in the event according 
+   //to this factor: 1==all particles are left
+   //                0==all are removed 
+   
+    fTrackRejectionFactor=trf;
+    fHBTprocessor->SetTrackRejectionFactor(trf);
+  }
+/*******************************************************************/
+
+void AliGenHBTprocessor::SetRefControl(Int_t rc) //default 2
+ {
+ //Sets the Refernce Control Switch
+ //switch wether read reference histograms from file =1
+ //              compute from input events =2 - default
+   fReferenceControl=rc;
+   fHBTprocessor->SetRefControl(rc);
+ }
+/*******************************************************************/
+
+void AliGenHBTprocessor::SetPIDs(Int_t pid1,Int_t pid2)
+  {
+   //default pi+ pi-
+   //Sets PDG Codes of particles to be processed, default \\Pi^{+} and \\Pi^{-}
+   //This method accepts PDG codes which is
+   //in opposite to THBBProcessor which accepts GEANT PID
+    if ( (pid1 == 0) && (pid2 == 0) )
+    {
+      Error("AliGenHBTprocessor::SetPIDs","Sensless Particle Types setting: 0 0, Ignoring\n");
+    }
+    
+    fPid[0]=pid1;
+    fPid[1]=pid2;
+    
+    if(pid1 == pid2)
+       {
+        fHBTprocessor->SetPIDs(IdFromPDG(pid1) ,0);
+        SetNPIDtypes(1);
+        SetSwitchType(1);
+       }
+    else
+       { 
+        fHBTprocessor->SetPIDs(IdFromPDG(pid1) ,IdFromPDG(pid2));
+        SetNPIDtypes(2);
+       }
+  }
+/*******************************************************************/
+
+void AliGenHBTprocessor::SetNPIDtypes(Int_t npidt)
+  {
+    //Number ofparticle types to be processed - default 2
+    //see AliGenHBTprocessor::SetPIDs
+    
+    fNPidTypes = npidt;
+    fHBTprocessor->SetNPIDtypes(npidt);
+  }
+/*******************************************************************/
+
+void AliGenHBTprocessor::SetDeltap(Float_t deltp) 
+  {
+  //default = 0.1 GeV
+  //maximum range for random momentum shifts in GeV/c;
+  //px,py,pz independent; Default = 0.1 GeV/c.
+    fDeltap=deltp;
+    fHBTprocessor->SetDeltap(deltp); 
+  }
+/*******************************************************************/
+
+void AliGenHBTprocessor::SetMaxIterations(Int_t maxiter) 
+  { 
+  //maximum # allowed iterations thru all the 
+  //tracks for each event. Default = 50.
+  //If maxit=0, then calculate the correlations 
+  //for the input set of events without doing the
+  //track adjustment procedure.
+    
+    fMaxit=maxiter;
+    fHBTprocessor->SetMaxIterations(maxiter);
+  }
+
+/*******************************************************************/
+void AliGenHBTprocessor::SetDelChi(Float_t dc)
+  {
+    //min % change in total chi-square for which
+    //the track adjustment iterations may stop,
+    //Default = 0.1%.
+    
+    fDelchi=dc;
+    fHBTprocessor->SetDelChi(dc);
+  }
+/*******************************************************************/
+
+void AliGenHBTprocessor::SetIRand(Int_t irnd) 
+  { 
+    //if fact dummy - used only for compatibility
+    //we are using AliRoot built in RNG
+    //dummy in fact since we are using aliroot build-in RNG
+    //Sets renaodom number generator
+    fIrand=irnd;
+    fHBTprocessor->SetIRand(irnd);
+  }
+/*******************************************************************/
+      
+void AliGenHBTprocessor::SetLambda(Float_t lam) 
+  { 
+  //default = 0.6
+  // Sets Chaoticity Parameter
+    fLambda = lam;
+    fHBTprocessor->SetLambda(lam);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetR1d(Float_t r) 
+  {
+    //default 7.0
+    //Sets Spherical source model radius (fm)
+    fR1d = r;
+    fHBTprocessor->SetR1d(r);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetRSide(Float_t rs) 
+  {
+   //default rs = 6.0
+   //Rside,Rout,Rlong  = Non-spherical Bertsch-Pratt source model (fm)
+   
+    fRside = rs;
+    fHBTprocessor->SetRSide(rs);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetROut(Float_t ro) 
+  {
+    //default ro = 7.0
+    //Rside,Rout,Rlong  = Non-spherical Bertsch-Pratt source model (fm)
+    fRout = ro;
+    fHBTprocessor->SetROut(ro);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetRLong(Float_t rl) 
+  {
+    //default rl = 4.0
+    //Rside,Rout,Rlong  = Non-spherical Bertsch-Pratt source model (fm)
+    fRlong = rl;
+    fHBTprocessor->SetRLong(rl);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetRPerp(Float_t rp) 
+  {
+   //default rp = 6.0
+   //Rperp,Rparallel,R0= Non-spherical Yano-Koonin-Podgoretski source model (fm).
+    fRperp = rp;
+    fHBTprocessor->SetRPerp(rp);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetRParallel(Float_t rprl) 
+  { 
+   //default rprl = 4.0
+   //Rperp,Rparallel,R0= Non-spherical Yano-Koonin-Podgoretski source model (fm).
+    fRparallel = rprl;
+    fHBTprocessor->SetRParallel(rprl);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetR0(Float_t r0) 
+  {
+  //default r0 = 4.0
+  //Rperp,Rparallel,R0= Non-spherical Yano-Koonin-Podgoretski source model (fm).
+    fR0 = r0;
+    fHBTprocessor->SetR0(r0);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetQ0(Float_t q0) 
+  { 
+  //default q0 = 9.0
+  //Sets Q0                = NA35 Coulomb parameter for finite source size in (GeV/c)
+  //                         if fSwitchCoulomb = 2
+  //                       = Spherical Coulomb source radius in (fm) 
+  //                         if switchCoulomb = 3, used to interpolate the
+  //                         input Pratt/Cramer discrete Coulomb source
+  //                         radii tables.
+    fQ0 = q0;
+    fHBTprocessor->SetQ0(q0);
+  }
+
+/*******************************************************************/
+void AliGenHBTprocessor::SetSwitch1D(Int_t s1d) 
+  {
+//default s1d = 3
+// Sets fSwitch1d   
+//                          = 0 to not compute the 1D two-body //orrelations.
+//                          = 1 to compute this using Q-invariant
+//                          = 2 to compute this using Q-total
+//                          = 3 to compute this using Q-3-ve//tor
+
+    fSwitch1d = s1d;
+    fHBTprocessor->SetSwitch1D(s1d);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetSwitch3D(Int_t s3d) 
+  {
+//default s3d = 0
+// Sets fSwitch3d
+//                         = 0 to not compute the 3D two-body correlations.
+//                         = 1 to compute this using the side-out-long form
+//                         = 2 to compute this using the Yanno-Koonin-Pogoredskij form.   
+     
+    fSwitch3d = s3d;
+    fHBTprocessor->SetSwitch3D(s3d);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetSwitchType(Int_t st)
+  {
+//default st = 3
+//  Sets  switch_type       = 1 to fit only the like pair correlations
+//                          = 2 to fit only the unlike pair correlations
+//                          = 3 to fit both the like and unlike pair correl.
+//See SetPIDs and Init
+//If only one particle type is set, unly==1 makes sens
+  
+    fSwitchType = st;
+    fHBTprocessor->SetSwitchType(st);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetSwitchCoherence(Int_t sc)
+  {
+// default  sc = 0
+//        switchCoherence  = 0 for purely incoherent source (but can have
+//                              lambda < 1.0)
+//                          = 1 for mixed incoherent and coherent source
+  
+    fSwitchCoherence = sc;
+    fHBTprocessor->SetSwitchCoherence(sc);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetSwitchCoulomb(Int_t scol) 
+  {
+//default scol = 2
+//        switchCoulomb    = 0 no Coulomb correction
+//                          = 1 Point source Gamow correction only
+//                          = 2 NA35 finite source size correction
+//                          = 3 Pratt/Cramer finite source size correction;
+//                              interpolated from input tables.
+    fSwitchCoulomb =scol;
+    fHBTprocessor->SetSwitchCoulomb(scol);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetSwitchFermiBose(Int_t sfb)
+  {
+//default sfb = 1
+//        switchFermiBose =  1 Boson pairs
+//                          = -1 Fermion pairs
+
+    fSwitchFermiBose = sfb;
+    fHBTprocessor->SetSwitchFermiBose(sfb);
+  }
+/*******************************************************************/
+void AliGenHBTprocessor::SetPtRange(Float_t ptmin, Float_t ptmax)
+ {
+// default ptmin = 0.1, ptmax = 0.98
+//Sets Pt range (GeV)
+   AliGenerator::SetPtRange(ptmin,ptmax);
+   fHBTprocessor->SetPtRange(ptmin,ptmax);
+ }
+
+/*******************************************************************/
+void AliGenHBTprocessor::SetPxRange(Float_t pxmin, Float_t pxmax)
+ {
+//default pxmin = -1.0, pxmax = 1.0
+//Sets Px range 
+  fPxMin =pxmin;
+  fPxMax =pxmax;
+  fHBTprocessor->SetPxRange(pxmin,pxmax);
+ }
+/*******************************************************************/
+void AliGenHBTprocessor::SetPyRange(Float_t pymin, Float_t pymax)
+ {
+//default  pymin = -1.0, pymax = 1.0
+//Sets Py range 
+  fPyMin =pymin;
+  fPyMax =pymax;
+   fHBTprocessor->SetPyRange(pymin,pymax);
+ }
+/*******************************************************************/
+void AliGenHBTprocessor::SetPzRange(Float_t pzmin, Float_t pzmax)
+ {
+//default pzmin = -3.6, pzmax = 3.6
+//Sets Py range
+   fPzMin =pzmin;
+   fPzMax =pzmax; 
+   fHBTprocessor->SetPzRange(pzmin,pzmax);
+ }
+void AliGenHBTprocessor::SetMomentumRange(Float_t pmin, Float_t pmax)
+ {
+ //default pmin=0, pmax=0
+ //Do not use this method!
+    MayNotUse("AliGenHBTprocessor::SetMomentumRange Method is Dummy");
+ }
+ /*******************************************************************/
+void AliGenHBTprocessor::SetPhiRange(Float_t phimin, Float_t phimax)
+ {
+//
+//Sets \\Phi range  
+  AliGenerator::SetPhiRange(phimin,phimax);
+  
+  fHBTprocessor->SetPhiRange(phimin,phimax);
+ }
+/*******************************************************************/
+void AliGenHBTprocessor::SetEtaRange(Float_t etamin, Float_t etamax)
+ {
+//default etamin = -1.5, etamax = 1.5
+//Sets \\Eta range   
+   fEtaMin= etamin;
+   fEtaMax =etamax;
+   fHBTprocessor->SetEtaRange(etamin,etamax);
+   
+   //set the azimothal angle range in the AliGeneraor - 
+   //to keep coherency between azimuthal angle and pseudorapidity
+   //DO NOT CALL this->SetThetaRange, because it calls this method (where we are) 
+   //which must cause INFINITE LOOP
+   AliGenerator::SetThetaRange(RadiansToDegrees(EtaToTheta(fEtaMin)), 
+                               RadiansToDegrees(EtaToTheta(fEtaMax)));
+   
+ }
+/*******************************************************************/
+void AliGenHBTprocessor::SetThetaRange(Float_t thetamin, Float_t thetamax)
+{
+  //default thetamin = 0, thetamax = 180
+  //Azimuthal angle, override AliGenerator method which uses widely (i.e. wrapper generators)
+  //core fortran HBTProcessor uses Eta (pseudorapidity)
+  //so these methods has to be synchronized
+  
+  AliGenerator::SetThetaRange(thetamin,thetamax);
+  
+  SetEtaRange( ThetaToEta(fThetaMin) , ThetaToEta(fThetaMax) );
+
+}
+  
+/*******************************************************************/
+void AliGenHBTprocessor::SetNPtBins(Int_t nptbin)
+ {
+  //default nptbin = 50
+  //set number of Pt bins  
+   fNPtBins= nptbin; 
+   fHBTprocessor->SetNPtBins(nptbin);
+ }
+/*******************************************************************/
+void AliGenHBTprocessor::SetNPhiBins(Int_t nphibin)
+{ 
+  //default nphibin = 50
+  //set number of Phi bins
+  fNPhiBins=nphibin;
+  fHBTprocessor->SetNPhiBins(nphibin);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetNEtaBins(Int_t netabin)
+{
+  //default netabin = 50
+  //set number of Eta bins
+  fNEtaBins = netabin;
+  fHBTprocessor->SetNEtaBins(netabin);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetNPxBins(Int_t npxbin)
+{
+  //default  npxbin = 20
+  //set number of Px bins
+  fNPxBins = npxbin; 
+  fHBTprocessor->SetNPxBins(npxbin);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetNPyBins(Int_t npybin)
+{
+  //default  npybin = 20
+  //set number of Py bins
+  fNPyBins = npybin;
+  fHBTprocessor->SetNPyBins(npybin);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetNPzBins(Int_t npzbin)
+{
+  //default npzbin = 70
+  //set number of Pz bins
+  fNPzBins = npzbin;
+  fHBTprocessor->SetNPzBins(npzbin);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetNBins1DFineMesh(Int_t n)
+{
+//default n = 10
+//Sets the number of bins in the 1D mesh
+   fN1dFine =n;
+   fHBTprocessor->SetNBins1DFineMesh(n);
+   
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetBinSize1DFineMesh(Float_t x)
+{
+//default x=0.01
+//Sets the bin size in the 1D mesh
+   fBinsize1dFine = x;
+   fHBTprocessor->SetBinSize1DFineMesh(x);
+}
+/*******************************************************************/
+      
+void AliGenHBTprocessor::SetNBins1DCoarseMesh(Int_t n)
+{
+//default n =2
+//Sets the number of bins in the coarse 1D mesh
+  fN1dCoarse =n;
+  fHBTprocessor->SetNBins1DCoarseMesh(n);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetBinSize1DCoarseMesh(Float_t x)
+{
+//default x=0.05
+//Sets the bin size in the coarse 1D mesh
+  fBinsize1dCoarse =x;
+  fHBTprocessor->SetBinSize1DCoarseMesh(x);
+}
+/*******************************************************************/
+      
+void AliGenHBTprocessor::SetNBins3DFineMesh(Int_t n)
+{
+//default n = 8
+//Sets the number of bins in the 3D mesh
+  fN3dFine =n;
+  fHBTprocessor->SetNBins3DFineMesh(n);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetBinSize3DFineMesh(Float_t x)
+{
+//default x=0.01
+//Sets the bin size in the 3D mesh
+  fBinsize3dFine =x;
+  fHBTprocessor->SetBinSize3DFineMesh(x);
+}
+/*******************************************************************/
+      
+void AliGenHBTprocessor::SetNBins3DCoarseMesh(Int_t n)
+{
+//default n = 2
+//Sets the number of bins in the coarse 3D mesh
+
+  fN3dCoarse = n;
+  fHBTprocessor->SetNBins3DCoarseMesh(n);
+}
+/*******************************************************************/
+void AliGenHBTprocessor::SetBinSize3DCoarseMesh(Float_t x)
+{
+//default x=0.08
+//Sets the bin size in the coarse 3D mesh
+  fBinsize3dCoarse = x;
+  fHBTprocessor->SetBinSize3DCoarseMesh(x);
+}
+/*******************************************************************/
+      
+void AliGenHBTprocessor::SetNBins3DFineProjectMesh(Int_t n )
+{
+//default n =3
+//Sets the number of bins in the fine project mesh
+  fN3dFineProject = n;
+  fHBTprocessor->SetNBins3DFineProjectMesh(n);
+}
+/*******************************************************************/
+
+
+/*******************************************************************/
+
+
+
+
+
+
+void AliGenHBTprocessor::DefineParticles()
+{
+  //
+  // Load standard numbers for GEANT particles and PDG conversion
+  fNPDGCodes = 0; //this is done in the constructor - but in any case ...
+  
+  fPDGCode[fNPDGCodes++]=-99;   //  0 = unused location
+  fPDGCode[fNPDGCodes++]=22;    //  1 = photon
+  fPDGCode[fNPDGCodes++]=-11;   //  2 = positron
+  fPDGCode[fNPDGCodes++]=11;    //  3 = electron
+  fPDGCode[fNPDGCodes++]=12;    //  4 = neutrino e
+  fPDGCode[fNPDGCodes++]=-13;   //  5 = muon +
+  fPDGCode[fNPDGCodes++]=13;    //  6 = muon -
+  fPDGCode[fNPDGCodes++]=111;   //  7 = pi0
+  fPDGCode[fNPDGCodes++]=211;   //  8 = pi+
+  fPDGCode[fNPDGCodes++]=-211;  //  9 = pi-
+  fPDGCode[fNPDGCodes++]=130;   // 10 = Kaon Long
+  fPDGCode[fNPDGCodes++]=321;   // 11 = Kaon +
+  fPDGCode[fNPDGCodes++]=-321;  // 12 = Kaon -
+  fPDGCode[fNPDGCodes++]=2112;  // 13 = Neutron
+  fPDGCode[fNPDGCodes++]=2212;  // 14 = Proton
+  fPDGCode[fNPDGCodes++]=-2212; // 15 = Anti Proton
+  fPDGCode[fNPDGCodes++]=310;   // 16 = Kaon Short
+  fPDGCode[fNPDGCodes++]=221;   // 17 = Eta
+  fPDGCode[fNPDGCodes++]=3122;  // 18 = Lambda
+  fPDGCode[fNPDGCodes++]=3222;  // 19 = Sigma +
+  fPDGCode[fNPDGCodes++]=3212;  // 20 = Sigma 0
+  fPDGCode[fNPDGCodes++]=3112;  // 21 = Sigma -
+  fPDGCode[fNPDGCodes++]=3322;  // 22 = Xi0
+  fPDGCode[fNPDGCodes++]=3312;  // 23 = Xi-
+  fPDGCode[fNPDGCodes++]=3334;  // 24 = Omega-
+  fPDGCode[fNPDGCodes++]=-2112; // 25 = Anti Proton
+  fPDGCode[fNPDGCodes++]=-3122; // 26 = Anti Proton
+  fPDGCode[fNPDGCodes++]=-3222; // 27 = Anti Sigma -
+  fPDGCode[fNPDGCodes++]=-3212; // 28 = Anti Sigma 0
+  fPDGCode[fNPDGCodes++]=-3112; // 29 = Anti Sigma 0
+  fPDGCode[fNPDGCodes++]=-3322; // 30 = Anti Xi 0
+  fPDGCode[fNPDGCodes++]=-3312; // 31 = Anti Xi +
+  fPDGCode[fNPDGCodes++]=-3334; // 32 = Anti Omega +
+}  
+
+/*******************************************************************/
+Int_t AliGenHBTprocessor::IdFromPDG(Int_t pdg) const 
+{
+  //
+  // Return Geant3 code from PDG and pseudo ENDF code
+  //
+  for(Int_t i=0;i<fNPDGCodes;++i)
+    if(pdg==fPDGCode[i]) return i;
+  return -1;
+}
+Int_t AliGenHBTprocessor::PDGFromId(Int_t id) const
+{
+  //
+  // Return PDG code and pseudo ENDF code from Geant3 code
+  //
+  if(id>0 && id<fNPDGCodes) return fPDGCode[id];
+  else return -1;
+}
+Double_t AliGenHBTprocessor::ThetaToEta(Double_t arg)
+ {
+  //converts etha(azimuthal angle) to Eta (pseudorapidity). Argument in radians
+   
+   if(arg>= TMath::Pi()) return  709.0;//This number is the biggest wich not crashes exp(x)
+   if(arg<= 0.0) return -709.0;//
+   
+   arg -= TMath::Pi()/2.;
+   if (arg > 0.0) 
+    { 
+      return -TMath::Log( TMath::Tan(arg/2.)) ;
+    }
+   else 
+    {
+      return TMath::Log( TMath::Tan(-arg/2.)) ;
+    }
+ }
+                                  
+/*******************************************************************/
+/******      ROUTINES    USED    FOR     COMMUNUCATION      ********/
+/********************     WITH      FORTRAN     ********************/
+/*******************************************************************/
+
+#ifndef WIN32
+  # define hbtpran hbtpran_  
+  # define alihbtp_puttrack alihbtp_puttrack_
+  # define alihbtp_gettrack alihbtp_gettrack_
+  # define alihbtp_getnumberevents alihbtp_getnumberevents_
+  # define alihbtp_getnumbertracks  alihbtp_getnumbertracks_
+  # define alihbtp_initialize alihbtp_initialize_
+  # define alihbtp_setactiveeventnumber alihbtp_setactiveeventnumber_
+  # define alihbtp_setparameters alihbtp_setparameters_
+  # define type_ofCall
+
+#else
+  # define hbtpran HBTPRAN
+  # define alihbtp_puttrack ALIHBTP_PUTTRACK
+  # define alihbtp_gettrack ALIHBTP_GETTRACK
+  # define alihbtp_getnumberevents ALIHBTP_GETNUMBEREVENTS
+  # define alihbtp_getnumbertracks  ALIHBTP_GETNUMBERTRACKS
+  # define alihbtp_initialize ALIHBTP_INITIALIZE
+  # define alihbtp_setactiveeventnumber ALIHBTP_SETACTIVEEVENTNUMBER
+  # define alihbtp_setparameters ALIHBTP_SETPARAMETERS
+  # define type_ofCall  _stdcall
+#endif    
+
+#include "AliGenCocktailAfterBurner.h"
+#include <string.h>
+/*******************************************************************/
+
+AliGenCocktailAfterBurner*  GetGenerator()
+ {
+   // This function has two tasks:
+   // Check if environment is OK (exist gAlice and generator)
+   // Returns pointer to genrator
+   //to be changed with TFolders
+
+   if(!gAlice)
+    {
+      cout<<endl<<"ERROR: There is NO gALICE! Check what you are doing!"<<endl;
+      gROOT->Fatal("AliGenHBTprocessor.cxx: GetGenerator()",
+                    "\nRunning HBT Processor without gAlice... Exiting \n");
+      return 0x0;
+    }
+   AliGenerator * gen = gAlice->Generator();
+   
+   if (!gen) 
+    {
+      gAlice->Fatal("AliGenHBTprocessor.cxx: GetGenerator()",
+                   "\nThere is no generator in gAlice, exiting\n");
+      return 0x0;
+    }
+
+   //we do not sure actual type of the genetator
+   //and simple casting is risky - we use ROOT machinery to do safe cast
+   
+   TClass* cabclass = AliGenCocktailAfterBurner::Class(); //get AliGenCocktailAfterBurner TClass
+   TClass* genclass = gen->IsA();//get TClass of the generator we got from galice 
+   //use casting implemented in TClass
+   //cast gen to cabclass
+   AliGenCocktailAfterBurner* cab=(AliGenCocktailAfterBurner*)genclass->DynamicCast(cabclass,gen);
+                                                                        
+   if (cab == 0x0)//if generator that we got is not AliGenCocktailAfterBurner or its descendant we get null
+   {              //then quit with error
+      gAlice->Fatal("AliGenHBTprocessor.cxx: GetGenerator()",
+                    "\nThe main Generator is not a AliGenCocktailAfterBurner, exiting\n");
+      return 0x0;
+   }
+   //   cout<<endl<<"Got generator"<<endl;
+   return cab;
+   
+ }
+/*******************************************************************/
+
+AliGenHBTprocessor* GetAliGenHBTprocessor()
+{
+//returns pointer to the current instance of AliGenHBTprocessor in
+//AliGenCocktailAfterBurner (can be more than one)
+//
+ AliGenCocktailAfterBurner* gen = GetGenerator();
+ AliGenerator* g = gen->GetCurrentGenerator();
+ if(g == 0x0)
+  {
+    gAlice->Fatal("AliGenHBTprocessor.cxx: GetAliGenHBTprocessor()",
+                  "Can not get the current generator. Exiting");
+    return 0x0;
+  }
+  
+ TClass* hbtpclass = AliGenHBTprocessor::Class(); //get AliGenCocktailAfterBurner TClass
+ TClass* gclass = g->IsA();//get TClass of the current generator we got from CAB
+ AliGenHBTprocessor* hbtp = (AliGenHBTprocessor*)gclass->DynamicCast(hbtpclass,g);//try to cast 
+ if (hbtp == 0x0)
+   {
+      gAlice->Fatal("AliGenHBTprocessor.cxx: GetAliGenHBTprocessor()",
+                    "\nCurrernt generator in AliGenCocktailAfterBurner is not a AliGenHBTprocessor, exiting\n");
+      return 0x0;
+   }
+ return hbtp;
+}
+
+/*******************************************************************/
+extern "C" void type_ofCall alihbtp_setparameters()
+ {
+   //dummy
+ }
+
+extern "C" void type_ofCall  alihbtp_initialize()
+ {
+   //dummy
+ }
+
+/*******************************************************************/
+
+extern "C" void type_ofCall alihbtp_getnumberevents(Int_t &nev)
+ {
+  //passes number of events to the fortran 
+   if(gDebug) cout<<"alihbtp_getnumberevents("<<nev<<") ....";
+   AliGenCocktailAfterBurner* gen = GetGenerator();
+   if(!gen)
+    {
+     nev = -1;
+     return;
+    } 
+     
+    nev = gen->GetNumberOfEvents();
+    
+   if(gDebug>5) cout<<"EXITED N Ev = "<<nev<<endl; 
+   
+ }
+
+/*******************************************************************/
+
+extern "C" void type_ofCall  alihbtp_setactiveeventnumber(Int_t & nev)
+ {
+//sets active event in generator (AliGenCocktailAfterBurner)
+
+   if(gDebug>5) cout<<"alihbtp_setactiveeventnumber("<<nev<<") ....";
+   if(gDebug>0) cout<<"Asked for event "<<nev-1<<endl;
+   AliGenCocktailAfterBurner* gen = GetGenerator();
+   if(!gen) return;
+   gen->SetActiveEventNumber(nev - 1); //fortran numerates events from 1 to N
+   
+   if(gDebug>5) cout<<"EXITED returned "<<nev<<endl; 
+   
+ }
+/*******************************************************************/
+extern "C" void type_ofCall  alihbtp_getnumbertracks(Int_t &ntracks)
+ {
+//passes number of particles in active event to the fortran  
+   if(gDebug>5) cout<<"alihbtp_getnumbertracks("<<ntracks<<") ....";
+
+   AliGenCocktailAfterBurner* gen = GetGenerator();
+   if (!gen) 
+    {
+     ntracks = -1;
+     return;
+    } 
+
+   ntracks = gen->GetActiveStack()->GetNprimary();
+   if(gDebug>5) cout<<"EXITED Ntracks = "<<ntracks<<endl; 
+ }
+/*******************************************************************/
+extern "C" void type_ofCall  
+   alihbtp_puttrack(Int_t & n,Int_t& flag, Float_t& px, 
+                    Float_t& py, Float_t& pz, Int_t& geantpid)
+ {
+//sets new parameters (momenta) in track number n
+// in the active event
+// n - number of the track in active event
+// flag - flag of the track
+// px,py,pz - momenta
+// geantpid - type of the particle - Geant Particle ID
+   if(gDebug>5) cout<<"alihbtp_puttrack("<<n<<") ....";
+
+   AliGenCocktailAfterBurner* gen = GetGenerator();
+   if(!gen) return;
+   
+   TParticle * track = gen->GetActiveStack()->Particle(n-1);
+   
+   AliGenHBTprocessor* g = GetAliGenHBTprocessor();
+       
+   //check to be deleted 
+   if (geantpid != (g->IdFromPDG( track->GetPdgCode() )))
+    {
+      cerr<<endl<<" AliGenHBTprocessor.cxx: alihbtp_puttrack: SOMETHING IS GOING BAD:\n   GEANTPIDS ARE NOT THE SAME"<<endl;
+    }
+   
+   if(gDebug>0)
+     if (px != track->Px()) 
+       cout<<"Px diff. = "<<px - track->Px()<<endl;
+   
+   if(gDebug>3) cout<<" track->GetPdgCode() --> "<<track->GetPdgCode()<<endl;
+   
+   
+   
+   Float_t m = track->GetMass();
+   Float_t e = TMath::Sqrt(m*m+px*px+py*py+pz*pz);
+   track->SetMomentum(px,py,pz,e);
+   
+   g->SetHbtPStatusCode(flag,n-1);
+   
+   if(gDebug>5) cout<<"EXITED "<<endl; 
+ }
+
+/*******************************************************************/
+
+extern "C" void type_ofCall  
+  alihbtp_gettrack(Int_t & n,Int_t & flag, Float_t & px, 
+                   Float_t & py, Float_t & pz, Int_t & geantpid)
+  
+ {
+//passes track parameters to the fortran
+// n - number of the track in active event
+// flag - flag of the track
+// px,py,pz - momenta
+// geantpid - type of the particle - Geant Particle ID
+   if(gDebug>5) cout<<"alihbtp_gettrack("<<n<<") ....";
+   AliGenCocktailAfterBurner* gen = GetGenerator();
+
+   if (!gen) 
+    {
+     n = -1;
+     flag =-1;
+     px = py = pz = -1;
+     geantpid = -1;
+     return;
+    } 
+
+   TParticle * track = gen->GetActiveStack()->Particle(n-1);
+   AliGenHBTprocessor* g = GetAliGenHBTprocessor();
+   
+   flag = g->GetHbtPStatusCode(n-1);
+
+   px = (Float_t)track->Px();
+   py = (Float_t)track->Py();
+   pz = (Float_t)track->Pz();
+  
+   geantpid = g->IdFromPDG( track->GetPdgCode() );
+  
+   if(gDebug>5) cout<<"EXITED "<<endl; 
+ }
+
+/*******************************************************************/
+extern "C" Float_t type_ofCall hbtpran(Int_t &)
+{
+//interface to the random number generator
+  return sRandom->Rndm();
+}        
+
+/*******************************************************************/
+
+
+/*******************************************************************/
diff --git a/THbtp/AliGenHBTprocessor.h b/THbtp/AliGenHBTprocessor.h
new file mode 100644 (file)
index 0000000..fd14265
--- /dev/null
@@ -0,0 +1,246 @@
+// Implementation of the interface for THBTprocessor
+// which is a wrapper itself to Fortran 
+// program "HBT processor" written by Lanny Ray
+// Author: Piotr Krzysztof Skowronski <Piotr.Skowronski@cern.ch>
+// 
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+#ifndef ALIGENHBTPROCESSOR_H
+#define ALIGENHBTPROCESSOR_H
+
+#include <TPDGCode.h>
+
+#include "AliGenerator.h"
+
+class THBTprocessor;
+class TClonesArray;
+
+enum {kHBTPMaxParticleTypes = 50};
+
+class AliGenHBTprocessor : public AliGenerator 
+{ 
+//Wrapper class for THBTProcessor 
+//which is a wrapper itself to Fortran 
+//program "HBT processor" written by Lanny Ray
+//
+//Piotr.Skowronski@cern.ch
+
+  public:
+    AliGenHBTprocessor();
+    AliGenHBTprocessor(const AliGenHBTprocessor& in);
+    virtual ~AliGenHBTprocessor();
+
+    virtual void Init();
+    virtual void Generate();
+    virtual void GetParticles(TClonesArray * particles) const;
+    Int_t        IdFromPDG(Int_t pdg) const;
+    Int_t        PDGFromId(Int_t id) const;
+
+    Int_t        GetHbtPStatusCode(Int_t part) const; 
+    void         SetHbtPStatusCode(Int_t hbtstatcode, Int_t part);
+/************* S E T T E R S ******************/  
+
+    virtual void SetTrackRejectionFactor(Float_t trf = 1.0);
+
+    virtual void SetRefControl(Int_t rc =2);
+    virtual void SetPIDs(Int_t pid1 = kPiPlus,Int_t pid2 = kPiMinus); //PDG Codes of particles to be processed, default \\Pi^{+} and \\Pi^{-}
+    virtual void SetNPIDtypes(Int_t npidt = 2); //Number ofparticle types to be processed
+    virtual void SetDeltap(Float_t deltp = 0.1); //maximum range for random momentum shifts in GeV/c;
+                                                 //px,py,pz independent; Default = 0.1 GeV/c.
+    virtual void SetMaxIterations(Int_t maxiter = 50);//
+    virtual void SetDelChi(Float_t dc = 0.1);
+    virtual void SetIRand(Int_t irnd = 76564) ;
+     
+    virtual void SetLambda(Float_t lam = 0.6);
+    virtual void SetR1d(Float_t r = 7.0) ;
+    virtual void SetRSide(Float_t rs = 6.0);
+    virtual void SetROut(Float_t ro = 7.0) ;
+    virtual void SetRLong(Float_t rl = 4.0) ;
+    virtual void SetRPerp(Float_t rp = 6.0);
+    virtual void SetRParallel(Float_t rprl = 4.0);
+    virtual void SetR0(Float_t r0 = 4.0) ;
+    virtual void SetQ0(Float_t q0 = 9.0) ;
+    virtual void SetSwitch1D(Int_t s1d = 3);
+    virtual void SetSwitch3D(Int_t s3d = 0) ;
+    virtual void SetSwitchType(Int_t st = 3);
+    virtual void SetSwitchCoherence(Int_t sc = 0);
+    virtual void SetSwitchCoulomb(Int_t scol = 2);
+    virtual void SetSwitchFermiBose(Int_t sfb = 1);
+    
+    virtual void SetMomentumRange(Float_t pmin=0, Float_t pmax=0); //Dummy method
+    virtual void SetPtRange(Float_t ptmin = 0.1, Float_t ptmax = 0.98);
+    virtual void SetPxRange(Float_t pxmin = -1.0, Float_t pxmax = 1.0);
+    virtual void SetPyRange(Float_t pymin = -1.0, Float_t pymax = 1.0);  
+    virtual void SetPzRange(Float_t pzmin = -3.6, Float_t pzmax = 3.6);
+    
+    virtual void SetPhiRange(Float_t phimin = 0.0, Float_t phimax = 360.0);//Angle in degrees
+                                                                                   //coherent with AliGenCocktail
+                                                                                   //incohernet with AliGenerator
+    virtual void SetEtaRange(Float_t etamin = -1.5, Float_t etamax = 1.5);//Pseudorapidity
+    void SetThetaRange(Float_t thetamin = 0, Float_t thetamax = 180); //Azimuthal angle, override AliGenerator method
+                                                                      //which uses this, core fortran HBTProcessor uses Eta (pseudorapidity)
+                                                                      //so these methods has to be synchronized         
+
+    virtual void SetNPtBins(Int_t nptbin = 50);
+    virtual void SetNPhiBins(Int_t nphibin = 50);
+    virtual void SetNEtaBins(Int_t netabin = 50);
+    virtual void SetNPxBins(Int_t npxbin = 20);
+    virtual void SetNPyBins(Int_t npybin = 20);
+    virtual void SetNPzBins(Int_t npzbin = 70);
+   
+    
+    virtual void SetNBins1DFineMesh(Int_t n = 10);
+    virtual void SetBinSize1DFineMesh(Float_t x=0.01);
+      
+    virtual void SetNBins1DCoarseMesh(Int_t n =2 );
+    virtual void SetBinSize1DCoarseMesh(Float_t x=0.05);
+      
+    virtual void SetNBins3DFineMesh(Int_t n = 8);
+    virtual void SetBinSize3DFineMesh(Float_t x=0.01);
+      
+    virtual void SetNBins3DCoarseMesh(Int_t n = 2);
+    virtual void SetBinSize3DCoarseMesh(Float_t x=0.08);
+      
+    virtual void SetNBins3DFineProjectMesh(Int_t n =3 );
+/***********************************************************************/
+/* * * * * * *    P R O T E C T E D   A R E A    * * * * * * * * * * * */ 
+/***********************************************************************/
+  protected:
+    
+    THBTprocessor * fHBTprocessor;       //pointer to generator (TGenerator)
+    Int_t         **fHbtPStatCodes;      //! hbtp status codes of particles
+    Int_t           fNPDGCodes;          //! Number of defined particles   
+    Int_t           fPDGCode[kHBTPMaxParticleTypes]; //! PDG codes (for conversion PDG<->Geant)
+    void            DefineParticles();   //initiates array with PDG codes
+    void            InitStatusCodes();   //Initiates status codes (allocates memory and sets everything to zero) 
+    void            CleanStatusCodes();   //deletes array with status codes
+    /**********   P A R A M E T E R S  OF THE GENERATOR****************/
+          
+    Float_t fTrackRejectionFactor; //variates in range 0.0 <-> 1.0
+                                   //Describes the factor of particles rejected from the output.
+                                   //Used only in case of low muliplicity particles e.g. lambdas.
+                                   //Processor generates addisional particles and builds the 
+                                   //correletions on such a statistics.
+                                   //At the end these particels are left in the event according 
+                                   //to this factor: 1==all particles are left
+                                   //                0==all are removed
+      Int_t fReferenceControl;     //switch wether read reference histograms from file =1
+                                   //              compute from input events =2 - default
+      Int_t fPrintFull;             // Full print out option - each event
+      Int_t fPrintSectorData;       // Print sector overflow diagnostics
+      Int_t fNPidTypes;             // # particle ID types to correlate
+      Int_t fPid[2];                // Geant particle ID #s, max of 2 types
+      Int_t fNevents ;              // # events in input event text file
+      Int_t fSwitch1d;              // Include 1D correlations
+      Int_t fSwitch3d;              // Include 3D correlations
+      Int_t fSwitchType ;           // For like, unlike or both PID pairs
+      Int_t fSwitchCoherence;       // To include incoh/coher mixed source
+      Int_t fSwitchCoulomb;         // Coulomb correction selection options
+      Int_t fSwitchFermiBose;      // For fermions or bosons
+
+//   Counters:
+
+      Int_t fEventLineCounter;     // Input event text file line counter
+      Int_t fMaxit;                  // Max # iterations in track adjustment
+      Int_t fIrand;                  // Random # starting seed (Def=12345)      
+//                                    //    line counter
+
+//   Correlation Model Parameters:
+
+      Float_t    fLambda;               // Chaoticity parameter
+      Float_t    fR1d;                   // Spherical source radius (fm)
+      Float_t    fRside;                  // 3D Bertsch-Pratt source 'side' R (fm)
+      Float_t    fRout;                   // 3D Bertsch-Pratt source 'out'  R (fm)
+      Float_t    fRlong;                  // 3D Bertsch-Pratt source 'long' R (fm)
+      Float_t    fRperp;                  // 3D YKP source transverse radius  (fm)
+      Float_t    fRparallel;              // 3D YKP source longitudinal radius(fm)
+      Float_t    fR0;                     // 3D YKP source emission time durat(fm)
+      Float_t    fQ0;                     // NA35 Coulomb parameter (GeV/c) or
+//                                    // Coul radius for Pratt finite src (fm)
+
+//   Search Control Parameters:
+
+
+      Float_t    fDeltap;                 // Max limit for x,y,z momt shifts(GeV/c)
+      Float_t    fDelchi;                 // Min% change in Chi-Sq to stop iterat.
+
+
+//   Particle Masses:
+
+
+  /**********   M E S H  ****************/      
+
+
+      Int_t fNPtBins;                  // # one-body pt bins
+      Int_t fNPhiBins;                 // # one-body phi bins
+      Int_t fNEtaBins;                 // # one-body eta bins
+     
+      Int_t fN1dFine;                  // # bins for 1D, Fine Mesh
+      Int_t fN1dCoarse;                // # bins for 1D, Coarse Mesh
+      Int_t fN1dTotal;                 // Total # bins for 1D
+      Int_t fN3dFine ;                 // # bins for 3D, Fine Mesh
+      Int_t fN3dCoarse;                // # bins for 3D, Coarse Mesh
+      Int_t fN3dTotal;                 // Total # bins for 3D
+      Int_t fN3dFineProject;          // # 3D fine mesh bins to sum over for
+
+//   Momentum Space Sectors for Track Sorting:
+
+      Int_t fNPxBins;                  // # sector bins in px
+      Int_t fNPyBins;                  // # sector bins in py
+      Int_t fNPzBins;                  // # sector bins in pz
+      Int_t fNSectors;                  // Total # sectors in 3D momentum space
+
+     
+      Float_t    fPtBinSize ;          // One-body pt bin size in (GeV/c)
+
+      
+      Float_t    fPhiBinSize;          // One-body phi bin size in (degrees)
+      
+      Float_t    fEtaBinSize ;         // One-body eta bin size
+      Float_t    fEtaMin;              // One-body eta min
+      Float_t    fEtaMax;              // One-body eta max
+//   Two-Body Histograms and Correlation Mesh for 1D and 3D distributions:
+//                                       // projections onto single axis.
+
+      Float_t    fBinsize1dFine;       // Bin Size - 1D, Fine Mesh in (GeV/c)
+      Float_t    fBinsize1dCoarse;     // Bin Size - 1D, Coarse Mesh in (GeV/c)
+      Float_t    fQmid1d;               // q (GeV/c) at fine-coarse mesh boundary
+      Float_t    fQmax1d;               // Max q (GeV/c) for 1D distributions
+      Float_t    fBinsize3dFine;       // Bin Size - 3D, Fine Mesh in (GeV/c)
+      Float_t    fBinsize3dCoarse;     // Bin Size - 3D, Coarse Mesh in (GeV/c)
+      Float_t    fQmid3d;               // q (GeV/c) at fine-coarse mesh boundary
+      Float_t    fQmax3d;               // Max q (GeV/c) for 3D distributions
+
+      Float_t    fPxMin;                // Sector range in px in GeV/c
+      Float_t    fPxMax;                //--//--
+      Float_t    fDelpx;                 // Mom. space sector cell size - px(GeV/c)     
+      
+      Float_t    fPyMin;                // Sector range in py in GeV/c 
+      Float_t    fPyMax;                // --//--
+      Float_t    fDelpy;                 // Mom. space sector cell size - py(GeV/c)     
+
+      Float_t    fPzMin;                // Sector range in pz in GeV/c min
+      Float_t    fPzMax;                // Sector range in pz in GeV/c max
+      Float_t    fDelpz;                 // Mom. space sector cell size - pz(GeV/c)
+
+
+      /******* P R O T E C T E D   M E T H O D S  *****/
+  private:
+  public:    
+        //conveerts Eta (pseudorapidity) to etha(azimuthal angle). Returns radians 
+    static Double_t EtaToTheta(Double_t arg){return 2.*TMath::ATan(TMath::Exp(-arg));}
+        //converts etha(azimuthal angle) to Eta (pseudorapidity). Argument in radians
+    static Double_t ThetaToEta(Double_t arg);
+        //converts Degrees To Radians
+    static Double_t DegreesToRadians(Double_t arg){return arg*TMath::Pi()/180.;}
+          //converts Radians To Degrees 
+    static Double_t RadiansToDegrees(Double_t arg){return arg*180./TMath::Pi();}
+    
+    ClassDef(AliGenHBTprocessor,1) // Interface class for AliMevsim
+    
+};
+#include <Riostream.h>
+#endif
index cbbb694..35ae621 100644 (file)
@@ -5,5 +5,5 @@
 #pragma link off all functions;
  
 #pragma link C++ class THBTprocessor+;
-
+#pragma link C++ class AliGenHBTprocessor+;
 #endif
index 75351fa..a22ef7e 100644 (file)
@@ -1,8 +1,8 @@
-SRCS= THBTprocessor.cxx 
+SRCS= THBTprocessor.cxx AliGenHBTprocessor.cxx
 
 HDRS= $(SRCS:.cxx=.h)
 
-DHDR= THBTprocessorLinkDef.h
+DHDR= THBTprocessorLinkDef.h 
 
 EXPORT:=HBTprocCOMMON.h THBTprocessor.h
 
diff --git a/THijing/AliGenHijing.cxx b/THijing/AliGenHijing.cxx
new file mode 100644 (file)
index 0000000..27a567d
--- /dev/null
@@ -0,0 +1,822 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*
+$Log$
+Revision 1.47  2003/01/14 10:50:18  alibrary
+Cleanup of STEER coding conventions
+
+Revision 1.46  2003/01/07 14:12:33  morsch
+Provides collision geometry.
+
+Revision 1.45  2002/12/16 09:44:49  morsch
+Default for fRadiation is 3.
+
+Revision 1.44  2002/10/14 14:55:35  hristov
+Merging the VirtualMC branch to the main development branch (HEAD)
+
+Revision 1.42.4.1  2002/08/28 15:06:50  alibrary
+Updating to v3-09-01
+
+Revision 1.43  2002/08/09 12:09:52  morsch
+Direct gamma trigger correctly included.
+
+Revision 1.42  2002/03/12 11:07:08  morsch
+Add status code of particle to SetTrack call.
+
+Revision 1.41  2002/03/05 11:25:33  morsch
+- New quenching options
+- Correction in CheckTrigger()
+
+Revision 1.40  2002/02/12 11:05:53  morsch
+Get daughter indices right.
+
+Revision 1.39  2002/02/12 09:16:39  morsch
+Correction in SelectFlavor()
+
+Revision 1.38  2002/02/12 08:53:21  morsch
+SetNoGammas can be used to inhibit writing of gammas and pi0.
+
+Revision 1.37  2002/02/08 16:50:50  morsch
+Add name and title in constructor.
+
+Revision 1.36  2002/01/31 20:17:55  morsch
+Allow for triggered jets with simplified topology: Exact pT, back-to-back
+
+Revision 1.35  2001/12/13 07:56:25  hristov
+Set pointers to zero in the default constructor
+
+Revision 1.34  2001/12/11 16:55:42  morsch
+Correct initialization for jet phi-range.
+
+Revision 1.33  2001/12/05 10:18:51  morsch
+Possibility of kinematic biasing of jet phi range. (J. Klay)
+
+Revision 1.32  2001/11/28 13:51:11  morsch
+Introduce kinematic biasing (etamin, etamax) of jet trigger. Bookkeeping
+(number of trials) done in AliGenHijingEventHeader.
+
+Revision 1.31  2001/11/06 12:30:34  morsch
+Add Boost() method to boost all particles to LHC lab frame. Needed for asymmetric collision systems.
+
+Revision 1.30  2001/10/21 18:35:56  hristov
+Several pointers were set to zero in the default constructors to avoid memory management problems
+
+Revision 1.29  2001/10/15 08:12:24  morsch
+- Vertex smearing with truncated gaussian.
+- Store triggered jet info before and after final state radiation into mc-heade
+
+Revision 1.28  2001/10/08 11:55:25  morsch
+Store 4-momenta of trigegred jets in event header.
+Possibility to switch of initial and final state radiation.
+
+Revision 1.27  2001/10/08 07:13:14  morsch
+Add setter for minimum transverse momentum of triggered jet.
+
+Revision 1.26  2001/10/04 08:12:24  morsch
+Redefinition of stable condition.
+
+Revision 1.25  2001/07/27 17:09:36  morsch
+Use local SetTrack, KeepTrack and SetHighWaterMark methods
+to delegate either to local stack or to stack owned by AliRun.
+(Piotr Skowronski, A.M.)
+
+Revision 1.24  2001/07/20 09:34:56  morsch
+Count the number of spectator neutrons and protons and add information
+to the event header. (Chiara Oppedisano)
+
+Revision 1.23  2001/07/13 17:30:22  morsch
+Derive from AliGenMC.
+
+Revision 1.22  2001/06/11 13:09:23  morsch
+- Store cross-Section and number of binary collisions as a function of impact parameter
+- Pass AliGenHijingEventHeader to gAlice.
+
+Revision 1.21  2001/02/28 17:35:24  morsch
+Consider elastic interactions (ks = 1 and ks = 11) as spectator (Chiara Oppedisano)
+
+Revision 1.20  2001/02/14 15:50:40  hristov
+The last particle in event marked using SetHighWaterMark
+
+Revision 1.19  2000/12/21 16:24:06  morsch
+Coding convention clean-up
+
+Revision 1.18  2000/12/06 17:46:30  morsch
+Avoid random numbers 1 and 0.
+
+Revision 1.17  2000/12/04 11:22:03  morsch
+Init of sRandom as in 1.15
+
+Revision 1.16  2000/12/02 11:41:39  morsch
+Use SetRandom() to initialize random number generator in constructor.
+
+Revision 1.15  2000/11/30 20:29:02  morsch
+Initialise static variable sRandom in constructor: sRandom = fRandom;
+
+Revision 1.14  2000/11/30 07:12:50  alibrary
+Introducing new Rndm and QA classes
+
+Revision 1.13  2000/11/09 17:40:27  morsch
+Possibility to select/unselect spectator protons and neutrons.
+Method SetSpectators(Int_t spect) added. (FCA, Ch. Oppedisano)
+
+Revision 1.12  2000/10/20 13:38:38  morsch
+Debug printouts commented.
+
+Revision 1.11  2000/10/20 13:22:26  morsch
+- skip particle type 92 (string)
+- Charmed and beauty baryions (5122, 4122) are considered as stable consistent with
+  mesons.
+
+Revision 1.10  2000/10/17 15:10:20  morsch
+Write first all the parent particles to the stack and then the final state particles.
+
+Revision 1.9  2000/10/17 13:38:59  morsch
+Protection against division by zero in EvaluateCrossSection() and KinematicSelection(..)     (FCA)
+
+Revision 1.8  2000/10/17 12:46:31  morsch
+Protect EvaluateCrossSections() against division by zero.
+
+Revision 1.7  2000/10/02 21:28:06  fca
+Removal of useless dependecies via forward declarations
+
+Revision 1.6  2000/09/11 13:23:37  morsch
+Write last seed to file (fortran lun 50) and reed back from same lun using calls to
+luget_hijing and luset_hijing.
+
+Revision 1.5  2000/09/07 16:55:40  morsch
+fHijing->Initialize(); after change of parameters. (Dmitri Yurevitch Peressounko)
+
+Revision 1.4  2000/07/11 18:24:56  fca
+Coding convention corrections + few minor bug fixes
+
+Revision 1.3  2000/06/30 12:08:36  morsch
+In member data: char* replaced by TString, Init takes care of resizing the strings to
+8 characters required by Hijing.
+
+Revision 1.2  2000/06/15 14:15:05  morsch
+Add possibility for heavy flavor selection: charm and beauty.
+
+Revision 1.1  2000/06/09 20:47:27  morsch
+AliGenerator interface class to HIJING using THijing (test version)
+
+*/
+
+
+
+// Generator using HIJING as an external generator
+// The main HIJING options are accessable for the user through this interface.
+// Uses the THijing implementation of TGenerator.
+//
+// andreas.morsch@cern.ch
+
+#include <TArrayI.h>
+#include <TGraph.h>
+#include <THijing.h>
+#include <TLorentzVector.h>
+#include <TPDGCode.h>
+#include <TParticle.h>
+
+#include "AliGenHijing.h"
+#include "AliGenHijingEventHeader.h"
+#include "AliRun.h"
+
+
+  ClassImp(AliGenHijing)
+
+AliGenHijing::AliGenHijing()
+                 :AliGenMC()
+{
+// Constructor
+    fParticles = 0;
+    fHijing    = 0;
+    fDsigmaDb  = 0;
+    fDnDb      = 0;
+}
+
+AliGenHijing::AliGenHijing(Int_t npart)
+    :AliGenMC(npart)
+{
+// Default PbPb collisions at 5. 5 TeV
+//
+    fName = "Hijing";
+    fTitle= "Particle Generator using HIJING";
+
+    SetEnergyCMS();
+    SetImpactParameterRange();
+    SetTarget();
+    SetProjectile();
+    SetBoostLHC();
+    SetJetEtaRange();
+    SetJetPhiRange();
+    
+    fKeep       =  0;
+    fQuench     =  1;
+    fShadowing  =  1;
+    fTrigger    =  0;
+    fDecaysOff  =  1;
+    fEvaluate   =  0;
+    fSelectAll  =  0;
+    fFlavor     =  0;
+    fSpectators =  1;
+    fDsigmaDb   =  0;
+    fDnDb       =  0;
+    fPtMinJet   = -2.5;        
+    fRadiation  =  3;
+    fEventVertex.Set(3);
+//
+    SetSimpleJets();
+    SetNoGammas();
+//
+    fParticles = new TClonesArray("TParticle",10000);    
+//
+// Set random number generator   
+    sRandom = fRandom;
+    fHijing = 0;
+
+}
+
+AliGenHijing::AliGenHijing(const AliGenHijing & Hijing)
+{
+// copy constructor
+}
+
+
+AliGenHijing::~AliGenHijing()
+{
+// Destructor
+    if ( fDsigmaDb) delete  fDsigmaDb;  
+    if ( fDnDb)     delete  fDnDb;  
+    delete fParticles;
+}
+
+void AliGenHijing::Init()
+{
+// Initialisation
+    fFrame.Resize(8);
+    fTarget.Resize(8);
+    fProjectile.Resize(8);
+    
+    SetMC(new THijing(fEnergyCMS, fFrame, fProjectile, fTarget, 
+                     fAProjectile, fZProjectile, fATarget, fZTarget, 
+                     fMinImpactParam, fMaxImpactParam));
+
+    fHijing=(THijing*) fgMCEvGen;
+    fHijing->SetIHPR2(2,  fRadiation);
+    fHijing->SetIHPR2(3,  fTrigger);
+    fHijing->SetIHPR2(6,  fShadowing);
+    fHijing->SetIHPR2(12, fDecaysOff);    
+    fHijing->SetIHPR2(21, fKeep);
+    fHijing->SetHIPR1(10, fPtMinJet);  
+    fHijing->SetHIPR1(50, fSimpleJet);
+//
+//  Quenching
+//
+//
+//  fQuench = 0:  no quenching
+//  fQuench = 1:  hijing default
+//  fQuench = 2:  new LHC  parameters for HIPR1(11) and HIPR1(14)
+//  fQuench = 3:  new RHIC parameters for HIPR1(11) and HIPR1(14)
+//  fQuench = 4:  new LHC  parameters with log(e) dependence
+//  fQuench = 5:  new RHIC parameters with log(e) dependence
+    fHijing->SetIHPR2(50, 0);
+    if (fQuench > 0) 
+       fHijing->SetIHPR2(4,  1);
+    else
+       fHijing->SetIHPR2(4,  0);
+// New LHC parameters from Xin-Nian Wang
+    if (fQuench == 2) {
+       fHijing->SetHIPR1(14, 1.1);
+       fHijing->SetHIPR1(11, 3.7);
+    } else if (fQuench == 3) {
+       fHijing->SetHIPR1(14, 0.20);
+       fHijing->SetHIPR1(11, 2.5);
+    } else if (fQuench == 4) {
+       fHijing->SetIHPR2(50, 1);
+       fHijing->SetHIPR1(14, 4.*0.34);
+       fHijing->SetHIPR1(11, 3.7);
+    } else if (fQuench == 5) {
+       fHijing->SetIHPR2(50, 1);
+       fHijing->SetHIPR1(14, 0.34);
+       fHijing->SetHIPR1(11, 2.5);
+    }
+    
+    
+    
+//
+//  Initialize Hijing  
+//    
+    fHijing->Initialize();
+//
+    if (fEvaluate) EvaluateCrossSections();
+//
+}
+
+void AliGenHijing::Generate()
+{
+// Generate one event
+
+  Float_t polar[3]    =   {0,0,0};
+  Float_t origin[3]   =   {0,0,0};
+  Float_t origin0[3]  =   {0,0,0};
+  Float_t p[3], random[6];
+  Float_t tof;
+
+//  converts from mm/c to s
+  const Float_t kconv = 0.001/2.999792458e8;
+//
+  Int_t nt  = 0;
+  Int_t jev = 0;
+  Int_t j, kf, ks, imo;
+  kf = 0;
+    
+
+    
+  fTrials = 0;
+  for (j = 0;j < 3; j++) origin0[j] = fOrigin[j];
+  if(fVertexSmear == kPerEvent) {
+      Float_t dv[3];
+      dv[2] = 1.e10;
+      while(TMath::Abs(dv[2]) > fCutVertexZ*fOsigma[2]) {
+         Rndm(random,6);
+         for (j=0; j < 3; j++) {
+             dv[j] = fOsigma[j]*TMath::Cos(2*random[2*j]*TMath::Pi())*
+                 TMath::Sqrt(-2*TMath::Log(random[2*j+1]));
+         }
+      }
+      for (j=0; j < 3; j++) origin0[j] += dv[j];
+  } else if (fVertexSmear == kPerTrack) {
+//         fHijing->SetMSTP(151,0);
+      for (j = 0; j < 3; j++) {
+//           fHijing->SetPARP(151+j, fOsigma[j]*10.);
+      }
+  }
+  while(1)
+  {
+//    Generate one event
+// --------------------------------------------------------------------------
+      fSpecn   = 0;  
+      fSpecp   = 0;
+// --------------------------------------------------------------------------
+      fHijing->GenerateEvent();
+      fTrials++;
+      fHijing->ImportParticles(fParticles,"All");
+      if (fTrigger != kNoTrigger) {
+         if (!CheckTrigger()) continue;
+      }
+      if (fLHC) Boost();
+      
+      
+      Int_t np = fParticles->GetEntriesFast();
+      printf("\n **************************************************%d\n",np);
+      Int_t nc = 0;
+      if (np == 0 ) continue;
+      Int_t i;
+      Int_t* newPos     = new Int_t[np];
+      Int_t* pSelected  = new Int_t[np];
+
+      for (i = 0; i < np; i++) {
+         newPos[i]    = i;
+         pSelected[i] = 0;
+      }
+      
+//      Get event vertex
+//
+      TParticle *  iparticle = (TParticle *) fParticles->At(0);
+      fEventVertex[0] = origin0[0];
+      fEventVertex[1] = origin0[1];    
+      fEventVertex[2] = origin0[2];
+      
+//
+//      First select parent particles
+//
+
+      for (i = 0; i < np; i++) {
+         iparticle = (TParticle *) fParticles->At(i);
+
+// Is this a parent particle ?
+         if (Stable(iparticle)) continue;
+//
+         Bool_t  selected             =  kTRUE;
+         Bool_t  hasSelectedDaughters =  kFALSE;
+         
+         
+         kf        = iparticle->GetPdgCode();
+         ks        = iparticle->GetStatusCode();
+         if (kf == 92) continue;
+           
+         if (!fSelectAll) selected = KinematicSelection(iparticle, 0) && 
+                              SelectFlavor(kf);
+         hasSelectedDaughters = DaughtersSelection(iparticle);
+//
+// Put particle on the stack if it is either selected or 
+// it is the mother of at least one seleted particle
+//
+         if (selected || hasSelectedDaughters) {
+             nc++;
+             pSelected[i] = 1;
+         } // selected
+      } // particle loop parents
+//
+// Now select the final state particles
+//
+
+      for (i = 0; i<np; i++) {
+         TParticle *  iparticle = (TParticle *) fParticles->At(i);
+// Is this a final state particle ?
+         if (!Stable(iparticle)) continue;
+      
+         Bool_t  selected             =  kTRUE;
+         kf        = iparticle->GetPdgCode();
+         ks        = iparticle->GetStatusCode();
+         
+// --------------------------------------------------------------------------
+// Count spectator neutrons and protons
+         if(ks == 0 || ks == 1 || ks == 10 || ks == 11){
+             if(kf == kNeutron) fSpecn += 1;
+             if(kf == kProton)  fSpecp += 1;
+         }
+// --------------------------------------------------------------------------
+//         
+         if (!fSelectAll) {
+             selected = KinematicSelection(iparticle,0)&&SelectFlavor(kf);
+             if (!fSpectators && selected) selected = (ks != 0 && ks != 1 && ks != 10
+                                                       && ks != 11);
+         }
+//
+// Put particle on the stack if selected
+//
+         if (selected) {
+             nc++;
+             pSelected[i] = 1;
+         } // selected
+      } // particle loop final state
+//
+// Write particles to stack
+//
+      for (i = 0; i<np; i++) {
+         TParticle *  iparticle = (TParticle *) fParticles->At(i);
+         Bool_t  hasMother   = (iparticle->GetFirstMother()     >=0);
+         Bool_t  hasDaughter = (iparticle->GetFirstDaughter()   >=0);
+
+         if (pSelected[i]) {
+             kf   = iparticle->GetPdgCode();
+             ks   = iparticle->GetStatusCode();
+             p[0] = iparticle->Px();
+             p[1] = iparticle->Py();
+             p[2] = iparticle->Pz();
+             origin[0] = origin0[0]+iparticle->Vx()/10;
+             origin[1] = origin0[1]+iparticle->Vy()/10;
+             origin[2] = origin0[2]+iparticle->Vz()/10;
+             tof = kconv*iparticle->T();
+             imo = -1;
+             TParticle* mother = 0;
+             if (hasMother) {
+                 imo = iparticle->GetFirstMother();
+                 mother = (TParticle *) fParticles->At(imo);
+                 imo = (mother->GetPdgCode() != 92) ? imo = newPos[imo] : -1;
+             } // if has mother   
+             Bool_t tFlag = (fTrackIt && !hasDaughter);
+             SetTrack(tFlag,imo,kf,p,origin,polar,
+                      tof,kPNoProcess,nt, 1., ks);
+             KeepTrack(nt);
+             newPos[i] = nt;
+         } // if selected
+      } // particle loop
+      delete[] newPos;
+      delete[] pSelected;
+      
+      printf("\n I've put %i particles on the stack \n",nc);
+      if (nc > 0) {
+         jev += nc;
+         if (jev >= fNpart || fNpart == -1) {
+             fKineBias = Float_t(fNpart)/Float_t(fTrials);
+             printf("\n Trials: %i %i %i\n",fTrials, fNpart, jev);
+             break;
+         }
+      }
+  } // event loop
+  MakeHeader();
+  SetHighWaterMark(nt);
+}
+
+void AliGenHijing::KeepFullEvent()
+{
+    fKeep=1;
+}
+
+void AliGenHijing::EvaluateCrossSections()
+{
+//     Glauber Calculation of geometrical x-section
+//
+    Float_t xTot       = 0.;          // barn
+    Float_t xTotHard   = 0.;          // barn 
+    Float_t xPart      = 0.;          // barn
+    Float_t xPartHard  = 0.;          // barn 
+    Float_t sigmaHard  = 0.1;         // mbarn
+    Float_t bMin       = 0.;
+    Float_t bMax       = fHijing->GetHIPR1(34)+fHijing->GetHIPR1(35);
+    const Float_t kdib = 0.2;
+    Int_t   kMax       = Int_t((bMax-bMin)/kdib)+1;
+
+
+    printf("\n Projectile Radius (fm): %f \n",fHijing->GetHIPR1(34));
+    printf("\n Target     Radius (fm): %f \n",fHijing->GetHIPR1(35));    
+    Int_t i;
+    Float_t oldvalue= 0.;
+
+    Float_t* b   = new Float_t[kMax];
+    Float_t* si1 = new Float_t[kMax];    
+    Float_t* si2 = new Float_t[kMax];    
+    
+    for (i = 0; i < kMax; i++)
+    {
+       Float_t xb  = bMin+i*kdib;
+       Float_t ov;
+       ov=fHijing->Profile(xb);
+       Float_t gb  =  2.*0.01*fHijing->GetHIPR1(40)*kdib*xb*(1.-TMath::Exp(-fHijing->GetHINT1(12)*ov));
+       Float_t gbh =  2.*0.01*fHijing->GetHIPR1(40)*kdib*xb*sigmaHard*ov;
+       xTot+=gb;
+       xTotHard += gbh;
+       printf("profile %f %f %f\n", xb, ov, fHijing->GetHINT1(12));
+       
+       if (xb > fMinImpactParam && xb < fMaxImpactParam)
+       {
+           xPart += gb;
+           xPartHard += gbh;
+       }
+       
+       if(oldvalue) if ((xTot-oldvalue)/oldvalue<0.0001) break;
+       oldvalue = xTot;
+       printf("\n Total cross section (barn): %d %f %f \n",i, xb, xTot);
+       printf("\n Hard  cross section (barn): %d %f %f \n\n",i, xb, xTotHard);
+       if (i>0) {
+           si1[i] = gb/kdib;
+           si2[i] = gbh/gb;
+           b[i]  = xb;
+       }
+    }
+
+    printf("\n Total cross section (barn): %f \n",xTot);
+    printf("\n Hard  cross section (barn): %f \n \n",xTotHard);
+    printf("\n Partial       cross section (barn): %f %f \n",xPart, xPart/xTot*100.);
+    printf("\n Partial  hard cross section (barn): %f %f \n",xPartHard, xPartHard/xTotHard*100.);
+
+//  Store result as a graph
+    b[0] = 0;
+    si1[0] = 0;
+    si2[0]=si2[1];
+    
+    fDsigmaDb  = new TGraph(i, b, si1);
+    fDnDb      = new TGraph(i, b, si2);
+}
+
+Bool_t AliGenHijing::DaughtersSelection(TParticle* iparticle)
+{
+//
+// Looks recursively if one of the daughters has been selected
+//
+//    printf("\n Consider daughters %d:",iparticle->GetPdgCode());
+    Int_t imin = -1;
+    Int_t imax = -1;
+    Int_t i;
+    Bool_t hasDaughters = (iparticle->GetFirstDaughter() >=0);
+    Bool_t selected = kFALSE;
+    if (hasDaughters) {
+       imin = iparticle->GetFirstDaughter();
+       imax = iparticle->GetLastDaughter();       
+       for (i = imin; i <= imax; i++){
+           TParticle *  jparticle = (TParticle *) fParticles->At(i);   
+           Int_t ip = jparticle->GetPdgCode();
+           if (KinematicSelection(jparticle,0)&&SelectFlavor(ip)) {
+               selected=kTRUE; break;
+           }
+           if (DaughtersSelection(jparticle)) {selected=kTRUE; break; }
+       }
+    } else {
+       return kFALSE;
+    }
+    return selected;
+}
+
+
+Bool_t AliGenHijing::SelectFlavor(Int_t pid)
+{
+// Select flavor of particle
+// 0: all
+// 4: charm and beauty
+// 5: beauty
+    Bool_t res = 0;
+    
+    if (fFlavor == 0) {
+       res = kTRUE;
+    } else {
+       Int_t ifl = TMath::Abs(pid/100);
+       if (ifl > 10) ifl/=10;
+       res = (fFlavor == ifl);
+    }
+//
+//  This part if gamma writing is inhibited
+    if (fNoGammas) 
+       res = res && (pid != kGamma && pid != kPi0);
+//
+    return res;
+}
+
+Bool_t AliGenHijing::Stable(TParticle*  particle)
+{
+// Return true for a stable particle
+//
+    
+    if (particle->GetFirstDaughter() < 0 )
+    {
+       return kTRUE;
+    } else {
+       return kFALSE;
+    }
+}
+
+
+void AliGenHijing::Boost()
+{
+//
+// Boost cms into LHC lab frame
+//
+    Double_t dy    = - 0.5 * TMath::Log(Double_t(fZProjectile) * Double_t(fATarget) / 
+                                     (Double_t(fZTarget)    * Double_t(fAProjectile)));
+    Double_t beta  = TMath::TanH(dy);
+    Double_t gamma = 1./TMath::Sqrt(1.-beta*beta);
+    Double_t gb    = gamma * beta;
+
+    printf("\n Boosting particles to lab frame %f %f %f", dy, beta, gamma);
+    
+    Int_t i;
+    Int_t np = fParticles->GetEntriesFast();
+    for (i = 0; i < np; i++) 
+    {
+       TParticle* iparticle = (TParticle*) fParticles->At(i);
+
+       Double_t e   = iparticle->Energy();
+       Double_t px  = iparticle->Px();
+       Double_t py  = iparticle->Py();
+       Double_t pz  = iparticle->Pz();
+
+       Double_t eb  = gamma * e -      gb * pz;
+       Double_t pzb =   -gb * e +   gamma * pz;
+
+       iparticle->SetMomentum(px, py, pzb, eb);
+    }
+}
+
+
+void AliGenHijing::MakeHeader()
+{
+// Builds the event header, to be called after each event
+    AliGenEventHeader* header = new AliGenHijingEventHeader("Hijing");
+    ((AliGenHijingEventHeader*) header)->SetNProduced(fHijing->GetNATT());
+    ((AliGenHijingEventHeader*) header)->SetImpactParameter(fHijing->GetHINT1(19));
+    ((AliGenHijingEventHeader*) header)->SetTotalEnergy(fHijing->GetEATT());
+    ((AliGenHijingEventHeader*) header)->SetHardScatters(fHijing->GetJATT());
+    ((AliGenHijingEventHeader*) header)->SetParticipants(fHijing->GetNP(), fHijing->GetNT());
+    ((AliGenHijingEventHeader*) header)->SetCollisions(fHijing->GetN0(),
+                                                      fHijing->GetN01(),
+                                                      fHijing->GetN10(),
+                                                      fHijing->GetN11());
+    ((AliGenHijingEventHeader*) header)->SetSpectators(fSpecn, fSpecp);
+
+// 4-momentum vectors of the triggered jets.
+//
+// Before final state gluon radiation.
+    TLorentzVector* jet1 = new TLorentzVector(fHijing->GetHINT1(21), 
+                                             fHijing->GetHINT1(22),
+                                             fHijing->GetHINT1(23),
+                                             fHijing->GetHINT1(24));
+
+    TLorentzVector* jet2 = new TLorentzVector(fHijing->GetHINT1(31), 
+                                             fHijing->GetHINT1(32),
+                                             fHijing->GetHINT1(33),
+                                             fHijing->GetHINT1(34));
+// After final state gluon radiation.
+    TLorentzVector* jet3 = new TLorentzVector(fHijing->GetHINT1(26), 
+                                             fHijing->GetHINT1(27),
+                                             fHijing->GetHINT1(28),
+                                             fHijing->GetHINT1(29));
+
+    TLorentzVector* jet4 = new TLorentzVector(fHijing->GetHINT1(36), 
+                                             fHijing->GetHINT1(37),
+                                             fHijing->GetHINT1(38),
+                                             fHijing->GetHINT1(39));
+    ((AliGenHijingEventHeader*) header)->SetJets(jet1, jet2, jet3, jet4);
+// Bookkeeping for kinematic bias
+    ((AliGenHijingEventHeader*) header)->SetTrials(fTrials);
+// Event Vertex
+    header->SetPrimaryVertex(fEventVertex);
+    gAlice->SetGenEventHeader(header);   
+    fCollisionGeometry = (AliGenHijingEventHeader*)  header;
+}
+
+Bool_t AliGenHijing::CheckTrigger()
+{
+// Check the kinematic trigger condition
+//
+    Bool_t   triggered = kFALSE;
+    if (fTrigger == 1) {
+//
+//  jet-jet Trigger    
+       
+       TLorentzVector* jet1 = new TLorentzVector(fHijing->GetHINT1(26), 
+                                                 fHijing->GetHINT1(27),
+                                                 fHijing->GetHINT1(28),
+                                                 fHijing->GetHINT1(29));
+       
+       TLorentzVector* jet2 = new TLorentzVector(fHijing->GetHINT1(36), 
+                                                 fHijing->GetHINT1(37),
+                                                 fHijing->GetHINT1(38),
+                                                 fHijing->GetHINT1(39));
+       Double_t eta1      = jet1->Eta();
+       Double_t eta2      = jet2->Eta();
+       Double_t phi1      = jet1->Phi();
+       Double_t phi2      = jet2->Phi();
+//    printf("\n Trigger: %f %f %f %f",
+//        fEtaMinJet, fEtaMaxJet, fPhiMinJet, fPhiMaxJet);
+       if (
+           (eta1 < fEtaMaxJet && eta1 > fEtaMinJet &&  
+            phi1 < fPhiMaxJet && phi1 > fPhiMinJet) 
+           ||
+           (eta2 < fEtaMaxJet && eta2 > fEtaMinJet &&  
+            phi2 < fPhiMaxJet && phi2 > fPhiMinJet)
+           ) 
+           triggered = kTRUE;
+    } else if (fTrigger == 2) {
+//  Gamma Jet
+//
+       Int_t np = fParticles->GetEntriesFast();
+       for (Int_t i = 0; i < np; i++) {
+           TParticle* part = (TParticle*) fParticles->At(i);
+           Int_t kf = part->GetPdgCode();
+           Int_t ks = part->GetStatusCode();
+           if (kf == 22 && ks == 40) {
+               Float_t phi = part->Phi();
+               Float_t eta = part->Eta();
+               if  (eta < fEtaMaxJet && 
+                    eta > fEtaMinJet &&
+                    phi < fPhiMaxJet && 
+                    phi > fPhiMinJet) {
+                   triggered = 1;
+                   break;
+               } // check phi,eta within limits
+           } // direct gamma ? 
+       } // particle loop
+    } // fTrigger == 2
+    return triggered;
+}
+
+
+
+
+AliGenHijing& AliGenHijing::operator=(const  AliGenHijing& rhs)
+{
+// Assignment operator
+    return *this;
+}
+
+#ifndef WIN32
+# define rluget_hijing rluget_hijing_
+# define rluset_hijing rluset_hijing_
+# define rlu_hijing    rlu_hijing_
+# define type_of_call
+#else
+# define rluget_hijing RLUGET_HIJING
+# define rluset_hijing RLUSET_HIJING
+# define rlu_hijing    RLU_HIJING
+# define type_of_call _stdcall
+#endif
+
+
+extern "C" {
+  void type_of_call rluget_hijing(Int_t & /*lfn*/, Int_t & /*move*/)
+  {printf("Dummy version of rluget_hijing reached\n");}
+
+  void type_of_call rluset_hijing(Int_t & /*lfn*/, Int_t & /*move*/)
+  {printf("Dummy version of rluset_hijing reached\n");}
+
+  Double_t type_of_call rlu_hijing(Int_t & /*idum*/) 
+  {
+      Float_t r;
+      do r=sRandom->Rndm(); while(0 >= r || r >= 1);
+      return r;
+  }
+}
diff --git a/THijing/AliGenHijing.h b/THijing/AliGenHijing.h
new file mode 100644 (file)
index 0000000..20bce42
--- /dev/null
@@ -0,0 +1,158 @@
+#ifndef ALIGENHIJING_H
+#define ALIGENHIJING_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+// Generator using HIJING as an external generator
+// The main HIJING options are accessable for the user through this interface.
+// andreas.morsch@cern.ch
+
+#include "AliGenMC.h"
+#include <TString.h>
+#include <TArrayI.h>
+
+class THijing;
+class TArrayI;
+class TParticle;
+class TClonesArray;
+class TGraph;
+
+class AliGenHijing : public AliGenMC
+{
+    enum {kNoTrigger, kHardProcesses, kDirectPhotons};
+
+ public:
+    AliGenHijing();
+    AliGenHijing(Int_t npart);
+    AliGenHijing(const AliGenHijing &Hijing);
+    virtual ~AliGenHijing();
+    virtual void    Generate();
+    virtual void    Init();
+    // set centre of mass energy
+    virtual void    SetEnergyCMS(Float_t energy=5500) {fEnergyCMS=energy;}
+    virtual void    SetReferenceFrame(TString frame="CMS")
+       {fFrame=frame;}
+    virtual void    SetProjectile(TString proj="A", Int_t a=208, Int_t z=82)
+       {fProjectile = proj; fAProjectile = a; fZProjectile = z;}    
+    virtual void    SetTarget(TString tar="A", Int_t a=208, Int_t z=82)
+       {fTarget = tar; fATarget = a; fZTarget = z;}    
+    virtual void    SetImpactParameterRange(Float_t bmin = 0, Float_t bmax = 15.)
+       {fMinImpactParam=bmin; fMaxImpactParam=bmax;}
+    virtual void    KeepFullEvent();
+    virtual void    SetJetQuenching(Int_t flag=1)     {fQuench     = flag;}
+    virtual void    SetShadowing(Int_t flag=1)        {fShadowing  = flag;}
+    virtual void    SetDecaysOff(Int_t flag=1)        {fDecaysOff  = flag;}
+    virtual void    SetTrigger(Int_t flag=kNoTrigger) {fTrigger    = flag;}
+    virtual void    SetFlavor(Int_t flag=0)           {fFlavor     = flag;}
+    virtual void    SetEvaluate(Int_t flag=0)         {fEvaluate   = flag;}
+    virtual void    SetSelectAll(Int_t flag=0)        {fSelectAll  = flag;}
+    virtual void    SetRadiation(Int_t flag=3)        {fRadiation  = flag;}    
+    virtual void    SetSpectators(Int_t spects=1)     {fSpectators = spects;}
+    virtual void    SetPtJet(Float_t ptmin)           {fPtMinJet   = ptmin;}
+    virtual void    SetSimpleJets(Int_t flag=0)       {fSimpleJet  = flag;}
+    virtual void    SetNoGammas(Int_t flag=0)         {fNoGammas   = flag;}
+           
+    virtual void    SetJetEtaRange(Float_t etamin = -20., Float_t etamax = 20.)
+       {fEtaMinJet = etamin; fEtaMaxJet = etamax;}
+    virtual void    SetJetPhiRange(Float_t phimin = -180., Float_t phimax = 180.)
+       {fPhiMinJet = TMath::Pi()*phimin/180.; fPhiMaxJet = TMath::Pi()*phimax/180.;}
+    virtual void    SetBoostLHC(Int_t flag = 0)         {fLHC        = flag;}
+// Getters
+    virtual Float_t GetEnergyCMS() {return fEnergyCMS;}
+    virtual TString GetReferenceFrame() {return fFrame;}
+    virtual void    GetProjectile(TString& tar, Int_t& a, Int_t& z)
+       {tar = fProjectile; a = fAProjectile; z = fZProjectile;}    
+    virtual void    GetTarget(TString& tar, Int_t& a, Int_t& z)
+       {tar = fTarget; a = fATarget; z = fZTarget;}    
+    virtual void    GetImpactParameterRange(Float_t& bmin, Float_t& bmax)
+       {bmin = fMinImpactParam; bmax = fMaxImpactParam;}
+    virtual Int_t   GetJetQuenching()                    {return fQuench;}
+    virtual Int_t   GetShadowing()                       {return fShadowing;}
+    virtual Int_t   GetTrigger(Int_t flag=kNoTrigger)    {return fTrigger;}
+    virtual Int_t   GetFlavor(Int_t flag=0)              {return fFlavor;}
+    virtual Int_t   GetRadiation(Int_t flag=3)           {return fRadiation;}    
+    virtual Int_t   GetSpectators(Int_t spects=1)        {return fSpectators;}
+    virtual Float_t GetPtJet(Float_t ptmin)              {return fPtMinJet;}
+    virtual void    GetJetEtaRange(Float_t& etamin, Float_t& etamax)
+       {etamin = fEtaMinJet; etamax = fEtaMaxJet;}
+    virtual void    GetJetPhiRange(Float_t& phimin, Float_t& phimax)
+       {phimin = fPhiMinJet*180./TMath::Pi(); phimax = fPhiMaxJet*180./TMath::Pi();}
+    
+
+// Physics Routines
+    virtual Bool_t ProvidesCollisionGeometry() {return kTRUE;}
+    virtual AliCollisionGeometry* CollisionGeometry() {return fCollisionGeometry;}
+    virtual void EvaluateCrossSections();
+    virtual void Boost();
+    virtual TGraph* CrossSection()     {return fDsigmaDb;}
+    virtual TGraph* BinaryCollisions() {return fDnDb;}
+    virtual Bool_t  CheckTrigger();
+//
+    AliGenHijing &  operator=(const AliGenHijing & rhs);
+ protected:
+    Bool_t SelectFlavor(Int_t pid);
+    void   MakeHeader();
+
+ protected:
+    TString     fFrame;         // Reference frame 
+    TString     fProjectile;    // Projectile
+    TString     fTarget;        // Target
+    Int_t       fAProjectile;    // Projectile A
+    Int_t       fZProjectile;    // Projectile Z
+    Int_t       fATarget;        // Target A
+    Int_t       fZTarget;        // Target Z
+    Float_t     fMinImpactParam; // minimum impact parameter
+    Float_t     fMaxImpactParam; // maximum impact parameter   
+    Int_t       fKeep;           // Flag to keep full event information
+    Int_t       fQuench;         // Flag to switch on jet quenching
+    Int_t       fShadowing;      // Flag to switch on nuclear effects on parton distribution function
+    Int_t       fDecaysOff;      // Flag to turn off decays of pi0, K_s, D, Lambda, sigma
+    Int_t       fTrigger;        // Trigger type
+    Int_t       fEvaluate;       // Evaluate total and partial cross-sections
+    Int_t       fSelectAll;      // Flag to write the full event
+    Int_t       fFlavor;         // Selected particle flavor 4: charm+beauty 5: beauty
+    Float_t     fEnergyCMS;      // Centre of mass energy
+    Float_t     fKineBias;       // Bias from kinematic selection
+    Int_t       fTrials;         // Number of trials
+    TArrayI     fParentSelect;   // Parent particles to be selected 
+    TArrayI     fChildSelect;    // Decay products to be selected
+    Float_t     fXsection;       // Cross-section
+    THijing    *fHijing;         // Hijing
+    Float_t     fPtHardMin;      // lower pT-hard cut 
+    Float_t     fPtHardMax;      // higher pT-hard cut
+    Int_t       fSpectators;     // put spectators on stack
+    TGraph*     fDsigmaDb;       // dSigma/db for the system
+    TGraph*     fDnDb;           // dNBinaryCollisions/db
+    Float_t     fPtMinJet;       // Minimum Pt of triggered Jet
+    Float_t     fEtaMinJet;      // Minimum eta of triggered Jet
+    Float_t     fEtaMaxJet;      // Maximum eta of triggered Jet
+    Float_t     fPhiMinJet;      // At least one of triggered Jets must be in this
+    Float_t     fPhiMaxJet;      // phi range
+    Int_t       fRadiation;      // Flag to switch on/off initial and final state radiation
+    Int_t       fSimpleJet;      // Flag to produce simple tiggered jet topology
+    Int_t       fNoGammas;       // Don't write gammas if flag "on"
+    
+// ZDC proposal (by Chiara) to store num. of SPECTATORS protons and neutrons
+    Int_t      fSpecn;          // Num. of spectator neutrons
+    Int_t      fSpecp;          // Num. of spectator protons
+    Int_t       fLHC;            // Assume LHC as lab frame
+    TClonesArray* fParticles;    // Particle List
+    
+ private:
+    // adjust the weight from kinematic cuts
+    void   AdjustWeights();
+    // check seleted daughters
+    Bool_t DaughtersSelection(TParticle* iparticle);
+    // check if stable
+    Bool_t Stable(TParticle*  particle);
+    
+    ClassDef(AliGenHijing,4) // AliGenerator interface to Hijing
+};
+#endif
+
+
+
+
+
diff --git a/THijing/AliGenHijingEventHeader.cxx b/THijing/AliGenHijingEventHeader.cxx
new file mode 100644 (file)
index 0000000..45497d3
--- /dev/null
@@ -0,0 +1,26 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*
+$Log$
+Revision 1.1  2000/06/15 15:47:48  morsch
+Proposal for an event header class for generated events.
+
+*/
+
+#include "AliGenHijingEventHeader.h"
+ClassImp(AliGenHijingEventHeader)
+
+
diff --git a/THijing/AliGenHijingEventHeader.h b/THijing/AliGenHijingEventHeader.h
new file mode 100644 (file)
index 0000000..0005471
--- /dev/null
@@ -0,0 +1,48 @@
+#ifndef ALIGENHIJINGEVENTHEADER_H
+#define ALIGENHIJINGEVENTHEADER_H
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+#include <TLorentzVector.h>
+
+#include "AliGenEventHeader.h"
+#include "AliCollisionGeometry.h"
+
+class AliGenHijingEventHeader : public AliGenEventHeader, public AliCollisionGeometry
+{
+ public:
+
+  AliGenHijingEventHeader(const char* name){;}
+  AliGenHijingEventHeader(){;}
+  virtual ~AliGenHijingEventHeader() {}
+  // Getters
+  Float_t TotalEnergy()  {return fTotalEnergy;} 
+  Int_t   Trials() {return fTrials;}
+  
+         
+  // Setters
+  void SetTotalEnergy(Float_t energy)  {fTotalEnergy=energy;}
+  void SetJets(TLorentzVector* jet1, TLorentzVector* jet2,
+              TLorentzVector* jet3, TLorentzVector* jet4)
+      {fJet1 = *jet1; fJet2 = *jet2; fJetFsr1 = *jet3; fJetFsr2 = *jet4;}
+  void GetJets(TLorentzVector& jet1, TLorentzVector& jet2,
+              TLorentzVector& jet3, TLorentzVector& jet4)  
+      {jet1 = fJet1; jet2 = fJet2; jet3 = fJetFsr1; jet4 = fJetFsr2;}
+  void SetTrials(Int_t trials) {fTrials = trials;}
+         
+protected:
+  Float_t fTotalEnergy;              // Total energy of produced particles
+  Int_t   fTrials;                   // Number of trials to fulfill trigger condition
+  
+  TLorentzVector  fJet1;             // 4-Momentum-Vector of first   triggered jet  
+  TLorentzVector  fJet2;             // 4-Momentum-Vector of second  triggered jet     
+  TLorentzVector  fJetFsr1;          // 4-Momentum-Vector of first   triggered jet  
+  TLorentzVector  fJetFsr2;          // 4-Momentum-Vector of second  triggered jet     
+  
+  ClassDef(AliGenHijingEventHeader,5) // Event header for hijing event
+};
+
+#endif
index 19c387c..2be7627 100644 (file)
@@ -5,5 +5,6 @@
 #pragma link off all functions;
  
 #pragma link C++ class THijing+;
-
+#pragma link C++ class AliGenHijing+;
+#pragma link C++ class AliGenHijingEventHeader+;
 #endif
index 4413491..b6bb2a8 100644 (file)
@@ -1,7 +1,7 @@
-SRCS= THijing.cxx
+SRCS= THijing.cxx AliGenHijing.cxx AliGenHijingEventHeader.cxx
 
-HDRS= THijing.h
+HDRS= $(SRCS:.cxx=.h) 
 
 DHDR:=THijingLinkDef.h
 
-EXPORT:=THijing.h
\ No newline at end of file
+EXPORT:=THijing.h  AliGenHijingEventHeader.h
\ No newline at end of file
diff --git a/TMEVSIM/AliGenMevSim.cxx b/TMEVSIM/AliGenMevSim.cxx
new file mode 100644 (file)
index 0000000..8f5fe4d
--- /dev/null
@@ -0,0 +1,197 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*
+$Log$
+Revision 1.3  2001/07/27 17:24:28  morsch
+Some bugs corrected. ( Piotr Skowronski)
+
+Revision 1.2  2001/03/28 07:31:48  hristov
+Loop variables declared only once (HP,Sun)
+
+Revision 1.1  2001/03/24 10:04:44  morsch
+MevSim interfaced through AliGenerator, first commit (Sylwester Radomski et al.)
+//Piotr Skowronski Line 104: fThetaMin-->fThetaMax
+*/
+
+//
+// Wrapper for MEVSIM generator.
+// It is using TMevSim to comunicate with fortarn code
+// 
+//     
+// Sylwester Radomski <radomski@if.pw.edu.pl>
+//
+
+//#include "TSystem.h"
+//#include "TUnixSystem.h"
+#include "TParticle.h"
+#include "TMevSim.h"
+
+#include "AliGenMevSim.h"
+#include "AliRun.h"
+
+
+ClassImp(AliGenMevSim)
+
+//____________________________________________________________________________
+AliGenMevSim::AliGenMevSim() : AliGenerator(-1) 
+{
+  //
+  // Standard creator
+  //
+  
+  fConfig = new AliMevSimConfig();
+  fMevSim = new TMevSim();
+  sRandom = fRandom;
+  
+}
+//____________________________________________________________________________
+AliGenMevSim::AliGenMevSim(AliMevSimConfig *config): AliGenerator(-1) {
+
+  fConfig = config;
+  fMevSim = new TMevSim(); 
+  sRandom = fRandom;
+}
+
+//____________________________________________________________________________
+AliGenMevSim::~AliGenMevSim() 
+{
+  //
+  // Standard destructor
+  //
+  if (fMevSim) delete fMevSim;
+}
+//____________________________________________________________________________
+void AliGenMevSim::SetConfig(AliMevSimConfig *config) {
+  
+  fConfig = config;
+}
+//____________________________________________________________________________
+void AliGenMevSim::AddParticleType(AliMevSimParticle *type) {
+
+  fMevSim->AddPartTypeParams((TMevSimPartTypeParams*)type);
+}
+//____________________________________________________________________________
+void AliGenMevSim::Init() 
+{
+  //
+  // Generator initialisation method
+  //
+
+  // fill data from AliMevSimConfig;
+
+  TMevSim *mevsim = fMevSim;
+
+  // geometry & momentum cut
+
+  if (TestBit(kPtRange))  mevsim->SetPtCutRange(fPtMin, fPtMax);
+  
+  if (TestBit(kPhiRange)) // from radians to 0 - 360 deg
+    mevsim->SetPhiCutRange( fPhiMin * 180 / TMath::Pi() , fPhiMax * 180 / TMath::Pi() );
+  
+  if (TestBit(kThetaRange)) // from theta to eta
+  {
+    mevsim->SetEtaCutRange( -TMath::Log( TMath::Tan(fThetaMax/2)) , -TMath::Log( TMath::Tan(fThetaMin/2)) );
+  }  
+
+  // mevsim specyfic parameters
+  
+  mevsim->SetModelType(fConfig->fModelType);
+  mevsim->SetReacPlaneCntrl(fConfig->fReacPlaneCntrl);
+  mevsim->SetPsiRParams(fConfig->fPsiRMean, fConfig->fPsiRStDev);
+  mevsim->SetMultFacParams(fConfig->fMultFacMean, fConfig->fMultFacStDev);
+  
+  // mevsim initialisation
+
+  mevsim->Initialize();
+}
+
+//____________________________________________________________________________
+void AliGenMevSim::Generate() 
+{
+  //
+  // Read the formatted output file and load the particles
+  // Temporary solution
+  //
+
+  Int_t i;
+
+  PDG_t pdg;
+  Float_t orgin[3] = {0,0,0};
+  Float_t polar[3] = {0,0,0};
+  Float_t p[3] = {1,1,1};
+  Float_t time = 0;
+  
+  const Int_t parent = -1;
+  Int_t id;
+
+  // vertexing 
+
+  VertexInternal();
+
+  orgin[0] = fVertex[0];
+  orgin[1] = fVertex[1];
+  orgin[2] = fVertex[2];
+
+  cout << "Vertex ";
+  for (i =0; i<3; i++)
+    cout << orgin[i] << "\t";
+  cout << endl;
+
+  Int_t nParticles = 0;
+
+  TClonesArray *particles = new TClonesArray("TParticle");
+  TParticle *particle;
+
+  fMevSim->GenerateEvent();
+  fNpart= fMevSim->ImportParticles(particles,"");
+
+  cout << "Found " << fNpart << " particles ..." << endl;
+  nParticles = fNpart;
+
+  for (i=0; i<nParticles; i++) {
+    
+    particle = (TParticle*) (*particles)[i];
+
+    pdg = (PDG_t) particle->GetPdgCode();
+    p[0] = particle->Px();
+    p[1] = particle->Py();
+    p[2] = particle->Pz();
+    
+    SetTrack(fTrackIt, parent, pdg, p, orgin, polar, time, kPPrimary, id);
+
+  }  
+  particles->Clear();
+  if (particles) delete particles;
+}
+//____________________________________________________________________________
+//____________________________________________________________________________
+
+
+#ifndef WIN32
+# define ran ran_
+# define type_of_call
+#else
+# define ran RAN
+# define type_of_call _stdcall
+#endif
+
+extern "C" Float_t type_of_call ran(Int_t &)
+{
+  return sRandom->Rndm(); 
+}
+
+//____________________________________________________________________________
diff --git a/TMEVSIM/AliGenMevSim.h b/TMEVSIM/AliGenMevSim.h
new file mode 100644 (file)
index 0000000..e5fdfc5
--- /dev/null
@@ -0,0 +1,44 @@
+#ifndef ALIGENMEVSIM_H
+#define ALIGENMEVSIM_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+// Implementation of the interface for TMevsim
+// Author: Sylwester Radomski <radomski@if.pw.edu.pl>
+
+
+#include "TString.h" 
+
+#include "AliGenerator.h"
+
+#include "AliMevSimConfig.h"
+#include "AliMevSimParticle.h"
+
+class TMevSim;
+
+class AliGenMevSim : public AliGenerator { 
+
+public:
+
+  AliGenMevSim();
+  AliGenMevSim(AliMevSimConfig *config);
+
+  virtual ~AliGenMevSim();
+
+  // 
+  virtual void SetConfig(AliMevSimConfig *config);
+  virtual void AddParticleType(AliMevSimParticle *type);
+
+  virtual void Init();
+  virtual void Generate();
+protected:
+  TMevSim * fMevSim;
+  AliMevSimConfig *fConfig;
+public:
+      
+  ClassDef(AliGenMevSim,1) // Interface class for AliMevsim
+    
+};
+#endif
diff --git a/TMEVSIM/AliMevSimConfig.cxx b/TMEVSIM/AliMevSimConfig.cxx
new file mode 100644 (file)
index 0000000..20879bd
--- /dev/null
@@ -0,0 +1,121 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*
+$Log$
+Revision 1.2  2001/03/24 10:08:30  morsch
+Log tag and copyright notice added.
+
+*/
+
+#include "AliMevSimConfig.h"
+
+ClassImp(AliMevSimConfig)
+
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+AliMevSimConfig::AliMevSimConfig() {
+
+  Init();
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+AliMevSimConfig::AliMevSimConfig(Int_t modelType) {
+
+  Init();
+  SetModelType(modelType);
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+AliMevSimConfig::~AliMevSimConfig() {
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+void AliMevSimConfig::Init() {
+
+  // Default Values
+
+  fModelType = 1;
+  fReacPlaneCntrl = 4;
+  fPsiRMean = fPsiRStDev = 0;
+
+  fMultFacMean  = 1.0;
+  fMultFacStDev = 0.0;
+
+  fNStDevMult = fNStDevTemp = fNStDevSigma = 3.0;
+  fNStDevExpVel = fNStdDevPSIr = fNStDevVn = fNStDevMultFac = 3.0;
+  
+  fNIntegPts = fNScanPts = 100;
+
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+void AliMevSimConfig::SetModelType(Int_t modelType) {
+
+  if (modelType < 0 || modelType > kMAX_MODEL)
+    Error("SetModelType","Wrog Model Type indentifier (%d)",modelType);
+
+  fModelType = modelType;
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+void AliMevSimConfig::SetRectPlane(Int_t ctrl, Float_t psiRMean, Float_t psiRStDev) {
+
+  if (ctrl < 0 || ctrl > kMAX_CTRL)
+    Error("SetReactPlane","Wrong Control Parameter (%d)", ctrl);
+
+  fReacPlaneCntrl = ctrl;
+  fPsiRMean = psiRMean;
+  fPsiRStDev = psiRStDev;
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+void AliMevSimConfig::SetMultFac(Float_t mean, Float_t stDev) {
+  
+  fMultFacMean = mean;
+  fMultFacStDev = stDev;
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+void AliMevSimConfig::SetStDev(Float_t mult, Float_t temp, Float_t sigma,
+  Float_t expVel, Float_t psiR, Float_t Vn, Float_t multFac) {
+
+  fNStDevMult = mult;
+  fNStDevTemp = temp;
+  fNStDevSigma = sigma;
+  fNStDevExpVel = expVel;
+  fNStdDevPSIr = psiR;
+  fNStDevVn = Vn;
+  fNStDevMultFac =multFac;
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
+void AliMevSimConfig::SetGrid(Int_t integr, Int_t scan) {
+
+  fNIntegPts = integr;
+  fNScanPts = scan;
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
diff --git a/TMEVSIM/AliMevSimConfig.h b/TMEVSIM/AliMevSimConfig.h
new file mode 100644 (file)
index 0000000..de01b2d
--- /dev/null
@@ -0,0 +1,59 @@
+#ifndef ALIMEVSIMCONFIG_H
+#define ALIMEVSIMCONFIG_H
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+#include "TObject.h"
+
+class AliGenMevSim;
+
+class AliMevSimConfig : public TObject {
+
+  friend class AliGenMevSim;
+
+ protected:
+
+  static const Int_t kMAX_MODEL = 4;
+  static const Int_t kMAX_CTRL = 4;
+
+  Int_t fModelType;
+
+  Int_t fReacPlaneCntrl;
+  Float_t fPsiRMean, fPsiRStDev;
+
+  Float_t fMultFacMean, fMultFacStDev;
+
+  Float_t fNStDevMult, fNStDevTemp, fNStDevSigma, fNStDevExpVel, fNStdDevPSIr, fNStDevVn, fNStDevMultFac;
+
+  Int_t fNIntegPts;
+  Int_t fNScanPts;
+
+  void Init();
+
+ public:
+
+  AliMevSimConfig();
+  AliMevSimConfig(Int_t modelType);
+
+  ~AliMevSimConfig();
+
+  void SetModelType(Int_t modelType);
+
+  void SetRectPlane(Int_t ctrl, Float_t psiRMean = 0, Float_t psiRStDev = 0);
+  void SetMultFac(Float_t mean, Float_t stDev);
+
+  void SetStDev(Float_t mult, Float_t temp, Float_t sigma,
+               Float_t expVel, Float_t psiR, Float_t Vn, Float_t multFac);
+
+  void SetGrid(Int_t integr, Int_t scan);
+
+
+  ClassDef(AliMevSimConfig,1)
+
+};
+
+
+#endif
diff --git a/TMEVSIM/AliMevSimParticle.cxx b/TMEVSIM/AliMevSimParticle.cxx
new file mode 100644 (file)
index 0000000..a924152
--- /dev/null
@@ -0,0 +1,81 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*
+$Log$
+Revision 1.3  2001/10/21 18:35:56  hristov
+Several pointers were set to zero in the default constructors to avoid memory management problems
+
+Revision 1.2  2001/03/24 10:08:30  morsch
+Log tag and copyright notice added.
+
+*/
+
+#include "AliMevSimParticle.h"
+
+
+ClassImp(AliMevSimParticle)
+
+   
+///////////////////////////////////////////////////////////////////////////////////////
+
+AliMevSimParticle::AliMevSimParticle()
+  : TMevSimPartTypeParams() {
+  fConv = 0;
+}
+
+///////////////////////////////////////////////////////////////////////////////////////
+
+AliMevSimParticle::AliMevSimParticle(PDG_t pdg, Int_t multmean, Int_t multvc, 
+                 Float_t tempmean, Float_t tempstdev, Float_t sigmamean,
+                 Float_t sigmastdev, Float_t expvelmean, Float_t expvelstdev)
+
+  : TMevSimPartTypeParams(0, multmean, multvc, tempmean, tempstdev, 
+                         sigmamean, sigmastdev, expvelmean, expvelstdev)  {
+
+
+  // Calculate geant ID from pdg
+  fConv = new TMevSimConverter();
+  fPdg = pdg;
+  if (fConv) fGPid = fConv->IdFromPDG(pdg);  
+
+}
+
+///////////////////////////////////////////////////////////////////////////////////////
+
+AliMevSimParticle::~AliMevSimParticle() {
+}
+
+///////////////////////////////////////////////////////////////////////////////////////
+
+void  AliMevSimParticle::SetPDG(PDG_t pdg) {
+
+  fPdg = pdg;
+  fGPid = fConv->IdFromPDG(pdg);
+}
+
+///////////////////////////////////////////////////////////////////////////////////////
+
+PDG_t AliMevSimParticle::GetPDG() {
+  
+  fPdg = (PDG_t)fConv->PDGFromId(fGPid);
+  return fPdg;
+}
+
+///////////////////////////////////////////////////////////////////////////////////////
+
+
+  
+
diff --git a/TMEVSIM/AliMevSimParticle.h b/TMEVSIM/AliMevSimParticle.h
new file mode 100644 (file)
index 0000000..40dd522
--- /dev/null
@@ -0,0 +1,47 @@
+
+#ifndef ALIMEVSIMPARTICLE_H
+#define ALIMEVSIMPARTICLE_H
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+#include <TPDGCode.h>
+
+#include "TMevSimConverter.h"
+#include "TMevSimPartTypeParams.h"
+
+class AliMevSimParticle :public TMevSimPartTypeParams {
+
+ protected:
+  
+  PDG_t fPdg;
+  TMevSimConverter *fConv;
+
+ public:
+  
+  ///////////////////////////////////////////////////////////////////////////////////////
+  
+  AliMevSimParticle();
+  AliMevSimParticle(PDG_t pdg, Int_t multmean, Int_t multvc, 
+                   Float_t tempmean, Float_t tempstdev, Float_t sigmamean,
+                   Float_t sigmastdev, Float_t expvelmean, Float_t expvelstdev);
+  
+  virtual ~AliMevSimParticle();
+  
+  ///////////////////////////////////////////////////////////////////////////////////////
+  
+  virtual void        SetPDG(PDG_t pdg);
+  virtual PDG_t       GetPDG();
+  
+  ///////////////////////////////////////////////////////////////////////////////////////
+  
+  ClassDef(AliMevSimParticle,1)  
+    
+  ///////////////////////////////////////////////////////////////////////////////////////
+    
+};
+
+#endif
+
index d02119e..fa51f51 100644 (file)
@@ -7,5 +7,7 @@
 #pragma link C++ class TMevSim+;
 #pragma link C++ class TMevSimPartTypeParams+;
 #pragma link C++ class TMevSimConverter+;
-
+#pragma link C++ class AliGenMevSim+;
+#pragma link C++ class AliMevSimConfig+;
+#pragma link C++ class AliMevSimParticle+;
 #endif
index bbb8466..03aaf86 100644 (file)
@@ -1,9 +1,13 @@
-SRCS= TMevSim.cxx TMevSimPartTypeParams.cxx TMevSimConverter.cxx
+SRCS= TMevSim.cxx TMevSimPartTypeParams.cxx TMevSimConverter.cxx \
+AliGenMevSim.cxx AliMevSimConfig.cxx AliMevSimParticle.cxx
 
 HDRS= \
 TMevSim.h \
 TMevSimPartTypeParams.h \
-TMevSimConverter.h
+TMevSimConverter.h \
+AliGenMevSim.h \
+AliMevSimConfig.h \
+AliMevSimParticle.h
 
 DHDR=TMevSimLinkDef.h