a9b7ed89c490a4d5c5885a1f5c0239ca2de0814e
[u/mrichter/AliRoot.git] / MONITOR / AliMonitorProcess.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 is the class for perfoming the monitoring process.                  //
21 //  It checks if a raw data file exists, loops over the events in the raw    //
22 //  data file, reconstructs TPC and ITS clusters and tracks, fills the       //
23 //  monitor histograms and sends the updated histograms to the clients.      //
24 //  Then the raw data file is deleted and it waits for a new file.           //
25 //                                                                           //
26 ///////////////////////////////////////////////////////////////////////////////
27
28
29 #include "AliMonitorProcess.h"
30 #include "AliMonitorTPC.h"
31 #include "AliMonitorITS.h"
32 #include "AliMonitorV0s.h"
33 #include "AliMonitorHLT.h"
34 #include "AliMonitorHLTHough.h"
35 #include "AliRawReaderRoot.h"
36 #include "AliLoader.h"
37 #include "AliRun.h"
38 #include "AliTPC.h"
39 #include "AliTPCclustererMI.h"
40 #include "AliTPCtrackerMI.h"
41 #include "AliITS.h"
42 #include "AliITSclustererV2.h"
43 #include "AliITStrackerV2.h"
44 #include "AliITSLoader.h"
45 #include "AliV0vertexer.h"
46 #include <TGrid.h>
47 #include <TSocket.h>
48 #include <TServerSocket.h>
49 #include <TMessage.h>
50 #include <TGridResult.h>
51 #include <TROOT.h>
52 #ifdef ALI_HLT
53 #include <AliLevel3.h>
54 #include <AliL3Transform.h>
55 #include <AliL3Track.h>
56 #include <AliL3TrackArray.h>
57 #include <AliL3StandardIncludes.h>
58 #include <AliL3HoughMaxFinder.h>
59 #include <AliL3HoughBaseTransformer.h>
60 #include <AliL3Hough.h>
61 #include <AliL3ClusterFitter.h>
62 #include <AliL3Vertex.h>
63 #include <AliL3Fitter.h>
64 #include <AliL3DDLDataFileHandler.h>
65 #endif
66
67 ClassImp(AliMonitorProcess) 
68
69
70 const Int_t AliMonitorProcess::fgkPort = 9327;
71
72
73 //_____________________________________________________________________________
74 AliMonitorProcess::AliMonitorProcess(
75 #if ROOT_VERSION_CODE <= 199169   // 3.10/01
76                                      const char* /*alienHost*/,
77 #else
78                                      const char* alienHost,
79 #endif
80                                      const char* alienDir,
81                                      const char* fileNameGalice)
82 {
83 // initialize the monitoring process and the monitor histograms
84
85 #if ROOT_VERSION_CODE <= 199169   // 3.10/01
86   fGrid = TGrid::Connect("alien", gSystem->Getenv("USER"));
87 #else
88   fGrid = TGrid::Connect(alienHost, gSystem->Getenv("USER"));
89 #endif
90   if (!fGrid || fGrid->IsZombie() || !fGrid->IsConnected()) {
91     delete fGrid;
92     Fatal("AliMonitorProcess", "could not connect to alien");
93   }
94 #if ROOT_VERSION_CODE <= 199169   // 3.10/01
95   fGrid->cd(alienDir);
96 #else
97   fAlienDir = alienDir;
98 #endif
99   fLogicalFileName = "";
100   fFileName = "";
101
102   fRunLoader = AliRunLoader::Open(fileNameGalice);
103   if (!fRunLoader) Fatal("AliMonitorProcess", 
104                          "could not get run loader from file %s", 
105                          fileNameGalice);
106
107   fRunLoader->CdGAFile();
108   fTPCParam = AliTPC::LoadTPCParam(gFile);
109   if (!fTPCParam) Fatal("AliMonitorProcess", "could not load TPC parameters");
110
111   fRunLoader->LoadgAlice();
112   gAlice = fRunLoader->GetAliRun();
113   if (!gAlice) Fatal("AliMonitorProcess", "no gAlice object found");
114   AliITS* its = (AliITS*) gAlice->GetModule("ITS");
115   if (!its) Fatal("AliMonitorProcess", "no ITS detector found");
116   fITSgeom = its->GetITSgeom();
117   if (!fITSgeom) Fatal("AliMonitorProcess", "could not load ITS geometry");
118
119 #ifdef ALI_HLT
120 // Init TPC parameters for HLT
121   Bool_t isinit=AliL3Transform::Init(const_cast<char*>(fileNameGalice),kTRUE);
122   if(!isinit){
123     cerr << "Could not create transform settings, please check log for error messages!" << endl;
124     return;
125   }
126 #endif
127
128   fRunNumber = 0;
129   fSubRunNumber = 0;
130   fNEvents = 0;
131   fNEventsMin = 1;
132   fWriteHistoList = kFALSE;
133
134   fTopFolder = new TFolder("Monitor", "monitor histograms");
135   fTopFolder->SetOwner(kTRUE);
136
137   fMonitors.Add(new AliMonitorTPC(fTPCParam));
138   fMonitors.Add(new AliMonitorITS(fITSgeom));
139   fMonitors.Add(new AliMonitorV0s);
140 #ifdef ALI_HLT
141   fMonitors.Add(new AliMonitorHLT(fTPCParam));
142   fMonitors.Add(new AliMonitorHLTHough(fTPCParam));
143 #endif
144
145   for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
146     ((AliMonitor*) fMonitors[iMonitor])->CreateHistos(fTopFolder);
147   }
148
149   TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
150   while (TFolder* folder = (TFolder*) iFolder->Next()) folder->SetOwner(kTRUE);
151   delete iFolder;
152
153   fFile = TFile::Open("monitor_tree.root", "RECREATE");
154   if (!fFile || !fFile->IsOpen()) {
155     Fatal("AliMonitorProcess", "could not open file for tree");
156   }
157   fTree = new TTree("MonitorTree", "tree for monitoring");
158   for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
159     ((AliMonitor*) fMonitors[iMonitor])->CreateBranches(fTree);
160   }
161   gROOT->cd();
162
163   fServerSocket = new TServerSocket(fgkPort, kTRUE);
164   fServerSocket->SetOption(kNoBlock, 1);
165   fDisplaySocket = NULL;
166   CheckForConnections();
167 #ifdef ALI_HLT
168   fHLT = NULL;
169 #endif
170
171   SetStatus(kStopped);
172   fStopping = kFALSE;
173
174   fInterruptHandler = new AliMonitorInterruptHandler(this);
175   gSystem->AddSignalHandler(fInterruptHandler);
176 }
177
178 //_____________________________________________________________________________
179 AliMonitorProcess::AliMonitorProcess(const AliMonitorProcess& process) :
180   TObject(process)
181 {
182   Fatal("AliMonitorProcess", "copy constructor not implemented");
183 }
184
185 //_____________________________________________________________________________
186 AliMonitorProcess& AliMonitorProcess::operator = (const AliMonitorProcess& 
187                                                   /*process*/)
188 {
189   Fatal("operator =", "assignment operator not implemented");
190   return *this;
191 }
192
193 //_____________________________________________________________________________
194 AliMonitorProcess::~AliMonitorProcess()
195 {
196 // clean up
197
198   fMonitors.Delete();
199   delete fTopFolder;
200   delete fRunLoader;
201
202   delete fServerSocket;
203   fSockets.Delete();
204   delete fDisplaySocket;
205
206 #if ROOT_VERSION_CODE <= 199169   // 3.10/01
207   fGrid->Close();
208 #endif
209   delete fGrid;
210
211   fFile->Close();
212   delete fFile;
213   gSystem->Unlink("monitor_tree.root");
214
215 #ifdef ALI_HLT
216   delete fHLT;
217   delete fHLTHough;
218 #endif
219
220   gSystem->RemoveSignalHandler(fInterruptHandler);
221   delete fInterruptHandler;
222 }
223
224
225 //_____________________________________________________________________________
226 const char* AliMonitorProcess::GetRevision()
227 {
228   return "$Revision$";
229 }
230
231
232 //_____________________________________________________________________________
233 void AliMonitorProcess::SetStatus(EStatus status)
234 {
235 // set the current status and process system events
236
237   fStatus = status;
238   gSystem->ProcessEvents();
239 }
240
241
242 //_____________________________________________________________________________
243 void AliMonitorProcess::Run()
244 {
245 // run the monitor process: 
246 //  check for a raw data file, process the raw data file and delete it
247
248   fStopping = kFALSE;
249
250   while (!fStopping) {
251     SetStatus(kWaiting);
252     while (!CheckForNewFile()) {
253       CheckForConnections();
254       SetStatus(kWaiting);
255       if (fStopping) break;
256       gSystem->Sleep(10);
257     }
258     if (fStopping) break;
259
260     ProcessFile();
261   }
262
263   WriteHistos();
264
265   fStopping = kFALSE;
266   SetStatus(kStopped);
267 }
268
269
270 //_____________________________________________________________________________
271 void AliMonitorProcess::Stop()
272 {
273 // set the fStopping flag to terminate the monitor process after the current
274 // event was processed
275
276   if (GetStatus() != kStopped) fStopping = kTRUE;
277 }
278
279
280 //_____________________________________________________________________________
281 void AliMonitorProcess::ProcessFile(const char* fileName)
282 {
283 // create a file with monitor histograms for a single file
284
285   if (GetStatus() != kStopped) {
286     Error("ProcessFile", "ProcessFile can not be called"
287           " while the monitor process is running");
288     return;
289   }
290
291   fFileName = fileName;
292   Int_t nEventMin = fNEventsMin;
293   fNEventsMin = 1;
294   ProcessFile();
295   WriteHistos();
296   fNEventsMin = nEventMin;
297   SetStatus(kStopped);
298 }
299
300
301 //_____________________________________________________________________________
302 Bool_t AliMonitorProcess::CheckForNewFile()
303 {
304 // check whether a new file was registered in alien
305
306 #if ROOT_VERSION_CODE <= 199169   // 3.10/01
307   TGridResult* result = fGrid->Ls();
308 #else
309 //  Grid_ResultHandle_t handle = fGrid->OpenDir(fAlienDir);
310   TDatime datime;
311   char findName[256];
312   sprintf(findName, "*_%d_*.root", datime.GetDate());
313   Grid_ResultHandle_t handle = fGrid->Find(fAlienDir, findName);
314   if (!handle) {
315     Error("CheckForNewFile", "could not open alien directory %s", 
316           fAlienDir.Data());
317     return kFALSE;
318   }
319   TGridResult* result = fGrid->CreateGridResult(handle);
320 #endif
321   Long_t maxDate = -1;
322   Long_t maxTime = -1;
323   TString fileName;
324
325 #if ROOT_VERSION_CODE <= 199169   // 3.10/01
326   while (const char* entry = result->Next()) {
327 #else
328   while (Grid_Result_t* resultEntry = result->Next()) {
329     const char* entry = resultEntry->name.c_str();
330 #endif
331     if (rindex(entry, '/')) entry = rindex(entry, '/')+1;
332     // entry = host_date_time.root
333     TString entryCopy(entry);
334     char* p = const_cast<char*>(entryCopy.Data());
335     if (!strtok(p, "_") || !p) continue;  // host name
336     char* dateStr = strtok(NULL, "_");
337     if (!dateStr || !p) continue;
338     char* timeStr = strtok(NULL, ".");
339     if (!timeStr || !p) continue;
340     Long_t date = atoi(dateStr);
341     Long_t time = atoi(timeStr);
342
343     if ((date > maxDate) || ((date == maxDate) && (time > maxTime))) {
344       maxDate = date;
345       maxTime = time;
346       fileName = entry;
347     }
348   }
349
350   delete result;
351   if (maxDate < 0) return kFALSE;  // no files found
352   if (fLogicalFileName.CompareTo(fileName) == 0) return kFALSE;  // no new file
353
354   fLogicalFileName = fileName;
355 #if ROOT_VERSION_CODE <= 199169   // 3.10/01
356   result = fGrid->GetPhysicalFileNames(fLogicalFileName.Data());
357   fFileName = result->Next();
358 #else
359   fileName = fAlienDir + "/" + fLogicalFileName;
360   handle = fGrid->GetPhysicalFileNames(fileName.Data());
361   if (!handle) {
362     Error("CheckForNewFile", "could not get physical file names for %s", 
363           fileName.Data());
364     return kFALSE;
365   }
366   result = fGrid->CreateGridResult(handle);
367   result->Reset();
368   Grid_Result_t* resultEntry = result->Next();
369   if (!resultEntry) {
370     Error("CheckForNewFile", "could not get physical file names for %s", 
371           fileName.Data());
372     return kFALSE;
373   }
374   fFileName = resultEntry->name2.c_str();
375   fFileName.ReplaceAll("castor:/", "rfio:/");
376 #endif
377   delete result;
378
379   return kTRUE;
380 }
381
382 //_____________________________________________________________________________
383 Bool_t AliMonitorProcess::ProcessFile()
384 {
385 // loop over all events in the raw data file, run the reconstruction
386 // and fill the monitor histograms
387
388   Int_t nEvents = GetNumberOfEvents(fFileName);
389   if (nEvents <= 0) return kFALSE;
390   Info("ProcessFile", "found %d event(s) in file %s", 
391        nEvents, fFileName.Data());
392 #ifdef ALI_HLT
393   CreateHLT(fFileName);
394   CreateHLTHough(fFileName);
395 #endif
396
397   // loop over the events
398   for (Int_t iEvent = 0; iEvent < nEvents; iEvent++) {
399     CheckForConnections();
400     SetStatus(kReading);
401     fRunLoader->SetEventNumber(0);
402     AliRawReaderRoot rawReader(fFileName, iEvent);
403     if (fStopping) break;
404     if (rawReader.GetRunNumber() != fRunNumber) {
405       WriteHistos();
406       StartNewRun();
407       fRunNumber = rawReader.GetRunNumber();
408       fEventNumber[0] = rawReader.GetEventId()[0];
409       fEventNumber[1] = rawReader.GetEventId()[1];
410       fSubRunNumber = 0;
411       if (fStopping) break;
412     }
413
414     // monitor only central physics events
415     if (rawReader.GetType() != 7) continue;
416     if ((rawReader.GetAttributes()[0] & 0x02) == 0) continue;
417     Info("ProcessFile", "run: %d  event: %d %d\n", rawReader.GetRunNumber(), 
418          rawReader.GetEventId()[0], rawReader.GetEventId()[1]);
419
420     CheckForConnections();
421     if (!ReconstructTPC(&rawReader)) return kFALSE;
422     if (fStopping) break;
423     CheckForConnections();
424     if (!ReconstructITS(&rawReader)) return kFALSE;
425     if (fStopping) break;
426     CheckForConnections();
427     if (!ReconstructV0s()) return kFALSE;
428     if (fStopping) break;
429     CheckForConnections();
430     if (!ReconstructHLT(iEvent)) return kFALSE;
431     if (fStopping) break;
432     CheckForConnections();
433     if (!ReconstructHLTHough(iEvent)) return kFALSE;
434     if (fStopping) break;
435
436     if (fDisplaySocket) fDisplaySocket->Send("new event");
437
438     Info("ProcessFile", "filling histograms...");
439     for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
440       CheckForConnections();
441       SetStatus(kFilling);
442       ((AliMonitor*) fMonitors[iMonitor])->FillHistos(fRunLoader, &rawReader);
443       if (fStopping) break;
444     }
445     if (fStopping) break;
446
447     Info("ProcessFile", "updating histograms...");
448     CheckForConnections();
449     SetStatus(kUpdating);
450     TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
451     while (TFolder* folder = (TFolder*) iFolder->Next()) {
452       TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
453       while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
454         histo->Update();
455       }
456       delete iHisto;
457     }
458     delete iFolder;
459     if (fStopping) break;
460
461     Info("ProcessFile", "filling the tree...");
462     fTree->Fill();
463
464     Info("ProcessFile", "broadcasting histograms...");
465     CheckForConnections();
466     BroadcastHistos();
467
468     fNEvents++;
469     if (fStopping) break;
470   }
471
472 #ifdef ALI_HLT
473   delete fHLT;
474   fHLT = NULL;
475 #endif
476
477   return kTRUE;
478 }
479
480 //_____________________________________________________________________________
481 void AliMonitorProcess::Reset()
482 {
483 // write the current histograms to a file and reset them
484
485   if (fSubRunNumber == 0) fSubRunNumber++;
486   WriteHistos();
487   StartNewRun();
488   fSubRunNumber++;
489 }
490
491
492 //_____________________________________________________________________________
493 UInt_t AliMonitorProcess::GetEventPeriodNumber() const
494 {
495 // get the period number from the event id
496
497   return (fEventNumber[1] >> 4);
498 }
499
500 //_____________________________________________________________________________
501 UInt_t AliMonitorProcess::GetEventOrbitNumber() const
502 {
503 // get the orbit number from the event id
504
505   return ((fEventNumber[1] & 0x000F) << 20) + (fEventNumber[0] >> 12);
506 }
507
508 //_____________________________________________________________________________
509 UInt_t AliMonitorProcess::GetEventBunchNumber() const
510 {
511 // get the bunch number from the event id
512
513   return (fEventNumber[0] % 0x0FFF);
514 }
515
516 //_____________________________________________________________________________
517 Int_t AliMonitorProcess::GetNumberOfEvents(const char* fileName) const
518 {
519 // determine the number of events in the given raw data file
520
521   Int_t nEvents = -1;
522
523   TFile* file = TFile::Open(fileName);
524   if (!file || !file->IsOpen()) {
525     Error("GetNumberOfEvents", "could not open file %s", fileName);
526     if (file) delete file;
527     return -1;
528   }
529
530   TTree* tree = (TTree*) file->Get("RAW");
531   if (!tree) {
532     Error("GetNumberOfEvents", "could not find tree with raw data");
533   } else {
534     nEvents = (Int_t) tree->GetEntries();
535   }
536   file->Close();
537   delete file;
538
539   return nEvents;
540 }
541
542 //_____________________________________________________________________________
543 Bool_t AliMonitorProcess::ReconstructTPC(AliRawReader* rawReader)
544 {
545 // find TPC clusters and tracks
546
547   SetStatus(kRecTPC);
548
549   AliLoader* tpcLoader = fRunLoader->GetLoader("TPCLoader");
550   if (!tpcLoader) {
551     Error("ReconstructTPC", "no TPC loader found");
552     return kFALSE;
553   }
554   gSystem->Unlink("TPC.RecPoints.root");
555   gSystem->Unlink("TPC.Tracks.root");
556
557   // cluster finder
558   Info("ReconstructTPC", "reconstructing clusters...");
559   tpcLoader->LoadRecPoints("recreate");
560   AliTPCclustererMI clusterer(fTPCParam);
561   tpcLoader->MakeRecPointsContainer();
562   clusterer.SetOutput(tpcLoader->TreeR());
563   clusterer.Digits2Clusters(rawReader);
564   tpcLoader->WriteRecPoints("OVERWRITE");
565
566   // track finder
567   Info("ReconstructTPC", "reconstructing tracks...");
568   tpcLoader->LoadTracks("recreate");
569   {
570     AliTPCtrackerMI tracker(fTPCParam);
571     tracker.Clusters2Tracks();
572     tracker.UnloadClusters();
573   }
574
575   tpcLoader->UnloadRecPoints();
576   tpcLoader->UnloadTracks();
577   return kTRUE;
578 }
579
580 //_____________________________________________________________________________
581 Bool_t AliMonitorProcess::ReconstructITS(AliRawReader* rawReader)
582 {
583 // find ITS clusters and tracks
584
585   SetStatus(kRecITS);
586
587   AliLoader* itsLoader = fRunLoader->GetLoader("ITSLoader");
588   if (!itsLoader) {
589     Error("ReconstructITS", "no ITS loader found");
590     return kFALSE;
591   }
592   AliLoader* tpcLoader = fRunLoader->GetLoader("TPCLoader");
593   if (!tpcLoader) {
594     Error("ReconstructITS", "no TPC loader found");
595     return kFALSE;
596   }
597   gSystem->Unlink("ITS.RecPoints.root");
598   gSystem->Unlink("ITS.Tracks.root");
599
600   // cluster finder
601   Info("ReconstructITS", "reconstructing clusters...");
602   itsLoader->LoadRecPoints("recreate");
603   AliITSclustererV2 clusterer(fITSgeom);
604   itsLoader->MakeRecPointsContainer();
605   clusterer.Digits2Clusters(rawReader);
606
607   // track finder
608   Info("ReconstructITS", "reconstructing tracks...");
609   itsLoader->LoadTracks("recreate");
610   itsLoader->MakeTracksContainer();
611   tpcLoader->LoadTracks();
612   AliITStrackerV2 tracker(fITSgeom);
613   tracker.LoadClusters(itsLoader->TreeR());
614   tracker.Clusters2Tracks(tpcLoader->TreeT(), itsLoader->TreeT());
615   tracker.UnloadClusters();
616   itsLoader->WriteTracks("OVERWRITE");
617
618   itsLoader->UnloadRecPoints();
619   itsLoader->UnloadTracks();
620   tpcLoader->UnloadTracks();
621   return kTRUE;
622 }
623
624 //_____________________________________________________________________________
625 Bool_t AliMonitorProcess::ReconstructV0s()
626 {
627 // find V0s
628
629   SetStatus(kRecV0s);
630
631   AliITSLoader* itsLoader = (AliITSLoader*) fRunLoader->GetLoader("ITSLoader");
632   if (!itsLoader) {
633     Error("ReconstructV0", "no ITS loader found");
634     return kFALSE;
635   }
636   gSystem->Unlink("ITS.V0s.root");
637
638   // V0 finder
639   Info("ReconstructV0s", "reconstructing V0s...");
640   itsLoader->LoadTracks("read");
641   itsLoader->LoadV0s("recreate");
642   AliV0vertexer vertexer;
643   TTree* tracks = itsLoader->TreeT();
644   if (!tracks) {
645     Error("ReconstructV0s", "no ITS tracks tree found");
646     return kFALSE;
647   }
648   if (!itsLoader->TreeV0()) itsLoader->MakeTree("V0");
649   TTree* v0s = itsLoader->TreeV0();
650   vertexer.Tracks2V0vertices(tracks, v0s);
651   itsLoader->WriteV0s("OVERWRITE");
652
653   itsLoader->UnloadTracks();
654   itsLoader->UnloadV0s();
655   return kTRUE;
656 }
657
658 //_____________________________________________________________________________
659 #ifdef ALI_HLT
660 void AliMonitorProcess::CreateHLT(const char* fileName)
661 {
662
663 // create the HLT (Level3) object
664
665   if (fHLT) delete fHLT;
666
667   char name[256];
668   strcpy(name, fileName);
669   fHLT = new AliLevel3(name);
670   fHLT->Init("./", AliLevel3::kRaw, 1);
671
672   fHLT->SetClusterFinderParam(-1, -1, kTRUE);
673   
674   Int_t phiSegments = 50;
675   Int_t etaSegments = 100;
676   Int_t trackletlength = 3;
677   Int_t tracklength = 20;//40 or 5
678   Int_t rowscopetracklet = 2;
679   Int_t rowscopetrack = 10;
680   Double_t minPtFit = 0;
681   Double_t maxangle = 0.1745;
682   Double_t goodDist = 5;
683   Double_t maxphi = 0.1;
684   Double_t maxeta = 0.1;
685   Double_t hitChi2Cut = 15;//100 or 15
686   Double_t goodHitChi2 = 5;//20 or 5
687   Double_t trackChi2Cut = 10;//50 or 10
688   fHLT->SetTrackerParam(phiSegments, etaSegments, 
689                         trackletlength, tracklength,
690                         rowscopetracklet, rowscopetrack,
691                         minPtFit, maxangle, goodDist, hitChi2Cut,
692                         goodHitChi2, trackChi2Cut, 50, maxphi, maxeta, kTRUE);
693   
694   fHLT->WriteFiles("./hlt/");  
695 }
696
697 //_____________________________________________________________________________
698 void AliMonitorProcess::CreateHLTHough(const char* fileName)
699 {
700
701 // create the HLT Hough transform (L3Hough) object
702
703   if (fHLTHough) delete fHLTHough;
704
705   char name[256];
706   strcpy(name, fileName);
707
708   fHLTHough = new AliL3Hough();
709   fHLTHough->SetThreshold(4);
710   fHLTHough->SetTransformerParams(140,150,0.5,-1);
711   fHLTHough->SetPeakThreshold(9000,-1);// or 6000
712   fHLTHough->Init("./", kFALSE, 50, kFALSE,0,name);
713   fHLTHough->SetAddHistograms();
714   //  fHLTHough->GetMaxFinder()->SetThreshold(14000);
715
716 }
717 #endif
718
719 //_____________________________________________________________________________
720 Bool_t AliMonitorProcess::ReconstructHLT(
721 #ifdef ALI_HLT
722   Int_t iEvent
723 #else
724   Int_t /* iEvent */
725 #endif
726 )
727 {
728 // run the HLT cluster and track finder
729
730   SetStatus(kRecHLT);
731
732 #ifndef ALI_HLT
733   Warning("ReconstructHLT", "the code was compiled without HLT support");
734   return kTRUE;
735
736 #else
737   gSystem->Exec("rm -rf hlt");
738   gSystem->MakeDirectory("hlt");
739   if (!fHLT) return kFALSE;
740
741   fHLT->ProcessEvent(0, 35, iEvent);
742
743   // remove the event number from the file names
744   char command[256];
745   sprintf(command, "rename points_%d points hlt/*.raw", iEvent);
746   gSystem->Exec(command);
747   sprintf(command, "rename tracks_tr_%d tracks_tr hlt/*.raw", iEvent);
748   gSystem->Exec(command);
749   sprintf(command, "rename tracks_gl_%d tracks_gl hlt/*.raw", iEvent);
750   gSystem->Exec(command);
751   sprintf(command, "rename tracks_%d tracks hlt/*.raw", iEvent);
752   gSystem->Exec(command);
753   return kTRUE;
754 #endif
755 }
756
757 //_____________________________________________________________________________
758 Bool_t AliMonitorProcess::ReconstructHLTHough(
759 #ifdef ALI_HLT
760   Int_t iEvent
761 #else
762   Int_t /* iEvent */
763 #endif
764 )
765 {
766 // run the HLT Hough transformer
767
768   SetStatus(kRecHLT);
769
770 #ifndef ALI_HLT
771   Warning("ReconstructHLTHough", "the code was compiled without HLT support");
772   return kTRUE;
773
774 #else
775   gSystem->Exec("rm -rf hlt/hough");
776   gSystem->MakeDirectory("hlt/hough");
777   gSystem->Exec("rm -rf hlt/fitter");
778   gSystem->MakeDirectory("hlt/fitter");
779   if (!fHLTHough) return kFALSE;
780
781   //  fHLTHough->Process(0, 35);
782   // Loop over TPC sectors and process the data
783   for(Int_t i=0; i<=35; i++)
784     {
785       fHLTHough->ReadData(i,iEvent);
786       fHLTHough->Transform();
787       //      if(fHLTHough->fAddHistograms)
788       fHLTHough->AddAllHistograms();
789       fHLTHough->FindTrackCandidates();
790       fHLTHough->AddTracks();
791     }
792   fHLTHough->WriteTracks("./hlt/hough");
793
794   // Run cluster fitter
795   AliL3ClusterFitter *fitter = new AliL3ClusterFitter("./hlt");
796
797   // Set debug flag for the cluster fitter
798   //  fitter->Debug();
799
800   // Setting fitter parameters
801   fitter->SetInnerWidthFactor(1,1.5);
802   fitter->SetOuterWidthFactor(1,1.5);
803   fitter->SetNmaxOverlaps(5);
804   
805   //fitter->SetChiSqMax(5,kFALSE); //isolated clusters
806   fitter->SetChiSqMax(5,kTRUE);  //overlapping clusters
807
808   Int_t rowrange[2] = {0,AliL3Transform::GetNRows()-1};
809
810   // Takes input from global hough tracks produced by HT
811   fitter->LoadSeeds(rowrange,kFALSE,iEvent);
812
813   UInt_t ndigits;
814
815   for(Int_t islice = 0; islice <= 35; islice++)
816     {
817       for(Int_t ipatch = 0; ipatch < AliL3Transform::GetNPatches(); ipatch++)
818         {
819           // Read digits
820           fHLTHough->GetMemHandler(ipatch)->Free();
821           fHLTHough->GetMemHandler(ipatch)->Init(islice,ipatch);
822           AliL3DigitRowData *digits = (AliL3DigitRowData *)fHLTHough->GetMemHandler(ipatch)->AliAltroDigits2Memory(ndigits,iEvent);
823
824           fitter->Init(islice,ipatch);
825           fitter->SetInputData(digits);
826           fitter->FindClusters();
827           fitter->WriteClusters();
828         }
829     }
830
831   // Refit of the clusters
832   AliL3Vertex vertex;
833   //The seeds are the input tracks from circle HT
834   AliL3TrackArray *tracks = fitter->GetSeeds();
835   AliL3Fitter *ft = new AliL3Fitter(&vertex,1);
836
837   ft->LoadClusters("./hlt/fitter/",iEvent,kFALSE);
838   for(Int_t i=0; i<tracks->GetNTracks(); i++)
839     {
840       AliL3Track *track = tracks->GetCheckedTrack(i);
841       if(!track) continue;
842       if(track->GetNHits() < 20) continue;
843       ft->SortTrackClusters(track);
844       ft->FitHelix(track);
845       ft->UpdateTrack(track);
846     }
847   delete ft;
848         
849   //Write the final tracks
850   fitter->WriteTracks(20);
851
852   delete fitter;
853
854   // remove the event number from the file names
855   char command[256];
856   sprintf(command, "rename tracks_%d tracks hlt/hough/*.raw", iEvent);
857   gSystem->Exec(command);
858   sprintf(command, "rename tracks_%d tracks hlt/fitter/*.raw", iEvent);
859   gSystem->Exec(command);
860   sprintf(command, "rename points_%d points hlt/fitter/*.raw", iEvent);
861   gSystem->Exec(command);
862   return kTRUE;
863 #endif
864 }
865
866 //_____________________________________________________________________________
867 Bool_t AliMonitorProcess::WriteHistos()
868 {
869 // write the monitor tree and the monitor histograms to the file 
870 // "monitor_<run number>[_<sub_run_number>].root"
871 // if at least fNEventsMin events were monitored
872
873   SetStatus(kWriting);
874
875   // rename tree file and create a new one
876   fFile->cd();
877   fTree->Write();
878   fFile->Close();
879   delete fFile;
880
881   char fileName[256];
882   sprintf(fileName, "monitor_tree_%d.root", fRunNumber);
883   if (fSubRunNumber > 0) {
884     sprintf(fileName, "monitor_tree_%d_%d.root", fRunNumber, fSubRunNumber);
885   }
886   if (fNEvents < fNEventsMin) {
887     gSystem->Unlink("monitor_tree.root");
888   } else {
889     gSystem->Rename("monitor_tree.root", fileName);
890   }
891
892   fFile = TFile::Open("monitor_tree.root", "RECREATE");
893   if (!fFile || !fFile->IsOpen()) {
894     Fatal("WriteHistos", "could not open file for tree");
895   }
896   fTree = new TTree("MonitorTree", "tree for monitoring");
897   for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
898     ((AliMonitor*) fMonitors[iMonitor])->CreateBranches(fTree);
899   }
900   gROOT->cd();
901
902   // write the histograms
903   if (fNEvents < fNEventsMin) return kTRUE;
904
905   if (!fWriteHistoList) {
906     TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
907     while (TFolder* folder = (TFolder*) iFolder->Next()) {
908       TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
909       while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
910         histo->ResetList();
911       }
912       delete iHisto;
913     }
914     delete iFolder;
915   }
916
917   Bool_t result = kTRUE;
918   sprintf(fileName, "monitor_%d.root", fRunNumber);
919   if (fSubRunNumber > 0) {
920     sprintf(fileName, "monitor_%d_%d.root", fRunNumber, fSubRunNumber);
921   }
922   TFile* file = TFile::Open(fileName, "recreate");
923   if (!file || !file->IsOpen()) {
924     Error("WriteHistos", "could not open file %s", fileName);
925     result = kFALSE;
926   } else {
927     fTopFolder->Write();
928     file->Close();
929   }
930   if (file) delete file;
931
932   return result;
933 }
934
935 //_____________________________________________________________________________
936 void AliMonitorProcess::StartNewRun()
937 {
938 // reset the histograms for a new run
939
940   SetStatus(kResetting);
941   TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
942   while (TFolder* folder = (TFolder*) iFolder->Next()) {
943     TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
944     while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
945       histo->Reset();
946     }
947     delete iHisto;
948   }
949   delete iFolder;
950
951   fNEvents = 0;
952 }
953
954
955 //_____________________________________________________________________________
956 void AliMonitorProcess::CheckForConnections()
957 {
958 // check if new clients want to connect and add them to the list of sockets
959
960   TSocket* socket;
961   while ((socket = fServerSocket->Accept()) != (TSocket*)-1) {
962     socket->SetOption(kNoBlock, 1);
963     char socketType[256];
964     if (socket->Recv(socketType, 255) <= 0) {
965       gSystem->Sleep(1000);
966       if (socket->Recv(socketType, 255) <= 0) {
967         TInetAddress adr = socket->GetInetAddress();
968         Error("CheckForConnections", "no socket type received - "
969               "disconnect client:\n %s (%s), port %d\n",
970               adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
971         delete socket;
972         continue;
973       }
974     }
975     if (strcmp(socketType, "client") == 0) {
976       fSockets.Add(socket);
977       TInetAddress adr = socket->GetInetAddress();
978       Info("CheckForConnections", "new client:\n %s (%s), port %d\n",
979            adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
980       if (fNEvents > 0) BroadcastHistos(socket);
981     } else if (strcmp(socketType, "display") == 0) {
982       if (fDisplaySocket) {
983         fDisplaySocket->Close();
984         delete fDisplaySocket;
985       }
986       fDisplaySocket = socket;
987       TInetAddress adr = socket->GetInetAddress();
988       Info("CheckForConnections", "new display:\n %s (%s), port %d\n",
989            adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
990     } else {
991       TInetAddress adr = socket->GetInetAddress();
992       Error("CheckForConnections", "unknown socket type %s - "
993             "disconnect client:\n %s (%s), port %d\n", socketType,
994             adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
995       delete socket;
996       continue;
997     }
998   }
999
1000   // remove finished or invalid clients
1001   for (Int_t iSocket = 0; iSocket < fSockets.GetEntriesFast(); iSocket++) {
1002     socket = (TSocket*) fSockets[iSocket];
1003     if (!socket) continue;
1004     char controlMessage[256];
1005     if (socket->Recv(controlMessage, 255)) {
1006       if (strcmp(controlMessage, "disconnect") == 0) {
1007         TInetAddress adr = socket->GetInetAddress();
1008         Info("CheckForConnections",
1009              "disconnect client:\n %s (%s), port %d\n",
1010              adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1011         delete fSockets.RemoveAt(iSocket);
1012         continue;
1013       }
1014     }
1015     if (!socket->IsValid()) {
1016       // remove invalid sockets from the list
1017       TInetAddress adr = socket->GetInetAddress();
1018       Error("CheckForConnections", 
1019             "disconnect invalid client:\n %s (%s), port %d\n",
1020             adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1021       delete fSockets.RemoveAt(iSocket);
1022     }
1023   }
1024   fSockets.Compress();
1025 }
1026
1027 //_____________________________________________________________________________
1028 void AliMonitorProcess::BroadcastHistos(TSocket* toSocket)
1029 {
1030 // send the monitor histograms to the clients
1031
1032   SetStatus(kBroadcasting);
1033   TMessage message(kMESS_OBJECT);
1034   message.WriteObject(fTopFolder); 
1035
1036   for (Int_t iSocket = 0; iSocket < fSockets.GetEntriesFast(); iSocket++) {
1037     TSocket* socket = (TSocket*) fSockets[iSocket];
1038     if (!socket) continue;
1039     if (toSocket && (socket != toSocket)) continue;
1040
1041     // send control message
1042     if (!socket->IsValid() || (socket->Send("histograms") <= 0)) {
1043       TInetAddress adr = socket->GetInetAddress();
1044       Error("BroadcastHistos", "connection to client failed - "
1045             "disconnect client:\n %s (%s), port %d\n",
1046             adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1047       delete fSockets.RemoveAt(iSocket);
1048     }
1049
1050     // receive control message
1051     char controlMessage[256];
1052     Int_t result = socket->Recv(controlMessage, 255);
1053     if (result <= 0) {
1054       gSystem->Sleep(1000);  // wait one second and try again
1055       result = socket->Recv(controlMessage, 255);
1056     }
1057     if (result <= 0) {
1058       TInetAddress adr = socket->GetInetAddress();
1059       Error("BroadcastHistos", "no response from client - "
1060             "disconnect client:\n %s (%s), port %d\n",
1061             adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1062       delete fSockets.RemoveAt(iSocket);
1063       continue;
1064     }
1065     if (strcmp(controlMessage, "ok") != 0) {
1066       TInetAddress adr = socket->GetInetAddress();
1067       Error("BroadcastHistos", "no \"ok\" message from client - "
1068             "disconnect client:\n %s (%s), port %d\n",
1069             adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1070       delete fSockets.RemoveAt(iSocket);
1071       continue;
1072     }
1073
1074     socket->SetOption(kNoBlock, 0);
1075     if (socket->Send(message) < 0) {
1076       // remove the socket from the list if there was an error
1077       TInetAddress adr = socket->GetInetAddress();
1078       Error("BroadcastHistos", "sending histograms failed - "
1079             "disconnect client:\n %s (%s), port %d\n",
1080             adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1081       delete fSockets.RemoveAt(iSocket);
1082     } else {
1083       gSystem->Sleep(100);
1084       socket->SetOption(kNoBlock, 1);
1085     }
1086   }
1087   fSockets.Compress();
1088 }
1089
1090
1091 //_____________________________________________________________________________
1092 AliMonitorProcess::AliMonitorInterruptHandler::AliMonitorInterruptHandler
1093   (AliMonitorProcess* process):
1094   TSignalHandler(kSigUser1, kFALSE), 
1095   fProcess(process) 
1096 {
1097 // constructor: set process
1098 }
1099
1100 //_____________________________________________________________________________
1101 AliMonitorProcess::AliMonitorInterruptHandler::AliMonitorInterruptHandler
1102   (const AliMonitorInterruptHandler& handler):
1103   TSignalHandler(handler)
1104 {
1105 // copy constructor
1106
1107   Fatal("AliMonitorInterruptHandler", "copy constructor not implemented");
1108 }
1109
1110 //_____________________________________________________________________________
1111 AliMonitorProcess::AliMonitorInterruptHandler& 
1112   AliMonitorProcess::AliMonitorInterruptHandler::operator = 
1113   (const AliMonitorInterruptHandler& /*handler*/) 
1114 {
1115 // assignment operator
1116
1117   Fatal("operator =", "assignment operator not implemented"); 
1118   return *this;
1119 }
1120
1121 //_____________________________________________________________________________
1122 Bool_t AliMonitorProcess::AliMonitorInterruptHandler::Notify() 
1123 {
1124 // interrupt signal -> stop process
1125
1126   Info("Notify", "the monitoring process will be stopped.");
1127   fProcess->Stop(); 
1128   return kTRUE;
1129 }