1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // this program performs local monitoring on a GDC by running the HLT code //
22 // If an argument is given, this is taken as the name of a date file which //
23 // is used instead of the local node. //
24 // The program can be stopped by pressing CTRL-C. //
26 ///////////////////////////////////////////////////////////////////////////////
29 #include <TSysEvtHandler.h>
34 #include "AliRawReaderDate.h"
37 #include <AliL3StandardIncludes.h>
38 #include <AliL3Transform.h>
39 #include <AliL3MemHandler.h>
40 #include <AliL3TrackArray.h>
41 #include <AliL3HoughMaxFinder.h>
42 #include <AliL3HoughBaseTransformer.h>
43 #include <AliL3Hough.h>
44 #include <AliL3Benchmark.h>
45 #include <AliKalmanTrack.h>
46 #include "AliITSgeom.h"
48 #include "AliMagFMaps.h"
49 #include <AliL3ITSclusterer.h>
50 #include <AliL3ITSVertexerZ.h>
51 #include <AliL3ITStracker.h>
55 // avoid loading pythia and pdf
60 //_____________________________________________________________________________
61 class AliGDCInterruptHandler : public TSignalHandler {
63 AliGDCInterruptHandler();
64 Bool_t Notify() {fStop = kTRUE; return kTRUE;};
65 Bool_t Stop() const {return fStop;};
67 Bool_t fStop; // CTRL-C pressed
70 //_____________________________________________________________________________
71 AliGDCInterruptHandler::AliGDCInterruptHandler() :
72 TSignalHandler(kSigInterrupt, kFALSE)
78 //_____________________________________________________________________________
80 int main(int argc, char** argv)
82 // set ROOT in batch mode
86 FILE* file = fopen("monitorGDC.log", "w");
90 // get data from a file or online from this node
93 status = monitorSetDataSource(argv[1]);
95 status = monitorSetDataSource(":");
97 if (status) ::Fatal("monitorSetDataSource", monitorDecodeError(status));
99 // monitor only a sample of physics events
100 char* table[] = {"Physics event", "yes", NULL};
101 status = monitorDeclareTable(table);
102 if (status) ::Fatal("monitorDeclareTable", monitorDecodeError(status));
104 // declare this monitoring program to DATE
105 status = monitorDeclareMp("GDC physics monitoring");
106 if (status) ::Fatal("monitorDeclareMp", monitorDecodeError(status));
108 // initialize HLT transformations
109 if (!AliL3Transform::Init("./", kFALSE)) {
110 ::Fatal("AliL3Transform::Init", "HLT initialization failed");
112 AliESD *esd = new AliESD;
113 AliKalmanTrack::SetConvConst(
114 1000/0.299792458/AliL3Transform::GetSolenoidField()
116 AliITSgeom *geom = new AliITSgeom();
117 geom->ReadNewFile("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymmFMD.det");
119 AliMagF* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k5kG);
120 AliTracker::SetFieldMap(field);
122 // create the signal handler
123 AliGDCInterruptHandler* handler = new AliGDCInterruptHandler;
124 gSystem->AddSignalHandler(handler);
128 while (!handler->Stop()) {
129 // get the next event
130 status = monitorGetEventDynamic(&ptr);
131 if (status == (Int_t)MON_ERR_EOF) break;
132 if (status) ::Fatal("monitorGetEventDynamic", monitorDecodeError(status));
136 gSystem->Sleep(1000); // sleep for 1 second
140 AliRawReaderDate rawReader(ptr);
141 // if ((rawReader.GetAttributes()[0] & 0x02) != 0) {
143 // Int_t errorCode = rawReader.CheckData();
145 if (errorCode && (errorCode != AliRawReader::kErrSize)) {
147 if (file) fprintf(file, "%s\n", time.AsString());
148 if (file) fprintf(file, "run: %d event: %d %d\n",
149 rawReader.GetRunNumber(),
150 rawReader.GetEventId()[0],
151 rawReader.GetEventId()[1]);
152 fprintf(file, "ERROR: %d\n\n", errorCode);
156 AliL3Benchmark *fBenchmark = new AliL3Benchmark();
157 fBenchmark->Start("Overall timing");
159 // ITS clusterer and vertexer
160 fBenchmark->Start("ITS Clusterer");
161 AliL3ITSclusterer clusterer(geom);
162 AliRawReader *itsrawreader=new AliRawReaderDate(ptr);
163 TTree* treeClusters = new TTree("TreeL3ITSclusters"," "); //make a tree
164 clusterer.Digits2Clusters(itsrawreader,treeClusters);
165 fBenchmark->Stop("ITS Clusterer");
167 AliL3ITSVertexerZ vertexer;
168 AliESDVertex *vertex = vertexer.FindVertexForCurrentEvent(geom,treeClusters);
170 Double_t vtxErr[3]={0.005,0.005,0.010};
171 vertex->GetXYZ(vtxPos);
172 // vertex->GetSigmaXYZ(vtxErr);
173 esd->SetVertex(vertex);
175 // TPC Hough reconstruction
176 Float_t ptmin = 0.1*AliL3Transform::GetSolenoidField();
177 Float_t zvertex = vtxPos[2];
179 // Run the Hough Transformer
180 fBenchmark->Start("Init");
181 AliL3Hough *hough1 = new AliL3Hough();
183 hough1->SetThreshold(4);
184 hough1->CalcTransformerParams(ptmin);
185 hough1->SetPeakThreshold(70,-1);
186 // printf("Pointer is %x\n",ptr);
187 hough1->Init("./", kFALSE, 100, kFALSE,4,0,(Char_t*)ptr,zvertex);
188 hough1->SetAddHistograms();
189 fBenchmark->Stop("Init");
191 fBenchmark->Start("Init");
192 AliL3Hough *hough2 = new AliL3Hough();
194 hough2->SetThreshold(4);
195 hough2->CalcTransformerParams(ptmin);
196 hough2->SetPeakThreshold(70,-1);
197 // printf("Pointer is %x\n",ptr);
198 hough2->Init("./", kFALSE, 100, kFALSE,4,0,(Char_t*)ptr,zvertex);
199 hough2->SetAddHistograms();
200 fBenchmark->Stop("Init");
202 Int_t nglobaltracks = 0;
204 hough1->StartProcessInThread(0,17);
205 hough2->StartProcessInThread(18,35);
207 // gSystem->Sleep(20000);
208 if(hough1->WaitForThreadFinish())
209 ::Fatal("AliL3Hough::WaitForThreadFinish"," Can not join the required thread! ");
210 if(hough2->WaitForThreadFinish())
211 ::Fatal("AliL3Hough::WaitForThreadFinish"," Can not join the required thread! ");
213 gSystem->MakeDirectory("hough1");
214 hough1->WriteTracks("./hough1");
215 gSystem->MakeDirectory("hough2");
216 hough2->WriteTracks("./hough2");
219 for(int slice=0; slice<=17; slice++)
221 // cout<<"Processing slice "<<slice<<endl;
222 fBenchmark->Start("ReadData");
223 hough1->ReadData(slice,0);
224 fBenchmark->Stop("ReadData");
225 fBenchmark->Start("Transform");
227 fBenchmark->Stop("Transform");
228 hough1->AddAllHistogramsRows();
229 hough1->FindTrackCandidatesRow();
230 fBenchmark->Start("AddTracks");
232 fBenchmark->Stop("AddTracks");
234 // AliL3TrackArray* tracks = (AliL3TrackArray*)hough1->GetTracks(0);
235 // nglobaltracks += tracks->GetNTracks();
237 for(int slice=18; slice<=35; slice++)
239 // cout<<"Processing slice "<<slice<<endl;
240 fBenchmark->Start("ReadData");
241 hough2->ReadData(slice,0);
242 fBenchmark->Stop("ReadData");
243 fBenchmark->Start("Transform");
245 fBenchmark->Stop("Transform");
246 hough2->AddAllHistogramsRows();
247 hough2->FindTrackCandidatesRow();
248 fBenchmark->Start("AddTracks");
250 fBenchmark->Stop("AddTracks");
252 // AliL3TrackArray* tracks = (AliL3TrackArray*)hough2->GetTracks(0);
253 // nglobaltracks += tracks->GetNTracks();
256 nglobaltracks += hough1->FillESD(esd);
257 nglobaltracks += hough2->FillESD(esd);
260 AliL3ITStracker itsTracker(geom);
261 itsTracker.SetVertex(vtxPos,vtxErr);
263 itsTracker.LoadClusters(treeClusters);
264 itsTracker.Clusters2Tracks(esd);
265 itsTracker.UnloadClusters();
267 fBenchmark->Stop("Overall timing");
269 if (file) fprintf(file, "%s\n", time.AsString());
270 if (file) fprintf(file, "run: %d event: %d %d\n",
271 rawReader.GetRunNumber(),
272 rawReader.GetEventId()[0],
273 rawReader.GetEventId()[1]);
274 if (errorCode) fprintf(file, "ERROR: %d\n", errorCode);
276 if (file) fprintf(file, "Hough Transformer found %d tracks\n",nglobaltracks);
278 hough1->DoBench("hough1");
279 hough2->DoBench("hough2");
280 fBenchmark->Analyze("overall");
282 FILE* bench = fopen("hough1.dat", "r");
283 while (bench && !feof(bench)) {
285 if (!fgets(buffer, 256, bench)) break;
286 fprintf(file, "%s", buffer);
291 FILE* bench = fopen("hough2.dat", "r");
292 while (bench && !feof(bench)) {
294 if (!fgets(buffer, 256, bench)) break;
295 fprintf(file, "%s", buffer);
300 FILE* bench = fopen("overall.dat", "r");
301 while (bench && !feof(bench)) {
303 if (!fgets(buffer, 256, bench)) break;
304 fprintf(file, "%s", buffer);
307 fprintf(file, "\n\n");
318 // read run, event, detector, DDL numbers and data size
319 AliRawReaderDate rawReader(ptr);
321 printf("\n%s\n", time.AsString());
322 printf("run: %d event: %d %d\n", rawReader.GetRunNumber(),
323 rawReader.GetEventId()[0], rawReader.GetEventId()[1]);
324 while (rawReader.ReadMiniHeader()) {
325 printf(" detector: %d DDL: %d size: %d\n",
326 rawReader.GetDetectorID(), rawReader.GetDDLID(),
327 rawReader.GetDataSize());
332 gSystem->Sleep(100); // sleep for 0.1 second
335 gSystem->ProcessEvents();
336 if (file) fflush(file);
340 gSystem->RemoveSignalHandler(handler);
341 if (file) fclose(file);
347 int main(int /*argc*/, char** /*argv*/)
349 ::Fatal("main", "this program was compiled without DATE");