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"
59 ClassImp(AliTPCDigitizer)
61 //___________________________________________
62 AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer(),fDebug(0)
65 // Default ctor - don't use it
70 //___________________________________________
71 AliTPCDigitizer::AliTPCDigitizer(AliDigitizationInput* digInput)
72 :AliDigitizer(digInput),fDebug(0)
75 // ctor which should be used
77 AliDebug(2,"(AliDigitizationInput* digInput) was processed");
80 //------------------------------------------------------------------------
81 AliTPCDigitizer::~AliTPCDigitizer()
88 //------------------------------------------------------------------------
89 Bool_t AliTPCDigitizer::Init()
97 //------------------------------------------------------------------------
98 void AliTPCDigitizer::Digitize(Option_t* option)
100 DigitizeFast(option);
102 //------------------------------------------------------------------------
103 void AliTPCDigitizer::DigitizeFast(Option_t* option)
106 // merge input tree's with summable digits
107 //output stored in TreeTPCD
110 TString optionString = option;
111 if (!strcmp(optionString.Data(),"deb")) {
112 cout<<"AliTPCDigitizer:::DigitizeFast called with option deb "<<endl;
115 //get detector and geometry
118 AliRunLoader *rl, *orl;
119 AliLoader *gime, *ogime;
123 Warning("DigitizeFast","gAlice is NULL. Loading from input 0");
124 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
127 Error("DigitizeFast","Can not find Run Loader for input 0. Can not proceed.");
133 AliTPC *pTPC = (AliTPC *) gAlice->GetModule("TPC");
134 AliTPCParam * param = pTPC->GetParam();
136 //sprintf(s,param->GetTitle());
137 snprintf(s,100,"%s",param->GetTitle());
138 //sprintf(ss,"75x40_100x60");
139 snprintf(ss,100,"75x40_100x60");
141 printf("2 pad-length geom hits with 3 pad-lenght geom digits...\n");
143 param=new AliTPCParamSR();
146 //sprintf(ss,"75x40_100x60_150x60");
147 snprintf(ss,100,"75x40_100x60_150x60");
148 if(strcmp(s,ss)!=0) {
149 printf("No TPC parameters found...\n");
154 pTPC->GenerNoise(500000); //create table with noise
156 Int_t nInputs = fDigInput->GetNinputs();
157 Int_t * masks = new Int_t[nInputs];
158 for (Int_t i=0; i<nInputs;i++)
159 masks[i]= fDigInput->GetMask(i);
160 Short_t **pdig= new Short_t*[nInputs]; //pointers to the expanded digits array
161 Int_t **ptr= new Int_t*[nInputs]; //pointers to the expanded tracks array
162 Bool_t *active= new Bool_t[nInputs]; //flag for active input segments
165 //create digits array for given sectors
168 //create branch's in TPC treeD
169 orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
170 ogime = orl->GetLoader("TPCLoader");
171 TTree * tree = ogime->TreeD();
172 AliSimDigits * digrow = new AliSimDigits;
176 ogime->MakeTree("D");
177 tree = ogime->TreeD();
179 tree->Branch("Segment","AliSimDigits",&digrow);
181 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
182 for (Int_t i1=0;i1<nInputs; i1++)
186 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i1));
187 gime = rl->GetLoader("TPCLoader");
188 gime->LoadSDigits("read");
189 TTree * treear = gime->TreeS();
193 cerr<<"AliTPCDigitizer: Input tree with SDigits not found in"
194 <<" input "<< i1<<endl;
195 for (Int_t i2=0;i2<i1+1; i2++){
197 if(digarr[i2]) delete digarr[i2];
207 //sprintf(phname,"lhcphase%d",i1);
208 snprintf(phname,100,"lhcphase%d",i1);
209 TParameter<float> *ph = (TParameter<float>*)treear->GetUserInfo()
210 ->FindObject("lhcphase0");
212 cerr<<"AliTPCDigitizer: LHC phase not found in"
213 <<" input "<< i1<<endl;
214 for (Int_t i2=0;i2<i1+1; i2++){
215 if(digarr[i2]) delete digarr[i2];
224 tree->GetUserInfo()->Add(new TParameter<float>(phname,ph->GetVal()));
226 if (treear->GetIndex()==0)
227 treear->BuildIndex("fSegmentID","fSegmentID");
228 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
236 param->SetZeroSup(2);
238 Int_t zerosup = param->GetZeroSup();
239 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
240 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
242 //Loop over segments of the TPC
244 for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++)
247 if (!param->AdjustSectorRow(segmentID,sec,row))
249 cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
252 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
253 AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec); // noise per given sector
254 digrow->SetID(segmentID);
259 Bool_t digitize = kFALSE;
260 for (Int_t i=0;i<nInputs; i++)
263 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i));
264 gime = rl->GetLoader("TPCLoader");
266 if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
267 digarr[i]->ExpandBuffer();
268 digarr[i]->ExpandTrackBuffer();
269 nrows = digarr[i]->GetNRows();
270 ncols = digarr[i]->GetNCols();
272 if (!GetRegionOfInterest() || (i == 0)) digitize = kTRUE;
276 if (GetRegionOfInterest() && !digitize) break;
278 if (!digitize) continue;
280 digrow->Allocate(nrows,ncols);
281 digrow->AllocateTrack(3);
284 Int_t label[1000]; //stack for 300 events
287 Int_t nElems = nrows*ncols;
289 for (Int_t i=0;i<nInputs; i++)
291 pdig[i] = digarr[i]->GetDigits();
292 ptr[i] = digarr[i]->GetTracks();
295 Short_t *pdig1= digrow->GetDigits();
296 Int_t *ptr1= digrow->GetTracks() ;
300 for (Int_t elem=0;elem<nElems; elem++)
306 for (Int_t i=0;i<nInputs; i++) if (active[i])
308 // q += digarr[i]->GetDigitFast(rows,col);
311 for (Int_t tr=0;tr<3;tr++)
313 // Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
314 Int_t lab = ptr[i][tr*nElems];
315 if ( (lab > 1) && *(pdig[i])>zerosup)
317 label[labptr]=lab+masks[i];
324 q/=16.; //conversion factor
325 Float_t gain = gainROC->GetValue(row,elem/nrows); // get gain for given - pad-row pad
327 //printf("problem\n");
330 Float_t noisePad = noiseROC->GetValue(row,elem/nrows);
331 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
332 Float_t noise = pTPC->GetNoise();
337 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
338 //digrow->SetDigitFast((Short_t)q,rows,col);
340 for (Int_t tr=0;tr<3;tr++)
343 ptr1[tr*nElems] = label[tr];
352 digrow->GlitchFilter();
354 digrow->CompresBuffer(1,zerosup);
355 digrow->CompresTrackBuffer(1);
357 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
358 } //for (Int_t n=0; n<param->GetNRowsTotal(); n++)
361 orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
362 ogime = orl->GetLoader("TPCLoader");
363 ogime->WriteDigits("OVERWRITE");
365 //fDigInput->GetTreeDTPC()->Write(0,TObject::kOverwrite);
368 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
378 //------------------------------------------------------------------------
379 void AliTPCDigitizer::DigitizeSave(Option_t* option)
382 // merge input tree's with summable digits
383 //output stored in TreeTPCD
385 TString optionString = option;
386 if (!strcmp(optionString.Data(),"deb")) {
387 cout<<"AliTPCDigitizer::Digitize: called with option deb "<<endl;
390 //get detector and geometry
391 AliRunLoader *rl, *orl;
392 AliLoader *gime, *ogime;
395 orl = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
396 ogime = orl->GetLoader("TPCLoader");
398 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
399 //gime = rl->GetLoader("TPCLoader");
400 rl->GetLoader("TPCLoader");
402 AliRun* alirun = rl->GetAliRun();
404 AliTPC *pTPC = (AliTPC *) alirun->GetModule("TPC");
405 AliTPCParam * param = pTPC->GetParam();
406 pTPC->GenerNoise(500000); //create teble with noise
407 printf("noise %f \n", param->GetNoise()*param->GetNoiseNormFac());
409 Int_t nInputs = fDigInput->GetNinputs();
410 // stupid protection...
411 if (nInputs <= 0) return;
413 Int_t * masks = new Int_t[nInputs];
414 for (Int_t i=0; i<nInputs;i++)
415 masks[i]= fDigInput->GetMask(i);
417 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
418 for(Int_t ii=0;ii<nInputs;ii++) digarr[ii]=0;
420 for (Int_t i1=0;i1<nInputs; i1++)
424 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i1));
425 gime = rl->GetLoader("TPCLoader");
427 TTree * treear = gime->TreeS();
430 cerr<<" TPC - not existing input = \n"<<i1<<" ";
432 for(Int_t i=0; i<nInputs; i++) delete digarr[i];
437 TBranch * br = treear->GetBranch("fSegmentID");
438 if (br) br->GetFile()->cd();
439 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
442 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
443 gime = rl->GetLoader("TPCLoader");
444 Stat_t nentries = gime->TreeS()->GetEntries();
447 //create branch's in TPC treeD
448 AliSimDigits * digrow = new AliSimDigits;
449 TTree * tree = ogime->TreeD();
451 tree->Branch("Segment","AliSimDigits",&digrow);
452 param->SetZeroSup(2);
454 Int_t zerosup = param->GetZeroSup();
455 //Loop over segments of the TPC
457 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
458 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
459 for (Int_t n=0; n<nentries; n++) {
460 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(0));
461 gime = rl->GetLoader("TPCLoader");
462 gime->TreeS()->GetEvent(n);
464 digarr[0]->ExpandBuffer();
465 digarr[0]->ExpandTrackBuffer();
468 for (Int_t i=1;i<nInputs; i++){
469 // fDigInput->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
470 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i));
471 gime = rl->GetLoader("TPCLoader");
472 gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
473 digarr[i]->ExpandBuffer();
474 digarr[i]->ExpandTrackBuffer();
475 if ((digarr[0]->GetID()-digarr[i]->GetID())>0)
481 if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
482 cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
486 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
487 AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec); // noise per given sector
488 digrow->SetID(digarr[0]->GetID());
490 Int_t nrows = digarr[0]->GetNRows();
491 Int_t ncols = digarr[0]->GetNCols();
492 digrow->Allocate(nrows,ncols);
493 digrow->AllocateTrack(3);
496 Int_t label[1000]; //stack for 300 events
501 for (Int_t rows=0;rows<nrows; rows++){
502 for (Int_t col=0;col<ncols; col++){
507 for (Int_t i=0;i<nInputs; i++){
508 q += digarr[i]->GetDigitFast(rows,col);
511 for (Int_t tr=0;tr<3;tr++) {
512 Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
513 //Int_t lab = ptr[i][tr*nElems];
515 label[labptr]=lab+masks[i];
523 q/=16.; //conversion factor
524 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
525 Float_t gain = gainROC->GetValue(row,col);
527 Float_t noisePad = noiseROC->GetValue(row, col);
529 Float_t noise = pTPC->GetNoise();
535 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
536 digrow->SetDigitFast((Short_t)q,rows,col);
537 // *pdig1 =Short_t(q);
538 for (Int_t tr=0;tr<3;tr++){
540 ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
541 //ptr1[tr*nElems] = label[tr];
543 // ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
544 // ptr1[tr*nElems] = 1;
552 digrow->CompresBuffer(1,zerosup);
553 digrow->CompresTrackBuffer(1);
555 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
557 // printf("end TPC merging - end -Tree %s\t%p\n",fDigInput->GetInputTreeH(0)->GetName(),fDigInput->GetInputTreeH(0)->GetListOfBranches()->At(3));
558 //fDigInput->GetTreeDTPC()->Write(0,TObject::kOverwrite);
559 ogime->WriteDigits("OVERWRITE");
561 for (Int_t i=1;i<nInputs; i++)
563 rl = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(i));
564 gime = rl->GetLoader("TPCLoader");
565 gime->UnloadSDigits();
567 ogime->UnloadDigits();
570 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];