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