KaliVeda  1.12/06
Heavy-Ion Analysis Toolkit
KVAlphaCalibration.h
Go to the documentation of this file.
1 /*
2  $Id: KVAlphaCalibration.h, v 1.0 2019/05/14 17:00:00 lemarie Exp $
3  $Revision: 1.0 $
4  $Date: 2019/05/14 17:00:00 $
5  $Author: lemarie $
6 */
7 
8 #ifndef _KVCALIBRATION_H_
9 #define _KVCALIBRATION_H_
10 
11 #include "TH1.h"
12 #include "TF1.h"
13 #include "TSpectrum.h"
14 #include "TGraph.h"
15 #include <string>
16 #include <TMath.h>
17 #include "TCanvas.h"
18 #include "TSystem.h"
19 #include "TFile.h"
20 #include "TStyle.h"
21 #include "TROOT.h"
22 #include <vector>
23 #include <cstdlib>
24 #include <iostream>
25 #include <limits>
26 
71 
72 protected :
80  std::vector<double> GaussianFitResults;
81  std::vector<double> GaussianFitResultsError;
82  double FunctionToFit(double* x, double* par); //the model we use to fit the peak
83  std::vector<double> InitializationPeak;
84 
85 
86 private :
91  std::vector<double> MeanOfPeak;
92  std::vector<double> IntensityOfPeak;
93  void HistoInit(TH1* h); //-
95  void SetFunction(std::string FunctionName = "x*[0]");
96 
97 public :
98 
99  KVAlphaCalibration(int NumberOfPeak_); //-
100  KVAlphaCalibration(int NumberOfPeak_, TH1* h); //-
102 
103  void Init(int);
104  void SetHisto(TH1* h);
105  void SetParameters(double SigmaOfTSpectrum_ = 1., double SigmaOfGaussian_ = 1., double ThresholdOfTSpectrum_ = 0.5, double IsOriginAtZero_ = false);
106  void AddPeak(double Energy_, double Intensity_); //-
107  void SetHistRange(double xmin, double xmax);
108  double GetGaussianFitParameter(int); //-
109  double GetGaussianFitParError(int);
110  double GetInitializationFitParameter(int);
111 
112 
113 
114  void FitInit(bool debug_ = false);
115  void FitSpectrum(bool debug_ = false);
116  void FitAll(bool debug_ = false);
117 
118  void DrawResult(bool WhatToDraw = true);
119  void PrintResult(void);
120 
121 };
122 
123 #endif
Set up and run the calibration of siliciums.
void AddPeak(double Energy_, double Intensity_)
void FitInit(bool debug_=false)
void SetHistRange(double xmin, double xmax)
double FunctionToFit(double *x, double *par)
double GetGaussianFitParError(int)
void FitSpectrum(bool debug_=false)
double InitializationFitResults[2]
Array that contain the results of the initialization fit.
~KVAlphaCalibration()
Default destructor.
std::vector< double > GaussianFitResults
Array that contains results of the final fit.
double SigmaOfTSpectrum
Width of the peaks that TSpectrum will search, chosen by user.
TF1 * GaussianFit
function that will fit the peaks
std::vector< double > IntensityOfPeak
array of peaks amlitude. Set by user
double GetGaussianFitParameter(int)
std::vector< double > InitializationPeak
TH1 * Histo
histogram that contains the data to fit. Set by user
void DrawResult(bool WhatToDraw=true)
double GetInitializationFitParameter(int)
TF1 * InitializationFit
function that will do the initial fit
TF1 GetFunction(void)
std::vector< double > MeanOfPeak
array peaks energy. Set by user
double ThresholdOfTSpectrum
Lowest peak that the TSpectrum has to search for. It is the ratio of the highest peak that the TSpect...
void FitAll(bool debug_=false)
This function calls the FitInit and FitSpectrum function.
TSpectrum * spec
the method that will search for peaks
std::vector< double > GaussianFitResultsError
Array that contains the error of the results of the final fit.
void SetFunction(std::string FunctionName="x*[0]")
void SetParameters(double SigmaOfTSpectrum_=1., double SigmaOfGaussian_=1., double ThresholdOfTSpectrum_=0.5, double IsOriginAtZero_=false)
KVAlphaCalibration(int NumberOfPeak_)
int NumberOfPeak
Number of Peak present in the histogram. Chose by user.
double SigmaOfGaussian
Width of the peak that the model will fit chosen by user. It is possible that it has to be different ...
TGraph * factorGraph
Graph where the initial fit is done.
void SetHisto(TH1 *h)
Set the histogram that contains the data.