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 Kr cluster class
21 // Origin: Adam Matyja, INP PAN, adam.matyja@ifj.edu.pl
22 //-----------------------------------------------------------------
25 Instruction - how to use that:
26 There are two macros prepared. One is for preparing clusters from MC
27 samples: FindKrClusters.C. The output is kept in TPC.RecPoints.root.
28 The other macro is prepared for data analysis: FindKrClustersRaw.C.
29 The output is created for each processed file in root file named adc.root.
30 For each data subsample the same named file is created. So be careful
31 do not overwrite them.
37 To run clusterizaton for MC type:
40 If you don't want to use the standard selection criteria then you
43 // load the standard setup
44 AliRunLoader* rl = AliRunLoader::Open("galice.root");
45 AliTPCLoader *tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");
48 gAlice=rl->GetAliRun();
49 TDirectory *cwd = gDirectory;
50 AliTPCv4 *tpc = (AliTPCv4*)gAlice->GetDetector("TPC");
51 Int_t ver = tpc->IsVersion();
53 AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");
54 AliTPCDigitsArray *digarr=new AliTPCDigitsArray;
59 Int_t nevmax=rl->GetNumberOfEvents();
60 for(Int_t nev=0;nev<nevmax ;nev++){
62 TTree* input_tree= tpcl->TreeD();//load tree with digits
63 digarr->ConnectTree(input_tree);
64 TTree *output_tree =tpcl->TreeR();//load output tree
66 AliTPCclustererKr *clusters = new AliTPCclustererKr();
67 clusters->SetParam(param);
68 clusters->SetInput(input_tree);
69 clusters->SetOutput(output_tree);
70 clusters->SetDigArr(digarr);
72 //If you want to change the cluster finder parameters for MC there are
75 //1. signal threshold (everything below the given number is treated as 0)
76 clusters->SetMinAdc(3);
78 //2. number of neighbouring timebins to be considered
79 clusters->SetMinTimeBins(2);
81 //3. distance of the cluster center to the center of a pad in pad-padrow plane
82 //(in cm). Remenber that this is still quantified by pad size.
83 clusters->SetMaxPadRangeCm(2.5);
85 //4. distance of the cluster center to the center of a padrow in pad-padrow
86 //plane (in cm). Remenber that this is still quantified by pad size.
87 clusters->SetMaxRowRangeCm(3.5);
89 //5. distance of the cluster center to the max time bin on a pad (in tackts)
90 //ie. fabs(centerT - time)<7
91 clusters->SetMaxTimeRange(7);
93 //6. cut reduce peak at 0. There are noises which appear mostly as two
94 //timebins on one pad.
95 clusters->SetValueToSize(3.5);
99 tpcl->WriteRecPoints("OVERWRITE");
101 delete rl;//cleans everything
104 ********* DATA *********
107 To run clusterizaton for DATA for file named raw_data.root type:
108 .x FindKrClustersRaw.C("raw_data.root")
110 If you want to change some criteria do the following:
113 // remove Altro warnings
115 AliLog::SetClassDebugLevel("AliTPCRawStream",-5);
116 AliLog::SetClassDebugLevel("AliRawReaderDate",-5);
117 AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5);
118 AliLog::SetModuleDebugLevel("RAW",-5);
121 // Get database with noises
123 // char *ocdbpath = gSystem->Getenv("OCDB_PATH");
124 char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";
126 ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";
128 AliCDBManager * man = AliCDBManager::Instance();
129 man->SetDefaultStorage(ocdbpath);
131 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
132 AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
135 TFile *hfile=new TFile("adc.root","RECREATE","ADC file");
136 // Create a ROOT Tree
137 TTree *mytree = new TTree("Kr","Krypton cluster tree");
140 const char *fileName="data.root";
141 AliRawReader *reader = new AliRawReaderRoot(fileName);
142 //AliRawReader *reader = new AliRawReaderDate(fileName);
144 AliAltroRawStreamFast* stream = new AliAltroRawStreamFast(reader);
145 stream->SelectRawData("TPC");
148 AliTPCclustererKr *clusters = new AliTPCclustererKr();
149 clusters->SetOutput(mytree);
150 clusters->SetRecoParam(0);//standard reco parameters
151 AliTPCParamSR *param=new AliTPCParamSR();
152 clusters->SetParam(param);//TPC parameters(sectors, timebins, etc.)
154 //set cluster finder parameters (from data):
155 //1. zero suppression parameter
156 clusters->SetZeroSup(param->GetZeroSup());
159 clusters->SetFirstBin(60);
162 clusters->SetLastBin(950);
165 clusters->SetMaxNoiseAbs(2);
167 //5. maximal amount of sigma of noise
168 clusters->SetMaxNoiseSigma(3);
170 //The remaining parameters are the same paramters as for MC (see MC section
172 clusters->SetMinAdc(3);
173 clusters->SetMinTimeBins(2);
174 clusters->SetMaxPadRangeCm(2.5);
175 clusters->SetMaxRowRangeCm(3.5);
176 clusters->SetMaxTimeRange(7);
177 clusters->SetValueToSize(3.5);
179 while (reader->NextEvent()) {
180 clusters->FinderIO(reader);
190 #include "AliTPCclustererKr.h"
191 #include "AliTPCclusterKr.h"
195 #include "AliPadMax.h"
196 #include "AliSimDigits.h"
197 #include "AliTPCv4.h"
198 #include "AliTPCParam.h"
199 #include "AliTPCDigitsArray.h"
200 #include "AliTPCvtpr.h"
201 #include "AliTPCClustersRow.h"
206 #include "AliTPCTransform.h"
208 //used in raw data finder
209 #include "AliTPCROC.h"
210 #include "AliTPCCalPad.h"
211 #include "AliTPCAltroMapping.h"
212 #include "AliTPCcalibDB.h"
213 #include "AliTPCRawStream.h"
214 #include "AliTPCRecoParam.h"
215 #include "AliTPCReconstructor.h"
216 #include "AliRawReader.h"
217 #include "AliTPCCalROC.h"
219 ClassImp(AliTPCclustererKr)
222 AliTPCclustererKr::AliTPCclustererKr()
251 // default constructor
255 AliTPCclustererKr::AliTPCclustererKr(const AliTPCclustererKr ¶m)
286 fParam = param.fParam;
287 fRecoParam = param.fRecoParam;
288 fRawData = param.fRawData;
289 fRowCl = param.fRowCl ;
290 fInput = param.fInput ;
291 fOutput = param.fOutput;
292 fDigarr = param.fDigarr;
293 fZeroSup = param.fZeroSup ;
294 fFirstBin = param.fFirstBin ;
295 fLastBin = param.fLastBin ;
296 fMaxNoiseAbs = param.fMaxNoiseAbs ;
297 fMaxNoiseSigma = param.fMaxNoiseSigma ;
298 fMinAdc = param.fMinAdc;
299 fMinTimeBins = param.fMinTimeBins;
300 // fMaxPadRange = param.fMaxPadRange ;
301 // fMaxRowRange = param.fMaxRowRange ;
302 fMaxTimeRange = param.fMaxTimeRange;
303 fValueToSize = param.fValueToSize;
304 fMaxPadRangeCm = param.fMaxPadRangeCm;
305 fMaxRowRangeCm = param.fMaxRowRangeCm;
306 fDebugLevel = param.fDebugLevel;
307 fHistoRow = param.fHistoRow ;
308 fHistoPad = param.fHistoPad ;
309 fHistoTime = param.fHistoTime;
310 fHistoRowPad = param.fHistoRowPad;
314 AliTPCclustererKr & AliTPCclustererKr::operator = (const AliTPCclustererKr & param)
317 // assignment operator
319 fParam = param.fParam;
320 fRecoParam = param.fRecoParam;
321 fRawData = param.fRawData;
322 fRowCl = param.fRowCl ;
323 fInput = param.fInput ;
324 fOutput = param.fOutput;
325 fDigarr = param.fDigarr;
326 fZeroSup = param.fZeroSup ;
327 fFirstBin = param.fFirstBin ;
328 fLastBin = param.fLastBin ;
329 fMaxNoiseAbs = param.fMaxNoiseAbs ;
330 fMaxNoiseSigma = param.fMaxNoiseSigma ;
331 fMinAdc = param.fMinAdc;
332 fMinTimeBins = param.fMinTimeBins;
333 // fMaxPadRange = param.fMaxPadRange ;
334 // fMaxRowRange = param.fMaxRowRange ;
335 fMaxTimeRange = param.fMaxTimeRange;
336 fValueToSize = param.fValueToSize;
337 fMaxPadRangeCm = param.fMaxPadRangeCm;
338 fMaxRowRangeCm = param.fMaxRowRangeCm;
339 fDebugLevel = param.fDebugLevel;
340 fHistoRow = param.fHistoRow ;
341 fHistoPad = param.fHistoPad ;
342 fHistoTime = param.fHistoTime;
343 fHistoRowPad = param.fHistoRowPad;
347 AliTPCclustererKr::~AliTPCclustererKr()
354 void AliTPCclustererKr::SetRecoParam(AliTPCRecoParam *recoParam)
357 // set reconstruction parameters
360 fRecoParam = recoParam;
362 //set default parameters if not specified
363 fRecoParam = AliTPCReconstructor::GetRecoParam();
364 if (!fRecoParam) fRecoParam = AliTPCRecoParam::GetLowFluxParam();
370 ////____________________________________________________________________________
372 void AliTPCclustererKr::SetInput(TTree * tree)
375 // set input tree with digits
378 if (!fInput->GetBranch("Segment")){
379 cerr<<"AliTPCclusterKr::FindClusterKr(): no proper input tree !\n";
385 void AliTPCclustererKr::SetOutput(TTree * tree)
391 AliTPCClustersRow clrow;
392 AliTPCClustersRow *pclrow=&clrow;
393 clrow.SetClass("AliTPCclusterKr");
394 clrow.SetArray(1); // to make Clones array
395 fOutput->Branch("Segment","AliTPCClustersRow",&pclrow,32000,200);
398 ////____________________________________________________________________________
400 Int_t AliTPCclustererKr::FinderIO()
402 // Krypton cluster finder for simulated events from MC
405 Error("Digits2Clusters", "input tree not initialised");
410 Error("Digits2Clusters", "output tree not initialised");
418 Int_t AliTPCclustererKr::FinderIO(AliRawReader* rawReader)
420 // Krypton cluster finder for the TPC raw data
422 // fParam must be defined before
424 if(rawReader)fRawData=kTRUE; //set flag to data
427 Error("Digits2Clusters", "output tree not initialised");
431 fParam->SetMaxTBin(fRecoParam->GetLastBin());//set number of timebins from reco -> param
432 // used later for memory allocation
434 Bool_t isAltro=kFALSE;
436 AliTPCROC * roc = AliTPCROC::Instance();
437 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
438 AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
440 AliTPCRawStream input(rawReader,(AliAltroMapping**)mapping);
442 const Int_t kNIS = fParam->GetNInnerSector();//number of inner sectors
443 const Int_t kNOS = fParam->GetNOuterSector();//number of outer sectors
444 const Int_t kNS = kNIS + kNOS;//all sectors
446 // //set cluster finder parameters
447 // SetZeroSup(fParam->GetZeroSup());//zero suppression parameter
450 // SetMaxNoiseAbs(2);
451 // SetMaxNoiseSigma(3);
454 AliTPCDigitsArray *digarr=new AliTPCDigitsArray(kFALSE);//data not sim
455 digarr->Setup(fParam);//as usually parameters
460 for(Int_t iSec = 0; iSec < kNS; iSec++) {
461 AliTPCCalROC * noiseROC;
463 noiseROC =new AliTPCCalROC(iSec);//noise=0
466 noiseROC = noiseTPC->GetCalROC(iSec); // noise per given sector
468 Int_t nRows = 0; //number of rows in sector
469 Int_t nDDLs = 0; //number of DDLs
470 Int_t indexDDL = 0; //DDL index
472 nRows = fParam->GetNRowLow();
476 nRows = fParam->GetNRowUp();
478 indexDDL = (iSec-kNIS) * 4 + kNIS * 2;
482 // Load the raw data for corresponding DDLs
485 rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
489 // Allocate memory for rows in sector (pads(depends on row) x timebins)
490 for(Int_t iRow = 0; iRow < nRows; iRow++) {
491 digarr->CreateRow(iSec,iRow);
492 }//end loop over rows
495 rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
498 // Begin loop over altro data
500 while (input.Next()) {
502 //check sector consistency
503 if (input.GetSector() != iSec)
504 AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSec,input.GetSector()));
506 Short_t iRow = input.GetRow();
507 Short_t iPad = input.GetPad();
508 Short_t iTimeBin = input.GetTime();
511 fHistoRow->Fill(iRow);
512 fHistoPad->Fill(iPad);
513 fHistoTime->Fill(iTimeBin);
514 fHistoRowPad->Fill(iPad,iRow);
515 }else if(fDebugLevel>=0&&fDebugLevel<72){
516 if(iSec==fDebugLevel){
517 fHistoRow->Fill(iRow);
518 fHistoPad->Fill(iPad);
519 fHistoTime->Fill(iTimeBin);
520 fHistoRowPad->Fill(iPad,iRow);
522 }else if(fDebugLevel==73){
524 fHistoRow->Fill(iRow);
525 fHistoPad->Fill(iPad);
526 fHistoTime->Fill(iTimeBin);
527 fHistoRowPad->Fill(iPad,iRow);
529 }else if(fDebugLevel==74){
531 fHistoRow->Fill(iRow);
532 fHistoPad->Fill(iPad);
533 fHistoTime->Fill(iTimeBin);
534 fHistoRowPad->Fill(iPad,iRow);
538 //check row consistency
539 if (iRow < 0 || iRow >= nRows){
540 AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",
545 //check pad consistency
546 if (iPad < 0 || iPad >= (Short_t)(roc->GetNPads(iSec,iRow))) {
547 AliError(Form("Pad index (%d) outside the range (%d -> %d) !",
548 iPad, 0, roc->GetNPads(iSec,iRow) ));
552 //check time consistency
553 if ( iTimeBin < fRecoParam->GetFirstBin() || iTimeBin >= fRecoParam->GetLastBin()){
554 //cout<<iTimeBin<<endl;
556 AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",
557 iTimeBin, 0, fRecoParam->GetLastBin() -1));
561 Int_t signal = input.GetSignal();
562 if (signal <= fZeroSup ||
563 iTimeBin < fFirstBin ||
566 digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
569 if (!noiseROC) continue;
570 Double_t noiseOnPad = noiseROC->GetValue(iRow,iPad);//noise on given pad and row in sector
571 if (noiseOnPad > fMaxNoiseAbs){
572 digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
573 continue; // consider noisy pad as dead
575 if(signal <= fMaxNoiseSigma * noiseOnPad){
576 digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
580 digarr->GetRow(iSec,iRow)->SetDigitFast(signal,iTimeBin,iPad);
582 }//end of loop over altro data
583 }//end of loop over sectors
586 if(isAltro) FindClusterKrIO();
592 ////____________________________________________________________________________
593 Int_t AliTPCclustererKr::FindClusterKrIO()
595 //fParam and fDigarr must be set to run this method
598 // SetMinAdc(3);//usually is 3
599 // SetMinTimeBins(2);//should be 2 - the best result of shape in MC
600 //// SetMaxPadRange(4);
601 //// SetMaxRowRange(3);
602 // SetMaxTimeRange(7);
603 // SetValueToSize(3.5);//3.5
604 // SetMaxPadRangeCm(2.5);
605 // SetMaxRowRangeCm(3.5);
607 Int_t clusterCounter=0;
608 const Short_t nTotalSector=fParam->GetNSector();//number of sectors
609 for(Short_t iSec=0; iSec<nTotalSector; ++iSec){
611 //vector of maxima for each sector
612 //std::vector<AliPadMax*> maximaInSector;
613 TObjArray *maximaInSector=new TObjArray();//to store AliPadMax*
616 // looking for the maxima on the pad
619 const Short_t kNRows=fParam->GetNRow(iSec);//number of rows in sector
620 for(Short_t iRow=0; iRow<kNRows; ++iRow){
621 AliSimDigits *digrow;
623 digrow = (AliSimDigits*)fDigarr->GetRow(iSec,iRow);//real data
625 digrow = (AliSimDigits*)fDigarr->LoadRow(iSec,iRow);//MC
627 if(digrow){//if pointer exist
628 digrow->ExpandBuffer(); //decrunch
629 const Short_t kNPads = digrow->GetNCols(); // number of pads
630 const Short_t kNTime = digrow->GetNRows(); // number of timebins
631 for(Short_t iPad=0;iPad<kNPads;iPad++){
633 Short_t timeBinMax=-1;//timebin of maximum
634 Short_t valueMaximum=-1;//value of maximum in adc
635 Short_t increaseBegin=-1;//timebin when increase starts
636 Short_t sumAdc=0;//sum of adc on the pad in maximum surrounding
637 bool ifIncreaseBegin=true;//flag - check if increasing started
638 bool ifMaximum=false;//flag - check if it could be maximum
640 for(Short_t iTimeBin=0;iTimeBin<kNTime;iTimeBin++){
641 Short_t adc = digrow->GetDigitFast(iTimeBin,iPad);
642 if(adc<fMinAdc){//standard was 3
644 if(iTimeBin-increaseBegin<fMinTimeBins){//at least 2 time bins
649 ifIncreaseBegin=true;
653 //insert maximum, default values and set flags
654 //Double_t xCord,yCord;
655 //GetXY(iSec,iRow,iPad,xCord,yCord);
656 Double_t x[]={iRow,iPad,iTimeBin};
658 AliTPCTransform trafo;
659 trafo.Transform(x,i,0,1);
661 AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,
671 //maximaInSector.push_back(oneMaximum);
672 maximaInSector->AddLast(oneMaximum);
678 ifIncreaseBegin=true;
685 ifIncreaseBegin=false;
686 increaseBegin=iTimeBin;
689 if(adc>valueMaximum){
695 if(iTimeBin==kNTime-1 && ifMaximum && kNTime-increaseBegin>fMinTimeBins){//on the edge
696 //at least 3 timebins
697 //insert maximum, default values and set flags
698 //Double_t xCord,yCord;
699 //GetXY(iSec,iRow,iPad,xCord,yCord);
700 Double_t x[]={iRow,iPad,iTimeBin};
702 AliTPCTransform trafo;
703 trafo.Transform(x,i,0,1);
704 AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,
714 //maximaInSector.push_back(oneMaximum);
715 maximaInSector->AddLast(oneMaximum);
721 ifIncreaseBegin=true;
726 }//end loop over timebins
727 }//end loop over pads
729 // cout<<"Pointer does not exist!!"<<endl;
730 }//end if poiner exists
731 }//end loop over rows
733 // cout<<"EF" <<maximaInSector->GetEntriesFast()<<" E "<<maximaInSector->GetEntries()<<" "<<maximaInSector->GetLast()<<endl;
734 maximaInSector->Compress();
735 // GetEntriesFast() - liczba wejsc w array of maxima
736 //cout<<"EF" <<maximaInSector->GetEntriesFast()<<" E"<<maximaInSector->GetEntries()<<endl;
744 Short_t maxTimeBin=0;
751 // for( std::vector<AliPadMax*>::iterator mp1 = maximaInSector.begin();
752 // mp1 != maximaInSector.end(); ++mp1 ) {
753 for(Int_t it1 = 0; it1 < maximaInSector->GetEntriesFast(); ++it1 ) {
755 AliPadMax *mp1=(AliPadMax *)maximaInSector->At(it1);
756 AliTPCclusterKr *tmp=new AliTPCclusterKr();
759 Int_t clusterValue=0;
760 clusterValue+=(mp1)->GetSum();
761 list<Short_t> nUsedRows;
762 nUsedRows.push_back((mp1)->GetRow());
764 maxDig =(mp1)->GetAdc() ;
765 maxSumAdc =(mp1)->GetSum() ;
766 maxTimeBin =(mp1)->GetTime();
767 maxPad =(mp1)->GetPad() ;
768 maxRow =(mp1)->GetRow() ;
773 AliSimDigits *digrowTmp;
775 digrowTmp = (AliSimDigits*)fDigarr->GetRow(iSec,(mp1)->GetRow());
777 digrowTmp = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp1)->GetRow());
780 digrowTmp->ExpandBuffer(); //decrunch
782 for(Short_t itb=(mp1)->GetBegin(); itb<((mp1)->GetEnd())+1; itb++){
783 Short_t adcTmp = digrowTmp->GetDigitFast(itb,(mp1)->GetPad());
784 AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp1)->GetPad(),(mp1)->GetRow(),(mp1)->GetX(),(mp1)->GetY(),(mp1)->GetT());
785 //tmp->fCluster.push_back(vtpr);
786 tmp->AddDigitToCluster(vtpr);
788 tmp->SetCenter();//set centr of the cluster
790 //maximaInSector.erase(mp1);
792 //cout<<maximaInSector->GetEntriesFast()<<" ";
793 maximaInSector->RemoveAt(it1);
794 maximaInSector->Compress();
796 // cout<<maximaInSector->GetEntriesFast()<<" "<<endl;
798 // for( std::vector<AliPadMax*>::iterator mp2 = maximaInSector.begin();
799 // mp2 != maximaInSector.end(); ++mp2 ) {
800 for(Int_t it2 = 0; it2 < maximaInSector->GetEntriesFast(); ++it2 ) {
801 AliPadMax *mp2=(AliPadMax *)maximaInSector->At(it2);
803 // if(abs(maxRow - (*mp2)->GetRow()) < fMaxRowRange && //3
804 // abs(maxPad - (*mp2)->GetPad()) < fMaxPadRange && //4
806 if(TMath::Abs(tmp->GetCenterX() - (mp2)->GetX()) < fMaxPadRangeCm &&
807 TMath::Abs(tmp->GetCenterY() - (mp2)->GetY()) < fMaxRowRangeCm &&
808 TMath::Abs(tmp->GetCenterT() - (mp2)->GetT()) < fMaxTimeRange){//7
810 clusterValue+=(mp2)->GetSum();
813 nUsedRows.push_back((mp2)->GetRow());
815 AliSimDigits *digrowTmp1;
817 digrowTmp1 = (AliSimDigits*)fDigarr->GetRow(iSec,(mp2)->GetRow());
819 digrowTmp1 = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp2)->GetRow());
821 digrowTmp1->ExpandBuffer(); //decrunch
823 for(Short_t itb=(mp2)->GetBegin(); itb<(mp2)->GetEnd()+1; itb++){
824 Short_t adcTmp = digrowTmp1->GetDigitFast(itb,(mp2)->GetPad());
825 AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp2)->GetPad(),(mp2)->GetRow(),(mp2)->GetX(),(mp2)->GetY(),(mp2)->GetT());
826 //tmp->fCluster.push_back(vtpr);
827 tmp->AddDigitToCluster(vtpr);
830 tmp->SetCenter();//set center of the cluster
832 //which one is bigger
833 if( (mp2)->GetAdc() > maxDig ){
834 maxDig =(mp2)->GetAdc() ;
835 maxSumAdc =(mp2)->GetSum() ;
836 maxTimeBin =(mp2)->GetTime();
837 maxPad =(mp2)->GetPad() ;
838 maxRow =(mp2)->GetRow() ;
839 maxX =(mp2)->GetX() ;
840 maxY =(mp2)->GetY() ;
841 maxT =(mp2)->GetT() ;
842 } else if ( (mp2)->GetAdc() == maxDig ){
843 if( (mp2)->GetSum() > maxSumAdc){
844 maxDig =(mp2)->GetAdc() ;
845 maxSumAdc =(mp2)->GetSum() ;
846 maxTimeBin =(mp2)->GetTime();
847 maxPad =(mp2)->GetPad() ;
848 maxRow =(mp2)->GetRow() ;
849 maxX =(mp2)->GetX() ;
850 maxY =(mp2)->GetY() ;
851 maxT =(mp2)->GetT() ;
854 maximaInSector->RemoveAt(it2);
855 maximaInSector->Compress();
858 //maximaInSector.erase(mp2);
864 //through out ADC=6,7 on 1 pad, 2 tb and ADC=12 on 2 pads,2 tb
865 //if(nUsedPads==1 && clusterValue/tmp->fCluster.size()<3.6)continue;
866 //if(nUsedPads==2 && clusterValue/tmp->fCluster.size()<3.1)continue;
868 //through out clusters on the edge and noise
869 //if(clusterValue/tmp->fCluster.size()<fValueToSize)continue;
870 if(clusterValue/(tmp->GetSize())<fValueToSize)continue;
872 tmp->SetADCcluster(clusterValue);
873 tmp->SetNPads(nUsedPads);
874 tmp->SetMax(AliTPCvtpr(maxDig,maxTimeBin,maxPad,maxRow,maxX,maxY,maxT));
880 tmp->SetNRows(nUsedRows.size());
886 //save each cluster into file
888 AliTPCClustersRow *clrow= new AliTPCClustersRow();
890 clrow->SetClass("AliTPCclusterKr");
892 fOutput->GetBranch("Segment")->SetAddress(&clrow);
895 TClonesArray * arr = fRowCl->GetArray();
896 AliTPCclusterKr * cl = new ((*arr)[tmpCluster]) AliTPCclusterKr(*tmp);
900 //end of save each cluster into file adc.root
904 cout<<"Number of clusters in event: "<<clusterCounter<<endl;
909 ////____________________________________________________________________________
912 void AliTPCclustererKr::GetXY(Short_t sec,Short_t row,Short_t pad,Double_t& xGlob,Double_t& yGlob){
914 //gives global XY coordinate of the pad
917 Double_t yLocal = fParam->GetPadRowRadii(sec,row);//radius of row in sector in cm
919 Short_t padmax = fParam->GetNPads(sec,row);//number of pads in a given row
921 if(sec<fParam->GetNInnerSector())padXSize=0.4;
923 Double_t xLocal=(pad+0.5-padmax/2.)*padXSize;//x-value of the center of pad
926 fParam->AdjustCosSin((Int_t)sec,cos,sin);//return sinus and cosinus of the sector
928 Double_t xGlob1 = xLocal * cos + yLocal * sin;
929 Double_t yGlob1 = -xLocal * sin + yLocal * cos;
935 xGlob = xGlob1 * TMath::Cos(rot) + yGlob1 * TMath::Sin(rot);
936 yGlob = -xGlob1 * TMath::Sin(rot) + yGlob1 * TMath::Cos(rot);
939 if(sec<18||(sec>=36&&sec<54)) xGlob =-1*xGlob;