MUON/MUONTRKda.cxx

   1 /*
   2
   3 Version 1 for MUONTRKda MUON tracking
   4 Working version for pedestal
   5 Framework for gain computing
   6 (Ch. Finck)
   7
   8
   9 Rem:  AliMUON2DMap stores all channels, even those which are not connected
  10       if pedMean == -1, channel not connected to a pad
  11
  12 Pb with static root libraries ??
  13
  14
  15 */
  16 extern "C" {
  17 #include <daqDA.h>
  18 }
  19
  20 #include "event.h"
  21 #include "monitor.h"
  22
  23 #include <Riostream.h>
  24 #include <stdio.h>
  25 #include <stdlib.h>
  26
  27 //AliRoot
  28 #include "AliMUONRawStreamTracker.h"
  29 #include "AliMUONDspHeader.h"
  30 #include "AliMUONBlockHeader.h"
  31 #include "AliMUONBusStruct.h"
  32 #include "AliMUONDDLTracker.h"
  33 #include "AliMUONV2DStore.h"
  34 #include "AliMUON2DMap.h"
  35 #include "AliMUONCalibParamNF.h"
  36 #include "AliMUONObjectPair.h"
  37 #include "AliMUONVDataIterator.h"
  38 #include "AliMpDDLStore.h"
  39 #include "AliMpIntPair.h"
  40
  41 #include "AliRawReaderDate.h"
  42
  43
  44 //ROOT
  45 #include "TFile.h"
  46 #include "TTree.h"
  47 #include "TH1F.h"
  48 #include "TString.h"
  49 #include "TStopwatch.h"
  50 #include "TMath.h"
  51 #include "TTimeStamp.h"
  52 #include "TGraphErrors.h"
  53 #include "TF1.h"
  54
  55
  56 // global variables
  57 AliMUONV2DStore* pedestalStore =  new AliMUON2DMap(kFALSE);
  58 const Int_t kNchannels = 64;
  59 Int_t nManu = 0;
  60 Int_t nChannel = 0;
  61 UInt_t runNumber = 0;
  62 Int_t nEvents = 0;
  63 TH1F* pedMeanHisto = 0x0;
  64 TH1F* pedSigmaHisto = 0x0;
  65 Char_t histoFileName[256];
  66 TString command("ped");
  67
  68 // funtions
  69
  70 //________________
  71 Double_t fitFunc(Double_t *x, Double_t *par)
  72 {
  73     //fit function for gains
  74
  75     Float_t xx = x[0];
  76     Double_t f;
  77
  78     if (xx < par[3])
  79         f = par[0] + par[1]*xx;
  80     else
  81         f= par[0] + par[1]*xx + par[2]*xx*xx;
  82
  83     return f;
  84 }
  85
  86
  87
  88 //__________
  89 void MakePed(Int_t busPatchId, Int_t manuId, Int_t channelId, Int_t charge)
  90 {
  91
  92     AliMUONVCalibParam* ped =
  93         static_cast<AliMUONVCalibParam*>(pedestalStore->Get(busPatchId, manuId));
  94
  95     if (!ped) {
  96       nManu++;
  97       ped = new AliMUONCalibParamNF(2, kNchannels, -1.); // put default wise -1, not connected channel
  98       pedestalStore->Set(busPatchId, manuId, ped, kFALSE);
  99     }
 100
 101     if (nEvents == 1) {
 102       ped->SetValueAsFloat(channelId, 0, 0.);
 103       ped->SetValueAsFloat(channelId, 1, 0.);
 104     }
 105
 106     Float_t pedMean  = ped->ValueAsFloat(channelId, 0) + charge;
 107     Float_t pedSigma = ped->ValueAsFloat(channelId, 1) + charge*charge;
 108
 109     ped->SetValueAsFloat(channelId, 0, pedMean);
 110     ped->SetValueAsFloat(channelId, 1, pedSigma);
 111
 112 }
 113
 114 //________________
 115 void MakePedStore(TString flatOutputFile = "")
 116 {
 117   TTimeStamp date;
 118   Float_t pedMean;
 119   Float_t pedSigma;
 120   ofstream fileout;
 121   Int_t busPatchId;
 122   Int_t manuId;
 123   Int_t channelId;
 124
 125  // histo
 126   sprintf(histoFileName,"mutrkped-%d.root",runNumber);
 127   TFile*  histoFile = new TFile(histoFileName,"RECREATE","MUON Tracking pedestals");
 128
 129   Char_t name[255];
 130   Char_t title[255];
 131   sprintf(name,"pedmean_allch");
 132   sprintf(title,"Pedestal mean all channels");
 133   Int_t nx = 1000;
 134   Int_t xmin = 0;
 135   Int_t xmax = 1000;
 136   pedMeanHisto = new TH1F(name,title,nx,xmin,xmax);
 137   pedMeanHisto->SetDirectory(histoFile);
 138
 139   sprintf(name,"pedsigma_allch");
 140   sprintf(title,"Pedestal sigma all channels");
 141   nx = 200;
 142   xmin = 0;
 143   xmax = 50;
 144   pedSigmaHisto = new TH1F(name,title,nx,xmin,xmax);
 145   pedSigmaHisto->SetDirectory(histoFile);
 146
 147   TTree* tree = new TTree("t","Pedestal tree");
 148   tree->Branch("bp",&busPatchId,"bp/I");
 149   tree->Branch("manu",&manuId,",manu/I");
 150   tree->Branch("channel",&channelId,",channel/I");
 151
 152   if (!flatOutputFile.IsNull()) {
 153     fileout.open(flatOutputFile.Data());
 154     fileout<<"//===========================================================================" << endl;
 155     fileout<<"//                       Pedestal file calculated by MUONTRKda"<<endl;
 156     fileout<<"//===========================================================================" << endl;
 157     fileout<<"//       * Run           : " << runNumber << endl;
 158     fileout<<"//       * Date          : " << date.AsString("l") <<endl;
 159     fileout<<"//       * Statictics    : " << nEvents << endl;
 160     fileout<<"//       * # of MANUS    : " << nManu << endl;
 161     fileout<<"//       * # of channels : " << nChannel << endl;
 162     fileout<<"//"<<endl;
 163     fileout<<"//---------------------------------------------------------------------------" << endl;
 164     fileout<<"//---------------------------------------------------------------------------" << endl;
 165     fileout<<"//format : BUS_PATCH MANU_ID CHANNEL MEAN SIGMA"<<endl;
 166     fileout<<"//---------------------------------------------------------------------------" << endl;
 167
 168   }
 169
 170   // iterator over pedestal
 171   AliMUONVDataIterator* it = pedestalStore->Iterator();
 172
 173   AliMUONObjectPair* p;
 174
 175   while ( ( p = dynamic_cast<AliMUONObjectPair*>(it->Next() ) ) )
 176   {
 177     AliMUONVCalibParam* ped = dynamic_cast<AliMUONVCalibParam*>(p->Value());
 178     AliMpIntPair* pair      = dynamic_cast<AliMpIntPair*>(p->Key());
 179     busPatchId              = pair->GetFirst();
 180     manuId                  = pair->GetSecond();
 181
 182     for (channelId = 0; channelId < ped->Size() ; ++channelId) {
 183
 184       pedMean  = ped->ValueAsFloat(channelId, 0);
 185
 186       if (pedMean > 0) { // connected channels
 187
 188         ped->SetValueAsFloat(channelId, 0, pedMean/(Float_t)nEvents);
 189
 190         pedMean  = ped->ValueAsFloat(channelId, 0);
 191         pedSigma = ped->ValueAsFloat(channelId, 1);
 192
 193         ped->SetValueAsFloat(channelId, 1, TMath::Sqrt(TMath::Abs(pedSigma/(Float_t)nEvents - pedMean*pedMean)));
 194
 195         pedMean  = ped->ValueAsFloat(channelId, 0) + 0.5 ;
 196         pedSigma = ped->ValueAsFloat(channelId, 1);
 197
 198
 199         if (!flatOutputFile.IsNull()) {
 200           fileout << "\t" << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t"
 201                   << pedMean <<"\t"<< pedSigma << endl;
 202         }
 203
 204         pedMeanHisto->Fill(pedMean);
 205         pedSigmaHisto->Fill(pedSigma);
 206
 207         tree->Fill();
 208       }
 209     }
 210
 211     delete p; // cos create a new object for each iteration
 212   }
 213
 214   // file outputs
 215   if (!flatOutputFile.IsNull())
 216       fileout.close();
 217
 218   histoFile->Write();
 219   histoFile->Close();
 220
 221 //  delete tree;
 222
 223
 224 // not usable need an update of DATE ->5.28, but it compiles
 225 // uncomment when DATE version ok.
 226 // setenv DATE_FES_PATH
 227 // setenv DATE_RUN_NUMBER
 228 // setenv DATE_ROLE_NAME
 229 // setenv DATE_DETECTOR_CODE
 230
 231 //   if (!flatOutputFile.IsNull()) {
 232
 233 //     flatOutputFile.Prepend("./");
 234
 235 //     status = daqDA_FES_storeFile(flatOutputFile.Data(),"MUONTRKda_results");
 236 //     if (status) {
 237 //       printf("Failed to export file : %d\n",status);
 238 //       return -1;
 239 //     }
 240 //  }
 241
 242 }
 243
 244 //________________
 245 void MakePedStoreForGain()
 246 {
 247     // store pedestal map in root file
 248
 249     Int_t injCharge = 200;
 250
 251     TTree* tree = 0x0;
 252
 253     // compute and store pedestals
 254     MakePedStore();
 255
 256     // store in root file
 257     sprintf(histoFileName,"mutrkgain.root");
 258
 259     TString mode("UPDATE");
 260
 261     if (command.Contains("cre")) {
 262         mode = "RECREATE";
 263     }
 264     TFile* histoFile = new TFile(histoFileName, mode.Data(), "MUON Tracking Gains");
 265
 266     // second argument should be the injected charge, taken from config crocus file
 267     // put also info about run number could be usefull
 268     AliMpIntPair* pair   = new AliMpIntPair(runNumber, injCharge);
 269
 270     if (mode.CompareTo("UPDATE") == 0) {
 271       tree = (TTree*)histoFile->Get("t");
 272       tree->SetBranchAddress("run",&pair);
 273       tree->SetBranchAddress("ped",&pedestalStore);
 274
 275     } else {
 276       tree = new TTree("t","Pedestal tree");
 277       tree->Branch("run", "AliMpIntPair",&pair);
 278       tree->Branch("ped", "AliMUON2DMap",&pedestalStore);
 279       tree->SetBranchAddress("run",&pair);
 280       tree->SetBranchAddress("ped",&pedestalStore);
 281
 282     }
 283
 284     tree->Fill();
 285     tree->Write("t", TObject::kOverwrite); // overwrite the tree
 286     histoFile->Close();
 287
 288     delete pair;
 289 }
 290
 291 //________________
 292 void MakeGainStore(TString flatOutputFile)
 293 {
 294
 295     // open file mutrkgain.root
 296     // read again the pedestal for the calibration runs (9 runs ?)
 297     // need the injection charge from config file (to be done)
 298     // For each channel make a TGraphErrors (mean, sigma) vs injected charge
 299     // Fit with a polynomial fct
 300     // store the result in a flat file.
 301
 302     Float_t pedMean[10];
 303     Float_t pedSigma[10];
 304     Float_t injCharge[10];
 305     Float_t injChargeErr[10];
 306
 307     ofstream fileout;
 308     TTimeStamp date;
 309
 310     if (!flatOutputFile.IsNull()) {
 311       fileout.open(flatOutputFile.Data());
 312       fileout<<"//===========================================================================" << endl;
 313       fileout<<"//                       Pedestal file calculated by MUONTRKda"<<endl;
 314       fileout<<"//===========================================================================" << endl;
 315       fileout<<"//       * Run           : " << runNumber << endl;
 316       fileout<<"//       * Date          : " << date.AsString("l") <<endl;
 317       fileout<<"//       * Statictics    : " << nEvents << endl;
 318       fileout<<"//       * # of MANUS    : " << nManu << endl;
 319       fileout<<"//       * # of channels : " << nChannel << endl;
 320       fileout<<"//"<<endl;
 321       fileout<<"//---------------------------------------------------------------------------" << endl;
 322       fileout<<"//---------------------------------------------------------------------------" << endl;
 323       fileout<<"//format : BUS_PATCH MANU_ID CHANNEL GAIN0 GAIN1 GAIN2 XLIM CHI2"<<endl;
 324       fileout<<"//---------------------------------------------------------------------------" << endl;
 325
 326     }
 327
 328
 329     sprintf(histoFileName,"mutrkgain.root");
 330     TFile*  histoFile = new TFile(histoFileName);
 331
 332     AliMUON2DMap* map[10];
 333     AliMUONVCalibParam* ped[10];
 334     AliMpIntPair* pair;
 335     AliMpIntPair* run[10];
 336
 337     //read back from root file
 338     TTree* tree = (TTree*)histoFile->Get("t");
 339     Int_t nEntries = tree->GetEntries();
 340
 341     //pb with static libraries of Root to read back info
 342     for (Int_t i = 0; i < nEntries; ++i) {
 343       printf("nEntry %d\n", i);
 344       tree->SetBranchAddress("ped",&map[i]);
 345       tree->SetBranchAddress("run",&run[i]);
 346       tree->GetEvent(i);
 347     }
 348
 349     // Q = f(ADC)
 350     TF1* func = new TF1("func",fitFunc, 0., 2500., 4);
 351
 352     // iterates over the first pedestal run
 353     AliMUONVDataIterator* it = map[0]->Iterator();
 354
 355     AliMUONObjectPair* p;
 356
 357     while ( ( p = dynamic_cast<AliMUONObjectPair*>(it->Next() ) ) )
 358     {
 359       ped[0]  = dynamic_cast<AliMUONVCalibParam*>(p->Value());
 360       pair    = dynamic_cast<AliMpIntPair*>(p->Key());
 361
 362       Int_t busPatchId = pair->GetFirst();
 363       Int_t manuId     = pair->GetSecond();
 364
 365       // read back pedestal from the other runs for the given (bupatch, manu)
 366       printf("nEntries %d\n", nEntries);
 367       for (Int_t i = 1; i < nEntries; ++i) {
 368         ped[i] = static_cast<AliMUONVCalibParam*>(map[i]->Get(busPatchId, manuId));
 369       }
 370
 371       // compute for each channel the gain parameters
 372       for ( Int_t channelId = 0; channelId < ped[0]->Size() ; ++channelId ) {
 373
 374         Int_t n = 0;
 375         for (Int_t i = 0; i < nEntries; ++i) {
 376
 377           if (!ped[i]) continue; //shouldn't happen.
 378           pedMean[i]      = ped[i]->ValueAsFloat(channelId, 0);
 379           pedSigma[i]     = ped[i]->ValueAsFloat(channelId, 1);
 380           injCharge[i]    = (Float_t)run[i]->GetSecond();
 381           injChargeErr[i] = 1.;
 382
 383           if (pedMean[i] < 0) n--; // not connected
 384
 385           if (pedSigma[i] <= 0) pedSigma[i] = 1.; // should not happen.
 386           n++;
 387         }
 388
 389         if (n > 4) {
 390           //fit if you can.... not working with present static version of Root
 391           TGraph *gain = new TGraphErrors(n, pedMean, pedSigma, injCharge, injChargeErr);
 392           //should set some initial parameters
 393           func->SetParameter(0,-300);  // a0
 394           func->SetParameter(1, 1.);   // a1
 395           func->SetParameter(2, 0.00001);// a2
 396           func->SetParameter(3, 1100.); // xlim in ADC
 397
 398           gain->Fit("func","q");
 399
 400           if (!flatOutputFile.IsNull()) {
 401             fileout << "\t" << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t"
 402                     << setw(8) << func->GetParameter(0)  <<"\t"<< setw(8) << func->GetParameter(1)
 403                     << setw(8) << func->GetParameter(2)  <<"\t"<< setw(5) << func->GetParameter(3)
 404                     << "\t" << func->GetChisquare() << endl;
 405           }
 406           delete gain;
 407         }
 408
 409       }
 410       delete p;
 411     }
 412
 413     // file outputs for gain
 414     if (!flatOutputFile.IsNull())
 415         fileout.close();
 416
 417     delete func;
 418
 419 }
 420
 421 //*************************************************************//
 422
 423 // main routine
 424 int main(Int_t argc, Char_t **argv)
 425 {
 426
 427     Int_t printLevel = 0;  // Global variable defined as extern in the others .cxx files
 428     Int_t skipEvents = 0;
 429     Int_t maxEvents  = 1000000;
 430     Float_t nSigma = 3;
 431     Int_t threshold = -1;
 432     Char_t inputFile[256];
 433     TString flatOutputFile;
 434     TString crocusOutputFile;
 435     TString crocusConfigFile;
 436
 437     Int_t totalParityErr = 0;
 438     Int_t totalGlitchErr = 0;
 439
 440 // option handler
 441
 442    // decode the input line
 443   for (Int_t i = 1; i < argc; i++) // argument 0 is the executable name
 444   {
 445       Char_t* arg;
 446
 447       arg = argv[i];
 448       if (arg[0] != '-') continue;
 449       switch (arg[1])
 450         {
 451         case 'f' :
 452           i++;
 453           sprintf(inputFile,argv[i]);
 454           break;
 455         case 'a' :
 456           i++;
 457           flatOutputFile = argv[i];
 458           break;
 459         case 'o' :
 460           i++;
 461           crocusOutputFile = argv[i];
 462           break;
 463         case 'c' :
 464           i++;
 465           crocusConfigFile = argv[i];
 466           break;
 467         case 'e' :
 468           i++;
 469           command = argv[i];
 470           break;
 471         case 'd' :
 472           i++;
 473           printLevel=atoi(argv[i]);
 474           break;
 475         case 's' :
 476           i++;
 477           skipEvents=atoi(argv[i]);
 478           break;
 479         case 'n' :
 480           i++;
 481           sscanf(argv[i],"%d",&maxEvents);
 482           break;
 483         case 'p' :
 484           i++;
 485           sscanf(argv[i],"%f",&nSigma);
 486           break;
 487         case 't' :
 488           i++;
 489           sscanf(argv[i],"%d",&threshold);
 490           break;
 491         case 'h' :
 492           i++;
 493           printf("\n******************* %s usage **********************",argv[0]);
 494           printf("\n%s -options, the available options are :",argv[0]);
 495           printf("\n-h help                   (this screen)");
 496           printf("\n");
 497           printf("\n Input");
 498           printf("\n-f <raw data file>        (default = %s)",inputFile);
 499           printf("\n-c <Crocus config. file>  (default = %s)",crocusConfigFile.Data());
 500           printf("\n");
 501           printf("\n Output");
 502           printf("\n-a <Flat ASCII file>      (default = %s)",flatOutputFile.Data());
 503           printf("\n-o <CROUCUS cmd file>     (default = %s)",crocusOutputFile.Data());
 504           printf("\n");
 505           printf("\n Options");
 506           printf("\n-d <print level>          (default = %d)",printLevel);
 507           printf("\n-s <skip events>          (default = %d)",skipEvents);
 508           printf("\n-n <max events>           (default = %d)",maxEvents);
 509           printf("\n-p <n sigmas>             (default = %f)",nSigma);
 510           printf("\n-t <threshold (-1 = no)>  (default = %d)",threshold);
 511           printf("\n-e <execute ped/gain>     (default = %s)",command.Data());
 512           printf("\n-e <gain create>           make gain & create a new root file");
 513           printf("\n-e <gain>                  make gain & update root file");
 514           printf("\n-e <gain compute>          make gain & compute gains");
 515
 516           printf("\n\n");
 517           exit(-1);
 518         default :
 519           printf("%s : bad argument %s (please check %s -h)\n",argv[0],argv[i],argv[0]);
 520           argc = 2; exit(-1); // exit if error
 521         } // end of switch
 522     } // end of for i
 523
 524   // set command to lower case
 525   command.ToLower();
 526
 527   // decoding the events
 528
 529   Int_t status;
 530   Int_t nDateEvents = 0;
 531
 532   void* event;
 533
 534    // containers
 535   AliMUONDDLTracker*       ddlTracker = 0x0;
 536   AliMUONBlockHeader*      blkHeader  = 0x0;
 537   AliMUONDspHeader*        dspHeader  = 0x0;
 538   AliMUONBusStruct*        busStruct  = 0x0;
 539
 540   pedMeanHisto = 0x0;
 541   pedSigmaHisto = 0x0;
 542
 543   TStopwatch timers;
 544
 545   timers.Start(kTRUE);
 546
 547   // once we have a configuration file in db
 548   // copy locally a file from daq detector config db
 549   // The current detector is identified by detector code in variable
 550   // DATE_DETECTOR_CODE. It must be defined.
 551   // If environment variable DAQDA_TEST_DIR is defined, files are copied from DAQDA_TEST_DIR
 552   // instead of the database. The usual environment variables are not needed.
 553   if (!crocusConfigFile.IsNull()) {
 554     status = daqDA_DB_getFile("myconfig", crocusConfigFile.Data());
 555     if (status) {
 556       printf("Failed to get config file : %d\n",status);
 557       return -1;
 558     }
 559   }
 560
 561
 562   status = monitorSetDataSource(inputFile);
 563   if (status) {
 564     cerr << "ERROR : monitorSetDataSource status (hex) = " << hex << status
 565               << " " << monitorDecodeError(status) << endl;
 566     return -1;
 567   }
 568   status = monitorDeclareMp("MUON Tracking monitoring");
 569   if (status) {
 570     cerr << "ERROR : monitorDeclareMp status (hex) = " << hex << status
 571               << " " << monitorDecodeError(status) << endl;
 572     return -1;
 573   }
 574
 575   cout << "MUONTRKda : Reading data from file " << inputFile <<endl;
 576
 577   while(1)
 578   {
 579     if (nEvents >= maxEvents) break;
 580     if (nEvents && nEvents % 100 == 0)
 581         cout<<"Cumulated events " << nEvents << endl;
 582
 583     // check shutdown condition
 584     if (daqDA_checkShutdown())
 585         break;
 586
 587     // Skip Events if needed
 588     while (skipEvents) {
 589       status = monitorGetEventDynamic(&event);
 590       skipEvents--;
 591     }
 592
 593     // starts reading
 594     status = monitorGetEventDynamic(&event);
 595     if (status < 0)  {
 596       cout<<"EOF found"<<endl;
 597       break;
 598     }
 599
 600     nDateEvents++;
 601
 602     // decoding rawdata headers
 603     AliRawReader *rawReader = new AliRawReaderDate(event);
 604
 605     Int_t eventType = rawReader->GetType();
 606     runNumber = rawReader->GetRunNumber();
 607
 608     if (eventType != PHYSICS_EVENT)
 609         continue; // for the moment
 610
 611     nEvents++;
 612
 613     // decoding MUON payload
 614     AliMUONRawStreamTracker* rawStream  = new AliMUONRawStreamTracker(rawReader);
 615
 616     // loops over DDL
 617     while((status = rawStream->NextDDL())) {
 618
 619       if (printLevel) printf("\niDDL %d\n", rawStream->GetDDL());
 620
 621       ddlTracker =  rawStream->GetDDLTracker();
 622
 623       // loop over block (CRT) structure
 624       Int_t nBlock = ddlTracker->GetBlkHeaderEntries();
 625       for(Int_t iBlock = 0; iBlock < nBlock ;iBlock++){
 626
 627         blkHeader = ddlTracker->GetBlkHeaderEntry(iBlock);
 628          if (printLevel) printf("Block Total length %d\n",blkHeader->GetTotalLength());
 629
 630         // loop over DSP (FRT) structure
 631         Int_t nDsp = blkHeader->GetDspHeaderEntries();
 632         for(Int_t iDsp = 0; iDsp < nDsp ;iDsp++){   //DSP loop
 633
 634           dspHeader =  blkHeader->GetDspHeaderEntry(iDsp);
 635            if (printLevel) printf("Dsp length %d\n",dspHeader->GetTotalLength());
 636
 637           // loop over BusPatch structure
 638           Int_t nBusPatch = dspHeader->GetBusPatchEntries();
 639           for(Int_t iBusPatch = 0; iBusPatch < nBusPatch; iBusPatch++) {
 640
 641             busStruct = dspHeader->GetBusPatchEntry(iBusPatch);
 642
 643              if (printLevel) printf("busPatchId %d", busStruct->GetBusPatchId());
 644              if (printLevel) printf(" BlockId %d", busStruct->GetBlockId());
 645              if (printLevel) printf(" DspId %d\n", busStruct->GetDspId());
 646
 647              Int_t busPatchId = busStruct->GetBusPatchId();
 648
 649             // loop over data
 650             Int_t dataSize = busStruct->GetLength();
 651             for (Int_t iData = 0; iData < dataSize; iData++) {
 652
 653               if (nEvents == 1)
 654                   nChannel++;
 655               Int_t manuId     = busStruct->GetManuId(iData);
 656               Int_t channelId  = busStruct->GetChannelId(iData);
 657               Int_t charge     = busStruct->GetCharge(iData);
 658
 659               if (printLevel) printf("manuId: %d, channelId: %d charge: %d\n", manuId,
 660                                       channelId, charge);
 661
 662               MakePed(busPatchId, manuId, channelId, charge);
 663
 664             } // iData
 665           } // iBusPatch
 666         } // iDsp
 667       } // iBlock
 668
 669       totalParityErr += rawStream->GetPayLoad()->GetParityErrors();
 670       totalGlitchErr += rawStream->GetPayLoad()->GetGlitchErrors();
 671     } // NextDDL
 672
 673     delete rawReader;
 674     delete rawStream;
 675
 676   } // while (1)
 677
 678   nEvents -= totalGlitchErr; // skipped event cos of Glitch errors
 679
 680
 681   cout << " Total number of parity errors in the Run: " << totalParityErr << endl;
 682   cout << " Total number of glitch errors in the Run: " << totalGlitchErr << endl;
 683
 684   if (command.CompareTo("ped") == 0)
 685       MakePedStore(flatOutputFile);
 686
 687   // option gain -> update root file with pedestal results
 688   // gain + create -> recreate root file
 689   // gain + comp -> update root file and compute gain parameters
 690
 691   if (command.Contains("gain"))
 692       MakePedStoreForGain();
 693
 694   if (command.Contains("comp"))
 695       MakeGainStore(flatOutputFile);
 696
 697
 698   delete pedestalStore;
 699
 700   timers.Stop();
 701
 702   if (!(crocusConfigFile.IsNull()))
 703       cout << "MUONTRKda : CROCUS command file generated : " << crocusOutputFile.Data() << endl;
 704   else
 705       cout << "MUONTRKda : WARNING no CROCUS command file generated" << endl;
 706
 707   cout << "MUONTRKda : Flat ASCII file generated     : " << flatOutputFile << endl;
 708   cout << "MUONTRKda : Histo file generated          : " << histoFileName  << endl;
 709   cout << "MUONTRKda : Nb of DATE events     = "         << nDateEvents    << endl;
 710   cout << "MUONTRKda : Nb of events used     = "         << nEvents        << endl;
 711
 712   printf("Execution time : R:%7.2fs C:%7.2fs\n", timers.RealTime(), timers.CpuTime());
 713
 714   return status;
 715 }