#include "TCanvas.h"
#include "TGraph.h"
#include "TFile.h"
#include "TH1.h"
#include "TStyle.h"
#include "TROOT.h"
#include "TStyle.h"
#include "TAxis.h"
#include "TMath.h"
#include "TAxis.h"
#include "TGraph.h"
#include "TPad.h"
#include "TLatex.h"
#include "TApplication.h"
#include "func.h"

int fillraw(char* data_dir, TH1F** hh)
{
  // 
  int ch_list[32] = {0,1,99,3,4,5,6,99,8,9,99,11,12,13,99,15,99,99,99,99,99,2,99,99,99,99,7,99,10,99,14,99};
  char fname[256];
  sprintf(fname,"%s/Raw_Data.txt",data_dir);
  FILE* data = fopen(fname,"r");
  if(!data)
    {
      cout<<"data file not found at "<<fname<<endl;
      return 0;
    }
  int nl = 0;
  char line[256];
  int tmpc, tmpadc;
  while(fgets(line,256,data))
    {
      if(sscanf(line, "%d %d", &tmpc,&tmpadc)==2)
	{
	  if(ch_list[tmpc] < 0 || ch_list[tmpc] >= 16) continue;
	  else hh[ch_list[tmpc]]->Fill(tmpadc);
	  nl++;
	  if(nl%100000 == 0) cout<<nl<<" lines processed"<<endl;
	}
    }
  //  hh[1]->Draw();
  return nl;
}

void replay()
{
  gStyle->SetOptFit(1);
  //  char data_dir[20] = "20110216_1137";
  char data_dir[512];
  int input_offset = 2;
  sprintf(data_dir,"./data/QDC/%s",gApplication->Argv(input_offset));
  int xtalk_flag = atoi(gApplication->Argv(1+input_offset));
  int save_flag = atoi(gApplication->Argv(2+input_offset));

  cout<<"Replaying data in "<<data_dir<<endl;
  cout<<"Cross-talk fit: "<<xtalk_flag<<endl;
  cout<<"Save result: "<<save_flag<<endl;

  TH1F* hh[32];
  char hname[256];
  char htitle[256];
  for(int i = 0; i<32; i++)
    {
      sprintf(hname,"histo_ch%02d",i);
      sprintf(htitle,"ADC Channel #%02d; ADC; Counts",i);
      hh[i] = new TH1F(hname,htitle,4096,0,4096);
    }
  int nn = fillraw(data_dir,hh);
  cout<<nn<<" lines were read."<<endl;
  if(!nn) return;

  int nch = 16;

  double par[7], par_err[7], xmin, xmax;
  double mu[16],gain[16],sigma[16],ped[16],sigma_p[16],xtalk[16],ch[16];
  double mu_e[16],gain_e[16],sig_e[16],ped_e[16],sig_p_e[16],xt_e[16];

  TCanvas* canvas[1];
  for(int i = 0; i<1; i++)
    {
      sprintf(hname,"c%d",i+1);
      sprintf(htitle,"ADC Channel #%02d - %02d",i*16,i*16+15);
      canvas[i] = new TCanvas(hname,htitle,1500,960);
      canvas[i]->Divide(4,4);
    }
  for(int i = 0; i<nch; i++)
    {
      int ci = (int)floor(i/16);
      canvas[ci]->cd(i%16+1);
      cout<<"Fitting "<<hh[i]->GetTitle()<<endl;
      fit(hh[i],par,par_err,xtalk_flag);
      xmin = par[5]-par[3];
      xmax = par[5]+10.*par[3];
      hh[i]->SetAxisRange(xmin,xmax,"X");
      canvas[ci]->Update();
      mu[i] = par[0];
      mu_e[i] = par_err[0];
      sigma[i] = fabs(par[1]);
      sig_e[i] = par_err[1];
      gain[i] = par[3];
      gain_e[i] = par_err[3];
      sigma_p[i] = fabs(par[4]);
      sig_p_e[i] = par_err[4];
      ped[i] = par[5];
      ped_e[i] = par_err[5];
      xtalk[i] = par[6];
      xt_e[i] = par_err[6];
      ch[i] = i;
    }

  // Plot uniformity
  TGraph* gr_mu = new TGraph(nch,ch,mu);
  gr_mu->SetTitle("Number of Photon Electrons; Channel; NPE");
  gr_mu->SetFillColor(kCyan+2);
  gr_mu->GetXaxis()->SetLimits(-1,nch);
  gr_mu->SetMinimum(0);
  gr_mu->SetMaximum(gr_mu->GetYaxis()->GetXmax()*1.2);
  //  gr_mu->SetFillStyle(3001);
  TGraph* gr_sigma = new TGraph(nch,ch,sigma);
  gr_sigma->SetTitle("Sigma of Single Photon Signal; Channel; Sigma");
  gr_sigma->SetFillColor(kMagenta-2);
  gr_sigma->GetXaxis()->SetLimits(-1,nch);
  gr_sigma->SetMinimum(0);
  gr_sigma->SetMaximum(gr_sigma->GetYaxis()->GetXmax()*1.2);
  //  gr_sigma->SetFillStyle(3001);
  TGraph* gr_gain = new TGraph(nch,ch,gain);
  gr_gain->SetTitle("Gain; Channel; Gain");
  gr_gain->SetFillColor(kRed-3);
  gr_gain->GetXaxis()->SetLimits(-1,nch);
  gr_gain->SetMinimum(0);
  gr_gain->SetMaximum(gr_gain->GetYaxis()->GetXmax()*1.2);
  //  gr_gain->SetFillStyle(3001);
  TGraph* gr_sigma_p = new TGraph(nch,ch,sigma_p);
  gr_sigma_p->SetTitle("Sigma of Pedestal; Channel; Sigma_{ped}");
  gr_sigma_p->SetFillColor(kTeal+3);
  gr_sigma_p->GetXaxis()->SetLimits(-1,nch);
  gr_sigma_p->SetMinimum(0);
  gr_sigma_p->SetMaximum(gr_sigma_p->GetYaxis()->GetXmax()*1.2);
  //  gr_sigma_p->SetFillStyle(3001);
  TGraph* gr_ped = new TGraph(nch,ch,ped);
  gr_ped->SetTitle("Pedestal; Channel; Pedestal");
  gr_ped->SetFillColor(kAzure+3);
  gr_ped->GetXaxis()->SetLimits(-1,nch);
  gr_ped->SetMinimum(0);
  gr_ped->SetMaximum(gr_ped->GetYaxis()->GetXmax()*1.2);
  //  gr_ped->SetFillStyle(3001);
  TGraph* gr_xtalk = new TGraph(nch,ch,xtalk);
  gr_xtalk->SetTitle("Cross-Talk; Channel; Cross-Talk");
  gr_xtalk->SetFillColor(kOrange+7);
  gr_xtalk->GetXaxis()->SetLimits(-1,nch);
  gr_xtalk->SetMinimum(0);
  gr_xtalk->SetMaximum(gr_xtalk->GetYaxis()->GetXmax()*1.2);
  //  gr_xtalk->SetFillStyle(3001);

  //  gr_mu->SetAxisRange(-1,4);

  TLatex ltext;
  char text[256];

  TCanvas* c2 = new TCanvas("c2","Uniformity",800,960);
  c2->Divide(2,3);
  c2->cd(1);
  gr_gain->Draw("AB");
  sprintf(text,"%.2f #pm %.1f%%",gr_gain->GetMean(2),gr_gain->GetRMS(2)*100./gr_gain->GetMean(2));
  ltext.DrawLatex(nch*0.3,gr_gain->GetMaximum()*0.88,text);

  c2->cd(2);
  gr_mu->Draw("AB");
  sprintf(text,"%.2f #pm %.1f%%",gr_mu->GetMean(2),gr_mu->GetRMS(2)*100./gr_mu->GetMean(2));
  ltext.DrawLatex(nch*0.3,gr_mu->GetMaximum()*0.88,text);

  c2->cd(3);
  gr_ped->Draw("AB");
  sprintf(text,"%.2f #pm %.1f%%",gr_ped->GetMean(2),gr_ped->GetRMS(2)*100./gr_ped->GetMean(2));
  ltext.DrawLatex(nch*0.3,gr_ped->GetMaximum()*0.88,text);

  c2->cd(4);
  gr_xtalk->Draw("AB");
  sprintf(text,"%.2f #pm %.1f%%",gr_xtalk->GetMean(2),gr_xtalk->GetRMS(2)*100./gr_xtalk->GetMean(2));
  ltext.DrawLatex(nch*0.3,gr_xtalk->GetMaximum()*0.88,text);

  c2->cd(5);
  gr_sigma->Draw("AB");
  sprintf(text,"%.2f #pm %.1f%%",gr_sigma->GetMean(2),gr_sigma->GetRMS(2)*100./gr_sigma->GetMean(2));
  ltext.DrawLatex(nch*0.3,gr_sigma->GetMaximum()*0.88,text);

  c2->cd(6);
  gr_sigma_p->Draw("AB");
  sprintf(text,"%.2f #pm %.1f%%",gr_sigma_p->GetMean(2),gr_sigma_p->GetRMS(2)*100./gr_sigma_p->GetMean(2));
  ltext.DrawLatex(nch*0.3,gr_sigma_p->GetMaximum()*0.88,text);

  // save results
  if(save_flag)
    {
      char fname[256];
      sprintf(fname,"%s/All_Histogram.png",data_dir);
      canvas[0]->Print(fname);
      sprintf(fname,"%s/Distribution.png",data_dir);
      c2->Print(fname);
      sprintf(fname,"%s/Results.txt",data_dir);
      FILE* data = fopen(fname,"w");
      fprintf(data,"channel\t npe\t npe_e\t gain\t gain_e\t pedestal\t ped_e\t xtalk\t xtalk_e\t sigma\t sig_e\t sigma_p\t sig_p_e\n");
      for(int i = 0; i<nch; i++)
	{
	  fprintf(data, "%d\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\t %.6f\n",i,mu[i],mu_e[i],gain[i],gain_e[i],ped[i],ped_e[i],xtalk[i],xt_e[i],sigma[i],sig_e[i],sigma_p[i],sig_p_e[i]);
	}
    }
}