8ce81660635d4526935301b865fd865eb75c104e
[u/mrichter/AliRoot.git] / TIsajet / AliGenIsajet.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 //
19 // Wrapper for ISAJET generator.
20 // It is using ISAJET to comunicate with fortarn code
21 // 
22 //      
23 // 
24 //
25
26 #include <Riostream.h>
27
28 //#include "TSystem.h"
29 //#include "TUnixSystem.h"
30 #include <TParticle.h>
31 #include "TIsajet.h"
32 #include "AliMC.h"
33 #include "AliIsajetRndm.h"
34 #include "AliGenIsajet.h"
35 #include "AliRun.h"
36
37
38
39 ClassImp(AliGenIsajet)
40
41  AliGenIsajet::AliGenIsajet() :
42    AliGenMC(),
43    fKineBias(1),
44    fTrials(0)
45  {  
46   
47   }
48 //____________________________________________________________________________
49 AliGenIsajet::AliGenIsajet(Int_t npart) 
50   :AliGenMC(npart),
51    fKineBias(1),
52    fTrials(0)
53 {
54   
55     fIsajet = new TIsajet();
56   AliIsajetRndm::SetIsajetRandom(GetRandom());
57   
58 }
59
60 AliGenIsajet::AliGenIsajet(const AliGenIsajet &  Isajet) 
61   :AliGenMC(Isajet),
62    fKineBias(1),
63    fTrials(0)
64 {
65   Isajet.Copy(*this);
66   
67 }
68 //____________________________________________________________________________
69 AliGenIsajet::~AliGenIsajet() 
70
71   //
72   // Standard destructor
73   //
74   if (fIsajet) delete fIsajet;
75 }
76
77 //____________________________________________________________________________
78 void AliGenIsajet::Init() 
79 {
80   //
81   // Generator initialisation method
82   //
83  cout<<"isajet4"<<endl;
84    fIsajet->SetECM(2000);
85  cout<<"isajet5"<<endl;
86   fIsajet->SetJobtype("MINBIAS");
87  cout<<"isajet6"<<endl;
88  fIsajet->SetIDIN(0,1);
89  cout<<"isajet7"<<endl;
90   fIsajet->SetIDIN(1,1);
91  cout<<"isajet8"<<endl;
92    fIsajet->Initialise();
93  cout<<"isajet9"<<endl;
94 }
95
96 //____________________________________________________________________________
97 void AliGenIsajet::Generate() 
98 {
99
100   Float_t polar[3] =   {0,0,0};
101   Float_t origin[3]=   {0,0,0};
102   Float_t origin0[3]=  {0,0,0};
103   Float_t p[4], random[6];
104   Int_t j, kf, ks, imo;
105   Int_t nt=0; 
106   Int_t jev=0;
107
108   fTrials=0;
109   
110  static TClonesArray *particles;
111  //cout<<"generate"<<endl;
112  if(!particles) particles=new TClonesArray("TParticle",10000);
113  const Float_t kconv=0.001/2.999792458e8;
114
115   VertexInternal();
116
117   origin[0] = fVertex[0];
118   origin[1] = fVertex[1];
119   origin[2] = fVertex[2];
120   //cout<<"generate1"<<endl;
121
122    while(1)
123           {
124       //cout<<"generate2"<<endl;
125   fIsajet->GenerateEvent();
126   //cout<<"generate3"<<endl;
127   fTrials++;
128   fIsajet->ImportParticles(particles,"All");
129   //cout<<"generate4"<<endl;
130   Int_t np = particles->GetEntriesFast()-1;
131   //cout<<"generate51"<<endl;
132    if (np == 0 ) continue;
133
134        Int_t nc=0;
135        Int_t * newPos = new Int_t[np];
136        //cout<<"generate52"<<endl;
137
138         for (Int_t i = 0; i<np; i++) *(newPos+i)=-1;
139         //cout<<"generate5"<<endl;
140   for (Int_t i = 0; i<np; i++) 
141  {
142      TParticle *  iparticle       = (TParticle *) particles->At(i);
143      //cout<<"generate6"<<endl;
144      imo = iparticle->GetFirstMother();
145      kf        = iparticle->GetPdgCode();
146      ks        = iparticle->GetStatusCode();
147              if (ks != 3 &&
148               KinematicSelection(iparticle,0))
149           {
150               nc++;
151      p[0]=iparticle->Px();
152      p[1]=iparticle->Py();
153      p[2]=iparticle->Pz();
154      p[3]=iparticle->Energy();
155      Float_t tof = kconv*iparticle->T();
156      Int_t   iparent = (imo > -1) ? newPos[imo] : -1;
157      Int_t   trackIt = (ks == 1) && fTrackIt;
158      PushTrack(trackIt, iparent, kf,
159                           p[0], p[1], p[2], p[3],
160                           origin[0], origin[1], origin[2],
161                           tof,
162                           polar[0], polar[1], polar[2],
163                           kPPrimary, nt,1.,ks);
164       KeepTrack(nt);
165         newPos[i]=nt;
166           }
167   }// end of for: particle loop
168       if (newPos) delete[] newPos;
169       // MakeHeader();
170         if (nc > 0) {
171             jev+=nc;
172             if (jev >= fNpart || fNpart == -1) {
173                 fKineBias=Float_t(fNpart)/Float_t(fTrials);
174                 break;
175             }}
176  
177   }
178  SetHighWaterMark(nt);
179 }
180 AliGenIsajet& AliGenIsajet::operator=(const  AliGenIsajet& rhs)
181 {
182 // Assignment operator
183     rhs.Copy(*this);
184     return (*this);
185 }
186
187 //____________________________________________________________________________
188 //____________________________________________________________________________
189