//
// Author:E.Gorini, Physics Department, Salento University, 9/2009
//
// Modifications by: Mosè Giordano
//
#include <TH1>
#include <TF1>
#include <TLegend>
#include <TCanvas>
#include <Riostream>
#include <vector>

void GenFit(string file){
  // initialization
  gROOT->Reset();
  gStyle->SetOptStat(1);
  gStyle->SetOptFit(111111);
  
  // here open the data file
  ifstream in;
  // uscire se file non esiste. 
  if(in.open(file.c_str()))exit; 
  // reads plot title and X and Y labels
  string Title,Xtitle,Ytitle;
  getline(in,Title);
  getline(in,Xtitle); 
  getline(in,Ytitle);
  cout << "Plot Title is: \""<< Title <<"\", ";
  cout << "x axis title is: \""<< Xtitle <<"\", ";
  cout << "y axis title is: \""<< Ytitle <<"\".\n";
  // reads function
  string Function;
  getline(in,Function);
  cout << "Function to fit is: "<< Function << "\n";
  // read parameter number and names
  int numpar; 
  in >> numpar;
  vector<string> parnames;
  for (int np=0;np<numpar;np++){
    string sbuf;
    in >> sbuf;
    parnames.push_back(sbuf);
  }
  cout << "Number of Parameters is "<< parnames.size() << " with value: ";
  // read parameter initial values
  TVector pvalue(numpar);
  for (int np=0;np<numpar;np++)in >> pvalue(np);
  for (int np=0;np<numpar;np++){
    cout << parnames[np]<< "="<< pvalue[np] << ", ";  
  }
  cout << "\n";
  // create the vectors, limited at 100 for the points, 20 for the parameters
  // vectors will be anyway resized if necessary.....
  Int_t n = 100, m = 20;
  TVectorF xx(n),yy(n),ex(n),ey(n),pull(n);
  Float_t xxx,yyy,exx,eyy;
  bool flag;
  Int_t nlines=0;
  // now reads data x,y, ex,ey 
  while (in.good()) {
    in >> flag >> xxx >> yyy >> exx >> eyy;
    if(flag){
      xx(nlines)=xxx;
      yy(nlines)=yyy;
      ex(nlines)=exx;
      ey(nlines)=eyy;
      nlines++;
      flag=false;
    } 
  }
  cout << "found " << nlines << " data points, ";
  //resize here TVectorF
  xx.ResizeTo(nlines);
  yy.ResizeTo(nlines);
  ex.ResizeTo(nlines);
  ey.ResizeTo(nlines);
  //   xx.Print();
  in.close();
    
  // create and set canvas style
  TCanvas *c1 = new TCanvas("c1",Title.c_str(),10,10,900,700);
  c1->SetFillColor(17);
  c1->SetFrameFillColor(40);
  c1->SetGrid(); 
  c1->Divide(1,2);
  c1->cd(1);  
  
  // plot the measurement 
  gr = new TGraphErrors(xx,yy,ex,ey);
  gr->SetTitle(Title.c_str());
  gr->SetMarkerColor(4);
  gr->SetMarkerStyle(21);
  gr->GetXaxis()->SetTitle(Xtitle.c_str());
  gr->GetYaxis()->SetTitle(Ytitle.c_str());
  gr->Draw("AP");
  c1->Update();
  
  // create the function
  TF1 *fff = new TF1("fff",Function.c_str(),1,1000);
  
  // set parameter names and initial values
  for (int np=0;np<numpar;np++){
    fff->SetParName(np,parnames[np].c_str());
    fff->SetParameter(np,pvalue(np));
  }	
  
  // fit   
  cout << "Start Fitting.....\n"; 
  gr->Fit("fff");
  
  c1->cd(2);
  // make plot of pulls
  for (int bin=0;bin<nlines;bin++)pull(bin) = (yy(bin)-fff->Eval(xx(bin),0.0,0.0))/ey(bin);
  
  re = new TGraph(xx,pull);
  re->SetTitle(Title.c_str());
  re->SetMarkerColor(4);
  re->SetMarkerStyle(21);
  re->SetMaximum(5.0);
  re->SetMinimum(-5.0);
  // draw a line 
  tl= new TLine(xx(0),0.0,xx(nlines-1),0.0);
  tl->Draw();
  re->Draw("AP");
  tl->Paint("AP");
  c1->Update();
  // save the plots in pdf
  
  string TitleOut=Title+".pdf";
  c1->Print(TitleOut.c_str());  
}