GamR  0.0.0
GammaROOT
Loading...
Searching...
No Matches
Fit.cc
Go to the documentation of this file.
1#include <iostream>
2#include <string.h>
3
4/* ROOT */
5#include <TArrow.h>
6#include <TH1.h>
7#include <TMarker.h>
8#include <TMath.h>
9#include <TMatrixD.h>
10#include <TString.h>
11#include <TText.h>
12#include <TPolyMarker.h>
13#include <TRatioPlot.h>
14#include <TFitResult.h>
15
16#include "Fit.hh"
17#include "FitGuesses.hh"
18#include <utils/Utilities.hh>
19#include <toolkit/Peak.hh>
20
21namespace GamR {
22 namespace Spect {
23 //extern PeakFitGuesses *gFitGuesses;
24
28
29 GamR::Spect::PeakFit::PeakFit(double Low, double High) : fLow(Low), fHigh(High) {
30 parameters = {};
31 fFitGuesses = gFitGuesses;
32 }
33
34 GamR::Spect::PeakFit::PeakFit(double Low, double High, PeakFitGuesses *FitGuesses) : fLow(Low), fHigh(High) {
35 parameters = {};
36 fFitGuesses = FitGuesses;
37 }
38
39 void GamR::Spect::PeakFit::SetOpts(Option_t *option) {
40 TString opts(option);
41 opts.ToLower();
42
43 /* iWidthsEnergy sets the sigma according to an empirical energy formula */
44 //int iWidthsEnergy = 0;
45
46 if (opts.Contains("z")) {
47 parameters.iQuiet = 1;
48 }
49 if (opts.Contains("f")) {
50 parameters.iFixWidths = 1;
51 }
52 if (opts.Contains("ff")) {
53 parameters.iFixWidthsFile = 1;
54 }
55 if (opts.Contains("q")) {
56 parameters.iQuadBack = 1;
57 }
58 if (opts.Contains("c")) {
59 parameters.iConstantBack = 1;
60 }
61 if (opts.Contains("t")) {
62 parameters.iTails1 = 1;
63 }
64 if (opts.Contains("s")) {
65 parameters.iStep = 1;
66 }
67 if (opts.Contains("g")) {
68 parameters.iTails2 = 1;
69 }
70 if (opts.Contains("e")) {
71 parameters.iFixEnergy = 1;
72 }
73 if (opts.Contains("n")) {
74 parameters.iNoFit = 1;
75 }
76
77 if (!(parameters.iStep || parameters.iTails1 || parameters.iTails2)) {
78 fPeakType = GamR::TK::Gaussian;
79 }
80 else if (parameters.iStep && !parameters.iTails1 && !parameters.iTails2) {
81 fPeakType = GamR::TK::StepGaussian;
82 }
83 else if (!parameters.iStep && parameters.iTails1) {
84 fPeakType = GamR::TK::OneTailGaussian;
85 }
86 else if (!parameters.iStep && parameters.iTails2) {
87 fPeakType = GamR::TK::TwoTailGaussian;
88 }
89 else if (parameters.iStep && parameters.iTails1) {
90 fPeakType = GamR::TK::OneTailStepGaussian;
91 }
92 else if (parameters.iStep && parameters.iTails2) {
93 fPeakType = GamR::TK::TwoTailStepGaussian;
94 }
95
96 }
97
98 void GamR::Spect::PeakFit::SetBackground() {
99 iParamCount = 0;
100 int iPeak = 0;
101 TString sFuncString;
102 if (parameters.iQuadBack) {
103 sFuncString.Form("[0] + [1]*(x-%f) + [2]*pow(x-%f, 2)", (fHigh+fLow)/2, (fHigh+fLow)/2);
104 iParamCount = iParamCount + 3;
105 } else if (parameters.iConstantBack) {
106 sFuncString.Form("[0]");
107 iParamCount = iParamCount + 1;
108 }
109 else {
110 sFuncString.Form("[0] + [1]*(x-%f)", (fHigh+fLow)/2);
111 iParamCount = iParamCount + 2;
112 }
113
114
115 fBackground = new TF1("fBackground", sFuncString.Data(), fLow, fHigh);
116 }
117
118 void GamR::Spect::PeakFit::SetPeaks(std::vector<double> Peaks) {
119 //make a map
120 std::map<std::string,double> mPeaks;
121 //loop over Peaks and generate names
122 for (int iPeak; iPeak<Peaks.size(); ++iPeak) {
123 std::string name = "Peak" + std::to_string(iPeak);
124 mPeaks[name] = Peaks[iPeak];
125 }
126 GamR::Spect::PeakFit::SetPeaks(mPeaks);
127 }
128
129 void GamR::Spect::PeakFit::SetPeaks(std::map<std::string,double> Peaks) {
130 for (auto &mPeak : Peaks) {
131 auto peakKey = mPeak.first;
132 auto peak = mPeak.second;
133 AddPeak(peakKey, fPeakType);
134 }
135 }
136
137 void GamR::Spect::PeakFit::SetIndices(std::map<std::string,GamR::TK::FPeak*> Peaks) {
138 std::string firstKey;
139 int i = 0;
140 for (auto &mPeak : Peaks) {
141 auto peakKey = mPeak.first;
142 auto peak = mPeak.second;
143
144 fPeakParamInds[peakKey]={};
145
146 int iAmpParam = iParamCount;
147 iParamCount++;
148 fPeakParamInds[peakKey].iAmplitude = iAmpParam;
149 int iCentParam = iParamCount;
150 iParamCount++;
151 fPeakParamInds[peakKey].iCentroid = iCentParam;
152 int iWidthParam;
153 int iSkewAmpParam;
154 int iSkewWidthParam;
155 int iStepAmpParam;
156 int iSkewAmpParam2;
157 int iSkewWidthParam2;
158 if (parameters.iFixWidthsFile) {
159 iWidthParam = -1;
160 }
161 else if (parameters.iFixWidths) {
162 if (i == 0) {
163 iWidthParam = iParamCount;
164 iParamCount++;
165 firstKey = peakKey;
166 } else {
167 iWidthParam = fPeakParamInds[firstKey].iWidth;
168 }
169 } else {
170 iWidthParam = iParamCount;
171 iParamCount++;
172 }
173 fPeakParamInds[peakKey].iWidth = iWidthParam;
174 if (fPeakType == GamR::TK::StepGaussian ||
175 fPeakType == GamR::TK::OneTailStepGaussian ||
176 fPeakType == GamR::TK::TwoTailStepGaussian) {
177 if (parameters.iFixWidthsFile) {
178 iStepAmpParam = -1;
179 }
180 else if (parameters.iFixWidths) {
181 if (i == 0) {
182 iStepAmpParam = iParamCount;
183 iParamCount++;
184 firstKey = peakKey;
185 } else {
186 iStepAmpParam = fPeakParamInds[firstKey].iStepAmplitude;
187 }
188 } else {
189 iStepAmpParam = iParamCount;
190 iParamCount++;
191 }
192 fPeakParamInds[peakKey].iStepAmplitude = iStepAmpParam;
193 }
194
195 if (fPeakType == GamR::TK::OneTailGaussian ||
196 fPeakType == GamR::TK::OneTailStepGaussian ||
197 fPeakType == GamR::TK::TwoTailGaussian ||
198 fPeakType == GamR::TK::TwoTailStepGaussian) {
199
200 if (parameters.iFixWidthsFile) {
201 iSkewAmpParam = -1;
202 iSkewWidthParam = -1;
203 }
204 else if (parameters.iFixWidths ) {
205 if (i == 0) {
206 iSkewAmpParam = iParamCount;
207 iParamCount++;
208 iSkewWidthParam = iParamCount;
209 iParamCount++;
210 } else {
211 iSkewAmpParam = fPeakParamInds[firstKey].iSkewAmplitude;
212 iSkewWidthParam = fPeakParamInds[firstKey].iSkewWidth;
213 }
214 } else {
215 iSkewAmpParam = iParamCount;
216 iParamCount++;
217 iSkewWidthParam = iParamCount;
218 iParamCount++;
219 }
220 fPeakParamInds[peakKey].iSkewAmplitude = iSkewAmpParam;
221 fPeakParamInds[peakKey].iSkewWidth = iSkewWidthParam;
222
223 if (fPeakType == GamR::TK::TwoTailGaussian ||
224 fPeakType == GamR::TK::TwoTailStepGaussian) {
225
226 if (parameters.iFixWidthsFile) {
227 iSkewAmpParam2 = -1;
228 iSkewWidthParam2 = -1;
229 }
230 else if (parameters.iFixWidths) {
231 if (i == 0) {
232 iSkewAmpParam2 = iParamCount;
233 iParamCount++;
234 iSkewWidthParam2 = iParamCount;
235 iParamCount++;
236 } else {
237 iSkewAmpParam2 = fPeakParamInds[firstKey].iSkewAmplitude2;
238 iSkewWidthParam2 = fPeakParamInds[firstKey].iSkewWidth2;
239 }
240 } else {
241 iSkewAmpParam2 = iParamCount;
242 iParamCount++;
243 iSkewWidthParam2 = iParamCount;
244 iParamCount++;
245 }
246 fPeakParamInds[peakKey].iSkewAmplitude2 = iSkewAmpParam2;
247 fPeakParamInds[peakKey].iSkewWidth2 = iSkewWidthParam2;
248 }
249 }
250
251 ++i;
252 }
253 /* this creates the total function from all the peaks + background */
254 SetFunction();
255 fTotal->SetNpx(500);
256
257 /* set parameter names */
258 fTotal->SetParName(0, "c");
259 if (!parameters.iConstantBack) {
260 fTotal->SetParName(1, "b");
261 }
262 if (parameters.iQuadBack) {
263 fTotal->SetParName(2, "a");
264 }
265
266 i=0;
267 for (auto &mpp : fPeakParamInds) {
268 auto peakKey = mpp.first;
269 auto pp = mpp.second;
270 fTotal->SetParName(pp.iAmplitude, ("Amp" + std::to_string(i)).c_str());
271 fTotal->SetParName(pp.iCentroid, ("Cent" + std::to_string(i)).c_str());
272 if (parameters.iFixWidths && !parameters.iFixWidthsFile) {
273 fTotal->SetParName(pp.iWidth, "Sigma");
274 } else {
275 fTotal->SetParName(pp.iWidth, ("Sigma" + std::to_string(i)).c_str());
276 }
277 if (fPeakType == GamR::TK::StepGaussian ||
278 fPeakType == GamR::TK::OneTailStepGaussian ||
279 fPeakType == GamR::TK::TwoTailStepGaussian) {
280 if (parameters.iFixWidths && !parameters.iFixWidthsFile) {
281 fTotal->SetParName(pp.iStepAmplitude, "H");
282 }
283 else {
284 fTotal->SetParName(pp.iStepAmplitude, ("H" + std::to_string(i)).c_str());
285 }
286 }
287
288 if (fPeakType == GamR::TK::OneTailGaussian ||
289 fPeakType == GamR::TK::OneTailStepGaussian ||
290 fPeakType == GamR::TK::TwoTailGaussian ||
291 fPeakType == GamR::TK::TwoTailStepGaussian) {
292
293 if (parameters.iFixWidths && !parameters.iFixWidthsFile) {
294 fTotal->SetParName(pp.iSkewAmplitude, "R");
295 fTotal->SetParName(pp.iSkewWidth, "Beta");
296 } else {
297 fTotal->SetParName(pp.iSkewAmplitude, ("R" + std::to_string(i)).c_str());
298 fTotal->SetParName(pp.iSkewWidth, ("Beta" + std::to_string(i)).c_str());
299 }
300
301 if (fPeakType == GamR::TK::TwoTailGaussian ||
302 fPeakType == GamR::TK::TwoTailStepGaussian) {
303 if (parameters.iFixWidths && !parameters.iFixWidthsFile) {
304 fTotal->SetParName(pp.iSkewAmplitude2, "R2");
305 fTotal->SetParName(pp.iSkewWidth2, "Beta2");
306 } else {
307 fTotal->SetParName(pp.iSkewAmplitude2, ("R2_" + std::to_string(i)).c_str());
308 fTotal->SetParName(pp.iSkewWidth2, ("Beta2_" + std::to_string(i)).c_str());
309 }
310 }
311 }
312 ++i;
313 }
314 }
315
316 void GamR::Spect::PeakFit::SetGuesses(TH1 *hist, std::map<std::string,GamR::TK::FPeak*> mPeaks, std::vector<std::string> FixPeaks /*= {}*/) {
317 std::map<std::string,double> mcentroids;
318 for (auto &mpeak : mPeaks) {
319 auto peakKey = mpeak.first;
320 auto peak = mpeak.second;
321 mcentroids[peakKey]=peak->GetCent();
322 }
323 GamR::Spect::PeakFit::SetGuesses(hist,mcentroids,FixPeaks);
324 }
325
326 void GamR::Spect::PeakFit::SetGuesses(TH1 *hist, std::map<std::string, double> Peaks, std::vector<std::string> FixPeaks /*= {}*/) {
327 std::cout << "Inside SetGuesses" << std::endl;
328 std::cout << "checking FitGuesses" << std::endl;
329 fFitGuesses->Print();
330 for (auto &peakKey : FixPeaks) {
331 FixPeakEnergy(peakKey);
332 }
333
334 /* set initial guesses */
335 fNData = hist->FindBin(fHigh) - hist->FindBin(fLow) + 1;
336 double slope = (hist->GetBinContent(hist->FindBin(fHigh)) - hist->GetBinContent(hist->FindBin(fLow))) / (fHigh - fLow);
337 if (!parameters.iConstantBack) {
338 fTotal->SetParameter(1, slope);
339 }
340 double offset = (hist->GetBinContent(hist->FindBin(fLow)) +
341 hist->GetBinContent(hist->FindBin(fHigh)))/2.0;
342 fTotal->SetParameter(0, offset);
343 if (parameters.iQuadBack) {
344 fTotal->SetParameter(2, 0.01);
345 }
346 //for (int i = 0; i < (int)Peaks.size(); i++) {
347 double histMax = hist->GetMaximum();
348 for (auto &peak : Peaks) {
349 auto peakKey = peak.first;
350 auto cent = peak.second;
351
352 //std::cout << peakKey << " " << cent << std::endl;
353 //int nBin = hist->FindBin(peaks[i]);
354 int nBin = hist->FindBin(cent);
355 //auto pp = fPeakParamInds[i];
356 auto pp = fPeakParamInds[peakKey];
357 fTotal->SetParameter(pp.iAmplitude, hist->GetBinContent(nBin) - (offset + slope*(cent - (fHigh + fLow)/2.0)));
358 //make sure peaks are not negative
359 fTotal->SetParLimits(pp.iAmplitude, 0, histMax*2.0);
360 //fTotal->SetParameter(pp.iCentroid, Peaks[i]);
361 fTotal->SetParameter(pp.iCentroid, cent);
362 //std::cout << pp.iCentroid << " " << cent << std::endl;
363 if (fPeakEnFix.find(peakKey) != fPeakEnFix.end()) {
364 if ( fPeakEnFix[peakKey] == 1 ) {
365 fTotal->FixParameter(pp.iCentroid, cent);
366 }
367 }
368
369
370 //sigma = 5*bin width (?)
371 //fTotal->SetParameter(pp.iWidth, 5.0*hist->GetBinWidth(i));
372 if (!parameters.iFixWidthsFile) {
373 fTotal->SetParameter(pp.iWidth, fFitGuesses->fWidth.val);
374 fTotal->SetParLimits(pp.iWidth, fFitGuesses->fWidth.low, fFitGuesses->fWidth.high);
375 }
376
377 if (fPeakType == GamR::TK::StepGaussian ||
378 fPeakType == GamR::TK::OneTailStepGaussian ||
379 fPeakType == GamR::TK::TwoTailStepGaussian) {
380 /*
381 double low;
382 double high;
383 fTotal->GetParLimits(pp.iStepAmplitude,low,high);
384 if (low==high) {
385 fTotal->FixParameter(pp.iStepAmplitude, 0.);
386 } else {
387 fTotal->SetParameter(pp.iStepAmplitude, (high-low)/2.);
388 }
389 */
390 if (!parameters.iFixWidthsFile) {
391 fTotal->FixParameter(pp.iStepAmplitude, 0.);
392 }
393 }
394
395 if (fPeakType == GamR::TK::OneTailGaussian ||
396 fPeakType == GamR::TK::OneTailStepGaussian ||
397 fPeakType == GamR::TK::TwoTailGaussian ||
398 fPeakType == GamR::TK::TwoTailStepGaussian) {
399
400 /*
401 fTotal->SetParameter(skewWidthParameters[i], 1.);
402 fTotal->SetParLimits(skewWidthParameters[i], 0, 10);
403 fTotal->SetParameter(skewAmplitudeParameters[i], 10.);
404 fTotal->SetParLimits(skewAmplitudeParameters[i], 0, 50);
405 fTotal->SetParameter(stepAmplitudeParameters[i], 1.);
406 fTotal->SetParLimits(stepAmplitudeParameters[i], 0, 20);
407 */
408 if (!parameters.iFixWidthsFile) {
409 fTotal->FixParameter(pp.iSkewWidth, 2.5);
410 fTotal->FixParameter(pp.iSkewAmplitude, 0.0);
411 }
412 if (fPeakType == GamR::TK::TwoTailGaussian ||
413 fPeakType == GamR::TK::TwoTailStepGaussian) {
414
415 if (!parameters.iFixWidthsFile) {
416 fTotal->FixParameter(pp.iSkewWidth2, 5);
417 fTotal->FixParameter(pp.iSkewAmplitude2, 0.0);
418 }
419 }
420 }
421 }
422 }
423
424 void GamR::Spect::PeakFit::SetLimit(std::string peakKey, std::string parameterKey, double lowerLimit, double upperLimit) {
425 auto pp = fPeakParamInds[peakKey];
426 auto peak = fPeaks[peakKey];
427 if (!peak) std::cout << "Error: Peak '" << peakKey << "' not found." << std::endl;
428 auto func = peak->GetFunc();
429 func->SetParLimits(func->GetParNumber((TString)parameterKey),lowerLimit,upperLimit);
430
431 }
432
433 void GamR::Spect::PeakFit::Fit(TH1 *hist, Option_t *foption) {
434 //if (!parameters.iQuiet) {
435 // fFitGuesses->Print();
436 //}
437 std::string fopts(foption);
438 fopts = fopts + " R";
439
440 fFitGuesses->fScale.val = 1./(hist->GetBinCenter(2)-hist->GetBinCenter(1));
441
442 //
443 if (parameters.iFixWidthsFile) {
444 FILE *file = fopen("FitWidths.dat", "ra");
445 if (!file) { std::cerr << "Error! file FitWidths.dat not opened" << std::endl; return;}
446 std::stringstream ss;
447 char cline[2048];
448
449 gWidth = new TGraph();
450 gStep = new TGraph();
451 gSkew = new TGraph();
452 gSkewAmp = new TGraph();
453 gSkew2 = new TGraph();
454 gSkewAmp2 = new TGraph();
455
456 while(std::fgets(cline, sizeof cline, file)!=NULL) {
457 std::string line(cline);
458 if (line.size() == 0) { continue; }
459 if (line[0] == '#') { continue; }
460 if (line[0] == ';') { continue; }
461 if (line[0] == '!') { continue; }
462
463 ss.clear();
464 ss.str(line);
465
466 double energy, width, step, skew, skewamp, skew2, skewamp2;
467 ss >> energy >> width >> step >> skew >> skewamp >> skew2 >> skewamp2;
468
469 gWidth->AddPoint(energy, width);
470 gStep->AddPoint(energy,step);
471 gSkew->AddPoint(energy,skew);
472 gSkewAmp->AddPoint(energy,skewamp);
473 gSkew2->AddPoint(energy,skew2);
474 gSkewAmp2->AddPoint(energy,skewamp2);
475 }
476 }
477
478 // first fit, only gaussian
479 // this seems to be broken (sometimes?) with tails when tail amplitude fixed at zero
480 if (fPeakType == GamR::TK::Gaussian) {
481 if (!parameters.iNoFit) {
482 hist->Fit(fTotal, fopts.c_str());
483 }
484 }
485
486 // with step now
487 if (fPeakType == GamR::TK::StepGaussian ||
488 fPeakType == GamR::TK::OneTailStepGaussian ||
489 fPeakType == GamR::TK::TwoTailStepGaussian) {
490 for (auto &mpp : fPeakParamInds) {
491 auto peakKey = mpp.first;
492 auto pp = mpp.second;
493 if (!parameters.iFixWidthsFile) {
494 fTotal->SetParameter(pp.iStepAmplitude, fFitGuesses->fStepAmp.val);
495 fTotal->SetParLimits(pp.iStepAmplitude, fFitGuesses->fStepAmp.low, fFitGuesses->fStepAmp.high);
496 }
497 /*
498 double low;
499 double high;
500 fTotal->GetParLimits(pp.iStepAmplitude,low,high);
501 fTotal->SetParameter(pp.iStepAmplitude, 1.);
502 if (low==high) {
503 fTotal->SetParLimits(pp.iStepAmplitude, 0, 20.);
504 }
505 */
506 }
507 if (!parameters.iNoFit) {
508 hist->Fit(fTotal, fopts.c_str());
509 }
510 }
511
512 // 2 more fits for tails, iterative
513 if (fPeakType == GamR::TK::OneTailGaussian ||
514 fPeakType == GamR::TK::OneTailStepGaussian ||
515 fPeakType == GamR::TK::TwoTailGaussian ||
516 fPeakType == GamR::TK::TwoTailStepGaussian) {
517 for (auto &mpp : fPeakParamInds) {
518 auto peakKey = mpp.first;
519 auto pp = mpp.second;
520 if (!parameters.iFixWidthsFile) {
521 fTotal->SetParameter(pp.iSkewWidth, fFitGuesses->fSkewWidth.val);
522 fTotal->SetParLimits(pp.iSkewWidth, fFitGuesses->fSkewWidth.low, fFitGuesses->fSkewWidth.high);
523 fTotal->SetParameter(pp.iSkewAmplitude, fFitGuesses->fSkewAmp.val);
524 fTotal->SetParLimits(pp.iSkewAmplitude, fFitGuesses->fSkewAmp.low, fFitGuesses->fSkewAmp.high);
525 }
526 }
527 if (!parameters.iNoFit) {
528 hist->Fit(fTotal, fopts.c_str());
529 }
530 }
531
532 if (fPeakType == GamR::TK::TwoTailGaussian ||
533 fPeakType == GamR::TK::TwoTailStepGaussian) {
534 for (auto &mpp : fPeakParamInds) {
535 auto peakKey = mpp.first;
536 auto pp = mpp.second;
537 if (!parameters.iFixWidthsFile) {
538 fTotal->SetParameter(pp.iSkewWidth2, fFitGuesses->fSkewWidth2.val);
539 fTotal->SetParLimits(pp.iSkewWidth2, fFitGuesses->fSkewWidth2.low, fFitGuesses->fSkewWidth2.high);
540 fTotal->SetParameter(pp.iSkewAmplitude2, fFitGuesses->fSkewAmp2.val);
541 fTotal->SetParLimits(pp.iSkewAmplitude2, fFitGuesses->fSkewAmp2.low, fFitGuesses->fSkewAmp2.high);
542 }
543 }
544 if (!parameters.iNoFit) {
545 hist->Fit(fTotal, fopts.c_str());
546 }
547 }
548
549 }
550
551 void GamR::Spect::PeakFit::SetResults(TH1* hist) {
552 /* set results in all the components of the results object from fTotal*/
553
554 GetBackground()->SetParameter(0, fTotal->GetParameter(0));
555 if (!parameters.iConstantBack) {
556 GetBackground()->SetParameter(1, fTotal->GetParameter(1));
557 }
558 if (parameters.iQuadBack) {
559 GetBackground()->SetParameter(2, fTotal->GetParameter(2));
560 }
561 SetChi2(hist->Chisquare(fTotal, "R"));
562 for (auto &mpeak : fPeaks) {
563 auto peakKey = mpeak.first;
564 auto peak = mpeak.second;
565 auto pp = fPeakParamInds[peakKey];
566 switch(fPeakType) {
567 case GamR::TK::Gaussian : {
568 if (parameters.iFixWidthsFile) {
569 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
570 fTotal->GetParameter(pp.iCentroid),
571 gWidth->Eval(fTotal->GetParameter(pp.iCentroid))/2.3548});
572 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
573 fTotal->GetParError(pp.iCentroid),
574 0});
575 }
576 else {
577 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
578 fTotal->GetParameter(pp.iCentroid),
579 fTotal->GetParameter(pp.iWidth)});
580 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
581 fTotal->GetParError(pp.iCentroid),
582 0.0});
583 }
584 break;
585 }
586 case GamR::TK::StepGaussian : {
587 if (parameters.iFixWidthsFile) {
588 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
589 fTotal->GetParameter(pp.iCentroid),
590 gWidth->Eval(fTotal->GetParameter(pp.iCentroid))/2.3548,
591 gStep->Eval(fTotal->GetParameter(pp.iCentroid))});
592 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
593 fTotal->GetParError(pp.iCentroid),
594 0.0,
595 0.0});
596 }
597 else {
598 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
599 fTotal->GetParameter(pp.iCentroid),
600 fTotal->GetParameter(pp.iWidth),
601 fTotal->GetParameter(pp.iStepAmplitude)});
602 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
603 fTotal->GetParError(pp.iCentroid),
604 fTotal->GetParError(pp.iWidth),
605 fTotal->GetParError(pp.iStepAmplitude)});
606 }
607 break;
608 }
609 case GamR::TK::OneTailGaussian : {
610 if (parameters.iFixWidthsFile) {
611 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
612 fTotal->GetParameter(pp.iCentroid),
613 gWidth->Eval(fTotal->GetParameter(pp.iCentroid))/2.3548,
614 gSkewAmp->Eval(fTotal->GetParameter(pp.iCentroid)),
615 gSkew->Eval(fTotal->GetParameter(pp.iCentroid))});
616 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
617 fTotal->GetParError(pp.iCentroid),
618 0.0,
619 0.0,
620 0.0});
621 }
622 else {
623 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
624 fTotal->GetParameter(pp.iCentroid),
625 fTotal->GetParameter(pp.iWidth),
626 fTotal->GetParameter(pp.iSkewAmplitude),
627 fTotal->GetParameter(pp.iSkewWidth)});
628 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
629 fTotal->GetParError(pp.iCentroid),
630 fTotal->GetParError(pp.iWidth),
631 fTotal->GetParError(pp.iSkewAmplitude),
632 fTotal->GetParError(pp.iSkewWidth)});
633 }
634 break;
635 }
636 case GamR::TK::TwoTailGaussian : {
637 if (parameters.iFixWidthsFile) {
638 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
639 fTotal->GetParameter(pp.iCentroid),
640 gWidth->Eval(fTotal->GetParameter(pp.iCentroid))/2.3548,
641 gSkewAmp->Eval(fTotal->GetParameter(pp.iCentroid)),
642 gSkew->Eval(fTotal->GetParameter(pp.iCentroid)),
643 gSkewAmp2->Eval(fTotal->GetParameter(pp.iCentroid)),
644 gSkew2->Eval(fTotal->GetParameter(pp.iCentroid))});
645 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
646 fTotal->GetParError(pp.iCentroid),
647 0.0,
648 0.0,
649 0.0,
650 0.0,
651 0.0});
652 }
653 else {
654 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
655 fTotal->GetParameter(pp.iCentroid),
656 fTotal->GetParameter(pp.iWidth),
657 fTotal->GetParameter(pp.iSkewAmplitude),
658 fTotal->GetParameter(pp.iSkewWidth),
659 fTotal->GetParameter(pp.iSkewAmplitude2),
660 fTotal->GetParameter(pp.iSkewWidth2)});
661 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
662 fTotal->GetParError(pp.iCentroid),
663 fTotal->GetParError(pp.iWidth),
664 fTotal->GetParError(pp.iSkewAmplitude),
665 fTotal->GetParError(pp.iSkewWidth),
666 fTotal->GetParError(pp.iSkewAmplitude2),
667 fTotal->GetParError(pp.iSkewWidth2)});
668 }
669 break;
670 }
671 case GamR::TK::OneTailStepGaussian : {
672 if (parameters.iFixWidthsFile) {
673 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
674 fTotal->GetParameter(pp.iCentroid),
675 gWidth->Eval(fTotal->GetParameter(pp.iCentroid))/2.3548,
676 gSkewAmp->Eval(fTotal->GetParameter(pp.iCentroid)),
677 gSkew->Eval(fTotal->GetParameter(pp.iCentroid)),
678 gStep->Eval(fTotal->GetParameter(pp.iCentroid))});
679 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
680 fTotal->GetParError(pp.iCentroid),
681 0.0,
682 0.0,
683 0.0,
684 0.0});
685
686 }
687 else {
688 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
689 fTotal->GetParameter(pp.iCentroid),
690 fTotal->GetParameter(pp.iWidth),
691 fTotal->GetParameter(pp.iSkewAmplitude),
692 fTotal->GetParameter(pp.iSkewWidth),
693 fTotal->GetParameter(pp.iStepAmplitude)});
694 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
695 fTotal->GetParError(pp.iCentroid),
696 fTotal->GetParError(pp.iWidth),
697 fTotal->GetParError(pp.iSkewAmplitude),
698 fTotal->GetParError(pp.iSkewWidth),
699 fTotal->GetParError(pp.iStepAmplitude)});
700 }
701 break;
702 }
703 case GamR::TK::TwoTailStepGaussian : {
704 if (parameters.iFixWidthsFile) {
705 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
706 fTotal->GetParameter(pp.iCentroid),
707 gWidth->Eval(fTotal->GetParameter(pp.iCentroid))/2.3548,
708 gSkewAmp->Eval(fTotal->GetParameter(pp.iCentroid)),
709 gSkew->Eval(fTotal->GetParameter(pp.iCentroid)),
710 gSkewAmp2->Eval(fTotal->GetParameter(pp.iCentroid)),
711 gSkew2->Eval(fTotal->GetParameter(pp.iCentroid)),
712 gStep->Eval(fTotal->GetParameter(pp.iCentroid))});
713 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
714 fTotal->GetParError(pp.iCentroid),
715 0.0,
716 0.0,
717 0.0,
718 0.0,
719 0.0,
720 0.0});
721 }
722 else {
723 peak->SetParameters({fTotal->GetParameter(pp.iAmplitude),
724 fTotal->GetParameter(pp.iCentroid),
725 fTotal->GetParameter(pp.iWidth),
726 fTotal->GetParameter(pp.iSkewAmplitude),
727 fTotal->GetParameter(pp.iSkewWidth),
728 fTotal->GetParameter(pp.iSkewAmplitude2),
729 fTotal->GetParameter(pp.iSkewWidth2),
730 fTotal->GetParameter(pp.iStepAmplitude)});
731 peak->SetParErrors({fTotal->GetParError(pp.iAmplitude),
732 fTotal->GetParError(pp.iCentroid),
733 fTotal->GetParError(pp.iWidth),
734 fTotal->GetParError(pp.iSkewAmplitude),
735 fTotal->GetParError(pp.iSkewWidth),
736 fTotal->GetParError(pp.iSkewAmplitude2),
737 fTotal->GetParError(pp.iSkewWidth2),
738 fTotal->GetParError(pp.iStepAmplitude)});
739 }
740 break;
741 }
742 }
743 }
744 }
745
746 void GamR::Spect::PeakFit::PrintResults() {
747 std::cout << "Chi2 = " << fChi2 << ", ndata = " << fNData << ", NDF = " << fNData - fNPars << ", Chi2/NDF = " << fChi2/(double)(fNData - fNPars) << std::endl;
748 TString sPrintString;
749 if (!parameters.iQuiet) {
750 sPrintString.Form("Peak Centroid Height FWHM Area Counts ");
751 std::printf("%s\n", sPrintString.Data());
752 }
753
754 int i=0;
755 for (auto &mPeak : fPeaks) {
756 auto peakKey = mPeak.first;
757 auto peak = mPeak.second;
758 char line[120];
759 *line = '\0';
760 GamR::Utils::wrresult(line, peak->GetCent(),
761 peak->GetCentError(), 15);
762 GamR::Utils::wrresult(line, peak->GetAmp(),
763 peak->GetAmpError(), 30);
764 GamR::Utils::wrresult(line, peak->GetWidth(),
765 peak->GetWidthError(), 45);
766 GamR::Utils::wrresult(line, peak->GetArea(),
767 peak->GetAreaError(), 60);
768 GamR::Utils::wrresult(line, peak->GetArea()*fFitGuesses->fScale.val,
769 peak->GetAreaError()*fFitGuesses->fScale.val, 60);
770 if (!parameters.iQuiet) {
771 std::printf("%-5d %s\n", i, line);
772 }
773 ++i;
774 }
775 }
776
787 GamR::Spect::PeakFit *FitPeaks(TH1 *hist, double Low, double High, std::vector<double> Peaks, Option_t *foption /*=""*/,
788 Option_t *option /*=""*/, std::vector<std::string> FixPeaks) {
789 std::map<std::string, double> mPeaks;
790 for (int iPeak=0; iPeak<Peaks.size(); ++iPeak) {
791 std::string peakKey;
792 peakKey = "Peak" + std::to_string(iPeak);
793 mPeaks[peakKey] = Peaks[iPeak];
794 }
795 return GamR::Spect::FitPeaks(hist, Low, High, mPeaks, foption, option, FixPeaks);
796 }
797
798 GamR::Spect::PeakFit *FitPeaks(TH1 *hist, double Low, double High, std::map<std::string,double> Peaks, Option_t *foption /*=""*/,
799 Option_t *option /*=""*/, std::vector<std::string> FixPeaks) {
800
801 GamR::Spect::PeakFit *retval = new GamR::Spect::PeakFit(Low, High);
802
803 std::cout << retval->GetLow() << " " << retval->GetHigh() << std::endl;
804 //setup the fit: make the relevant GamR::TK::FPeak objects, set the background, make the function
805 //to fit, and set initial guesses in a sensible way
806
807 retval->SetUp(hist, Peaks, option, FixPeaks);
808 //execute the fit to the histogram
809
810 retval->Fit(hist, foption);
811
812 //draw total function after fitting
813 retval->GetTotal()->Draw("same");
814
815 //set the results of the fitted total function to all the sub-functions:
816 //background, each peak's individual parameters
817 retval->SetResults(hist);
818
819 //print results
820 retval->PrintResults();
821
822 return retval;
823 }
824
836 GamR::Spect::PeakFit *FitPeaks(TCanvas *canvas /*=gPad->GetCanvas()*/, Option_t *foption /*=""*/, Option_t *option /*=""*/)
837 {
838 TH1D *hist = GamR::Utils::GetHist1D(canvas);
839 if (!hist) { return NULL; }
840
841 TString opts(option);
842 opts.ToLower();
843
844 int iQuiet = 0;
845
846 if (opts.Contains("z")) {
847 iQuiet = 1;
848 }
849
850 /* get clicks for background and peaks */
851
852 TMarker *marker;
853 TObject *obj;
854
855 int peaks = -2;
856
857 double low = 0;
858 double high = 0;
859 std::vector<double> centroids;
860 std::vector<TArrow *> arrows;
861 std::vector<TText *> labels;
862 std::vector<TLine *> bounds;
863
864 canvas->SetCrosshair();
865
866 canvas->Update();
867
868 std::string canvasname = canvas->GetName();
869
870 GamR::Utils::Clicker click;
871 canvas->Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", "GamR::Utils::Clicker", &click, "GetClick(Int_t,Int_t,Int_t,TObject*)");
872
873 while (true) {
874 if (peaks == -2) {
875 if (!iQuiet){
876 std::cout << "Click for lower bound of region..." << std::endl;
877 }
878 } else if (peaks == -1) {
879 if (!iQuiet){
880 std::cout << "Click for upper bound of region..." << std::endl;
881 }
882 } else {
883 if (!iQuiet){
884 std::cout << "Click for peak centroid, press any key to fit..." << std::endl;
885 }
886 }
887 click.waiting=false;
888 obj = canvas->WaitPrimitive();
889 if (!obj)
890 break;
891 if (strncmp(obj->ClassName(), "TMarker", 7) == 0) {
892 marker = (TMarker *)obj;
893 if (peaks == -2) {
894 low = marker->GetX();
895 if (!iQuiet){
896 std::cout << "Lower bound: " << low << std::endl;
897 }
898 double low_y = canvas->GetUymin();
899 if (canvas->GetLogy()) {
900 low_y = std::pow(10, low_y);
901 }
902 TLine *line = new TLine(marker->GetX(), hist->GetBinContent(hist->FindBin(marker->GetX())), marker->GetX(),
903 low_y );
904 line->Draw();
905 bounds.push_back(line);
906 } else if (peaks == -1) {
907 high = marker->GetX();
908 if (!iQuiet){
909 std::cout << "Upper bound: " << high << std::endl;
910 }
911 double low_y = canvas->GetUymin();
912 if (canvas->GetLogy()) {
913 low_y = std::pow(10, low_y);
914 }
915 TLine *line = new TLine(marker->GetX(), hist->GetBinContent(hist->FindBin(marker->GetX())), marker->GetX(),
916 low_y);
917 line->Draw();
918 bounds.push_back(line);
919 } else {
920 centroids.push_back(marker->GetX());
921 if (!iQuiet){
922 std::cout << "Centroid: " << marker->GetX() << std::endl;
923 }
924 TArrow *arrow = new TArrow(marker->GetX(), 0.8 * hist->GetBinContent(hist->FindBin(marker->GetX())),
925 marker->GetX(), 0.95 * hist->GetBinContent(hist->FindBin(marker->GetX())), 0.01);
926 TString s;
927 s.Form("%d", peaks);
928 TText *text =
929 new TText(marker->GetX(), 0.75 * hist->GetBinContent(hist->FindBin(marker->GetX())), s.Data());
930 text->SetTextSize(0.03);
931 text->SetTextAlign(21);
932
933 arrow->Draw();
934 text->Draw();
935
936 arrows.push_back(arrow);
937 labels.push_back(text);
938 }
939 peaks = peaks + 1;
940 delete marker;
941 }
942 }
943
944 canvas->Disconnect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", &click, "GetClick(Int_t,Int_t,Int_t,TObject*)");
945
946 if (peaks < 1) {
947 if (!iQuiet){
948 std::cout << "defined no peaks!" << std::endl;
949 }
950 canvas->SetCrosshair(0);
951 hist->Draw("hist");
952 return NULL;
953 }
954
955 std::vector<std::string> FixEnergies;
956 if (opts.Contains("e")) {
957 std::cout << "Fixing peak energies: " << std::endl;
958 double energy;
959 int index = 0;
960 for (auto &cent : centroids) {
961 energy = 0;
962 std::cout << "Peak " << index << " [ " << cent << " ] : ";
963 std::string line;
964 std::getline(std::cin, line);
965
966 if (line.empty()) { index += 1; continue; }
967 std::stringstream ss;
968 ss << line;
969 ss >> energy;
970
971 if (energy != 0) {
972 centroids[index] = energy;
973 FixEnergies.push_back("Peak"+std::to_string(index));
974 }
975 index += 1;
976 }
977 }
978
979 GamR::Spect::PeakFit *retval = FitPeaks(hist, low, high, centroids, foption, option, FixEnergies);
980
981 TF1 *func = retval->GetTotal();
982
983 // draw background, individual peaks, and arrows
984 TF1 *background = retval->GetBackground();
985
986 background->Draw("same");
987
988 if ((int)centroids.size() >= 1) {
989 int i=0;
990 //for (int i = 0; i < (int)centroids.size(); i++) {
991 for (auto &mpeak : retval->GetPeaks()) {
992 auto peakKey = mpeak.first;
993 auto peak = mpeak.second;
994 //GamR::TK::FPeak *peak = retval->Get(i);
995 //TF1 *single = retval->GetPeakBG(i);
996 TF1 *single = retval->GetPeakBG(peakKey);
997 single->SetNpx(500);
998 single->SetLineColor(kBlue+1);
999 single->SetLineWidth(2);
1000 single->Draw("same");
1001
1002 TF1 *step = retval->GetPeakStepBG(peakKey);
1003 if (step != NULL) {
1004 step->SetNpx(500);
1005 step->SetLineColor(kBlack);
1006 step->SetLineWidth(1);
1007 step->Draw("same");
1008 }
1009
1010 if ((retval->GetPeakType() != GamR::TK::Gaussian) &&
1011 (retval->GetPeakType() != GamR::TK::StepGaussian)) {
1012 TF1 *gauss = retval->GetPeakGaussBG(peakKey);
1013 if (gauss != NULL) {
1014 gauss->SetNpx(500);
1015 gauss->SetLineWidth(1);
1016 gauss->SetLineColor(kGreen+2);
1017 gauss->Draw("same");
1018 }
1019 }
1020
1021 TArrow *arrow = arrows[i];
1022 canvas->GetListOfPrimitives()->Remove(arrow);
1023 arrow->SetX1(peak->GetCent());
1024 arrow->SetY1(0.8 *
1025 single->Eval(peak->GetCent()));
1026 arrow->SetX2(peak->GetCent());
1027 arrow->SetY2(0.95 *
1028 single->Eval(peak->GetCent()));
1029
1030 TString s;
1031 s.Form("%d", i);
1032 TText *text = labels[i];
1033 canvas->GetListOfPrimitives()->Remove(text);
1034 text->SetX(peak->GetCent());
1035 text->SetY(0.75 *
1036 single->Eval(peak->GetCent()));
1037 text->SetTextSize(0.03);
1038 text->SetTextAlign(21);
1039
1040 canvas->Modified();
1041 canvas->Update();
1042 arrow->Draw();
1043 text->Draw();
1044 canvas->Update();
1045 ++i;
1046 }
1047 }
1048
1049 func->Draw("same");
1050
1051 //TRatioPlot *rp = new TRatioPlot(hist);
1052 //rp->Draw();
1053 // draw residuals
1054 // calculate y offset
1055
1056 int iLowBin = hist->FindBin(low);
1057 int iHighBin = hist->FindBin(high);
1058 double min = canvas->GetUymin();
1059 double offset = min + (std::min(hist->GetBinContent(iLowBin), hist->GetBinContent(iHighBin)) - min) / 2.0;
1060
1061 TH1D *residuals = new TH1D("residuals", "residuals", iHighBin - iLowBin + 1,
1062 hist->GetBinCenter(iLowBin) - hist->GetBinWidth(iLowBin) / 2.0,
1063 hist->GetBinCenter(iHighBin) + hist->GetBinWidth(iHighBin) / 2.0);
1064 for (int i = 1; i <= iHighBin - iLowBin + 1; i++) {
1065 residuals->SetBinContent(i, hist->GetBinContent(i + iLowBin - 1) -
1066 func->Eval(hist->GetBinCenter(i + iLowBin - 1)) + offset);
1067 residuals->SetBinError(i, hist->GetBinError(i + iLowBin - 1));
1068 }
1069
1070 TF1 *zeroFunc = new TF1("zero", "[0]", low, high);
1071 zeroFunc->SetParameter(0, offset);
1072 residuals->SetMarkerSize(0.5);
1073 residuals->Draw("hist E1 same");
1074 zeroFunc->Draw("same");
1075
1076
1077 // update and wait for key press
1078 canvas->Update();
1079 if (!iQuiet){
1080 std::cout << "Press any key to exit" << std::endl;
1081 }
1082 click.waiting=false;
1083 canvas->WaitPrimitive();
1084
1085
1086 // clean up
1087 residuals->Delete();
1088 hist->Draw("hist");
1089
1090 canvas->SetCrosshair(0);
1091
1092 return retval;
1093 } /* function FitPeaks */
1094
1095 void PeakFit::AddPeak(std::string peakKey, GamR::TK::FPeak* peak) {
1096 fPeakKeys[(int)fPeaks.size()] = peakKey;
1097 fPeaks[peakKey] = peak;
1098 }
1099
1100 void PeakFit::AddPeak(std::string peakKey, GamR::TK::PeakType peaktype) {
1101 if (GetNPeaks()==0){
1102 fPeakType = peaktype;
1103 } else {
1104 if (fPeakType!=peaktype) {
1105 std::cout << "Error: peak types must be the same.";
1106 return;
1107 }
1108 }
1109 switch(fPeakType) {
1110 case GamR::TK::Gaussian: {
1111 GamR::TK::GaussianPeak *gpeak = new GamR::TK::GaussianPeak(fLow, fHigh);
1112 AddPeak(peakKey, gpeak);
1113 break;
1114 }
1115 case GamR::TK::StepGaussian: {
1116 GamR::TK::StepGaussianPeak *sgpeak = new GamR::TK::StepGaussianPeak(fLow, fHigh);
1117 AddPeak(peakKey, sgpeak);
1118 break;
1119 }
1120 case GamR::TK::OneTailGaussian: {
1121 GamR::TK::OneTailGaussianPeak *otgpeak = new GamR::TK::OneTailGaussianPeak(fLow, fHigh);
1122 AddPeak(peakKey, otgpeak);
1123 break;
1124 }
1125 case GamR::TK::TwoTailGaussian: {
1126 GamR::TK::TwoTailGaussianPeak *ttgpeak = new GamR::TK::TwoTailGaussianPeak(fLow, fHigh);
1127 AddPeak(peakKey, ttgpeak);
1128 break;
1129 }
1130 case GamR::TK::OneTailStepGaussian: {
1131 GamR::TK::OneTailStepGaussianPeak *otsgpeak = new GamR::TK::OneTailStepGaussianPeak(fLow, fHigh);
1132 AddPeak(peakKey, otsgpeak);
1133 break;
1134 }
1135 case GamR::TK::TwoTailStepGaussian: {
1136 GamR::TK::TwoTailStepGaussianPeak *ttsgpeak = new GamR::TK::TwoTailStepGaussianPeak(fLow, fHigh);
1137 AddPeak(peakKey, ttsgpeak);
1138 break;
1139 }
1140 }
1141 }
1142
1143 void PeakFit::AddPeak(std::string peakKey, GamR::TK::PeakType peaktype, double centroid) {
1144 //make correct type of peak and add, this will break if a new type of peak is added.
1145 AddPeak(peakKey, peaktype);
1146 SetCent(peakKey,centroid);
1147 }
1148
1149 double PeakFit::GetAmp(std::string peakKey) {
1150 return fPeaks[peakKey]->GetAmp();
1151 }
1152
1153 double PeakFit::GetAmpError(std::string peakKey) {
1154 return fPeaks[peakKey]->GetAmpError();
1155 }
1156
1157 double PeakFit::GetCent(std::string peakKey) {
1158 return fPeaks[peakKey]->GetCent();
1159 }
1160
1161 double PeakFit::GetCentError(std::string peakKey) {
1162 return fPeaks[peakKey]->GetCentError();
1163 }
1164
1165 double PeakFit::GetWidth(std::string peakKey) {
1166 return fPeaks[peakKey]->GetWidth();
1167 }
1168
1169 double PeakFit::GetWidthError(std::string peakKey) {
1170 return fPeaks[peakKey]->GetWidthError();
1171 }
1172
1173
1174 void PeakFit::SetAmp(std::string peakKey, double amp) {
1175 fPeaks[peakKey]->SetAmp(amp);
1176 }
1177
1178 void PeakFit::SetAmpError(std::string peakKey, double error) {
1179 fPeaks[peakKey]->SetAmpError(error);
1180 }
1181
1182 void PeakFit::SetCent(std::string peakKey, double centroid) {
1183 fPeaks[peakKey]->SetCent(centroid);
1184 }
1185
1186 void PeakFit::SetCentError(std::string peakKey, double error) {
1187 fPeaks[peakKey]->SetCentError(error);
1188 }
1189
1190 void PeakFit::SetWidth(std::string peakKey, double width) {
1191 fPeaks[peakKey]->SetWidth(width);
1192 }
1193
1194 void PeakFit::SetWidthError(std::string peakKey, double error) {
1195 fPeaks[peakKey]->SetWidthError(error);
1196 }
1197
1198 void PeakFit::SetFunction() {
1199 fTotal = new TF1("func",
1200 [&](double *x, double *p){
1201 double ret = 0;
1202 fBackground->SetParameter(0, p[0]);
1203 if (!parameters.iConstantBack) {
1204 fBackground->SetParameter(1, p[1]);
1205 }
1206 if (parameters.iQuadBack) {
1207 fBackground->SetParameter(2, p[2]);
1208 }
1209 ret += fBackground->Eval(x[0]);
1210
1211 for (auto &mpp : fPeakParamInds) {
1212 auto peakKey = mpp.first;
1213 auto pp = mpp.second;
1214 switch(fPeakType) {
1215 case GamR::TK::Gaussian: {
1216 if (parameters.iFixWidthsFile) {
1217 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1218 p[pp.iCentroid],
1219 gWidth->Eval(p[pp.iCentroid])/2.3548});
1220 }
1221 else {
1222 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1223 p[pp.iCentroid],
1224 p[pp.iWidth]});
1225 }
1226 break;
1227 }
1228 case GamR::TK::StepGaussian: {
1229 if (parameters.iFixWidthsFile) {
1230 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1231 p[pp.iCentroid],
1232 gWidth->Eval(p[pp.iCentroid])/2.3548,
1233 gStep->Eval(p[pp.iCentroid])});
1234 }
1235 else {
1236 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1237 p[pp.iCentroid],
1238 p[pp.iWidth],
1239 p[pp.iStepAmplitude]});
1240 }
1241 break;
1242 }
1243 case GamR::TK::OneTailGaussian: {
1244 if (parameters.iFixWidthsFile) {
1245 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1246 p[pp.iCentroid],
1247 gWidth->Eval(p[pp.iCentroid])/2.3548,
1248 gSkewAmp->Eval(p[pp.iCentroid]),
1249 gSkew->Eval(p[pp.iCentroid])});
1250 }
1251 else {
1252 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1253 p[pp.iCentroid],
1254 p[pp.iWidth],
1255 p[pp.iSkewAmplitude],
1256 p[pp.iSkewWidth]});
1257 }
1258 break;
1259 }
1260 case GamR::TK::TwoTailGaussian: {
1261 if (parameters.iFixWidthsFile) {
1262 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1263 p[pp.iCentroid],
1264 gWidth->Eval(p[pp.iCentroid])/2.3548,
1265 gSkewAmp->Eval(p[pp.iCentroid]),
1266 gSkew->Eval(p[pp.iCentroid]),
1267 gSkewAmp2->Eval(p[pp.iCentroid]),
1268 gSkew2->Eval(p[pp.iCentroid])});
1269 }
1270 else {
1271 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1272 p[pp.iCentroid],
1273 p[pp.iWidth],
1274 p[pp.iSkewAmplitude],
1275 p[pp.iSkewWidth],
1276 p[pp.iSkewAmplitude2],
1277 p[pp.iSkewWidth2]});
1278 }
1279 break;
1280 }
1281 case GamR::TK::OneTailStepGaussian: {
1282 if (parameters.iFixWidthsFile) {
1283 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1284 p[pp.iCentroid],
1285 gWidth->Eval(p[pp.iCentroid])/2.3548,
1286 gSkewAmp->Eval(p[pp.iCentroid]),
1287 gSkew->Eval(p[pp.iCentroid]),
1288 gStep->Eval(p[pp.iCentroid])});
1289 }
1290 else {
1291 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1292 p[pp.iCentroid],
1293 p[pp.iWidth],
1294 p[pp.iSkewAmplitude],
1295 p[pp.iSkewWidth],
1296 p[pp.iStepAmplitude]});
1297 }
1298 break;
1299 }
1300 case GamR::TK::TwoTailStepGaussian: {
1301 if (parameters.iFixWidthsFile) {
1302 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1303 p[pp.iCentroid],
1304 gWidth->Eval(p[pp.iCentroid])/2.3548,
1305 gSkewAmp->Eval(p[pp.iCentroid]),
1306 gSkew->Eval(p[pp.iCentroid]),
1307 gSkewAmp2->Eval(p[pp.iCentroid]),
1308 gSkew2->Eval(p[pp.iCentroid]),
1309 gStep->Eval(p[pp.iCentroid])});
1310 }
1311 else {
1312 fPeaks[peakKey]->SetParameters({p[pp.iAmplitude],
1313 p[pp.iCentroid],
1314 p[pp.iWidth],
1315 p[pp.iSkewAmplitude],
1316 p[pp.iSkewWidth],
1317 p[pp.iSkewAmplitude2],
1318 p[pp.iSkewWidth2],
1319 p[pp.iStepAmplitude]});
1320 }
1321 break;
1322 }
1323 }
1324
1325 ret += fPeaks[peakKey]->GetFunc()->Eval(x[0]);
1326 }
1327 return ret;
1328 }, fLow, fHigh, iParamCount);
1329
1330 for (auto & mPeak : fPeaks) {
1331 auto peakKey = mPeak.first;
1332 auto peak = mPeak.second;
1333 auto func = peak->GetFunc();
1334 auto pp = fPeakParamInds[peakKey];
1335 switch(fPeakType) {
1336 case GamR::TK::Gaussian: {
1337 double low;
1338 double high;
1339 func->GetParLimits(0,low,high);
1340 if (low!=high) {
1341 fTotal->SetParLimits(pp.iAmplitude, low, high);
1342 }
1343 func->GetParLimits(1,low,high);
1344 if (low!=high) {
1345 fTotal->SetParLimits(pp.iCentroid, low, high);
1346 }
1347 func->GetParLimits(2,low,high);
1348 if (low!=high) {
1349 fTotal->SetParLimits(pp.iWidth, low, high);
1350 }
1351 break;
1352 }
1353 case GamR::TK::StepGaussian: {
1354 double low;
1355 double high;
1356 func->GetParLimits(0,low,high);
1357 if (low!=high) {
1358 fTotal->SetParLimits(pp.iAmplitude, low, high);
1359 }
1360 func->GetParLimits(1,low,high);
1361 if (low!=high) {
1362 fTotal->SetParLimits(pp.iCentroid, low, high);
1363 }
1364 if (!parameters.iFixWidthsFile) {
1365 func->GetParLimits(2,low,high);
1366 if (low!=high) {
1367 fTotal->SetParLimits(pp.iWidth, low, high);
1368 }
1369 func->GetParLimits(3,low,high);
1370 if (low!=high) {
1371 fTotal->SetParLimits(pp.iStepAmplitude, low, high);
1372 }
1373 }
1374 break;
1375 }
1376 case GamR::TK::OneTailGaussian: {
1377 double low;
1378 double high;
1379 func->GetParLimits(0,low,high);
1380 if (low!=high) {
1381 fTotal->SetParLimits(pp.iAmplitude, low, high);
1382 }
1383 func->GetParLimits(1,low,high);
1384 if (low!=high) {
1385 fTotal->SetParLimits(pp.iCentroid, low, high);
1386 }
1387 if (!parameters.iFixWidthsFile) {
1388 func->GetParLimits(2,low,high);
1389 if (low!=high) {
1390 fTotal->SetParLimits(pp.iWidth, low, high);
1391 }
1392 func->GetParLimits(3,low,high);
1393 if (low!=high) {
1394 fTotal->SetParLimits(pp.iSkewAmplitude, low, high);
1395 }
1396 func->GetParLimits(4,low,high);
1397 if (low!=high) {
1398 fTotal->SetParLimits(pp.iSkewWidth, low, high);
1399 }
1400 }
1401 break;
1402 }
1403 case GamR::TK::TwoTailGaussian: {
1404 double low;
1405 double high;
1406 func->GetParLimits(0,low,high);
1407 if (low!=high) {
1408 fTotal->SetParLimits(pp.iAmplitude, low, high);
1409 }
1410 func->GetParLimits(1,low,high);
1411 if (low!=high) {
1412 fTotal->SetParLimits(pp.iCentroid, low, high);
1413 }
1414 if (!parameters.iFixWidthsFile) {
1415 func->GetParLimits(2,low,high);
1416 if (low!=high) {
1417 fTotal->SetParLimits(pp.iWidth, low, high);
1418 }
1419 func->GetParLimits(3,low,high);
1420 if (low!=high) {
1421 fTotal->SetParLimits(pp.iSkewAmplitude, low, high);
1422 }
1423 func->GetParLimits(4,low,high);
1424 if (low!=high) {
1425 fTotal->SetParLimits(pp.iSkewWidth, low, high);
1426 }
1427 func->GetParLimits(5,low,high);
1428 if (low!=high) {
1429 fTotal->SetParLimits(pp.iSkewAmplitude2, low, high);
1430 }
1431 func->GetParLimits(6,low,high);
1432 if (low!=high) {
1433 fTotal->SetParLimits(pp.iSkewWidth2, low, high);
1434 }
1435 }
1436 break;
1437 }
1438 case GamR::TK::OneTailStepGaussian: {
1439 double low;
1440 double high;
1441 func->GetParLimits(0,low,high);
1442 if (low!=high) {
1443 fTotal->SetParLimits(pp.iAmplitude, low, high);
1444 }
1445 func->GetParLimits(1,low,high);
1446 if (low!=high) {
1447 fTotal->SetParLimits(pp.iCentroid, low, high);
1448 }
1449 if (!parameters.iFixWidthsFile) {
1450 func->GetParLimits(2,low,high);
1451 if (low!=high) {
1452 fTotal->SetParLimits(pp.iWidth, low, high);
1453 }
1454 func->GetParLimits(3,low,high);
1455 if (low!=high) {
1456 fTotal->SetParLimits(pp.iSkewAmplitude, low, high);
1457 }
1458 func->GetParLimits(4,low,high);
1459 if (low!=high) {
1460 fTotal->SetParLimits(pp.iSkewWidth, low, high);
1461 }
1462 func->GetParLimits(5,low,high);
1463 if (low!=high) {
1464 fTotal->SetParLimits(pp.iStepAmplitude, low, high);
1465 }
1466 }
1467 break;
1468 }
1469 case GamR::TK::TwoTailStepGaussian: {
1470 double low;
1471 double high;
1472 func->GetParLimits(0,low,high);
1473 if (low!=high) {
1474 fTotal->SetParLimits(pp.iAmplitude, low, high);
1475 }
1476 func->GetParLimits(1,low,high);
1477 if (low!=high) {
1478 fTotal->SetParLimits(pp.iCentroid, low, high);
1479 }
1480 if (!parameters.iFixWidthsFile) {
1481 func->GetParLimits(2,low,high);
1482 if (low!=high) {
1483 fTotal->SetParLimits(pp.iWidth, low, high);
1484 }
1485 func->GetParLimits(3,low,high);
1486 if (low!=high) {
1487 fTotal->SetParLimits(pp.iSkewAmplitude, low, high);
1488 }
1489 func->GetParLimits(4,low,high);
1490 if (low!=high) {
1491 fTotal->SetParLimits(pp.iSkewWidth, low, high);
1492 }
1493 func->GetParLimits(5,low,high);
1494 if (low!=high) {
1495 fTotal->SetParLimits(pp.iSkewAmplitude2, low, high);
1496 }
1497 func->GetParLimits(6,low,high);
1498 if (low!=high) {
1499 fTotal->SetParLimits(pp.iSkewWidth2, low, high);
1500 }
1501 func->GetParLimits(7,low,high);
1502 if (low!=high) {
1503 fTotal->SetParLimits(pp.iStepAmplitude, low, high);
1504 }
1505 }
1506 break;
1507 }
1508 }
1509 }
1510
1511 fNPars = fTotal->GetNpar();
1512 }
1513
1514 TF1 *PeakFit::GetPeakBG(std::string peakKey) {
1515 TF1 *PeakBG = new TF1("PeakBG", [&](double *x, double *p){ return (FindPeak((int)p[0]))->GetFunc()->Eval(x[0]) + fBackground->Eval(x[0]);}, fLow, fHigh, 1);
1516 int i=0;
1517 for (auto &peak : fPeaks) {
1518 auto key = peak.first;
1519 if (key==peakKey) {
1520 PeakBG->SetParameter(0, i);
1521 break;
1522 }
1523 ++i;
1524 }
1525 return PeakBG;
1526 }
1527
1528 TF1 *PeakFit::GetPeakGaussBG(std::string peakKey) {
1529 if (FindPeak(0)->GetGaussFunc() == NULL) { return NULL; }
1530 TF1 *GaussPeakBG = new TF1("GaussPeakBG", [&](double *x, double *p){ return (FindPeak((int)p[0]))->GetGaussFunc()->Eval(x[0]) + fBackground->Eval(x[0]);}, fLow, fHigh, 1);
1531 int i=0;
1532 for (auto &peak : fPeaks) {
1533 auto key = peak.first;
1534 if (key==peakKey) {
1535 GaussPeakBG->SetParameter(0, i);
1536 break;
1537 }
1538 ++i;
1539 }
1540 return GaussPeakBG;
1541 }
1542
1543 TF1 *PeakFit::GetPeakStepBG(std::string peakKey) {
1544 if (FindPeak(0)->GetStepFunc() == NULL) { return NULL; }
1545 TF1 *PeakStepBG = new TF1("PeakStepBG", [&](double *x, double *p){ return (FindPeak((int)p[0]))->GetStepFunc()->Eval(x[0]) + fBackground->Eval(x[0]);}, fLow, fHigh, 1);
1546 int i=0;
1547 for (auto &peak : fPeaks) {
1548 auto key = peak.first;
1549 if (key==peakKey) {
1550 PeakStepBG->SetParameter(0, i);
1551 break;
1552 }
1553 ++i;
1554 }
1555 return PeakStepBG;
1556 }
1557
1558 void PeakFit::ToText(std::string filename, std::string delimiter, int nPoints) {
1559 std::ofstream outfile(filename);
1560
1561 std::vector<TF1*> peaks;
1562
1563 for (int i=0; i<GetNPeaks(); ++i) {
1564 peaks.push_back(GetPeakBG(i));
1565 }
1566
1567 for (int i=0; i<nPoints; ++i) {
1568 double x = fLow + (fHigh-fLow)*(double)(i)/(double)(nPoints-1);
1569 outfile << x << delimiter << fTotal->Eval(x) << delimiter;
1570 outfile << fBackground->Eval(x) << delimiter;
1571 for (int j=0; j<GetNPeaks(); ++j) {
1572 outfile << peaks[j]->Eval(x) << delimiter;
1573 }
1574 outfile << std::endl;
1575 }
1576
1577 outfile.close();
1578 }
1579
1580 TF1 *FitGaussians(TH1 *hist, double Low, double High, std::vector<double> Peaks, Option_t *foption /*=""*/,
1581 Option_t *option /*=""*/) {
1582 auto retval = FitPeaks(hist, Low, High, Peaks, foption, option);
1583 TF1 *out = retval->GetTotal();
1584 return out;
1585 }
1586
1587 TF1 *FitGaussians(TCanvas *canvas /*= gPad->GetCanvas()*/, Option_t *foption /*= ""*/, Option_t *option /*= ""*/) {
1588 auto retval = FitPeaks(canvas, foption, option);
1589 TF1 *out = retval->GetTotal();
1590 return out;
1591 }
1592
1593 TSpectrum *FindPeaks(TCanvas *canvas /* = gPad->GetCanvas()*/,
1594 double sigma /* = 2*/,
1595 Option_t *option /*=""*/,
1596 double threshold /* = 0.05*/) {
1597 TH1* hist = GamR::Utils::GetHist1D(canvas);
1598
1599 TSpectrum *spectrum = new TSpectrum();
1600 spectrum->Search(hist, sigma, option, threshold);
1601
1602 //eliminate old labels
1603 while (true) {
1604 TText *text = (TText*)hist->GetListOfFunctions()->FindObject("PeakLabel");
1605 if (text) {
1606 hist->GetListOfFunctions()->Remove(text);
1607 delete text;
1608 }
1609 else {
1610 break;
1611 }
1612 }
1613
1614 for (int i=0; i<spectrum->GetNPeaks(); ++i) {
1615 double x = spectrum->GetPositionX()[i];
1616 double y = spectrum->GetPositionY()[i];
1617 char t[100];
1618 sprintf(t,"%.1f", x);
1619 TText *text = new TText(x, y, t);
1620 text->SetName("PeakLabel");
1621 text->SetTextAngle(90);
1622 text->SetTextSize(0.025);
1623 hist->GetListOfFunctions()->Add(text);
1624 }
1625 spectrum->Print();
1626 hist->Draw(option);
1627 return spectrum;
1628 }
1629
1630 // TSpectrum *FindPeaks(TCanvas *canvas /* = gPad->GetCanvas()*/,
1631 // double sigma /* = 2*/,
1632 // Option_t *option /*=""*/,
1633 // double threshold /* = 0.05*/) {
1634 // TH1 *hist = GamR::Utils::GetHist1D(canvas);
1635
1636 // const int nBins = hist->GetNbinsX();
1637
1638 // double in_hist[nBins];
1639 // double out_hist[nBins];
1640
1641 // for (int i=0; i<nBins; ++i) {
1642 // in_hist[i] = hist->GetBinContent(i+1);
1643 // }
1644
1645 // TSpectrum *spectrum = new TSpectrum();
1646 // spectrum->SearchHighRes(&in_hist[0], &out_hist[0], nBins, sigma, threshold*100, true, 3, true, 3);
1647
1648 // while (true) {
1649 // TText *text = (TText*)hist->GetListOfFunctions()->FindObject("PeakLabel");
1650 // if (text) {
1651 // hist->GetListOfFunctions()->Remove(text);
1652 // delete text;
1653 // }
1654 // else {
1655 // break;
1656 // }
1657 // }
1658 // TPolyMarker *marks = new TPolyMarker();
1659 // for (int i=0; i<spectrum->GetNPeaks(); ++i) {
1660 // double x = spectrum->GetPositionX()[i];
1661 // int bin = 1 + Int_t(x + 0.5);
1662 // double le = hist->GetBinLowEdge(bin);
1663 // double he = hist->GetBinLowEdge(bin+1);
1664 // double frac = x - floor(x);
1665 // x = le + frac*(he-le);
1666
1667 // double y = hist->GetBinContent(bin);
1668 // marks->SetNextPoint(x, y);
1669 // char t[100];
1670 // sprintf(t,"%.1f", x);
1671 // TText *text = new TText(x, y, t);
1672 // text->SetName("PeakLabel");
1673 // text->SetTextAngle(90);
1674 // text->SetTextSize(0.025);
1675 // hist->GetListOfFunctions()->Add(text);
1676 // }
1677
1678 // hist->GetListOfFunctions()->Add(marks);
1679 // marks->SetMarkerStyle(23);
1680 // marks->SetMarkerColor(kRed);
1681 // marks->SetMarkerSize(1.3);
1682 // hist->Draw(option);
1683 // return spectrum;
1684 // }
1685
1686 } /* namespace Spect */
1687} /* namespace GamR */
ss
void ss(std::vector< int > indexes, TCanvas *canvas, Option_t *option)
Definition Abbreviations.cc:205
GamR::Spect::PeakFit::SetWidth
void SetWidth(std::string peakKey, double width)
Definition Fit.cc:1190
GamR::Spect::PeakFit::ToText
void ToText(std::string filename, std::string delimiter=" ", int nPoints=1000)
Definition Fit.cc:1558
GamR::Spect::PeakFit::GetPeakBG
TF1 * GetPeakBG(std::string peakKey)
Definition Fit.cc:1514
GamR::Spect::PeakFit::GetPeakStepBG
TF1 * GetPeakStepBG(std::string peakKey)
Definition Fit.cc:1543
GamR::Spect::PeakFit::GetLow
double GetLow()
Definition Fit.hh:96
GamR::Spect::PeakFit::GetPeakType
GamR::TK::PeakType GetPeakType()
Definition Fit.hh:92
GamR::Spect::PeakFit::SetPeaks
void SetPeaks(std::vector< double > Peaks)
Definition Fit.cc:118
GamR::Spect::PeakFit::GetBackground
TF1 * GetBackground()
Definition Fit.hh:93
GamR::Spect::PeakFit::SetOpts
void SetOpts(Option_t *option)
Definition Fit.cc:39
GamR::Spect::PeakFit::PrintResults
void PrintResults()
Definition Fit.cc:746
GamR::Spect::PeakFit::fFitGuesses
PeakFitGuesses * fFitGuesses
Definition Fit.hh:86
GamR::Spect::PeakFit::SetCentError
void SetCentError(std::string peakKey, double error)
Definition Fit.cc:1186
GamR::Spect::PeakFit::SetLimit
void SetLimit(std::string peakKey, std::string parameterKey, double lowerLimit, double upperLimit)
Definition Fit.cc:424
GamR::Spect::PeakFit::GetAmpError
double GetAmpError(std::string peakKey)
Definition Fit.cc:1153
GamR::Spect::PeakFit::GetCent
double GetCent(std::string peakKey)
Definition Fit.cc:1157
GamR::Spect::PeakFit::GetPeaks
std::map< std::string, GamR::TK::FPeak * > GetPeaks()
Definition Fit.hh:190
GamR::Spect::PeakFit::SetWidthError
void SetWidthError(std::string peakKey, double error)
Definition Fit.cc:1194
GamR::Spect::PeakFit::SetCent
void SetCent(std::string peakKey, double centroid)
Definition Fit.cc:1182
GamR::Spect::PeakFit::SetAmp
void SetAmp(std::string peakKey, double amp)
Definition Fit.cc:1174
GamR::Spect::PeakFit::AddPeak
void AddPeak(std::string peakKey, GamR::TK::FPeak *peak)
Definition Fit.cc:1095
GamR::Spect::PeakFit::SetUp
void SetUp(TH1 *hist, std::map< std::string, double > Peaks, Option_t *option, std::vector< std::string > FixPeaks={})
Definition Fit.hh:129
GamR::Spect::PeakFit::GetNPeaks
size_t GetNPeaks()
Definition Fit.hh:98
GamR::Spect::PeakFit::FixPeakEnergy
void FixPeakEnergy(std::string peakKey)
Definition Fit.hh:173
GamR::Spect::PeakFit::GetWidthError
double GetWidthError(std::string peakKey)
Definition Fit.cc:1169
GamR::Spect::PeakFit::SetResults
void SetResults(TH1 *hist)
Definition Fit.cc:551
GamR::Spect::PeakFit::GetTotal
TF1 * GetTotal()
Definition Fit.hh:94
GamR::Spect::PeakFit::GetPeakGaussBG
TF1 * GetPeakGaussBG(std::string peakKey)
Definition Fit.cc:1528
GamR::Spect::PeakFit::GetHigh
double GetHigh()
Definition Fit.hh:97
GamR::Spect::PeakFit::SetIndices
void SetIndices(std::map< std::string, GamR::TK::FPeak * > Peaks)
Definition Fit.cc:137
GamR::Spect::PeakFit::SetChi2
void SetChi2(double Chi2)
Definition Fit.hh:128
GamR::Spect::PeakFit::GetCentError
double GetCentError(std::string peakKey)
Definition Fit.cc:1161
GamR::Spect::PeakFit::GetAmp
double GetAmp(std::string peakKey)
Definition Fit.cc:1149
GamR::Spect::PeakFit::Fit
void Fit(TH1 *hist, Option_t *foption)
Definition Fit.cc:433
GamR::Spect::PeakFit::SetBackground
void SetBackground()
Definition Fit.cc:98
GamR::Spect::PeakFit::SetAmpError
void SetAmpError(std::string peakKey, double error)
Definition Fit.cc:1178
GamR::Spect::PeakFit::GetWidth
double GetWidth(std::string peakKey)
Definition Fit.cc:1165
GamR::TK::PeakType
PeakType
Definition Peak.hh:33
GamR::TK::OneTailStepGaussian
@ OneTailStepGaussian
Definition Peak.hh:33
GamR::TK::StepGaussian
@ StepGaussian
Definition Peak.hh:33
GamR::TK::Gaussian
@ Gaussian
Definition Peak.hh:33
GamR::TK::TwoTailGaussian
@ TwoTailGaussian
Definition Peak.hh:33
GamR::TK::OneTailGaussian
@ OneTailGaussian
Definition Peak.hh:33
GamR::TK::TwoTailStepGaussian
@ TwoTailStepGaussian
Definition Peak.hh:33
GamR::Spect::FitGaussians
TF1 * FitGaussians(TH1 *hist, double Low, double High, std::vector< double > Peaks, Option_t *foption, Option_t *option)
Definition Fit.cc:1580
GamR::Spect::FindPeaks
TSpectrum * FindPeaks(TCanvas *canvas, double sigma, Option_t *option, double threshold)
Definition Fit.cc:1593
GamR::Spect::gFitGuesses
PeakFitGuesses * gFitGuesses
GamR::Spect::FitPeaks
GamR::Spect::PeakFit * FitPeaks(TH1 *hist, double Low, double High, std::vector< double > Peaks, Option_t *foption, Option_t *option, std::vector< std::string > FixPeaks)
Definition Fit.cc:787
GamR::Utils::GetHist1D
TH1D * GetHist1D(TVirtualPad *canvas)
Definition Utilities.cc:182
GamR::Utils::wrresult
int wrresult(char *out, float value, float err, int minlen)
Definition Utilities.cc:241
GamR
Definition Gain.cc:19