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 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. //
26 ///////////////////////////////////////////////////////////////////////////////
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"
39 #include "AliTPCclustererMI.h"
40 #include "AliTPCtrackerMI.h"
42 #include "AliITSclustererV2.h"
43 #include "AliITStrackerV2.h"
44 #include "AliITSLoader.h"
45 #include "AliV0vertexer.h"
48 #include <TServerSocket.h>
50 #include <TGridResult.h>
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>
67 ClassImp(AliMonitorProcess)
70 const Int_t AliMonitorProcess::fgkPort = 9327;
73 //_____________________________________________________________________________
74 AliMonitorProcess::AliMonitorProcess(
75 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
76 const char* /*alienHost*/,
78 const char* alienHost,
81 const char* fileNameGalice)
83 // initialize the monitoring process and the monitor histograms
85 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
86 fGrid = TGrid::Connect("alien", gSystem->Getenv("USER"));
88 fGrid = TGrid::Connect(alienHost, gSystem->Getenv("USER"));
90 if (!fGrid || fGrid->IsZombie() || !fGrid->IsConnected()) {
92 Fatal("AliMonitorProcess", "could not connect to alien");
94 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
99 fLogicalFileName = "";
102 fRunLoader = AliRunLoader::Open(fileNameGalice);
103 if (!fRunLoader) Fatal("AliMonitorProcess",
104 "could not get run loader from file %s",
107 fRunLoader->CdGAFile();
108 fTPCParam = AliTPC::LoadTPCParam(gFile);
109 if (!fTPCParam) Fatal("AliMonitorProcess", "could not load TPC parameters");
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");
120 // Init TPC parameters for HLT
121 Bool_t isinit=AliL3Transform::Init(const_cast<char*>(fileNameGalice),kTRUE);
123 cerr << "Could not create transform settings, please check log for error messages!" << endl;
132 fWriteHistoList = kFALSE;
134 fTopFolder = new TFolder("Monitor", "monitor histograms");
135 fTopFolder->SetOwner(kTRUE);
137 fMonitors.Add(new AliMonitorTPC(fTPCParam));
138 fMonitors.Add(new AliMonitorITS(fITSgeom));
139 fMonitors.Add(new AliMonitorV0s);
141 fMonitors.Add(new AliMonitorHLT(fTPCParam));
142 fMonitors.Add(new AliMonitorHLTHough(fTPCParam));
145 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
146 ((AliMonitor*) fMonitors[iMonitor])->CreateHistos(fTopFolder);
149 TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
150 while (TFolder* folder = (TFolder*) iFolder->Next()) folder->SetOwner(kTRUE);
153 fFile = TFile::Open("monitor_tree.root", "RECREATE");
154 if (!fFile || !fFile->IsOpen()) {
155 Fatal("AliMonitorProcess", "could not open file for tree");
157 fTree = new TTree("MonitorTree", "tree for monitoring");
158 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
159 ((AliMonitor*) fMonitors[iMonitor])->CreateBranches(fTree);
163 fServerSocket = new TServerSocket(fgkPort, kTRUE);
164 fServerSocket->SetOption(kNoBlock, 1);
165 fDisplaySocket = NULL;
166 CheckForConnections();
174 fInterruptHandler = new AliMonitorInterruptHandler(this);
175 gSystem->AddSignalHandler(fInterruptHandler);
178 //_____________________________________________________________________________
179 AliMonitorProcess::AliMonitorProcess(const AliMonitorProcess& process) :
182 Fatal("AliMonitorProcess", "copy constructor not implemented");
185 //_____________________________________________________________________________
186 AliMonitorProcess& AliMonitorProcess::operator = (const AliMonitorProcess&
189 Fatal("operator =", "assignment operator not implemented");
193 //_____________________________________________________________________________
194 AliMonitorProcess::~AliMonitorProcess()
202 delete fServerSocket;
204 delete fDisplaySocket;
206 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
213 gSystem->Unlink("monitor_tree.root");
220 gSystem->RemoveSignalHandler(fInterruptHandler);
221 delete fInterruptHandler;
225 //_____________________________________________________________________________
226 const char* AliMonitorProcess::GetRevision()
232 //_____________________________________________________________________________
233 void AliMonitorProcess::SetStatus(EStatus status)
235 // set the current status and process system events
238 gSystem->ProcessEvents();
242 //_____________________________________________________________________________
243 void AliMonitorProcess::Run()
245 // run the monitor process:
246 // check for a raw data file, process the raw data file and delete it
252 while (!CheckForNewFile()) {
253 CheckForConnections();
255 if (fStopping) break;
258 if (fStopping) break;
270 //_____________________________________________________________________________
271 void AliMonitorProcess::Stop()
273 // set the fStopping flag to terminate the monitor process after the current
274 // event was processed
276 if (GetStatus() != kStopped) fStopping = kTRUE;
280 //_____________________________________________________________________________
281 void AliMonitorProcess::ProcessFile(const char* fileName)
283 // create a file with monitor histograms for a single file
285 if (GetStatus() != kStopped) {
286 Error("ProcessFile", "ProcessFile can not be called"
287 " while the monitor process is running");
291 fFileName = fileName;
292 Int_t nEventMin = fNEventsMin;
296 fNEventsMin = nEventMin;
301 //_____________________________________________________________________________
302 Bool_t AliMonitorProcess::CheckForNewFile()
304 // check whether a new file was registered in alien
306 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
307 TGridResult* result = fGrid->Ls();
309 // Grid_ResultHandle_t handle = fGrid->OpenDir(fAlienDir);
312 sprintf(findName, "*_%d_*.root", datime.GetDate());
313 Grid_ResultHandle_t handle = fGrid->Find(fAlienDir, findName);
315 Error("CheckForNewFile", "could not open alien directory %s",
319 TGridResult* result = fGrid->CreateGridResult(handle);
325 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
326 while (const char* entry = result->Next()) {
328 while (Grid_Result_t* resultEntry = result->Next()) {
329 const char* entry = resultEntry->name.c_str();
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);
343 if ((date > maxDate) || ((date == maxDate) && (time > maxTime))) {
351 if (maxDate < 0) return kFALSE; // no files found
352 if (fLogicalFileName.CompareTo(fileName) == 0) return kFALSE; // no new file
354 fLogicalFileName = fileName;
355 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
356 result = fGrid->GetPhysicalFileNames(fLogicalFileName.Data());
357 fFileName = result->Next();
359 fileName = fAlienDir + "/" + fLogicalFileName;
360 handle = fGrid->GetPhysicalFileNames(fileName.Data());
362 Error("CheckForNewFile", "could not get physical file names for %s",
366 result = fGrid->CreateGridResult(handle);
368 Grid_Result_t* resultEntry = result->Next();
370 Error("CheckForNewFile", "could not get physical file names for %s",
374 fFileName = resultEntry->name2.c_str();
375 fFileName.ReplaceAll("castor:/", "rfio:/");
382 //_____________________________________________________________________________
383 Bool_t AliMonitorProcess::ProcessFile()
385 // loop over all events in the raw data file, run the reconstruction
386 // and fill the monitor histograms
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());
393 CreateHLT(fFileName);
394 CreateHLTHough(fFileName);
397 // loop over the events
398 for (Int_t iEvent = 0; iEvent < nEvents; iEvent++) {
399 CheckForConnections();
401 fRunLoader->SetEventNumber(0);
402 AliRawReaderRoot rawReader(fFileName, iEvent);
403 if (fStopping) break;
404 if (rawReader.GetRunNumber() != fRunNumber) {
407 fRunNumber = rawReader.GetRunNumber();
408 fEventNumber[0] = rawReader.GetEventId()[0];
409 fEventNumber[1] = rawReader.GetEventId()[1];
411 if (fStopping) break;
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]);
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;
436 if (fDisplaySocket) fDisplaySocket->Send("new event");
438 Info("ProcessFile", "filling histograms...");
439 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
440 CheckForConnections();
442 ((AliMonitor*) fMonitors[iMonitor])->FillHistos(fRunLoader, &rawReader);
443 if (fStopping) break;
445 if (fStopping) break;
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()) {
459 if (fStopping) break;
461 Info("ProcessFile", "filling the tree...");
464 Info("ProcessFile", "broadcasting histograms...");
465 CheckForConnections();
469 if (fStopping) break;
480 //_____________________________________________________________________________
481 void AliMonitorProcess::Reset()
483 // write the current histograms to a file and reset them
485 if (fSubRunNumber == 0) fSubRunNumber++;
492 //_____________________________________________________________________________
493 UInt_t AliMonitorProcess::GetEventPeriodNumber() const
495 // get the period number from the event id
497 return (fEventNumber[1] >> 4);
500 //_____________________________________________________________________________
501 UInt_t AliMonitorProcess::GetEventOrbitNumber() const
503 // get the orbit number from the event id
505 return ((fEventNumber[1] & 0x000F) << 20) + (fEventNumber[0] >> 12);
508 //_____________________________________________________________________________
509 UInt_t AliMonitorProcess::GetEventBunchNumber() const
511 // get the bunch number from the event id
513 return (fEventNumber[0] % 0x0FFF);
516 //_____________________________________________________________________________
517 Int_t AliMonitorProcess::GetNumberOfEvents(const char* fileName) const
519 // determine the number of events in the given raw data file
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;
530 TTree* tree = (TTree*) file->Get("RAW");
532 Error("GetNumberOfEvents", "could not find tree with raw data");
534 nEvents = (Int_t) tree->GetEntries();
542 //_____________________________________________________________________________
543 Bool_t AliMonitorProcess::ReconstructTPC(AliRawReader* rawReader)
545 // find TPC clusters and tracks
549 AliLoader* tpcLoader = fRunLoader->GetLoader("TPCLoader");
551 Error("ReconstructTPC", "no TPC loader found");
554 gSystem->Unlink("TPC.RecPoints.root");
555 gSystem->Unlink("TPC.Tracks.root");
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");
567 Info("ReconstructTPC", "reconstructing tracks...");
568 tpcLoader->LoadTracks("recreate");
570 AliTPCtrackerMI tracker(fTPCParam);
572 tracker.LoadClusters();
573 tracker.Clusters2Tracks();
574 tracker.WriteTracks();
575 tracker.UnloadClusters();
578 tpcLoader->UnloadRecPoints();
579 tpcLoader->UnloadTracks();
583 //_____________________________________________________________________________
584 Bool_t AliMonitorProcess::ReconstructITS(AliRawReader* rawReader)
586 // find ITS clusters and tracks
590 AliLoader* itsLoader = fRunLoader->GetLoader("ITSLoader");
592 Error("ReconstructITS", "no ITS loader found");
595 AliLoader* tpcLoader = fRunLoader->GetLoader("TPCLoader");
597 Error("ReconstructITS", "no TPC loader found");
600 gSystem->Unlink("ITS.RecPoints.root");
601 gSystem->Unlink("ITS.Tracks.root");
604 Info("ReconstructITS", "reconstructing clusters...");
605 itsLoader->LoadRecPoints("recreate");
606 AliITSclustererV2 clusterer(fITSgeom);
607 itsLoader->MakeRecPointsContainer();
608 clusterer.Digits2Clusters(rawReader);
611 Info("ReconstructITS", "reconstructing tracks...");
612 itsLoader->LoadTracks("recreate");
613 itsLoader->MakeTracksContainer();
614 tpcLoader->LoadTracks();
615 AliITStrackerV2 tracker(fITSgeom);
616 tracker.LoadClusters(itsLoader->TreeR());
617 tracker.Clusters2Tracks(tpcLoader->TreeT(), itsLoader->TreeT());
618 tracker.UnloadClusters();
619 itsLoader->WriteTracks("OVERWRITE");
621 itsLoader->UnloadRecPoints();
622 itsLoader->UnloadTracks();
623 tpcLoader->UnloadTracks();
627 //_____________________________________________________________________________
628 Bool_t AliMonitorProcess::ReconstructV0s()
634 AliITSLoader* itsLoader = (AliITSLoader*) fRunLoader->GetLoader("ITSLoader");
636 Error("ReconstructV0", "no ITS loader found");
639 gSystem->Unlink("ITS.V0s.root");
642 Info("ReconstructV0s", "reconstructing V0s...");
643 itsLoader->LoadTracks("read");
644 itsLoader->LoadV0s("recreate");
645 AliV0vertexer vertexer;
646 TTree* tracks = itsLoader->TreeT();
648 Error("ReconstructV0s", "no ITS tracks tree found");
651 if (!itsLoader->TreeV0()) itsLoader->MakeTree("V0");
652 TTree* v0s = itsLoader->TreeV0();
653 vertexer.Tracks2V0vertices(tracks, v0s);
654 itsLoader->WriteV0s("OVERWRITE");
656 itsLoader->UnloadTracks();
657 itsLoader->UnloadV0s();
661 //_____________________________________________________________________________
663 void AliMonitorProcess::CreateHLT(const char* fileName)
666 // create the HLT (Level3) object
668 if (fHLT) delete fHLT;
671 strcpy(name, fileName);
672 fHLT = new AliLevel3(name);
673 fHLT->Init("./", AliLevel3::kRaw, 1);
675 fHLT->SetClusterFinderParam(-1, -1, kTRUE);
677 Int_t phiSegments = 50;
678 Int_t etaSegments = 100;
679 Int_t trackletlength = 3;
680 Int_t tracklength = 20;//40 or 5
681 Int_t rowscopetracklet = 2;
682 Int_t rowscopetrack = 10;
683 Double_t minPtFit = 0;
684 Double_t maxangle = 0.1745;
685 Double_t goodDist = 5;
686 Double_t maxphi = 0.1;
687 Double_t maxeta = 0.1;
688 Double_t hitChi2Cut = 15;//100 or 15
689 Double_t goodHitChi2 = 5;//20 or 5
690 Double_t trackChi2Cut = 10;//50 or 10
691 fHLT->SetTrackerParam(phiSegments, etaSegments,
692 trackletlength, tracklength,
693 rowscopetracklet, rowscopetrack,
694 minPtFit, maxangle, goodDist, hitChi2Cut,
695 goodHitChi2, trackChi2Cut, 50, maxphi, maxeta, kTRUE);
697 fHLT->WriteFiles("./hlt/");
700 //_____________________________________________________________________________
701 void AliMonitorProcess::CreateHLTHough(const char* fileName)
704 // create the HLT Hough transform (L3Hough) object
706 if (fHLTHough) delete fHLTHough;
709 strcpy(name, fileName);
711 fHLTHough = new AliL3Hough();
712 fHLTHough->SetThreshold(4);
713 fHLTHough->SetTransformerParams(140,150,0.5,-1);
714 fHLTHough->SetPeakThreshold(9000,-1);// or 6000
715 fHLTHough->Init("./", kFALSE, 50, kFALSE,0,name);
716 fHLTHough->SetAddHistograms();
717 // fHLTHough->GetMaxFinder()->SetThreshold(14000);
722 //_____________________________________________________________________________
723 Bool_t AliMonitorProcess::ReconstructHLT(
731 // run the HLT cluster and track finder
736 Warning("ReconstructHLT", "the code was compiled without HLT support");
740 gSystem->Exec("rm -rf hlt");
741 gSystem->MakeDirectory("hlt");
742 if (!fHLT) return kFALSE;
744 fHLT->ProcessEvent(0, 35, iEvent);
746 // remove the event number from the file names
748 sprintf(command, "rename points_%d points hlt/*.raw", iEvent);
749 gSystem->Exec(command);
750 sprintf(command, "rename tracks_tr_%d tracks_tr hlt/*.raw", iEvent);
751 gSystem->Exec(command);
752 sprintf(command, "rename tracks_gl_%d tracks_gl hlt/*.raw", iEvent);
753 gSystem->Exec(command);
754 sprintf(command, "rename tracks_%d tracks hlt/*.raw", iEvent);
755 gSystem->Exec(command);
760 //_____________________________________________________________________________
761 Bool_t AliMonitorProcess::ReconstructHLTHough(
769 // run the HLT Hough transformer
774 Warning("ReconstructHLTHough", "the code was compiled without HLT support");
778 gSystem->Exec("rm -rf hlt/hough");
779 gSystem->MakeDirectory("hlt/hough");
780 gSystem->Exec("rm -rf hlt/fitter");
781 gSystem->MakeDirectory("hlt/fitter");
782 if (!fHLTHough) return kFALSE;
784 // fHLTHough->Process(0, 35);
785 // Loop over TPC sectors and process the data
786 for(Int_t i=0; i<=35; i++)
788 fHLTHough->ReadData(i,iEvent);
789 fHLTHough->Transform();
790 // if(fHLTHough->fAddHistograms)
791 fHLTHough->AddAllHistograms();
792 fHLTHough->FindTrackCandidates();
793 fHLTHough->AddTracks();
795 fHLTHough->WriteTracks("./hlt/hough");
797 // Run cluster fitter
798 AliL3ClusterFitter *fitter = new AliL3ClusterFitter("./hlt");
800 // Set debug flag for the cluster fitter
803 // Setting fitter parameters
804 fitter->SetInnerWidthFactor(1,1.5);
805 fitter->SetOuterWidthFactor(1,1.5);
806 fitter->SetNmaxOverlaps(5);
808 //fitter->SetChiSqMax(5,kFALSE); //isolated clusters
809 fitter->SetChiSqMax(5,kTRUE); //overlapping clusters
811 Int_t rowrange[2] = {0,AliL3Transform::GetNRows()-1};
813 // Takes input from global hough tracks produced by HT
814 fitter->LoadSeeds(rowrange,kFALSE,iEvent);
818 for(Int_t islice = 0; islice <= 35; islice++)
820 for(Int_t ipatch = 0; ipatch < AliL3Transform::GetNPatches(); ipatch++)
823 fHLTHough->GetMemHandler(ipatch)->Free();
824 fHLTHough->GetMemHandler(ipatch)->Init(islice,ipatch);
825 AliL3DigitRowData *digits = (AliL3DigitRowData *)fHLTHough->GetMemHandler(ipatch)->AliAltroDigits2Memory(ndigits,iEvent);
827 fitter->Init(islice,ipatch);
828 fitter->SetInputData(digits);
829 fitter->FindClusters();
830 fitter->WriteClusters();
834 // Refit of the clusters
836 //The seeds are the input tracks from circle HT
837 AliL3TrackArray *tracks = fitter->GetSeeds();
838 AliL3Fitter *ft = new AliL3Fitter(&vertex,1);
840 ft->LoadClusters("./hlt/fitter/",iEvent,kFALSE);
841 for(Int_t i=0; i<tracks->GetNTracks(); i++)
843 AliL3Track *track = tracks->GetCheckedTrack(i);
845 if(track->GetNHits() < 20) continue;
846 ft->SortTrackClusters(track);
848 ft->UpdateTrack(track);
852 //Write the final tracks
853 fitter->WriteTracks(20);
857 // remove the event number from the file names
859 sprintf(command, "rename tracks_%d tracks hlt/hough/*.raw", iEvent);
860 gSystem->Exec(command);
861 sprintf(command, "rename tracks_%d tracks hlt/fitter/*.raw", iEvent);
862 gSystem->Exec(command);
863 sprintf(command, "rename points_%d points hlt/fitter/*.raw", iEvent);
864 gSystem->Exec(command);
869 //_____________________________________________________________________________
870 Bool_t AliMonitorProcess::WriteHistos()
872 // write the monitor tree and the monitor histograms to the file
873 // "monitor_<run number>[_<sub_run_number>].root"
874 // if at least fNEventsMin events were monitored
878 // rename tree file and create a new one
885 sprintf(fileName, "monitor_tree_%d.root", fRunNumber);
886 if (fSubRunNumber > 0) {
887 sprintf(fileName, "monitor_tree_%d_%d.root", fRunNumber, fSubRunNumber);
889 if (fNEvents < fNEventsMin) {
890 gSystem->Unlink("monitor_tree.root");
892 gSystem->Rename("monitor_tree.root", fileName);
895 fFile = TFile::Open("monitor_tree.root", "RECREATE");
896 if (!fFile || !fFile->IsOpen()) {
897 Fatal("WriteHistos", "could not open file for tree");
899 fTree = new TTree("MonitorTree", "tree for monitoring");
900 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
901 ((AliMonitor*) fMonitors[iMonitor])->CreateBranches(fTree);
905 // write the histograms
906 if (fNEvents < fNEventsMin) return kTRUE;
908 if (!fWriteHistoList) {
909 TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
910 while (TFolder* folder = (TFolder*) iFolder->Next()) {
911 TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
912 while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
920 Bool_t result = kTRUE;
921 sprintf(fileName, "monitor_%d.root", fRunNumber);
922 if (fSubRunNumber > 0) {
923 sprintf(fileName, "monitor_%d_%d.root", fRunNumber, fSubRunNumber);
925 TFile* file = TFile::Open(fileName, "recreate");
926 if (!file || !file->IsOpen()) {
927 Error("WriteHistos", "could not open file %s", fileName);
933 if (file) delete file;
938 //_____________________________________________________________________________
939 void AliMonitorProcess::StartNewRun()
941 // reset the histograms for a new run
943 SetStatus(kResetting);
944 TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
945 while (TFolder* folder = (TFolder*) iFolder->Next()) {
946 TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
947 while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
958 //_____________________________________________________________________________
959 void AliMonitorProcess::CheckForConnections()
961 // check if new clients want to connect and add them to the list of sockets
964 while ((socket = fServerSocket->Accept()) != (TSocket*)-1) {
965 socket->SetOption(kNoBlock, 1);
966 char socketType[256];
967 if (socket->Recv(socketType, 255) <= 0) {
968 gSystem->Sleep(1000);
969 if (socket->Recv(socketType, 255) <= 0) {
970 TInetAddress adr = socket->GetInetAddress();
971 Error("CheckForConnections", "no socket type received - "
972 "disconnect client:\n %s (%s), port %d\n",
973 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
978 if (strcmp(socketType, "client") == 0) {
979 fSockets.Add(socket);
980 TInetAddress adr = socket->GetInetAddress();
981 Info("CheckForConnections", "new client:\n %s (%s), port %d\n",
982 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
983 if (fNEvents > 0) BroadcastHistos(socket);
984 } else if (strcmp(socketType, "display") == 0) {
985 if (fDisplaySocket) {
986 fDisplaySocket->Close();
987 delete fDisplaySocket;
989 fDisplaySocket = socket;
990 TInetAddress adr = socket->GetInetAddress();
991 Info("CheckForConnections", "new display:\n %s (%s), port %d\n",
992 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
994 TInetAddress adr = socket->GetInetAddress();
995 Error("CheckForConnections", "unknown socket type %s - "
996 "disconnect client:\n %s (%s), port %d\n", socketType,
997 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1003 // remove finished or invalid clients
1004 for (Int_t iSocket = 0; iSocket < fSockets.GetEntriesFast(); iSocket++) {
1005 socket = (TSocket*) fSockets[iSocket];
1006 if (!socket) continue;
1007 char controlMessage[256];
1008 if (socket->Recv(controlMessage, 255)) {
1009 if (strcmp(controlMessage, "disconnect") == 0) {
1010 TInetAddress adr = socket->GetInetAddress();
1011 Info("CheckForConnections",
1012 "disconnect client:\n %s (%s), port %d\n",
1013 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1014 delete fSockets.RemoveAt(iSocket);
1018 if (!socket->IsValid()) {
1019 // remove invalid sockets from the list
1020 TInetAddress adr = socket->GetInetAddress();
1021 Error("CheckForConnections",
1022 "disconnect invalid client:\n %s (%s), port %d\n",
1023 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1024 delete fSockets.RemoveAt(iSocket);
1027 fSockets.Compress();
1030 //_____________________________________________________________________________
1031 void AliMonitorProcess::BroadcastHistos(TSocket* toSocket)
1033 // send the monitor histograms to the clients
1035 SetStatus(kBroadcasting);
1036 TMessage message(kMESS_OBJECT);
1037 message.WriteObject(fTopFolder);
1039 for (Int_t iSocket = 0; iSocket < fSockets.GetEntriesFast(); iSocket++) {
1040 TSocket* socket = (TSocket*) fSockets[iSocket];
1041 if (!socket) continue;
1042 if (toSocket && (socket != toSocket)) continue;
1044 // send control message
1045 if (!socket->IsValid() || (socket->Send("histograms") <= 0)) {
1046 TInetAddress adr = socket->GetInetAddress();
1047 Error("BroadcastHistos", "connection to client failed - "
1048 "disconnect client:\n %s (%s), port %d\n",
1049 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1050 delete fSockets.RemoveAt(iSocket);
1053 // receive control message
1054 char controlMessage[256];
1055 Int_t result = socket->Recv(controlMessage, 255);
1057 gSystem->Sleep(1000); // wait one second and try again
1058 result = socket->Recv(controlMessage, 255);
1061 TInetAddress adr = socket->GetInetAddress();
1062 Error("BroadcastHistos", "no response from client - "
1063 "disconnect client:\n %s (%s), port %d\n",
1064 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1065 delete fSockets.RemoveAt(iSocket);
1068 if (strcmp(controlMessage, "ok") != 0) {
1069 TInetAddress adr = socket->GetInetAddress();
1070 Error("BroadcastHistos", "no \"ok\" message from client - "
1071 "disconnect client:\n %s (%s), port %d\n",
1072 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1073 delete fSockets.RemoveAt(iSocket);
1077 socket->SetOption(kNoBlock, 0);
1078 if (socket->Send(message) < 0) {
1079 // remove the socket from the list if there was an error
1080 TInetAddress adr = socket->GetInetAddress();
1081 Error("BroadcastHistos", "sending histograms failed - "
1082 "disconnect client:\n %s (%s), port %d\n",
1083 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort());
1084 delete fSockets.RemoveAt(iSocket);
1086 gSystem->Sleep(100);
1087 socket->SetOption(kNoBlock, 1);
1090 fSockets.Compress();
1094 //_____________________________________________________________________________
1095 AliMonitorProcess::AliMonitorInterruptHandler::AliMonitorInterruptHandler
1096 (AliMonitorProcess* process):
1097 TSignalHandler(kSigUser1, kFALSE),
1100 // constructor: set process
1103 //_____________________________________________________________________________
1104 AliMonitorProcess::AliMonitorInterruptHandler::AliMonitorInterruptHandler
1105 (const AliMonitorInterruptHandler& handler):
1106 TSignalHandler(handler)
1110 Fatal("AliMonitorInterruptHandler", "copy constructor not implemented");
1113 //_____________________________________________________________________________
1114 AliMonitorProcess::AliMonitorInterruptHandler&
1115 AliMonitorProcess::AliMonitorInterruptHandler::operator =
1116 (const AliMonitorInterruptHandler& /*handler*/)
1118 // assignment operator
1120 Fatal("operator =", "assignment operator not implemented");
1124 //_____________________________________________________________________________
1125 Bool_t AliMonitorProcess::AliMonitorInterruptHandler::Notify()
1127 // interrupt signal -> stop process
1129 Info("Notify", "the monitoring process will be stopped.");