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