GammaWare  Head Version for release 0.9
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
H1Calibrator.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  * Copyright (C) 2004-2006 by Olivier Stezowski & Christophe Theisen *
3  * stezow(AT)ipnl.in2p3.fr, christophe.theisen(AT)cea.fr *
4  * *
5  * This program is free software; you can redistribute it and/or modify *
6  * it under the terms of the GNU General Public License as published by *
7  * the Free Software Foundation; either version 2 of the License, or *
8  * (at your option) any later version. *
9  * *
10  * This program is distributed in the hope that it will be useful, *
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
13  * GNU General Public License for more details. *
14  * *
15  * You should have received a copy of the GNU General Public License *
16  * along with this program; if not, write to the *
17  * Free Software Foundation, Inc., *
18  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
19  ***************************************************************************/
20 
23 #ifndef GW_H1CALIBRATOR_H
24 #include "H1Calibrator.h"
25 #endif
26 
27 #ifndef ROOT_TROOT
28 #include "TROOT.h"
29 #endif
30 #ifndef ROOT_TPad
31 #include "TPad.h"
32 #endif
33 #ifndef ROOT_TCanvas
34 #include "TCanvas.h"
35 #endif
36 #ifndef ROOT_TArrayD
37 #include "TArrayD.h"
38 #endif
39 #ifndef ROOT_TString
40 #include "TString.h"
41 #endif
42 
43 #include <iostream>
44 #include <fstream>
45 
46 using namespace std;
47 using namespace Gw;
48 
49 // ROOT dictionnary
50 ClassImp(H1Calibrator);
51 
52 //__________________________________________________________________________________________________
53 H1Calibrator::H1Calibrator()
54  : fXRawAxis(0x0),
55  fYRawAxis(0x0),
56  fZRawAxis(0x0),
57  fHistoCal(0x0)
58 {
59 
60 }
61 
62 //__________________________________________________________________________________________________
63 H1Calibrator::~H1Calibrator() //destructor
64 {
65  if (fXRawAxis != 0) delete fXRawAxis;
66  if (fYRawAxis != 0) delete fYRawAxis;
67  if (fZRawAxis != 0) delete fZRawAxis;
68  if (fHistoCal != 0) delete fHistoCal;
69 }
70 
71 //__________________________________________________________________________________________________
72 void H1Calibrator::Calibrate(const TH1 *histo, const TF1 *function, const Option_t *option)
73 {
74  if (!histo) {
75  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
76  return;
77  }
78 
79  if (!function) {
80  printf("Info in <H1Calibrator::Calibrate> function NOT defined !!\n");
81  return;
82  }
83 
84  TString opt = option;
85  opt.ToLower();
86  if ( (!opt.Contains("x")) && (!opt.Contains("y")) && (!opt.Contains("z")) ) opt = "x";
87  if (histo->GetDimension() == 1) opt = "x";
88  if ( (histo->GetDimension() == 2) && (!opt.Contains("x")) && (!opt.Contains("y")) ) opt = "x";
89 
90  //Set the calibrator
91  Set(histo);
92 
93  if (opt.Contains("x")) {
94  CalibOneAxis(function,"x");
95  }
96  if (opt.Contains("y")) {
97  CalibOneAxis(function,"y");
98  }
99  if (opt.Contains("z")) {
100  CalibOneAxis(function,"z");
101  }
102 
103 
104  DrawResult();
105 }
106 
107 //__________________________________________________________________________________________________
108 void H1Calibrator::Calibrate(const TH1 *histo, const Double_t a, const Double_t b,
109  const Double_t c, const Option_t *option)
110 {
111  //Calibration with an order 2 polynome (a + bX + c*X*X)
112  if ((b==0.0) && (c==0.0)) {
113  printf("Info in <H1Calibrator::Calibrate> Invalid parameters !!\n");
114  return;
115  }
116 
117  if (!histo) {
118  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
119  return;
120  }
121 
122  TString opt = option;
123  opt.ToLower();
124  if ( (!opt.Contains("x")) && (!opt.Contains("y")) && (!opt.Contains("z")) ) opt = "x";
125  if (histo->GetDimension() == 1) opt = "x";
126  if ( (histo->GetDimension() == 2) && (!opt.Contains("x")) && (!opt.Contains("y")) ) opt = "x";
127 
128  //Set the calibrator
129  Set(histo);
130 
131  if (opt.Contains("x")) {
132  CalibOneAxis(a,b,c,"x");
133  }
134  if (opt.Contains("y")) {
135  CalibOneAxis(a,b,c,"y");
136  }
137  if (opt.Contains("z")) {
138  CalibOneAxis(a,b,c,"z");
139  }
140 
141  DrawResult();
142 }
143 
144 //__________________________________________________________________________________________________
145 void H1Calibrator::Calibrate(const TH1 *histo, const Double_t a, const Double_t b, const Option_t *option)
146 {
147  //Calibration with an order 1 polynome (a + bX )
148  Calibrate(histo,a,b,0.0,option);
149 }
150 
151 //__________________________________________________________________________________________________
152 void H1Calibrator::Calibrate(const TH2 *histo, const TF1 *functionX, const TF1 *functionY)
153 {
154  if (!histo) {
155  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
156  return;
157  }
158  if (!functionX) {
159  printf("Info in <H1Calibrator::Calibrate> functionX NOT defined !!\n");
160  return;
161  }
162  if (!functionY) {
163  printf("Info in <H1Calibrator::Calibrate> functionY NOT defined !!\n");
164  return;
165  }
166 
167  //Set the calibrator
168  Set(histo);
169 
170  CalibOneAxis(functionX,"x");
171  CalibOneAxis(functionY,"y");
172 
173  DrawResult();
174 }
175 
176 //__________________________________________________________________________________________________
177 void H1Calibrator::Calibrate(const TH2 *histo, const Double_t ax, const Double_t bx, const Double_t cx,
178  const Double_t ay, const Double_t by, const Double_t cy)
179 {
180  if (!histo) {
181  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
182  return;
183  }
184 
185  //Set the calibrator
186  Set(histo);
187 
188  CalibOneAxis(ax,bx,cx,"x");
189  CalibOneAxis(ay,by,cy,"y");
190 
191  DrawResult();
192 }
193 
194 //__________________________________________________________________________________________________
195 void H1Calibrator::Calibrate(const TH2 *histo, const Double_t ax, const Double_t bx, const Double_t ay, const Double_t by)
196 {
197  Calibrate(histo,ax,bx,0.0,ay,by,0.0);
198 }
199 
200 //__________________________________________________________________________________________________
201 void H1Calibrator::Calibrate(const TH3 *histo, const TF1 *functionX, const TF1 *functionY, const TF1 *functionZ)
202 {
203  if (!histo) {
204  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
205  return;
206  }
207  if (!functionX) {
208  printf("Info in <H1Calibrator::Calibrate> functionX NOT defined !!\n");
209  return;
210  }
211  if (!functionY) {
212  printf("Info in <H1Calibrator::Calibrate> functionY NOT defined !!\n");
213  return;
214  }
215  if (!functionZ) {
216  printf("Info in <H1Calibrator::Calibrate> functionZ NOT defined !!\n");
217  return;
218  }
219 
220  //Set the calibrator
221  Set(histo);
222 
223  CalibOneAxis(functionX,"x");
224  CalibOneAxis(functionY,"y");
225  CalibOneAxis(functionZ,"z");
226 
227  DrawResult();
228 }
229 
230 //__________________________________________________________________________________________________
231 void H1Calibrator::Calibrate(const TH3 *histo, const TF1 *function1, const TF1 *function2, const Option_t *option)
232 {
233  if (!histo) {
234  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
235  return;
236  }
237  if (!function1) {
238  printf("Info in <H1Calibrator::Calibrate> function1 NOT defined !!\n");
239  return;
240  }
241  if (!function2) {
242  printf("Info in <H1Calibrator::Calibrate> function2 NOT defined !!\n");
243  return;
244  }
245 
246  TString opt = option;
247  opt.ToLower();
248  if ( (!opt.Contains("xy")) && (!opt.Contains("xz")) && (!opt.Contains("yz")) ) opt = "xy";
249 
250  //Set the calibrator
251  Set(histo);
252 
253  if (!opt.CompareTo("xy")) {
254  CalibOneAxis(function1,"x");
255  CalibOneAxis(function2,"y");
256  }
257  if (!opt.CompareTo("xz")) {
258  CalibOneAxis(function1,"x");
259  CalibOneAxis(function2,"z");
260  }
261  if (!opt.CompareTo("yz")) {
262  CalibOneAxis(function1,"y");
263  CalibOneAxis(function2,"z");
264  }
265 
266  DrawResult();
267 }
268 
269 //__________________________________________________________________________________________________
270 void H1Calibrator::Calibrate(const TH3 *histo, const Double_t ax, const Double_t bx, const Double_t cx,
271  const Double_t ay, const Double_t by, const Double_t cy, const Double_t az,
272  const Double_t bz, const Double_t cz)
273 {
274  if (!histo) {
275  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
276  return;
277  }
278 
279  //Set the calibrator
280  Set(histo);
281 
282  CalibOneAxis(ax,bx,cx,"x");
283  CalibOneAxis(ay,by,cy,"y");
284  CalibOneAxis(az,bz,cz,"z");
285 
286  DrawResult();
287 }
288 
289 //__________________________________________________________________________________________________
290 void H1Calibrator::Calibrate(const TH3 *histo, const Double_t ax, const Double_t bx, const Double_t ay,
291  const Double_t by, const Double_t az, const Double_t bz)
292 {
293  Calibrate(histo,ax,bx,0.0,ay,by,0.0,az,bz,0.0);
294 }
295 
296 //__________________________________________________________________________________________________
297 void H1Calibrator::Calibrate(const TH3 *histo, const Double_t a1, const Double_t b1, const Double_t c1,
298  const Double_t a2, const Double_t b2, const Double_t c2, const Option_t *option)
299 {
300  if (!histo) {
301  printf("Info in <H1Calibrator::Calibrate> histo NOT defined !!\n");
302  return;
303  }
304 
305  TString opt = option;
306  opt.ToLower();
307  if ( (!opt.Contains("xy")) && (!opt.Contains("xz")) && (!opt.Contains("yz")) ) opt = "xy";
308 
309  //Set the calibrator
310  Set(histo);
311 
312  if (!opt.CompareTo("xy")) {
313  CalibOneAxis(a1,b1,c1,"x");
314  CalibOneAxis(a2,b2,c2,"y");
315  }
316  if (!opt.CompareTo("xz")) {
317  CalibOneAxis(a1,b1,c1,"x");
318  CalibOneAxis(a2,b2,c2,"z");
319  }
320  if (!opt.CompareTo("yz")) {
321  CalibOneAxis(a1,b1,c1,"y");
322  CalibOneAxis(a2,b2,c2,"z");
323  }
324 
325  DrawResult();
326 }
327 
328 //__________________________________________________________________________________________________
329 void H1Calibrator::Calibrate(const TH3 *histo, const Double_t a1, const Double_t b1,
330  const Double_t a2, const Double_t b2, const Option_t *option)
331 {
332  Calibrate(histo,a1,b1,0.0,a2,b2,0.0,option);
333 }
334 
335 //__________________________________________________________________________________________________
336 TH1* H1Calibrator::GetHistoCal() const {
337  if (!fHistoCal) {
338  printf("Info in <H1Calibrator::GetHistoCal> fHistoCal not defined, return 0 !!\n");
339  return 0;
340  }
341  return fHistoCal;
342 }
343 
344 //__________________________________________________________________________________________________
345 void H1Calibrator::DrawResult(const Option_t *opt) const
346 { // draw the calibrated histogram in the current pad with same range of the Raw Histogram
347  // set the range of the calibrated histogram
348  fHistoCal->GetXaxis()->SetRange(fXRawAxis->GetFirst(),fXRawAxis->GetLast());
349  fHistoCal->GetYaxis()->SetRange(fYRawAxis->GetFirst(),fYRawAxis->GetLast());
350  fHistoCal->GetZaxis()->SetRange(fZRawAxis->GetFirst(),fZRawAxis->GetLast());
351  fHistoCal->Draw(opt);
352 }
353 
354 //__________________________________________________________________________________________________
355 void H1Calibrator::Set(const TH1 *histo)
356 { // Set the H1Calibrator with a new raw histogram - the old histogram and the 3 axis are deleted
357 
358  if (!histo) {
359  printf("Info in <H1Calibrator::Set> histo NOT defined !!\n");
360  return;
361  }
362 
363  // deleting old calibrated histogram and old axis
364  if (fHistoCal != 0) {delete fHistoCal; fHistoCal=0;}
365  if (fXRawAxis != 0) {delete fXRawAxis; fXRawAxis=0;}
366  if (fYRawAxis != 0) {delete fYRawAxis; fYRawAxis=0;}
367  if (fZRawAxis != 0) {delete fZRawAxis; fZRawAxis=0;}
368 
369  //Clone the histogram and the axis
370  fHistoCal = (TH1 *)histo->Clone("HistoCal");
371  fXRawAxis = (TAxis *)histo->GetXaxis()->Clone("XAxisCal");
372  fYRawAxis = (TAxis *)histo->GetYaxis()->Clone("YAxisCal");
373  fZRawAxis = (TAxis *)histo->GetZaxis()->Clone("ZAxisCal");
374  fHistoCal->SetTitle("Calibrated Histogram");
375 }
376 
377 //__________________________________________________________________________________________________
378 void H1Calibrator::CalibOneAxis(const TF1 *function, const Option_t *option)
379 {
380  if (!function) {
381  printf("Info in <H1Calibrator::CalibOneAxis> function NOT defined !!\n");
382  return;
383  }
384 
385  TString opt = option;
386  opt.ToLower();
387 
388  Double_t xmin = 0.0;
389  Double_t xmax = 0.0;
390 
391  if (!opt.CompareTo("x")) {
392  xmin = function->Eval(fXRawAxis->GetXmin());
393  xmax = function->Eval(fXRawAxis->GetXmax());
394  if (!fXRawAxis->IsVariableBinSize()) { //fixed
395  if (xmin < xmax) {
396  fHistoCal->GetXaxis()->SetLimits(xmin,xmax);
397  } else {
398  Invert();
399  fHistoCal->GetXaxis()->SetLimits(xmax,xmin);
400  }
401  } else { // variable bin size
402  Double_t *xbins = new Double_t[fXRawAxis->GetNbins()+1];
403  const TArrayD *array = fXRawAxis->GetXbins();
404  xmin = function->Eval(fXRawAxis->GetXmin());
405  if (xmin < xmax) {
406  for (Int_t i = 0 ; i < fXRawAxis->GetNbins()+1; i++) {
407  xbins[i] = function->Eval((*array)[i]);
408  }
409  } else {
410  Invert();
411  for (Int_t i=0 ; i<fXRawAxis->GetNbins()+1; i++) {
412  xbins[i] = function->Eval((*array)[fXRawAxis->GetNbins()-i]);
413  }
414  }
415  fHistoCal->GetXaxis()->Set(fXRawAxis->GetNbins(),xbins);
416  delete [] xbins;
417  }
418  }
419 
420  if (!opt.CompareTo("y")) {
421  xmin = function->Eval(fYRawAxis->GetXmin());
422  xmax = function->Eval(fYRawAxis->GetXmax());
423  if (!fYRawAxis->IsVariableBinSize()) { //fixed
424  if (xmin < xmax) {
425  fHistoCal->GetYaxis()->SetLimits(xmin,xmax);
426  } else {
427  Invert();
428  fHistoCal->GetYaxis()->SetLimits(xmax,xmin);
429  }
430  } else { // variable bin size
431  Double_t *xbins = new Double_t[fYRawAxis->GetNbins()+1];
432  const TArrayD *array = fYRawAxis->GetXbins();
433  if (xmin < xmax) {
434  for (Int_t i = 0 ; i < fYRawAxis->GetNbins()+1; i++) {
435  xbins[i] = function->Eval((*array)[i]);
436  }
437  } else {
438  Invert();
439  for (Int_t i = 0 ; i < fYRawAxis->GetNbins()+1; i++) {
440  xbins[i] = function->Eval((*array)[fYRawAxis->GetNbins()-i]);
441  }
442  }
443  fHistoCal->GetYaxis()->Set(fYRawAxis->GetNbins(),xbins);
444  delete [] xbins;
445  }
446  }
447 
448  if (!opt.CompareTo("z")) {
449  xmin = function->Eval(fZRawAxis->GetXmin());
450  xmax = function->Eval(fZRawAxis->GetXmax());
451  if (!fZRawAxis->IsVariableBinSize()) { //fixed
452  if (xmin < xmax) {
453  fHistoCal->GetZaxis()->SetLimits(xmin,xmax);
454  } else {
455  Invert();
456  fHistoCal->GetZaxis()->SetLimits(xmax,xmin);
457  }
458  } else { // variable bin size
459  Double_t *xbins = new Double_t[fZRawAxis->GetNbins()+1];
460  const TArrayD *array = fZRawAxis->GetXbins();
461  if (xmin < xmax) {
462  for (Int_t i=0 ; i<fZRawAxis->GetNbins()+1; i++) {
463  xbins[i] = function->Eval((*array)[i]);
464  }
465  } else {
466  Invert();
467  for (Int_t i = 0 ; i < fZRawAxis->GetNbins()+1; i++) {
468  xbins[i] = function->Eval((*array)[fZRawAxis->GetNbins()-i]);
469  }
470  }
471  fHistoCal->GetZaxis()->Set(fZRawAxis->GetNbins(),xbins);
472  delete [] xbins;
473  }
474  }
475 }
476 
477 //__________________________________________________________________________________________________
478 void H1Calibrator::CalibOneAxis(const Double_t a, const Double_t b, const Double_t c, const Option_t *option)
479 {
480  if ((b==0.0) && (c==0.0)) {
481  printf("Info in <H1Calibrator::CalibOneAxis> invalid parameters !!\n");
482  return;
483  }
484 
485  TF1 *function = new TF1("ftemp","pol2(0)");
486  function->SetParameter(0,a);
487  function->SetParameter(1,b);
488  function->SetParameter(2,c);
489 
490  CalibOneAxis(function,option);
491 
492  if (function != 0) delete function;
493 
494 }
495 
496 //__________________________________________________________________________________________________
497 void H1Calibrator::CalibOneAxis(const Double_t a, const Double_t b, const Option_t *option)
498 {
499  if ( b==0.0 ) {
500  printf("Info in <H1Calibrator::CalibOneAxis> invalid parameters !!\n");
501  return;
502  }
503 
504  TF1 *function = new TF1("ftemp","pol1(0)");
505  function->SetParameter(0,a);
506  function->SetParameter(1,b);
507 
508  CalibOneAxis(function,option);
509 
510  if (function != 0) delete function;
511 
512 }
513 
514 //__________________________________________________________________________________________________
515 void H1Calibrator::Invert()
516 {
517  if (fHistoCal == 0) {
518  printf("Info in <H1Calibrator::Invert> fHistoCal not defined !!\n");
519  return;
520  }
521  TH1 *HistoCalTemp = (TH1 *)fHistoCal->Clone("HistoCalTemp");
522  if (HistoCalTemp == 0) {
523  printf("Info in <H1Calibrator::Invert> HistoCalTemp not defined !!\n");
524  return;
525  }
526  if (fHistoCal->GetDimension() > 1) {
527  printf("Info in <H1Calibrator::Invert> Method not already defined for Histogram Dimension >= 2 !!\n");
528  return;
529  }
530 
531  fHistoCal->Reset();
532  for (Int_t i=0 ; i <= fHistoCal->GetNbinsX()+1 ; i++) {
533  fHistoCal->SetBinContent(i,HistoCalTemp->GetBinContent(fHistoCal->GetNbinsX()+1-i));
534  }
535  if (HistoCalTemp != 0) delete HistoCalTemp;
536 }
537 
538 
539 
540 
541 
542 
543 
544 
b
TBrowser * b
Definition: adf/Macros/Deprecated/Load.C:128
printf
printf("******************************************************************** \n")
Gw::H1Calibrator::~H1Calibrator
virtual ~H1Calibrator()
Definition: H1Calibrator.cpp:63
ClassImp
ClassImp(H1Calibrator)
Gw::H1Calibrator::Calibrate
void Calibrate(const TH1 *, const TF1 *, const Option_t *AxisOption="x")
Definition: H1Calibrator.cpp:72
histo
TH1F * histo[MaxValue]
Definition: ReadDaqAlone.C:31
Gw::H1Calibrator::DrawResult
void DrawResult(const Option_t *opt="") const
draw the calibrated histogram in the current pad with same range of the raw histogram ...
Definition: H1Calibrator.cpp:345
Gw::H1Calibrator::GetHistoCal
TH1 * GetHistoCal() const
return the calibrated histogram
Definition: H1Calibrator.cpp:336
Gw::H1Calibrator
H1Calibrator is a service class in order to calibrate a TH1 The raw histogram must never be changed T...
Definition: H1Calibrator.h:58
H1Calibrator.h
header file for the calibration facility