correct usage of output loader
[u/mrichter/AliRoot.git] / MUON / AliMUONDigitizerv2.cxx
CommitLineData
d1775029 1
2// Do the Digitization (Digit) from summable Digits (SDigit)
3// Allow the merging of signal file with background file(s).
4
5#include <Riostream.h>
6#include <TDirectory.h>
7#include <TFile.h>
8#include <TObjArray.h>
9#include <TPDGCode.h>
10#include <TTree.h>
11#include <TMath.h>
12
13#include "AliMUON.h"
14#include "AliMUONChamber.h"
15#include "AliMUONConstants.h"
16#include "AliMUONDigit.h"
17#include "AliMUONDigitizerv2.h"
18#include "AliMUONHit.h"
19#include "AliMUONHitMapA1.h"
20#include "AliMUONPadHit.h"
21#include "AliMUONTransientDigit.h"
22#include "AliRun.h"
23#include "AliRunDigitizer.h"
24#include "AliRunLoader.h"
25#include "AliLoader.h"
26
27ClassImp(AliMUONDigitizerv2)
28
29//___________________________________________
30AliMUONDigitizerv2::AliMUONDigitizerv2() :
31 AliDigitizer(),
32 fHitMap(0),
33 fTDList(0),
34 fTDCounter(0),
35 fDebug(0),
36 fMask(0),
37 fSignal(0)
38{
39// Default ctor - don't use it
40 if (GetDebug()>2)
41 cerr<<"AliMUONDigitizerv2::AliMUONDigitizerv2"
42 <<"(AliRunDigitizer* manager) was processed"<<endl;
43}
44
45//___________________________________________
46AliMUONDigitizerv2::AliMUONDigitizerv2(AliRunDigitizer* manager):
47 AliDigitizer(manager),
48 fHitMap(0),
49 fTDList(0),
50 fTDCounter(0),
51 fDebug(0),
52 fMask(0),
53 fSignal(0)
54{
55// ctor which should be used
56 if (GetDebug()>2)
57 cerr<<"AliMUONDigitizerv2::AliMUONDigitizerv2"
58 <<"(AliRunDigitizer* manager) was processed"<<endl;
59}
60
61//------------------------------------------------------------------------
62AliMUONDigitizerv2::~AliMUONDigitizerv2()
63{
64// Destructor
65}
66
67//------------------------------------------------------------------------
68void AliMUONDigitizerv2::AddTransientDigit(AliMUONTransientDigit * mTD)
69{
70 // Choosing the maping of the cathode plane of the chamber:
71 Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
72 fTDList->AddAtAndExpand(mTD, fTDCounter);
73 fHitMap[iNchCpl]->SetHit( mTD->PadX(), mTD->PadY(), fTDCounter);
74 fTDCounter++;
75}
76
77//------------------------------------------------------------------------
78Bool_t AliMUONDigitizerv2::ExistTransientDigit(AliMUONTransientDigit * mTD)
79{
80 // Choosing the maping of the cathode plane of the chamber:
81 Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
82 return( fHitMap[iNchCpl]->TestHit(mTD->PadX(), mTD->PadY()) );
83}
84
85//------------------------------------------------------------------------
86Bool_t AliMUONDigitizerv2::Init()
87{
88
89// Initialization
90 if (GetDebug()>2) Info("Init","AliMUONDigitizerv2::Init() starts");
91
92 //Loaders (We assume input0 to be the output too)
93 AliRunLoader * runloader; // Input loader
94 AliLoader * gime;
95
96 // Getting runloader
97 runloader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
98 if (runloader == 0x0) {
99 Error("Init","RunLoader is not in input file 0");
100 return kFALSE; // RunDigitizer is not working.
101 }
102 // Getting MUONloader
103 gime = runloader->GetLoader("MUONLoader");
104 gime->LoadSDigits("READ");
105 gime->LoadDigits("RECREATE");
106
107 return kTRUE;
108}
109
110//------------------------------------------------------------------------
111void AliMUONDigitizerv2::UpdateTransientDigit(Int_t /*track*/, AliMUONTransientDigit * mTD)
112{
113 // Choosing the maping of the cathode plane of the chamber:
114 Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
115 AliMUONTransientDigit *pdigit =
116 static_cast<AliMUONTransientDigit*>(fHitMap[iNchCpl]->GetHit(mTD->PadX(),mTD->PadY()));
117 // update charge
118 //
119 Int_t iqpad = mTD->Signal(); // charge per pad
120 pdigit->AddSignal(iqpad);
121 pdigit->AddPhysicsSignal(iqpad);
122 // update list of tracks
123 //
124 // List of tracks and trackcharge
125 Int_t itrack;
126 for(itrack=0;itrack<kMAXTRACKS;itrack++) {
127 pdigit->UpdateTrackList(mTD->Track(itrack), mTD->TrackCharge(itrack));
128 }
129}
130
131//--------------------------------------------------------------------------
132void AliMUONDigitizerv2::MakeTransientDigitFromSDigit(Int_t iChamber, AliMUONDigit * sDigit)
133{
134 Int_t digits[6];
135
136
137 //
138 // Creating a new TransientDigits from SDigit
139 digits[0] = sDigit->PadX(); // Padx of the Digit
140 digits[1] = sDigit->PadY(); // Pady of the Digit
141 digits[2] = sDigit->Cathode()+1; // Cathode plane
142 digits[3] = sDigit->Signal(); // Induced charge in the Pad
143 if (fSignal) digits[4] = sDigit->Signal();
144 else digits[4] = 0;
145 digits[5] = sDigit->Hit(); // Hit number in the list
146 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::MakeTransientDigitFromSDigit " <<
147 "PadX "<< digits[0] << " " <<
148 "PadY "<< digits[1] << " " <<
149 "Plane " << digits[2] << " " <<
150 "Charge " << digits[3] <<" " <<
151 "Hit " << digits[5] << endl;
152
153 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::MakeTransientDigitFromSDigit Creating AliMUONTransientDigit"<<endl;
154 AliMUONTransientDigit * mTD = new AliMUONTransientDigit(iChamber, digits);
155
156 // List of tracks and trackcharge
157 Int_t itrack;
158 for(itrack=0;itrack<kMAXTRACKS;itrack++) {
159 mTD->AddToTrackList(fMask+sDigit->Track(itrack), sDigit->TrackCharge(itrack));
160 }
161
162 if (!ExistTransientDigit(mTD)) {
163 AddTransientDigit(mTD);
164 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::MakeTransientDigitFromSDigit Adding TransientDigit"<<endl;
165 }
166 else {
167 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::MakeTransientDigitFromSDigit updating TransientDigit"<<endl;
168 UpdateTransientDigit(0, mTD);
169 delete mTD;
170 }
171}
172
173//-----------------------------------------------------------------------
174void AliMUONDigitizerv2::Exec(Option_t* option)
175{
176 TString optionString = option;
177 if (optionString.Data() == "deb") {
178 Info("Digitize","Called with option deb ");
179 fDebug = 3;
180 }
181
182 AliMUONChamber* chamber;
183 AliSegmentation* c1Segmentation; //Cathode plane c1 of the chamber
184 AliSegmentation* c2Segmentation; //Cathode place c2 of the chamber
185
186 if (GetDebug()>2) Info("Exec","AliMUONDigitizerv2::Exec() starts");
187 fTDList = new TObjArray;
188
189 //Loaders (We assume input0 to be the output too)
be42eafc 190 AliRunLoader *runloader, *runloaderOut; // Input / Output loaders
191 AliLoader *gime, *gimeOut;
d1775029 192
193 // Getting runloader
194 runloader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
195 if (runloader == 0x0) {
196 Error("Exec","RunLoader is not in input file 0");
197 return; // RunDigitizer is not working.
198 }
199 // Getting MUONloader
200 gime = runloader->GetLoader("MUONLoader");
201 if (gime->TreeS()==0x0) {
202 if (GetDebug()>2) Info("Exec","TreeS is not loaded yet. Loading...");
203 gime->LoadSDigits("READ");
be42eafc 204 if (GetDebug()>2) Info("Exec","Now treeS is %#x. MUONLoader is %#x",gime->TreeS(),gime);
d1775029 205 }
206
207 if (GetDebug()>2) Info("Exec","Loaders ready");
208
209 if (runloader->GetAliRun() == 0x0) runloader->LoadgAlice();
210 gAlice = runloader->GetAliRun();
211
212 // Getting Module MUON
213 AliMUON *pMUON = (AliMUON *) gAlice->GetDetector("MUON");
214 if (!pMUON) {
215 Error("Digitize","Module MUON not found in the input file");
216 return;
217 }
218 // Getting Muon data
219 AliMUONData * muondata = pMUON->GetMUONData();
d1775029 220
221 // Loading Event
222 Int_t currentevent = fManager->GetOutputEventNr();
223
224 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::Exec() Event Number is "<<currentevent <<endl;
225 if ( (currentevent<10) ||
226 (Int_t(TMath::Log10(currentevent)) == TMath::Log10(currentevent) ) )
227 cout <<"ALiMUONDigitizerv2::Exec() Event Number is "<< currentevent <<endl;
228
be42eafc 229 // Output runloader
230 runloaderOut = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
231 runloaderOut->SetNumberOfEventsPerFile(runloader->GetNumberOfEvents());
232 gimeOut = runloaderOut->GetLoader("MUONLoader");
233
234
235 if (gimeOut->TreeD() == 0x0) {
236 gimeOut->MakeDigitsContainer();
d1775029 237 }
be42eafc 238 TTree* treeD = gimeOut->TreeD();
239 muondata->SetLoader(gimeOut);
d1775029 240 muondata->MakeBranch("D");
241 muondata->SetTreeAddress("D");
242
243 // Array of pointer of the AliMUONHitMapA1:
244 // two HitMaps per chamber, or one HitMap per cahtode plane
245 fHitMap= new AliMUONHitMapA1* [2*AliMUONConstants::NCh()];
246
247 //Loop over chambers for the definition AliMUONHitMap
248 for (Int_t i=0; i<AliMUONConstants::NCh(); i++) {
249 chamber = &(pMUON->Chamber(i));
250 c1Segmentation = chamber->SegmentationModel(1); // Cathode plane 1
251 fHitMap[i] = new AliMUONHitMapA1(c1Segmentation, fTDList);
252 c2Segmentation = chamber->SegmentationModel(2); // Cathode plane 2
253 fHitMap[i+AliMUONConstants::NCh()] = new AliMUONHitMapA1(c2Segmentation, fTDList);
254 }
255
256// Loop over files to merge and to digitize
257 fSignal = kTRUE;
258 for (Int_t inputFile=0; inputFile<fManager->GetNinputs(); inputFile++) {
259 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::Exec() Input File is "<<inputFile<<endl;
be42eafc 260 if (inputFile == 0) {
d1775029 261 muondata->SetLoader(gime);
be42eafc 262 muondata->SetTreeAddress("S");
263 } else {
264 // Connect MUON Hit branch
265 fSignal = kFALSE;
266 runloader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(inputFile));
267 if (runloader == 0x0) {
268 cerr<<"AliMUONDigitizerv2::Digitize() RunLoader for inputFile "<<inputFile<< " not found !!! "<<endl;
269 }
270 gime = runloader->GetLoader("MUONLoader");
271 if (gime->TreeS() == 0x0) gime->LoadSDigits("READ");
272 muondata->SetLoader(gime);
273 muondata->SetTreeAddress("S");
d1775029 274 }
be42eafc 275
d1775029 276 // Setting the address
277 TTree *treeS = gime->TreeS();
278 if (treeS == 0x0) {
279 Error("Digitize","Can not get TreeS from input %d",inputFile);
280 Info("Digitize","Now treeS is %#x. MUONLoader is %#x",gime->TreeS(),gime);
281 return;
282 }
283 if (GetDebug()>2) {
284 cerr<<"AliMUONDigitizerv2::Exec inputFile is "<<inputFile<<" "<<endl;
285 cerr<<"AliMUONDigitizerv2::Exec treeS" << treeS <<endl;
286 }
be42eafc 287
d1775029 288 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::Exec Setting tree addresses"<<endl;
be42eafc 289
d1775029 290 fMask = fManager->GetMask(inputFile);
291 //
292 // Loop over SDigits
293 Int_t ndig,k;
294 AliMUONDigit *sDigit;
295 TClonesArray * muonSDigits;
296 for (Int_t ich = 0; ich < AliMUONConstants::NCh(); ich++) { // loop over chamber
297 muondata->ResetSDigits();
298 muondata->GetCathodeS(0);
299 //TClonesArray *
300 muonSDigits = muondata->SDigits(ich);
301 ndig = muonSDigits->GetEntriesFast();
be42eafc 302 // printf("\n 1 Found %d Sdigits in %p chamber %d", ndig, muonSDigits,ich);
d1775029 303 Int_t n = 0;
304 for (k = 0; k < ndig; k++) {
be42eafc 305 sDigit = (AliMUONDigit*) muonSDigits->UncheckedAt(k);
306 MakeTransientDigitFromSDigit(ich,sDigit);
d1775029 307 }
308 muondata->ResetSDigits();
309 muondata->GetCathodeS(1);
310 //TClonesArray *
311 muonSDigits = muondata->SDigits(ich);
312 ndig=muonSDigits->GetEntriesFast();
be42eafc 313 // printf("\n 2 Found %d Sdigits in %p chamber %d", ndig, muonSDigits,ich);
d1775029 314 n = 0;
315 for (k = 0; k < ndig; k++) {
be42eafc 316 sDigit = (AliMUONDigit*) muonSDigits->UncheckedAt(k);
317 MakeTransientDigitFromSDigit(ich,sDigit);
d1775029 318 }
319
320 } // SDigits loop end loop over chamber
321 } // end file loop
322 if (GetDebug()>2) cerr<<"AliMUONDigitizerv2::Exec End of Sdigits, file loops"<<endl;
be42eafc 323
d1775029 324 // Loop on cathodes
325 Int_t icat;
326 for(icat=0; icat<2; icat++) {
327 //
328 // Filling Digit List
329 Int_t tracks[kMAXTRACKS];
330 Int_t charges[kMAXTRACKS];
331 Int_t nentries = fTDList->GetEntriesFast();
332 Int_t digits[6];
333 for (Int_t nent = 0; nent < nentries; nent++) {
334 AliMUONTransientDigit *address = (AliMUONTransientDigit*)fTDList->At(nent);
335 if (address == 0) continue;
336 Int_t ich = address->Chamber();
337 Int_t q = address->Signal();
338 chamber = &(pMUON->Chamber(ich));
339 //
340 // Digit Response (noise, threshold, saturation, ...)
341 AliMUONResponse * response = chamber->ResponseModel();
342 q = response->DigitResponse(q,address);
343
344 if (!q) continue;
345
346 digits[0] = address->PadX();
347 digits[1] = address->PadY();
348 digits[2] = address->Cathode()-1;
349 digits[3] = q;
350 digits[4] = address->Physics();
351 digits[5] = address->Hit();
352
353 Int_t nptracks = address->GetNTracks();
354
355 if (nptracks > kMAXTRACKS) {
356 if (GetDebug() >0) {
357 cerr<<"AliMUONDigitizer:Exec nptracks > 10 "<<nptracks;
358 cerr<<"reset to max value "<<kMAXTRACKS<<endl;
359 }
360 nptracks = kMAXTRACKS;
361 }
362 if (nptracks > 2 && GetDebug() >2) {
363 cerr<<"AliMUONDigitizer::Exec nptracks > 2 "<<nptracks<<endl;
364 // printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q);
365 }
366 for (Int_t tr = 0; tr < nptracks; tr++) {
367 tracks[tr] = address->GetTrack(tr);
368 charges[tr] = address->GetCharge(tr);
369 } //end loop over list of tracks for one pad
370 // Sort list of tracks according to charge
371 if (nptracks > 1) {
372 SortTracks(tracks,charges,nptracks);
373 }
374 if (nptracks < kMAXTRACKS ) {
375 for (Int_t i = nptracks; i < kMAXTRACKS; i++) {
376 tracks[i] = -1;
377 charges[i] = 0;
378 }
379 }
380
381 // Add digits
382 if (GetDebug()>3) cerr<<"AliMUONDigitzerv2::Exec TransientDigit to Digit"<<endl;
383 if ( digits[2] == icat ) muondata->AddDigit(ich,tracks,charges,digits);
384// printf("test rm ich %d padX %d padY %d \n",ich, digits[0], digits[1]);
385 }
be42eafc 386 muondata->SetLoader(gimeOut);
d1775029 387 muondata->Fill("D");
388 muondata->ResetDigits();
389 } // end loop cathode
390 fTDList->Delete();
391
392 for(Int_t ii = 0; ii < 2*AliMUONConstants::NCh(); ++ii) {
393 if (fHitMap[ii]) {
394 delete fHitMap[ii];
395 fHitMap[ii] = 0;
396 }
397 }
398
399 if (GetDebug()>2)
400 cerr<<"AliMUONDigitizer::Exec: writing the TreeD: "
401 <<treeD->GetName()<<endl;
be42eafc 402
403 gimeOut->WriteDigits("OVERWRITE");
d1775029 404 delete [] fHitMap;
405 delete fTDList;
406 muondata->ResetSDigits();
407 gime->UnloadSDigits();
be42eafc 408 gimeOut->UnloadDigits();
d1775029 409}
410//------------------------------------------------------------------------
411void AliMUONDigitizerv2::SortTracks(Int_t *tracks,Int_t *charges,Int_t ntr)
412{
413 //
414 // Sort the list of tracks contributing to a given digit
415 // Only the 3 most significant tracks are acctually sorted
416 //
417
418 //
419 // Loop over signals, only 3 times
420 //
421
422 Int_t qmax;
423 Int_t jmax;
424 Int_t idx[3] = {-2,-2,-2};
425 Int_t jch[3] = {-2,-2,-2};
426 Int_t jtr[3] = {-2,-2,-2};
427 Int_t i,j,imax;
428
429 if (ntr<3) imax=ntr;
430 else imax=3;
431 for(i=0;i<imax;i++){
432 qmax=0;
433 jmax=0;
434
435 for(j=0;j<ntr;j++){
436
437 if((i == 1 && j == idx[i-1])
438 ||(i == 2 && (j == idx[i-1] || j == idx[i-2]))) continue;
439
440 if(charges[j] > qmax) {
441 qmax = charges[j];
442 jmax=j;
443 }
444 }
445
446 if(qmax > 0) {
447 idx[i]=jmax;
448 jch[i]=charges[jmax];
449 jtr[i]=tracks[jmax];
450 }
451
452 }
453
454 for(i=0;i<3;i++){
455 if (jtr[i] == -2) {
456 charges[i]=0;
457 tracks[i]=0;
458 } else {
459 charges[i]=jch[i];
460 tracks[i]=jtr[i];
461 }
462 }
463}