]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ZDC/ZDCLASERda.cxx
Script to create a random bad channel map.
[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*/
24
25#include <stdio.h>
26#include <stdlib.h>
27#include <Riostream.h>
28
29// DATE
30#include <event.h>
31#include <monitor.h>
32#include <daqDA.h>
33
34//ROOT
35#include <TRandom.h>
36#include <TH1F.h>
37#include <TF1.h>
38#include <TFile.h>
442e1b18 39#include <TFitter.h>
ead118d8 40
41//AliRoot
42#include <AliRawReaderDate.h>
442e1b18 43#include <AliRawEventHeaderBase.h>
ead118d8 44#include <AliZDCRawStream.h>
45
46
47/* Main routine
48 Arguments: list of DATE raw data files
49*/
50int main(int argc, char **argv) {
442e1b18 51
52 TFitter *minuitFit = new TFitter(4);
53 TVirtualFitter::SetFitter(minuitFit);
ead118d8 54
55 int status = 0;
56
57 /* log start of process */
58 printf("ZDC LASER program started\n");
59
60 /* check that we got some arguments = list of files */
61 if (argc<2) {
62 printf("Wrong number of arguments\n");
63 return -1;
64 }
65
66 // --- Histograms for LASER runs
67 // 20 signal channels + 2 reference PTMs
68 //
69 TH1F::AddDirectory(0);
70 TH1F *hZDCsideC[10], *hZDCsideA[10];
71 char nhistZDCC[50], nhistZDCA[50];
72 for(Int_t j=0; j<10; j++){
73 if(j<5){ // ZNs
74 sprintf(nhistZDCC,"ZNCtow%d",j);
75 sprintf(nhistZDCA,"ZNAtow%d",j);
76 }
77 else if(j>=5 && j<10){ // ZPs
78 sprintf(nhistZDCC,"ZPCtow%d",j);
79 sprintf(nhistZDCA,"ZPAtow%d",j);
80 }
81 hZDCsideC[j] = new TH1F(nhistZDCC, nhistZDCC, 100, 0., 1000.);
82 hZDCsideA[j] = new TH1F(nhistZDCA, nhistZDCA, 100, 0., 1000.);
83 }
84 TH1F *hPMRefsideC = new TH1F("hPMRefsideC","hPMRefsideC", 100,0.,1000.);
85 TH1F *hPMRefsideA = new TH1F("hPMRefsideA","hPMRefsideA", 100,0.,1000.);
86
87
88 /* open result file */
89 FILE *fp=NULL;
90 fp=fopen("./result.txt","a");
91 if (fp==NULL) {
92 printf("Failed to open file\n");
93 return -1;
94 }
95
442e1b18 96 FILE *mapFile4Shuttle;
97 const char *mapfName = "ZDCChMapping.dat";
ead118d8 98
99 /* report progress */
100 daqDA_progressReport(10);
101
102
103 /* init some counters */
104 int nevents_physics=0;
105 int nevents_total=0;
106
107 /* read the data files */
108 int n;
109 for (n=1;n<argc;n++) {
110
111 status=monitorSetDataSource( argv[n] );
112 if (status!=0) {
113 printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
114 return -1;
115 }
116
117 /* report progress */
118 /* in this example, indexed on the number of files */
119 daqDA_progressReport(10+80*n/argc);
120
121 /* read the file */
122 for(;;) {
123 struct eventHeaderStruct *event;
124 eventTypeType eventT;
125
126 /* get next event */
127 status=monitorGetEventDynamic((void **)&event);
128 if (status==MON_ERR_EOF) break; /* end of monitoring file has been reached */
129 if (status!=0) {
130 printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
131 return -1;
132 }
133
134 /* retry if got no event */
135 if (event==NULL) {
136 break;
137 }
138
442e1b18 139 // Initalize raw-data reading and decoding
140 AliRawReader *reader = new AliRawReaderDate((void*)event);
141 reader->Select("ZDC");
142 // --- Reading event header
143 UInt_t evtype = reader->GetType();
144 //printf("\n\t ZDCPEDESTALda -> ev. type %d\n",evtype);
145 //printf("\t ZDCPEDESTALda -> run # %d\n",reader->GetRunNumber());
146 //
147 AliZDCRawStream *rawStreamZDC = new AliZDCRawStream(reader);
148
149
150 /* use event - here, just write event id to result file */
151 eventT=event->eventType;
152
153 Int_t ich=0, adcMod[48], adcCh[48], sigCode[48], det[48], sec[48];
154 if(eventT==START_OF_DATA){
155
156 if(!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
157 else{
158 while(rawStreamZDC->Next()){
159 if(rawStreamZDC->IsChMapping()){
160 adcMod[ich] = rawStreamZDC->GetADCModFromMap(ich);
161 adcCh[ich] = rawStreamZDC->GetADCChFromMap(ich);
162 sigCode[ich] = rawStreamZDC->GetADCSignFromMap(ich);
163 det[ich] = rawStreamZDC->GetDetectorFromMap(ich);
164 sec[ich] = rawStreamZDC->GetTowerFromMap(ich);
165 ich++;
166 }
167 }
168 }
169 // --------------------------------------------------------
170 // --- Writing ascii data file for the Shuttle preprocessor
171 mapFile4Shuttle = fopen(mapfName,"w");
172 for(Int_t i=0; i<ich; i++){
173 fprintf(mapFile4Shuttle,"\t%d\t%d\t%d\t%d\t%d\t%d\n",i,
174 adcMod[i],adcCh[i],sigCode[i],det[i],sec[i]);
175 //
176 //printf("ZDCPEDESTALDA.cxx -> ch.%d mod %d, ch %d, code %d det %d, sec %d\n",
177 // i,adcMod[i],adcCh[i],sigCode[i],det[i],sec[i]);
178 }
179 fclose(mapFile4Shuttle);
180 }
ead118d8 181
182 /* use event - here, just write event id to result file */
183 eventT=event->eventType;
184
185 if(eventT==PHYSICS_EVENT){
186 //
442e1b18 187 // *** To analyze LASER events you MUST have a pedestal data file!!!
ead118d8 188 // *** -> check if a pedestal run has been analyzied
189 FILE *filePed=NULL;
190 filePed=fopen("./ZDCPedestal.dat","r");
191 if (filePed==NULL) {
192 printf("\t ERROR!!! You MUST have a ZDCPedestal.dat file!!!\n");
193 return -1;
194 }
195 // 144 = 48 in-time + 48 out-of-time + 48 correlations
196 Float_t readValues[2][144], MeanPed[44], MeanPedWidth[44],
197 MeanPedOOT[44], MeanPedWidthOOT[44],
198 CorrCoeff0[44], CorrCoeff1[44];
199 //
200 for(int jj=0; jj<144; jj++){
201 for(int ii=0; ii<2; ii++){
202 fscanf(filePed,"%f",&readValues[ii][jj]);
203 }
204 if(jj<48){
205 MeanPed[jj] = readValues[0][jj];
206 MeanPedWidth[jj] = readValues[1][jj];
207 }
208 else if(jj>48 && jj<96){
209 MeanPedOOT[jj-48] = readValues[0][jj];
210 MeanPedWidthOOT[jj-48] = readValues[1][jj];
211 }
212 else if(jj>144){
213 CorrCoeff0[jj-96] = readValues[0][jj];
214 CorrCoeff1[jj-96] = readValues[1][jj];;
215 }
216 }
217 //
442e1b18 218 // --- Reading data header
219 reader->ReadHeader();
ead118d8 220 const AliRawDataHeader* header = reader->GetDataHeader();
221 if(header) {
222 UChar_t message = header->GetAttributes();
223 if(message & 0x20){ // DEDICATED LASER RUN
442e1b18 224 //printf("\t STANDALONE_LASER_RUN raw data found\n");
ead118d8 225 continue;
226 }
227 else{
228 printf("\t NO STANDALONE_LASER_RUN raw data found\n");
229 return -1;
230 }
231 }
442e1b18 232 else{
233 printf("\t ATTENTION! No Raw Data Header found!!!\n");
234 return -1;
235 }
236
ead118d8 237 if (!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
442e1b18 238 //
239 // ----- Setting ch. mapping -----
240 for(Int_t jk=0; jk<48; jk++){
241 rawStreamZDC->SetMapADCMod(jk, adcMod[jk]);
242 rawStreamZDC->SetMapADCCh(jk, adcCh[jk]);
243 rawStreamZDC->SetMapADCSig(jk, sigCode[jk]);
244 rawStreamZDC->SetMapDet(jk, det[jk]);
245 rawStreamZDC->SetMapTow(jk, sec[jk]);
246 }
ead118d8 247 //
248 while(rawStreamZDC->Next()){
249 Int_t index=-1;
250 // Implemented only for HIGH gain chain
251 if((rawStreamZDC->IsADCDataWord()) && (rawStreamZDC->GetADCGain()==0)){
252 index = rawStreamZDC->GetADCChannel();
253 Float_t Pedestal = MeanPed[index];
254 Float_t CorrADC = rawStreamZDC->GetADCValue() - Pedestal;
255 if(rawStreamZDC->GetSector(0)==1){
256 if(rawStreamZDC->GetSector(1)==5) hPMRefsideC->Fill(CorrADC);
257 else hZDCsideC[rawStreamZDC->GetSector(1)]->Fill(CorrADC);
258 }
259 else if(rawStreamZDC->GetSector(0)==2){
260 hZDCsideC[rawStreamZDC->GetSector(1)+5]->Fill(CorrADC);
261 }
262 else if(rawStreamZDC->GetSector(0)==4){
263 if(rawStreamZDC->GetSector(1)==5) hPMRefsideA->Fill(CorrADC);
264 else hZDCsideA[rawStreamZDC->GetSector(1)]->Fill(CorrADC);
265 }
266 else if(rawStreamZDC->GetSector(0)==5){
267 hZDCsideA[rawStreamZDC->GetSector(1)+5]->Fill(CorrADC);
268 }
269 }//IsADCDataWord()
270 //
271 }
272 //
273 nevents_physics++;
274 //
275 delete reader;
276 delete rawStreamZDC;
277
278 }//(if PHYSICS_EVENT)
279 nevents_total++;
280
281 /* free resources */
282 free(event);
283
284 }
285 }
286
287 /* Analysis of the histograms */
288 //
289 Int_t maxBinC[10], maxBinA[10], maxBinRef[2];
290 Int_t nBinC[10], nBinA[10], nBinRef[2];
291 Float_t xMaxC[10], xMaxA[10], xMaxRef[2];
292 Float_t maxXvalC[10], maxXvalA[10], maxXvalRef[2];
293 Float_t xlowC[10], xlowA[10], xlowRef[10];
294 TF1 *funA[10], *funC[10], *funRef[2];
295 //
296 Float_t meanC[10], meanA[10], meanRef[2];
297 Float_t sigmaA[10], sigmaC[10], sigmaRef[10];
298 //
299 for(Int_t k=0; k<10; k++){
300 maxBinC[k] = hZDCsideC[k]->GetMaximumBin();
301 nBinC[k] = (hZDCsideC[k]->GetXaxis())->GetNbins();
302 xMaxC[k] = (hZDCsideC[k]->GetXaxis())->GetXmax();
303 maxXvalC[k] = maxBinC[k]*xMaxC[k]/nBinC[k];
304 //
305 if(maxXvalC[k]-100.<0.) {xlowC[k]=0.;}
306 else xlowC[k] = maxXvalC[k];
307 hZDCsideC[k]->Fit("gaus","Q","",xlowC[k],maxXvalC[k]+100.);
308 funC[k] = hZDCsideC[k]->GetFunction("gaus");
309 meanC[k] = (Float_t) (funC[k]->GetParameter(1));
310 sigmaC[k] = (Float_t) (funC[k]->GetParameter(2));
311 //
312 maxBinA[k] = hZDCsideA[k]->GetMaximumBin();
313 nBinA[k] = (hZDCsideA[k]->GetXaxis())->GetNbins();
314 xMaxA[k] = (hZDCsideA[k]->GetXaxis())->GetXmax();
315 maxXvalA[k] = maxBinA[k]*xMaxA[k]/nBinA[k];
316 //
317 if(maxXvalA[k]-100.<0.) {xlowA[k]=0.;}
318 else xlowA[k] = maxXvalA[k];
319 hZDCsideA[k]->Fit("gaus","Q","",xlowA[k],maxXvalA[k]+100.);
320 funA[k] = hZDCsideC[k]->GetFunction("gaus");
321 meanA[k] = (Float_t) (funA[k]->GetParameter(1));
322 sigmaA[k] = (Float_t) (funA[k]->GetParameter(2));
323 }
324 //
325 maxBinRef[0] = hPMRefsideC->GetMaximumBin();
326 nBinRef[0] = (hPMRefsideC->GetXaxis())->GetNbins();
327 xMaxRef[0] = (hPMRefsideC->GetXaxis())->GetXmax();
328 maxXvalRef[0] = maxBinRef[0]*xMaxRef[0]/nBinRef[0];
329 //
330 if(maxXvalRef[0]-100.<0.) {xlowRef[0]=0.;}
331 else xlowRef[0] = maxXvalRef[0];
332 hPMRefsideC->Fit("gaus","Q","",xlowRef[0],maxXvalRef[0]+100.);
333 funRef[0] = hPMRefsideC->GetFunction("gaus");
334 meanRef[0] = (Float_t) (funRef[0]->GetParameter(1));
335 sigmaRef[0] = (Float_t) (funRef[0]->GetParameter(2));
336 //
337 maxBinRef[1] = hPMRefsideA->GetMaximumBin();
338 nBinRef[1] = (hPMRefsideA->GetXaxis())->GetNbins();
339 xMaxRef[1] = (hPMRefsideA->GetXaxis())->GetXmax();
340 maxXvalRef[1] = maxBinRef[1]*xMaxRef[1]/nBinRef[1];
341 //
342 if(maxXvalRef[1]-100.<0.) {xlowRef[1]=0.;}
343 else xlowRef[1] = maxXvalRef[1];
344 hPMRefsideA->Fit("gaus","Q","",xlowRef[1],maxXvalRef[1]+100.);
345 funRef[1] = hPMRefsideA->GetFunction("gaus");
346 meanRef[1] = (Float_t) (funRef[1]->GetParameter(1));
347 sigmaRef[1] = (Float_t) (funRef[1]->GetParameter(2));
348 //
349 FILE *fileShuttle;
350 const char *fName = "ZDCLaser.dat";
351 fileShuttle = fopen(fName,"w");
352 for(Int_t i=0; i<10; i++) fprintf(fileShuttle,"\t%f\t%f\n",meanC[i], sigmaC[i]);
353 for(Int_t i=0; i<10; i++) fprintf(fileShuttle,"\t%f\t%f\n",meanA[i], sigmaA[i]);
354 for(Int_t i=0; i<2; i++) fprintf(fileShuttle,"\t%f\t%f\n",meanRef[i], sigmaRef[i]);
355 //
356 fclose(fileShuttle);
357 //
358 for(Int_t j=0; j<10; j++){
359 delete hZDCsideC[j];
360 delete hZDCsideA[j];
361 delete hPMRefsideC;
362 delete hPMRefsideA;
363 }
364
442e1b18 365 //delete minuitFit;
366 TVirtualFitter::SetFitter(0);
367
ead118d8 368 /* write report */
369 fprintf(fp,"Run #%s, received %d physics events out of %d\n",getenv("DATE_RUN_NUMBER"),nevents_physics,nevents_total);
370
371 /* close result file */
372 fclose(fp);
373
374 /* report progress */
375 daqDA_progressReport(90);
376
377 /* store the result file on FES */
442e1b18 378 status = daqDA_FES_storeFile(mapfName,"ZDCCHMAPPING_data");
379 if(status){
380 printf("Failed to export file : %d\n",status);
381 return -1;
382 }
383 //
ead118d8 384 status = daqDA_FES_storeFile(fName,"ZDCLASER_data");
385 if(status){
386 printf("Failed to export file : %d\n",status);
387 return -1;
388 }
389
390 /* report progress */
391 daqDA_progressReport(100);
392
393
394 return status;
395}