#include "TROOT.h"
 
#include "TFile.h"
 
#include "TBenchmark.h"
 
#include "TH1F.h"
 
#include "TH2F.h"
 
#include "TCanvas.h"
 
#include "TStyle.h"
 
#include "TPad.h"
 
#include "TF1.h"
 
#include "TGraph.h"
 
#include "TSpectrum.h"
 
 
 
#include "RTUtil.h"
 
 
 
#define TDC_BIN (25./1000.) //1 TDC bin in ns
 
#define MIKRO_BIN 0.49609/1000. //1 mikro step in mm
 
 
 
#define HFILL_COLOR 18
 
 
 
void mppc(char *fname, char *plopt="t", int nChannels = 1, double fitw=1.0, char *fitf="g")
 
{
 
        char fullname[256];
 
        
 
        //get ROOT file with histograms
 
        char fnameroot[1024];
 
        TFile * rootfile;
 
        TDirectory *dir;
 
        
 
        sprintf(fnameroot, "root/%s.root", fname); 
 
        rootfile = (TFile *) gROOT->FindObject(fname); 
 
        if(rootfile==NULL) rootfile = new TFile(fnameroot);
 
        if(rootfile==NULL) {
 
          printf("Cannot open root file %s!!!\n",fnameroot); 
 
          return;
 
        }
 
        dir = (TDirectory*) rootfile;
 
        
 
        
 
        // set draw style
 
        gStyle->SetOptStat("ne");
 
        gStyle->SetPalette(1, 0);
 
        
 
        gStyle->SetPaperSize(TStyle::kA4);
 
        gStyle->SetStatBorderSize(1);
 
        gStyle->SetFrameBorderMode(0);
 
        gStyle->SetFrameFillColor(0);
 
        gStyle->SetCanvasBorderMode(0);
 
        gStyle->SetPadBorderMode(0);
 
        gStyle->SetPadColor(0);
 
        gStyle->SetCanvasColor(0);
 
        gStyle->SetStatColor(0);
 
        gStyle->SetOptFit(11);
 
        gStyle->SetOptStat();
 
        gStyle->SetPadRightMargin(0.15);
 
        gStyle->SetPadLeftMargin(0.12);
 
        //gStyle->SetTitleYOffset(1.4);
 
        
 
  TCanvas *c[16];
 
  int cc=-1;
 
  char hname[256];
 
  TH1F *hp1d; TH2F *hp2d;
 
  
 
  
 
// 2d Scan --------------------------------------------------------------------------------------------
 
 
 
  if( strchr(plopt, '2') != NULL ) {                    
 
    c[++cc] = new TCanvas("2D Scan", fname, 650, 400, 600, 600);
 
                        
 
  int i=0;
 
                        gPad->SetLogz();
 
                                sprintf(hname, "h2d%d",i); hp2d = (TH2F *) dir->Get(hname);
 
                                
 
                                (hp2d->GetYaxis())->SetLabelOffset(0.01);
 
                                hp2d->SetTitle("; X [mm]; Y[mm]");
 
                                hp2d->DrawCopy("COLZ");
 
                
 
                sprintf(fullname, "ps/%s_2d.eps", fname);
 
                c[cc]->SaveAs(fullname);
 
  }
 
 
 
// 1d Scan Y ------------------------------------------------------------------------------------------
 
 
 
  if( strchr(plopt, 'y') != NULL ) {                    
 
    c[++cc] = new TCanvas("Y Scan", fname, 250, 200, 600, 600);
 
                        
 
  int i=0;
 
                        //gPad->SetLogy();
 
                                sprintf(hname, "hnhitsy%d",i); hp1d = (TH1F *) dir->Get(hname);
 
                                
 
                                hp1d->SetTitle("Y Scan; Y; Count");
 
                                hp1d->DrawCopy();
 
                
 
                sprintf(fullname, "ps/%s_y.eps", fname);
 
                c[cc]->SaveAs(fullname);
 
  }  
 
  
 
// 1d Scan X ------------------------------------------------------------------------------------------
 
 
 
  if( strchr(plopt, 'x') != NULL ) {                    
 
    c[++cc] = new TCanvas("X Scan", fname, 250, 200, 600, 600);
 
                        
 
  int i=0;
 
                        //gPad->SetLogy();
 
                                sprintf(hname, "hnhitsx%d",i); hp1d = (TH1F *) dir->Get(hname);
 
                                
 
                                hp1d->SetTitle("X Scan; X; Count");
 
                                hp1d->DrawCopy();
 
                
 
                sprintf(fullname, "ps/%s_x.eps", fname);
 
                c[cc]->SaveAs(fullname);
 
  }  
 
  
 
// Fits of corrected TDCs -----------------------------------------------------------------------------
 
 
 
  TF1 *fg = new TF1("fg", "gaus");
 
        TF1 *fgg = new TF1("fgg", "gaus(0)+gaus(3)");
 
        fgg->SetParNames("Constant","Mean","Sigma","Constant2","Mean2","Sigma2");
 
  if( strchr(plopt, 'f') != NULL ) {
 
        c[++cc] = new TCanvas("Corrected TDCs", fname, 50*(cc+1), 0, 600, 850);
 
                //c[cc]->Divide(2,4);
 
                c[cc]->Divide(nChannels);
 
 
 
                for(int i=0;i<nChannels;i++) {
 
                  c[cc]->cd(i+1); gPad->SetLogy();
 
                        sprintf(hname, "hctdc%d",i); hp1d = (TH1F *) dir->Get(hname);
 
                        hp1d->SetTitle("cTDC; TDC [ns]; N");
 
                        (hp1d->GetXaxis())->SetRangeUser(-40*TDC_BIN,160*TDC_BIN);
 
                        if( hp1d->GetMaximum() < 50 ) continue;
 
                        //(hp1d->GetXaxis())->SetRangeUser(-0.25,0.75);
 
                        
 
                        if( strcmp(fitf, "g")==0 ) {
 
                        fg->SetParameters(hp1d->GetMaximum(), hp1d->GetBinCenter(hp1d->GetMaximumBin()), 0.05);
 
                        hp1d->Fit(fg, "0ql", "", fg->GetParameter(1)-0.1*fitw, fg->GetParameter(1)+0.05*fitw);
 
                        hp1d->Fit(fg, "ql", "", fg->GetParameter(1)-0.1*fitw, fg->GetParameter(1)+0.05*fitw);
 
                        
 
                        printf("Ch[%d] Sigma =%6.1lfps\n",i,fg->GetParameter(2)*1000.);
 
      } 
 
      else if( strcmp(fitf, "gg")==0 ) {
 
                        //fgg->SetParameters(2600.,0*TDC_BIN,2*TDC_BIN, 100.,0*TDC_BIN,10*TDC_BIN);
 
                        fgg->SetParameters(hp1d->GetMaximum(), hp1d->GetBinCenter(hp1d->GetMaximumBin()), 0.03,
 
                                          hp1d->GetMaximum(), hp1d->GetBinCenter(hp1d->GetMaximumBin())/10., 0.1);
 
                        hp1d->Fit(fgg, "0ql", "", -50*TDC_BIN, 50*TDC_BIN);
 
                        hp1d->Fit(fgg, "ql", "", fgg->GetParameter(1)-2.5*fitw*fgg->GetParameter(2), 
 
                                                 fgg->GetParameter(1)+3.5*fitw*fgg->GetParameter(2));
 
                                                 
 
                        printf("Ch[%d] Sigma =%6.1lfps\n",i,fgg->GetParameter(2)*1000.);
 
/*
 
                  fg->SetRange(-50., 200.);
 
                        fg->SetParameters(fgg->GetParameter(0), fgg->GetParameter(1), fgg->GetParameter(2));
 
                        fg->SetLineColor(2);
 
                        fg->DrawCopy("LSAME");
 
                        fg->SetParameters(fgg->GetParameter(3), fgg->GetParameter(4), fgg->GetParameter(5));
 
                        fg->SetLineColor(3);
 
                        fg->DrawCopy("LSAME");
 
*/
 
      } else {printf("Wrong fit function (parameter 4)!!!\n"); return;}
 
                }       
 
                sprintf(fullname, "ps/%s_cTDC.eps", fname); 
 
                c[cc]->SaveAs(fullname);
 
  }  
 
 
 
// TDCs -----------------------------------------------------------------------------------------------
 
 
 
  if( strchr(plopt, 't') != NULL ) {    
 
        c[++cc] = new TCanvas("Raw TDC", fname, 0, 500, 600, 430);
 
                        
 
                int i=0;
 
                        gPad->SetLogy();
 
                                sprintf(hname, "htdc%d",i); hp1d = (TH1F *) dir->Get(hname);
 
                                hp1d->SetTitle("TDC; TDC [ns]; N");
 
                                (hp1d->GetXaxis())->SetRangeUser(-10,20);
 
                                hp1d->DrawCopy();
 
                                //if( hp1d->GetMaximum() < 10 ) continue;
 
                                
 
                          fg->SetParameters(hp1d->GetMaximum(), hp1d->GetBinCenter(hp1d->GetMaximumBin()), 0.25);         
 
                          fg->SetLineColor(kRed);
 
                          hp1d->Fit(fg,"0QL", "", fg->GetParameter(1)-0.5, fg->GetParameter(1)+0.5);
 
                          hp1d->Fit(fg,"QL", "", fg->GetParameter(1)-0.8, fg->GetParameter(1)+0.8);
 
 
 
                    printf("Ch[%d] RAW Sigma =%6.1lfps\n",i,fg->GetParameter(2)*1000.);
 
                
 
                sprintf(fullname, "ps/%s_TDC.eps", fname); 
 
                c[cc]->SaveAs(fullname);
 
  } 
 
 
 
  
 
        //rootfile->Close();
 
}