#include "TFile.h"
#include "TTree.h"
#include "TH1.h"
#include "TH2F.h"
#include "TF1.h"
#include "TROOT.h"
#include "TString.h"
#include <stdio.h>
#include <stdlib.h>
#include <bitset>

Double_t fitFunction(Double_t *x, Double_t *par) {
    // P(phi)=phi_1+P_0*Exp(-phi/phi_0)  determine phi_0
    Double_t fitf=par[0]+par[1]*TMath::Exp(-par[2]*x[0]);
    return (fitf);
}


int plot2(){
   const int n = 17051;
   Double_t xx[n]={0},y[n]={0},ex[n]={0},ey[n]={0};
   TCanvas *c1 = new TCanvas("c1","multipads",900,700);
   c1->SetGrid();
   
   //gStyle->SetOptStat(0);c1->Divide(2,2);
   Double_t DoseDeposited;
   Double_t Polarization;
   Int_t ExperimentNum;
   
    TFile *f = new TFile("PDPout.root");

   f->Get("PDPvar");
   // stuff
   PDPvar->SetBranchAddress("DoseDeposited",&DoseDeposited);
   PDPvar->SetBranchAddress("Polarization",&Polarization);  
   PDPvar->SetBranchAddress("ExperimentNum",&ExperimentNum);
   
   /*
   ExperimentNum coding
   # - Material - Temperature - Beam Voltage - Magnetic Field - Experiment Name
   1 - LiD - 1K - 48 GeV - N/A - Nucl. instr. and meth res. manual p35
   2 - NH3 - N/A - N/A - N/A - Annealing, Proceeding of the 4th... p28
   3 - NH3 - 1K - 2.5 GeV - 2.5 T - Proceeding of the 4th... p29 graph 3
   4 - NH3 - 1K - 2.5 GeV - 2.5 T - Proceeding of the 4th... p30 graph 1
   5 - pre-irradiated NH3 - 1K - 2.5 GeV - 2.5 T - Proceeding of the 4th... p30 graph 1   
   6 - pre-irradiated ND3 - 1K - 2.5 GeV - 2.5 T - Proceeding of the 4th... p30 graph 2
   7 - ND3 - N/A - N/A - 2.5 T - Annealing - Proceeding of the 4th... p31 graph 2
   http://p25ext.lanl.gov/~ming/SantaFe-DY/talks/Monday/Maxwell_DY.pdf
   8 - NH3 - 1K - 6.0 GeV - 5 T - SANE 72417 -  negative pol.
   9 - NH3 - 1K - 6.0 GeV - 5 T - SANE 72417 - positive pol.
   10- NH3 - 1K - 6.0 GeV - 5 T - SANE 72428 - negative pol.
   11- NH3 - 1K - 6.0 GeV - 5 T - SANE 72428 - positive pol.
   12- NH3 - 1K - 6.0 GeV - 5 T - SANE 72986 - negative pol.
   13- NH3 - 1K - 6.0 GeV - 5 T - SANE 72986 - positive pol.
   14- run 72162 positive pol.
   15- run 72162 negative pol.
   16- run 72164 negative pol.
   17- run 72164 positive pol.
   18- run 72657 negative pol.
   19- run 72657 positive pol.
   20- run 72669 negative pol.
   21- run 72669 posiitve pol.
   22- run 72824 negative pol.
   23- run 72824 positive pol.
   24- run 72837 negative pol.
   25- run 72837 positive pol.
   26- run 72984 negative pol.
   27- run 72984 positive pol.
   28- run 72986 negative pol.
   29- run 72986 positive pol.
   30- jlab2001 proton polarization
   31- eg1 dvcs dose top negative pol.
   32- eg1 dvcs dose top positive pol.
   33- eg1 dvcs dose bottom negative pol.
   34- eg1 dvcs dose bottom positive pol.
    */
   Int_t nentries = (Int_t)PDPvar->GetEntries();   
   //cout<<nentries<<endl;
   Int_t j;
   for (j=0;j<nentries;j++) { PDPvar->GetEntry(j);
       if (ExperimentNum == 25)
       {   //cut on dose to isolate points in fit
	   //if (DoseDeposited<7.1 && DoseDeposited>3.4) {			   
	   xx[j]  = DoseDeposited;
           y[j]  = -Polarization;
	   //}
       }     
   }
   
   TGraphErrors *gra = new TGraphErrors(j,xx,y,ex,ey);
   gra->SetTitle("TGraphErrors Example");
   gra->SetTitle("TGraphErrors Example");
   gra->SetMarkerColor(4);
   gra->SetMarkerStyle(20);
   gra->SetMarkerSize(0.8);
   gra->Draw("AP");     
   c1->Update();
   
   /*
   //Fit Function through range
   TF1 *fitFun = new TF1("fitFun",fitFunction,3.4,7.1,3);
   //set parameters of exponential fit if needed 
   //(probably have to play with these for each fit)
   fitFun->SetParLimits(0,0,100);   
   fitFun->SetParLimits(1,0,1000);
   fitFun->SetParLimits(2,0,100);   
   gra->Fit("fitFun","RV+","ep");	
   */

   return 0;
}