Adding configuration file for TOF da running during Physics runs.
[u/mrichter/AliRoot.git] / TOF / TOFda.cxx
1 /*
2
3 TOF DA for online calibration
4
5 Contact: Chiara.Zampolli@bo.infn.it
6 Link: www.bo.infn.it/~zampolli
7 Run Type: PHYSICS
8 DA Type: MON
9 Number of events needed: depending on the run, being run-level
10 Input Files: TOFdaTotal.root, to be updated if existing
11 Output Files: TOFdaRun.root, TOFdaTotal.root, both to be exported to the DAQ FXS
12 Trigger types used: PHYSICS_EVENT
13
14 */
15
16 #define FILE_TOTAL "TOFdaTotal.root"
17 #define FILE_RUN "TOFdaRun.root"
18
19 // DATE
20 #include "event.h"
21 #include "monitor.h"
22 #include "daqDA.h"
23
24 #include <stdio.h>
25 #include <stdlib.h>
26
27 //AliRoot
28 #include <AliTOFRawStream.h>
29 #include <AliRawReaderDate.h>
30 #include <AliRawReader.h>
31 #include <AliTOFGeometry.h>
32 #include <AliT0RawReader.h>
33 #include <AliDAQ.h>
34 #include <AliTOFHitData.h>
35 #include <AliTOFHitDataBuffer.h>
36 #include <AliTOFDaConfigHandler.h>
37
38 //ROOT
39 #include <TFile.h>
40 #include <TKey.h>
41 #include <TH2S.h>
42 #include <TObject.h>
43 #include <TMath.h>
44 #include <TSystem.h>
45 #include "TROOT.h"
46 #include "TPluginManager.h"
47 #include "TSAXParser.h"
48
49 /* Main routine
50       Arguments: 
51       1- monitoring data source
52 */
53 int main(int argc, char **argv) {
54
55   /* magic line from Rene */
56   gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
57                     "*",
58                     "TStreamerInfo",
59                     "RIO",
60                     "TStreamerInfo()");
61
62   AliTOFGeometry * geom = new AliTOFGeometry();
63
64   static const Int_t size = AliTOFGeometry::NPadXSector()*AliTOFGeometry::NSectors();
65   static const Int_t nbins = 500;
66   static const Int_t binmin = -20;
67   const Float_t c = 2.99792458E10; //speed of light [cm/s]
68   TH1F::AddDirectory(0);
69   TH2S * htofPartial = new TH2S("htof","histo with delays",
70                                 size,-0.5,size*1.-0.5,
71                                 nbins,binmin-0.5,nbins*1.+binmin-0.5);
72   
73   int status;
74
75   /* log start of process */
76   printf("TOF DA started\n");  
77
78   /* check that we got some arguments = list of files */
79   if (argc!=2) {
80     printf("Wrong number of arguments\n");
81     return -1;
82   }
83
84   /* retrieve config file */
85   int getConfigFile = daqDA_DB_getFile("TOFPhysicsConfig.xml","TOFPhysicsConfig.xml");
86   if (getConfigFile != 0){
87     printf("Failed to retrieve config file from DB! returning...\n");
88     return -1;
89   }
90
91   AliTOFDaConfigHandler* tofHandler = new AliTOFDaConfigHandler();
92   TSAXParser *parser = new TSAXParser();
93   parser->ConnectToHandler("AliTOFDaConfigHandler", tofHandler);  
94   if (parser->ParseFile("./TOFPhysicsConfig.xml") != 0) {
95     printf("Failed parsing config file! retunring... \n");
96     return -1;
97   }
98
99   Int_t debugFlag = tofHandler->GetDebugFlag();
100   printf("the debug flag is %i\n",debugFlag);
101   Int_t t0Flag = tofHandler->GetT0Flag();
102   printf("the debug flag is %i\n",t0Flag);
103
104   /* init some counters */
105   int nevents_physics=0;
106   int nevents_total=0;
107
108   Int_t iev=0;
109
110   Int_t nPDBEntriesToT = 0;
111   Int_t nDBEntriesToT = 0;
112   AliTOFHitData *HitData;
113   Int_t dummy = -1;
114   Int_t Volume[5];
115   for (Int_t i=0;i<5;i++) Volume[i]=-1;
116   AliTOFRawStream *rawStreamTOF = new AliTOFRawStream();
117   AliTOFHitDataBuffer DataBuffer;
118   AliTOFHitDataBuffer PackedDataBuffer;
119   Int_t nDBEntries = 0;
120   Int_t nPDBEntries = 0;
121   
122   struct eventHeaderStruct *event;
123   eventTypeType eventT;
124
125   /* define data source : this is argument 1 */  
126   status=monitorSetDataSource( argv[1] );
127   if (status!=0) {
128     printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
129     return -1;
130   }
131
132   /* declare monitoring program */
133   status=monitorDeclareMp( __FILE__ );
134   if (status!=0) {
135     printf("monitorDeclareMp() failed : %s\n",monitorDecodeError(status));
136     return -1;
137   }
138
139   /* define wait event timeout - 1s max */
140   monitorSetNowait();
141   monitorSetNoWaitNetworkTimeout(1000);
142   
143   /* main loop (infinite) */
144   for(;;) {
145     
146     /* check shutdown condition */
147     if (daqDA_checkShutdown()) break;
148     
149     /* get next event (blocking call until timeout) */
150     status=monitorGetEventDynamic((void **)&event);
151     if (status==MON_ERR_EOF) {
152       printf ("End of File detected\n");
153       break; /* end of monitoring file has been reached */
154     }
155     
156     if (status!=0) {
157       printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
158       break;
159     }
160     
161     /* retry if got no event */
162     if (event==NULL) continue;
163
164     iev++; 
165
166    /* use event - here, just write event id to result file */
167     nevents_total++;
168     eventT=event->eventType;
169     switch (event->eventType) {
170       
171       /* START OF RUN */
172     case START_OF_RUN:
173       break;
174       /* END START OF RUN */
175       
176       /* END OF RUN */
177     case END_OF_RUN:
178       break;
179       /* END END OF RUN */
180       
181     case PHYSICS_EVENT:
182       nevents_physics++;
183       AliRawReader *rawReader = new AliRawReaderDate((void*)event);
184       //rawReader->RequireHeader(kFALSE);
185
186       //T0 event
187       Int_t meantime = 0;     
188       AliT0RawReader *rawReaderT0 = new AliT0RawReader(rawReader,kTRUE);
189       if (!rawReaderT0->Next()) {
190         printf("T0: no raw data found!\n");
191       } 
192       else {
193         /*
194         Int_t allData[105][5];
195         for (Int_t i=0; i<105; i++) {
196           allData[i][0]=rawReaderT0->GetData(i,0);
197         }
198         meantime = allData[49][0];
199         */
200         //meantime = rawReaderT0->GetData(49,0); //OLD
201         meantime = (Int_t)((rawReaderT0->GetData(51,0)+rawReaderT0->GetData(52,0))/2.); //Alla
202         if (debugFlag > 0) {
203                 printf("time zero = %i (TDC bin)\n", meantime);   // debugging purpose
204                 printf("time zero (ns) = %i (%f) \n", meantime, (meantime*24.4-200)*1E-3);   // debugging purpose
205         }
206       }
207       
208       delete rawReaderT0;
209       rawReaderT0 = 0x0;
210       rawReader->Reset();
211       
212       //TOF event
213       dummy = -1;
214       for (Int_t ii=0; ii<5; ii++) Volume[ii]=-1;
215       rawStreamTOF->SetRawReader(rawReader);
216       //rawReader->ReadHeader();
217       rawStreamTOF->ResetBuffers();
218       rawStreamTOF->DecodeDDL(0, AliDAQ::NumberOfDdls("TOF") - 1,0);
219       nPDBEntriesToT = 0;
220       nDBEntriesToT = 0;
221       for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls("TOF"); iDDL++) {
222
223         /* read decoded data */
224         DataBuffer = rawStreamTOF->GetDataBuffer(iDDL);
225         PackedDataBuffer = rawStreamTOF->GetPackedDataBuffer(iDDL);
226         
227         /* get buffer entries */
228         nDBEntries = DataBuffer.GetEntries();
229         nPDBEntries = PackedDataBuffer.GetEntries();
230         nPDBEntriesToT+=nPDBEntries;
231         nDBEntriesToT+=nDBEntries;
232
233         //for (Int_t iHit = 0; iHit < nDBEntries; iHit++) {
234         // HitData = DataBuffer->GetHit(iHit);
235           /* store volume information */
236         // rawStreamTOF->EquipmentId2VolumeId(HitData, Volume);
237         //}
238         /* reset buffer */
239         DataBuffer.Reset();
240
241         /* read data buffer hits */
242         for (Int_t iHit = 0; iHit < nPDBEntries; iHit++) {
243           HitData = PackedDataBuffer.GetHit(iHit);
244           /* add volume information */
245           HitData->SetDDLID(iDDL);
246           rawStreamTOF->EquipmentId2VolumeId(HitData, Volume);
247           if (Volume[0]==-1 ||
248               Volume[1]==-1 ||
249               Volume[2]==-1 ||
250               Volume[3]==-1 ||
251               Volume[4]==-1) continue;
252           else {
253             dummy = Volume[3];
254             Volume[3] = Volume[4];
255             Volume[4] = dummy;
256             Int_t tofRaw = (Int_t)((Double_t)HitData->GetTime()*1E3/AliTOFGeometry::TdcBinWidth());
257             Int_t tof;
258             if (!t0Flag) tof = tofRaw;
259             else tof = tofRaw - meantime;
260             Int_t index = geom->GetIndex(Volume);
261             Float_t pos[3];
262             geom->GetPosPar(Volume,pos);
263             Float_t texp = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2])/c*1E9; //expected time in ns
264             Float_t texpBin = texp*1E3/AliTOFGeometry::TdcBinWidth(); //expected time in number of TDC bin
265             Int_t deltabin = tof-TMath::Nint(texpBin);   //to be used with real data; rounding expected time to Int_t
266             htofPartial->Fill(index,deltabin); //channel index start from 0, bin index from 1
267             //debugging printings
268             if (debugFlag > 1) {
269                     printf("tofRaw = %i, tof = %i \n",tofRaw,tof);
270             }
271             if (debugFlag > 2) {
272                     printf("sector %2d, plate %1d, strip %2d, padz %1d, padx %2d \n",Volume[0],Volume[1],Volume[2],Volume[3],Volume[4]); // too verbose
273                     printf("pos x = %f, pos y = %f, pos z = %f \n",pos[0],pos[1],pos[2]); // too verbose
274                     printf("expected time = %f (ns)\n",texp); // too verbose
275                     printf("expected time bin = %f (TDC bin)\n",texpBin); // too verbose
276                     printf("measured time bin = %i (TDC bin) with %f (ns) and ACQ bit = %i \n",tof, HitData->GetTime(), HitData->GetACQ()); // too verbose
277                     printf("index = %6d, deltabin = %d , filling index = %6d, and bin = %d\n",index, deltabin, index, deltabin); // too verbose
278             }
279
280           }
281           /* reset buffer */
282           PackedDataBuffer.Reset();
283         }
284       }
285       //if (debugFlag) {
286       //  printf(" Packed Hit Buffer Entries = %i \n",nPDBEntriesToT); // too verbose
287       //  printf(" Hit Buffer Entries = %i \n",nDBEntriesToT); // too verbose
288       //}
289
290       delete rawReader;
291       rawReader = 0x0;
292     }
293     
294     /* free resources */
295     free(event);
296     
297     /* exit when last event received, no need to wait for TERM signal */
298     if (eventT==END_OF_RUN) {
299       printf("EOR event detected\n");
300       break;
301     }
302
303   }
304
305   delete rawStreamTOF;
306   rawStreamTOF = 0x0;
307
308   delete geom;
309   geom = 0x0;
310
311   //write the Run level file   
312   TFile * fileRun = new TFile (FILE_RUN,"RECREATE"); 
313   htofPartial->Write();
314   fileRun->Close();
315
316   //write the Total file
317   TH2S *htoftot = 0x0;
318   TFile * filetot = 0x0;
319   Bool_t isThere=kFALSE;
320   const char *dirname = "./";
321   TString filename = FILE_TOTAL;
322   if((gSystem->FindFile(dirname,filename))!=NULL){
323     isThere=kTRUE;
324     printf("%s found \n",FILE_TOTAL);
325   }
326   if (isThere) {
327
328     TFile * filetot1 = new TFile (FILE_TOTAL,"READ"); 
329     //look for the file
330     if (!filetot1->IsZombie()){
331       printf("updating file %s \n",FILE_TOTAL);
332       TIter next(filetot1->GetListOfKeys());
333       TKey *key;
334       //look for the histogram
335       while ((key=(TKey*)next())){
336         const char * namekey = key->GetName();
337         if (strcmp(namekey,"htoftot")==0) {
338           printf(" histo found \n");
339           htoftot = (TH2S*) filetot1->Get("htoftot");
340           htoftot->AddDirectory(0);
341           htoftot->Add(htofPartial);
342           break;
343         }
344       }
345     }
346     filetot1->Close();
347     delete filetot1;
348     filetot1=0x0;
349   }
350   else {
351     printf(" no %s file found \n",FILE_TOTAL);
352     htoftot = new TH2S(*htofPartial);
353     htoftot->SetName("htoftot");
354     htoftot->AddDirectory(0);
355   }
356
357   filetot = new TFile (FILE_TOTAL,"RECREATE");
358   filetot->cd();
359   htoftot->Write();
360   filetot->Close();
361   
362   delete fileRun;
363   delete filetot;
364   delete htofPartial;
365   delete htoftot;
366
367   fileRun = 0x0;
368   filetot = 0x0;
369   htofPartial = 0x0;
370   htoftot = 0x0;
371
372   /* write report */
373   printf("Run #%s, received %d physics events out of %d\n",
374          getenv("DATE_RUN_NUMBER"),nevents_physics,nevents_total);
375
376   status = 0;
377
378   /* export file to FXS */
379   if (daqDA_FES_storeFile(FILE_RUN, "RUNLevel"))
380     status=-2;
381   if (daqDA_FES_storeFile(FILE_TOTAL, "DELAYS"))
382     status=-2;
383
384   return status;
385 }