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