void linear_fit(void)
{
// plot data and fit to linear function.

// #include <TMath.h>
#include <TRandom.h>

gROOT->Reset();
//TTree *Bfield = (TTree *) gROOT->FindObject("Bfield");
gStyle->SetPalette(1,0);
gStyle->SetOptStat(kFALSE);
gStyle->SetOptFit(kFALSE);
// gStyle->SetOptFit(kTRUE);
// gStyle->SetOptFit(1111);
gStyle->SetPadRightMargin(0.15);
gStyle->SetPadLeftMargin(0.15);
gStyle->SetPadBottomMargin(0.15);
gStyle->SetFillColor(0);
//

   char string[256];
   Int_t j,jj;

   // input data 
 
   Int_t npts=11;
   // add dummy point on end to be able to plot beyond last data
   Float_t x[12]={-1,1,2,3,4,5,6,7,8,9,10,20};
   Float_t y[12]={-1,1.1,2.2,3.3,4.4,5.5,6.6,7.7,8.8,9.9,11,-1};
   // Float_t y[12]={-1,1,1,1,1,1,1,1,1,1,1,-1};
   Float_t ex[12]={0.01,0.01,0.01,0.01,0.01,0.01,0.01,0.01,0.01,0.01,0.01,0.01};
   Float_t ey[12]={0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5};
   Float_t ydif[12]={1,-1,0,1,0,-1,1,-1,0,-1,0,-1};

  Float_t yoffset=5;
  Float_t yslope=0;
   for (j=1;j<npts+1;j++) {
        y[j] = y[j] + x[j]*yslope + yoffset + ydif[j]*ey[j];
       }

   // compute parameters by hand

   Double_t yi=0;
   Double_t xi=0;
   Double_t xi2=0;
   Double_t xiyi=0;
   Double_t n=0;
   
   for (j=1;j<npts;j++) {
         n++;
         yi = yi+y[j];
         xi = xi+x[j];
         xi2 = xi2 +x[j]*x[j];
         xiyi = xiyi + x[j]*y[j];
         }

     //  compute parameters

    Double_t delta = n*xi2-xi*xi;
    Double_t a = (yi*xi2-xi*xiyi)/delta;
    Double_t b = (n*xiyi - xi*yi)/delta;
    // Double_t siga = ey[0]/sqrt(n);
    Double_t siga = ey[0]*sqrt(xi2/delta);
    // Double_t sigb = ey[0]/sqrt(xi2);
    Double_t sigb = ey[0]*sqrt(n/delta);
    printf ("a = %.2f, b=%.2f, siga=%.2f, sigb=%.2f \n",a,b, siga, sigb);
 

   TGraphErrors *gr = new TGraphErrors (npts,x,y,ex,ey);
     
   //
   TCanvas *c1 = new TCanvas("c1","Plot linear function",200,10,700,500);
   Float_t xmin = 0;
   Float_t xmax = 15;
   Float_t ymin = 0;
   Float_t ymax = 20;
   c1->SetGridx();
   c1->SetGridy();
   c1->SetBorderMode(0);
   c1->SetFillColor(0);

   gr->GetXaxis()->SetRangeUser(xmin,xmax);
   gr->GetYaxis()->SetRangeUser(ymin,ymax);
   gr->GetXaxis()->SetTitleSize(0.07);
   gr->GetYaxis()->SetTitleSize(0.07);
   gr->GetXaxis()->SetTitle("x_{i}");
   gr->GetYaxis()->SetTitle("y_{i}");
   gr->SetMarkerColor(2);
   gr->SetMarkerStyle(20);
   gr->SetTitle("");
   gr->Fit("pol1","","",xmin,xmax);
   TF1 *fit = gr->GetFunction("pol1");
   Float_t p0 = fit->GetParameter(0); 
   Float_t p1 = fit->GetParameter(1); 
   Float_t p0err = fit->GetParError(0); 
   Float_t p1err = fit->GetParError(1); 
   Float_t chi2= fit->GetChisquare(); 
   Int_t ndf= fit->GetNDF();  
   Float_t prob= fit->GetProb(); 
   gr->Draw("AP");

   // TF1 *calc = new TF1("calc","x-53.2",53,75);
   // calc->SetLineStyle(2);
   // calc->Draw("same");

   /*sprintf (string,"Generated slope =%.2f\n",yslope);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.85,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();*/

   sprintf (string,"Fit function y(x) = a + bx");
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.85,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();

   sprintf (string,"a = %.2f #pm%.2f\n",a,siga);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.78,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();
   sprintf (string,"b = %.2f #pm%.2f\n",b,sigb);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.72,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();

   /*sprintf (string,"a = %.2f #pm%.2f\n",p0,p0err);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.78,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();
   sprintf (string,"b = %.2f #pm%.2f\n",p1,p1err);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.72,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();
   sprintf (string,"#chi^{2}/NDF = %.2f/%d\n",chi2,ndf);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.66,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();
   sprintf (string,"Probability= %.2f\n",prob);
   printf("string=%s",string);
   t1 = new TLatex(0.2,0.60,string);
   t1->SetNDC();
   t1->SetTextSize(0.05);
   t1->Draw();*/

   // c1->SaveAs("linear_fit_c1.pdf");
   c1->SaveAs("linear_fit_c1.png");
}