ErrAnalLab15 Fitter

From New IAC Wiki
Jump to navigation Jump to search

//MyLinFit(10,-1.2,-0.85,.042,0.045);

  1. ifndef __CINT__
  2. include "TMatrixD.h"
  3. include "TVectorD.h"
  4. endif

// data arrays are global variables

Int_t NumPoints=0;

Float_t x[100],y[100],errory[100];



void LoadDataArrays(Char_t *filename, Int_t Weight) //void LoadDataArrays(string filename, Int_t Weight) {

 ifstream in;


 NumPoints=0;


 in.open(filename);
 while(in.good() && NumPoints < 100){
   if(Weight==1.0)
     {

in >> x[NumPoints] >> y[NumPoints] ;

errory[NumPoints]=0.05;

     }
   else
     {

in >> x[NumPoints] >> y[NumPoints] ;

if(y[NumPoints]>0) errory[NumPoints]=sqrt(y[NumPoints]);

else errory[NumPoints]=abs(y[NumPoints]);

if(errory[NumPoints]==0) errory[NumPoints]=1.0;


     }
   NumPoints++;
 }
 in.close();

} void DataArrayPrint(void) {

 if(NumPoints)
   {
     printf("__X__\t__Y__\tEror Y\n");
     for(Int_t i=0;i<NumPoints;i++)
     	 printf("%g\t%g\t%g\n",x[i],y[i],errory[i]);
   }
 else
   {
     printf("Data arrays are not filled\n execute function \n\tLoadDataArrays\n");
   }

}


void LinContour(Int_t steps, Float_t a1min, Float_t a1max, Float_t a2min, Float_t a2max)) {

 //  LinContour(100,-1.2,-0.9,0.042,.048) 
 Int_t nRowsCols=2;
 TH2D *ChiMin=new TH2D("ChiMin","ChiMin",steps,a1min,a1max,steps,a2min,a2max);


 if(NumPoints==0)
   {
     printf("Data arrays are not filled\n execute function \n\tLoadDataArrays\n");
     return;
   }
   
 // ALl the data is loaded in 
 //     Alpha()=-1.00503 +- 0.0210647   0.0431444 +- 0.000360375    
 //ChiSqrd = 43.4671
 
 Double_t ChiSqrd,Yfit,Alpha[100];
 Double_t a1Step=(a1max-a1min)/steps;
 Double_t a2Step=(a2max-a2min)/steps;


 //  printf("a2Step=%g\t a2Step=%g\n",a1Step,a2Step);


 //  Alpha[0]=-1.00503;
 //  Alpha[1]=0.043144;


 Alpha[1]=a2min;
 for(Int_t k=0;k<steps;k++)
   {
     Alpha[1]+=a2Step;
     Alpha[0]=a1min;
     for(Int_t l=0;l<steps;l++)

{

     Alpha[0]+=a1Step;


     ChiSqrd=0;


     for(Int_t j=0;j<NumPoints;j++)

{ Yfit=0;

for(Int_t m=0;m<nRowsCols;m++) { Yfit+=Alpha[m]*pow(x[j],m);

} // printf("y[j]-Yfit=%g-%g=%g\n",y[j],Yfit,y[j]-Yfit);

ChiSqrd+= (y[j]-Yfit)*(y[j]-Yfit)/errory[j]/errory[j];

}

     //      printf("A1=%g\tA2=%g\tChiSqrd=%g\tnu=%d\n",Alpha[0],Alpha[1],ChiSqrd,NumPoints-nRowsCols);
     //      TH2D *ChiMin=new TH2D("ChiMin","ChiMin",steps,a1min,a1max,steps,a2min,a2max);
     ChiMin->Fill(Alpha[0],Alpha[1],ChiSqrd);


}

   }  
 //  ChiMin->Draw();
 ChiMin->Draw("lego");


}