#include <string>
#include <stdio.h>
#include <map>

//Function to set the display colors on the 2D Histogram
void SetColorGradientPalette() {
  const Int_t NRGBs = 5;
  const Int_t NCont = 255;
  
  Double_t stops[NRGBs] = { 0.00, 0.25, 0.50, 0.75, 1.00 };
  Double_t red[NRGBs]   = { 0.00, 0.00, 0.17, 0.87, 0.9 };
  Double_t green[NRGBs] = { 0.00, 0.81, 0.80, 1.00, 0.00 };
  Double_t blue[NRGBs]  = { 0.51, 0.60, 0.12, 0.00, 0.00 };
  TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
  gStyle->SetNumberContours(NCont);
}

void P0DRes(char* inputFile, bool option = 0) {
  //Setting graphing style parameters
  gROOT->SetStyle("Plain");
  
  //Creating a Chain associated with Tree
  TChain* xLayers = new TChain("ResidualX");
  TChain* yLayers = new TChain("ResidualY");
  
  //Adding files to Chain
  xLayers->AddFile(inputFile);
  yLayers->AddFile(inputFile);
  
  //Setting variable names
  double residual;
  double p0dule;
  
  //Getting Branches
  xLayers->SetBranchAddress("p0dule", &p0dule);
  xLayers->SetBranchAddress("residual", &residual);
  yLayers->SetBranchAddress("p0dule", &p0dule);
  yLayers->SetBranchAddress("residual", &residual);
  

  //Creating 2D Histograms
  xHist  = new TH2D("xHist","Residual Distribution X Layers",40,0,40,80,-20,20);
  yHist  = new TH2D("yHist","Residual Distribution Y Layers",40,0,40,80,-20,20);
  
  //Creating profiles
  xProfile  = new TProfile("xProfile","Residual Distribution X Layers",40,0,40);
  yProfile  = new TProfile("yProfile","Residual Distribution Y Layers",40,0,40);
  
  //Labeling axes
  xHist->GetXaxis()->SetTitle("P0Dule Number");
  xHist->GetYaxis()->SetTitle("Residual (mm)");
  yHist->GetXaxis()->SetTitle("P0Dule Number");
  yHist->GetYaxis()->SetTitle("Residual (mm)");
  
  //Labeling axes
  xProfile->GetXaxis()->SetTitle("P0Dule Number");
  xProfile->GetYaxis()->SetTitle("Residual (mm)");
  yProfile->GetXaxis()->SetTitle("P0Dule Number");
  yProfile->GetYaxis()->SetTitle("Residual (mm)");
  
  //Setting Style features
  xHist->SetLineWidth(3);
  yHist->SetLineWidth(3);
  xHist->SetStats(kFALSE);
  yHist->SetStats(kFALSE);
  
  //Setting Style features
  xProfile->SetLineWidth(3);
  yProfile->SetLineWidth(3);
  xProfile->SetStats(kFALSE);
  yProfile->SetStats(kFALSE);
  
  //Filling 2D Histograms
  Int_t nentriesx = (Int_t)xLayers->GetEntries();
  for ( int i=0; i<nentriesx; i++) {
    xLayers->GetEntry(i);
    xHist->Fill(p0dule,residual);
    xProfile->Fill(p0dule,residual);
  }
  
  Int_t nentriesy = (Int_t)yLayers->GetEntries();
  for ( int i=0; i<nentriesy; i++) {
    yLayers->GetEntry(i);
    yHist->Fill(p0dule,residual);
    yProfile->Fill(p0dule,residual);
  }
  
  //Filling lists with alignment parameters                                                          
  double x[40]={}; double y[40]={};
  double xErr[40]={}; double yErr[40]={};
  double layer[40]={};

  //Pulling mean and mean error from profiles
  for (int i=0;i<40;i++) {
    x[i] = xProfile->GetBinContent(i+1);
    y[i] = yProfile->GetBinContent(i+1);
    xErr[i] = xProfile->GetBinError(i+1);
    yErr[i] = yProfile->GetBinError(i+1);
  }
  
  //SuperP0Dule division definitions
  double highBound = 20; double lowBound = -20;
  TLine *sP0Dule1 = new TLine(7,lowBound,7,highBound);
  TLine *sP0Dule2 = new TLine(20,lowBound,20,highBound);
  TLine *sP0Dule3 = new TLine(33, lowBound, 33, highBound);
  sP0Dule1->SetLineWidth(2); sP0Dule2->SetLineWidth(2); sP0Dule3->SetLineWidth(2);
  sP0Dule1->SetLineColor(kGray); sP0Dule2->SetLineColor(kGray); sP0Dule3->SetLineColor(kGray);

  //create canvas
  TCanvas *c1 = new TCanvas("c1","Parameters of Alignment in X",200,10,700,500);
  c1->Divide(1,1);
  TCanvas *c2 = new TCanvas("c2","Parameters of Alignment in Y",200,10,700,500);
  c2->Divide(1,1);
  TCanvas *c3 = new TCanvas("c3","Profile of Alignment in X",200,10,700,500);
  c3->Divide(1,1);
  TCanvas *c4 = new TCanvas("c4","Profile of Alignment in Y",200,10,700,500);
  c4->Divide(1,1);
  
  c1->SetFrameLineWidth(3);
  c2->SetFrameLineWidth(3);
  c3->SetFrameLineWidth(3);
  c4->SetFrameLineWidth(3);

  //Draw Profiles
  c1->cd(1);
  xProfile->Draw();
  sP0Dule1->Draw(); sP0Dule2->Draw(); sP0Dule3->Draw();
  c1->Modified();
  c1->Update();
  
  c2->cd(1);
  yProfile->Draw();
  sP0Dule1->Draw(); sP0Dule2->Draw(); sP0Dule3->Draw();c2->Modified();
  c2->Update();
  
  //Draw 2D Histograms
  SetColorGradientPalette();
  c3->cd(1);
  xHist->Draw("COLZ");
  sP0Dule1->Draw(); sP0Dule2->Draw(); sP0Dule3->Draw();
  c3->Modified();
  c3->Update();
  
  c4->cd(1);
  SetColorGradientPalette();
  yHist->Draw("COLZ");
  sP0Dule1->Draw(); sP0Dule2->Draw(); sP0Dule3->Draw();
  c4->Modified();
  c4->Update();

  //To save output file option must equal 1
  if (option == 0) return;

  //Check to see if P0DInt.txt exists and read contents   
  //Contents of P0DInt.txt formatted like:
  // P0DuleNumber XParameter XParameterError YParameter YParameterError
  int inP0Dule=0; float inX=0; float inXErr=0; float inY=0; float inYErr=0;
  float xList[40]={}; float yList[40]={}; float xErrList[40]={}; float yErrList[40]={};
  FILE * P0DInt1;
  FILE * P0DInt;
  P0DInt1 = fopen ("P0DInt.txt","r");
  if (P0DInt1 == NULL) {
    perror ("File does not exist! Creating new file");
  }
  else {
    for (int i=0; i<40; i++) {
      fscanf (P0DInt1, "%i    %f      %f      %f      %f", &inP0Dule, &inX, &inXErr, &inY, &inYErr);
      xList[i] = inX;
      xErrList[i] = inXErr;
      yList[i] = inY;
      yErrList[i] = inYErr;
    }
    fclose (P0DInt1);
  }
  
  //Record new values for alignment
  P0DInt = fopen ("P0DInt.txt","w");
  for (int i=0; i<40; i++) {
    fprintf (P0DInt,"%i     %f      %f      %f      %f\n",i,
       xList[i]-x[i],sqrt(pow(xErrList[i],2)+pow(xErr[i],2)),
       yList[i]-y[i],sqrt(pow(yErrList[i],2)+pow(yErr[i],2)));
  }
}