For each data subsample the same named file is created. So be careful \r
do not overwrite them. \r
\r
+Additional selection criteria to select the GOLD cluster\r
+Example:\r
+// open file with clusters\r
+TFile f("Krypton.root");\r
+TTree * tree = (TTree*)f.Get("Kr")\r
+TCut cutR0("cutR0","fADCcluster/fSize<100"); // adjust it according v seetings - \r
+TCut cutR1("cutR1","fADCcluster/fSize>7"); // cosmic tracks and noise removal\r
+TCut cutR2("cutR2","fMax.fAdc/fADCcluster<0.2"); // digital noise removal\r
+TCut cutR3("cutR3","fMax.fAdc/fADCcluster>0.01"); // noise removal\r
+TCut cutS1("cutS1","fSize<200"); // adjust it according v seetings - cosmic tracks\r
+TCut cutAll = cutR0+cutR1+cutR2+cutR3+cutS1;\r
+This values are typical values to be applied in selectors\r
+\r
+\r
*\r
**** MC ****\r
*\r
#include "TTree.h"\r
#include "TH1F.h"\r
#include "TH2F.h"\r
+#include "TTreeStream.h"\r
+\r
+#include "AliTPCTransform.h"\r
\r
//used in raw data finder\r
#include "AliTPCROC.h"\r
#include "AliTPCAltroMapping.h"\r
#include "AliTPCcalibDB.h"\r
#include "AliTPCRawStream.h"\r
+#include "AliTPCRawStreamV3.h"\r
#include "AliTPCRecoParam.h"\r
#include "AliTPCReconstructor.h"\r
#include "AliRawReader.h"\r
#include "AliTPCCalROC.h"\r
+#include "AliRawEventHeaderBase.h"\r
\r
ClassImp(AliTPCclustererKr)\r
\r
AliTPCclustererKr::AliTPCclustererKr()\r
:TObject(),\r
fRawData(kFALSE),\r
- fRowCl(0),\r
fInput(0),\r
fOutput(0),\r
fParam(0),\r
fValueToSize(3.5),\r
fMaxPadRangeCm(2.5),\r
fMaxRowRangeCm(3.5),\r
+ fIsolCut(3),\r
fDebugLevel(-1),\r
fHistoRow(0),\r
fHistoPad(0),\r
fHistoTime(0),\r
- fHistoRowPad(0)\r
+ fHistoRowPad(0),\r
+ fTimeStamp(0),\r
+ fRun(0)\r
{\r
//\r
// default constructor\r
AliTPCclustererKr::AliTPCclustererKr(const AliTPCclustererKr ¶m)\r
:TObject(),\r
fRawData(kFALSE),\r
- fRowCl(0),\r
fInput(0),\r
fOutput(0),\r
fParam(0),\r
fValueToSize(3.5),\r
fMaxPadRangeCm(2.5),\r
fMaxRowRangeCm(3.5),\r
+ fIsolCut(3),\r
fDebugLevel(-1),\r
fHistoRow(0),\r
fHistoPad(0),\r
fHistoTime(0),\r
- fHistoRowPad(0)\r
+ fHistoRowPad(0),\r
+ fTimeStamp(0),\r
+ fRun(0)\r
{\r
//\r
// copy constructor\r
fParam = param.fParam;\r
fRecoParam = param.fRecoParam;\r
fRawData = param.fRawData;\r
- fRowCl = param.fRowCl ;\r
fInput = param.fInput ;\r
fOutput = param.fOutput;\r
fDigarr = param.fDigarr;\r
fValueToSize = param.fValueToSize;\r
fMaxPadRangeCm = param.fMaxPadRangeCm;\r
fMaxRowRangeCm = param.fMaxRowRangeCm;\r
+ fIsolCut = param.fIsolCut;\r
fDebugLevel = param.fDebugLevel;\r
fHistoRow = param.fHistoRow ;\r
fHistoPad = param.fHistoPad ;\r
fHistoTime = param.fHistoTime;\r
fHistoRowPad = param.fHistoRowPad;\r
+ fTimeStamp = param.fTimeStamp;\r
+ fRun = param.fRun;\r
\r
} \r
\r
fParam = param.fParam;\r
fRecoParam = param.fRecoParam;\r
fRawData = param.fRawData;\r
- fRowCl = param.fRowCl ;\r
fInput = param.fInput ;\r
fOutput = param.fOutput;\r
fDigarr = param.fDigarr;\r
fValueToSize = param.fValueToSize;\r
fMaxPadRangeCm = param.fMaxPadRangeCm;\r
fMaxRowRangeCm = param.fMaxRowRangeCm;\r
+ fIsolCut = param.fIsolCut;\r
fDebugLevel = param.fDebugLevel;\r
fHistoRow = param.fHistoRow ;\r
fHistoPad = param.fHistoPad ;\r
fHistoTime = param.fHistoTime;\r
fHistoRowPad = param.fHistoRowPad;\r
+ fTimeStamp = param.fTimeStamp;\r
+ fRun = param.fRun;\r
return (*this);\r
}\r
\r
//\r
// destructor\r
//\r
+ delete fOutput;\r
}\r
\r
void AliTPCclustererKr::SetRecoParam(AliTPCRecoParam *recoParam)\r
}\r
}\r
\r
-void AliTPCclustererKr::SetOutput(TTree * tree) \r
+void AliTPCclustererKr::SetOutput(TTree * /*tree*/) \r
{\r
//\r
- // set output\r
+ //dummy\r
//\r
- fOutput= tree; \r
- AliTPCClustersRow clrow;\r
- AliTPCClustersRow *pclrow=&clrow; \r
- clrow.SetClass("AliTPCclusterKr");\r
- clrow.SetArray(1); // to make Clones array\r
- fOutput->Branch("Segment","AliTPCClustersRow",&pclrow,32000,200); \r
+ fOutput = new TTreeSRedirector("Krypton.root");\r
}\r
\r
////____________________________________________________________________________\r
return 0;\r
}\r
\r
+\r
+\r
Int_t AliTPCclustererKr::FinderIO(AliRawReader* rawReader)\r
+{\r
+ // Krypton cluster finder for the TPC raw data\r
+ // this method is unsing AliAltroRawStreamV3\r
+ // fParam must be defined before\r
+ if (!rawReader) return 1;\r
+ //\r
+ fRawData=kTRUE; //set flag to data\r
+ \r
+ if (!fOutput) {\r
+ Error("Digits2Clusters", "output tree not initialised");\r
+ return 11;\r
+ }\r
+ \r
+ fParam->SetMaxTBin(fRecoParam->GetLastBin());//set number of timebins from reco -> param\r
+ // used later for memory allocation\r
+\r
+ AliRawEventHeaderBase* eventHeader = (AliRawEventHeaderBase*)rawReader->GetEventHeader();\r
+ if (eventHeader){\r
+ fTimeStamp = eventHeader->Get("Timestamp");\r
+ fRun = rawReader->GetRunNumber();\r
+ }\r
+\r
+\r
+ Bool_t isAltro=kFALSE;\r
+ \r
+ AliTPCROC * roc = AliTPCROC::Instance();\r
+ AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();\r
+ AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();\r
+ //\r
+ AliTPCRawStreamV3 input(rawReader,(AliAltroMapping**)mapping);\r
+ \r
+ const Int_t kNIS = fParam->GetNInnerSector();//number of inner sectors\r
+ const Int_t kNOS = fParam->GetNOuterSector();//number of outer sectors\r
+ const Int_t kNS = kNIS + kNOS;//all sectors\r
+ \r
+ \r
+ //crate TPC view\r
+ AliTPCDigitsArray *digarr=new AliTPCDigitsArray(kFALSE);//data not sim\r
+ digarr->Setup(fParam);//as usually parameters\r
+ \r
+ for(Int_t iSec = 0; iSec < kNS; iSec++) {\r
+ AliTPCCalROC * noiseROC;\r
+ AliTPCCalROC noiseDummy(iSec);\r
+ if(noiseTPC==0x0){\r
+ noiseROC = &noiseDummy;//noise=0\r
+ }else{\r
+ noiseROC = noiseTPC->GetCalROC(iSec); // noise per given sector\r
+ }\r
+ Int_t nRows = 0; //number of rows in sector\r
+ Int_t nDDLs = 0; //number of DDLs\r
+ Int_t indexDDL = 0; //DDL index\r
+ if (iSec < kNIS) {\r
+ nRows = fParam->GetNRowLow();\r
+ nDDLs = 2;\r
+ indexDDL = iSec * 2;\r
+ }else {\r
+ nRows = fParam->GetNRowUp();\r
+ nDDLs = 4;\r
+ indexDDL = (iSec-kNIS) * 4 + kNIS * 2;\r
+ }\r
+ \r
+ //\r
+ // Load the raw data for corresponding DDLs\r
+ //\r
+ rawReader->Reset();\r
+ rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);\r
+ \r
+ \r
+ while (input.NextDDL()){\r
+ // Allocate memory for rows in sector (pads(depends on row) x timebins)\r
+ if (!digarr->GetRow(iSec,0)){\r
+ for(Int_t iRow = 0; iRow < nRows; iRow++) {\r
+ digarr->CreateRow(iSec,iRow);\r
+ }//end loop over rows\r
+ }\r
+ //loop over pads\r
+ while ( input.NextChannel() ) {\r
+ Int_t iRow = input.GetRow();\r
+ Int_t iPad = input.GetPad();\r
+ //check row consistency\r
+ if (iRow < 0 ) continue;\r
+ if (iRow < 0 || iRow >= nRows){\r
+ AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",\r
+ iRow, 0, nRows -1));\r
+ continue;\r
+ }\r
+ \r
+ //check pad consistency\r
+ if (iPad < 0 || iPad >= (Int_t)(roc->GetNPads(iSec,iRow))) {\r
+ AliError(Form("Pad index (%d) outside the range (%d -> %d) !",\r
+ iPad, 0, roc->GetNPads(iSec,iRow) ));\r
+ continue;\r
+ }\r
+ \r
+ //loop over bunches\r
+ while ( input.NextBunch() ){\r
+ Int_t startTbin = (Int_t)input.GetStartTimeBin();\r
+ Int_t bunchlength = (Int_t)input.GetBunchLength();\r
+ const UShort_t *sig = input.GetSignals();\r
+ isAltro=kTRUE;\r
+ for (Int_t iTime = 0; iTime<bunchlength; iTime++){\r
+ Int_t iTimeBin=startTbin-iTime;\r
+ //\r
+ if(fDebugLevel==72){\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
+ }else if(fDebugLevel>=0&&fDebugLevel<72){\r
+ if(iSec==fDebugLevel){\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
+ }\r
+ }else if(fDebugLevel==73){\r
+ if(iSec<36){\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
+ }\r
+ }else if(fDebugLevel==74){\r
+ if(iSec>=36){\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
+ }\r
+ }\r
+ \r
+ //check time consistency\r
+ if ( iTimeBin < fRecoParam->GetFirstBin() || iTimeBin >= fRecoParam->GetLastBin()){\r
+ //cout<<iTimeBin<<endl;\r
+ continue;\r
+ AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",\r
+ iTimeBin, 0, fRecoParam->GetLastBin() -1));\r
+ }\r
+ //signal\r
+ Float_t signal=(Float_t)sig[iTime];\r
+ if (signal <= fZeroSup ||\r
+ iTimeBin < fFirstBin ||\r
+ iTimeBin > fLastBin\r
+ ) {\r
+ digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+ continue;\r
+ }\r
+ if (!noiseROC) continue;\r
+ Double_t noiseOnPad = noiseROC->GetValue(iRow,iPad);//noise on given pad and row in sector\r
+ if (noiseOnPad > fMaxNoiseAbs){\r
+ digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+ continue; // consider noisy pad as dead\r
+ }\r
+ if(signal <= fMaxNoiseSigma * noiseOnPad){\r
+ digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+ continue;\r
+ }\r
+ digarr->GetRow(iSec,iRow)->SetDigitFast(TMath::Nint(signal),iTimeBin,iPad);\r
+ }// end loop signals in bunch\r
+ }// end loop bunches\r
+ } // end loop pads\r
+ }// end ddl loop\r
+ }// end sector loop\r
+ SetDigArr(digarr);\r
+ if(isAltro) FindClusterKrIO();\r
+ delete digarr;\r
+ \r
+ return 0;\r
+}\r
+\r
+\r
+\r
+\r
+\r
+Int_t AliTPCclustererKr::FinderIOold(AliRawReader* rawReader)\r
{\r
// Krypton cluster finder for the TPC raw data\r
//\r
// fParam must be defined before\r
+ if (!rawReader) return 1;\r
\r
if(rawReader)fRawData=kTRUE; //set flag to data\r
-\r
+ \r
if (!fOutput) {\r
Error("Digits2Clusters", "output tree not initialised");\r
return 11;\r
}\r
-\r
+ \r
fParam->SetMaxTBin(fRecoParam->GetLastBin());//set number of timebins from reco -> param\r
// used later for memory allocation\r
\r
+ AliRawEventHeaderBase* eventHeader = (AliRawEventHeaderBase*)rawReader->GetEventHeader();\r
+ if (eventHeader){\r
+ fTimeStamp = eventHeader->Get("Timestamp");\r
+ fRun = rawReader->GetRunNumber();\r
+ }\r
+ \r
Bool_t isAltro=kFALSE;\r
-\r
+ \r
AliTPCROC * roc = AliTPCROC::Instance();\r
AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();\r
AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();\r
//\r
AliTPCRawStream input(rawReader,(AliAltroMapping**)mapping);\r
-\r
+ \r
const Int_t kNIS = fParam->GetNInnerSector();//number of inner sectors\r
const Int_t kNOS = fParam->GetNOuterSector();//number of outer sectors\r
const Int_t kNS = kNIS + kNOS;//all sectors\r
-\r
-// //set cluster finder parameters\r
-// SetZeroSup(fParam->GetZeroSup());//zero suppression parameter\r
-// SetFirstBin(60);\r
-// SetLastBin(950);\r
-// SetMaxNoiseAbs(2);\r
-// SetMaxNoiseSigma(3);\r
-\r
+ \r
+ \r
//crate TPC view\r
AliTPCDigitsArray *digarr=new AliTPCDigitsArray(kFALSE);//data not sim\r
digarr->Setup(fParam);//as usually parameters\r
-\r
+ \r
//\r
// Loop over sectors\r
//\r
nDDLs = 4;\r
indexDDL = (iSec-kNIS) * 4 + kNIS * 2;\r
}\r
-\r
+ \r
//\r
// Load the raw data for corresponding DDLs\r
//\r
rawReader->Reset();\r
rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);\r
-\r
+ \r
if(input.Next()) {\r
isAltro=kTRUE;\r
// Allocate memory for rows in sector (pads(depends on row) x timebins)\r
for(Int_t iRow = 0; iRow < nRows; iRow++) {\r
- digarr->CreateRow(iSec,iRow);\r
+ digarr->CreateRow(iSec,iRow);\r
}//end loop over rows\r
}\r
rawReader->Reset();\r
rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);\r
-\r
+ \r
//\r
// Begin loop over altro data\r
//\r
while (input.Next()) {\r
-\r
+ \r
//check sector consistency\r
if (input.GetSector() != iSec)\r
- AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSec,input.GetSector()));\r
+ AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSec,input.GetSector()));\r
\r
- Short_t iRow = input.GetRow();\r
- Short_t iPad = input.GetPad();\r
- Short_t iTimeBin = input.GetTime();\r
-\r
+ Int_t iRow = input.GetRow();\r
+ Int_t iPad = input.GetPad();\r
+ Int_t iTimeBin = input.GetTime();\r
+ \r
+ //\r
if(fDebugLevel==72){\r
- fHistoRow->Fill(iRow);\r
- fHistoPad->Fill(iPad);\r
- fHistoTime->Fill(iTimeBin);\r
- fHistoRowPad->Fill(iPad,iRow);\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
}else if(fDebugLevel>=0&&fDebugLevel<72){\r
- if(iSec==fDebugLevel){\r
- fHistoRow->Fill(iRow);\r
- fHistoPad->Fill(iPad);\r
- fHistoTime->Fill(iTimeBin);\r
- fHistoRowPad->Fill(iPad,iRow);\r
- }\r
+ if(iSec==fDebugLevel){\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
+ }\r
}else if(fDebugLevel==73){\r
- if(iSec<36){\r
- fHistoRow->Fill(iRow);\r
- fHistoPad->Fill(iPad);\r
- fHistoTime->Fill(iTimeBin);\r
- fHistoRowPad->Fill(iPad,iRow);\r
- }\r
+ if(iSec<36){\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
+ }\r
}else if(fDebugLevel==74){\r
- if(iSec>=36){\r
- fHistoRow->Fill(iRow);\r
- fHistoPad->Fill(iPad);\r
- fHistoTime->Fill(iTimeBin);\r
- fHistoRowPad->Fill(iPad,iRow);\r
- }\r
+ if(iSec>=36){\r
+ fHistoRow->Fill(iRow);\r
+ fHistoPad->Fill(iPad);\r
+ fHistoTime->Fill(iTimeBin);\r
+ fHistoRowPad->Fill(iPad,iRow);\r
+ }\r
}\r
-\r
+ \r
//check row consistency\r
+ if (iRow < 0 ) continue;\r
if (iRow < 0 || iRow >= nRows){\r
- AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",\r
- iRow, 0, nRows -1));\r
- continue;\r
+ AliError(Form("Pad-row index (%d) outside the range (%d -> %d) !",\r
+ iRow, 0, nRows -1));\r
+ continue;\r
}\r
-\r
+ \r
//check pad consistency\r
- if (iPad < 0 || iPad >= (Short_t)(roc->GetNPads(iSec,iRow))) {\r
- AliError(Form("Pad index (%d) outside the range (%d -> %d) !",\r
- iPad, 0, roc->GetNPads(iSec,iRow) ));\r
- continue;\r
+ if (iPad < 0 || iPad >= (Int_t)(roc->GetNPads(iSec,iRow))) {\r
+ AliError(Form("Pad index (%d) outside the range (%d -> %d) !",\r
+ iPad, 0, roc->GetNPads(iSec,iRow) ));\r
+ continue;\r
}\r
-\r
+ \r
//check time consistency\r
if ( iTimeBin < fRecoParam->GetFirstBin() || iTimeBin >= fRecoParam->GetLastBin()){\r
- //cout<<iTimeBin<<endl;\r
- continue;\r
- AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",\r
- iTimeBin, 0, fRecoParam->GetLastBin() -1));\r
+ //cout<<iTimeBin<<endl;\r
+ continue;\r
+ AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",\r
+ iTimeBin, 0, fRecoParam->GetLastBin() -1));\r
}\r
-\r
+ \r
//signal\r
Int_t signal = input.GetSignal();\r
- if (signal <= fZeroSup || \r
- iTimeBin < fFirstBin ||\r
- iTimeBin > fLastBin\r
- ) {\r
- digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
- continue;\r
- }\r
+ if (signal <= fZeroSup ||\r
+ iTimeBin < fFirstBin ||\r
+ iTimeBin > fLastBin\r
+ ) {\r
+ digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+ continue;\r
+ }\r
if (!noiseROC) continue;\r
Double_t noiseOnPad = noiseROC->GetValue(iRow,iPad);//noise on given pad and row in sector\r
if (noiseOnPad > fMaxNoiseAbs){\r
- digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
- continue; // consider noisy pad as dead\r
+ digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+ continue; // consider noisy pad as dead\r
}\r
if(signal <= fMaxNoiseSigma * noiseOnPad){\r
- digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
- continue;\r
+ digarr->GetRow(iSec,iRow)->SetDigitFast(0,iTimeBin,iPad);\r
+ continue;\r
}\r
-\r
digarr->GetRow(iSec,iRow)->SetDigitFast(signal,iTimeBin,iPad);\r
-\r
}//end of loop over altro data\r
}//end of loop over sectors\r
\r
SetDigArr(digarr);\r
if(isAltro) FindClusterKrIO();\r
delete digarr;\r
-\r
+ \r
return 0;\r
}\r
\r
+void AliTPCclustererKr::CleanSector(Int_t sector){\r
+ //\r
+ // clean isolated digits\r
+ // \r
+ const Int_t kNRows=fParam->GetNRow(sector);//number of rows in sector\r
+ for(Int_t iRow=0; iRow<kNRows; ++iRow){\r
+ AliSimDigits *digrow;\r
+ if(fRawData){\r
+ digrow = (AliSimDigits*)fDigarr->GetRow(sector,iRow);//real data\r
+ }else{\r
+ digrow = (AliSimDigits*)fDigarr->LoadRow(sector,iRow);//MC\r
+ }\r
+ if(!digrow) continue;\r
+ digrow->ExpandBuffer(); //decrunch\r
+ const Int_t kNPads = digrow->GetNCols(); // number of pads\r
+ const Int_t kNTime = digrow->GetNRows(); // number of timebins\r
+ for(Int_t iPad=1;iPad<kNPads-1;iPad++){\r
+ Short_t* val = digrow->GetDigitsColumn(iPad);\r
+\r
+ for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){\r
+ if (val[iTimeBin]<=0) continue;\r
+ if (val[iTimeBin-1]+val[iTimeBin+1]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+ if (val[iTimeBin-kNTime]+val[iTimeBin+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+ //\r
+ if (val[iTimeBin-1-kNTime]+val[iTimeBin+1+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+ if (val[iTimeBin+1-kNTime]+val[iTimeBin-1+kNTime]<fIsolCut) {val[iTimeBin]=0; continue;}\r
+\r
+ }\r
+ }\r
+ }\r
+}\r
+\r
+\r
////____________________________________________________________________________\r
Int_t AliTPCclustererKr::FindClusterKrIO()\r
{\r
- //fParam and fDigarr must be set to run this method\r
\r
-// //set parameters \r
-// SetMinAdc(3);//usually is 3 \r
-// SetMinTimeBins(2);//should be 2 - the best result of shape in MC\r
-//// SetMaxPadRange(4);\r
-//// SetMaxRowRange(3);\r
-// SetMaxTimeRange(7);\r
-// SetValueToSize(3.5);//3.5\r
-// SetMaxPadRangeCm(2.5);\r
-// SetMaxRowRangeCm(3.5);\r
+ //\r
+ //fParam and fDigarr must be set to run this method\r
+ //\r
\r
Int_t clusterCounter=0;\r
- const Short_t nTotalSector=fParam->GetNSector();//number of sectors\r
- for(Short_t iSec=0; iSec<nTotalSector; ++iSec){\r
- \r
+ const Int_t nTotalSector=fParam->GetNSector();//number of sectors\r
+ for(Int_t iSec=0; iSec<nTotalSector; ++iSec){\r
+ CleanSector(iSec);\r
+\r
//vector of maxima for each sector\r
//std::vector<AliPadMax*> maximaInSector;\r
TObjArray *maximaInSector=new TObjArray();//to store AliPadMax*\r
// looking for the maxima on the pad\r
//\r
\r
- const Short_t kNRows=fParam->GetNRow(iSec);//number of rows in sector\r
- for(Short_t iRow=0; iRow<kNRows; ++iRow){\r
+ const Int_t kNRows=fParam->GetNRow(iSec);//number of rows in sector\r
+ for(Int_t iRow=0; iRow<kNRows; ++iRow){\r
AliSimDigits *digrow;\r
if(fRawData){\r
digrow = (AliSimDigits*)fDigarr->GetRow(iSec,iRow);//real data\r
}\r
if(digrow){//if pointer exist\r
digrow->ExpandBuffer(); //decrunch\r
- const Short_t kNPads = digrow->GetNCols(); // number of pads\r
- const Short_t kNTime = digrow->GetNRows(); // number of timebins\r
- for(Short_t iPad=0;iPad<kNPads;iPad++){\r
+ const Int_t kNPads = digrow->GetNCols(); // number of pads\r
+ const Int_t kNTime = digrow->GetNRows(); // number of timebins\r
+ for(Int_t iPad=0;iPad<kNPads;iPad++){\r
\r
- Short_t timeBinMax=-1;//timebin of maximum \r
- Short_t valueMaximum=-1;//value of maximum in adc\r
- Short_t increaseBegin=-1;//timebin when increase starts\r
- Short_t sumAdc=0;//sum of adc on the pad in maximum surrounding\r
+ Int_t timeBinMax=-1;//timebin of maximum \r
+ Int_t valueMaximum=-1;//value of maximum in adc\r
+ Int_t increaseBegin=-1;//timebin when increase starts\r
+ Int_t sumAdc=0;//sum of adc on the pad in maximum surrounding\r
bool ifIncreaseBegin=true;//flag - check if increasing started\r
bool ifMaximum=false;//flag - check if it could be maximum\r
+ Short_t* val = digrow->GetDigitsColumn(iPad);\r
+ for(Int_t iTimeBin=1;iTimeBin<kNTime-1;iTimeBin++){\r
+ if (!ifMaximum) {\r
+ if (val[iTimeBin]==-1) break; // 0 until the end\r
+ for( ; iTimeBin<kNTime-2&&val[iTimeBin]<fMinAdc ;iTimeBin++) {}\r
+ }\r
+ //\r
+ Short_t adc = val[iTimeBin];\r
\r
- for(Short_t iTimeBin=0;iTimeBin<kNTime;iTimeBin++){\r
- Short_t adc = digrow->GetDigitFast(iTimeBin,iPad);\r
- if(adc<fMinAdc){//standard was 3\r
+ if(adc<fMinAdc){//standard was 3 for fMinAdc\r
if(ifMaximum){\r
if(iTimeBin-increaseBegin<fMinTimeBins){//at least 2 time bins\r
timeBinMax=-1;\r
continue;\r
}\r
//insert maximum, default values and set flags\r
- Double_t xCord,yCord;\r
- GetXY(iSec,iRow,iPad,xCord,yCord);\r
+ //Double_t xCord,yCord;\r
+ //GetXY(iSec,iRow,iPad,xCord,yCord);\r
+ Double_t x[]={iRow,iPad,iTimeBin};\r
+ Int_t i[]={iSec};\r
+ AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;\r
+\r
+ transform->Transform(x,i,0,1);\r
+ \r
AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,\r
timeBinMax,\r
iPad,\r
iRow,\r
- xCord,\r
- yCord,\r
- timeBinMax),\r
+ x[0],//xCord,\r
+ x[1],//yCord,\r
+ x[2]/*timeBinMax*/),\r
increaseBegin,\r
iTimeBin-1,\r
sumAdc);\r
- //maximaInSector.push_back(oneMaximum);\r
maximaInSector->AddLast(oneMaximum);\r
\r
timeBinMax=-1;\r
}\r
continue;\r
}\r
- \r
+\r
+\r
+\r
+\r
+\r
+\r
if(ifIncreaseBegin){\r
ifIncreaseBegin=false;\r
increaseBegin=iTimeBin;\r
if(iTimeBin==kNTime-1 && ifMaximum && kNTime-increaseBegin>fMinTimeBins){//on the edge\r
//at least 3 timebins\r
//insert maximum, default values and set flags\r
- Double_t xCord,yCord;\r
- GetXY(iSec,iRow,iPad,xCord,yCord);\r
+ //Double_t xCord,yCord;\r
+ //GetXY(iSec,iRow,iPad,xCord,yCord);\r
+ Double_t x[]={iRow,iPad,iTimeBin};\r
+ Int_t i[]={iSec};\r
+ //AliTPCTransform trafo;\r
+ //trafo.Transform(x,i,0,1);\r
+\r
+ AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;\r
+\r
+ transform->Transform(x,i,0,1);\r
+\r
AliPadMax *oneMaximum = new AliPadMax(AliTPCvtpr(valueMaximum,\r
timeBinMax,\r
iPad,\r
iRow,\r
- xCord,\r
- yCord,\r
- timeBinMax),\r
+ x[0],//xCord,\r
+ x[1],//yCord,\r
+ x[2]/*timeBinMax*/),\r
increaseBegin,\r
iTimeBin-1,\r
sumAdc);\r
- //maximaInSector.push_back(oneMaximum);\r
maximaInSector->AddLast(oneMaximum);\r
\r
timeBinMax=-1;\r
}//end if poiner exists\r
}//end loop over rows\r
\r
- // cout<<"EF" <<maximaInSector->GetEntriesFast()<<" E "<<maximaInSector->GetEntries()<<" "<<maximaInSector->GetLast()<<endl;\r
- maximaInSector->Compress();\r
- // GetEntriesFast() - liczba wejsc w array of maxima\r
- //cout<<"EF" <<maximaInSector->GetEntriesFast()<<" E"<<maximaInSector->GetEntries()<<endl;\r
-\r
- //\r
- // Making clusters\r
+ MakeClusters(maximaInSector,iSec,clusterCounter);\r
//\r
+ maximaInSector->SetOwner(kTRUE);\r
+ maximaInSector->Delete();\r
+ delete maximaInSector;\r
+ }//end sector for\r
+ cout<<"Number of clusters in event: "<<clusterCounter<<endl;\r
+ return 0;\r
+}\r
\r
- Short_t maxDig=0;\r
- Short_t maxSumAdc=0;\r
- Short_t maxTimeBin=0;\r
- Short_t maxPad=0;\r
- Short_t maxRow=0;\r
- Double_t maxX=0;\r
- Double_t maxY=0;\r
- \r
-// for( std::vector<AliPadMax*>::iterator mp1 = maximaInSector.begin();\r
-// mp1 != maximaInSector.end(); ++mp1 ) {\r
- for(Int_t it1 = 0; it1 < maximaInSector->GetEntriesFast(); ++it1 ) {\r
-\r
- AliPadMax *mp1=(AliPadMax *)maximaInSector->At(it1);\r
- AliTPCclusterKr *tmp=new AliTPCclusterKr();\r
- \r
- Short_t nUsedPads=1;\r
- Int_t clusterValue=0;\r
- clusterValue+=(mp1)->GetSum();\r
- list<Short_t> nUsedRows;\r
- nUsedRows.push_back((mp1)->GetRow());\r
-\r
- maxDig =(mp1)->GetAdc() ;\r
- maxSumAdc =(mp1)->GetSum() ;\r
- maxTimeBin =(mp1)->GetTime();\r
- maxPad =(mp1)->GetPad() ;\r
- maxRow =(mp1)->GetRow() ;\r
- maxX =(mp1)->GetX();\r
- maxY =(mp1)->GetY();\r
-\r
- AliSimDigits *digrowTmp;\r
- if(fRawData){\r
- digrowTmp = (AliSimDigits*)fDigarr->GetRow(iSec,(mp1)->GetRow());\r
- }else{\r
- digrowTmp = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp1)->GetRow());\r
- }\r
-\r
- digrowTmp->ExpandBuffer(); //decrunch\r
-\r
- for(Short_t itb=(mp1)->GetBegin(); itb<((mp1)->GetEnd())+1; itb++){\r
- Short_t adcTmp = digrowTmp->GetDigitFast(itb,(mp1)->GetPad());\r
- AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp1)->GetPad(),(mp1)->GetRow(),(mp1)->GetX(),(mp1)->GetY(),itb);\r
- //tmp->fCluster.push_back(vtpr);\r
- tmp->AddDigitToCluster(vtpr);\r
- }\r
- tmp->SetCenter();//set centr of the cluster\r
- \r
- //maximaInSector.erase(mp1);\r
- //mp1--;\r
- //cout<<maximaInSector->GetEntriesFast()<<" ";\r
- maximaInSector->RemoveAt(it1);\r
- maximaInSector->Compress();\r
- it1--;\r
- // cout<<maximaInSector->GetEntriesFast()<<" "<<endl;\r
-\r
-// for( std::vector<AliPadMax*>::iterator mp2 = maximaInSector.begin();\r
-// mp2 != maximaInSector.end(); ++mp2 ) {\r
- for(Int_t it2 = 0; it2 < maximaInSector->GetEntriesFast(); ++it2 ) {\r
- AliPadMax *mp2=(AliPadMax *)maximaInSector->At(it2);\r
-\r
-// if(abs(maxRow - (*mp2)->GetRow()) < fMaxRowRange && //3\r
-// abs(maxPad - (*mp2)->GetPad()) < fMaxPadRange && //4\r
- \r
- if(TMath::Abs(tmp->GetCenterX() - (mp2)->GetX()) < fMaxPadRangeCm &&\r
- TMath::Abs(tmp->GetCenterY() - (mp2)->GetY()) < fMaxRowRangeCm &&\r
- TMath::Abs(tmp->GetCenterT() - (mp2)->GetT()) < fMaxTimeRange){//7\r
- \r
- clusterValue+=(mp2)->GetSum();\r
-\r
- nUsedPads++;\r
- nUsedRows.push_back((mp2)->GetRow());\r
-\r
- AliSimDigits *digrowTmp1;\r
- if(fRawData){\r
- digrowTmp1 = (AliSimDigits*)fDigarr->GetRow(iSec,(mp2)->GetRow());\r
- }else{\r
- digrowTmp1 = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp2)->GetRow());\r
- }\r
- digrowTmp1->ExpandBuffer(); //decrunch\r
- \r
- for(Short_t itb=(mp2)->GetBegin(); itb<(mp2)->GetEnd()+1; itb++){\r
- Short_t adcTmp = digrowTmp1->GetDigitFast(itb,(mp2)->GetPad());\r
- AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp2)->GetPad(),(mp2)->GetRow(),(mp2)->GetX(),(mp2)->GetY(),itb);\r
- //tmp->fCluster.push_back(vtpr);\r
- tmp->AddDigitToCluster(vtpr);\r
- }\r
- \r
- tmp->SetCenter();//set center of the cluster\r
+void AliTPCclustererKr::MakeClusters(TObjArray * maximaInSector, Int_t iSec, Int_t &clusterCounter){\r
+ //\r
+ // Make clusters\r
+ //\r
\r
- //which one is bigger\r
- if( (mp2)->GetAdc() > maxDig ){\r
+ Int_t maxDig=0;\r
+ Int_t maxSumAdc=0;\r
+ Int_t maxTimeBin=0;\r
+ Int_t maxPad=0;\r
+ Int_t maxRow=0;\r
+ Double_t maxX=0;\r
+ Double_t maxY=0;\r
+ Double_t maxT=0;\r
+ Int_t entriesArr = maximaInSector->GetEntriesFast();\r
+ for(Int_t it1 = 0; it1 < entriesArr; ++it1 ) {\r
+ \r
+ AliPadMax *mp1=(AliPadMax *)maximaInSector->UncheckedAt(it1);\r
+ if (!mp1) continue;\r
+ AliTPCclusterKr clusterKr;\r
+ \r
+ Int_t nUsedPads=1;\r
+ Int_t clusterValue=0;\r
+ clusterValue+=(mp1)->GetSum();\r
+ list<Int_t> nUsedRows;\r
+ nUsedRows.push_back((mp1)->GetRow());\r
+ \r
+ maxDig =(mp1)->GetAdc() ;\r
+ maxSumAdc =(mp1)->GetSum() ;\r
+ maxTimeBin =(mp1)->GetTime();\r
+ maxPad =(mp1)->GetPad() ;\r
+ maxRow =(mp1)->GetRow() ;\r
+ maxX =(mp1)->GetX();\r
+ maxY =(mp1)->GetY();\r
+ maxT =(mp1)->GetT();\r
+ \r
+ AliSimDigits *digrowTmp;\r
+ if(fRawData){\r
+ digrowTmp = (AliSimDigits*)fDigarr->GetRow(iSec,(mp1)->GetRow());\r
+ }else{\r
+ digrowTmp = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp1)->GetRow());\r
+ }\r
+ \r
+ digrowTmp->ExpandBuffer(); //decrunch\r
+ \r
+ for(Int_t itb=(mp1)->GetBegin(); itb<((mp1)->GetEnd())+1; itb++){\r
+ Int_t adcTmp = digrowTmp->GetDigitUnchecked(itb,(mp1)->GetPad());\r
+ AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp1)->GetPad(),(mp1)->GetRow(),(mp1)->GetX(),(mp1)->GetY(),(mp1)->GetT());\r
+ clusterKr.AddDigitToCluster(vtpr);\r
+ }\r
+ clusterKr.SetCenter();//set centr of the cluster\r
+ \r
+ for(Int_t it2 = it1+1; it2 < entriesArr; ++it2 ) {\r
+ AliPadMax *mp2=(AliPadMax *)maximaInSector->UncheckedAt(it2);\r
+ if (!mp2) continue;\r
+ if (TMath::Abs(clusterKr.GetCenterX() - (mp2)->GetX()) > fMaxPadRangeCm) continue;\r
+ if (TMath::Abs(clusterKr.GetCenterY() - (mp2)->GetY()) > fMaxRowRangeCm) continue; \r
+ if (TMath::Abs(clusterKr.GetCenterT() - (mp2)->GetT()) > fMaxTimeRange) continue;\r
+\r
+ {\r
+ clusterValue+=(mp2)->GetSum();\r
+ \r
+ nUsedPads++;\r
+ nUsedRows.push_back((mp2)->GetRow());\r
+ \r
+ AliSimDigits *digrowTmp1;\r
+ if(fRawData){\r
+ digrowTmp1 = (AliSimDigits*)fDigarr->GetRow(iSec,(mp2)->GetRow());\r
+ }else{\r
+ digrowTmp1 = (AliSimDigits*)fDigarr->LoadRow(iSec,(mp2)->GetRow());\r
+ }\r
+ digrowTmp1->ExpandBuffer(); //decrunch\r
+ \r
+ for(Int_t itb=(mp2)->GetBegin(); itb<(mp2)->GetEnd()+1; itb++){\r
+ Int_t adcTmp = digrowTmp1->GetDigitUnchecked(itb,(mp2)->GetPad());\r
+ AliTPCvtpr *vtpr=new AliTPCvtpr(adcTmp,itb,(mp2)->GetPad(),(mp2)->GetRow(),(mp2)->GetX(),(mp2)->GetY(),(mp2)->GetT());\r
+ clusterKr.AddDigitToCluster(vtpr);\r
+ }\r
+ \r
+ clusterKr.SetCenter();//set center of the cluster\r
+ \r
+ //which one is bigger\r
+ if( (mp2)->GetAdc() > maxDig ){\r
+ maxDig =(mp2)->GetAdc() ;\r
+ maxSumAdc =(mp2)->GetSum() ;\r
+ maxTimeBin =(mp2)->GetTime();\r
+ maxPad =(mp2)->GetPad() ;\r
+ maxRow =(mp2)->GetRow() ;\r
+ maxX =(mp2)->GetX() ;\r
+ maxY =(mp2)->GetY() ;\r
+ maxT =(mp2)->GetT() ;\r
+ } else if ( (mp2)->GetAdc() == maxDig ){\r
+ if( (mp2)->GetSum() > maxSumAdc){\r
maxDig =(mp2)->GetAdc() ;\r
maxSumAdc =(mp2)->GetSum() ;\r
maxTimeBin =(mp2)->GetTime();\r
maxRow =(mp2)->GetRow() ;\r
maxX =(mp2)->GetX() ;\r
maxY =(mp2)->GetY() ;\r
- } else if ( (mp2)->GetAdc() == maxDig ){\r
- if( (mp2)->GetSum() > maxSumAdc){\r
- maxDig =(mp2)->GetAdc() ;\r
- maxSumAdc =(mp2)->GetSum() ;\r
- maxTimeBin =(mp2)->GetTime();\r
- maxPad =(mp2)->GetPad() ;\r
- maxRow =(mp2)->GetRow() ;\r
- maxX =(mp2)->GetX() ;\r
- maxY =(mp2)->GetY() ;\r
- }\r
+ maxT =(mp2)->GetT() ;\r
}\r
- maximaInSector->RemoveAt(it2);\r
- maximaInSector->Compress();\r
- it2--;\r
-\r
- //maximaInSector.erase(mp2);\r
- //mp2--;\r
}\r
- }//inside loop\r
-\r
- tmp->SetSize();\r
- //through out ADC=6,7 on 1 pad, 2 tb and ADC=12 on 2 pads,2 tb\r
- //if(nUsedPads==1 && clusterValue/tmp->fCluster.size()<3.6)continue;\r
- //if(nUsedPads==2 && clusterValue/tmp->fCluster.size()<3.1)continue;\r
-\r
- //through out clusters on the edge and noise\r
- //if(clusterValue/tmp->fCluster.size()<fValueToSize)continue;\r
- if(clusterValue/(tmp->GetSize())<fValueToSize)continue;\r
-\r
- tmp->SetADCcluster(clusterValue);\r
- tmp->SetNPads(nUsedPads);\r
- tmp->SetMax(AliTPCvtpr(maxDig,maxTimeBin,maxPad,maxRow,maxX,maxY,maxTimeBin));\r
- tmp->SetSec(iSec);\r
- tmp->SetSize();\r
-\r
- nUsedRows.sort();\r
- nUsedRows.unique();\r
- tmp->SetNRows(nUsedRows.size());\r
- tmp->SetCenter();\r
-\r
- clusterCounter++;\r
- \r
-\r
- //save each cluster into file\r
+ delete maximaInSector->RemoveAt(it2);\r
+ }\r
+ }//inside loop\r
+ delete maximaInSector->RemoveAt(it1); \r
+ clusterKr.SetSize();\r
+ //through out clusters on the edge and noise\r
+ //if(clusterValue/clusterKr.fCluster.size()<fValueToSize)continue;\r
+ if(clusterValue/(clusterKr.GetSize())<fValueToSize)continue;\r
+ \r
+ clusterKr.SetADCcluster(clusterValue);\r
+ clusterKr.SetNPads(nUsedPads);\r
+ clusterKr.SetMax(AliTPCvtpr(maxDig,maxTimeBin,maxPad,maxRow,maxX,maxY,maxT));\r
+ clusterKr.SetSec(iSec);\r
+ clusterKr.SetSize();\r
+ \r
+ nUsedRows.sort();\r
+ nUsedRows.unique();\r
+ clusterKr.SetNRows(nUsedRows.size());\r
+ clusterKr.SetCenter();\r
+ \r
+ clusterKr.SetRMS();//Set pad,row,timebin RMS\r
+ clusterKr.Set1D();//Set size in pads and timebins\r
\r
- AliTPCClustersRow *clrow= new AliTPCClustersRow();\r
- fRowCl=clrow;\r
- clrow->SetClass("AliTPCclusterKr");\r
- clrow->SetArray(1);\r
- fOutput->GetBranch("Segment")->SetAddress(&clrow);\r
+ clusterKr.SetTimeStamp(fTimeStamp);\r
+ clusterKr.SetRun(fRun);\r
\r
- Int_t tmpCluster=0;\r
- TClonesArray * arr = fRowCl->GetArray();\r
- AliTPCclusterKr * cl = new ((*arr)[tmpCluster]) AliTPCclusterKr(*tmp);\r
- \r
- fOutput->Fill(); \r
- delete clrow;\r
- //end of save each cluster into file adc.root\r
- }//outer loop\r
+ clusterCounter++;\r
+ \r
+ \r
+ //save each cluster into file\r
+ if (fOutput){\r
+ (*fOutput)<<"Kr"<<\r
+ "Cl.="<<&clusterKr<<\r
+ "\n";\r
+ }\r
+ //end of save each cluster into file adc.root\r
+ }//outer loop\r
+}\r
\r
- }//end sector for\r
- cout<<"Number of clusters in event: "<<clusterCounter<<endl;\r
\r
- return 0;\r
-}\r
\r
////____________________________________________________________________________\r
\r
\r
-void AliTPCclustererKr::GetXY(Short_t sec,Short_t row,Short_t pad,Double_t& xGlob,Double_t& yGlob){\r
+void AliTPCclustererKr::GetXY(Int_t sec,Int_t row,Int_t pad,Double_t& xGlob,Double_t& yGlob){\r
//\r
//gives global XY coordinate of the pad\r
//\r
\r
Double_t yLocal = fParam->GetPadRowRadii(sec,row);//radius of row in sector in cm\r
\r
- Short_t padmax = fParam->GetNPads(sec,row);//number of pads in a given row\r
+ Int_t padmax = fParam->GetNPads(sec,row);//number of pads in a given row\r
Float_t padXSize;\r
if(sec<fParam->GetNInnerSector())padXSize=0.4;\r
else padXSize=0.6;\r
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff