First version of kdtree (Alexander, Marian)
[u/mrichter/AliRoot.git] / MONITOR / monitorGDC.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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  **************************************************************************/
15
16 /* $Id$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 // this program performs local monitoring on a GDC by running the HLT code   //
21 //                                                                           //
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.                            //
25 //                                                                           //
26 ///////////////////////////////////////////////////////////////////////////////
27
28 #include <TError.h>
29 #include <TSysEvtHandler.h>
30 #ifdef ALI_DATE
31 #include <TROOT.h>
32 #include <TSystem.h>
33 #include <TDatime.h>
34 #include "AliRawReaderDate.h"
35 #include "event.h"
36 #include "monitor.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 "AliESDVertex.h"
46 #include <AliKalmanTrack.h>
47 #include "AliITSgeomTGeo.h"
48 #include "AliITSgeom.h"
49 #include "AliMagF.h"
50 #include "AliMagFMaps.h"
51 #include <AliHLTITSclusterer.h>
52 #include <AliHLTITSVertexerZ.h>
53 #include <AliHLTITStracker.h>
54 #endif
55
56 //_____________________________________________________________________________
57 class AliGDCInterruptHandler : public TSignalHandler {
58 public:
59   AliGDCInterruptHandler();
60   Bool_t Notify() {fStop = kTRUE; return kTRUE;};
61   Bool_t Stop() const {return fStop;};
62 private:
63   Bool_t fStop;  // CTRL-C pressed
64 };
65
66 //_____________________________________________________________________________
67 AliGDCInterruptHandler::AliGDCInterruptHandler() : 
68   TSignalHandler(kSigInterrupt, kFALSE) 
69 {
70   fStop = kFALSE;
71 };
72
73
74 //_____________________________________________________________________________
75 #ifdef ALI_DATE
76 int main(int argc, char** argv)
77 {
78   // set ROOT in batch mode
79   gROOT->SetBatch();   
80
81   // open a log file
82   FILE* file = fopen("monitorGDC.log", "w");
83
84   TDatime time;
85
86   // get data from a file or online from this node
87   Int_t status = 0;
88   if (argc > 1) {
89     status = monitorSetDataSource(argv[1]);
90   } else {
91     status = monitorSetDataSource(":");
92   }
93   if (status) ::Fatal("monitorSetDataSource", monitorDecodeError(status));
94
95   // monitor only a sample of physics events
96   char* table[] = {"Physics event", "yes", NULL};
97   status = monitorDeclareTable(table);
98   if (status) ::Fatal("monitorDeclareTable", monitorDecodeError(status));
99
100   // declare this monitoring program to DATE
101   status = monitorDeclareMp("GDC physics monitoring");
102   if (status) ::Fatal("monitorDeclareMp", monitorDecodeError(status));
103
104   // initialize HLT transformations
105   if (!AliHLTTransform::Init("./", kFALSE)) {
106     ::Fatal("AliHLTTransform::Init", "HLT initialization failed");
107   }
108   AliESDEvent *esd = new AliESDEvent;
109   esd->CreateStdContent();
110   //  AliKalmanTrack::SetConvConst(
111   //     1000/0.299792458/AliHLTTransform::GetSolenoidField()
112   //  );
113   AliITSgeom *geom = new AliITSgeom();
114   geom->ReadNewFile("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymmFMD.det");
115   if (!geom) return 1;
116   Int_t sfield = 0;
117   switch ((Int_t)(AliHLTTransform::GetSolenoidField()+0.5)) {
118   case 2:
119     sfield = AliMagFMaps::k2kG;
120     break;
121   case 4:
122     sfield = AliMagFMaps::k4kG;
123     break;
124   case 5:
125     sfield = AliMagFMaps::k5kG;
126     break;
127   default:
128     ::Fatal("AliHLTTransform::GetSolenoidField", "Incorrect magnetic field");
129   }
130   AliMagF* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., sfield);
131   AliTracker::SetFieldMap(field,kTRUE);
132
133   // Init PID
134   AliPID pid;
135
136   // create the signal handler
137   AliGDCInterruptHandler* handler = new AliGDCInterruptHandler;
138   gSystem->AddSignalHandler(handler);
139
140   // endless loop
141   void* ptr = NULL;
142   while (!handler->Stop()) {
143     // get the next event
144     status = monitorGetEventDynamic(&ptr);
145     if (status == (Int_t)MON_ERR_EOF) break;
146     if (status) ::Fatal("monitorGetEventDynamic", monitorDecodeError(status));
147
148     // if no new event
149     if (!ptr) {
150       gSystem->Sleep(1000);   // sleep for 1 second
151       continue;
152     }
153
154     AliRawReaderDate rawReader(ptr);
155     //    if ((rawReader.GetAttributes()[0] & 0x02) != 0) {
156
157     //      Int_t errorCode = rawReader.CheckData();
158     Int_t errorCode = 0;
159       if (errorCode && (errorCode != AliRawReader::kErrSize)) {
160         time.Set();
161         if (file) fprintf(file, "%s\n", time.AsString());
162         if (file) fprintf(file, "run: %d  event: %d %d\n", 
163                           rawReader.GetRunNumber(), 
164                           rawReader.GetEventId()[0], 
165                           rawReader.GetEventId()[1]);
166         fprintf(file, "ERROR: %d\n\n", errorCode);
167
168       } else {
169
170         AliHLTBenchmark *fBenchmark = new AliHLTBenchmark();
171         fBenchmark->Start("Overall timing");
172
173         // ITS clusterer and vertexer
174         fBenchmark->Start("ITS Clusterer");
175         AliHLTITSclusterer clusterer(0);
176         AliRawReader *itsrawreader=new AliRawReaderDate(ptr);
177         TTree* treeClusters = new TTree("TreeL3ITSclusters"," "); //make a tree
178         clusterer.Digits2Clusters(itsrawreader,treeClusters);
179         fBenchmark->Stop("ITS Clusterer");
180         
181         AliHLTITSVertexerZ vertexer;
182         AliESDVertex *vertex = vertexer.FindVertexForCurrentEvent(geom,treeClusters);
183         Double_t vtxPos[3];
184         Double_t vtxErr[3]={0.005,0.005,0.010};
185         vertex->GetXYZ(vtxPos);
186         //      vertex->GetSigmaXYZ(vtxErr);
187         esd->SetVertex(vertex);
188
189         // TPC Hough reconstruction
190         Float_t ptmin = 0.1*AliHLTTransform::GetSolenoidField();
191         Float_t zvertex = vtxPos[2];
192
193         // Run the Hough Transformer
194         fBenchmark->Start("Init");
195         AliHLTHough *hough1 = new AliHLTHough();
196
197         hough1->SetThreshold(4);
198         hough1->CalcTransformerParams(ptmin);
199         hough1->SetPeakThreshold(70,-1);
200         //      printf("Pointer is %x\n",ptr);
201         hough1->Init("./", kFALSE, 100, kFALSE,4,0,(Char_t*)ptr,zvertex);
202         hough1->SetAddHistograms();
203         fBenchmark->Stop("Init");
204
205         fBenchmark->Start("Init");
206         AliHLTHough *hough2 = new AliHLTHough();
207
208         hough2->SetThreshold(4);
209         hough2->CalcTransformerParams(ptmin);
210         hough2->SetPeakThreshold(70,-1);
211         //      printf("Pointer is %x\n",ptr);
212         hough2->Init("./", kFALSE, 100, kFALSE,4,0,(Char_t*)ptr,zvertex);
213         hough2->SetAddHistograms();
214         fBenchmark->Stop("Init");
215
216         Int_t nglobaltracks = 0;
217         /*
218         hough1->StartProcessInThread(0,17);
219         hough2->StartProcessInThread(18,35);
220
221         //      gSystem->Sleep(20000);
222         if(hough1->WaitForThreadFinish())
223           ::Fatal("AliHLTHough::WaitForThreadFinish"," Can not join the required thread! ");
224         if(hough2->WaitForThreadFinish())
225           ::Fatal("AliHLTHough::WaitForThreadFinish"," Can not join the required thread! ");
226
227         gSystem->MakeDirectory("hough1");
228         hough1->WriteTracks("./hough1");
229         gSystem->MakeDirectory("hough2");
230         hough2->WriteTracks("./hough2");
231         */
232
233         for(int slice=0; slice<=17; slice++)
234         {
235           //      cout<<"Processing slice "<<slice<<endl;
236           fBenchmark->Start("ReadData");
237           hough1->ReadData(slice,0);
238           fBenchmark->Stop("ReadData");
239           fBenchmark->Start("Transform");
240           hough1->Transform();
241           fBenchmark->Stop("Transform");
242           hough1->AddAllHistogramsRows();
243           hough1->FindTrackCandidatesRow();
244           fBenchmark->Start("AddTracks");
245           hough1->AddTracks();
246           fBenchmark->Stop("AddTracks");
247
248           //      AliHLTTrackArray* tracks = (AliHLTTrackArray*)hough1->GetTracks(0);
249           //      nglobaltracks += tracks->GetNTracks();
250         }
251         for(int slice=18; slice<=35; slice++)
252         {
253           //      cout<<"Processing slice "<<slice<<endl;
254           fBenchmark->Start("ReadData");
255           hough2->ReadData(slice,0);
256           fBenchmark->Stop("ReadData");
257           fBenchmark->Start("Transform");
258           hough2->Transform();
259           fBenchmark->Stop("Transform");
260           hough2->AddAllHistogramsRows();
261           hough2->FindTrackCandidatesRow();
262           fBenchmark->Start("AddTracks");
263           hough2->AddTracks();
264           fBenchmark->Stop("AddTracks");
265
266           //      AliHLTTrackArray* tracks = (AliHLTTrackArray*)hough2->GetTracks(0);
267           //      nglobaltracks += tracks->GetNTracks();
268         }
269
270         nglobaltracks += hough1->FillESD(esd);
271         nglobaltracks += hough2->FillESD(esd);
272
273         // ITS tracker
274         AliHLTITStracker itsTracker(0);
275         itsTracker.SetVertex(vtxPos,vtxErr);
276
277         itsTracker.LoadClusters(treeClusters);
278         itsTracker.Clusters2Tracks(esd);
279         itsTracker.UnloadClusters();
280
281         fBenchmark->Stop("Overall timing");
282         time.Set();
283         if (file) fprintf(file, "%s\n", time.AsString());
284         if (file) fprintf(file, "run: %d  event: %d %d\n", 
285                           rawReader.GetRunNumber(), 
286                           rawReader.GetEventId()[0], 
287                           rawReader.GetEventId()[1]);
288         if (errorCode) fprintf(file, "ERROR: %d\n", errorCode);
289
290         if (file) fprintf(file, "Hough Transformer found %d tracks\n",nglobaltracks);
291
292         hough1->DoBench("hough1");
293         hough2->DoBench("hough2");
294         fBenchmark->Analyze("overall");
295         if (file) {
296           FILE* bench = fopen("hough1.dat", "r");
297           while (bench && !feof(bench)) {
298             char buffer[256];
299             if (!fgets(buffer, 256, bench)) break;
300             fprintf(file, "%s", buffer);
301           }
302           fclose(bench);
303         }
304         if (file) {
305           FILE* bench = fopen("hough2.dat", "r");
306           while (bench && !feof(bench)) {
307             char buffer[256];
308             if (!fgets(buffer, 256, bench)) break;
309             fprintf(file, "%s", buffer);
310           }
311           fclose(bench);
312         }
313         if (file) {
314           FILE* bench = fopen("overall.dat", "r");
315           while (bench && !feof(bench)) {
316             char buffer[256];
317             if (!fgets(buffer, 256, bench)) break;
318             fprintf(file, "%s", buffer);
319           }
320           fclose(bench);
321           fprintf(file, "\n\n");
322         }
323
324         delete hough1;
325         delete hough2;
326
327         esd->Reset();
328       }
329     //    }
330
331     /*
332     // read run, event, detector, DDL numbers and data size
333     AliRawReaderDate rawReader(ptr);
334     time.Set();
335     printf("\n%s\n", time.AsString());
336     printf("run: %d  event: %d %d\n", rawReader.GetRunNumber(), 
337            rawReader.GetEventId()[0], rawReader.GetEventId()[1]);
338     while (rawReader.ReadMiniHeader()) {
339       printf(" detector: %d   DDL: %d  size: %d\n", 
340              rawReader.GetDetectorID(), rawReader.GetDDLID(), 
341              rawReader.GetDataSize());
342     }
343
344     */
345
346     gSystem->Sleep(100);   // sleep for 0.1 second
347     free(ptr);
348
349     gSystem->ProcessEvents();
350     if (file) fflush(file);
351
352   }
353
354   gSystem->RemoveSignalHandler(handler);
355   if (file) fclose(file);
356
357   return 0;
358 }
359
360 #else
361 int main(int /*argc*/, char** /*argv*/)
362 {
363   ::Fatal("main", "this program was compiled without DATE");
364
365   return 1;
366 }
367 #endif