Line data Source code
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 : /* $Id$ */
17 :
18 : ////////////////////////////////////////////////////////////////////////////
19 : // //
20 : // Transition Radiation Detector version 1 -- slow simulator //
21 : // //
22 : ////////////////////////////////////////////////////////////////////////////
23 :
24 : #include <TLorentzVector.h>
25 : #include <TMath.h>
26 : #include <TRandom.h>
27 : #include <TVirtualMC.h>
28 : #include <TGeoManager.h>
29 : #include <TGeoMatrix.h>
30 : #include <TGeoPhysicalNode.h>
31 :
32 : #include "AliTrackReference.h"
33 : #include "AliMC.h"
34 : #include "AliRun.h"
35 : #include "AliGeomManager.h"
36 :
37 : #include "AliTRDgeometry.h"
38 : #include "AliTRDCommonParam.h"
39 : #include "AliTRDsimTR.h"
40 : #include "AliTRDv1.h"
41 : #include "AliLog.h"
42 :
43 12 : ClassImp(AliTRDv1)
44 :
45 : //_____________________________________________________________________________
46 : AliTRDv1::AliTRDv1()
47 12 : :AliTRD()
48 12 : ,fTRon(kTRUE)
49 12 : ,fTR(NULL)
50 12 : ,fStepSize(0)
51 12 : ,fWion(0)
52 60 : {
53 : //
54 : // Default constructor
55 : //
56 :
57 24 : }
58 :
59 : //_____________________________________________________________________________
60 : AliTRDv1::AliTRDv1(const char *name, const char *title)
61 1 : :AliTRD(name,title)
62 1 : ,fTRon(kTRUE)
63 1 : ,fTR(NULL)
64 1 : ,fStepSize(0.1)
65 1 : ,fWion(0)
66 5 : {
67 : //
68 : // Standard constructor for Transition Radiation Detector version 1
69 : //
70 :
71 1 : SetBufferSize(128000);
72 :
73 2 : if (AliTRDCommonParam::Instance()->IsXenon()) {
74 1 : fWion = 23.53; // Ionization energy XeCO2 (85/15)
75 1 : }
76 0 : else if (AliTRDCommonParam::Instance()->IsArgon()) {
77 0 : fWion = 27.21; // Ionization energy ArCO2 (82/18)
78 : }
79 : else {
80 0 : AliFatal("Wrong gas mixture");
81 0 : exit(1);
82 : }
83 :
84 2 : }
85 :
86 : //_____________________________________________________________________________
87 : AliTRDv1::~AliTRDv1()
88 78 : {
89 : //
90 : // AliTRDv1 destructor
91 : //
92 :
93 13 : if (fTR) {
94 26 : delete fTR;
95 13 : fTR = 0;
96 13 : }
97 :
98 39 : }
99 :
100 : //_____________________________________________________________________________
101 : void AliTRDv1::AddAlignableVolumes() const
102 : {
103 : //
104 : // Create entries for alignable volumes associating the symbolic volume
105 : // name with the corresponding volume path. Needs to be syncronized with
106 : // eventual changes in the geometry.
107 : //
108 :
109 2 : TString volPath;
110 1 : TString symName;
111 :
112 1 : TString vpStr = "ALIC_1/B077_1/BSEGMO";
113 1 : TString vpApp1 = "_1/BTRD";
114 1 : TString vpApp2 = "_1";
115 1 : TString vpApp3a = "/UTR1_1/UTS1_1/UTI1_1/UT";
116 1 : TString vpApp3b = "/UTR2_1/UTS2_1/UTI2_1/UT";
117 1 : TString vpApp3c = "/UTR3_1/UTS3_1/UTI3_1/UT";
118 1 : TString vpApp3d = "/UTR4_1/UTS4_1/UTI4_1/UT";
119 :
120 1 : TString snStr = "TRD/sm";
121 1 : TString snApp1 = "/st";
122 1 : TString snApp2 = "/pl";
123 :
124 : //
125 : // The super modules
126 : // The symbolic names are: TRD/sm00
127 : // ...
128 : // TRD/sm17
129 : //
130 38 : for (Int_t isector = 0; isector < AliTRDgeometry::Nsector(); isector++) {
131 :
132 18 : volPath = vpStr;
133 18 : volPath += isector;
134 18 : volPath += vpApp1;
135 18 : volPath += isector;
136 18 : volPath += vpApp2;
137 :
138 18 : symName = snStr;
139 36 : symName += Form("%02d",isector);
140 :
141 54 : gGeoManager->SetAlignableEntry(symName.Data(),volPath.Data());
142 :
143 : }
144 :
145 : //
146 : // The readout chambers
147 : // The symbolic names are: TRD/sm00/st0/pl0
148 : // ...
149 : // TRD/sm17/st4/pl5
150 : //
151 : AliGeomManager::ELayerID idTRD1 = AliGeomManager::kTRD1;
152 : Int_t layer, modUID;
153 :
154 38 : for (Int_t isector = 0; isector < AliTRDgeometry::Nsector(); isector++) {
155 :
156 18 : if (fGeometry->GetSMstatus(isector) == 0) continue;
157 :
158 84 : for (Int_t istack = 0; istack < AliTRDgeometry::Nstack(); istack++) {
159 490 : for (Int_t ilayer = 0; ilayer < AliTRDgeometry::Nlayer(); ilayer++) {
160 :
161 210 : layer = idTRD1 + ilayer;
162 420 : modUID = AliGeomManager::LayerToVolUIDSafe(layer,isector*5+istack);
163 :
164 210 : Int_t idet = AliTRDgeometry::GetDetectorSec(ilayer,istack);
165 :
166 210 : volPath = vpStr;
167 210 : volPath += isector;
168 210 : volPath += vpApp1;
169 210 : volPath += isector;
170 210 : volPath += vpApp2;
171 210 : switch (isector) {
172 : case 17:
173 30 : if ((istack == 4) && (ilayer == 4)) {
174 1 : continue;
175 : }
176 29 : volPath += vpApp3d;
177 : break;
178 : case 13:
179 : case 14:
180 : case 15:
181 0 : if (istack == 2) {
182 0 : continue;
183 : }
184 0 : volPath += vpApp3c;
185 : break;
186 : case 11:
187 : case 12:
188 0 : volPath += vpApp3b;
189 : break;
190 : default:
191 180 : volPath += vpApp3a;
192 : };
193 418 : volPath += Form("%02d",idet);
194 209 : volPath += vpApp2;
195 :
196 209 : symName = snStr;
197 418 : symName += Form("%02d",isector);
198 209 : symName += snApp1;
199 209 : symName += istack;
200 209 : symName += snApp2;
201 209 : symName += ilayer;
202 :
203 : TGeoPNEntry *alignableEntry =
204 627 : gGeoManager->SetAlignableEntry(symName.Data(),volPath.Data(),modUID);
205 :
206 : // Add the tracking to local matrix following the TPC example
207 209 : if (alignableEntry) {
208 209 : TGeoHMatrix *globMatrix = alignableEntry->GetGlobalOrig();
209 209 : Double_t sectorAngle = 20.0 * (isector % 18) + 10.0;
210 418 : TGeoHMatrix *t2lMatrix = new TGeoHMatrix();
211 209 : t2lMatrix->RotateZ(sectorAngle);
212 418 : t2lMatrix->MultiplyLeft(&(globMatrix->Inverse()));
213 209 : alignableEntry->SetMatrix(t2lMatrix);
214 209 : }
215 : else {
216 0 : AliError(Form("Alignable entry %s is not valid!",symName.Data()));
217 : }
218 :
219 209 : }
220 : }
221 7 : }
222 :
223 1 : }
224 :
225 : //_____________________________________________________________________________
226 : void AliTRDv1::CreateGeometry()
227 : {
228 : //
229 : // Create the GEANT geometry for the Transition Radiation Detector - Version 1
230 : // This version covers the full azimuth.
231 : //
232 :
233 : // Check that FRAME is there otherwise we have no place where to put the TRD
234 2 : AliModule* frame = gAlice->GetModule("FRAME");
235 1 : if (!frame) {
236 0 : AliError("TRD needs FRAME to be present\n");
237 0 : return;
238 : }
239 :
240 : // Define the chambers
241 1 : AliTRD::CreateGeometry();
242 :
243 2 : }
244 :
245 : //_____________________________________________________________________________
246 : void AliTRDv1::CreateMaterials()
247 : {
248 : //
249 : // Create materials for the Transition Radiation Detector version 1
250 : //
251 :
252 2 : AliTRD::CreateMaterials();
253 :
254 1 : }
255 :
256 : //_____________________________________________________________________________
257 : void AliTRDv1::CreateTRhit(Int_t det)
258 : {
259 : //
260 : // Creates an electron cluster from a TR photon.
261 : // The photon is assumed to be created a the end of the radiator. The
262 : // distance after which it deposits its energy takes into account the
263 : // absorbtion of the entrance window and of the gas mixture in drift
264 : // volume.
265 : //
266 :
267 : // Maximum number of TR photons per track
268 : const Int_t kNTR = 50;
269 :
270 792 : TLorentzVector mom;
271 396 : TLorentzVector pos;
272 :
273 396 : Float_t eTR[kNTR];
274 396 : Int_t nTR;
275 :
276 : // Create TR photons
277 792 : TVirtualMC::GetMC()->TrackMomentum(mom);
278 792 : Float_t pTot = mom.Rho();
279 396 : fTR->CreatePhotons(11,pTot,nTR,eTR);
280 396 : if (nTR > kNTR) {
281 0 : AliFatal(Form("Boundary error: nTR = %d, kNTR = %d",nTR,kNTR));
282 : }
283 :
284 : // Loop through the TR photons
285 832 : for (Int_t iTR = 0; iTR < nTR; iTR++) {
286 :
287 20 : Float_t energyMeV = eTR[iTR] * 0.001;
288 20 : Float_t energyeV = eTR[iTR] * 1000.0;
289 : Float_t absLength = 0.0;
290 : Float_t sigma = 0.0;
291 :
292 : // Take the absorbtion in the entrance window into account
293 20 : Double_t muMy = fTR->GetMuMy(energyMeV);
294 20 : sigma = muMy * fFoilDensity;
295 20 : if (sigma > 0.0) {
296 40 : absLength = gRandom->Exp(1.0/sigma);
297 20 : if (absLength < AliTRDgeometry::MyThick()) {
298 1 : continue;
299 : }
300 : }
301 : else {
302 0 : continue;
303 : }
304 :
305 : // The absorbtion cross sections in the drift gas
306 : // Gas-mixture (Xe/CO2)
307 : Double_t muNo = 0.0;
308 38 : if (AliTRDCommonParam::Instance()->IsXenon()) {
309 19 : muNo = fTR->GetMuXe(energyMeV);
310 19 : }
311 0 : else if (AliTRDCommonParam::Instance()->IsArgon()) {
312 0 : muNo = fTR->GetMuAr(energyMeV);
313 0 : }
314 19 : Double_t muCO = fTR->GetMuCO(energyMeV);
315 38 : sigma = (fGasNobleFraction * muNo + (1.0 - fGasNobleFraction) * muCO)
316 19 : * fGasDensity
317 19 : * fTR->GetTemp();
318 :
319 : // The distance after which the energy of the TR photon
320 : // is deposited.
321 19 : if (sigma > 0.0) {
322 38 : absLength = gRandom->Exp(1.0/sigma);
323 38 : if (absLength > (AliTRDgeometry::DrThick()
324 19 : + AliTRDgeometry::AmThick())) {
325 0 : continue;
326 : }
327 : }
328 : else {
329 0 : continue;
330 : }
331 :
332 : // The position of the absorbtion
333 19 : Float_t posHit[3];
334 38 : TVirtualMC::GetMC()->TrackPosition(pos);
335 57 : posHit[0] = pos[0] + mom[0] / pTot * absLength;
336 57 : posHit[1] = pos[1] + mom[1] / pTot * absLength;
337 57 : posHit[2] = pos[2] + mom[2] / pTot * absLength;
338 :
339 : // Create the charge
340 19 : Int_t q = ((Int_t) (energyeV / fWion));
341 :
342 : // Add the hit to the array. TR photon hits are marked
343 : // by negative charge
344 38 : AddHit(gAlice->GetMCApp()->GetCurrentTrackNumber()
345 : ,det
346 19 : ,posHit
347 19 : ,-q
348 57 : ,TVirtualMC::GetMC()->TrackTime()*1.0e06
349 : ,kTRUE);
350 :
351 19 : }
352 :
353 396 : }
354 :
355 : //_____________________________________________________________________________
356 : void AliTRDv1::Init()
357 : {
358 : //
359 : // Initialise Transition Radiation Detector after geometry has been built.
360 : //
361 :
362 2 : AliTRD::Init();
363 :
364 3 : AliDebug(1,"Slow simulator\n");
365 :
366 : // Switch on TR simulation as default
367 1 : if (!fTRon) {
368 0 : AliInfo("TR simulation off");
369 0 : }
370 : else {
371 2 : fTR = new AliTRDsimTR();
372 : }
373 :
374 3 : AliDebug(1,"+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++");
375 :
376 1 : }
377 :
378 : //_____________________________________________________________________________
379 : void AliTRDv1::StepManager()
380 : {
381 : //
382 : // Slow simulator. Every charged track produces electron cluster as hits
383 : // along its path across the drift volume. The step size is fixed in
384 : // this version of the step manager.
385 : //
386 : // Works for Xe/CO2 as well as Ar/CO2
387 : //
388 :
389 : // PDG code electron
390 : const Int_t kPdgElectron = 11;
391 :
392 : Int_t layer = 0;
393 : Int_t stack = 0;
394 : Int_t sector = 0;
395 : Int_t det = 0;
396 : Int_t qTot;
397 :
398 448210 : Float_t hits[3];
399 : Double_t eDep;
400 :
401 : Bool_t drRegion = kFALSE;
402 : Bool_t amRegion = kFALSE;
403 :
404 224105 : TString cIdPath;
405 224105 : Char_t cIdSector[3];
406 224105 : cIdSector[2] = 0;
407 :
408 224105 : TString cIdCurrent;
409 224105 : TString cIdSensDr = "J";
410 224105 : TString cIdSensAm = "K";
411 224105 : Char_t cIdChamber[3];
412 224105 : cIdChamber[2] = 0;
413 :
414 224105 : TLorentzVector pos;
415 224105 : TLorentzVector mom;
416 :
417 224105 : const Int_t kNlayer = AliTRDgeometry::Nlayer();
418 224105 : const Int_t kNstack = AliTRDgeometry::Nstack();
419 224105 : const Int_t kNdetsec = kNlayer * kNstack;
420 :
421 : const Double_t kBig = 1.0e+12;
422 : const Float_t kEkinMinStep = 1.0e-5; // Minimum energy for the step size adjustment
423 :
424 : // Set the maximum step size to a very large number for all
425 : // neutral particles and those outside the driftvolume
426 672315 : if (!fPrimaryIonisation) TVirtualMC::GetMC()->SetMaxStep(kBig);
427 :
428 : // If not charged track or already stopped or disappeared, just return.
429 791524 : if ((!TVirtualMC::GetMC()->TrackCharge()) ||
430 238418 : TVirtualMC::GetMC()->IsTrackDisappeared()) {
431 104916 : return;
432 : }
433 :
434 : // Inside a sensitive volume?
435 357567 : cIdCurrent = TVirtualMC::GetMC()->CurrentVolName();
436 :
437 476756 : if (cIdSensDr == cIdCurrent[1]) {
438 : drRegion = kTRUE;
439 24282 : }
440 476756 : if (cIdSensAm == cIdCurrent[1]) {
441 : amRegion = kTRUE;
442 8516 : }
443 :
444 214096 : if ((!drRegion) &&
445 94907 : (!amRegion)) {
446 86391 : return;
447 : }
448 :
449 : // The hit coordinates and charge
450 65596 : TVirtualMC::GetMC()->TrackPosition(pos);
451 65596 : hits[0] = pos[0];
452 65596 : hits[1] = pos[1];
453 65596 : hits[2] = pos[2];
454 :
455 : // The sector number (0 - 17), according to standard coordinate system
456 65596 : cIdPath = gGeoManager->GetPath();
457 65596 : cIdSector[0] = cIdPath[21];
458 65596 : cIdSector[1] = cIdPath[22];
459 32798 : sector = atoi(cIdSector);
460 :
461 : // The plane and chamber number
462 65596 : cIdChamber[0] = cIdCurrent[2];
463 65596 : cIdChamber[1] = cIdCurrent[3];
464 65596 : Int_t idChamber = (atoi(cIdChamber) % kNdetsec);
465 32798 : stack = ((Int_t) idChamber / kNlayer);
466 32798 : layer = ((Int_t) idChamber % kNlayer);
467 :
468 : // The detector number
469 32798 : det = fGeometry->GetDetector(layer,stack,sector);
470 :
471 : // 0: InFlight 1:Entering 2:Exiting
472 : Int_t trkStat = 0;
473 :
474 : // Special hits only in the drift region
475 57080 : if ((drRegion) &&
476 48564 : (TVirtualMC::GetMC()->IsTrackEntering())) {
477 :
478 : // Create a track reference at the entrance of each
479 : // chamber that contains the momentum components of the particle
480 1040 : TVirtualMC::GetMC()->TrackMomentum(mom);
481 1040 : AddTrackReference(gAlice->GetMCApp()->GetCurrentTrackNumber(), AliTrackReference::kTRD);
482 : trkStat = 1;
483 :
484 : // Create the hits from TR photons if electron/positron is
485 : // entering the drift volume
486 1040 : if ((fTR) &&
487 520 : (fTRon) &&
488 1560 : (TMath::Abs(TVirtualMC::GetMC()->TrackPid()) == kPdgElectron)) {
489 396 : CreateTRhit(det);
490 : }
491 :
492 : }
493 40794 : else if ((amRegion) &&
494 17032 : (TVirtualMC::GetMC()->IsTrackExiting())) {
495 :
496 : // Create a track reference at the exit of each
497 : // chamber that contains the momentum components of the particle
498 970 : TVirtualMC::GetMC()->TrackMomentum(mom);
499 970 : AddTrackReference(gAlice->GetMCApp()->GetCurrentTrackNumber(), AliTrackReference::kTRD);
500 : trkStat = 2;
501 :
502 485 : }
503 :
504 : // Calculate the charge according to GEANT Edep
505 : // Create a new dEdx hit
506 98394 : eDep = TMath::Max(TVirtualMC::GetMC()->Edep(),0.0) * 1.0e+09;
507 32798 : qTot = (Int_t) (eDep / fWion);
508 65596 : if ((qTot) ||
509 32798 : (trkStat)) {
510 46412 : AddHit(gAlice->GetMCApp()->GetCurrentTrackNumber()
511 : ,det
512 23206 : ,hits
513 : ,qTot
514 69618 : ,TVirtualMC::GetMC()->TrackTime()*1.0e06
515 : ,drRegion);
516 : }
517 :
518 : // Set Maximum Step Size
519 : // Produce only one hit if Ekin is below cutoff
520 163990 : if ((TVirtualMC::GetMC()->Etot() - TVirtualMC::GetMC()->TrackMass()) < kEkinMinStep) {
521 93 : return;
522 : }
523 98115 : if (!fPrimaryIonisation) TVirtualMC::GetMC()->SetMaxStep(fStepSize);
524 :
525 256810 : }
|