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 // this program performs local monitoring on a GDC by running the HLT code //
22 // If an argument is given, this is taken as the name of a date file which //
23 // is used instead of the local node. //
24 // The program can be stopped by pressing CTRL-C. //
26 ///////////////////////////////////////////////////////////////////////////////
29 #include <TSysEvtHandler.h>
34 #include "AliRawReaderDate.h"
37 #include <AliHLTStandardIncludes.h>
38 #include <AliHLTTransform.h>
39 #include <AliHLTMemHandler.h>
40 #include <AliHLTTrackArray.h>
41 #include <AliHLTHoughMaxFinder.h>
42 #include <AliHLTHoughBaseTransformer.h>
43 #include <AliHLTHough.h>
44 #include <AliHLTBenchmark.h>
45 #include <AliKalmanTrack.h>
46 #include "AliITSgeomTGeo.h"
47 #include "AliITSgeom.h"
49 #include "AliMagFMaps.h"
50 #include <AliHLTITSclusterer.h>
51 #include <AliHLTITSVertexerZ.h>
52 #include <AliHLTITStracker.h>
55 //_____________________________________________________________________________
56 class AliGDCInterruptHandler : public TSignalHandler {
58 AliGDCInterruptHandler();
59 Bool_t Notify() {fStop = kTRUE; return kTRUE;};
60 Bool_t Stop() const {return fStop;};
62 Bool_t fStop; // CTRL-C pressed
65 //_____________________________________________________________________________
66 AliGDCInterruptHandler::AliGDCInterruptHandler() :
67 TSignalHandler(kSigInterrupt, kFALSE)
73 //_____________________________________________________________________________
75 int main(int argc, char** argv)
77 // set ROOT in batch mode
81 FILE* file = fopen("monitorGDC.log", "w");
85 // get data from a file or online from this node
88 status = monitorSetDataSource(argv[1]);
90 status = monitorSetDataSource(":");
92 if (status) ::Fatal("monitorSetDataSource", monitorDecodeError(status));
94 // monitor only a sample of physics events
95 char* table[] = {"Physics event", "yes", NULL};
96 status = monitorDeclareTable(table);
97 if (status) ::Fatal("monitorDeclareTable", monitorDecodeError(status));
99 // declare this monitoring program to DATE
100 status = monitorDeclareMp("GDC physics monitoring");
101 if (status) ::Fatal("monitorDeclareMp", monitorDecodeError(status));
103 // initialize HLT transformations
104 if (!AliHLTTransform::Init("./", kFALSE)) {
105 ::Fatal("AliHLTTransform::Init", "HLT initialization failed");
107 AliESDEvent *esd = new AliESDEvent;
108 esd->CreateStdContent();
109 // AliKalmanTrack::SetConvConst(
110 // 1000/0.299792458/AliHLTTransform::GetSolenoidField()
112 AliITSgeom *geom = new AliITSgeom();
113 geom->ReadNewFile("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymmFMD.det");
116 switch ((Int_t)(AliHLTTransform::GetSolenoidField()+0.5)) {
118 sfield = AliMagFMaps::k2kG;
121 sfield = AliMagFMaps::k4kG;
124 sfield = AliMagFMaps::k5kG;
127 ::Fatal("AliHLTTransform::GetSolenoidField", "Incorrect magnetic field");
129 AliMagF* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., sfield);
130 AliTracker::SetFieldMap(field,kTRUE);
135 // create the signal handler
136 AliGDCInterruptHandler* handler = new AliGDCInterruptHandler;
137 gSystem->AddSignalHandler(handler);
141 while (!handler->Stop()) {
142 // get the next event
143 status = monitorGetEventDynamic(&ptr);
144 if (status == (Int_t)MON_ERR_EOF) break;
145 if (status) ::Fatal("monitorGetEventDynamic", monitorDecodeError(status));
149 gSystem->Sleep(1000); // sleep for 1 second
153 AliRawReaderDate rawReader(ptr);
154 // if ((rawReader.GetAttributes()[0] & 0x02) != 0) {
156 // Int_t errorCode = rawReader.CheckData();
158 if (errorCode && (errorCode != AliRawReader::kErrSize)) {
160 if (file) fprintf(file, "%s\n", time.AsString());
161 if (file) fprintf(file, "run: %d event: %d %d\n",
162 rawReader.GetRunNumber(),
163 rawReader.GetEventId()[0],
164 rawReader.GetEventId()[1]);
165 fprintf(file, "ERROR: %d\n\n", errorCode);
169 AliHLTBenchmark *fBenchmark = new AliHLTBenchmark();
170 fBenchmark->Start("Overall timing");
172 // ITS clusterer and vertexer
173 fBenchmark->Start("ITS Clusterer");
174 AliHLTITSclusterer clusterer(0);
175 AliRawReader *itsrawreader=new AliRawReaderDate(ptr);
176 TTree* treeClusters = new TTree("TreeL3ITSclusters"," "); //make a tree
177 clusterer.Digits2Clusters(itsrawreader,treeClusters);
178 fBenchmark->Stop("ITS Clusterer");
180 AliHLTITSVertexerZ vertexer;
181 AliESDVertex *vertex = vertexer.FindVertexForCurrentEvent(geom,treeClusters);
183 Double_t vtxErr[3]={0.005,0.005,0.010};
184 vertex->GetXYZ(vtxPos);
185 // vertex->GetSigmaXYZ(vtxErr);
186 esd->SetVertex(vertex);
188 // TPC Hough reconstruction
189 Float_t ptmin = 0.1*AliHLTTransform::GetSolenoidField();
190 Float_t zvertex = vtxPos[2];
192 // Run the Hough Transformer
193 fBenchmark->Start("Init");
194 AliHLTHough *hough1 = new AliHLTHough();
196 hough1->SetThreshold(4);
197 hough1->CalcTransformerParams(ptmin);
198 hough1->SetPeakThreshold(70,-1);
199 // printf("Pointer is %x\n",ptr);
200 hough1->Init("./", kFALSE, 100, kFALSE,4,0,(Char_t*)ptr,zvertex);
201 hough1->SetAddHistograms();
202 fBenchmark->Stop("Init");
204 fBenchmark->Start("Init");
205 AliHLTHough *hough2 = new AliHLTHough();
207 hough2->SetThreshold(4);
208 hough2->CalcTransformerParams(ptmin);
209 hough2->SetPeakThreshold(70,-1);
210 // printf("Pointer is %x\n",ptr);
211 hough2->Init("./", kFALSE, 100, kFALSE,4,0,(Char_t*)ptr,zvertex);
212 hough2->SetAddHistograms();
213 fBenchmark->Stop("Init");
215 Int_t nglobaltracks = 0;
217 hough1->StartProcessInThread(0,17);
218 hough2->StartProcessInThread(18,35);
220 // gSystem->Sleep(20000);
221 if(hough1->WaitForThreadFinish())
222 ::Fatal("AliHLTHough::WaitForThreadFinish"," Can not join the required thread! ");
223 if(hough2->WaitForThreadFinish())
224 ::Fatal("AliHLTHough::WaitForThreadFinish"," Can not join the required thread! ");
226 gSystem->MakeDirectory("hough1");
227 hough1->WriteTracks("./hough1");
228 gSystem->MakeDirectory("hough2");
229 hough2->WriteTracks("./hough2");
232 for(int slice=0; slice<=17; slice++)
234 // cout<<"Processing slice "<<slice<<endl;
235 fBenchmark->Start("ReadData");
236 hough1->ReadData(slice,0);
237 fBenchmark->Stop("ReadData");
238 fBenchmark->Start("Transform");
240 fBenchmark->Stop("Transform");
241 hough1->AddAllHistogramsRows();
242 hough1->FindTrackCandidatesRow();
243 fBenchmark->Start("AddTracks");
245 fBenchmark->Stop("AddTracks");
247 // AliHLTTrackArray* tracks = (AliHLTTrackArray*)hough1->GetTracks(0);
248 // nglobaltracks += tracks->GetNTracks();
250 for(int slice=18; slice<=35; slice++)
252 // cout<<"Processing slice "<<slice<<endl;
253 fBenchmark->Start("ReadData");
254 hough2->ReadData(slice,0);
255 fBenchmark->Stop("ReadData");
256 fBenchmark->Start("Transform");
258 fBenchmark->Stop("Transform");
259 hough2->AddAllHistogramsRows();
260 hough2->FindTrackCandidatesRow();
261 fBenchmark->Start("AddTracks");
263 fBenchmark->Stop("AddTracks");
265 // AliHLTTrackArray* tracks = (AliHLTTrackArray*)hough2->GetTracks(0);
266 // nglobaltracks += tracks->GetNTracks();
269 nglobaltracks += hough1->FillESD(esd);
270 nglobaltracks += hough2->FillESD(esd);
273 AliHLTITStracker itsTracker(0);
274 itsTracker.SetVertex(vtxPos,vtxErr);
276 itsTracker.LoadClusters(treeClusters);
277 itsTracker.Clusters2Tracks(esd);
278 itsTracker.UnloadClusters();
280 fBenchmark->Stop("Overall timing");
282 if (file) fprintf(file, "%s\n", time.AsString());
283 if (file) fprintf(file, "run: %d event: %d %d\n",
284 rawReader.GetRunNumber(),
285 rawReader.GetEventId()[0],
286 rawReader.GetEventId()[1]);
287 if (errorCode) fprintf(file, "ERROR: %d\n", errorCode);
289 if (file) fprintf(file, "Hough Transformer found %d tracks\n",nglobaltracks);
291 hough1->DoBench("hough1");
292 hough2->DoBench("hough2");
293 fBenchmark->Analyze("overall");
295 FILE* bench = fopen("hough1.dat", "r");
296 while (bench && !feof(bench)) {
298 if (!fgets(buffer, 256, bench)) break;
299 fprintf(file, "%s", buffer);
304 FILE* bench = fopen("hough2.dat", "r");
305 while (bench && !feof(bench)) {
307 if (!fgets(buffer, 256, bench)) break;
308 fprintf(file, "%s", buffer);
313 FILE* bench = fopen("overall.dat", "r");
314 while (bench && !feof(bench)) {
316 if (!fgets(buffer, 256, bench)) break;
317 fprintf(file, "%s", buffer);
320 fprintf(file, "\n\n");
331 // read run, event, detector, DDL numbers and data size
332 AliRawReaderDate rawReader(ptr);
334 printf("\n%s\n", time.AsString());
335 printf("run: %d event: %d %d\n", rawReader.GetRunNumber(),
336 rawReader.GetEventId()[0], rawReader.GetEventId()[1]);
337 while (rawReader.ReadMiniHeader()) {
338 printf(" detector: %d DDL: %d size: %d\n",
339 rawReader.GetDetectorID(), rawReader.GetDDLID(),
340 rawReader.GetDataSize());
345 gSystem->Sleep(100); // sleep for 0.1 second
348 gSystem->ProcessEvents();
349 if (file) fflush(file);
353 gSystem->RemoveSignalHandler(handler);
354 if (file) fclose(file);
360 int main(int /*argc*/, char** /*argv*/)
362 ::Fatal("main", "this program was compiled without DATE");