source: trunk/MagicSoft/Mars/mhflux/MAlphaFitter.h@ 7384

Last change on this file since 7384 was 7384, checked in by tbretz, 19 years ago
*** empty log message ***
File size: 8.3 KB
Line 
1#ifndef MARS_MAlphaFitter
2#define MARS_MAlphaFitter
3
4#ifndef MARS_MParContainer
5#include "MParContainer.h"
6#endif
7
8#ifndef ROOT_TArrayD
9#include <TArrayD.h>
10#endif
11
12#ifndef ROOT_TF1
13#include <TF1.h>
14#endif
15
16class TH1D;
17class TH3D;
18
19class MAlphaFitter : public MParContainer
20{
21public:
22 enum ScaleMode_t {
23 kNone, // No scaling
24 kEntries, // scale by the number of entries in on and off
25 kIntegral, // scale by the integral in on and off
26 kOffRegion, // scale by the integral between fScaleMin, fScaleMax in on and off
27 kBackground, // scale by the integral between fBgMin, fBgMax in on and off
28 kLeastSquare, // not yet implemented
29 kUserScale // scale by fixed factor set by SetScaleUser
30 };
31 enum Strategy_t {
32 kSignificance,
33 kSignificanceChi2,
34 kSignificanceLogExcess,
35 kSignificanceExcess,
36 kExcess,
37 kGaussSigma
38 };
39 enum SignalFunc_t {
40 kGauss, kThetaSq
41 };
42
43private:
44 // Fitting Setup
45 Float_t fSigInt; // minimum of range to fit the signal
46 Float_t fSigMax; // maximum of range to fit the signal
47 Float_t fBgMin; // minimum of range to fit the background
48 Float_t fBgMax; // minimum of range to fit the background
49 Float_t fScaleMin; // minimum of range to determin the scale factor of the background
50 Float_t fScaleMax; // maximum of range to determin the scale factor of the background
51 Int_t fPolynomOrder; // order of polyom to be fitted to the background
52 Bool_t fFitBackground; // Backround fit: yes/no
53 SignalFunc_t fSignalFunc; // Type of signal function
54 // Result
55 Double_t fSignificance; // significance of signal
56 Double_t fSignificanceExc; // significance of excess
57 Double_t fEventsExcess; // calculated number of excess events (signal-bg)
58 Double_t fEventsSignal; // calculated number of signal events
59 Double_t fEventsBackground; // calculated number of bg events (fScaleFactor already applied)
60
61 Double_t fChiSqSignal; // Reduced (chi^2/NDF) chisq of signal fit
62 Double_t fChiSqBg; // Reduced (chi^2/NDF) chisq of bg fit
63 Double_t fIntegralMax; // Calculated bin border to which it was integrated
64 Double_t fScaleFactor; // Scale factor determined for off-data
65
66 TArrayD fCoefficients; // Fit result
67
68 // Function
69 TF1 *fFunc; // fit function (gauss + polynom)
70
71 // Scaling setup
72 ScaleMode_t fScaleMode; // scaling mode
73 Double_t fScaleUser; // user scale factor
74
75 // Minimization strategy
76 Strategy_t fStrategy; // How to calc minimization value
77
78public:
79 // Implementing the function yourself is only about 5% faster
80 MAlphaFitter(const char *name=0, const char *title=0) : fSigInt(15),
81 fSigMax(75), fBgMin(45), fBgMax(85), fScaleMin(40), fScaleMax(80),
82 fPolynomOrder(2), fFitBackground(kTRUE), fSignalFunc(kGauss),
83 fCoefficients(3+fPolynomOrder+1),
84 fFunc(new TF1("", Form("gaus(0) + pol%d(3)", fPolynomOrder), 0, 90)),
85 fScaleMode(kOffRegion), fScaleUser(1), fStrategy(kSignificance)
86 {
87 fName = name ? name : "MAlphaFitter";
88 fTitle = title ? title : "Fit alpha";
89
90 gROOT->GetListOfFunctions()->Remove(fFunc);
91 fFunc->SetName("Dummy");
92
93 Clear();
94 }
95
96 MAlphaFitter(const MAlphaFitter &f) : fFunc(0)
97 {
98 f.Copy(*this);
99 }
100 ~MAlphaFitter()
101 {
102 delete fFunc;
103 }
104
105 // TObject
106 void Clear(Option_t *o="");
107 void Print(Option_t *o="") const; //*MENU*
108 void Copy(TObject &o) const;
109
110 // Setter
111 void SetScaleUser(Float_t scale) { fScaleUser = scale; fScaleMode=kUserScale; }
112 void SetScaleMode(ScaleMode_t mode) { fScaleMode = mode; }
113 void SetMinimizationStrategy(Strategy_t mode) { fStrategy = mode; }
114 void SetSignalIntegralMax(Float_t s) { fSigInt = s; }
115 void SetSignalFitMax(Float_t s) { fSigMax = s; }
116 void SetBackgroundFitMin(Float_t s) { fBgMin = s; }
117 void SetBackgroundFitMax(Float_t s) { fBgMax = s; }
118 void SetScaleMin(Float_t s) { fScaleMin = s; }
119 void SetScaleMax(Float_t s) { fScaleMax = s; }
120 void SetPolynomOrder(Int_t s)
121 {
122 if (s==fPolynomOrder)
123 return;
124
125 fPolynomOrder = s;
126
127 SetSignalFunction(fSignalFunc);
128 }
129 void SetSignalFunction(SignalFunc_t func)
130 {
131 delete fFunc;
132 switch (func)
133 {
134 case kGauss:
135 fFunc=new TF1 ("", Form("gaus(0) + pol%d(3)", fPolynomOrder));
136 break;
137 case kThetaSq:
138 fFunc=new TF1 ("", Form("[0]*exp(-0.5*((sqrt(x)-[1])/[2])^2) + pol%d(3)", fPolynomOrder));
139 break;
140 }
141 fSignalFunc=func;
142 gROOT->GetListOfFunctions()->Remove(fFunc);
143 fFunc->SetName("Dummy");
144 fCoefficients.Set(3+fPolynomOrder+1);
145 fCoefficients.Reset();
146 }
147 void EnableBackgroundFit(Bool_t b=kTRUE) { fFitBackground=b; }
148
149 // Getter
150 Double_t GetSignalIntegralMax() const { return fSigInt; }
151
152 Double_t GetEventsExcess() const { return fEventsExcess; }
153 Double_t GetEventsSignal() const { return fEventsSignal; }
154 Double_t GetEventsBackground() const { return fEventsBackground; }
155
156 Double_t GetSignificance() const { return fSignificance; }
157 Double_t GetSignificanceExc() const { return fSignificanceExc; }
158 Double_t GetChiSqSignal() const { return fChiSqSignal; }
159 Double_t GetChiSqBg() const { return fChiSqBg; }
160 Double_t GetScaleFactor() const { return fScaleFactor; }
161 Double_t GetMinimizationValue() const;
162
163 Double_t GetGausSigma() const { return fCoefficients[2]; }
164 Double_t GetGausMu() const { return fCoefficients[1]; }
165 Double_t GetGausA() const { return fCoefficients[0]; }
166 Double_t GetCoefficient(Int_t i) const { return fCoefficients[i]; }
167 const TArrayD &GetCoefficients() const { return fCoefficients; }
168 Double_t Eval(Double_t d) const { return fFunc ? fFunc->Eval(d) : 0; }
169
170 // Interface to fit
171 Bool_t Fit(TH1D &h, Bool_t paint=kFALSE);
172 Bool_t Fit(const TH1D &on, const TH1D &off, Double_t alpha, Bool_t paint=kFALSE);
173 Bool_t Fit(TH1D &on, TH1D *off, Double_t alpha, Bool_t paint=kFALSE)
174 {
175 return off ? Fit(on, *off, alpha, paint) : Fit(on, paint);
176 }
177 Bool_t Fit(TH1D &on, TH1D *off, Bool_t paint=kFALSE)
178 {
179 return off ? Fit(on, *off, 1, paint) : Fit(on, paint);
180 }
181 Bool_t ScaleAndFit(TH1D &on, TH1D *off, Bool_t paint=kFALSE)
182 {
183 const Double_t alpha = off ? Scale(*off, on) : 1;
184 return off ? Fit(on, *off, alpha, paint) : Fit(on, paint);
185 }
186
187 Bool_t FitAlpha(const TH3D &h, Bool_t paint=kFALSE);
188 Bool_t FitEnergy(const TH3D &h, UInt_t bin, Bool_t paint=kFALSE);
189 Bool_t FitTheta(const TH3D &h, UInt_t bin, Bool_t paint=kFALSE);
190 //Bool_t FitTime(const TH3D &h, UInt_t bin, Bool_t paint=kFALSE);
191
192 Bool_t FitAlpha(const TH3D &on, const TH3D &off, Bool_t paint=kFALSE);
193 Bool_t FitEnergy(const TH3D &on, const TH3D &off, UInt_t bin, Bool_t paint=kFALSE);
194 Bool_t FitTheta(const TH3D &on, const TH3D &off, UInt_t bin, Bool_t paint=kFALSE);
195 //Bool_t FitTime(const TH3D &on, const TH3D &off, UInt_t bin, Bool_t paint=kFALSE);
196
197 Bool_t FitAlpha(const TH3D &on, const TH3D *off, Bool_t paint=kFALSE)
198 {
199 return off ? FitAlpha(on, *off, paint) : FitAlpha(on, paint);
200 }
201 Bool_t FitEnergy(const TH3D &on, const TH3D *off, UInt_t bin, Bool_t paint=kFALSE)
202 {
203 return off ? FitEnergy(on, *off, bin, paint) : FitEnergy(on, bin, paint);
204 }
205 Bool_t FitTheta(const TH3D &on, const TH3D *off, UInt_t bin, Bool_t paint=kFALSE)
206 {
207 return off ? FitTheta(on, *off, bin, paint) : FitTheta(on, bin, paint);
208 }/*
209 Bool_t FitTime(const TH3D &on, const TH3D *off, UInt_t bin, Bool_t paint=kFALSE)
210 {
211 return off ? FitTime(on, *off, bin, paint) : FitTime(on, bin, paint);
212 }*/
213
214 Double_t Scale(TH1D &off, const TH1D &on) const;
215
216 // Interface to result
217 void PaintResult(Float_t x=0.04, Float_t y=0.94, Float_t size=0.035) const;
218
219 // MTask
220 Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE);
221
222 ClassDef(MAlphaFitter, 1)
223};
224
225#endif
Note: See TracBrowser for help on using the repository browser.