Obsolete material removed.
[u/mrichter/AliRoot.git] / MUON / AliMUONv1.cxx
CommitLineData
a9e2aefa 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
d7c4fbc4 3 * SigmaEffect_thetadegrees *
a9e2aefa 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 *
2c799aa2 12 * about the suitability of this software for any purpeateose. It is *
a9e2aefa 13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
88cb7938 16/* $Id$ */
a9e2aefa 17
18/////////////////////////////////////////////////////////
acc838fd 19// Manager and hits classes for set:MUON version 1 //
a9e2aefa 20/////////////////////////////////////////////////////////
acc838fd 21
abaf7c9d 22#include <TRandom.h>
23#include <TF1.h>
116cbefd 24#include <TClonesArray.h>
acc838fd 25#include <TRandom.h>
acc838fd 26#include <TGeoMatrix.h>
88cb7938 27#include <TVirtualMC.h>
a9e2aefa 28
30178c30 29#include "AliMUONv1.h"
a9e2aefa 30#include "AliConst.h"
31#include "AliMUONChamber.h"
88cb7938 32#include "AliMUONConstants.h"
33#include "AliMUONFactory.h"
a9e2aefa 34#include "AliMUONHit.h"
8c449e83 35#include "AliMUONTriggerCircuit.h"
acc838fd 36#include "AliMUONVGeometryBuilder.h"
37#include "AliMUONChamberGeometry.h"
38#include "AliMUONGeometryEnvelope.h"
39#include "AliMUONGeometryConstituent.h"
88cb7938 40#include "AliMagF.h"
41#include "AliRun.h"
5d12ce38 42#include "AliMC.h"
a9e2aefa 43
44ClassImp(AliMUONv1)
45
46//___________________________________________
30178c30 47AliMUONv1::AliMUONv1()
48 : AliMUON(),
49 fTrackMomentum(), fTrackPosition(),fGlobalTransformation(0)
a9e2aefa 50{
51// Constructor
1391e633 52 fChambers = 0;
1391e633 53 fStepManagerVersionOld = kFALSE;
374ebd7d 54 fAngleEffect = kTRUE;
1391e633 55 fStepMaxInActiveGas = 0.6;
56 fStepSum = 0x0;
57 fDestepSum = 0x0;
58 fElossRatio = 0x0;
59 fAngleEffect10 = 0x0;
60 fAngleEffectNorma= 0x0;
61}
a9e2aefa 62//___________________________________________
63AliMUONv1::AliMUONv1(const char *name, const char *title)
1391e633 64 : AliMUON(name,title), fTrackMomentum(), fTrackPosition()
a9e2aefa 65{
66// Constructor
ba030c0e 67 // By default include all stations
ba030c0e 68 AliMUONFactory factory;
69 factory.Build(this, title);
c33d9661 70
71 fStepManagerVersionOld = kFALSE;
374ebd7d 72 fAngleEffect = kTRUE;
1391e633 73 fStepMaxInActiveGas = 0.6;
74
75 fStepSum = new Float_t [AliMUONConstants::NCh()];
76 fDestepSum = new Float_t [AliMUONConstants::NCh()];
77 for (Int_t i=0; i<AliMUONConstants::NCh(); i++) {
78 fStepSum[i] =0.0;
79 fDestepSum[i]=0.0;
80 }
81 // Ratio of particle mean eloss with respect MIP's Khalil Boudjemline, sep 2003, PhD.Thesis and Particle Data Book
82 fElossRatio = new TF1("ElossRatio","[0]+[1]*x+[2]*x*x+[3]*x*x*x+[4]*x*x*x*x",0.5,5.);
83 fElossRatio->SetParameter(0,1.02138);
84 fElossRatio->SetParameter(1,-9.54149e-02);
85 fElossRatio->SetParameter(2,+7.83433e-02);
86 fElossRatio->SetParameter(3,-9.98208e-03);
87 fElossRatio->SetParameter(4,+3.83279e-04);
88
89 // Angle effect in tracking chambers at theta =10 degres as a function of ElossRatio (Khalil BOUDJEMLINE sep 2003 Ph.D Thesis) (in micrometers)
90 fAngleEffect10 = new TF1("AngleEffect10","[0]+[1]*x+[2]*x*x",0.5,3.0);
91 fAngleEffect10->SetParameter(0, 1.90691e+02);
92 fAngleEffect10->SetParameter(1,-6.62258e+01);
93 fAngleEffect10->SetParameter(2,+1.28247e+01);
94 // Angle effect: Normalisation form theta=10 degres to theta between 0 and 10 (Khalil BOUDJEMLINE sep 2003 Ph.D Thesis)
95 // Angle with respect to the wires assuming that chambers are perpendicular to the z axis.
96 fAngleEffectNorma = new TF1("AngleEffectNorma","[0]+[1]*x+[2]*x*x+[3]*x*x*x",0.0,10.0);
97 fAngleEffectNorma->SetParameter(0,4.148);
98 fAngleEffectNorma->SetParameter(1,-6.809e-01);
99 fAngleEffectNorma->SetParameter(2,5.151e-02);
100 fAngleEffectNorma->SetParameter(3,-1.490e-03);
acc838fd 101
102 // Define the global transformation:
103 // Transformation from the old ALICE coordinate system to a new one:
104 // x->-x, z->-z
105 TGeoRotation* rotGlobal
106 = new TGeoRotation("rotGlobal", 90., 180., 90., 90., 180., 0.);
107 fGlobalTransformation = new TGeoCombiTrans(0., 0., 0., rotGlobal);
a9e2aefa 108}
109
fb1bf5c0 110//_____________________________________________________________________________
111AliMUONv1::AliMUONv1(const AliMUONv1& right)
112 : AliMUON(right)
113{
114 // copy constructor (not implemented)
115
116 Fatal("AliMUONv1", "Copy constructor not provided.");
117}
118
a9e2aefa 119//___________________________________________
acc838fd 120AliMUONv1::~AliMUONv1()
a9e2aefa 121{
acc838fd 122// Destructor
a9e2aefa 123
acc838fd 124 delete fGlobalTransformation;
125}
b17c0c87 126
fb1bf5c0 127//_____________________________________________________________________________
128AliMUONv1& AliMUONv1::operator=(const AliMUONv1& right)
129{
130 // assignement operator (not implemented)
131
132 // check assignement to self
133 if (this == &right) return *this;
134
135 Fatal("operator =", "Assignement operator not provided.");
136
137 return *this;
138}
139
acc838fd 140//__________________________________________________
141void AliMUONv1::CreateGeometry()
142{
b17c0c87 143//
acc838fd 144// Construct geometry using geometry builders.
a9e2aefa 145//
a9e2aefa 146
acc838fd 147 for (Int_t i=0; i<fGeometryBuilders->GetEntriesFast(); i++) {
3c084d9f 148
acc838fd 149 // Get the builder
150 AliMUONVGeometryBuilder* builder
151 = (AliMUONVGeometryBuilder*)fGeometryBuilders->At(i);
1e8fff9c 152
acc838fd 153 // Create geometry with each builder
154 if (builder) {
155 builder->CreateGeometry();
156 builder->SetTransformations();
157 }
158 }
3c084d9f 159
acc838fd 160 for (Int_t j=0; j<AliMUONConstants::NCh(); j++) {
3c084d9f 161
acc838fd 162 AliMUONChamberGeometry* geometry = Chamber(j).GetGeometry();
3c084d9f 163
acc838fd 164 if (!geometry) continue;
165 // Skip chambers with not defined geometry
166
167 // Loop over envelopes
168 const TObjArray* kEnvelopes = geometry->GetEnvelopes();
169 for (Int_t k=0; k<kEnvelopes->GetEntriesFast(); k++) {
170
171 // Get envelope
172 AliMUONGeometryEnvelope* env = (AliMUONGeometryEnvelope*)kEnvelopes->At(k);
173 const TGeoCombiTrans* kEnvTrans = env->GetTransformation();
6b82c1f0 174 const char* only = "ONLY";
175 if (env->IsMANY()) only = "MANY";
acc838fd 176
177 if (env->IsVirtual() && env->GetConstituents()->GetEntriesFast() == 0 ) {
178 // virtual envelope + nof constituents = 0
179 // => not allowed;
180 // empty virtual envelope has no sense
181 Fatal("CreateGeometry", "Virtual envelope must have constituents.");
182 return;
183 }
3c084d9f 184
acc838fd 185 if (!env->IsVirtual() && env->GetConstituents()->GetEntriesFast() > 0 ) {
186 // non virtual envelope + nof constituents > 0
187 // => not allowed;
188 // use VMC to place constituents
189 Fatal("CreateGeometry", "Non virtual envelope cannot have constituents.");
190 return;
191 }
3c084d9f 192
acc838fd 193 if (!env->IsVirtual() && env->GetConstituents()->GetEntriesFast() == 0 ) {
194 // non virtual envelope + nof constituents = 0
195 // => place envelope in ALICE by composed transformation:
196 // Tglobal * Tch * Tenv
197
198 // Compound chamber transformation with the envelope one
199 TGeoHMatrix total
200 = (*fGlobalTransformation) *
201 (*geometry->GetTransformation()) *
202 (*kEnvTrans);
203 PlaceVolume(env->GetName(), geometry->GetMotherVolume(),
6b82c1f0 204 env->GetCopyNo(), total, 0, 0, only);
acc838fd 205 }
3c084d9f 206
acc838fd 207 if (env->IsVirtual() && env->GetConstituents()->GetEntriesFast() > 0 ) {
208 // virtual envelope + nof constituents > 0
209 // => do not place envelope and place constituents
210 // in ALICE by composed transformation:
211 // Tglobal * Tch * Tenv * Tconst
212
213 for (Int_t l=0; l<env->GetConstituents()->GetEntriesFast(); l++) {
214 AliMUONGeometryConstituent* constituent
215 = (AliMUONGeometryConstituent*)env->GetConstituents()->At(l);
216
217 // Compound chamber transformation with the envelope one + the constituent one
218 TGeoHMatrix total
219 = (*fGlobalTransformation) *
220 (*geometry->GetTransformation()) *
221 (*kEnvTrans) *
222 (*constituent->GetTransformation());
223
224 PlaceVolume(constituent->GetName(), geometry->GetMotherVolume(),
225 constituent->GetCopyNo(), total,
6b82c1f0 226 constituent->GetNpar(), constituent->GetParam(), only);
acc838fd 227 }
228 }
229 }
230 }
231}
232//__________________________________________________________________
233Int_t AliMUONv1::GetChamberId(Int_t volId) const
234{
235// Check if the volume with specified volId is a sensitive volume (gas)
236// of some chamber and returns the chamber number;
237// if not sensitive volume - return 0.
238// ---
a9e2aefa 239
acc838fd 240/*
241 for (Int_t i = 1; i <= AliMUONConstants::NCh(); i++)
242 if (volId==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()) return i;
243*/
244 for (Int_t i = 1; i <= AliMUONConstants::NCh(); i++)
245 if ( ((AliMUONChamber*)(*fChambers)[i-1])->IsSensId(volId) ) return i;
a9e2aefa 246
acc838fd 247 return 0;
a9e2aefa 248}
acc838fd 249//________________________________________________________________
a9e2aefa 250void AliMUONv1::CreateMaterials()
251{
acc838fd 252
a9e2aefa 253 // *** DEFINITION OF AVAILABLE MUON MATERIALS ***
254 //
b64652f5 255 // Ar-CO2 gas (80%+20%)
a9e2aefa 256 Float_t ag1[3] = { 39.95,12.01,16. };
257 Float_t zg1[3] = { 18.,6.,8. };
258 Float_t wg1[3] = { .8,.0667,.13333 };
259 Float_t dg1 = .001821;
260 //
261 // Ar-buthane-freon gas -- trigger chambers
262 Float_t atr1[4] = { 39.95,12.01,1.01,19. };
263 Float_t ztr1[4] = { 18.,6.,1.,9. };
264 Float_t wtr1[4] = { .56,.1262857,.2857143,.028 };
265 Float_t dtr1 = .002599;
266 //
267 // Ar-CO2 gas
268 Float_t agas[3] = { 39.95,12.01,16. };
269 Float_t zgas[3] = { 18.,6.,8. };
270 Float_t wgas[3] = { .74,.086684,.173316 };
271 Float_t dgas = .0018327;
272 //
273 // Ar-Isobutane gas (80%+20%) -- tracking
274 Float_t ag[3] = { 39.95,12.01,1.01 };
275 Float_t zg[3] = { 18.,6.,1. };
276 Float_t wg[3] = { .8,.057,.143 };
277 Float_t dg = .0019596;
278 //
279 // Ar-Isobutane-Forane-SF6 gas (49%+7%+40%+4%) -- trigger
280 Float_t atrig[5] = { 39.95,12.01,1.01,19.,32.066 };
281 Float_t ztrig[5] = { 18.,6.,1.,9.,16. };
282 Float_t wtrig[5] = { .49,1.08,1.5,1.84,0.04 };
283 Float_t dtrig = .0031463;
284 //
285 // bakelite
286
287 Float_t abak[3] = {12.01 , 1.01 , 16.};
288 Float_t zbak[3] = {6. , 1. , 8.};
289 Float_t wbak[3] = {6. , 6. , 1.};
290 Float_t dbak = 1.4;
291
292 Float_t epsil, stmin, deemax, tmaxfd, stemax;
293
294 Int_t iSXFLD = gAlice->Field()->Integ();
295 Float_t sXMGMX = gAlice->Field()->Max();
296 //
297 // --- Define the various materials for GEANT ---
298 AliMaterial(9, "ALUMINIUM$", 26.98, 13., 2.7, 8.9, 37.2);
299 AliMaterial(10, "ALUMINIUM$", 26.98, 13., 2.7, 8.9, 37.2);
97ddb489 300 // Air
301 Float_t aAir[4]={12.0107,14.0067,15.9994,39.948};
302 Float_t zAir[4]={6.,7.,8.,18.};
303 Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};
304 Float_t dAir = 1.20479E-3;
305 AliMixture(15, "AIR$ ", aAir, zAir, dAir,4, wAir);
306 // AliMaterial(15, "AIR$ ", 14.61, 7.3, .001205, 30423.24, 67500);
a9e2aefa 307 AliMixture(19, "Bakelite$", abak, zbak, dbak, -3, wbak);
308 AliMixture(20, "ArC4H10 GAS$", ag, zg, dg, 3, wg);
309 AliMixture(21, "TRIG GAS$", atrig, ztrig, dtrig, -5, wtrig);
310 AliMixture(22, "ArCO2 80%$", ag1, zg1, dg1, 3, wg1);
311 AliMixture(23, "Ar-freon $", atr1, ztr1, dtr1, 4, wtr1);
312 AliMixture(24, "ArCO2 GAS$", agas, zgas, dgas, 3, wgas);
1e8fff9c 313 // materials for slat:
314 // Sensitive area: gas (already defined)
315 // PCB: copper
316 // insulating material and frame: vetronite
317 // walls: carbon, rohacell, carbon
318 Float_t aglass[5]={12.01, 28.09, 16., 10.8, 23.};
319 Float_t zglass[5]={ 6., 14., 8., 5., 11.};
320 Float_t wglass[5]={ 0.5, 0.105, 0.355, 0.03, 0.01};
321 Float_t dglass=1.74;
322
323 // rohacell: C9 H13 N1 O2
324 Float_t arohac[4] = {12.01, 1.01, 14.010, 16.};
325 Float_t zrohac[4] = { 6., 1., 7., 8.};
326 Float_t wrohac[4] = { 9., 13., 1., 2.};
327 Float_t drohac = 0.03;
328
329 AliMaterial(31, "COPPER$", 63.54, 29., 8.96, 1.4, 0.);
330 AliMixture(32, "Vetronite$",aglass, zglass, dglass, 5, wglass);
331 AliMaterial(33, "Carbon$", 12.01, 6., 2.265, 18.8, 49.9);
332 AliMixture(34, "Rohacell$", arohac, zrohac, drohac, -4, wrohac);
333
a9e2aefa 334
335 epsil = .001; // Tracking precision,
336 stemax = -1.; // Maximum displacement for multiple scat
337 tmaxfd = -20.; // Maximum angle due to field deflection
338 deemax = -.3; // Maximum fractional energy loss, DLS
339 stmin = -.8;
340 //
341 // Air
342 AliMedium(1, "AIR_CH_US ", 15, 1, iSXFLD, sXMGMX, tmaxfd, stemax, deemax, epsil, stmin);
343 //
344 // Aluminum
345
346 AliMedium(4, "ALU_CH_US ", 9, 0, iSXFLD, sXMGMX, tmaxfd, fMaxStepAlu,
347 fMaxDestepAlu, epsil, stmin);
348 AliMedium(5, "ALU_CH_US ", 10, 0, iSXFLD, sXMGMX, tmaxfd, fMaxStepAlu,
349 fMaxDestepAlu, epsil, stmin);
350 //
351 // Ar-isoC4H10 gas
352
353 AliMedium(6, "AR_CH_US ", 20, 1, iSXFLD, sXMGMX, tmaxfd, fMaxStepGas,
354 fMaxDestepGas, epsil, stmin);
355//
356 // Ar-Isobuthane-Forane-SF6 gas
357
358 AliMedium(7, "GAS_CH_TRIGGER ", 21, 1, iSXFLD, sXMGMX, tmaxfd, stemax, deemax, epsil, stmin);
359
360 AliMedium(8, "BAKE_CH_TRIGGER ", 19, 0, iSXFLD, sXMGMX, tmaxfd, fMaxStepAlu,
361 fMaxDestepAlu, epsil, stmin);
362
363 AliMedium(9, "ARG_CO2 ", 22, 1, iSXFLD, sXMGMX, tmaxfd, fMaxStepGas,
364 fMaxDestepAlu, epsil, stmin);
1e8fff9c 365 // tracking media for slats: check the parameters!!
366 AliMedium(11, "PCB_COPPER ", 31, 0, iSXFLD, sXMGMX, tmaxfd,
367 fMaxStepAlu, fMaxDestepAlu, epsil, stmin);
368 AliMedium(12, "VETRONITE ", 32, 0, iSXFLD, sXMGMX, tmaxfd,
369 fMaxStepAlu, fMaxDestepAlu, epsil, stmin);
370 AliMedium(13, "CARBON ", 33, 0, iSXFLD, sXMGMX, tmaxfd,
371 fMaxStepAlu, fMaxDestepAlu, epsil, stmin);
372 AliMedium(14, "Rohacell ", 34, 0, iSXFLD, sXMGMX, tmaxfd,
373 fMaxStepAlu, fMaxDestepAlu, epsil, stmin);
acc838fd 374
375
376
377 //.Materials specific to stations
378 // created via builders
379
380 for (Int_t i=0; i<fGeometryBuilders->GetEntriesFast(); i++) {
381
382 // Get the builder
383 AliMUONVGeometryBuilder* builder
384 = (AliMUONVGeometryBuilder*)fGeometryBuilders->At(i);
385
386 // Create materials with each builder
387 if (builder) builder->CreateMaterials();
388 }
a9e2aefa 389}
390
acc838fd 391//______________________________________________________________________________
392void AliMUONv1::PlaceVolume(const TString& name, const TString& mName,
393 Int_t copyNo, const TGeoHMatrix& matrix,
6b82c1f0 394 Int_t npar, Double_t* param, const char* only) const
acc838fd 395{
396// Place the volume specified by name with the given transformation matrix
397// ---
398
399 // Do not apply global transformation
400 // if mother volume == DDIP
401 // (as it is applied on this volume)
402 TGeoHMatrix transform(matrix);
403 if (mName == TString("DDIP")) {
404 transform = (*fGlobalTransformation) * transform;
405 // To be changed to (*fGlobalTransformation).inverse()
406 // when available in TGeo
407 // To make this correct also for a general case when
408 // (*fGlobalTransformation) * *fGlobalTransformation) != 1
409 }
410
411 // Decompose transformation
412 const Double_t* xyz = transform.GetTranslation();
413 const Double_t* rm = transform.GetRotationMatrix();
414
415 //cout << "Got translation: "
416 // << xyz[0] << " " << xyz[1] << " " << xyz[2] << endl;
417
418 //cout << "Got rotation: "
419 // << rm[0] << " " << rm[1] << " " << rm[2] << endl
420 // << rm[3] << " " << rm[4] << " " << rm[5] << endl
421 // << rm[6] << " " << rm[7] << " " << rm[8] << endl;
422
423 // Check for presence of rotation
424 // (will be nice to be available in TGeo)
425 const Double_t kTolerance = 1e-04;
426 Bool_t isRotation = true;
427 if (TMath::Abs(rm[0] - 1.) < kTolerance &&
428 TMath::Abs(rm[1] - 0.) < kTolerance &&
429 TMath::Abs(rm[2] - 0.) < kTolerance &&
430 TMath::Abs(rm[3] - 0.) < kTolerance &&
431 TMath::Abs(rm[4] - 1.) < kTolerance &&
432 TMath::Abs(rm[5] - 0.) < kTolerance &&
433 TMath::Abs(rm[6] - 0.) < kTolerance &&
434 TMath::Abs(rm[7] - 0.) < kTolerance &&
435 TMath::Abs(rm[8] - 1.) < kTolerance) isRotation = false;
436
437 Int_t krot = 0;
438 if (isRotation) {
439 TGeoRotation rot;
440 rot.SetMatrix(const_cast<Double_t*>(transform.GetRotationMatrix()));
441 Double_t theta1, phi1, theta2, phi2, theta3, phi3;
442 rot.GetAngles(theta1, phi1, theta2, phi2, theta3, phi3);
443
444 //cout << "angles: "
445 // << theta1 << " " << phi1 << " "
446 // << theta2 << " " << phi2 << " "
447 // << theta3 << " " << phi3 << endl;
448
449 AliMatrix(krot, theta1, phi1, theta2, phi2, theta3, phi3);
450 }
451
452 // Place the volume in ALIC
6b82c1f0 453 if (npar == 0)
454 gMC->Gspos(name, copyNo, mName, xyz[0], xyz[1], xyz[2] , krot, only);
acc838fd 455 else
6b82c1f0 456 gMC->Gsposp(name, copyNo, mName, xyz[0], xyz[1], xyz[2] , krot, only,
acc838fd 457 param, npar);
458
459}
a9e2aefa 460
acc838fd 461//___________________________________________
a9e2aefa 462void AliMUONv1::Init()
463{
464 //
465 // Initialize Tracking Chambers
466 //
467
9e1a0ddb 468 if(fDebug) printf("\n%s: Start Init for version 1 - CPC chamber type\n\n",ClassName());
e17592e9 469 Int_t i;
f665c1ea 470 for (i=0; i<AliMUONConstants::NCh(); i++) {
a9e2aefa 471 ( (AliMUONChamber*) (*fChambers)[i])->Init();
472 }
473
474 //
475 // Set the chamber (sensitive region) GEANT identifier
acc838fd 476 //
477 for (Int_t i=0; i<fGeometryBuilders->GetEntriesFast(); i++) {
478
479 // Get the builder
480 AliMUONVGeometryBuilder* builder
481 = (AliMUONVGeometryBuilder*)fGeometryBuilders->At(i);
482
483 // Set sesitive volumes with each builder
484 if (builder) builder->SetSensitiveVolumes();
485 }
486
487/*
488 //
489 // Set the chamber (sensitive region) GEANT identifier
b74f1c6a 490 ((AliMUONChamber*)(*fChambers)[0])->SetGid(gMC->VolId("S01G"));
491 ((AliMUONChamber*)(*fChambers)[1])->SetGid(gMC->VolId("S02G"));
b17c0c87 492
b74f1c6a 493 ((AliMUONChamber*)(*fChambers)[2])->SetGid(gMC->VolId("S03G"));
494 ((AliMUONChamber*)(*fChambers)[3])->SetGid(gMC->VolId("S04G"));
b17c0c87 495
1e8fff9c 496 ((AliMUONChamber*)(*fChambers)[4])->SetGid(gMC->VolId("S05G"));
497 ((AliMUONChamber*)(*fChambers)[5])->SetGid(gMC->VolId("S06G"));
b17c0c87 498
1e8fff9c 499 ((AliMUONChamber*)(*fChambers)[6])->SetGid(gMC->VolId("S07G"));
500 ((AliMUONChamber*)(*fChambers)[7])->SetGid(gMC->VolId("S08G"));
b17c0c87 501
1e8fff9c 502 ((AliMUONChamber*)(*fChambers)[8])->SetGid(gMC->VolId("S09G"));
503 ((AliMUONChamber*)(*fChambers)[9])->SetGid(gMC->VolId("S10G"));
b17c0c87 504
b74f1c6a 505 ((AliMUONChamber*)(*fChambers)[10])->SetGid(gMC->VolId("SG1A"));
506 ((AliMUONChamber*)(*fChambers)[11])->SetGid(gMC->VolId("SG2A"));
507 ((AliMUONChamber*)(*fChambers)[12])->SetGid(gMC->VolId("SG3A"));
508 ((AliMUONChamber*)(*fChambers)[13])->SetGid(gMC->VolId("SG4A"));
acc838fd 509*/
9e1a0ddb 510 if(fDebug) printf("\n%s: Finished Init for version 1 - CPC chamber type\n",ClassName());
a9e2aefa 511
512 //cp
9e1a0ddb 513 if(fDebug) printf("\n%s: Start Init for Trigger Circuits\n",ClassName());
f665c1ea 514 for (i=0; i<AliMUONConstants::NTriggerCircuit(); i++) {
a9e2aefa 515 ( (AliMUONTriggerCircuit*) (*fTriggerCircuits)[i])->Init(i);
516 }
9e1a0ddb 517 if(fDebug) printf("%s: Finished Init for Trigger Circuits\n",ClassName());
a9e2aefa 518 //cp
a9e2aefa 519}
1391e633 520
521//_______________________________________________________________________________
c33d9661 522void AliMUONv1::StepManager()
523{
d7c4fbc4 524 // Stepmanager for the chambers
525
1391e633 526 if (fStepManagerVersionOld) {
c33d9661 527 StepManagerOld();
528 return;
529 }
c33d9661 530
531 // Only charged tracks
532 if( !(gMC->TrackCharge()) ) return;
1391e633 533 // Only charged tracks
534
c33d9661 535 // Only gas gap inside chamber
536 // Tag chambers and record hits when track enters
acc838fd 537 static Int_t idvol=-1;
1391e633 538 Int_t iChamber=0;
539 Int_t id=0;
540 Int_t copy;
541 const Float_t kBig = 1.e10;
542
acc838fd 543
544 //
545 // Only gas gap inside chamber
546 // Tag chambers and record hits when track enters
c33d9661 547 id=gMC->CurrentVolID(copy);
acc838fd 548 iChamber = GetChamberId(id);
549 idvol = iChamber -1;
550
551 if (idvol == -1) return;
abaf7c9d 552
967a2a25 553 // Filling TrackRefs file for MUON. Our Track references are the active volume of the chambers
554 if ( (gMC->IsTrackEntering() || gMC->IsTrackExiting() ) )
555 AddTrackReference(gAlice->GetMCApp()->GetCurrentTrackNumber());
556
1391e633 557 if( gMC->IsTrackEntering() ) {
558 Float_t theta = fTrackMomentum.Theta();
4ac9d21e 559 if ((TMath::Pi()-theta)*kRaddeg>=15.) gMC->SetMaxStep(fStepMaxInActiveGas); // We use Pi-theta because z is negative
560 }
abaf7c9d 561
1391e633 562// if (GetDebug()) {
563// Float_t z = ( (AliMUONChamber*)(*fChambers)[idvol])->Z() ;
564// Info("StepManager Step","Active volume found %d chamber %d Z chamber is %f ",idvol,iChamber, z);
565// }
566 // Particule id and mass,
567 Int_t ipart = gMC->TrackPid();
568 Float_t mass = gMC->TrackMass();
569
570 fDestepSum[idvol]+=gMC->Edep();
571 // Get current particle id (ipart), track position (pos) and momentum (mom)
572 if ( fStepSum[idvol]==0.0 ) gMC->TrackMomentum(fTrackMomentum);
573 fStepSum[idvol]+=gMC->TrackStep();
abaf7c9d 574
1391e633 575// if (GetDebug()) {
576// Info("StepManager Step","iChamber %d, Particle %d, theta %f phi %f mass %f StepSum %f eloss %g",
acc838fd 577// iChamber,ipart, fTrackMomentum.Theta()*kRaddeg, fTrackMomentum.Phi()*kRaddeg, mass, fStepSum[idvol], gMC->Edep());
1391e633 578// Info("StepManager Step","Track Momentum %f %f %f", fTrackMomentum.X(), fTrackMomentum.Y(), fTrackMomentum.Z()) ;
579// gMC->TrackPosition(fTrackPosition);
580// Info("StepManager Step","Track Position %f %f %f",fTrackPosition.X(),fTrackPosition.Y(),fTrackPosition.Z()) ;
d08aff2d 581// }
1391e633 582
583 // Track left chamber or StepSum larger than fStepMaxInActiveGas
584 if ( gMC->IsTrackExiting() ||
585 gMC->IsTrackStop() ||
586 gMC->IsTrackDisappeared()||
587 (fStepSum[idvol]>fStepMaxInActiveGas) ) {
588
589 if ( gMC->IsTrackExiting() ||
acc838fd 590 gMC->IsTrackStop() ||
591 gMC->IsTrackDisappeared() ) gMC->SetMaxStep(kBig);
1391e633 592
593 gMC->TrackPosition(fTrackPosition);
594 Float_t theta = fTrackMomentum.Theta();
595 Float_t phi = fTrackMomentum.Phi();
596
d7c4fbc4 597 TLorentzVector backToWire( fStepSum[idvol]/2.*sin(theta)*cos(phi),
acc838fd 598 fStepSum[idvol]/2.*sin(theta)*sin(phi),
599 fStepSum[idvol]/2.*cos(theta),0.0 );
1391e633 600 // if (GetDebug())
601 // Info("StepManager Exit","Track Position %f %f %f",fTrackPosition.X(),fTrackPosition.Y(),fTrackPosition.Z()) ;
602 // if (GetDebug())
d7c4fbc4 603 // Info("StepManager Exit ","Track backToWire %f %f %f",backToWire.X(),backToWire.Y(),backToWire.Z()) ;
604 fTrackPosition-=backToWire;
1391e633 605
606 //-------------- Angle effect
607 // Ratio between energy loss of particle and Mip as a function of BetaGamma of particle (Energy/Mass)
608
30178c30 609 Float_t betaxGamma = fTrackMomentum.P()/mass;// pc/mc2
d7c4fbc4 610 Float_t sigmaEffect10degrees;
611 Float_t sigmaEffectThetadegrees;
612 Float_t eLossParticleELossMip;
613 Float_t yAngleEffect=0.;
614 Float_t thetawires = TMath::Abs( TMath::ASin( TMath::Sin(TMath::Pi()-theta) * TMath::Sin(phi) ) );// We use Pi-theta because z is negative
615
374ebd7d 616
617 if (fAngleEffect){
30178c30 618 if ( (betaxGamma >3.2) && (thetawires*kRaddeg<=15.) ) {
619 betaxGamma=TMath::Log(betaxGamma);
620 eLossParticleELossMip = fElossRatio->Eval(betaxGamma);
1391e633 621 // 10 degrees is a reference for a model (arbitrary)
d7c4fbc4 622 sigmaEffect10degrees=fAngleEffect10->Eval(eLossParticleELossMip);// in micrometers
1391e633 623 // Angle with respect to the wires assuming that chambers are perpendicular to the z axis.
d7c4fbc4 624 sigmaEffectThetadegrees = sigmaEffect10degrees/fAngleEffectNorma->Eval(thetawires*kRaddeg); // For 5mm gap
1391e633 625 if ( (iChamber==1) || (iChamber==2) )
acc838fd 626 sigmaEffectThetadegrees/=(1.09833e+00+1.70000e-02*(thetawires*kRaddeg)); // The gap is different (4mm)
d7c4fbc4 627 yAngleEffect=1.e-04*gRandom->Gaus(0,sigmaEffectThetadegrees); // Error due to the angle effect in cm
1391e633 628 }
374ebd7d 629 }
1391e633 630
631 // One hit per chamber
47464038 632 GetMUONData()->AddHit(fIshunt,
633 gAlice->GetMCApp()->GetCurrentTrackNumber(),
634 iChamber, ipart,
635 fTrackPosition.X(),
636 fTrackPosition.Y()+yAngleEffect,
637 fTrackPosition.Z(),
638 gMC->TrackTime(),
639 fTrackMomentum.P(),
640 theta,
641 phi,
642 fStepSum[idvol],
643 fDestepSum[idvol],
644 fTrackPosition.X(),
645 fTrackPosition.Y(),
646 fTrackPosition.Z());
647
1391e633 648// if (GetDebug()){
649// Info("StepManager Exit","Particle exiting from chamber %d",iChamber);
650// Info("StepManager Exit","StepSum %f eloss geant %g ",fStepSum[idvol],fDestepSum[idvol]);
651// Info("StepManager Exit","Track Position %f %f %f",fTrackPosition.X(),fTrackPosition.Y(),fTrackPosition.Z()) ;
652// }
653 fStepSum[idvol] =0; // Reset for the next event
654 fDestepSum[idvol]=0; // Reset for the next event
655 }
abaf7c9d 656}
5f91c9e8 657
4ac9d21e 658//__________________________________________
c33d9661 659void AliMUONv1::StepManagerOld()
a9e2aefa 660{
d7c4fbc4 661 // Old Stepmanager for the chambers
a9e2aefa 662 Int_t copy, id;
acc838fd 663 static Int_t idvol =-1;
a9e2aefa 664 static Int_t vol[2];
665 Int_t ipart;
666 TLorentzVector pos;
667 TLorentzVector mom;
668 Float_t theta,phi;
669 Float_t destep, step;
abaf7c9d 670
d7c4fbc4 671 static Float_t sstep;
1e8fff9c 672 static Float_t eloss, eloss2, xhit, yhit, zhit, tof, tlength;
2eb55fab 673 const Float_t kBig = 1.e10;
a9e2aefa 674 static Float_t hits[15];
675
676 TClonesArray &lhits = *fHits;
677
678 //
a9e2aefa 679 //
680 // Only charged tracks
681 if( !(gMC->TrackCharge()) ) return;
682 //
683 // Only gas gap inside chamber
684 // Tag chambers and record hits when track enters
a9e2aefa 685 id=gMC->CurrentVolID(copy);
5f91c9e8 686 vol[0] = GetChamberId(id);
687 idvol = vol[0] -1;
688
689 if (idvol == -1) return;
690
a9e2aefa 691 //
692 // Get current particle id (ipart), track position (pos) and momentum (mom)
693 gMC->TrackPosition(pos);
694 gMC->TrackMomentum(mom);
695
696 ipart = gMC->TrackPid();
a9e2aefa 697
698 //
699 // momentum loss and steplength in last step
700 destep = gMC->Edep();
701 step = gMC->TrackStep();
abaf7c9d 702 // cout<<"------------"<<step<<endl;
a9e2aefa 703 //
704 // record hits when track enters ...
705 if( gMC->IsTrackEntering()) {
abaf7c9d 706
a9e2aefa 707 gMC->SetMaxStep(fMaxStepGas);
708 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
709 Double_t rt = TMath::Sqrt(tc);
710 Double_t pmom = TMath::Sqrt(tc+mom[2]*mom[2]);
2eb55fab 711 Double_t tx = mom[0]/pmom;
712 Double_t ty = mom[1]/pmom;
713 Double_t tz = mom[2]/pmom;
714 Double_t s = ((AliMUONChamber*)(*fChambers)[idvol])
acc838fd 715 ->ResponseModel()
716 ->Pitch()/tz;
a9e2aefa 717 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
718 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
719 hits[0] = Float_t(ipart); // Geant3 particle type
2eb55fab 720 hits[1] = pos[0]+s*tx; // X-position for hit
721 hits[2] = pos[1]+s*ty; // Y-position for hit
722 hits[3] = pos[2]+s*tz; // Z-position for hit
a9e2aefa 723 hits[4] = theta; // theta angle of incidence
724 hits[5] = phi; // phi angle of incidence
ce3f5e87 725 hits[8] = 0;//PadHits does not exist anymore (Float_t) fNPadHits; // first padhit
a9e2aefa 726 hits[9] = -1; // last pad hit
2eb55fab 727 hits[10] = mom[3]; // hit momentum P
728 hits[11] = mom[0]; // Px
729 hits[12] = mom[1]; // Py
730 hits[13] = mom[2]; // Pz
a9e2aefa 731 tof=gMC->TrackTime();
2eb55fab 732 hits[14] = tof; // Time of flight
733 tlength = 0;
734 eloss = 0;
735 eloss2 = 0;
d7c4fbc4 736 sstep=0;
2eb55fab 737 xhit = pos[0];
738 yhit = pos[1];
739 zhit = pos[2];
681d067b 740 Chamber(idvol).ChargeCorrelationInit();
a9e2aefa 741 // Only if not trigger chamber
1e8fff9c 742
abaf7c9d 743// printf("---------------------------\n");
744// printf(">>>> Y = %f \n",hits[2]);
745// printf("---------------------------\n");
746
1e8fff9c 747
748
abaf7c9d 749 // if(idvol < AliMUONConstants::NTrackingCh()) {
acc838fd 750// //
751// // Initialize hit position (cursor) in the segmentation model
752// ((AliMUONChamber*) (*fChambers)[idvol])
753// ->SigGenInit(pos[0], pos[1], pos[2]);
abaf7c9d 754// } else {
acc838fd 755// //geant3->Gpcxyz();
756// //printf("In the Trigger Chamber #%d\n",idvol-9);
abaf7c9d 757// }
a9e2aefa 758 }
759 eloss2+=destep;
d7c4fbc4 760 sstep+=step;
abaf7c9d 761
d7c4fbc4 762 // cout<<sstep<<endl;
abaf7c9d 763
a9e2aefa 764 //
765 // Calculate the charge induced on a pad (disintegration) in case
766 //
767 // Mip left chamber ...
768 if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
769 gMC->SetMaxStep(kBig);
770 eloss += destep;
771 tlength += step;
772
802a864d 773 Float_t x0,y0,z0;
774 Float_t localPos[3];
775 Float_t globalPos[3] = {pos[0], pos[1], pos[2]};
802a864d 776 gMC->Gmtod(globalPos,localPos,1);
777
2eb55fab 778 if(idvol < AliMUONConstants::NTrackingCh()) {
a9e2aefa 779// tracking chambers
acc838fd 780 x0 = 0.5*(xhit+pos[0]);
781 y0 = 0.5*(yhit+pos[1]);
782 z0 = 0.5*(zhit+pos[2]);
a9e2aefa 783 } else {
784// trigger chambers
acc838fd 785 x0 = xhit;
786 y0 = yhit;
787 z0 = 0.;
a9e2aefa 788 }
789
1e8fff9c 790
ce3f5e87 791 // if (eloss >0) MakePadHits(x0,y0,z0,eloss,tof,idvol);
a9e2aefa 792
acc838fd 793
2eb55fab 794 hits[6] = tlength; // track length
795 hits[7] = eloss2; // de/dx energy loss
796
abaf7c9d 797
ce3f5e87 798 // if (fNPadHits > (Int_t)hits[8]) {
acc838fd 799 // hits[8] = hits[8]+1;
800 // hits[9] = 0: // PadHits does not exist anymore (Float_t) fNPadHits;
ce3f5e87 801 //}
2eb55fab 802//
803// new hit
804
a9e2aefa 805 new(lhits[fNhits++])
acc838fd 806 AliMUONHit(fIshunt, gAlice->GetMCApp()->GetCurrentTrackNumber(), vol,hits);
a9e2aefa 807 eloss = 0;
808 //
809 // Check additional signal generation conditions
810 // defined by the segmentation
a75f073c 811 // model (boundary crossing conditions)
812 // only for tracking chambers
a9e2aefa 813 } else if
a75f073c 814 ((idvol < AliMUONConstants::NTrackingCh()) &&
815 ((AliMUONChamber*) (*fChambers)[idvol])->SigGenCond(pos[0], pos[1], pos[2]))
a9e2aefa 816 {
817 ((AliMUONChamber*) (*fChambers)[idvol])
acc838fd 818 ->SigGenInit(pos[0], pos[1], pos[2]);
802a864d 819
820 Float_t localPos[3];
821 Float_t globalPos[3] = {pos[0], pos[1], pos[2]};
822 gMC->Gmtod(globalPos,localPos,1);
823
e0f71fb7 824 eloss += destep;
802a864d 825
ce3f5e87 826 // if (eloss > 0 && idvol < AliMUONConstants::NTrackingCh())
acc838fd 827 // MakePadHits(0.5*(xhit+pos[0]),0.5*(yhit+pos[1]),pos[2],eloss,tof,idvol);
a9e2aefa 828 xhit = pos[0];
829 yhit = pos[1];
e0f71fb7 830 zhit = pos[2];
831 eloss = 0;
a9e2aefa 832 tlength += step ;
833 //
834 // nothing special happened, add up energy loss
835 } else {
836 eloss += destep;
837 tlength += step ;
838 }
839}