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 **************************************************************************/
20 #include <TObjArray.h>
22 #include <TDirectory.h>
23 #include <Riostream.h>
25 #include "AliTPCDigitizer.h"
28 #include "AliTPCParam.h"
29 #include "AliTPCParamSR.h"
32 #include "AliRunDigitizer.h"
33 #include "AliSimDigits.h"
36 #include "AliTPCcalibDB.h"
37 #include "AliTPCCalPad.h"
38 #include "AliTPCCalROC.h"
40 ClassImp(AliTPCDigitizer)
42 //___________________________________________
43 AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer()
45 // Default ctor - don't use it
49 //___________________________________________
50 AliTPCDigitizer::AliTPCDigitizer(AliRunDigitizer* manager)
51 :AliDigitizer(manager)
53 // ctor which should be used
55 AliDebug(2,"(AliRunDigitizer* manager) was processed");
58 //------------------------------------------------------------------------
59 AliTPCDigitizer::~AliTPCDigitizer()
66 //------------------------------------------------------------------------
67 Bool_t AliTPCDigitizer::Init()
75 //------------------------------------------------------------------------
76 void AliTPCDigitizer::Exec(Option_t* option)
80 //------------------------------------------------------------------------
81 void AliTPCDigitizer::ExecFast(Option_t* option)
84 // merge input tree's with summable digits
85 //output stored in TreeTPCD
88 TString optionString = option;
89 if (optionString.Data() == "deb") {
90 cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
93 //get detector and geometry
96 AliRunLoader *rl, *orl;
97 AliLoader *gime, *ogime;
101 Warning("ExecFast","gAlice is NULL. Loading from input 0");
102 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
105 Error("ExecFast","Can not find Run Loader for input 0. Can not proceed.");
111 AliTPC *pTPC = (AliTPC *) gAlice->GetModule("TPC");
112 AliTPCParam * param = pTPC->GetParam();
114 sprintf(s,param->GetTitle());
115 sprintf(ss,"75x40_100x60");
117 printf("2 pad-length geom hits with 3 pad-lenght geom digits...\n");
119 param=new AliTPCParamSR();
122 sprintf(ss,"75x40_100x60_150x60");
123 if(strcmp(s,ss)!=0) {
124 printf("No TPC parameters found...\n");
129 pTPC->GenerNoise(500000); //create teble with noise
131 Int_t nInputs = fManager->GetNinputs();
132 Int_t * masks = new Int_t[nInputs];
133 for (Int_t i=0; i<nInputs;i++)
134 masks[i]= fManager->GetMask(i);
135 Short_t **pdig= new Short_t*[nInputs]; //pointers to the expanded digits array
136 Int_t **ptr= new Int_t*[nInputs]; //pointers to teh expanded tracks array
137 Bool_t *active= new Bool_t[nInputs]; //flag for active input segments
140 //create digits array for given sectors
143 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
144 for (Int_t i1=0;i1<nInputs; i1++)
148 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
149 gime = rl->GetLoader("TPCLoader");
150 gime->LoadSDigits("read");
151 TTree * treear = gime->TreeS();
155 cerr<<"AliTPCDigitizer: Input tree with SDigits not found in"
156 <<" input "<< i1<<endl;
160 if (treear->GetIndex()==0)
161 treear->BuildIndex("fSegmentID","fSegmentID");
162 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
166 //create branch's in TPC treeD
167 AliSimDigits * digrow = new AliSimDigits;
169 orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
170 ogime = orl->GetLoader("TPCLoader");
172 TTree * tree = ogime->TreeD();
175 ogime->MakeTree("D");
176 tree = ogime->TreeD();
178 tree->Branch("Segment","AliSimDigits",&digrow);
181 param->SetZeroSup(2);
183 Int_t zerosup = param->GetZeroSup();
184 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetPadGainFactor();
186 //Loop over segments of the TPC
188 for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++)
191 if (!param->AdjustSectorRow(segmentID,sec,row))
193 cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
196 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
197 digrow->SetID(segmentID);
202 Bool_t digitize = kFALSE;
203 for (Int_t i=0;i<nInputs; i++)
206 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
207 gime = rl->GetLoader("TPCLoader");
209 if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
210 digarr[i]->ExpandBuffer();
211 digarr[i]->ExpandTrackBuffer();
212 nrows = digarr[i]->GetNRows();
213 ncols = digarr[i]->GetNCols();
215 if (!fRegionOfInterest || (i == 0)) digitize = kTRUE;
219 if (fRegionOfInterest && !digitize) break;
221 if (!digitize) continue;
223 digrow->Allocate(nrows,ncols);
224 digrow->AllocateTrack(3);
227 Int_t label[1000]; //stack for 300 events
230 Int_t nElems = nrows*ncols;
232 for (Int_t i=0;i<nInputs; i++)
234 pdig[i] = digarr[i]->GetDigits();
235 ptr[i] = digarr[i]->GetTracks();
238 Short_t *pdig1= digrow->GetDigits();
239 Int_t *ptr1= digrow->GetTracks() ;
243 for (Int_t elem=0;elem<nElems; elem++)
249 for (Int_t i=0;i<nInputs; i++) if (active[i])
251 // q += digarr[i]->GetDigitFast(rows,col);
254 for (Int_t tr=0;tr<3;tr++)
256 // Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
257 Int_t lab = ptr[i][tr*nElems];
258 if ( (lab > 1) && *(pdig[i])>zerosup)
260 label[labptr]=lab+masks[i];
267 q/=16.; //conversion factor
268 Float_t gain = gainROC->GetValue(row,elem/nrows); // get gain for given - pad-row pad
273 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
274 Float_t noise = pTPC->GetNoise();
279 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
280 //digrow->SetDigitFast((Short_t)q,rows,col);
282 for (Int_t tr=0;tr<3;tr++)
285 ptr1[tr*nElems] = label[tr];
292 digrow->CompresBuffer(1,zerosup);
293 digrow->CompresTrackBuffer(1);
295 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
296 } //for (Int_t n=0; n<param->GetNRowsTotal(); n++)
299 orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
300 ogime = orl->GetLoader("TPCLoader");
301 ogime->WriteDigits("OVERWRITE");
303 //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
306 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
316 //------------------------------------------------------------------------
317 void AliTPCDigitizer::ExecSave(Option_t* option)
320 // merge input tree's with summable digits
321 //output stored in TreeTPCD
323 TString optionString = option;
324 if (optionString.Data() == "deb") {
325 cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
328 //get detector and geometry
329 AliRunLoader *rl, *orl;
330 AliLoader *gime, *ogime;
333 orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
334 ogime = orl->GetLoader("TPCLoader");
336 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
337 gime = rl->GetLoader("TPCLoader");
340 AliRun* alirun = rl->GetAliRun();
342 AliTPC *pTPC = (AliTPC *) alirun->GetModule("TPC");
343 AliTPCParam * param = pTPC->GetParam();
344 pTPC->GenerNoise(500000); //create teble with noise
345 printf("noise %f \n", param->GetNoise()*param->GetNoiseNormFac());
347 Int_t nInputs = fManager->GetNinputs();
348 Int_t * masks = new Int_t[nInputs];
349 for (Int_t i=0; i<nInputs;i++)
350 masks[i]= fManager->GetMask(i);
352 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
353 for (Int_t i1=0;i1<nInputs; i1++)
357 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
358 gime = rl->GetLoader("TPCLoader");
360 TTree * treear = gime->TreeS();
361 TBranch * br = treear->GetBranch("fSegmentID");
362 if (br) br->GetFile()->cd();
364 cerr<<" TPC - not existing input = \n"<<i1<<" ";
366 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
369 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
370 gime = rl->GetLoader("TPCLoader");
371 Stat_t nentries = gime->TreeS()->GetEntries();
374 //create branch's in TPC treeD
375 AliSimDigits * digrow = new AliSimDigits;
376 TTree * tree = ogime->TreeD();
378 tree->Branch("Segment","AliSimDigits",&digrow);
379 param->SetZeroSup(2);
381 Int_t zerosup = param->GetZeroSup();
382 //Loop over segments of the TPC
384 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetPadGainFactor();
385 for (Int_t n=0; n<nentries; n++) {
386 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
387 gime = rl->GetLoader("TPCLoader");
388 gime->TreeS()->GetEvent(n);
390 digarr[0]->ExpandBuffer();
391 digarr[0]->ExpandTrackBuffer();
394 for (Int_t i=1;i<nInputs; i++){
395 // fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
396 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
397 gime = rl->GetLoader("TPCLoader");
398 gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
399 digarr[i]->ExpandBuffer();
400 digarr[i]->ExpandTrackBuffer();
401 if ((digarr[0]->GetID()-digarr[i]->GetID())>0)
407 if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
408 cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
412 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
413 digrow->SetID(digarr[0]->GetID());
415 Int_t nrows = digarr[0]->GetNRows();
416 Int_t ncols = digarr[0]->GetNCols();
417 digrow->Allocate(nrows,ncols);
418 digrow->AllocateTrack(3);
421 Int_t label[1000]; //stack for 300 events
426 for (Int_t rows=0;rows<nrows; rows++){
427 for (Int_t col=0;col<ncols; col++){
432 for (Int_t i=0;i<nInputs; i++){
433 q += digarr[i]->GetDigitFast(rows,col);
436 for (Int_t tr=0;tr<3;tr++) {
437 Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
438 //Int_t lab = ptr[i][tr*nElems];
440 label[labptr]=lab+masks[i];
448 q/=16.; //conversion factor
449 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
450 Float_t gain = gainROC->GetValue(row,col);
452 Float_t noise = pTPC->GetNoise();
457 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
458 digrow->SetDigitFast((Short_t)q,rows,col);
459 // *pdig1 =Short_t(q);
460 for (Int_t tr=0;tr<3;tr++){
462 ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
463 //ptr1[tr*nElems] = label[tr];
465 // ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
466 // ptr1[tr*nElems] = 1;
474 digrow->CompresBuffer(1,zerosup);
475 digrow->CompresTrackBuffer(1);
477 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
479 // printf("end TPC merging - end -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
480 //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
481 ogime->WriteDigits("OVERWRITE");
483 for (Int_t i=1;i<nInputs; i++)
485 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
486 gime = rl->GetLoader("TPCLoader");
487 gime->UnloadSDigits();
489 ogime->UnloadDigits();
492 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];