#include <inttypes.h>
#include <iostream>
#include <fstream>
#include <utility>
#include <vector>
#include <map>
#include <stdio.h>
#include <stdlib.h>
#include <string>
#include <cstring>
#include <sstream>

#include "TLatex.h"
#include "TPaveStats.h"
#include "TGraphPainter.h"
#include "TString.h"
#include "TCollection.h"
#include "TCanvas.h"
#include "TFile.h"
#include "TH2D.h"
#include "TH1F.h"
#include "TF1.h"
#include "TGraph.h"
#include "TGraphErrors.h"
#include "TMinuit.h"
#include "TKey.h"
#include "TDatime.h"
#include "TAxis.h"
#include "TLine.h"
#include "TTree.h"
#include "TBranch.h"
#include "TStyle.h"
#include "TProfile.h"

Double_t sqrt_func (Double_t *x, Double_t *par) {

   Double_t sigma1=par[0];
   Double_t x1=x[0];

   Double_t func;

   func = sigma1/sqrt(x1);
   //  cout << "x1=" << x1 << " sigma1=" << sigma1 << " func=" << func << endl;

   return func;
}

void pedestal_studiesA_sigma4(Int_t crate_select=37, Int_t slot_select=3, Int_t channel_select=0, Int_t sample_pick=0)
{

  // pedestal_studies script modified to compute max and min rms from data 

gStyle->SetPalette(1,0);
gStyle->SetOptStat(kTRUE);;
gStyle->SetOptStat(111111);
// gStyle->SetOptStat(kFALSE);
gStyle->SetOptFit(kTRUE);
// gStyle->SetOptFit(kFALSE);
// gStyle->SetOptFit(1111);
gStyle->SetPadRightMargin(0.15);
gStyle->SetPadLeftMargin(0.15);
gStyle->SetPadBottomMargin(0.15); 

gStyle->SetLabelSize(0.05,"XYZ");
gStyle->SetTitleSize(0.05,"XYZ");
gStyle->SetTitleOffset(1.5,"Y");
gStyle->SetNdivisions(505,"XYZ");

// gStyle->SetFillColor(0);
//

  // Read input file generated by DAQ and converted to ROOT using DAQTree plugin
  
  // TFile* f = new TFile("bcal_s4_3281.0.root");
  TFile* f = new TFile("hd_rawdata_003620_000.root");
  //TTree* tr = (TTree*)f->Get("channeltree");
  TTree* tr = (TTree*)f->Get("Df250WindowRawData");
  
  uint32_t channelnum;         /// Arbitrary global channel number
  uint32_t eventnum;           /// Event number
  uint32_t rocid;              /// Crate number
  uint32_t slot;               /// Slot number in crate
  uint32_t channel;            /// Channel number in slot
  uint32_t itrigger; ///trigger number for cases when this hit was read in a multi-event block (DDAQAddress)
  uint32_t nsamples;           /// Number of samples in the event
  std::vector<uint32_t>* waveform = 0;   /// STL vector of samples in the waveform for the event
  uint32_t w_integral;         /// Sum of all samples in the waveform
  uint32_t w_min;              /// Minimum sample in the waveform
  uint32_t w_max;              /// Maximum sample in the waveform
  uint32_t w_samp1;            /// First sample in the waveform
  
  tr->SetBranchAddress("channelnum",&channelnum);
  tr->SetBranchAddress("eventnum",&eventnum);
  tr->SetBranchAddress("rocid",&rocid);
  tr->SetBranchAddress("slot",&slot);
  tr->SetBranchAddress("channel",&channel);
  tr->SetBranchAddress("itrigger",&itrigger);
  tr->SetBranchAddress("nsamples",&nsamples);
  tr->SetBranchAddress("waveform",&waveform);
  tr->SetBranchAddress("w_integral",&w_integral);
  tr->SetBranchAddress("w_min",&w_min);
  tr->SetBranchAddress("w_max",&w_max);
  tr->SetBranchAddress("w_samp1",&w_samp1);
  
// Combos is 1, 2, 5, 10
//  int combo[] = {1,2,5,10,50};
//  const int _crate = 6, _slot = 20, _ch = 16, _co = sizeof(combo)/sizeof(*combo), _samp = 100;
//  TH1F* histPed[_crate][_slot][_ch][_co];
//  TH1F* histSample[_crate][_slot][_ch][_samp];

  Int_t nbins=100;
  
  // Entries over events
  // const int max = 10000;
  const int max = 200;
  const int maxids = 768;

  /*Int_t crate_select=37;
  Int_t slot_select=3;
  Int_t channel_select=8;*/
  Int_t crate = crate_select;

  Int_t nsize = (Int_t) tr->GetEntries();
  Int_t maxevents = maxids*max > nsize? nsize : maxids*max;

  const int _samp = 100;

   char string[256];

   // Double_t pedestal=100;
   // Double_t ped_sigma=1.2;

   Double_t xmin=0;
   Double_t xmax=nbins;

   // open list output file

   ofstream outlist;
   outlist.open ("pedestal_studiesA_sigma4.list");
             

   TProfile *hpedave[maxids];
   TH2D *hpedmin = new TH2D("hpedmin","Min ped",48,0,48,16,0,16);
   TH2D *hpedmax = new TH2D("hpedmax","Max ped",48,0,48,16,0,16);
   TH2D *hpeddelta = new TH2D("hpeddelta","Delta ped",48,0,48,16,0,16);
   TH2D *hpedrms = new TH2D("hpedrms","RMS ped",48,0,48,16,0,16);

  for (int j=0; j<maxids; j++) {

    // cout << " j=" << j << " nhist=" << nhist << endl;

       // define histogram to keep results

      sprintf (string,"hpedave_h%03d\n",j);
      // printf("string=%s",string);
      hpedave[j] = new TProfile(string,string,max,0,max);
  }

  Int_t evnt_ndx=0;
  Double_t Window_sample[max];

  for (Int_t ievent=0; ievent<maxevents; ievent++){
    tr->GetEntry(ievent);
    if (ievent%1000==0)  cout << " event counter=" << ievent << endl;
    

    Int_t crateid=-1;
    /*if (rocid <39)     // use for south
      crateid = (rocid-37)*12 + slot -3;
    else
      crateid = (rocid-38)*12 + slot -3;
      if (slot > 10) crateid = crateid -2;*/

    if (rocid <33)       // use for north
      crateid = (rocid-31)*12 + slot -3;
    else
      crateid = (rocid-32)*12 + slot -3;
      if (slot > 10) crateid = crateid -2;

    crateid = crateid*16 + channel;

    cout << " ievent=" << ievent << " rocid=" << rocid << " slot=" << slot << " channel=" << channel << " crateid=" << crateid << endl; 
    outlist  << " ievent=" << ievent << " rocid=" << rocid << " slot=" << slot << " channel=" << channel << " crateid=" << crateid << endl;

    /*if (crate != rocid) continue;
    if (slot != slot_select) continue;
    if (channel != channel_select) continue;*/

    for (Int_t jj=0; jj<_samp; jj++){
      
      // Get adc value
      int adc_value = (*waveform)[jj];
      hpedave[crateid]->Fill(ievent/maxids,adc_value);
    }
   
    // outlist << "ievent=" << ievent << " crateid=" << crateid << " ievent/maxids=" << ievent/maxids << endl;
  } // end over event loop

  // define function to fit
  TF1 *sine = new TF1("sine","[0]*sin([1]*x+[3])+[2]",0,200);
  sine->SetParNames("amplitude","frequency","offset","phase");
  sine->SetParameters(1,3.8,98,0);


   TCanvas *c1 = new TCanvas("c1", "c1",200,10,1000,1000);
   c1->SetBorderMode(0);
   c1->SetFillColor(0);

   c1->SetGridx();
   c1->SetGridy();
   c1->Divide(4,4);

   for (Int_t j=0; j<16; j++) {
      c1->cd(j+1);
      // hpedave[j]->Fit("sine");
      hpedave[j]->Draw();
   }

   // find max and min values of pedestal

   for (Int_t j=0; j<maxids; j++) {
     Int_t cr=j/16;
     Int_t sl=j%16;
     Double_t tempmin=hpedave[j]->GetMinimum();
     Double_t tempmax=hpedave[j]->GetMaximum();
     Double_t temprms=hpedave[j]->GetRMS(2);

     outlist << " j=" << j << " cr=" << cr << " sl=" << sl << " tempmin=" << tempmin << " tempmax=" << tempmax << endl;
     hpedmin->Fill(cr,sl,tempmin);
     hpedmax->Fill(cr,sl,tempmax);
     hpeddelta->Fill(cr,sl,tempmax-tempmin);
     hpedrms->Fill(cr,sl,temprms);
   }

   Double_t Ped0=102;  // for R3620
   // Double_t Ped0=98;  // for R3281
   Double_t delta=4; 
   
   TCanvas *c2 = new TCanvas("c2", "c2",200,10,1000,1000);
   c2->SetBorderMode(0);
   c2->SetFillColor(0);
   c2->Divide(2,2);

   c2->cd(1);
   hpedmin->Draw("colz");
   hpedmin->GetZaxis()->SetRangeUser(Ped0-delta,Ped0+delta);
   c2->cd(2);
   hpedmax->Draw("colz");
   hpedmax->GetZaxis()->SetRangeUser(Ped0-delta,Ped0+delta);
   c2->cd(3);
   hpeddelta->Draw("colz");
   hpeddelta->GetZaxis()->SetRangeUser(0,+delta);;
   c2->cd(4);
   hpedrms->Draw("colz");
   hpedrms->GetZaxis()->SetRangeUser(0,+delta);
   
   TCanvas *c3 = new TCanvas("c3", "c3",200,10,1000,1000);
   c3->SetBorderMode(0);
   c3->SetFillColor(0);
   hpedave[0]->Draw();

   // ouput legend

   /*TLegend *leg = new TLegend(0.4,0.7,0.7,0.9);
   leg->AddEntry(gavepedDI,"Real Ave, Int Ped","p");
   leg->Draw();*/


   /*sprintf (string,"Mean Pedestal=%.2f\n",pedestal);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.92,string);
   t1->SetNDC();
   t1->SetTextSize(0.07);
   t1->Draw();*/

   sprintf (string,"pedestal_studiesA_sigma4.pdf(");
   c1->SaveAs(string);
   sprintf (string,"pedestal_studiesA_sigma4.pdf)");
   c2->SaveAs(string);
   sprintf (string,"pedestal_studiesA_sigma4_c3.C");
   c3->SaveAs(string); 

   outlist.close();




}