which included commits to RCS files with non-trunk default branches.
--- /dev/null
+#ifndef _ALIL3DEFS_H_
+#define _ALIL3DEFS_H_
+
+#include <Rtypes.h>
+
+const Int_t NRowsSlice = 174;
+const Int_t NRows[5][2] = {{ 0, 45},{46,77},{78,109},{110,141},{142,173}};
+const Double_t Pi = 3.14159265358979323846;
+const Double_t ToRad = Pi/180.;
+
+
+#endif /* _ALIL3DEFS_H_ */
--- /dev/null
+
+#include "AliL3Histogram.h"
+
+ClassImp(AliL3Histogram)
+
+
+AliL3Histogram::AliL3Histogram()
+{
+
+}
+
+
+AliL3Histogram::~AliL3Histogram()
+{
+ //Destructor
+
+}
--- /dev/null
+#ifndef ALIL3_HOUGHHISTOGRAM
+#define ALIL3_HOUGHHISTOGRAM
+
+#include "AliL3RootTypes.h"
+
+class AliL3HoughTransformer : public TH2F {
+
+ private:
+
+
+
+ public:
+ AliL3Histogram();
+
+
+ ClassDef(AliL3Histogram,1)
+
+};
+
+#endif
--- /dev/null
+#include <TClonesArray.h>
+#include <TH2.h>
+#include <TH1.h>
+#include <TFile.h>
+#include <AliRun.h>
+#include <TParticle.h>
+
+#include "AliL3TrackArray.h"
+#include "AliL3Transform.h"
+#include "AliL3HoughTransformer.h"
+#include "AliL3Defs.h"
+#include "AliL3HoughTrack.h"
+#include "AliL3HoughEval.h"
+
+ClassImp(AliL3HoughEval)
+
+
+AliL3HoughEval::AliL3HoughEval()
+{
+ //Default constructor
+ fTransform = new AliL3Transform();
+ fHoughTransformer = NULL;
+}
+
+
+AliL3HoughEval::AliL3HoughEval(AliL3HoughTransformer *transformer)
+{
+ //Constructor
+ fHoughTransformer = transformer;
+ fTransform = new AliL3Transform();
+}
+
+
+AliL3HoughEval::~AliL3HoughEval()
+{
+ //Destructor
+ if(fTransform)
+ delete fTransform;
+
+}
+
+
+void AliL3HoughEval::LookInsideRoad(AliL3TrackArray *tracks,TH2F *hist,TH2F *fake)
+{
+ //Look at rawdata along the road specified by the track candidates.
+
+ if(!fHoughTransformer)
+ {
+ printf("AliL3HoughEval: No transformer object\n");
+ return;
+ }
+
+ AliL3Digits *pixel;
+ Int_t patch = fHoughTransformer->fPatch;
+ Int_t slice = fHoughTransformer->fSlice;
+ Int_t num_of_pads_to_look = 1;
+ Int_t rows_to_miss = 2;
+
+ for(Int_t i=0; i<tracks->GetNTracks(); i++)
+ {
+ AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(i);
+ if(!track) {printf("No track\n"); return;}
+
+ Int_t nrow = 0;
+
+ for(Int_t padrow = NRows[patch][0]; padrow <= NRows[patch][1]; padrow++)
+ {
+
+ Float_t xyz[3];
+ //Get the crossing point of the track and current padrow:
+ //if(!track->GetCrossingPoint(slice,padrow,xyz))
+ if(!track->GetCrossingPoint(padrow,xyz))
+ {
+ printf("AliL3HoughEval::LookInsideRoad : Track does not cross line!!\n");
+ continue;
+ }
+
+ if(fake)
+ fake->Fill(xyz[0],xyz[1],1);
+ Int_t npixs = 0;
+
+ //Get the pixels along the track candidate
+ for(pixel=(AliL3Digits*)fHoughTransformer->fRowContainer[padrow].first; pixel!=0; pixel=(AliL3Digits*)pixel->nextRowPixel)
+ {
+
+ Int_t sector,row;
+ Float_t xyz_pix[3];
+ fTransform->Slice2Sector(slice,padrow,sector,row);
+ fTransform->Raw2Local(xyz_pix,sector,row,pixel->fPad,pixel->fTime); //y alone pad
+
+
+ //check if we are inside road
+ if(fabs(xyz_pix[1] - xyz[1]) > num_of_pads_to_look*fTransform->GetPadPitchWidthLow()) continue;
+ npixs++;
+
+ if(hist)
+ {
+ //fTransform->Local2Global(xyz_pix,slice);
+ hist->Fill(xyz_pix[0],xyz_pix[1],pixel->fCharge);
+ }
+
+ }
+ if(npixs > 0)
+ {
+ nrow++;
+ }
+ }
+
+ if(nrow < NRows[patch][1]-NRows[patch][0]-rows_to_miss)
+ tracks->Remove(i);
+
+ }
+
+ tracks->Compress();
+
+
+}
+
+void AliL3HoughEval::DisplaySlice(TH2F *hist)
+{
+
+ AliL3Digits *pixel;
+
+ for(Int_t padrow = 0; padrow < 174; padrow++)
+ {
+
+ for(pixel=(AliL3Digits*)fHoughTransformer->fRowContainer[padrow].first; pixel!=0; pixel=(AliL3Digits*)pixel->nextRowPixel)
+ {
+
+ Int_t sector,row;
+ Float_t xyz_pix[3];
+ fTransform->Slice2Sector(fHoughTransformer->fSlice,padrow,sector,row);
+ fTransform->Raw2Local(xyz_pix,sector,row,pixel->fPad,pixel->fTime); //y alone pad
+
+ if(hist)
+ {
+ hist->Fill(xyz_pix[0],xyz_pix[1],pixel->fCharge);
+ }
+
+ }
+ }
+
+}
+
+
+void AliL3HoughEval::CompareMC(Char_t *rootfile,AliL3TrackArray *merged_tracks,Float_t *eta)
+{
+
+ Int_t slice = fSlice;
+
+ TFile *file = new TFile(rootfile);
+ file->cd();
+
+ AliRun *gAlice = (AliRun*)file->Get("gAlice");
+ gAlice->GetEvent(0);
+
+ TClonesArray *particles=gAlice->Particles();
+ Int_t n=particles->GetEntriesFast();
+ Float_t torad=TMath::Pi()/180;
+ Float_t phi_min = slice*20 - 10;
+ Float_t phi_max = slice*20 + 10;
+
+ for (Int_t j=0; j<n; j++) {
+ TParticle *p=(TParticle*)particles->UncheckedAt(j);
+ if (p->GetFirstMother()>=0) continue; //secondary particle
+ if(p->Eta() < eta[0] || p->Eta() > eta[1]) continue;
+ Double_t ptg=p->Pt(),pxg=p->Px(),pyg=p->Py();//,pzg=p->Pz();
+ Double_t phi_part = TMath::ATan2(pyg,pxg);
+ if (phi_part < 0) phi_part += 2*TMath::Pi();
+
+ if(phi_part < phi_min*torad || phi_part > phi_max*torad) {continue;}
+ if(ptg<0.100) continue;
+
+ AliL3HoughTrack *sel_track=0;
+
+ Double_t min_dist = 10000;
+
+ for(Int_t t=0; t<merged_tracks->GetNTracks(); t++)
+ {
+ AliL3HoughTrack *track = (AliL3HoughTrack*)merged_tracks->GetCheckedTrack(t);
+ if(!track) {printf("AliL3HoughEval: NO TRACK\n"); continue;}
+ Float_t phi0[2] = {track->GetPhi0(),0};
+ fTransform->Local2GlobalAngle(phi0,slice);
+ Double_t dPt = ptg - track->GetPt();
+ Double_t dPhi0 = phi_part - phi0[0];
+ Double_t distance = pow(dPt,2) + pow(dPhi0,2);
+ if(distance < min_dist)
+ {
+ min_dist = distance;
+ sel_track = track;
+ }
+
+ }
+ if(sel_track)
+ sel_track->SetBestMCid(j,min_dist);
+
+ Double_t dpt = fabs(ptg-sel_track->GetPt());
+ Double_t dphi0 = fabs(phi_part-sel_track->GetPhi0());
+ //printf("Found match, min_dist %f dPt %f dPhi0 %f\n",min_dist,dpt,dphi0);
+ }
+ file->Close();
+ delete file;
+
+
+}
--- /dev/null
+#ifndef ALIL3_HOUGH_Eval
+#define ALIL3_HOUGH_Eval
+
+#include "AliL3RootTypes.h"
+
+class TClonesArray;
+class TObjArray;
+class TH2F;
+class TH1F;
+class AliL3Transform;
+class AliL3TrackArray;
+class AliL3HoughTransformer;
+
+class AliL3HoughEval : public TObject {
+
+ private:
+
+ Int_t fSlice;
+ AliL3HoughTransformer *fHoughTransformer;
+ AliL3Transform *fTransform; //!
+
+ public:
+ AliL3HoughEval();
+ AliL3HoughEval(AliL3HoughTransformer *transformer);
+ virtual ~AliL3HoughEval();
+
+ void LookInsideRoad(AliL3TrackArray *tracks,TH2F *hist=0,TH2F *fake=0);
+ void DisplaySlice(TH2F *hist);
+ void CompareMC(Char_t *rootfile,AliL3TrackArray *merged_tracks,Float_t *eta);
+
+ void SetTransformer(AliL3HoughTransformer *t) {fHoughTransformer=t;}
+
+ ClassDef(AliL3HoughEval,1)
+
+};
+
+#endif
--- /dev/null
+#ifdef __CINT__
+
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+
+#pragma link C++ class AliL3HoughTransformer;
+#pragma link C++ class AliL3HoughPixel;
+//#pragma link C++ class AliL3HoughTrack;
+#pragma link C++ class AliL3HoughMaxFinder;
+#pragma link C++ class AliL3HoughEval;
+#pragma link C++ class AliL3HoughMerge;
+
+#endif
+
--- /dev/null
+#include <string.h>
+
+#include <TH2.h>
+
+#include "AliL3TrackArray.h"
+#include "AliL3HoughTrack.h"
+#include "AliL3HoughMaxFinder.h"
+
+ClassImp(AliL3HoughMaxFinder)
+
+
+AliL3HoughMaxFinder::AliL3HoughMaxFinder()
+{
+ //Default constructor
+ fThreshold = 0;
+ //fTracks = 0;
+ fHistoType=0;
+}
+
+
+AliL3HoughMaxFinder::AliL3HoughMaxFinder(Char_t *histotype)
+{
+ //Constructor
+
+ //fTracks = new AliL3TrackArray("AliL3HoughTrack");
+ if(strcmp(histotype,"KappaPhi")==0) fHistoType='c';
+ if(strcmp(histotype,"DPsi")==0) fHistoType='l';
+
+}
+
+
+AliL3HoughMaxFinder::~AliL3HoughMaxFinder()
+{
+ //Destructor
+
+}
+
+
+AliL3TrackArray *AliL3HoughMaxFinder::FindMaxima(TH2F *hist,Int_t *rowrange,Int_t ref_row)
+{
+ //Locate all the maxima in input histogram.
+ //Maxima is defined as bins with more entries than the
+ //immediately neighbouring bins.
+
+ Int_t xmin = hist->GetXaxis()->GetFirst();
+ Int_t xmax = hist->GetXaxis()->GetLast();
+ Int_t ymin = hist->GetYaxis()->GetFirst();
+ Int_t ymax = hist->GetYaxis()->GetLast();
+ Int_t bin[9],track_counter=0;
+ Stat_t value[9];
+
+ AliL3TrackArray *tracks = new AliL3TrackArray("AliL3HoughTrack");
+ AliL3HoughTrack *track;
+
+ for(Int_t xbin=xmin+1; xbin<xmax-1; xbin++)
+ {
+ for(Int_t ybin=ymin+1; ybin<ymax-1; ybin++)
+ {
+ bin[0] = hist->GetBin(xbin-1,ybin-1);
+ bin[1] = hist->GetBin(xbin,ybin-1);
+ bin[2] = hist->GetBin(xbin+1,ybin-1);
+ bin[3] = hist->GetBin(xbin-1,ybin);
+ bin[4] = hist->GetBin(xbin,ybin);
+ bin[5] = hist->GetBin(xbin+1,ybin);
+ bin[6] = hist->GetBin(xbin-1,ybin+1);
+ bin[7] = hist->GetBin(xbin,ybin+1);
+ bin[8] = hist->GetBin(xbin+1,ybin+1);
+ value[0] = hist->GetBinContent(bin[0]);
+ value[1] = hist->GetBinContent(bin[1]);
+ value[2] = hist->GetBinContent(bin[2]);
+ value[3] = hist->GetBinContent(bin[3]);
+ value[4] = hist->GetBinContent(bin[4]);
+ value[5] = hist->GetBinContent(bin[5]);
+ value[6] = hist->GetBinContent(bin[6]);
+ value[7] = hist->GetBinContent(bin[7]);
+ value[8] = hist->GetBinContent(bin[8]);
+
+ if(value[4] <= fThreshold) continue;//central bin below threshold
+
+ if(value[4]>value[0] && value[4]>value[1] && value[4]>value[2]
+ && value[4]>value[3] && value[4]>value[5] && value[4]>value[6]
+ && value[4]>value[7] && value[4]>value[8])
+ {
+ //Found a local maxima
+ Float_t max_x = hist->GetXaxis()->GetBinCenter(xbin);
+ Float_t max_y = hist->GetYaxis()->GetBinCenter(ybin);
+
+ track = (AliL3HoughTrack*)tracks->NextTrack();
+
+
+ switch(fHistoType)
+ {
+ case 'c':
+ track->SetTrackParameters(max_x,max_y,(Int_t)value[4]);
+ break;
+ case 'l':
+ track->SetLineParameters(max_x,max_y,(Int_t)value[4],rowrange,ref_row);
+ break;
+ default:
+ printf("AliL3HoughMaxFinder: Error in tracktype\n");
+ }
+
+ track_counter++;
+
+
+ }
+ else
+ continue; //not a maxima
+ }
+ }
+ tracks->QSort();
+ return tracks;
+}
--- /dev/null
+#ifndef ALIL3_HOUGH_MaxFinder
+#define ALIL3_HOUGH_MaxFinder
+
+#include "AliL3RootTypes.h"
+
+class AliL3TrackArray;
+
+class AliL3HoughMaxFinder : public TObject {
+
+ private:
+
+ Int_t fThreshold;
+ //AliL3TrackArray *fTracks; //!
+
+ Char_t fHistoType;
+
+ public:
+ AliL3HoughMaxFinder();
+ AliL3HoughMaxFinder(Char_t *histotype);
+ virtual ~AliL3HoughMaxFinder();
+
+ AliL3TrackArray *FindMaxima(TH2F *hist,Int_t *rowrange=0,Int_t ref_row=0);
+ void SetThreshold(Int_t f) {fThreshold = f;}
+
+ //AliL3TrackArray *GetTracks() {return fTracks;}
+
+ ClassDef(AliL3HoughMaxFinder,1)
+
+};
+
+#endif
--- /dev/null
+#include <TH2.h>
+
+#include "AliL3TrackArray.h"
+#include "AliL3HoughTrack.h"
+#include "AliL3HoughMerge.h"
+#include "AliL3HoughTransformer.h"
+
+ClassImp(AliL3HoughMerge)
+
+
+AliL3HoughMerge::AliL3HoughMerge()
+{
+ //Default constructor
+
+ fInTracks = NULL;
+ fOutTracks = NULL;
+ fNPatches = 5;
+}
+
+
+AliL3HoughMerge::AliL3HoughMerge(Int_t slice,Int_t row_patches)
+{
+ //Constructor
+
+ //fInTracks = (AliL3TrackArray**)new Byte_t[row_patches*sizeof(AliL3TrackArray*)];
+ fInTracks = new AliL3TrackArray*[row_patches];
+ fNPatches = row_patches;
+ for(Int_t i=0; i<row_patches; i++)
+ fInTracks[i] = new AliL3TrackArray("AliL3HoughTrack");
+
+ fOutTracks = new AliL3TrackArray("AliL3HoughTrack");
+}
+
+
+AliL3HoughMerge::~AliL3HoughMerge()
+{
+ //Destructor
+ if(fInTracks)
+ delete fInTracks;
+ if(fOutTracks)
+ delete fOutTracks;
+
+}
+
+void AliL3HoughMerge::FillTracks(AliL3TrackArray *tracks,Int_t patch)
+{
+ fInTracks[patch]->AddTracks(tracks); //copies tracks to new trackarray. Does not delete the track objects.
+}
+
+void AliL3HoughMerge::FillHisto(TH2F *merge_hist)
+{
+
+ for(Int_t pat=0; pat < fNPatches; pat++)
+ {
+ for(Int_t t=0; t<fInTracks[pat]->GetNTracks(); t++)
+ {
+ AliL3HoughTrack *tr = (AliL3HoughTrack*)fInTracks[pat]->GetCheckedTrack(t);
+ if(!tr) {printf("AliL3HoughMerge NO TRACK\n"); continue;}
+ merge_hist->Fill(tr->GetKappa(),tr->GetPhi0(),tr->GetNHits());
+ }
+ }
+}
+
+/*
+void AliL3HoughMerge::MergeLines()
+{
+
+
+
+}
+*/
--- /dev/null
+#ifndef ALIL3_HOUGH_Merge
+#define ALIL3_HOUGH_Merge
+
+#include "AliL3RootTypes.h"
+
+class TClonesArray;
+class TObjArray;
+class TH2F;
+class AliL3TrackArray;
+
+class AliL3HoughMerge : public TObject {
+
+ private:
+
+ AliL3TrackArray **fInTracks; //!
+ AliL3TrackArray *fOutTracks; //!
+
+ Int_t fNPatches;
+ public:
+ AliL3HoughMerge();
+ AliL3HoughMerge(Int_t slice,Int_t row_patches=5);
+ virtual ~AliL3HoughMerge();
+
+ void FillTracks(AliL3TrackArray *tracks,Int_t patch);
+ void FillHisto(TH2F *merge_hist);
+
+ ClassDef(AliL3HoughMerge,1)
+
+};
+
+#endif
--- /dev/null
+
+#include "AliL3Transform.h"
+#include "AliL3HoughPixel.h"
+
+ClassImp(AliL3HoughPixel)
+
+AliL3HoughPixel::AliL3HoughPixel()
+{
+
+}
+
+AliL3HoughPixel::AliL3HoughPixel(Int_t pad,Int_t time,Short_t signal)
+ // Int_t slice,Int_t padrow)
+{
+ //Constructor
+ fPad = pad;
+ fTime = time;
+ fSignal = signal;
+ /*fSlice = slice;
+ fPadrow = padrow;
+ */
+}
+
+AliL3HoughPixel::~AliL3HoughPixel()
+{
+
+}
+/*
+void AliL3HoughPixel::SetCoordinates()
+{
+
+ AliL3Transform *transform = new AliL3Transform();
+
+ Float_t xyz[3];
+ Int_t sector,row;
+ transform->Slice2Sector(fSlice,fPadrow,sector,row);
+ transform->Raw2Global(xyz,sector,row,fPad,fTime);
+ fX = xyz[0];
+ fY = xyz[1];
+ fZ = xyz[2];
+ fPhi = transform->GetPhi(xyz);
+ fEta = transform->GetEta(xyz);
+
+}
+*/
--- /dev/null
+#ifndef ALI_PIXEL
+#define ALI_PIXEL
+
+#include <TObject.h>
+
+class AliL3HoughPixel : public TObject {
+
+ private:
+ Short_t fSignal;
+ Int_t fPad;
+ Int_t fTime;
+ /*Int_t fSlice;
+ Int_t fPadrow;
+
+ Float_t fX;
+ Float_t fY;
+ Float_t fZ;
+ Double_t fPhi;
+ Double_t fEta;
+ */
+ public:
+
+ AliL3HoughPixel();
+ //AliL3HoughPixel(Int_t pad,Int_t time,Short_t signal,Int_t slice,Int_t padrow);
+ AliL3HoughPixel(Int_t pad,Int_t time,Short_t signal);
+ virtual ~AliL3HoughPixel();
+
+ /*Int_t GetPad() {return fPad;}
+ Int_t GetTime() {return fTime;}
+ Short_t GetSignal() {return fSignal;}
+ Double_t GetPhi() {return fPhi;}
+ Double_t GetEta() {return fEta;}
+ Float_t GetX() {return fX;}
+ Float_t GetY() {return fY;}
+ Float_t GetZ() {return fZ;}
+
+ void SetCoordinates();
+ */
+ AliL3HoughPixel *nextRowPixel;
+
+ ClassDef(AliL3HoughPixel,1)
+};
+
+
+#endif
--- /dev/null
+#include <math.h>
+
+#include <TFile.h>
+#include <TTree.h>
+#include <TH2.h>
+
+#include "AliTPCParam.h"
+#include "AliSimDigits.h"
+
+#include "AliL3Defs.h"
+#include "AliL3Transform.h"
+#include "AliL3HoughPixel.h"
+#include "AliL3HoughTransformer.h"
+#include "AliL3HoughTrack.h"
+#include "AliL3TrackArray.h"
+
+ClassImp(AliL3HoughTransformer)
+
+AliL3HoughTransformer::AliL3HoughTransformer()
+{
+ //Default constructor
+ fTransform = 0;
+ fEtaMin = 0;
+ fEtaMax = 0;
+ fSlice = 0;
+ fPatch = 0;
+ fRowContainer = 0;
+ fPhiRowContainer = 0;
+ fNumOfPadRows=0;
+ fContainerBounds=0;
+ fNDigits=0;
+ fIndex = 0;
+}
+
+
+AliL3HoughTransformer::AliL3HoughTransformer(Int_t slice,Int_t patch,Float_t *etarange,Int_t phi_segments)
+{
+ //Constructor
+
+ fTransform = new AliL3Transform();
+
+
+ fEtaMin = etarange[0];
+ fEtaMax = etarange[1];
+ fSlice = slice;
+ fPatch = patch;
+ fNPhiSegments = phi_segments;
+ fNumOfPadRows=NRowsSlice;
+}
+
+
+AliL3HoughTransformer::~AliL3HoughTransformer()
+{
+ //Destructor
+ if(fRowContainer)
+ delete [] fRowContainer;
+ if(fTransform)
+ delete fTransform;
+ if(fPhiRowContainer)
+ delete [] fPhiRowContainer;
+ if(fIndex)
+ {
+ for(Int_t i=0; i<fNDigits; i++)
+ delete [] fIndex[i];
+ delete [] fIndex;
+ }
+}
+
+void AliL3HoughTransformer::InitTemplates(TH2F *hist)
+{
+
+ AliL3Digits *pixel;
+
+ Int_t ymin = hist->GetYaxis()->GetFirst();
+ Int_t ymax = hist->GetYaxis()->GetLast();
+ Int_t nbinsy = hist->GetNbinsY();
+
+ fIndex = new Int_t*[fNDigits];
+ for(Int_t i=0; i<fNDigits; i++)
+ fIndex[i] = new Int_t[nbinsy+1];
+
+ Int_t sector,row;
+ for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
+ {
+
+ for(pixel=(AliL3Digits*)fRowContainer[padrow].first; pixel!=0; pixel=(AliL3Digits*)pixel->nextRowPixel)
+ {
+ Float_t xyz[3];
+ fTransform->Slice2Sector(fSlice,padrow,sector,row);
+ fTransform->Raw2Local(xyz,sector,row,pixel->fPad,pixel->fTime);
+
+ Double_t r_pix = sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
+
+ Double_t phi_pix = fTransform->GetPhi(xyz);
+ Short_t signal = pixel->fCharge;
+ Int_t index = pixel->fIndex;
+ if(index >= fNDigits)
+ printf("AliL3HoughTransformer::InitTemplates : Index error! %d\n",index);
+ for(Int_t p=ymin; p<=ymax; p++)
+ {
+
+ Double_t phi0 = hist->GetYaxis()->GetBinCenter(p);
+ Double_t kappa = 2*sin(phi_pix-phi0)/r_pix;
+
+ Int_t bin = hist->FindBin(kappa,phi0);
+ if(fIndex[index][p]!=0)
+ printf("AliL3HoughTransformer::InitTemplates : Overlapping indexes\n");
+ fIndex[index][p] = bin;
+ }
+ }
+ }
+
+}
+
+
+void AliL3HoughTransformer::CountBins()
+{
+
+ Int_t middle_row = 87; //middle of the slice
+
+ Double_t r_in_bundle = fTransform->Row2X(middle_row);
+ // Double_t phi_min = (fSlice*20 - 10)*ToRad;
+ //Double_t phi_max = (fSlice*20 + 10)*ToRad;
+ Double_t phi_min = -15*ToRad;
+ Double_t phi_max = 15*ToRad;
+
+ Double_t phi_slice = (phi_max - phi_min)/fNPhiSegments;
+ Double_t min_phi0 = 10000;
+ Double_t max_phi0 = 0;
+ Double_t min_kappa = 100000;
+ Double_t max_kappa = 0;
+
+ Int_t xbin = 60;
+ Int_t ybin = 60;
+ Float_t xrange[2] = {-0.006 , 0.006}; //Pt 0.2->
+ Float_t yrange[2] = {-0.26 , 0.26}; //slice 2 0.55->0.88
+
+ TH2F *histo = new TH2F("histo","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+
+ for(Int_t padrow=NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
+ {
+ for(Int_t pad=0; pad < fTransform->GetNPads(padrow); pad++)
+ {
+ for(Int_t time = 0; time < fTransform->GetNTimeBins(); time++)
+ {
+ Float_t xyz[3];
+ Int_t sector,row;
+ fTransform->Slice2Sector(fSlice,padrow,sector,row);
+ fTransform->Raw2Global(xyz,sector,row,pad,time);
+ Double_t eta = fTransform->GetEta(xyz);
+ if(eta < fEtaMin || eta > fEtaMax) continue;
+ fTransform->Global2Local(xyz,sector);
+ Double_t r_pix = sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
+ Double_t phi_pix = fTransform->GetPhi(xyz);
+
+ for(Int_t p=0; p<=fNPhiSegments; p++)
+ {
+ Double_t phi_in_bundle = phi_min + p*phi_slice;
+
+ Double_t tanPhi0 = (r_pix*sin(phi_in_bundle)-r_in_bundle*sin(phi_pix))/(r_pix*cos(phi_in_bundle)-r_in_bundle*cos(phi_pix));
+
+ Double_t phi0 = atan(tanPhi0);
+ // if(phi0 < 0.55 || phi0 > 0.88) continue;
+
+ //if(phi0 < 0) phi0 = phi0 +2*Pi;
+ //Double_t kappa = sin(phi_in_bundle - phi0)*2/r_in_bundle;
+
+ Double_t angle = phi_pix - phi0;
+ Double_t kappa = 2*sin(angle)/r_pix;
+ histo->Fill(kappa,phi0,1);
+ if(phi0 < min_phi0)
+ min_phi0 = phi0;
+ if(phi0 > max_phi0)
+ max_phi0 = phi0;
+ if(kappa < min_kappa)
+ min_kappa = kappa;
+ if(kappa > max_kappa)
+ max_kappa = kappa;
+
+ }
+
+ }
+
+ }
+
+ }
+ Int_t count=0,bi=0;
+
+ Int_t xmin = histo->GetXaxis()->GetFirst();
+ Int_t xmax = histo->GetXaxis()->GetLast();
+ Int_t ymin = histo->GetYaxis()->GetFirst();
+ Int_t ymax = histo->GetYaxis()->GetLast();
+
+
+ for(Int_t xbin=xmin+1; xbin<xmax; xbin++)
+ {
+ for(Int_t ybin=ymin+1; ybin<ymax; ybin++)
+ {
+ bi++;
+ Int_t bin = histo->GetBin(xbin,ybin);
+ if(histo->GetBinContent(bin)>0)
+ count++;
+ }
+ }
+
+
+ printf("Number of possible tracks in this region %d, bins %d\n",count,bi);
+ printf("Phi, min %f max %f\n",min_phi0,max_phi0);
+ printf("Kappa, min %f max %f\n",min_kappa,max_kappa);
+ histo->Draw("box");
+}
+
+
+void AliL3HoughTransformer::Transform2Circle(TH2F *hist,Int_t middle_row)
+{
+ //Transformation is done with respect to local coordinates in slice.
+
+
+ AliL3Digits *pix1;
+ Int_t sector,row;
+
+ //Define a common point
+
+ /*
+ Double_t rowdist1 = fTransform->Row2X(middle_row-1);
+ Double_t rowdist2 = fTransform->Row2X(middle_row);
+ Double_t r_in_bundle = rowdist1 + (rowdist1-rowdist2)/2;
+ */
+
+ Double_t r_in_bundle = fTransform->Row2X(middle_row);
+ //Make overlap between slices
+ Double_t phi_min = -15*ToRad;
+ Double_t phi_max = 15*ToRad;
+
+ Double_t phi_slice = (phi_max - phi_min)/fNPhiSegments;
+
+ for(Int_t p=0; p <= fNPhiSegments; p++)
+ {
+ Double_t phi_in_bundle = phi_min + p*phi_slice;
+ //printf("phi %f in slice %d patch %d middle row %f\n",phi_in_bundle/ToRad,fSlice,fPatch,r_in_bundle);
+
+ for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
+ //for(Int_t padrow = middle_row; padrow <= 173; padrow++)
+ //for(Int_t padrow = 0; padrow <= middle_row; padrow++)
+ {
+
+ for(pix1=(AliL3Digits*)fRowContainer[padrow].first; pix1!=0; pix1=(AliL3Digits*)pix1->nextRowPixel)
+ {
+
+ Float_t xyz[3];
+ fTransform->Slice2Sector(fSlice,padrow,sector,row);
+ fTransform->Raw2Local(xyz,sector,row,pix1->fPad,pix1->fTime);
+ //fTransform->Raw2Global(xyz,sector,row,pix1->fPad,pix1->fTime);
+
+ Double_t r_pix = sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
+
+ Double_t phi_pix = fTransform->GetPhi(xyz);
+
+ //Double_t tanPhi0 = (r_pix*sin(phi_in_bundle)-r_in_bundle*sin(phi_pix))/(r_pix*cos(phi_in_bundle)-r_in_bundle*cos(phi_pix));
+ Double_t tanPhi0 = (r_in_bundle*sin(phi_pix)-r_pix*sin(phi_in_bundle))/(r_in_bundle*cos(phi_pix)-r_pix*cos(phi_in_bundle));
+
+ Double_t phi0 = atan(tanPhi0);
+ //if(padrow > middle_row)
+ //phi0 = -phi0;
+ //if(phi0 < 0.55 || phi0 > 0.88) continue;
+
+ Double_t angle = phi_pix - phi0;
+ Double_t kappa = 2*sin(angle)/r_pix;
+
+ //Double_t angle = phi_in_bundle - phi0;
+ //Double_t kappa = 2*sin(angle)/r_in_bundle;
+
+ //if(kappa < -0.006 || kappa > 0.006) continue;
+
+ Short_t signal = pix1->fCharge;
+
+ hist->Fill(kappa,phi0,signal);
+
+
+ }
+
+ }
+ }
+
+}
+
+
+void AliL3HoughTransformer::Transform2Circle(TH2F *hist)
+{
+ //Transformation is done with respect to local coordinates in slice.
+ //Transform every pixel into whole phirange, using parametrisation:
+ //kappa = 2*sin(phi-phi0)/R
+
+ printf("Transforming 1 pixel only\n");
+
+ AliL3Digits *pix1;
+ Int_t sector,row;
+
+ Int_t ymin = hist->GetYaxis()->GetFirst();
+ Int_t ymax = hist->GetYaxis()->GetLast();
+ Int_t nbinsy = hist->GetNbinsY();
+
+ for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
+ {
+
+ for(pix1=(AliL3Digits*)fRowContainer[padrow].first; pix1!=0; pix1=(AliL3Digits*)pix1->nextRowPixel)
+ {
+
+ Short_t signal = pix1->fCharge;
+ Int_t index = pix1->fIndex;
+
+ for(Int_t p=0; p <= nbinsy; p++)
+ hist->AddBinContent(fIndex[index][p],signal);
+
+ }
+
+ }
+
+
+}
+
+/*
+ void AliL3HoughTransformer::Transform2Circle(TH2F *hist)
+ {
+ //Transformation is done with respect to local coordinates in slice.
+ //Transform every pixel into whole phirange, using parametrisation:
+ //kappa = 2*sin(phi-phi0)/R
+
+ printf("Transforming 1 pixel only\n");
+
+ AliL3Digits *pix1;
+ Int_t sector,row;
+
+ Int_t ymin = hist->GetYaxis()->GetFirst();
+ Int_t ymax = hist->GetYaxis()->GetLast();
+
+ for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
+ {
+
+ for(pix1=(AliL3Digits*)fRowContainer[padrow].first; pix1!=0; pix1=(AliL3Digits*)pix1->nextRowPixel)
+ {
+
+ Float_t xyz[3];
+ fTransform->Slice2Sector(fSlice,padrow,sector,row);
+ fTransform->Raw2Local(xyz,sector,row,pix1->fPad,pix1->fTime);
+
+ Double_t r_pix = sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
+
+ Double_t phi_pix = fTransform->GetPhi(xyz);
+ Short_t signal = pix1->fCharge;
+
+ for(Int_t p=ymin+1; p<=ymax; p++)
+ {
+ Double_t phi0 = hist->GetYaxis()->GetBinCenter(p);
+ Double_t kappa = 2*sin(phi_pix-phi0)/r_pix;
+ hist->Fill(kappa,phi0,signal);
+ }
+
+ }
+
+ }
+
+
+ }
+*/
+
+void AliL3HoughTransformer::TransformLines2Circle(TH2F *hist,AliL3TrackArray *tracks)
+{
+
+ for(Int_t i=0; i<tracks->GetNTracks(); i++)
+ {
+ AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(i);
+ if(!track) {printf("AliL3HoughTransformer::TransformLines2Circle : NO TRACK IN ARRAY\n"); continue;}
+
+ Double_t xmin = fTransform->Row2X(track->GetFirstRow());
+ Double_t xmax = fTransform->Row2X(track->GetLastRow());
+
+ Double_t a = -1./tan(track->GetPsiLine());
+ Double_t b = track->GetDLine()/sin(track->GetPsiLine());
+
+ Double_t ymin = a*xmin + b;
+ Double_t ymax = a*xmax + b;
+
+ Double_t middle_x = xmin + (xmax-xmin)/2;
+ Double_t middle_y = ymin + (ymax-ymin)/2;
+
+ Double_t r_middle = sqrt(middle_x*middle_x + middle_y*middle_y);
+ Double_t phi = atan2(middle_y,middle_x);
+ Double_t phi0 = 2*phi - track->GetPsiLine();
+ Double_t kappa = 2*sin(phi-phi0)/r_middle;
+ hist->Fill(kappa,phi0,track->GetWeight());
+
+ }
+
+
+}
+
+void AliL3HoughTransformer::Transform2Line(TH2F *hist,Int_t ref_row,Int_t *rowrange,Double_t *phirange,TH2F *raw)
+{
+ //Transform every pixel into histogram, using parametrisation:
+ //D = x*cos(psi) + y*sin(psi)
+
+ // printf("In Transform; rowrange %d %d ref_row %d phirange %f %f\n",rowrange[0],rowrange[1],ref_row,phirange[0],phirange[1]);
+
+ AliL3Digits *pix1;
+
+ Int_t xmin = hist->GetXaxis()->GetFirst();
+ Int_t xmax = hist->GetXaxis()->GetLast();
+
+ Double_t x0 = fTransform->Row2X(ref_row);
+ Double_t y0 = 0;
+
+ Int_t sector,row;
+ //for(Int_t padrow = NRows[fPatch][0]; padrow <= NRows[fPatch][1]; padrow++)
+
+ Double_t phi_min = -10*ToRad;
+ Double_t phi_max = 10*ToRad;
+ Double_t delta_phi = (phi_max-phi_min)/fNPhiSegments;
+
+
+ Int_t phi_min_index = (Int_t)((phirange[0]*ToRad-phi_min)/delta_phi);
+ Int_t phi_max_index = (Int_t)((phirange[1]*ToRad-phi_min)/delta_phi);
+
+
+ Int_t index;
+
+ for(Int_t phi=phi_min_index; phi <= phi_max_index; phi++)
+ {
+ for(Int_t padrow = rowrange[0]; padrow <= rowrange[1]; padrow++)
+ {
+
+ index = phi*fNumOfPadRows + padrow;
+ //printf("Looping index %d\n",index);
+ if(index > fContainerBounds || index < 0)
+ {
+ printf("AliL3HoughTransformer::Transform2Line : index %d out of range \n",index);
+ return;
+ }
+ //for(pix1=(AliL3Digits*)fRowContainer[padrow].first; pix1!=0; pix1=(AliL3Digits*)pix1->nextRowPixel)
+ for(pix1=(AliL3Digits*)fPhiRowContainer[index].first; pix1!=0; pix1=(AliL3Digits*)pix1->nextPhiRowPixel)
+ {
+ //printf("Transforming pixel in index %d pad %d time %d padrow %d\n",index,pix1->fPad,pix1->fTime,padrow);
+ Float_t xyz[3];
+ fTransform->Slice2Sector(fSlice,padrow,sector,row);
+ fTransform->Raw2Local(xyz,sector,row,pix1->fPad,pix1->fTime);
+
+ if(raw)
+ raw->Fill(xyz[0],xyz[1],pix1->fCharge);
+
+ xyz[0] = xyz[0]-x0;
+ xyz[1] = xyz[1]-y0;
+
+ //printf("calculating...");
+ for(Int_t d=xmin+1; d<=xmax; d++)
+ {
+ Double_t psi = hist->GetXaxis()->GetBinCenter(d);
+ Double_t D = xyz[0]*cos(psi) + xyz[1]*sin(psi);
+
+ Short_t signal = pix1->fCharge;
+ hist->Fill(psi,D,signal);
+ //printf("Filling histo, psi %f D %f\n",psi,D);
+ }
+ //printf("done\n");
+ }
+ //printf(" done\n");
+ }
+
+ }
+
+}
+
+void AliL3HoughTransformer::GetPixels(Char_t *rootfile,TH2F *hist)
+{
+
+ TFile *file = new TFile(rootfile);
+ file->cd();
+
+ AliTPCParam *param=(AliTPCParam *)file->Get("75x40_100x60");
+ TTree *t=(TTree*)file->Get("TreeD_75x40_100x60");
+
+ AliSimDigits da, *digarr=&da;
+ t->GetBranch("Segment")->SetAddress(&digarr);
+ Stat_t num_of_entries=t->GetEntries();
+
+ Int_t digit_counter=0;
+ Float_t xyz[3];
+ Double_t eta,phi;
+
+ Int_t nrows = NRows[fPatch][1] - NRows[fPatch][0] + 1;
+ printf("nrows %d slice %d patch %d\n",nrows,fSlice,fPatch);
+
+ if(fRowContainer)
+ delete [] fRowContainer;
+ fRowContainer = new AliL3HoughContainer[fNumOfPadRows];
+ memset(fRowContainer,0,fNumOfPadRows*sizeof(AliL3HoughContainer));
+
+
+ fContainerBounds = (fNPhiSegments+1)*(fNumOfPadRows+1);
+ printf("Allocating %d bytes to container of size %d\n",fContainerBounds*sizeof(AliL3HoughContainer),fContainerBounds);
+ if(fPhiRowContainer)
+ delete [] fPhiRowContainer;
+ fPhiRowContainer = new AliL3HoughContainer[fContainerBounds];
+ memset(fPhiRowContainer,0,fContainerBounds*sizeof(AliL3HoughContainer));
+
+ Double_t phi_min = -10*ToRad;
+ Double_t phi_max = 10*ToRad;
+ Double_t delta_phi = (phi_max-phi_min)/fNPhiSegments;
+ Int_t index;
+ digit_counter=0;
+
+ for (Int_t i=0; i<num_of_entries; i++)
+ {
+ t->GetEvent(i);
+ Int_t sector;
+ Int_t row;
+ param->AdjustSectorRow(digarr->GetID(),sector,row);
+ Int_t slice,padrow;
+ fTransform->Sector2Slice(slice,padrow,sector,row);
+ if(slice != fSlice) continue;
+ if(padrow < NRows[fPatch][0]) continue;
+ if(padrow > NRows[fPatch][1]) break;
+ digarr->First();
+ do {
+ Int_t time=digarr->CurrentRow();
+ Int_t pad=digarr->CurrentColumn();
+ Short_t signal=digarr->CurrentDigit();
+ if(time < param->GetMaxTBin()-1 && time > 0)
+ if(digarr->GetDigit(time-1,pad) <= param->GetZeroSup()
+ && digarr->GetDigit(time+1,pad) <= param->GetZeroSup())
+ continue;
+
+
+ fTransform->Raw2Global(xyz,sector,row,pad,time);
+ eta = fTransform->GetEta(xyz);
+ if(eta < fEtaMin || eta > fEtaMax) continue;
+ fTransform->Global2Local(xyz,sector);
+
+
+ phi = fTransform->GetPhi(xyz);
+ if(hist)
+ hist->Fill(xyz[0],xyz[1],signal);
+
+ AliL3Digits *dig = new AliL3Digits;
+ dig->fIndex = digit_counter;
+ digit_counter++;
+ dig->fCharge = signal;
+ dig->fPad = pad;
+ dig->fTime = time;
+
+ if(fRowContainer[padrow].first == NULL)
+ fRowContainer[padrow].first = (void*)dig;
+ else
+ ((AliL3Digits*)(fRowContainer[padrow].last))->nextRowPixel=dig;
+ fRowContainer[padrow].last = (void*)dig;
+
+ Int_t phi_index = (Int_t)((phi-phi_min)/delta_phi);
+ index = phi_index*fNumOfPadRows + padrow;
+ if(phi_index > fContainerBounds || phi_index < 0)
+ {
+ printf("AliL3HoughTransform::GetPixels : index out of range %d\n",phi_index);
+ continue;
+ }
+
+ if(fPhiRowContainer[index].first == NULL)
+ fPhiRowContainer[index].first = (void*)dig;
+ else
+ ((AliL3Digits*)(fPhiRowContainer[index].last))->nextPhiRowPixel = dig;
+ fPhiRowContainer[index].last=(void*)dig;
+
+ }while (digarr->Next());
+
+ }
+
+ fNDigits = digit_counter;
+ printf("digitcounter %d\n",digit_counter);
+ printf("Allocated %d bytes to pixels\n",digit_counter*sizeof(AliL3Digits));
+ file->Close();
+ delete file;
+
+}
--- /dev/null
+#ifndef ALIL3_HOUGHTRANSFORMER
+#define ALIL3_HOUGHTRANSFORMER
+
+#include "AliL3RootTypes.h"
+
+struct AliL3Digits
+{
+ UShort_t fCharge;
+ UChar_t fPad;
+ UShort_t fTime;
+ Int_t fIndex;
+ AliL3Digits *nextPhiRowPixel;
+ AliL3Digits *nextRowPixel;
+};
+
+struct AliL3HoughContainer
+{
+ void *first;
+ void *last;
+};
+
+class TH2F;
+class AliL3Transform;
+class AliL3HoughEvaluate;
+class AliL3TrackArray;
+
+class AliL3HoughTransformer : public TObject {
+
+ private:
+
+ friend class AliL3HoughEval;
+
+ AliL3Transform *fTransform; //!
+
+ Int_t fNPhiSegments; //Number of patches in phi.
+ Float_t fEtaMin;
+ Float_t fEtaMax;
+ Int_t fNumEtaSegments;
+ Int_t fNumOfPadRows;
+
+ AliL3HoughContainer *fRowContainer; //!
+ AliL3HoughContainer *fPhiRowContainer; //!
+ Int_t fContainerBounds;
+ Int_t fNDigits;
+ Int_t **fIndex; //!
+
+
+ Int_t fSlice;
+ Int_t fPatch;
+
+ public:
+ AliL3HoughTransformer();
+ AliL3HoughTransformer(Int_t slice,Int_t patch,Float_t *etarange,Int_t phi_segments);
+ virtual ~AliL3HoughTransformer();
+
+ void Transform2Circle(TH2F *hist,Int_t middle_row);
+ void Transform2Circle(TH2F *hist);
+ void Transform2Line(TH2F *hist,Int_t ref_row,Int_t *rowrange,Double_t *phirange,TH2F *raw=0);
+ void TransformLines2Circle(TH2F *hist,AliL3TrackArray *tracks);
+ void GetPixels(Char_t *rootfile,TH2F *hist=0);
+ void InitTemplates(TH2F *hist);
+ void CountBins();
+
+ ClassDef(AliL3HoughTransformer,1)
+
+};
+
+#endif
--- /dev/null
+
+//Author: Anders Strand Vestbo
+//Author: Uli Frankenfeld
+//Last Modified: 05.12.2000
+
+#include "AliL3Logging.h"
+#include "AliL3Transform.h"
+//#include <TFile.h>
+#include <math.h>
+//___________________________
+// AliL3Transform
+//
+// Transformation class for ALICE TPC.
+
+ClassImp(AliL3Transform);
+
+
+AliL3Transform::AliL3Transform(){
+ //constructor
+ Init();
+}
+
+
+AliL3Transform::~AliL3Transform(){
+}
+
+void AliL3Transform::Init(){
+ //sector:
+ fNTimeBins = 512; //?uli
+ fNRowLow = 55;
+ fNRowUp = 119;
+ fNSectorLow = 36;
+ fNSectorUp = 36;
+ fNSector = 72;
+ fPadPitchWidthLow = 0.400000;
+ fPadPitchWidthUp = 0.600000;
+ fZWidth = 0.56599998474121093750;
+ fZSigma = 0.22880849748219134199;
+
+ //slices:
+ fNSlice = 36;
+ fNRow = 174;
+ fPi = 3.14159265358979323846;
+ for(Int_t i=0;i<36;i++){
+ fCos[i] = cos(fPi/9*i);
+ fSin[i] = sin(fPi/9*i);
+ }
+
+ fX[0] = 88.850006103515625;
+ fX[1] = 89.600006103515625;
+ fX[2] = 90.350006103515625;
+ fX[3] = 91.100006103515625;
+ fX[4] = 91.850006103515625;
+ fX[5] = 92.600006103515625;
+ fX[6] = 93.350006103515625;
+ fX[7] = 94.100006103515625;
+ fX[8] = 94.850006103515625;
+ fX[9] = 95.600006103515625;
+ fX[10] = 96.350006103515625;
+ fX[11] = 97.100006103515625;
+ fX[12] = 97.850006103515625;
+ fX[13] = 98.600006103515625;
+ fX[14] = 99.350006103515625;
+ fX[15] = 100.100006103515625;
+ fX[16] = 100.850006103515625;
+ fX[17] = 101.600006103515625;
+ fX[18] = 102.350006103515625;
+ fX[19] = 103.100006103515625;
+ fX[20] = 103.850006103515625;
+ fX[21] = 104.600006103515625;
+ fX[22] = 105.350006103515625;
+ fX[23] = 106.100006103515625;
+ fX[24] = 106.850006103515625;
+ fX[25] = 107.600006103515625;
+ fX[26] = 108.350006103515625;
+ fX[27] = 109.100006103515625;
+ fX[28] = 109.850006103515625;
+ fX[29] = 110.600006103515625;
+ fX[30] = 111.350006103515625;
+ fX[31] = 112.100006103515625;
+ fX[32] = 112.850006103515625;
+ fX[33] = 113.600006103515625;
+ fX[34] = 114.350006103515625;
+ fX[35] = 115.100006103515625;
+ fX[36] = 115.850006103515625;
+ fX[37] = 116.600006103515625;
+ fX[38] = 117.350006103515625;
+ fX[39] = 118.100006103515625;
+ fX[40] = 118.850006103515625;
+ fX[41] = 119.600006103515625;
+ fX[42] = 120.350006103515625;
+ fX[43] = 121.100006103515625;
+ fX[44] = 121.850006103515625;
+ fX[45] = 122.600006103515625;
+ fX[46] = 123.350006103515625;
+ fX[47] = 124.100006103515625;
+ fX[48] = 124.850006103515625;
+ fX[49] = 125.600006103515625;
+ fX[50] = 126.350006103515625;
+ fX[51] = 127.100006103515625;
+ fX[52] = 127.850006103515625;
+ fX[53] = 128.600006103515625;
+ fX[54] = 129.350006103515625;
+ fX[55] = 132.574996948242188;
+ fX[56] = 133.574996948242188;
+ fX[57] = 134.574996948242188;
+ fX[58] = 135.574996948242188;
+ fX[59] = 136.574996948242188;
+ fX[60] = 137.574996948242188;
+ fX[61] = 138.574996948242188;
+ fX[62] = 139.574996948242188;
+ fX[63] = 140.574996948242188;
+ fX[64] = 141.574996948242188;
+ fX[65] = 142.574996948242188;
+ fX[66] = 143.574996948242188;
+ fX[67] = 144.574996948242188;
+ fX[68] = 145.574996948242188;
+ fX[69] = 146.574996948242188;
+ fX[70] = 147.574996948242188;
+ fX[71] = 148.574996948242188;
+ fX[72] = 149.574996948242188;
+ fX[73] = 150.574996948242188;
+ fX[74] = 151.574996948242188;
+ fX[75] = 152.574996948242188;
+ fX[76] = 153.574996948242188;
+ fX[77] = 154.574996948242188;
+ fX[78] = 155.574996948242188;
+ fX[79] = 156.574996948242188;
+ fX[80] = 157.574996948242188;
+ fX[81] = 158.574996948242188;
+ fX[82] = 159.574996948242188;
+ fX[83] = 160.574996948242188;
+ fX[84] = 161.574996948242188;
+ fX[85] = 162.574996948242188;
+ fX[86] = 163.574996948242188;
+ fX[87] = 164.574996948242188;
+ fX[88] = 165.574996948242188;
+ fX[89] = 166.574996948242188;
+ fX[90] = 167.574996948242188;
+ fX[91] = 168.574996948242188;
+ fX[92] = 169.574996948242188;
+ fX[93] = 170.574996948242188;
+ fX[94] = 171.574996948242188;
+ fX[95] = 172.574996948242188;
+ fX[96] = 173.574996948242188;
+ fX[97] = 174.574996948242188;
+ fX[98] = 175.574996948242188;
+ fX[99] = 176.574996948242188;
+ fX[100] = 177.574996948242188;
+ fX[101] = 178.574996948242188;
+ fX[102] = 179.574996948242188;
+ fX[103] = 180.574996948242188;
+ fX[104] = 181.574996948242188;
+ fX[105] = 182.574996948242188;
+ fX[106] = 183.574996948242188;
+ fX[107] = 184.574996948242188;
+ fX[108] = 185.574996948242188;
+ fX[109] = 186.574996948242188;
+ fX[110] = 187.574996948242188;
+ fX[111] = 188.574996948242188;
+ fX[112] = 189.574996948242188;
+ fX[113] = 190.574996948242188;
+ fX[114] = 191.574996948242188;
+ fX[115] = 192.574996948242188;
+ fX[116] = 193.574996948242188;
+ fX[117] = 194.574996948242188;
+ fX[118] = 195.574996948242188;
+ fX[119] = 196.574996948242188;
+ fX[120] = 197.574996948242188;
+ fX[121] = 198.574996948242188;
+ fX[122] = 199.574996948242188;
+ fX[123] = 200.574996948242188;
+ fX[124] = 201.574996948242188;
+ fX[125] = 202.574996948242188;
+ fX[126] = 203.574996948242188;
+ fX[127] = 204.574996948242188;
+ fX[128] = 205.574996948242188;
+ fX[129] = 206.574996948242188;
+ fX[130] = 207.574996948242188;
+ fX[131] = 208.574996948242188;
+ fX[132] = 209.574996948242188;
+ fX[133] = 210.574996948242188;
+ fX[134] = 211.574996948242188;
+ fX[135] = 212.574996948242188;
+ fX[136] = 213.574996948242188;
+ fX[137] = 214.574996948242188;
+ fX[138] = 215.574996948242188;
+ fX[139] = 216.574996948242188;
+ fX[140] = 217.574996948242188;
+ fX[141] = 218.574996948242188;
+ fX[142] = 219.574996948242188;
+ fX[143] = 220.574996948242188;
+ fX[144] = 221.574996948242188;
+ fX[145] = 222.574996948242188;
+ fX[146] = 223.574996948242188;
+ fX[147] = 224.574996948242188;
+ fX[148] = 225.574996948242188;
+ fX[149] = 226.574996948242188;
+ fX[150] = 227.574996948242188;
+ fX[151] = 228.574996948242188;
+ fX[152] = 229.574996948242188;
+ fX[153] = 230.574996948242188;
+ fX[154] = 231.574996948242188;
+ fX[155] = 232.574996948242188;
+ fX[156] = 233.574996948242188;
+ fX[157] = 234.574996948242188;
+ fX[158] = 235.574996948242188;
+ fX[159] = 236.574996948242188;
+ fX[160] = 237.574996948242188;
+ fX[161] = 238.574996948242188;
+ fX[162] = 239.574996948242188;
+ fX[163] = 240.574996948242188;
+ fX[164] = 241.574996948242188;
+ fX[165] = 242.574996948242188;
+ fX[166] = 243.574996948242188;
+ fX[167] = 244.574996948242188;
+ fX[168] = 245.574996948242188;
+ fX[169] = 246.574996948242188;
+ fX[170] = 247.574996948242188;
+ fX[171] = 248.574996948242188;
+ fX[172] = 249.574996948242188;
+ fX[173] = 250.574996948242188;
+ fNPads[0] = 71;
+ fNPads[1] = 71;
+ fNPads[2] = 71;
+ fNPads[3] = 73;
+ fNPads[4] = 73;
+ fNPads[5] = 73;
+ fNPads[6] = 75;
+ fNPads[7] = 75;
+ fNPads[8] = 75;
+ fNPads[9] = 77;
+ fNPads[10] = 77;
+ fNPads[11] = 77;
+ fNPads[12] = 79;
+ fNPads[13] = 79;
+ fNPads[14] = 79;
+ fNPads[15] = 81;
+ fNPads[16] = 81;
+ fNPads[17] = 81;
+ fNPads[18] = 83;
+ fNPads[19] = 83;
+ fNPads[20] = 83;
+ fNPads[21] = 85;
+ fNPads[22] = 85;
+ fNPads[23] = 85;
+ fNPads[24] = 87;
+ fNPads[25] = 87;
+ fNPads[26] = 87;
+ fNPads[27] = 89;
+ fNPads[28] = 89;
+ fNPads[29] = 89;
+ fNPads[30] = 89;
+ fNPads[31] = 91;
+ fNPads[32] = 91;
+ fNPads[33] = 91;
+ fNPads[34] = 93;
+ fNPads[35] = 93;
+ fNPads[36] = 93;
+ fNPads[37] = 95;
+ fNPads[38] = 95;
+ fNPads[39] = 95;
+ fNPads[40] = 97;
+ fNPads[41] = 97;
+ fNPads[42] = 97;
+ fNPads[43] = 99;
+ fNPads[44] = 99;
+ fNPads[45] = 99;
+ fNPads[46] = 101;
+ fNPads[47] = 101;
+ fNPads[48] = 101;
+ fNPads[49] = 103;
+ fNPads[50] = 103;
+ fNPads[51] = 103;
+ fNPads[52] = 105;
+ fNPads[53] = 105;
+ fNPads[54] = 105;
+ fNPads[55] = 73;
+ fNPads[56] = 73;
+ fNPads[57] = 73;
+ fNPads[58] = 75;
+ fNPads[59] = 75;
+ fNPads[60] = 75;
+ fNPads[61] = 75;
+ fNPads[62] = 77;
+ fNPads[63] = 77;
+ fNPads[64] = 77;
+ fNPads[65] = 79;
+ fNPads[66] = 79;
+ fNPads[67] = 79;
+ fNPads[68] = 81;
+ fNPads[69] = 81;
+ fNPads[70] = 81;
+ fNPads[71] = 81;
+ fNPads[72] = 83;
+ fNPads[73] = 83;
+ fNPads[74] = 83;
+ fNPads[75] = 85;
+ fNPads[76] = 85;
+ fNPads[77] = 85;
+ fNPads[78] = 85;
+ fNPads[79] = 87;
+ fNPads[80] = 87;
+ fNPads[81] = 87;
+ fNPads[82] = 89;
+ fNPads[83] = 89;
+ fNPads[84] = 89;
+ fNPads[85] = 91;
+ fNPads[86] = 91;
+ fNPads[87] = 91;
+ fNPads[88] = 91;
+ fNPads[89] = 93;
+ fNPads[90] = 93;
+ fNPads[91] = 93;
+ fNPads[92] = 95;
+ fNPads[93] = 95;
+ fNPads[94] = 95;
+ fNPads[95] = 95;
+ fNPads[96] = 97;
+ fNPads[97] = 97;
+ fNPads[98] = 97;
+ fNPads[99] = 99;
+ fNPads[100] = 99;
+ fNPads[101] = 99;
+ fNPads[102] = 101;
+ fNPads[103] = 101;
+ fNPads[104] = 101;
+ fNPads[105] = 101;
+ fNPads[106] = 103;
+ fNPads[107] = 103;
+ fNPads[108] = 103;
+ fNPads[109] = 105;
+ fNPads[110] = 105;
+ fNPads[111] = 105;
+ fNPads[112] = 105;
+ fNPads[113] = 107;
+ fNPads[114] = 107;
+ fNPads[115] = 107;
+ fNPads[116] = 109;
+ fNPads[117] = 109;
+ fNPads[118] = 109;
+ fNPads[119] = 111;
+ fNPads[120] = 111;
+ fNPads[121] = 111;
+ fNPads[122] = 111;
+ fNPads[123] = 113;
+ fNPads[124] = 113;
+ fNPads[125] = 113;
+ fNPads[126] = 115;
+ fNPads[127] = 115;
+ fNPads[128] = 115;
+ fNPads[129] = 115;
+ fNPads[130] = 117;
+ fNPads[131] = 117;
+ fNPads[132] = 117;
+ fNPads[133] = 119;
+ fNPads[134] = 119;
+ fNPads[135] = 119;
+ fNPads[136] = 121;
+ fNPads[137] = 121;
+ fNPads[138] = 121;
+ fNPads[139] = 121;
+ fNPads[140] = 123;
+ fNPads[141] = 123;
+ fNPads[142] = 123;
+ fNPads[143] = 125;
+ fNPads[144] = 125;
+ fNPads[145] = 125;
+ fNPads[146] = 125;
+ fNPads[147] = 127;
+ fNPads[148] = 127;
+ fNPads[149] = 127;
+ fNPads[150] = 129;
+ fNPads[151] = 129;
+ fNPads[152] = 129;
+ fNPads[153] = 129;
+ fNPads[154] = 131;
+ fNPads[155] = 131;
+ fNPads[156] = 131;
+ fNPads[157] = 133;
+ fNPads[158] = 133;
+ fNPads[159] = 133;
+ fNPads[160] = 135;
+ fNPads[161] = 135;
+ fNPads[162] = 135;
+ fNPads[163] = 135;
+ fNPads[164] = 137;
+ fNPads[165] = 137;
+ fNPads[166] = 137;
+ fNPads[167] = 139;
+ fNPads[168] = 139;
+ fNPads[169] = 139;
+ fNPads[170] = 139;
+ fNPads[171] = 141;
+ fNPads[172] = 141;
+ fNPads[173] = 141;
+}
+
+
+Double_t AliL3Transform::GetEta(Float_t *xyz)
+{
+ Double_t r3 = sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]+xyz[2]*xyz[2]);
+ Double_t eta = 0.5 * log((r3+xyz[2])/(r3-xyz[2]));
+ return eta;
+}
+
+Double_t AliL3Transform::GetPhi(Float_t *xyz)
+{
+
+ Double_t phi = atan2(xyz[1],xyz[0]);
+ //if(phi<0) phi=phi+2*TMath::Pi();
+ return phi;
+}
+
+
+Bool_t AliL3Transform::Slice2Sector(Int_t slice, Int_t slicerow, Int_t & sector, Int_t &row) const{
+ if(slicerow<0&&slicerow>=fNRow) return kFALSE;
+ if(slice<0||slice>=fNSlice) return kFALSE;
+
+ if(slicerow<fNRowLow){
+ sector = slice;
+ row = slicerow;
+ }
+ else {
+ sector = slice+fNSlice;
+ row = slicerow-fNRowLow;
+ }
+ return kTRUE;
+}
+
+Bool_t AliL3Transform::Sector2Slice(Int_t & slice, Int_t sector) const{
+ if(sector<0||sector>=fNSector) return kFALSE;
+ if(sector<fNSectorLow) slice = sector;
+ else slice = sector - fNSectorLow;
+ return kTRUE;
+}
+
+Bool_t AliL3Transform::Sector2Slice(Int_t & slice, Int_t & slicerow,Int_t sector, Int_t row) const{
+ if(sector<0||sector>=fNSector||row<0) return kFALSE;
+ if(sector<fNSectorLow){
+ if(row>=fNRowLow) return kFALSE;
+ slice = sector;
+ slicerow = row;
+ }
+ else{
+ if(row>=fNRowUp) return kFALSE;
+ slice = sector - fNSectorLow;
+ slicerow = row + fNRowLow;
+ }
+ return kTRUE;
+}
+
+Double_t AliL3Transform::Row2X(Int_t slicerow){
+ if(slicerow<0||slicerow>=fNRow) return 0;
+ return fX[slicerow];
+}
+
+void AliL3Transform::Local2Global(Float_t *xyz,Int_t slice)
+{
+ //Transformation to global coordinate system
+ Float_t x0 = xyz[0];
+ Float_t y0 = xyz[1];
+ Float_t cs,sn;
+ cs = fCos[slice];
+ sn = fSin[slice];
+ xyz[0]=x0*cs-y0*sn;
+ xyz[1]=x0*sn+y0*cs;
+ xyz[2]=xyz[2];//global z=local z
+}
+
+void AliL3Transform::Local2GlobalAngle(Float_t *angle,Int_t slice){
+ angle[0] = fmod(angle[0]+slice*fPi/9,2*fPi);
+}
+
+void AliL3Transform::Global2LocalAngle(Float_t *angle,Int_t slice){
+ angle[0] = angle[0]-slice*fPi/9;
+ if(angle[0]<0) angle[0]+=2*fPi;
+}
+
+void AliL3Transform::Raw2Local(Float_t *xyz,Int_t sector,Int_t row,Float_t pad,Float_t time)
+{
+ //Transformation from rawdata to local coordinate system
+
+ Int_t slice,slicerow;
+ Sector2Slice(slice, slicerow, sector, row);
+
+ xyz[0]=Row2X(slicerow);
+ Int_t npads= fNPads[slicerow];
+ if(sector<fNSectorLow)
+ xyz[1]=(pad-0.5*(npads-1))*fPadPitchWidthLow;
+ else
+ xyz[1]=(pad-0.5*(npads-1))*fPadPitchWidthUp;
+ xyz[2]=fZWidth*time-3.*fZSigma;
+ Int_t sign=-1;
+ Int_t nis=fNSectorLow;
+ Int_t nos=fNSectorUp;
+
+ if((sector<nis)/2 || ((sector-nis)<nos/2)) sign=1;
+ xyz[2]=sign*(250.-xyz[2]);
+
+}
+
+void AliL3Transform::Local2Global(Float_t *xyz,Int_t sector,Int_t row)
+{
+ //Transformation to global coordinate system
+ Int_t slice,slicerow;
+ Sector2Slice(slice, slicerow, sector, row);
+ Float_t r=Row2X(slicerow);
+ Float_t cs = fCos[slice];
+ Float_t sn = fSin[slice];
+
+ xyz[0]=r*cs-xyz[1]*sn;
+ xyz[1]=r*sn+xyz[1]*cs;
+ xyz[2]=xyz[2];//global z=local z
+}
+
+Double_t AliL3Transform::GetMaxY(Int_t slicerow)
+{
+
+ if(slicerow < fNRowLow)
+ return fPadPitchWidthLow*fNPads[slicerow]/2;
+
+ else
+ return fPadPitchWidthUp*fNPads[slicerow]/2;
+}
+
+void AliL3Transform::Global2Local(Float_t *xyz,Int_t sector)
+{
+ Int_t slice;
+ Sector2Slice(slice, sector);
+ Float_t cs = fCos[slice];
+ Float_t sn = fSin[slice];
+ Float_t x1 = xyz[0]*cs + xyz[1]*sn;
+ Float_t y1 = -xyz[0]*sn + xyz[1]*cs;
+ xyz[0] = x1;
+ xyz[1] = y1;
+}
+
+void AliL3Transform::Raw2Global(Float_t *xyz,Int_t sector,Int_t row,Float_t pad,Float_t time)
+{
+ //Transformation from raw to global coordinates
+
+ Raw2Local(xyz,sector,row,pad,time);
+ Local2Global(xyz,sector,row);
+}
+
+void AliL3Transform::Local2Raw(Float_t *xyz,Int_t sector,Int_t row)
+{
+ //Transformation from local coordinates to raw
+
+ Int_t slice,slicerow;
+ Sector2Slice(slice, slicerow, sector, row);
+
+ if(sector<fNSectorLow)
+ xyz[1]=xyz[1]/fPadPitchWidthLow+0.5*(fNPads[slicerow]-1);
+ else
+ xyz[1]=xyz[1]/fPadPitchWidthUp+0.5*(fNPads[slicerow]-1);
+ Int_t sign=-1;
+ Int_t nis=fNSectorLow;
+ Int_t nos=fNSectorUp;
+
+ if ((sector<nis/2) || ((sector-nis)<nos/2)) sign=1;
+ xyz[2]=250-sign*xyz[2];
+ xyz[2]=(xyz[2]+3.*fZSigma)/fZWidth;
+}
+
+void AliL3Transform::Global2Raw(Float_t *xyz,Int_t sector,Int_t row)
+{
+ //Transformation from global coordinates to raw.
+
+ Global2Local(xyz,sector);
+ Local2Raw(xyz,sector,row);
+
+}
--- /dev/null
+#ifndef ALIL3TRANSFORM_H
+#define ALIL3TRANSFORM_H
+
+
+
+//#include "AliTPCParam.h"
+#include "AliL3RootTypes.h"
+class TFile;
+
+class AliL3Transform {
+ private:
+// AliTPCParam *fParam;
+ Int_t fNTimeBins;
+ Int_t fNRowLow;
+ Int_t fNRowUp;
+ Int_t fNSectorLow;
+ Int_t fNSectorUp;
+ Double_t fPadPitchWidthLow;
+ Double_t fPadPitchWidthUp;
+ Double_t fZWidth;
+ Double_t fZSigma;
+ Int_t fNSector;
+ Int_t fNSlice;
+ Int_t fNRow;
+ Double_t fPi;
+ Double_t fCos[36];
+ Double_t fSin[36];
+ Double_t fX[174];
+ Int_t fNPads[174];
+ public:
+ AliL3Transform();
+ virtual ~AliL3Transform();
+ void Init();
+
+ Double_t GetPadPitchWidthLow() {return fPadPitchWidthLow;}
+ Double_t GetPadPitchWidthUp() {return fPadPitchWidthUp;}
+
+ Bool_t Slice2Sector(Int_t slice, Int_t slicerow, Int_t & sector, Int_t &row) const;
+
+ Bool_t Sector2Slice(Int_t & slice, Int_t sector) const;
+ Bool_t Sector2Slice(Int_t & slice, Int_t & slicerow,Int_t sector, Int_t row) const;
+
+ Double_t Row2X(Int_t slicerow);
+ Int_t GetNPads(Int_t row){return (row<174)?fNPads[row]:0;}
+ Int_t GetNTimeBins(){return fNTimeBins;}
+
+ Double_t GetEta(Float_t *xyz);
+ Double_t GetPhi(Float_t *xyz);
+ Double_t GetMaxY(Int_t slicerow);
+ void Local2Global(Float_t *xyz,Int_t slice);
+ void Local2GlobalAngle(Float_t *angle,Int_t slice);
+ void Global2LocalAngle(Float_t *angle,Int_t slice);
+
+ void Raw2Local(Float_t *xyz,Int_t sector,Int_t row,Float_t pad,Float_t time);
+ void Local2Global(Float_t *xyz,Int_t sector,Int_t row);
+ void Global2Local(Float_t *xyz,Int_t sector);
+ void Raw2Global(Float_t *xyz,Int_t sector,Int_t row,Float_t pad,Float_t time);
+ void Local2Raw(Float_t *xyz,Int_t sector,Int_t row);
+ void Global2Raw(Float_t *xyz,Int_t sector,Int_t row);
+
+ ClassDef(AliL3Transform,1)
+};
+
+
+#endif
--- /dev/null
+############################### TPC Makefile ##################################
+
+# Include machine specific definitions
+
+include $(ALICE_ROOT)/conf/GeneralDef
+include $(ALICE_ROOT)/conf/MachineDef.$(ALICE_TARGET)
+
+PACKAGE = AliL3Hough
+
+# C++ sources
+
+
+SRCS = AliL3HoughTransformer.cxx AliL3HoughPixel.cxx \
+ AliL3HoughMaxFinder.cxx AliL3HoughEval.cxx AliL3HoughMerge.cxx
+
+
+# C++ Headers
+
+HDRS = $(SRCS:.cxx=.h) AliL3Defs.h AliL3HoughLinkDef.h
+# Library dictionary
+
+DICT = AliL3HoughCint.cxx
+DICTH = $(DICT:.cxx=.h)
+DICTO = $(patsubst %.cxx,tgt_$(ALICE_TARGET)/%.o,$(DICT))
+
+# FORTRAN Objectrs
+
+FOBJS = $(FSRCS:.f=.o)
+
+# C Objects
+
+COBJS = $(CSRCS:.c=.o)
+
+# C++ Objects
+
+OBJS = $(patsubst %.cxx,tgt_$(ALICE_TARGET)/%.o,$(SRCS)) $(DICTO)
+
+# C++ compilation flags
+
+CXXFLAGS = $(CXXOPTS) -g -I$(ROOTSYS)/include -I. -I$(ALICE_ROOT)/include/ -I$(ALICE_ROOT)/TPC \
+ -I$(ALICE_ROOT)/CONTAINERS -I$(ALICE)/level3code/src
+
+#CXXFLAGS = $(CXXOPTS) -g -Wall -I$(ROOTSYS)/include -I. -I $(ALICE_ROOT)/TPC -I$(ALICE_ROOT)/include/ -DCOMPILING
+
+# FORTRAN compilation flags
+
+FFLAGS = $(FOPT)
+
+##### TARGETS #####
+
+# Target
+
+SLIBRARY = $(ALICE)/mylibs/libAliL3Hough.$(SL)
+
+default: $(SLIBRARY)
+
+$(ALICE)/mylibs/libAliL3Hough.$(SL): $(OBJS)
+
+$(DICT): $(HDRS)
+
+depend: $(SRCS)
+
+TOCLEAN = $(OBJS) *Cint.h *Cint.cxx
+
+############################### General Macros ################################
+
+include $(ALICE_ROOT)/conf/GeneralMacros
+
+############################ Dependencies #####################################
+
+include tgt_$(ALICE_TARGET)/Make-depend
+
+###########
+so:
+ rm -fr $(ALICE)/mylibs/libAliL3Hough.so
+clean:
+ rm -fr tgt_Linux/*.o
+ rm -fr $(ALICE)/mylibs/libAliL3Hough.so
+ rm -fr $(DICT) $(DICTH) $(DICTO)
--- /dev/null
+void hough(char *rootfile,int patch,Bool_t MC=false)
+{
+ gStyle->SetOptStat(0);
+
+ Int_t xbin = 60;
+ Int_t ybin = 60;
+ Float_t xrange[2] = {-0.006 , 0.006}; //Pt 0.2->
+ //Float_t yrange[2] = {-0.17 , 0.17}; //slice 2 0.55->0.88
+ Float_t yrange[2] = {-0.26 , 0.26};// -15 - 15
+ //Float_t yrange[2] = {0.55 , 0.88}; //slice 2 0.55->0.88
+ TH2F *hist = new TH2F("hist","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ hist->GetXaxis()->SetTitle("#kappa [cm^{-1}]");
+ hist->GetYaxis()->SetTitle("#Phi [rad]");
+ SetTH1Options(hist);
+
+ Int_t xr[2] = {0,250};
+ Int_t yr[2] = {-125,125};
+ TH1F *ntracks = new TH1F("ntracks","",100,0,200);
+ TH2F *raw = new TH2F("raw","",250,xr[0],xr[1],250,yr[0],yr[1]);
+ TH2F *road = new TH2F("road","",250,0,250,250,yr[0],yr[1]);
+ TH2F *fake = new TH2F("fake","",250,0,250,250,-125,125);
+ TH2F *peaks = new TH2F("peaks","",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ peaks->SetMarkerStyle(3);
+ peaks->SetMarkerColor(2);
+
+ int slice = 2,n_phi_segments=30;
+ //float eta[2] = {0.3,0.4};
+ float eta[2] = {0.03,0.04};
+
+ AliL3HoughTransformer *a = new AliL3HoughTransformer(slice,patch,eta,n_phi_segments);
+ a->GetPixels(rootfile,raw);
+ a->InitTemplates(hist);
+
+ //a->Transform2Circle(hist,78);
+ TStopwatch sw;
+ a->Transform2Circle(hist);
+ sw.Stop();
+ printf("Transformation done in %f ms\n",sw.CpuTime()*1000);
+
+ AliL3HoughMaxFinder *b = new AliL3HoughMaxFinder("KappaPhi");
+ //b->SetThreshold(10000);
+ AliL3TrackArray *tracks = b->FindMaxima(hist);
+
+ AliL3HoughEval *c = new AliL3HoughEval(a);
+ printf("number of peaks before %d\n",tracks->GetNTracks());
+
+ /*
+ for(int i=0; i<tracks->GetNTracks(); i++)
+ {
+ AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(i);
+ if(!track) printf("no track\n");
+ printf("Pt %f Phi %f kappa %f weigth %d charge %d\n",track->GetPt(),track->GetPhi0(),track->GetKappa(),track->GetNHits(),track->GetCharge());
+ }
+ */
+ //c->LookInsideRoad(a,tracks,road,fake,ntracks);
+ c->LookInsideRoad(tracks);
+ printf("Number of trackcandidates %d\n",tracks->GetNTracks());
+
+ AliL3Transform *transform = new AliL3Transform();
+ Int_t NRows[5][2] = {{ 0, 45},{46,77},{78,109},{110,141},{142,173}};
+ for(Int_t i=0; i<tracks->GetNTracks(); i++)
+ {
+ AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(i);
+ if(!track) {printf("NO TRACK\n"); continue;}
+ peaks->Fill(track->GetKappa(),track->GetPhi0(),1);
+ //for(int row=0; row<174; row++)
+ for(int row=NRows[patch][0]; row<=NRows[patch][1]; row++)
+ {
+ Float_t xyz[3];
+ track->GetCrossingPoint(row,xyz);
+ //track->GetCrossingPoint(slice,row,xyz);
+ //printf("Track does not cross line at row %d\n",row);
+ //transform->Local2Global(xyz,slice);
+ road->Fill(xyz[0],xyz[1],1);
+ }
+ //float an[2] = {track->GetPhi0(),0};
+ //transform->Local2GlobalAngle(an,slice);
+ printf("pt %f phi0 %f weight %d\n",track->GetPt(),track->GetPhi0(),track->GetWeight());
+ }
+
+
+ TCanvas *c1 = new TCanvas("c1","",900,900);
+ SetCanvasOptions(c1);
+ c1->Divide(2,2);
+ c1->cd(1);
+ hist->Draw("lego1");
+ // peaks->Draw("same");
+
+ // TCanvas *c2 = new TCanvas("c2","",2);
+ c1->cd(3);
+ raw->Draw("");
+
+ //TCanvas *c3 = new TCanvas("c3","",2);
+ //road->SetMarkerStyle(5);
+ //road->SetMarkerColor(3);
+ c1->cd(4);
+ road->SetMarkerStyle(9);
+ road->Draw();
+ /*
+ TCanvas *c4 = new TCanvas("c4","",2);
+ fake->Draw("cont1");
+
+ TCanvas *c5 = new TCanvas("c5","",2);
+ ntracks->Draw();
+ */
+ delete tracks;
+ delete a;
+ delete b;
+ delete c;
+
+ //----------------------------------------------------------------------------
+ if(MC)
+ {
+ TFile *file = new TFile(rootfile);
+ file->cd();
+
+ gAlice = (AliRun*)file->Get("gAlice");
+ gAlice->GetEvent(0);
+
+ TClonesArray *particles=gAlice->Particles();
+ Int_t n=particles->GetEntriesFast();
+ Float_t torad=TMath::Pi()/180;
+ Float_t phi_min = slice*20 - 10;
+ Float_t phi_max = slice*20 + 10;
+
+ for (Int_t j=0; j<n; j++) {
+ TParticle *p=(TParticle*)particles->UncheckedAt(j);
+ if (p->GetFirstMother()>=0) continue; //secondary particle
+ if(p->Eta() < eta[0] || p->Eta() > eta[1]) continue;
+ Double_t ptg=p->Pt(),pxg=p->Px(),pyg=p->Py(),pzg=p->Pz();
+ Double_t phi_part = TMath::ATan2(pyg,pxg);
+ if (phi_part < 0) phi_part += 2*TMath::Pi();
+
+ if(phi_part < phi_min*torad || phi_part > phi_max*torad) {continue;}
+ if(ptg<0.100) continue;
+
+ printf("Generated particle %d, with pt %f phi %f\n",p->GetPdgCode(),ptg,phi_part);
+ }
+ file->Close();
+ delete file;
+ }
+
+}
--- /dev/null
+void hough_line(char *rootfile,int patch,Bool_t MC=false)
+{
+ gStyle->SetOptStat(0);
+
+
+ Double_t pi = TMath::Pi();
+ Int_t xbin = 60;
+ Int_t ybin = 60;
+ Float_t xrange[2] = {-pi/2,pi/2};
+ Float_t yrange[2] = {-40,40};
+
+ TH2F *hist = new TH2F("hist","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ hist->GetXaxis()->SetTitle("#Psi");
+ hist->GetYaxis()->SetTitle("D");
+ SetTH1Options(hist);
+
+ Int_t xr[2] = {0,250};
+ Int_t yr[2] = {-125,125};
+ TH1F *ntracks = new TH1F("ntracks","",100,0,200);
+ TH2F *raw = new TH2F("raw","",250,xr[0],xr[1],250,yr[0],yr[1]);
+ TH2F *road = new TH2F("road","",250,0,250,250,yr[0],yr[1]);
+ TH2F *container = new TH2F("container","",250,0,250,250,yr[0],yr[1]);
+ TH2F *peaks = new TH2F("peaks","",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ peaks->SetMarkerStyle(3);
+ peaks->SetMarkerColor(2);
+
+ int slice = 2,n_phi_segments=20;
+ float eta[2] = {0.3,0.4};
+ //float eta[2] = {0.05,0.06};
+
+ AliL3HoughTransformer *a = new AliL3HoughTransformer(slice,patch,eta,n_phi_segments);
+
+ const Int_t NRows[5][2] = {{ 0, 45},{46,77},{78,109},{110,141},{142,173}};
+ const int ref_row[5]={22,60,93,125,157};
+ int rowrange[2] = {6,12};
+ double phirange[2] = {-10.,10.};
+ a->GetPixels(rootfile,raw);
+ a->Transform2Line(hist,9,rowrange,phirange);
+
+ AliL3HoughMaxFinder *b = new AliL3HoughMaxFinder("DPsi");
+ AliL3TrackArray *tracks = b->FindMaxima(hist,rowrange,5);
+
+ TH2F *cross = new TH2F("cross","",250,xr[0],xr[1],250,yr[0],yr[1]);
+ for(int i=0; i<2; i++)
+ //for(int i=0; i<tracks->GetNTracks(); i++)
+ {
+ AliL3HoughTrack *track = (AliL3HoughTrack*)tracks->GetCheckedTrack(i);
+ if(!track) continue;
+ Double_t xy[2];
+ for(int j=0; j<174; j++)
+ {
+ track->GetLineCrossingPoint(j,xy);
+ cross->Fill(xy[0],xy[1],1);
+ }
+ printf("Segment %f %f weight %d\n",track->GetPsiLine(),track->GetDLine(),track->GetWeight());
+ }
+ cross->SetMarkerStyle(6);
+ cross->SetMarkerColor(2);
+
+
+ TCanvas *c1 = new TCanvas("c1","",900,900);
+ SetCanvasOptions(c1);
+ c1->Divide(2,2);
+ c1->cd(1);
+ hist->Draw("box");
+ peaks->SetMarkerStyle(3);
+ peaks->SetMarkerColor(2);
+ //peaks->Draw("same");
+
+ c1->cd(2);
+ container->Draw();
+
+ c1->cd(3);
+ raw->Draw("");
+ cross->Draw("same");
+ return;
+ //TCanvas *c3 = new TCanvas("c3","",2);
+ //road->SetMarkerStyle(5);
+ //road->SetMarkerColor(3);
+ c1->cd(4);
+ road->SetMarkerStyle(9);
+ road->Draw();
+ /*
+ TCanvas *c4 = new TCanvas("c4","",2);
+ fake->Draw("cont1");
+
+ TCanvas *c5 = new TCanvas("c5","",2);
+ ntracks->Draw();
+ */
+ delete a;
+ delete b;
+ delete c;
+
+ //----------------------------------------------------------------------------
+ if(MC)
+ {
+ TFile *file = new TFile(rootfile);
+ file->cd();
+
+ gAlice = (AliRun*)file->Get("gAlice");
+ gAlice->GetEvent(0);
+
+ TClonesArray *particles=gAlice->Particles();
+ Int_t n=particles->GetEntriesFast();
+ Float_t torad=TMath::Pi()/180;
+ Float_t phi_min = slice*20 - 10;
+ Float_t phi_max = slice*20 + 10;
+
+ for (Int_t j=0; j<n; j++) {
+ TParticle *p=(TParticle*)particles->UncheckedAt(j);
+ if (p->GetFirstMother()>=0) continue; //secondary particle
+ if(p->Eta() < eta[0] || p->Eta() > eta[1]) continue;
+ Double_t ptg=p->Pt(),pxg=p->Px(),pyg=p->Py(),pzg=p->Pz();
+ Double_t phi_part = TMath::ATan2(pyg,pxg);
+ if (phi_part < 0) phi_part += 2*TMath::Pi();
+
+ if(phi_part < phi_min*torad || phi_part > phi_max*torad) {continue;}
+ if(ptg<0.100) continue;
+
+ printf("Generated particle %d, with pt %f phi %f\n",p->GetPdgCode(),ptg,phi_part);
+ }
+ file->Close();
+ delete file;
+ }
+
+}
--- /dev/null
+void hough_line_merge(char *rootfile,int patch,Bool_t MC=false)
+{
+ gStyle->SetOptStat(0);
+
+
+ Double_t pi = TMath::Pi();
+ Int_t xbin = 60;
+ Int_t ybin = 60;
+ Float_t xrange[2] = {-pi/2,pi/2};
+ Float_t yrange[2] = {-40,40};
+
+ TH2F *hist;
+
+ Int_t xr[2] = {0,250};
+ Int_t yr[2] = {-125,125};
+ TH1F *ntracks = new TH1F("ntracks","",100,0,200);
+ TH2F *raw = new TH2F("raw","",250,xr[0],xr[1],250,yr[0],yr[1]);
+ TH2F *road = new TH2F("road","",250,0,250,250,yr[0],yr[1]);
+ TH2F *fake = new TH2F("fake","",250,0,250,250,-125,125);
+ TH2F *peaks = new TH2F("peaks","",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ peaks->SetMarkerStyle(3);
+ peaks->SetMarkerColor(2);
+
+ int slice = 2,n_phi_segments=20;
+ float eta[2] = {0.3,0.4};
+ //float eta[2] = {0.05,0.06};
+
+ AliL3HoughTransformer *a = new AliL3HoughTransformer(slice,patch,eta,n_phi_segments);
+ a->GetPixels(rootfile,raw);
+
+ AliL3HoughMaxFinder *b = new AliL3HoughMaxFinder("DPsi");
+
+ const Int_t NRows[5][2] = {{ 0, 45},{46,77},{78,109},{110,141},{142,173}};
+ //const int ref_row[4]={5,15,25,35};
+ //int rowrange[4][2]={{0,10},{10,20},{20,30},{30,45}};
+ int ref_row[7] = {3,9,15,21,27,33,39};
+ int rowrange[7][2] = {{0,6},{6,12},{12,18},{18,24},{24,30},{30,36},{36,42}};
+
+ double phirange[7][2] = {{-10.,-5.},{-7.5,-2.5},{-5.,0.},{-2.5,2.5},{0.,5.},{2.5,7.5},{5.,10.}};
+ int npatches = 7;
+ int count=0;
+
+ //AliL3TrackArray **track_pieces = new AliL3TrackArray*[npatches];
+ AliL3TrackArray *track_pieces = new AliL3TrackArray("AliL3HoughTrack");
+
+ double pran[2] = {-10,10};
+
+ // int pa=3;
+ for(int row_pat=0; row_pat < 7; row_pat++)
+ {
+
+ for(int pat=0; pat<7; pat++)
+ {
+
+ TH2F *histo = new TH2F("histo","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+
+ printf("transforming ref_row %d rowrange %d %d phirange %f %f\n",ref_row[row_pat],rowrange[row_pat][0],rowrange[row_pat][1],phirange[pat][0],phirange[pat][1]);
+
+ a->Transform2Line(histo,ref_row[row_pat],rowrange[row_pat],phirange[pat]);
+ //a->Transform2Line(histo,ref_row[row_pat],rowrange[row_pat],pran);
+
+ AliL3TrackArray *tmp = (AliL3TrackArray*)b->FindMaxima(histo,rowrange[row_pat],ref_row[row_pat]);
+
+ track_pieces->AddTracks(tmp);
+
+ //AliL3HoughTrack *tra = (AliL3HoughTrack*)b->GetTracks()->GetCheckedTrack(0);
+ //printf("psi %f rowrange %d %d\n",tra->GetPsiLine(),tra->GetFirstRow(),tra->GetLastRow());
+
+ delete histo;
+ delete tmp;
+ }
+
+ }
+
+
+ printf("Transformation complete, genereted %d tracks\n",track_pieces->GetNTracks());
+
+
+ //track_pieces->Sort();
+ TH2F *cross = new TH2F("cross","",250,xr[0],xr[1],250,yr[0],yr[1]);
+ Double_t xy[2];
+ //for(int i=0; i<npatches; i++)
+ for(int i=0; i<4; i++)
+ {
+ if(track_pieces->GetNTracks()==0) continue;
+ AliL3HoughTrack *track = (AliL3HoughTrack*)track_pieces->GetCheckedTrack(i);
+ if(!track) {printf("NO TRACK\n"); continue;}
+ //for(int j=rowrange[pa][0]; j<=rowrange[pa][1]; j++)
+ /*
+ for(int j=0; j<174; j++)
+ {
+ track->GetLineCrossingPoint(j,xy);
+ cross->Fill(xy[0],xy[1],1);
+ }
+ */
+ //printf("weight %d psi %f D %f rowrange %d %d\n",track->GetWeight(),track->GetPsiLine(),track->GetDLine(),track->GetFirstRow(),track->GetLastRow());
+
+ }
+ cross->SetMarkerStyle(6);
+ cross->SetMarkerColor(2);
+ //cross->Draw();
+
+
+
+ AliL3HoughTransformer *a2 = new AliL3HoughTransformer(slice,pat,eta,n_phi_segments);
+
+ xbin = 60 ;
+ ybin = 60;
+ xrange[0] = -0.006;
+ xrange[1] = 0.006;
+ yrange[0] = -0.26;
+ yrange[1] = 0.26;
+ hist = new TH2F("hist","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+
+ TH2F *cross_helix = new TH2F("cross_helix","",250,xr[0],xr[1],250,yr[0],yr[1]);
+
+ a2->TransformLines2Circle(hist,(AliL3TrackArray*)track_pieces);
+
+ AliL3HoughMaxFinder *b2 = new AliL3HoughMaxFinder("KappaPhi");
+ AliL3TrackArray *final_tracks = b2->FindMaxima(hist);
+
+ printf("Number of track candidates before checking %d\n",final_tracks->GetNTracks());
+ AliL3HoughEval *eval = new AliL3HoughEval(slice);
+ eval->SetTransformer(a);
+ eval->LookInsideRoad(final_tracks);
+
+ printf("Total number of tracks after checking %d\n",final_tracks->GetNTracks());
+ for(int i=0; i<3; i++)//final_tracks->GetNTracks(); i++)
+ {
+
+ AliL3HoughTrack *track = (AliL3HoughTrack*)final_tracks->GetCheckedTrack(i);
+ if(!track) continue;
+
+ for(int r=0; r<174; r++)
+ {
+ Float_t xyz_helix[3];
+ track->GetCrossingPoint(r,xyz_helix);
+ cross_helix->Fill(xyz_helix[0],xyz_helix[1],1);
+ }
+
+ peaks->Fill(track->GetKappa(),track->GetPhi0(),1);
+ printf("Found track, pt %f phi0 %f\n",track->GetPt(),track->GetPhi0());
+ }
+ cross_helix->SetMarkerStyle(6);
+ cross_helix->SetMarkerColor(2);
+ TCanvas *c1 = new TCanvas("c1","",900,900);
+ SetCanvasOptions(c1);
+ c1->Divide(2,2);
+ c1->cd(1);
+ hist->Draw("box");
+
+ peaks->SetMarkerStyle(3);
+ peaks->SetMarkerColor(2);
+ //peaks->Draw("same");
+
+ c1->cd(2);
+ cross_helix->Draw();
+
+ c1->cd(3);
+ raw->Draw("");
+ cross_helix->DrawCopy("same");
+ //cross->Draw("same");
+
+ c1->cd(4);
+ //hist->Draw("box");
+ raw->DrawCopy();
+
+ delete track_pieces;
+ delete a;
+ delete b;
+ delete a2;
+ delete eval;
+ //----------------------------------------------------------------------------
+ if(MC)
+ {
+ TFile *file = new TFile(rootfile);
+ file->cd();
+
+ gAlice = (AliRun*)file->Get("gAlice");
+ gAlice->GetEvent(0);
+
+ TClonesArray *particles=gAlice->Particles();
+ Int_t n=particles->GetEntriesFast();
+ Float_t torad=TMath::Pi()/180;
+ Float_t phi_min = slice*20 - 10;
+ Float_t phi_max = slice*20 + 10;
+
+ for (Int_t j=0; j<n; j++) {
+ TParticle *p=(TParticle*)particles->UncheckedAt(j);
+ if (p->GetFirstMother()>=0) continue; //secondary particle
+ if(p->Eta() < eta[0] || p->Eta() > eta[1]) continue;
+ Double_t ptg=p->Pt(),pxg=p->Px(),pyg=p->Py(),pzg=p->Pz();
+ Double_t phi_part = TMath::ATan2(pyg,pxg);
+ if (phi_part < 0) phi_part += 2*TMath::Pi();
+
+ if(phi_part < phi_min*torad || phi_part > phi_max*torad) {continue;}
+ if(ptg<0.100) continue;
+
+ printf("Generated particle %d, with pt %f phi %f\n",p->GetPdgCode(),ptg,phi_part);
+ }
+ file->Close();
+ delete file;
+ }
+
+}
--- /dev/null
+void hough_merge(char *rootfile,Bool_t MC=false)
+{
+ gStyle->SetOptStat(0);
+
+ Int_t xbin = 60;
+ Int_t ybin = 60;
+ Float_t xrange[2] = {-0.006 , 0.006}; //Pt 0.1->4GeV 0.006-0.006
+ //Float_t yrange[2] = {-0.17 , 0.17}; //slice 2 0.55->0.88
+ Float_t yrange[2] = {-0.26 , 0.26}; //slice 2 0.55->0.88
+ //Float_t yrange[2] = {0.55 , 0.88}; //slice 2 0.55->0.88
+
+ Int_t xr[2] = {0,250};
+ Int_t yr[2] = {-125,125};
+ TH1F *ntracks = new TH1F("ntracks","",100,0,200);
+ TH2F *raw = new TH2F("raw","",250,xr[0],xr[1],250,yr[0],yr[1]);
+ TH2F *road = new TH2F("road","",250,0,250,250,yr[0],yr[1]);
+ TH2F *fake = new TH2F("fake","",300,0,300,250,-125,125);
+ TH2F *peaks = new TH2F("peaks","",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ peaks->SetMarkerStyle(3);
+ peaks->SetMarkerColor(2);
+
+ int slice = 2,patch=0,n_phi_segments=40;
+ float eta[2] = {0.3,0.4};
+ // float eta[2] = {0.04,0.05};
+
+ AliL3HoughTransformer *a;
+ AliL3HoughMaxFinder *b;
+ AliL3HoughEval *c = new AliL3HoughEval(slice);
+ AliL3HoughMerge *d = new AliL3HoughMerge(slice);
+ TH2F *hist;
+ for(int pat=0; pat<5; pat++)
+ {
+ a = new AliL3HoughTransformer(slice,pat,eta,n_phi_segments);
+ hist = new TH2F("hist","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ a->GetPixels(rootfile,raw);
+ a->Transform2Circle(hist,87);
+
+ b = new AliL3HoughMaxFinder("KappaPsi");
+ AliL3TrackArray *tracks = b->FindMaxima(hist);
+ c->SetTransformer(a);
+ c->LookInsideRoad(tracks);
+ d->FillTracks(tracks,pat);
+ delete a;
+ delete b;
+ }
+ return;
+
+
+ TH2F *merge_hist = new TH2F("merge_hist","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+
+ d->FillHisto(merge_hist);
+
+
+ b = new AliL3HoughMaxFinder(slice,pat);
+ b->FindMaxima(merge_hist);
+
+ AliL3TrackArray *mtrs = b->GetTracks();
+
+ c->CompareMC(rootfile,mtrs,eta);
+
+ AliL3Transform *transform = new AliL3Transform();
+ for(int tr=0; tr<mtrs->GetNTracks(); tr++)
+ {
+ AliL3HoughTrack *t = (AliL3HoughTrack*)mtrs->GetCheckedTrack(tr);
+ if(!t) {printf("NO TRACK11\n"); continue;}
+ if(t->GetMCid() < 0) {printf("Track %d was fake, weigth %d\n",tr,t->GetNHits()); continue;}
+ printf("Found pt %f weigth %d\n",t->GetPt(),t->GetNHits());
+ //printf("charge %f\n",t->GetCharge());
+ for(Int_t j=0; j<174; j++)
+ {
+ Float_t xyz[3];
+ if(!t->GetCrossingPoint(j,xyz))
+ {
+ printf("Track does not cross line\n");
+ continue;
+ }
+ road->Fill(xyz[0],xyz[1],1);
+
+ }
+
+
+ }
+
+
+ TCanvas *c1 = new TCanvas("c1","",800,800);
+ SetCanvasOptions(c1);
+ c1->Divide(2,2);
+ c1->cd(1);
+ merge_hist->Draw("box");
+
+ //TCanvas *c2 = new TCanvas("c2","",2);
+ c1->cd(3);
+ raw->Draw();
+
+ //TCanvas *c3 = new TCanvas("c3","",2);
+ c1->cd(4);
+ road->SetMarkerStyle(7);
+ road->Draw();
+
+ delete b;
+ delete d;
+ delete c;
+ //----------------------------------------------------------------------------
+ if(MC)
+ {
+ TFile *file = new TFile(rootfile);
+ file->cd();
+
+ gAlice = (AliRun*)file->Get("gAlice");
+ gAlice->GetEvent(0);
+
+ TClonesArray *particles=gAlice->Particles();
+ Int_t n=particles->GetEntriesFast();
+ Float_t torad=TMath::Pi()/180;
+ Float_t phi_min = slice*20 - 10;
+ Float_t phi_max = slice*20 + 10;
+
+ for (Int_t j=0; j<n; j++) {
+ TParticle *p=(TParticle*)particles->UncheckedAt(j);
+ if (p->GetFirstMother()>=0) continue; //secondary particle
+ if(p->Eta() < eta[0] || p->Eta() > eta[1]) continue;
+ Double_t ptg=p->Pt(),pxg=p->Px(),pyg=p->Py(),pzg=p->Pz();
+ Double_t phi_part = TMath::ATan2(pyg,pxg);
+ if (phi_part < 0) phi_part += 2*TMath::Pi();
+
+ if(phi_part < phi_min*torad || phi_part > phi_max*torad) {continue;}
+ if(ptg<0.100) continue;
+ int found = 0;
+ for(int m=0; m<mtrs->GetNTracks(); m++)
+ {
+ AliL3HoughTrack *tra = (AliL3HoughTrack*)mtrs->GetCheckedTrack(m);
+ if(tra->GetMCid() != j) continue;
+ found = 1;
+ double dPt = fabs(ptg-tra->GetPt());
+ float phi0[2] = {tra->GetPhi0(),0};
+ transform->Local2GlobalAngle(phi0,slice);
+ double dPhi0 = fabs(phi_part-phi0[0]);
+ printf("Difference Pt %f Phi0 %f\n",dPt,dPhi0);
+
+ }
+ if(!found) printf("Track %d was not found\n",j);
+ printf("Generated particle %d, with pt %f phi %f\n",p->GetPdgCode(),ptg,phi_part);
+ }
+ file->Close();
+ delete file;
+ }
+
+}
--- /dev/null
+void hough_mergehistos(char *rootfile)
+{
+ gStyle->SetOptStat(0);
+
+ Int_t xbin = 60;
+ Int_t ybin = 60;
+ Float_t xrange[2] = {-0.006 , 0.006}; //Pt 0.1->4GeV 0.006-0.006
+ //Float_t yrange[2] = {-0.17 , 0.17}; //slice 2 0.55->0.88
+ Float_t yrange[2] = {-0.26 , 0.26}; //slice 2 0.55->0.88
+ //Float_t yrange[2] = {0.55 , 0.88}; //slice 2 0.55->0.88
+
+ Int_t xr[2] = {0,250};
+ Int_t yr[2] = {-125,125};
+ TH1F *ntracks = new TH1F("ntracks","",100,0,200);
+ TH2F *raw = new TH2F("raw","",250,xr[0],xr[1],250,yr[0],yr[1]);
+ TH2F *road = new TH2F("road","",250,0,250,250,yr[0],yr[1]);
+ TH2F *fake = new TH2F("fake","",300,0,300,250,-125,125);
+ TH2F *peaks = new TH2F("peaks","",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ peaks->SetMarkerStyle(3);
+ peaks->SetMarkerColor(2);
+
+ int slice = 2,patch=0,n_phi_segments=40;
+ float eta[2] = {0.3,0.4};
+ // float eta[2] = {0.04,0.05};
+
+ AliL3HoughTransformer *a;
+ AliL3HoughMaxFinder *b = new AliL3HoughMaxFinder("KappaPsi");
+
+ TH2F **hist = new TH2F*[5];
+
+ for(int pat=0; pat<5; pat++)
+ {
+ a = new AliL3HoughTransformer(slice,pat,eta,n_phi_segments);
+ hist[pat] = new TH2F("hist","Parameter space",xbin,xrange[0],xrange[1],ybin,yrange[0],yrange[1]);
+ a->GetPixels(rootfile,raw);
+ a->InitTemplates(hist[pat]);
+ a->Transform2Circle(hist[pat]);
+
+ //AliL3TrackArray *tracks = b->FindMaxima(hist[pat]);
+
+ delete a;
+
+ }
+
+
+ TCanvas *c1 = new TCanvas("c1","",800,800);
+ c1->Divide(2,2);
+ c1->cd(1);
+ hist[0]->DrawCopy("box");
+ c1->cd(2);
+ hist[4]->DrawCopy("box");
+ c1->cd(3);
+ hist[2]->DrawCopy("box");
+
+ c1->cd(4);
+ hist[4]->Add(hist[0],hist[2]);
+ hist[4]->Draw("box");
+
+ delete b;
+
+
+ delete [] hist;
+}
--- /dev/null
+
+{
+ printf("\nWELCOME to the magic world of Level3\n\n");
+
+ //load library neccesary for TPC alone
+ gSystem->Load("$(ROOTSYS)/lib/libPhysics");
+ gSystem->Load("$(ROOTSYS)/lib/libEG");
+ gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libSTEER");
+ gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libCONTAINERS");
+ gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libTPC");
+ gSystem->Load("/usr/local/franken/lib/MLUC/lib/linux-i386/libMLUC.so");
+
+ gSystem->Load("$(ALICE)/mylibs/libAliL3");
+ gSystem->Load("$(ALICE)/mylibs/libAliL3Hough");
+
+ gROOT->LoadMacro("$(ALICE)/dev/XFunct.C");
+ gStyle->SetStatBorderSize(1);
+ gStyle->SetTitleBorderSize(0);
+ printf("TPC libraries loaded\n");
+ printf("Level3 libraries loaded\n");
+}