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.
33 Additional selection criteria to select the GOLD cluster
35 // open file with clusters
36 TFile f("Krypton.root");
37 TTree * tree = (TTree*)f.Get("Kr")
38 TCut cutR0("cutR0","fADCcluster/fSize<100"); // adjust it according v seetings -
39 TCut cutR1("cutR1","fADCcluster/fSize>7"); // cosmic tracks and noise removal
40 TCut cutR2("cutR2","fMax.fAdc/fADCcluster<0.2"); // digital noise removal
41 TCut cutR3("cutR3","fMax.fAdc/fADCcluster>0.01"); // noise removal
42 TCut cutS1("cutS1","fSize<200"); // adjust it according v seetings - cosmic tracks
43 TCut cutAll = cutR0+cutR1+cutR2+cutR3+cutS1;
44 This values are typical values to be applied in selectors
51 To run clusterizaton for MC type:
54 If you don't want to use the standard selection criteria then you
57 // load the standard setup
58 AliRunLoader* rl = AliRunLoader::Open("galice.root");
59 AliTPCLoader *tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");
62 gAlice=rl->GetAliRun();
63 TDirectory *cwd = gDirectory;
64 AliTPCv4 *tpc = (AliTPCv4*)gAlice->GetDetector("TPC");
65 Int_t ver = tpc->IsVersion();
67 AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");
68 AliTPCDigitsArray *digarr=new AliTPCDigitsArray;
73 Int_t nevmax=rl->GetNumberOfEvents();
74 for(Int_t nev=0;nev<nevmax ;nev++){
76 TTree* input_tree= tpcl->TreeD();//load tree with digits
77 digarr->ConnectTree(input_tree);
78 TTree *output_tree =tpcl->TreeR();//load output tree
80 AliTPCclustererKr *clusters = new AliTPCclustererKr();
81 clusters->SetParam(param);
82 clusters->SetInput(input_tree);
83 clusters->SetOutput(output_tree);
84 clusters->SetDigArr(digarr);
86 //If you want to change the cluster finder parameters for MC there are
89 //1. signal threshold (everything below the given number is treated as 0)
90 clusters->SetMinAdc(3);
92 //2. number of neighbouring timebins to be considered
93 clusters->SetMinTimeBins(2);
95 //3. distance of the cluster center to the center of a pad in pad-padrow plane
96 //(in cm). Remenber that this is still quantified by pad size.
97 clusters->SetMaxPadRangeCm(2.5);
99 //4. distance of the cluster center to the center of a padrow in pad-padrow
100 //plane (in cm). Remenber that this is still quantified by pad size.
101 clusters->SetMaxRowRangeCm(3.5);
103 //5. distance of the cluster center to the max time bin on a pad (in tackts)
104 //ie. fabs(centerT - time)<7
105 clusters->SetMaxTimeRange(7);
107 //6. cut reduce peak at 0. There are noises which appear mostly as two
108 //timebins on one pad.
109 clusters->SetValueToSize(3.5);
112 clusters->finderIO();
113 tpcl->WriteRecPoints("OVERWRITE");
115 delete rl;//cleans everything
118 ********* DATA *********
121 To run clusterizaton for DATA for file named raw_data.root type:
122 .x FindKrClustersRaw.C("raw_data.root")
124 If you want to change some criteria do the following:
127 // remove Altro warnings
129 AliLog::SetClassDebugLevel("AliRawReaderDate",-5);
130 AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5);
131 AliLog::SetModuleDebugLevel("RAW",-5);
134 // Get database with noises
136 // char *ocdbpath = gSystem->Getenv("OCDB_PATH");
137 char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";
139 ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";
141 AliCDBManager * man = AliCDBManager::Instance();
142 man->SetDefaultStorage(ocdbpath);
144 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
145 AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
148 TFile *hfile=new TFile("adc.root","RECREATE","ADC file");
149 // Create a ROOT Tree
150 TTree *mytree = new TTree("Kr","Krypton cluster tree");
153 const char *fileName="data.root";
154 AliRawReader *reader = new AliRawReaderRoot(fileName);
155 //AliRawReader *reader = new AliRawReaderDate(fileName);
157 AliAltroRawStreamFast* stream = new AliAltroRawStreamFast(reader);
158 stream->SelectRawData("TPC");
161 AliTPCclustererKr *clusters = new AliTPCclustererKr();
162 clusters->SetOutput(mytree);
163 clusters->SetRecoParam(0);//standard reco parameters
164 AliTPCParamSR *param=new AliTPCParamSR();
165 clusters->SetParam(param);//TPC parameters(sectors, timebins, etc.)
167 //set cluster finder parameters (from data):
168 //1. zero suppression parameter
169 clusters->SetZeroSup(param->GetZeroSup());
172 clusters->SetFirstBin(60);
175 clusters->SetLastBin(950);
178 clusters->SetMaxNoiseAbs(2);
180 //5. maximal amount of sigma of noise
181 clusters->SetMaxNoiseSigma(3);
183 //The remaining parameters are the same paramters as for MC (see MC section
185 clusters->SetMinAdc(3);
186 clusters->SetMinTimeBins(2);
187 clusters->SetMaxPadRangeCm(2.5);
188 clusters->SetMaxRowRangeCm(3.5);
189 clusters->SetMaxTimeRange(7);
190 clusters->SetValueToSize(3.5);
192 while (reader->NextEvent()) {
193 clusters->FinderIO(reader);
203 #include "AliTPCclustererKr.h"
204 #include "AliTPCclusterKr.h"
208 #include "AliPadMax.h"
209 #include "AliSimDigits.h"
210 //#include "AliTPCv4.h"
211 #include "AliTPCParam.h"
212 #include "AliTPCDigitsArray.h"
213 #include "AliTPCvtpr.h"
214 #include "AliTPCClustersRow.h"
218 #include "TTreeStream.h"
220 #include "AliTPCTransform.h"
222 //used in raw data finder
223 #include "AliTPCROC.h"
224 #include "AliTPCCalPad.h"
225 #include "AliTPCAltroMapping.h"
226 #include "AliTPCcalibDB.h"
227 #include "AliTPCRawStreamV3.h"
228 #include "AliTPCRecoParam.h"
229 #include "AliTPCReconstructor.h"
230 #include "AliRawReader.h"
231 #include "AliTPCCalROC.h"
232 #include "AliRawEventHeaderBase.h"
238 ClassImp(AliTPCclustererKr)
241 AliTPCclustererKr::AliTPCclustererKr()
272 // default constructor
276 AliTPCclustererKr::AliTPCclustererKr(const AliTPCclustererKr ¶m)
309 fParam = param.fParam;
310 fRecoParam = param.fRecoParam;
311 fRawData = param.fRawData;
312 fInput = param.fInput ;
313 fOutput = param.fOutput;
314 fDigarr = param.fDigarr;
315 fZeroSup = param.fZeroSup ;
316 fFirstBin = param.fFirstBin ;
317 fLastBin = param.fLastBin ;
318 fMaxNoiseAbs = param.fMaxNoiseAbs ;
319 fMaxNoiseSigma = param.fMaxNoiseSigma ;
320 fMinAdc = param.fMinAdc;
321 fMinTimeBins = param.fMinTimeBins;
322 // fMaxPadRange = param.fMaxPadRange ;
323 // fMaxRowRange = param.fMaxRowRange ;
324 fMaxTimeRange = param.fMaxTimeRange;
325 fValueToSize = param.fValueToSize;
326 fMaxPadRangeCm = param.fMaxPadRangeCm;
327 fMaxRowRangeCm = param.fMaxRowRangeCm;
328 fIsolCut = param.fIsolCut;
329 fDebugLevel = param.fDebugLevel;
330 fHistoRow = param.fHistoRow ;
331 fHistoPad = param.fHistoPad ;
332 fHistoTime = param.fHistoTime;
333 fHistoRowPad = param.fHistoRowPad;
334 fTimeStamp = param.fTimeStamp;
339 AliTPCclustererKr & AliTPCclustererKr::operator = (const AliTPCclustererKr & param)
342 // assignment operator
344 if (this == ¶m) return (*this);
346 fParam = param.fParam;
347 fRecoParam = param.fRecoParam;
348 fRawData = param.fRawData;
349 fInput = param.fInput ;
350 fOutput = param.fOutput;
351 fDigarr = param.fDigarr;
352 fZeroSup = param.fZeroSup ;
353 fFirstBin = param.fFirstBin ;
354 fLastBin = param.fLastBin ;
355 fMaxNoiseAbs = param.fMaxNoiseAbs ;
356 fMaxNoiseSigma = param.fMaxNoiseSigma ;
357 fMinAdc = param.fMinAdc;
358 fMinTimeBins = param.fMinTimeBins;
359 // fMaxPadRange = param.fMaxPadRange ;
360 // fMaxRowRange = param.fMaxRowRange ;
361 fMaxTimeRange = param.fMaxTimeRange;
362 fValueToSize = param.fValueToSize;
363 fMaxPadRangeCm = param.fMaxPadRangeCm;
364 fMaxRowRangeCm = param.fMaxRowRangeCm;
365 fIsolCut = param.fIsolCut;
366 fDebugLevel = param.fDebugLevel;
367 fHistoRow = param.fHistoRow ;
368 fHistoPad = param.fHistoPad ;
369 fHistoTime = param.fHistoTime;
370 fHistoRowPad = param.fHistoRowPad;
371 fTimeStamp = param.fTimeStamp;
376 AliTPCclustererKr::~AliTPCclustererKr()
384 void AliTPCclustererKr::SetRecoParam(AliTPCRecoParam *recoParam)
387 // set reconstruction parameters
390 fRecoParam = recoParam;
392 //set default parameters if not specified
393 fRecoParam = AliTPCReconstructor::GetRecoParam();
394 if (!fRecoParam) fRecoParam = AliTPCRecoParam::GetLowFluxParam();
400 ////____________________________________________________________________________
402 void AliTPCclustererKr::SetInput(TTree * tree)
405 // set input tree with digits
408 if (!fInput->GetBranch("Segment")){
409 cerr<<"AliTPCclusterKr::FindClusterKr(): no proper input tree !\n";
415 void AliTPCclustererKr::SetOutput(TTree * /*tree*/)
420 fOutput = new TTreeSRedirector("Krypton.root");
423 ////____________________________________________________________________________
425 Int_t AliTPCclustererKr::FinderIO()
427 // Krypton cluster finder for simulated events from MC
430 Error("Digits2Clusters", "input tree not initialised");
435 Error("Digits2Clusters", "output tree not initialised");
445 Int_t AliTPCclustererKr::FinderIO(AliRawReader* rawReader)
447 // Krypton cluster finder for the TPC raw data
448 // this method is unsing AliAltroRawStreamV3
449 // fParam must be defined before
450 if (!rawReader) return 1;
452 fRawData=kTRUE; //set flag to data
455 Error("Digits2Clusters", "output tree not initialised");
459 fParam->SetMaxTBin(fRecoParam->GetLastBin());//set number of timebins from reco -> param
460 // used later for memory allocation
462 AliRawEventHeaderBase* eventHeader = (AliRawEventHeaderBase*)rawReader->GetEventHeader();
464 fTimeStamp = eventHeader->Get("Timestamp");
465 fRun = rawReader->GetRunNumber();
469 Bool_t isAltro=kFALSE;
471 AliTPCROC * roc = AliTPCROC::Instance();
472 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
473 AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
475 AliTPCRawStreamV3 input(rawReader,(AliAltroMapping**)mapping);
477 const Int_t kNIS = fParam->GetNInnerSector();//number of inner sectors
478 const Int_t kNOS = fParam->GetNOuterSector();//number of outer sectors
479 const Int_t kNS = kNIS + kNOS;//all sectors
483 AliTPCDigitsArray *digarr=new AliTPCDigitsArray(kFALSE);//data not sim
484 digarr->Setup(fParam);//as usually parameters
486 for(Int_t iSec = 0; iSec < kNS; iSec++) {
487 AliTPCCalROC * noiseROC;
488 AliTPCCalROC noiseDummy(iSec);
490 noiseROC = &noiseDummy;//noise=0
492 noiseROC = noiseTPC->GetCalROC(iSec); // noise per given sector
494 Int_t nRows = 0; //number of rows in sector
495 Int_t nDDLs = 0; //number of DDLs
496 Int_t indexDDL = 0; //DDL index
498 nRows = fParam->GetNRowLow();
502 nRows = fParam->GetNRowUp();
504 indexDDL = (iSec-kNIS) * 4 + kNIS * 2;
508 // Load the raw data for corresponding DDLs
511 rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
514 while (input.NextDDL()){
515 // Allocate memory for rows in sector (pads(depends on row) x timebins)
516 if (!digarr->GetRow(iSec,0)){
517 for(Int_t iRow = 0; iRow < nRows; iRow++) {
518 digarr->CreateRow(iSec,iRow);
519 }//end loop over rows
522 while ( input.NextChannel() ) {
523 Int_t iRow = input.GetRow();
524 Int_t iPad = input.GetPad();
525 //check row consistency
526 if (iRow < 0 ) continue;
527 if (iRow < 0 || iRow >= nRows){
528 AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",
533 //check pad consistency
534 if (iPad < 0 || iPad >= (Int_t)(roc->GetNPads(iSec,iRow))) {
535 AliError(Form("Pad index (%d) outside the range (%d -> %d) !",
536 iPad, 0, roc->GetNPads(iSec,iRow) ));
541 while ( input.NextBunch() ){
542 Int_t startTbin = (Int_t)input.GetStartTimeBin();
543 Int_t bunchlength = (Int_t)input.GetBunchLength();
544 const UShort_t *sig = input.GetSignals();
546 for (Int_t iTime = 0; iTime<bunchlength; iTime++){
547 Int_t iTimeBin=startTbin-iTime;
550 fHistoRow->Fill(iRow);
551 fHistoPad->Fill(iPad);
552 fHistoTime->Fill(iTimeBin);
553 fHistoRowPad->Fill(iPad,iRow);
554 }else if(fDebugLevel>=0&&fDebugLevel<72){
555 if(iSec==fDebugLevel){
556 fHistoRow->Fill(iRow);
557 fHistoPad->Fill(iPad);
558 fHistoTime->Fill(iTimeBin);
559 fHistoRowPad->Fill(iPad,iRow);
561 }else if(fDebugLevel==73){
563 fHistoRow->Fill(iRow);
564 fHistoPad->Fill(iPad);
565 fHistoTime->Fill(iTimeBin);
566 fHistoRowPad->Fill(iPad,iRow);
568 }else if(fDebugLevel==74){
570 fHistoRow->Fill(iRow);
571 fHistoPad->Fill(iPad);
572 fHistoTime->Fill(iTimeBin);
573 fHistoRowPad->Fill(iPad,iRow);
577 //check time consistency
578 if ( iTimeBin < fRecoParam->GetFirstBin() || iTimeBin >= fRecoParam->GetLastBin()){
579 //cout<<iTimeBin<<endl;
581 AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",
582 iTimeBin, 0, fRecoParam->GetLastBin() -1));
585 Float_t signal=(Float_t)sig[iTime];
586 if (signal <= fZeroSup ||
587 iTimeBin < fFirstBin ||
590 digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
593 if (!noiseROC) continue;
594 Double_t noiseOnPad = noiseROC->GetValue(iRow,iPad);//noise on given pad and row in sector
595 if (noiseOnPad > fMaxNoiseAbs){
596 digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
597 continue; // consider noisy pad as dead
599 if(signal <= fMaxNoiseSigma * noiseOnPad){
600 digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);
603 digarr->GetRow(iSec,iRow)->SetDigitFast(TMath::Nint(signal),iTimeBin,iPad);
604 }// end loop signals in bunch
610 if(isAltro) FindClusterKrIO();
616 void AliTPCclustererKr::CleanSector(Int_t sector){
618 // clean isolated digits
620 const Int_t kNRows=fParam->GetNRow(sector);//number of rows in sector
621 for(Int_t iRow=0; iRow<kNRows; ++iRow){
622 AliSimDigits *digrow;
624 digrow = (AliSimDigits*)fDigarr->GetRow(sector,iRow);//real data
626 digrow = (AliSimDigits*)fDigarr->LoadRow(sector,iRow);//MC
628 if(!digrow) continue;
629 digrow->ExpandBuffer(); //decrunch
630 const Int_t kNPads = digrow->GetNCols(); // number of pads
631 const Int_t kNTime = digrow->GetNRows(); // number of timebins
632 for(Int_t iPad=1;iPad<kNPads-1;iPad++){
633 Short_t* val = digrow->GetDigitsColumn(iPad);
635 for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){
636 if (val[iTimeBin]<=0) continue;
637 if (val[iTimeBin-1]+val[iTimeBin+1]<fIsolCut) {val[iTimeBin]=0; continue;}
638 if (val[iTimeBin-kNTime]+val[iTimeBin+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}
640 if (val[iTimeBin-1-kNTime]+val[iTimeBin+1+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}
641 if (val[iTimeBin+1-kNTime]+val[iTimeBin-1+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}
649 ////____________________________________________________________________________
650 Int_t AliTPCclustererKr::FindClusterKrIO()
654 //fParam and fDigarr must be set to run this method
657 Int_t clusterCounter=0;
658 const Int_t nTotalSector=fParam->GetNSector();//number of sectors
659 for(Int_t iSec=0; iSec<nTotalSector; ++iSec){
662 //vector of maxima for each sector
663 //std::vector<AliPadMax*> maximaInSector;
664 TObjArray *maximaInSector=new TObjArray();//to store AliPadMax*
667 // looking for the maxima on the pad
670 const Int_t kNRows=fParam->GetNRow(iSec);//number of rows in sector
671 for(Int_t iRow=0; iRow<kNRows; ++iRow){
672 AliSimDigits *digrow;
674 digrow = (AliSimDigits*)fDigarr->GetRow(iSec,iRow);//real data
676 digrow = (AliSimDigits*)fDigarr->LoadRow(iSec,iRow);//MC
678 if(digrow){//if pointer exist
679 digrow->ExpandBuffer(); //decrunch
680 const Int_t kNPads = digrow->GetNCols(); // number of pads
681 const Int_t kNTime = digrow->GetNRows(); // number of timebins
682 for(Int_t iPad=0;iPad<kNPads;iPad++){
684 Int_t timeBinMax=-1;//timebin of maximum
685 Int_t valueMaximum=-1;//value of maximum in adc
686 Int_t increaseBegin=-1;//timebin when increase starts
687 Int_t sumAdc=0;//sum of adc on the pad in maximum surrounding
688 bool ifIncreaseBegin=true;//flag - check if increasing started
689 bool ifMaximum=false;//flag - check if it could be maximum
690 Short_t* val = digrow->GetDigitsColumn(iPad);
691 for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){
693 if (val[iTimeBin]==-1) break; // 0 until the end
694 for( ; iTimeBin<kNTime-2&&val[iTimeBin]<fMinAdc ;iTimeBin++) {}
697 Short_t adc = val[iTimeBin];
699 if(adc<fMinAdc){//standard was 3 for fMinAdc
701 if(iTimeBin-increaseBegin<fMinTimeBins){//at least 2 time bins
706 ifIncreaseBegin=true;
710 //insert maximum, default values and set flags
711 //Double_t xCord,yCord;
712 //GetXY(iSec,iRow,iPad,xCord,yCord);
713 Double_t x[]={static_cast<Double_t>(iRow),static_cast<Double_t>(iPad),static_cast<Double_t>(iTimeBin)};
715 AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
717 transform->Transform(x,i,0,1);
719 AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,
729 maximaInSector->AddLast(oneMaximum);
735 ifIncreaseBegin=true;
747 ifIncreaseBegin=false;
748 increaseBegin=iTimeBin;
751 if(adc>valueMaximum){
757 if(iTimeBin==kNTime-1 && ifMaximum && kNTime-increaseBegin>fMinTimeBins){//on the edge
758 //at least 3 timebins
759 //insert maximum, default values and set flags
760 //Double_t xCord,yCord;
761 //GetXY(iSec,iRow,iPad,xCord,yCord);
762 Double_t x[]={static_cast<Double_t>(iRow),static_cast<Double_t>(iPad),static_cast<Double_t>(iTimeBin)};
764 //AliTPCTransform trafo;
765 //trafo.Transform(x,i,0,1);
767 AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
769 transform->Transform(x,i,0,1);
771 AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,
781 maximaInSector->AddLast(oneMaximum);
787 ifIncreaseBegin=true;
792 }//end loop over timebins
793 }//end loop over pads
795 // cout<<"Pointer does not exist!!"<<endl;
796 }//end if poiner exists
797 }//end loop over rows
799 MakeClusters(maximaInSector,iSec,clusterCounter);
801 maximaInSector->SetOwner(kTRUE);
802 maximaInSector->Delete();
803 delete maximaInSector;
805 cout<<"Number of clusters in event: "<<clusterCounter<<endl;
809 void AliTPCclustererKr::MakeClusters(TObjArray * maximaInSector, Int_t iSec, Int_t &clusterCounter){
822 Int_t entriesArr = maximaInSector->GetEntriesFast();
823 for(Int_t it1 = 0; it1 < entriesArr; ++it1 ) {
825 AliPadMax *mp1=(AliPadMax *)maximaInSector->UncheckedAt(it1);
827 AliTPCclusterKr clusterKr;
830 Int_t clusterValue=0;
831 clusterValue+=(mp1)->GetSum();
832 list<Int_t> nUsedRows;
833 nUsedRows.push_back((mp1)->GetRow());
835 maxDig =(mp1)->GetAdc() ;
836 maxSumAdc =(mp1)->GetSum() ;
837 maxTimeBin =(mp1)->GetTime();
838 maxPad =(mp1)->GetPad() ;
839 maxRow =(mp1)->GetRow() ;
844 AliSimDigits *digrowTmp;
846 digrowTmp = (AliSimDigits*)fDigarr->GetRow(iSec,(mp1)->GetRow());
848 digrowTmp = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp1)->GetRow());
851 digrowTmp->ExpandBuffer(); //decrunch
853 for(Int_t itb=(mp1)->GetBegin(); itb<((mp1)->GetEnd())+1; itb++){
854 Int_t adcTmp = digrowTmp->GetDigitUnchecked(itb,(mp1)->GetPad());
855 AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp1)->GetPad(),(mp1)->GetRow(),(mp1)->GetX(),(mp1)->GetY(),(mp1)->GetT());
856 clusterKr.AddDigitToCluster(vtpr);
858 clusterKr.SetCenter();//set centr of the cluster
860 for(Int_t it2 = it1+1; it2 < entriesArr; ++it2 ) {
861 AliPadMax *mp2=(AliPadMax *)maximaInSector->UncheckedAt(it2);
863 if (TMath::Abs(clusterKr.GetCenterX() - (mp2)->GetX()) > fMaxPadRangeCm) continue;
864 if (TMath::Abs(clusterKr.GetCenterY() - (mp2)->GetY()) > fMaxRowRangeCm) continue;
865 if (TMath::Abs(clusterKr.GetCenterT() - (mp2)->GetT()) > fMaxTimeRange) continue;
868 clusterValue+=(mp2)->GetSum();
871 nUsedRows.push_back((mp2)->GetRow());
873 AliSimDigits *digrowTmp1;
875 digrowTmp1 = (AliSimDigits*)fDigarr->GetRow(iSec,(mp2)->GetRow());
877 digrowTmp1 = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp2)->GetRow());
879 digrowTmp1->ExpandBuffer(); //decrunch
881 for(Int_t itb=(mp2)->GetBegin(); itb<(mp2)->GetEnd()+1; itb++){
882 Int_t adcTmp = digrowTmp1->GetDigitUnchecked(itb,(mp2)->GetPad());
883 AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp2)->GetPad(),(mp2)->GetRow(),(mp2)->GetX(),(mp2)->GetY(),(mp2)->GetT());
884 clusterKr.AddDigitToCluster(vtpr);
887 clusterKr.SetCenter();//set center of the cluster
889 //which one is bigger
890 if( (mp2)->GetAdc() > maxDig ){
891 maxDig =(mp2)->GetAdc() ;
892 maxSumAdc =(mp2)->GetSum() ;
893 maxTimeBin =(mp2)->GetTime();
894 maxPad =(mp2)->GetPad() ;
895 maxRow =(mp2)->GetRow() ;
896 maxX =(mp2)->GetX() ;
897 maxY =(mp2)->GetY() ;
898 maxT =(mp2)->GetT() ;
899 } else if ( (mp2)->GetAdc() == maxDig ){
900 if( (mp2)->GetSum() > maxSumAdc){
901 maxDig =(mp2)->GetAdc() ;
902 maxSumAdc =(mp2)->GetSum() ;
903 maxTimeBin =(mp2)->GetTime();
904 maxPad =(mp2)->GetPad() ;
905 maxRow =(mp2)->GetRow() ;
906 maxX =(mp2)->GetX() ;
907 maxY =(mp2)->GetY() ;
908 maxT =(mp2)->GetT() ;
911 delete maximaInSector->RemoveAt(it2);
914 delete maximaInSector->RemoveAt(it1);
916 //through out clusters on the edge and noise
917 //if(clusterValue/clusterKr.fCluster.size()<fValueToSize)continue;
918 if(clusterValue/(clusterKr.GetSize())<fValueToSize)continue;
920 clusterKr.SetADCcluster(clusterValue);
921 clusterKr.SetNPads(nUsedPads);
922 clusterKr.SetMax(AliTPCvtpr(maxDig,maxTimeBin,maxPad,maxRow,maxX,maxY,maxT));
923 clusterKr.SetSec(iSec);
928 clusterKr.SetNRows(nUsedRows.size());
929 clusterKr.SetCenter();
931 clusterKr.SetRMS();//Set pad,row,timebin RMS
932 clusterKr.Set1D();//Set size in pads and timebins
934 clusterKr.SetTimeStamp(fTimeStamp);
935 clusterKr.SetRun(fRun);
940 //save each cluster into file
946 //end of save each cluster into file adc.root
952 ////____________________________________________________________________________
955 void AliTPCclustererKr::GetXY(Int_t sec,Int_t row,Int_t pad,Double_t& xGlob,Double_t& yGlob){
957 //gives global XY coordinate of the pad
960 Double_t yLocal = fParam->GetPadRowRadii(sec,row);//radius of row in sector in cm
962 Int_t padmax = fParam->GetNPads(sec,row);//number of pads in a given row
964 if(sec<fParam->GetNInnerSector())padXSize=0.4;
966 Double_t xLocal=(pad+0.5-padmax/2.)*padXSize;//x-value of the center of pad
969 fParam->AdjustCosSin((Int_t)sec,cos,sin);//return sinus and cosinus of the sector
971 Double_t xGlob1 = xLocal * cos + yLocal * sin;
972 Double_t yGlob1 = -xLocal * sin + yLocal * cos;
978 xGlob = xGlob1 * TMath::Cos(rot) + yGlob1 * TMath::Sin(rot);
979 yGlob = -xGlob1 * TMath::Sin(rot) + yGlob1 * TMath::Cos(rot);
982 if(sec<18||(sec>=36&&sec<54)) xGlob =-1*xGlob;