/// ecal.C (based on SCA_data.C)
/// ----------
/// This is a C++ root macro written by Jeff Martin (2012),
/// commented by Blair Jamieson (2013).
///
/// This macro shows how to create a graph with error bars,
/// from a text file of data in two columns.
///
/// Note: C++ comments start with // (or ///)
// The following function takes a TF1 root function
// object, and from that prints all of the fit
// information from the last root Fit that was done
void PrintMoreFitInfo( TF1* tf1 ){
int npars = tf1->GetNpar();
// extract the fit information:
// get the fit values and uncertainties:
double * apars = tf1->GetParameters();
double * aparerrs = tf1->GetParErrors();
// create storage for correlation and covariance matrices
double *covmat;
double *cormat;
covmat = malloc( npars * npars * sizeof( double ) );
cormat = malloc( npars * npars * sizeof( double ) );
// get the covariance matrix:
gMinuit->mnemat( covmat,npars);
// print the covariance matrix:
std::cout<<" Covariance Matrix: "<GetParName(i)<<" ";
for(int j=0; jGetChisquare();
double andof = tf1->GetNDF();
double apval = tf1->GetProb();
printf(" Fit chi2/NDOF = %5.2lf / % 5d, p-value=%5.2lf\n",achi2,int(andof),apval);
// calculate the correlation matrix from covariance matrix:
for (int i=0; iGetParName(i)<<" ";
for(int j=0; jSetOptStat(1);
gStyle->SetOptFit(1111);
/// Create a canvas Note that TCanvas is a C++ class that is from
/// the ROOT library. This class describes an object that is a
/// canvas (window) to draw graphics in.) c1 is now an instance of
/// a TCanvas, with the title "A Simple Graph with error bars"
/// Note that the "new" command in C++ is used to allocate new
/// memory for an object, and return a pointer to the object.
/// A pointer is just a fancy way of saying the address in memory
/// where the object is stored.
/// Note that the canvas is created using the new operator, rather
/// than as an object, because of the way persistence works in C++.
/// Anything created using the new operator stays in memory until
/// the delete operator is called. This is useful because we want
/// the canvas to stay there even when this script ends, so that we
/// can see our graph.
c1 = new TCanvas("c1","Energy Calibration");
/// Define some variables used to make the graph define MAX to be a
/// constant integer for the maximum number of points to read from
/// the input file (the maximum number of points in graph)
const Int_t MAX=1000;
/// Define a variable to use as an index into our arrays as we read
/// in the data as i (C++ index always start counting from zero)
Int_t i=0;
/// Define the x, y, error in x, and error in y arrays of floating
/// point values (Double_t is a double length floating point number)
/// The [MAX] defines the maximum length of the array.
Double_t x[MAX],y[MAX],ex[MAX],ey[MAX];
std::string labels[MAX];
/// ifstream is a standard C++ class that defines objects which read
/// in from a file. fin is an ifstream object (an instance of the
/// ifstream class. Note that in this case we create an instance
/// without using the new command, so the ifstream object will
/// automatically be deleted when the script is done
ifstream fin;
/// open is a method (function) of the ifstream class, which is used
/// to define which file to read from.
fin.open("ecal.dat");
/// In C++ the while loop is used to loop until some condition
/// inside the brackets is met. while(1) will continue to loop
/// over the comands inside the curly brackets, and while(0) will
/// break out of the loop and proceed to the next command after
/// the end of the curly brackts.
///
/// fin>> is the ifstream operator to read a value from the file
/// opened by fin.open. Note that it tries to read in a value of
/// the same type as x[i] and y[i] (floating point numbers)
/// if the fin>> operator doesn't get any data it returns zero,
/// ending the loop
while(fin>>labels[i]>>x[i]>>ex[i]>>y[i]>>ey[i]) {
ey[i] += 4.0; /// add a fixed systematic uncertainty to the fit channel
i++; /// ++ increments by 1 the value in i
}
/// TGraphErrors is a ROOT class which defines a graph with
/// error bars. The arguments to create a new instance are:
/// i the number of points in the graph
/// x the array of x values
/// y the array of y values
/// ex the array of errors in x values
/// ey the array of errors in the y values
/// Below gr is a pointer to a newly created TGraphErrors object
gr = new TGraphErrors(i,x,y,ex,ey);
/// Several of the methods (functions) which operate on TGraphErrors
/// to set the title, xaxis and yaxis labels are below
gr->SetTitle("Beta spectroscopy energy calibration");
gr->GetXaxis()->SetTitle("True Energy (keV)");
gr->GetYaxis()->SetTitle("Fit Channel (ADC)");
/// The command below sets the shape of the data point markers
gr->SetMarkerStyle(21);
/// The command below draws the graph with axes "A", and points "P"
/// it is drawn to the currently active cnvas (only c1 exists now,
/// so that is where it is drawn)
gr->Draw("AP");
// Update is a method of the TCanvas class, which updates the canvas
c1->Update();
/// Define a line for fitting
///
/// To simplify calculation of energies, we want to invert the parameters being fit.
/// Ie. the fit is to Channel = A*Energy +B
/// but we would perfer to have constants m, b in : Energy = m*Channel + b
/// using C=Channel, E=Energy:
/// E = m*C + b
/// E-b = m*C
/// C = (1/m)*E + (-b/m)
///
/// Therefore fit to the function "1.0/[0] * ( x - [1])"
/// The TF1 object has several different ways of being initialized.
/// In the version below, the arguments in the constructor of the TF1 are:
/// 1) "fline" is the arbitrarily chosen unique name
/// 2) "1.0/[0]*x+[1]/[0]" is the function definition as a character string,
/// the [0], [1] are treated as parameters, x,y,z are variables
/// 3) 0. is the minimum range for the function
/// 4) 1000. is the maximum range for the function
TF1 * fline = new TF1( "fline", "1.0/[0]*(x-[1])", 0.0, 1200.0 );
/// Set labels for parameters
fline->SetParName(0, "1/Slope ");
fline->SetParName(1, "-Intercept/Slope ");
/// Set the initial guess for the parameters of the line
fline->SetParameters( 1.0, 0.0 );
/// do the fit
gr->Fit( fline );
/// Make a residual graph
Double_t resid[MAX];
Double_t reside[MAX];
for (int j=0; j*Eval( x[j] );
}
c2 = new TCanvas("c2","Energy Calibration Residual");
gresid = new TGraphErrors(i,x,resid,ex,ey);
gresid->SetTitle("Energy calibration residual");
gresid->GetXaxis()->SetTitle("True Energy (keV)");
gresid->GetYaxis()->SetTitle("Fit Residual (ADC)");
/// The command below sets the shape of the data point markers
gresid->SetMarkerStyle(21);
/// The command below draws the graph with axes "A", and points "P"
/// it is drawn to the currently active cnvas (only c1 exists now,
/// so that is where it is drawn)
gresid->Draw("AP");
/// Print the covariance matrix:
PrintMoreFitInfo( fline );
}
*