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 ///////////////////////////////////////////////////////////////////////////////
30 #include <TGridResult.h>
33 #include <TServerSocket.h>
39 #include "AliITSclustererV2.h"
40 #include "AliITStrackerV2.h"
41 #include "AliLoader.h"
42 #include "AliMonitorHLT.h"
43 #include "AliMonitorHLTHough.h"
44 #include "AliMonitorITS.h"
45 #include "AliMonitorProcess.h"
46 #include "AliMonitorTPC.h"
47 #include "AliMonitorV0s.h"
48 #include "AliRawReaderRoot.h"
51 #include "AliTPCclustererMI.h"
52 #include "AliTPCtrackerMI.h"
53 #include "AliV0vertexer.h"
55 #include <AliL3StandardIncludes.h>
56 #include <AliL3MemHandler.h>
57 #include <AliL3ClusterFitter.h>
58 #include <AliL3Fitter.h>
59 #include <AliL3Hough.h>
60 #include <AliL3HoughBaseTransformer.h>
61 #include <AliL3StandardIncludes.h>
62 #include <AliL3Track.h>
63 #include <AliL3TrackArray.h>
64 #include <AliL3Transform.h>
65 #include <AliL3Vertex.h>
66 #include <AliLevel3.h>
68 ClassImp(AliMonitorProcess)
71 const Int_t AliMonitorProcess::fgkPort = 9327;
74 //_____________________________________________________________________________
75 AliMonitorProcess::AliMonitorProcess(
76 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
77 const char* /*alienHost*/,
79 const char* alienHost,
82 const char* selection,
83 const char* fileNameGalice):
84 fSelection(selection),
99 fWriteHistoList(kFALSE),
108 fDisplaySocket(NULL),
113 fInterruptHandler(NULL)
115 // initialize the monitoring process and the monitor histograms
117 fSelection = selection;
119 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
120 fGrid = TGrid::Connect("alien", gSystem->Getenv("USER"));
122 fGrid = TGrid::Connect(alienHost, gSystem->Getenv("USER"));
124 if (!fGrid || fGrid->IsZombie() || !fGrid->IsConnected()) {
126 AliFatal("could not connect to alien");
128 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
132 fRunLoader = AliRunLoader::Open(fileNameGalice);
133 if (!fRunLoader) AliFatal(Form("could not get run loader from file %s",
136 fRunLoader->CdGAFile();
137 fTPCParam = AliTPC::LoadTPCParam(gFile);
138 if (!fTPCParam) AliFatal("could not load TPC parameters");
140 fRunLoader->LoadgAlice();
141 gAlice = fRunLoader->GetAliRun();
142 if (!gAlice) AliFatal("no gAlice object found");
143 AliITS* its = (AliITS*) gAlice->GetModule("ITS");
144 if (!its) AliFatal("no ITS detector found");
145 fITSgeom = its->GetITSgeom();
146 if (!fITSgeom) AliFatal("could not load ITS geometry");
148 // Init TPC parameters for HLT
149 Bool_t isinit=AliL3Transform::Init(const_cast<char*>(fileNameGalice),kTRUE);
151 AliFatal("Could not create transform settings, please check log for error messages!");
154 fTopFolder = new TFolder("Monitor", "monitor histograms");
155 fTopFolder->SetOwner(kTRUE);
157 if (IsSelected("TPC")) fMonitors.Add(new AliMonitorTPC(fTPCParam));
158 if (IsSelected("ITS")) fMonitors.Add(new AliMonitorITS(fITSgeom));
159 if (IsSelected("V0s")) fMonitors.Add(new AliMonitorV0s);
160 if (IsSelected("HLTConfMap")) fMonitors.Add(new AliMonitorHLT(fTPCParam));
161 if (IsSelected("HLTHough")) fMonitors.Add(new AliMonitorHLTHough(fTPCParam));
163 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
164 ((AliMonitor*) fMonitors[iMonitor])->CreateHistos(fTopFolder);
167 TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
168 while (TFolder* folder = (TFolder*) iFolder->Next()) folder->SetOwner(kTRUE);
171 fFile = TFile::Open("monitor_tree.root", "RECREATE");
172 if (!fFile || !fFile->IsOpen()) {
173 AliFatal("could not open file for tree");
175 fTree = new TTree("MonitorTree", "tree for monitoring");
176 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
177 ((AliMonitor*) fMonitors[iMonitor])->CreateBranches(fTree);
181 fServerSocket = new TServerSocket(fgkPort, kTRUE);
182 fServerSocket->SetOption(kNoBlock, 1);
183 CheckForConnections();
185 fInterruptHandler = new AliMonitorInterruptHandler(this);
186 gSystem->AddSignalHandler(fInterruptHandler);
189 //_____________________________________________________________________________
190 AliMonitorProcess::AliMonitorProcess(const AliMonitorProcess& process) :
199 fLogicalFileName(""),
208 fWriteHistoList(kFALSE),
217 fDisplaySocket(NULL),
222 fInterruptHandler(NULL)
225 AliFatal("copy constructor not implemented");
228 //_____________________________________________________________________________
229 AliMonitorProcess& AliMonitorProcess::operator = (const AliMonitorProcess&
232 AliFatal("assignment operator not implemented");
236 //_____________________________________________________________________________
237 AliMonitorProcess::~AliMonitorProcess()
245 delete fServerSocket;
247 delete fDisplaySocket;
249 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
256 gSystem->Unlink("monitor_tree.root");
261 gSystem->RemoveSignalHandler(fInterruptHandler);
262 delete fInterruptHandler;
266 //_____________________________________________________________________________
267 const char* AliMonitorProcess::GetRevision()
273 //_____________________________________________________________________________
274 void AliMonitorProcess::SetStatus(EStatus status)
276 // set the current status and process system events
279 gSystem->ProcessEvents();
283 //_____________________________________________________________________________
284 void AliMonitorProcess::Run()
286 // run the monitor process:
287 // check for a raw data file, process the raw data file and delete it
293 while (!CheckForNewFile()) {
294 CheckForConnections();
296 if (fStopping) break;
299 if (fStopping) break;
311 //_____________________________________________________________________________
312 void AliMonitorProcess::Stop()
314 // set the fStopping flag to terminate the monitor process after the current
315 // event was processed
317 if (GetStatus() != kStopped) fStopping = kTRUE;
321 //_____________________________________________________________________________
322 void AliMonitorProcess::ProcessFile(const char* fileName)
324 // create a file with monitor histograms for a single file
326 if (GetStatus() != kStopped) {
327 AliError("ProcessFile can not be called"
328 " while the monitor process is running");
332 fFileName = fileName;
333 Int_t nEventMin = fNEventsMin;
337 fNEventsMin = nEventMin;
342 //_____________________________________________________________________________
343 Bool_t AliMonitorProcess::CheckForNewFile()
345 // check whether a new file was registered in alien
347 #if ROOT_VERSION_CODE < ROOT_VERSION(5,0,0)
348 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
349 TGridResult* result = fGrid->Ls();
353 sprintf(dirName, "%s/adc-%d", fAlienDir.Data(), datime.GetDate());
355 sprintf(findName, "*.root");
356 Grid_ResultHandle_t handle = fGrid->Find(dirName, findName);
358 AliError(Form("could not open alien directory %s",
362 TGridResult* result = fGrid->CreateGridResult(handle);
368 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
369 while (const char* entry = result->Next()) {
371 while (Grid_Result_t* resultEntry = result->Next()) {
372 const char* entry = resultEntry->name.c_str();
374 if (strrchr(entry, '/')) entry = strrchr(entry, '/')+1;
375 // entry = host_date_time.root
376 TString entryCopy(entry);
377 char* p = const_cast<char*>(entryCopy.Data());
378 if (!strtok(p, "_") || !p) continue; // host name
379 char* dateStr = strtok(NULL, "_");
380 if (!dateStr || !p) continue;
381 char* timeStr = strtok(NULL, ".");
382 if (!timeStr || !p) continue;
383 Long_t date = atoi(dateStr);
384 Long_t time = atoi(timeStr);
386 if ((date > maxDate) || ((date == maxDate) && (time > maxTime))) {
394 if (maxDate < 0) return kFALSE; // no files found
395 if (fLogicalFileName.CompareTo(fileName) == 0) return kFALSE; // no new file
397 fLogicalFileName = fileName;
398 #if ROOT_VERSION_CODE <= 199169 // 3.10/01
399 result = fGrid->GetPhysicalFileNames(fLogicalFileName.Data());
400 fFileName = result->Next();
402 fileName = dirName + ("/" + fLogicalFileName);
403 handle = fGrid->GetPhysicalFileNames(fileName.Data());
405 AliError(Form("could not get physical file names for %s",
409 result = fGrid->CreateGridResult(handle);
411 Grid_Result_t* resultEntry = result->Next();
413 AliError(Form("could not get physical file names for %s",
417 fFileName = resultEntry->name2.c_str();
418 fFileName.ReplaceAll("castor:/", "rfio:/");
422 Error("CheckForNewFile", "needs to be ported to new TGrid");
427 //_____________________________________________________________________________
428 Bool_t AliMonitorProcess::ProcessFile()
430 // loop over all events in the raw data file, run the reconstruction
431 // and fill the monitor histograms
433 Int_t nEvents = GetNumberOfEvents(fFileName);
434 if (nEvents <= 0) return kFALSE;
435 AliDebug(1, Form("found %d event(s) in file %s",
436 nEvents, fFileName.Data()));
437 if (IsSelected("HLTConfMap")) CreateHLT(fFileName);
438 if (IsSelected("HLTHough")) CreateHLTHough(fFileName);
440 // loop over the events
441 for (Int_t iEvent = 0; iEvent < nEvents; iEvent++) {
442 CheckForConnections();
444 fRunLoader->SetEventNumber(0);
445 AliRawReaderRoot rawReader(fFileName, iEvent);
446 if (fStopping) break;
447 if (rawReader.GetRunNumber() != fRunNumber) {
450 fRunNumber = rawReader.GetRunNumber();
451 fEventNumber[0] = rawReader.GetEventId()[0];
452 fEventNumber[1] = rawReader.GetEventId()[1];
454 if (fStopping) break;
457 // monitor only central physics events
458 if (rawReader.GetType() != 7) continue;
459 if ((rawReader.GetAttributes()[0] & 0x02) == 0) continue;
460 AliInfo(Form("run: %d event: %d %d\n", rawReader.GetRunNumber(),
461 rawReader.GetEventId()[0], rawReader.GetEventId()[1]));
464 if (IsSelected("TPC")) {
465 CheckForConnections();
466 if (!ReconstructTPC(&rawReader, &esd)) return kFALSE;
467 if (fStopping) break;
469 if (IsSelected("ITS")) {
470 CheckForConnections();
471 if (!ReconstructITS(&rawReader, &esd)) return kFALSE;
472 if (fStopping) break;
474 if (IsSelected("V0s")) {
475 CheckForConnections();
476 if (!ReconstructV0s(&esd)) return kFALSE;
477 if (fStopping) break;
479 if (IsSelected("HLTConfMap")) {
480 CheckForConnections();
481 if (!ReconstructHLT(iEvent)) return kFALSE;
482 if (fStopping) break;
484 if (IsSelected("HLTHough")) {
485 CheckForConnections();
486 if (!ReconstructHLTHough(iEvent)) return kFALSE;
487 if (fStopping) break;
490 if (fDisplaySocket) fDisplaySocket->Send("new event");
492 AliDebug(1, "filling histograms...");
493 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
494 CheckForConnections();
496 ((AliMonitor*) fMonitors[iMonitor])->FillHistos(fRunLoader, &rawReader,
498 if (fStopping) break;
500 if (fStopping) break;
502 AliDebug(1, "updating histograms...");
503 CheckForConnections();
504 SetStatus(kUpdating);
505 TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
506 while (TFolder* folder = (TFolder*) iFolder->Next()) {
507 TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
508 while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
514 if (fStopping) break;
516 AliDebug(1, "filling the tree...");
519 AliDebug(1, "broadcasting histograms...");
520 CheckForConnections();
524 if (fStopping) break;
539 //_____________________________________________________________________________
540 void AliMonitorProcess::Reset()
542 // write the current histograms to a file and reset them
544 if (fSubRunNumber == 0) fSubRunNumber++;
551 //_____________________________________________________________________________
552 UInt_t AliMonitorProcess::GetEventPeriodNumber() const
554 // get the period number from the event id
556 return (fEventNumber[1] >> 4);
559 //_____________________________________________________________________________
560 UInt_t AliMonitorProcess::GetEventOrbitNumber() const
562 // get the orbit number from the event id
564 return ((fEventNumber[1] & 0x000F) << 20) + (fEventNumber[0] >> 12);
567 //_____________________________________________________________________________
568 UInt_t AliMonitorProcess::GetEventBunchNumber() const
570 // get the bunch number from the event id
572 return (fEventNumber[0] % 0x0FFF);
575 //_____________________________________________________________________________
576 Int_t AliMonitorProcess::GetNumberOfEvents(const char* fileName) const
578 // determine the number of events in the given raw data file
582 TFile* file = TFile::Open(fileName);
583 if (!file || !file->IsOpen()) {
584 AliError(Form("could not open file %s", fileName));
585 if (file) delete file;
589 TTree* tree = (TTree*) file->Get("RAW");
591 AliError("could not find tree with raw data");
593 nEvents = (Int_t) tree->GetEntries();
601 //_____________________________________________________________________________
602 Bool_t AliMonitorProcess::ReconstructTPC(AliRawReader* rawReader, AliESD* esd)
604 // find TPC clusters and tracks
608 AliLoader* tpcLoader = fRunLoader->GetLoader("TPCLoader");
610 AliError("no TPC loader found");
613 gSystem->Unlink("TPC.RecPoints.root");
616 AliDebug(1, "reconstructing clusters...");
617 tpcLoader->LoadRecPoints("recreate");
618 AliTPCclustererMI clusterer(fTPCParam);
619 tpcLoader->MakeRecPointsContainer();
620 clusterer.SetOutput(tpcLoader->TreeR());
621 clusterer.Digits2Clusters(rawReader);
622 tpcLoader->WriteRecPoints("OVERWRITE");
625 AliDebug(1, "reconstructing tracks...");
626 AliTPCtrackerMI tracker(fTPCParam);
627 tracker.LoadClusters(tpcLoader->TreeR());
628 tracker.Clusters2Tracks(esd);
629 tracker.UnloadClusters();
630 tpcLoader->UnloadRecPoints();
635 //_____________________________________________________________________________
636 Bool_t AliMonitorProcess::ReconstructITS(AliRawReader* rawReader, AliESD* esd)
638 // find ITS clusters and tracks
642 AliLoader* itsLoader = fRunLoader->GetLoader("ITSLoader");
644 AliError("no ITS loader found");
647 gSystem->Unlink("ITS.RecPoints.root");
650 AliDebug(1, "reconstructing clusters...");
651 itsLoader->LoadRecPoints("recreate");
652 AliITSclustererV2 clusterer(fITSgeom);
653 itsLoader->MakeRecPointsContainer();
654 clusterer.Digits2Clusters(rawReader);
657 AliDebug(1, "reconstructing tracks...");
658 AliITStrackerV2 tracker(fITSgeom);
659 tracker.LoadClusters(itsLoader->TreeR());
660 tracker.Clusters2Tracks(esd);
661 tracker.UnloadClusters();
663 itsLoader->UnloadRecPoints();
667 //_____________________________________________________________________________
668 Bool_t AliMonitorProcess::ReconstructV0s(AliESD* esd)
675 AliDebug(1, "reconstructing V0s...");
676 AliV0vertexer vertexer;
678 esd->GetVertex()->GetXYZ(vtx);
679 vertexer.SetVertex(vtx);
680 vertexer.Tracks2V0vertices(esd);
685 //_____________________________________________________________________________
686 void AliMonitorProcess::CreateHLT(const char* fileName)
689 // create the HLT (Level3) object
691 if (fHLT) delete fHLT;
694 strcpy(name, fileName);
695 fHLT = new AliLevel3(name);
696 fHLT->Init("./", AliLevel3::kRaw, 1);
698 fHLT->SetClusterFinderParam(-1, -1, kTRUE);
700 Int_t phiSegments = 50;
701 Int_t etaSegments = 100;
702 Int_t trackletlength = 3;
703 Int_t tracklength = 20;//40 or 5
704 Int_t rowscopetracklet = 2;
705 Int_t rowscopetrack = 10;
706 Double_t minPtFit = 0;
707 Double_t maxangle = 0.1745;
708 Double_t goodDist = 5;
709 Double_t maxphi = 0.1;
710 Double_t maxeta = 0.1;
711 Double_t hitChi2Cut = 15;//100 or 15
712 Double_t goodHitChi2 = 5;//20 or 5
713 Double_t trackChi2Cut = 10;//50 or 10
714 fHLT->SetTrackerParam(phiSegments, etaSegments,
715 trackletlength, tracklength,
716 rowscopetracklet, rowscopetrack,
717 minPtFit, maxangle, goodDist, hitChi2Cut,
718 goodHitChi2, trackChi2Cut, 50, maxphi, maxeta, kTRUE);
720 fHLT->WriteFiles("./hlt/");
723 //_____________________________________________________________________________
724 void AliMonitorProcess::CreateHLTHough(const char* fileName)
727 // create the HLT Hough transform (L3Hough) object
729 if (fHLTHough) delete fHLTHough;
732 strcpy(name, fileName);
734 fHLTHough = new AliL3Hough();
735 fHLTHough->SetThreshold(4);
736 fHLTHough->SetTransformerParams(140,150,0.5,-1);
737 fHLTHough->SetPeakThreshold(9000,-1);// or 6000
738 fHLTHough->Init("./", kFALSE, 50, kFALSE,0,name);
739 fHLTHough->SetAddHistograms();
740 // fHLTHough->GetMaxFinder()->SetThreshold(14000);
744 //_____________________________________________________________________________
745 Bool_t AliMonitorProcess::ReconstructHLT(Int_t iEvent)
747 // run the HLT cluster and track finder
751 gSystem->Exec("rm -rf hlt");
752 gSystem->MakeDirectory("hlt");
753 if (!fHLT) return kFALSE;
755 fHLT->ProcessEvent(0, 35, iEvent);
757 // remove the event number from the file names
759 sprintf(command, "rename points_%d points hlt/*.raw", iEvent);
760 gSystem->Exec(command);
761 sprintf(command, "rename tracks_tr_%d tracks_tr hlt/*.raw", iEvent);
762 gSystem->Exec(command);
763 sprintf(command, "rename tracks_gl_%d tracks_gl hlt/*.raw", iEvent);
764 gSystem->Exec(command);
765 sprintf(command, "rename tracks_%d tracks hlt/*.raw", iEvent);
766 gSystem->Exec(command);
770 //_____________________________________________________________________________
771 Bool_t AliMonitorProcess::ReconstructHLTHough(Int_t iEvent)
773 // run the HLT Hough transformer
777 gSystem->Exec("rm -rf hlt/hough");
778 gSystem->MakeDirectory("hlt/hough");
779 gSystem->Exec("rm -rf hlt/fitter");
780 gSystem->MakeDirectory("hlt/fitter");
781 if (!fHLTHough) return kFALSE;
783 // fHLTHough->Process(0, 35);
784 // Loop over TPC sectors and process the data
785 for(Int_t i=0; i<=35; i++)
787 fHLTHough->ReadData(i,iEvent);
788 fHLTHough->Transform();
789 // if(fHLTHough->fAddHistograms)
790 fHLTHough->AddAllHistograms();
791 fHLTHough->FindTrackCandidates();
792 fHLTHough->AddTracks();
794 fHLTHough->WriteTracks("./hlt/hough");
796 // Run cluster fitter
797 AliL3ClusterFitter *fitter = new AliL3ClusterFitter("./hlt");
799 // Set debug flag for the cluster fitter
802 // Setting fitter parameters
803 fitter->SetInnerWidthFactor(1,1.5);
804 fitter->SetOuterWidthFactor(1,1.5);
805 fitter->SetNmaxOverlaps(5);
807 //fitter->SetChiSqMax(5,kFALSE); //isolated clusters
808 fitter->SetChiSqMax(5,kTRUE); //overlapping clusters
810 Int_t rowrange[2] = {0,AliL3Transform::GetNRows()-1};
812 // Takes input from global hough tracks produced by HT
813 fitter->LoadSeeds(rowrange,kFALSE,iEvent);
817 for(Int_t islice = 0; islice <= 35; islice++)
819 for(Int_t ipatch = 0; ipatch < AliL3Transform::GetNPatches(); ipatch++)
822 fHLTHough->GetMemHandler(ipatch)->Free();
823 fHLTHough->GetMemHandler(ipatch)->Init(islice,ipatch);
824 AliL3DigitRowData *digits = (AliL3DigitRowData *)fHLTHough->GetMemHandler(ipatch)->AliAltroDigits2Memory(ndigits,iEvent);
826 fitter->Init(islice,ipatch);
827 fitter->SetInputData(digits);
828 fitter->FindClusters();
829 fitter->WriteClusters();
833 // Refit of the clusters
835 //The seeds are the input tracks from circle HT
836 AliL3TrackArray *tracks = fitter->GetSeeds();
837 AliL3Fitter *ft = new AliL3Fitter(&vertex,1);
839 ft->LoadClusters("./hlt/fitter/",iEvent,kFALSE);
840 for(Int_t i=0; i<tracks->GetNTracks(); i++)
842 AliL3Track *track = tracks->GetCheckedTrack(i);
844 if(track->GetNHits() < 20) continue;
845 ft->SortTrackClusters(track);
847 track->UpdateToFirstPoint();
851 //Write the final tracks
852 fitter->WriteTracks(20);
856 // remove the event number from the file names
858 sprintf(command, "rename tracks_%d tracks hlt/hough/*.raw", iEvent);
859 gSystem->Exec(command);
860 sprintf(command, "rename tracks_%d tracks hlt/fitter/*.raw", iEvent);
861 gSystem->Exec(command);
862 sprintf(command, "rename points_%d points hlt/fitter/*.raw", iEvent);
863 gSystem->Exec(command);
867 //_____________________________________________________________________________
868 Bool_t AliMonitorProcess::WriteHistos()
870 // write the monitor tree and the monitor histograms to the file
871 // "monitor_<run number>[_<sub_run_number>].root"
872 // if at least fNEventsMin events were monitored
876 // rename tree file and create a new one
883 sprintf(fileName, "monitor_tree_%d.root", fRunNumber);
884 if (fSubRunNumber > 0) {
885 sprintf(fileName, "monitor_tree_%d_%d.root", fRunNumber, fSubRunNumber);
887 if (fNEvents < fNEventsMin) {
888 gSystem->Unlink("monitor_tree.root");
890 gSystem->Rename("monitor_tree.root", fileName);
893 fFile = TFile::Open("monitor_tree.root", "RECREATE");
894 if (!fFile || !fFile->IsOpen()) {
895 AliFatal("could not open file for tree");
897 fTree = new TTree("MonitorTree", "tree for monitoring");
898 for (Int_t iMonitor = 0; iMonitor < fMonitors.GetEntriesFast(); iMonitor++) {
899 ((AliMonitor*) fMonitors[iMonitor])->CreateBranches(fTree);
903 // write the histograms
904 if (fNEvents < fNEventsMin) return kTRUE;
906 if (!fWriteHistoList) {
907 TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
908 while (TFolder* folder = (TFolder*) iFolder->Next()) {
909 TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
910 while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
918 Bool_t result = kTRUE;
919 sprintf(fileName, "monitor_%d.root", fRunNumber);
920 if (fSubRunNumber > 0) {
921 sprintf(fileName, "monitor_%d_%d.root", fRunNumber, fSubRunNumber);
923 TFile* file = TFile::Open(fileName, "recreate");
924 if (!file || !file->IsOpen()) {
925 AliError(Form("could not open file %s", fileName));
931 if (file) delete file;
936 //_____________________________________________________________________________
937 void AliMonitorProcess::StartNewRun()
939 // reset the histograms for a new run
941 SetStatus(kResetting);
942 TIterator* iFolder = fTopFolder->GetListOfFolders()->MakeIterator();
943 while (TFolder* folder = (TFolder*) iFolder->Next()) {
944 TIterator* iHisto = folder->GetListOfFolders()->MakeIterator();
945 while (AliMonitorPlot* histo = (AliMonitorPlot*) iHisto->Next()) {
956 //_____________________________________________________________________________
957 void AliMonitorProcess::CheckForConnections()
959 // check if new clients want to connect and add them to the list of sockets
962 while ((socket = fServerSocket->Accept()) != (TSocket*)-1) {
963 socket->SetOption(kNoBlock, 1);
964 char socketType[256];
965 if (socket->Recv(socketType, 255) <= 0) {
966 gSystem->Sleep(1000);
967 if (socket->Recv(socketType, 255) <= 0) {
968 TInetAddress adr = socket->GetInetAddress();
969 AliError(Form("no socket type received - "
970 "disconnect client: %s (%s), port %d",
971 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
976 if (strcmp(socketType, "client") == 0) {
977 fSockets.Add(socket);
978 TInetAddress adr = socket->GetInetAddress();
979 AliInfo(Form("new client: %s (%s), port %d",
980 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
981 if (fNEvents > 0) BroadcastHistos(socket);
982 } else if (strcmp(socketType, "display") == 0) {
983 if (fDisplaySocket) {
984 fDisplaySocket->Close();
985 delete fDisplaySocket;
987 fDisplaySocket = socket;
988 TInetAddress adr = socket->GetInetAddress();
989 AliInfo(Form("new display: %s (%s), port %d",
990 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
992 TInetAddress adr = socket->GetInetAddress();
993 AliError(Form("unknown socket type %s - "
994 "disconnect client: %s (%s), port %d", socketType,
995 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
1001 // remove finished or invalid clients
1002 for (Int_t iSocket = 0; iSocket < fSockets.GetEntriesFast(); iSocket++) {
1003 socket = (TSocket*) fSockets[iSocket];
1004 if (!socket) continue;
1005 char controlMessage[256];
1006 if (socket->Recv(controlMessage, 255)) {
1007 if (strcmp(controlMessage, "disconnect") == 0) {
1008 TInetAddress adr = socket->GetInetAddress();
1009 AliInfo(Form("disconnect client: %s (%s), port %d",
1010 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
1011 delete fSockets.RemoveAt(iSocket);
1015 if (!socket->IsValid()) {
1016 // remove invalid sockets from the list
1017 TInetAddress adr = socket->GetInetAddress();
1018 AliError(Form("disconnect invalid client: %s (%s), port %d",
1019 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
1020 delete fSockets.RemoveAt(iSocket);
1023 fSockets.Compress();
1026 //_____________________________________________________________________________
1027 void AliMonitorProcess::BroadcastHistos(TSocket* toSocket)
1029 // send the monitor histograms to the clients
1031 SetStatus(kBroadcasting);
1032 TMessage message(kMESS_OBJECT);
1033 message.WriteObject(fTopFolder);
1035 for (Int_t iSocket = 0; iSocket < fSockets.GetEntriesFast(); iSocket++) {
1036 TSocket* socket = (TSocket*) fSockets[iSocket];
1037 if (!socket) continue;
1038 if (toSocket && (socket != toSocket)) continue;
1040 // send control message
1041 if (!socket->IsValid() || (socket->Send("histograms") <= 0)) {
1042 TInetAddress adr = socket->GetInetAddress();
1043 AliError(Form("connection to client failed - "
1044 "disconnect client: %s (%s), port %d",
1045 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
1046 delete fSockets.RemoveAt(iSocket);
1049 // receive control message
1050 char controlMessage[256];
1051 Int_t result = socket->Recv(controlMessage, 255);
1053 gSystem->Sleep(1000); // wait one second and try again
1054 result = socket->Recv(controlMessage, 255);
1057 TInetAddress adr = socket->GetInetAddress();
1058 AliError(Form("no response from client - "
1059 "disconnect client: %s (%s), port %d",
1060 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
1061 delete fSockets.RemoveAt(iSocket);
1064 if (strcmp(controlMessage, "ok") != 0) {
1065 TInetAddress adr = socket->GetInetAddress();
1066 AliError(Form("no \"ok\" message from client - "
1067 "disconnect client: %s (%s), port %d",
1068 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
1069 delete fSockets.RemoveAt(iSocket);
1073 socket->SetOption(kNoBlock, 0);
1074 if (socket->Send(message) < 0) {
1075 // remove the socket from the list if there was an error
1076 TInetAddress adr = socket->GetInetAddress();
1077 AliError(Form("sending histograms failed - "
1078 "disconnect client: %s (%s), port %d",
1079 adr.GetHostName(), adr.GetHostAddress(), adr.GetPort()));
1080 delete fSockets.RemoveAt(iSocket);
1082 gSystem->Sleep(100);
1083 socket->SetOption(kNoBlock, 1);
1086 fSockets.Compress();
1090 //_____________________________________________________________________________
1091 AliMonitorProcess::AliMonitorInterruptHandler::AliMonitorInterruptHandler
1092 (AliMonitorProcess* process):
1093 TSignalHandler(kSigUser1, kFALSE),
1096 // constructor: set process
1099 //_____________________________________________________________________________
1100 AliMonitorProcess::AliMonitorInterruptHandler::AliMonitorInterruptHandler
1101 (const AliMonitorInterruptHandler& handler):
1102 TSignalHandler(handler)
1106 AliFatal("copy constructor not implemented");
1109 //_____________________________________________________________________________
1110 AliMonitorProcess::AliMonitorInterruptHandler&
1111 AliMonitorProcess::AliMonitorInterruptHandler::operator =
1112 (const AliMonitorInterruptHandler& /*handler*/)
1114 // assignment operator
1116 AliFatal("assignment operator not implemented");
1120 //_____________________________________________________________________________
1121 Bool_t AliMonitorProcess::AliMonitorInterruptHandler::Notify()
1123 // interrupt signal -> stop process
1125 AliInfo("the monitoring process will be stopped.");