Line data Source code
1 : /**************************************************************************
2 : * Copyright(c) 1998-2007, 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 : // Generator for particles according generic functions
17 : // For high pt performance studies realistic distribution of the particles - local density in jets to be approximated
18 : //
19 : //
20 : // TF1 * fF1Momentum; // 1/momentum distribution function inGeV
21 : // TF1 * fFPhi; // phi distribution function in rad - if not set flat 0-2pi used
22 : // TF1 * fFTheta; // theta distribution function in rad - if not set flat pi/4-3pi/4 used
23 : // TF3 * fFPosition; // position distribution function in cm - TO FIX if not set 0 used ()
24 : // TF1 * fFPdg; // pdg distribution function - if not set flat jet pdg used
25 : // We assume that the moment, postion and PDG code of particles are independent
26 : // Only tracks/particle crossing the reference radius at given z range
27 : //
28 : // Origin: marian.ivanov@cern.ch
29 :
30 :
31 : /*
32 : To test generator for particular setting run
33 : AliGenPerformance::TestAliGenPerformance(Int_t nEvents, TF1 *f1pt, TF1 *fpdg){
34 : For distribution of
35 : fF1Momentum=new TF1("f1pt","1-10*x",0,0.1);
36 : AliGenPerformance::TestAliGenPerformance(1000, fF1Momentum,0)
37 : pt distribution of charged primary particles described by powerlaw with slope -1.7
38 :
39 : gSystem->Load("libpythia6");
40 : gSystem->Load("libEGPythia6");
41 : gSystem->Load("liblhapdf");
42 : gSystem->Load("libAliPythia6");
43 : //
44 : AliGenPerformance::TestAliGenPerformance(5000,0);
45 : TFile * f = TFile::Open("testAliGenPerformance.root");
46 : testGener.Draw("pt>>his(100,10,100)","charge!=0&&abs(fKF)<2000","");
47 : f1pt = new TF1("f1pt","1/(0.01+x)",0,0.1);
48 : his->Fit(f1);
49 :
50 : */
51 :
52 :
53 :
54 :
55 : #include <TParticle.h>
56 : #include <TF1.h>
57 : #include <TF3.h>
58 : #include <TDatabasePDG.h>
59 :
60 : #include "AliRun.h"
61 : #include "AliLog.h"
62 : #include "AliESDtrack.h"
63 : #include "AliESDVertex.h"
64 : #include "AliGenPerformance.h"
65 : #include "AliGenEventHeader.h"
66 : #include "TMCParticle.h"
67 : #include <AliPythia.h>
68 : #include <TTreeStream.h>
69 :
70 6 : ClassImp(AliGenPerformance)
71 :
72 : //-----------------------------------------------------------------------------
73 : AliGenPerformance::AliGenPerformance():
74 0 : AliGenerator(),
75 0 : fNJets(1), // mean number of jets per event
76 0 : fF1Momentum(0), // momentum distribution function
77 0 : fFPhi(0), // phi distribution function
78 0 : fFTheta(0), // theta distribution function
79 0 : fFPosition(0), // position distribution function
80 0 : fFPdg(0), // pdg distribution function
81 0 : fTestStream(0) // test stream - used for tuning of parameters of generator
82 0 : {
83 : //
84 : // Default constructor
85 : //
86 0 : SetNumberParticles(1);
87 0 : }
88 :
89 : AliGenPerformance::AliGenPerformance(const AliGenPerformance& func):
90 0 : AliGenerator(),
91 0 : fNJets(func.fNJets), //
92 0 : fF1Momentum(func.fF1Momentum), // momentum distribution function
93 0 : fFPhi(func.fFPhi), // phi distribution function
94 0 : fFTheta(func.fFTheta), // theta distribution function
95 0 : fFPosition(func.fFPosition), // position distribution function
96 0 : fFPdg(func.fFPdg) // pdg distribution function
97 0 : {
98 : // Copy constructor
99 0 : SetNumberParticles(1);
100 0 : }
101 :
102 : AliGenPerformance & AliGenPerformance::operator=(const AliGenPerformance& func)
103 : {
104 : // Assigment operator
105 0 : if(&func == this) return *this;
106 0 : fNJets = func.fNJets;
107 0 : fF1Momentum = func.fF1Momentum;
108 0 : fFPhi = func.fFPhi;
109 0 : fFTheta = func.fFTheta;
110 0 : fFPosition = func.fFPosition;
111 0 : fFPdg = func.fFPdg;
112 0 : return *this;
113 0 : }
114 :
115 :
116 : //-----------------------------------------------------------------------------
117 : void AliGenPerformance::Generate()
118 : {
119 : //
120 : // Generate one muon
121 : //
122 : Int_t naccepted =0;
123 0 : Int_t njets=gRandom->Poisson(fNJets);
124 0 : TDatabasePDG *databasePDG = TDatabasePDG::Instance();
125 :
126 0 : for (Int_t iparton=0; iparton<njets; iparton++){
127 : //
128 : //
129 : //
130 0 : Float_t mom[3];
131 0 : Float_t posf[3];
132 0 : Double_t pos[3];
133 0 : Int_t pdg;
134 0 : Double_t ptot, pt, phi, theta;
135 : //
136 0 : if (fF1Momentum){
137 0 : ptot = 1./fF1Momentum->GetRandom();
138 0 : }else{
139 0 : ptot = 0.001+fF1Momentum->GetRandom()*0.2;
140 0 : ptot/=ptot;
141 : }
142 0 : if (fFPhi){
143 0 : phi = fFPhi->GetRandom();
144 0 : }else{
145 0 : phi = gRandom->Rndm()*TMath::TwoPi();
146 : }
147 0 : if (fFTheta){
148 0 : theta = fFTheta->GetRandom();
149 0 : }else{
150 0 : theta = (gRandom->Rndm()-0.5)*TMath::Pi()*0.5 +TMath::Pi()/2;
151 : }
152 0 : pt = ptot*TMath::Sin(theta);
153 0 : mom[0] = pt*TMath::Cos(phi);
154 0 : mom[1] = pt*TMath::Sin(phi);
155 0 : mom[2] = ptot*TMath::Cos(theta);
156 0 : pos[0]=fOrigin[0]; pos[1]=fOrigin[1]; pos[2]=fOrigin[2];
157 : //
158 0 : if (fFPosition) fFPosition->GetRandom3(pos[0],pos[1],pos[2]); // ????
159 : //
160 0 : posf[0]=pos[0];
161 0 : posf[1]=pos[1];
162 0 : posf[2]=pos[2];
163 0 : if (fFPdg){
164 0 : pdg = TMath::Nint(fFPdg->GetRandom());
165 0 : }else{
166 0 : pdg = 1+TMath::Nint(gRandom->Rndm()*5.);
167 : }
168 0 : Float_t polarization[3]= {0,0,0};
169 0 : Int_t nt;
170 : //
171 0 : AliPythia *py=AliPythia::Instance();
172 0 : py->Py1ent(-1, -pdg, ptot, theta, phi);
173 0 : py->Py1ent( 2, pdg, ptot, theta, phi+TMath::Pi());
174 0 : py->Pyexec();
175 0 : TObjArray * array = py->GetPrimaries();
176 0 : Int_t nParticles=array->GetEntries();
177 : //array->Print();
178 0 : for (Int_t iparticle=2; iparticle<nParticles; iparticle++){
179 0 : TMCParticle * mcParticle= (TMCParticle*)array->At(iparticle);
180 0 : Int_t flavour = mcParticle->GetKF();
181 0 : mom[0]=mcParticle->GetPx();
182 0 : mom[1]=mcParticle->GetPy();
183 0 : mom[2]=mcParticle->GetPz();
184 : //
185 0 : if (!fTestStream) PushTrack(fTrackIt,-1,flavour,mom, posf, polarization,0,kPPrimary,nt);
186 0 : if (fTestStream){
187 0 : TParticlePDG * pdgParticle=databasePDG->GetParticle(mcParticle->GetKF());
188 0 : if (pdgParticle){
189 0 : Double_t charge=pdgParticle->Charge();
190 0 : Double_t mass=pdgParticle->Mass();
191 0 : Double_t pt=TMath::Sqrt(mcParticle->GetPx()*mcParticle->GetPx()+mcParticle->GetPy()*mcParticle->GetPy());
192 0 : (*fTestStream)<<"testGener"<<
193 0 : "njets="<<njets<<
194 0 : "ptot="<<ptot<<
195 0 : "theta="<<theta<<
196 0 : "phi="<<phi<<
197 0 : "pdg="<<pdg<<
198 0 : "charge="<<charge<<
199 0 : "mass="<<mass<<
200 0 : "pt="<<pt<<
201 0 : "nParticles="<<nParticles<<
202 0 : "ipart="<<iparticle<<
203 0 : "mcParticle.="<<mcParticle<<
204 : "\n";
205 0 : }
206 0 : }
207 0 : naccepted++;
208 : }
209 0 : }
210 : // AliGenEventHeader* header = new AliGenEventHeader("THn");
211 : //gAlice->SetGenEventHeader(header);
212 :
213 : return;
214 0 : }
215 :
216 :
217 : //-----------------------------------------------------------------------------
218 : void AliGenPerformance::Init()
219 : {
220 : //
221 : // Initialisation, check consistency of selected ranges
222 : //
223 :
224 :
225 0 : printf("************ AliGenPerformance ****************\n");
226 0 : printf("************************************************\n");
227 0 : if (!fF1Momentum){
228 0 : AliInfo("Momentum distribution function not specified");
229 0 : }
230 0 : if (!fFPhi){
231 0 : AliInfo("phi distribution function not specified");
232 0 : }
233 0 : if (!fFTheta){
234 0 : AliInfo("Theta distribution function not specified");
235 0 : }
236 0 : if (!fFPosition){
237 0 : AliInfo("Position distribution function not specified");
238 0 : }
239 0 : if (!fFPdg){
240 0 : AliInfo("PDG distribution function not specified");
241 0 : }
242 :
243 0 : return;
244 : }
245 :
246 : void AliGenPerformance::SetFunctions(TF1 * momentum, TF1 *fphi, TF1 *ftheta,TF3 * position, TF1* pdg){
247 : //
248 : // Set the function
249 : //
250 0 : fF1Momentum = momentum;
251 0 : fFPhi = fphi;
252 0 : fFTheta = ftheta;
253 0 : fFPosition = position;
254 0 : fFPdg = pdg;
255 0 : }
256 :
257 : void AliGenPerformance::TestAliGenPerformance(Int_t nEvents, TF1 *f1pt, TF1 *fpdg){
258 : //
259 : // test the genrator class - write particle to the tree
260 : //
261 0 : AliGenPerformance *genPerformance= new AliGenPerformance;
262 0 : genPerformance->SetNJets(1);
263 0 : if (!f1pt) f1pt = new TF1("f1pt","1-10*x",0,0.1);
264 0 : if (!fpdg) fpdg = new TF1("f1pt","x",1,6);
265 0 : genPerformance->SetFunctions(f1pt,0,0,0,fpdg);
266 0 : TTreeSRedirector*pcstream = new TTreeSRedirector("testAliGenPerformance.root","recreate");
267 0 : genPerformance->SetStreamer(pcstream);
268 0 : for (Int_t i=0; i<nEvents;i++){
269 0 : genPerformance->Generate();
270 : }
271 0 : delete pcstream;
272 0 : }
273 :
|
