]>
Commit | Line | Data |
---|---|---|
e3b21c34 | 1 | /**************************************************************************\r |
2 | * Copyright(c) 2009-2010, ALICE Experiment at CERN, All rights reserved. *\r | |
3 | * *\r | |
4 | * Author: The ALICE Off-line Project. *\r | |
5 | * Contributors are mentioned in the code where appropriate. *\r | |
6 | * *\r | |
7 | * Permission to use, copy, modify and distribute this software and its *\r | |
8 | * documentation strictly for non-commercial purposes is hereby granted *\r | |
9 | * without fee, provided that the above copyright notice appears in all *\r | |
10 | * copies and that both the copyright notice and this permission notice *\r | |
11 | * appear in the supporting documentation. The authors make no claims *\r | |
12 | * about the suitability of this software for any purpose. It is *\r | |
13 | * provided "as is" without express or implied warranty. *\r | |
14 | **************************************************************************/\r | |
15 | \r | |
16 | //\r | |
17 | // A very simple model based on the article "Physics Letters B325 (1994) 7-12".\r | |
18 | // The only addition is to model the freeze-out at which the\r | |
19 | // light nuclei can form for the energies of the LHC, since PYTHIA does not give\r | |
20 | // any detailed spatial description of the generated particles at the level of a\r | |
21 | // few fermis.\r | |
22 | //\r | |
23 | // There are two options to place the nucleons: a thermal picture where all\r | |
24 | // nucleons are placed randomly and homogeneously in an spherical volume of a\r | |
25 | // few fermis, and expansion one where they are projected on its surface.\r | |
26 | //\r | |
27 | // A (anti)deuteron will form if there is a pair of (anti)proton-(anti)neutron\r | |
28 | // with momentum difference less than ~ 300MeV and relative distance less than a\r | |
29 | // ~ 2.1fm. Only 3/4 of this clusters are accepted due to the deuteron spin.\r | |
30 | //\r | |
31 | // Author: Eulogio Serradilla <eulogio.serradilla@ciemat.es>,\r | |
32 | // Arturo Menchaca <menchaca@fisica.unam.mx>\r | |
33 | //\r | |
34 | \r | |
539d3694 | 35 | #include "Riostream.h"\r |
36 | #include "TParticle.h"\r | |
37 | #include "AliStack.h"\r | |
38 | #include "AliGenDeuteron.h"\r | |
39 | #include "AliGenCocktailAfterBurner.h"\r | |
40 | #include "TMath.h"\r | |
41 | #include "TMCProcess.h"\r | |
42 | #include "TList.h"\r | |
43 | #include "TVector3.h"\r | |
539d3694 | 44 | #include "AliMC.h"\r |
4a33c50d | 45 | #include "AliRun.h"\r |
e3b21c34 | 46 | #include "TArrayF.h"\r |
47 | #include "AliGenCocktailEventHeader.h"\r | |
539d3694 | 48 | \r |
49 | ClassImp(AliGenDeuteron)\r | |
50 | \r | |
539d3694 | 51 | AliGenDeuteron::AliGenDeuteron(Int_t sign, Double_t pmax, Double_t rmax, Double_t rsrc, Int_t model):\r |
52 | fDeuteronMass(1.87561282), //pdg\r | |
53 | fPmax(pmax), // GeV/c\r | |
54 | fRmax(rmax*1.e-13), //cm\r | |
55 | fRsrc(rsrc*1.e-13), // cm\r | |
56 | fSpinProb(0.75),\r | |
57 | fModel(model),\r | |
58 | fSign(1)\r | |
59 | {\r | |
60 | fSign = sign > 0 ? 1:-1;\r | |
61 | }\r | |
62 | \r | |
63 | AliGenDeuteron::~AliGenDeuteron()\r | |
64 | {\r | |
65 | }\r | |
66 | \r | |
67 | void AliGenDeuteron::Init()\r | |
68 | {\r | |
69 | // Standard AliGenerator Initializer\r | |
70 | }\r | |
71 | \r | |
72 | void AliGenDeuteron::Generate()\r | |
73 | {\r | |
74 | // Modify the stack of each event, adding the new particles at the end and\r | |
75 | // not transporting the parents.\r | |
76 | \r | |
77 | Info("Generate","freeze-out model : %d (0 expansion, 1 thermal)",fModel);\r | |
78 | Info("Generate","relative momentum: %g GeV/c",fPmax);\r | |
79 | Info("Generate","relative distance: %g cm",fRmax);\r | |
80 | Info("Generate","source radius : %g cm",fRsrc);\r | |
81 | Info("Generate","spin probability : %g ",fSpinProb);\r | |
82 | Info("Generate","sign : %d ",fSign);\r | |
83 | \r | |
539d3694 | 84 | Int_t ntr;\r |
85 | \r | |
86 | // Get the cocktail generator\r | |
87 | AliGenCocktailAfterBurner* gener = (AliGenCocktailAfterBurner*)gAlice->GetMCApp()->Generator();\r | |
88 | \r | |
89 | // Loop over events\r | |
90 | for(Int_t ns = 0; ns < gener->GetNumberOfEvents(); ++ns)\r | |
91 | {\r | |
92 | gener->SetActiveEventNumber(ns);\r | |
93 | \r | |
94 | AliStack* stack = gener->GetStack(ns); // Stack of event ns\r | |
95 | if(!stack)\r | |
96 | {\r | |
97 | Info("Generate", "no stack, exiting");\r | |
98 | return;\r | |
99 | }\r | |
100 | \r | |
101 | // find emerging nucleons from the collision of current event\r | |
102 | \r | |
103 | TList* protons = new TList();\r | |
104 | protons->SetOwner(kFALSE);\r | |
105 | TList* neutrons = new TList();\r | |
106 | neutrons->SetOwner(kFALSE);\r | |
107 | \r | |
e3b21c34 | 108 | // primary vertex of current event\r |
109 | TArrayF primVtx;\r | |
110 | gener->GetActiveEventHeader()->PrimaryVertex(primVtx);\r | |
111 | TVector3 r0(primVtx[0],primVtx[1],primVtx[2]);\r | |
539d3694 | 112 | \r |
113 | for (Int_t i=0; i<stack->GetNtrack(); ++i)\r | |
114 | {\r | |
115 | TParticle* iParticle = stack->Particle(i);\r | |
116 | \r | |
117 | if(iParticle->TestBit(kDoneBit)) continue;\r | |
118 | // select only nucleons within the freeze-out volume\r | |
119 | TVector3 r(iParticle->Vx(),iParticle->Vy(),iParticle->Vz());\r | |
120 | if((r-r0).Mag()>fRsrc) continue;\r | |
121 | \r | |
122 | Int_t pdgCode = iParticle->GetPdgCode();\r | |
123 | if(pdgCode == fSign*2212)// (anti)proton\r | |
124 | {\r | |
e3b21c34 | 125 | FixProductionVertex(iParticle);\r |
539d3694 | 126 | protons->Add(iParticle);\r |
127 | }\r | |
128 | else if(pdgCode == fSign*2112) // (anti)neutron\r | |
129 | {\r | |
e3b21c34 | 130 | FixProductionVertex(iParticle);\r |
539d3694 | 131 | neutrons->Add(iParticle);\r |
132 | }\r | |
133 | }\r | |
134 | \r | |
135 | // coalescence conditions\r | |
136 | // A.J. Baltz et al., Phys. lett B 325(1994)7\r | |
137 | \r | |
138 | TIter p_next(protons);\r | |
139 | while(TParticle* n0 = (TParticle*) p_next())\r | |
140 | {\r | |
141 | TParticle* iProton = 0;\r | |
142 | TParticle* jNeutron = 0;\r | |
143 | \r | |
144 | // Select the first nucleon\r | |
145 | iProton = n0;\r | |
146 | \r | |
147 | TVector3 v0(n0->Vx(), n0->Vy(), n0->Vz());\r | |
148 | TVector3 p0(n0->Px(), n0->Py(), n0->Pz()), p1(0,0,0);\r | |
149 | \r | |
150 | // See if there is a neibourgh...\r | |
151 | TIter n_next(neutrons);\r | |
152 | while(TParticle* n1 = (TParticle*) n_next() )\r | |
153 | {\r | |
154 | TVector3 v(n1->Vx(), n1->Vy(), n1->Vz());\r | |
155 | TVector3 p(n1->Px(), n1->Py(), n1->Pz());\r | |
156 | \r | |
157 | if((p-p0).Mag() > fPmax) continue;\r | |
158 | if((v-v0).Mag() > fRmax) continue;\r | |
159 | \r | |
160 | jNeutron = n1;\r | |
161 | p1 = p;\r | |
162 | break;\r | |
163 | }\r | |
164 | \r | |
165 | if(jNeutron == 0) continue; // with next proton\r | |
166 | \r | |
e3b21c34 | 167 | if(Rndm() > fSpinProb) continue;\r |
539d3694 | 168 | \r |
169 | // neutron captured!\r | |
170 | \r | |
171 | TVector3 pN = p0+p1;\r | |
172 | TVector3 vN(iProton->Vx(), iProton->Vy(), iProton->Vz()); // production vertex = p = n = collision vertex\r | |
173 | \r | |
174 | // E^2 = p^2 + m^2\r | |
175 | Double_t EN = TMath::Sqrt(pN.Mag2()+fDeuteronMass*fDeuteronMass);\r | |
176 | \r | |
177 | // Add a new (anti)deuteron to the current event stack\r | |
178 | stack->PushTrack(1, stack->Particles()->IndexOf(iProton), fSign*1000010020,\r | |
179 | pN.X(), pN.Y(), pN.Z(), EN,\r | |
180 | vN.X(), vN.Y(), vN.Z(), iProton->T(),\r | |
181 | 0., 0., 0., kPNCapture, ntr,1.,0);\r | |
182 | \r | |
183 | //Info("Generate","neutron capture (NEW DEUTERON)");\r | |
184 | \r | |
185 | // Set bit not to transport for the nucleons\r | |
186 | iProton->SetBit(kDoneBit);\r | |
187 | jNeutron->SetBit(kDoneBit);\r | |
188 | \r | |
189 | // Remove from the list\r | |
190 | neutrons->Remove(jNeutron);\r | |
191 | }\r | |
192 | \r | |
193 | protons->Clear("nodelete");\r | |
194 | neutrons->Clear("nodelete");\r | |
195 | \r | |
196 | delete protons;\r | |
197 | delete neutrons;\r | |
198 | }\r | |
199 | \r | |
200 | Info("Generate","Coalescence afterburner: DONE");\r | |
201 | }\r | |
202 | \r | |
203 | // create the freeze-out nucleon distribution around the collision vertex\r | |
e3b21c34 | 204 | void AliGenDeuteron::FixProductionVertex(TParticle* i)\r |
539d3694 | 205 | {\r |
4a33c50d | 206 | // Fix for the production vertex\r |
539d3694 | 207 | Double_t x,y,z;\r |
208 | \r | |
209 | if(fModel == kThermal) // homogeneous volume\r | |
210 | {\r | |
211 | // random (r,theta,phi)\r | |
e3b21c34 | 212 | Double_t r = fRsrc*TMath::Power(Rndm(),1./3.);\r |
213 | Double_t theta = TMath::ACos(2.*Rndm()-1.);\r | |
214 | Double_t phi = 2.*TMath::Pi()*Rndm();\r | |
539d3694 | 215 | \r |
216 | // transform coordenates\r | |
217 | x = r*TMath::Sin(theta)*TMath::Cos(phi);\r | |
218 | y = r*TMath::Sin(theta)*TMath::Sin(phi);\r | |
219 | z = r*TMath::Cos(theta);\r | |
220 | }\r | |
221 | else // projection into the surface of an sphere\r | |
222 | {\r | |
223 | x = fRsrc*TMath::Sin(i->Theta())*TMath::Cos(i->Phi());\r | |
224 | y = fRsrc*TMath::Sin(i->Theta())*TMath::Sin(i->Phi());\r | |
225 | z = fRsrc*TMath::Cos(i->Theta());\r | |
226 | }\r | |
227 | \r | |
e3b21c34 | 228 | // assume nucleons at the freeze-out have the same production vertex\r |
539d3694 | 229 | i->SetProductionVertex(i->Vx()+x, i->Vy()+y, i->Vz()+z, i->T());\r |
230 | }\r |