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