From 917e711b48434c7f14546db30a7a057948bb63db Mon Sep 17 00:00:00 2001 From: cvetan Date: Fri, 11 Jun 2004 12:13:58 +0000 Subject: [PATCH] Coding violations... --- HLT/hough/AliL3Histogram.cxx | 24 + HLT/hough/AliL3Histogram.h | 57 ++- HLT/hough/AliL3Hough.cxx | 66 +-- HLT/hough/AliL3Hough.h | 113 +++-- HLT/hough/AliL3HoughMaxFinder.cxx | 606 +++++++++++++------------ HLT/hough/AliL3HoughMaxFinder.h | 104 ++--- HLT/hough/AliL3HoughTransformer.cxx | 258 ++++++----- HLT/hough/AliL3HoughTransformer.h | 36 +- HLT/hough/AliL3HoughTransformerRow.cxx | 543 +++++++++++----------- HLT/hough/AliL3HoughTransformerRow.h | 74 +-- 10 files changed, 972 insertions(+), 909 deletions(-) diff --git a/HLT/hough/AliL3Histogram.cxx b/HLT/hough/AliL3Histogram.cxx index 327769824bb..06db92748ec 100644 --- a/HLT/hough/AliL3Histogram.cxx +++ b/HLT/hough/AliL3Histogram.cxx @@ -29,6 +29,7 @@ ClassImp(AliL3Histogram) AliL3Histogram::AliL3Histogram() { + // Default constructor fNxbins = 0; fNybins = 0; fNcells = 0; @@ -54,6 +55,7 @@ AliL3Histogram::AliL3Histogram(Char_t *name,Char_t */*id*/, Int_t nxbin,Double_t xmin,Double_t xmax, Int_t nybin,Double_t ymin,Double_t ymax) { + // Normal constructor strcpy(fName,name); fNxbins = nxbin; @@ -93,6 +95,7 @@ AliL3Histogram::~AliL3Histogram() void AliL3Histogram::Reset() { + // Reset histogram contents if(fContent) for(Int_t i=0; i fXmax) return 0; @@ -179,6 +189,7 @@ Int_t AliL3Histogram::FindXbin(Double_t x) const Int_t AliL3Histogram::FindYbin(Double_t y) const { + // Finds the bin which correspond to y if(y < fYmin || y > fYmax) return 0; @@ -187,6 +198,7 @@ Int_t AliL3Histogram::FindYbin(Double_t y) const Int_t AliL3Histogram::GetBin(Int_t xbin,Int_t ybin) const { + // Returns the bin which correspond to xbin and ybin if(xbin < fFirstXbin || xbin > fLastXbin) return 0; if(ybin < fFirstYbin || ybin > fLastYbin) @@ -197,6 +209,7 @@ Int_t AliL3Histogram::GetBin(Int_t xbin,Int_t ybin) const Int_t AliL3Histogram::GetLabelBin(Int_t xbin,Int_t ybin) const { + // Returns the corresponding bin with the mc labels if(xbin < fFirstXbin || xbin > fLastXbin) return -1; if(ybin < fFirstYbin || ybin > fLastYbin) @@ -207,6 +220,7 @@ Int_t AliL3Histogram::GetLabelBin(Int_t xbin,Int_t ybin) const Int_t AliL3Histogram::GetBinContent(Int_t bin) const { + // Return the bin content if(bin >= fNcells) { LOG(AliL3Log::kError,"AliL3Histogram::GetBinContent","array")<= fNcells) { @@ -247,6 +263,7 @@ void AliL3Histogram::SetBinContent(Int_t bin,Int_t value) void AliL3Histogram::AddBinContent(Int_t xbin,Int_t ybin,Int_t weight) { + // Adds weight to bin content Int_t bin = GetBin(xbin,ybin); #ifdef _IFON_ if(bin == 0) @@ -258,6 +275,7 @@ void AliL3Histogram::AddBinContent(Int_t xbin,Int_t ybin,Int_t weight) void AliL3Histogram::AddBinContent(Int_t bin,Int_t weight) { + // Adds weight to bin content if(bin < 0 || bin > fNcells) { LOG(AliL3Log::kError,"AliL3Histogram::AddBinContent","array")< fLastXbin) { LOG(AliL3Log::kError,"AliL3Histogram::GetBinCenterX","xbin") @@ -316,6 +335,7 @@ Double_t AliL3Histogram::GetBinCenterX(Int_t xbin) const Double_t AliL3Histogram::GetBinCenterY(Int_t ybin) const { + // Returns the position of the center of a bin if(ybin < fFirstYbin || ybin > fLastYbin) { LOG(AliL3Log::kError,"AliL3Histogram::GetBinCenterY","ybin") @@ -328,6 +348,7 @@ Double_t AliL3Histogram::GetBinCenterY(Int_t ybin) const Double_t AliL3Histogram::GetPreciseBinCenterX(Float_t xbin) const { + // Returns the position of the center of a bin using precise values inside the bin if(xbin < (fFirstXbin-0.5) || xbin > (fLastXbin+0.5)) { LOG(AliL3Log::kError,"AliL3Histogram::GetBinCenterX","xbin") @@ -340,6 +361,7 @@ Double_t AliL3Histogram::GetPreciseBinCenterX(Float_t xbin) const Double_t AliL3Histogram::GetPreciseBinCenterY(Float_t ybin) const { + // Returns the position of the center of a bin using precise values inside the bin if(ybin < (fFirstYbin-0.5) || ybin > (fLastYbin+0.5)) { LOG(AliL3Log::kError,"AliL3Histogram::GetBinCenterY","ybin") @@ -352,6 +374,7 @@ Double_t AliL3Histogram::GetPreciseBinCenterY(Float_t ybin) const void AliL3Histogram::Draw(Char_t *option) { + // Fill the contents of the corresponding ROOT histogram and draws it #ifdef use_root if(!fRootHisto) CreateRootHisto(); @@ -375,6 +398,7 @@ void AliL3Histogram::Draw(Char_t *option) void AliL3Histogram::CreateRootHisto() { + // Create ROOT histogram out of AliL3Histogram #ifdef use_root fRootHisto = new TH2F(fName,"",fNxbins,fXmin,fXmax,fNybins,fYmin,fYmax); return; diff --git a/HLT/hough/AliL3Histogram.h b/HLT/hough/AliL3Histogram.h index 003ba80a375..9145e08d515 100644 --- a/HLT/hough/AliL3Histogram.h +++ b/HLT/hough/AliL3Histogram.h @@ -1,7 +1,7 @@ // @(#) $Id$ -#ifndef ALIL3_HISTOGRAM -#define ALIL3_HISTOGRAM +#ifndef ALIL3HISTOGRAM_H +#define ALIL3HISTOGRAM_H #include "AliL3StandardIncludes.h" #include "AliL3RootTypes.h" @@ -12,33 +12,6 @@ #endif class AliL3Histogram { - - private: - Double_t fBinwidthX; - Double_t fBinwidthY; - - protected: - Int_t *fContent; //! - Char_t fName[100]; - Int_t fNxbins; - Int_t fNybins; - Int_t fNcells; - Int_t fEntries; - Int_t fFirstXbin; - Int_t fFirstYbin; - Int_t fLastXbin; - Int_t fLastYbin; - Int_t fThreshold; - - Double_t fXmin; - Double_t fYmin; - Double_t fXmax; - Double_t fYmax; - -#ifdef use_root - TH2F *fRootHisto; -#endif - public: AliL3Histogram(); @@ -92,6 +65,32 @@ class AliL3Histogram { Int_t GetNbinsX() const {return fNxbins;} Int_t GetNbinsY() const {return fNybins;} Int_t GetNEntries() const {return fEntries;} + + protected: + Int_t *fContent; //! + Char_t fName[100]; // Name of the histogram + Int_t fNxbins; // Number of bins in the histogram + Int_t fNybins; // Number of bins in the histogram + Int_t fNcells; // Overall number of bins in the histogram + Int_t fEntries; // Number of entries in the histogram + Int_t fFirstXbin; // First active bin + Int_t fFirstYbin; // First active bin + Int_t fLastXbin; // Last active bin + Int_t fLastYbin; // Last active bin + Int_t fThreshold; // Bin content threshold + + Double_t fXmin; // Lower limit in X + Double_t fYmin; // Lower limit in Y + Double_t fXmax; // Upper limit in X + Double_t fYmax; // Upper limit in Y + +#ifdef use_root + TH2F *fRootHisto; // Corresponding ROOT histogram +#endif + + private: + Double_t fBinwidthX; // Bin width of the Hough space + Double_t fBinwidthY; // Bin width of the Hough space ClassDef(AliL3Histogram,1) //2D histogram class diff --git a/HLT/hough/AliL3Hough.cxx b/HLT/hough/AliL3Hough.cxx index 933a818ecdc..4ca24880200 100644 --- a/HLT/hough/AliL3Hough.cxx +++ b/HLT/hough/AliL3Hough.cxx @@ -102,11 +102,12 @@ AliL3Hough::AliL3Hough() #endif } -AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr) +AliL3Hough::AliL3Hough(Char_t *path,Bool_t binary,Int_t netasegments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr) { + //Normal constructor fBinary = binary; strcpy(fPath,path); - fNEtaSegments = n_eta_segments; + fNEtaSegments = netasegments; fAddHistograms = kFALSE; fDoIterative = kFALSE; fWriteDigits = kFALSE; @@ -188,11 +189,12 @@ void AliL3Hough::CleanUp() //cout << "Cleaned class mem " << endl; } -void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr,Float_t zvertex) +void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t netasegments,Bool_t bit8,Int_t tv,Char_t *infile,Char_t *ptr,Float_t zvertex) { + //Normal init of the AliL3Hough fBinary = binary; strcpy(fPath,path); - fNEtaSegments = n_eta_segments; + fNEtaSegments = netasegments; fWriteDigits = kFALSE; fUse8bits = bit8; fVersion = tv; @@ -217,6 +219,7 @@ void AliL3Hough::Init(Char_t *path,Bool_t binary,Int_t n_eta_segments,Bool_t bit void AliL3Hough::Init(Bool_t doit, Bool_t addhists) { + // Init fDoIterative = doit; fAddHistograms = addhists; @@ -311,7 +314,7 @@ void AliL3Hough::Init(Bool_t doit, Bool_t addhists) void AliL3Hough::SetTransformerParams(Float_t ptres,Float_t ptmin,Float_t ptmax,Int_t ny,Int_t patch) { - + // Setup the parameters for the Hough Transformer Int_t mrow; Float_t psi=0; if(patch==-1) @@ -351,6 +354,7 @@ void AliL3Hough::SetTransformerParams(Float_t ptres,Float_t ptmin,Float_t ptmax, /* void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t patch) { + // Setup the parameters for the Hough Transformer Int_t mrow=80; Double_t lineradius = sqrt(pow(AliL3Transform::Row2X(mrow),2) + pow(AliL3Transform::GetMaxY(mrow),2)); @@ -371,6 +375,7 @@ void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t patc */ void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t ptmin,Int_t /*patch*/) { + // Setup the parameters for the Hough Transformer Int_t mrow=79; @@ -414,6 +419,7 @@ void AliL3Hough::SetTransformerParams(Int_t nx,Int_t ny,Float_t lpt,Float_t phi) void AliL3Hough::SetThreshold(Int_t t3,Int_t patch) { + // Set digits threshold if(patch==-1) { Int_t i=0; @@ -426,6 +432,7 @@ void AliL3Hough::SetThreshold(Int_t t3,Int_t patch) void AliL3Hough::SetPeakThreshold(Int_t threshold,Int_t patch) { + // Set Peak Finder threshold if(patch==-1) { Int_t i=0; @@ -509,7 +516,7 @@ void AliL3Hough::ReadData(Int_t slice,Int_t eventnr) fEvent=eventnr; } -void AliL3Hough::Transform(Int_t *row_range) +void AliL3Hough::Transform(Int_t *rowrange) { //Transform all data given to the transformer within the given slice //(after ReadData(slice)) @@ -522,10 +529,10 @@ void AliL3Hough::Transform(Int_t *row_range) if((fVersion != 4) || (i == 0)) fHoughTransformer[i]->Reset();//Reset the histograms fBenchmark->Start("Hough Transform"); - if(!row_range) + if(!rowrange) fHoughTransformer[i]->TransformCircle(); else - fHoughTransformer[i]->TransformCircleC(row_range,1); + fHoughTransformer[i]->TransformCircleC(rowrange,1); fBenchmark->Stop("Hough Transform"); } cpuTime = GetCpuTime() - initTime; @@ -535,6 +542,7 @@ void AliL3Hough::Transform(Int_t *row_range) void AliL3Hough::MergePatches() { + // Merge patches if they are not summed if(fAddHistograms) //Nothing to merge here return; fMerger->MergePatches(kTRUE); @@ -542,6 +550,7 @@ void AliL3Hough::MergePatches() void AliL3Hough::MergeInternally() { + // Merge patches internally if(fAddHistograms) fInterMerger->FillTracks(fTracks[0]); else @@ -607,7 +616,7 @@ void AliL3Hough::ProcessPatchIter(Int_t patch) //Process patch in a iterative way. //transform + peakfinding + evaluation + transform +... - Int_t num_of_tries = 5; + Int_t numoftries = 5; AliL3HoughBaseTransformer *tr = fHoughTransformer[patch]; AliL3TrackArray *tracks = fTracks[patch]; tracks->Reset(); @@ -617,7 +626,7 @@ void AliL3Hough::ProcessPatchIter(Int_t patch) ev->SetNumOfRowsToMiss(3); ev->SetNumOfPadsToLook(2); AliL3Histogram *hist; - for(Int_t t=0; tReset(); tr->TransformCircle(); @@ -721,6 +730,7 @@ void AliL3Hough::AddAllHistogramsRows() void AliL3Hough::AddTracks() { + // Add current slice slice tracks to the global list of found tracks if(!fTracks[0]) { cerr<<"AliL3Hough::AddTracks : No tracks"<Start("Find Maxima"); - for(Int_t i=0; iReset(); @@ -814,22 +825,23 @@ void AliL3Hough::FindTrackCandidatesRow() void AliL3Hough::FindTrackCandidates() { + // Find AliL3HoughTransformer track candidates if(fVersion == 4) { LOG(AliL3Log::kError,"AliL3Hough::FindTrackCandidatesRow()","") <<"Incompatible Peak Finder version!"<Start("Find Maxima"); - for(Int_t i=0; iReset(); @@ -874,7 +886,7 @@ void AliL3Hough::InitEvaluate() fEval[i]->InitTransformer(fHoughTransformer[i]); } -Int_t AliL3Hough::Evaluate(Int_t road_width,Int_t nrowstomiss) +Int_t AliL3Hough::Evaluate(Int_t roadwidth,Int_t nrowstomiss) { //Evaluate the tracks, by looking along the road in the raw data. //If track does not cross all padrows - rows2miss, it is removed from the arrray. @@ -894,7 +906,7 @@ Int_t AliL3Hough::Evaluate(Int_t road_width,Int_t nrowstomiss) return 0; } - Int_t removed_tracks=0; + Int_t removedtracks=0; AliL3TrackArray *tracks=0; if(fAddHistograms) @@ -909,7 +921,7 @@ Int_t AliL3Hough::Evaluate(Int_t road_width,Int_t nrowstomiss) } for(Int_t i=0; iGetNHits() < AliL3Transform::GetNRows() - nrowstomiss) { tracks->Remove(j); - removed_tracks++; + removedtracks++; } } tracks->Compress(); tracks->QSort(); } - return removed_tracks; + return removedtracks; } -void AliL3Hough::EvaluatePatch(Int_t i,Int_t road_width,Int_t nrowstomiss) +void AliL3Hough::EvaluatePatch(Int_t i,Int_t roadwidth,Int_t nrowstomiss) { //Evaluate patch i. fEval[i]->InitTransformer(fHoughTransformer[i]); - fEval[i]->SetNumOfPadsToLook(road_width); + fEval[i]->SetNumOfPadsToLook(roadwidth); fEval[i]->SetNumOfRowsToMiss(nrowstomiss); //fEval[i]->RemoveFoundTracks(); @@ -1016,6 +1028,7 @@ void AliL3Hough::MergeEtaSlices() void AliL3Hough::WriteTracks(Char_t *path) { + // Write found tracks into file //cout<<"AliL3Hough::WriteTracks : Sorting the tracsk"<QSort(); @@ -1030,6 +1043,7 @@ void AliL3Hough::WriteTracks(Char_t *path) void AliL3Hough::WriteTracks(Int_t slice,Char_t *path) { + // Write found tracks slice by slice into file AliL3MemHandler mem; Char_t fname[100]; @@ -1054,8 +1068,8 @@ void AliL3Hough::WriteTracks(Int_t slice,Char_t *path) void AliL3Hough::WriteDigits(Char_t *outfile) { -#ifdef use_aliroot //Write the current data to a new rootfile. +#ifdef use_aliroot for(Int_t i=0; iGetFirstYbin(); Int_t ymax = hist->GetLastYbin(); Int_t bin; - Double_t value,max_value=0; + Double_t value,maxvalue=0; - Int_t max_xbin=0,max_ybin=0; + Int_t maxxbin=0,maxybin=0; for(Int_t xbin=xmin; xbin<=xmax; xbin++) { for(Int_t ybin=ymin; ybin<=ymax; ybin++) { bin = hist->GetBin(xbin,ybin); value = hist->GetBinContent(bin); - if(value>max_value) + if(value>maxvalue) { - max_value = value; - max_xbin = xbin; - max_ybin = ybin; + maxvalue = value; + maxxbin = xbin; + maxybin = ybin; } } } - if(max_value == 0) + if(maxvalue == 0) return; if(fNPeaks > fNMax) @@ -207,29 +210,29 @@ void AliL3HoughMaxFinder::FindAbsMaxima() return; } - Double_t max_x = hist->GetBinCenterX(max_xbin); - Double_t max_y = hist->GetBinCenterY(max_ybin); - fXPeaks[fNPeaks] = max_x; - fYPeaks[fNPeaks] = max_y; - fWeight[fNPeaks] = (Int_t)max_value; + Double_t maxx = hist->GetBinCenterX(maxxbin); + Double_t maxy = hist->GetBinCenterY(maxybin); + fXPeaks[fNPeaks] = maxx; + fYPeaks[fNPeaks] = maxy; + fWeight[fNPeaks] = (Int_t)maxvalue; fNPeaks++; #ifndef no_root if(fNtuppel) { - Int_t bin3 = hist->GetBin(max_xbin-1,max_ybin); - Int_t bin5 = hist->GetBin(max_xbin+1,max_ybin); - Int_t bin1 = hist->GetBin(max_xbin,max_ybin-1); - Int_t bin7 = hist->GetBin(max_xbin,max_ybin+1); + Int_t bin3 = hist->GetBin(maxxbin-1,maxybin); + Int_t bin5 = hist->GetBin(maxxbin+1,maxybin); + Int_t bin1 = hist->GetBin(maxxbin,maxybin-1); + Int_t bin7 = hist->GetBin(maxxbin,maxybin+1); - fNtuppel->Fill(max_x,max_y,max_value,hist->GetBinContent(bin3),hist->GetBinContent(bin5),hist->GetBinContent(bin1),hist->GetBinContent(bin7)); + fNtuppel->Fill(maxx,maxy,maxvalue,hist->GetBinContent(bin3),hist->GetBinContent(bin5),hist->GetBinContent(bin1),hist->GetBinContent(bin7)); } #endif } void AliL3HoughMaxFinder::FindBigMaxima() { - + // Another Peak finder AliL3Histogram *hist = fCurrentHisto; if(hist->GetNEntries() == 0) @@ -239,32 +242,32 @@ void AliL3HoughMaxFinder::FindBigMaxima() Int_t xmax = hist->GetLastXbin(); Int_t ymin = hist->GetFirstYbin(); Int_t ymax = hist->GetLastYbin(); - Int_t bin[25],bin_index; + Int_t bin[25],binindex; Double_t value[25]; for(Int_t xbin=xmin+2; xbinGetBin(xb,yb); - value[bin_index]=hist->GetBinContent(bin[bin_index]); - bin_index++; + bin[binindex]=hist->GetBin(xb,yb); + value[binindex]=hist->GetBinContent(bin[binindex]); + binindex++; } } if(value[12]==0) continue; Int_t b=0; while(1) { - if(value[b] > value[12] || b==bin_index) break; + if(value[b] > value[12] || b==binindex) break; b++; //printf("b %d\n",b); } - if(b == bin_index) + if(b == binindex) { //Found maxima if(fNPeaks > fNMax) @@ -273,10 +276,10 @@ void AliL3HoughMaxFinder::FindBigMaxima() return; } - Double_t max_x = hist->GetBinCenterX(xbin); - Double_t max_y = hist->GetBinCenterY(ybin); - fXPeaks[fNPeaks] = max_x; - fYPeaks[fNPeaks] = max_y; + Double_t maxx = hist->GetBinCenterX(xbin); + Double_t maxy = hist->GetBinCenterY(ybin); + fXPeaks[fNPeaks] = maxx; + fYPeaks[fNPeaks] = maxy; fNPeaks++; } } @@ -331,8 +334,8 @@ void AliL3HoughMaxFinder::FindMaxima(Int_t threshold) && value[4]>value[7] && value[4]>value[8]) { //Found a local maxima - Float_t max_x = fCurrentHisto->GetBinCenterX(xbin); - Float_t max_y = fCurrentHisto->GetBinCenterY(ybin); + Float_t maxx = fCurrentHisto->GetBinCenterX(xbin); + Float_t maxy = fCurrentHisto->GetBinCenterY(ybin); if((Int_t)value[4] <= threshold) continue;//central bin below threshold if(fNPeaks >= fNMax) @@ -348,8 +351,8 @@ void AliL3HoughMaxFinder::FindMaxima(Int_t threshold) if(value[1]/value[4] > fGradY && value[7]/value[4] > fGradY) continue; - fXPeaks[fNPeaks] = max_x; - fYPeaks[fNPeaks] = max_y; + fXPeaks[fNPeaks] = maxx; + fYPeaks[fNPeaks] = maxy; fWeight[fNPeaks] = (Int_t)value[4]; fNPeaks++; @@ -358,13 +361,13 @@ void AliL3HoughMaxFinder::FindMaxima(Int_t threshold) Bool_t bigger = kFALSE; for(Int_t p=0; p fWeight[p]) //this peak is bigger. { - fXPeaks[p] = max_x; - fYPeaks[p] = max_y; + fXPeaks[p] = maxx; + fYPeaks[p] = maxy; fWeight[p] = (Int_t)value[4]; } else @@ -373,8 +376,8 @@ void AliL3HoughMaxFinder::FindMaxima(Int_t threshold) } if(!bigger) //there were no overlapping peaks. { - fXPeaks[fNPeaks] = max_x; - fYPeaks[fNPeaks] = max_y; + fXPeaks[fNPeaks] = maxx; + fYPeaks[fNPeaks] = maxy; fWeight[fNPeaks] = (Int_t)value[4]; fNPeaks++; } @@ -385,13 +388,13 @@ void AliL3HoughMaxFinder::FindMaxima(Int_t threshold) } -struct Window +struct AliL3Window { - Int_t start; - Int_t sum; + Int_t fStart; // Start + Int_t fSum; // Sum }; -void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) +void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cutratio) { //Peak finder which looks for peaks with a certain shape. //The first step involves a pre-peak finder, which looks for peaks @@ -421,17 +424,17 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) //Start by looking for pre-peaks: - Window **local_maxima = new Window*[hist->GetNbinsY()]; + AliL3Window **localmaxima = new AliL3Window*[hist->GetNbinsY()]; Short_t *nmaxs = new Short_t[hist->GetNbinsY()]; - Int_t n,last_sum,sum; - Bool_t sum_was_rising; + Int_t n,lastsum,sum; + Bool_t sumwasrising; for(Int_t ybin=ymin; ybin<=ymax; ybin++) { - local_maxima[ybin-ymin] = new Window[hist->GetNbinsX()]; + localmaxima[ybin-ymin] = new AliL3Window[hist->GetNbinsX()]; nmaxs[ybin-ymin] = 0; - sum_was_rising=0; - last_sum=0; + sumwasrising=0; + lastsum=0; n=0; for(Int_t xbin=xmin; xbin<=xmax-kappawindow; xbin++) { @@ -439,21 +442,21 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) for(Int_t lbin=xbin; lbinGetBinContent(hist->GetBin(lbin,ybin)); - if(sum < last_sum) + if(sum < lastsum) { if(sum > fThreshold) - if(sum_was_rising)//Previous sum was a local maxima + if(sumwasrising)//Previous sum was a local maxima { - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].start = xbin-1; - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].sum = last_sum; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fStart = xbin-1; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fSum = lastsum; nmaxs[ybin-ymin]++; } - sum_was_rising=0; + sumwasrising=0; } else if(sum > 0) - sum_was_rising=1; - last_sum=sum; + sumwasrising=1; + lastsum=sum; } } @@ -465,13 +468,13 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) { for(Int_t i=0; iGetBinContent(hist->GetBin(k,ybin)) > maxvalue) { maxvalue = hist->GetBinContent(hist->GetBin(k,ybin)); @@ -481,27 +484,27 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) //start expanding in the psi-direction: - Int_t lb = local_maxima[ybin-ymin][i].start; + Int_t lb = localmaxima[ybin-ymin][i].fStart; //Int_t ystart=ybin; - starts[ybin] = local_maxima[ybin-ymin][i].start; + starts[ybin] = localmaxima[ybin-ymin][i].fStart; maxs[ybin] = maxxbin; Int_t yl=ybin-1,nybins=1; - //cout<<"Starting search at ybin "<= ymin) { Bool_t found=0; for(Int_t j=0; j= lb && local_maxima[yl-ymin][j].sum > 0) + if( localmaxima[yl-ymin][j].fStart - lb < 0) continue; + if( localmaxima[yl-ymin][j].fStart < lb + kappawindow + match && + localmaxima[yl-ymin][j].fStart >= lb && localmaxima[yl-ymin][j].fSum > 0) { - //cout<<"match at ybin "<GetBinCenterY(yl)<<" start "<GetBinCenterY(yl)<<" start "<GetBinContent(hist->GetBin(k,yl)) > maxvalue) { @@ -517,11 +520,11 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) } } nybins++; - starts[yl] = local_maxima[yl-ymin][j].start; + starts[yl] = localmaxima[yl-ymin][j].fStart; maxs[yl] = lmaxxbin; - local_maxima[yl-ymin][j].sum=-1; //Mark as used + localmaxima[yl-ymin][j].fSum=-1; //Mark as used found=1; - lb = local_maxima[yl-ymin][j].start; + lb = localmaxima[yl-ymin][j].fStart; break;//Since we found a match in this bin, we dont have to search it anymore, goto next bin. } } @@ -557,9 +560,9 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) //cout<<"ratio "< cut_ratio || lower_ratio > cut_ratio) + if(upperratio > cutratio || lowerratio > cutratio) truepeak=kFALSE; if(truepeak) @@ -654,41 +657,41 @@ void AliL3HoughMaxFinder::FindAdaptedPeaks(Int_t kappawindow,Float_t cut_ratio) } for(Int_t i=0; iGetNbinsY(); i++) - delete local_maxima[i]; + delete localmaxima[i]; - delete [] local_maxima; + delete [] localmaxima; delete [] nmaxs; delete [] starts; delete [] maxs; } -struct PreYPeak +struct AliL3PreYPeak { - Int_t start_position; - Int_t end_position; - Int_t min_value; - Int_t max_value; - Int_t prev_value; - Int_t left_value; - Int_t right_value; + Int_t fStartPosition; // Start position in X + Int_t fEndPosition; // End position in X + Int_t fMinValue; // Minimum value inside the prepeak + Int_t fMaxValue; // Maximum value inside the prepeak + Int_t fPrevValue; // Neighbour values + Int_t fLeftValue; // Neighbour values + Int_t fRightValue; // Neighbour values }; -struct Pre2DPeak +struct AliL3Pre2DPeak { - Float_t x; - Float_t y; - Float_t size_x; - Float_t size_y; - Int_t start_x; - Int_t start_y; - Int_t end_x; - Int_t end_y; - Float_t weight; + Float_t fX; // X coordinate of the preak + Float_t fY; // Y coordinate of the preak + Float_t fSizeX; // Size of the peak + Float_t fSizeY; // Size of the peak + Int_t fStartX; // Start position of the peak + Int_t fStartY; // Start position of the peak + Int_t fEndX; // End position of the peak + Int_t fEndY; // End position of the peak + Float_t fWeight; // Weight assigned to the peak }; void AliL3HoughMaxFinder::FindAdaptedRowPeaks(Int_t kappawindow,Int_t xsize,Int_t ysize) { - + // Peak finder which is working over the Hough Space provided by the AliL3HoughTransformerRow class AliL3Histogram *hist = fCurrentHisto; if(!hist) @@ -707,140 +710,140 @@ void AliL3HoughMaxFinder::FindAdaptedRowPeaks(Int_t kappawindow,Int_t xsize,Int_ //Start by looking for pre-peaks: - PreYPeak **local_maxima = new PreYPeak*[hist->GetNbinsY()]; + AliL3PreYPeak **localmaxima = new AliL3PreYPeak*[hist->GetNbinsY()]; Short_t *nmaxs = new Short_t[hist->GetNbinsY()]; - Int_t last_value=0,value=0; + Int_t lastvalue=0,value=0; for(Int_t ybin=ymin; ybin<=ymax; ybin++) { - local_maxima[ybin-ymin] = new PreYPeak[hist->GetNbinsX()]; + localmaxima[ybin-ymin] = new AliL3PreYPeak[hist->GetNbinsX()]; nmaxs[ybin-ymin] = 0; - last_value = 0; + lastvalue = 0; Bool_t found = 0; for(Int_t xbin=xmin; xbin<=xmax; xbin++) { value = hist->GetBinContent(hist->GetBin(xbin,ybin)); if(value > 0) { - if(abs(value - last_value) > 1) + if(abs(value - lastvalue) > 1) { if(found) { - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].right_value = value; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fRightValue = value; nmaxs[ybin-ymin]++; } - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].start_position = xbin; - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].end_position = xbin; - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].min_value = value; - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].max_value = value; - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].prev_value = 0; - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].left_value = last_value; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fStartPosition = xbin; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fEndPosition = xbin; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fMinValue = value; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fMaxValue = value; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fPrevValue = 0; + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fLeftValue = lastvalue; found = 1; } - if(abs(value - last_value) <= 1) + if(abs(value - lastvalue) <= 1) { - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].end_position = xbin; - if(value>local_maxima[ybin-ymin][nmaxs[ybin-ymin]].max_value) - local_maxima[ybin-ymin][nmaxs[ybin-ymin]].max_value = value; - if(valuelocalmaxima[ybin-ymin][nmaxs[ybin-ymin]].fMaxValue) + localmaxima[ybin-ymin][nmaxs[ybin-ymin]].fMaxValue = value; + if(value= ymin; ybin--) { for(Int_t i=0; i= ymin) + while(localy >= ymin) { Bool_t found=0; - for(Int_t j=0; j= (temp_x_start - kappawindow))) + if( (localmaxima[localy-ymin][j].fStartPosition <= (tempxend + kappawindow)) && (localmaxima[localy-ymin][j].fEndPosition >= (tempxstart - kappawindow))) { - if(((local_maxima[local_y-ymin][j].min_value <= local_max_value) && (local_maxima[local_y-ymin][j].min_value >= local_min_value)) || - ((local_maxima[local_y-ymin][j].max_value >= local_min_value) && (local_maxima[local_y-ymin][j].max_value <= local_max_value))) + if(((localmaxima[localy-ymin][j].fMinValue <= localmaxvalue) && (localmaxima[localy-ymin][j].fMinValue >= localminvalue)) || + ((localmaxima[localy-ymin][j].fMaxValue >= localminvalue) && (localmaxima[localy-ymin][j].fMaxValue <= localmaxvalue))) { - if(local_maxima[local_y-ymin][j].end_position > local_x_end) - local_x_end = local_maxima[local_y-ymin][j].end_position; - if(local_maxima[local_y-ymin][j].start_position < local_x_start) - local_x_start = local_maxima[local_y-ymin][j].start_position; - temp_x_start = local_maxima[local_y-ymin][j].start_position; - temp_x_end = local_maxima[local_y-ymin][j].end_position; - if(local_maxima[local_y-ymin][j].min_value < local_min_value) - local_min_value = local_maxima[local_y-ymin][j].min_value; - if(local_maxima[local_y-ymin][j].max_value > local_max_value) - local_max_value = local_maxima[local_y-ymin][j].max_value; - if(local_maxima[local_y-ymin][j].right_value > local_right_value) - local_right_value = local_maxima[local_y-ymin][j].right_value; - if(local_maxima[local_y-ymin][j].left_value > local_left_value) - local_left_value = local_maxima[local_y-ymin][j].left_value; - local_maxima[local_y-ymin][j].min_value = -1; + if(localmaxima[localy-ymin][j].fEndPosition > localxend) + localxend = localmaxima[localy-ymin][j].fEndPosition; + if(localmaxima[localy-ymin][j].fStartPosition < localxstart) + localxstart = localmaxima[localy-ymin][j].fStartPosition; + tempxstart = localmaxima[localy-ymin][j].fStartPosition; + tempxend = localmaxima[localy-ymin][j].fEndPosition; + if(localmaxima[localy-ymin][j].fMinValue < localminvalue) + localminvalue = localmaxima[localy-ymin][j].fMinValue; + if(localmaxima[localy-ymin][j].fMaxValue > localmaxvalue) + localmaxvalue = localmaxima[localy-ymin][j].fMaxValue; + if(localmaxima[localy-ymin][j].fRightValue > localrightvalue) + localrightvalue = localmaxima[localy-ymin][j].fRightValue; + if(localmaxima[localy-ymin][j].fLeftValue > localleftvalue) + localleftvalue = localmaxima[localy-ymin][j].fLeftValue; + localmaxima[localy-ymin][j].fMinValue = -1; found = 1; nybins++; break; } else { - if(local_max_value > local_maxima[local_y-ymin][j].prev_value) - local_maxima[local_y-ymin][j].prev_value = local_max_value; - if(local_maxima[local_y-ymin][j].max_value > local_next_value) - local_next_value = local_maxima[local_y-ymin][j].max_value; + if(localmaxvalue > localmaxima[localy-ymin][j].fPrevValue) + localmaxima[localy-ymin][j].fPrevValue = localmaxvalue; + if(localmaxima[localy-ymin][j].fMaxValue > localnextvalue) + localnextvalue = localmaxima[localy-ymin][j].fMaxValue; } } } - if(!found || local_y == ymin)//no more local maximas to be matched, so write the final peak and break the expansion: + if(!found || localy == ymin)//no more local maximas to be matched, so write the final peak and break the expansion: { - if((nybins > ysize) && ((local_x_end-local_x_start+1) > xsize) && (local_max_value > local_prev_value) && (local_max_value > local_next_value) && (local_max_value > local_left_value) && (local_max_value > local_right_value)) - // if((nybins > ysize) && ((local_x_end-local_x_start+1) > xsize)) + if((nybins > ysize) && ((localxend-localxstart+1) > xsize) && (localmaxvalue > localprevvalue) && (localmaxvalue > localnextvalue) && (localmaxvalue > localleftvalue) && (localmaxvalue > localrightvalue)) + // if((nybins > ysize) && ((localxend-localxstart+1) > xsize)) { - maxima[nmaxima].x = ((Float_t)local_x_start+(Float_t)local_x_end)/2.0; - maxima[nmaxima].y = ((Float_t)ybin+(Float_t)(local_y+1))/2.0; - maxima[nmaxima].size_x = local_x_end-local_x_start+1; - maxima[nmaxima].size_y = nybins; - maxima[nmaxima].weight = (local_min_value+local_max_value)/2; - maxima[nmaxima].start_x = local_x_start; - maxima[nmaxima].end_x = local_x_end; - maxima[nmaxima].start_y = local_y +1; - maxima[nmaxima].end_y = ybin; + maxima[nmaxima].fX = ((Float_t)localxstart+(Float_t)localxend)/2.0; + maxima[nmaxima].fY = ((Float_t)ybin+(Float_t)(localy+1))/2.0; + maxima[nmaxima].fSizeX = localxend-localxstart+1; + maxima[nmaxima].fSizeY = nybins; + maxima[nmaxima].fWeight = (localminvalue+localmaxvalue)/2; + maxima[nmaxima].fStartX = localxstart; + maxima[nmaxima].fEndX = localxend; + maxima[nmaxima].fStartY = localy +1; + maxima[nmaxima].fEndY = ybin; #ifdef do_mc - // cout<<"Peak found at: "<<((Float_t)local_x_start+(Float_t)local_x_end)/2.0<<" "<<((Float_t)ybin+(Float_t)(local_y+1))/2.0<<" "< maxima[j].size_x*maxima[j].size_y) - maxima[j].weight = -maxima[j].weight; - if(maxima[i].size_x*maxima[i].size_y < maxima[j].size_x*maxima[j].size_y) - maxima[i].weight = -maxima[i].weight; + if(maxima[i].fSizeX*maxima[i].fSizeY > maxima[j].fSizeX*maxima[j].fSizeY) + maxima[j].fWeight = -maxima[j].fWeight; + if(maxima[i].fSizeX*maxima[i].fSizeY < maxima[j].fSizeX*maxima[j].fSizeY) + maxima[i].fWeight = -maxima[i].fWeight; #ifdef do_mc - // cout<<"Merge peaks "< 0) { - fXPeaks[fNPeaks] = hist->GetPreciseBinCenterX(maxima[i].x); - fYPeaks[fNPeaks] = hist->GetPreciseBinCenterY(maxima[i].y); - fWeight[fNPeaks] = (Int_t)maxima[i].weight; + if(maxima[i].fWeight > 0) { + fXPeaks[fNPeaks] = hist->GetPreciseBinCenterX(maxima[i].fX); + fYPeaks[fNPeaks] = hist->GetPreciseBinCenterY(maxima[i].fY); + fWeight[fNPeaks] = (Int_t)maxima[i].fWeight; #ifdef do_mc - cout<<"Final Peak found at: "<= 1) continue; - if((maxima[i].start_x <= fENDXPeaks[j]+1) && (maxima[i].end_x >= fSTARTXPeaks[j]-1)) { - if((maxima[i].start_y <= fENDYPeaks[j]+1) && (maxima[i].end_y >= fSTARTYPeaks[j]-1)) { + if((maxima[i].fStartX <= fENDXPeaks[j]+1) && (maxima[i].fEndX >= fSTARTXPeaks[j]-1)) { + if((maxima[i].fStartY <= fENDYPeaks[j]+1) && (maxima[i].fEndY >= fSTARTYPeaks[j]-1)) { //merge merged = kTRUE; - fXPeaks[fNPeaks] = (hist->GetPreciseBinCenterX(maxima[i].x)+(fENDETAPeaks[j]-fSTARTETAPeaks[j]+1)*fXPeaks[j])/(fENDETAPeaks[j]-fSTARTETAPeaks[j]+2); - fYPeaks[fNPeaks] = (hist->GetPreciseBinCenterY(maxima[i].y)+(fENDETAPeaks[j]-fSTARTETAPeaks[j]+1)*fYPeaks[j])/(fENDETAPeaks[j]-fSTARTETAPeaks[j]+2); - fWeight[fNPeaks] = (Int_t)maxima[i].weight + fWeight[j]; - fSTARTXPeaks[fNPeaks] = maxima[i].start_x; - fSTARTYPeaks[fNPeaks] = maxima[i].start_y; - fENDXPeaks[fNPeaks] = maxima[i].end_x; - fENDYPeaks[fNPeaks] = maxima[i].end_y; + fXPeaks[fNPeaks] = (hist->GetPreciseBinCenterX(maxima[i].fX)+(fENDETAPeaks[j]-fSTARTETAPeaks[j]+1)*fXPeaks[j])/(fENDETAPeaks[j]-fSTARTETAPeaks[j]+2); + fYPeaks[fNPeaks] = (hist->GetPreciseBinCenterY(maxima[i].fY)+(fENDETAPeaks[j]-fSTARTETAPeaks[j]+1)*fYPeaks[j])/(fENDETAPeaks[j]-fSTARTETAPeaks[j]+2); + fWeight[fNPeaks] = (Int_t)maxima[i].fWeight + fWeight[j]; + fSTARTXPeaks[fNPeaks] = maxima[i].fStartX; + fSTARTYPeaks[fNPeaks] = maxima[i].fStartY; + fENDXPeaks[fNPeaks] = maxima[i].fEndX; + fENDYPeaks[fNPeaks] = maxima[i].fEndY; fSTARTETAPeaks[fNPeaks] = fSTARTETAPeaks[j]; fENDETAPeaks[fNPeaks] = fCurrentEtaSlice; fNPeaks++; @@ -971,16 +974,16 @@ void AliL3HoughMaxFinder::FindAdaptedRowPeaks(Int_t kappawindow,Int_t xsize,Int_ } } } - fXPeaks[fNPeaks] = hist->GetPreciseBinCenterX(maxima[i].x); - fYPeaks[fNPeaks] = hist->GetPreciseBinCenterY(maxima[i].y); + fXPeaks[fNPeaks] = hist->GetPreciseBinCenterX(maxima[i].fX); + fYPeaks[fNPeaks] = hist->GetPreciseBinCenterY(maxima[i].fY); if(!merged) - fWeight[fNPeaks] = (Int_t)maxima[i].weight; + fWeight[fNPeaks] = (Int_t)maxima[i].fWeight; else - fWeight[fNPeaks] = -(Int_t)maxima[i].weight; - fSTARTXPeaks[fNPeaks] = maxima[i].start_x; - fSTARTYPeaks[fNPeaks] = maxima[i].start_y; - fENDXPeaks[fNPeaks] = maxima[i].end_x; - fENDYPeaks[fNPeaks] = maxima[i].end_y; + fWeight[fNPeaks] = -(Int_t)maxima[i].fWeight; + fSTARTXPeaks[fNPeaks] = maxima[i].fStartX; + fSTARTYPeaks[fNPeaks] = maxima[i].fStartY; + fENDXPeaks[fNPeaks] = maxima[i].fEndX; + fENDYPeaks[fNPeaks] = maxima[i].fEndY; fSTARTETAPeaks[fNPeaks] = fCurrentEtaSlice; fENDETAPeaks[fNPeaks] = fCurrentEtaSlice; fNPeaks++; @@ -989,20 +992,20 @@ void AliL3HoughMaxFinder::FindAdaptedRowPeaks(Int_t kappawindow,Int_t xsize,Int_ fN2PeaksPrevEtaSlice = fNPeaks; for(Int_t i=0; iGetNbinsY(); i++) - delete local_maxima[i]; + delete localmaxima[i]; - delete [] local_maxima; + delete [] localmaxima; delete [] nmaxs; } -void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides) +void AliL3HoughMaxFinder::FindPeak1(Int_t ywindow,Int_t xbinsides) { //Testing mutliple peakfinding. //The algorithm searches the histogram for prepreaks by looking in windows - //for each bin on the xaxis. The size of these windows is controlled by y_window. + //for each bin on the xaxis. The size of these windows is controlled by ywindow. //Then the prepreaks are sorted according to their weight (sum inside window), //and the peak positions are calculated by taking the weighted mean in both - //x and y direction. The size of the peak in x-direction is controlled by x_bin_sides. + //x and y direction. The size of the peak in x-direction is controlled by xbinsides. if(!fCurrentHisto) { @@ -1012,13 +1015,13 @@ void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides) if(fCurrentHisto->GetNEntries()==0) return; - //Int_t y_window=2; - //Int_t x_bin_sides=1; + //Int_t ywindow=2; + //Int_t xbinsides=1; //Float_t max_kappa = 0.001; //Float_t max_phi0 = 0.08; - Int_t max_sum=0; + Int_t maxsum=0; Int_t xmin = fCurrentHisto->GetFirstXbin(); Int_t xmax = fCurrentHisto->GetLastXbin(); @@ -1026,40 +1029,40 @@ void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides) Int_t ymax = fCurrentHisto->GetLastYbin(); Int_t nbinsx = fCurrentHisto->GetNbinsX()+1; - AxisWindow **windowPt = new AxisWindow*[nbinsx]; - AxisWindow **anotherPt = new AxisWindow*[nbinsx]; + AliL3AxisWindow **windowPt = new AliL3AxisWindow*[nbinsx]; + AliL3AxisWindow **anotherPt = new AliL3AxisWindow*[nbinsx]; for(Int_t i=0; iGetBin(xbin,b); - sum_in_window += (Int_t)fCurrentHisto->GetBinContent(bin); + suminwindow += (Int_t)fCurrentHisto->GetBinContent(bin); } - if(sum_in_window > max_sum) + if(suminwindow > maxsum) { - max_sum = sum_in_window; - windowPt[xbin]->ymin = ybin; - windowPt[xbin]->ymax = ybin + y_window; - windowPt[xbin]->weight = sum_in_window; - windowPt[xbin]->xbin = xbin; + maxsum = suminwindow; + windowPt[xbin]->fYmin = ybin; + windowPt[xbin]->fYmax = ybin + ywindow; + windowPt[xbin]->fWeight = suminwindow; + windowPt[xbin]->fXbin = xbin; } } } @@ -1071,19 +1074,19 @@ void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides) for(Int_t i=0; ixbin; + Int_t xbin = windowPt[i]->fXbin; if(xbin xmax-1) continue; //Check if this is really a local maxima - if(anotherPt[xbin-1]->weight > anotherPt[xbin]->weight || - anotherPt[xbin+1]->weight > anotherPt[xbin]->weight) + if(anotherPt[xbin-1]->fWeight > anotherPt[xbin]->fWeight || + anotherPt[xbin+1]->fWeight > anotherPt[xbin]->fWeight) continue; - for(Int_t j=windowPt[i]->ymin; jymax; j++) + for(Int_t j=windowPt[i]->fYmin; jfYmax; j++) { //Calculate the mean in y direction: - Int_t bin = fCurrentHisto->GetBin(windowPt[i]->xbin,j); + Int_t bin = fCurrentHisto->GetBin(windowPt[i]->fXbin,j); top += (fCurrentHisto->GetBinCenterY(j))*(fCurrentHisto->GetBinContent(bin)); butt += fCurrentHisto->GetBinContent(bin); } @@ -1091,7 +1094,7 @@ void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides) if(butt < fThreshold) continue; - fXPeaks[fNPeaks] = fCurrentHisto->GetBinCenterX(windowPt[i]->xbin); + fXPeaks[fNPeaks] = fCurrentHisto->GetBinCenterX(windowPt[i]->fXbin); fYPeaks[fNPeaks] = top/butt; fWeight[fNPeaks] = (Int_t)butt; //cout<<"mean in y "<xbin<<" content "<FindXbin(fXPeaks[i]); - if(xbin - x_bin_sides < xmin || xbin + x_bin_sides > xmax) continue; + if(xbin - xbinsides < xmin || xbin + xbinsides > xmax) continue; top=butt=0; ytop=0,ybutt=0; w=0; - prev = xbin - x_bin_sides+1; - for(Int_t j=xbin-x_bin_sides; j<=xbin+x_bin_sides; j++) + prev = xbin - xbinsides+1; + for(Int_t j=xbin-xbinsides; j<=xbin+xbinsides; j++) { /* //Check if the windows are overlapping: @@ -1125,10 +1128,10 @@ void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides) prev = j; */ - top += fCurrentHisto->GetBinCenterX(j)*anotherPt[j]->weight; - butt += anotherPt[j]->weight; + top += fCurrentHisto->GetBinCenterX(j)*anotherPt[j]->fWeight; + butt += anotherPt[j]->fWeight; - for(Int_t k=anotherPt[j]->ymin; kymax; k++) + for(Int_t k=anotherPt[j]->fYmin; kfYmax; k++) { Int_t bin = fCurrentHisto->GetBin(j,k); ytop += (fCurrentHisto->GetBinCenterY(k))*(fCurrentHisto->GetBinContent(bin)); @@ -1163,12 +1166,12 @@ void AliL3HoughMaxFinder::FindPeak1(Int_t y_window,Int_t x_bin_sides) delete [] anotherPt; } -void AliL3HoughMaxFinder::SortPeaks(struct AxisWindow **a,Int_t first,Int_t last) +void AliL3HoughMaxFinder::SortPeaks(struct AliL3AxisWindow **a,Int_t first,Int_t last) { //General sorting routine //Sort according to PeakCompare() - static struct AxisWindow *tmp; + static struct AliL3AxisWindow *tmp; static int i; // "static" to save stack space int j; @@ -1205,10 +1208,11 @@ void AliL3HoughMaxFinder::SortPeaks(struct AxisWindow **a,Int_t first,Int_t last } -Int_t AliL3HoughMaxFinder::PeakCompare(struct AxisWindow *a,struct AxisWindow *b) +Int_t AliL3HoughMaxFinder::PeakCompare(struct AliL3AxisWindow *a,struct AliL3AxisWindow *b) const { - if(a->weight < b->weight) return 1; - if(a->weight > b->weight) return -1; + // Peak comparison based on peaks weight + if(a->fWeight < b->fWeight) return 1; + if(a->fWeight > b->fWeight) return -1; return 0; } @@ -1234,8 +1238,8 @@ void AliL3HoughMaxFinder::FindPeak(Int_t t1,Double_t t2,Int_t t3) Int_t nbinsx = hist->GetNbinsX()+1; Int_t *m = new Int_t[nbinsx]; - Int_t *m_low = new Int_t[nbinsx]; - Int_t *m_up = new Int_t[nbinsx]; + Int_t *mlow = new Int_t[nbinsx]; + Int_t *mup = new Int_t[nbinsx]; recompute: //this is a goto. @@ -1243,11 +1247,11 @@ void AliL3HoughMaxFinder::FindPeak(Int_t t1,Double_t t2,Int_t t3) for(Int_t i=0; i m[xbin]) //Max value locally in this xbin { m[xbin]=sum; - m_low[xbin]=ybin; - m_up[xbin]=ybin + t1; + mlow[xbin]=ybin; + mup[xbin]=ybin + t1; } } - if(m[xbin] > max_x) //Max value globally in x-direction + if(m[xbin] > maxx) //Max value globally in x-direction { - max_xbin = xbin; - max_x = m[xbin];//sum; + maxxbin = xbin; + maxx = m[xbin];//sum; } } - //printf("max_xbin %d max_x %d m_low %d m_up %d\n",max_xbin,max_x,m_low[max_xbin],m_up[max_xbin]); - //printf("ylow %f yup %f\n",hist->GetBinCenterY(m_low[max_xbin]),hist->GetBinCenterY(m_up[max_xbin])); + //printf("maxxbin %d maxx %d mlow %d mup %d\n",maxxbin,maxx,mlow[maxxbin],mup[maxxbin]); + //printf("ylow %f yup %f\n",hist->GetBinCenterY(mlow[maxxbin]),hist->GetBinCenterY(mup[maxxbin])); //Determine a width in the x-direction - Int_t x_low=0,x_up=0; + Int_t xlow=0,xup=0; - for(Int_t xbin=max_xbin-1; xbin >= xmin; xbin--) + for(Int_t xbin=maxxbin-1; xbin >= xmin; xbin--) { - if(m[xbin] < max_x*t2) + if(m[xbin] < maxx*t2) { - x_low = xbin+1; + xlow = xbin+1; break; } } - for(Int_t xbin = max_xbin+1; xbin <=xmax; xbin++) + for(Int_t xbin = maxxbin+1; xbin <=xmax; xbin++) { - if(m[xbin] < max_x*t2) + if(m[xbin] < maxx*t2) { - x_up = xbin-1; + xup = xbin-1; break; } } - Double_t top=0,butt=0,value,x_peak; - if(x_up - x_low + 1 > t3) + Double_t top=0,butt=0,value,xpeak; + if(xup - xlow + 1 > t3) { t1 -= 1; - printf("\nxrange out if limit x_up %d x_low %d t1 %d\n\n",x_low,x_up,t1); + printf("\nxrange out if limit xup %d xlow %d t1 %d\n\n",xlow,xup,t1); if(t1 > 1) goto recompute; else { - x_peak = hist->GetBinCenterX(max_xbin); + xpeak = hist->GetBinCenterX(maxxbin); goto moveon; } } - //printf("xlow %f xup %f\n",hist->GetBinCenterX(x_low),hist->GetBinCenterX(x_up)); - //printf("Spread in x %d\n",x_up-x_low +1); + //printf("xlow %f xup %f\n",hist->GetBinCenterX(xlow),hist->GetBinCenterX(xup)); + //printf("Spread in x %d\n",xup-xlow +1); //Now, calculate the center of mass in x-direction - for(Int_t xbin=x_low; xbin <= x_up; xbin++) + for(Int_t xbin=xlow; xbin <= xup; xbin++) { value = hist->GetBinCenterX(xbin); top += value*m[xbin]; butt += m[xbin]; } - x_peak = top/butt; + xpeak = top/butt; moveon: //Find the peak in y direction: - Int_t x_l = hist->FindXbin(x_peak); - if(hist->GetBinCenterX(x_l) > x_peak) - x_l--; + Int_t xl = hist->FindXbin(xpeak); + if(hist->GetBinCenterX(xl) > xpeak) + xl--; - Int_t x_u = x_l + 1; + Int_t xu = xl + 1; - if(hist->GetBinCenterX(x_l) > x_peak || hist->GetBinCenterX(x_u) <= x_peak) - printf("\nAliL3HoughMaxFinder::FindPeak : Wrong xrange %f %f %f\n\n",hist->GetBinCenterX(x_l),x_peak,hist->GetBinCenterX(x_u)); + if(hist->GetBinCenterX(xl) > xpeak || hist->GetBinCenterX(xu) <= xpeak) + printf("\nAliL3HoughMaxFinder::FindPeak : Wrong xrange %f %f %f\n\n",hist->GetBinCenterX(xl),xpeak,hist->GetBinCenterX(xu)); - //printf("\nxlow %f xup %f\n",hist->GetBinCenterX(x_l),hist->GetBinCenterX(x_u)); + //printf("\nxlow %f xup %f\n",hist->GetBinCenterX(xl),hist->GetBinCenterX(xu)); value=top=butt=0; - //printf("ylow %f yup %f\n",hist->GetBinCenterY(m_low[x_l]),hist->GetBinCenterY(m_up[x_l])); - //printf("ylow %f yup %f\n",hist->GetBinCenterY(m_low[x_u]),hist->GetBinCenterY(m_up[x_u])); + //printf("ylow %f yup %f\n",hist->GetBinCenterY(mlow[xl]),hist->GetBinCenterY(mup[xl])); + //printf("ylow %f yup %f\n",hist->GetBinCenterY(mlow[xu]),hist->GetBinCenterY(mup[xu])); - for(Int_t ybin=m_low[x_l]; ybin <= m_up[x_l]; ybin++) + for(Int_t ybin=mlow[xl]; ybin <= mup[xl]; ybin++) { value = hist->GetBinCenterY(ybin); - bin = hist->GetBin(x_l,ybin); + bin = hist->GetBin(xl,ybin); top += value*hist->GetBinContent(bin); butt += hist->GetBinContent(bin); } - Double_t y_peak_low = top/butt; + Double_t ypeaklow = top/butt; - //printf("y_peak_low %f\n",y_peak_low); + //printf("ypeaklow %f\n",ypeaklow); value=top=butt=0; - for(Int_t ybin=m_low[x_u]; ybin <= m_up[x_u]; ybin++) + for(Int_t ybin=mlow[xu]; ybin <= mup[xu]; ybin++) { value = hist->GetBinCenterY(ybin); - bin = hist->GetBin(x_u,ybin); + bin = hist->GetBin(xu,ybin); top += value*hist->GetBinContent(bin); butt += hist->GetBinContent(bin); } - Double_t y_peak_up = top/butt; + Double_t ypeakup = top/butt; - //printf("y_peak_up %f\n",y_peak_up); + //printf("ypeakup %f\n",ypeakup); - Double_t x_value_up = hist->GetBinCenterX(x_u); - Double_t x_value_low = hist->GetBinCenterX(x_l); + Double_t xvalueup = hist->GetBinCenterX(xu); + Double_t xvaluelow = hist->GetBinCenterX(xl); - Double_t y_peak = (y_peak_low*(x_value_up - x_peak) + y_peak_up*(x_peak - x_value_low))/(x_value_up - x_value_low); + Double_t ypeak = (ypeaklow*(xvalueup - xpeak) + ypeakup*(xpeak - xvaluelow))/(xvalueup - xvaluelow); //Find the weight: - //bin = hist->FindBin(x_peak,y_peak); + //bin = hist->FindBin(xpeak,ypeak); //Int_t weight = (Int_t)hist->GetBinContent(bin); //AliL3HoughTrack *track = new AliL3HoughTrack(); - //track->SetTrackParameters(x_peak,y_peak,weight); - fXPeaks[fNPeaks]=x_peak; - fYPeaks[fNPeaks]=y_peak; + //track->SetTrackParameters(xpeak,ypeak,weight); + fXPeaks[fNPeaks]=xpeak; + fYPeaks[fNPeaks]=ypeak; fWeight[fNPeaks]=(Int_t)hist->GetBinContent(bin); fNPeaks++; delete [] m; - delete [] m_low; - delete [] m_up; + delete [] mlow; + delete [] mup; //return track; } -Float_t AliL3HoughMaxFinder::GetXPeakSize(Int_t i) +Float_t AliL3HoughMaxFinder::GetXPeakSize(Int_t i) const { + // Get X size of a peak if(i<0 || i>fNMax) { STDCERR<<"AliL3HoughMaxFinder::GetXPeakSize : Invalid index "<GetBinWidthX(); - return BinWidth*(fENDXPeaks[i]-fSTARTXPeaks[i]+1); + Float_t binwidth = fCurrentHisto->GetBinWidthX(); + return binwidth*(fENDXPeaks[i]-fSTARTXPeaks[i]+1); } -Float_t AliL3HoughMaxFinder::GetYPeakSize(Int_t i) +Float_t AliL3HoughMaxFinder::GetYPeakSize(Int_t i) const { + // Get Y size of a peak if(i<0 || i>fNMax) { STDCERR<<"AliL3HoughMaxFinder::GetYPeak : Invalid index "<GetBinWidthY(); - return BinWidth*(fENDYPeaks[i]-fSTARTYPeaks[i]+1); + Float_t binwidth = fCurrentHisto->GetBinWidthY(); + return binwidth*(fENDYPeaks[i]-fSTARTYPeaks[i]+1); } diff --git a/HLT/hough/AliL3HoughMaxFinder.h b/HLT/hough/AliL3HoughMaxFinder.h index b1ff64eface..64ec6cf40dc 100644 --- a/HLT/hough/AliL3HoughMaxFinder.h +++ b/HLT/hough/AliL3HoughMaxFinder.h @@ -1,7 +1,7 @@ // @(#) $Id$ -#ifndef ALIL3_HOUGH_MaxFinder -#define ALIL3_HOUGH_MaxFinder +#ifndef ALIL3HOUGHMAXFINDER_H +#define ALIL3HOUGHMAXFINDER_H #include "AliL3RootTypes.h" #include "AliL3StandardIncludes.h" @@ -11,43 +11,15 @@ class AliL3TrackArray; class AliL3HoughTrack; class TNtuple; -struct AxisWindow +struct AliL3AxisWindow { - Int_t ymin; - Int_t ymax; - Int_t xbin; - Int_t weight; + Int_t fYmin; // min Y + Int_t fYmax; // max Y + Int_t fXbin; // X bin + Int_t fWeight; // weight }; class AliL3HoughMaxFinder { - - private: - - Int_t fThreshold; - Int_t fCurrentEtaSlice; - AliL3Histogram *fCurrentHisto; //! - - Float_t fGradX; - Float_t fGradY; - Float_t *fXPeaks; //! - Float_t *fYPeaks; //! - Int_t *fSTARTXPeaks; //! - Int_t *fSTARTYPeaks; //! - Int_t *fENDXPeaks; //! - Int_t *fENDYPeaks; //! - Int_t *fSTARTETAPeaks; //! - Int_t *fENDETAPeaks; //! - Int_t *fWeight; //! - Int_t fN1PeaksPrevEtaSlice; - Int_t fN2PeaksPrevEtaSlice; - Int_t fNPeaks; - Int_t fNMax; - - Char_t fHistoType; - -#ifndef no_root - TNtuple *fNtuppel; //! -#endif public: AliL3HoughMaxFinder(); @@ -62,14 +34,14 @@ class AliL3HoughMaxFinder { void FindAbsMaxima(); void FindBigMaxima(); void FindMaxima(Int_t threshold=0); - void FindAdaptedPeaks(Int_t nkappawindow,Float_t cut_ratio); + void FindAdaptedPeaks(Int_t nkappawindow,Float_t cutratio); //Peak finder for HoughTransformerRow void FindAdaptedRowPeaks(Int_t kappawindow,Int_t xsize,Int_t ysize); //More sophisticated peak finders: void FindPeak(Int_t t1,Double_t t2,Int_t t3); - void FindPeak1(Int_t y_window=2,Int_t x_bin_sides=1); - void SortPeaks(struct AxisWindow **a,Int_t first,Int_t last); - Int_t PeakCompare(struct AxisWindow *a,struct AxisWindow *b); + void FindPeak1(Int_t ywindow=2,Int_t xbinsides=1); + void SortPeaks(struct AliL3AxisWindow **a,Int_t first,Int_t last); + Int_t PeakCompare(struct AliL3AxisWindow *a,struct AliL3AxisWindow *b) const; //Setters: void SetGradient(Float_t x,Float_t y) {fGradX=x; fGradY=y;} @@ -78,20 +50,48 @@ class AliL3HoughMaxFinder { void SetEtaSlice(Int_t etaslice) {fCurrentEtaSlice = etaslice;} //Getters: - Float_t GetXPeak(Int_t i); - Float_t GetYPeak(Int_t i); - Float_t GetXPeakSize(Int_t i); - Float_t GetYPeakSize(Int_t i); - Int_t GetWeight(Int_t i); - Int_t GetStartEta(Int_t i); - Int_t GetEndEta(Int_t i); - Int_t GetEntries() {return fNPeaks;} + Float_t GetXPeak(Int_t i) const; + Float_t GetYPeak(Int_t i) const; + Float_t GetXPeakSize(Int_t i) const; + Float_t GetYPeakSize(Int_t i) const; + Int_t GetWeight(Int_t i) const; + Int_t GetStartEta(Int_t i) const; + Int_t GetEndEta(Int_t i) const; + Int_t GetEntries() const {return fNPeaks;} + + private: + + Int_t fThreshold; // Threshold for Peak Finder + Int_t fCurrentEtaSlice; // Current eta slice being processed + AliL3Histogram *fCurrentHisto; //! + + Float_t fGradX; // Gradient threshold inside Peak Finder + Float_t fGradY; // Gradient threshold inside Peak Finder + Float_t *fXPeaks; //! + Float_t *fYPeaks; //! + Int_t *fSTARTXPeaks; //! + Int_t *fSTARTYPeaks; //! + Int_t *fENDXPeaks; //! + Int_t *fENDYPeaks; //! + Int_t *fSTARTETAPeaks; //! + Int_t *fENDETAPeaks; //! + Int_t *fWeight; //! + Int_t fN1PeaksPrevEtaSlice; // Index of the first peak in the previous eta slice + Int_t fN2PeaksPrevEtaSlice; // Index of the last peak in the previous eta slice + Int_t fNPeaks; // Index of the last accumulated peak + Int_t fNMax; // Maximum allowed number of peaks + + Char_t fHistoType; // Histogram type + +#ifndef no_root + TNtuple *fNtuppel; //! +#endif ClassDef(AliL3HoughMaxFinder,1) //Maximum finder class }; -inline Float_t AliL3HoughMaxFinder::GetXPeak(Int_t i) +inline Float_t AliL3HoughMaxFinder::GetXPeak(Int_t i) const { if(i<0 || i>fNMax) { @@ -101,7 +101,7 @@ inline Float_t AliL3HoughMaxFinder::GetXPeak(Int_t i) return fXPeaks[i]; } -inline Float_t AliL3HoughMaxFinder::GetYPeak(Int_t i) +inline Float_t AliL3HoughMaxFinder::GetYPeak(Int_t i) const { if(i<0 || i>fNMax) { @@ -112,7 +112,7 @@ inline Float_t AliL3HoughMaxFinder::GetYPeak(Int_t i) } -inline Int_t AliL3HoughMaxFinder::GetWeight(Int_t i) +inline Int_t AliL3HoughMaxFinder::GetWeight(Int_t i) const { if(i<0 || i>fNMax) { @@ -122,7 +122,7 @@ inline Int_t AliL3HoughMaxFinder::GetWeight(Int_t i) return fWeight[i]; } -inline Int_t AliL3HoughMaxFinder::GetStartEta(Int_t i) +inline Int_t AliL3HoughMaxFinder::GetStartEta(Int_t i) const { if(i<0 || i>fNMax) { @@ -132,7 +132,7 @@ inline Int_t AliL3HoughMaxFinder::GetStartEta(Int_t i) return fSTARTETAPeaks[i]; } -inline Int_t AliL3HoughMaxFinder::GetEndEta(Int_t i) +inline Int_t AliL3HoughMaxFinder::GetEndEta(Int_t i) const { if(i<0 || i>fNMax) { diff --git a/HLT/hough/AliL3HoughTransformer.cxx b/HLT/hough/AliL3HoughTransformer.cxx index cca5afcd32a..883b05e68be 100644 --- a/HLT/hough/AliL3HoughTransformer.cxx +++ b/HLT/hough/AliL3HoughTransformer.cxx @@ -39,7 +39,7 @@ AliL3HoughTransformer::AliL3HoughTransformer() #endif } -AliL3HoughTransformer::AliL3HoughTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments,Bool_t DoEtaOverlap,Bool_t /*DoMC*/) : AliL3HoughBaseTransformer(slice,patch,n_eta_segments) +AliL3HoughTransformer::AliL3HoughTransformer(Int_t slice,Int_t patch,Int_t netasegments,Bool_t DoEtaOverlap,Bool_t /*DoMC*/) : AliL3HoughBaseTransformer(slice,patch,netasegments) { //Normal constructor fParamSpace = 0; @@ -53,6 +53,7 @@ AliL3HoughTransformer::AliL3HoughTransformer(Int_t slice,Int_t patch,Int_t n_eta AliL3HoughTransformer::~AliL3HoughTransformer() { + // Dtor DeleteHistograms(); #ifdef do_mc if(fTrackID) @@ -69,6 +70,7 @@ AliL3HoughTransformer::~AliL3HoughTransformer() void AliL3HoughTransformer::DeleteHistograms() { + // Clean up if(!fParamSpace) return; for(Int_t i=0; i= GetNEtaSegments() || eta_index < 0) + // Return a pointer to the histogram which contains etaindex eta slice + if(!fParamSpace || etaindex >= GetNEtaSegments() || etaindex < 0) return 0; - if(!fParamSpace[eta_index]) + if(!fParamSpace[etaindex]) return 0; - return fParamSpace[eta_index]; + return fParamSpace[etaindex]; } -Double_t AliL3HoughTransformer::GetEta(Int_t eta_index,Int_t /*slice*/) +Double_t AliL3HoughTransformer::GetEta(Int_t etaindex,Int_t /*slice*/) { - Double_t eta_slice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments(); + // Return eta calculated in the middle of the eta slice + Double_t etaslice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments(); Double_t eta=0; if(fEtaOverlap) { - Int_t index = eta_index + 1; - eta=(Double_t)((index)*eta_slice); + Int_t index = etaindex + 1; + eta=(Double_t)((index)*etaslice); } else - eta=(Double_t)((eta_index+0.5)*eta_slice); + eta=(Double_t)((etaindex+0.5)*etaslice); return eta; } @@ -232,9 +240,9 @@ void AliL3HoughTransformer::TransformCircle() //Transform the input data with a circle HT. //The function loops over all the data, and transforms each pixel with the equations: // - //kappa = 2/R*sin(phi - phi0) + //kappa = 2/r*sin(phi - phi0) // - //where R = sqrt(x*x +y*y), and phi = arctan(y/x) + //where r = sqrt(x*x +y*y), and phi = arctan(y/x) // //Each pixel then transforms into a curve in the (kappa,phi0)-space. In order to find //which histogram in which the pixel should be transformed, the eta-value is calculated @@ -283,21 +291,21 @@ void AliL3HoughTransformer::TransformCircle() Double_t eta = AliL3Transform::GetEta(xyz); //Get the corresponding index, which determines which histogram to fill: - Int_t eta_index = GetEtaIndex(eta); + Int_t etaindex = GetEtaIndex(eta); - if(eta_index < 0 || eta_index >= GetNEtaSegments()) + if(etaindex < 0 || etaindex >= GetNEtaSegments()) continue; //Get the correct histogrampointer: - AliL3Histogram *hist = fParamSpace[eta_index]; + AliL3Histogram *hist = fParamSpace[etaindex]; if(!hist) { - cerr<<"AliL3HoughTransformer::TransformCircle : Error getting histogram in index "<GetFirstYbin(); b<=hist->GetLastYbin(); b++) { Float_t phi0 = hist->GetBinCenterY(b); - Float_t kappa = 2*sin(phi - phi0)/R; + Float_t kappa = 2*sin(phi - phi0)/r; //hist->Fill(kappa,phi0,(int)rint(log((Float_t)charge))); hist->Fill(kappa,phi0,charge); //hist->Fill(kappa,phi0,1); @@ -319,11 +327,11 @@ void AliL3HoughTransformer::TransformCircle() if(label < 0) break; UInt_t c; for(c=0; c= AliL3Transform::GetLastRow(GetPatch())) minrow = AliL3Transform::GetFirstRow(GetPatch()); if(maxrow < AliL3Transform::GetFirstRow(GetPatch()) || maxrow >= AliL3Transform::GetLastRow(GetPatch())) @@ -388,15 +396,15 @@ void AliL3HoughTransformer::TransformCircleC(Int_t *row_range,Int_t every) } Int_t bound = (GetNEtaSegments()+1)*(AliL3Transform::GetNRows(GetPatch())+1); - EtaContainer *etaPt = new EtaContainer[bound]; - memset(etaPt,0,bound*sizeof(EtaContainer)); + AliL3EtaContainer *etaPt = new AliL3EtaContainer[bound]; + memset(etaPt,0,bound*sizeof(AliL3EtaContainer)); - Digit *digits = new Digit[counter]; - cout<<"Allocating "< 0 && index < bound) { - if(etaPt[index].first == 0) - etaPt[index].first = &digits[counter]; + if(etaPt[index].fFirst == 0) + etaPt[index].fFirst = &digits[counter]; else - (etaPt[index].last)->next = &digits[counter]; - etaPt[index].last = &digits[counter]; + (etaPt[index].fLast)->fNext = &digits[counter]; + etaPt[index].fLast = &digits[counter]; } } else { - Int_t eta_index[2]; - GetEtaIndexes(eta,eta_index); + Int_t etaindex[2]; + GetEtaIndexes(eta,etaindex); Int_t index[2]; - index[0] = (GetNEtaSegments()+1)*(i-AliL3Transform::GetFirstRow(GetPatch())) + eta_index[0]; - index[1] = (GetNEtaSegments()+1)*(i-AliL3Transform::GetFirstRow(GetPatch())) + eta_index[1]; + index[0] = (GetNEtaSegments()+1)*(i-AliL3Transform::GetFirstRow(GetPatch())) + etaindex[0]; + index[1] = (GetNEtaSegments()+1)*(i-AliL3Transform::GetFirstRow(GetPatch())) + etaindex[1]; if(index[0] == index[1]) { cerr<<"Same etaindexes "< 0 && ind < bound) { - if(etaPt[ind].first == 0) - etaPt[ind].first = &digits[counter]; + if(etaPt[ind].fFirst == 0) + etaPt[ind].fFirst = &digits[counter]; else - (etaPt[ind].last)->next = &digits[counter]; - etaPt[ind].last = &digits[counter]; + (etaPt[ind].fLast)->fNext = &digits[counter]; + etaPt[ind].fLast = &digits[counter]; } ind = index[1]; if(ind > 0 && ind < bound) { - if(etaPt[ind].first == 0) - etaPt[ind].first = &digits[counter]; + if(etaPt[ind].fFirst == 0) + etaPt[ind].fFirst = &digits[counter]; else - (etaPt[ind].last)->next = &digits[counter]; - etaPt[ind].last = &digits[counter]; + (etaPt[ind].fLast)->fNext = &digits[counter]; + etaPt[ind].fLast = &digits[counter]; } } @@ -476,7 +484,7 @@ void AliL3HoughTransformer::TransformCircleC(Int_t *row_range,Int_t every) cout<<"Doing the combinatorics"<next) + for(dPt1 = (AliL3Digit*)etaPt[index1].fFirst; dPt1 != 0; dPt1 = (AliL3Digit*)dPt1->fNext) { for(Int_t j=i+every; j<=maxrow; j+=every) { Int_t index2 = (GetNEtaSegments()+1)*(j-AliL3Transform::GetFirstRow(GetPatch())) + e; - for(dPt2 = (Digit*)etaPt[index2].first; dPt2 != 0; dPt2 = (Digit*)dPt2->next) + for(dPt2 = (AliL3Digit*)etaPt[index2].fFirst; dPt2 != 0; dPt2 = (AliL3Digit*)dPt2->fNext) { - if(dPt1->row == dPt2->row) + if(dPt1->fRow == dPt2->fRow) { cerr<<"same row; indexes "<r; - phi1 = dPt1->phi; - r2 = dPt2->r; - phi2 = dPt2->phi; - phi_0 = atan( (r2*sin(phi1)-r1*sin(phi2))/(r2*cos(phi1)-r1*cos(phi2)) ); - kappa = 2*sin(phi2-phi_0)/r2; - tot_charge = dPt1->charge + dPt2->charge; - hist->Fill(kappa,phi_0,tot_charge); + r1 = dPt1->fR; + phi1 = dPt1->fPhi; + r2 = dPt2->fR; + phi2 = dPt2->fPhi; + phi0 = atan( (r2*sin(phi1)-r1*sin(phi2))/(r2*cos(phi1)-r1*cos(phi2)) ); + kappa = 2*sin(phi2-phi0)/r2; + totcharge = dPt1->fCharge + dPt2->fCharge; + hist->Fill(kappa,phi0,totcharge); } } @@ -528,7 +536,7 @@ void AliL3HoughTransformer::TransformCircleC(Int_t *row_range,Int_t every) } -void AliL3HoughTransformer::TransformLine(Int_t *row_range,Float_t *phirange) +void AliL3HoughTransformer::TransformLine(Int_t *rowrange,Float_t *phirange) { //Do a line transform on the data. @@ -543,10 +551,10 @@ void AliL3HoughTransformer::TransformLine(Int_t *row_range,Float_t *phirange) Int_t minrow = AliL3Transform::GetFirstRow(GetPatch()); Int_t maxrow = AliL3Transform::GetLastRow(GetPatch()); - if(row_range) + if(rowrange) { - minrow = row_range[0]; - maxrow = row_range[1]; + minrow = rowrange[0]; + maxrow = rowrange[1]; if(minrow < AliL3Transform::GetFirstRow(GetPatch()) || minrow >= AliL3Transform::GetLastRow(GetPatch())) minrow = AliL3Transform::GetFirstRow(GetPatch()); if(maxrow < AliL3Transform::GetFirstRow(GetPatch()) || maxrow >= AliL3Transform::GetLastRow(GetPatch())) @@ -585,17 +593,17 @@ void AliL3HoughTransformer::TransformLine(Int_t *row_range,Float_t *phirange) continue; } Float_t eta = AliL3Transform::GetEta(xyz); - Int_t eta_index = GetEtaIndex(eta);//(Int_t)(eta/etaslice); - if(eta_index < 0 || eta_index >= GetNEtaSegments()) + Int_t etaindex = GetEtaIndex(eta);//(Int_t)(eta/etaslice); + if(etaindex < 0 || etaindex >= GetNEtaSegments()) continue; xyz[0] = xyz[0] - AliL3Transform::Row2X(minrow); //Get the correct histogram: - AliL3Histogram *hist = fParamSpace[eta_index]; + AliL3Histogram *hist = fParamSpace[etaindex]; if(!hist) { - printf("AliL3HoughTransformer::TransformLine : Error getting histogram in index %d\n",eta_index); + printf("AliL3HoughTransformer::TransformLine : Error getting histogram in index %d\n",etaindex); continue; } for(Int_t xbin=hist->GetFirstXbin(); xbinGetLastXbin(); xbin++) @@ -610,18 +618,19 @@ void AliL3HoughTransformer::TransformLine(Int_t *row_range,Float_t *phirange) } -struct LDigit { - Int_t row; - Int_t charge; - Float_t y; - LDigit *next; +struct AliL3LDigit { + Int_t fRow; // Digit rowpad + Int_t fCharge; // Digit charge + Float_t fY; // Y position of the digit in the local coor system + AliL3LDigit *fNext; // Next digit }; -struct LEtaContainer { - LDigit *first; - LDigit *last; +struct AliL3LEtaContainer { + AliL3LDigit *fFirst; //First digit + AliL3LDigit *fLast; //Last digit }; void AliL3HoughTransformer::TransformLineC(Int_t *rowrange,Float_t *phirange) { + //Circle transform ?? AliL3DigitRowData *tempPt = GetDataPointer(); if(!tempPt) LOG(AliL3Log::kError,"AliL3HoughTransformer::TransformCircleC","Data") @@ -636,12 +645,12 @@ void AliL3HoughTransformer::TransformLineC(Int_t *rowrange,Float_t *phirange) } Int_t bound = (GetNEtaSegments()+1)*(AliL3Transform::GetNRows(GetPatch())+1); - LEtaContainer *etaPt = new LEtaContainer[bound]; - memset(etaPt,0,bound*sizeof(LEtaContainer)); + AliL3LEtaContainer *etaPt = new AliL3LEtaContainer[bound]; + memset(etaPt,0,bound*sizeof(AliL3LEtaContainer)); - LDigit *digits = new LDigit[counter]; - cout<<"Allocating "< phirange[1]) continue; - digits[counter].row = i; - digits[counter].y = xyz[1]; - digits[counter].charge = charge; + digits[counter].fRow = i; + digits[counter].fY = xyz[1]; + digits[counter].fCharge = charge; - Int_t eta_index = GetEtaIndex(eta); - Int_t index = (GetNEtaSegments()+1)*(i-AliL3Transform::GetFirstRow(GetPatch())) + eta_index; + Int_t etaindex = GetEtaIndex(eta); + Int_t index = (GetNEtaSegments()+1)*(i-AliL3Transform::GetFirstRow(GetPatch())) + etaindex; if(index > 0 && index < bound) { - if(etaPt[index].first == 0) - etaPt[index].first = &digits[counter]; + if(etaPt[index].fFirst == 0) + etaPt[index].fFirst = &digits[counter]; else - (etaPt[index].last)->next = &digits[counter]; - etaPt[index].last = &digits[counter]; + (etaPt[index].fLast)->fNext = &digits[counter]; + etaPt[index].fLast = &digits[counter]; } counter++; } @@ -687,7 +696,7 @@ void AliL3HoughTransformer::TransformLineC(Int_t *rowrange,Float_t *phirange) cout<<"Doing the combinatorics"<next) + for(dPt1 = (AliL3LDigit*)etaPt[index1].fFirst; dPt1 != 0; dPt1 = (AliL3LDigit*)dPt1->fNext) { for(Int_t j=i+1; j<=rowrange[1]; j++) { Int_t index2 = (GetNEtaSegments()+1)*(j-AliL3Transform::GetFirstRow(GetPatch())) + e; - for(dPt2 = (LDigit*)etaPt[index2].first; dPt2 != 0; dPt2 = (LDigit*)dPt2->next) + for(dPt2 = (AliL3LDigit*)etaPt[index2].fFirst; dPt2 != 0; dPt2 = (AliL3LDigit*)dPt2->fNext) { - if(dPt1->row == dPt2->row) + if(dPt1->fRow == dPt2->fRow) { cerr<<"same row; indexes "<row) - AliL3Transform::Row2X(rowrange[0]); - float x2 = AliL3Transform::Row2X(dPt2->row) - AliL3Transform::Row2X(rowrange[0]); - float y1 = dPt1->y; - float y2 = dPt2->y; + float x1 = AliL3Transform::Row2X(dPt1->fRow) - AliL3Transform::Row2X(rowrange[0]); + float x2 = AliL3Transform::Row2X(dPt2->fRow) - AliL3Transform::Row2X(rowrange[0]); + float y1 = dPt1->fY; + float y2 = dPt2->fY; float theta = atan2(x2-x1,y1-y2); float rho = x1*cos(theta)+y1*sin(theta); hist->Fill(theta,rho,1);//dPt1->charge+dPt2->charge); @@ -736,8 +745,9 @@ void AliL3HoughTransformer::TransformLineC(Int_t *rowrange,Float_t *phirange) delete [] digits; } -Int_t AliL3HoughTransformer::GetTrackID(Int_t eta_index,Double_t kappa,Double_t psi) +Int_t AliL3HoughTransformer::GetTrackID(Int_t etaindex,Double_t kappa,Double_t psi) { + // Returns the MC label for a given peak found in the Hough space if(!fDoMC) { cerr<<"AliL3HoughTransformer::GetTrackID : Flag switched off"< GetNEtaSegments()) + if(etaindex < 0 || etaindex > GetNEtaSegments()) { - cerr<<"AliL3HoughTransformer::GetTrackID : Wrong etaindex "<FindBin(kappa,psi); Int_t label=-1; Int_t max=0; for(UInt_t i=0; i max) { max = nhits; - label = fTrackID[eta_index][bin].fLabel[i]; + label = fTrackID[etaindex][bin].fLabel[i]; } } //nhits = max; diff --git a/HLT/hough/AliL3HoughTransformer.h b/HLT/hough/AliL3HoughTransformer.h index d6596be889d..0c9ca856668 100644 --- a/HLT/hough/AliL3HoughTransformer.h +++ b/HLT/hough/AliL3HoughTransformer.h @@ -1,7 +1,7 @@ // @(#) $Id$ -#ifndef ALIL3_HOUGHTRANSFORMER -#define ALIL3_HOUGHTRANSFORMER +#ifndef ALIL3HOUGHTRANSFORMER_H +#define ALIL3HOUGHTRANSFORMER_H #include "AliL3RootTypes.h" #include "AliL3HoughBaseTransformer.h" @@ -9,21 +9,10 @@ class AliL3Histogram; class AliL3HoughTransformer : public AliL3HoughBaseTransformer { - - private: - - AliL3Histogram **fParamSpace; //! -#ifdef do_mc - TrackIndex **fTrackID; //! -#endif - Bool_t fDoMC; - Bool_t fEtaOverlap; - - void DeleteHistograms(); public: AliL3HoughTransformer(); - AliL3HoughTransformer(Int_t slice,Int_t patch,Int_t n_eta_segments,Bool_t DoEtaOverlap=kFALSE,Bool_t DoMC=kFALSE); + AliL3HoughTransformer(Int_t slice,Int_t patch,Int_t netasegments,Bool_t DoEtaOverlap=kFALSE,Bool_t DoMC=kFALSE); virtual ~AliL3HoughTransformer(); void CreateHistograms(Float_t ptmin,Float_t ptmax,Float_t ptres,Int_t nybin,Float_t psi); @@ -32,16 +21,27 @@ class AliL3HoughTransformer : public AliL3HoughBaseTransformer { Int_t nybin,Float_t ymin,Float_t ymax); void Reset(); void TransformCircle(); - void TransformCircleC(Int_t *row_range,Int_t every=1); + void TransformCircleC(Int_t *rowrange,Int_t every=1); void TransformLine(Int_t *rowrange=0,Float_t *phirange=0); void TransformLineC(Int_t *rowrange,Float_t *phirange); Int_t GetEtaIndex(Double_t eta); void GetEtaIndexes(Double_t eta,Int_t *indexes); - AliL3Histogram *GetHistogram(Int_t eta_index); - Double_t GetEta(Int_t eta_index,Int_t slice); - Int_t GetTrackID(Int_t eta_index,Double_t kappa,Double_t psi); + AliL3Histogram *GetHistogram(Int_t etaindex); + Double_t GetEta(Int_t etaindex,Int_t slice); + Int_t GetTrackID(Int_t etaindex,Double_t kappa,Double_t psi); + private: + + AliL3Histogram **fParamSpace; //! +#ifdef do_mc + TrackIndex **fTrackID; //! +#endif + Bool_t fDoMC; // Calculate mc labels or not + Bool_t fEtaOverlap; // Allow overlapping of eta slice or not + + void DeleteHistograms(); + ClassDef(AliL3HoughTransformer,1) //Normal Hough transformation class }; diff --git a/HLT/hough/AliL3HoughTransformerRow.cxx b/HLT/hough/AliL3HoughTransformerRow.cxx index 80b3ff759c0..3fc92943e49 100644 --- a/HLT/hough/AliL3HoughTransformerRow.cxx +++ b/HLT/hough/AliL3HoughTransformerRow.cxx @@ -26,18 +26,18 @@ using namespace std; ClassImp(AliL3HoughTransformerRow) -UChar_t **AliL3HoughTransformerRow::fRowCount = 0; -UChar_t **AliL3HoughTransformerRow::fGapCount = 0; -UChar_t **AliL3HoughTransformerRow::fCurrentRowCount = 0; +UChar_t **AliL3HoughTransformerRow::fgRowCount = 0; +UChar_t **AliL3HoughTransformerRow::fgGapCount = 0; +UChar_t **AliL3HoughTransformerRow::fgCurrentRowCount = 0; #ifdef do_mc -TrackIndex **AliL3HoughTransformerRow::fTrackID = 0; +TrackIndex **AliL3HoughTransformerRow::fgTrackID = 0; #endif -UChar_t *AliL3HoughTransformerRow::fTrackNRows = 0; -UChar_t *AliL3HoughTransformerRow::fTrackFirstRow = 0; -UChar_t *AliL3HoughTransformerRow::fTrackLastRow = 0; +UChar_t *AliL3HoughTransformerRow::fgTrackNRows = 0; +UChar_t *AliL3HoughTransformerRow::fgTrackFirstRow = 0; +UChar_t *AliL3HoughTransformerRow::fgTrackLastRow = 0; -Float_t AliL3HoughTransformerRow::fBeta1 = AliL3Transform::Row2X(79)/pow(AliL3Transform::Row2X(79),2); -Float_t AliL3HoughTransformerRow::fBeta2 = (AliL3Transform::Row2X(158)+6.0)/pow((AliL3Transform::Row2X(158)+6.0),2); +Float_t AliL3HoughTransformerRow::fgBeta1 = AliL3Transform::Row2X(79)/pow(AliL3Transform::Row2X(79),2); +Float_t AliL3HoughTransformerRow::fgBeta2 = (AliL3Transform::Row2X(158)+6.0)/pow((AliL3Transform::Row2X(158)+6.0),2); AliL3HoughTransformerRow::AliL3HoughTransformerRow() { @@ -51,7 +51,7 @@ AliL3HoughTransformerRow::AliL3HoughTransformerRow() fLUTbackwardZ2=0; } -AliL3HoughTransformerRow::AliL3HoughTransformerRow(Int_t slice,Int_t patch,Int_t n_eta_segments,Bool_t /*DoMC*/,Float_t zvertex) : AliL3HoughBaseTransformer(slice,patch,n_eta_segments,zvertex) +AliL3HoughTransformerRow::AliL3HoughTransformerRow(Int_t slice,Int_t patch,Int_t netasegments,Bool_t /*DoMC*/,Float_t zvertex) : AliL3HoughBaseTransformer(slice,patch,netasegments,zvertex) { //Normal constructor fParamSpace = 0; @@ -69,69 +69,71 @@ AliL3HoughTransformerRow::AliL3HoughTransformerRow(Int_t slice,Int_t patch,Int_t AliL3HoughTransformerRow::~AliL3HoughTransformerRow() { + //Destructor DeleteHistograms(); #ifdef do_mc - if(fTrackID) + if(fgTrackID) { for(Int_t i=0; iGetNbinsX()+3)/2; Int_t ncellsy = (hist->GetNbinsY()+3)/2; Int_t ncells = ncellsx*ncellsy; - if(!fTrackID) + if(!fgTrackID) { LOG(AliL3Log::kInformational,"AliL3HoughTransformerRow::CreateHistograms()","") - <<"Transformer: Allocating "<GetNbinsX()+2)*(hist->GetNbinsY()+2); - if(!fRowCount) + if(!fgRowCount) { LOG(AliL3Log::kInformational,"AliL3HoughTransformerRow::CreateHistograms()","") - <<"Transformer: Allocating "<GetFirstXbin(); @@ -251,17 +256,17 @@ void AliL3HoughTransformerRow::CreateHistograms(Int_t nxbin,Float_t xmin,Float_t { //cvetan: we get strange warning on gcc-2.95 //warning: large integer implicitly truncated to unsigned type - fTrackNRows[xbin + ybin*nxbins] = 99999; + fgTrackNRows[xbin + ybin*nxbins] = 99999; for(Int_t deltay = -1; deltay <= 1; deltay += 2) { for(Int_t deltax = -1; deltax <= 1; deltax += 2) { Float_t xtrack = hist->GetPreciseBinCenterX((Float_t)xbin+0.5*(Float_t)deltax); Float_t ytrack = hist->GetPreciseBinCenterY((Float_t)ybin+0.5*(Float_t)deltay); - Float_t psi = atan((xtrack-ytrack)/(fBeta1-fBeta2)); - Float_t kappa = 2.0*(xtrack*cos(psi)-fBeta1*sin(psi)); + Float_t psi = atan((xtrack-ytrack)/(fgBeta1-fgBeta2)); + Float_t kappa = 2.0*(xtrack*cos(psi)-fgBeta1*sin(psi)); track.SetTrackParameters(kappa,psi,1); - Bool_t first_row = kFALSE; + Bool_t firstrow = kFALSE; UInt_t maxfirstrow = 0; UInt_t maxlastrow = 0; UInt_t curfirstrow = 0; @@ -273,15 +278,15 @@ void AliL3HoughTransformerRow::CreateHistograms(Int_t nxbin,Float_t xmin,Float_t AliL3Transform::LocHLT2Raw(hit,0,j); if(hit[1]>=0 && hit[1]= (maxlastrow-maxfirstrow)) { maxfirstrow = curfirstrow; maxlastrow = curlastrow; @@ -293,14 +298,14 @@ void AliL3HoughTransformerRow::CreateHistograms(Int_t nxbin,Float_t xmin,Float_t maxfirstrow = curfirstrow; maxlastrow = curlastrow; } - if((maxlastrow-maxfirstrow) < fTrackNRows[xbin + ybin*nxbins]) { - fTrackNRows[xbin + ybin*nxbins] = maxlastrow-maxfirstrow; - fTrackFirstRow[xbin + ybin*nxbins] = maxfirstrow; - fTrackLastRow[xbin + ybin*nxbins] = maxlastrow; + if((maxlastrow-maxfirstrow) < fgTrackNRows[xbin + ybin*nxbins]) { + fgTrackNRows[xbin + ybin*nxbins] = maxlastrow-maxfirstrow; + fgTrackFirstRow[xbin + ybin*nxbins] = maxfirstrow; + fgTrackLastRow[xbin + ybin*nxbins] = maxlastrow; } } } - // cout<<" fTrackNRows "<GetNbinsY()+3)/2; Int_t ncells = ncellsx*ncellsy; for(Int_t i=0; iGetNbinsX()+2)*(hist->GetNbinsY()+2); for(Int_t i=0; i= GetNEtaSegments() || eta_index < 0) + // Return a pointer to the histogram which contains etaindex eta slice + if(!fParamSpace || etaindex >= GetNEtaSegments() || etaindex < 0) return 0; - if(!fParamSpace[eta_index]) + if(!fParamSpace[etaindex]) return 0; - return fParamSpace[eta_index]; + return fParamSpace[etaindex]; } -Double_t AliL3HoughTransformerRow::GetEta(Int_t eta_index,Int_t /*slice*/) +Double_t AliL3HoughTransformerRow::GetEta(Int_t etaindex,Int_t /*slice*/) { - Double_t eta_slice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments(); + // Return eta calculated in the middle of the eta slice + Double_t etaslice = (GetEtaMax()-GetEtaMin())/GetNEtaSegments(); Double_t eta=0; - eta=(Double_t)((eta_index+0.5)*eta_slice); + eta=(Double_t)((etaindex+0.5)*etaslice); return eta; } -struct EtaRow { - UChar_t start_pad; - UChar_t end_pad; - Bool_t found; - Float_t start_y; +struct AliL3EtaRow { + UChar_t fStartPad; //First pad in the cluster + UChar_t fEndPad; //Last pad in the cluster + Bool_t fIsFound; //Is the cluster already found + Float_t fStartY; //Y position of the first pad in the cluster #ifdef do_mc - Int_t mc_labels[MaxTrack]; + Int_t fMcLabels[MaxTrack]; //Array to store mc labels inside cluster #endif }; void AliL3HoughTransformerRow::TransformCircle() { - Float_t beta1 = fBeta1; - Float_t beta2 = fBeta2; - Float_t beta1_minus_beta2 = fBeta1 - fBeta2; + //Do the Hough Transform + Float_t beta1 = fgBeta1; + Float_t beta2 = fgBeta2; + Float_t beta1minusbeta2 = fgBeta1 - fgBeta2; - Int_t n_eta_segments = GetNEtaSegments(); - Double_t eta_min = GetEtaMin(); - Double_t eta_slice = (GetEtaMax() - eta_min)/n_eta_segments; + Int_t netasegments = GetNEtaSegments(); + Double_t etamin = GetEtaMin(); + Double_t etaslice = (GetEtaMax() - etamin)/netasegments; - Int_t lower_threshold = GetLowerThreshold(); + Int_t lowerthreshold = GetLowerThreshold(); //Assumes that all the etaslice histos are the same! AliL3Histogram *h = fParamSpace[0]; - Float_t y_min = h->GetYmin(); + Float_t ymin = h->GetYmin(); //Float_t y_max = h->GetYmax(); - Float_t hist_bin = h->GetBinWidthY(); - Int_t first_bin = h->GetFirstYbin(); - Int_t last_bin = h->GetLastYbin(); - Float_t x_min = h->GetXmin(); - Float_t x_max = h->GetXmax(); - Float_t x_bin = (x_max-x_min)/h->GetNbinsX(); - Int_t first_binx = h->GetFirstXbin()-1; - Int_t last_binx = h->GetLastXbin()+1; + Float_t histbin = h->GetBinWidthY(); + Int_t firstbin = h->GetFirstYbin(); + Int_t lastbin = h->GetLastYbin(); + Float_t xmin = h->GetXmin(); + Float_t xmax = h->GetXmax(); + Float_t xbin = (xmax-xmin)/h->GetNbinsX(); + Int_t firstbinx = h->GetFirstXbin()-1; + Int_t lastbinx = h->GetLastXbin()+1; Int_t nbinx = h->GetNbinsX()+2; - UChar_t last_pad; - Int_t last_eta_index; - EtaRow *eta_clust = new EtaRow[n_eta_segments]; + UChar_t lastpad; + Int_t lastetaindex; + AliL3EtaRow *etaclust = new AliL3EtaRow[netasegments]; AliL3DigitRowData *tempPt = GetDataPointer(); if(!tempPt) @@ -434,7 +442,7 @@ void AliL3HoughTransformerRow::TransformCircle() } Int_t ipatch = GetPatch(); - Double_t pad_pitch = AliL3Transform::GetPadPitchWidth(ipatch); + Double_t padpitch = AliL3Transform::GetPadPitchWidth(ipatch); Int_t islice = GetSlice(); Float_t *fLUTz; Float_t *fLUTz2; @@ -452,9 +460,9 @@ void AliL3HoughTransformerRow::TransformCircle() { Int_t npads = AliL3Transform::GetNPads((Int_t)i)-1; - last_pad = 255; + lastpad = 255; //Flush eta clusters array - memset(eta_clust,0,n_eta_segments*sizeof(EtaRow)); + memset(etaclust,0,netasegments*sizeof(AliL3EtaRow)); Float_t x = AliL3Transform::Row2X((Int_t)i); Float_t x2 = x*x; @@ -473,17 +481,17 @@ void AliL3HoughTransformerRow::TransformCircle() for(UInt_t j=0; jfNDigit; j++) { UShort_t charge = digPt[j].fCharge; - if((Int_t)charge <= lower_threshold) + if((Int_t)charge <= lowerthreshold) continue; UChar_t pad = digPt[j].fPad; UShort_t time = digPt[j].fTime; - if(pad != last_pad) + if(pad != lastpad) { - y = (pad-0.5*npads)*pad_pitch; + y = (pad-0.5*npads)*padpitch; Float_t y2 = y*y; r2 = x2 + y2; - last_eta_index = -1; + lastetaindex = -1; } //Transform data to local cartesian coordinates: @@ -493,22 +501,22 @@ void AliL3HoughTransformerRow::TransformCircle() Double_t r = sqrt(r2+z2); Double_t eta = 0.5 * log((r+z)/(r-z)); //Get the corresponding index, which determines which histogram to fill: - Int_t eta_index = (Int_t)((eta-eta_min)/eta_slice); + Int_t etaindex = (Int_t)((eta-etamin)/etaslice); #ifndef do_mc - if(eta_index == last_eta_index) continue; + if(etaindex == lastetaindex) continue; #endif - // cout<<" Digit at patch "<= n_eta_segments) + if(etaindex < 0 || etaindex >= netasegments) continue; - if(!eta_clust[eta_index].found) + if(!etaclust[etaindex].fIsFound) { - eta_clust[eta_index].start_pad = pad; - eta_clust[eta_index].end_pad = pad; - eta_clust[eta_index].start_y = y - pad_pitch/2.0; - eta_clust[eta_index].found = 1; + etaclust[etaindex].fStartPad = pad; + etaclust[etaindex].fEndPad = pad; + etaclust[etaindex].fStartY = y - padpitch/2.0; + etaclust[etaindex].fIsFound = 1; #ifdef do_mc if(fDoMC) { @@ -518,10 +526,10 @@ void AliL3HoughTransformerRow::TransformCircle() if(label < 0) break; UInt_t c; for(c=0; c<(MaxTrack-1); c++) - if(eta_clust[eta_index].mc_labels[c] == label || eta_clust[eta_index].mc_labels[c] == 0) + if(etaclust[etaindex].fMcLabels[c] == label || etaclust[etaindex].fMcLabels[c] == 0) break; // if(c == MaxTrack-1) cerr<<"AliL3HoughTransformer::TransformCircle : Cluster array reached maximum!! "<=0) - kappa2 += delta_kappa2; + Float_t starty = etaclust[etaindex].fStartY; + if(etaclust[etaindex].fStartPad == 0) + starty -= 0.0*padpitch; + Float_t r1 = x2 + starty*starty; + Float_t xoverr1 = x/r1; + Float_t startyoverr1 = starty/r1; + Float_t endy = (etaclust[etaindex].fEndPad-0.5*(npads-1))*padpitch; + if(etaclust[etaindex].fEndPad == npads) + endy += 0.0*padpitch; + Float_t r2 = x2 + endy*endy; + Float_t xoverr2 = x/r2; + Float_t endyoverr2 = endy/r2; + Float_t a1 = beta1minusbeta2/(xoverr1-beta2); + Float_t b1 = (xoverr1-beta1)/(xoverr1-beta2); + Float_t a2 = beta1minusbeta2/(xoverr2-beta2); + Float_t b2 = (xoverr2-beta1)/(xoverr2-beta2); + + Float_t kappa1 = (a1*startyoverr1+b1*ymin-xmin)/xbin; + Float_t deltaalpha1 = b1*histbin/xbin; + if(b1<0) + kappa1 += deltaalpha1; + Float_t kappa2 = (a2*endyoverr2+b2*ymin-xmin)/xbin; + Float_t deltaalpha2 = b2*histbin/xbin; + if(b2>=0) + kappa2 += deltaalpha2; //Fill the histogram along the phirange - for(Int_t b=first_bin; b<=last_bin; b++, kappa1 += delta_kappa1, kappa2 += delta_kappa2) + for(Int_t b=firstbin; b<=lastbin; b++, kappa1 += deltaalpha1, kappa2 += deltaalpha2) { Int_t binx1 = 1 + (Int_t)kappa1; - if(binx1>last_binx) continue; - if(binx1lastbinx) continue; + if(binx1last_binx) binx2 = last_binx; + if(binx2lastbinx) binx2 = lastbinx; #ifdef do_mc if(binx2=0) - kappa2 += delta_kappa2; + Float_t starty = etaclust[etaindex].fStartY; + if(etaclust[etaindex].fStartPad == 0) + starty -= 0.0*padpitch; + Float_t r1 = x2 + starty*starty; + Float_t xoverr1 = x/r1; + Float_t startyoverr1 = starty/r1; + Float_t endy = (etaclust[etaindex].fEndPad-0.5*(npads-1))*padpitch; + if(etaclust[etaindex].fEndPad == npads) + endy += 0.0*padpitch; + Float_t r2 = x2 + endy*endy; + Float_t xoverr2 = x/r2; + Float_t endyoverr2 = endy/r2; + Float_t a1 = beta1minusbeta2/(xoverr1-beta2); + Float_t b1 = (xoverr1-beta1)/(xoverr1-beta2); + Float_t a2 = beta1minusbeta2/(xoverr2-beta2); + Float_t b2 = (xoverr2-beta1)/(xoverr2-beta2); + + Float_t kappa1 = (a1*startyoverr1+b1*ymin-xmin)/xbin; + Float_t deltaalpha1 = b1*histbin/xbin; + if(b1<0) + kappa1 += deltaalpha1; + Float_t kappa2 = (a2*endyoverr2+b2*ymin-xmin)/xbin; + Float_t deltaalpha2 = b2*histbin/xbin; + if(b2>=0) + kappa2 += deltaalpha2; //Fill the histogram along the phirange - for(Int_t b=first_bin; b<=last_bin; b++, kappa1 += delta_kappa1, kappa2 += delta_kappa2) + for(Int_t b=firstbin; b<=lastbin; b++, kappa1 += deltaalpha1, kappa2 += deltaalpha2) { Int_t binx1 = 1 + (Int_t)kappa1; - if(binx1>last_binx) continue; - if(binx1lastbinx) continue; + if(binx1last_binx) binx2 = last_binx; + if(binx2lastbinx) binx2 = lastbinx; #ifdef do_mc if(binx2 GetNEtaSegments()) + if(etaindex < 0 || etaindex > GetNEtaSegments()) { LOG(AliL3Log::kWarning,"AliL3HoughTransformerRow::GetTrackID","Data") - <<"Wrong etaindex "<FindLabelBin(kappa,psi); if(bin==-1) { LOG(AliL3Log::kWarning,"AliL3HoughTransformerRow::GetTrackID()","") @@ -808,30 +817,32 @@ Int_t AliL3HoughTransformerRow::GetTrackID(Int_t eta_index,Double_t kappa,Double Int_t max=0; for(UInt_t i=0; i<(MaxTrack-1); i++) { - Int_t nhits=fTrackID[eta_index][bin].fNHits[i]; + Int_t nhits=fgTrackID[etaindex][bin].fNHits[i]; if(nhits == 0) break; if(nhits > max) { max = nhits; - label = fTrackID[eta_index][bin].fLabel[i]; + label = fgTrackID[etaindex][bin].fLabel[i]; } } Int_t label2=-1; Int_t max2=0; for(UInt_t i=0; i<(MaxTrack-1); i++) { - Int_t nhits=fTrackID[eta_index][bin].fNHits[i]; + Int_t nhits=fgTrackID[etaindex][bin].fNHits[i]; if(nhits == 0) break; if(nhits > max2) { - if(fTrackID[eta_index][bin].fLabel[i]!=label) { + if(fgTrackID[etaindex][bin].fLabel[i]!=label) { max2 = nhits; - label2 = fTrackID[eta_index][bin].fLabel[i]; + label2 = fgTrackID[etaindex][bin].fLabel[i]; } } } - LOG(AliL3Log::kDebug,"AliL3HoughTransformerRow::GetTrackID()","") - <<" TrackID"<<" label "<