]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - EMCAL/AliEMCALDigit.cxx
update from Marco on v2 task
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALDigit.cxx
... / ...
CommitLineData
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/* $Id$ */
17
18//_________________________________________________________________________
19// EMCAL digit:
20// A Digit is the sum of the energy lost in an EMCAL Tower
21// It also stores information on Primary, and enterring particle
22// tracknumbers Digits are created using AliEMCALSDigitizer, followed
23// by AliEMCALDigitizer
24//
25//*-- Author: Sahal Yacoob (LBL)
26// based on : AliPHOSDigit
27//__________________________________________________________________________
28
29// --- ROOT system ---
30#include <Riostream.h>
31#include <TMath.h>
32
33// --- Standard library ---
34
35// --- AliRoot header files ---
36
37#include "AliEMCALDigit.h"
38#include "AliEMCALGeometry.h"
39#include "AliLog.h"
40
41using std::cout;
42using std::endl;
43ClassImp(AliEMCALDigit)
44
45//____________________________________________________________________________
46AliEMCALDigit::AliEMCALDigit() :
47AliDigitNew(),
48 fAmpFloat(0.),
49 fNSamples(0),
50 fSamples(0x0),
51 fNSamplesHG(0),
52 fSamplesHG(0x0),
53 fNprimary(0),
54 fNMaxPrimary(5),
55 fPrimary(0x0),
56 fDEPrimary(0x0),
57 fNiparent(0),
58 fNMaxiparent(5),
59 fIparent(0x0),
60 fDEParent(0x0),
61 fMaxIter(0),
62 fTime(0.),
63 fTimeR(0.),
64 fChi2(0.),
65 fNDF(0),
66 fDigitType(kUnknown),
67 fAmpCalib(-1)
68{
69 // default ctor
70
71 // Need to initialise for reading old files
72 fPrimary = new Int_t[fNMaxPrimary] ;
73 fDEPrimary = new Float_t[fNMaxPrimary] ;
74 fIparent = new Int_t[fNMaxiparent] ;
75 fDEParent = new Float_t[fNMaxiparent] ;
76 for ( Int_t i = 0; i < fNMaxPrimary ; i++) {
77 fPrimary[i] = -1 ;
78 fDEPrimary[i] = 0 ;
79 }
80
81 for ( Int_t i = 0; i < fNMaxiparent ; i++) {
82 fIparent[i] = -1 ;
83 fDEParent[i] = 0 ;
84 }
85}
86
87//____________________________________________________________________________
88AliEMCALDigit::AliEMCALDigit(Int_t primary, Int_t iparent, Int_t id, Float_t digEnergy, Float_t time, Int_t type, Int_t index, Float_t chi2, Int_t ndf, Float_t dE)
89 : AliDigitNew(),
90 fAmpFloat(digEnergy),
91 fNSamples(0),
92 fSamples(0x0),
93 fNSamplesHG(0),
94 fSamplesHG(0x0),
95 fNprimary(0),
96 fNMaxPrimary(25),
97 fPrimary(0x0),
98 fDEPrimary(0x0),
99 fNiparent(0),
100 fNMaxiparent(150),
101 fIparent(0x0),
102 fDEParent(0x0),
103 fMaxIter(5),
104 fTime(time),
105 fTimeR(time),
106 fChi2(chi2),
107 fNDF(ndf),
108 fDigitType(type),
109 fAmpCalib(-1)
110{
111 // ctor with all data
112
113 // data memebrs of the base class (AliNewDigit)
114 fAmp = 0;
115 fId = id ;
116 fIndexInList = index ;
117
118 // data member
119 fPrimary = new Int_t[fNMaxPrimary] ;
120 fDEPrimary = new Float_t[fNMaxPrimary] ;
121 fIparent = new Int_t[fNMaxiparent] ;
122 fDEParent = new Float_t[fNMaxiparent] ;
123 if( primary != -1){
124 fNprimary = 1 ;
125 fPrimary[0] = primary ;
126 fDEPrimary[0] = dE ;
127 fNiparent = 1 ;
128 fIparent[0] = iparent ;
129 fDEParent[0] = dE ;
130 }
131 else{ //If the contribution of this primary smaller than fDigitThreshold (AliEMCALv1)
132 fNprimary = 0 ;
133 fPrimary[0] = -1 ;
134 fDEPrimary[0] = 0 ;
135 fNiparent = 0 ;
136 fIparent[0] = -1 ;
137 fDEParent[0] = 0 ;
138 }
139 Int_t i ;
140 for ( i = 1; i < fNMaxPrimary ; i++) {
141 fPrimary[i] = -1 ;
142 fDEPrimary[i] = 0 ;
143 }
144
145 for ( i = 1; i< fNMaxiparent ; i++) {
146 fIparent[i] = -1 ;
147 fDEParent[i] = 0 ;
148 }
149}
150
151//____________________________________________________________________________
152AliEMCALDigit::AliEMCALDigit(const AliEMCALDigit & digit)
153 : AliDigitNew(digit),
154 fAmpFloat(digit.fAmpFloat),
155 fNSamples(digit.fNSamples),
156 fSamples(),
157 fNSamplesHG(digit.fNSamplesHG),
158 fSamplesHG(),
159 fNprimary(digit.fNprimary),
160 fNMaxPrimary(digit.fNMaxPrimary),
161 fPrimary(),
162 fDEPrimary(),
163 fNiparent(digit.fNiparent),
164 fNMaxiparent(digit.fNMaxiparent),
165 fIparent(),
166 fDEParent(),
167 fMaxIter(digit.fMaxIter),
168 fTime(digit.fTime),
169 fTimeR(digit.fTimeR),
170 fChi2(digit.fChi2),
171 fNDF(digit.fNDF),
172 fDigitType(digit.fDigitType),
173 fAmpCalib(digit.fAmpCalib)
174{
175 // copy ctor
176 // data memebrs of the base class (AliNewDigit)
177 fAmp = digit.fAmp ;
178 fId = digit.fId;
179 fIndexInList = digit.fIndexInList ;
180
181 // data members
182 if (fSamples ) delete [] fSamples ; fSamples = NULL ;
183 if (fSamplesHG ) delete [] fSamplesHG ; fSamplesHG = NULL ;
184 if (fPrimary ) delete [] fPrimary ; fPrimary = NULL ;
185 if (fDEPrimary) delete [] fDEPrimary ; fDEPrimary = NULL ;
186 if (fIparent ) delete [] fIparent ; fIparent = NULL ;
187 if (fDEParent) delete [] fDEParent ; fDEParent = NULL ;
188
189 if (fNSamples){
190 fSamples = new Int_t[fNSamples];
191 for (Int_t i=0; i < digit.fNSamples; i++) fSamples[i] = digit.fSamples[i];
192 }
193
194 if (fNSamplesHG){
195 fSamplesHG = new Int_t[fNSamplesHG];
196 for (Int_t i=0; i < digit.fNSamplesHG; i++) fSamplesHG[i] = digit.fSamplesHG[i];
197 }
198
199
200 if (fNMaxPrimary){
201 fPrimary = new Int_t [fNMaxPrimary] ;
202 fDEPrimary = new Float_t[fNMaxPrimary] ;
203 for ( Int_t i = 0; i < fNMaxPrimary ; i++) {
204 fPrimary[i] = digit.fPrimary[i] ;
205 fDEPrimary[i] = digit.fDEPrimary[i] ;
206 }
207 }
208
209 if (fNMaxiparent){
210 fIparent = new Int_t [fNMaxiparent] ;
211 fDEParent = new Float_t[fNMaxiparent] ;
212 for (Int_t j = 0; j< fNMaxiparent ; j++) {
213 fIparent[j] = digit.fIparent[j] ;
214 fDEParent[j] = digit.fDEParent[j] ;
215 }
216 }
217
218}
219
220//____________________________________________________________________________
221AliEMCALDigit::~AliEMCALDigit()
222{
223 // Delete array of primaries if any
224
225 if (fSamples ) delete [] fSamples ; fSamples = NULL ;
226 if (fSamplesHG ) delete [] fSamplesHG ; fSamplesHG = NULL ;
227 if (fPrimary ) delete [] fPrimary ; fPrimary = NULL ;
228 if (fDEPrimary) delete [] fDEPrimary ; fDEPrimary = NULL ;
229 if (fIparent ) delete [] fIparent ; fIparent = NULL ;
230 if (fDEParent) delete [] fDEParent ; fDEParent = NULL ;
231
232}
233
234
235//____________________________________________________________________________
236void AliEMCALDigit::Clear(const Option_t*)
237{
238 // Delete array of primaries if any
239 if (fSamples ) delete [] fSamples ; fSamples = NULL ;
240 if (fSamplesHG ) delete [] fSamplesHG ; fSamplesHG = NULL ;
241 if (fPrimary ) delete [] fPrimary ; fPrimary = NULL ;
242 if (fDEPrimary) delete [] fDEPrimary ; fDEPrimary = NULL ;
243 if (fIparent ) delete [] fIparent ; fIparent = NULL ;
244 if (fDEParent) delete [] fDEParent ; fDEParent = NULL ;
245
246}
247
248//____________________________________________________________________________
249Int_t AliEMCALDigit::Compare(const TObject * obj) const
250{
251 // Compares two digits with respect to its Id
252 // to sort according increasing Id
253
254 Int_t rv = 2 ;
255
256 AliEMCALDigit * digit = (AliEMCALDigit *)obj ;
257
258 Int_t iddiff = fId - digit->GetId() ;
259
260 if ( iddiff > 0 )
261 rv = 1 ;
262 else if ( iddiff < 0 )
263 rv = -1 ;
264 else
265 rv = 0 ;
266
267 return rv ;
268
269}
270
271//____________________________________________________________________________
272Float_t AliEMCALDigit::GetEta() const
273{
274 //return pseudorapidity for this digit
275 // should be change in EMCALGeometry - 19-nov-04
276 Float_t eta=-10., phi=-10.;
277 Int_t id = GetId();
278 const AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
279 g->EtaPhiFromIndex(id,eta,phi);
280 return eta ;
281}
282
283//____________________________________________________________________________
284Float_t AliEMCALDigit::GetPhi() const
285{
286 //return phi coordinate of digit
287 // should be change in EMCALGeometry - 19-nov-04
288 Float_t eta=-10., phi=-10.;
289 Int_t id = GetId();
290 const AliEMCALGeometry *g = AliEMCALGeometry::GetInstance();
291 g->EtaPhiFromIndex(id,eta,phi);
292 return phi ;
293}
294
295//____________________________________________________________________________
296Bool_t AliEMCALDigit::GetFALTROSample(const Int_t iSample, Int_t& timeBin, Int_t& amp) const
297{
298 //Get FALTRO sample in time bin iSample
299 if (iSample >= fNSamples || iSample < 0 || fDigitType==kTrigger) return kFALSE;
300
301 amp = fSamples[iSample] & 0xFFF;
302 timeBin = (fSamples[iSample] >> 12) & 0xFF;
303
304 return kTRUE;
305}
306//____________________________________________________________________________
307Bool_t AliEMCALDigit::GetALTROSampleLG(const Int_t iSample, Int_t& timeBin, Int_t& amp) const
308{
309 //Get Low Gain ALTRO sample in time bin iSample
310 if (iSample >= fNSamples || iSample < 0 || fDigitType==kLG) return kFALSE;
311
312 amp = fSamples[iSample] & 0xFFF;
313 timeBin = (fSamples[iSample] >> 12) & 0xFF;
314
315 return kTRUE;
316}
317
318//____________________________________________________________________________
319void AliEMCALDigit::SetALTROSamplesLG(const Int_t nSamples, Int_t *samples)
320{
321 //Set array of ALTRO samples, Low Gain or FALTRO
322 fNSamples = nSamples;
323 fSamples = new Int_t[fNSamples];
324 for (Int_t i=0; i < fNSamples; i++) fSamples[i] = samples[i];
325}
326
327//____________________________________________________________________________
328void AliEMCALDigit::SetALTROSamplesHG(const Int_t nSamples, Int_t *samples)
329{
330 //Set array of ALTRO samples, High Gain.
331 fNSamplesHG = nSamples;
332 fSamplesHG = new Int_t[fNSamplesHG];
333 for (Int_t i=0; i < fNSamplesHG; i++) fSamplesHG[i] = samples[i];
334}
335
336
337//____________________________________________________________________________
338Bool_t AliEMCALDigit::GetALTROSampleHG(const Int_t iSample, Int_t& timeBin, Int_t& amp) const
339{
340 //Get High Gain ALTRO sample in time bin iSample
341 if (iSample >= fNSamplesHG || iSample < 0 || fDigitType==kHG) return kFALSE;
342
343 amp = fSamplesHG[iSample] & 0xFFF;
344 timeBin = (fSamplesHG[iSample] >> 12) & 0xFF;
345
346 return kTRUE;
347}
348
349
350//____________________________________________________________________________
351Int_t AliEMCALDigit::GetPrimary(Int_t index) const
352{
353 // retrieves the primary particle number given its index in the list
354 if ( (index <= fNprimary) && (index > 0)){
355 return fPrimary[index-1] ;
356 }
357
358 return -1 ;
359}
360
361//____________________________________________________________________________
362Float_t AliEMCALDigit::GetDEPrimary(Int_t index) const
363{
364 // retrieves the primary particle energy contribution
365 // given its index in the list
366 if ( (index <= fNprimary) && (index > 0)){
367 return fDEPrimary[index-1] ;
368 }
369
370 return 0 ;
371
372}
373
374//____________________________________________________________________________
375Int_t AliEMCALDigit::GetIparent(Int_t index) const
376{
377 // retrieves the primary particle number given its index in the list
378 if ( index <= fNiparent && index > 0){
379 return fIparent[index-1] ;
380 }
381
382 return -1 ;
383
384}
385
386//____________________________________________________________________________
387Float_t AliEMCALDigit::GetDEParent(Int_t index) const
388{
389 // retrieves the parent particle energy contribution
390 // given its index in the list
391 if ( (index <= fNiparent) && (index > 0)){
392 return fDEParent[index-1] ;
393 }
394
395 return 0;
396}
397
398//____________________________________________________________________________
399void AliEMCALDigit::ShiftPrimary(Int_t shift){
400 //shifts primary number to BIG offset, to separate primary in different TreeK
401 Int_t index ;
402 for(index = 0; index <fNprimary; index++ ){
403 fPrimary[index] = fPrimary[index]+ shift * 10000000 ;}
404 for(index =0; index <fNiparent; index++){
405 fIparent[index] = fIparent[index] + shift * 10000000 ;}
406}
407
408
409//____________________________________________________________________________
410AliEMCALDigit& AliEMCALDigit::operator= (const AliEMCALDigit &digit)
411{
412 // assignment operator
413
414 if(&digit == this) return *this;
415
416 fAmpFloat = digit.fAmpFloat;
417 fNSamples = digit.fNSamples;
418 fNSamplesHG = digit.fNSamplesHG;
419 fNprimary = digit.fNprimary;
420 fNMaxPrimary = digit.fNMaxPrimary;
421 fNiparent = digit.fNiparent;
422 fNMaxiparent = digit.fNMaxiparent;
423 fMaxIter = digit.fMaxIter;
424 fTime = digit.fTime;
425 fTimeR = digit.fTimeR;
426 fChi2 = digit.fChi2;
427 fNDF = digit.fNDF;
428 fDigitType = digit.fDigitType;
429 fAmpCalib = digit.fAmpCalib;
430 fAmp = digit.fAmp ;
431 fId = digit.fId;
432 fIndexInList = digit.fIndexInList ;
433
434 // data members
435 if (fSamples ) delete [] fSamples ; fSamples = NULL ;
436 if (fSamplesHG ) delete [] fSamplesHG ; fSamplesHG = NULL ;
437 if (fPrimary ) delete [] fPrimary ; fPrimary = NULL ;
438 if (fDEPrimary) delete [] fDEPrimary ; fDEPrimary = NULL ;
439 if (fIparent ) delete [] fIparent ; fIparent = NULL ;
440 if (fDEParent) delete [] fDEParent ; fDEParent = NULL ;
441
442 if (fNSamples){
443 fSamples = new Int_t[fNSamples];
444 for (Int_t i=0; i < digit.fNSamples; i++) fSamples[i] = digit.fSamples[i];
445 }
446
447 if (fNSamplesHG){
448 fSamplesHG = new Int_t[fNSamplesHG];
449 for (Int_t i=0; i < digit.fNSamplesHG; i++) fSamplesHG[i] = digit.fSamplesHG[i];
450 }
451
452
453 if (fNMaxPrimary){
454 fPrimary = new Int_t [fNMaxPrimary] ;
455 fDEPrimary = new Float_t[fNMaxPrimary] ;
456 for ( Int_t i = 0; i < fNMaxPrimary ; i++) {
457 fPrimary[i] = digit.fPrimary[i] ;
458 fDEPrimary[i] = digit.fDEPrimary[i] ;
459 }
460 }
461
462 if (fNMaxiparent){
463 fIparent = new Int_t [fNMaxiparent] ;
464 fDEParent = new Float_t[fNMaxiparent] ;
465 for (Int_t j = 0; j< fNMaxiparent ; j++) {
466 fIparent[j] = digit.fIparent[j] ;
467 fDEParent[j] = digit.fDEParent[j] ;
468 }
469 }
470
471 return *this;
472
473}
474
475//____________________________________________________________________________
476Bool_t AliEMCALDigit::operator==(AliEMCALDigit const & digit) const
477{
478 // Two digits are equal if they have the same Id
479
480 if ( fId == digit.fId )
481 return kTRUE ;
482 else
483 return kFALSE ;
484}
485
486//____________________________________________________________________________
487AliEMCALDigit AliEMCALDigit::operator+(const AliEMCALDigit &digit)
488{
489 // Adds the amplitude of digits and completes the list of primary particles
490 // if amplitude is larger than
491
492 fAmpFloat += digit.fAmpFloat ;
493 for (Int_t i=0; i < fNSamples ; i++) fSamples[i] += digit.fSamples[i];
494 for (Int_t i=0; i < fNSamplesHG; i++) fSamplesHG[i] += digit.fSamplesHG[i];
495
496 fAmp += digit.fAmp ;
497 if(fTime > digit.fTime)
498 fTime = digit.fTime ;
499 if (digit.fTimeR < fTimeR)
500 fTimeR = digit.fTimeR ;
501
502 Int_t max1 = fNprimary ;
503 Int_t max2 = fNiparent ;
504 Int_t index ;
505 for (index = 0 ; index < digit.fNprimary ; index++){
506 Bool_t newPrim = kTRUE ;
507 Int_t old ;
508 for ( old = 0 ; (old < max1) && newPrim; old++) { //already have this primary?
509 if(fPrimary[old] == digit.fPrimary[index]) {
510 newPrim = kFALSE;
511 fDEPrimary[old] += digit.fDEPrimary[index];
512 }
513 }
514 if (newPrim) {
515 if(max1<fNMaxPrimary){
516 fPrimary[max1] = digit.fPrimary[index] ;
517 fDEPrimary[max1] = digit.fDEPrimary[index] ;
518 fNprimary++ ;
519 max1++;
520 }
521 if(fNprimary==fNMaxPrimary) {
522
523 TString mess = " NMaxPrimary = " ;
524 mess += fNMaxPrimary ;
525 mess += " is too small" ;
526 AliFatal(mess.Data()) ;
527
528 }
529 }
530 }
531
532 for (index = 0 ; index < digit.fNiparent ; index++){
533 Bool_t newParent = kTRUE ;
534 Int_t old ;
535 for ( old = 0 ; (old < max2) && newParent; old++) { //already have this primary?
536 if(fIparent[old] == digit.fIparent[index]) {
537 newParent = kFALSE;
538 fDEParent[old] += digit.fDEParent[index];
539 }
540 }
541 if(newParent){
542 if(max2<fNMaxiparent) {
543 fIparent[max2] = digit.fIparent[index] ;
544 fDEParent[max2] = digit.fDEParent[index] ;
545 fNiparent++ ;
546 max2++;
547 }
548 if(fNiparent==fNMaxiparent) {
549
550 TString mess = " NMaxiparent = " ;
551 mess += fNMaxiparent ;
552 mess += " is too small" ;
553 AliFatal(mess.Data()) ;
554
555 }
556 }
557 }
558
559 return *this ;
560}
561
562//____________________________________________________________________________
563AliEMCALDigit AliEMCALDigit::operator*(Float_t factor)
564{
565 // Multiplies the amplitude by a factor
566
567 //Float_t tempo = static_cast<Float_t>(fAmp) ;
568 //tempo *= factor ;
569 //fAmp = static_cast<Int_t>(TMath::Floor(tempo)) ;
570
571 fAmpFloat *= factor;
572 for (Int_t i=0; i < fNSamples ; i++) fSamples[i] = Int_t(factor*fSamples[i]);
573 for (Int_t i=0; i < fNSamplesHG; i++) fSamplesHG[i] = Int_t(factor*fSamplesHG[i]);
574
575 for(Int_t i=0; i < fNprimary; i++)
576 fDEPrimary[i] *= factor;
577 for(Int_t i=0; i < fNiparent; i++)
578 fDEParent[i] *= factor;
579
580 return *this ;
581}
582
583//____________________________________________________________________________
584ostream& operator << ( ostream& out , const AliEMCALDigit & digit)
585{
586 // Prints the data of the digit
587
588 out << "ID " << digit.fId << " Energy = " << digit.fAmp << " Time = " << digit.fTime << endl ;
589 for(Int_t i=0;i<digit.fNprimary;i++)
590 out << "Primary " << i+1 << " = " << digit.fPrimary[i]
591 << " : DE " << digit.fDEPrimary[i] << endl ;
592
593 for(Int_t j=0;j<digit.fNiparent;j++)
594 out << "Iparent " << j+1 << " = " << digit.fIparent[j]
595 << " : DE " << digit.fDEParent[j] << endl ;
596 out << "Position in list = " << digit.fIndexInList << endl ;
597 return out ;
598}
599
600//____________________________________________________________________________
601void AliEMCALDigit::Print(const Option_t* /*opt*/) const
602{
603 //Print
604 printf("===\nDigit id: %4d / Energy %2.3f ; Time %e ; Time samples %d ; Chi2 %2.3f, NDF %d, Type? %d \n",
605 fId, fAmpFloat,fTime, fNSamples, fChi2, fNDF, fDigitType);
606 if(fDigitType==kTrigger){
607 printf("FALTRO: ");
608 for (Int_t i=0; i < GetNFALTROSamples(); i++)
609 {
610 Int_t timeBin, amp;
611 GetFALTROSample(i, timeBin, amp);
612 printf(" (%d,%d) ",timeBin,amp);
613 }
614 printf("\n");
615 }//trigger
616 else{
617 printf("ALTRO, Low Gain: ");
618 for (Int_t i=0; i < GetNALTROSamplesLG(); i++)
619 {
620 Int_t timeBin, amp;
621 GetALTROSampleLG(i, timeBin, amp);
622 printf(" (%d,%d) ",timeBin,amp);
623 }
624 printf("\n");
625 printf("ALTRO, High Gain: ");
626 for (Int_t i=0; i < GetNALTROSamplesHG(); i++)
627 {
628 Int_t timeBin, amp;
629 GetALTROSampleHG(i, timeBin, amp);
630 printf(" (%d,%d) ",timeBin,amp);
631 }
632 printf("\n");
633 }//trigger
634
635}
636
637
638