Line data Source code
1 : /**************************************************************************
2 : * Copyright(c) 2007, 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 : //This class produces PHOS digits of one event
19 : //using AliPHOSRawFitter.
20 : //
21 : // For example:
22 : // TClonesArray *digits = new TClonesArray("AliPHOSDigit",100);
23 : // AliRawReader* rawReader = new AliRawReaderDate("2006run2211.raw");
24 : // AliPHOSRawDecoder dc(rawReader);
25 : // while (rawReader->NextEvent()) {
26 : // AliPHOSRawDigiProducer producer;
27 : // producer.MakeDigits(digits,&dc);
28 : // }
29 :
30 : // Author: Boris Polichtchouk
31 :
32 : // --- ROOT system ---
33 : #include "TClonesArray.h"
34 :
35 : // --- AliRoot header files ---
36 : #include "AliPHOSRawDigiProducer.h"
37 : #include "AliPHOSRawFitterv0.h"
38 : #include "AliPHOSGeometry.h"
39 : #include "AliPHOSDigit.h"
40 : #include "AliPHOSCalibData.h"
41 : #include "AliPHOSPulseGenerator.h"
42 : #include "AliCaloRawStreamV3.h"
43 : #include "AliDAQ.h"
44 : #include "AliRawReader.h"
45 : #include "AliLog.h"
46 :
47 22 : ClassImp(AliPHOSRawDigiProducer)
48 :
49 : AliPHOSCalibData * AliPHOSRawDigiProducer::fgCalibData = 0 ;
50 :
51 : //--------------------------------------------------------------------------------------
52 : AliPHOSRawDigiProducer::AliPHOSRawDigiProducer():
53 0 : TObject(),
54 0 : fSubtractL1phase(kTRUE),
55 0 : fEmcMinE(0.),
56 0 : fCpvMinE(0.),
57 0 : fSampleQualityCut(1.),
58 0 : fSampleToSec(0.),
59 0 : fEmcCrystals(0),
60 0 : fGeom(0),
61 0 : fPulseGenerator(0),
62 0 : fRawReader(0),
63 0 : fRawStream(0),
64 0 : fADCValuesLG(0),
65 0 : fADCValuesHG(0)
66 0 : {
67 : // Default constructor
68 :
69 0 : fGeom=AliPHOSGeometry::GetInstance() ;
70 0 : if(!fGeom) fGeom = AliPHOSGeometry::GetInstance("IHEP");
71 :
72 0 : fEmcCrystals=fGeom->GetNCristalsInModule()*fGeom->GetNModules() ;
73 0 : fPulseGenerator = new AliPHOSPulseGenerator();
74 0 : GetCalibrationParameters() ;
75 :
76 0 : }
77 : //--------------------------------------------------------------------------------------
78 : AliPHOSRawDigiProducer::AliPHOSRawDigiProducer(AliRawReader *rawReader,
79 : AliAltroMapping **mapping):
80 4 : TObject(),
81 4 : fSubtractL1phase(kTRUE),
82 4 : fEmcMinE(0.),
83 4 : fCpvMinE(0.),
84 4 : fSampleQualityCut(1.),
85 4 : fSampleToSec(0.),
86 4 : fEmcCrystals(0),
87 4 : fGeom(0),
88 4 : fPulseGenerator(0),
89 4 : fRawReader(rawReader),
90 4 : fRawStream(0),
91 4 : fADCValuesLG(0),
92 4 : fADCValuesHG(0)
93 20 : {
94 : // Default constructor
95 :
96 8 : fGeom=AliPHOSGeometry::GetInstance() ;
97 4 : if(!fGeom) fGeom = AliPHOSGeometry::GetInstance("IHEP");
98 :
99 8 : fEmcCrystals=fGeom->GetNCristalsInModule()*fGeom->GetNModules() ;
100 12 : fPulseGenerator = new AliPHOSPulseGenerator();
101 4 : GetCalibrationParameters() ;
102 :
103 16 : fRawStream = new AliCaloRawStreamV3(rawReader,"PHOS",mapping);
104 : // Select only data in ALTRO format and skip STU, the last PHOS DDL
105 8 : rawReader->Select("PHOS",0,AliDAQ::NumberOfDdls("PHOS")-2);
106 :
107 8 : }
108 : //--------------------------------------------------------------------------------------
109 : AliPHOSRawDigiProducer::AliPHOSRawDigiProducer(const AliPHOSRawDigiProducer &dp):
110 0 : TObject(),
111 0 : fSubtractL1phase(kTRUE),
112 0 : fEmcMinE(0.),
113 0 : fCpvMinE(0.),
114 0 : fSampleQualityCut(1.),
115 0 : fSampleToSec(0.),
116 0 : fEmcCrystals(0),
117 0 : fGeom(0),
118 0 : fPulseGenerator(0),
119 0 : fRawReader(0),
120 0 : fRawStream(0),
121 0 : fADCValuesLG(0),
122 0 : fADCValuesHG(0)
123 :
124 0 : {
125 : // Copy constructor
126 :
127 0 : fEmcMinE = dp.fEmcMinE ;
128 0 : fCpvMinE = dp.fCpvMinE ;
129 0 : fSampleQualityCut = dp.fSampleQualityCut;
130 0 : fSampleToSec = dp.fSampleToSec ;
131 0 : fEmcCrystals = dp.fEmcCrystals ;
132 0 : fPulseGenerator = new AliPHOSPulseGenerator();
133 0 : fGeom = dp.fGeom ;
134 0 : }
135 : //--------------------------------------------------------------------------------------
136 : AliPHOSRawDigiProducer& AliPHOSRawDigiProducer::operator= (const AliPHOSRawDigiProducer &dp)
137 : {
138 : // Assign operator
139 :
140 0 : if(&dp == this) return *this;
141 :
142 0 : fSubtractL1phase = dp.fSubtractL1phase ;
143 0 : fEmcMinE = dp.fEmcMinE ;
144 0 : fCpvMinE = dp.fCpvMinE ;
145 0 : fSampleQualityCut = dp.fSampleQualityCut ;
146 0 : fSampleToSec = dp.fSampleToSec ;
147 0 : fEmcCrystals = dp.fEmcCrystals ;
148 0 : fGeom = dp.fGeom ;
149 0 : if(fPulseGenerator) delete fPulseGenerator ;
150 0 : fPulseGenerator = new AliPHOSPulseGenerator();
151 0 : return *this;
152 0 : }
153 : //--------------------------------------------------------------------------------------
154 : AliPHOSRawDigiProducer::~AliPHOSRawDigiProducer()
155 16 : {
156 : // Destructor
157 12 : if(fPulseGenerator) delete fPulseGenerator ;
158 4 : fPulseGenerator=0 ;
159 8 : delete fRawStream;
160 8 : delete [] fADCValuesLG;
161 8 : delete [] fADCValuesHG;
162 8 : }
163 : //--------------------------------------------------------------------------------------
164 : void AliPHOSRawDigiProducer::MakeDigits(TClonesArray *digits, AliPHOSRawFitterv0* fitter)
165 : {
166 : // Create a temporary array of LG digits and then make digits from raw data
167 :
168 0 : TClonesArray *tmpLG = new TClonesArray("AliPHOSDigit",10000) ;
169 0 : MakeDigits(digits, tmpLG, fitter);
170 0 : tmpLG->Delete();
171 0 : delete tmpLG;
172 0 : }
173 : //--------------------------------------------------------------------------------------
174 : void AliPHOSRawDigiProducer::MakeDigits(TClonesArray *digits, TClonesArray *tmpDigLG, AliPHOSRawFitterv0* fitter)
175 : {
176 : //Makes the job.
177 : //TClonesArray *digits, *tmpDigLG and raw data fitter should be provided by calling function.
178 : Int_t iDigit=0 ;
179 8 : Int_t relId[4], absId=-1, caloFlag=-1;
180 :
181 : const Double_t baseLine=1. ; //Minimal energy of digit in ADC ch.
182 :
183 : //Calculate conversion coeff. from Sample time step to seconds
184 : //If OCDB contains negative or zero value - use one from RCU trailer
185 : //Negative value in OCDB is used only for simulation of raw digits
186 4 : if(fgCalibData->GetSampleTimeStep()>0.)
187 4 : fSampleToSec=fgCalibData->GetSampleTimeStep() ;
188 : else
189 0 : fSampleToSec=fRawStream->GetTSample() ;
190 :
191 : // Clear a temporary array for LowGain digits
192 4 : tmpDigLG->Clear();
193 : Int_t ilgDigit=0 ;
194 :
195 : //Let fitter subtract pedestals in case of ZS
196 4 : fitter->SetCalibData(fgCalibData) ;
197 :
198 62 : while (fRawStream->NextDDL()) {
199 : // Skip STU DDL
200 27 : if (fRawStream->GetDDLNumber() == fgkSTUDDL) continue;
201 181 : while (fRawStream->NextChannel()) {
202 154 : relId[0] = 5 - fRawStream->GetModule() ; // counts from 1 to 5
203 154 : relId[1] = 0; // 0=EMC
204 154 : relId[2] = fRawStream->GetCellX() + 1; // counts from 1 to 64
205 154 : relId[3] = fRawStream->GetCellZ() + 1; // counts from 1 to 56
206 154 : caloFlag = fRawStream->GetCaloFlag(); // 0=LG, 1=HG, 2=TRU
207 :
208 154 : if(caloFlag!=0 && caloFlag!=1) continue; //TRU data!
209 :
210 154 : fitter->SetChannelGeo(relId[0],relId[2],relId[3],caloFlag);
211 :
212 308 : if(fitter->GetAmpOffset()==0 && fitter->GetAmpThreshold()==0) {
213 154 : short value = fRawStream->GetAltroCFG1();
214 154 : bool ZeroSuppressionEnabled = (value >> 15) & 0x1;
215 154 : if(ZeroSuppressionEnabled) {
216 0 : short offset = (value >> 10) & 0xf;
217 0 : short threshold = value & 0x3ff;
218 0 : fitter->SubtractPedestals(kFALSE);
219 0 : fitter->SetAmpOffset(offset);
220 0 : fitter->SetAmpThreshold(threshold);
221 0 : }
222 154 : }
223 :
224 154 : fGeom->RelToAbsNumbering(relId, absId);
225 :
226 154 : fitter->SetNBunches(0);
227 : Int_t sigStart =0 ;
228 : Int_t sigLength=0 ;
229 616 : while (fRawStream->NextBunch()) { //Take the first in time bunch
230 154 : const UShort_t *sig = fRawStream->GetSignals();
231 154 : sigStart = fRawStream->GetStartTimeBin();
232 154 : sigLength = fRawStream->GetBunchLength();
233 154 : fitter->Eval(sig,sigStart,sigLength);
234 154 : if (caloFlag == AliCaloRawStreamV3::kLowGain) {
235 88 : delete [] fADCValuesLG;
236 46 : fADCValuesLG = new Int_t[sigLength];
237 3204 : for (Int_t i=0; i<sigLength; i++)
238 1556 : fADCValuesLG[sigLength-i-1] = sig[i];
239 46 : }
240 108 : else if (caloFlag == AliCaloRawStreamV3::kHighGain) {
241 212 : delete [] fADCValuesHG;
242 108 : fADCValuesHG = new Int_t[sigLength];
243 14134 : for (Int_t i=0; i<sigLength; i++)
244 6959 : fADCValuesHG[sigLength-i-1] = sig[i];
245 108 : }
246 : } // End of NextBunch()
247 :
248 :
249 154 : Double_t energy = fitter->GetEnergy() ;
250 154 : Double_t time = fitter->GetTime() ;
251 154 : if(energy<=baseLine) //in ADC channels
252 27 : continue ;
253 :
254 : //remove digits with bad shape. Fitter should calculate quality so that
255 : //in default case quality [0,1], while larger values of quality mean somehow
256 : //corrupted samples, 999 means obviously corrupted sample.
257 : //It is difficult to fit samples with overflow (even setting cut on overflow values)
258 : //because too few points are left to fit. So we do not evaluate samples with overflow
259 :
260 127 : if(fitter->GetSignalQuality() > fSampleQualityCut && !(fitter->IsOverflow()))
261 0 : continue ;
262 :
263 127 : energy = CalibrateE(energy,relId,!caloFlag) ;
264 :
265 : //convert time from sample bin units to s
266 127 : time*=fSampleToSec ;
267 : //CalibrateT moved to Clusterizer
268 : // time = CalibrateT(time,relId,!caloFlag) ;
269 : // subtract RCU L1 phase (L1Phase is in seconds) w.r.t. L0:
270 : //Very strange behaviour of electronics, but cross-checkes several times...
271 127 : if(fSubtractL1phase){
272 127 : if( fRawStream->GetL1Phase()<55.*1.e-9 ) //for phase=0,25,50
273 127 : time -= fRawStream->GetL1Phase();
274 : else //for phase 75
275 0 : time += 25.*1.e-9 ;
276 : }
277 :
278 127 : if(energy <= 0.)
279 0 : continue;
280 :
281 127 : if (caloFlag == AliCaloRawStreamV3::kLowGain) {
282 19 : new((*tmpDigLG)[ilgDigit]) AliPHOSDigit(-1,absId,(Float_t)energy,(Float_t)time);
283 38 : if (sigLength>0 && fADCValuesLG!=0)
284 19 : static_cast<AliPHOSDigit*>(tmpDigLG->At(ilgDigit))->SetALTROSamplesLG(sigLength,fADCValuesLG);
285 19 : ilgDigit++ ;
286 19 : }
287 108 : else if (caloFlag == AliCaloRawStreamV3::kHighGain) {
288 216 : if(fitter->IsOverflow()) //Keep this digit to replace it by Low Gain later.
289 : //If there is no LogGain it wil be removed by cut on Min E
290 110 : new((*digits)[iDigit]) AliPHOSDigit(-1,absId,-1.f,(Float_t)time);
291 : else
292 106 : new((*digits)[iDigit]) AliPHOSDigit(-1,absId,(Float_t)energy,(Float_t)time);
293 216 : if (sigLength>0 && fADCValuesHG!=0)
294 108 : static_cast<AliPHOSDigit*>(digits->At(iDigit))->SetALTROSamplesHG(sigLength,fADCValuesHG);
295 108 : iDigit++;
296 108 : }
297 127 : } // End of NextChannel()
298 :
299 : //Now scan created LG and HG digits and keep only those which appeared in both lists
300 : //replace energy of HighGain digits only if there is overflow
301 : //negative energy (overflow)
302 27 : digits->Sort() ;
303 27 : tmpDigLG->Sort() ;
304 : Int_t iLG = 0;
305 27 : Int_t nLG1 = tmpDigLG->GetEntriesFast()-1 ;
306 :
307 1284 : for(Int_t iDig=0 ; iDig < digits->GetEntriesFast() ; iDig++) {
308 1845 : AliPHOSDigit * digHG = dynamic_cast<AliPHOSDigit*>(digits->At(iDig)) ;
309 615 : if (!digHG) continue;
310 1835 : AliPHOSDigit * digLG = dynamic_cast<AliPHOSDigit*>(tmpDigLG->At(iLG)) ;
311 2649 : while(digLG && iLG<nLG1 && digHG->GetId()> digLG->GetId()){
312 90 : iLG++ ;
313 270 : digLG = dynamic_cast<AliPHOSDigit*>(tmpDigLG->At(iLG)) ;
314 : }
315 615 : absId=digHG->GetId() ;
316 615 : fGeom->AbsToRelNumbering(absId,relId) ;
317 :
318 1220 : if(digLG && digHG->GetId() == digLG->GetId()){ //we found pair
319 113 : if(digHG->GetEnergy()<0.){ //This is overflow in HG
320 2 : digHG->SetTime(digLG->GetTime()) ;
321 2 : digHG->SetEnergy(digLG->GetEnergy()) ;
322 2 : digHG->SetLG(kTRUE) ;
323 2 : }
324 : }
325 : else{ //no pair - remove
326 502 : if(digHG->GetEnergy()<0.) //no pair, in saturation
327 0 : digits->RemoveAt(iDig) ;
328 : }
329 615 : }
330 : } // End of NextDDL()
331 :
332 4 : CleanDigits(digits) ;
333 :
334 4 : }
335 : //____________________________________________________________________________
336 : Double_t AliPHOSRawDigiProducer::CalibrateE(Double_t amp, Int_t* relId, Bool_t isLowGain)
337 : {
338 : // Convert EMC LG amplitude to HG (multipli by ~16)
339 : // Calibration parameters are taken from calibration data base
340 127 : if(fgCalibData){
341 127 : Int_t module = relId[0];
342 127 : Int_t column = relId[3];
343 127 : Int_t row = relId[2];
344 127 : if(relId[1]==0) { // this is EMC
345 127 : if(isLowGain){
346 19 : amp*= fgCalibData->GetHighLowRatioEmc(module,column,row);
347 19 : }
348 127 : return amp ;
349 : }
350 0 : }
351 0 : return 0;
352 127 : }
353 : //____________________________________________________________________________
354 : Double_t AliPHOSRawDigiProducer::CalibrateT(Double_t time, Int_t * relId, Bool_t /* isLowGain */)
355 : {
356 : //Calibrate time
357 0 : if(fgCalibData){
358 0 : Int_t module = relId[0];
359 0 : Int_t column = relId[3];
360 0 : Int_t row = relId[2];
361 0 : if(relId[1]==0) { // this is EMC
362 0 : time += fgCalibData->GetTimeShiftEmc(module,column,row);
363 0 : return time ;
364 : }
365 0 : }
366 :
367 0 : return -999.;
368 0 : }
369 : //____________________________________________________________________________
370 : void AliPHOSRawDigiProducer::CleanDigits(TClonesArray * digits)
371 : {
372 : // remove digits with amplitudes below threshold.
373 : // remove digits in bad channels
374 : // sort digits with icreasing AbsId
375 :
376 : //remove digits in bad map and below threshold
377 : Bool_t isBadMap = 0 ;
378 8 : if(fgCalibData->GetNumOfEmcBadChannels()){
379 : isBadMap=1 ;
380 0 : }
381 :
382 224 : for(Int_t i=0; i<digits->GetEntriesFast(); i++){
383 108 : AliPHOSDigit * digit = static_cast<AliPHOSDigit*>(digits->At(i)) ;
384 108 : if(!digit)
385 0 : continue ;
386 216 : if ( (IsInEMC(digit) && digit->GetEnergy() < fEmcMinE) ||
387 108 : (IsInCPV(digit) && digit->GetEnergy() < fCpvMinE) ){
388 0 : digits->RemoveAt(i) ;
389 0 : continue ;
390 : }
391 108 : if(isBadMap){ //check bad map now
392 0 : Int_t relid[4] ;
393 0 : fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
394 0 : if(fgCalibData->IsBadChannelEmc(relid[0],relid[3],relid[2])){
395 0 : digits->RemoveAt(i) ;
396 0 : }
397 0 : }
398 108 : }
399 :
400 : //Compress, sort and set indexes
401 4 : digits->Compress() ;
402 : // digits->Sort(); already sorted earlier
403 224 : for (Int_t i = 0 ; i < digits->GetEntriesFast() ; i++) {
404 108 : AliPHOSDigit *digit = static_cast<AliPHOSDigit*>( digits->At(i) ) ;
405 108 : digit->SetIndexInList(i) ;
406 : }
407 4 : }
408 : //____________________________________________________________________________
409 : Bool_t AliPHOSRawDigiProducer::IsInEMC(AliPHOSDigit * digit) const
410 : {
411 : // Tells if (true) or not (false) the digit is in a PHOS-EMC module
412 216 : return digit->GetId() <= fEmcCrystals ;
413 :
414 : }
415 :
416 : //____________________________________________________________________________
417 : Bool_t AliPHOSRawDigiProducer::IsInCPV(AliPHOSDigit * digit) const
418 : {
419 : // Tells if (true) or not (false) the digit is in a PHOS-CPV module
420 216 : return digit->GetId() > fEmcCrystals ;
421 : }
422 : //____________________________________________________________________________
423 : void AliPHOSRawDigiProducer::GetCalibrationParameters()
424 : {
425 : // Set calibration parameters:
426 : // if calibration database exists, they are read from database,
427 : // otherwise, reconstruction stops in the constructor of AliPHOSCalibData
428 : //
429 : // It is a user responsilibity to open CDB before reconstruction, for example:
430 : // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");
431 :
432 8 : if (!fgCalibData){
433 2 : fgCalibData = new AliPHOSCalibData(-1); //use AliCDBManager's run number
434 1 : }
435 4 : if (fgCalibData->GetCalibDataEmc() == 0)
436 0 : AliFatal("Calibration parameters for PHOS EMC not found. Stop reconstruction.\n");
437 4 : if (fgCalibData->GetCalibDataCpv() == 0)
438 0 : AliFatal("Calibration parameters for PHOS CPV not found. Stop reconstruction.\n");
439 4 : }
|
