Updated DA from Brigitte. Steerable clock range for pedestal evaluation - default...
[u/mrichter/AliRoot.git] / VZERO / VZEROda.cxx
1 /*********************************************************************************
2 - Contact:    Brigitte Cheynis     b.cheynis@ipnl.in2p3.fr
3 - Link:       http
4 - Raw data test file :          
5 - Reference run number :              
6 - Run Type:   PHYSICS
7 - DA Type:    MON
8 - Number of events needed: >=500
9 - Input Files:  argument list
10 - Output Files: local files  VZERO_Histos.root, V0_Pedestals.dat (Online mapping)
11                 FXS file     V0_Ped_Width_Gain.dat (Offline mapping)
12 - Trigger types used: PHYSICS_EVENT
13 **********************************************************************************/
14
15 /**********************************************************************************
16 *                                                                                 *
17 * VZERO Detector Algorithm used for extracting calibration parameters             *
18 *                                                                                 *
19 * This program connects to the DAQ data source passed as argument.                *
20 * It computes calibration parameters, populates local "./V0_Ped_Width_Gain.dat"   *            
21 * file, exports it to the FES, and stores it into DAQ DB                          *
22 * The program exits when being asked to shut down (daqDA_checkshutdown)           *
23 * or on End of Run event.                                                         *
24 * We have 128 channels instead of 64 as expected for V0 due to the two sets of    *
25 * charge integrators which are used by the FEE ...                                *
26 *                                                                                 *
27 ***********************************************************************************/
28
29 // DATE
30 #include "event.h"
31 #include "monitor.h"
32 #include "daqDA.h"
33
34 //AliRoot
35 #include <AliVZERORawStream.h>
36 #include <AliRawReaderDate.h>
37 #include <AliRawReader.h>
38 #include <AliDAQ.h>
39
40 // standard
41 #include <stdio.h>
42 #include <stdlib.h>
43
44 //ROOT
45 #include "TROOT.h"
46 #include "TPluginManager.h"
47 #include <TFile.h>
48 #include <TH1F.h>
49 #include <TMath.h>
50
51 Int_t GetOfflineChannel(Int_t channel);
52
53 /* Main routine --- Arguments: monitoring data source */
54       
55 int main(int argc, char **argv) {
56
57 /* magic line from Cvetan */
58   gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
59                     "*",
60                     "TStreamerInfo",
61                     "RIO",
62                     "TStreamerInfo()");
63   int status;
64   if (argc!=2) {
65      printf("Wrong number of arguments\n");
66      return -1;
67   }
68   
69   // Online values (using FEE channel numbering), 
70   // stored into local V0_Pedestals.dat:
71   Double_t ADCmean[128];
72   Double_t ADCsigma[128];
73   Double_t PEDmean[128];
74   Double_t PEDsigma[128];
75   // Offline values(same but ordered as in aliroot for offliners)
76   // stored into V0_Ped_Width_Gain.dat: 
77   Double_t ADCmean_Off[128];
78   Double_t ADCsigma_Off[128];
79   Double_t PEDmean_Off[128];
80   Double_t PEDsigma_Off[128];
81        
82 //___________________________________________________
83 // Get cuts from V00DA.config file
84
85   Int_t    kClockMin;   // = 16;   LHC Clock Min for pedestal calculation
86   Int_t    kClockMax;   // = 19;   LHC Clock Max for pedestal calculation
87   Int_t    kLowCut;     // = 60;   low cut on signal distribution - to be tuned
88   Int_t    kHighCut;    // = 50;   high cut on pedestal distribution - to be tuned
89   
90   status = daqDA_DB_getFile("V00DA.config","./V00DA.config");
91   if (status) {
92       printf("Failed to get Config file (V00DA.config) from DAQ DB, status=%d\n", status);
93       printf("Take default values of parameters for pedestal calculation \n");
94       kClockMin  =  16; 
95       kClockMax  =  19; 
96       kLowCut    =  60;   
97       kHighCut   =  50;  
98   } else {
99       /* open the config file and retrieve cuts */
100       FILE *fpConfig = fopen("V00DA.config","r");
101       int res = fscanf(fpConfig,"%d %d %d %d ",&kClockMin,&kClockMax,&kLowCut,&kHighCut);
102       if(res!=4) {
103             printf("Failed to get values from Config file (V00DA.config): wrong file format - 4 integers are expected - \n");
104             kClockMin  =  16; 
105             kClockMax  =  19; 
106             kLowCut    =  60;   
107             kHighCut   =  50; 
108       }
109       fclose(fpConfig);
110   }
111   
112   printf("LHC Clock Min for pedestal calculation = %d; LHC Clock Max for pedestal calculation = %d; LowCut on signal = %d ; HighCut on pedestal = %d\n",
113           kClockMin, kClockMax, kLowCut, kHighCut);
114
115 //___________________________________________________
116 // Book HISTOGRAMS - dynamics of p-p collisions -
117       
118   char     ADCname[6]; 
119   char     PEDname[6]; 
120   TH1F     *hADCname[128];
121   TH1F     *hPEDname[128];  
122   
123   char  texte[12];
124   for (Int_t i=0; i<128; i++) {
125        sprintf(ADCname,"hADC%d",i);
126        sprintf(texte,"ADC cell%d",i);
127        hADCname[i]  = new TH1F(ADCname,texte,1024,-0.5, 1023.5);
128        sprintf(PEDname,"hPED%d",i);
129        sprintf(texte,"PED cell%d",i);
130        hPEDname[i]  = new TH1F(PEDname,texte,1024,-0.5, 1023.5);
131   }
132 //___________________________________________________ 
133
134   /* open result file to be exported to FES */
135   FILE *fpLocal=NULL;
136   fpLocal=fopen("./V0_Pedestals.dat","w");
137   if (fpLocal==NULL) {
138       printf("Failed to open local result file\n");
139       return -1;}
140    
141   /* open result file to be exported to FES */
142   FILE *fp=NULL;
143   fp=fopen("./V0_Ped_Width_Gain.dat","w");
144   if (fp==NULL) {
145       printf("Failed to open result file\n");
146       return -1;}
147
148   /* define data source : this is argument 1 */  
149   status=monitorSetDataSource( argv[1] );
150   if (status!=0) {
151      printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
152      return -1;
153   }
154
155   /* declare monitoring program */
156   status=monitorDeclareMp( __FILE__ );
157   if (status!=0) {
158     printf("monitorDeclareMp() failed : %s\n",monitorDecodeError(status));
159     return -1;
160   }
161
162   /* define wait event timeout - 1s max */
163   monitorSetNowait();
164   monitorSetNoWaitNetworkTimeout(1000);
165   
166   /* init counters on events */
167   int nevents_physics=0;
168   int nevents_total=0;
169
170   /* loop on events (infinite) */
171   for(;;) {
172       struct eventHeaderStruct *event;
173       eventTypeType eventT;
174
175       /* check shutdown condition */
176       if (daqDA_checkShutdown()) {break;}   
177
178       /* get next event (blocking call until timeout) */
179       status=monitorGetEventDynamic((void **)&event);
180       if (status==MON_ERR_EOF) {
181           printf ("End of File detected\n");
182       break; /* end of monitoring file has been reached */
183       }
184     
185       if (status!=0) {
186            printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
187       break;
188       }
189
190       /* retry if got no event */
191       if (event==NULL) continue;
192                
193       /* decode event */
194       eventT=event->eventType;
195         
196       switch (event->eventType){
197       
198       case START_OF_RUN:
199            break;
200       
201       case END_OF_RUN:
202            printf("End Of Run detected\n");
203            break;
204       
205       case PHYSICS_EVENT:
206            nevents_physics++;
207                          
208            AliRawReader *rawReader = new AliRawReaderDate((void*)event);
209   
210            AliVZERORawStream* rawStream  = new AliVZERORawStream(rawReader); 
211            rawStream->Next();   
212            for(Int_t i=0; i<64; i++) {
213                 Int_t nFlag = 0;
214                 for(Int_t j=kClockMin; j <= kClockMax; j++) {  // Check flags on clock range used for pedestal calculation
215                    if((rawStream->GetBBFlag(i,j)) || (rawStream->GetBGFlag(i,j))) nFlag++; 
216                 }
217                 if(nFlag == 0){       // Fill 64*2 pedestal histograms  - 2 integrators -
218                    for(Int_t j=kClockMin;j <= kClockMax;j++){
219                        Int_t Integrator = rawStream->GetIntegratorFlag(i,j);
220                        Float_t pedestal = (float)(rawStream->GetPedestal(i,j));
221                        hPEDname[i + 64 * Integrator]->Fill(pedestal);
222                    }    
223                 } 
224                 if((rawStream->GetBBFlag(i,10)) || (rawStream->GetBGFlag(i,10))){ // Charge
225                     Int_t Integrator = rawStream->GetIntegratorFlag(i,10);
226                     Float_t charge = (float)(rawStream->GetADC(i));   // Fill 64*2 ADCmax histograms 
227                     hADCname[i + 64 * Integrator]->Fill(charge);
228                 }                          
229            }    
230            delete rawStream;
231            rawStream = 0x0;      
232            delete rawReader;
233            rawReader = 0x0;                                                                      
234       } // end of switch on event type 
235         
236       nevents_total++;
237       /* free resources */
238       free(event);
239         
240       /* exit when last event received, no need to wait for TERM signal */
241       if (eventT==END_OF_RUN) {
242            printf("End Of Run event detected\n");
243       break;
244     }
245
246   }  // loop over events
247   
248   printf("%d physics events processed\n",nevents_physics);
249     
250 //___________________________________________________________________________
251 //  Computes mean values, converts FEE channels into Offline AliRoot channels
252 //  and dumps the ordered values into the output text file for SHUTTLE
253         
254   for(Int_t i=0; i<128; i++) {
255       hPEDname[i]->GetXaxis()->SetRange(0,kHighCut);
256       PEDmean[i]  = hPEDname[i]->GetMean(); 
257       PEDsigma[i] = hPEDname[i]->GetRMS(); 
258       hADCname[i]->GetXaxis()->SetRange(kLowCut,1024);
259       ADCmean[i]  = hADCname[i]->GetMean();
260       ADCsigma[i] = hADCname[i]->GetRMS(); 
261 //      printf(" i = %d, %.3f %.3f %.3f %.3f\n",i,PEDmean[i],PEDsigma[i],ADCmean[i],ADCsigma[i]);
262       fprintf(fpLocal," %.3f %.3f %.3f %.3f\n",PEDmean[i],PEDsigma[i],
263                                                ADCmean[i],ADCsigma[i]);
264       if (i < 64) {
265           Int_t j = GetOfflineChannel(i);     
266           PEDmean_Off[j]  = PEDmean[i];
267           PEDsigma_Off[j] = PEDsigma[i];
268           ADCmean_Off[j]  = ADCmean[i];
269           ADCsigma_Off[j] = ADCsigma[i]; }
270       else{
271           Int_t j = GetOfflineChannel(i-64);     
272           PEDmean_Off[j+64]  = PEDmean[i];
273           PEDsigma_Off[j+64] = PEDsigma[i];
274           ADCmean_Off[j+64]  = ADCmean[i];
275           ADCsigma_Off[j+64] = ADCsigma[i];
276       }
277   }
278   
279   for(Int_t j=0; j<128; j++) {
280 //      printf(" j = %d, %.3f %.3f %.3f %.3f\n",j,PEDmean_Off[j],PEDsigma_Off[j],
281 //                                                ADCmean_Off[j],ADCsigma_Off[j]);
282       fprintf(fp," %.3f %.3f %.3f %.3f\n",PEDmean_Off[j],PEDsigma_Off[j],
283                                           ADCmean_Off[j],ADCsigma_Off[j]);                                     
284   }
285    
286 //________________________________________________________________________
287 // Write root file with histos for users further check - just in case - 
288
289   TFile *histoFile = new TFile("VZERO_histos.root","RECREATE");
290     
291   for (Int_t i=0; i<128; i++) {
292        hADCname[i]->GetXaxis()->SetRange(0,1024);
293        hADCname[i]->Write(); 
294        hPEDname[i]->Write(); }
295
296   histoFile->Close(); 
297   delete histoFile;
298   
299 //________________________________________________________________________
300    
301   /* close local result file and FXS result file*/
302   fclose(fpLocal);
303   fclose(fp);
304  
305   /* export result file to FES */
306   status=daqDA_FES_storeFile("./V0_Ped_Width_Gain.dat","V00da_results");
307   if (status)    {
308       printf("Failed to export file : %d\n",status);
309       return -1; }
310
311   /* store result file into Online DB */
312   status=daqDA_DB_storeFile("./V0_Pedestals.dat","V00da_results");
313   if (status)    {
314       printf("Failed to store file into Online DB: %d\n",status);
315       return -1; }
316       
317   return status;
318 }
319
320  Int_t GetOfflineChannel(Int_t channel) {
321
322 // Channel mapping Online - Offline:
323  
324  Int_t fOfflineChannel[64] = {39, 38, 37, 36, 35, 34, 33, 32, 
325                               47, 46, 45, 44, 43, 42, 41, 40, 
326                               55, 54, 53, 52, 51, 50, 49, 48, 
327                               63, 62, 61, 60, 59, 58, 57, 56,
328                                7,  6,  5,  4,  3,  2,  1,  0, 
329                               15, 14, 13, 12, 11, 10,  9,  8,
330                               23, 22, 21, 20, 19, 18, 17, 16, 
331                               31, 30, 29, 28, 27, 26, 25, 24};
332  return fOfflineChannel[channel];                     
333 }