--- /dev/null
+/**************************************************************************
+ * 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());
+ }
+}
+
--- /dev/null
+#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
#pragma link off all functions;
#pragma link C++ class TEpEmGen+;
+#pragma link C++ class AliGenEpEmv1+;
#endif
-SRCS= TEpEmGen.cxx
+SRCS= TEpEmGen.cxx AliGenEpEmv1.cxx
-HDRS= TEpEmGen.h
+HDRS= TEpEmGen.h AliGenEpEmv1.h
DHDR:=TEPEMGENLinkDef.h
--- /dev/null
+/* 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();
+}
+
+/*******************************************************************/
+
+
+/*******************************************************************/
--- /dev/null
+// 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
#pragma link off all functions;
#pragma link C++ class THBTprocessor+;
-
+#pragma link C++ class AliGenHBTprocessor+;
#endif
-SRCS= THBTprocessor.cxx
+SRCS= THBTprocessor.cxx AliGenHBTprocessor.cxx
HDRS= $(SRCS:.cxx=.h)
-DHDR= THBTprocessorLinkDef.h
+DHDR= THBTprocessorLinkDef.h
EXPORT:=HBTprocCOMMON.h THBTprocessor.h
--- /dev/null
+/**************************************************************************
+ * 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;
+ }
+}
--- /dev/null
+#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
+
+
+
+
+
--- /dev/null
+/**************************************************************************
+ * 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)
+
+
--- /dev/null
+#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
#pragma link off all functions;
#pragma link C++ class THijing+;
-
+#pragma link C++ class AliGenHijing+;
+#pragma link C++ class AliGenHijingEventHeader+;
#endif
-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
--- /dev/null
+/**************************************************************************
+ * 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();
+}
+
+//____________________________________________________________________________
--- /dev/null
+#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
--- /dev/null
+/**************************************************************************
+ * 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;
+}
+
+//////////////////////////////////////////////////////////////////////////////////////////////////
+
--- /dev/null
+#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
--- /dev/null
+/**************************************************************************
+ * 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;
+}
+
+///////////////////////////////////////////////////////////////////////////////////////
+
+
+
+
--- /dev/null
+
+#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
+
#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
-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