Bug on digit output. Now treeD event per cathode
[u/mrichter/AliRoot.git] / MUON / AliMUONDigitizerv1.cxx
CommitLineData
8dbbc4e3 1
2
3#include <Riostream.h>
4#include <TDirectory.h>
5#include <TFile.h>
6#include <TObjArray.h>
7#include <TPDGCode.h>
8#include <TTree.h>
9
10#include "AliMUON.h"
11#include "AliMUONChamber.h"
12#include "AliMUONConstants.h"
13#include "AliMUONDigit.h"
14#include "AliMUONDigitizerv1.h"
15#include "AliMUONHit.h"
16#include "AliMUONHitMapA1.h"
17#include "AliMUONPadHit.h"
18#include "AliMUONTransientDigit.h"
19#include "AliRun.h"
20#include "AliRunDigitizer.h"
21
22ClassImp(AliMUONDigitizerv1)
23
24//___________________________________________
25AliMUONDigitizerv1::AliMUONDigitizerv1() :AliDigitizer()
26{
27// Default ctor - don't use it
28 fHitMap = 0;
29 fTDList = 0;
30 if (GetDebug()>2)
31 cerr<<"AliMUONDigitizerv1::AliMUONDigitizerv1"
32 <<"(AliRunDigitizer* manager) was processed"<<endl;
33}
34
35//___________________________________________
36AliMUONDigitizerv1::AliMUONDigitizerv1(AliRunDigitizer* manager)
37 :AliDigitizer(manager)
38{
39// ctor which should be used
40 fHitMap = 0;
41 fTDList = 0;
42 fDebug = 0;
43 fHits = new TClonesArray("AliMUONHit",1000);
44 if (GetDebug()>2)
45 cerr<<"AliMUONDigitizerv1::AliMUONDigitizerv1"
46 <<"(AliRunDigitizer* manager) was processed"<<endl;
47}
48
49//------------------------------------------------------------------------
50AliMUONDigitizerv1::~AliMUONDigitizerv1()
51{
52// Destructor
53 delete fHits;
54}
55
56//------------------------------------------------------------------------
57void AliMUONDigitizerv1::AddTransientDigit(AliMUONTransientDigit * mTD)
58{
59 // Choosing the maping of the cathode plane of the chamber:
60 Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
61 fTDList->AddAtAndExpand(mTD, fTDCounter);
62 fHitMap[iNchCpl]->SetHit( mTD->PadX(), mTD->PadY(), fTDCounter);
63 fTDCounter++;
64}
65
66//------------------------------------------------------------------------
67Bool_t AliMUONDigitizerv1::ExistTransientDigit(AliMUONTransientDigit * mTD)
68{
69 // Choosing the maping of the cathode plane of the chamber:
70 Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
71 return( fHitMap[iNchCpl]->TestHit(mTD->PadX(), mTD->PadY()) );
72}
73
74//------------------------------------------------------------------------
75Bool_t AliMUONDigitizerv1::Init()
76{
77// Initialization
78 return kTRUE;
79}
80
81//------------------------------------------------------------------------
82void AliMUONDigitizerv1::UpdateTransientDigit(Int_t track, AliMUONTransientDigit * mTD)
83{
84 // Choosing the maping of the cathode plane of the chamber:
85 Int_t iNchCpl= mTD->Chamber() + (mTD->Cathode()-1) * AliMUONConstants::NCh();
86 AliMUONTransientDigit *pdigit =
87 static_cast<AliMUONTransientDigit*>(fHitMap[iNchCpl]->GetHit(mTD->PadX(),mTD->PadY()));
88 // update charge
89 //
90 Int_t iqpad = mTD->Signal(); // charge per pad
91 pdigit->AddSignal(iqpad);
92 pdigit->AddPhysicsSignal(iqpad);
93 // update list of tracks
94 //
95 Int_t charge;
96 track=+ fMask;
97 if (fSignal) charge = iqpad;
98 else charge = kBgTag;
99 pdigit->UpdateTrackList(track,charge);
100}
101
102
103//--------------------------------------------------------------------------
104void AliMUONDigitizerv1::MakeTransientDigit(Int_t track, Int_t iHit, AliMUONHit * mHit)
105{
106 AliMUON *pMUON = (AliMUON *) gAlice->GetModule("MUON");
107 if (!pMUON) {
108 cerr<<"AliMUONDigitizerv1::Digitize Error:"
109 <<" module MUON not found in the input file"<<endl;
110 }
111 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit starts"<<endl;
112 Int_t ichamber = mHit->Chamber()-1;
113 AliMUONChamber & chamber = pMUON->Chamber(ichamber);
114 Float_t xhit = mHit->X();
115 Float_t yhit = mHit->Y();
116 Float_t zhit = mHit->Z();
117 Float_t eloss= mHit->Eloss();
118 Float_t tof = mHit->Age();
119 // Variables for chamber response from AliMUONChamber::DisIntegration
120 Float_t newdigit[6][500]; // Pad information
121 Int_t nnew=0; // Number of touched Pads per hit
122 Int_t digits[6];
123
124 //
125 // Calls the charge disintegration method of the current chamber
126 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit calling AliMUONChamber::DisIngtegration starts"<<endl;
127 chamber.DisIntegration(eloss, tof, xhit, yhit, zhit, nnew, newdigit);
128 // Creating a new TransientDigits from hit
129 for(Int_t iTD=0; iTD<nnew; iTD++) {
130 digits[0] = Int_t(newdigit[1][iTD]); // Padx of the Digit
131 digits[1] = Int_t(newdigit[2][iTD]); // Pady of the Digit
132 digits[2] = Int_t(newdigit[5][iTD]); // Cathode plane
133 digits[3] = Int_t(newdigit[3][iTD]); // Induced charge in the Pad
134 if (fSignal) digits[4] = Int_t(newdigit[3][iTD]);
135 else digits[4] = 0;
136 digits[5] = iHit+fMask; // Hit number in the list
137 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit " <<
138 "PadX "<< digits[0] << " " <<
139 "PadY "<< digits[1] << " " <<
140 "Plane " << digits[2] << " " <<
141 "Charge " << digits[3] <<" " <<
142 "Hit " << digits[5] << endl;
143 // list of tracks
144 Int_t charge;
145 track += fMask;
146 if (fSignal) charge = digits[3];
147 else charge = kBgTag;
148 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit Creating AliMUONTransientDigit"<<endl;
149 AliMUONTransientDigit * mTD = new AliMUONTransientDigit(ichamber, digits);
150 mTD->AddToTrackList(track,charge);
151
152 if (!ExistTransientDigit(mTD)) {
153 AddTransientDigit(mTD);
154 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit Adding TransientDigit"<<endl;
155 }
156 else {
157 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::MakeTransientDigit updating TransientDigit"<<endl;
158 UpdateTransientDigit(track, mTD);
159 delete mTD;
160 }
161 }
162}
163//-----------------------------------------------------------------------
164void AliMUONDigitizerv1::Exec(Option_t* option)
165{
166 TString optionString = option;
167 if (optionString.Data() == "deb") {
168 cout<<"AliMUONDigitizerv1::Exec: called with option deb "<<endl;
169 fDebug = 3;
170 }
171
172 AliMUONChamber* chamber;
173 AliSegmentation* c1Segmentation; //Cathode plane c1 of the chamber
174 AliSegmentation* c2Segmentation; //Cathode place c2 of the chamber
175
176 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Digitize() starts"<<endl;
177 fTDList = new TObjArray;
178
179 // Getting Module MUON
180 AliMUON *pMUON = (AliMUON *) gAlice->GetModule("MUON");
181 if (!pMUON) {
182 cerr<<"AliMUONDigitizerv1::Digitize Error:"
183 <<" module MUON not found in the input file"<<endl;
184 return;
185 }
186 // New branch for MUON digit in the tree of digits
187 pMUON->MakeBranchInTreeD(fManager->GetTreeD());
188
189 // Array of pointer of the AliMUONHitMapA1:
190 // two HitMaps per chamber, or one HitMap per cahtode plane
191 fHitMap= new AliMUONHitMapA1* [2*AliMUONConstants::NCh()];
192
193 //Loop over chambers for the definition AliMUONHitMap
194 for (Int_t i=0; i<AliMUONConstants::NCh(); i++) {
195 chamber = &(pMUON->Chamber(i));
196 c1Segmentation = chamber->SegmentationModel(1); // Cathode plane 1
197 fHitMap[i] = new AliMUONHitMapA1(c1Segmentation, fTDList);
198 c2Segmentation = chamber->SegmentationModel(2); // Cathode plane 2
199 fHitMap[i+AliMUONConstants::NCh()] = new AliMUONHitMapA1(c2Segmentation, fTDList);
200 }
201
202// Loop over files to merge and to digitize
203 fSignal = kTRUE;
204 for (Int_t inputFile=0; inputFile<fManager->GetNinputs(); inputFile++) {
205 // Connect MUON Hit branch
206 if (inputFile > 0 ) fSignal = kFALSE;
207 TBranch *branchHits = 0;
208 TTree *treeH = fManager->GetInputTreeH(inputFile);
209 if (GetDebug()>2) {
210 cerr<<"AliMUONDigitizerv1::Exec inputFile is "<<inputFile<<" "<<endl;
211 cerr<<"AliMUONDigitizerv1::Exec treeH, fHits "<<treeH<<" "<<fHits<<endl;
212 }
213 if (treeH && fHits) {
214 branchHits = treeH->GetBranch("MUON");
215 if (branchHits) {
391cd891 216 fHits->Clear();
8dbbc4e3 217 branchHits->SetAddress(&fHits);
218 }
219 else
220 Error("Exec","branch MUON was not found");
221 }
222 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Exec branchHits = "<<branchHits<<endl;
223
224 fMask = fManager->GetMask(inputFile);
225 //
226 // Loop over tracks
227 Int_t itrack;
228 Int_t ntracks = (Int_t) treeH->GetEntries();
229 for (itrack = 0; itrack < ntracks; itrack++) {
230 if (GetDebug()>2) cerr<<"AliMUONDigitizerv1::Exec itrack = "<<itrack<<endl;
5989ad6a 231 fHits->Clear();
8dbbc4e3 232 branchHits->GetEntry(itrack);
233 //
234 // Loop over hits
235 Int_t ihit, ichamber;
236 AliMUONHit* mHit;
237 for(ihit = 0; ihit < fHits->GetEntriesFast(); ihit++) {
238 mHit = static_cast<AliMUONHit*>(fHits->At(ihit));
239 ichamber = mHit->Chamber()-1; // chamber number
240 if (ichamber > AliMUONConstants::NCh()-1) {
241 cerr<<"AliMUONDigitizer: ERROR: "
242 <<"fNch > AliMUONConstants::NCh()-1, fNch, NCh(): "
243 <<ichamber<<", "<< AliMUONConstants::NCh()<<endl;
244 return;
245 }
246 chamber = &(pMUON->Chamber(ichamber));
247 //
248 //Dumping Hit content:
249 if (GetDebug()>2) {
250 cerr<<"AliMuonDigitizerv1::Exec ihit, ichamber, x, y, z, eloss " <<
251 ihit << " " <<
252 mHit->Chamber() << " " <<
253 mHit->X() << " " <<
254 mHit->Y() << " " <<
255 mHit->Z() << " " <<
256 mHit->Eloss() << " " << endl;
257 }
258 //
259 // Inititializing Correlation
260 chamber->ChargeCorrelationInit();
261 if (ichamber < AliMUONConstants::NTrackingCh()) {
262 // Tracking Chamber
263 // Initialize hit position (cursor) in the segmentation model
264 chamber->SigGenInit(mHit->X(), mHit->Y(), mHit->Z());
265 } else {
266 // Trigger Chamber
267 }
268 MakeTransientDigit(itrack, ihit, mHit);
269 } // hit loop
270 } // track loop
271 } // end file loop
272 if (GetDebug()>2) cerr<<"AliMUONDigitizer::Exec End of hits, track and file loops"<<endl;
273
391cd891 274
275 // Loop on cathodes
276 Int_t icat;
277 for(icat=0; icat<2; icat++) {
8dbbc4e3 278 //
391cd891 279 // Filling Digit List
280 Int_t tracks[kMAXTRACKS];
281 Int_t charges[kMAXTRACKS];
282 Int_t nentries = fTDList->GetEntriesFast();
283 Int_t digits[6];
284 for (Int_t nent = 0; nent < nentries; nent++) {
285 AliMUONTransientDigit *address = (AliMUONTransientDigit*)fTDList->At(nent);
286 if (address == 0) continue;
287 Int_t ich = address->Chamber();
288 Int_t q = address->Signal();
289 chamber = &(pMUON->Chamber(ich));
290 //
291 // Digit Response (noise, threshold, saturation, ...)
292 AliMUONResponse * response = chamber->ResponseModel();
293 q = response->DigitResponse(q,address);
294
295 if (!q) continue;
296
297 digits[0] = address->PadX();
298 digits[1] = address->PadY();
299 digits[2] = address->Cathode()-1;
300 digits[3] = q;
301 digits[4] = address->Physics();
302 digits[5] = address->Hit();
303
304 Int_t nptracks = address->GetNTracks();
305
306 if (nptracks > kMAXTRACKS) {
307 if (GetDebug() >0) {
308 cerr<<"AliMUONDigitizer:Exec nptracks > 10 "<<nptracks;
309 cerr<<"reset to max value "<<kMAXTRACKS<<endl;
310 }
311 nptracks = kMAXTRACKS;
8dbbc4e3 312 }
391cd891 313 if (nptracks > 2 && GetDebug() >2) {
314 cerr<<"AliMUONDigitizer::Exec nptracks > 2 "<<nptracks<<endl;
315 // printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q);
316 }
317 for (Int_t tr = 0; tr < nptracks; tr++) {
318 tracks[tr] = address->GetTrack(tr);
319 charges[tr] = address->GetCharge(tr);
320 } //end loop over list of tracks for one pad
321 // Sort list of tracks according to charge
322 if (nptracks > 1) {
323 SortTracks(tracks,charges,nptracks);
8dbbc4e3 324 }
391cd891 325 if (nptracks < kMAXTRACKS ) {
326 for (Int_t i = nptracks; i < kMAXTRACKS; i++) {
327 tracks[i] = 0;
328 charges[i] = 0;
329 }
330 }
331
332 // fill digits
333 if (GetDebug()>2) cerr<<"AliMUONDigitzerv1::Exex TransientDigit to Digit"<<endl;
334 if ( digits[2] == icat ) pMUON->AddDigits(ich,tracks,charges,digits);
8dbbc4e3 335 }
391cd891 336 fManager->GetTreeD()->Fill();
337 pMUON->ResetDigits(); //
338 } // end loop cathode
339 fTDList->Delete();
8dbbc4e3 340
341 for(Int_t ii = 0; ii < 2*AliMUONConstants::NCh(); ++ii) {
342 if (fHitMap[ii]) {
343 delete fHitMap[ii];
344 fHitMap[ii] = 0;
345 }
346 }
347
348 if (GetDebug()>2)
349 cerr<<"AliMUONDigitizer::Exec: writing the TreeD: "
350 <<fManager->GetTreeD()->GetName()<<endl;
351 fManager->GetTreeD()->Write(0,TObject::kOverwrite);
352 delete [] fHitMap;
353 delete fTDList;
354
355 if (fHits) fHits->Clear();
356}
357
358//------------------------------------------------------------------------
359void AliMUONDigitizerv1::SortTracks(Int_t *tracks,Int_t *charges,Int_t ntr)
360{
361 //
362 // Sort the list of tracks contributing to a given digit
363 // Only the 3 most significant tracks are acctually sorted
364 //
365
366 //
367 // Loop over signals, only 3 times
368 //
369
370 Int_t qmax;
371 Int_t jmax;
372 Int_t idx[3] = {-2,-2,-2};
373 Int_t jch[3] = {-2,-2,-2};
374 Int_t jtr[3] = {-2,-2,-2};
375 Int_t i,j,imax;
376
377 if (ntr<3) imax=ntr;
378 else imax=3;
379 for(i=0;i<imax;i++){
380 qmax=0;
381 jmax=0;
382
383 for(j=0;j<ntr;j++){
384
385 if((i == 1 && j == idx[i-1])
386 ||(i == 2 && (j == idx[i-1] || j == idx[i-2]))) continue;
387
388 if(charges[j] > qmax) {
389 qmax = charges[j];
390 jmax=j;
391 }
392 }
393
394 if(qmax > 0) {
395 idx[i]=jmax;
396 jch[i]=charges[jmax];
397 jtr[i]=tracks[jmax];
398 }
399
400 }
401
402 for(i=0;i<3;i++){
403 if (jtr[i] == -2) {
404 charges[i]=0;
405 tracks[i]=0;
406 } else {
407 charges[i]=jch[i];
408 tracks[i]=jtr[i];
409 }
410 }
411}