ZDC/ZDCLASERda.cxx

   1 /*
   2
   3 This program reads the DAQ data files passed as argument using the monitoring library.
   4
   5 It computes the average event size and populates local "./result.txt" file with the
   6 result.
   7
   8 The program reports about its processing progress.
   9
  10 Messages on stdout are exported to DAQ log system.
  11
  12 DA for ZDC standalone pedestal runs
  13
  14 Contact: Chiara.Oppedisano@to.infn.it
  15 Link:
  16 Run Type: STANDALONE_LASER_RUN
  17 DA Type: LDC
  18 Number of events needed: no constraint (tipically ~10^3)
  19 Input Files:
  20 Output Files: ZDCLaser.dat
  21 Trigger Types Used: Standalone Trigger
  22
  23 */
  24 #define PEDDATA_FILE  "ZDCPedestal.dat"
  25 #define MAPDATA_FILE  "ZDCChMapping.dat"
  26 #define LASHISTO_FILE "ZDCLaserHisto.root"
  27 #define LASDATA_FILE  "ZDCLaserCalib.dat"
  28
  29 #include <stdio.h>
  30 #include <stdlib.h>
  31 #include <Riostream.h>
  32
  33 // DATE
  34 #include <event.h>
  35 #include <monitor.h>
  36 #include <daqDA.h>
  37
  38 //ROOT
  39 #include <TROOT.h>
  40 #include <TPluginManager.h>
  41 #include <TH1F.h>
  42 #include <TF1.h>
  43 #include <TFile.h>
  44 #include <TFitter.h>
  45 #include "TMinuitMinimizer.h"
  46
  47 //AliRoot
  48 #include <AliRawReaderDate.h>
  49 #include <AliRawEventHeaderBase.h>
  50 #include <AliZDCRawStream.h>
  51
  52
  53 /* Main routine
  54       Arguments: list of DATE raw data files
  55 */
  56 int main(int argc, char **argv) {
  57
  58   gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
  59                                         "*",
  60                                         "TStreamerInfo",
  61                                         "RIO",
  62                                         "TStreamerInfo()");
  63
  64   TMinuitMinimizer m;
  65   gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer", "Minuit","TMinuitMinimizer",
  66       "Minuit", "TMinuitMinimizer(const char *)");
  67   TVirtualFitter::SetDefaultFitter("Minuit");
  68
  69   int status = 0;
  70   int const kNChannels = 24;
  71
  72   /* log start of process */
  73   printf("\n ZDC LASER program started\n");
  74
  75   /* check that we got some arguments = list of files */
  76   if (argc<2) {
  77     printf("Wrong number of arguments\n");
  78     return -1;
  79   }
  80
  81   // --- Histograms for LASER runs
  82   //     20 signal channels + 2 reference PTMs
  83   //
  84   TH1F::AddDirectory(0);
  85   // --- Histos for reference PMTs (high gain chains)
  86   TH1F *hPMRefChg = new TH1F("hPMRefChg","hPMRefChg", 100,0.,1400.);
  87   TH1F *hPMRefAhg = new TH1F("hPMRefAhg","hPMRefAhg", 100,0.,1400.);
  88   TH1F *hPMRefClg = new TH1F("hPMRefClg","hPMRefClg", 100,0.,4000.);
  89   TH1F *hPMRefAlg = new TH1F("hPMRefAlg","hPMRefAlg", 100,0.,4000.);
  90   //
  91   // --- Histos for detector PMTs
  92   TH1F *hZNChg[5], *hZPChg[5], *hZNAhg[5], *hZPAhg[5], *hZEMhg[2];
  93   TH1F *hZNClg[5], *hZPClg[5], *hZNAlg[5], *hZPAlg[5], *hZEMlg[2];
  94   char hnamZNChg[20], hnamZPChg[20], hnamZNAhg[20], hnamZPAhg[20];
  95   char hnamZNClg[20], hnamZPClg[20], hnamZNAlg[20], hnamZPAlg[20];
  96   char hnamZEMhg[20], hnamZEMlg[20];
  97   for(Int_t j=0; j<5; j++){
  98     sprintf(hnamZNChg,"ZNChg-tow%d",j);
  99     sprintf(hnamZPChg,"ZPChg-tow%d",j);
 100     sprintf(hnamZNAhg,"ZNAhg-tow%d",j);
 101     sprintf(hnamZPAhg,"ZPAhg-tow%d",j);
 102     //
 103     hZNChg[j] = new TH1F(hnamZNChg, hnamZNChg, 100, 0., 1400.);
 104     hZPChg[j] = new TH1F(hnamZPChg, hnamZPChg, 100, 0., 1400.);
 105     hZNAhg[j] = new TH1F(hnamZNAhg, hnamZNAhg, 100, 0., 1400.);
 106     hZPAhg[j] = new TH1F(hnamZPAhg, hnamZPAhg, 100, 0., 1400.);
 107     //
 108     sprintf(hnamZNClg,"ZNClg-tow%d",j);
 109     sprintf(hnamZPClg,"ZPClg-tow%d",j);
 110     sprintf(hnamZNAlg,"ZNAlg-tow%d",j);
 111     sprintf(hnamZPAlg,"ZPAlg-tow%d",j);
 112     //
 113     hZNClg[j] = new TH1F(hnamZNClg, hnamZNClg, 100, 0., 4000.);
 114     hZPClg[j] = new TH1F(hnamZPClg, hnamZPClg, 100, 0., 4000.);
 115     hZNAlg[j] = new TH1F(hnamZNAlg, hnamZNAlg, 100, 0., 4000.);
 116     hZPAlg[j] = new TH1F(hnamZPAlg, hnamZPAlg, 100, 0., 4000.);
 117     //
 118     if(j<2){
 119       sprintf(hnamZEMhg,"ZEM%dhg",j);
 120       sprintf(hnamZEMlg,"ZEM%dlg",j);
 121       //
 122       hZEMhg[j] = new TH1F(hnamZEMhg, hnamZEMhg, 100, 0., 1400.);
 123       hZEMlg[j] = new TH1F(hnamZEMlg, hnamZEMlg, 100, 0., 4000.);
 124     }
 125   }
 126
 127   /* open result file */
 128   FILE *fp=NULL;
 129   fp=fopen("./result.txt","a");
 130   if (fp==NULL) {
 131     printf("Failed to open file\n");
 132     return -1;
 133   }
 134   /* report progress */
 135   daqDA_progressReport(10);
 136
 137   // *** To analyze LASER events you MUST have a pedestal data file!!!
 138   // *** -> check if a pedestal run has been analyzed
 139   int read = 0;
 140   read = daqDA_DB_getFile(PEDDATA_FILE, PEDDATA_FILE);
 141   if(read){
 142     printf("\t ERROR!!! ZDCPedestal.dat file NOT FOUND in DAQ db!!!\n");
 143     return -1;
 144   }
 145   else printf("\t ZDCPedestal.dat file retrieved from DAQ db\n");
 146
 147   FILE *filePed = fopen(PEDDATA_FILE,"r");
 148   if (filePed==NULL) {
 149     printf("\t ERROR!!! Can't open ZDCPedestal.dat file!!!\n");
 150     return -1;
 151   }
 152
 153   // 144 = 48 in-time + 48 out-of-time + 48 correlations
 154   Float_t readValues[2][6*kNChannels];
 155   Float_t MeanPedhg[kNChannels], MeanPedlg[kNChannels];
 156   Float_t CorrCoeff0[2*kNChannels], CorrCoeff1[2*kNChannels];
 157   // ***************************************************
 158   //   Unless we have a narrow correlation to fit we
 159   //    don't fit and store in-time vs. out-of-time
 160   //    histograms -> mean pedstal subtracted!!!!!!
 161   // ***************************************************
 162   //
 163   for(int jj=0; jj<6*kNChannels; jj++){
 164     for(int ii=0; ii<2; ii++){
 165        fscanf(filePed,"%f",&readValues[ii][jj]);
 166     }
 167     if(jj<kNChannels){
 168       MeanPedhg[jj] = readValues[0][jj];
 169       //printf("\t MeanPedhg[%d] = %1.1f\n",jj, MeanPedhg[jj]);
 170     }
 171     else if(jj>=kNChannels && jj<2*kNChannels){
 172       MeanPedlg[jj-kNChannels] = readValues[0][jj];
 173       //printf("\t MeanPedlg[%d] = %1.1f\n",jj-kNChannels, MeanPedlg[jj-kNChannels]);
 174     }
 175     else if(jj>4*kNChannels){
 176       CorrCoeff0[jj-4*kNChannels] = readValues[0][jj];
 177       CorrCoeff1[jj-4*kNChannels] = readValues[1][jj];;
 178     }
 179   }
 180
 181   FILE *mapFile4Shuttle;
 182
 183   /* report progress */
 184   daqDA_progressReport(20);
 185
 186
 187   /* init some counters */
 188   int nevents_physics=0;
 189   int nevents_total=0;
 190
 191   /* read the data files */
 192   int n;
 193   for(n=1;n<argc;n++) {
 194
 195     status=monitorSetDataSource( argv[n] );
 196     if (status!=0) {
 197       printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
 198       return -1;
 199     }
 200
 201     /* report progress */
 202     /* in this example, indexed on the number of files */
 203     daqDA_progressReport(20+70*n/argc);
 204
 205     /* read the file */
 206     for(;;) {
 207       struct eventHeaderStruct *event;
 208       eventTypeType eventT;
 209
 210       /* get next event */
 211       status=monitorGetEventDynamic((void **)&event);
 212       if (status==MON_ERR_EOF) break; /* end of monitoring file has been reached */
 213       if (status!=0) {
 214         printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
 215         return -1;
 216       }
 217
 218       /* retry if got no event */
 219       if (event==NULL) {
 220         break;
 221       }
 222
 223       // Initalize raw-data reading and decoding
 224       AliRawReader *reader = new AliRawReaderDate((void*)event);
 225       reader->Select("ZDC");
 226       // --- Reading event header
 227       //UInt_t evtype = reader->GetType();
 228       //printf("\n\t ZDCLASERda -> ev. type %d\n",evtype);
 229       //printf("\t ZDCLASERda -> run # %d\n",reader->GetRunNumber());
 230       //
 231       AliZDCRawStream *rawStreamZDC = new AliZDCRawStream(reader);
 232
 233
 234       /* use event - here, just write event id to result file */
 235       eventT=event->eventType;
 236
 237       Int_t ich=0;
 238       Int_t adcMod[2*kNChannels], adcCh[2*kNChannels];
 239       Int_t sigCode[2*kNChannels], det[2*kNChannels], sec[2*kNChannels];
 240
 241       if(eventT==START_OF_DATA){
 242
 243         rawStreamZDC->SetSODReading(kTRUE);
 244
 245         // --------------------------------------------------------
 246         // --- Writing ascii data file for the Shuttle preprocessor
 247         mapFile4Shuttle = fopen(MAPDATA_FILE,"w");
 248         if(!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
 249         else{
 250           while((rawStreamZDC->Next()) && (ich<2*kNChannels)){
 251             if(rawStreamZDC->IsChMapping()){
 252               adcMod[ich] = rawStreamZDC->GetADCModFromMap(ich);
 253               adcCh[ich] = rawStreamZDC->GetADCChFromMap(ich);
 254               sigCode[ich] = rawStreamZDC->GetADCSignFromMap(ich);
 255               det[ich] = rawStreamZDC->GetDetectorFromMap(ich);
 256               sec[ich] = rawStreamZDC->GetTowerFromMap(ich);
 257               //
 258               fprintf(mapFile4Shuttle,"\t%d\t%d\t%d\t%d\t%d\t%d\n",
 259                 ich,adcMod[ich],adcCh[ich],sigCode[ich],det[ich],sec[ich]);
 260               //
 261               //printf("ZDCPEDESTALda.cxx -> %d mod %d ch %d, code %d det %d, sec %d\n",
 262               //   ich,adcMod[ich],adcCh[ich],sigCode[ich],det[ich],sec[ich]);
 263               //
 264               ich++;
 265             }
 266           }
 267         }
 268         fclose(mapFile4Shuttle);
 269       }// SOD event
 270
 271       if(eventT==PHYSICS_EVENT){
 272         // --- Reading data header
 273         reader->ReadHeader();
 274         const AliRawDataHeader* header = reader->GetDataHeader();
 275         if(header) {
 276          UChar_t message = header->GetAttributes();
 277          if(message & 0x30){ // DEDICATED LASER RUN
 278             //printf("\t STANDALONE_LASER_RUN raw data found\n");
 279             continue;
 280          }
 281          else{
 282             printf("ZDCLASERda.cxx -> NO STANDALONE_LASER_RUN raw data found\n");
 283             return -1;
 284          }
 285         }
 286         else{
 287            printf("\t ATTENTION! No Raw Data Header found!!!\n");
 288            return -1;
 289         }
 290
 291         rawStreamZDC->SetSODReading(kTRUE);
 292
 293         if (!rawStreamZDC->Next()) printf(" \t No raw data found!! \n");
 294         //
 295         // ----- Setting ch. mapping -----
 296         for(Int_t jk=0; jk<2*kNChannels; jk++){
 297           rawStreamZDC->SetMapADCMod(jk, adcMod[jk]);
 298           rawStreamZDC->SetMapADCCh(jk, adcCh[jk]);
 299           rawStreamZDC->SetMapADCSig(jk, sigCode[jk]);
 300           rawStreamZDC->SetMapDet(jk, det[jk]);
 301           rawStreamZDC->SetMapTow(jk, sec[jk]);
 302         }
 303         //
 304         while(rawStreamZDC->Next()){
 305           Int_t index=-1;
 306           Int_t detector = rawStreamZDC->GetSector(0);
 307           Int_t sector = rawStreamZDC->GetSector(1);
 308
 309           if(rawStreamZDC->IsADCDataWord() && !(rawStreamZDC->IsUnderflow())
 310              && !(rawStreamZDC->IsOverflow()) && detector!=-1){
 311
 312             //printf("  IsADCWord %d, IsUnderflow %d, IsOverflow %d\n",
 313             //  rawStreamZDC->IsADCDataWord(),rawStreamZDC->IsUnderflow(),rawStreamZDC->IsOverflow());
 314
 315             if(sector!=5){ // Physics signals
 316               if(detector==1) index = sector;        // *** ZNC
 317               else if(detector==2) index = sector+5; // *** ZPC
 318               else if(detector==3) index = sector+9; // *** ZEM
 319               else if(detector==4) index = sector+12;// *** ZNA
 320               else if(detector==5) index = sector+17;// *** ZPA
 321             }
 322             else{ // Reference PMs
 323               index = (detector-1)/3+22;
 324             }
 325             //
 326             if(index==-1) printf("ERROR in ZDCLASERda.cxx -> det %d quad %d res %d index %d ADC %d\n",
 327               detector, sector, rawStreamZDC->GetADCGain(), index, rawStreamZDC->GetADCValue());
 328
 329             Float_t Pedestal=0.;
 330             if(rawStreamZDC->GetADCGain()==0) Pedestal = MeanPedhg[index];
 331             else if(rawStreamZDC->GetADCGain()==1) Pedestal = MeanPedlg[index];
 332             //
 333             Float_t CorrADC = rawStreamZDC->GetADCValue() - Pedestal;
 334             //
 335             //printf("\tdet %d sec %d res %d index %d ped %1.0f ADCcorr %1.0f\n",
 336             //  detector, sector, rawStreamZDC->GetADCGain(), index, Pedestal,CorrADC);
 337
 338             // **** Detector PMs
 339             if(sector!=5){
 340               if(rawStreamZDC->GetADCGain()==0){ // --- High gain chain ---
 341                 // ---- side C
 342                 if(detector==1) hZNChg[sector]->Fill(CorrADC);
 343                 else if(detector==2) hZPChg[sector]->Fill(CorrADC);
 344                 // ---- side A
 345                 else if(detector==4) hZNAhg[sector]->Fill(CorrADC);
 346                 else if(detector==5) hZPAhg[sector]->Fill(CorrADC);
 347                 // ---- ZEM
 348                 else if(detector==3) hZEMhg[sector-1]->Fill(CorrADC);
 349               }
 350               else if(rawStreamZDC->GetADCGain()==1){ // --- Low gain chain ---
 351                 // ---- side C
 352                 if(detector==1) hZNClg[sector]->Fill(CorrADC);
 353                 else if(detector==2) hZPClg[sector]->Fill(CorrADC);
 354                 // ---- side A
 355                 else if(detector==4) hZNAlg[sector]->Fill(CorrADC);
 356                 else if(detector==5) hZPAlg[sector]->Fill(CorrADC);
 357                 // ---- ZEM
 358                 else if(detector==3) hZEMlg[sector-1]->Fill(CorrADC);
 359               }
 360             }
 361             // **** Reference PMs
 362             else if(sector==5){
 363               if(rawStreamZDC->GetADCGain()==0){ // --- High gain chain ---
 364                 // ---- PMRef chain side C
 365                 if(detector==1) hPMRefChg->Fill(CorrADC);
 366                 // ---- PMRef side A
 367                 else if(detector==4) hPMRefAhg->Fill(CorrADC);
 368               }
 369               else if(rawStreamZDC->GetADCGain()==1){ // --- Low gain chain ---
 370                 // ---- PMRef chain side C
 371                 if(detector==1) hPMRefClg->Fill(CorrADC);
 372                 // ---- PMRef side A
 373                 else if(detector==4) hPMRefAlg->Fill(CorrADC);
 374               }
 375             }
 376           }//IsADCDataWord()+NOunderflow+NOoverflow
 377           //
 378          }
 379          //
 380          nevents_physics++;
 381          //
 382          delete reader;
 383          delete rawStreamZDC;
 384
 385       }//(if PHYSICS_EVENT)
 386       nevents_total++;
 387
 388       /* free resources */
 389       free(event);
 390
 391     }
 392   }
 393
 394   /* Analysis of the histograms */
 395   //
 396   Int_t det[2*kNChannels], quad[2*kNChannels];
 397   Int_t maxBin[2*kNChannels], nBin[2*kNChannels];
 398   Float_t xMax[2*kNChannels], maxXval[2*kNChannels], xlow[2*kNChannels];
 399   Float_t mean[2*kNChannels], sigma[2*kNChannels];
 400   TF1 *fun[2*kNChannels];
 401
 402   // ******** High gain chain ********
 403   for(Int_t k=0; k<5; k++){
 404     // --- ZNC
 405     det[k] = 1;
 406     quad[k] = k;
 407     maxBin[k] = hZNChg[k]->GetMaximumBin();
 408     nBin[k] = (hZNChg[k]->GetXaxis())->GetNbins();
 409     xMax[k] = (hZNChg[k]->GetXaxis())->GetXmax();
 410     if(nBin[k]!=0) maxXval[k] = maxBin[k]*xMax[k]/nBin[k];
 411     if(maxXval[k]-150.<0.) xlow[k]=0.;
 412     else xlow[k] = maxXval[k]-150.;
 413     // checking if histos are empty
 414     if(hZNChg[k]->GetEntries() == 0){
 415       printf("\n WARNING! Empty LASER histos -> ending DA WITHOUT writing output\n\n");
 416       return -1;
 417     }
 418     //
 419     hZNChg[k]->Fit("gaus","Q","",xlow[k],maxXval[k]+150.);
 420     fun[k] = hZNChg[k]->GetFunction("gaus");
 421     mean[k]  = (Float_t) (fun[k]->GetParameter(1));
 422     sigma[k] = (Float_t) (fun[k]->GetParameter(2));
 423     // --- ZPC
 424     det[k+5] = 2;
 425     quad[k+5] = k;
 426     maxBin[k+5] = hZPChg[k]->GetMaximumBin();
 427     nBin[k+5] = (hZPChg[k]->GetXaxis())->GetNbins();
 428     xMax[k+5] = (hZPChg[k]->GetXaxis())->GetXmax();
 429     if(nBin[k+5]!=0) maxXval[k+5] = maxBin[k+5]*xMax[k+5]/nBin[k+5];
 430     if(maxXval[k+5]-150.<0.) xlow[k+5]=0.;
 431     else xlow[k+5] = maxXval[k+5]-150.;
 432     hZPChg[k]->Fit("gaus","Q","",xlow[k+5],maxXval[k+5]+150.);
 433     fun[k+5] = hZPChg[k]->GetFunction("gaus");
 434     mean[k+5]  = (Float_t) (fun[k+5]->GetParameter(1));
 435     sigma[k+5] = (Float_t) (fun[k+5]->GetParameter(2));
 436     // --- ZEM1
 437     if(k<2){
 438       det[k+10] = 3;
 439       quad[k+10] = k+1;
 440       maxBin[k+10] = hZEMhg[k]->GetMaximumBin();
 441       nBin[k+10] = (hZEMhg[k]->GetXaxis())->GetNbins();
 442       xMax[k+10] = (hZEMhg[k]->GetXaxis())->GetXmax();
 443       if(nBin[k+10]!=0) maxXval[k+10] = maxBin[k+10]*xMax[k+10]/nBin[k+10];
 444       if(maxXval[k+10]-150.<0.) xlow[k+10]=0.;
 445       else xlow[k+10] = maxXval[k+10]-150.;
 446       hZEMhg[k]->Fit("gaus","Q","",xlow[k+10],maxXval[k+10]+150.);
 447       fun[k+10] = hZEMhg[k]->GetFunction("gaus");
 448       mean[k+10]  = (Float_t) (fun[k+10]->GetParameter(1));
 449       sigma[k+10] = (Float_t) (fun[k+10]->GetParameter(2));
 450     }
 451     // --- ZNA
 452     det[k+12] = 4;
 453     quad[k+12] = k;
 454     maxBin[k+12] = hZNAhg[k]->GetMaximumBin();
 455     nBin[k+12] = (hZNAhg[k]->GetXaxis())->GetNbins();
 456     xMax[k+12] = (hZNAhg[k]->GetXaxis())->GetXmax();
 457     if(nBin[k+12]!=0) maxXval[k+12] = maxBin[k+12]*xMax[k+12]/nBin[k+12];
 458     if(maxXval[k+12]-150.<0.) xlow[k+12]=0.;
 459     else xlow[k+12] = maxXval[k+12]-150.;
 460     hZNAhg[k]->Fit("gaus","Q","",xlow[k+12],maxXval[k+12]+150.);
 461     fun[k+12] = hZNAhg[k]->GetFunction("gaus");
 462     mean[k+12]  = (Float_t) (fun[k+12]->GetParameter(1));
 463     sigma[k+12] = (Float_t) (fun[k+12]->GetParameter(2));
 464     // --- ZPA
 465     det[k+17] = 4;
 466     quad[k+17] = 5;
 467     maxBin[k+17] = hZPAhg[k]->GetMaximumBin();
 468     nBin[k+17] = (hZPAhg[k]->GetXaxis())->GetNbins();
 469     xMax[k+17] = (hZPAhg[k]->GetXaxis())->GetXmax();
 470     if(nBin[k+17]!=0) maxXval[k+17] = maxBin[k+17]*xMax[k+17]/nBin[k+17];
 471     if(maxXval[k+17]-150.<0.) xlow[k+17]=0.;
 472     else xlow[k+17] = maxXval[k+17]-150.;
 473     hZPAhg[k]->Fit("gaus","Q","",xlow[k+17],maxXval[k+17]+150.);
 474     fun[k+17] = hZPAhg[k]->GetFunction("gaus");
 475     mean[k+17]  = (Float_t) (fun[k+17]->GetParameter(1));
 476     sigma[k+17] = (Float_t) (fun[k+17]->GetParameter(2));
 477   }
 478   // ~~~~~~~~ PM Ref side C ~~~~~~~~
 479   det[22] = 1;
 480   quad[22] = 5;
 481   maxBin[22] = hPMRefChg->GetMaximumBin();
 482   nBin[22] = (hPMRefChg->GetXaxis())->GetNbins();
 483   xMax[22] = (hPMRefChg->GetXaxis())->GetXmax();
 484   if(nBin[22]!=0) maxXval[22] = maxBin[22]*xMax[22]/nBin[22];
 485   if(maxXval[22]-150.<0.) xlow[22]=0.;
 486   else xlow[22] = maxXval[22];
 487   hPMRefChg->Fit("gaus","Q","",xlow[22],maxXval[22]+150.);
 488   fun[22] = hPMRefChg->GetFunction("gaus");
 489   mean[22]  = (Float_t) (fun[22]->GetParameter(1));
 490   sigma[22] = (Float_t) (fun[22]->GetParameter(2));
 491   // ~~~~~~~~ PM Ref side A ~~~~~~~~
 492   det[23] = 4;
 493   quad[23] = 5;
 494   maxBin[23] = hPMRefAhg->GetMaximumBin();
 495   nBin[23] = (hPMRefAhg->GetXaxis())->GetNbins();
 496   xMax[23] = (hPMRefAhg->GetXaxis())->GetXmax();
 497   if(nBin[23]!=0) maxXval[23] = maxBin[23]*xMax[23]/nBin[23];
 498   if(maxXval[23]-100.<0.) xlow[23]=0.;
 499   else xlow[23] = maxXval[23];
 500   hPMRefAhg->Fit("gaus","Q","",xlow[23],maxXval[23]+100.);
 501   fun[23] = hPMRefAhg->GetFunction("gaus");
 502   mean[23]  = (Float_t) (fun[23]->GetParameter(1));
 503   sigma[23] = (Float_t) (fun[23]->GetParameter(2));
 504
 505   // ******** Low gain chain ********
 506   Int_t kOffset = 24;
 507   for(Int_t k=0; k<5; k++){
 508     // --- ZNC
 509     det[k+kOffset] = 1;
 510     quad[k+kOffset] = k;
 511     maxBin[k+kOffset] = hZNClg[k]->GetMaximumBin();
 512     nBin[k+kOffset] = (hZNClg[k]->GetXaxis())->GetNbins();
 513     xMax[k+kOffset] = (hZNClg[k]->GetXaxis())->GetXmax();
 514     if(nBin[k+kOffset]!=0) maxXval[k+kOffset] = maxBin[k+kOffset]*xMax[k+kOffset]/nBin[k+kOffset];
 515     if(maxXval[k+kOffset]-150.<0.) xlow[k+kOffset]=0.;
 516     else xlow[k+kOffset] = maxXval[k+kOffset]-150.;
 517     hZNClg[k]->Fit("gaus","Q","",xlow[k+kOffset],maxXval[k+kOffset]+150.);
 518     fun[k+kOffset] = hZNClg[k]->GetFunction("gaus");
 519     mean[k+kOffset]  = (Float_t) (fun[k+kOffset]->GetParameter(1));
 520     sigma[k+kOffset] = (Float_t) (fun[k+kOffset]->GetParameter(2));
 521     // --- ZPC
 522     det[k+kOffset+5] = 2;
 523     quad[k+kOffset+5] = k;
 524     maxBin[k+kOffset+5] = hZPClg[k]->GetMaximumBin();
 525     nBin[k+kOffset+5] = (hZPClg[k]->GetXaxis())->GetNbins();
 526     xMax[k+kOffset+5] = (hZPClg[k]->GetXaxis())->GetXmax();
 527     if(nBin[k+kOffset+5]!=0) maxXval[k+kOffset+5] = maxBin[k+kOffset+5]*xMax[k+kOffset+5]/nBin[k+kOffset+5];
 528     if(maxXval[k+kOffset+5]-150.<0.) xlow[k+kOffset+5]=0.;
 529     else xlow[k+kOffset+5] = maxXval[k+kOffset+5]-150.;
 530     hZPClg[k]->Fit("gaus","Q","",xlow[k+kOffset+5],maxXval[k+kOffset+5]+150.);
 531     fun[k+kOffset+5] = hZPClg[k]->GetFunction("gaus");
 532     mean[k+kOffset+5]  = (Float_t) (fun[k+kOffset+5]->GetParameter(1));
 533     sigma[k+kOffset+5] = (Float_t) (fun[k+kOffset+5]->GetParameter(2));
 534     // --- ZEM1
 535     if(k+kOffset<2){
 536       det[k+kOffset+10] = 3;
 537       quad[k+kOffset+10] = k+1;
 538       maxBin[k+kOffset+10] = hZEMlg[k]->GetMaximumBin();
 539       nBin[k+kOffset+10] = (hZEMlg[k]->GetXaxis())->GetNbins();
 540       xMax[k+kOffset+10] = (hZEMlg[k]->GetXaxis())->GetXmax();
 541       if(nBin[k+kOffset+10]!=0) maxXval[k+kOffset+10] = maxBin[k+kOffset+10]*xMax[k+kOffset+10]/nBin[k+kOffset+10];
 542       if(maxXval[k+kOffset+10]-150.<0.) xlow[k+kOffset+10]=0.;
 543       else xlow[k+kOffset+10] = maxXval[k+kOffset+10]-150.;
 544       hZEMlg[k]->Fit("gaus","Q","",xlow[k+kOffset+10],maxXval[k+kOffset+10]+150.);
 545       fun[k+kOffset+10] = hZEMlg[k]->GetFunction("gaus");
 546       mean[k+kOffset+10]  = (Float_t) (fun[k+kOffset+10]->GetParameter(1));
 547       sigma[k+kOffset+10] = (Float_t) (fun[k+kOffset+10]->GetParameter(2));
 548     }
 549     // --- ZNA
 550     det[k+kOffset+12] = 4;
 551     quad[k+kOffset+12] = k;
 552     maxBin[k+kOffset+12] = hZNAlg[k]->GetMaximumBin();
 553     nBin[k+kOffset+12] = (hZNAlg[k]->GetXaxis())->GetNbins();
 554     xMax[k+kOffset+12] = (hZNAlg[k]->GetXaxis())->GetXmax();
 555     if(nBin[k+kOffset+12]!=0) maxXval[k+kOffset+12] = maxBin[k+kOffset+12]*xMax[k+kOffset+12]/nBin[k+kOffset+12];
 556     if(maxXval[k+kOffset+12]-150.<0.) xlow[k+kOffset+12]=0.;
 557     else xlow[k+kOffset+12] = maxXval[k+kOffset+12]-150.;
 558     hZNAlg[k]->Fit("gaus","Q","",xlow[k+kOffset+12],maxXval[k+kOffset+12]+150.);
 559     fun[k+kOffset+12] = hZNAlg[k]->GetFunction("gaus");
 560     mean[k+kOffset+12]  = (Float_t) (fun[k+kOffset+12]->GetParameter(1));
 561     sigma[k+kOffset+12] = (Float_t) (fun[k+kOffset+12]->GetParameter(2));
 562     // --- ZPA
 563     det[k+kOffset+17] = 5;
 564     quad[k+kOffset+17] = k;
 565     maxBin[k+kOffset+17] = hZPAlg[k]->GetMaximumBin();
 566     nBin[k+kOffset+17] = (hZPAlg[k]->GetXaxis())->GetNbins();
 567     xMax[k+kOffset+17] = (hZPAlg[k]->GetXaxis())->GetXmax();
 568     if(nBin[k+kOffset+17]!=0) maxXval[k+kOffset+17] = maxBin[k+kOffset+17]*xMax[k+kOffset+17]/nBin[k+kOffset+17];
 569     if(maxXval[k+kOffset+17]-150.<0.) xlow[k+kOffset+17]=0.;
 570     else xlow[k+kOffset+17] = maxXval[k+kOffset+17]-150.;
 571     hZPAlg[k]->Fit("gaus","Q","",xlow[k+kOffset+17],maxXval[k+kOffset+17]+150.);
 572     fun[k+kOffset+17] = hZPAlg[k]->GetFunction("gaus");
 573     mean[k+kOffset+17]  = (Float_t) (fun[k+kOffset+17]->GetParameter(1));
 574     sigma[k+kOffset+17] = (Float_t) (fun[k+kOffset+17]->GetParameter(2));
 575   }
 576   // ~~~~~~~~ PM Ref side C ~~~~~~~~
 577   det[46] = 1;
 578   quad[46] = 5;
 579   maxBin[46] = hPMRefClg->GetMaximumBin();
 580   nBin[46] = (hPMRefClg->GetXaxis())->GetNbins();
 581   xMax[46] = (hPMRefClg->GetXaxis())->GetXmax();
 582   if(nBin[46]!=0) maxXval[46] = maxBin[46]*xMax[46]/nBin[46];
 583   if(maxXval[46]-150.<0.) xlow[46]=0.;
 584   else xlow[46] = maxXval[46];
 585   hPMRefClg->Fit("gaus","Q","",xlow[46],maxXval[46]+150.);
 586   fun[46] = hPMRefClg->GetFunction("gaus");
 587   mean[46]  = (Float_t) (fun[46]->GetParameter(1));
 588   sigma[46] = (Float_t) (fun[46]->GetParameter(2));
 589   // ~~~~~~~~ PM Ref side A ~~~~~~~~
 590   det[47] = 4;
 591   quad[47] = 5;
 592   maxBin[47] = hPMRefAlg->GetMaximumBin();
 593   nBin[47] = (hPMRefAlg->GetXaxis())->GetNbins();
 594   xMax[47] = (hPMRefAlg->GetXaxis())->GetXmax();
 595   if(nBin[47]!=0) maxXval[47] = maxBin[47]*xMax[47]/nBin[47];
 596   if(maxXval[47]-100.<0.) xlow[47]=0.;
 597   else xlow[47] = maxXval[47];
 598   hPMRefAlg->Fit("gaus","Q","",xlow[47],maxXval[47]+100.);
 599   fun[47] = hPMRefAlg->GetFunction("gaus");
 600   mean[47]  = (Float_t) (fun[47]->GetParameter(1));
 601   sigma[47] = (Float_t) (fun[47]->GetParameter(2));
 602
 603   FILE *fileShuttle;
 604   fileShuttle = fopen(LASDATA_FILE,"w");
 605   for(Int_t i=0; i<2*kNChannels; i++){
 606     fprintf(fileShuttle,"\t%d\t%d\t%f\t%f\n",det[i],quad[i],mean[i], sigma[i]);
 607   }
 608   //
 609   fclose(fileShuttle);
 610   /* report progress */
 611   daqDA_progressReport(80);
 612   //
 613   TFile *histofile = new TFile(LASHISTO_FILE,"RECREATE");
 614   histofile->cd();
 615   for(int j=0; j<5; j++){
 616      hZNChg[j]->Write();
 617      hZPChg[j]->Write();
 618      hZNAhg[j]->Write();
 619      hZPAhg[j]->Write();
 620      hZNClg[j]->Write();
 621      hZPClg[j]->Write();
 622      hZNAlg[j]->Write();
 623      hZPAlg[j]->Write();
 624      if(j<2){
 625        hZEMhg[j]->Write();
 626        hZEMlg[j]->Write();
 627     }
 628   }
 629   hPMRefChg->Write();
 630   hPMRefAhg->Write();
 631   hPMRefClg->Write();
 632   hPMRefAlg->Write();
 633   //
 634   histofile->Close();
 635   //
 636   for(Int_t j=0; j<5; j++){
 637     delete hZNChg[j];
 638     delete hZPChg[j];
 639     delete hZNAhg[j];
 640     delete hZPAhg[j];
 641     delete hZNClg[j];
 642     delete hZPClg[j];
 643     delete hZNAlg[j];
 644     delete hZPAlg[j];
 645     if(j<2){
 646       delete hZEMhg[j];
 647       delete hZEMlg[j];
 648     }
 649   }
 650   delete hPMRefChg;
 651   delete hPMRefAhg;
 652   delete hPMRefClg;
 653   delete hPMRefAlg;
 654
 655   /* write report */
 656   fprintf(fp,"Run #%s, received %d physics events out of %d\n",getenv("DATE_RUN_NUMBER"),nevents_physics,nevents_total);
 657
 658   /* close result file */
 659   fclose(fp);
 660
 661   /* report progress */
 662   daqDA_progressReport(90);
 663
 664   /* store the result file on FES */
 665   // [1] File with mapping
 666   status = daqDA_FES_storeFile(MAPDATA_FILE,MAPDATA_FILE);
 667   if(status){
 668     printf("Failed to export file : %d\n",status);
 669     return -1;
 670   }
 671   //
 672   // [2] File with laser data
 673   status = daqDA_FES_storeFile(LASDATA_FILE,LASDATA_FILE);
 674   if(status){
 675     printf("Failed to export file : %d\n",status);
 676     return -1;
 677   }
 678   // [3] File with laser histos
 679   status = daqDA_FES_storeFile(LASHISTO_FILE,LASHISTO_FILE);
 680   if(status){
 681     printf("Failed to export pedestal histos file to DAQ FES\n");
 682     return -1;
 683   }
 684
 685   /* report progress */
 686   daqDA_progressReport(100);
 687
 688   return status;
 689 }