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