2c8bc630a2256e2ef850bf41750195f235ae236d
[u/mrichter/AliRoot.git] / TGeant3 / AliGeant3.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 /*
17 $Log$
18 Revision 1.17  2001/06/15 09:31:23  morsch
19 In gudcay: write only first generation decay products to stack to respect the possibility of secondary, tertiary, ... vertices during tracking.
20
21 Revision 1.16  2001/05/16 14:57:23  alibrary
22 New files for folders and Stack
23
24 Revision 1.15  2001/03/20 06:28:49  alibrary
25 New detector loop split in 2
26
27 Revision 1.14  2000/12/20 08:39:39  fca
28 Support for Cerenkov and process list in Virtual MC
29
30 Revision 1.13  2000/11/30 07:12:54  alibrary
31 Introducing new Rndm and QA classes
32
33 Revision 1.12  2000/11/06 11:35:46  morsch
34 Call BuildGeometry() after Init() to be able to share common detector parameters.
35
36 Revision 1.11  2000/10/04 16:30:22  fca
37 Add include for exit()
38
39 Revision 1.10  2000/10/02 21:28:16  fca
40 Removal of useless dependecies via forward declarations
41
42 Revision 1.9  2000/09/12 14:36:17  morsch
43 In gudcay(): call ForceDecay() before Decay()
44
45 Revision 1.8  2000/09/06 14:56:34  morsch
46 gudcay() method implemented.
47 Decays are performed using  the AliDecayer interface. The pointer to the instance of AliDecayer
48 is obtained from geant3 (will be gMC once it has been implemented for Geant4).
49
50 Revision 1.7  2000/07/13 16:19:10  fca
51 Mainly coding conventions + some small bug fixes
52
53 Revision 1.6  2000/07/11 18:24:59  fca
54 Coding convention corrections + few minor bug fixes
55
56 Revision 1.5  2000/05/20 14:49:48  fca
57 Call gdebug at the end of gustep
58
59 Revision 1.4  2000/04/26 10:17:32  fca
60 Changes in Lego for G4 compatibility
61
62 Revision 1.3  2000/04/07 11:12:35  fca
63 G4 compatibility changes
64
65 Revision 1.2  2000/02/29 19:11:17  fca
66 Move gucode into AliGeant3.cxx
67
68 Revision 1.1  2000/02/23 16:25:25  fca
69 AliVMC and AliGeant3 classes introduced
70 ReadEuclid moved from AliRun to AliModule
71
72 */
73
74 #include <stdlib.h>
75
76 #include <TParticle.h>
77 #include <TStopwatch.h>
78
79 #include "AliDecayer.h"
80 #include "AliGeant3.h"
81 #include "AliRun.h"
82 #include "TGeant3.h"
83 #include "AliCallf77.h" 
84 #include "AliModule.h"
85 #include "AliMagF.h"
86 #include "AliGenerator.h"
87
88 #ifndef WIN32
89
90 # define rxgtrak rxgtrak_
91 # define rxouth  rxouth_
92 # define rxinh   rxinh_
93 #else
94
95 # define rxgtrak RXGTRAK 
96 # define rxouth  RXOUTH
97 # define rxinh   RXINH
98 #endif
99
100 ClassImp(AliGeant3)
101
102 AliGeant3::AliGeant3(const char *title) : 
103   TGeant3(title) {}
104
105 //____________________________________________________________________________
106 void AliGeant3::FinishGeometry()
107 {
108   //
109   // Finalise geometry construction
110   //
111   TGeant3::FinishGeometry();
112   //Create the color table
113   SetColors();
114 }
115
116 //____________________________________________________________________________
117 void AliGeant3::Init()
118 {
119     //
120     //=================Create Materials and geometry
121     //
122     TStopwatch stw;
123     TObjArray *modules = gAlice->Modules();
124     TIter next(modules);
125     AliModule *detector;
126     printf("Geometry creation:\n");
127     while((detector = (AliModule*)next())) {
128       stw.Start();
129       // Initialise detector materials and geometry
130       detector->CreateMaterials();
131       detector->CreateGeometry();
132       printf("%10s R:%.2fs C:%.2fs\n",
133              detector->GetName(),stw.RealTime(),stw.CpuTime());
134     }
135     //Terminate building of geometry
136     FinishGeometry();
137     
138     printf("Initialisation:\n");
139     next.Reset();
140     while((detector = (AliModule*)next())) {
141       stw.Start();
142       // Initialise detector and display geometry
143       detector->Init();
144       detector->BuildGeometry();
145       printf("%10s R:%.2fs C:%.2fs\n",
146              detector->GetName(),stw.RealTime(),stw.CpuTime());
147     }
148
149 }
150
151 //____________________________________________________________________________
152 void AliGeant3::ProcessRun(Int_t nevent)
153 {
154   //
155   // Process the run
156   //
157   Int_t todo = TMath::Abs(nevent);
158   for (Int_t i=0; i<todo; i++) {
159   // Process one run (one run = one event)
160      gAlice->BeginEvent();
161      ProcessEvent();
162      gAlice->FinishEvent();
163   }
164 }
165
166 //_____________________________________________________________________________
167 void AliGeant3::ProcessEvent()
168 {
169   //
170   // Process one event
171   //
172   Gtrigi();
173   Gtrigc();
174   Gtrig();
175 }
176
177 //_____________________________________________________________________________
178 void AliGeant3::SetColors()
179 {
180   //
181   // Set the colors for all the volumes
182   // this is done sequentially for all volumes
183   // based on the number of their medium
184   //
185   Int_t kv, icol;
186   Int_t jvolum=fGclink->jvolum;
187   //Int_t jtmed=fGclink->jtmed;
188   //Int_t jmate=fGclink->jmate;
189   Int_t nvolum=fGcnum->nvolum;
190   char name[5];
191   //
192   //    Now for all the volumes
193   for(kv=1;kv<=nvolum;kv++) {
194     //     Get the tracking medium
195     Int_t itm=Int_t (fZq[fZlq[jvolum-kv]+4]);
196     //     Get the material
197     //Int_t ima=Int_t (fZq[fZlq[jtmed-itm]+6]);
198     //     Get z
199     //Float_t z=fZq[fZlq[jmate-ima]+7];
200     //     Find color number
201     //icol = Int_t(z)%6+2;
202     //icol = 17+Int_t(z*150./92.);
203     //icol = kv%6+2;
204     icol = itm%6+2;
205     strncpy(name,(char*)&fZiq[jvolum+kv],4);
206     name[4]='\0';
207     Gsatt(name,"COLO",icol);
208   }
209 }
210
211 //_____________________________________________________________________________
212 //
213 //                 Interfaces to Fortran
214 //
215 //_____________________________________________________________________________
216
217 extern "C" void type_of_call  rxgtrak (Int_t &mtrack, Int_t &ipart, Float_t *pmom, 
218                                        Float_t &e, Float_t *vpos, Float_t *polar,
219                                        Float_t &tof)
220 {
221   //
222   //     Fetches next track from the ROOT stack for transport. Called by the
223   //     modified version of GTREVE.
224   //
225   //              Track number in the ROOT stack. If MTRACK=0 no
226   //      mtrack  more tracks are left in the stack to be
227   //              transported.
228   //      ipart   Particle code in the GEANT conventions.
229   //      pmom[3] Particle momentum in GeV/c
230   //      e       Particle energy in GeV
231   //      vpos[3] Particle position
232   //      tof     Particle time of flight in seconds
233   //
234   Int_t pdg;
235   gAlice->GetNextTrack(mtrack, pdg, pmom, e, vpos, polar, tof);
236   ipart = gMC->IdFromPDG(pdg);
237   mtrack++;
238 }
239
240 //_____________________________________________________________________________
241 extern "C" void type_of_call  rxouth ()
242 {
243   //
244   // Called by Gtreve at the end of each primary track
245   //
246   gAlice->FinishPrimary();
247 }
248
249 //_____________________________________________________________________________
250 extern "C" void type_of_call  rxinh ()
251 {
252   //
253   // Called by Gtreve at the beginning of each primary track
254   //
255   gAlice->BeginPrimary();
256 }
257
258 #ifndef WIN32
259 #  define gudigi gudigi_
260 #  define guhadr guhadr_
261 #  define guout  guout_
262 #  define guphad guphad_
263 #  define gudcay gudcay_
264 #  define guiget guiget_
265 #  define guinme guinme_
266 #  define guinti guinti_
267 #  define gunear gunear_
268 #  define guskip guskip_
269 #  define guview guview_
270 #  define gupara gupara_
271 #  define gudtim gudtim_
272 #  define guplsh guplsh_
273 #  define gutrev gutrev_
274 #  define gutrak gutrak_
275 #  define guswim guswim_
276 #  define gufld  gufld_
277 #  define gustep gustep_
278 #  define gukine gukine_
279 #  define uglast uglast_
280
281 #  define gheish gheish_
282 #  define flufin flufin_
283 #  define gfmfin gfmfin_
284 #  define gpghei gpghei_
285 #  define fldist fldist_
286 #  define gfmdis gfmdis_
287 #  define ghelx3 ghelx3_
288 #  define ghelix ghelix_
289 #  define grkuta grkuta_
290 #  define gtrack gtrack_
291 #  define gtreveroot gtreveroot_
292 #  define glast  glast_
293
294 #else
295 #  define gudigi GUDIGI
296 #  define guhadr GUHADR
297 #  define guout  GUOUT
298 #  define guphad GUPHAD
299 #  define gudcay GUDCAY
300 #  define guiget GUIGET
301 #  define guinme GUINME
302 #  define guinti GUINTI
303 #  define gunear GUNEAR
304 #  define guskip GUSKIP
305 #  define guview GUVIEW
306 #  define gupara GUPARA
307 #  define gudtim GUDTIM
308 #  define guplsh GUPLSH
309 #  define gutrev GUTREV
310 #  define gutrak GUTRAK
311 #  define guswim GUSWIM
312 #  define gufld  GUFLD
313 #  define gustep GUSTEP
314 #  define gukine GUKINE
315 #  define uglast UGLAST
316
317 #  define gheish GHEISH
318 #  define flufin FLUFIN
319 #  define gfmfin GFMFIN
320 #  define gpghei GPGHEI
321 #  define fldist FLDIST
322 #  define gfmdis GFMDIS
323 #  define ghelx3 GHELX3
324 #  define ghelix GHELIX
325 #  define grkuta GRKUTA
326 #  define gtrack GTRACK
327 #  define gtreveroot GTREVEROOT
328 #  define glast  GLAST
329
330 #endif
331
332 extern "C" type_of_call void gheish();
333 extern "C" type_of_call void flufin();
334 extern "C" type_of_call void gfmfin();
335 extern "C" type_of_call void gpghei();
336 extern "C" type_of_call void fldist();
337 extern "C" type_of_call void gfmdis();
338 extern "C" type_of_call void ghelx3(Float_t&, Float_t&, Float_t*, Float_t*);
339 extern "C" type_of_call void ghelix(Float_t&, Float_t&, Float_t*, Float_t*);
340 extern "C" type_of_call void grkuta(Float_t&, Float_t&, Float_t*, Float_t*);
341 extern "C" type_of_call void gtrack();
342 extern "C" type_of_call void gtreveroot();
343 extern "C" type_of_call void glast();
344
345 extern "C" type_of_call {
346
347 //______________________________________________________________________
348 void gudigi() 
349 {
350 //
351 //    ******************************************************************
352 //    *                                                                *
353 //    *       User routine to digitize one event                       *
354 //    *                                                                *
355 //    *    ==>Called by : GTRIG                                        *
356 //    *                                                                *
357 //    ******************************************************************
358
359 }
360
361
362 //______________________________________________________________________
363 void guhadr()
364 {
365 //
366 //    ******************************************************************
367 //    *                                                                *
368 //    *       User routine to generate one hadronic interaction        *
369 //    *                                                                *
370 //    *    ==>Called by : GTHADR,GTNEUT                                *
371 //    *                                                                *
372 //    ******************************************************************
373 //
374 //
375 //    ------------------------------------------------------------------
376 //
377       TGeant3* geant3 = (TGeant3*) gMC;
378       Int_t ihadr=geant3->Gcphys()->ihadr;
379       if (ihadr<4)       gheish();
380       else if (ihadr==4) flufin();
381       else               gfmfin();
382 }
383
384 //______________________________________________________________________
385 void guout()
386 {
387 //
388 //    ******************************************************************
389 //    *                                                                *
390 //    *       User routine called at the end of each event             *
391 //    *                                                                *
392 //    *    ==>Called by : GTRIG                                        *
393 //    *                                                                *
394 //    ******************************************************************
395 //
396 //
397 //    ------------------------------------------------------------------
398 //
399 }
400
401 //______________________________________________________________________
402 void guphad()
403 {
404 //
405 //    ******************************************************************
406 //    *                                                                *
407 //    *       User routine to compute Hadron. inter. probabilities     *
408 //    *                                                                *
409 //    *    ==>Called by : GTHADR,GTNEUT                                *
410 //    *                                                                *
411 //    ******************************************************************
412 //
413 //
414 //    ------------------------------------------------------------------
415 //
416       TGeant3* geant3 = (TGeant3*) gMC;
417       Int_t ihadr=geant3->Gcphys()->ihadr;
418       if (ihadr<4)       gpghei();
419       else if (ihadr==4) fldist();
420       else               gfmdis();
421 }
422
423 //______________________________________________________________________
424 void gudcay()
425 {
426 //
427 //    ******************************************************************
428 //    *                                                                *
429 //    *       User routine to decay particles                          *
430 //    *                                                                *
431 //    *    ==>Called by : GDECAY                                       *
432 //    *                                                                *
433 //    ******************************************************************
434 //
435 //
436 //    ------------------------------------------------------------------
437 //
438     
439     TGeant3* geant3=(TGeant3*) gMC;
440     // set decay table
441     gMC->Decayer()->ForceDecay();
442
443 // Initialize 4-momentum vector    
444     Int_t ipart = geant3->Gckine()->ipart;
445     TLorentzVector p;
446     
447     p[0]=geant3->Gctrak()->vect[3];
448     p[1]=geant3->Gctrak()->vect[4];
449     p[2]=geant3->Gctrak()->vect[5];
450     p[3]=geant3->Gctrak()->vect[6];    
451     
452 // Convert from geant to lund particle code
453     Int_t iplund=gMC->PDGFromId(ipart);
454 // Particle list
455     static TClonesArray *particles;
456     if(!particles) particles=new TClonesArray("TParticle",1000);
457 // Decay
458     gMC->Decayer()->Decay(iplund, &p);
459     
460 // Fetch Particles
461     Int_t np = geant3->Decayer()->ImportParticles(particles);
462     if (np <=1) return;
463
464     TParticle *  iparticle = (TParticle *) particles->At(0);
465     Int_t ipF = 0, ipL = 0 ;
466     Int_t i,j;
467
468 // Array to flag deselected particles
469     Int_t* pFlag = new Int_t[np];
470     for (i=0; i<np; i++) pFlag[i]=0;
471 // Particle loop
472     for (i=1; i < np; i++) 
473     {
474         iparticle = (TParticle *) particles->At(i);
475         ipF = iparticle->GetFirstDaughter();
476         ipL = iparticle->GetLastDaughter();     
477         Int_t kf = iparticle->GetPdgCode();
478         Int_t ks = iparticle->GetStatusCode();
479 //
480 // Deselect daughters of deselected particles
481 // and jump skip the current particle
482         if (pFlag[i] == 1) {
483             if (ipF >= 0) for (j=ipF-1; j<ipL; j++) pFlag[j]=1;
484             continue;
485         } // deselected ??
486 // Particles with long life-time are put on the stack for further tracking
487 // Decay products are deselected
488 //      
489         if (ks != 1) { 
490             Double_t lifeTime = gMC->Decayer()->GetLifetime(kf);
491             if (lifeTime > (Double_t) 1.e-15) {
492                 if (ipF >= 0) for (j=ipF-1; j<ipL; j++) pFlag[j]=1;
493             } else{
494                 continue;
495             }
496         } // ks==1 ?
497 // Skip neutrinos
498         if (kf==12 || kf ==-12) continue;
499         if (kf==14 || kf ==-14) continue;
500         if (kf==16 || kf ==-16) continue;
501         
502         Int_t index=geant3->Gcking()->ngkine;
503 // Put particle on geant stack
504 // momentum vector
505         
506         (geant3->Gcking()->gkin[index][0]) = iparticle->Px();
507         (geant3->Gcking()->gkin[index][1]) = iparticle->Py();
508         (geant3->Gcking()->gkin[index][2]) = iparticle->Pz();
509         (geant3->Gcking()->gkin[index][3]) = iparticle->Energy();
510         Int_t ilu = gMC->IdFromPDG(kf);
511
512 // particle type        
513         (geant3->Gcking()->gkin[index][4]) = Float_t(ilu);
514 // position
515         (geant3->Gckin3()->gpos[index][0]) = geant3->Gctrak()->vect[0];
516         (geant3->Gckin3()->gpos[index][1]) = geant3->Gctrak()->vect[1];
517         (geant3->Gckin3()->gpos[index][2]) = geant3->Gctrak()->vect[2];
518 // time of flight offset (mm)
519         (geant3->Gcking()->tofd[index])    = 0.;
520 // increase stack counter
521         (geant3->Gcking()->ngkine)=index+1;
522     }
523 //
524     if (pFlag) delete[] pFlag;
525 }
526
527 //______________________________________________________________________
528 void guiget(Int_t&, Int_t&, Int_t&)
529 {
530 //
531 //    ******************************************************************
532 //    *                                                                *
533 //    *       User routine for interactive control of GEANT            *
534 //    *                                                                *
535 //    *    ==>Called by : <GXINT>, GINCOM                              *
536 //    *                                                                *
537 //    ******************************************************************
538 //
539 //
540 //    ------------------------------------------------------------------
541 //
542 }
543
544 //______________________________________________________________________
545 void guinme(Float_t*, Int_t&, Float_t*, Int_t& IYES)
546 {
547 //
548 //    **********************************************
549 //    *                                            *
550 //    *    USER ROUTINE TO PROVIDE GINME FUNCTION  *
551 //    *    FOR ALL USER SHAPES IDENTIFIED BY THE   *
552 //    *    SHAPE NUMBER SH. POINT IS GIVEN IN X    *
553 //    *    THE PARAMETERS ARE GIVEN IN P. IYES IS  *
554 //    *    RETURNED 1 IF POINT IS IN, 0 IF POINT   *
555 //    *    IS OUT AND LESS THAN ZERO IF SHAPE      *
556 //    *    NUMBER IS NOT SUPPORTED.                *
557 //    *                                            *
558 //    *    ==>Called by : GINME                    *
559 //    *                                            *
560 //    **********************************************
561 //
562       IYES=-1;
563 }
564
565 //______________________________________________________________________
566 void guinti()
567 {
568 //
569 //    ******************************************************************
570 //    *                                                                *
571 //    *       User routine for interactive version                     *
572 //    *                                                                *
573 //    *    ==>Called by : <GXINT>,  GINTRI                             *
574 //    *                                                                *
575 //    ******************************************************************
576 //
577 //
578 //    ------------------------------------------------------------------
579 //
580 }
581
582 //______________________________________________________________________
583 void gunear(Int_t&, Int_t&, Float_t*, Int_t&)
584 {
585 //
586 //    ******************************************************************
587 //    *                                                                *
588 //    *    User search                                                 *
589 //    *       ISEARC to identify the given volume                      *
590 //    *       ICALL  to identify the calling routine                   *
591 //    *              1 GMEDIA like                                     *
592 //    *              2 GNEXT like                                      *
593 //    *       X      coordinates (+direction for ICALL=2)              *
594 //    *       JNEAR  address of default list of neighbours             *
595 //    *              (list to be overwriten by user)                   *
596 //    *                                                                *
597 //    *    Called by : GFTRAC, GINVOL, GTMEDI, GTNEXT, GNEXT, GMEDIA   *
598 //    *                                                                *
599 //    ******************************************************************
600 //
601 //
602 //    ------------------------------------------------------------------
603 //
604 }
605
606 //______________________________________________________________________
607 void guskip(Int_t& ISKIP)
608 {
609 //
610 //    ******************************************************************
611 //    *                                                                *
612 //    *   User routine to skip unwanted tracks                         *
613 //    *                                                                *
614 //    *   Called by : GSSTAK                                           *
615 //    *   Author    : F.Bruyant                                        *
616 //    *                                                                *
617 //    ******************************************************************
618 //
619 //
620 //    ------------------------------------------------------------------
621 //
622       ISKIP = 0;
623 }
624
625 //______________________________________________________________________
626 void guswim(Float_t& CHARGE, Float_t& STEP, Float_t* VECT, Float_t* VOUT)
627 {
628 //
629 //    ******************************************************************
630 //    *                                                                *
631 //    *       User routine to control tracking of one track            *
632 //    *       in a magnetic field                                      *
633 //    *                                                                *
634 //    *    ==>Called by : GTELEC,GTHADR,GTMUON                         *
635 //    *                                                                *
636 //    ******************************************************************
637 //
638 //
639 //    ------------------------------------------------------------------
640 //
641   TGeant3* geant3 = (TGeant3*) gMC;
642   Int_t ifield=geant3->Gctmed()->ifield;
643   Float_t fieldm=geant3->Gctmed()->fieldm;
644
645   if (ifield==3) {
646     Float_t fldcharge = fieldm*CHARGE;
647     ghelx3(fldcharge,STEP,VECT,VOUT);
648   }
649   else if (ifield==2) ghelix(CHARGE,STEP,VECT,VOUT);
650   else                grkuta(CHARGE,STEP,VECT,VOUT);
651 }
652
653 //______________________________________________________________________
654 void guview(Int_t&, Int_t&, DEFCHARD, Int_t& DEFCHARL)
655 {
656 //
657 //    ******************************************************************
658 //    *                                                                *
659 //    *       User routine for interactive version                     *
660 //    *                                                                *
661 //    *    ==>Called by : <GXINT>, GINC1                               *
662 //    *                                                                *
663 //    ******************************************************************
664 //
665 //
666 //    ------------------------------------------------------------------
667 //
668 }
669
670 //______________________________________________________________________
671 void gupara()
672 {
673 //
674 //    ******************************************************************
675 //    *                                                                *
676 //    *       User routine called every time a particle falls below    *
677 //    *       parametrization threshold. This routine should create    *
678 //    *       the parametrization stack, and, when this is full,       *
679 //    *       parametrize the shower and track the geantinos.          *
680 //    *                                                                *
681 //    *    ==>Called by : GTRACK                                       *
682 //    *                                                                *
683 //    ******************************************************************
684 //
685 //
686 //    ------------------------------------------------------------------
687 //
688 }
689
690 //______________________________________________________________________
691 Float_t gudtim(Float_t&, Float_t&, Int_t&, Int_t&)
692 {
693 //
694 //    ******************************************************************
695 //    *                                                                *
696 //    *       User function called by GCDRIF to return drift time      *
697 //    *                                                                *
698 //    *    ==>Called by : GCDRIF                                       *
699 //    *                                                                *
700 //    ******************************************************************
701 //
702 //
703 //    ------------------------------------------------------------------
704 //
705       return 0;
706 }
707
708
709 //______________________________________________________________________
710 Float_t guplsh(Int_t&, Int_t&)
711 {
712 //
713 //    ******************************************************************
714 //    *                                                                *
715 //    *                                                                *
716 //    *    ==>Called by : GLISUR                                       *
717 //    *                                                                *
718 //    ******************************************************************
719 //
720 //
721 //    ------------------------------------------------------------------
722 //
723 //
724 //*** By default this defines perfect smoothness
725       return 1;
726 }
727
728 //______________________________________________________________________
729 void gutrak()
730 {
731 //
732 //    ******************************************************************
733 //    *                                                                *
734 //    *       User routine to control tracking of one track            *
735 //    *                                                                *
736 //    *    ==>Called by : GTREVE                                       *
737 //    *                                                                *
738 //    ******************************************************************
739 //
740 //
741 //    ------------------------------------------------------------------
742 //
743      gAlice->PreTrack();
744
745      gtrack();
746
747      gAlice->PostTrack();
748 }
749
750 //______________________________________________________________________
751 void gutrev()
752 {
753 //
754 //    ******************************************************************
755 //    *                                                                *
756 //    *       User routine to control tracking of one event            *
757 //    *                                                                *
758 //    *    ==>Called by : GTRIG                                        *
759 //    *                                                                *
760 //    ******************************************************************
761 //
762 //
763 //    ------------------------------------------------------------------
764 //
765   gtreveroot();
766 }
767
768
769 //______________________________________________________________________
770 void gufld(Float_t *x, Float_t *b)
771 {
772       if(gAlice->Field()) {
773          gAlice->Field()->Field(x,b);
774       } else {
775          printf("No mag field defined!\n");
776          b[0]=b[1]=b[2]=0;
777       }
778 }
779
780 //______________________________________________________________________
781 void gustep()
782 {
783 //
784 //    ******************************************************************
785 //    *                                                                *
786 //    *       User routine called at the end of each tracking step     *
787 //    *       INWVOL is different from 0 when the track has reached    *
788 //    *              a volume boundary                                 *
789 //    *       ISTOP is different from 0 if the track has stopped       *
790 //    *                                                                *
791 //    *    ==>Called by : GTRACK                                       *
792 //    *                                                                *
793 //    ******************************************************************
794 //
795
796
797   TLorentzVector x;
798   Float_t r;
799   Int_t ipp, jk, id, nt;
800   Float_t polar[3]={0,0,0};
801   Float_t mom[3];
802   AliMCProcess pProc;
803
804   
805   TGeant3* geant3 = (TGeant3*) gMC;
806   //     Stop particle if outside user defined tracking region 
807   gMC->TrackPosition(x);
808   r=TMath::Sqrt(x[0]*x[0]+x[1]*x[1]);
809   if (r > gAlice->TrackingRmax() || TMath::Abs(x[2]) > gAlice->TrackingZmax()) {
810         gMC->StopTrack();
811   }
812
813   // --- Add new created particles 
814   if (gMC->NSecondaries() > 0) {
815     pProc=gMC->ProdProcess(0);
816     for (jk = 0; jk < geant3->Gcking()->ngkine; ++jk) {
817       ipp = Int_t (geant3->Gcking()->gkin[jk][4]+0.5);
818       // --- Skip neutrinos! 
819       if (ipp != 4) {
820         gAlice->SetTrack(1,gAlice->CurrentTrack(),gMC->PDGFromId(ipp), geant3->Gcking()->gkin[jk], 
821                          geant3->Gckin3()->gpos[jk], polar,geant3->Gctrak()->tofg, pProc, nt);
822       }
823     }
824   }
825   // Cherenkov photons here
826   if ( geant3->Gckin2()->ngphot ) {
827     for (jk = 0; jk < geant3->Gckin2()->ngphot; ++jk) {
828       mom[0]=geant3->Gckin2()->xphot[jk][3]*geant3->Gckin2()->xphot[jk][6];
829       mom[1]=geant3->Gckin2()->xphot[jk][4]*geant3->Gckin2()->xphot[jk][6];
830       mom[2]=geant3->Gckin2()->xphot[jk][5]*geant3->Gckin2()->xphot[jk][6];
831       gAlice->SetTrack(1, gAlice->CurrentTrack(), gMC->PDGFromId(50),
832                        mom,                             //momentum
833                        geant3->Gckin2()->xphot[jk],     //position
834                        &geant3->Gckin2()->xphot[jk][7], //polarisation
835                        geant3->Gckin2()->xphot[jk][10], //time of flight
836                        kPCerenkov, nt);
837       }
838   }
839   // --- Particle leaving the setup ?
840   if (!gMC->IsTrackOut()) 
841     if ((id=gAlice->DetFromMate(geant3->Gctmed()->numed)) >= 0) gAlice->StepManager(id);
842
843   // --- Standard GEANT debug routine 
844   if(geant3->Gcflag()->idebug) geant3->Gdebug();
845 }
846
847 //______________________________________________________________________
848 void gukine ()
849 {
850 //
851 //    ******************************************************************
852 //    *                                                                *
853 //    *       Read or Generates Kinematics for primary tracks          *
854 //    *                                                                *
855 //    *    ==>Called by : GTRIG                                        *
856 //    *                                                                *
857 //    ******************************************************************
858 //
859 //
860 //    ------------------------------------------------------------------
861 //
862   gAlice->Generator()->Generate();
863 }
864
865
866 //______________________________________________________________________
867 void uglast()
868 {
869 //
870 //    ******************************************************************
871 //    *                                                                *
872 //    *       User routine called at the end of the run                *
873 //    *                                                                *
874 //    *    ==>Called by : GRUN                                         *
875 //    *                                                                *
876 //    ******************************************************************
877 //
878 //
879 }
880 }
881