2dcfed76923703a0d08061f72a643428e940cf21
[u/mrichter/AliRoot.git] / TRD / AliTRDv2.cxx
1 ///////////////////////////////////////////////////////////////////////////////
2 //                                                                           //
3 //  Transition Radiation Detector version 2 -- detailed simulation           //
4 //                                                                           //
5 //Begin_Html
6 /*
7 <img src="gif/AliTRDv2Class.gif">
8 */
9 //End_Html
10 //                                                                           //
11 //                                                                           //
12 ///////////////////////////////////////////////////////////////////////////////
13
14 #include <TMath.h>
15 #include <TRandom.h>
16 #include <TVector.h>
17 #include <TGeometry.h>
18 #include <TNode.h>
19 #include <TPGON.h> 
20
21 #include "GParticle.h"
22 #include "AliTRDv2.h"
23 #include "AliRun.h"
24 #include "AliConst.h"
25 #include "AliMC.h"
26
27 ClassImp(AliTRDv2)
28
29 //_____________________________________________________________________________
30 AliTRDv2::AliTRDv2(const char *name, const char *title) 
31          :AliTRD(name, title) 
32 {
33   //
34   // Standard constructor for Transition Radiation Detector version 2
35   //
36   fIdSenO1 = fIdSenO2 = fIdSenO3 = fIdSenO4 = fIdSenO5 = fIdSenO6 = 0;
37   fIdSenI1 = fIdSenI2 = fIdSenI3 = fIdSenI4 = fIdSenI5 = fIdSenI6 = 0;
38   SetBufferSize(128000);
39 }
40  
41 //_____________________________________________________________________________
42 void AliTRDv2::CreateGeometry()
43 {
44   //
45   // Create geometry for the Transition Radiation Detector version 2
46   // This version covers the full azimuth. 
47   //
48   //Begin_Html
49   /*
50     <img src="gif/AliTRDv2.gif">
51   */
52   //End_Html
53   //Begin_Html
54   /*
55     <img src="gif/AliTRDv2Tree.gif">
56   */
57   //End_Html
58
59
60   // --- Name Conventions : 
61   
62   //     The mother volume and the support structure 
63   //        TRD     --> Mother TRD volume                          (Air) 
64   //        UTRS    --> The sectors of the detector                (Air) 
65   //        UTSP    --> The aluminum part of the support structure (Al) 
66   //        UTII(O) --> The inner parts of the support structure   (Air) 
67   
68   //     The chambers 
69   //        UCI1-6  --> The frame of the inner chambers            (C) 
70   //        UCO1-6  --> The frame of the outer chambers            (C) 
71   //        UII1-6  --> The inner part of the inner chambers       (Air) 
72   //        UIO1-6  --> The inner part of the outer chambers       (Air) 
73   
74   //     The layers inside a chamber 
75   //        UT1I(O) --> Radiator layer                             (CO2) 
76   //        UT2I(O) --> Polyethylene layer                         (PE) 
77   //        UT3I(O) --> Mylar layer                                (Mylar) 
78   //        UXI1-6  --> Xe/C02 layer in the inner chambers         (Xe/C02) 
79   //        UXO1-6  --> Xe/C02 layer in the outer chambers         (Xe/C02) 
80   //        UT5I(O) --> Cu layer (pads/sensitive)                  (Cu) 
81   //        UT6I(O) --> Kapton layer                               (Kapton) 
82   //        UT7I(O) --> NOMEX layer                                (C) 
83   //        UT8I(O) --> Readout layer                              (Al) 
84   
85   
86   // --- Contains geometry information 
87   
88   // --- Number of sectors in the full detector 
89   // --- Number of modules in each sector 
90   // --- z-Coordinates of the TRD-frame 
91   // --- r-Coordinates of the TRD-frame 
92   // --- Thickness of the aluminium of the support frame 
93   // --- Thickness of the interior of the support frame 
94   // --- Thickness of the carbon chamber frame 
95   // --- Thickness and z-position of the PE-layer in the radiator 
96   // --- Thickness and z-position of the radiator 
97   // --- Thickness and z-position of the mylar-layer 
98   // --- Thickness and z-position of the Xe/C02-layer 
99   // --- Thickness and z-position of the Cu-layer (Pads) 
100   // --- Thickness and z-position of the kapton-layer 
101   // --- Thickness and z-position of the NOMEX-layer 
102   //     Simple C-layer for the time being 
103   // --- Thickness and z-position of the readout-layer 
104   // --- Parameter for the arrays 
105   // --- Number of different chambers 
106   // --- Number of rotation matrices 
107
108   AliMC* pMC = AliMC::GetMC();
109   
110   Float_t xpos, ypos, zpos;
111   Int_t icham;
112   Int_t idmat[2];
113   Float_t widma, widmi;
114   Float_t lendifc, widdifc, heightc;
115   Float_t par_ic[3], par_oc[3], par_mo[10], par_sp[4], par_ch[3];
116   
117   Int_t *idtmed = gAlice->Idtmed();
118   
119   //************************************************************************
120   
121   //     Definition of Volumes 
122   
123   //************************************************************************
124   
125   const Int_t ncham = 6; //Number of different chambers
126   
127   widmi  = rmin*TMath::Sin(kPI/nsect);
128   widma  = rmax*TMath::Sin(kPI/nsect);
129   // --- Some parameter for the chambers 
130   lendifc = (zmax1-zmax2)/nmodul;
131   heightc = (rmax-rmin)/nmodul;
132   widdifc = (widma - widmi)/nmodul;
133   // --- Definition of the Mother volume for the TRD (Air) 
134   par_mo[0] = 0.;
135   par_mo[1] = 360.;
136   par_mo[2] = nsect;
137   par_mo[3] = 2.;
138   par_mo[4] = -zmax1;
139   par_mo[5] = rmin;
140   par_mo[6] = rmax;
141   par_mo[7] = zmax1;
142   par_mo[8] = rmin;
143   par_mo[9] = rmax;
144   pMC->Gsvolu("TRD ", "PGON", idtmed[1301], par_mo, 10);
145   // --- Divide the mother volume into sectors 
146   pMC->Gsdvn("UTRS", "TRD ", nsect, 2);
147   // --- Definition of the aluminum part of the support structure (Al) 
148   par_sp[0] = widmi - inframe/2;
149   par_sp[1] = widma - inframe/2;
150   par_sp[2] = zmax1/2;
151   par_sp[3] = (rmax-rmin)/2;
152   pMC->Gsvolu("UTSP", "TRD1", idtmed[1300], par_sp, 4);
153   // --- Definition of the inner part of the support structure (Air) 
154   par_sp[0] = widmi - inframe/2 - alfram1/2;
155   par_sp[1] = widma - inframe/2 - alfram1/2;
156   par_sp[2] = zmax1/4 -alfram2/2;
157   par_sp[3] = (rmax-rmin)/2;
158   pMC->Gsvolu("UTII", "TRD1", idtmed[1301], par_sp, 4);
159   pMC->Gsvolu("UTIO", "TRD1", idtmed[1301], par_sp, 4);
160   // --- Definition of the chambers 
161   char ctagc[5], ctagi[5], ctagx[5];
162   for (icham = 1; icham <= ncham; ++icham) {
163     // --- Carbon frame of the inner chambers (C) 
164     par_ch[0] = widmi + (icham-1) * widdifc - (inframe+alfram1)/2;
165     par_ch[1] = zmax1/4 -alfram2/2;
166     par_ch[2] = heightc/2.;
167     sprintf(ctagc,"UCI%1d",icham);
168     pMC->Gsvolu(ctagc, "BOX ", idtmed[1306], par_ch, 3);
169     // --- Inner part of the inner chambers (Air) 
170     par_ch[0] -= ccframe;
171     par_ch[1] -= ccframe;
172     sprintf(ctagi,"UII%1d",icham);
173     pMC->Gsvolu(ctagi, "BOX ", idtmed[1301], par_ch, 3);
174     // --- Carbon frame of the outer chambers (C) 
175     par_ch[0] = widmi + (icham - 1) * widdifc - 2.;
176     par_ch[1] = (icham - 6) * lendifc / 2. + zmax1/4 -alfram2/2;
177     par_ch[2] = heightc / 2.;
178     sprintf(ctagc,"UCO%1d",icham);
179     pMC->Gsvolu(ctagc, "BOX ", idtmed[1306], par_ch, 3);
180     // --- Inner part of the outer chambers (Air) 
181     par_ch[0] -= ccframe;
182     par_ch[1] -= ccframe;
183     sprintf(ctagi,"UIO%1d",icham);
184     pMC->Gsvolu(ctagi, "BOX ", idtmed[1301], par_ch, 3);
185   }
186   // --- Definition of the layers in each inner chamber 
187   par_ic[0] = -1.;
188   par_ic[1] = -1.;
189   // --- Radiator layer 
190   par_ic[2] = rathick/2;;
191   pMC->Gsvolu("UT1I", "BOX ", idtmed[1311], par_ic, 3);
192   // --- Polyethylene layer 
193   par_ic[2] = pethick/2;
194   pMC->Gsvolu("UT2I", "BOX ", idtmed[1302], par_ic, 3);
195   // --- Mylar layer 
196   par_ic[2] = mythick/2;
197   pMC->Gsvolu("UT3I", "BOX ", idtmed[1307], par_ic, 3);
198   // --- Xe/CO2 layer 
199   par_ic[2] = 1.8;
200   for (icham = 1; icham <= 6; ++icham) {
201     sprintf(ctagx,"UXI%1d",icham);
202     pMC->Gsvolu(ctagx, "BOX ", idtmed[1308], par_ic, 3);
203   }
204   // --- Cu layer 
205   par_ic[2] = cuthick/2;
206   pMC->Gsvolu("UT5I", "BOX ", idtmed[1304], par_ic, 3);
207   // --- Kapton layer 
208   par_ic[2] = kathick/2;
209   pMC->Gsvolu("UT6I", "BOX ", idtmed[1310], par_ic, 3);
210   // --- NOMEX layer 
211   par_ic[2] = nothick/2;
212   pMC->Gsvolu("UT7I", "BOX ", idtmed[1309], par_ic, 3);
213   // --- Read out layer 
214   par_ic[2] = rothick/2;
215   pMC->Gsvolu("UT8I", "BOX ", idtmed[1305], par_ic, 3);
216   // --- Definition of the layers in each outer chamber 
217   par_oc[0] = -1.;
218   par_oc[1] = -1.;
219   // --- Radiator layer 
220   par_oc[2] = rathick/2;;
221   pMC->Gsvolu("UT1O", "BOX ", idtmed[1311], par_oc, 3);
222   // --- Polyethylene layer 
223   par_oc[2] = pethick/2;
224   pMC->Gsvolu("UT2O", "BOX ", idtmed[1302], par_oc, 3);
225   // --- Mylar layer 
226   par_oc[2] = mythick/2;
227   pMC->Gsvolu("UT3O", "BOX ", idtmed[1307], par_oc, 3);
228   // --- Xe/CO2 layer 
229   par_oc[2] = xethick/2;
230   for (icham = 1; icham <= ncham; ++icham) {
231     sprintf(ctagx,"UXO%1d",icham);
232     pMC->Gsvolu(ctagx, "BOX ", idtmed[1308], par_oc, 3);
233   }
234   // --- Cu layer 
235   par_oc[2] = cuthick/2;
236   pMC->Gsvolu("UT5O", "BOX ", idtmed[1304], par_oc, 3);
237   // --- Kapton layer 
238   par_oc[2] = kathick/2;
239   pMC->Gsvolu("UT6O", "BOX ", idtmed[1310], par_oc, 3);
240   // --- NOMEX layer 
241   par_oc[2] = nothick/2;
242   pMC->Gsvolu("UT7O", "BOX ", idtmed[1309], par_oc, 3);
243   // --- Read out layer 
244   par_oc[2] = rothick/2;
245   pMC->Gsvolu("UT8O", "BOX ", idtmed[1305], par_oc, 3);
246   //*************************************************************************
247   
248   //     Positioning of Volumes 
249   
250   //************************************************************************
251   // --- The rotation matrices 
252   AliMatrix(idmat[0], 90., 90., 180., 0., 90., 0.);
253   AliMatrix(idmat[1], 90., 90.,   0., 0., 90., 0.);
254   // --- Position of the layers in a chamber 
255   for (icham = 1; icham <= ncham; ++icham) {
256     // --- The inner chambers 
257     sprintf(ctagi,"UII%1d",icham);
258     sprintf(ctagx,"UXI%1d",icham);
259     pMC->Gspos("UT8I", icham, ctagi, 0., 0., rozpos, 0, "ONLY");
260     pMC->Gspos("UT7I", icham, ctagi, 0., 0., nozpos, 0, "ONLY");
261     pMC->Gspos("UT6I", icham, ctagi, 0., 0., kazpos, 0, "ONLY");
262     pMC->Gspos("UT5I", icham, ctagi, 0., 0., cuzpos, 0, "ONLY");
263     pMC->Gspos(ctagx,  1,     ctagi, 0., 0., xezpos, 0, "ONLY");
264     pMC->Gspos("UT3I", icham, ctagi, 0., 0., myzpos, 0, "ONLY");
265     pMC->Gspos("UT1I", icham, ctagi, 0., 0., razpos, 0, "ONLY");
266     // --- The outer chambers 
267     sprintf(ctagi,"UIO%d",icham);
268     sprintf(ctagx,"UXO%d",icham);
269     pMC->Gspos("UT8O", icham, ctagi, 0., 0., rozpos, 0, "ONLY");
270     pMC->Gspos("UT7O", icham, ctagi, 0., 0., nozpos, 0, "ONLY");
271     pMC->Gspos("UT6O", icham, ctagi, 0., 0., kazpos, 0, "ONLY");
272     pMC->Gspos("UT5O", icham, ctagi, 0., 0., cuzpos, 0, "ONLY");
273     pMC->Gspos(ctagx,  1,     ctagi, 0., 0., xezpos, 0, "ONLY");
274     pMC->Gspos("UT3O", icham, ctagi, 0., 0., myzpos, 0, "ONLY");
275     pMC->Gspos("UT1O", icham, ctagi, 0., 0., razpos, 0, "ONLY");
276   }
277   pMC->Gspos("UT2I", 1, "UT1I", 0., 0., 0., 0, "ONLY");
278   pMC->Gspos("UT2O", 1, "UT1O", 0., 0., 0., 0, "ONLY");
279   // --- Position of the inner part of the chambers in the carbon-frames 
280   for (icham = 1; icham <= ncham; ++icham) {
281     // --- The inner chambers 
282     sprintf(ctagi,"UII%1d",icham);
283     sprintf(ctagc,"UCI%1d",icham);
284     pMC->Gspos(ctagi, 1, ctagc, 0., 0., 0., 0, "ONLY");
285     // --- The outer chambers 
286     sprintf(ctagi,"UIO%1d",icham);
287     sprintf(ctagc,"UCO%1d",icham);
288     pMC->Gspos(ctagi, 1, ctagc, 0., 0., 0., 0, "ONLY");
289   }
290   // --- Position of the chambers in the full TRD-setup 
291   for (icham = 1; icham <= ncham; ++icham) {
292     // --- The inner chambers 
293     xpos = 0.;
294     ypos = 0.;
295     zpos = (icham - .5) * heightc - (rmax-rmin)/2;
296     sprintf(ctagc,"UCI%1d",icham);
297     pMC->Gspos(ctagc, 1, "UTII", xpos, ypos, zpos, 0, "ONLY");
298     // --- The outer chambers 
299     xpos = 0.;
300     ypos = 0. - (icham - 6) * lendifc / 2.;
301     zpos = (icham - .5) * heightc - (rmax-rmin)/2;
302     sprintf(ctagc,"UCO%1d",icham);
303     pMC->Gspos(ctagc, 1, "UTIO", xpos, ypos, zpos, 0, "ONLY");
304   }
305   // --- Position of the inner parts of the support frame 
306   xpos = 0.;
307   ypos = zmax1/4;
308   zpos = 0.;
309   pMC->Gspos("UTIO", 1, "UTSP", xpos,-ypos, zpos, 0, "ONLY");
310   pMC->Gspos("UTII", 1, "UTSP", xpos, ypos, zpos, 0, "ONLY");
311   // --- Position of the support frame in the TRD-sectors 
312   xpos = (rmax+rmin)/2;
313   ypos = 0.;
314   zpos = zmax1/2;
315   pMC->Gspos("UTSP", 1, "UTRS", xpos, ypos, zpos, idmat[0], "ONLY");
316   pMC->Gspos("UTSP", 2, "UTRS", xpos, ypos,-zpos, idmat[1], "ONLY");
317   // --- Position of TRD mother volume in ALICE experiment 
318   pMC->Gspos("TRD ", 1, "ALIC", 0., 0., 0., 0, "ONLY");
319 }
320
321 //_____________________________________________________________________________
322 void AliTRDv2::DrawModule()
323 {
324   //
325   // Draw a shaded view of the Transition Radiation Detector version 1
326   //
327
328   AliMC* pMC = AliMC::GetMC();
329   
330   // Set everything unseen
331   pMC->Gsatt("*", "seen", -1);
332   // 
333   // Set ALIC mother transparent
334   pMC->Gsatt("ALIC","SEEN",0);
335   //
336   // Set the volumes visible
337   pMC->Gsatt("TRD","SEEN",0);
338   pMC->Gsatt("UTRS","SEEN",0);
339   pMC->Gsatt("UTSP","SEEN",0);
340   pMC->Gsatt("UTII","SEEN",0);
341   pMC->Gsatt("UTIO","SEEN",0);
342   pMC->Gsatt("UCI1","SEEN",0);
343   pMC->Gsatt("UII1","SEEN",0);
344   pMC->Gsatt("UCO1","SEEN",0);
345   pMC->Gsatt("UIO1","SEEN",0);
346   pMC->Gsatt("UCI2","SEEN",0);
347   pMC->Gsatt("UII2","SEEN",0);
348   pMC->Gsatt("UCO2","SEEN",0);
349   pMC->Gsatt("UIO2","SEEN",0);
350   pMC->Gsatt("UCI3","SEEN",0);
351   pMC->Gsatt("UII3","SEEN",0);
352   pMC->Gsatt("UCO3","SEEN",0);
353   pMC->Gsatt("UIO3","SEEN",0);
354   pMC->Gsatt("UCI4","SEEN",0);
355   pMC->Gsatt("UII4","SEEN",0);
356   pMC->Gsatt("UCO4","SEEN",0);
357   pMC->Gsatt("UIO4","SEEN",0);
358   pMC->Gsatt("UCI5","SEEN",0);
359   pMC->Gsatt("UII5","SEEN",0);
360   pMC->Gsatt("UCO5","SEEN",0);
361   pMC->Gsatt("UIO5","SEEN",0);
362   pMC->Gsatt("UCI6","SEEN",0);
363   pMC->Gsatt("UII6","SEEN",0);
364   pMC->Gsatt("UCO6","SEEN",0);
365   pMC->Gsatt("UIO6","SEEN",0);
366   pMC->Gsatt("UT1I","SEEN",1);
367   pMC->Gsatt("UXI1","SEEN",1);
368   pMC->Gsatt("UXI2","SEEN",1);
369   pMC->Gsatt("UXI3","SEEN",1);
370   pMC->Gsatt("UXI4","SEEN",1);
371   pMC->Gsatt("UXI5","SEEN",1);
372   pMC->Gsatt("UXI6","SEEN",1);
373   pMC->Gsatt("UT1O","SEEN",1);
374   pMC->Gsatt("UXO1","SEEN",1);
375   pMC->Gsatt("UXO2","SEEN",1);
376   pMC->Gsatt("UXO3","SEEN",1);
377   pMC->Gsatt("UXO4","SEEN",1);
378   pMC->Gsatt("UXO5","SEEN",1);
379   pMC->Gsatt("UXO6","SEEN",1);
380   //
381   pMC->Gdopt("hide", "on");
382   pMC->Gdopt("shad", "on");
383   pMC->Gsatt("*", "fill", 7);
384   pMC->SetClipBox(".");
385   pMC->SetClipBox("*", 0, 2000, -2000, 2000, -2000, 2000);
386   pMC->DefaultRange();
387   pMC->Gdraw("alic", 40, 30, 0, 12, 9.4, .021, .021);
388   pMC->Gdhead(1111, "Transition Radiation Detector Version 2");
389   pMC->Gdman(18, 4, "MAN");
390   pMC->Gdopt("hide", "off");
391 }
392
393 //_____________________________________________________________________________
394 void AliTRDv2::CreateMaterials()
395 {
396   //
397   // Create materials for the Transition Radiation Detector version 2
398   //
399   printf("TRD: Slow simulation with fine geometry\n");
400   AliTRD::CreateMaterials();
401 }
402
403 //_____________________________________________________________________________
404 void AliTRDv2::Init() 
405 {
406   //
407   // Initialise Transition Radiation Detector after geometry has been built
408   //
409   AliTRD::Init();
410   AliMC* pMC = AliMC::GetMC();
411   fIdSenI1 = pMC->VolId("UXI1");
412   fIdSenI2 = pMC->VolId("UXI2");
413   fIdSenI3 = pMC->VolId("UXI3");
414   fIdSenI4 = pMC->VolId("UXI4");
415   fIdSenI5 = pMC->VolId("UXI5");
416   fIdSenI6 = pMC->VolId("UXI6");
417   
418   fIdSenO1 = pMC->VolId("UXO1");
419   fIdSenO2 = pMC->VolId("UXO2");
420   fIdSenO3 = pMC->VolId("UXO3");
421   fIdSenO4 = pMC->VolId("UXO4");
422   fIdSenO5 = pMC->VolId("UXO5");
423   fIdSenO6 = pMC->VolId("UXO6");
424 }
425
426 //_____________________________________________________________________________
427 void AliTRDv2::StepManager()
428 {
429   //
430   // Called at every step in the Transition Radiation Detector version 2
431   //
432   Int_t          idSens, icSens;
433   Int_t          iPla, iCha, iSec;
434   Int_t          iOut;
435   Int_t          vol[3]; 
436   Int_t          iPid;
437   
438   const Float_t  kBig   = 1.0E+12;
439   
440   Float_t        random[1];
441   Float_t        charge;
442   Float_t        betaGamma, pp;
443   Float_t        aMass;
444   Float_t        eLos, qTot;
445   Float_t        hits[4];
446   Float_t        mom[4];
447   Float_t        pTot;
448   
449   TClonesArray  &lhits = *fHits;
450   
451   AliMC* pMC = AliMC::GetMC();
452   
453   // Ionization energy
454   // taken from: Ionization Measurements in High Energy Physics, Springer
455   const Double_t kWion    = 23.0E-9;
456   // Maximum energy for e+ e- g for the step-size calculation
457   const Float_t  kPTotMax = 0.002;
458   // Plateau value of the energy-loss for electron in xenon
459   // taken from: Allison + Comb, Ann. Rev. Nucl. Sci. (1980), 30, 253
460   const Float_t  kPlateau = 1.70;
461   // dN1/dx|min for the gas mixture (90% Xe + 10% CO2)
462   // taken from: Ionization Measurements in High Energy Physics, Springer
463   const Float_t  kPrim    = 43.68;
464   
465   // Set the maximum step size to a very large number for all 
466   // neutral particles and those outside the driftvolume
467   pMC->SetMaxStep(kBig); 
468   
469   // Use only charged tracks 
470   if (( pMC->TrackCharge()   ) &&
471       (!pMC->TrackStop()     ) && 
472       (!pMC->TrackDisappear())) {
473     
474     // Find the sensitive volume
475     idSens = pMC->CurrentVol(0,icSens);
476     iPla   = 0;
477     iOut   = 0;
478     if      (idSens == fIdSenI1) iPla = 1; 
479     else if (idSens == fIdSenI2) iPla = 2;
480     else if (idSens == fIdSenI3) iPla = 3;
481     else if (idSens == fIdSenI4) iPla = 4;
482     else if (idSens == fIdSenI5) iPla = 5;
483     else if (idSens == fIdSenI6) iPla = 6;
484     else if (idSens == fIdSenO1) {
485       iPla = 1;
486       iOut = 1; }
487     else if (idSens == fIdSenO2) {
488       iPla = 2;
489       iOut = 1; }
490     else if (idSens == fIdSenO3) {
491       iPla = 3;
492       iOut = 1; }
493     else if (idSens == fIdSenO4) {
494       iPla = 4;
495       iOut = 1; }
496     else if (idSens == fIdSenO5) {
497       iPla = 5;
498       iOut = 1; }
499     else if (idSens == fIdSenO6) {
500       iPla = 6;
501       iOut = 1; }
502     
503     // Inside a sensitive volume?
504     if (iPla) {
505       
506       // Calculate the energy loss 
507       // ((1/E^2.0) distribution for the fluctuations)
508       pMC->Rndm(random,1);
509       eLos = Eloss(random[0]);
510       
511       // The amount of charge created
512       qTot = (Float_t) ((Int_t) (eLos / kWion) + 1);
513       
514       // The sector number
515       pMC->CurrentVolOff(5,0,iSec);
516       
517       // The chamber number
518       pMC->CurrentVolOff(4,0,iCha);
519       if (iOut) 
520         iCha = (iCha - 1) * 3 + 1;
521       else      
522         iCha =  iCha + 1;
523       
524       vol[0]  = iSec;
525       vol[1]  = iCha;
526       vol[2]  = iPla;
527
528       pMC->TrackPosition(hits);
529       hits[3] = qTot;
530       
531       // Add the hit
532       new(lhits[fNhits++]) AliTRDhit(fIshunt,gAlice->CurrentTrack(),vol,hits);
533       
534       pMC->TrackMomentum(mom);
535       pTot = mom[3];
536       
537       // New step size for electrons only if momentum is small enough
538       iPid = pMC->TrackPid();
539       if ( (iPid >  3) ||
540            ((iPid <= 3) && (pTot < kPTotMax))) {
541         aMass     = pMC->TrackMass();
542         betaGamma = pTot / aMass;
543         pp        = kPrim * BetheBloch(betaGamma);
544         // Take charge > 1 into account
545         charge = pMC->TrackCharge();
546         if (TMath::Abs(charge) > 1) pp = pp * charge*charge;
547       }
548       // Electrons above 20 Mev/c are at the plateau
549       else {
550         pp = kPrim * kPlateau;
551       }
552       
553       // Calculate the maximum step size for the next tracking step
554       do 
555         pMC->Rndm(random,1);
556       while ((random[0] == 1.) || (random[0] == 0.));
557       pMC->SetMaxStep( - TMath::Log(random[0]) / pp);
558       
559     }
560   }
561 }