]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ZDC/ZDCLASERda.cxx
recMC.C submitMC.sh
[u/mrichter/AliRoot.git] / ZDC / ZDCLASERda.cxx
CommitLineData
ead118d8 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
442e1b18 17DA Type: LDC
ead118d8 18Number of events needed: no constraint (tipically ~10^3)
19Input Files:
20Output Files: ZDCLaser.dat
21Trigger Types Used: Standalone Trigger
22
23*/
78beff0d 24#define PEDDATA_FILE "ZDCPedestal.dat"
218f916a 25#define MAPDATA_FILE "ZDCChMapping.dat"
26#define LASDATA_FILE "ZDCLaserCalib.dat"
ead118d8 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>
442e1b18 42#include <TFitter.h>
ead118d8 43
44//AliRoot
45#include <AliRawReaderDate.h>
442e1b18 46#include <AliRawEventHeaderBase.h>
ead118d8 47#include <AliZDCRawStream.h>
48
49
50/* Main routine
51 Arguments: list of DATE raw data files
52*/
53int main(int argc, char **argv) {
442e1b18 54
55 TFitter *minuitFit = new TFitter(4);
56 TVirtualFitter::SetFitter(minuitFit);
ead118d8 57
58 int status = 0;
59
60 /* log start of process */
78beff0d 61 printf("\nZDC LASER program started\n");
ead118d8 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);
78beff0d 73 //
bd94dbdd 74 // --- Histos for reference PMTs (high gain chains)
75 TH1F *hPMRefC = new TH1F("hPMRefC","hPMRefC", 100,0.,1400.);
76 TH1F *hPMRefA = new TH1F("hPMRefA","hPMRefA", 100,0.,1400.);
78beff0d 77 //
bd94dbdd 78 // --- Histos for detector PMTs (just high gain chain)
79 TH1F *hZNC[5], *hZPC[5], *hZNA[5], *hZPA[5];
80 char hnamZNC[20], hnamZPC[20], hnamZNA[20], hnamZPA[20];
81 for(Int_t j=0; j<5; j++){
82 sprintf(hnamZNC,"ZNC-tow%d",j);
83 sprintf(hnamZPC,"ZPC-tow%d",j);
84 sprintf(hnamZNA,"ZNA-tow%d",j);
85 sprintf(hnamZPA,"ZPA-tow%d",j);
86 //
87 hZNC[j] = new TH1F(hnamZNC, hnamZNC, 100, 0., 1400.);
88 hZPC[j] = new TH1F(hnamZPC, hnamZPC, 100, 0., 1400.);
89 hZNA[j] = new TH1F(hnamZNA, hnamZNA, 100, 0., 1400.);
90 hZPA[j] = new TH1F(hnamZPA, hnamZPA, 100, 0., 1400.);
91 }
ead118d8 92
93 /* open result file */
94 FILE *fp=NULL;
95 fp=fopen("./result.txt","a");
96 if (fp==NULL) {
97 printf("Failed to open file\n");
98 return -1;
99 }
100
442e1b18 101 FILE *mapFile4Shuttle;
78beff0d 102
103 // *** To analyze LASER events you MUST have a pedestal data file!!!
104 // *** -> check if a pedestal run has been analyzied
105 int read = 0;
218f916a 106 read = daqDA_DB_getFile(PEDDATA_FILE,PEDDATA_FILE);
78beff0d 107 if(read){
108 printf("\t ERROR!!! ZDCPedestal.dat file NOT FOUND in DAQ db!!!\n");
109 return -1;
110 }
111 else printf("\t ZDCPedestal.dat file retrieved from DAQ db\n");
112
113 FILE *filePed = fopen(PEDDATA_FILE,"r");
114 if (filePed==NULL) {
115 printf("\t ERROR!!! Can't open ZDCPedestal.dat file!!!\n");
116 return -1;
117 }
118
119 // 144 = 48 in-time + 48 out-of-time + 48 correlations
120 Float_t readValues[2][144], MeanPed[44], MeanPedWidth[44],
121 MeanPedOOT[44], MeanPedWidthOOT[44];
122 // ***************************************************
123 // Unless we have a narrow correlation to fit we
124 // don't fit and store in-time vs. out-of-time
125 // histograms -> mean pedstal subtracted!!!!!!
126 // ***************************************************
127 //Float_t CorrCoeff0[44], CorrCoeff1[44];
128 //
129 for(int jj=0; jj<144; jj++){
130 for(int ii=0; ii<2; ii++){
131 fscanf(filePed,"%f",&readValues[ii][jj]);
132 }
133 if(jj<48){
134 MeanPed[jj] = readValues[0][jj];
135 MeanPedWidth[jj] = readValues[1][jj];
136 //printf("\t MeanPed[%d] = %1.1f\n",jj, MeanPed[jj]);
137 }
138 else if(jj>48 && jj<96){
139 MeanPedOOT[jj-48] = readValues[0][jj];
140 MeanPedWidthOOT[jj-48] = readValues[1][jj];
141 }
142 /*else if(jj>144){
143 CorrCoeff0[jj-96] = readValues[0][jj];
144 CorrCoeff1[jj-96] = readValues[1][jj];;
145 }
146 */
147 }
ead118d8 148
149 /* report progress */
150 daqDA_progressReport(10);
151
152
153 /* init some counters */
154 int nevents_physics=0;
155 int nevents_total=0;
156
157 /* read the data files */
158 int n;
159 for (n=1;n<argc;n++) {
160
161 status=monitorSetDataSource( argv[n] );
162 if (status!=0) {
163 printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
164 return -1;
165 }
166
167 /* report progress */
168 /* in this example, indexed on the number of files */
169 daqDA_progressReport(10+80*n/argc);
170
171 /* read the file */
172 for(;;) {
173 struct eventHeaderStruct *event;
174 eventTypeType eventT;
175
176 /* get next event */
177 status=monitorGetEventDynamic((void **)&event);
178 if (status==MON_ERR_EOF) break; /* end of monitoring file has been reached */
179 if (status!=0) {
180 printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
181 return -1;
182 }
183
184 /* retry if got no event */
185 if (event==NULL) {
186 break;
187 }
188
442e1b18 189 // Initalize raw-data reading and decoding
190 AliRawReader *reader = new AliRawReaderDate((void*)event);
191 reader->Select("ZDC");
192 // --- Reading event header
78beff0d 193 //UInt_t evtype = reader->GetType();
194 //printf("\n\t ZDCLASERda -> ev. type %d\n",evtype);
195 //printf("\t ZDCLASERda -> run # %d\n",reader->GetRunNumber());
442e1b18 196 //
197 AliZDCRawStream *rawStreamZDC = new AliZDCRawStream(reader);
198
199
200 /* use event - here, just write event id to result file */
201 eventT=event->eventType;
202
203 Int_t ich=0, adcMod[48], adcCh[48], sigCode[48], det[48], sec[48];
204 if(eventT==START_OF_DATA){
7f4bde92 205
206 rawStreamZDC->SetSODReading(kTRUE);
442e1b18 207
208 if(!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
209 else{
210 while(rawStreamZDC->Next()){
211 if(rawStreamZDC->IsChMapping()){
212 adcMod[ich] = rawStreamZDC->GetADCModFromMap(ich);
213 adcCh[ich] = rawStreamZDC->GetADCChFromMap(ich);
214 sigCode[ich] = rawStreamZDC->GetADCSignFromMap(ich);
215 det[ich] = rawStreamZDC->GetDetectorFromMap(ich);
216 sec[ich] = rawStreamZDC->GetTowerFromMap(ich);
217 ich++;
218 }
219 }
220 }
221 // --------------------------------------------------------
222 // --- Writing ascii data file for the Shuttle preprocessor
218f916a 223 mapFile4Shuttle = fopen(MAPDATA_FILE,"w");
442e1b18 224 for(Int_t i=0; i<ich; i++){
225 fprintf(mapFile4Shuttle,"\t%d\t%d\t%d\t%d\t%d\t%d\n",i,
226 adcMod[i],adcCh[i],sigCode[i],det[i],sec[i]);
227 //
78beff0d 228 //printf("ZDCLASERDA.cxx -> ch.%d mod %d, ch %d, code %d det %d, sec %d\n",
442e1b18 229 // i,adcMod[i],adcCh[i],sigCode[i],det[i],sec[i]);
230 }
231 fclose(mapFile4Shuttle);
232 }
ead118d8 233
234 /* use event - here, just write event id to result file */
235 eventT=event->eventType;
236
237 if(eventT==PHYSICS_EVENT){
ead118d8 238 //
442e1b18 239 // --- Reading data header
240 reader->ReadHeader();
ead118d8 241 const AliRawDataHeader* header = reader->GetDataHeader();
242 if(header) {
243 UChar_t message = header->GetAttributes();
244 if(message & 0x20){ // DEDICATED LASER RUN
442e1b18 245 //printf("\t STANDALONE_LASER_RUN raw data found\n");
ead118d8 246 continue;
247 }
248 else{
249 printf("\t NO STANDALONE_LASER_RUN raw data found\n");
250 return -1;
251 }
252 }
442e1b18 253 else{
254 printf("\t ATTENTION! No Raw Data Header found!!!\n");
255 return -1;
256 }
257
7f4bde92 258 rawStreamZDC->SetSODReading(kTRUE);
259
ead118d8 260 if (!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
442e1b18 261 //
262 // ----- Setting ch. mapping -----
263 for(Int_t jk=0; jk<48; jk++){
264 rawStreamZDC->SetMapADCMod(jk, adcMod[jk]);
265 rawStreamZDC->SetMapADCCh(jk, adcCh[jk]);
266 rawStreamZDC->SetMapADCSig(jk, sigCode[jk]);
267 rawStreamZDC->SetMapDet(jk, det[jk]);
268 rawStreamZDC->SetMapTow(jk, sec[jk]);
269 }
ead118d8 270 //
271 while(rawStreamZDC->Next()){
272 Int_t index=-1;
bd94dbdd 273 Int_t detector = rawStreamZDC->GetSector(0);
274
275 if(rawStreamZDC->IsADCDataWord() && !(rawStreamZDC->IsUnderflow())
276 && !(rawStreamZDC->IsOverflow()) && detector!=-1){
277
278 printf(" IsADCWord %d, IsUnderflow %d, IsOverflow %d\n",
279 rawStreamZDC->IsADCDataWord(),rawStreamZDC->IsUnderflow(),rawStreamZDC->IsOverflow());
280
281 if(rawStreamZDC->GetSector(1)!=5){ // Physics signals
282 if(detector==1) index = rawStreamZDC->GetSector(1); // *** ZNC
283 else if(detector==2) index = rawStreamZDC->GetSector(1)+5; // *** ZPC
284 else if(detector==4) index = rawStreamZDC->GetSector(1)+12;// *** ZNA
285 else if(detector==5) index = rawStreamZDC->GetSector(1)+17;// *** ZPA
286 }
287 else{ // Reference PMs
288 index = (detector-1)/3+22;
289 }
290
ead118d8 291 Float_t Pedestal = MeanPed[index];
292 Float_t CorrADC = rawStreamZDC->GetADCValue() - Pedestal;
78beff0d 293
bd94dbdd 294 // **** Detector PMs
295 if(rawStreamZDC->GetSector(1)!=5 && rawStreamZDC->GetADCGain()==0){
296 // ---- side C
297 hZNC[rawStreamZDC->GetSector(1)]->Fill(CorrADC);
298 hZPC[rawStreamZDC->GetSector(1)]->Fill(CorrADC);
299 // ---- side A
300 hZNA[rawStreamZDC->GetSector(1)]->Fill(CorrADC);
301 hZPA[rawStreamZDC->GetSector(1)]->Fill(CorrADC);
ead118d8 302 }
bd94dbdd 303 // **** Reference PMs
304 if(rawStreamZDC->GetSector(1)==5 && rawStreamZDC->GetADCGain()==0){
305 // ---- PMRef chain side C
306 if(detector==1) hPMRefC->Fill(CorrADC);
307 // ---- PMRef side A
308 else if(detector==4) hPMRefA->Fill(CorrADC);
309 }
310 }//IsADCDataWord()+NOunderflow+NOoverflow
ead118d8 311 //
312 }
313 //
314 nevents_physics++;
315 //
316 delete reader;
317 delete rawStreamZDC;
318
319 }//(if PHYSICS_EVENT)
320 nevents_total++;
321
322 /* free resources */
323 free(event);
324
325 }
326 }
327
328 /* Analysis of the histograms */
329 //
bd94dbdd 330 Int_t maxBin[22], nBin[22];
331 Float_t xMax[22], maxXval[22], xlow[22];
332 Float_t mean[22], sigma[22];
333 TF1 *fun[4];
6f427255 334
bd94dbdd 335 for(Int_t k=0; k<5; k++){
336 // --- ZNC
337 maxBin[k] = hZNC[k]->GetMaximumBin();
338 nBin[k] = (hZNC[k]->GetXaxis())->GetNbins();
339 xMax[k] = (hZNC[k]->GetXaxis())->GetXmax();
340 if(nBin[k]!=0) maxXval[k] = maxBin[k]*xMax[k]/nBin[k];
341 //
342 if(maxXval[k]-150.<0.) xlow[k]=0.;
343 else xlow[k] = maxXval[k]-150.;
344 hZNC[k]->Fit("gaus","Q","",xlow[k],maxXval[k]+150.);
345 fun[k] = hZNC[k]->GetFunction("gaus");
346 mean[k] = (Float_t) (fun[k]->GetParameter(1));
347 sigma[k] = (Float_t) (fun[k]->GetParameter(2));
348 // --- ZPC
349 maxBin[k+5] = hZPC[k]->GetMaximumBin();
350 nBin[k+5] = (hZPC[k]->GetXaxis())->GetNbins();
351 xMax[k+5] = (hZPC[k]->GetXaxis())->GetXmax();
352 if(nBin[k+5]!=0) maxXval[k+5] = maxBin[k+5]*xMax[k+5]/nBin[k+5];
353 //
354 if(maxXval[k+5]-150.<0.) xlow[k+5]=0.;
355 else xlow[k+5] = maxXval[k+5]-150.;
356 hZPC[k]->Fit("gaus","Q","",xlow[k+5],maxXval[k+5]+150.);
357 fun[k+5] = hZPC[k]->GetFunction("gaus");
358 mean[k+5] = (Float_t) (fun[k+5]->GetParameter(1));
359 sigma[k+5] = (Float_t) (fun[k+5]->GetParameter(2));
360 // --- ZNA
361 maxBin[k+10] = hZNA[k]->GetMaximumBin();
362 nBin[k+10] = (hZNA[k]->GetXaxis())->GetNbins();
363 xMax[k+10] = (hZNA[k]->GetXaxis())->GetXmax();
364 if(nBin[k+10]!=0) maxXval[k+10] = maxBin[k+10]*xMax[k+10]/nBin[k+10];
365 //
366 if(maxXval[k+10]-150.<0.) xlow[k+10]=0.;
367 else xlow[k+10] = maxXval[k+10]-150.;
368 hZNA[k]->Fit("gaus","Q","",xlow[k+10],maxXval[k+10]+150.);
369 fun[k+10] = hZNA[k]->GetFunction("gaus");
370 mean[k+10] = (Float_t) (fun[k+10]->GetParameter(1));
371 sigma[k+10] = (Float_t) (fun[k+10]->GetParameter(2));
372 // --- ZPA
373 maxBin[k+15] = hZPA[k]->GetMaximumBin();
374 nBin[k+15] = (hZPA[k]->GetXaxis())->GetNbins();
375 xMax[k+15] = (hZPA[k]->GetXaxis())->GetXmax();
376 if(nBin[k+15]!=0) maxXval[k+15] = maxBin[k+15]*xMax[k+15]/nBin[k+15];
377 //
378 if(maxXval[k+15]-150.<0.) xlow[k+15]=0.;
379 else xlow[k+15] = maxXval[k+15]-150.;
380 hZPA[k]->Fit("gaus","Q","",xlow[k+15],maxXval[k+15]+150.);
381 fun[k+15] = hZPA[k]->GetFunction("gaus");
382 mean[k+15] = (Float_t) (fun[k+15]->GetParameter(1));
383 sigma[k+15] = (Float_t) (fun[k+15]->GetParameter(2));
384
385 }
6f427255 386
bd94dbdd 387 // ~~~~~~~~ PM Ref side C ~~~~~~~~
388 maxBin[20] = hPMRefC->GetMaximumBin();
389 nBin[20] = (hPMRefC->GetXaxis())->GetNbins();
390 xMax[20] = (hPMRefC->GetXaxis())->GetXmax();
391 if(nBin[20]!=0) maxXval[20] = maxBin[20]*xMax[20]/nBin[20];
392 //
393 if(maxXval[20]-150.<0.) xlow[20]=0.;
394 else xlow[20] = maxXval[20];
395 hPMRefC->Fit("gaus","Q","",xlow[20],maxXval[20]+150.);
396 fun[20] = hPMRefC->GetFunction("gaus");
397 mean[20] = (Float_t) (fun[20]->GetParameter(1));
398 sigma[20] = (Float_t) (fun[20]->GetParameter(2));
6f427255 399
bd94dbdd 400 // ~~~~~~~~ PM Ref side A ~~~~~~~~
401 maxBin[21] = hPMRefA->GetMaximumBin();
402 nBin[21] = (hPMRefA->GetXaxis())->GetNbins();
403 xMax[21] = (hPMRefA->GetXaxis())->GetXmax();
404 if(nBin[21]!=0) maxXval[21] = maxBin[21]*xMax[21]/nBin[21];
78beff0d 405 //
bd94dbdd 406 if(maxXval[21]-100.<0.) xlow[21]=0.;
407 else xlow[21] = maxXval[21];
408 hPMRefA->Fit("gaus","Q","",xlow[21],maxXval[21]+100.);
409 fun[21] = hPMRefA->GetFunction("gaus");
410 mean[21] = (Float_t) (fun[21]->GetParameter(1));
411 sigma[21] = (Float_t) (fun[21]->GetParameter(2));
412
ead118d8 413 FILE *fileShuttle;
218f916a 414 fileShuttle = fopen(LASDATA_FILE,"w");
2e33be6a 415 Int_t det[22] = {1,1,1,1,1,2,2,2,2,2,4,4,4,4,4,5,5,5,5,5,1,4};
416 Int_t quad[22] = {0,1,2,3,4,0,1,2,3,4,0,1,2,3,4,0,1,2,3,4,5,5};
bd94dbdd 417 for(Int_t i=0; i<22; i++){
2e33be6a 418 fprintf(fileShuttle,"\t%d\t%d\t%f\t%f\n",det[i],quad[i],mean[i], sigma[i]);
bd94dbdd 419 }
ead118d8 420 //
421 fclose(fileShuttle);
422 //
bd94dbdd 423 for(Int_t j=0; j<5; j++){
424 delete hZNC[j];
425 delete hZPC[j];
426 delete hZNA[j];
427 delete hZPA[j];
428 }
429 delete hPMRefC;
430 delete hPMRefA;
ead118d8 431
442e1b18 432 //delete minuitFit;
433 TVirtualFitter::SetFitter(0);
434
ead118d8 435 /* write report */
436 fprintf(fp,"Run #%s, received %d physics events out of %d\n",getenv("DATE_RUN_NUMBER"),nevents_physics,nevents_total);
437
438 /* close result file */
439 fclose(fp);
440
441 /* report progress */
442 daqDA_progressReport(90);
443
444 /* store the result file on FES */
218f916a 445 status = daqDA_FES_storeFile(MAPDATA_FILE,MAPDATA_FILE);
442e1b18 446 if(status){
447 printf("Failed to export file : %d\n",status);
448 return -1;
449 }
450 //
218f916a 451 status = daqDA_FES_storeFile(LASDATA_FILE,LASDATA_FILE);
ead118d8 452 if(status){
453 printf("Failed to export file : %d\n",status);
454 return -1;
455 }
456
457 /* report progress */
458 daqDA_progressReport(100);
459
460
461 return status;
462}