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 | |
22 | ClassImp(AliMUONDigitizerv1) |
23 | |
24 | //___________________________________________ |
25 | AliMUONDigitizerv1::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 | //___________________________________________ |
36 | AliMUONDigitizerv1::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 | //------------------------------------------------------------------------ |
50 | AliMUONDigitizerv1::~AliMUONDigitizerv1() |
51 | { |
52 | // Destructor |
53 | delete fHits; |
54 | } |
55 | |
56 | //------------------------------------------------------------------------ |
57 | void 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 | //------------------------------------------------------------------------ |
67 | Bool_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 | //------------------------------------------------------------------------ |
75 | Bool_t AliMUONDigitizerv1::Init() |
76 | { |
77 | // Initialization |
78 | return kTRUE; |
79 | } |
80 | |
81 | //------------------------------------------------------------------------ |
82 | void 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 | //-------------------------------------------------------------------------- |
104 | void 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 | //----------------------------------------------------------------------- |
164 | void 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) { |
216 | fHits->Delete(); |
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; |
231 | fHits->Delete(); |
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 | |
274 | // |
275 | // Filling Digit List |
276 | Int_t tracks[kMAXTRACKS]; |
277 | Int_t charges[kMAXTRACKS]; |
278 | Int_t nentries = fTDList->GetEntriesFast(); |
279 | Int_t digits[6]; |
280 | for (Int_t nent = 0; nent < nentries; nent++) { |
281 | AliMUONTransientDigit *address = (AliMUONTransientDigit*)fTDList->At(nent); |
282 | if (address == 0) continue; |
283 | Int_t ich = address->Chamber(); |
284 | Int_t q = address->Signal(); |
285 | chamber = &(pMUON->Chamber(ich)); |
286 | // |
287 | // Digit Response (noise, threshold, saturation, ...) |
288 | AliMUONResponse * response = chamber->ResponseModel(); |
289 | q = response->DigitResponse(q,address); |
290 | |
291 | if (!q) continue; |
292 | |
293 | digits[0] = address->PadX(); |
294 | digits[1] = address->PadY(); |
295 | digits[2] = address->Cathode(); |
296 | digits[3] = q; |
297 | digits[4] = address->Physics(); |
298 | digits[5] = address->Hit(); |
299 | |
300 | Int_t nptracks = address->GetNTracks(); |
301 | |
302 | if (nptracks > kMAXTRACKS) { |
303 | if (GetDebug() >0) { |
304 | cerr<<"AliMUONDigitizer:Exec nptracks > 10 "<<nptracks; |
305 | cerr<<"reset to max value "<<kMAXTRACKS<<endl; |
306 | } |
307 | nptracks = kMAXTRACKS; |
308 | } |
309 | if (nptracks > 2 && GetDebug() >2) { |
310 | cerr<<"AliMUONDigitizer::Exec nptracks > 2 "<<nptracks<<endl; |
311 | // printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q); |
312 | } |
313 | for (Int_t tr = 0; tr < nptracks; tr++) { |
314 | tracks[tr] = address->GetTrack(tr); |
315 | charges[tr] = address->GetCharge(tr); |
316 | } //end loop over list of tracks for one pad |
317 | // Sort list of tracks according to charge |
318 | if (nptracks > 1) { |
319 | SortTracks(tracks,charges,nptracks); |
320 | } |
321 | if (nptracks < kMAXTRACKS ) { |
322 | for (Int_t i = nptracks; i < kMAXTRACKS; i++) { |
323 | tracks[i] = 0; |
324 | charges[i] = 0; |
325 | } |
326 | } |
327 | |
328 | // fill digits |
329 | if (GetDebug()>2) cerr<<"AliMUONDigitzerv1::Exex TransientDigit to Digit"<<endl; |
330 | pMUON->AddDigits(ich,tracks,charges,digits); |
331 | } |
332 | fManager->GetTreeD()->Fill(); |
333 | pMUON->ResetDigits(); // |
334 | fTDList->Delete(); |
335 | |
336 | for(Int_t ii = 0; ii < 2*AliMUONConstants::NCh(); ++ii) { |
337 | if (fHitMap[ii]) { |
338 | delete fHitMap[ii]; |
339 | fHitMap[ii] = 0; |
340 | } |
341 | } |
342 | |
343 | if (GetDebug()>2) |
344 | cerr<<"AliMUONDigitizer::Exec: writing the TreeD: " |
345 | <<fManager->GetTreeD()->GetName()<<endl; |
346 | fManager->GetTreeD()->Write(0,TObject::kOverwrite); |
347 | delete [] fHitMap; |
348 | delete fTDList; |
349 | |
350 | if (fHits) fHits->Clear(); |
351 | } |
352 | |
353 | //------------------------------------------------------------------------ |
354 | void AliMUONDigitizerv1::SortTracks(Int_t *tracks,Int_t *charges,Int_t ntr) |
355 | { |
356 | // |
357 | // Sort the list of tracks contributing to a given digit |
358 | // Only the 3 most significant tracks are acctually sorted |
359 | // |
360 | |
361 | // |
362 | // Loop over signals, only 3 times |
363 | // |
364 | |
365 | Int_t qmax; |
366 | Int_t jmax; |
367 | Int_t idx[3] = {-2,-2,-2}; |
368 | Int_t jch[3] = {-2,-2,-2}; |
369 | Int_t jtr[3] = {-2,-2,-2}; |
370 | Int_t i,j,imax; |
371 | |
372 | if (ntr<3) imax=ntr; |
373 | else imax=3; |
374 | for(i=0;i<imax;i++){ |
375 | qmax=0; |
376 | jmax=0; |
377 | |
378 | for(j=0;j<ntr;j++){ |
379 | |
380 | if((i == 1 && j == idx[i-1]) |
381 | ||(i == 2 && (j == idx[i-1] || j == idx[i-2]))) continue; |
382 | |
383 | if(charges[j] > qmax) { |
384 | qmax = charges[j]; |
385 | jmax=j; |
386 | } |
387 | } |
388 | |
389 | if(qmax > 0) { |
390 | idx[i]=jmax; |
391 | jch[i]=charges[jmax]; |
392 | jtr[i]=tracks[jmax]; |
393 | } |
394 | |
395 | } |
396 | |
397 | for(i=0;i<3;i++){ |
398 | if (jtr[i] == -2) { |
399 | charges[i]=0; |
400 | tracks[i]=0; |
401 | } else { |
402 | charges[i]=jch[i]; |
403 | tracks[i]=jtr[i]; |
404 | } |
405 | } |
406 | } |