Macros for creating OCDB
[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),
58b93f36 56 fNSlices(1),
820b4d9e 57 fSDigitsPerPlane(0),
58 fDigitsPerPlane(0),
59 fRecPointsPerPlane(0),
60 fSideDigits(0),
61 fSegmentation(0),
62 fNameGeomFile(0),
58b93f36 63 fChargeDispersion(25.e-4),
820b4d9e 64 fSingleStepForChargeDispersion(0),
58b93f36 65 fNStepForChargeDispersion(4),
66 fDensitySupportOverSi(0.15)
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),
58b93f36 79 fNSlices(1),
820b4d9e 80 fSDigitsPerPlane(0),
81 fDigitsPerPlane(0),
82 fRecPointsPerPlane(0),
83 fSideDigits(0),
84 fSegmentation(0),
85 fNameGeomFile(0),
58b93f36 86 fChargeDispersion(25.e-4),
820b4d9e 87 fSingleStepForChargeDispersion(0),
58b93f36 88 fNStepForChargeDispersion(4),
89 fDensitySupportOverSi(0.15)
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),
58b93f36 106 fNSlices(1),
820b4d9e 107 fSDigitsPerPlane(0),
108 fDigitsPerPlane(0),
109 fRecPointsPerPlane(0),
110 fSideDigits(0),
111 fSegmentation(0),
112 fNameGeomFile(0),
58b93f36 113 fChargeDispersion(25.e-4),
820b4d9e 114 fSingleStepForChargeDispersion(0),
58b93f36 115 fNStepForChargeDispersion(4),
116 fDensitySupportOverSi(0.15)
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
a2b7dc2a 204void AliMFT::AddAlignableVolumes() {
205
206 // Create entries for alignable volumes associating the symbolic volume
207 // name with the corresponding volume path. Needs to be syncronized with
208 // eventual changes in the geometry.
209
210 TString sysName = "MFT";
211 TString volPath = "/ALIC_1/MFT_0";
212
213 if (!gGeoManager->SetAlignableEntry(sysName.Data(),volPath.Data())) {
214 AliFatal(Form("Alignable entry %s not created. Volume path %s not valid", sysName.Data(), volPath.Data()));
215 }
216
217}
218
219//====================================================================================================================================================
220
820b4d9e 221void AliMFT::StepManager() {
222
223 // Full Step Manager
224
d4643a10 225 AliDebug(2, Form("Entering StepManager: gMC->CurrentVolName() = %s", gMC->CurrentVolName()));
226
820b4d9e 227 if (!fSegmentation) AliFatal("No segmentation available"); // DO WE HAVE A SEGMENTATION???
228
229 if (!(this->IsActive())) return;
230 if (!(gMC->TrackCharge())) return;
231
232 TString planeNumber = gMC->CurrentVolName();
233 TString detElemNumber = gMC->CurrentVolName();
234 if (planeNumber.Contains("support")) return;
235 if (planeNumber.Contains("readout")) return;
236 planeNumber.Remove(0,9);
237 detElemNumber.Remove(0,19);
238 planeNumber.Remove(2);
239 detElemNumber.Remove(3);
240 Int_t detElemID = fSegmentation->GetDetElemID(planeNumber.Atoi(), detElemNumber.Atoi());
241
242 if (gMC->IsTrackExiting()) {
243 AddTrackReference(gAlice->GetMCApp()->GetCurrentTrackNumber(), AliTrackReference::kMFT);
244 }
245
246 static TLorentzVector position, momentum;
247 static AliMFTHit hit;
248
249 Int_t status = 0;
250
251 // Track status
252 if (gMC->IsTrackInside()) status += 1;
253 if (gMC->IsTrackEntering()) status += 2;
254 if (gMC->IsTrackExiting()) status += 4;
255 if (gMC->IsTrackOut()) status += 8;
256 if (gMC->IsTrackDisappeared()) status += 16;
257 if (gMC->IsTrackStop()) status += 32;
258 if (gMC->IsTrackAlive()) status += 64;
259
260 // ---------- Fill hit structure
261
262 hit.SetDetElemID(detElemID);
263 hit.SetPlane(planeNumber.Atoi());
264 hit.SetTrack(gAlice->GetMCApp()->GetCurrentTrackNumber());
265
266 gMC->TrackPosition(position);
267 gMC->TrackMomentum(momentum);
268
d4643a10 269 AliDebug(1, Form("AliMFT::StepManager()->%s Hit #%06d (x=%f, y=%f, z=%f) belongs to track %02d\n",
270 gMC->CurrentVolName(), fNhits, position.X(), position.Y(), position.Z(), gAlice->GetMCApp()->GetCurrentTrackNumber()));
820b4d9e 271
272 hit.SetPosition(position);
273 hit.SetTOF(gMC->TrackTime());
274 hit.SetMomentum(momentum);
275 hit.SetStatus(status);
276 hit.SetEloss(gMC->Edep());
277 // hit.SetShunt(GetIshunt());
278// if (gMC->IsTrackEntering()) {
279// hit.SetStartPosition(position);
280// hit.SetStartTime(gMC->TrackTime());
281// hit.SetStartStatus(status);
282// return; // don't save entering hit.
283// }
284
285 // Fill hit structure with this new hit.
286 new ((*fHits)[fNhits++]) AliMFTHit(hit);
287
288 // Save old position... for next hit.
289// hit.SetStartPosition(position);
290// hit.SetStartTime(gMC->TrackTime());
291// hit.SetStartStatus(status);
292
293 return;
294
295}
296
297//====================================================================================================================================================
298
299TGeoVolumeAssembly* AliMFT::CreateVol() {
300
301 // method to create the MFT Geometry (silicon circular planes)
302
303 if (!fSegmentation) CreateGeometry();
304
305 TGeoVolumeAssembly *vol = new TGeoVolumeAssembly("MFT");
306 TGeoMedium *silicon = gGeoManager->GetMedium("MFT_Si");
307 TGeoMedium *readout = gGeoManager->GetMedium("MFT_Readout");
308 TGeoMedium *support = gGeoManager->GetMedium("MFT_Support");
309 for (Int_t iPar=0; iPar<8; iPar++) AliInfo(Form("silicon->GetParam(%d) = %f", iPar, silicon->GetParam(iPar)));
310 for (Int_t iPar=0; iPar<8; iPar++) AliInfo(Form("readout->GetParam(%d) = %f", iPar, readout->GetParam(iPar)));
311 for (Int_t iPar=0; iPar<8; iPar++) AliInfo(Form("support->GetParam(%d) = %f", iPar, support->GetParam(iPar)));
312
313 Double_t origin[3] = {0};
314
315 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
316
317 AliDebug(1, Form("Creating volumes for MFT plane %02d",iPlane));
318
319 AliMFTPlane *plane = fSegmentation->GetPlane(iPlane);
320
321 // --------- support element(s)
322
323 origin[0] = 0.5*(plane->GetSupportElement(0)->GetAxis(0)->GetXmax() + plane->GetSupportElement(0)->GetAxis(0)->GetXmin());
324 origin[1] = 0.5*(plane->GetSupportElement(0)->GetAxis(1)->GetXmax() + plane->GetSupportElement(0)->GetAxis(1)->GetXmin());
325 origin[2] = -0.5*(plane->GetSupportElement(0)->GetAxis(2)->GetXmax() + plane->GetSupportElement(0)->GetAxis(2)->GetXmin());
326 TGeoVolume *supportElem = gGeoManager->MakeTube(Form("MFT_plane%02d_support", iPlane), support,
327 plane->GetRMinSupport(),
328 plane->GetRMaxSupport(),
329 0.5*(plane->GetSupportElement(0)->GetAxis(2)->GetXmax() -
330 plane->GetSupportElement(0)->GetAxis(2)->GetXmin()) );
d4643a10 331 AliDebug(2, Form("Created vol %s", supportElem->GetName()));
332 supportElem->SetLineColor(kCyan-9);
820b4d9e 333 vol -> AddNode(supportElem, 0, new TGeoTranslation(origin[0], origin[1], origin[2]));
334
335 AliDebug(1, "support elements created!");
336
337 // --------- active elements
338
339 for (Int_t iActive=0; iActive<plane->GetNActiveElements(); iActive++) {
340
341 Double_t dx = 0.5*TMath::Abs(plane->GetActiveElement(iActive)->GetAxis(0)->GetXmax() - plane->GetActiveElement(iActive)->GetAxis(0)->GetXmin());
342 Double_t dy = 0.5*TMath::Abs(plane->GetActiveElement(iActive)->GetAxis(1)->GetXmax() - plane->GetActiveElement(iActive)->GetAxis(1)->GetXmin());
343 Double_t dz = 0.5*TMath::Abs(plane->GetActiveElement(iActive)->GetAxis(2)->GetXmax() - plane->GetActiveElement(iActive)->GetAxis(2)->GetXmin());
344 dz /= Double_t(fNSlices);
345
346 origin[0] = 0.5*(plane->GetActiveElement(iActive)->GetAxis(0)->GetXmax() + plane->GetActiveElement(iActive)->GetAxis(0)->GetXmin());
347 origin[1] = 0.5*(plane->GetActiveElement(iActive)->GetAxis(1)->GetXmax() + plane->GetActiveElement(iActive)->GetAxis(1)->GetXmin());
348
349 for (Int_t iSlice=0; iSlice<fNSlices; iSlice++) {
350 origin[2] = -0.5*(plane->GetActiveElement(iActive)->GetAxis(2)->GetXmin() + 2*dz*(iSlice+1) + plane->GetActiveElement(iActive)->GetAxis(2)->GetXmin() + 2*dz*(iSlice) );
351 TGeoVolume *activeElem = gGeoManager->MakeBox(Form("MFT_plane%02d_active%03d_slice%02d", iPlane, iActive, iSlice), silicon, dx, dy, dz);
d4643a10 352 AliDebug(2, Form("Created vol %s", activeElem->GetName()));
353 activeElem->SetLineColor(kGreen);
820b4d9e 354 vol -> AddNode(activeElem, 0, new TGeoTranslation(origin[0], origin[1], origin[2]));
355 }
356
357 }
358
359 AliDebug(1, "active elements created!");
360
361 // --------- readout elements
362
363 for (Int_t iReadout=0; iReadout<plane->GetNReadoutElements(); iReadout++) {
364
365 Double_t dx = 0.5*TMath::Abs(plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmax() - plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmin());
366 Double_t dy = 0.5*TMath::Abs(plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmax() - plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmin());
367 Double_t dz = 0.5*TMath::Abs(plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmax() - plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmin());
368
369 origin[0] = 0.5*(plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmax() + plane->GetReadoutElement(iReadout)->GetAxis(0)->GetXmin());
370 origin[1] = 0.5*(plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmax() + plane->GetReadoutElement(iReadout)->GetAxis(1)->GetXmin());
371 origin[2] = -0.5*(plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmax() + plane->GetReadoutElement(iReadout)->GetAxis(2)->GetXmin());
372
373 TGeoVolume *readoutElem = gGeoManager->MakeBox(Form("MFT_plane%02d_readout%03d", iPlane, iReadout), readout, dx, dy, dz);
d4643a10 374 AliDebug(2, Form("Created vol %s", readoutElem->GetName()));
375 readoutElem->SetLineColor(kRed);
820b4d9e 376 vol -> AddNode(readoutElem, 0, new TGeoTranslation(origin[0], origin[1], origin[2]));
377
378 }
379
380 AliDebug(1, "readout elements created!");
381
382 }
383
384 return vol;
385
386}
387
388//====================================================================================================================================================
389
390void AliMFT::Hits2SDigits(){
391
392 // Interface method invoked from AliSimulation to create a list of sdigits corresponding to list of hits. Every hit generates one sdigit.
393
394 AliDebug(1,"Start Hits2SDigits.");
395
396 if (!fSegmentation) CreateGeometry();
397
398 if (!fLoader->TreeH()) fLoader->LoadHits();
399
400 if (!fLoader->TreeS()) {
401
402 for (Int_t iEvt=0;iEvt<fLoader->GetRunLoader()->GetNumberOfEvents(); iEvt++) {
403
404 fLoader->GetRunLoader()->GetEvent(iEvt);
405 fLoader->MakeTree("S");
406 MakeBranch("S");
407 SetTreeAddress();
408
409 AliDebug(1, Form("Event %03d: fLoader->TreeH()->GetEntries() = %2d", iEvt, Int_t(fLoader->TreeH()->GetEntries())));
410
411 for (Int_t iTrack=0; iTrack<fLoader->TreeH()->GetEntries(); iTrack++) {
412 fLoader->TreeH()->GetEntry(iTrack);
413 Hits2SDigitsLocal(Hits(), GetSDigitsList(), iTrack); // convert these hits to a list of sdigits
414 }
415
416 fLoader->TreeS()->Fill();
417 fLoader->WriteSDigits("OVERWRITE");
418 ResetSDigits();
419
420 }
421 }
422
423 fLoader->UnloadHits();
424 fLoader->UnloadSDigits();
425
426 AliDebug(1,"Stop Hits2SDigits.");
427
428}
429
430//====================================================================================================================================================
431
432void AliMFT::Hits2SDigitsLocal(TClonesArray *hits, const TObjArray *pSDig, Int_t track) {
433
434 // Add sdigits of these hits to the list
435
58b93f36 436 AliDebug(1, "Entering Hits2SDigitsLocal");
820b4d9e 437
438 if (!fSegmentation) CreateGeometry();
439
440 TClonesArray *pSDigList[fNMaxPlanes];
441 for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) pSDigList[iPlane] = NULL;
442 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
443 pSDigList[iPlane] = (TClonesArray*) (*pSDig)[iPlane];
444 AliDebug(1,Form("Entries of pSDigList %3d; plane: %02d,",pSDigList[iPlane]->GetEntries(),iPlane));
445 if (!track && pSDigList[iPlane]->GetEntries()!=0) AliErrorClass("Some of sdigits lists is not empty");
446 }
447
448 for (Int_t iHit=0; iHit<hits->GetEntries(); iHit++) {
449
450 AliMFTHit *hit = (AliMFTHit*) hits->At(iHit);
451
58b93f36 452 // Creating "main digit"
453
454 AliMFTDigit *mainSDigit = new AliMFTDigit();
455 mainSDigit->SetEloss(hit->GetEloss());
456 mainSDigit->SetDetElemID(hit->GetDetElemID());
457 mainSDigit->SetPlane(hit->GetPlane());
458 mainSDigit->AddMCLabel(hit->GetTrack());
820b4d9e 459
460 Int_t xPixel = -1;
461 Int_t yPixel = -1;
58b93f36 462 if (fSegmentation->Hit2PixelID(hit->X(), hit->Y(), mainSDigit->GetDetElemID(), xPixel, yPixel)) {
463 mainSDigit->SetPixID(xPixel, yPixel, 0);
464 mainSDigit->SetPixWidth(fSegmentation->GetPixelSizeX(mainSDigit->GetDetElemID()),
465 fSegmentation->GetPixelSizeY(mainSDigit->GetDetElemID()),
466 fSegmentation->GetPixelSizeZ(mainSDigit->GetDetElemID()));
467 mainSDigit->SetPixCenter(fSegmentation->GetPixelCenterX(mainSDigit->GetDetElemID(), xPixel),
468 fSegmentation->GetPixelCenterY(mainSDigit->GetDetElemID(), yPixel),
469 fSegmentation->GetPixelCenterZ(mainSDigit->GetDetElemID(), 0));
470 new ((*fSideDigits)[fSideDigits->GetEntries()]) AliMFTDigit(*mainSDigit);
820b4d9e 471 AliDebug(1, Form("Created new sdigit (%f, %f, %f) from hit (%f, %f, %f)",
58b93f36 472 mainSDigit->GetPixelCenterX(), mainSDigit->GetPixelCenterY(), mainSDigit->GetPixelCenterZ(), hit->X(), hit->Y(), hit->Z()));
820b4d9e 473 }
474
58b93f36 475 // creating "side digits" to simulate the effect of charge dispersion
820b4d9e 476
820b4d9e 477 Double_t pi4 = TMath::Pi()/4.;
478 for (Int_t iStep=0; iStep<fNStepForChargeDispersion; iStep++) {
479 Double_t shift = (iStep+1) * fSingleStepForChargeDispersion;
480 for (Int_t iAngle=0; iAngle<8; iAngle++) {
481 Double_t shiftX = shift*TMath::Cos(iAngle*pi4);
482 Double_t shiftY = shift*TMath::Sin(iAngle*pi4);
58b93f36 483 if (fSegmentation->Hit2PixelID(hit->X()+shiftX, hit->Y()+shiftY, hit->GetDetElemID(), xPixel, yPixel)) {
820b4d9e 484 Bool_t digitExists = kFALSE;
58b93f36 485 for (Int_t iSideDigit=0; iSideDigit<fSideDigits->GetEntries(); iSideDigit++) {
486 if (xPixel==((AliMFTDigit*) fSideDigits->At(iSideDigit))->GetPixelX() &&
487 yPixel==((AliMFTDigit*) fSideDigits->At(iSideDigit))->GetPixelY()) {
488 digitExists = kTRUE;
489 break;
820b4d9e 490 }
491 }
492 if (!digitExists) {
58b93f36 493 AliMFTDigit *sideSDigit = new AliMFTDigit();
494 sideSDigit->SetEloss(0.);
495 sideSDigit->SetDetElemID(hit->GetDetElemID());
496 sideSDigit->SetPlane(hit->GetPlane());
497 sideSDigit->AddMCLabel(hit->GetTrack());
498 sideSDigit->SetPixID(xPixel, yPixel, 0);
499 sideSDigit->SetPixWidth(fSegmentation->GetPixelSizeX(sideSDigit->GetDetElemID()),
500 fSegmentation->GetPixelSizeY(sideSDigit->GetDetElemID()),
501 fSegmentation->GetPixelSizeZ(sideSDigit->GetDetElemID()));
502 sideSDigit->SetPixCenter(fSegmentation->GetPixelCenterX(sideSDigit->GetDetElemID(), xPixel),
503 fSegmentation->GetPixelCenterY(sideSDigit->GetDetElemID(), yPixel),
504 fSegmentation->GetPixelCenterZ(sideSDigit->GetDetElemID(), 0));
505 new ((*fSideDigits)[fSideDigits->GetEntries()]) AliMFTDigit(*sideSDigit);
820b4d9e 506 }
507 }
508 }
509 }
58b93f36 510
511 // ------------ In case we should simulate a rectangular pattern of pixel...
512
513 if (fSegmentation->GetPlane(mainSDigit->GetPlane())->HasPixelRectangularPatternAlongY()) {
514 for (Int_t iSDigit=0; iSDigit<fSideDigits->GetEntries(); iSDigit++) {
515 AliMFTDigit *mySDig = (AliMFTDigit*) (fSideDigits->At(iSDigit));
bcaf50eb 516 if (mySDig->GetPixelX()%2 == mySDig->GetPixelY()%2) { // both pair or both odd
58b93f36 517 xPixel = mySDig->GetPixelX();
518 yPixel = mySDig->GetPixelY()+1;
519 if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixel, yPixel)) {
520 AliMFTDigit *newSDigit = new AliMFTDigit();
521 newSDigit->SetEloss(0.);
522 newSDigit->SetDetElemID(mySDig->GetDetElemID());
523 newSDigit->SetPlane(mySDig->GetDetElemID());
524 newSDigit->SetPixID(xPixel, yPixel, 0);
525 newSDigit->SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit->GetDetElemID()),
526 fSegmentation->GetPixelSizeY(newSDigit->GetDetElemID()),
527 fSegmentation->GetPixelSizeZ(newSDigit->GetDetElemID()));
528 newSDigit->SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit->GetDetElemID(), xPixel),
529 fSegmentation->GetPixelCenterY(newSDigit->GetDetElemID(), yPixel),
530 fSegmentation->GetPixelCenterZ(newSDigit->GetDetElemID(), 0));
531 new ((*fSideDigits)[fSideDigits->GetEntries()]) AliMFTDigit(*newSDigit);
bcaf50eb 532 }
533 }
534 else { // pair-odd
58b93f36 535 xPixel = mySDig->GetPixelX();
536 yPixel = mySDig->GetPixelY()-1;
537 if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixel, yPixel)) {
538 AliMFTDigit *newSDigit = new AliMFTDigit();
539 newSDigit->SetEloss(0.);
540 newSDigit->SetDetElemID(mySDig->GetDetElemID());
541 newSDigit->SetPlane(mySDig->GetPlane());
542 newSDigit->SetPixID(xPixel, yPixel, 0);
543 newSDigit->SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit->GetDetElemID()),
544 fSegmentation->GetPixelSizeY(newSDigit->GetDetElemID()),
545 fSegmentation->GetPixelSizeZ(newSDigit->GetDetElemID()));
546 newSDigit->SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit->GetDetElemID(), xPixel),
547 fSegmentation->GetPixelCenterY(newSDigit->GetDetElemID(), yPixel),
548 fSegmentation->GetPixelCenterZ(newSDigit->GetDetElemID(), 0));
549 new ((*fSideDigits)[fSideDigits->GetEntries()]) AliMFTDigit(*newSDigit);
bcaf50eb 550 }
551 }
552 }
553 }
58b93f36 554
555 // -------- checking which pixels switched on have their diode actually within the charge dispersion radius
556
557 for (Int_t iSDigit=0; iSDigit<fSideDigits->GetEntries(); iSDigit++) {
558 AliMFTDigit *mySDig = (AliMFTDigit*) (fSideDigits->At(iSDigit));
559 Double_t distance = TMath::Sqrt(TMath::Power(mySDig->GetPixelCenterX()-hit->X(),2) + TMath::Power(mySDig->GetPixelCenterY()-hit->Y(),2));
560 if (fSegmentation->GetPlane(mySDig->GetPlane())->HasPixelRectangularPatternAlongY()) {
561 if (mySDig->GetPixelX()%2 == mySDig->GetPixelY()%2) { // both pair or both odd
562 if (distance<fChargeDispersion) {
563 AliDebug(1, Form("Created new sdigit (%f, %f, %f) from hit (%f, %f, %f)",
564 mySDig->GetPixelCenterX(), mySDig->GetPixelCenterY(), mySDig->GetPixelCenterZ(), hit->X(), hit->Y(), hit->Z()));
565 new ((*pSDigList[mySDig->GetPlane()])[pSDigList[mySDig->GetPlane()]->GetEntries()]) AliMFTDigit(*mySDig);
566 xPixel = mySDig->GetPixelX();
567 yPixel = mySDig->GetPixelY()+1;
568 if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixel, yPixel)) {
569 AliMFTDigit *newSDigit = new AliMFTDigit();
570 newSDigit->SetEloss(0.);
571 newSDigit->SetDetElemID(mySDig->GetDetElemID());
572 newSDigit->SetPlane(mySDig->GetPlane());
573 newSDigit->SetPixID(xPixel, yPixel, 0);
574 newSDigit->SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit->GetDetElemID()),
575 fSegmentation->GetPixelSizeY(newSDigit->GetDetElemID()),
576 fSegmentation->GetPixelSizeZ(newSDigit->GetDetElemID()));
577 newSDigit->SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit->GetDetElemID(), xPixel),
578 fSegmentation->GetPixelCenterY(newSDigit->GetDetElemID(), yPixel),
579 fSegmentation->GetPixelCenterZ(newSDigit->GetDetElemID(), 0));
580 AliDebug(1, Form("Created new sdigit (%f, %f, %f) from hit (%f, %f, %f)",
581 newSDigit->GetPixelCenterX(), newSDigit->GetPixelCenterY(), newSDigit->GetPixelCenterZ(), hit->X(), hit->Y(), hit->Z()));
582 new ((*pSDigList[newSDigit->GetPlane()])[pSDigList[newSDigit->GetPlane()]->GetEntries()]) AliMFTDigit(*newSDigit);
583 }
584 }
585 }
586 }
587 else {
588 if (distance<fChargeDispersion) {
589 AliDebug(1, Form("Created new sdigit (%f, %f, %f) from hit (%f, %f, %f)",
590 mySDig->GetPixelCenterX(), mySDig->GetPixelCenterY(), mySDig->GetPixelCenterZ(), hit->X(), hit->Y(), hit->Z()));
591 new ((*pSDigList[mySDig->GetPlane()])[pSDigList[mySDig->GetPlane()]->GetEntries()]) AliMFTDigit(*mySDig);
592 }
593 }
594 }
595
596 fSideDigits->Delete();
597
bcaf50eb 598 }
599
58b93f36 600 AliDebug(1,"Exiting Hits2SDigitsLocal");
820b4d9e 601
602}
603
604//====================================================================================================================================================
605
606void AliMFT::MakeBranch(Option_t *option) {
607
608 printf("AliMFT::MakeBranch(...)\n");
609
610 // Create Tree branches
611 AliDebug(1, Form("Start with option= %s.",option));
612
613 const Int_t kBufSize = 4000;
614
615 const Char_t *cH = strstr(option,"H");
616 const Char_t *cD = strstr(option,"D");
617 const Char_t *cS = strstr(option,"S");
618
619 if (cH && fLoader->TreeH()) {
620 CreateHits();
621 MakeBranchInTree(fLoader->TreeH(), "MFT", &fHits, kBufSize, 0);
622 }
623
624 if (cS && fLoader->TreeS()) {
625 CreateSDigits();
626 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) MakeBranchInTree(fLoader->TreeS(),
627 Form("Plane_%02d",iPlane),
628 &((*fSDigitsPerPlane)[iPlane]),
629 kBufSize, 0);
630 }
631
632 if (cD && fLoader->TreeD()) {
633 CreateDigits();
634 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) MakeBranchInTree(fLoader->TreeD(),
635 Form("Plane_%02d",iPlane),
636 &((*fDigitsPerPlane)[iPlane]),
637 kBufSize, 0);
638 }
639
640 AliDebug(1,"Stop.");
641
642}
643
644//====================================================================================================================================================
645
646void AliMFT::SetTreeAddress() {
647
648 AliDebug(1, "AliMFT::SetTreeAddress()");
649
650 //Set branch address for the Hits and Digits Tree.
651 AliDebug(1, "Start.");
652
653 AliDebug(1, Form("AliMFT::SetTreeAddress Hits fLoader->TreeH() = %p\n", fLoader->TreeH()));
654 if (fLoader->TreeH() && fLoader->TreeH()->GetBranch("MFT")) {
655 CreateHits();
656 fLoader->TreeH()->SetBranchAddress("MFT", &fHits);
657 }
658
659 AliDebug(1, Form("AliMFT::SetTreeAddress SDigits fLoader->TreeS() = %p\n", fLoader->TreeS()));
660 if (fLoader->TreeS() && fLoader->TreeS()->GetBranch("Plane_00")) {
661 CreateSDigits();
662 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
663 fLoader->TreeS()->SetBranchAddress(Form("Plane_%02d",iPlane), &((*fSDigitsPerPlane)[iPlane]));
664 }
665 }
666
667 AliDebug(1, Form("AliMFT::SetTreeAddress Digits fLoader->TreeD() = %p\n", fLoader->TreeD()));
668 if (fLoader->TreeD() && fLoader->TreeD()->GetBranch("Plane_00")) {
669 CreateDigits();
670 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
671 fLoader->TreeD()->SetBranchAddress(Form("Plane_%02d",iPlane), &((*fDigitsPerPlane)[iPlane]));
672 }
673 }
674
675 AliDebug(1, Form("AliMFT::SetTreeAddress RecPoints fLoader->TreeR() = %p\n", fLoader->TreeR()));
676 if (fLoader->TreeR() && fLoader->TreeR()->GetBranch("Plane_00")) {
677 CreateRecPoints();
678 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
679 fLoader->TreeR()->SetBranchAddress(Form("Plane_%02d",iPlane), &((*fRecPointsPerPlane)[iPlane]));
680 }
681 }
682
683 AliDebug(1,"Stop.");
684
685}
686
687//====================================================================================================================================================
688
689void AliMFT::SetGeometry() {
690
691 printf("AliMFT::SetGeometry\n");
692
693 fSegmentation = new AliMFTSegmentation(fNameGeomFile.Data());
694
695 fNPlanes = fSegmentation->GetNPlanes();
696
697}
698
699//====================================================================================================================================================
700
701void AliMFT::CreateHits() {
702
703 // create array of hits
704
705 AliDebug(1, "AliMFT::CreateHits()");
706
707 if (fHits) return;
708 fHits = new TClonesArray("AliMFTHit");
709
710}
711
712//====================================================================================================================================================
713
714void AliMFT::CreateSDigits() {
715
716 // create sdigits list
717
718 AliDebug(1, "AliMFT::CreateSDigits()");
719
720 if (fSDigitsPerPlane) return;
721 fSDigitsPerPlane = new TObjArray(fNPlanes);
722 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) fSDigitsPerPlane->AddAt(new TClonesArray("AliMFTDigit"), iPlane);
723
724 fSideDigits = new TClonesArray("AliMFTDigit");
725
726}
727
728//====================================================================================================================================================
729
730void AliMFT::CreateDigits() {
731
732 // create digits list
733
734 AliDebug(1, "AliMFT::CreateDigits()");
735
736 if (fDigitsPerPlane) return;
737 fDigitsPerPlane = new TObjArray(fNPlanes);
738 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) fDigitsPerPlane->AddAt(new TClonesArray("AliMFTDigit"), iPlane);
739
740}
741
742//====================================================================================================================================================
743
744void AliMFT::CreateRecPoints() {
745
746 // create recPoints list
747
748 AliDebug(1, "AliMFT::CreateRecPoints()");
749
750 if (fRecPointsPerPlane) return;
751 fRecPointsPerPlane = new TObjArray(fNPlanes);
752 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) fRecPointsPerPlane->AddAt(new TClonesArray("AliMFTCluster"), iPlane);
753
754}
755
756//====================================================================================================================================================