#include <math.h>

#include <pthread.h>
#include <stdio.h>
#include <stdlib.h>
#include <assert.h>
#include <string.h> 
#include <unistd.h>

#include "RootHeader.h"

#include "TTree.h"
#include "TProfile.h"
#include "TH1F.h"
#include "TFile.h"



void calib_ped_pol(int crate, int slot, int range = 0)
{


  //  TFile *f = new TFile("tagh_4.root");  //   DAQ 3300
  //  TFile *f = new TFile("tagh_5.root");      //   DAQ 3200
  //  TFile *f = new TFile("tagh_6.root");      //   Check DAC Settings

  //  TFile *f = new TFile("ps_7.root");      //   Check DAC Settings
  //  TFile *f = new TFile("ps_9.root");      //   DAC 3200
  //  TFile *f = new TFile("ps_10.root");     //   DAC 3300
  //  TFile *f = new TFile("ps_11.root");     //   DAC Check 

  //  PSC
  //  TFile *f = new TFile("psc_5.root");         //   DAC  3200 
  //  TFile *f = new TFile("psc_6.root");         //   DAC 3300
  //  TFile *f = new TFile("psc_7.root");             //   DAC 3300

  //  TFile *f = new TFile("ps_calib_58.root");          //   DAC  3200 
  //  TFile *f = new TFile("ps_calib_59.root");          //   DAC  3300
  //   TFile *f = new TFile("ps_calib_60.root");         //   DAC  3300

  //   TFile *f = new TFile("fcal_cosmic_126.root");     //   DAC  3200
  
  //  TFile *f = new TFile("tagm_calib_169.root");       //   DAC  3200
  //  TFile *f = new TFile("tagm_calib_170.root");       //   DAC  3300

  //  TFile *f = new TFile("tagm_calib_171.root");       //   DAC  3300

  //  TFile *f = new TFile("tagm_calib_203.root");       //   DAC  3200
  //  TFile *f = new TFile("tagm_calib_204.root");       //   DAC  3200

  // PS check
  //  TFile *f = new TFile("tagm_calib_302.root"); 

  // New PS calibration for 2015
  //  TFile *f = new TFile("ps_play_back_505.root"); 
  //  TFile *f = new TFile("ps_play_back_506.root");
  
  //  TFile *f = new TFile("ps_play_back_507.root");
  
  //  TFile *f = new TFile("calib_15.root");
  
  // TFile *f = new TFile("hd_rawdata_004224_000.root");

  // FALL RUN PS2, SLOTS 3,4,5
  //  TFile *f = new TFile("calib_4427.root");
  //  TFile *f = new TFile("calib_4428.root"); bad

  //  TFile *f = new TFile("calib_4431.root");
  
  //  TFile *f = new TFile("calib_4434.root");

  // Check PS1
  //  TFile *f = new TFile("hd_rawdata_004436_000.root");

  //  TFile *f = new TFile("hd_rawdata_004436_000.root");

  // TAGM

  //  TFile *f = new TFile("hd_rawdata_004459_000.root");

  //  TFile *f = new TFile("hd_rawdata_004727_000.root");

  // DAC 3200
  // TFile *f = new TFile("hd_rawdata_010059_000.root");

  // DAC 3300
  TFile *f = new TFile("hd_rawdata_010060_000.root");

  gStyle->SetOptStat(11111111);
  

  TH1F *hped = new TH1F("PS ped mean", "PS ped mean",1100,3000.,4100.);
  TH1F *hped_sig = new TH1F("PS ped sigma", "PS ped sigma",200,0.,4.);


  int write_file = 1;

  Int_t x_min = 3500;
  Int_t x_max = 4100;

  //  Int_t x_min = 90;
  //  Int_t x_max = 120;

  Float_t  ped[21][16];
  Float_t  ped_w[21][16];

  //  Int_t slot_min  = 3;
  //  Int_t slot_max  = 7;

  Int_t slot_min  = 13;
  Int_t slot_max  = 16;

  FILE *out_file;    
  
  if(write_file == 1){

    out_file = fopen("tpol_3300.txt","w");

    //    out_file = fopen("tagm_3300.txt","w");

    //    out_file = fopen("3300.txt","w");
    //    out_file = fopen("3200.txt","w");
    
    //    out_file = fopen("3200_cr0.txt","w");
    //    out_file = fopen("3300_cr0.txt","w");
  }
  
  TCanvas *c1[21];


  if(slot != 0){  
    slot_min = slot;
    slot_max = slot+1;
  } 


  cout << " SLOT_MIN = " << slot_min << endl;
  cout << " SLOT_MAX = " << slot_max << endl;
 

  memset(ped,0,sizeof(ped));
  memset(ped_w,0,sizeof(ped_w));


  for(Int_t sl = slot_min; sl <= slot_max; sl++){

    char name[30];
    sprintf(name," Crate %d, Slot %d", crate, sl);


    c1[sl] = new TCanvas(name,name,800,1000);
    c1[sl]->Divide(4,4);

    for(Int_t ch = 0; ch < 16; ch++){

      c1[sl]->cd(ch+1);

      char title[30];
      sprintf(title,"ped_cr%d_bd%d_ch%d",crate,sl,ch);
      
      TH1F * h1 = (TH1F*)f->Get(title); 

      
      cout << " Read histogram " << title <<  "  "  <<  h1  <<  endl;
      

      if(h1 == NULL){
	cout <<  " Channel don't exist " <<  title << endl;
	continue;
      }

      if(range == 1){
	h1->GetXaxis()->SetRangeUser(x_min,x_max);
      }

      h1->Fit("gaus");

      Double_t mean  = 0;
      Double_t sigma = 0;

      TF1 *fit = h1->GetFunction("gaus");
      mean   =  fit->GetParameter(1);
      sigma  =  fit->GetParameter(2);
      
      cout << sl << endl;

      ped[sl][ch]   = mean;
      ped_w[sl][ch] = sigma;

      if(sigma > 3){
	cout << " CRATE =  "  <<  crate << "   SLOT = "  <<  sl  <<  " Channel = " << ch  << endl;
	cout << " Mean =  "   <<  mean <<  " Sigma = "  <<  sigma << endl;
	getchar();
      }      

      h1->Draw();
      
      hped->Fill(mean);
      hped_sig->Fill(sigma);

    }   // Loop over channels

    if(write_file == 1){
      
      fprintf(out_file,"CRATE %d    SLOT  %d \n",crate,sl);
      
      fprintf(out_file,"PED    %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f\n",
	      ped[sl][0], ped[sl][1], ped[sl][2],ped[sl][3],ped[sl][4],ped[sl][5],ped[sl][6],ped[sl][7],
	      ped[sl][8], ped[sl][9], ped[sl][10],ped[sl][11],ped[sl][12],ped[sl][13],ped[sl][14],ped[sl][15]);
    
      fprintf(out_file,"SIG    %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f  %7.2f\n",
	      ped_w[sl][0], ped_w[sl][1], ped_w[sl][2],ped_w[sl][3],ped_w[sl][4],ped_w[sl][5],ped_w[sl][6],ped_w[sl][7],
	      ped_w[sl][8], ped_w[sl][9], ped_w[sl][10],ped_w[sl][11],ped_w[sl][12],ped_w[sl][13],ped_w[sl][14],ped_w[sl][15]);
      
    }
    
  }     // Loop over boards
  
  TCanvas *pp = new TCanvas("pp","pp",800,1000);

  pp->Divide(1,2);
  
  pp->cd(1);
  hped->Draw();

  pp->cd(2);
  hped_sig->Draw();



  if(write_file == 1) fclose(out_file);

}