#include <cometEventLoop.hxx>
#include "ICOMETRawEvent.hxx"
#include "IMidasFile.hxx"
#include <TApplication.h>
#include <memory>
#include <IEventFolder.hxx>
#include <IHandle.hxx> 
#include "IMidasBank.hxx"
#include <IRealDatum.hxx>
#include <IEventFolder.hxx>
#include <IOADatabase.hxx>

#include <TFile.h>
#include "ITPCPedestal.hxx"
#include <iostream>
#include <fstream> // needed for filehandling

using namespace COMET; 

#define MAXNDCCS 18 
#define MAXNFEMS  4
#define NFECS     6
#define NASICS    4
#define NCHANNELS 79



ITPCPedestal::ITPCPedestal( std::string fileName,std::string textfileName, std::string badtextfileName )
  : fFileName( fileName ),
    fTextFileName( textfileName ),
    fBadChanFileName( badtextfileName ),
    fFile(0)
 // setting member variables in the constructor 
{
  COMETVerbose("ITPCPedestal::ITPCPedestal()");
}

ITPCPedestal::~ITPCPedestal()
{
  
}

// Initialization of task (at start of first file)
int ITPCPedestal::Initialize( void ) {

  //  tpcpeddecoder = new ITpcPedDecoder(3.,10.,500.); // Fit to +/- 3sigma, range of 10-500.
  tpcpeddecoder = new ITpcPedDecoder(2.,10.,500.); // Fit to +/- 3sigma, range of 10-500.

// output statement if verbosity is enabled

  COMETVerbose("ITPCdEdxCalib::Initialize( void )"); //use -v -v option to enable

  //  fFile->cd();

  bzero(chanisbad,NDCC*NFEM*NFEC*NASIC*NCHAN*sizeof(int));
  bzero(fitbad,NDCC*NFEM*sizeof(int));
  bzero(chanmap,NDCC*NFEM*sizeof(int));

  gPedMean = new TGraph();  // Graph of mean of pedestal vs. DCC/FEM/FEC/ASIC/Channel
  gPedMean->SetName("gPedMean"); // if Name is not set, will not appear in output file
  gPedMean->SetMarkerStyle(24);

  gPedLocalMean = new TGraph();  // Graph of mean of pedestal vs. DCC/FEM/FEC/ASIC/Channel
  gPedLocalMean->SetName("gPedLocalMean"); // if Name is not set, will not appear in output file
  gPedLocalMean->SetMarkerStyle(24);


  gPedLocalThres = new TGraph();  // Graph of mean of pedestal vs. DCC/FEM/FEC/ASIC/Channel
  gPedLocalThres->SetName("gPedLocalThres"); // if Name is not set, will not appear in output file
  gPedLocalThres->SetMarkerStyle(24);

  gPedSigma = new TGraph();  // Graph of RMS of pedestal vs. DCC/FEM/FEC/ASIC/Channel
  gPedSigma->SetName("gPedSigma"); // if Name is not set, will not appear in output file
  gPedSigma->SetMarkerStyle(23);

  gPedEntries = new TGraph();  // Graph of number of entries vs. DCC/FEM/FEC/ASIC/Channel
  gPedEntries->SetName("gPedEntries"); // if Name is not set, will not appear in output file
  gPedEntries->SetMarkerStyle(20);

  gPedChi = new TGraph();  // Graph of goodness of fit vs. DCC/FEM/FEC/ASIC/Channel
  gPedChi->SetName("gPedChi"); // if Name is not set, will not appear in output file
  gPedChi->SetMarkerStyle(20);


  return 0;
  
}
 

int ITPCPedestal::Process(COMET::ICOMETEvent& event){

    tpcpeddecoder->Process(event);

    return 0;
}



int ITPCPedestal::Finalize( void ){

  COMETVerbose("ITPCdEdxCalib::Finalize( void )");

  COMETVerbose( "Check for any bad channels in the " << fBadChanFileName.c_str() << "file" ); 
  // Read in bad channels from file
  //  char name[256];
  // sprintf(name,"badchannels.dat");
  std::fstream BadFile(fBadChanFileName.c_str(),std::ios::in);
  if (BadFile.fail()){
    COMETError("Failed to open bad channel text file (default badchannels.dat). Code will assuming all channels are good, can cause problems for particular channels");
  }

  do{
    int DCC,Fem,Fec,Asic,Ch;
    BadFile>>DCC>>Fem>>Fec>>Asic>>Ch;
    if ( !BadFile.fail() ){
      if( (Fem >= 0 && Fem < MAXNFEMS) && (Fec >= 0 && Fec < NFECS) && 
    (Asic >= 0 && Asic < NASICS) && (Ch >= 0 && Ch < NCHANNELS) ) {
  chanisbad[DCC][Fem][Fec][Asic][Ch] = 1;
      }
    }
  } while ( !BadFile.fail() );
  BadFile.close();

  // Following calibration text file code
  // FILE *pedFile;
  //  const char pedestalfile = fTextFileName;
  //pedFile = fopen(fTextFileName,"w");

  COMETVerbose("Open Pedestal file");
  std::fstream pedFile(fTextFileName.c_str(), std::ios::out);
  //  std::fstream pedFile("pedestals.dat", std::ios::out);
  //  if (pedFile.fail()) {
  //  std::cout << "Could not open " << fTextFileName << std::endl;
  //  }
  
  int ctr =0;
  double nsigmas = 4.5;
  int mean = 250;
  
  int entries = 0; // counter for histograms under 30 events
  int entries0 = 0; // counter for histogram entries == 0 events
  int histo = 0; // counter for histograms existance
  int badchannel = 0; // counter for number of bad channels
  int widerms = 0; // counter for wide RMS
  int badfit1 = 0; // counter for bad fit according to root
  int badfit2 = 0; // counter for bad fit via testing mean vs. fitted mean
  int totalchannels = 0; // counter to make sure total number of channels is correct
  
  
  COMETVerbose("Begin reading histograms");
  
  for( int idcc = 0; idcc < NDCC; idcc++ ) 
    for( int ifem = 0; ifem < NFEM; ifem++ ) 
      for( int ifec = 0; ifec < NFEC; ifec++ ) 
  for( int iasic = 0; iasic < NASIC; iasic++ ) 
    for( int ichan = 0; ichan < NCHAN; ichan++ ) {
      
      TH1F *h =tpcpeddecoder->getHistogram(idcc,ifem,ifec,iasic,ichan);
      // entries, mean, sigma
      gPedMean->SetPoint(ctr,ctr,-1.);      
      gPedSigma->SetPoint(ctr,ctr,-1.);      
      gPedChi->SetPoint(ctr,ctr,-1.);      
      
      if(!h){
        COMETVerbose("Histogram not found for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan);
        histo++;
      }
        
      if( h ){
        

        Int_t numentries = static_cast<Int_t>(h->GetEntries());
        Double_t hmean = h->GetMean();
        Double_t hrms = h->GetRMS();
        
        gPedEntries->SetPoint(ctr,ctr,numentries);
    
        //	      COMETVerbose("Have histogram for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan << " with entries=" << numentries);
        if(ichan>2&&ichan!=15&&ichan!=28&&ichan!=53&&ichan!=66){
    if(numentries<30){
      COMETVerbose("Less than 30 entries per channel, run more events; for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan << " with entries=" << numentries);
      entries++;
      if(numentries==0){
        COMETVerbose("Channel has 0 events unexpectedly, possible problem with DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan << " with entries=" << numentries);
        entries0++;
      }
    }
        }		  
        
        if(numentries>10){
        Double_t fmean = tpcpeddecoder->getMean(idcc,ifem,ifec,iasic,ichan);
        gPedMean->SetPoint(ctr,ctr,fmean);      
        COMETVerbose("Found Mean ");		
        
        Double_t frms = tpcpeddecoder->getSigma(idcc,ifem,ifec,iasic,ichan);
        gPedSigma->SetPoint(ctr,ctr,frms);
    
        Double_t prob = tpcpeddecoder->getProb(idcc,ifem,ifec,iasic,ichan);
        gPedChi->SetPoint(ctr,ctr,prob);      
        
        COMETVerbose("Mean: " << hmean << " fitted mean: " << fmean << " RMS: " << hrms << " fitted RMS " << frms << " entries: " << numentries);
        // Test if Gaussian fit had any problems
        // Removed as a check, is apparently meaningless. Means are consitent with fitted means, same for RMS.
        /*
    if(prob<0.01||prob>99.99){
    
    COMETVerbose("Low probability of TFit of " << prob << " for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan);
    COMETVerbose("Mean: " << hmean << " fitted mean: " << fmean << " RMS: " << hrms << " fitted RMS " << frms << " entries: " << numentries);
    badfit2++;		    
        }
        */	
    
        // Check if mean is far away from fitted mean
        
        if((fmean-hmean)>hrms){
    COMETVerbose("Fitted mean is far from histogram mean for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan);
    COMETVerbose("Mean: " << hmean << " fitted mean: " << fmean << " RMS: " << hrms << " fitted RMS " << frms << " entries: " << numentries);
    badfit1++;		    
        }


        // Check if mean 0, indicates a bad fit or bad histogram not otherwise caught
        if(ichan>2&&ichan!=15&&ichan!=28&&ichan!=53&&ichan!=66){
    if(fmean<0.1){
      COMETVerbose("Fitted mean is very small for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan);
      COMETVerbose("Mean: " << hmean << " fitted mean: " << fmean << " RMS: " << hrms << " fitted RMS " << frms << " entries: " << numentries);
      fitbad[idcc][ifem]++;
      badfit2++;		    
    }
        }
        int localmean  = (int)(mean-fmean); 
        
        int localthresh  = (int)(mean + nsigmas*frms);
    
        if(frms>8.||frms<1.){
    // Assign a large threshold to any wide RMS channels, so called "bad channels"
    localthresh = 350;
    widerms++;
    // Write to bad channel file
    if(!chanisbad[idcc][ifem][ifec][iasic][ichan]){
      
      std::fstream BadFileWrite(fBadChanFileName.c_str(),std::ios::out | std::ios::app);
      BadFileWrite << idcc << " " << ifem << " " << ifec << " " << iasic << " " << ichan << std::endl;
      BadFileWrite.close();
    }
        }
    
        if ( chanisbad[idcc][ifem][ifec][iasic][ichan] ){
    COMETVerbose("Known bad channel for DCC will be given a high threshold (350) for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan);
    localthresh = 350;
    badchannel++;
        }
        
        gPedLocalMean->SetPoint(ctr,ctr,localmean);      
        gPedLocalThres->SetPoint(ctr,ctr,localthresh);      
        
        //		COMETVerbose("Write to pedestal file for DCC " << idcc << " FEM " << ifem << " FEC " << ifec << " ASIC " << iasic << " CHAN" << ichan);
        pedFile << idcc << " " << ifem << " " << ifec << " " << iasic << " " << ichan << " " << localmean << " " << localthresh << std::endl;	
        totalchannels++;
        chanmap[idcc][ifem]++;
        } // number of entries >10 
      } // if histo gram exists
      ctr++;
    }
  

      
  std::cout << "== Summary of pedestal file generation ==" << std::endl;
  std::cout << "Investigate any of following errors before uploading pedestal:" << std::endl;
  std::cout << "# of missing histograms (i.e. DCC not present): " << histo << std::endl;
  std::cout << "# of histograms with events <30 (can create poor fits): " << entries << std::endl;
  std::cout << "# of histograms with events == 0 (issue with channel or data): " << entries0 << std::endl;
  std::cout << "# of fit failures (mean inconsistent with fitted mean): " << badfit1 << std::endl;
  std::cout << "# fits with mean~0 (possible pedestal bank corruption or channel issue): " << badfit2 << std::endl;
  std::cout << "# of wide RMS channels: " << widerms << std::endl;
  std::cout << "# of bad channels from text file: " << badchannel << std::endl;
  if(totalchannels!=124416){
    std::cout << "Total channels written to file incorrect: " << totalchannels << " expect 124416" << std::endl;
    std::cout << "Expect 72 * 24 = 1728 per each FEM " << std::endl;
    std::cout << "         FEM0   FEM1  FEM2   FEM3" << std::endl;
    for( int idcc = 0; idcc < NDCC; idcc++ ) 
      std::cout << std::setw(6) <<  "DCC" << idcc << " " << std::setw(6) <<  chanmap[idcc][0] << " " << std::setw(6) <<   chanmap[idcc][1] << " " << std::setw(6) <<   chanmap[idcc][2] << " " << std::setw(6) <<   chanmap[idcc][3] <<  std::endl;
    
  }

  if(badfit2>0){
    std::cout << "         FEM0   FEM1  FEM2   FEM3" << std::endl;
    for( int idcc = 0; idcc < NDCC; idcc++ ) 
      //      std::cout << std::setw(4) <<  "DCC" << idcc << "   " <<  fitbad[idcc][0] << "   " <<   fitbad[idcc][1] << "   " <<   fitbad[idcc][2] << "   " <<   fitbad[idcc][3] << "   " <<  std::endl;
      std::cout << std::setw(6) <<  "DCC" << idcc << " " << std::setw(6) <<  fitbad[idcc][0] << " " << std::setw(6) <<   fitbad[idcc][1] << " " << std::setw(6) <<   fitbad[idcc][2] << " " << std::setw(6) <<   fitbad[idcc][3] <<  std::endl;
  }

    

  TList *list = gDirectory->GetList();  
  fFile = new TFile( fFileName.c_str(), "RECREATE" );
  //    TFile f("output_ped.root","RECREATE");
  list->Write("allhists",TObject::kSingleKey);
  
  //f.Close();
  // fFile->cd();
  
    
  fFile->WriteTObject(gPedMean);
  fFile->WriteTObject(gPedEntries);
  fFile->WriteTObject(gPedSigma);
  fFile->WriteTObject(gPedChi);
  fFile->WriteTObject(gPedLocalMean);
  fFile->WriteTObject(gPedLocalThres);

  
  fFile->Close(); // causes a seg fault
    
  
    // fclose(pedFile);
  pedFile.close();
  return 0;
  

};
  

void ITPCPedestal::SetFileName( const std::string& fileName )
{
  COMETVerbose("ITPCPedestal::SetFileName(" << fileName << ")" );
  fFileName = fileName;
}

void ITPCPedestal::SetTextFileName( const std::string& textfileName )
{
  COMETVerbose("ITPCPedestal::SetTextFileName(" << textfileName << ")" );
  fTextFileName = textfileName;

}

void ITPCPedestal::SetBadChanFileName( const std::string& badtextfileName )
{
  COMETVerbose("ITPCPedestal::SetBadChanFileName(" << badtextfileName << ")" );
  fBadChanFileName = badtextfileName;

}

bool ITPCPedestal::SetOption( std::string option, std::string value )
{
  if ( option == "tpc-ped-file" ) {
    SetTextFileName( value );
    return true;
  }

  if ( option == "bad-chan-file" ) {
    SetBadChanFileName( value );
    return true;
  }

  if ( option == "tpc-dn-root-file" ) {
    SetFileName( value );
    return true;
  }
  return false;
}