1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // high level filter algorithm for TPC using a hough transformation //
22 ///////////////////////////////////////////////////////////////////////////////
24 #include <TStopwatch.h>
26 #include "AliHLTStandardIncludes.h"
27 #include "AliHLTLogging.h"
28 #include "AliHLTTransform.h"
29 #include "AliHLTHough.h"
31 #include <AliHLTITSclusterer.h>
32 #include <AliHLTITSVertexerZ.h>
33 #include <AliHLTITStracker.h>
35 #include "AliHoughFilter.h"
37 #include "AliRawReaderRoot.h"
39 #include <AliMagFMaps.h>
40 #include <AliKalmanTrack.h>
41 #include <AliITSgeom.h>
42 #include <AliESDVertex.h>
44 ClassImp(AliHoughFilter)
46 //_____________________________________________________________________________
47 AliHoughFilter::AliHoughFilter():
51 // default constructor
54 AliHLTLog::fgLevel = AliHLTLog::kError;
55 if (AliDebugLevel() > 0) AliHLTLog::fgLevel = AliHLTLog::kWarning;
56 if (AliDebugLevel() > 1) AliHLTLog::fgLevel = AliHLTLog::kInformational;
57 if (AliDebugLevel() > 2) AliHLTLog::fgLevel = AliHLTLog::kDebug;
59 // Init TPC HLT geometry
60 const char *path = gSystem->Getenv("ALICE_ROOT");
61 Char_t pathname[1024];
62 strcpy(pathname,path);
63 strcat(pathname,"/HLT/src");
64 if (!AliHLTTransform::Init(pathname, kFALSE))
65 AliError("HLT initialization failed!");
67 // Init magnetic field
68 AliMagF* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k5kG);
69 AliTracker::SetFieldMap(field,kTRUE);
70 fPtmin = 0.1*AliHLTTransform::GetSolenoidField();
73 fITSgeom = new AliITSgeom();
74 fITSgeom->ReadNewFile("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymmFMD.det");
76 AliError("ITS geometry not found!");
82 //_____________________________________________________________________________
83 Bool_t AliHoughFilter::Filter(AliRawEvent* event, AliESDEvent* esd)
85 // Run fast online reconstruction
86 // based on the HLT tracking algorithms
88 TStopwatch globaltimer;
94 TTree *treeITSclusters = new TTree("TreeL3ITSclusters"," "); //make a tree
95 treeITSclusters->SetDirectory(0);
97 RunITSclusterer(event,treeITSclusters);
98 RunITSvertexer(esd,treeITSclusters);
99 RunTPCtracking(event,esd);
100 RunITStracking(esd,treeITSclusters);
101 delete treeITSclusters;
103 AliInfo(Form("Event filter has finished in %f seconds\n\n\n\n\n\n",globaltimer.RealTime()));
108 //_____________________________________________________________________________
109 void AliHoughFilter::RunITSclusterer(AliRawEvent* event, TTree *treeClusters)
112 // The clusters are stored in a tree
117 AliHLTITSclusterer clusterer(0);
118 AliRawReader *itsrawreader=new AliRawReaderRoot(event);
119 clusterer.Digits2Clusters(itsrawreader,treeClusters);
121 AliInfo(Form("ITS clusterer has finished in %f seconds\n",timer.RealTime()));
126 //_____________________________________________________________________________
127 void AliHoughFilter::RunITSvertexer(AliESDEvent* esd, TTree *treeClusters)
130 // Store the result in the ESD
135 AliHLTITSVertexerZ vertexer;
136 AliESDVertex *vertex = vertexer.FindVertexForCurrentEvent(fITSgeom,treeClusters);
137 esd->SetVertex(vertex);
138 AliInfo(Form("ITS vertexer has finished in %f seconds\n",timer.RealTime()));
142 //_____________________________________________________________________________
143 void AliHoughFilter::RunTPCtracking(AliRawEvent* event, AliESDEvent* esd)
145 // Run hough transform tracking in TPC
146 // The z of the vertex is taken from the ESD
147 // The result of the tracking is stored in the ESD
151 const AliESDVertex *vertex = esd->GetVertex();
152 Float_t zvertex = vertex->GetZv();
154 AliHLTHough *hough1 = new AliHLTHough();
156 hough1->SetThreshold(4);
157 hough1->CalcTransformerParams(fPtmin);
158 hough1->SetPeakThreshold(70,-1);
159 hough1->Init(100,4,event,zvertex);
160 hough1->SetAddHistograms();
162 AliHLTHough *hough2 = new AliHLTHough();
164 hough2->SetThreshold(4);
165 hough2->CalcTransformerParams(fPtmin);
166 hough2->SetPeakThreshold(70,-1);
167 hough2->Init(100,4,event,zvertex);
168 hough2->SetAddHistograms();
170 Int_t nglobaltracks = 0;
171 /* In case we run HLT code in 2 threads */
172 hough1->StartProcessInThread(0,17);
173 hough2->StartProcessInThread(18,35);
175 if(hough1->WaitForThreadFinish())
176 ::Fatal("AliHLTHough::WaitForThreadFinish"," Can not join the required thread! ");
177 if(hough2->WaitForThreadFinish())
178 ::Fatal("AliHLTHough::WaitForThreadFinish"," Can not join the required thread! ");
180 /* In case we run HLT code in the main thread
181 for(Int_t slice=0; slice<=17; slice++)
183 hough1->ReadData(slice,0);
185 hough1->AddAllHistogramsRows();
186 hough1->FindTrackCandidatesRow();
189 for(Int_t slice=18; slice<=35; slice++)
191 hough2->ReadData(slice,0);
193 hough2->AddAllHistogramsRows();
194 hough2->FindTrackCandidatesRow();
199 nglobaltracks += hough1->FillESD(esd);
200 nglobaltracks += hough2->FillESD(esd);
202 /* In case we want to debug the ESD
203 gSystem->MakeDirectory("hough1");
204 hough1->WriteTracks("./hough1");
205 gSystem->MakeDirectory("hough2");
206 hough2->WriteTracks("./hough2");
212 AliInfo(Form("\nHough Transformer has found %d TPC tracks in %f seconds\n\n", nglobaltracks,timer.RealTime()));
216 //_____________________________________________________________________________
217 void AliHoughFilter::RunITStracking(AliESDEvent* esd, TTree *treeClusters)
219 // Run the ITS tracker
220 // The tracks from the HT TPC tracking are used as seeds
221 // The prologated tracks are updated in the ESD
226 Double_t vtxErr[3]={0.005,0.005,0.010};
227 const AliESDVertex *vertex = esd->GetVertex();
228 vertex->GetXYZ(vtxPos);
230 AliHLTITStracker itsTracker(0);
231 itsTracker.SetVertex(vtxPos,vtxErr);
233 itsTracker.LoadClusters(treeClusters);
234 itsTracker.Clusters2Tracks(esd);
235 itsTracker.UnloadClusters();
236 AliInfo(Form("ITS tracker has finished in %f seconds\n",timer.RealTime()));