1 /**************************************************************************
2 * Copyright(c) 1998-2000, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 Class for creating of the sumable digits and digits from MC data
21 The input : ideal signals (Hits->Diffusion->Attachment -Ideal signal)
22 The output: raw digits
25 1. Pad by pad gain map
27 3. The dead channels identified - zerro noise for corresponding pads
28 In this case the outpu equal zerro
37 #include <TObjArray.h>
39 #include <TDirectory.h>
40 #include <Riostream.h>
41 #include <TParameter.h>
43 #include "AliTPCDigitizer.h"
46 #include "AliTPCParam.h"
47 #include "AliTPCParamSR.h"
49 #include "AliLoader.h"
51 #include "AliDigitizationInput.h"
52 #include "AliSimDigits.h"
55 #include "AliTPCcalibDB.h"
56 #include "AliTPCCalPad.h"
57 #include "AliTPCCalROC.h"
62 ClassImp(AliTPCDigitizer)
64 //___________________________________________
65 AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer(),fDebug(0)
68 // Default ctor - don't use it
73 //___________________________________________
74 AliTPCDigitizer::AliTPCDigitizer(AliDigitizationInput* digInput)
75 :AliDigitizer(digInput),fDebug(0)
78 // ctor which should be used
80 AliDebug(2,"(AliDigitizationInput* digInput) was processed");
83 //------------------------------------------------------------------------
84 AliTPCDigitizer::~AliTPCDigitizer()
91 //------------------------------------------------------------------------
92 Bool_t AliTPCDigitizer::Init()
100 //------------------------------------------------------------------------
101 void AliTPCDigitizer::Digitize(Option_t* option)
103 DigitizeFast(option);
105 //------------------------------------------------------------------------
106 void AliTPCDigitizer::DigitizeFast(Option_t* option)
109 // merge input tree's with summable digits
110 //output stored in TreeTPCD
113 TString optionString = option;
114 if (!strcmp(optionString.Data(),"deb")) {
115 cout<<"AliTPCDigitizer:::DigitizeFast called with option deb "<<endl;
118 //get detector and geometry
121 AliRunLoader *rl, *orl;
122 AliLoader *gime, *ogime;
126 Warning("DigitizeFast","gAlice is NULL. Loading from input 0");
127 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
130 Error("DigitizeFast","Can not find Run Loader for input 0. Can not proceed.");
136 AliTPC *pTPC = (AliTPC *) gAlice->GetModule("TPC");
137 AliTPCParam * param = pTPC->GetParam();
139 //sprintf(s,param->GetTitle());
140 snprintf(s,100,"%s",param->GetTitle());
141 //sprintf(ss,"75x40_100x60");
142 snprintf(ss,100,"75x40_100x60");
144 printf("2 pad-length geom hits with 3 pad-lenght geom digits...\n");
146 param=new AliTPCParamSR();
149 //sprintf(ss,"75x40_100x60_150x60");
150 snprintf(ss,100,"75x40_100x60_150x60");
151 if(strcmp(s,ss)!=0) {
152 printf("No TPC parameters found...\n");
157 pTPC->GenerNoise(500000); //create table with noise
159 Int_t nInputs = fDigInput->GetNinputs();
160 Int_t * masks = new Int_t[nInputs];
161 for (Int_t i=0; i<nInputs;i++)
162 masks[i]= fDigInput->GetMask(i);
163 Short_t **pdig= new Short_t*[nInputs]; //pointers to the expanded digits array
164 Int_t **ptr= new Int_t*[nInputs]; //pointers to the expanded tracks array
165 Bool_t *active= new Bool_t[nInputs]; //flag for active input segments
168 //create digits array for given sectors
171 //create branch's in TPC treeD
172 orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
173 ogime = orl->GetLoader("TPCLoader");
174 TTree * tree = ogime->TreeD();
175 AliSimDigits * digrow = new AliSimDigits;
179 ogime->MakeTree("D");
180 tree = ogime->TreeD();
182 tree->Branch("Segment","AliSimDigits",&digrow);
184 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
185 for (Int_t i1=0;i1<nInputs; i1++)
189 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i1));
190 gime = rl->GetLoader("TPCLoader");
191 gime->LoadSDigits("read");
192 TTree * treear = gime->TreeS();
196 cerr<<"AliTPCDigitizer: Input tree with SDigits not found in"
197 <<" input "<< i1<<endl;
198 for (Int_t i2=0;i2<i1+1; i2++){
200 if(digarr[i2]) delete digarr[i2];
210 //sprintf(phname,"lhcphase%d",i1);
211 snprintf(phname,100,"lhcphase%d",i1);
212 TParameter<float> *ph = (TParameter<float>*)treear->GetUserInfo()
213 ->FindObject("lhcphase0");
215 cerr<<"AliTPCDigitizer: LHC phase not found in"
216 <<" input "<< i1<<endl;
217 for (Int_t i2=0;i2<i1+1; i2++){
218 if(digarr[i2]) delete digarr[i2];
227 tree->GetUserInfo()->Add(new TParameter<float>(phname,ph->GetVal()));
229 if (treear->GetIndex()==0)
230 treear->BuildIndex("fSegmentID","fSegmentID");
231 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
239 param->SetZeroSup(2);
241 Int_t zerosup = param->GetZeroSup();
242 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
243 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
245 //Loop over segments of the TPC
247 for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++)
250 if (!param->AdjustSectorRow(segmentID,sec,row))
252 cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
255 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
256 AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec); // noise per given sector
257 digrow->SetID(segmentID);
262 Bool_t digitize = kFALSE;
263 for (Int_t i=0;i<nInputs; i++)
266 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i));
267 gime = rl->GetLoader("TPCLoader");
269 if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
270 digarr[i]->ExpandBuffer();
271 digarr[i]->ExpandTrackBuffer();
272 nrows = digarr[i]->GetNRows();
273 ncols = digarr[i]->GetNCols();
275 if (!GetRegionOfInterest() || (i == 0)) digitize = kTRUE;
279 if (GetRegionOfInterest() && !digitize) break;
281 if (!digitize) continue;
283 digrow->Allocate(nrows,ncols);
284 digrow->AllocateTrack(3);
287 Int_t label[1000]; //stack for 300 events
290 Int_t nElems = nrows*ncols;
292 for (Int_t i=0;i<nInputs; i++)
294 pdig[i] = digarr[i]->GetDigits();
295 ptr[i] = digarr[i]->GetTracks();
298 Short_t *pdig1= digrow->GetDigits();
299 Int_t *ptr1= digrow->GetTracks() ;
303 for (Int_t elem=0;elem<nElems; elem++)
309 for (Int_t i=0;i<nInputs; i++) if (active[i])
311 // q += digarr[i]->GetDigitFast(rows,col);
314 for (Int_t tr=0;tr<3;tr++)
316 // Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
317 Int_t lab = ptr[i][tr*nElems];
318 if ( (lab > 1) && *(pdig[i])>zerosup)
320 label[labptr]=lab+masks[i];
327 q/=16.; //conversion factor
328 Float_t gain = gainROC->GetValue(row,elem/nrows); // get gain for given - pad-row pad
330 //printf("problem\n");
333 Float_t noisePad = noiseROC->GetValue(row,elem/nrows);
334 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
335 Float_t noise = pTPC->GetNoise();
340 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
341 //digrow->SetDigitFast((Short_t)q,rows,col);
343 for (Int_t tr=0;tr<3;tr++)
346 ptr1[tr*nElems] = label[tr];
355 digrow->GlitchFilter();
357 digrow->CompresBuffer(1,zerosup);
358 digrow->CompresTrackBuffer(1);
360 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
361 } //for (Int_t n=0; n<param->GetNRowsTotal(); n++)
364 orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
365 ogime = orl->GetLoader("TPCLoader");
366 ogime->WriteDigits("OVERWRITE");
368 //fDigInput->GetTreeDTPC()->Write(0,TObject::kOverwrite);
371 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
381 //------------------------------------------------------------------------
382 void AliTPCDigitizer::DigitizeSave(Option_t* option)
385 // merge input tree's with summable digits
386 //output stored in TreeTPCD
388 TString optionString = option;
389 if (!strcmp(optionString.Data(),"deb")) {
390 cout<<"AliTPCDigitizer::Digitize: called with option deb "<<endl;
393 //get detector and geometry
394 AliRunLoader *rl, *orl;
395 AliLoader *gime, *ogime;
398 orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
399 ogime = orl->GetLoader("TPCLoader");
401 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
402 //gime = rl->GetLoader("TPCLoader");
403 rl->GetLoader("TPCLoader");
405 AliRun* alirun = rl->GetAliRun();
407 AliTPC *pTPC = (AliTPC *) alirun->GetModule("TPC");
408 AliTPCParam * param = pTPC->GetParam();
409 pTPC->GenerNoise(500000); //create teble with noise
410 printf("noise %f \n", param->GetNoise()*param->GetNoiseNormFac());
412 Int_t nInputs = fDigInput->GetNinputs();
413 // stupid protection...
414 if (nInputs <= 0) return;
416 Int_t * masks = new Int_t[nInputs];
417 for (Int_t i=0; i<nInputs;i++)
418 masks[i]= fDigInput->GetMask(i);
420 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
421 for(Int_t ii=0;ii<nInputs;ii++) digarr[ii]=0;
423 for (Int_t i1=0;i1<nInputs; i1++)
427 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i1));
428 gime = rl->GetLoader("TPCLoader");
430 TTree * treear = gime->TreeS();
433 cerr<<" TPC - not existing input = \n"<<i1<<" ";
435 for(Int_t i=0; i<nInputs; i++) delete digarr[i];
440 TBranch * br = treear->GetBranch("fSegmentID");
441 if (br) br->GetFile()->cd();
442 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
445 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
446 gime = rl->GetLoader("TPCLoader");
447 Stat_t nentries = gime->TreeS()->GetEntries();
450 //create branch's in TPC treeD
451 AliSimDigits * digrow = new AliSimDigits;
452 TTree * tree = ogime->TreeD();
454 tree->Branch("Segment","AliSimDigits",&digrow);
455 param->SetZeroSup(2);
457 Int_t zerosup = param->GetZeroSup();
458 //Loop over segments of the TPC
460 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
461 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
462 for (Int_t n=0; n<nentries; n++) {
463 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
464 gime = rl->GetLoader("TPCLoader");
465 gime->TreeS()->GetEvent(n);
467 digarr[0]->ExpandBuffer();
468 digarr[0]->ExpandTrackBuffer();
471 for (Int_t i=1;i<nInputs; i++){
472 // fDigInput->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
473 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i));
474 gime = rl->GetLoader("TPCLoader");
475 gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
476 digarr[i]->ExpandBuffer();
477 digarr[i]->ExpandTrackBuffer();
478 if ((digarr[0]->GetID()-digarr[i]->GetID())>0)
484 if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
485 cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
489 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
490 AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec); // noise per given sector
491 digrow->SetID(digarr[0]->GetID());
493 Int_t nrows = digarr[0]->GetNRows();
494 Int_t ncols = digarr[0]->GetNCols();
495 digrow->Allocate(nrows,ncols);
496 digrow->AllocateTrack(3);
499 Int_t label[1000]; //stack for 300 events
504 for (Int_t rows=0;rows<nrows; rows++){
505 for (Int_t col=0;col<ncols; col++){
510 for (Int_t i=0;i<nInputs; i++){
511 q += digarr[i]->GetDigitFast(rows,col);
514 for (Int_t tr=0;tr<3;tr++) {
515 Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
516 //Int_t lab = ptr[i][tr*nElems];
518 label[labptr]=lab+masks[i];
526 q/=16.; //conversion factor
527 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
528 Float_t gain = gainROC->GetValue(row,col);
530 Float_t noisePad = noiseROC->GetValue(row, col);
532 Float_t noise = pTPC->GetNoise();
538 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
539 digrow->SetDigitFast((Short_t)q,rows,col);
540 // *pdig1 =Short_t(q);
541 for (Int_t tr=0;tr<3;tr++){
543 ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
544 //ptr1[tr*nElems] = label[tr];
546 // ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
547 // ptr1[tr*nElems] = 1;
555 digrow->CompresBuffer(1,zerosup);
556 digrow->CompresTrackBuffer(1);
558 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
560 // printf("end TPC merging - end -Tree %s\t%p\n",fDigInput->GetInputTreeH(0)->GetName(),fDigInput->GetInputTreeH(0)->GetListOfBranches()->At(3));
561 //fDigInput->GetTreeDTPC()->Write(0,TObject::kOverwrite);
562 ogime->WriteDigits("OVERWRITE");
564 for (Int_t i=1;i<nInputs; i++)
566 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i));
567 gime = rl->GetLoader("TPCLoader");
568 gime->UnloadSDigits();
570 ogime->UnloadDigits();
573 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];