Adding printouts (R. Preghenella)
[u/mrichter/AliRoot.git] / TOF / DA / TOFphysicsda.cxx
1 /*
2  
3  TOF DA for online calibration
4  
5  Version: 1.0
6  Contact: Roberto.Preghenella@bo.infn.it
7  Reference Run: 115401
8  Run Type: PHYSICS
9  DA Type: MON
10  Number of events needed:
11  Input Files: no input
12  Output Files: TOFdaHits.root TOFdaReadout.root
13  Event types used: PHYSICS_EVENT
14  
15  */
16
17 #define FILE_HITS "TOFdaHits.root"
18 #define FILE_READOUT "TOFdaReadout.root"
19
20 #define EVENT_DETECTOR_TOF 5
21 #define READOUT_INFO_HISTO 1
22
23 // DATE
24 #include "event.h"
25 #include "monitor.h"
26 #include "daqDA.h"
27
28 #include <stdio.h>
29 #include <stdlib.h>
30 #include <errno.h>
31
32 //AliRoot
33 #include "TROOT.h"
34 #include "AliTOFRawStream.h"
35 #include "AliRawReaderDate.h"
36 #include "AliRawReader.h"
37 #include "AliDAQ.h"
38 #include "AliTOFHitData.h"
39 #include "AliTOFHitDataBuffer.h"
40 #include "AliTOFDaConfigHandler.h"
41 #include "AliTOFHitField.h"
42 #include "AliLog.h"
43 #include "AliTOFGeometry.h"
44 #include "AliTOFDecoderV2.h"
45 #include "AliTOFDecoderSummaryData.h"
46 #include "AliTOFDRMSummaryData.h"
47 #include "AliTOFTRMSummaryData.h"
48 #include "AliTOFChainSummaryData.h"
49 #include "AliTOFTDCHitBuffer.h"
50 #include "AliTOFTDCHit.h"
51 #include "AliTOFTDCErrorBuffer.h"
52 #include "AliTOFTDCError.h"
53
54 //ROOT
55 #include "TFile.h"
56 #include "TKey.h"
57 #include "TH2S.h"
58 #include "TObject.h"
59 #include "TMath.h"
60 #include "TSystem.h"
61 #include "TROOT.h"
62 #include "TPluginManager.h"
63 #include "TSAXParser.h"
64 #include "TTree.h"
65
66 /* Main routine
67  Arguments:
68  1- monitoring data source
69  */
70 int
71 main(int argc, char **argv)
72 {
73   
74   printf("in main\n"); // remove me
75     
76     /* magic line from Rene */
77     gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
78                                           "*",
79                                           "TStreamerInfo",
80                                           "RIO",
81                                           "TStreamerInfo()");
82     
83     
84     /* log start of process */
85     printf("TOF DA started\n");
86     
87     /* check that we got some arguments = list of files */
88     if (argc!=2) {
89         printf("Wrong number of arguments\n");
90         return -1;
91     }
92     
93     /*
94      * CONFIG
95      */
96     
97     printf("retrieve config file\n"); // remove me
98
99
100     /* retrieve config file */
101     int getConfigFile = daqDA_DB_getFile("TOFPhysicsConfig.xml","TOFPhysicsConfig.xml");
102     if (getConfigFile != 0){
103         printf("Failed to retrieve config file from DB! returning...\n");
104         return -1;
105     }
106
107     printf("parse config file\n"); // remove me
108
109     /* parse config file */
110     AliTOFDaConfigHandler* tofHandler = new AliTOFDaConfigHandler();
111     TSAXParser *parser = new TSAXParser();
112     parser->ConnectToHandler("AliTOFDaConfigHandler", tofHandler);
113     if (parser->ParseFile("./TOFPhysicsConfig.xml") != 0) {
114         printf("Failed parsing config file! retunring... \n");
115         return -1;
116     }
117     /* setup config params */
118     Int_t meanMultiplicity = tofHandler->GetMeanMultiplicity(); /* average expected TOF multiplicity */
119     Int_t maxHits = tofHandler->GetMaxHits(); /* max number of hits to be collected */
120     printf("current settings:\n");
121     printf(" - meanMultiplicity = %d\n", meanMultiplicity);
122     printf(" - maxHits          = %d\n", maxHits);
123     /* constants */
124     const Int_t nChannels = 157248;
125     Int_t noiseCheckTrigger = 10; /* first noise check after 10 events */
126     Float_t meanChannelRate = (Float_t)meanMultiplicity / (Float_t)nChannels; /* average expected channel rate (hits/event) */
127     Float_t noiseThreshold = 10. * meanChannelRate; /* noise threshold (hits/event) */
128     Int_t minNoiseHits = 10; /* min number of channel hits to check noise */
129     /* counters and flags */
130     Int_t nPhysicsEvents, nCollectedPhysicsEvents, totHits;
131     Int_t nChHits[nChannels];
132     Bool_t inhibitCollection;
133     Bool_t noiseFlag[nChannels];
134     /* variables */
135     Int_t ddl, slot, trm, chain, tdc, channel, index, timebin, totbin, deltaBC, l0l1latency, det[5], dummy;
136     Float_t noiseHitThreshold;
137     
138     /*
139      * INIT
140      */
141     
142     printf("init\n"); // remove me
143
144     /* init counters and flags */
145     nPhysicsEvents = 0;
146     nCollectedPhysicsEvents = 0;
147     totHits = 0;
148     inhibitCollection = kFALSE;
149     for (Int_t ich = 0; ich < nChannels; ich++) {
150         nChHits[ich] = 0;
151         noiseFlag[ich] = kFALSE;
152     }
153     
154     /* TOF raw data handling */
155     AliTOFRawStream *rawStream = new AliTOFRawStream();
156     AliTOFDecoderV2 *decoder = rawStream->GetDecoderV2();
157     AliTOFDecoderSummaryData *decodersd;
158     AliTOFDRMSummaryData *drmsd;
159     //  AliTOFLTMSummaryData *ltmsd;
160     AliTOFTRMSummaryData *trmsd;
161     AliTOFChainSummaryData *chainsd;
162     AliTOFTDCHitBuffer *hitBuffer;
163     AliTOFTDCHit *hit;
164     AliTOFTDCErrorBuffer *errorBuffer;
165     AliTOFTDCError *error;
166     UShort_t errorFlags;
167     UChar_t *data = 0x0;
168     Int_t dataSize;
169     Int_t dataWords;
170     Int_t currentDDL;
171     const AliRawDataHeader *currentCDH;
172     const AliRawDataHeaderV3 *currentCDHV3;
173     Int_t currentMiniEventID;
174     Int_t currentEventID1;
175     Int_t currentL0BCID ;
176     Int_t currentBunchID;
177     Bool_t skipTRM, skipChain;
178     Int_t trmIndex, chainIndex, tdcIndex;
179     Double_t chainEff;
180     
181     printf("open files and create histos\n"); // remove me
182
183     /* open HITS output file */
184     TFile *fileOutHits = new TFile(FILE_HITS, "RECREATE");
185     /* create hit field data structure */
186     AliTOFHitField *hitField = new AliTOFHitField();
187     /* create temporary tree */
188     TTree *tempTree = new TTree("tempTree", "temporary tree");
189     tempTree->Branch("hit", "AliTOFHitField", &hitField);
190     /* create output tree */
191     TTree *outTree = new TTree("hitTree", "hit tree");
192     outTree->Branch("hit", "AliTOFHitField", &hitField);
193     
194     /* open READOUT output file */
195     TFile *fileOutReadout = new TFile(FILE_READOUT, "RECREATE");
196     /* create chain readout efficiency histo */
197     TH1F *hChainEfficiency = new TH1F("hChainEfficiency", "Chain efficiency;chain;efficiency", 1440, 0., 1440.);
198     
199 #if READOUT_INFO_HISTO
200     /* create TRM data histo */
201     TH1F *hTRMData = new TH1F("hTRMData", "TRM data;TRM;frequency", 720, 0., 720.);
202     /* create TRM empty event frequency histo */
203     TH1F *hTRMEmptyEvent = new TH1F("hTRMEmptyEvent", "TRM empty event error;TRM;frequency", 720, 0., 720.);
204     /* create TRM bad event counter frequency histo */
205     TH1F *hTRMBadEventCounter = new TH1F("hTRMBadEventCounter", "TRM bad event counter;TRM;frequency", 720, 0., 720.);
206     /* create TRM bad CRC frequency histo */
207     TH1F *hTRMBadCRC = new TH1F("hTRMBadCRC", "TRM bad CRC;TRM;frequency", 720, 0., 720.);
208     
209     /* create chain data histo */
210     TH1F *hChainData = new TH1F("hChainData", "Chain data;chain;frequency", 1440, 0., 1440.);
211     /* create chain bad status frequency histo */
212     TH1F *hChainBadStatus = new TH1F("hChainBadStatus", "Chain bad status;chain;frequency", 1440, 0., 1440.);
213     /* create chain bad event counter frequency histo */
214     TH1F *hChainBadEventCounter = new TH1F("hChainBadEventCounter", "Chain bad event counter;chain;status;frequency", 1440, 0., 1440.);
215     
216     /* create TDC error frequency histo */
217     TH1F *hTDCError = new TH1F("hTDCError", "TDC error;TDC;frequency", 21600, 0., 21600.);
218     /* create TDC error flags frequency histo */
219     TH2F *hTDCErrorFlags = new TH2F("hTDCErrorFlags", "TDC error flags;TDC;error flag;frequency", 21600, 0., 21600., 15, 0., 15);
220 #endif
221     
222     
223     /*
224      * ONLINE MONITOR
225      */
226     
227     printf("monitor declare table\n"); // remove me
228
229     AliLog::SetGlobalLogLevel(AliLog::kFatal);
230     struct eventHeaderStruct *event;
231     int ret;
232     /* define monitoring table */
233     char *monTable[5] = {
234       (char *)"ALL", (char *)"no",
235       (char *)"PHY", (char *)"yes",
236       (char *)NULL
237     };
238     ret = monitorDeclareTable(monTable);
239     if (ret != 0) {
240         printf("monitorDeclareTable() failed: %s\n", monitorDecodeError(ret));
241         return -1;
242     }
243
244     printf("monitor set data source\n"); // remove me 
245
246     /* define data source : this is argument 1 */
247     ret = monitorSetDataSource(argv[1]);
248     if (ret != 0) {
249         printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(ret));
250         return -1;
251     }
252
253     printf("monitor declare mp\n"); // remove me
254
255     /* declare monitoring program */
256     ret = monitorDeclareMp("TOFdaPhysics");
257     if (ret != 0) {
258         printf("monitorDeclareMp() failed : %s\n",monitorDecodeError(ret));
259         return -1;
260     }
261
262     printf("monitor set no wait\n"); // remove me
263
264     /* define wait event timeout - 1s max */
265     monitorSetNowait();
266
267     printf("monitor set no wait network timeout\n"); // remove me
268
269     monitorSetNoWaitNetworkTimeout(1000);
270     
271     /* variables */
272
273     printf("endless loop\n"); // remove me
274
275     /* loop over events */
276     while (1) {
277         
278         /* check shutdown condition */
279         if (daqDA_checkShutdown()) break;
280         
281         /*
282          * NOISE CHECK
283          */
284         
285         /* check inhibit collection */
286         if (!inhibitCollection) {
287             /* check number of events and check noise */
288             if (nCollectedPhysicsEvents >= noiseCheckTrigger || totHits >= maxHits) {
289                 noiseHitThreshold = noiseThreshold * nCollectedPhysicsEvents;
290                 printf("noise check triggered after %d events: threshold is %f hits\n", nCollectedPhysicsEvents, noiseHitThreshold);
291                 /* loop over all channels */
292                 for (Int_t ich = 0; ich < nChannels; ich++) {
293                     /* check */
294                     if (nChHits[ich] < minNoiseHits || noiseFlag[ich] || nChHits[ich] < noiseHitThreshold) continue;
295                     printf("channel %d tagged as noisy (%d hits): disabled\n", ich, nChHits[ich]);
296                     noiseFlag[ich] = kTRUE;
297                     totHits -= nChHits[ich];
298                 } /* end of loop over all channels */
299                 /* set new noise check trigger value */
300                 noiseCheckTrigger *= 10;
301             } /* end of noise check */
302             
303             /* inhibit hit collection when maximum number of hits exceeded */
304             if (totHits >= maxHits) {
305                 printf("maximum number of hits exceeded (%d): inhibit hit collection\n", maxHits);
306                 inhibitCollection = kTRUE;
307             }
308         }
309         
310         /*
311          * GET EVENT
312          */
313         
314         /* get next event (blocking call until timeout) */
315         ret = monitorGetEventDynamic((void **)&event);
316         if (ret == MON_ERR_EOF) {
317             printf ("End of File detected\n");
318             break; /* end of monitoring file has been reached */
319         }
320         if (ret != 0) {
321             printf("monitorGetEventDynamic() failed (ret=%d errno=%d): %s\n", ret, errno, monitorDecodeError(ret));
322             break;
323         }
324         /* retry if got no event */
325         if (event==NULL) continue;
326         /* check TOF in partecipating detectors */
327         if (!TEST_DETECTOR_IN_PATTERN(event->eventDetectorPattern, EVENT_DETECTOR_TOF)) {
328             free(event);
329             continue;
330         }
331         /* check event type */
332         if (event->eventType != PHYSICS_EVENT) {
333             printf("not a physics event: %d\n", event->eventType);
334             free(event);
335             continue;
336         }
337         /* increment number of physics events */
338         nPhysicsEvents++;
339         if (!inhibitCollection) nCollectedPhysicsEvents++;
340         
341         /*
342          * DECODE EVENT
343          */
344         
345         /* create and setup raw reader */
346         AliRawReader *rawReader = new AliRawReaderDate((void *)event);
347         rawReader->Reset();
348         rawReader->Select("TOF", 0, AliDAQ::NumberOfDdls("TOF") - 1);
349         /* setup raw stream */
350         rawStream->SetRawReader(rawReader);
351         
352         /* loop over DDLs - rawReader->ReadHeader() */
353         while (rawReader->ReadHeader()) {
354             
355             /* read equipment data */
356             dataSize = rawReader->GetDataSize();
357             data = new UChar_t[dataSize];
358             if (!rawReader->ReadNext(data, dataSize)){
359                 delete [] data;
360                 continue;
361             }
362             
363             /* decode data */
364             dataWords = dataSize / 4;
365             decoder->Decode((UInt_t *)data, dataWords);
366             delete [] data;
367             
368             /* read equipment info */
369             currentDDL = rawReader->GetDDLID();
370             currentCDH = rawReader->GetDataHeader();
371             currentCDHV3 = rawReader->GetDataHeaderV3();
372             currentMiniEventID = currentCDH ? currentCDH->GetMiniEventID() : -1;
373             currentMiniEventID = currentCDHV3 ? currentCDHV3->GetMiniEventID(): currentMiniEventID;
374             currentEventID1 = currentCDH ? currentCDH->GetEventID1() : -1;
375             currentEventID1 = currentCDHV3? currentCDHV3->GetEventID1() : currentEventID1;
376             
377             /* read decoder summary data */
378             decodersd = decoder->GetDecoderSummaryData();
379             
380             /* check DRM header/trailer */
381             drmsd = decodersd->GetDRMSummaryData();
382             if (!drmsd->GetHeader() || !drmsd->GetTrailer()) continue;
383             
384             /* get DRM data */
385             currentL0BCID = drmsd->GetL0BCID();
386             
387             /* loop over TRM to get hits */
388             for (Int_t itrm = 0; itrm < 10; itrm++) {
389                 trmsd = drmsd->GetTRMSummaryData(itrm);
390                 trmIndex = itrm + 10 * currentDDL;
391                 skipTRM = kFALSE;
392                 
393                 /* check header/trailer */
394                 if (!trmsd->GetHeader() || !trmsd->GetTrailer()) continue;
395                 
396 #if READOUT_INFO_HISTO
397                 /* fill TRM data */
398                 hTRMData->Fill(trmIndex);
399                 /* fill TRM empty event */
400                 if (trmsd->GetEBit() != 0) {
401                     hTRMEmptyEvent->Fill(trmIndex);
402                     skipTRM = kTRUE;
403                 }
404                 /* fill TRM bad event counter */
405                 if (trmsd->GetEventCounter() != drmsd->GetLocalEventCounter()) {
406                     hTRMBadEventCounter->Fill(trmIndex);
407                     skipTRM = kTRUE;
408                 }
409                 /* fill TRM bad CRC */
410                 if (trmsd->GetEventCRC() != trmsd->GetDecoderCRC()) {
411                     hTRMBadCRC->Fill(trmIndex);
412                     skipTRM = kTRUE;
413                 }
414 #else
415                 /* check bad condition and skip TRM */
416                 if ( trmsd->GetEBit() != 0 ||
417                     trmsd->GetEventCounter() != drmsd->GetLocalEventCounter() ||
418                     trmsd->GetEventCRC() != trmsd->GetDecoderCRC() ) continue;
419 #endif
420                 
421                 /* loop over chains */
422                 for (Int_t ichain = 0; ichain < 2; ichain++) {
423                     chainsd = trmsd->GetChainSummaryData(ichain);
424                     chainIndex = ichain + 2 * itrm + 20 * currentDDL;
425                     skipChain = kFALSE;
426                     
427                     /* check header/trailer */
428                     if (!chainsd->GetHeader() || !chainsd->GetTrailer()) continue;
429                     currentBunchID = chainsd->GetBunchID();
430                     hitBuffer = chainsd->GetTDCPackedHitBuffer();
431                     errorBuffer = chainsd->GetTDCErrorBuffer();
432                     
433 #if READOUT_INFO_HISTO
434                     /* fill chain data */
435                     hChainData->Fill(chainIndex);
436                     /* check chain bad status */
437                     if (chainsd->GetStatus() != 0) {
438                         hChainBadStatus->Fill(chainIndex);
439                         skipChain = kTRUE;
440                     }
441                     /* check chain bad event counter */
442                     if (chainsd->GetEventCounter() != drmsd->GetLocalEventCounter()) {
443                         hChainBadEventCounter->Fill(chainIndex);
444                         skipChain = kTRUE;
445                     }
446                     /* fill TDC error frequency histo */
447                     for (Int_t ierr = 0; ierr < errorBuffer->GetEntries(); ierr++) {
448                         error = errorBuffer->GetError(ierr);
449                         tdc = error->GetTDCID();
450                         tdcIndex = tdc + 15 * ichain + 30 * itrm + 300 * currentDDL;
451                         hTDCError->Fill(tdcIndex);
452                         errorFlags = error->GetErrorFlags();
453                         for (Int_t ierflg = 0; ierflg < 15; ierflg++)
454                             if (errorFlags & (1 << ierflg))
455                                 hTDCErrorFlags->Fill(tdcIndex, ierflg);
456                     }
457 #else
458                     /* check bad condition and skip chain */
459                     if ( chainsd->GetStatus() != 0 ||
460                         chainsd->GetEventCounter() != drmsd->GetLocalEventCounter() ) continue;
461 #endif
462                     
463                     /*
464                      * CHAIN READOUT EFFICIENCY
465                      */
466                     
467                     /* compute number of available channels removing TDCs in error */
468                     chainEff = (120. - 8. * errorBuffer->GetEntries()) / 120.;
469                     /* fill chain readout efficiency histo */
470                     if (!skipTRM && !skipChain)
471                         hChainEfficiency->Fill(chainIndex, chainEff);
472                     
473                     /*
474                      * HIT MANIPULATION
475                      */
476                     
477                     /* check inhibit collection */
478                     if (inhibitCollection) continue;
479                     
480                     /* loop over hits in buffer */
481                     for (Int_t ihit = 0; ihit < hitBuffer->GetEntries(); ihit++) {
482                         
483                         /* get hit */
484                         hit = hitBuffer->GetHit(ihit);
485                         
486                         /* get channel info */
487                         ddl = currentDDL;
488                         slot = trmsd->GetSlotID();
489                         trm = slot - 3;
490                         chain = chainsd->GetChain();
491                         tdc = hit->GetTDCID();
492                         channel = hit->GetChan();
493                         /* get index */
494                         rawStream->EquipmentId2VolumeId(ddl, slot, chain, tdc, channel, det);
495                         dummy = det[4];
496                         det[4] = det[3];
497                         det[3] = dummy;
498                         /* check valid index */
499                         if (det[0] < 0 || det[0] > 17 ||
500                             det[1] < 0 || det[1] > 5 ||
501                             det[2] < 0 || det[2] > 18 ||
502                             det[3] < 0 || det[3] > 1 ||
503                             det[4] < 0 || det[4] > 47) continue;
504                         index = AliTOFGeometry::GetIndex(det);
505                         
506                         /* check noise flag */
507                         if (noiseFlag[index]) continue;
508                         /* increment number of channel hits and total hits */
509                         nChHits[index]++;
510                         totHits++;
511                         /* get signal info */
512                         timebin = hit->GetHitTime();
513                         totbin = hit->GetTOTWidth();
514                         deltaBC = currentBunchID - currentEventID1;
515                         l0l1latency = currentMiniEventID - currentL0BCID;
516                         /* set hit field data */
517                         hitField->SetIndex(index);
518                         hitField->SetTimeBin(timebin);
519                         hitField->SetTOTBin(totbin);
520                         hitField->SetDeltaBC(deltaBC);
521                         hitField->SetL0L1Latency(l0l1latency);
522                         /* fill temp tree */
523                         tempTree->Fill();
524                         
525                     } /* end of loop over hits in buffer */
526                 } /* end of loop over chains */
527             } /* end of loop over TRMs */
528         } /* end of loop over DDLs - rawReader->ReadHeader() */
529         
530         /* delete raw reader */
531         delete rawReader;
532         /* free event */
533         free(event);
534         
535     } /* end of loop over events */
536     
537     /* final noise check */
538     noiseHitThreshold = noiseThreshold * nCollectedPhysicsEvents;
539     printf("final noise check after collectiong %d events: threshold is %f hits\n", nCollectedPhysicsEvents, noiseHitThreshold);
540     /* loop over all channels */
541     for (Int_t ich = 0; ich < nChannels; ich++) {
542         /* check */
543         if (nChHits[ich] < minNoiseHits || noiseFlag[ich] || nChHits[ich] < noiseHitThreshold) continue;
544         printf("channel %d tagged as noisy (%d hits): disabled\n", ich, nChHits[ich]);
545         noiseFlag[ich] = kTRUE;
546         totHits -= nChHits[ich];
547     } /* end of loop over all channels */
548     
549     /* copy hits into output tree from temp tree */
550     printf("copy hits from temporary tree into output tree\n");
551     printf("temporary tree contains %d hits\n", (Int_t)tempTree->GetEntries());
552     for (Int_t ihit = 0; ihit < tempTree->GetEntries(); ihit++) {
553         /* get entry */
554         tempTree->GetEntry(ihit);
555         /* check noise flag */
556         if (noiseFlag[hitField->GetIndex()]) continue;
557         /* fill output tree */
558         outTree->Fill();
559     } /* end of copy hits into output tree from temp tree */
560     printf("output tree contains %d hits\n", (Int_t)outTree->GetEntries());
561     
562     /* write output tree on HITS file */
563     fileOutHits->cd();
564     outTree->Write();
565     fileOutHits->Close();
566     /* export file to FXS */
567     if (daqDA_FES_storeFile(FILE_HITS, "HITS"))
568         return -2;
569     
570 #if READOUT_INFO_HISTO
571     hTRMData->Sumw2();
572     hTRMEmptyEvent->Sumw2();
573     hTRMBadEventCounter->Sumw2();
574     hTRMBadCRC->Sumw2();
575     
576     hChainData->Sumw2();
577     hChainBadStatus->Sumw2();
578     hChainBadEventCounter->Sumw2();
579     
580     hTDCError->Sumw2();
581     hTDCErrorFlags->Sumw2();
582     
583     /* divide histos */
584     if (nPhysicsEvents > 0) {
585         hTRMEmptyEvent->Divide(hTRMData);
586         hTRMBadEventCounter->Divide(hTRMData);
587         hTRMBadCRC->Divide(hTRMData);
588         hTRMData->Scale(1. / nPhysicsEvents);
589         
590         hChainBadStatus->Divide(hChainData);
591         hChainBadEventCounter->Divide(hChainData);
592         hChainData->Scale(1. / nPhysicsEvents);
593         
594         hTDCError->Scale(1. / nPhysicsEvents);
595         hTDCErrorFlags->Scale(1. / nPhysicsEvents);
596     }
597 #endif
598     
599     /* scale chain efficiency by number of physics events */
600     printf("found %d physics events\n", nPhysicsEvents);
601     hChainEfficiency->Sumw2();
602     if (nPhysicsEvents > 0) hChainEfficiency->Scale(1. / (nPhysicsEvents));
603     
604     /* write efficiency histo on READOUT file */
605     fileOutReadout->cd();
606 #if READOUT_INFO_HISTO
607     hTRMData->Write();
608     hTRMEmptyEvent->Write();
609     hTRMBadEventCounter->Write();
610     hTRMBadCRC->Write();
611     
612     hChainData->Write();
613     hChainBadStatus->Write();
614     hChainBadEventCounter->Write();
615     
616     hTDCError->Write();
617     hTDCErrorFlags->Write();
618 #endif
619     hChainEfficiency->Write();
620     fileOutReadout->Close();
621     /* export file to FXS */
622     if (daqDA_FES_storeFile(FILE_READOUT, "READOUT"))
623         return -2;
624     
625     return 0;
626 }