]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - ZDC/ZDCLASERda.cxx
DAs upgraded to AliZDCRawStream class
[u/mrichter/AliRoot.git] / ZDC / ZDCLASERda.cxx
... / ...
CommitLineData
1/*
2
3This program reads the DAQ data files passed as argument using the monitoring library.
4
5It computes the average event size and populates local "./result.txt" file with the
6result.
7
8The program reports about its processing progress.
9
10Messages on stdout are exported to DAQ log system.
11
12DA for ZDC standalone pedestal runs
13
14Contact: Chiara.Oppedisano@to.infn.it
15Link:
16Run Type: STANDALONE_LASER_RUN
17DA Type: LDC
18Number of events needed: no constraint (tipically ~10^3)
19Input Files:
20Output Files: ZDCLaser.dat
21Trigger Types Used: Standalone Trigger
22
23*/
24#define PEDDATA_FILE "ZDCPedestal.dat"
25#define MAPDATA_FILE "ZDCChMapping.dat"
26#define LASDATA_FILE "ZDCLaserCalib.dat"
27
28#include <stdio.h>
29#include <stdlib.h>
30#include <Riostream.h>
31
32// DATE
33#include <event.h>
34#include <monitor.h>
35#include <daqDA.h>
36
37//ROOT
38#include <TRandom.h>
39#include <TH1F.h>
40#include <TF1.h>
41#include <TFile.h>
42#include <TFitter.h>
43
44//AliRoot
45#include <AliRawReaderDate.h>
46#include <AliRawEventHeaderBase.h>
47#include <AliZDCRawStream.h>
48
49
50/* Main routine
51 Arguments: list of DATE raw data files
52*/
53int main(int argc, char **argv) {
54
55 TFitter *minuitFit = new TFitter(4);
56 TVirtualFitter::SetFitter(minuitFit);
57
58 int status = 0;
59
60 /* log start of process */
61 printf("\nZDC LASER program started\n");
62
63 /* check that we got some arguments = list of files */
64 if (argc<2) {
65 printf("Wrong number of arguments\n");
66 return -1;
67 }
68
69 // --- Histograms for LASER runs
70 // 20 signal channels + 2 reference PTMs
71 //
72 TH1F::AddDirectory(0);
73 //
74 TH1F *hPMRefChg = new TH1F("hPMRefChg","hPMRefChg", 100,0.,1000.);
75 TH1F *hPMRefAhg = new TH1F("hPMRefAhg","hPMRefAhg", 100,0.,1000.);
76 //
77 TH1F *hPMRefClg = new TH1F("hPMRefClg","hPMRefClg", 100,0.,4000.);
78 TH1F *hPMRefAlg = new TH1F("hPMRefAlg","hPMRefAlg", 100,0.,4000.);
79
80
81 /* open result file */
82 FILE *fp=NULL;
83 fp=fopen("./result.txt","a");
84 if (fp==NULL) {
85 printf("Failed to open file\n");
86 return -1;
87 }
88
89 FILE *mapFile4Shuttle;
90
91 // *** To analyze LASER events you MUST have a pedestal data file!!!
92 // *** -> check if a pedestal run has been analyzied
93 int read = 0;
94 read = daqDA_DB_getFile(PEDDATA_FILE,PEDDATA_FILE);
95 if(read){
96 printf("\t ERROR!!! ZDCPedestal.dat file NOT FOUND in DAQ db!!!\n");
97 return -1;
98 }
99 else printf("\t ZDCPedestal.dat file retrieved from DAQ db\n");
100
101 FILE *filePed = fopen(PEDDATA_FILE,"r");
102 if (filePed==NULL) {
103 printf("\t ERROR!!! Can't open ZDCPedestal.dat file!!!\n");
104 return -1;
105 }
106
107 // 144 = 48 in-time + 48 out-of-time + 48 correlations
108 Float_t readValues[2][144], MeanPed[44], MeanPedWidth[44],
109 MeanPedOOT[44], MeanPedWidthOOT[44];
110 // ***************************************************
111 // Unless we have a narrow correlation to fit we
112 // don't fit and store in-time vs. out-of-time
113 // histograms -> mean pedstal subtracted!!!!!!
114 // ***************************************************
115 //Float_t CorrCoeff0[44], CorrCoeff1[44];
116 //
117 for(int jj=0; jj<144; jj++){
118 for(int ii=0; ii<2; ii++){
119 fscanf(filePed,"%f",&readValues[ii][jj]);
120 }
121 if(jj<48){
122 MeanPed[jj] = readValues[0][jj];
123 MeanPedWidth[jj] = readValues[1][jj];
124 //printf("\t MeanPed[%d] = %1.1f\n",jj, MeanPed[jj]);
125 }
126 else if(jj>48 && jj<96){
127 MeanPedOOT[jj-48] = readValues[0][jj];
128 MeanPedWidthOOT[jj-48] = readValues[1][jj];
129 }
130 /*else if(jj>144){
131 CorrCoeff0[jj-96] = readValues[0][jj];
132 CorrCoeff1[jj-96] = readValues[1][jj];;
133 }
134 */
135 }
136
137 /* report progress */
138 daqDA_progressReport(10);
139
140
141 /* init some counters */
142 int nevents_physics=0;
143 int nevents_total=0;
144
145 /* read the data files */
146 int n;
147 for (n=1;n<argc;n++) {
148
149 status=monitorSetDataSource( argv[n] );
150 if (status!=0) {
151 printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
152 return -1;
153 }
154
155 /* report progress */
156 /* in this example, indexed on the number of files */
157 daqDA_progressReport(10+80*n/argc);
158
159 /* read the file */
160 for(;;) {
161 struct eventHeaderStruct *event;
162 eventTypeType eventT;
163
164 /* get next event */
165 status=monitorGetEventDynamic((void **)&event);
166 if (status==MON_ERR_EOF) break; /* end of monitoring file has been reached */
167 if (status!=0) {
168 printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
169 return -1;
170 }
171
172 /* retry if got no event */
173 if (event==NULL) {
174 break;
175 }
176
177 // Initalize raw-data reading and decoding
178 AliRawReader *reader = new AliRawReaderDate((void*)event);
179 reader->Select("ZDC");
180 // --- Reading event header
181 //UInt_t evtype = reader->GetType();
182 //printf("\n\t ZDCLASERda -> ev. type %d\n",evtype);
183 //printf("\t ZDCLASERda -> run # %d\n",reader->GetRunNumber());
184 //
185 AliZDCRawStream *rawStreamZDC = new AliZDCRawStream(reader);
186
187
188 /* use event - here, just write event id to result file */
189 eventT=event->eventType;
190
191 Int_t ich=0, adcMod[48], adcCh[48], sigCode[48], det[48], sec[48];
192 if(eventT==START_OF_DATA){
193
194 rawStreamZDC->SetSODReading(kTRUE);
195
196 if(!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
197 else{
198 while(rawStreamZDC->Next()){
199 if(rawStreamZDC->IsChMapping()){
200 adcMod[ich] = rawStreamZDC->GetADCModFromMap(ich);
201 adcCh[ich] = rawStreamZDC->GetADCChFromMap(ich);
202 sigCode[ich] = rawStreamZDC->GetADCSignFromMap(ich);
203 det[ich] = rawStreamZDC->GetDetectorFromMap(ich);
204 sec[ich] = rawStreamZDC->GetTowerFromMap(ich);
205 ich++;
206 }
207 }
208 }
209 // --------------------------------------------------------
210 // --- Writing ascii data file for the Shuttle preprocessor
211 mapFile4Shuttle = fopen(MAPDATA_FILE,"w");
212 for(Int_t i=0; i<ich; i++){
213 fprintf(mapFile4Shuttle,"\t%d\t%d\t%d\t%d\t%d\t%d\n",i,
214 adcMod[i],adcCh[i],sigCode[i],det[i],sec[i]);
215 //
216 //printf("ZDCLASERDA.cxx -> ch.%d mod %d, ch %d, code %d det %d, sec %d\n",
217 // i,adcMod[i],adcCh[i],sigCode[i],det[i],sec[i]);
218 }
219 fclose(mapFile4Shuttle);
220 }
221
222 /* use event - here, just write event id to result file */
223 eventT=event->eventType;
224
225 if(eventT==PHYSICS_EVENT){
226 //
227 // --- Reading data header
228 reader->ReadHeader();
229 const AliRawDataHeader* header = reader->GetDataHeader();
230 if(header) {
231 UChar_t message = header->GetAttributes();
232 if(message & 0x20){ // DEDICATED LASER RUN
233 //printf("\t STANDALONE_LASER_RUN raw data found\n");
234 continue;
235 }
236 else{
237 printf("\t NO STANDALONE_LASER_RUN raw data found\n");
238 return -1;
239 }
240 }
241 else{
242 printf("\t ATTENTION! No Raw Data Header found!!!\n");
243 return -1;
244 }
245
246 rawStreamZDC->SetSODReading(kTRUE);
247
248 if (!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
249 //
250 // ----- Setting ch. mapping -----
251 for(Int_t jk=0; jk<48; jk++){
252 rawStreamZDC->SetMapADCMod(jk, adcMod[jk]);
253 rawStreamZDC->SetMapADCCh(jk, adcCh[jk]);
254 rawStreamZDC->SetMapADCSig(jk, sigCode[jk]);
255 rawStreamZDC->SetMapDet(jk, det[jk]);
256 rawStreamZDC->SetMapTow(jk, sec[jk]);
257 }
258 //
259 while(rawStreamZDC->Next()){
260 Int_t index=-1;
261 // Getting data only for reference PMTs (sector[1]=5)
262 if((rawStreamZDC->IsADCDataWord()) && (rawStreamZDC->GetSector(1)==5)){
263 index = rawStreamZDC->GetADCChannel();
264 Float_t Pedestal = MeanPed[index];
265 Float_t CorrADC = rawStreamZDC->GetADCValue() - Pedestal;
266
267 // ==== HIGH GAIN CHAIN
268 if(rawStreamZDC->GetADCGain() == 0){
269 // %%%%% PMRef chain side C
270 if(rawStreamZDC->GetSector(0)==1) hPMRefChg->Fill(CorrADC);
271 // %%%%% PMRef side A
272 else if(rawStreamZDC->GetSector(0)==4) hPMRefAhg->Fill(CorrADC);
273 }
274 // ==== LOW GAIN CHAIN
275 else{
276 // %%%%% PMRef chain side C
277 if(rawStreamZDC->GetSector(0)==1) hPMRefClg->Fill(CorrADC);
278 // %%%%% PMRef side A
279 else if(rawStreamZDC->GetSector(0)==4) hPMRefAlg->Fill(CorrADC);
280 }
281 }//IsADCDataWord()
282 //
283 }
284 //
285 nevents_physics++;
286 //
287 delete reader;
288 delete rawStreamZDC;
289
290 }//(if PHYSICS_EVENT)
291 nevents_total++;
292
293 /* free resources */
294 free(event);
295
296 }
297 }
298
299 /* Analysis of the histograms */
300 //
301 Int_t maxBinRef[4], nBinRef[4];
302 Float_t xMaxRef[4], maxXvalRef[4], xlowRef[4];
303 Float_t meanRef[2], sigmaRef[2];
304 TF1 *funRef[4];
305
306 // ~~~~~~~~ PM Ref side C high gain chain ~~~~~~~~
307 maxBinRef[0] = hPMRefChg->GetMaximumBin();
308 nBinRef[0] = (hPMRefChg->GetXaxis())->GetNbins();
309 xMaxRef[0] = (hPMRefChg->GetXaxis())->GetXmax();
310 maxXvalRef[0] = maxBinRef[0]*xMaxRef[0]/nBinRef[0];
311 //
312 if(maxXvalRef[0]-100.<0.) {xlowRef[0]=0.;}
313 else xlowRef[0] = maxXvalRef[0];
314 hPMRefChg->Fit("gaus","Q","",xlowRef[0],maxXvalRef[0]+100.);
315 funRef[0] = hPMRefChg->GetFunction("gaus");
316 meanRef[0] = (Float_t) (funRef[0]->GetParameter(1));
317 sigmaRef[0] = (Float_t) (funRef[0]->GetParameter(2));
318
319 // ~~~~~~~~ PM Ref side A high gain chain ~~~~~~~~
320 maxBinRef[1] = hPMRefAhg->GetMaximumBin();
321 nBinRef[1] = (hPMRefAhg->GetXaxis())->GetNbins();
322 xMaxRef[1] = (hPMRefAhg->GetXaxis())->GetXmax();
323 maxXvalRef[1] = maxBinRef[1]*xMaxRef[1]/nBinRef[1];
324 //
325 if(maxXvalRef[1]-100.<0.) {xlowRef[1]=0.;}
326 else xlowRef[1] = maxXvalRef[1];
327 hPMRefAhg->Fit("gaus","Q","",xlowRef[1],maxXvalRef[1]+100.);
328 funRef[1] = hPMRefAhg->GetFunction("gaus");
329 meanRef[1] = (Float_t) (funRef[1]->GetParameter(1));
330 sigmaRef[1] = (Float_t) (funRef[1]->GetParameter(2));
331
332 // ~~~~~~~~ PM Ref side C low gain chain ~~~~~~~~
333 maxBinRef[2] = hPMRefClg->GetMaximumBin();
334 nBinRef[2] = (hPMRefClg->GetXaxis())->GetNbins();
335 xMaxRef[2] = (hPMRefClg->GetXaxis())->GetXmax();
336 maxXvalRef[2] = maxBinRef[2]*xMaxRef[2]/nBinRef[2];
337 //
338 if(maxXvalRef[2]-100.<0.) {xlowRef[2]=0.;}
339 else xlowRef[2] = maxXvalRef[2];
340 hPMRefClg->Fit("gaus","Q","",xlowRef[2],maxXvalRef[2]+100.);
341 funRef[2] = hPMRefClg->GetFunction("gaus");
342 meanRef[2] = (Float_t) (funRef[2]->GetParameter(1));
343 sigmaRef[2] = (Float_t) (funRef[2]->GetParameter(2));
344
345 // ~~~~~~~~ PM Ref side A low gain chain ~~~~~~~~
346 maxBinRef[3] = hPMRefAlg->GetMaximumBin();
347 nBinRef[3] = (hPMRefAlg->GetXaxis())->GetNbins();
348 xMaxRef[3] = (hPMRefAlg->GetXaxis())->GetXmax();
349 maxXvalRef[3] = maxBinRef[3]*xMaxRef[3]/nBinRef[3];
350 //
351 if(maxXvalRef[3]-100.<0.) {xlowRef[3]=0.;}
352 else xlowRef[3] = maxXvalRef[3];
353 hPMRefAlg->Fit("gaus","Q","",xlowRef[3],maxXvalRef[3]+100.);
354 funRef[3] = hPMRefAlg->GetFunction("gaus");
355 meanRef[3] = (Float_t) (funRef[3]->GetParameter(1));
356 sigmaRef[3] = (Float_t) (funRef[3]->GetParameter(2));
357 //
358 FILE *fileShuttle;
359 fileShuttle = fopen(LASDATA_FILE,"w");
360 for(Int_t i=0; i<4; i++) fprintf(fileShuttle,"\t%f\t%f\n",meanRef[i], sigmaRef[i]);
361 //
362 fclose(fileShuttle);
363 //
364 delete hPMRefChg;
365 delete hPMRefAhg;
366 delete hPMRefClg;
367 delete hPMRefAlg;
368
369 //delete minuitFit;
370 TVirtualFitter::SetFitter(0);
371
372 /* write report */
373 fprintf(fp,"Run #%s, received %d physics events out of %d\n",getenv("DATE_RUN_NUMBER"),nevents_physics,nevents_total);
374
375 /* close result file */
376 fclose(fp);
377
378 /* report progress */
379 daqDA_progressReport(90);
380
381 /* store the result file on FES */
382 status = daqDA_FES_storeFile(MAPDATA_FILE,MAPDATA_FILE);
383 if(status){
384 printf("Failed to export file : %d\n",status);
385 return -1;
386 }
387 //
388 status = daqDA_FES_storeFile(LASDATA_FILE,LASDATA_FILE);
389 if(status){
390 printf("Failed to export file : %d\n",status);
391 return -1;
392 }
393
394 /* report progress */
395 daqDA_progressReport(100);
396
397
398 return status;
399}