#include <iostream>
#include <fstream>
#include <sstream>
#include <string>
#include <iomanip>
#include <limits>
#include <math.h>
#include <TArc.h>
#include <TArrow.h>
#include <TBox.h>
#include <TBrowser.h>
#include <TCanvas.h>
#include <TCondition.h>
#include <TDatime.h>
#include <TDirectory.h>
#include <TEllipse.h>
#include <TGraphErrors.h>
#include <TGButton.h>
#include <TGButtonGroup.h>
#include <TGClient.h>
#include <TGComboBox.h>
#include <TGProgressBar.h>
#include <TGLabel.h>
#include <TGTextEntry.h>
#include <TGFrame.h>
#include <TGraph.h>
#include <TGraphErrors.h>
#include <TF1.h>
#include <TFile.h>
#include <TH1D.h>
#include <TH2D.h>
#include <TLatex.h>
#include <TLeaf.h>
#include <TLine.h>
#include "TMath.h"
#include <TObjectTable.h>
#include <TPaveText.h>
#include <TPolyLine.h>
#include <TProfile.h>
#include <TRandom.h>
#include <TRootEmbeddedCanvas.h>
#include <TStyle.h>
#include <TTree.h>
#include "TSystem.h"
#include <TVectorD.h>
#include <TDatime.h>
#include <sys/types.h>
#include <time.h>
#include <TROOT.h>
#include <TError.h>
#include "TVirtualFFT.h"
#define SizeOfArray(a) (sizeof(a)/sizeof(a[0]))
#include <classes2016/myCovariance.C>
#include <classes2016/MIPcounter2012.C>
#include <classes2016/myRCCR.C>
#include <classes2016/parabola2012.C>
#include <classes2016/interpolate2012.C>
#include <classes2016/RingBuffer.C>
#include <classes2016/EventHeader.C>
#include <classes2016/miniBurstHeader.C>
#include <classes2016/miniEventHeader.C>
#include <classes2016/eventEssentials.C>
#include <classes2016/BurstHeader.C>
#include <classes2016/RunHeader.C>
#include <classes2016/ProcessingStatus.C>
#include <classes2016/CherenkovBuffer.C>
#include <classes2016/CherenkovCalibration.C>
#include <classes2016/TimingCuts.C>
#include <classes2016/TReference.C>
#include <classes2016/TSignal.C>
#include <classes2016/TChannel.C>
#include <classes2016/TModule.C>
#include <classes2016/SetMIPrunHVgroup.C>
#include <classes2016/DisplaySpike.C>
#include <classes2016/DisplayModules.C>


using std::ios;
using std::cout;
using std::cin;
using std::endl;
using std::istringstream;
using std::ostringstream;

std::ifstream indata;
std::ifstream RFdata;
std::ofstream spilldata;

std::ifstream fitdata;


miniBurstHeader miniBH;
miniEventHeader miniEH;
eventEssentials     evtEssentials;

const int nMod=7;
int myRun;
const bool saveSpillProfile=true;
char myText[200];
char fName[200];
char treeName[200];
char branchName[200];

char version_out[100];

//1019,1060,1016,1008
//----------------------------- Global Objects ------------------------
MIPcounter *myMIPcounter;




//------------------------------------------------------- global histograms for results ----------------------//
TH1D    *SignalFactor1;

TLatex *myLatex;
TLine  *myLine;

//double T0noise=210.; // interval for noise measurements
//double T1noise=85.;
// in 2012 Tbunch shifted left by 50 ns!
double T0noise=160.; // interval for noise measurements
double T1noise=35.;



// histograms for RF study //
TH1D *RFperiod;
TH1D *RFamplitude;
TH1D *RFratio;
TH1D *RFt0;
TH1D *RFt0diff;
TH1D *RFresidual;
// end of RF histograms;


TTree *runTree;
TTree *eventsTree;
TTree *eventsTreeModules;
TTree *eventsTreeCounters;
TTree *spillsTree;
TFile *treeFile;
TFile *treeFileModules;
TFile *treeFileCounters;



// end of histograms for MIP study;

int nPulses;
double logRateBin;

TH2D     *CoGXY;
TH1D     *hGlobalPhase;
TH1D     *hBunchPhase;

TCanvas  *PhaseCanvas;

//****************************** canvas and plots for rate history Cuts ***************

TCanvas *PilotCanvas; 

//************************************ RF correction for M2 and Cherenkov *************
double  sum0rf[2][2]; // M2,Cherenkov and Straight,Delayed;
double  sum1rf[2][2];
double  sum2rf[2][2];
double  sum3rf[2][2];
//*************************************************************************************
// FFT noise subtraction stuff
int emptyCount;
int eventsDumped;
const double FFTfreqScale = 0.15789; 
TH2D *noiseReduction;
//*************************************************************************************
// FFT tailed signal size for 12.5ns
int nSampleTail125 = 2000;
//*************************************************************************************

double HAmax0[4], HAmax1[4], HAavr0[4], HAavr1[4]; // high activity parameter per module
double RCtime[100];
bool   flagRCtime[100];

int    nRC;
int    irc1; // where the modules drift time start
double *rateIntegrals;
double *slowIntegrals;
double *beamiMAX;
double *beamiMIN;
bool   *minSET;
bool   *maxSET;

//TH1D *MaxFuzz;



int blockCounter;

//------------------------------------------------------- Classes etc ----------------------------------------//

CherenkovCalibration *CherCal;
struct scintMIPs {
  double scintMIP[3]; //scintillator MIPs for the Run in MPI FFT units (integral)
};

const char *module_name[nMod] = { 
  "HEC", 
  "EMEC", 
  "FCal025", 
  "FCal010", 
  "Beam",
  "Monitor",
  "Integrator"
};

int NADC, kMd[96], kCh[96], kDl[96], kGn[96]; // detector to ADC mapping -> channel, module, gain and delay
int TDCMaxMin[2][2]; 
int nBBlock=0;
bool isOnline=false ; //if TRUE, supresses some advanced options (like logbook use) 
                      // for online express analysis
bool showMode = false; // shows spill profiles, events prehistroy, events themselves, Cherenkov plots too 

unsigned short GetWord();
void fFADC(unsigned short BWC, unsigned short BID, TModule* Setup, EventHeader* EH);
void fTDC(unsigned short BWC, unsigned short BID, EventHeader* EH);
void CalculateIntegrator(TModule* Setup, ProcessingStatus* ProcStatus);
int  LoadDetectorMap(Int_t Run, TModule* Setup);
bool LoadPMTdelays(Int_t Run, TModule* Setup);
bool LoadMIP2data(Int_t Run, TModule* Setup);
bool LoadMIP2dataGroup(Int_t Run, TModule* Setup);
bool LoadChannelStatus(Int_t Run, TModule* Setup);
bool PrepareTimeSlices(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus);
bool ProcessOneRun(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
void RunReport(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer*CherBuff);
void SetModuleTiming(Int_t Run, TModule* Setup);
void SavePulseShapes(Int_t Run, TModule* Setup);
unsigned int Unpack(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
bool AnalyseEvent(TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
void FitModule(TModule* Setup, int Condition, ProcessingStatus* ProcStat);
int  FillHistograms(TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
unsigned int EventPreProcessing(TModule *Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
int DisplayAllPulsesDetails(TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
int DisplayAllPulses2012(TModule* Setup,    ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
int DisplayAllCounters2012(TModule* Setup,  ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
int DisplayForMaslennikov(TModule* Setup,   ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff);
double GlobalPhase2011(TModule* Setup, ProcessingStatus* ProcStat);
bool   GlobalPhase2012(TModule* Setup, ProcessingStatus* ProcStat);
void MyExecuteEvent(Int_t event, Int_t px, Int_t py, TObject *sel);
void ResetHistograms(TModule* Setup);
//bool SetChannelHiLoCuts2013(Int_t myRun, TModule* Setup, ProcessingStatus* ProcStatus);
bool TuneStreamCherenkovTiming(TModule* Setup, CherenkovBuffer* CherBuff, ProcessingStatus* ProcStatus);
bool LoadRFcorrection(Int_t myRun, TModule* Setup);
bool LoadFFTcorrection(Int_t myRun, TModule* Setup);
bool LoadTransferFunction(Int_t iOpt, TModule* Setup);

int iqq;

TCanvas *SummaryCanvas;
TCanvas *SpectraCanvas;
TCanvas *SignalCanvas;
TCanvas *HiLoCanvas;

TCanvas *EventCanvasSpecial;
TCanvas *EventCanvasNew;
TCanvas *EventCanvas;
TCanvas *EventCanvasModule;

///////////////////////////////////////////////////////////////////////////////////////////////////
void HiLum2016_0222_tree4() {
  bool Dbg=false;
  showMode=false;


  TGClient *TC = new TGClient();
  int displayWidth  = TC->GetDisplayWidth();
  int displayHeight = TC->GetDisplayHeight();
  //printf("w=%d h=%d\n", displayWidth, displayHeight);

  cout << "***********************************************************************" << endl;
  cout << "*Analysis for 2012-2013     HiLum runs                                *" << endl;
  cout << "*Dynamic pedestal for every MPI ADC channel;                          *" << endl;
  cout << "*3sigma pedestal cleaning                                             *" << endl;
  cout << "*Bunch timing done by RF record                                       *" << endl;
  cout << "*Uses results from HiLum20120617computeRFcorrection.C                 *" << endl;
  cout << "*Uses files:      reports/pickupReductionFFTcorr00run****binary.dat   *" << endl;
  cout << "*                 reports/pickupReductionFFTcorr00run****map.dat      *" << endl;
  cout << "*   for FFT pickup reduction/filtering                                *" << endl;
  cout << "*                                                                     *" << endl;
  cout << "* Uses resut of HiLum2013mip0607.C and Summary2013mip0607.C           *" << endl;
  cout << "*  to fill low intensity spectra                                      *" << endl;
  cout << "* Uses results of 0613step2 to select single-proton events and file   *" << endl;
  cout << "*  single-proton spectra for S1, S2, S3 and Cherenkov                 *" << endl;
  cout << "*                                                                     *" << endl;
  cout << "*  Uses LikeliHood MIP counting from Summary2013mip1115step3_100.C;   *" << endl;
  cout << "*  Computes integral for cluster also                                 *" << endl;
  cout << "*  Position using integral over all reconstructed filled bunches      *" << endl;
  cout << "*                                                                     *" << endl;
  cout << "*  !!! does al lthe job, but saves only S123 vs Rate study !!!!!      *" << endl;
  cout << "*                                                                     *" << endl;
  cout << "*  tries 0.2/0.05 GADC leak correction and its impact on the rate     *" << endl;
  cout << "*  estimation                                                         *" << endl;
  cout << "*                                                                     *" << endl;
  cout << "* 2015.08.13                                                          *" << endl;
  cout << "***********************************************************************" << endl;


  if (Dbg) {
    printf ("\n************************************************\n");
    printf ("    Attention! Running in junky debug mode !    \n");
    printf ("************************************************\n");
  }

  sprintf(version_out, "%s", "2016_0222_tree2");

  myLine  = new TLine();
  myLatex = new TLatex();

  char TitlC[200];
  char NameC[200];
  RunHeader *RH;

  myMIPcounter = new MIPcounter(100);
/*
  evtEssentials.tInt[0]=1.0e6; // ns, 1ms integration
  evtEssentials.tInt[1]=1.0e7; // ns, 10ms integration
  evtEssentials.tInt[2]=3.0e7; // ns, 30ms integration
  evtEssentials.tInt[3]=1.0e8; // ns, 100ms integration
*/
  int dWidth  = displayWidth;
  if (dWidth>1800) dWidth=1800;
  PilotCanvas = new TCanvas("PilotCanvas", "Pilot Canvas", 
                            displayWidth-int(dWidth/2), 0, int(dWidth/2), int(dWidth/2)); // used for channels significance 
  PilotCanvas -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                           "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");
  PilotCanvas->Iconify();
//  noiseReduction = new TH2D("noiseReduction", "HEC Sum Noise Reduction", 200, 0., 25., 200, 0., 25.);
  noiseReduction = new TH2D("noiseReduction", "FCal Sum Noise Reduction", 200, 0., 50., 200, 0., 50.);
  noiseReduction->GetXaxis()->SetTitle("Raw noise, HG ADC");
  noiseReduction->GetYaxis()->SetTitle("Reduced noise, HG ADC");
  noiseReduction->SetStats(kFALSE);

  int dw=int(displayWidth/2.);
//  if (dw>1000) {
//    dw=1000;
//  }
  EventCanvasSpecial = new TCanvas("EventCanvasSpecial", "Event Canvas Special", 
                   dw/2,         0, 
                   dw/2, int(0.9*displayHeight-10.)); // used for event display
  EventCanvasNew = new TCanvas("EventCanvasNew", "Event Canvas", 
                   0,            0, 
                   dw/2, int(0.9*displayHeight-10.)); // used for event display
  EventCanvas = new TCanvas("EventCanvas", "Event Canvas", 
                   3*dw/2,        0, 
                   dw/2,        int(0.9*displayHeight-10.)); // used for event display
  EventCanvasModule = new TCanvas("EventCanvasModule", "Event Canvas for Modules", 
                   4*dw/2,        0, 
                   dw/2,        int(0.9*displayHeight-10.)); // used for event display


  EventCanvasSpecial -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                           "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");
  EventCanvasNew -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                           "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");
  EventCanvas    -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                           "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");

  EventCanvas->Iconify();
  EventCanvasNew->Iconify();
  EventCanvasSpecial->Iconify();


  int canvasWidth  = displayWidth/4-20;
  int canvasHeight = canvasWidth*6/4;
//  canvasWidth=int(0.5*canvasWidth);
  canvasHeight=int(2*canvasWidth);
  if (canvasHeight>displayHeight-10) canvasHeight=displayHeight-10;

  SummaryCanvas = new TCanvas("SummaryCanvas", " ", 3*(canvasWidth+5), 0, canvasWidth, canvasHeight);
  SummaryCanvas -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                           "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");
  SummaryCanvas->Iconify();

  SpectraCanvas = new TCanvas("SpectraCanvas", " ", 0, 0, canvasWidth, canvasHeight);
  SpectraCanvas -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                           "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");
  SpectraCanvas -> Iconify();

  SignalCanvas = new TCanvas("SignalCanvas", " " , canvasWidth+5, 0, canvasWidth, canvasHeight);
  SignalCanvas -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                           "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");
  SignalCanvas -> Iconify();

  HiLoCanvas = new TCanvas("HiLoCanvas", " ", 2*(canvasWidth+5), 0, canvasWidth, canvasHeight);
  HiLoCanvas -> Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", 0, 0, 
                        "MyExecuteEvent(Int_t,Int_t,Int_t,TObject*)");
  HiLoCanvas -> Iconify();


  ProcessingStatus *ProcStatus = new ProcessingStatus(false, true); // do not change those settings!
  ProcStatus->SetOnline(isOnline);
  ProcStatus->GetEventHeader()->EmptyCut=0.4;
  ProcStatus->GetEventHeader()->NotEmptyCut=0.9;
  RH = ProcStatus->GetRunHeader();


  //ProcStatus -> SetReportMode(2);
  ProcStatus -> SetReportMode(5);
  ProcStatus -> SetReportMode(1);
  if (showMode) ProcStatus -> SetReportMode(1);
  if (!ProcStatus -> LoadLocation()) {
    exit(-1);
  }

  if (isOnline) {
    ProcStatus -> SetLogbookFile("../doc/runbook2012dummy.elog"); // defective logbook file
    ProcStatus -> SetRunListBasic(); // prepares runs list to process - basic version for Online!
  } else {
    ProcStatus -> SetLogbookFile("../doc/runbook2013.elog"); // use normal logbook file
    ProcStatus -> SetRunList(); // prepares runs list to process
  }
  ProcStatus -> SetDisplay(displayWidth, displayHeight);    
   
  if (ProcStatus -> GetLocation()==1) {  
    sprintf(NameC,  "../cherenkov/data/");
    sprintf(TitlC, "../cherenkov/spill_lists/");
  } else if (ProcStatus -> GetLocation()==0) {
    sprintf(NameC,  "/raid04/HiLum/cherenkov/data/");
    sprintf(TitlC, "/raid04/HiLum/cherenkov/spill_lists/");
  }
  CherenkovBuffer  *CherBuffer = new CherenkovBuffer(NameC, TitlC, ProcStatus);
  CherBuffer->SetCutMode(1);


  CherBuffer->SetdT(164.97);
  CherBuffer->SetCherenkovTimeLine(25.);
  if (!isOnline) CherBuffer->PrepareHistograms();


 // Create Cherenkov MPI Calibration object // 
  CherCal = new CherenkovCalibration(ProcStatus);

  // Prepare the modules channel status // 
  TModule *Setup = new TModule[nMod];

  

  // prepare modules and channels //
  for (int iMod=0; iMod<nMod; iMod++) {
    Setup[iMod].SetID(iMod);
    Setup[iMod].SetName(module_name[iMod]);
    Setup[iMod].SetSlices(252);
    Setup[iMod].ClearEmptyEventPedestal();

    if (iMod==0) { // HEC
      Setup[iMod].GetChannel(0)->SetChannelPosition(-1.0,  1.0);
      Setup[iMod].GetChannel(1)->SetChannelPosition( 1.0,  1.0);
      Setup[iMod].GetChannel(2)->SetChannelPosition(-1.0, -1.0);
      Setup[iMod].GetChannel(3)->SetChannelPosition( 1.0, -1.0);
    } else if (iMod==1) { // EMEC
      Setup[iMod].GetChannel(0)->SetChannelPosition(-1.0,  1.0);
      Setup[iMod].GetChannel(1)->SetChannelPosition(-1.0, -1.0);
      Setup[iMod].GetChannel(2)->SetChannelPosition( 1.0,  1.0);
      Setup[iMod].GetChannel(3)->SetChannelPosition( 1.0, -1.0);
    } else if (iMod<4) { // FCals
      Setup[iMod].GetChannel(0)->SetChannelPosition( 1.0,  1.0);
      Setup[iMod].GetChannel(1)->SetChannelPosition(-1.0,  1.0);
      Setup[iMod].GetChannel(2)->SetChannelPosition( 1.0, -1.0);
      Setup[iMod].GetChannel(3)->SetChannelPosition(-1.0, -1.0);
    } 

    if (iMod==4) {
      Setup[iMod].GetChannel(0)->SetName("S1");
      Setup[iMod].GetChannel(1)->SetName("S2");
      Setup[iMod].GetChannel(2)->SetName("S3");
      Setup[iMod].GetChannel(3)->SetName("Ch");
    } else if (iMod==5) {
      Setup[iMod].GetChannel(0)->SetName("S4");
      Setup[iMod].GetChannel(1)->SetName("S5");
      Setup[iMod].GetChannel(2)->SetName("RF");
      Setup[iMod].GetChannel(3)->Off();
    } else if (iMod==6) {
      Setup[iMod].GetChannel(0)->SetName("2.60ms");
      Setup[iMod].GetChannel(1)->SetName("18.4ms");
      Setup[iMod].GetChannel(2)->SetName("33.0ms");
      Setup[iMod].GetChannel(3)->SetName("65.0ms");
    }
    int iCh0=0;
    int iCh1=5;
    if (iMod>3) iCh1=4; // sums (#4) for modules only!
    for (int iChan=iCh0; iChan<iCh1; iChan++) {
      Setup[iMod].GetChannel(iChan)->On();
    }
    Setup[iMod].SetProcStatus(ProcStatus);
  }

  PhaseCanvas = new TCanvas("PhaseCanvas", "Phase test ", 0, 0, 500, 800);
  hGlobalPhase = new TH1D("GlobalPhase", "Global Phase", 1000, 0., 1000.);
  hBunchPhase  = new TH1D("BunchPhase",  "Bunch Phase",  1000, 0., 1000.);
  PhaseCanvas->Iconify();

  nRC=35;
  minSET        = (bool*)   calloc(nRC*20, sizeof(bool));
  maxSET        = (bool*)   calloc(nRC*20, sizeof(bool));

  // RF plots ---->
  RFperiod    = new TH1D[2];
  RFamplitude = new TH1D[2];
  RFt0        = new TH1D[2];
  RFresidual  = new TH1D[2];
  RFt0diff    = new TH1D("RFt0diff",  "RF T0 difference",       250,  -0.5,  0.5);
  RFratio     = new TH1D("RFratio",   "Delayed/Straight amplitude ratio", 1000, 0.9, 1.0);

  for (int iDelay=0; iDelay<2; iDelay++) {
    sprintf(NameC, "RFperiod%d", iDelay);
    sprintf(TitlC, "RFperiod %d", iDelay);
    RFperiod[iDelay].SetBins(250, 164.75, 165.25);
    RFperiod[iDelay].SetNameTitle(NameC, TitlC);
    RFperiod[iDelay].SetLineColor(iDelay+1);

    sprintf(NameC, "RFamplitude%d", iDelay);
    sprintf(TitlC, "RFamplitude %d", iDelay);
    RFamplitude[iDelay].SetBins(600, 330., 630.);
    RFamplitude[iDelay].SetNameTitle(NameC, TitlC);
    RFamplitude[iDelay].SetLineColor(iDelay+1);
 
    sprintf(NameC, "RFt0%d", iDelay);
    sprintf(TitlC, "RFt0 %d", iDelay);
    RFt0[iDelay].SetBins(50, -25., 25.);
    RFt0[iDelay].SetNameTitle(NameC, TitlC);
    RFt0[iDelay].SetLineColor(iDelay+1);

    sprintf(NameC, "RFresidual%d", iDelay);
    sprintf(TitlC, "RFresidual RMS %d", iDelay);
    RFresidual[iDelay].SetBins(100, 0., 10.);
    RFresidual[iDelay].SetNameTitle(NameC, TitlC);
    RFresidual[iDelay].SetLineColor(iDelay+1);
  }
  

// <---- end of RF plots
  
  SignalFactor1       = new TH1D("SignalFactorHist", "Signal Factor for Global Phase", 55, -0.05, 1.05); 
  SignalFactor1 -> SetStats(0);
// <---- end of Signal Factor
 

  EventHeader* EH = ProcStatus->GetEventHeader();

  while ((myRun=ProcStatus->GetNextRunToProcess())>0) {

    printf("[main] Processing run#%d\n", myRun);
    bool RunReady = true;

    //if (!isOnline) {
    RunReady = RunReady && ProcStatus->SetRun(myRun, isOnline); // load logbook data 
                                                                  //(more to be moved here)
    //}
    printf("[main] run was set, RunReady=%d\n", RunReady);
    printf("[main] beam counters in = %d\n", RH->S123inBeam);
//    if (RH->S123inBeam) continue;


    int iAtt = int(0.5+RH->AttnSet);
    printf("[main] AttnSet=%f, iAtt=%d\n", RH->AttnSet, iAtt);
    if (!LoadTransferFunction(-iAtt, Setup)) {
      printf("[main] Transfer Function not loaded\n");
      exit(-1);
    }

    if (RH->S123inBeam) {
      if (!LoadPMTdelays(myRun, Setup))      continue;
      if (!LoadMIP2dataGroup(myRun, Setup))  continue; // uses HV groups to get MIP cuts
    }

    NADC=LoadDetectorMap(myRun, Setup);              // load MPI ADC mapping
    
    printf("[main] detector map loaded NADC=%d\n", NADC);
    if (!NADC) {
      printf("[main] Run#%d --> ADC Mapping load failed!\n", myRun);
      exit(-1);
    }

    if (myRun) {
      Setup[6].SetIntegrator(myRun);                  // prepare parameters for integrators
    }
    
    if (!LoadChannelStatus(myRun, Setup)) {  // channel status loaded?
      printf("[main] Run#%d channel status load failed\n", myRun);
      continue;
    }

    if (RunReady 
//                && ((SetChannelHiLoCuts2011(myRun, Setup, ProcStatus)  // load prehistory rules, 2011 version
               //&& SetChannelHiLoCuts(myRun, Setup, ProcStatus) // load prehistory rules 
//               && LoadCoG(myRun, Setup)) 
//                 || isOnline)) 
              )                 {                        // for High/Low crisis resolution
      if (ProcStatus->RunIsGood(myRun)>0) {
        SetModuleTiming(myRun, Setup);                 // load timing correction to allign all signals   
        CherCal->SetCorrection();
        double ChHVcorrMPIF=CherCal->GetCorrectionMPIF(); // find HV correction
        double ChHVcorrGADC=CherCal->GetCorrectionGADC(); // find HV correction
        Setup[4].GetChannel(3)->SetGainCorrection(ChHVcorrMPIF); // set gain correction for Cherenkov counter
        CherBuffer->SetGainCorrection(ChHVcorrGADC);

        ResetHistograms(Setup);
        emptyCount=0;
        eventsDumped=0;
        RH->PedestalEvents=0;

        for (int iMod=0; iMod<4; iMod++) {
          int irc=nRC-4+iMod;
          RCtime[irc]=ProcStatus->GetTdrift(iMod);
        }

        RH->ClearLowGainOverflows();
        RH->notEmpty=0;

//        continue;
        sprintf(fName, "/raid04/HiLum/analysis2012/ntuples/v4/run%04d_v4.root", myRun);
        treeFile = new TFile(fName, "RECREATE");

        sprintf(treeName, "run%04d", myRun);
        runTree = new TTree(treeName, treeName);
        runTree -> Branch("S123inBeam",     RH->GetS123pointer(),                 "S123inBeam/O");
        runTree -> Branch("HVmod",          RH->GetHVsetPointer(),                "HVmod[4]/D");
        runTree -> Branch("Attenuator",     RH->GetAttenuatorPointer(),           "Attenuator/D");
        runTree -> Branch("AntennaModule",  RH->GetAntennaPointer(),              "AntennaModule/I");
        runTree -> Branch("CryoX",          RH->GetCryoXPointer(),                "CryoX[3]/D");
        runTree -> Branch("ScintMIPs",      myMIPcounter->GetMIPaveragePointer(), "ScintMIPs[3]/D");
       

        sprintf(treeName, "run%04dspills", myRun);
        spillsTree = new TTree(treeName, treeName);
        spillsTree -> Branch("miniBirstHeader", &miniBH.burst,                   "burst/I:ChMPIF/F:ChGADC/F");
        spillsTree -> Branch("nLaps",      CherBuffer->GetNlapsPointer(),        "nLaps/I");
        spillsTree -> Branch("LapSum",     CherBuffer->GetLapsSumPointer(),      "lapSum[nLaps]/F"); // all filled bunches, summed over laps
        spillsTree -> Branch("LapFill",    CherBuffer->GetFilledSumPointer(),    "lapFill[30]/D"); // accumulated amplitudes for 30 bunches
        spillsTree -> Branch("LapBool",    CherBuffer->GetFilledSamplePointer(), "lapBool[30]/O"); // true for filled, false for empty
        spillsTree -> Branch("LapActn",    ProcStatus->GetBunchActionPointer(),  "lapActn[30]/I"); // 0 - no action 1 - signal 2 - ped
        spillsTree -> Branch("HVbeam",     RH->GetHVbeamPointer(),               "HVbeam[4]/D");   // S123 from GADC record
        spillsTree -> Branch("HVbeamRaw",  RH->GetHVbeamRawPointer(),            "HVbeamRaw[4]/D");   // S123 from GADC record


        sprintf(treeName, "run%04devents", myRun);
        eventsTree = new TTree(treeName, treeName);
        eventsTree -> Branch("EventHeaderBranch", &miniEH.RFperiod, "RFperiod/D:RFfreq:RFamp:TimeGlobal:SignalFactor:burst/I:event:burstEvent:scaler10mhz:scaler6mhz:firstGated:GADCoverflows:isEmpty/O:notEmpty:filledSegment:inHotZone:isSpiky");// branch for event header;
        eventsTree -> Branch("EventGADC",       &evtEssentials.GADCevent[0],      "GADCevent[40]/F");
        eventsTree -> Branch("EventGADC",       &evtEssentials.GADCeventTail[0],  "GADCeventTail[40]/F");
        eventsTree -> Branch("MPICevent",       &evtEssentials.MPICeventHG[0],    "MPICeventHG[40]/F");
        eventsTree -> Branch("MPICevent",       &evtEssentials.MPICeventLG[0],    "MPICeventLG[40]/F");
        eventsTree -> Branch("firstGADCsample", &evtEssentials.firstGADCsample,   "firstGADCsample/I");
        eventsTree -> Branch("bunchAction",     &evtEssentials.bunchAction[0],    "bunchAction[40]/I");

        eventsTree -> Branch("scintFFT",        &evtEssentials.scintFFT[0],     "scintFFT[120]/F");

        eventsTree -> Branch("isIsolated",      &evtEssentials.isIsolated[0],   "isIsolated[40]/O");
        eventsTree -> Branch("isCluster0",      &evtEssentials.isCluster0[0],   "isCluster0[40]/O");
        eventsTree -> Branch("isCluster1",      &evtEssentials.isCluster1[0],   "isCluster1[40]/O");
        eventsTree -> Branch("inScope",         &evtEssentials.inScope[0],      "inScope[40]/O");
       

        eventsTree -> Branch("HGoverflow",      &evtEssentials.HGoverflow[0],   "HGoverflow[4]/O");
        eventsTree -> Branch("LGoverflow",      &evtEssentials.LGoverflow[0],   "LGoverflow[4]/O");

        eventsTree -> Branch("RFclust0",        &evtEssentials.RFclust0[0],     "RFclust0[5]/I");
        eventsTree -> Branch("RFclust1",        &evtEssentials.RFclust1[0],     "RFclust1[5]/I");
  
        eventsTree -> Branch("miniModule",      &evtEssentials.miniModule[0],   "miniModule[200]/F");

        for (int iChan=0; iChan<4; iChan++) {
          TModule *Module = &Setup[4];
          TChannel *Channel = Module->GetChannel(iChan);
          for (int iGain=0; iGain<2; iGain++) {
            TSignal *Signal = Channel->GetSignalLong(iGain);
            sprintf(branchName, "Beam_bunch_integral_ch%d_gain%d", iChan, iGain);
            eventsTree -> Branch(branchName, Signal->GetSum2013Pointer(), "sum2013[50]/F");
            sprintf(branchName, "Beam_bunch_pedestal_ch%d_gain%d", iChan, iGain);
            eventsTree -> Branch(branchName, Signal->GetPed2013Pointer(), "ped2013[50]/F");
            sprintf(branchName, "Beam_bunch_Thalf0_ch%d_gain%d", iChan, iGain);
            eventsTree -> Branch(branchName, Signal->GetThalf0Pointer(), "Thalf0[50]/F");
            sprintf(branchName, "Beam_bunch_Thalf1_ch%d_gain%d", iChan, iGain);
            eventsTree -> Branch(branchName, Signal->GetThalf1Pointer(), "Thalf1[50]/F");
          }
        }

//==========================================================
        sprintf(fName, "/raid04/HiLum/analysis2012/ntuples/v4/run%04dmodules_v4.root", myRun);
        treeFileModules = new TFile(fName, "RECREATE");

        sprintf(treeName, "run%04deventModules", myRun);
        eventsTreeModules = new TTree(treeName, treeName);

        for (int iMod=0; iMod<4; iMod++) {
          TModule *Module = &Setup[iMod];
          int iChan1=5;
          for (int iChan=0; iChan<iChan1; iChan++) {
            TChannel *Channel = Module->GetChannel(iChan);
            for (int iGain=0; iGain<2; iGain++) {
              TSignal *Signal = Channel->GetSignalLong(iGain);
              if (iGain==0) { // Pedestal Windows are same for both gains, do it once
                sprintf(branchName, "nPedWindows_mod%dchan%d", iMod, iChan);
                eventsTreeModules -> Branch(branchName, Signal->GetNPedWindowsPointer(), "nPedWindows/I");
                sprintf(branchName, "firstPedWindowSample_mod%dchan%d", iMod, iChan);
                eventsTreeModules -> Branch(branchName, Signal->GetFirstSamplePointer(), "firstSample[25]/I");
                sprintf(branchName, "lastPedWindowSample_mod%dchan%d",  iMod, iChan);
                eventsTreeModules -> Branch(branchName, Signal->GetLastSamplePointer(),  "lastSample[25]/I");
              }
              sprintf(branchName, "signal_mod%dchan%dgain%d", iMod, iChan, iGain); // MPI ADC signals, mV
              eventsTreeModules -> Branch(branchName, Signal->GetAmpIFpoint(), "ampIF[512]/F"); 
              sprintf(branchName, "current_mod%dchan%dgain%d", iMod, iChan, iGain); // Module currents, mA
              eventsTreeModules -> Branch(branchName, Signal->GetAmpFinalFpoint(), "signalFinalF[512]/F");
            }
          }
        }
//==========================================================
        sprintf(fName, "/raid04/HiLum/analysis2012/ntuples/v4/run%04dcounters_v4.root", myRun);
        treeFileCounters = new TFile(fName, "RECREATE");
        sprintf(treeName, "run%04deventCounters", myRun);
        eventsTreeCounters = new TTree(treeName, treeName);

        for (int iMod=4; iMod<5; iMod++) {
          TModule *Module = &Setup[iMod];
          int iChan1=4;
          for (int iChan=0; iChan<iChan1; iChan++) {
            TChannel *Channel = Module->GetChannel(iChan);
            for (int iGain=0; iGain<2; iGain++) {
              TSignal *Signal = Channel->GetSignalLong(iGain);
              if (iGain==0) { // Pedestal Windows are same for both gains, do it once
                sprintf(branchName, "nPedWindows_mod%dchan%d", iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetNPedWindowsPointer(), "nPedWindows/I");
                sprintf(branchName, "firstPedWindowSample_mod%dchan%d", iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetFirstSamplePointer(), "firstSample[25]/I");
                sprintf(branchName, "lastPedWindowSample_mod%dchan%d",  iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetLastSamplePointer(),  "lastSample[25]/I");
                
                sprintf(branchName, "t0pmtREC_mod%dchan%d",  iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetT0pmtRecPointer(),  "t0pmtREC/D");
                sprintf(branchName, "t1pmtREC_mod%dchan%d",  iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetT1pmtRecPointer(),  "t1pmtREC/D");
                sprintf(branchName, "time2013start_mod%dchan%d",  iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetT02013pointer(),    "time2013start/D");
                sprintf(branchName, "time2013step_mod%dchan%d",  iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetDT2013pointer(),    "time2013step/D");
                sprintf(branchName, "ns2013_mod%dchan%d",  iMod, iChan);
                eventsTreeCounters -> Branch(branchName, Signal->GetNS2013pointer(),    "ns2013/I");
               
              }
              sprintf(branchName, "signal_mod%dchan%dgain%d", iMod, iChan, iGain); // MPI ADC signals, mV
              eventsTreeCounters -> Branch(branchName, Signal->GetAmpIFpoint(),       "ampIF[512]/F"); 
              sprintf(branchName, "amp2013_mod%dchan%dgain%d",  iMod, iChan, iGain); // scintillators, FFT
              eventsTreeCounters -> Branch(branchName, Signal->GetAmp2013pointer(),   "amp2013[ns2013]/F");
            }
          }
        }
///////////////////////////////////////////////////////////////////////
        if (ProcessOneRun(myRun, Setup, ProcStatus, CherBuffer)) {
          RunReport(myRun, Setup, ProcStatus, CherBuffer);  
        } else if (ProcStatus->GetReportMode()>0) {
          printf("[main] **** Run#%d analysis failed! ****\n", myRun);
        }
///////////////////////////////////////////////////////////////////////
        eventsTree->Print();
        printf("[Main] writing general root file ->\n");
        treeFile->cd();
        runTree->Write();
        spillsTree->Write();
        eventsTree->Write();
        printf("[Main] general root file written\n");
        treeFile->Close();
        printf("[Main] general root file closed\n");

        treeFileModules->cd();
        eventsTreeModules->Write();
        printf("[Main] modules root file written\n");
        treeFileModules->Close();
        printf("[Main] modules root file closed\n");


        printf("[Main] writing Counters root file ->\n");
        treeFileCounters->cd();
        eventsTreeCounters->Write();
        printf("[Main] Counters root file written\n");
        treeFileCounters->Close();
        printf("[Main] Counters root file closed\n");

       
       
       

//        eventsTree->Delete("all");
//        printf("[Main] tree deleted\n");

        if (ProcStatus->GetReportMode()<0) {
          printf("[main] Run#%d ended with ReportMode=%d\n", myRun, ProcStatus->GetReportMode());
          break;
	}
      }

    } else {
      printf ("[main] Run#%d --> Logbook error\n", myRun);
    }
  } // end of run loop 
  printf ("[main] requested runs were processed, bye! \n");
  // plot global historeams 
  return;
}

/////////////////////////////////////////////////////////////////////////////////////////////////

void ResetHistograms(TModule* Setup) {

for (int iMod=0; iMod<6; iMod++) {

  Setup[iMod].ClearPulseShapes();
  Setup[iMod].ClearPulseShapesSplit();
  Setup[iMod].ResetSignalNoise();
  int iCh0=0;
  int iCh1=5;
  if (iMod>3) iCh1=4;
  
  for (int iChan=iCh0; iChan<iCh1; iChan++) {
    Setup[iMod].GetChannel(iChan)->GetHighToLow()->Reset();
    Setup[iMod].GetChannel(iChan)->GetHighVsLow()->Reset();
    Setup[iMod].GetChannel(iChan)->GetFullEventSpectrum()->Reset();
    Setup[iMod].GetChannel(iChan)->GetPedestalSpectrum() ->Reset();
  }

}
}

/////////////////////////////////////////////////////////////////////////////////////////////////
bool ProcessOneRun(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff) {

  RunHeader*   RH=ProcStatus->GetRunHeader();
  unsigned short Block;
  char OnLineDataFile[100]; // reserved for DAQ data file name
  char CherenkovCompressedDataFile[100]; // reserved for Chrerenkov data file name
  sprintf(CherenkovCompressedDataFile, "../cherenkov/compressed/run_%d.dat", Run);

  sprintf(OnLineDataFile, "/raid04/HiLum/data/run_%d.dat", Run);
  indata.open(OnLineDataFile, ios::binary | ios::in);

  if(!indata.is_open()) {
    printf("[ProcessOneRun] Cannot open File %s\n", OnLineDataFile);
    indata.clear();
    sprintf(OnLineDataFile, "../data/run_%d.dat", Run);
    indata.open(OnLineDataFile, ios::binary | ios::in);
  }

  if(!indata.is_open()) {
    printf("[ProcessOneRun] Cannot open File %s\n", OnLineDataFile);
    indata.clear();
    sprintf(OnLineDataFile, "../data0/run_%d.dat", Run);
    indata.open(OnLineDataFile, ios::binary | ios::in);
  }



  if(!indata.is_open()) {
    printf("[ProcessOneRun] Cannot open File %s\n", OnLineDataFile);
    indata.clear();
    sprintf(OnLineDataFile, "../data/run_%d.dat", Run);
    indata.open(OnLineDataFile, ios::binary | ios::in);
    if(!indata.is_open()) {
      printf("[ProcessOneRun] Cannot open File %s\n", OnLineDataFile);
      indata.clear();
      return false;
    }
  }

 
  for (int iMod=0; iMod<6; iMod++) {
    Setup[iMod].ClearTrains();
    //Setup[iMod].ClearPulseShapesSplit();
    Setup[iMod].ClearPulseShapes();
  }

  if (!LoadRFcorrection(Run, Setup)) {
    printf("[ProcessOneRun] failed to load RF correction data for Run#%d\n", Run);
    indata.close();
    indata.clear();
//    exit(-1);
    return false;
  }

  if (!LoadFFTcorrection(Run, Setup)) {
    printf("[ProcessOneRun] failed to load FFT pickup correction data for Run#%d\n", Run);
    indata.close();
    indata.clear();
    return false;
  }
 
  printf("[ProcessOneRun] Processing File %s\n", OnLineDataFile);
  blockCounter=0;
 
  char BeamFitResultsFile[100]; // File with beam counters fit results 
  sprintf(BeamFitResultsFile, "BeamFitResults/beamFit00_%d.dat", Run);
  if (ProcStatus->GetUseFit()) fitdata.open(BeamFitResultsFile, ios::in);


  if(!fitdata.is_open() && ProcStatus->GetUseFit()) {
    printf("[ProcessOneRun] Cannot open File %s\n", BeamFitResultsFile);
    indata.close();
    fitdata.clear();
    return false;
  }

  printf("***********************\nProcessing run#%d attn=%f\n************************\n", myRun, RH->AttnSet);
  runTree->Fill();
  printf("run tree filled\n");
  do {
    Block=Unpack(myRun, Setup, ProcStatus, CherBuff);
    if (ProcStatus->GetReportMode()==0) break; 
  } while (Block != 10 && ProcStatus->GetReportMode()>0  &&            // Continue until RunSummary block found
           ProcStatus->GetReportMode() != 7              &&            // next run
           (ProcStatus->GetEventHeader()->event<5000 || !isOnline) );  // Or limited number of 
                                                                         //events in online mode
  indata.close();
  indata.clear();
  if (ProcStatus->GetUseFit()) fitdata.close();
  printf("[ProcessOneRun] ReporMode=%d\n", ProcStatus->GetReportMode());
  if (ProcStatus->GetReportMode() == 7) ProcStatus->SetReportMode(6);
  if (ProcStatus->GetReportMode()>0) return true; else return false;
}



unsigned int Unpack(int Run, TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff) {

/// unpack the MPI data
/// take action accordingly to record type

  bool Dbg=false;

  unsigned short word;
  unsigned short TWC,BWC,BID;
  bool BBlock=true;

  RunHeader*   RH = ProcStatus->GetRunHeader();
  EventHeader* EH = ProcStatus->GetEventHeader();
  BurstHeader* BH = ProcStatus->GetBurstHeader();

  Dbg=false;
//  Dbg=true;

  if (Dbg) printf("[Unpack] start for Run#%d!\n", Run);

  word=GetWord();
  
  if (Dbg) printf("[Unpack->Entry] got word %d, firstRunEvent=%d\n", word, ProcStatus->IsFirstRunEvent());
//==================================================================================
// Is word a Beginning-of-File marker?
  if ((word==1234) && ProcStatus->IsFirstRunEvent()) {
    if (Dbg) cout <<"[Unpack 1234] Word=" << word << " Found beginning of file"<< endl;
    TWC=GetWord();
    if (Dbg) cout <<"[Unpack 1234] Found beginning of Logical Record with TWC="<<TWC<<" words"<<endl;
    BBlock=false;
    return (0);
  }

// Assume we're at the beginning of a Block
  BWC=word;
  BID=GetWord();

  if(Dbg) cout << "[Unpack] BID=" << BID << " BWC=" << BWC << endl;

//====================================================================================
// Is the Block a RunHeader?
  if(BID==9) {
    if (Dbg) cout <<"[Unpack] Found RunHeader Block";
    RH->run=GetWord();
    if (Dbg) cout << " for runno " << RH->run << endl;
    if (Run != RH->run) {
      printf("[Unpack] Wrong run number in a Run Header\n");
      exit(-2);
    }

    for (int i=0;i<BWC-3;i++) {
      word=GetWord(); 
      printf("[Unpack] RunHeader#%2d %d\n", i, word); 
    }
    word=GetWord();  
    for (int iMod=0; iMod<nMod; iMod++) {
      Setup[iMod].SetSlices(word);
    }
    printf("[Unpack] calling PrepareTimeSlices\n");
    PrepareTimeSlices(Run, Setup, ProcStatus);
    printf("[Unpack] PrepareTimeSlices done\n");


    BBlock=false;

    // load TDC max-min values for the run

    TDCMaxMin[0][0]=235;
    TDCMaxMin[0][1]=1576;
    TDCMaxMin[1][0]=257;
    TDCMaxMin[1][1]=1574;

    nBBlock=0;

    RH->LFcount=0;
    RH->cloneCount=0;
    RH->trainCount=0;
    RH->spikyCount=0;
    RH->grayEvents=0;
    BH->cloneCount=0;
    BH->trainCount=0;
    BH->spikyCount=0;
    
   
  } 
    
//====================================================================================
// Is the Block an EventHeader?
  if (BID==1) {
    Dbg=false;
    if (Dbg) cout <<"[Unpack] Found EventHeader Block";
    EH->run=GetWord();
    EH->Time  = GetWord()*65536+GetWord();
    EH->cTime = ctime(&EH->Time);
    EH->event=GetWord();
    if (Dbg) cout << "[Unpack]  Evt#" << EH->event;
    EH->burstEvent=GetWord();
    EH->burst=GetWord();

    if (EH->burstEvent==0) {
      ProcStatus->SetBeginningOfBurst(true);
      RH->nSpills++;
    }
    EH->pattern[0]=GetWord();
    EH->pattern[1]=GetWord();
    EH->pattern[2]=GetWord();
    EH->pattern[3]=GetWord();
    word=GetWord(); // alignment word
    EH->pattern[4]=GetWord();
    EH->pattern[5]=GetWord();

    if (ProcStatus->GetReportMode()==1) {
      printf("[Unpack]           Run#%d Spill#%04d Event%05d(%03d)\n", 
             EH->run, EH->burst, EH->event, EH->burstEvent);
    }
    //if (EH->run==678 && EH->burst==50 && EH->event==50) ProcStatus->SetReportMode(1);
    //if (EH->run==820 && EH->event>310) {ProcStatus->SetReportMode(1); Dbg=true;}
    //if (EH->run==747 && EH->burst==56) {ProcStatus->SetReportMode(1); Dbg=true;}
    CherBuff->SetCherenkovEventAvailable(true);


    ProcStatus->SetFirstBurstEvent(EH->burstEvent==0);  // sets flag marking First Event in a Spill (just if needed) 

    if (EH->burstEvent==0) {

      //ProcStatus->SetReportMode(1);
      EH->scaler6mhz=0; // will be reloaded, needed here for proper spill integration 
      BH->EventShown=false;
      BH->Time  = EH->Time;
      BH->cTime = EH->cTime;
      BH->burst = EH->burst;
      BH->run   = EH->run;
      printf("\n[Unpack] Run%d Burst#%2d Event#%5d(%d)  TimeBurst= %s", 
               EH->run, EH->burst, EH->event, EH->burstEvent, EH->cTime);
      tm * ptm;
      ptm = gmtime(&EH->Time);
      int dst = 0;
      if (ptm->tm_mon==10) dst=1;
      //printf("Start Spill Moscow Time = %4d-%02d-%02d %2d:%2d:%2d\n", 
      //        1900+ptm->tm_year, ptm->tm_mon+1, ptm->tm_mday+(ptm->tm_hour+4-dst)/24, 
      //        (ptm->tm_hour+4-dst)%24, ptm->tm_min, ptm->tm_sec);

      ProcStatus->SetBunchPatternSet(false);// set to false, will be set to true in LoadSpillSummary;
      if (ProcStatus->LoadSpillSummary(1900+ptm->tm_year, ptm->tm_mon+1, 
                                       ptm->tm_mday+(ptm->tm_hour+4-dst)/24, 
                                       (ptm->tm_hour+4-dst)%24, ptm->tm_min, ptm->tm_sec)) {
        printf("[Unpack] Burst File loaded\n");

        PrepareTimeSlices(Run, Setup, ProcStatus); ///?!

        if (saveSpillProfile) {
          sprintf(myText, "/raid04/HiLum/analysis2012/GADCspillShapes/Run%04dspill%04da.dat", EH->run, EH->burst);
          spilldata.open(myText, ios::out);
          int    nLaps=CherBuff->GetNlaps();
          double tLap0=CherBuff->GetLapTime(0);
          double tLap1=CherBuff->GetLapTime(nLaps-1);
          spilldata << nLaps << " " << tLap0 << " " << tLap1 << " " << CherBuff->GetHVspill() << endl;
          double t0=tLap0; 
          double t1=t0+100000.; // integrate in 0.1mSec intervals
          double w0=0.;
          double wt=0.;
          double wGADC=0.;
          for (int iLap=0; iLap<nLaps; iLap++) {
            double tLap=CherBuff->GetLapTime(iLap);
            if (tLap>t1 && w0>0.0) { // save average rate, reset sums, move to next interval
              wt/=w0; 
              wGADC/=w0;
              spilldata << wt/1.0e6 << " " << wGADC << endl; // save time in msec
              w0=0.;
              wt=0.;
              wGADC=0.;
              t0=t1;
              t1=t0+100000.;
            }
            w0+=1.0; // count laps
            wt+=tLap;
            wGADC+=CherBuff->GetLap(iLap);
          }
          if (w0>0.0) { // save integral, reset sums, move to next interval
            wt/=w0; 
            wGADC/=w0;
            spilldata << wt << " " << wGADC << endl;
          }
          spilldata << tLap0 << " " << 0. << endl;

          spilldata.close();
          spilldata.clear();
        }

        BH->SetGoodSpill(true);
        //exit(-1);
      } else {
        printf("[Unpack] No Burst File for burst#%d run#%d\n", BH->burst, BH->run);
        BH->SetGoodSpill(false);
        //exit(-1);
      }

      if (!isOnline) {//  use Gated Cherenkov and other data in offline only! 


        if (CherBuff->GetSpillFile(1900+ptm->tm_year, 
                                  ptm->tm_mon+1, 
                                  ptm->tm_mday+(ptm->tm_hour+4-dst)/24, 
                                 (ptm->tm_hour+4-dst)%24, ptm->tm_min, ptm->tm_sec)) {
          printf("Gated Cherenkov data loaded\n");

          if (ProcStatus->GetReportMode()==4) {
            ProcStatus->SetReportMode(1);
          }

          BH->goodGatedCherenkov=true;
          BH->CherenkovLoaded=true;
//          double HVnew=CherBuff->GetCherenkovHV(); //now done inside of the Cherenkov Buffer
//          RH->HVbeam[3]=HVnew;
          CherCal->SetCorrection(); // it will use HV from RH->HVbeam[3] !
          double ChHVcorrMPIF=CherCal->GetCorrectionMPIF(); // find HV correction
          double ChHVcorrGADC=CherCal->GetCorrectionGADC(); // find HV correction
          Setup[4].GetChannel(3)->SetGainCorrection(ChHVcorrMPIF); // set gain correction for Cherenkov counter
          CherBuff->SetGainCorrection(ChHVcorrGADC);

          miniBH.burst  = EH->burst;
          miniBH.ChMPIF = (float) ChHVcorrMPIF;
          miniBH.ChGADC = (float) ChHVcorrGADC;
          spillsTree->Fill();

          printf("[Unpack] Cherenkov HV from Spill Summary = %.3f, new Correction=%.3fGADC/track  %.3fMPIF/track\n", 
                 RH->HVbeam[3], 1./ChHVcorrGADC, 1./ChHVcorrMPIF); 
        } else {
          BH->goodGatedCherenkov=false;
          printf("[Unpack] no Gated Cherenkov data available for burst#%d run#%d\n\n", BH->burst, BH->run);
        }


        if (CherBuff->LoadAlignment(BH->run,BH->burst)) { // load MPI vs Gated event alignment data

        } else {
          BH->goodGatedCherenkov=false;
          printf("[Unpack] no Gated Alignment data available for burst#%d run#%d\n\n", BH->burst, BH->run);
        }

        if (//ProcStatus->GetIntensityStatus(EH->burst) != 0 ||
          ProcStatus->GetSummaryStatus(EH->burst)   != 0 ||
          ProcStatus->GetSpillStatusHV2(EH->burst)  != 0 || 
          ProcStatus->GetSpillStatusHC2(EH->burst)  != 0) {
          BH->goodGatedCherenkov=false;
          printf("[Unpack] Spill rejected by STATUS: Intensity=%d Summary=%d HV=%d HC=%d \n\n", 
               ProcStatus->GetIntensityStatus(EH->burst), 
               ProcStatus->GetSummaryStatus(EH->burst), 
               ProcStatus->GetSpillStatusHV2(EH->burst),
               ProcStatus->GetSpillStatusHC2(EH->burst));
        } else if (BH->goodGatedCherenkov) {
        //CherBuff->DisplaySpill();
        //int iqq;
        //cout << "Push!" ;
        //cin >> iqq;
        //if (iqq==0) exit(-1);
        }

        if (showMode) CherBuff->DisplaySpill();
      }

      if (Dbg) cout << " burst time = " << BH->cTime; 
      //Dbg=false;
      EH->scaler10mhzOld=0;
      EH->scaler6mhzOld=0;
    } else {
      EH->scaler10mhzOld=EH->scaler10mhz;
      EH->scaler6mhzOld=EH->scaler6mhz;
    }

    if (EH->pattern[0] != 9 && EH->pattern[0] != 5 && Dbg)cout << " non-beam Pattern  -> "<< EH->pattern[0] << endl;
    if (Dbg) cout << endl;

    if (myRun>300) {
      EH->scaler10mhz=GetWord()*65536+GetWord();
      EH->tms = (double) EH->scaler10mhz / 10000.; // event clock, seconds
      EH->scaler6mhz=GetWord()*65536+GetWord();
      for(int i=0;i<BWC-18;i++) {
        word=GetWord();
      }
    }

    if (Dbg) printf("[Unpack] CherenkovSpillAvailable=%d\n", CherBuff->GetCherenkovSpillAvailable());
    if (CherBuff->GetCherenkovSpillAvailable()) CherBuff->SetEvent6mhz();

    BBlock=false;
    if (Dbg || EH->burstEvent==0) printf("[Unpack] Event Header loaded 6MHz=%d 10MHz=%d\n\n", EH->scaler6mhz, EH->scaler10mhz);

  }

//====================================================================================
// Is the Block a FADC?
  if (BID==36) {
    fFADC(BWC, BID, Setup, EH);
    BBlock=false;
  }

//====================================================================================
// Is the Block a TDC?

  if (BID==39) {
    fTDC(BWC, BID, EH);
    BBlock=false;
  }

//====================================================================================
//====================================================================================
// Is the Block an Event Summary?
//====================================================================================
//====================================================================================
  if (BID==32) {

    //Dbg=true;



    EH->event=GetWord();    
    bool rep1=(ProcStatus->GetReportMode()==1);
    int iAtt=ProcStatus->GetAttenuatorStatus();

    if (Dbg || rep1) cout <<"\n\n[Unpack] Found Event Summary Block  for evtno = " << EH->event << endl;
    // In run 22 the BWC is 5, not 2
    // Change BWC to 3 here
    // to leave the end-of-record markers
    BWC=3;
    for (int i=0;i<BWC-2;i++) {
      word=GetWord();
    }
    BBlock=false;


    // do the Event Analysis here //
    if (Dbg || (ProcStatus->GetReportMode()<10100 && ProcStatus->GetReportMode()>10000)) {
      printf("-------------------------------------\n");
      printf(" [Unpack->Debug] Processing Evt#%d Trig#%d iAtt=%d\n", EH->event, EH->pattern[0], iAtt);
      printf("------------------------------------\n");
    }
    
    bool result=AnalyseEvent(Setup, ProcStatus, CherBuff);

    if (!result && ProcStatus->GetReportMode()==1) {
      ProcStatus->SetReportMode(ProcStatus->WhatToDo("[Main->Bad Event]", -100));
    }
    
    if(ProcStatus->IsFirstRunEvent()) { // do what you need with the first event in the run, drop flag
          /*
          .....
          */
      ProcStatus->SetFirstRunEvent(false);
    }
    if(ProcStatus->IsFirstBurstEvent()) { // do what you need with the first event in the spill, drop flag
          /*
          .....
          */
      ProcStatus->SetFirstBurstEvent(false);
    }


     // now time to save raw amplitudes for a history and do Clone Catching!
    if (Dbg) cout <<"[Unpack] Saving Raw for evt#" << EH->event << "  Good=" << result << endl;
    for (int iMod=0; iMod<6; iMod++) {
      Setup[iMod].SaveRaw();
    }
    if (Dbg) cout <<"[Unpack] Done! evt#" << EH->event << "  Good=" << result << endl;
    Dbg=false;
     
  }
//====================================================================================
// Is the block a Slow Control?
  if (BID==70) {
    if (Dbg) cout <<"Found SlowControl Block"<<endl;
    //if (Dbg) cout <<endl;
    for (int i=0;i<BWC-1;i++) {
      word=GetWord();
    }
    BBlock=false;
  }

//====================================================================================
// Is the block a Temp Probes?
  if (BID==72) {
    if (Dbg) cout <<"Found TempProbes Block" << endl;
    for(int i=0;i<BWC-1;i++) {
      word=GetWord();
    }
    BBlock=false;
  }

//====================================================================================
// Is the block a Purity Monitor?
  if (BID==75) {
    if (Dbg) cout <<"Found PurityMonitor Block"<<endl;
    //if (Dbg) cout <<endl;
    for (int i=0;i<BWC-1;i++) {
      word=GetWord();
    }
    BBlock=false;
  }

//====================================================================================
// Is the block a Burst Summary?
  if (BID==65) {
    printf("[Unpack->Burst Summary]    All Samples -> %1.2e protons\n", CherBuff->GetSpillIntegral());
    printf("[Unpack->Burst Summary] Filled Samples -> %1.2e protons\n", CherBuff->GetSpillIntegralFilled());
    if (Dbg || (ProcStatus->GetReportMode()<10201 && ProcStatus->GetReportMode()>10000)) {
      cout << "[Unpack] Found Burst Summary Block BWC=" << BWC << endl;
      for (int i=0;i<BWC-1;i++) {
        word=GetWord();
        printf(" #%3d %d\n", i, word);
      }
    } else {
      for (int i=0;i<BWC-1;i++) {
        word=GetWord();
      }
    }
    BBlock=false;
    //if (CherBuff->GetHistogramsPrepared() && ProcStatus->GetReportMode()>10000) {
    printf("[Unpack->Burst Summary] Gated Cherenkov HistogramsPrepared=%d\n", CherBuff->GetHistogramsPrepared());
    if (CherBuff->GetHistogramsPrepared() != 0 && showMode) {
      CherBuff->ShowStreamPlots();
      //printf("Push!");
      //int i0=0;
      //cin >> i0;
      //int ReportMode=ProcStatus->WhatToDo("[Unpack]", -200); 
      //ProcStatus->SetReportMode(ReportMode);
    }


    if(ProcStatus->IsFirstRunBurst()) { // do what you need with the first burst (AKA spill) in the run, drop flag
        /*
        .....
        */
        ProcStatus->SetFirstRunBurst(false);
    }
    if (BH->cloneCount>0) printf("[Unpack] *** %d clone events in this spill!\n", BH->cloneCount);
    BH->cloneCount=0;
    if (BH->trainCount>0) printf("[Unpack] *** %d train events in this spill!\n", BH->trainCount);
    BH->trainCount=0;
    if (BH->spikyCount>0) printf("[Unpack] *** %d spiky events in this spill!\n", BH->spikyCount);
    BH->spikyCount=0;
  }

//====================================================================================
// Is the block a Beam Intensity?
  if (BID==66) {
    Dbg=false;
    printf("[Unpack->Beam Intensity]    All Samples -> %1.2e protons\n", CherBuff->GetSpillIntegral());
    printf("[Unpack->Beam Intensity] Filled Samples -> %1.2e protons\n", CherBuff->GetSpillIntegralFilled());
    if (Dbg) cout <<"Found BeamIntensity Block BWC=" << BWC << endl;

// code from JR ->//
    
    union mix_t {double D; 
                 unsigned short sc[4]; 
                 unsigned int ic[2]; 
                 unsigned char c[8];} mix;

    short size_sc;
    short id_sc;
    short size_bd;
    short id_bd;

    int spill;
    int spillBC;
    int scalers[8];
    //int hodoX[16];
    //int hodoY[16];
    int calib;
    int xSEC[16];
    int ySEC[16];
    int xSECt=0;
    int ySECt=0;
    double ionC;

    size_sc=GetWord(); 
    id_sc=GetWord(); 
    spill=(int)GetWord()+65536*(int)GetWord();

    for (int i=0;i<8;i++) {
      scalers[i]=(int)GetWord()+65536*(int)GetWord();
    }
    if(Dbg) {
      cout << "size=" << size_sc << " id=" << id_sc << " spill=" << spill << endl;
      cout <<"scalers="; 
      for (int i=0; i<8; i++) {
        cout << scalers[i] << " ";
      } 
      cout << endl;
    }

    for (int i=0;i<64;i++) { // skip over hodoX and hodoY
      word=GetWord();
    } 
    size_bd=GetWord(); 
    id_bd=GetWord(); 
    spillBC=(int)GetWord()+65536*(int)GetWord();

    if (Dbg) {cout << "size_bd=" << size_bd << " id_bd=" << id_bd << " spillBC=" << spillBC << endl;}
    for (int i=0;i<4;i++) {
      mix.c[2*i+1]=indata.get(); 
      mix.c[2*i]=indata.get();
    }
    ionC=mix.D;

    if (Dbg) cout << "ionC=" << ionC << endl;
    calib=(int)GetWord()+65536*(int)GetWord();

    if (Dbg) cout << "Ion Chamber calibration constant=" << calib << endl;
    for (int i=0; i<16; i++) {
      xSEC[i]=(int)GetWord()+65536*(int)GetWord(); 
      xSECt+=xSEC[i];
    }
    for (int i=0; i<16; i++) {
      ySEC[i]=(int)GetWord()+65536*(int)GetWord(); 
      ySECt+=ySEC[i];
    }
    if(Dbg) {cout << "xSEC="; for(int i=0;i<16;i++) {cout << xSEC[i] << " ";} cout << xSECt << endl;}
    if(Dbg) {cout << "ySEC="; for(int i=0;i<16;i++) {cout << ySEC[i] << " ";} cout << ySECt << endl;}
    int nS=1+4+16+64+4+sizeof(double)/2+2+64;
    for(int i=0;i<BWC-nS;i++) {
      word=GetWord(); // cout << word<<" "; if(((i+1)/16)*16==(i+1)) {cout <<endl;}
    }
    cout << "[Intensity Record]  SpillNo. " << spill << " IonCham " << ionC << " xSEC "<< xSECt << " ySEC "<< ySECt << endl;

// <- code from JR end //

/*
    cout << "Push!";
    int iqq=0;
    cin >> iqq;
    if (iqq==0) exit(-1);
*/
    BBlock=false;
    Dbg=false;
  }

//====================================================================================
// Is the block a Beam Intensity at Run Start?
  if (BID==67) {
    if (Dbg) cout <<"Found Run Start BeamIntensity Block BID="<< BID << " size=" << BWC << endl;
    //cout <<"Found Run Start BeamIntensity Block BID="<< BID << " size=" << BWC << endl;
    for (int i=0;i<BWC-1;i++) {
      word=GetWord();
    }
    BBlock=false;
  }

//====================================================================================
// Is the Block a RunSummary?
  if (BID==10) {
    Dbg=true;
    if (Dbg) cout <<"Found RunSummary Block";
    ProcStatus->GetRunHeader()->run=GetWord(); 
    ProcStatus->GetRunHeader()->eventLast=GetWord();
    if (Dbg) cout << " for runno = " << ProcStatus->GetRunHeader()->run 
                  << " and " << ProcStatus->GetRunHeader()->eventLast << " events" << endl;
    for (int i=0;i<BWC-3;i++) {
      word=GetWord();
    }
    BBlock=false;
    Dbg=false;  
  }

  if ((BWC==65535)&&(BID==65535)) {
    //Dbg=true;
    if (Dbg) cout << "Found end of record " << endl << endl;
    TWC=GetWord();
    if (Dbg) cout << "Found beginning of Logical Record with TWC=" << TWC << " words" << endl;
    BBlock=false;
  }

  // Is the Block unidentified?
  if (BBlock) {
    if (Dbg) cout << "Found unidentified Block No. " << BID << " with BWC = " << BWC << endl;
    for (int i=0;i<BWC-1;i++) {
      word=GetWord();
      //if (Dbg) cout << "   " << word << endl;
    }
    BBlock=false;

    nBBlock++;

    printf(">>> unidentified Block BID= %d  BWC= %d, evtno= %d\n", BID, BWC, EH->event);
    if (nBBlock>10) {
      exit(-2);
    }
  }
  //  printf("just before return(BID), evtno= %d\n",evtno);

  Dbg=false;

  return (BID);
}

bool AnalyseEvent(TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff) {

  bool result=true;
  bool Dbg=false;
  bool rep1=(ProcStatus->GetReportMode()==1);
  int GoodEvent=1;
  //int iAtt=ProcStatus->GetAttenuatorStatus();
  //Dbg=true;

  RunHeader*   RH = ProcStatus->GetRunHeader();
  EventHeader* EH = ProcStatus->GetEventHeader();
  BurstHeader* BH = ProcStatus->GetBurstHeader();

  if (Dbg || rep1) {
    printf("[AnalyseEvent] Run#%d Evt#%d(%d) GoodSpill=%d BunchPattern=%d\n",
           EH->run, EH->event, EH->burstEvent, BH->GetGoodSpill(), ProcStatus->CheckBunchPattern());
  }

  if (!BH->GetGoodSpill())              return false;
  if (!ProcStatus->CheckBunchPattern()) return false;
  if ((EH->scaler6mhz > 0.61*EH->scaler10mhz || EH->scaler6mhz < 0.59*EH->scaler10mhz) && 
      (EH->scaler6mhz!=0 || EH->scaler10mhz!=0) && !ProcStatus->GetOnline()) { // common sense scaler check: fails sometimes!!
    printf("[AnalyseEvent] Burst#%d Evt#%d(%d) corrupt scaler data 6MHz=%d 10MHz=%d\n",  
            EH->burst, EH->event, EH->burstEvent, EH->scaler6mhz, EH->scaler10mhz);
    //int PM=ProcStatus->WhatToDo("Unpack] corrupt scaler data", 2);
    //ProcStatus->SetReportMode(PM); 
    return false;
  }

  GoodEvent = EventPreProcessing(Setup, ProcStatus, CherBuff); // global phase, pedestals, 
                                                     // FFT pickup subtraction
                                                     // Sums for modules
  if (GoodEvent != 1) return false;

  CherBuff->SetdT(EH->GetRFperiod());
  if (ProcStatus->GetReportMode()==1 || Dbg) printf("[AnalyseEvent] EventPreProcessing result=%d\n", GoodEvent);

// check if event is empty or not .//
  bool EventNotEmpty = EH->notEmpty;// based on signal/noise 
 
  if (!EH->notEmpty && !EH->isEmpty) {
    printf("\n**************************************************************\n");
    printf("[AnalyseEvent] -> Event is not full and not empty\n");
    for (int iMod=0; iMod<5; iMod++) {
      if (iMod==RH->AntennaModule) continue;
      printf("[AnalyseEvent] %10s ", Setup[iMod].GetName());
      for (int iChan=0; iChan<4; iChan++) {
        printf("%6.2f ", Setup[iMod].GetChannel(iChan)->GetSignalFactor());
      }
      printf("\n");
    }
    printf("[AnalyseEvent] -> Evt#%d(%d) SignalFactor=%.2f nTrig=%d\n", EH->event, EH->burstEvent, EH->SignalFactor, EH->nTrig);
    printf("**************************************************************\n");
    
    if (ProcStatus->GetReportMode() != 2 && 
        ProcStatus->GetReportMode() != 5) ProcStatus->SetReportMode(1);
    //if (ProcStatus->GetReportMode() !=1 ) GoodEvent=0;
    //if (!RH->S123inBeam) GoodEvent=0;
    RH->grayEvents++;
  }

  if (ProcStatus->GetReportMode() == 1) {  
    int nw=ProcStatus->GetPedWindows();
    printf("[AnalyseEvent] pedestal windows -> %d\n", nw);
    for (int iw=0; iw<nw; iw++) {
      printf("[AnalyseEvent] %d ", iw);
      for (int iMod=0; iMod<4; iMod++) {
        double t0=Setup[iMod].GetChannel(4)->GetSignalLong(1)->GetFirstSample(iw);
        double t1=Setup[iMod].GetChannel(4)->GetSignalLong(1)->GetLastSample(iw);
        printf("| %s->(%6.0f...%6.0f)", Setup[iMod].GetName(), t0, t1); 
      }
      printf("\n");
    }
  }

  if (EventNotEmpty) {
    ProcStatus->SetBeginningOfBurst(false);
  }

  if (GoodEvent!=1)      return false;
  if (EH->burstEvent==0) return false; // avoid some problems with first event GADC and phase //

//  if (CherBuff->hotDisplayRange() && EventNotEmpty) {
//    printf("\n\n[AnalyseEvent] **********> Event is in HOT Gated Display!\n");
//    ProcStatus->SetReportMode(1);
//  }

  if (CherBuff->hotEvent() && !ProcStatus->GetOnline()) {
    printf("[AnalyseEvent] **********> Event#%d is in HOT zone!\n", EH->event);
    //GoodEvent=0;
    //ProcStatus->SetReportMode(1);
    EH->SetHotZone(true);
  } else {
    EH->SetHotZone(false);
  }



  if (EH->isClone) {
    printf("[AnalyseEvent]    Clone!\n");
    GoodEvent=0;
  }
  if (EH->isTrain) {
    printf("[AnalyseEvent]    Train!\n");
    GoodEvent=0;
  }
  if (EH->isSpiky) {
    printf("[AnalyseEvent]    Spiky=%f\n", EH->Spiky);
    //GoodEvent=0;
  }
  if (!BH->goodGatedCherenkov) {
    if (rep1 || Dbg) {
      printf("[AnalyseEvent] -> Skipping event %d of spill%d - no Gated Cherenkov Data found \n", EH->event, EH->burst);
    }
    GoodEvent=0;
  }

  if (GoodEvent==0) {
    if (rep1 || Dbg) printf("[AnalyseEvent] -> Skipping event %d of spill%d as BAD EVENT \n", EH->event, EH->burst);
    return false;
  }

  int CherOver = 0;
  if (BH->goodGatedCherenkov) {
    CherOver = CherBuff->CountOverflows(4.0e6); // count overflows in 4 previous milliseconds
  }
  if (ProcStatus->GetReportMode()==1 || Dbg) printf("[AnalyseEvent] Cherenkov Overflow fraction = %d \n", CherOver);
  //if (CherOver>10) GoodEvent=0;
  miniEH.GADCoverflows=CherOver;
  if (Setup[4].GetChannel(3)->LowGainOverflowDetected()) {
    printf("[AnalyseEvent] **** Cherenkov MPI low gain overflow ******\n");
  }

  if (EventNotEmpty) {

    RH->notEmpty++;

    for (int iMod=0; iMod<6; iMod++) {
      int iCh0=0;
      int iCh1=5;
      if (iMod==4) iCh1=4; // no 'module sums' for scintillators
      if (iMod==5) iCh1=2; // there are only 2 monitor counters
      for (int iCh=iCh0; iCh<iCh1; iCh++) {
        if (Setup[iMod].GetChannel(iCh)->LowGainOverflowDetected()) {
          RH->addLowGainOverflow(iMod, iCh);
        }
      }
    }
  }

  if (Dbg) printf("[AnalyseEvent] Good Event=%d\n", GoodEvent);
  if (!GoodEvent==1) {  // work with good events only! //
    return false;
  }

  if (rep1 || Dbg) printf("[AnalyseEvent] Processing good event : run#%d burst#%d evt#%d(%d)\n",
                               EH->run, EH->burst, EH->event, EH->burstEvent);

  int Nbunches = ProcStatus->GetNbunches();
  int Nfilled  = ProcStatus->GetNfilled();

  if (Dbg || rep1) printf("[AnalyseEvent] %d bunches detected\n", Nbunches);
  if (Nfilled==2) {
    printf("[AnalyseEvent] %d filled bunches detected, %d good for signal\n", Nfilled, Nbunches);
    //ProcStatus->SetReportMode(1);
    //rep1=true;
  }
  if (Nfilled==1) { // do not use "late" bunches in Single Bunch Analysis;
    if (ProcStatus->GetBunchTime(0)>4000.) {
      printf("[AnalyseEvent] late bunch detected, %.0fns\n", ProcStatus->GetBunchTime(0));
      return false;
    }
  }


 // pre-compute times and signals //
  for (int iMod=0; iMod<5; iMod++) {
    Setup[iMod].PreCompute2013();
  }  
/*-------------------------------------------------------------------------------------*/

  if (!CherBuff->CheckFilledSegmet()) {
    printf("[AnalyseEvent] Event#%d(%d) out of filled spill range\n", EH->event, EH->burstEvent);
//    return false;
  }

  miniEH.filledSegment = CherBuff->CheckFilledSegmet();
  miniEH.burst         = EH->burst;
  miniEH.event         = EH->event;
  miniEH.burstEvent    = EH->burstEvent;
  miniEH.scaler10mhz   = EH->scaler10mhz;
  miniEH.scaler6mhz    = EH->scaler6mhz;
  miniEH.isEmpty       = EH->isEmpty; 
  miniEH.notEmpty      = EH->notEmpty;
  miniEH.inHotZone     = EH->inHotZone;
  miniEH.isSpiky       = EH->isSpiky;
  miniEH.RFperiod      = EH->RFperiod;  
  miniEH.RFfreq        = EH->RFfreq ;
  miniEH.RFamp         = EH->RFamp ;
  miniEH.TimeGlobal    = EH->TimeGlobal ;
  miniEH.SignalFactor  = EH->SignalFactor ;




  FillHistograms(Setup, ProcStatus, CherBuff);


  //Dbg=true;
  //if (EH->notEmpty && ProcStatus->GetReportMode()==1 || Dbg) {
  //if (EH->isEmpty && ProcStatus->GetReportMode()==1 || Dbg) {
  if (ProcStatus->GetReportMode()==1 || Dbg) {
    //ProcStatus->SetReportMode(1);
    printf("[AnalyseEvent] Calling DisplayAllPulses2012     GoodEvent=%d\n", GoodEvent);
    DisplayModules(Setup, ProcStatus, EventCanvasModule);
    iqq=DisplayAllPulses2012(Setup, ProcStatus, CherBuff);
    if (RH->S123inBeam) {
      iqq=DisplayAllCounters2012(Setup, ProcStatus, CherBuff);
    }
    //iqq=DisplayForMaslennikov(Setup, ProcStatus, CherBuff);
    //iqq=DisplayAllPulsesDetails(Setup, ProcStatus, CherBuff);
    iqq=ProcStatus->WhatToDo("[AnalyseEvent]", -100);
    ProcStatus->SetReportMode(iqq);
  }
  Dbg=false;
  return result;

}

//------------------------------------------------------------------------------------------------
int FillHistograms (TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff) {
// 2013 version, uses spectrum hists from 'inside' of TChannel objects
//
EventHeader* EH = ProcStatus->GetEventHeader();
RunHeader*   RH = ProcStatus->GetRunHeader();
double T00 = EH->GetTimeGlobal();
double Trf = EH->GetRFperiod();
//int    iOff= ProcStatus->GetEventHeader()->GetFirstBunch();
bool notEmpty = EH->notEmpty;// based on signal/noise for signals
bool isEmpty  = EH->isEmpty;
bool DumpIt   = false;
int  iModX    = RH->AntennaModule;
int  FillResult=0;
//int Nbunches=ProcStatus->GetNbunches(); // get number of bunches visible in the event

evtEssentials.Clear();

//double Rate =0.;
//double Rate1=0.;
//double Rate10=0.;
//double Rate01=0.;
int iGC=CherBuff->GetSampleIndex(0);
evtEssentials.firstGADCsample=iGC;
miniEH.firstGated=iGC;

double T1=EH->scaler10mhz*100; // ns
//double T01=T1-1.0e5;            // 0.1ms
//double T0 =T1-1.0e6;            // 1ms
//double T10=T1-1.0e7;            // 10ms
double T01=T1-1.0e6;            // 1ms
double T0 =T1-1.0e7;            // 10ms
double T10=T1-1.0e8;            // 100ms



//-> retrieve GADC amplitudes for the event 

double tRF=0.;
for (int ir=0; ir<40; ir++) {
  evtEssentials.inScope[ir]   =false; 
  evtEssentials.isIsolated[ir]=false;
  evtEssentials.isCluster0[ir]=false;
  evtEssentials.isCluster1[ir]=false;
  evtEssentials.bunchAction[ir]  =0;
  evtEssentials.GADCevent[ir]    =0.0;
  evtEssentials.GADCeventTail[ir]=0.0;
  evtEssentials.MPICeventHG[ir]    =0.0;
  evtEssentials.MPICeventLG[ir]    =0.0;
}

for (int ir=0; tRF<6200. && ir<40; ir++) {

  tRF=T00+ir*Trf;
  evtEssentials.tRFb[ir]=tRF;
  evtEssentials.GADCevent[ir]=CherBuff->GetSampleAdjustD(ir);         // dynamic pedestal
  evtEssentials.GADCeventTail[ir]=CherBuff->GetSampleAdjustDtail(ir); // dynamic pedestal with tail subtraction
  evtEssentials.inScope[ir]=true;
  if (ProcStatus->GetReportMode()==1) {
    printf("[FillHistograms] ir=%2d tRF=%5.1f ", ir, tRF);
    for (int iCh=0; iCh<4; iCh++) {
      printf(" %2s->%d", Setup[4].GetChannel(iCh)->GetName(), 
                         Setup[4].GetChannel(iCh)->GetClusterPMTbunchSET(ir));
    }
  }
  for (int iCh=0; iCh<4; iCh++) {
    TChannel *Channel=Setup[4].GetChannel(iCh);
    evtEssentials.inScope[ir]=evtEssentials.inScope[ir] && Channel->GetClusterPMTbunchSET(ir);
  }
  if (ProcStatus->GetReportMode()==1) {
    printf(" inScope=%d", evtEssentials.inScope[ir]);
  }
  if (!evtEssentials.inScope[ir]) {
    if (ProcStatus->GetReportMode()==1) printf("\n");
    continue;
  }

  evtEssentials.MPICeventHG[ir]=Setup[4].GetChannel(3)->GetClusterPMTbunchAmp2(ir, 1); // MPI cherenkov
  evtEssentials.MPICeventLG[ir]=Setup[4].GetChannel(3)->GetClusterPMTbunchAmp2(ir, 0); // MPI cherenkov

  int iGC=CherBuff->GetSampleIndex(ir);   //      fing it's position in Gated record
  int iAct=ProcStatus->GetBunchAction(iGC%30); // find if it is full or empty
  int iAc00=ProcStatus->GetBunchAction((iGC-1)%30);
  int iAc01=ProcStatus->GetBunchAction((iGC-2)%30);
  int iAc10=ProcStatus->GetBunchAction((iGC+1)%30);
  int iAc11=ProcStatus->GetBunchAction((iGC+2)%30);
  evtEssentials.isIsolated[ir]=(iAct==1 && iAc00!=1 && iAc10!=1 && iAc01!=1 && iAc10!=1);
  evtEssentials.isCluster0[ir]=(iAct==1 && iAc00!=1 && iAc01!=1);
  evtEssentials.isCluster1[ir]=(iAct==1 && iAc10!=1 && iAc11!=1);
  evtEssentials.bunchAction[ir]=iAct;
  EH->SetBunchAction(ir, iAct);

  if (ProcStatus->GetReportMode()==1) {
    printf(" pattern=%d%d>%d<%d%d isol=%d clust0=%d clust1=%d\n", 
            iAc01, iAc00, iAct, iAc10, iAc11, 
            evtEssentials.isIsolated[ir],
            evtEssentials.isCluster0[ir],
            evtEssentials.isCluster1[ir]);
  }
  

}



int iClust=0; // filled cluster, single bunch or a group of 
double t0, t1;

for (int ir=0; evtEssentials.tRFb[ir]<6200. && ir<40; ir++) {
  
//printf("ir=%2d inScope=%d\n", ir, evtEssentials.inScope[ir]);

  if (!evtEssentials.inScope[ir]) continue;

  if (iClust==5) {
    printf("\n\n[FillHistograms] too many clusters!\n\n");
    ProcStatus->SetReportMode(1);
    exit(-1);
    return 1;
  }

  if (ProcStatus->GetReportMode() != 2  ||  DumpIt) {

    printf("[FillHistograms] %2d isol=%d%d%d", ir, 
            evtEssentials.isIsolated[ir], 
            evtEssentials.isCluster0[ir], 
            evtEssentials.isCluster1[ir]);
    for (int iCh=0; iCh<4; iCh++) {
      double beamMIP=1.0;
      if (iCh<3) beamMIP=myMIPcounter->GetMIPaverage(iCh);
      printf("|%2s A=%7.1f(%d) ", 
                  Setup[4].GetChannel(iCh)->GetName(),
                  Setup[4].GetChannel(iCh)->GetClusterPMTbunchAmp(ir)/beamMIP,
                  Setup[4].GetChannel(iCh)->GetClusterPMTgain(ir));
    }
    printf("|GADC=%5.1f", evtEssentials.GADCevent[ir]);

    t0=evtEssentials.tRFb[ir];
    if (evtEssentials.isIsolated[ir]) {
      for (int iMod=0; iMod<4; iMod++) {
        if (iMod==iModX) continue;
        bool good=true;
        double tInt=Setup[iMod].GetIntegralTime();
        double modIntegral=Setup[iMod].GetIntegral(4, 1, t0, t0+tInt+12.5, good);
        if (good) printf("|%5s A=%7.1ffC ", Setup[iMod].GetName(), modIntegral);
      }
    } else if (evtEssentials.isCluster0[ir]) {
      t0=evtEssentials.tRFb[ir];
      t1=-1000.;
      for (int ir1=ir+1; evtEssentials.tRFb[ir1]<6200. && ir1<40; ir1++) {
        if (evtEssentials.isCluster1[ir1]) {
          t1=evtEssentials.tRFb[ir1];
          break;
        }
      }
      if (t1>-1000.) { // end of cluster found
        for (int iMod=0; iMod<4; iMod++) {
          if (iMod==iModX) continue;
          bool good=true;
          double tInt=Setup[iMod].GetIntegralTime();
          double modIntegral=Setup[iMod].GetIntegral(4, 1, t0, t1+tInt+12.5, good);
          if (good) printf("|%5s A=%7.1ffC ", Setup[iMod].GetName(), modIntegral);
        }
        
      } 
    }
    printf("\n");
  } // end of dump

  t0=evtEssentials.tRFb[ir];
  t1=-1000.;
  int ir0=ir;

  if (evtEssentials.isCluster0[ir]) { // signal cluster found
    evtEssentials.RFclust0[iClust]=ir;
    while(!evtEssentials.isCluster1[ir0] && ir0<40) {ir0++;}
    evtEssentials.RFclust1[iClust]=ir0;
    if (ir0<40) t1=evtEssentials.tRFb[ir0];
  } else {
    continue;
  }

  if (t1==-1000.) {
    printf("[Fill Histogram] Closing bunch not found for cluster starting at bunch#%d t=%.1f\n",
           ir, t0);
    ProcStatus->SetReportMode(1);
    continue;
  }

  int id=0;
  bool good=true;

  for (int iMod=0; iMod<4; iMod++) {
    if (iMod==iModX) continue; // skip the antena module
    double tInt=Setup[iMod].GetIntegralTime();

    evtEssentials.HGoverflow[iMod]=Setup[iMod].CheckHGoverflow();
    evtEssentials.LGoverflow[iMod]=Setup[iMod].CheckLGoverflow();

    t1=evtEssentials.tRFb[evtEssentials.RFclust1[iClust]]+Setup[iMod].GetIntegralTime()+12.5;
 
    for (int iChan=0; iChan<5; iChan++) {
      for (int iGain=0; iGain<2; iGain++) {
        id=iGain+2*iChan+10*iMod+40*iClust;
        evtEssentials.miniModule[id]=Setup[iMod].GetIntegral(iChan, iGain, t0, t1, good);
      }
    }
  } // end of module loop
  iClust++;
} // end of bunch loop


eventsTree->Fill();// fill branch for event header;
eventsTreeModules->Fill();
eventsTreeCounters->Fill();
//if (ProcStatus->GetReportMode()==1) eventsTree -> Print();



FillResult=1;
return FillResult;
}// end of FillHistogram




//===========================================================================================
unsigned int EventPreProcessing (TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer* CherBuff) {
// compute event phase, trains, pedestals, empty/full  where applicable
EventHeader* EH=ProcStatus->GetEventHeader();
BurstHeader* BH=ProcStatus->GetBurstHeader();
RunHeader*   RH=ProcStatus->GetRunHeader();
int trigger=EH->pattern[0];
int goodEvent=1;
bool Dbg=false;
bool res1=ProcStatus->GetReportMode()==1;
int iModX=RH->AntennaModule;

EH->isClone=false;
EH->isTrain=false;
EH->isSpiky=false;
EH->S123inBeam = RH->S123inBeam;

if (Dbg) printf("[EventPreProcessing] Evt#%d pattern=%d was loaded; \n", EH->event, trigger);

if (trigger != 9 && trigger != 5) {
  if (Dbg || res1) printf("[EventPreProcessing] Evt#%d pattern=%d was loaded; \n", EH->event, trigger);
  goodEvent=0;
}


for (int iMod=0;iMod<nMod;iMod++) {
  Setup[iMod].CountTrains();
  if (Setup[iMod].IsClone()) {
    if (ProcStatus->GetReportMode()<10100 && ProcStatus->GetReportMode()>10000) {
      printf("[main->EventPreProcessing] Evt#%d found Clone signal in %s\t with %d reprtitions, skipping the analysis\n", 
              EH->event, Setup[iMod].GetName(),  Setup[iMod].GetIdentical());
    }
    EH->isClone=true;
    goodEvent=-2;// throw away the Clones!!
      //break;
  } else if (Setup[iMod].IsTrain()) {
    printf("------Evt#%d\t %s\t total matches: %3d Train!\n", 
           EH->event, Setup[iMod].GetName(), Setup[iMod].GetIdentical());      
    EH->isTrain=true;
    goodEvent=-2;// throw away the Trains!!
  }
}

if (EH->isClone) { // clones invading !!!
  BH->cloneCount++;
  RH->cloneCount++;
  return goodEvent;
} else if (EH->isTrain) { // trains found !!!
  BH->trainCount++;
  RH->trainCount++;
  return goodEvent;
}  
  
bool GoodPhase = false;   double Tg=0.0;

if (EH->run<852) {
  Tg = GlobalPhase2011(Setup, ProcStatus)-50.; // old
  while (Tg>0.0) Tg-=165.00;
  EH->SetTimeGlobal(Tg);
  EH->SetRFperiod(165.0);
  GoodPhase=true;
} else {
  GoodPhase=GlobalPhase2012(Setup, ProcStatus); // finder of an RF phase
}

if (!GoodPhase) {
  ProcStatus->SetReportMode(1);
  printf("[main->EventPreProcessing] Run#%d Event#%d Bad Phase detected GoodPhase=%d\n", 
          EH->run, EH->event, GoodPhase);
//  exit(-1);
//  ProcStatus->SetReportMode(1);
  return 0;
}

Tg = EH->GetTimeGlobal();
 
for (int iMod=0; iMod<6; iMod++) {
  if (iMod==4 || iMod==5) { // apply RF correction for the AZ shaper channels// 
    Setup[iMod].SubtractRFcross();
  }
  Setup[iMod].ComputePrePed();
} 


int iGC=CherBuff->GetSampleIndex(0); // position of event's sample #0 in Gated Cherenkov buffer
//int iGC30=iGC%30;                    // bunch alignment in the event
// compute SpillPattern noises and signals for all detectors
for (int iMod=0; iMod<6; iMod++) {
  if (iMod==iModX) continue; // do not spend time in disabled module
  int iCh0=0;
  int iCh1=4;
  if (iMod==5) {
    iCh1=2;
  }
  for (int iCh=iCh0; iCh<iCh1; iCh++) {
    TChannel *Channel=Setup[iMod].GetChannel(iCh);
    bool zz=Channel->ComputeVarianceByBunch();
    if (!zz) {
      printf("[Main->EventPreProcessing] ComputeVarianceByBunch FAILED\n");
      exit(-1);
    }
    for (int ib=0; ib<30; ib++) {
      Channel->ComputeRFsig(ib);
    }
  }
}

//printf("[Main->EventPreProcessing] scintillators in beam = %d\n",  RH->S123inBeam);

double SignalPresent[6][4];
for (int i=0; i<6; i++) {
  for (int j=0; j<4; j++) SignalPresent[i][j]=0.0;
}
  //for (int iMod=0; iMod<6; iMod++) {
int iMod0=0;
int iMod1=5;
//if (RH->S123inBeam) iMod0=4; // no modules for self-trigger at low intensity
int nTrig=0;
double SumAvr=0.0;
double SigSum1=0.0;
double SigMax=-1.0e10;
double SigMin=1.0e10;
//double trigCut[5]={0.5, 0.5, 0.5, 0.5, 0.5};
double trigCut[5]={0.825, 0.825, 0.825, 0.825, 0.825};

for (int iMod=iMod0; iMod<iMod1; iMod++) { //use FCal-in-beam, Ch and S123 (if in beam)//
  if (iMod==iModX) continue; // do not spend time in disabled module
  int iCh0=0;
  int iCh1=4;
  if (iMod==5) {
    iCh1=2; // 
  }
  for (int iCh=iCh0; iCh<iCh1; iCh++) {
      
    TChannel *Channel=Setup[iMod].GetChannel(iCh);
    int    iMax=0;
    SumAvr=0.0;
    SigMax=-1.0e10;
    SigMin=1.0e10;
    for (int ib=0; ib<30; ib++) {
      SigSum1=0.0;
      for (int iGain=1; iGain<2; iGain++) {
        for (int iDelay=0; iDelay<2; iDelay++) {
          TSignal *Signal=Channel->GetSignalShort(2*iGain+iDelay);
          SigSum1+=Signal->GetRFdet(ib)/Signal->GetRFsig(ib); 
        }
        //if (res1 && (iMod==4 || iMod==2) && iCh==3) printf("\n");
      }           
      if (SigSum1>SigMax) {iMax=ib; SigMax=SigSum1;}
      if (SigSum1<SigMin) SigMin=SigSum1;
      SumAvr+=SigSum1;
    }
    SumAvr=(SumAvr-SigMax)/29.;
    double signalFactor=(SigMax-SigMin)/(SigMax+SigMin);
    Channel->FillSignalFactorTH1D(signalFactor);
    SignalPresent[iMod][iCh]=signalFactor;
    Channel->SetSignalFactor(signalFactor);

    if (signalFactor>trigCut[iMod]) {
      if (RH->S123inBeam) {
        if (iMod==4) nTrig++; // only scintillators & cherenkov
      } else {
        if (iMod<4)  nTrig++;                 // modules and
        else if (iMod==4 && iCh==3) nTrig++; // cherenkov
      }
    }
    if (Dbg || res1) {
      printf("[Main->EventPreProcessing] %10s ", 
             Setup[iMod].GetChannel(iCh)->GetName());
      printf(" iMax=%02d Smin=%6.1f Smax=%6.1f Savr=%6.1f Factor=%.2f signal=%d", 
              iMax, SigMin, SigMax, SumAvr, signalFactor, signalFactor>trigCut[iMod]);
      if (iMod<4 && RH->HVset[iMod]!=0. && !RH->S123inBeam) {
        printf(" used\n");
      } else if (iMod==4 && RH->S123inBeam && iCh<3) {
        printf(" used\n");
      } else if (iMod==4 && iCh==3) {
        printf(" used\n");
      } else {
        printf("\n");
      }
    } 
  }
}

double SumMax=-1.0e10;
double SumMin=1.0e10;
SumAvr=0.0;
iMod0=0;
iMod1=5;
if (RH->S123inBeam) iMod0=4; // only Ch and S123 if in beam
int    iMaxSum=-1;
//int ib=iGC30;
//if (ib<0) ib=0;
if (Dbg || res1)printf("[Main->EventPreProcessing] Sig/Noise-> ");
for (int ib=0; ib<30; ib++) {
  double SumDet=0.0;
  for (int iMod=iMod0; iMod<5; iMod++) { //use FCal-in-beam, //
    if (iMod==iModX) continue; // do not spend time in disabled module
    if (!ProcStatus->GetOnline()) {
      if (iMod<4 && RH->HVset[iMod]==0.) continue;
    }
    int iCh0=0;
    int iCh1=4;
    if (iMod==4 && !RH->S123inBeam) { // S123 not in the beam, cherenkov only
      iCh0=3; 
    }
    for (int iCh=iCh0; iCh<iCh1; iCh++) {
      TChannel *Channel=Setup[iMod].GetChannel(iCh);
      for (int iGain=1; iGain<2; iGain++) {
        for (int iDelay=0; iDelay<2; iDelay++) {
          TSignal *Signal=Channel->GetSignalShort(2*iGain+iDelay);
          SumDet+=Signal->GetRFdet(ib)/Signal->GetRFsig(ib);
        }
      }
    }
  }
  SumAvr+=SumDet;
  if (SumDet>SumMax) {
    SumMax=SumDet;
    iMaxSum=ib;
  }
  if (SumDet<SumMin) {
    SumMin=SumDet;
  }     
  if (Dbg || res1)printf(" %6.1f ", SumDet);
}
SumAvr-=SumMax;
SumAvr/=29.;
if (Dbg || res1)printf("\n[Main->EventPreProcessing] Min=%6.1f Max=%6.1f Avr=%6.1f\n", SumMin, SumMax, SumAvr);
  
double signalFactor=0.0;
if (SumMax+SumMin>0.0) {
  signalFactor=(SumMax-SumMin)/(SumMax+SumMin);
}
if (!RH->S123inBeam && nTrig>4 && signalFactor<SignalPresent[4][3]) signalFactor=SignalPresent[4][3];// cherenkov 

if (Dbg || res1) {
  printf("[Main->EventPreProcessing] RawNoise Global iMax=%d  Sf=%.2f nTrig=%d\n", iMaxSum, signalFactor, nTrig);
} 

EH->SetNtrig(nTrig);
EH->SetSignalFactor(signalFactor);

if (iMaxSum==-1) {
  printf("/n/n************/n[Main->EventPreProcessing] Hottest Bunch not found\n************\n"); 
  ProcStatus->SetReportMode(1);
  iMaxSum=0; 
}

int iMold=iMaxSum;
iMaxSum=iGC%30; 

if (EH->notEmpty && iGC%30 != iMaxSum && !ProcStatus->GetOnline()) {
  printf("[Main->EventPreProcessing] evt#%d(%d) GADC offset failed iMaxSum=%d iGated=%d\n", 
          EH->event, EH->burstEvent, iMold, iGC%30); 
  EH->Print();  
  if (EH->burstEvent>0) ProcStatus->SetReportMode(1);
  //ProcStatus->SetReportMode(1);
  //return 0
}

if (iMaxSum<0) {
  printf("[Main->EventPreProcessing] evt#%d(%d)  GADC offset failed iMaxSum=%d iGated=%d(%d) SignalFactor=%.2f Ntrig=%d\n", 
          EH->event, EH->burstEvent, iMaxSum, iGC, iGC%30, signalFactor, nTrig); 
  if (!EH->notEmpty) iMaxSum=0; // do not care
  else if (EH->burstEvent>0) return 0; // skip bad event
  else               iMaxSum=0; // event#0...            
}

EH->SetFirstBunch(iMaxSum);

//printf("[EventPreProcessing] burst#%d burstEvent=%d notEmpty=%d ReportMode=%d  iMaxSum=%d iGated=%d\n",
//       EH->burst, EH->burstEvent, EH->notEmpty, ProcStatus->GetReportMode(), iMaxSum, iGC%30);

if (EH->burstEvent==0 && ProcStatus->GetReportMode()==6) {
  ProcStatus->SetReportMode(5);
}

if (EH->notEmpty && ProcStatus->GetReportMode()==5) {
  ProcStatus->SetReportMode(1);
}



ProcStatus->PrepareBunches();
SignalFactor1->Fill(signalFactor);

//for (int iMod=0; iMod<6; iMod++) {
for (int iMod=0; iMod<5; iMod++) { // skip the Monitor, who cares?
  Setup[iMod].SetNbunchs(ProcStatus->GetNbunches());
  Setup[iMod].SetPedestalWindowsNew();
  if (Setup[iMod].GetNPedWindows()<2) goodEvent=-10;
}

if (goodEvent==-10) {
  printf("[EventPreProcessing] not enough pedestal windows in event!\n");
  for (int iMod=0; iMod<6; iMod++) {
    printf("  %10s->%d", Setup[iMod].GetName(), Setup[iMod].GetNPedWindows());
  }
  printf("\n-------------------------------------------------------\n");
  return goodEvent;
}

if (EH->isEmpty) { // use empty event for static pedestal estimation
  if (EH->burstEvent==0) {
    //printf("[Main->EventPreProcessing] clearing empty event pedestals\n");
    for (int iMod=0; iMod<6; iMod++) {
      Setup[iMod].ClearEmptyEventPedestal();
    }
  }
  if (ProcStatus->GetBeginningOfBurst()) {
    if (BH->GetGoodSpill()) RH->PedestalEvents++; // count pedestals in good spills only!
    if (ProcStatus->GetReportMode()==1) {
      printf("[Main->EventPreProcessing] updating empty event pedestals\n");
    }
    for (int iMod=0; iMod<6; iMod++) {
      Setup[iMod].UpdateEmptyEventPedestal();
    }
  }
}  // end of empty event handling //

for (int iMod=0;iMod<6;iMod++) {
  if (!Setup[iMod].CheckEmptyEventPedestals()) {
    if (Dbg || res1) {
      printf("[Main->EventPreProcessing] No empty event pedestals for %s\n", Setup[iMod].GetName());
    }
    goodEvent=-4;
    return goodEvent;
  }
  if (Dbg) {
    printf("[Main->EventPreProcessing] Subtracting empty event pedestals in %s %d\n", 
              Setup[iMod].GetName(), iMod);
  }
  Setup[iMod].SubtractPedestals();
  if (Dbg) { 
    printf("[Main->EventPreProcessing] Setting LONG pulses in %s\n", Setup[iMod].GetName());
  }
  Setup[iMod].SetPulses();  
}


TModule  *ModuleX = &Setup[iModX];
//if (ProcStatus->GetReportMode() != 2) {
//  printf("[Main->EventPreProcessing] using %s as Antena\n", ModuleX->GetName());
//}
ModuleX->ComputeSumPulse(); // compute Sums for Antenna signal
TChannel *ChannelX = ModuleX->GetChannel(4); // get the #4 Sum channel of Anetnna module
TSignal  *SignalX  = ChannelX->GetSignalLong(1); // access HG Sum signal of the Antenna 
SignalX->StraightFFTransform();
SignalX->ComputeFFTAmpPhase();
//SignalX->FilterFFT(500., 0.);
double LF_Frac=SignalX->GetLowFrequencyFFT(5); // checks fraction of a^2 in bins 1,2,3,4

if (LF_Frac>.35) {
  printf("[Main->EventPreProcessing] ************************************************\n");
  printf("[Main->EventPreProcessing] Evt#%d(%d) %s  LF excess =%f.3\n\t", 
        EH->event, EH->burstEvent, ModuleX->GetName(), LF_Frac);
  for (int ib=1; ib<10; ib++) {
    printf(" %6.1f ",  SignalX->GetFFTamp(ib));
  }
  printf("\n");
  printf("[Main->EventPreProcessing] ************************************************\n");
  RH->LFcount++;
  return 0;
}

iMod0=0; // MIP stage 2

for (int iMod=iMod0; iMod<6; iMod++) {
  if (iMod==iModX) continue; // skip the "antena" FCal Module
  int iCh0=0;
  int iCh1=4;
  if (iMod==4) iCh1=4; // beam
  if (iMod==5) iCh1=2; // only 2 monitors are there
  if (ProcStatus->GetReportMode() == 1) {
    printf("[Main->EventPreProcessing] AttnSet=%.1f Module->%s Antena->%s\n", 
            RH->AttnSet, Setup[iMod].GetName(), Setup[iModX].GetName());
  }
  if (RH->AttnSet < 30.0 || iMod>3) Setup[iMod].SubtractCoherentNoiseFFT(SignalX, iCh0, iCh1, 1.5);
  Setup[iMod].PrepareAmpI();    // prepare interpolated long signal ampI
  if (iMod<4) Setup[iMod].ComputeSumPulse();
  Setup[iMod].ComputeFFTAmpPhase(); // prepare amplitude and pfase for following analysis
}   // module
  


if (EH->isEmpty) { // use empty event for static pedestal estimation in LONG signals
  if (ProcStatus->GetBeginningOfBurst()) {
    //printf("[Main->EventPreProcessing] updating empty event LONG pedestals\n");
    for (int iMod=0; iMod<6; iMod++) {
      Setup[iMod].UpdateEmptyEventPedestalLong();
    }
  }
}  // end of empty event handling //double rawNoise=Setup[iMod].ComputeRawNoise();


//CalculateIntegrator(Setup, ProcStatus);

// seek for spikes in signals
double maxFuzz=0.0;
for (int iMod=0; iMod<6; iMod++) {
  TModule *Module = &Setup[iMod];
  for (int iChan=0; iChan<4; iChan++) {
    for (int iGain=0; iGain<2; iGain++) {
      TChannel *Channel = Module->GetChannel(iChan);
      TSignal *Signal = Channel->GetSignalLong(iGain);
      if (Signal->GetStatus()) {
        maxFuzz+=Signal->GetFuzz();  
      }
    }
  }
}

EH->Spiky   =  maxFuzz;
if (maxFuzz>1.1) {
  //goodEvent=-3;// throw away the Spikes !!
  EH->isSpiky=true;
  BH->spikyCount++;
  RH->spikyCount++;
} else {
  EH->isSpiky=false;
}
return goodEvent;
}

//===================================================================================================
void FitModule(TModule* Setup, int Condition, ProcessingStatus* ProcStat) {
  int i1=Condition/5;
  int iChan=Condition-i1*5-1;
  int bunch=i1/10;
  int iMod=i1-bunch*10;
  printf("[FitModule] Calling Module Fit routine for Module#%d Channel#%d bunch#%d\n",
	 iMod, iChan, bunch);
  Setup[iMod].Fit(ProcStat, bunch, iChan);
  return;
}


void  CalculateIntegrator(TModule *Setup, ProcessingStatus* ProcStatus){
  bool Dbg=false;
  TModule *Module = &Setup[6];
  EventHeader* EH = ProcStatus->GetEventHeader();
  double gc=Setup[4].GetChannel(3)->GetGainCorrection(); // conversion of High Gain ADC to Single Particles
  double gc1=gc*Setup[4].GetChannel(3)->GetHighLow(); // conversion of Low Gain ADC to Single Particles
  if (Dbg) printf("[CalculateIntegrator] calculating for %s\n", Module->GetName());
  for (int iChan=0; iChan<4; iChan++) {
    TSignal *Signal = Module->GetChannel(iChan)->GetSignalShort(0);
    if (Dbg) printf("[CalculateIntegrator] trying %s -> ", Signal->GetName());  
    if (Module->GetChannel(iChan)->GetStatus()) {
      if (Signal->CalculateIntegrator()) {
        if (EH->pattern[0]==18 || 
             (EH->run>329 && 
              EH->burstEvent<2 && 
              EH->isEmpty)) { // or early event with no signal in runs
                              // with extended integrator gate
          Signal->SetIntegratorPedestal();
          if (Dbg || (ProcStatus->GetReportMode()<10101 && ProcStatus->GetReportMode()>10000)) 
             printf("[CalculateIntegrator] Integrator Pedestals Set!\n"); 
        }
	Signal->SubtractIntegratorPedestal();
        double integrator=Signal->GetIntegratorLevel();
        EH->integ[iChan]=integrator;     // raw integrator, low gain ADC
        EH->integC[iChan]=integrator*gc1; // integrator now in protons
        if (integrator*gc1>0.01) {
          EH->integL[iChan]=log(integrator*gc1)/log(10.0); // Log10 HV corr integrator
        } else {
          EH->integL[iChan]=-2.0;
        }
        if (Dbg) printf("  integrator#%d I=%.2f\n", iChan, Signal->GetIntegratorLevel());
      } else {
	printf("[CalculateIntegrator]  integrator calculation failed!\n Chan#%d\n", iChan);
	exit(-1);
      }
    } else {
      printf(" channel#%d %s has status 0\n", iChan, Signal->GetName());
    }
  }
  if (Dbg) {
    for (int iChan=0; iChan<4; iChan++) {
      printf("  integrator#%d I=%.2f\n", iChan, EH->integ[iChan]);
    }
  } 
  Dbg=false;
  return;
}



void RunReport(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus, CherenkovBuffer *CherBuff) {
  
  if (ProcStatus->GetReportMode()!=1) {
    EventCanvasNew->Update();
    EventCanvas->Update();
    EventCanvasNew->Iconify();
    EventCanvas->Iconify();
  }
   // plot and save run historams ETC

  RunHeader   *RH = ProcStatus->GetRunHeader();
  EventHeader *EH = ProcStatus->GetEventHeader();
  int iModX=RH->AntennaModule; // save Sums for the FFT correction;
  //int dWidth  = ProcStatus->GetDisplayWidth();
  char titleC[180];
  //char nameC[180];

  TH1D    *dummy   = new TH1D("dummy",   "dummy",   2500, -50., 8000.);
  TH2D    *dummy2D = new TH2D("dummy2D", "dummy2D", 200, -1., 1., 200, -1., 1.);
  TLatex  *TL = new TLatex; 
  dummy  ->SetStats(0);
  dummy2D->SetStats(0);
  dummy2D->GetXaxis()->Delete();
  dummy2D->GetYaxis()->Delete();
  printf("[RunReport] doing run report for Run#%d\n", Run);
  char padName[40];
  char reportLine[100];


//============================================================================== //
  PilotCanvas->Clear();
  PilotCanvas->Divide(4,4);
 
  int iZone=1;
  for (int iMod=0; iMod<5; iMod++) {
    if (iMod==iModX) continue;
    for (int iChan=0; iChan<4; iChan++) {
      PilotCanvas->cd(iZone++);
      Setup[iMod].GetChannel(iChan)->GetSignalFactorTH1D()->DrawCopy();
    }
  } 
  PilotCanvas->Update();
  PilotCanvas->Show();
//============================================================================== //
 
  PhaseCanvas->Clear();
  PhaseCanvas->Divide(2,3);
  //PhaseCanvas->Iconify();

  PhaseCanvas->cd(1);
  gStyle->SetOptStat(111110);
  RFt0[0].DrawCopy();
  RFt0[1].DrawCopy("same");

  PhaseCanvas->cd(2);
  gStyle->SetOptStat(111110);
  RFt0diff->DrawCopy();

  PhaseCanvas->cd(3);
  gStyle->SetOptStat(111110);
  double amin=1.0e11;
  double amax=-1.0e11;
  double a1=RFamplitude[0].GetMean()+5.0*RFamplitude[0].GetRMS();
  if (a1>amax) amax=a1;
  a1=RFamplitude[1].GetMean()+5.0*RFamplitude[1].GetRMS();
  if (a1>amax) amax=a1;
  a1=RFamplitude[0].GetMean()-5.0*RFamplitude[0].GetRMS();
  if (a1<amin) amin=a1;
  a1=RFamplitude[1].GetMean()-5.0*RFamplitude[1].GetRMS();
  if (a1<amin) amin=a1;
  RFamplitude[0].GetXaxis()->SetRangeUser(amin, amax);
  RFamplitude[0].DrawCopy();
  RFamplitude[1].DrawCopy("same");

  PhaseCanvas->cd(4);
  gStyle->SetOptStat(111110);
  amin=RFratio->GetMean()-5.0*RFratio->GetRMS();
  amax=RFratio->GetMean()+5.0*RFratio->GetRMS();
  RFratio->GetXaxis()->SetRangeUser(amin, amax);
  RFratio->DrawCopy();

  PhaseCanvas->cd(5);
  gStyle->SetOptStat(111110);
  RFperiod[0].DrawCopy();
  RFperiod[1].DrawCopy("same");

  PhaseCanvas->cd(6);
  gStyle->SetOptStat(111110);
  RFresidual[0].DrawCopy();
  RFresidual[1].DrawCopy("same");

  PhaseCanvas->Update();
  PhaseCanvas->Show();
//=============================================================================


////
  sprintf(titleC,"Summary for run#%d", Run);
  SummaryCanvas -> SetEditable(true);
  SummaryCanvas->SetTitle(titleC);
  SummaryCanvas->Clear();
  SummaryCanvas->Divide(2,3);

  SummaryCanvas -> Update();
  SummaryCanvas -> Show();
  SummaryCanvas -> SetEditable(false);
  
  
/************************************************************************/
/*


/////////////////////////////////////////////////////////
// pulse shapes here!
/////////////////////////////////////////////////////////

  for (int iMod=0; iMod<6; iMod++) {
    printf("[RunReport] Showing shapes for %s\n", Setup[iMod].GetName()); 
    Setup[iMod].SavePulseShapeSplit();
    Setup[iMod].ShowPulseShapeSplit();
  }

*/

// save profiles for calibration purposes



/*******************************************************************************************/
/*  End of RunReport */
  printf("Run#%d finished successfully\n", myRun);
  int ReportMode=ProcStatus->WhatToDo("[RunReport]", -300);
  ProcStatus->SetReportMode(ReportMode); 

  TL->Delete();
  dummy->Delete();
  dummy2D->Delete();

  return;
}




Double_t FitP3(Double_t *x, Double_t *par) {
  return par[3]+x[0]*(par[2]+x[0]*(par[1]+x[0]*par[0]));
}

//====================================================================================
/*
bool SetChannelHiLoCuts(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus) {
// checks for Hi/Low crysis cuts. not needed in 2012 run;
  char nameC[100];
  sprintf(nameC, "HiLoCuts/CutsRun%d.dat", Run);
  indata.open(nameC, ios::in);

  if (Run != ProcStatus->GetRun()) {
    printf("[SetChannelHiLoCuts] wrong run number-> Run=%d, ProcStatus->GetRun()=%d\n",
           Run, ProcStatus->GetRun());
  }

  for (int i=0; i<nRC-4; i++)    flagRCtime[i]=false; // clear integral computation flags except for last four - module drift times
  for (int i=nRC-4; i<nRC; i++)  flagRCtime[i]=true;
  
  if (!indata.is_open()) {
    printf("[SetChannelHiLoCuts] No High/Low cuts file %s for run#%d\n", nameC, Run);
    indata.clear();
    return false;
  }

  if (!indata.good()) {
    printf("[SetChannelHiLoCuts] High/Low cuts file %s for run#%d is bad\n", nameC, Run);
    indata.close();
    return false;
  }

  int nErr=0;
  for (int iMod=0; iMod<4; iMod++) {
    for (int iChan=0; iChan<4; iChan++) {
      TChannel *Channel =Setup[iMod].GetChannel(iChan);
      int iMd, iCh, nSt; 
      double hiLo, sigma, trc, rate, eff;
      indata >> iMd >> iCh >> nSt >> hiLo >> sigma >> trc >> rate >> eff;
      if (iMd != iMod || iCh != iChan) {
        printf("[SetChannelHiLoCuts] XXXXX Error in file %s, Module->(%d, %d), Channel->(%d, %d)\n",
               nameC, iMod, iMd, iChan, iCh);
        cout << iMd  << " " << iCh   << " " << nSt << " " <<
                hiLo << " " << sigma << " " << trc << " " << 
                rate << " " << eff   << endl;
        nErr++;
      }
      if (nSt>1000 && hiLo!=0.0 && sigma!=0.0 && sigma<0.6 && eff>0.3) Channel->SetHighLow(hiLo);

      
      int irc=-1;
      for (int i=0; i<nRC; i++) {
        if (RCtime[i]==trc) {
          irc=i;
          break;
        }
      }
      Channel->SetHLcorrection(hiLo, sigma, irc, trc, rate, eff);
    
      if (irc<0 && trc!=0.0) {
        printf("[SetChannelHiLoCuts] RC=%f not found in RCtime array ->\n", trc);
        for (int i=0; i<nRC; i++) {
          printf(" %6.2e", RCtime[i]);
          if ((i+1)%5==0) printf("\n");
        }
        //exit(-1);
        nErr++;
      } else {
        flagRCtime[irc]=true; // set computation flag for this integral
      }
    }
  }

  indata.close();
  return (nErr==0);
}

bool SetChannelHiLoCuts2011(Int_t Run, TModule* Setup, ProcessingStatus* ProcStatus) {
  char nameC[100];
  sprintf(nameC, "HiLoCuts/Cuts2011.06_Run%d.dat", Run);
  if (Run != ProcStatus->GetRun()) {
    printf("[Main->SetChannelHiLoCuts] wrong run number-> Run=%d, ProcStatus->GetRun()=%d\n",
           Run, ProcStatus->GetRun());
  }
  indata.open(nameC, ios::in);

  for (int i=0; i<nRC-4; i++)    flagRCtime[i]=false; // clear integral computation flags except for last four - module drift times
  for (int i=nRC-4; i<nRC; i++)  flagRCtime[i]=true;
  
  if (!indata.is_open()) {
    printf("[Main->SetChannelHiLoCuts2011] No High/Low cuts for run#%d\n", Run);
    indata.clear();
    return false;
  }

  if (!indata.good()) {
    printf("[Main->SetChannelHiLoCuts2011] High/Low cuts file %s for run#%d is bad\n", nameC, Run);
    indata.close();
    return false;
  }

  int nErr=0;
  for (int iMod=0; iMod<4; iMod++) {
    for (int iChan=0; iChan<4; iChan++) {
      TChannel *Channel =Setup[iMod].GetChannel(iChan);
      int iMd, iCh, nSt; 
      double hiLo, sigma, trc, rate, eff;
      double bin1, bin2;
      int nbins;
      indata >> iMd >> iCh >> nSt >> hiLo >> sigma >> trc >> rate >> eff >> eff >> bin1 >> bin2 >> nbins;
      //cout << iMd  << " " << iCh   << " " << nSt << " " <<
      //          hiLo << " " << sigma << " " << trc << " " << 
      //          rate << " " << eff <<  " " << bin1 << " " << bin2 << " " << nbins << endl;
      if (iMd != iMod || iCh != iChan) {
        printf("[SetChannelHiLoCuts2011] XXXXX Error in file %s, Module->(%d, %d), Channel->(%d, %d)\n",
               nameC, iMod, iMd, iChan, iCh);
        nErr++;
      }
      if (nSt>100 && hiLo!=0.0 && sigma!=0.0 && sigma<1.0) Channel->SetHighLow(hiLo);

      
      int irc=-1;
      for (int i=0; i<nRC; i++) {
        if (RCtime[i]==trc) {
          irc=i;
          break;
        }
      }
      Channel->SetHLcorrection2011(hiLo, sigma, nSt, irc, trc, rate, eff, bin1, bin2, nbins);
      if (nbins>200) {
        printf("[SetChannelHiLoCuts2011] too many HiLi rate bins requested!\n");
        nErr++;
      }
      for (int ib=1; ib<nbins; ib++) {
        int ib1;
        double xb;
        int nstat; 
        double mean;
        double rmsUp;
        double rmsDw;
        indata >> ib1 >> xb >> nstat >> mean >> rmsUp >> rmsDw;
        Channel->SetHLcorrection2011bin(ib1, xb, nstat, mean, rmsUp, rmsDw);
        //printf( "     ix=%3d x=%5.3f n=%4d Median=%5.2f rms+=%.1f rms-=%.1f\n",
        //         ib1, xb, nstat, mean, rmsUp, rmsDw);
        if (ib != ib1) {
          nErr++;
          printf("[SetChannelHiLoCuts2011] bin index mismatch! %d %d\n", ib, ib1);
        }
      }
    
      if (irc<0 && trc!=0.0) {
        printf("[SetChannelHiLoCuts2011] RC=%f not found in RCtime array ->\n", trc);
        for (int i=0; i<nRC; i++) {
          printf(" %6.2e", RCtime[i]);
          if ((i+1)%5==0) printf("\n");
        }
        //exit(-1);
        nErr++;
      } else {
        flagRCtime[irc]=true; // set computation flag for this integral
      }
    }
  }

  indata.close();
  return (nErr==0);
}
*/



bool TuneStreamCherenkovTiming(TModule* Setup, CherenkovBuffer* CherBuff, ProcessingStatus* ProcStatus) {
//  check synchronization of MPI data with CherenkovGatedStream //
  int    iOffbest=0;
  double covBest=-1.;        
  double corBest=-1.; 
  double rmsBest=0.;
  bool Dbg=false;

  int Nbunches=ProcStatus->GetNbunches();
  if (Dbg && ProcStatus->GetReportMode()==1) {
    for (int bunch=0; bunch<Nbunches; bunch++) {
      double Tbunch=ProcStatus->GetBunchTime(bunch);
      double CherenkovADC=Setup[4].GetChannel(3)->GetResponseADC(bunch);
      double CherenkovStream1=CherBuff->GetSampleAdjustD(Tbunch);
      double CherenkovStream2=CherBuff->GetSampleAdjustD(int(Tbunch/CherBuff->GetDT()));
      printf(" bunch#%d T=%5.1f A=%5.1f GatedSample=%d Gated1=%5.1f Gated2=%5.1f\n", 
             bunch, Tbunch, CherenkovADC, int(Tbunch/CherBuff->GetDT()),
             CherenkovStream1, CherenkovStream2);
    }  
  }

  double stream0=15.; // simulate 15 events at zero for proper correlation leverage
  double stream1=0.;
  double stream2=0.;
  for (int is=0;is<39;is++) {
    double CherenkovStream=CherBuff->GetSampleAdjustD(is);
    stream0+=1.0;
    stream1+=CherenkovStream;
    stream2+=CherenkovStream*CherenkovStream;
  }
  stream1/=stream0;
  stream2/=stream0;
  stream2-=stream1*stream1;
  stream2=sqrt(stream2);

  for (int iOff=-4; iOff<5; iOff++) {
    double s0=Nbunches/2; // to simulate point (0,0) 
    double a1=0.; 
    double b1=0.; 
    double a2=0.; 
    double b2=0.; 
    double ab=0.; 
    for (int bunch=0; bunch<Nbunches; bunch++) {
      double Tbunch=ProcStatus->GetBunchTime(bunch);
      double CherenkovADC=Setup[4].GetChannel(3)->GetResponseADC(bunch);
      double CherenkovStream=CherBuff->GetSampleAdjustD(Tbunch+iOff*CherBuff->GetDT());
      s0+=1.0;
      a1+=CherenkovADC;
      b1+=CherenkovStream;
      a2+=CherenkovADC*CherenkovADC;
      b2+=CherenkovStream*CherenkovStream;
      ab+=CherenkovADC*CherenkovStream;
    }
    a1/=s0;
    b1/=s0;
    a2/=s0;
    b2/=s0;
    ab/=s0;
    a2-=a1*a1;
    b2-=b1*b1;
    ab-=a1*b1;

    if (Dbg && ProcStatus->GetReportMode()==1) {
      printf("  iOff=%2d corr=%5.2f cov=%.1f\n", iOff, ab/sqrt(a2*b2), ab);
    }
    if (ab>covBest) {
      iOffbest=iOff;
      covBest=ab;
      corBest=ab/sqrt(a2*b2);
      rmsBest=sqrt(b2);
    }
  } // offset looper

  if (iOffbest != 0 && rmsBest>stream2 && corBest>0.6 && rmsBest>6.0/corBest*CherBuff->GetSpillNoise() )  {
    printf("------- best Offset=%d cov=%.1f corr=%.2f S-factor=%3.1f bestN=%5.1f spillN=%5.1f streamN=%5.1f\n", 
              iOffbest, covBest, corBest, ProcStatus->GetEventHeader()->SignalFactor,
              rmsBest, CherBuff->GetSpillNoise(), stream2);
    printf("      changing  Gated Stream offset from %d to %d\n",  
              CherBuff->GetEventOffset(), CherBuff->GetEventOffset()-iOffbest);
    CherBuff->SetEventOffset(CherBuff->GetEventOffset()-iOffbest);
    printf("----------------------------------------------------\n");
    //ProcStatus->SetReportMode(1);
    return false;
  }
  
  return true;
}

//--------------------------------------------------------------------------------//
void MyExecuteEvent(Int_t event, Int_t px, Int_t py, TObject *sel)
{
  if(sel) {
    //
  }
  TCanvas *c = (TCanvas *) gTQSender;

  if ( (event == kButton1Double) ) {
    char cname[256];
    TPad *pad = c->Pick(px, py, 0);
    //printf("d-clicked on pad %s #%d\n", pad->GetName(), pad->GetNumber());
    if(pad == 0) {
      std::cout << "Can;t get pad at " << px << " " << py << std::endl;
      return;
    }
    if(pad == c) {
      std::cout << "You have selected parent object" << std::endl;
      return;
    }
    float xyratio =  float(pad->GetAbsWNDC()*pad->GetCanvas()->GetWindowWidth()) / float(pad->GetAbsHNDC()*pad->GetCanvas()->GetWindowHeight());
    //std::cout << " xyratio: " << xyratio << std::endl;
    TPad *pad_new = (TPad *) pad->Clone();
    pad_new->SetPad(0.0, 0.0, 1.0, 1.0);
    //pad_new -> SetEditable(true);

    sprintf(cname,"%s_%d",pad->GetTitle(), (int)time(0));
    int iY=600;
    int iX=(int)(600*xyratio);
    if (iX<600) {
      iY=iY*600/iX;
      iX=600;
    }
    TCanvas *c1 = new TCanvas(cname, cname, iX, iY);
    c1->cd();
    if (cname[0]=='S' && cname[1]=='p' && cname[2]=='e' && cname[3]=='c') {
      gStyle->SetOptStat(111110);
    } 
    pad_new->SetEditable(kTRUE);
    pad_new->Draw();    
    pad_new->Modified();
    pad_new->Update();
//    pad_new->SetEditable(kFALSE);
//     c1->RaiseWindow();
//     c1->Modified();
//     c1->Update();
    c1->SetWindowSize(iX+20, iY+40);
    c1->Update();
    gSystem->ProcessEvents();
  }
}

#include <classes/UnpackFunctions.C>
#include <classes/DisplayAllPulses2012.C>
#include <classes/DisplayAllCounters2012.C>
#include <classes/DisplayAllPulsesDetails.C>
#include <classes/DisplayForMaslennikov.C>
#include <classes/GlobalPhase2011.C>
#include <classes/GlobalPhase2012.C>
#include <classes/LoadTransferFunction.C>
#include <classes/LoadChannelStatus.C>
#include <classes/LoadFFTcorrection.C>
#include <classes/LoadDetectorMap.C>
#include <classes/LoadPMTdelays.C>
#include <classes/LoadMIP2data.C>
#include <classes/LoadMIP2dataGroup.C>
#include <classes/LoadRFcorrection.C>
#include <classes/SavePulseShapes.C>
#include <classes/SetModuleTiming.C>
#include <classes/PrepareTimeSlices.C>