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(),fDebug(0)
46 // Default ctor - don't use it
51 //___________________________________________
52 AliTPCDigitizer::AliTPCDigitizer(AliRunDigitizer* manager)
53 :AliDigitizer(manager),fDebug(0)
56 // ctor which should be used
58 AliDebug(2,"(AliRunDigitizer* manager) was processed");
61 //------------------------------------------------------------------------
62 AliTPCDigitizer::~AliTPCDigitizer()
69 //------------------------------------------------------------------------
70 Bool_t AliTPCDigitizer::Init()
78 //------------------------------------------------------------------------
79 void AliTPCDigitizer::Exec(Option_t* option)
83 //------------------------------------------------------------------------
84 void AliTPCDigitizer::ExecFast(Option_t* option)
87 // merge input tree's with summable digits
88 //output stored in TreeTPCD
91 TString optionString = option;
92 if (optionString.Data() == "deb") {
93 cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
96 //get detector and geometry
99 AliRunLoader *rl, *orl;
100 AliLoader *gime, *ogime;
104 Warning("ExecFast","gAlice is NULL. Loading from input 0");
105 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
108 Error("ExecFast","Can not find Run Loader for input 0. Can not proceed.");
114 AliTPC *pTPC = (AliTPC *) gAlice->GetModule("TPC");
115 AliTPCParam * param = pTPC->GetParam();
117 sprintf(s,param->GetTitle());
118 sprintf(ss,"75x40_100x60");
120 printf("2 pad-length geom hits with 3 pad-lenght geom digits...\n");
122 param=new AliTPCParamSR();
125 sprintf(ss,"75x40_100x60_150x60");
126 if(strcmp(s,ss)!=0) {
127 printf("No TPC parameters found...\n");
132 pTPC->GenerNoise(500000); //create teble with noise
134 Int_t nInputs = fManager->GetNinputs();
135 Int_t * masks = new Int_t[nInputs];
136 for (Int_t i=0; i<nInputs;i++)
137 masks[i]= fManager->GetMask(i);
138 Short_t **pdig= new Short_t*[nInputs]; //pointers to the expanded digits array
139 Int_t **ptr= new Int_t*[nInputs]; //pointers to teh expanded tracks array
140 Bool_t *active= new Bool_t[nInputs]; //flag for active input segments
143 //create digits array for given sectors
146 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
147 for (Int_t i1=0;i1<nInputs; i1++)
151 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
152 gime = rl->GetLoader("TPCLoader");
153 gime->LoadSDigits("read");
154 TTree * treear = gime->TreeS();
158 cerr<<"AliTPCDigitizer: Input tree with SDigits not found in"
159 <<" input "<< i1<<endl;
163 if (treear->GetIndex()==0)
164 treear->BuildIndex("fSegmentID","fSegmentID");
165 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
169 //create branch's in TPC treeD
170 AliSimDigits * digrow = new AliSimDigits;
172 orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
173 ogime = orl->GetLoader("TPCLoader");
175 TTree * tree = ogime->TreeD();
178 ogime->MakeTree("D");
179 tree = ogime->TreeD();
181 tree->Branch("Segment","AliSimDigits",&digrow);
184 param->SetZeroSup(2);
186 Int_t zerosup = param->GetZeroSup();
187 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetPadGainFactor();
189 //Loop over segments of the TPC
191 for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++)
194 if (!param->AdjustSectorRow(segmentID,sec,row))
196 cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
199 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
200 digrow->SetID(segmentID);
205 Bool_t digitize = kFALSE;
206 for (Int_t i=0;i<nInputs; i++)
209 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
210 gime = rl->GetLoader("TPCLoader");
212 if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
213 digarr[i]->ExpandBuffer();
214 digarr[i]->ExpandTrackBuffer();
215 nrows = digarr[i]->GetNRows();
216 ncols = digarr[i]->GetNCols();
218 if (!fRegionOfInterest || (i == 0)) digitize = kTRUE;
222 if (fRegionOfInterest && !digitize) break;
224 if (!digitize) continue;
226 digrow->Allocate(nrows,ncols);
227 digrow->AllocateTrack(3);
230 Int_t label[1000]; //stack for 300 events
233 Int_t nElems = nrows*ncols;
235 for (Int_t i=0;i<nInputs; i++)
237 pdig[i] = digarr[i]->GetDigits();
238 ptr[i] = digarr[i]->GetTracks();
241 Short_t *pdig1= digrow->GetDigits();
242 Int_t *ptr1= digrow->GetTracks() ;
246 for (Int_t elem=0;elem<nElems; elem++)
252 for (Int_t i=0;i<nInputs; i++) if (active[i])
254 // q += digarr[i]->GetDigitFast(rows,col);
257 for (Int_t tr=0;tr<3;tr++)
259 // Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
260 Int_t lab = ptr[i][tr*nElems];
261 if ( (lab > 1) && *(pdig[i])>zerosup)
263 label[labptr]=lab+masks[i];
270 q/=16.; //conversion factor
271 Float_t gain = gainROC->GetValue(row,elem/nrows); // get gain for given - pad-row pad
276 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
277 Float_t noise = pTPC->GetNoise();
282 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
283 //digrow->SetDigitFast((Short_t)q,rows,col);
285 for (Int_t tr=0;tr<3;tr++)
288 ptr1[tr*nElems] = label[tr];
295 digrow->CompresBuffer(1,zerosup);
296 digrow->CompresTrackBuffer(1);
298 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
299 } //for (Int_t n=0; n<param->GetNRowsTotal(); n++)
302 orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
303 ogime = orl->GetLoader("TPCLoader");
304 ogime->WriteDigits("OVERWRITE");
306 //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
309 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
319 //------------------------------------------------------------------------
320 void AliTPCDigitizer::ExecSave(Option_t* option)
323 // merge input tree's with summable digits
324 //output stored in TreeTPCD
326 TString optionString = option;
327 if (optionString.Data() == "deb") {
328 cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
331 //get detector and geometry
332 AliRunLoader *rl, *orl;
333 AliLoader *gime, *ogime;
336 orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
337 ogime = orl->GetLoader("TPCLoader");
339 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
340 gime = rl->GetLoader("TPCLoader");
343 AliRun* alirun = rl->GetAliRun();
345 AliTPC *pTPC = (AliTPC *) alirun->GetModule("TPC");
346 AliTPCParam * param = pTPC->GetParam();
347 pTPC->GenerNoise(500000); //create teble with noise
348 printf("noise %f \n", param->GetNoise()*param->GetNoiseNormFac());
350 Int_t nInputs = fManager->GetNinputs();
351 Int_t * masks = new Int_t[nInputs];
352 for (Int_t i=0; i<nInputs;i++)
353 masks[i]= fManager->GetMask(i);
355 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
356 for (Int_t i1=0;i1<nInputs; i1++)
360 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
361 gime = rl->GetLoader("TPCLoader");
363 TTree * treear = gime->TreeS();
364 TBranch * br = treear->GetBranch("fSegmentID");
365 if (br) br->GetFile()->cd();
367 cerr<<" TPC - not existing input = \n"<<i1<<" ";
369 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
372 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
373 gime = rl->GetLoader("TPCLoader");
374 Stat_t nentries = gime->TreeS()->GetEntries();
377 //create branch's in TPC treeD
378 AliSimDigits * digrow = new AliSimDigits;
379 TTree * tree = ogime->TreeD();
381 tree->Branch("Segment","AliSimDigits",&digrow);
382 param->SetZeroSup(2);
384 Int_t zerosup = param->GetZeroSup();
385 //Loop over segments of the TPC
387 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetPadGainFactor();
388 for (Int_t n=0; n<nentries; n++) {
389 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
390 gime = rl->GetLoader("TPCLoader");
391 gime->TreeS()->GetEvent(n);
393 digarr[0]->ExpandBuffer();
394 digarr[0]->ExpandTrackBuffer();
397 for (Int_t i=1;i<nInputs; i++){
398 // fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
399 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
400 gime = rl->GetLoader("TPCLoader");
401 gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
402 digarr[i]->ExpandBuffer();
403 digarr[i]->ExpandTrackBuffer();
404 if ((digarr[0]->GetID()-digarr[i]->GetID())>0)
410 if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
411 cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
415 AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec); // pad gains per given sector
416 digrow->SetID(digarr[0]->GetID());
418 Int_t nrows = digarr[0]->GetNRows();
419 Int_t ncols = digarr[0]->GetNCols();
420 digrow->Allocate(nrows,ncols);
421 digrow->AllocateTrack(3);
424 Int_t label[1000]; //stack for 300 events
429 for (Int_t rows=0;rows<nrows; rows++){
430 for (Int_t col=0;col<ncols; col++){
435 for (Int_t i=0;i<nInputs; i++){
436 q += digarr[i]->GetDigitFast(rows,col);
439 for (Int_t tr=0;tr<3;tr++) {
440 Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
441 //Int_t lab = ptr[i][tr*nElems];
443 label[labptr]=lab+masks[i];
451 q/=16.; //conversion factor
452 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
453 Float_t gain = gainROC->GetValue(row,col);
455 Float_t noise = pTPC->GetNoise();
460 if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
461 digrow->SetDigitFast((Short_t)q,rows,col);
462 // *pdig1 =Short_t(q);
463 for (Int_t tr=0;tr<3;tr++){
465 ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
466 //ptr1[tr*nElems] = label[tr];
468 // ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
469 // ptr1[tr*nElems] = 1;
477 digrow->CompresBuffer(1,zerosup);
478 digrow->CompresTrackBuffer(1);
480 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
482 // printf("end TPC merging - end -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
483 //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
484 ogime->WriteDigits("OVERWRITE");
486 for (Int_t i=1;i<nInputs; i++)
488 rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
489 gime = rl->GetLoader("TPCLoader");
490 gime->UnloadSDigits();
492 ogime->UnloadDigits();
495 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];