]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MFT/AliMFT.cxx
Analysis files updated
[u/mrichter/AliRoot.git] / MFT / AliMFT.cxx
CommitLineData
820b4d9e 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//
18// Main class of the ALICE Muon Forward Tracker
19//
20// Contact author: antonio.uras@cern.ch
21//
22//====================================================================================================================================================
23
24#include "AliLog.h"
25#include "TFile.h"
26#include "TGeoManager.h"
27#include "TGeoVolume.h"
28#include "TGeoMatrix.h"
29#include "TVirtualMC.h"
30#include "TClonesArray.h"
31#include "TGeoGlobalMagField.h"
32#include "AliRun.h"
33#include "AliLoader.h"
34#include "AliDetector.h"
35#include "AliMC.h"
36#include "AliMagF.h"
37#include "AliMFT.h"
38#include "AliMFTHit.h"
39#include "AliMFTDigit.h"
40#include "AliMFTCluster.h"
41#include "AliTrackReference.h"
42#include "AliMFTSegmentation.h"
43#include "AliMFTDigitizer.h"
44#include "AliMFTPlane.h"
45#include "TString.h"
46#include "TObjArray.h"
47
48ClassImp(AliMFT)
49
50//====================================================================================================================================================
51
52AliMFT::AliMFT():
53 AliDetector(),
54 fVersion(1),
55 fNPlanes(0),
56 fNSlices(0),
57 fSDigitsPerPlane(0),
58 fDigitsPerPlane(0),
59 fRecPointsPerPlane(0),
60 fSideDigits(0),
61 fSegmentation(0),
62 fNameGeomFile(0),
63 fChargeDispersion(0),
64 fSingleStepForChargeDispersion(0),
65 fNStepForChargeDispersion(0),
d4643a10 66 fDensitySupportOverSi(0.1)
820b4d9e 67{
68
69 // default constructor
70
71}
72
73//====================================================================================================================================================
74
75AliMFT::AliMFT(const Char_t *name, const Char_t *title):
76 AliDetector(name, title),
77 fVersion(1),
78 fNPlanes(0),
79 fNSlices(0),
80 fSDigitsPerPlane(0),
81 fDigitsPerPlane(0),
82 fRecPointsPerPlane(0),
83 fSideDigits(0),
84 fSegmentation(0),
85 fNameGeomFile(0),
86 fChargeDispersion(0),
87 fSingleStepForChargeDispersion(0),
88 fNStepForChargeDispersion(0),
d4643a10 89 fDensitySupportOverSi(0.1)
820b4d9e 90{
91
92 fNameGeomFile = "AliMFTGeometry.root";
93
94 SetGeometry();
95
96 Init();
97
98}
99
100//====================================================================================================================================================
101
102AliMFT::AliMFT(const Char_t *name, const Char_t *title, Char_t *nameGeomFile):
103 AliDetector(name, title),
104 fVersion(1),
105 fNPlanes(0),
106 fNSlices(0),
107 fSDigitsPerPlane(0),
108 fDigitsPerPlane(0),
109 fRecPointsPerPlane(0),
110 fSideDigits(0),
111 fSegmentation(0),
112 fNameGeomFile(0),
113 fChargeDispersion(0),
114 fSingleStepForChargeDispersion(0),
115 fNStepForChargeDispersion(0),
d4643a10 116 fDensitySupportOverSi(0.1)
820b4d9e 117{
118
119 fNameGeomFile = nameGeomFile;
120
121 SetGeometry();
122
123 Init();
124
125}
126
127//====================================================================================================================================================
128
129AliMFT::~AliMFT() {
130
131 if (fSDigitsPerPlane) { fSDigitsPerPlane->Delete(); delete fSDigitsPerPlane; }
132 if (fDigitsPerPlane) { fDigitsPerPlane->Delete(); delete fDigitsPerPlane; }
133 if (fRecPointsPerPlane) { fRecPointsPerPlane->Delete(); delete fRecPointsPerPlane; }
134
135}
136
137//====================================================================================================================================================
138
139void AliMFT::CreateMaterials() {
140
141 // Definition of MFT materials - to be updated to the most recent values
142
143 AliInfo("Start MFT materials");
144
145 // data from PDG booklet 2002 density [gr/cm^3] rad len [cm] abs len [cm]
146 Float_t aAir[4]={12,14,16,36} , zAir[4]={6,7,8,18} , wAir[4]={0.000124,0.755267,0.231781,0.012827} , dAir=0.00120479; Int_t nAir=4; // Air mixture
147 Float_t aSi = 28.085 , zSi = 14 , dSi = 2.329 , radSi = 21.82/dSi , absSi = 108.4/dSi ; // Silicon
148
149 Int_t matId = 0; // tmp material id number
150 Int_t unsens = 0, sens=1; // sensitive or unsensitive medium
151 Int_t itgfld = 3; // type of field intergration 0 no field -1 user in guswim 1 Runge Kutta 2 helix 3 const field along z
152 Float_t maxfld = 5.; // max field value
153
154 Float_t tmaxfd = -10.0; // max deflection angle due to magnetic field in one step
155 Float_t stemax = 0.001; // max step allowed [cm]
156 Float_t deemax = -0.2; // maximum fractional energy loss in one step 0<deemax<=1
157 Float_t epsil = 0.001; // tracking precision [cm]
158 Float_t stmin = -0.001; // minimum step due to continuous processes [cm] (negative value: choose it automatically)
159
160 Float_t tmaxfdSi = 0.1; // max deflection angle due to magnetic field in one step
161 Float_t stemaxSi = 5.0e-4; // maximum step allowed [cm]
162 Float_t deemaxSi = 0.1; // maximum fractional energy loss in one step 0<deemax<=1
163 Float_t epsilSi = 0.5e-4; // tracking precision [cm]
164 Float_t stminSi = -0.001; // minimum step due to continuous processes [cm] (negative value: choose it automatically)
165
166 Int_t isxfld = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Integ(); // from CreateMaterials in STRUCT/AliPIPEv3.cxx
167 Float_t sxmgmx = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Max(); // from CreateMaterials in STRUCT/AliPIPEv3.cxx
168
169 AliMixture(++matId,"Air", aAir, zAir, dAir, nAir, wAir);
170 AliMedium(kAir, "Air", matId, unsens, itgfld, maxfld, tmaxfd, stemax, deemax, epsil, stmin);
171
172 AliMaterial(++matId, "Si", aSi, zSi, dSi, radSi, absSi );
173 AliMedium(kSi, "Si", matId, sens, isxfld, sxmgmx, tmaxfdSi, stemaxSi, deemaxSi, epsilSi, stminSi);
174
175 AliMaterial(++matId, "Readout", aSi, zSi, dSi, radSi, absSi );
176 AliMedium(kReadout, "Readout", matId, unsens, isxfld, sxmgmx, tmaxfdSi, stemaxSi, deemaxSi, epsilSi, stminSi);
177
d4643a10 178 AliMaterial(++matId, "Support", aSi, zSi, dSi*fDensitySupportOverSi, radSi/fDensitySupportOverSi, absSi/fDensitySupportOverSi);
820b4d9e 179 AliMedium(kSupport, "Support", matId, unsens, isxfld, sxmgmx, tmaxfdSi, stemaxSi, deemaxSi, epsilSi, stminSi);
180
181 AliInfo("End MFT materials");
182
183}
184
185//====================================================================================================================================================
186
187void AliMFT::CreateGeometry() {
188
189 // Creates detailed geometry simulation (currently GEANT volumes tree)
190
191 AliInfo("Start MFT preliminary version building");
192 if(!gMC->IsRootGeometrySupported()) return;
193 TGeoVolumeAssembly *vol = CreateVol();
194 AliInfo("TGeoVolumeAssembly created!");
195 gGeoManager->GetVolume("ALIC")->AddNode(vol,0);
196 AliInfo("Stop MFT preliminary version building");
197
198 if (fNStepForChargeDispersion) fSingleStepForChargeDispersion = fChargeDispersion/Double_t(fNStepForChargeDispersion);
199
200}
201
202//====================================================================================================================================================
203
204void AliMFT::StepManager() {
205
206 // Full Step Manager
207
d4643a10 208 AliDebug(2, Form("Entering StepManager: gMC->CurrentVolName() = %s", gMC->CurrentVolName()));
209
820b4d9e 210 if (!fSegmentation) AliFatal("No segmentation available"); // DO WE HAVE A SEGMENTATION???
211
212 if (!(this->IsActive())) return;
213 if (!(gMC->TrackCharge())) return;
214
215 TString planeNumber = gMC->CurrentVolName();
216 TString detElemNumber = gMC->CurrentVolName();
217 if (planeNumber.Contains("support")) return;
218 if (planeNumber.Contains("readout")) return;
219 planeNumber.Remove(0,9);
220 detElemNumber.Remove(0,19);
221 planeNumber.Remove(2);
222 detElemNumber.Remove(3);
223 Int_t detElemID = fSegmentation->GetDetElemID(planeNumber.Atoi(), detElemNumber.Atoi());
224
225 if (gMC->IsTrackExiting()) {
226 AddTrackReference(gAlice->GetMCApp()->GetCurrentTrackNumber(), AliTrackReference::kMFT);
227 }
228
229 static TLorentzVector position, momentum;
230 static AliMFTHit hit;
231
232 Int_t status = 0;
233
234 // Track status
235 if (gMC->IsTrackInside()) status += 1;
236 if (gMC->IsTrackEntering()) status += 2;
237 if (gMC->IsTrackExiting()) status += 4;
238 if (gMC->IsTrackOut()) status += 8;
239 if (gMC->IsTrackDisappeared()) status += 16;
240 if (gMC->IsTrackStop()) status += 32;
241 if (gMC->IsTrackAlive()) status += 64;
242
243 // ---------- Fill hit structure
244
245 hit.SetDetElemID(detElemID);
246 hit.SetPlane(planeNumber.Atoi());
247 hit.SetTrack(gAlice->GetMCApp()->GetCurrentTrackNumber());
248
249 gMC->TrackPosition(position);
250 gMC->TrackMomentum(momentum);
251
d4643a10 252 AliDebug(1, Form("AliMFT::StepManager()->%s Hit #%06d (x=%f, y=%f, z=%f) belongs to track %02d\n",
253 gMC->CurrentVolName(), fNhits, position.X(), position.Y(), position.Z(), gAlice->GetMCApp()->GetCurrentTrackNumber()));
820b4d9e 254
255 hit.SetPosition(position);
256 hit.SetTOF(gMC->TrackTime());
257 hit.SetMomentum(momentum);
258 hit.SetStatus(status);
259 hit.SetEloss(gMC->Edep());
260 // hit.SetShunt(GetIshunt());
261// if (gMC->IsTrackEntering()) {
262// hit.SetStartPosition(position);
263// hit.SetStartTime(gMC->TrackTime());
264// hit.SetStartStatus(status);
265// return; // don't save entering hit.
266// }
267
268 // Fill hit structure with this new hit.
269 new ((*fHits)[fNhits++]) AliMFTHit(hit);
270
271 // Save old position... for next hit.
272// hit.SetStartPosition(position);
273// hit.SetStartTime(gMC->TrackTime());
274// hit.SetStartStatus(status);
275
276 return;
277
278}
279
280//====================================================================================================================================================
281
282TGeoVolumeAssembly* AliMFT::CreateVol() {
283
284 // method to create the MFT Geometry (silicon circular planes)
285
286 if (!fSegmentation) CreateGeometry();
287
288 TGeoVolumeAssembly *vol = new TGeoVolumeAssembly("MFT");
289 TGeoMedium *silicon = gGeoManager->GetMedium("MFT_Si");
290 TGeoMedium *readout = gGeoManager->GetMedium("MFT_Readout");
291 TGeoMedium *support = gGeoManager->GetMedium("MFT_Support");
292 for (Int_t iPar=0; iPar<8; iPar++) AliInfo(Form("silicon->GetParam(%d) = %f", iPar, silicon->GetParam(iPar)));
293 for (Int_t iPar=0; iPar<8; iPar++) AliInfo(Form("readout->GetParam(%d) = %f", iPar, readout->GetParam(iPar)));
294 for (Int_t iPar=0; iPar<8; iPar++) AliInfo(Form("support->GetParam(%d) = %f", iPar, support->GetParam(iPar)));
295
296 Double_t origin[3] = {0};
297
298 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
299
300 AliDebug(1, Form("Creating volumes for MFT plane %02d",iPlane));
301
302 AliMFTPlane *plane = fSegmentation->GetPlane(iPlane);
303
304 // --------- support element(s)
305
306 origin[0] = 0.5*(plane->GetSupportElement(0)->GetAxis(0)->GetXmax() + plane->GetSupportElement(0)->GetAxis(0)->GetXmin());
307 origin[1] = 0.5*(plane->GetSupportElement(0)->GetAxis(1)->GetXmax() + plane->GetSupportElement(0)->GetAxis(1)->GetXmin());
308 origin[2] = -0.5*(plane->GetSupportElement(0)->GetAxis(2)->GetXmax() + plane->GetSupportElement(0)->GetAxis(2)->GetXmin());
309 TGeoVolume *supportElem = gGeoManager->MakeTube(Form("MFT_plane%02d_support", iPlane), support,
310 plane->GetRMinSupport(),
311 plane->GetRMaxSupport(),
312 0.5*(plane->GetSupportElement(0)->GetAxis(2)->GetXmax() -
313 plane->GetSupportElement(0)->GetAxis(2)->GetXmin()) );
d4643a10 314 AliDebug(2, Form("Created vol %s", supportElem->GetName()));
315 supportElem->SetLineColor(kCyan-9);
820b4d9e 316 vol -> AddNode(supportElem, 0, new TGeoTranslation(origin[0], origin[1], origin[2]));
317
318 AliDebug(1, "support elements created!");
319
320 // --------- active elements
321
322 for (Int_t iActive=0; iActive<plane->GetNActiveElements(); iActive++) {
323
324 Double_t dx = 0.5*TMath::Abs(plane->GetActiveElement(iActive)->GetAxis(0)->GetXmax() - plane->GetActiveElement(iActive)->GetAxis(0)->GetXmin());
325 Double_t dy = 0.5*TMath::Abs(plane->GetActiveElement(iActive)->GetAxis(1)->GetXmax() - plane->GetActiveElement(iActive)->GetAxis(1)->GetXmin());
326 Double_t dz = 0.5*TMath::Abs(plane->GetActiveElement(iActive)->GetAxis(2)->GetXmax() - plane->GetActiveElement(iActive)->GetAxis(2)->GetXmin());
327 dz /= Double_t(fNSlices);
328
329 origin[0] = 0.5*(plane->GetActiveElement(iActive)->GetAxis(0)->GetXmax() + plane->GetActiveElement(iActive)->GetAxis(0)->GetXmin());
330 origin[1] = 0.5*(plane->GetActiveElement(iActive)->GetAxis(1)->GetXmax() + plane->GetActiveElement(iActive)->GetAxis(1)->GetXmin());
331
332 for (Int_t iSlice=0; iSlice<fNSlices; iSlice++) {
333 origin[2] = -0.5*(plane->GetActiveElement(iActive)->GetAxis(2)->GetXmin() + 2*dz*(iSlice+1) + plane->GetActiveElement(iActive)->GetAxis(2)->GetXmin() + 2*dz*(iSlice) );
334 TGeoVolume *activeElem = gGeoManager->MakeBox(Form("MFT_plane%02d_active%03d_slice%02d", iPlane, iActive, iSlice), silicon, dx, dy, dz);
d4643a10 335 AliDebug(2, Form("Created vol %s", activeElem->GetName()));
336 activeElem->SetLineColor(kGreen);
820b4d9e 337 vol -> AddNode(activeElem, 0, new TGeoTranslation(origin[0], origin[1], origin[2]));
338 }
339
340 }
341
342 AliDebug(1, "active elements created!");
343
344 // --------- readout elements
345
346 for (Int_t iReadout=0; iReadout<plane->GetNReadoutElements(); iReadout++) {
347
348 Double_t dx = 0.5*TMath::Abs(plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmax() - plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmin());
349 Double_t dy = 0.5*TMath::Abs(plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmax() - plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmin());
350 Double_t dz = 0.5*TMath::Abs(plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmax() - plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmin());
351
352 origin[0] = 0.5*(plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmax() + plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmin());
353 origin[1] = 0.5*(plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmax() + plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmin());
354 origin[2] = -0.5*(plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmax() + plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmin());
355
356 TGeoVolume *readoutElem = gGeoManager->MakeBox(Form("MFT_plane%02d_readout%03d", iPlane, iReadout), readout, dx, dy, dz);
d4643a10 357 AliDebug(2, Form("Created vol %s", readoutElem->GetName()));
358 readoutElem->SetLineColor(kRed);
820b4d9e 359 vol -> AddNode(readoutElem, 0, new TGeoTranslation(origin[0], origin[1], origin[2]));
360
361 }
362
363 AliDebug(1, "readout elements created!");
364
365 }
366
367 return vol;
368
369}
370
371//====================================================================================================================================================
372
373void AliMFT::Hits2SDigits(){
374
375 // Interface method invoked from AliSimulation to create a list of sdigits corresponding to list of hits. Every hit generates one sdigit.
376
377 AliDebug(1,"Start Hits2SDigits.");
378
379 if (!fSegmentation) CreateGeometry();
380
381 if (!fLoader->TreeH()) fLoader->LoadHits();
382
383 if (!fLoader->TreeS()) {
384
385 for (Int_t iEvt=0;iEvt<fLoader->GetRunLoader()->GetNumberOfEvents(); iEvt++) {
386
387 fLoader->GetRunLoader()->GetEvent(iEvt);
388 fLoader->MakeTree("S");
389 MakeBranch("S");
390 SetTreeAddress();
391
392 AliDebug(1, Form("Event %03d: fLoader->TreeH()->GetEntries() = %2d", iEvt, Int_t(fLoader->TreeH()->GetEntries())));
393
394 for (Int_t iTrack=0; iTrack<fLoader->TreeH()->GetEntries(); iTrack++) {
395 fLoader->TreeH()->GetEntry(iTrack);
396 Hits2SDigitsLocal(Hits(), GetSDigitsList(), iTrack); // convert these hits to a list of sdigits
397 }
398
399 fLoader->TreeS()->Fill();
400 fLoader->WriteSDigits("OVERWRITE");
401 ResetSDigits();
402
403 }
404 }
405
406 fLoader->UnloadHits();
407 fLoader->UnloadSDigits();
408
409 AliDebug(1,"Stop Hits2SDigits.");
410
411}
412
413//====================================================================================================================================================
414
415void AliMFT::Hits2SDigitsLocal(TClonesArray *hits, const TObjArray *pSDig, Int_t track) {
416
417 // Add sdigits of these hits to the list
418
419 AliDebug(1, "Start Hits2SDigitsLocal");
420
421 if (!fSegmentation) CreateGeometry();
422
423 TClonesArray *pSDigList[fNMaxPlanes];
424 for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) pSDigList[iPlane] = NULL;
425 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
426 pSDigList[iPlane] = (TClonesArray*) (*pSDig)[iPlane];
427 AliDebug(1,Form("Entries of pSDigList %3d; plane: %02d,",pSDigList[iPlane]->GetEntries(),iPlane));
428 if (!track && pSDigList[iPlane]->GetEntries()!=0) AliErrorClass("Some of sdigits lists is not empty");
429 }
430
431 for (Int_t iHit=0; iHit<hits->GetEntries(); iHit++) {
432
433 AliMFTHit *hit = (AliMFTHit*) hits->At(iHit);
434
435 AliMFTDigit sDigit;
436 sDigit.SetEloss(hit->GetEloss());
437 sDigit.SetDetElemID(hit->GetDetElemID());
438 sDigit.SetPlane(hit->GetPlane());
439 sDigit.AddMCLabel(hit->GetTrack());
440
441 Int_t xPixel = -1;
442 Int_t yPixel = -1;
443 if (fSegmentation->Hit2PixelID(hit->X(), hit->Y(), sDigit.GetDetElemID(), xPixel, yPixel)) {
444 sDigit.SetPixID(xPixel, yPixel, 0);
445 sDigit.SetPixWidth(fSegmentation->GetPixelSizeX(sDigit.GetDetElemID()),
446 fSegmentation->GetPixelSizeY(sDigit.GetDetElemID()),
447 fSegmentation->GetPixelSizeZ(sDigit.GetDetElemID()));
448 sDigit.SetPixCenter(fSegmentation->GetPixelCenterX(sDigit.GetDetElemID(), xPixel),
449 fSegmentation->GetPixelCenterY(sDigit.GetDetElemID(), yPixel),
450 fSegmentation->GetPixelCenterZ(sDigit.GetDetElemID(), 0));
451 new ((*pSDigList[sDigit.GetPlane()])[pSDigList[sDigit.GetPlane()]->GetEntries()]) AliMFTDigit(sDigit);
452 AliDebug(1, Form("Created new sdigit (%f, %f, %f) from hit (%f, %f, %f)",
453 sDigit.GetPixelCenterX(), sDigit.GetPixelCenterY(), sDigit.GetPixelCenterZ(), hit->X(), hit->Y(), hit->Z()));
454// AliDebug(1, Form("Created new sdigit from hit: residual is (%f, %f, %f)",
455// sDigit.GetPixelCenterX()-hit->X(), sDigit.GetPixelCenterY()-hit->Y(), sDigit.GetPixelCenterZ()-hit->Z()));
456 }
457
458 // creating "side hits" to simulate the effect of charge dispersion
459
460 Int_t xPixelNew = -1;
461 Int_t yPixelNew = -1;
462 Double_t x0 = hit->X();
463 Double_t y0 = hit->Y();
464 Double_t pi4 = TMath::Pi()/4.;
465 for (Int_t iStep=0; iStep<fNStepForChargeDispersion; iStep++) {
466 Double_t shift = (iStep+1) * fSingleStepForChargeDispersion;
467 for (Int_t iAngle=0; iAngle<8; iAngle++) {
468 Double_t shiftX = shift*TMath::Cos(iAngle*pi4);
469 Double_t shiftY = shift*TMath::Sin(iAngle*pi4);
470 if (fSegmentation->Hit2PixelID(x0+shiftX, y0+shiftY, hit->GetDetElemID(), xPixelNew, yPixelNew)) {
471 Bool_t digitExists = kFALSE;
472 if (xPixelNew==xPixel && yPixelNew==yPixel) digitExists = kTRUE;
473 if (!digitExists) {
474 for (Int_t iSideDigit=0; iSideDigit<fSideDigits->GetEntries(); iSideDigit++) {
475 if (xPixelNew==((AliMFTDigit*) fSideDigits->At(iSideDigit))->GetPixelX() &&
476 yPixelNew==((AliMFTDigit*) fSideDigits->At(iSideDigit))->GetPixelY()) {
477 digitExists = kTRUE;
478 break;
479 }
480 }
481 }
482 if (!digitExists) {
483 AliMFTDigit newSDigit;
484 newSDigit.SetEloss(0.);
485 newSDigit.SetDetElemID(hit->GetDetElemID());
486 newSDigit.SetPlane(hit->GetPlane());
487 newSDigit.AddMCLabel(hit->GetTrack());
488 newSDigit.SetPixID(xPixelNew, yPixelNew, 0);
489 newSDigit.SetPixWidth(fSegmentation->GetPixelSizeX(sDigit.GetDetElemID()),
490 fSegmentation->GetPixelSizeY(sDigit.GetDetElemID()),
491 fSegmentation->GetPixelSizeZ(sDigit.GetDetElemID()));
492 newSDigit.SetPixCenter(fSegmentation->GetPixelCenterX(sDigit.GetDetElemID(), xPixelNew),
493 fSegmentation->GetPixelCenterY(sDigit.GetDetElemID(), yPixelNew),
494 fSegmentation->GetPixelCenterZ(sDigit.GetDetElemID(), 0));
495 new ((*fSideDigits)[fSideDigits->GetEntries()]) AliMFTDigit(newSDigit);
496 }
497 }
498 }
499 }
500
501 for (Int_t iSideDigit=0; iSideDigit<fSideDigits->GetEntries(); iSideDigit++) {
bcaf50eb 502 AliMFTDigit *newSDigit = (AliMFTDigit*) fSideDigits->At(iSideDigit);
503 new ((*pSDigList[sDigit.GetPlane()])[pSDigList[sDigit.GetPlane()]->GetEntries()]) AliMFTDigit(*newSDigit);
820b4d9e 504 }
505
bcaf50eb 506 fSideDigits->Delete();
820b4d9e 507
508 }
bcaf50eb 509
510 // ------------ In case we should simulate a rectangular pattern of pixel...
511
512 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
513 if (fSegmentation->GetPlane(iPlane)->HasPixelRectangularPatternAlongY()) {
514 Int_t nSDigits = pSDigList[iPlane]->GetEntries();
515 for (Int_t iSDigit=0; iSDigit<nSDigits; iSDigit++) {
516 AliMFTDigit *mySDig = (AliMFTDigit*) (pSDigList[iPlane]->At(iSDigit));
517 if (mySDig->GetPixelX()%2 == mySDig->GetPixelY()%2) { // both pair or both odd
518 Int_t xPixelNew = mySDig->GetPixelX();
519 Int_t yPixelNew = mySDig->GetPixelY()+1;
520 if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixelNew, yPixelNew)) {
521 AliMFTDigit newSDigit;
522 newSDigit.SetEloss(0.);
523 newSDigit.SetDetElemID(mySDig->GetDetElemID());
524 newSDigit.SetPlane(iPlane);
525 newSDigit.SetPixID(xPixelNew, yPixelNew, 0);
526 newSDigit.SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit.GetDetElemID()),
527 fSegmentation->GetPixelSizeY(newSDigit.GetDetElemID()),
528 fSegmentation->GetPixelSizeZ(newSDigit.GetDetElemID()));
529 newSDigit.SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit.GetDetElemID(), xPixelNew),
530 fSegmentation->GetPixelCenterY(newSDigit.GetDetElemID(), yPixelNew),
531 fSegmentation->GetPixelCenterZ(newSDigit.GetDetElemID(), 0));
532 new ((*pSDigList[iPlane])[pSDigList[iPlane]->GetEntries()]) AliMFTDigit(newSDigit);
533 }
534 }
535 else { // pair-odd
536 Int_t xPixelNew = mySDig->GetPixelX();
537 Int_t yPixelNew = mySDig->GetPixelY()-1;
538 if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixelNew, yPixelNew)) {
539 AliMFTDigit newSDigit;
540 newSDigit.SetEloss(0.);
541 newSDigit.SetDetElemID(mySDig->GetDetElemID());
542 newSDigit.SetPlane(iPlane);
543 newSDigit.SetPixID(xPixelNew, yPixelNew, 0);
544 newSDigit.SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit.GetDetElemID()),
545 fSegmentation->GetPixelSizeY(newSDigit.GetDetElemID()),
546 fSegmentation->GetPixelSizeZ(newSDigit.GetDetElemID()));
547 newSDigit.SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit.GetDetElemID(), xPixelNew),
548 fSegmentation->GetPixelCenterY(newSDigit.GetDetElemID(), yPixelNew),
549 fSegmentation->GetPixelCenterZ(newSDigit.GetDetElemID(), 0));
550 new ((*pSDigList[iPlane])[pSDigList[iPlane]->GetEntries()]) AliMFTDigit(newSDigit);
551 }
552 }
553 }
554 }
555 }
556
557 //------------------------------------------------------------------------
820b4d9e 558
559 AliDebug(1,"Stop Hits2SDigitsLocal");
560
561}
562
563//====================================================================================================================================================
564
565void AliMFT::MakeBranch(Option_t *option) {
566
567 printf("AliMFT::MakeBranch(...)\n");
568
569 // Create Tree branches
570 AliDebug(1, Form("Start with option= %s.",option));
571
572 const Int_t kBufSize = 4000;
573
574 const Char_t *cH = strstr(option,"H");
575 const Char_t *cD = strstr(option,"D");
576 const Char_t *cS = strstr(option,"S");
577
578 if (cH && fLoader->TreeH()) {
579 CreateHits();
580 MakeBranchInTree(fLoader->TreeH(), "MFT", &fHits, kBufSize, 0);
581 }
582
583 if (cS && fLoader->TreeS()) {
584 CreateSDigits();
585 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) MakeBranchInTree(fLoader->TreeS(),
586 Form("Plane_%02d",iPlane),
587 &((*fSDigitsPerPlane)[iPlane]),
588 kBufSize, 0);
589 }
590
591 if (cD && fLoader->TreeD()) {
592 CreateDigits();
593 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) MakeBranchInTree(fLoader->TreeD(),
594 Form("Plane_%02d",iPlane),
595 &((*fDigitsPerPlane)[iPlane]),
596 kBufSize, 0);
597 }
598
599 AliDebug(1,"Stop.");
600
601}
602
603//====================================================================================================================================================
604
605void AliMFT::SetTreeAddress() {
606
607 AliDebug(1, "AliMFT::SetTreeAddress()");
608
609 //Set branch address for the Hits and Digits Tree.
610 AliDebug(1, "Start.");
611
612 AliDebug(1, Form("AliMFT::SetTreeAddress Hits fLoader->TreeH() = %p\n", fLoader->TreeH()));
613 if (fLoader->TreeH() && fLoader->TreeH()->GetBranch("MFT")) {
614 CreateHits();
615 fLoader->TreeH()->SetBranchAddress("MFT", &fHits);
616 }
617
618 AliDebug(1, Form("AliMFT::SetTreeAddress SDigits fLoader->TreeS() = %p\n", fLoader->TreeS()));
619 if (fLoader->TreeS() && fLoader->TreeS()->GetBranch("Plane_00")) {
620 CreateSDigits();
621 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
622 fLoader->TreeS()->SetBranchAddress(Form("Plane_%02d",iPlane), &((*fSDigitsPerPlane)[iPlane]));
623 }
624 }
625
626 AliDebug(1, Form("AliMFT::SetTreeAddress Digits fLoader->TreeD() = %p\n", fLoader->TreeD()));
627 if (fLoader->TreeD() && fLoader->TreeD()->GetBranch("Plane_00")) {
628 CreateDigits();
629 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
630 fLoader->TreeD()->SetBranchAddress(Form("Plane_%02d",iPlane), &((*fDigitsPerPlane)[iPlane]));
631 }
632 }
633
634 AliDebug(1, Form("AliMFT::SetTreeAddress RecPoints fLoader->TreeR() = %p\n", fLoader->TreeR()));
635 if (fLoader->TreeR() && fLoader->TreeR()->GetBranch("Plane_00")) {
636 CreateRecPoints();
637 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
638 fLoader->TreeR()->SetBranchAddress(Form("Plane_%02d",iPlane), &((*fRecPointsPerPlane)[iPlane]));
639 }
640 }
641
642 AliDebug(1,"Stop.");
643
644}
645
646//====================================================================================================================================================
647
648void AliMFT::SetGeometry() {
649
650 printf("AliMFT::SetGeometry\n");
651
652 fSegmentation = new AliMFTSegmentation(fNameGeomFile.Data());
653
654 fNPlanes = fSegmentation->GetNPlanes();
655
656}
657
658//====================================================================================================================================================
659
660void AliMFT::CreateHits() {
661
662 // create array of hits
663
664 AliDebug(1, "AliMFT::CreateHits()");
665
666 if (fHits) return;
667 fHits = new TClonesArray("AliMFTHit");
668
669}
670
671//====================================================================================================================================================
672
673void AliMFT::CreateSDigits() {
674
675 // create sdigits list
676
677 AliDebug(1, "AliMFT::CreateSDigits()");
678
679 if (fSDigitsPerPlane) return;
680 fSDigitsPerPlane = new TObjArray(fNPlanes);
681 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) fSDigitsPerPlane->AddAt(new TClonesArray("AliMFTDigit"), iPlane);
682
683 fSideDigits = new TClonesArray("AliMFTDigit");
684
685}
686
687//====================================================================================================================================================
688
689void AliMFT::CreateDigits() {
690
691 // create digits list
692
693 AliDebug(1, "AliMFT::CreateDigits()");
694
695 if (fDigitsPerPlane) return;
696 fDigitsPerPlane = new TObjArray(fNPlanes);
697 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) fDigitsPerPlane->AddAt(new TClonesArray("AliMFTDigit"), iPlane);
698
699}
700
701//====================================================================================================================================================
702
703void AliMFT::CreateRecPoints() {
704
705 // create recPoints list
706
707 AliDebug(1, "AliMFT::CreateRecPoints()");
708
709 if (fRecPointsPerPlane) return;
710 fRecPointsPerPlane = new TObjArray(fNPlanes);
711 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) fRecPointsPerPlane->AddAt(new TClonesArray("AliMFTCluster"), iPlane);
712
713}
714
715//====================================================================================================================================================