1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 **************************************************************************/
16 /* $Id: AliTPCclustererKr.cxx,v 1.7 2008/02/06 17:24:53 matyja Exp $ */
18 //-----------------------------------------------------------------
19 // Implementation of the TPC cluster class
21 // Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl
22 //-----------------------------------------------------------------
24 #include "AliTPCclustererKr.h"
25 #include "AliTPCclusterKr.h"
28 #include "AliPadMax.h"
29 #include "AliSimDigits.h"
31 #include "AliTPCParam.h"
32 #include "AliTPCDigitsArray.h"
33 #include "AliTPCvtpr.h"
34 #include "AliTPCClustersRow.h"
37 //used in raw data finder
38 #include "AliTPCROC.h"
39 #include "AliTPCCalPad.h"
40 #include "AliTPCAltroMapping.h"
41 #include "AliTPCcalibDB.h"
42 #include "AliTPCRawStream.h"
43 #include "AliTPCRecoParam.h"
44 #include "AliTPCReconstructor.h"
45 #include "AliRawReader.h"
46 #include "AliTPCCalROC.h"
48 ClassImp(AliTPCclustererKr)
51 AliTPCclustererKr::AliTPCclustererKr()
60 fIsOldRCUFormat(kFALSE)
63 // default constructor
67 AliTPCclustererKr::AliTPCclustererKr(const AliTPCclustererKr ¶m)
76 fIsOldRCUFormat(kFALSE)
81 fParam = param.fParam;
82 fRecoParam=param.fRecoParam;
83 fIsOldRCUFormat=param.fIsOldRCUFormat;
84 fRawData=param.fRawData;
85 fRowCl =param.fRowCl ;
86 fInput =param.fInput ;
87 fOutput=param.fOutput;
88 fDigarr=param.fDigarr;
91 AliTPCclustererKr & AliTPCclustererKr::operator = (const AliTPCclustererKr & param)
93 fParam = param.fParam;
94 fRecoParam=param.fRecoParam;
95 fIsOldRCUFormat=param.fIsOldRCUFormat;
96 fRawData=param.fRawData;
97 fRowCl =param.fRowCl ;
98 fInput =param.fInput ;
99 fOutput=param.fOutput;
100 fDigarr=param.fDigarr;
104 AliTPCclustererKr::~AliTPCclustererKr()
111 void AliTPCclustererKr::SetRecoParam(AliTPCRecoParam *recoParam)
114 fRecoParam = recoParam;
116 //set default parameters if not specified
117 fRecoParam = AliTPCReconstructor::GetRecoParam();
118 if (!fRecoParam) fRecoParam = AliTPCRecoParam::GetLowFluxParam();
124 ////____________________________________________________________________________
126 void AliTPCclustererKr::SetInput(TTree * tree)
129 // set input tree with digits
132 if (!fInput->GetBranch("Segment")){
133 cerr<<"AliTPCclusterKr::findClusterKr(): no proper input tree !\n";
139 void AliTPCclustererKr::SetOutput(TTree * tree)
145 AliTPCClustersRow clrow;
146 AliTPCClustersRow *pclrow=&clrow;
147 clrow.SetClass("AliTPCclusterKr");
148 clrow.SetArray(1); // to make Clones array
149 fOutput->Branch("Segment","AliTPCClustersRow",&pclrow,32000,200);
152 ////____________________________________________________________________________
154 Int_t AliTPCclustererKr::finderIO()
156 // Krypton cluster finder for simulated events from MC
159 Error("Digits2Clusters", "input tree not initialised");
164 Error("Digits2Clusters", "output tree not initialised");
172 Int_t AliTPCclustererKr::finderIO(AliRawReader* rawReader)
174 // Krypton cluster finder for the TPC raw data
176 // fParam must be defined before
178 // consider noiceROC or not
180 if(rawReader)fRawData=kTRUE; //set flag to data
183 Error("Digits2Clusters", "output tree not initialised");
187 fParam->SetMaxTBin(fRecoParam->GetLastBin());//set number of timebins from reco -> param
188 // used later for memory allocation
190 Bool_t isAltro=kFALSE;
192 AliTPCROC * roc = AliTPCROC::Instance();
193 //AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
194 AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
196 AliTPCRawStream input(rawReader,(AliAltroMapping**)mapping);
198 const Int_t kNIS = fParam->GetNInnerSector();//number of inner sectors
199 const Int_t kNOS = fParam->GetNOuterSector();//number of outer sectors
200 const Int_t kNS = kNIS + kNOS;//all sectors
201 Int_t zeroSup = fParam->GetZeroSup();//zero suppression parameter
204 AliTPCDigitsArray *digarr=new AliTPCDigitsArray(kFALSE);//data not sim
205 digarr->Setup(fParam);//as usually parameters
210 for(Int_t sec = 0; sec < kNS; sec++) {
211 //AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec); // noise per given sector
213 Int_t nRows = 0; //number of rows in sector
214 Int_t nDDLs = 0; //number of DDLs
215 Int_t indexDDL = 0; //DDL index
217 nRows = fParam->GetNRowLow();
221 nRows = fParam->GetNRowUp();
223 indexDDL = (sec-kNIS) * 4 + kNIS * 2;
227 // Load the raw data for corresponding DDLs
230 input.SetOldRCUFormat(fIsOldRCUFormat);
231 rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
235 // Allocate memory for rows in sector (pads(depends on row) x timebins)
236 for(Int_t row = 0; row < nRows; row++) {
237 digarr->CreateRow(sec,row);
238 }//end loop over rows
240 rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
243 // Begin loop over altro data
245 while (input.Next()) {
247 //check sector consistency
248 if (input.GetSector() != sec)
249 AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",sec,input.GetSector()));
251 //check row consistency
252 Short_t iRow = input.GetRow();
253 if (iRow < 0 || iRow >= nRows){
254 AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",
259 //check pad consistency
260 Short_t iPad = input.GetPad();
261 if (iPad < 0 || iPad >= (Short_t)(roc->GetNPads(sec,iRow))) {
262 AliError(Form("Pad index (%d) outside the range (%d -> %d) !",
263 iPad, 0, roc->GetNPads(sec,iRow) ));
267 //check time consistency
268 Short_t iTimeBin = input.GetTime();
269 if ( iTimeBin < fRecoParam->GetFirstBin() || iTimeBin >= fRecoParam->GetLastBin()){
270 //cout<<iTimeBin<<endl;
272 AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",
273 iTimeBin, 0, fRecoParam->GetLastBin() -1));
277 Int_t signal = input.GetSignal();
278 if (signal <= zeroSup) {
280 digarr->GetRow(sec,iRow)->SetDigitFast(0,iTimeBin,iPad);
283 digarr->GetRow(sec,iRow))->SetDigitFast(signal,iTimeBin,iPad);
284 }//end of loop over altro data
285 }//end of loop over sectors
288 if(isAltro) findClusterKrIO();
294 ////____________________________________________________________________________
295 Int_t AliTPCclustererKr::findClusterKrIO()
297 //fParam and fDigarr must be set to run this method
299 Int_t cluster_counter=0;
300 const Short_t NTotalSector=fParam->GetNSector();//number of sectors
301 for(Short_t sec=0; sec<NTotalSector; sec++){
303 //vector of maxima for each sector
304 std::vector<AliPadMax*> maxima_in_sector;
306 const Short_t Nrows=fParam->GetNRow(sec);//number of rows in sector
307 for(Short_t row=0; row<Nrows; row++){
308 AliSimDigits *digrow;
310 digrow = (AliSimDigits*)fDigarr->GetRow(sec,row);//real data
312 digrow = (AliSimDigits*)fDigarr->LoadRow(sec,row);//MC
314 if(digrow){//if pointer exist
315 digrow->ExpandBuffer(); //decrunch
316 const Short_t npads = digrow->GetNCols(); // number of pads
317 const Short_t ntime = digrow->GetNRows(); // number of timebins
318 for(Short_t np=0;np<npads;np++){
320 Short_t tb_max=-1;//timebin of maximum
321 Short_t value_maximum=-1;//value of maximum in adc
322 Short_t increase_begin=-1;//timebin when increase starts
323 Short_t sum_adc=0;//sum of adc on the pad in maximum surrounding
324 bool if_increase_begin=true;//flag - check if increasing start
325 bool if_maximum=false;//flag - check if it could be maximum
327 for(Short_t nt=0;nt<ntime;nt++){
328 Short_t adc = digrow->GetDigitFast(nt,np);
331 if(nt-1-increase_begin<1){//at least 2 time bins
336 if_increase_begin=true;
340 //insert maximum, default values and set flags
341 AliPadMax *one_maximum = new AliPadMax(AliTPCvtpr(value_maximum,
348 maxima_in_sector.push_back(one_maximum);
354 if_increase_begin=true;
360 if(if_increase_begin){
361 if_increase_begin=false;
365 if(adc>value_maximum){
371 if(nt==ntime-1 && if_maximum && ntime-1-increase_begin>2){//on the edge
372 //insert maximum, default values and set flags
373 AliPadMax *one_maximum = new AliPadMax(AliTPCvtpr(value_maximum,
380 maxima_in_sector.push_back(one_maximum);
386 if_increase_begin=true;
391 }//end loop over timebins
392 }//end loop over pads
394 // cout<<"Pointer does not exist!!"<<endl;
395 }//end if poiner exists
396 }//end loop over rows
400 Short_t max_sum_adc=0;
405 for( std::vector<AliPadMax*>::iterator mp1 = maxima_in_sector.begin();
406 mp1 != maxima_in_sector.end(); ++mp1 ) {
408 AliTPCclusterKr *tmp=new AliTPCclusterKr();
410 Short_t n_used_pads=1;
411 Short_t cluster_value=0;
412 cluster_value+=(*mp1)->GetSum();
414 max_dig =(*mp1)->GetAdc() ;
415 max_sum_adc =(*mp1)->GetSum() ;
416 max_nt =(*mp1)->GetTime();
417 max_np =(*mp1)->GetPad() ;
418 max_row =(*mp1)->GetRow() ;
420 AliSimDigits *digrow_tmp;
422 digrow_tmp = (AliSimDigits*)fDigarr->GetRow(sec,(*mp1)->GetRow());
424 digrow_tmp = (AliSimDigits*)fDigarr->LoadRow(sec,(*mp1)->GetRow());
427 digrow_tmp->ExpandBuffer(); //decrunch
428 for(Short_t tb=(*mp1)->GetBegin(); tb<((*mp1)->GetEnd())+1; tb++){
429 Short_t adc_tmp = digrow_tmp->GetDigitFast(tb,(*mp1)->GetPad());
430 AliTPCvtpr *vtpr=new AliTPCvtpr(adc_tmp,tb,(*mp1)->GetPad(),(*mp1)->GetRow());
431 tmp->fCluster.push_back(vtpr);
434 maxima_in_sector.erase(mp1);
437 for( std::vector<AliPadMax*>::iterator mp2 = maxima_in_sector.begin();
438 mp2 != maxima_in_sector.end(); ++mp2 ) {
440 Short_t two_pad_max=4;
441 Short_t two_row_max=3;
443 if(abs(max_row - (*mp2)->GetRow())<two_row_max &&
444 abs(max_np - (*mp2)->GetPad())<two_pad_max &&
445 abs(max_nt - (*mp2)->GetTime())<7){
447 cluster_value+=(*mp2)->GetSum();
450 AliSimDigits *digrow_tmp1;
452 digrow_tmp1 = (AliSimDigits*)fDigarr->GetRow(sec,(*mp2)->GetRow());
454 digrow_tmp1 = (AliSimDigits*)fDigarr->LoadRow(sec,(*mp2)->GetRow());
456 digrow_tmp1->ExpandBuffer(); //decrunch
458 for(Short_t tb=(*mp2)->GetBegin(); tb<(*mp2)->GetEnd()+1; tb++){
459 Short_t adc_tmp = digrow_tmp1->GetDigitFast(tb,(*mp2)->GetPad());
460 AliTPCvtpr *vtpr=new AliTPCvtpr(adc_tmp,tb,(*mp2)->GetPad(),(*mp2)->GetRow());
461 tmp->fCluster.push_back(vtpr);
464 //which one is bigger
465 if( (*mp2)->GetAdc() > max_dig ){
466 max_dig =(*mp2)->GetAdc() ;
467 max_sum_adc =(*mp2)->GetSum() ;
468 max_nt =(*mp2)->GetTime();
469 max_np =(*mp2)->GetPad() ;
470 max_row =(*mp2)->GetRow() ;
471 } else if ( (*mp2)->GetAdc() == max_dig ){
472 if( (*mp2)->GetSum() > max_sum_adc){
473 max_dig =(*mp2)->GetAdc() ;
474 max_sum_adc =(*mp2)->GetSum() ;
475 max_nt =(*mp2)->GetTime();
476 max_np =(*mp2)->GetPad() ;
477 max_row =(*mp2)->GetRow() ;
480 maxima_in_sector.erase(mp2);
485 //through out ADC=6,7 on 1 pad, 2 tb and ADC=12 on 2 pads,2 tb
486 if(n_used_pads==1 && cluster_value/tmp->fCluster.size()<3.6)continue;
487 if(n_used_pads==2 && cluster_value/tmp->fCluster.size()<3.1)continue;
489 tmp->SetADCcluster(cluster_value);
490 tmp->SetNpads(n_used_pads);
491 tmp->SetMax(AliTPCvtpr(max_dig,max_nt,max_np,max_row));
497 //save each cluster into file
499 AliTPCClustersRow *clrow= new AliTPCClustersRow();
501 clrow->SetClass("AliTPCclusterKr");
503 fOutput->GetBranch("Segment")->SetAddress(&clrow);
506 TClonesArray * arr = fRowCl->GetArray();
507 AliTPCclusterKr * cl = new ((*arr)[tmp_cluster]) AliTPCclusterKr(*tmp);
512 //end of save each cluster into file adc.root
516 cout<<"Number of clusters in event: "<<cluster_counter<<endl;
521 ////____________________________________________________________________________