GamR  0.0.0
GammaROOT
Loading...
Searching...
No Matches
GaussianDistr.hh
Go to the documentation of this file.
1#ifndef GAMROOT_ANGDIST_GAUSSIANDISTR_HH
2#define GAMROOT_ANGDIST_GAUSSIANDISTR_HH
3
4/* STD */
5#include <string>
6#include <vector>
7
8/* ROOT */
9#include <TF1.h>
10#include <TH1.h>
11
12/*GamR */
13#include <utils/Utilities.hh>
14#include <toolkit/Gate.hh>
15#include <nucleus/Transition.hh>
16
17namespace GamR {
18 namespace AngDist {
20 private:
21 std::vector<double> Bhfs;
22 std::vector<double> RelAmps;
23 int nSamples;
24 int Limit;
25
26 std::vector<double> ParVals;
27 std::vector<double> ParLimsLow;
28 std::vector<double> ParLimsHigh;
29 std::vector<double> ParErrs;
30 std::vector<int> ParFix;
31
32 class GlobalChiSquare;
33
34 public:
35 double mean;
36 double FWHM;
37 double amp;
38 double offset;
39 double g;
40
41 double Low;
42 double High;
43
44 TH1 *Hist;
45 TF1 *FitFunc;
46
48
49 GaussianDistr(double L, double H, std::vector<double> pars) : nSamples(50), Limit(20), Low(L), High(H), ParVals(pars) {
50 ParLimsLow = {-1.0, -1.5, -100.0, 0.0, -100.0};
51 ParLimsHigh = {1.0, 1.5, 100.0, 50.0, 100.0};
52 ParFix = {0, 1, 0, 0, 0};
53 ParErrs = {0.01, 0.01, 0.01, 0.01, 0.01};
54 //set initial guesses into member variables
55 SetPars();
56
57 FixedField = 0;
58 }
59
60 void Normalise();
61 std::vector<double> GetBhfs() { return Bhfs; }
62 std::vector<double> GetRelAmps() { return RelAmps; }
63 void SetSamples(int i) { nSamples = i; }
64 void SetLimit(int i) { Limit = i; }
65
66 void SetAmp(double a) { amp = a; ParVals[0] = amp;}
67 void SetG(double gfac) { g = gfac; ParVals[1] = g;}
68 void SetMean(double m) { mean = m; ParVals[2] = mean;}
69 void SetFWHM(double w) { FWHM = w; ParVals[3] = FWHM;}
70 void SetOffset(double o) { offset = o; ParVals[4] = offset;}
71 void SetHist(TH1* hist) { Hist = hist; }
72
73 void SetField();
74
75 double GetG() { return g; }
76 double GetGErr() { return ParErrs[1]; }
77 double GetAmp() { return amp; }
78 double GetAmpErr() { return ParErrs[0]; }
79 double GetMean() { return mean; }
80 double GetMeanErr() { return ParErrs[2]; }
81 double GetFWHM() { return FWHM; }
82 double GetFWHMErr() { return ParErrs[3]; }
83 double GetOffset() { return offset; }
84 double GetOffsetErr() { return ParErrs[4]; }
85
86 void FixField() { ParFix[2] = 1; ParFix[3] = 1; FixedField = 1;}
87 void UnFixField() { ParFix[2] = 0; ParFix[3] = 0; FixedField = 0;}
88 void FixG() { ParFix[1] = 1; }
89 void UnFixG() { ParFix[1] = 0; }
90 void FixAmp() { ParFix[0] = 1; }
91 void UnFixAmp() { ParFix[0] = 1; }
92 void FixOffset() { ParFix[4] = 1; }
93 void UnFixOffset() { ParFix[4] = 0; }
94 std::vector<double> GetPars() { return ParVals; }
95 void SetPars() { SetPars(ParVals); }
96 void SetPars(std::vector<double> pars) { ParVals = pars; amp = pars[0]; g = pars[1]; mean = pars[2]; FWHM = pars[3]; offset = pars[4]; }
97
98 double RatioFunc(double x, double a, double gfac, double off);
99 double operator()(double *x, double *p);
100 TF1* GetFunc();
101
102 void Fit(int quiet);
103
104 void Print() {
105 std::cout << " g amp offset Bhf FWHM " << std::endl;
106 std::cout << GamR::Utils::wrresult(GetG(), GetGErr()) << " "
110 << GamR::Utils::wrresult(GetFWHM(), GetFWHMErr()) << " " << std::endl;
111 }
112
113
114 };
115
116 class DoubleGaussianDistr { // this is for two precessions overlayed, same field
117 private:
118 std::vector<double> Bhfs;
119 std::vector<double> RelAmps;
120 int nSamples;
121 int Limit;
122
123 std::vector<double> ParVals;
124 std::vector<double> ParLimsLow;
125 std::vector<double> ParLimsHigh;
126 std::vector<double> ParErrs;
127 std::vector<int> ParFix;
128
129 class GlobalChiSquare;
130
131 public:
132 double amp1;
133 double amp2;
134 double g1;
135 double g2;
136 double tau1;
137 double tau2;
138 double n1n2;
139 double mean;
140 double FWHM;
141 double offset;
142
143 double Low;
144 double High;
145
146 TH1 *Hist;
148
150
151 DoubleGaussianDistr(double L, double H, std::vector<double> pars) : nSamples(50), Limit(20), Low(L), High(H), ParVals(pars) {
152 ParLimsLow = {-1.0, -1.0, -1.5, -1.5, 0.0, 0.0, 0.0, -100.0, 0.0, -100.0};
153 ParLimsHigh = {1.0, 1.0, 1.5, 1.5, 1000.0, 1000.0, 2.0, 100.0, 50.0, 100.0};
154 ParFix = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
155 ParErrs = {0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01};
156 SetPars();
157 FixedField = 0;
158 }
159
160 void Normalise();
161 std::vector<double> GetBhfs() { return Bhfs; }
162 std::vector<double> GetRelAmps() { return RelAmps; }
163 void SetSamples(int i) { nSamples = i; }
164 void SetLimit(int i) { Limit = i; }
165 void SetAmp1(double a) { amp1 = a; ParVals[0] = amp1; }
166 void SetAmp2(double a) { amp2 = a; ParVals[1] = amp2; }
167 void SetG1(double gfac) { g1 = gfac; ParVals[2] = g1; }
168 void SetG2(double gfac) { g2 = gfac; ParVals[3] = g2;}
169 void SetTau1(double tau) { tau1 = tau; ParVals[4] = tau1;}
170 void SetTau2(double tau) { tau2 = tau; ParVals[5] = tau2;}
171 void SetN1N2(double n) { n1n2 = n; ParVals[6] = n1n2;}
172 void SetMean(double m) { mean = m; ParVals[7] = mean;}
173 void SetFWHM(double w) { FWHM = w; ParVals[8] = FWHM;}
174 void SetOffset(double o) { offset = o; ParVals[9] = offset;}
175 void SetHist(TH1* hist) { Hist = hist; }
176 void SetField();
177
178 double GetG1() { return g1; }
179 double GetG1Err() { return ParErrs[2]; }
180 double GetG2() { return g2; }
181 double GetG2Err() { return ParErrs[3]; }
182 double GetAmp1() { return amp1; }
183 double GetAmp1Err() { return ParErrs[0]; }
184 double GetAmp2() { return amp2; }
185 double GetAmp2Err() { return ParErrs[1]; }
186 double GetTau1() { return tau1; }
187 double GetTau1Err() { return ParErrs[4]; }
188 double GetTau2() { return tau2; }
189 double GetTau2Err() { return ParErrs[5]; }
190 double GetN1N2() { return n1n2; }
191 double GetN1N2Err() { return ParErrs[6]; }
192 double GetMean() { return mean; }
193 double GetMeanErr() { return ParErrs[7]; }
194 double GetFWHM() { return FWHM; }
195 double GetFWHMErr() { return ParErrs[8]; }
196 double GetOffset() { return offset; }
197 double GetOffsetErr() { return ParErrs[9]; }
198
199 void FixField() { ParFix[7] = 1; ParFix[8] = 1; FixedField = 1;}
200 void UnFixField() { ParFix[7] = 0; ParFix[8] = 0; FixedField = 0;}
201 void FixG1() { ParFix[2] = 1; }
202 void UnFixG1() { ParFix[2] = 0; }
203 void FixG2() { ParFix[3] = 1; }
204 void UnFixG2() { ParFix[3] = 0; }
205 void FixAmp1() { ParFix[0] = 1; }
206 void UnFixAmp1() { ParFix[0] = 1; }
207 void FixAmp2() { ParFix[1] = 1; }
208 void UnFixAmp2() { ParFix[1] = 1; }
209 void FixTau1() { ParFix[4] = 1; }
210 void UnFixTau1() { ParFix[4] = 1; }
211 void FixTau2() { ParFix[5] = 1; }
212 void UnFixTau2() { ParFix[5] = 1; }
213 void FixN1N2() { ParFix[6] = 1; }
214 void UnFixN1N2() { ParFix[7] = 1; }
215 void FixOffset() { ParFix[9] = 1; }
216 void UnFixOffset() { ParFix[9] = 0; }
217 std::vector<double> GetPars() { std::vector<double> pars = { amp1, amp2, g1, g2, tau1, tau2, n1n2, mean, FWHM, offset }; return pars; }
218
219 void SetPars() { SetPars(ParVals); }
220 void SetPars(std::vector<double> pars) { ParVals = pars; amp1 = pars[0]; amp2 = pars[1]; g1 = pars[2]; g2 = pars[3]; tau1 = pars[4]; tau2 = pars[5]; n1n2 = pars[6]; mean = pars[7]; FWHM = pars[8]; offset = pars[9]; }
221
222 void SetParLimits(int i, double low, double high) { ParLimsLow[i] = low; ParLimsHigh[i] = high; }
223
224 double RatioFunc(double x, double a1, double a2, double gfac1, double gfac2, double t1, double t2, double n, double off);
225 double operator()(double *x, double *p);
226 TF1* GetFunc();
227
228 void Fit(int quiet);
229
230 void Print() {
231 std::cout << "g1 " << GamR::Utils::wrresult(GetG1(), GetG1Err()) << std::endl;
232 std::cout << "g2 " << GamR::Utils::wrresult(GetG2(), GetG2Err()) << std::endl;
233 std::cout << "amp1 " << GamR::Utils::wrresult(GetAmp1(), GetAmp1Err()) << std::endl;
234 std::cout << "amp2 " << GamR::Utils::wrresult(GetAmp2(), GetAmp2Err()) << std::endl;
235 std::cout << "tau1 " << GamR::Utils::wrresult(GetTau1(), GetTau1Err()) << std::endl;
236 std::cout << "tau2 " << GamR::Utils::wrresult(GetTau2(), GetTau2Err()) << std::endl;
237 std::cout << "n1/n2 " << GamR::Utils::wrresult(GetN1N2(), GetN1N2Err()) << std::endl;
238 std::cout << "Bhf " << GamR::Utils::wrresult(GetMean(), GetMeanErr()) << std::endl;
239 std::cout << "FWHM " << GamR::Utils::wrresult(GetFWHM(), GetFWHMErr()) << std::endl;
240 std::cout << "Offset " << GamR::Utils::wrresult(GetOffset(), GetOffsetErr()) << std::endl;
241 }
242
243 };
244
245 } // namespace AngDist
246} // namespace GamR
247
248#endif
void SetPars(std::vector< double > pars)
double operator()(double *x, double *p)
void SetParLimits(int i, double low, double high)
std::vector< double > GetBhfs()
std::vector< double > GetPars()
std::vector< double > GetRelAmps()
DoubleGaussianDistr(double L, double H, std::vector< double > pars)
double RatioFunc(double x, double a1, double a2, double gfac1, double gfac2, double t1, double t2, double n, double off)
void SetPars(std::vector< double > pars)
double RatioFunc(double x, double a, double gfac, double off)
std::vector< double > GetPars()
GaussianDistr(double L, double H, std::vector< double > pars)
std::vector< double > GetRelAmps()
double operator()(double *x, double *p)
std::vector< double > GetBhfs()
int wrresult(char *out, float value, float err, int minlen)
Definition Utilities.cc:241
Definition Gain.cc:19