/[PAMELA software]/PamCut/CollectionActions/Histo2DActions/Histo2DAction/Histo2DAction.h
ViewVC logotype

Contents of /PamCut/CollectionActions/Histo2DActions/Histo2DAction/Histo2DAction.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.2 - (show annotations) (download)
Tue Oct 27 10:17:48 2009 UTC (15 years, 1 month ago) by pam-fi
Branch: MAIN
Changes since 1.1: +23 -12 lines
File MIME type: text/plain
Updates in the Doxygen documentation.

1 /*
2 * Histo2DAction.h
3 *
4 * Created on: 01/set/2009
5 * Author: Nicola Mori
6 */
7
8 /*! @file Histo2DAction.h The Histo2DAction class declaration file. */
9
10 #ifndef HISTO2DACTION_H_
11 #define HISTO2DACTION_H_
12
13 #include "../../CollectionAction/CollectionAction.h"
14
15 #include <TH2I.h>
16 #include <TH2F.h>
17 #include <TH2D.h>
18
19 /*! @brief An abstract action that fills a 2-dimensional histogram.
20 *
21 * This abstract class provides all the common methods needed to produce 2-dimensional histograms,
22 * like reading the bins from a file, setting up the root histogram and so on. It handles histograms
23 * both as a ROOT histogram and as a matrix. The OnGood() method has to be implemented in children
24 * classes, with the specific physical variables to fill the histogram. Eg., a dE/dx vs. rigidity
25 * histogram class derived from Histo2DAction could have this OnGood() implementation:
26 *
27 * Fill(event->GetTrack(0)->GetTrkTrack()->GetRigidity(), event->GetTrack(0)->GetTrkTrack()->GetDEDX());
28 *
29 * Fill() automatically fills the ROOT and the matrix histogram.
30 * The template argument HistoType is the type of histogram (Int_t, Double_t,...). In current
31 * implementation, ROOT histograms are only supported for Int_t, Float_t and Double_t. For all the
32 * other types, only the the matrix histogram with text file output will be produced. If you want to
33 * implement some new types of ROOT histograms, please specialize the template definitions of _CreateHisto()
34 * and ~Histo2DAction() (see the specializations in the .cpp file).
35 *
36 * */
37 template<class HistoType>
38 class Histo2DAction: public CollectionAction {
39 public:
40
41 /*! @brief Constructor.
42 *
43 * The histogram binning is defined by the content of the bins vector.
44 *
45 * @param actionName The action's name.
46 * @param title The ROOT histogram title.
47 * @param outFileBase The file base name for the histogram output. If "", no output file will be produced.
48 * @param mode The mode of ROOT file creation (see documentation of TFile constructor
49 * in ROOT's reference guide).
50 * @param outRoot If true, an output ROOT file named outFileBase + ".root" will be produced.
51 * @param outText If true, an output text file named outFileBase + ".txt" will be produced. It will overwrite an
52 * eventually existing file with the same name.
53 */
54 Histo2DAction(const char *actionName, TString title = "", TString outFileBase = "", TString mode = "UPDATE",
55 bool outRoot = false, bool outText = false);
56
57 /*! @brief Destructor
58 *
59 * This has to be specialized for each kind of ROOT histogram.
60 *
61 * */
62 ~Histo2DAction();
63
64 /*! @brief Sets the X axis using a binning read from a a vector.
65 *
66 * @param label The axis' label.
67 * @param bins A vector containing the histogram binning (with sign, eventually; from lowest to highest).
68 */
69 void SetXAxis(TString label, vector<float> &bins);
70 /*! @brief Sets the X axis reading the binning from a file.
71 *
72 * @param label The axis' label.
73 * @param binsFile The file containing the bins (with sign, eventually; from lowest to highest).
74 */
75 void SetXAxis(TString label, TString binsFile);
76 /*! @brief Sets the X axis specifying the binning parameters.
77 *
78 * @param label The axis' label.
79 * @param nBins The number of bins.
80 * @param min The lower limit of the axis.
81 * @param max The upper limit of the axis.
82 * @param logBinning If true, bins will be logarithmically spaced.
83 */
84 void SetXAxis(TString label, unsigned int nBins, float min, float max, bool logBinning = false);
85
86 /*! @brief Sets the Y axis using a binning read from a a vector.
87 *
88 * @param label The axis' label.
89 * @param bins A vector containing the histogram binning (with sign, eventually; from lowest to highest).
90 */
91 void SetYAxis(TString label, vector<float> &bins);
92 /*! @brief Sets the Y axis reading the binning from a file.
93 *
94 * @param label The axis' label.
95 * @param binsFile The file containing the bins (with sign, eventually; from lowest to highest).
96 */
97 void SetYAxis(TString label, TString binsFile);
98
99 /*! @brief Sets the Y axis specifying the binning parameters.
100 *
101 * @param label The axis' label.
102 * @param nBins The number of bins.
103 * @param min The lower limit of the axis.
104 * @param max The upper limit of the axis.
105 * @param logBinning If true, bins will be logarithmically spaced.
106 */
107 void SetYAxis(TString label, unsigned int nBins, float min, float max, bool logBinning = false);
108
109 /*! @brief Sets the ROOT histogram's title.
110 *
111 * @param title The histogram title as it will appear on the histogram itself.
112 */
113 void SetTitle(TString &title) {
114 _title = title;
115 }
116 /*! @brief Sets up the histogram
117 *
118 * This routine effectively prepares the histogram, after the desired parameters has been set by #SetXAxis() and #SetYAxis().
119 *
120 * @param events Pointer to PamLevel2 events (unused).
121 */
122 void Setup(PamLevel2 *events) {
123 CollectionAction::Setup(events);
124 _InitHistos();
125
126 }
127
128 /*! @ brief Post-analysis tasks.
129 *
130 * This method writes the output on a ROOT file if the proper option was set in the constructor.
131 */
132 void Finalize();
133
134 /*! @brief Returns the histogram.
135 *
136 * Element [i][j] is the value of the i-th bin on the Y axis and on the j-th bin on the X axis. This sort
137 * of index inversion is due to logically maintain the Y values on the vertical axis of the matrix, which
138 * is indexed by the row number (which conventionally is the first index of a matrix element).
139 *
140 * @return A reference to a matrix containing the values of the bins of the histogram.
141 */
142 SimpleMatrix<HistoType> &GetHisto() {
143 return _histo;
144 }
145
146 /*! Fills the ROOT and the vector histogram.
147 *
148 * @param value The value of the X coordinate associated to the event.
149 * @param value The value of the Y coordinate associated to the event.
150 * @param weight The weight which will be applied to the event.
151 */
152 void Fill(double xValue, double yValue, double weight = 1.);
153
154 /*! @brief Gets the X overflow histogram.
155 *
156 * The returned vector contains the histogram filled with the X overflow events, eg.,
157 * those events whose X variable falls above the upper limit of the X axis. This histogram
158 * is binned like the Y axis.
159 *
160 * @return The X overflow histogram.
161 */
162 vector<HistoType>& GetXOverflow() {
163 return _xOverflow;
164 }
165
166 /*! @brief Gets the X underflow histogram.
167 *
168 * The returned vector contains the histogram filled with the X underflow events, eg.,
169 * those events whose X variable falls below the lower limit of the X axis. This histogram
170 * is binned like the Y axis.
171 *
172 * @return The X underflow histogram.
173 */
174 vector<HistoType>& GetXUnderflow() {
175 return _xUnderflow;
176 }
177
178 /*! @brief Gets the Y overflow histogram.
179 *
180 * The returned vector contains the histogram filled with the Y overflow events, eg.,
181 * those events whose Y variable falls above the upper limit of the Y axis. This histogram
182 * is binned like the X axis.
183 *
184 * @return The Y overflow histogram.
185 */
186 vector<HistoType>& GetYOverflow() {
187 return _yOverflow;
188 }
189
190 /*! @brief Gets the Y underflow histogram.
191 *
192 * The returned vector contains the histogram filled with the Y underflow events, eg.,
193 * those events whose Y variable falls below the lower limit of the Y axis. This histogram
194 * is binned like the X axis.
195 *
196 * @return The Y underflow histogram.
197 */
198 vector<HistoType>& GetYUnderflow() {
199 return _yUnderflow;
200 }
201
202 /*! @brief Gets the X-Y overflow histogram.
203 *
204 * This single-bin histogram is built with events which fall above the lower limits of both the
205 * X and Y axis.
206 *
207 * @return The X-Y overflow histogram.
208 */
209 HistoType GetXOverYOverflow() {
210 return _xOverYOverflow;
211 }
212
213 /*! @brief Gets the X-Y underflow histogram.
214 *
215 * This single-bin histogram is built with events fall below the lower limits of both the X and Y axis.
216 *
217 * @return The X-Y overflow histogram.
218 */
219 HistoType GetXUnderYUnderflow() {
220 return _xUnderYUnderflow;
221 }
222 /*! @brief Gets the (X over - Y under)flow histogram.
223 *
224 * This single-bin histogram is built with events fall above the lower limit of the Xaxis and below
225 * that of the Y axis.
226 *
227 * @return The (X over -Y under)flow histogram.
228 */
229 HistoType GetXOverYUnderflow() {
230 return _xOverYUnderflow;
231 }
232
233 /*! @brief Gets the (X under - Y over)flow histogram.
234 *
235 * This single-bin histogram is built with events fall below the lower limit of the X axis and above
236 * that of the Y axis.
237 *
238 * @return The (X under -Y over)flow histogram.
239 */
240 HistoType GetXUnderYOverflow() {
241 return _xUnderYOverflow;
242 }
243
244 protected:
245
246 /*! @brief The vector containing the limits of the X bins(from lower to higher). */
247 std::vector<float> _xBins;
248 /*! @brief The vector containing the limits of the Y bins(from lower to higher). */
249 std::vector<float> _yBins;
250 /*! @brief A matrix containing the value of the histogram for each X-Y bin. */
251 SimpleMatrix<HistoType> _histo;
252 /*! @brief The ROOT histogram. */
253 TH2 *_rootHisto;
254
255 private:
256
257 vector<HistoType> _xUnderflow, _xOverflow, _yUnderflow, _yOverflow;
258 HistoType _xUnderYUnderflow, _xOverYOverflow, _xUnderYOverflow, _xOverYUnderflow;
259 TString _outFileBase;
260 TString _mode;
261 TString _title, _xLabel, _yLabel;
262 bool _outRoot;
263 bool _outText;
264
265 void _CreateHisto();
266 void _InitHistos();
267 };
268
269 template<class HistoType>
270 void Histo2DAction<HistoType>::_CreateHisto() {
271
272 // No ROOT histogram for generic type; see template specializations in .cpp file.
273 _rootHisto = NULL;
274 }
275
276 template<class HistoType>
277 void Histo2DAction<HistoType>::_InitHistos() {
278
279 _CreateHisto();
280
281 if (_xBins.size() < 2) {
282 _xBins.resize(2);
283 _xBins[0] = 0.;
284 _xBins[1] = 1.;
285 }
286
287 if (_yBins.size() < 2) {
288 _yBins.resize(2);
289 _yBins[0] = 0.;
290 _yBins[1] = 1.;
291 }
292
293 if (_rootHisto) {
294 Double_t *auxXArray = new Double_t[_xBins.size()];
295
296 for (unsigned int i = 0; i < _xBins.size(); i++) {
297 auxXArray[i] = _xBins[i];
298 }
299
300 Double_t *auxYArray = new Double_t[_yBins.size()];
301
302 for (unsigned int i = 0; i < _yBins.size(); i++) {
303 auxYArray[i] = _yBins[i];
304 }
305
306 _rootHisto->SetBins(_xBins.size() - 1, auxXArray, _yBins.size() - 1, auxYArray);
307 _rootHisto->SetName(GetName());
308 _rootHisto->SetTitle(_title);
309 _rootHisto->SetXTitle(_xLabel);
310 _rootHisto->SetYTitle(_yLabel);
311
312 delete[] auxXArray;
313 delete[] auxYArray;
314 }
315
316 /* The row index (first) corresponds to the position on the vertical (Y) axis. */
317 _histo.Resize(_yBins.size() - 1, _xBins.size() - 1);
318
319 _xUnderflow.resize(_yBins.size());
320 _xOverflow.resize(_yBins.size());
321 _yUnderflow.resize(_xBins.size());
322 _yOverflow.resize(_xBins.size());
323
324 }
325
326 template<class HistoType>
327 Histo2DAction<HistoType>::~Histo2DAction() {
328 }
329
330 template<class HistoType>
331 Histo2DAction<HistoType>::Histo2DAction(const char *actionName, TString title, TString outFileBase, TString mode,
332 bool outRoot, bool outText) :
333 CollectionAction(actionName), _xBins(0), _yBins(0), _histo(0, 0), _rootHisto(NULL), _xUnderflow(0), _xOverflow(0),
334 _yUnderflow(0), _yOverflow(0), _outFileBase(outFileBase), _mode(mode), _title(title), _xLabel(""), _yLabel(""),
335 _outRoot(outRoot), _outText(outText) {
336
337 }
338
339 template<class HistoType>
340 void Histo2DAction<HistoType>::SetXAxis(TString label, vector<float> &bins) {
341
342 _xBins = bins;
343 _xLabel = label;
344
345 }
346
347 template<class HistoType>
348 void Histo2DAction<HistoType>::SetYAxis(TString label, vector<float> &bins) {
349
350 _yBins = bins;
351 _yLabel = label;
352
353 }
354
355 template<class HistoType>
356 void Histo2DAction<HistoType>::SetXAxis(TString label, TString binsFile) {
357
358 // Reading the bins from file
359 ifstream binListFile;
360 binListFile.open(binsFile);
361
362 TString auxString;
363 _xBins.resize(0);
364 while (!binListFile.eof()) {
365 binListFile >> auxString;
366 if (auxString != "") {
367 _xBins.push_back(auxString.Atof());
368 }
369 }
370 binListFile.close();
371
372 _xLabel = label;
373
374 }
375
376 template<class HistoType>
377 void Histo2DAction<HistoType>::SetYAxis(TString label, TString binsFile) {
378
379 // Reading the bins from file
380 ifstream binListFile;
381 binListFile.open(binsFile);
382
383 TString auxString;
384 _yBins.resize(0);
385 while (!binListFile.eof()) {
386 binListFile >> auxString;
387 if (auxString != "") {
388 _yBins.push_back(auxString.Atof());
389 }
390 }
391 binListFile.close();
392
393 _yLabel = label;
394
395 }
396
397 template<class HistoType>
398 void Histo2DAction<HistoType>::SetXAxis(TString label, unsigned int nBins, float min, float max, bool logBinning) {
399
400 _xBins.resize(nBins + 1);
401
402 if (!logBinning || (logBinning && min <= 0.)) {
403
404 #ifdef DEBUGPAMCUT
405 if (logbinning && rigMin <= 0.)
406 cout << "Warning: logarithmic binning was chosen for X axis but min <= 0. Using linear binning."
407 #endif
408
409 float step = (max - min) / nBins;
410 for (unsigned int i = 0; i < nBins + 1; i++) {
411 _xBins[i] = min + i * step;
412 }
413
414 }
415 else {
416
417 double maxExp = log10(max / min);
418 for (unsigned int i = 0; i < nBins + 1; i++) {
419 _xBins[i] = min * pow(10., (double) i / ((double) nBins) * maxExp);
420 }
421
422 }
423
424 _xLabel = label;
425 }
426
427 template<class HistoType>
428 void Histo2DAction<HistoType>::SetYAxis(TString label, unsigned int nBins, float min, float max, bool logBinning) {
429
430 _yBins.resize(nBins + 1);
431
432 if (!logBinning || (logBinning && min <= 0.)) {
433
434 #ifdef DEBUGPAMCUT
435 if (logbinning && rigMin <= 0.)
436 cout << "Warning: logarithmic binning was chosen for Y axis but min <= 0. Using linear binning."
437 #endif
438
439 float step = (max - min) / nBins;
440 for (unsigned int i = 0; i < nBins + 1; i++) {
441 _yBins[i] = min + i * step;
442 }
443
444 }
445 else {
446
447 double maxExp = log10(max / min);
448 for (unsigned int i = 0; i < nBins + 1; i++) {
449 _yBins[i] = min * pow(10., (double) i / ((double) nBins) * maxExp);
450 }
451
452 }
453
454 _yLabel = label;
455 }
456
457 template<class HistoType>
458 inline void Histo2DAction<HistoType>::Fill(double xValue, double yValue, double weight) {
459
460 _rootHisto->Fill(xValue, yValue, weight);
461
462 int xBin = 0, yBin = 0;
463
464 bool UOflow = false;
465
466 if (xValue < _xBins.front()) {
467 xBin = -1;
468 UOflow = true;
469 }
470
471 if (xValue > _xBins.back()) {
472 xBin = _xBins.size();
473 UOflow = true;
474 }
475
476 if (yValue < _yBins.front()) {
477 yBin = -1;
478 UOflow = true;
479 }
480
481 if (yValue > _yBins.back()) {
482 yBin = _yBins.size();
483 UOflow = true;
484 }
485
486 if (xBin == 0) {
487
488 while (xValue >= _xBins[xBin]) {
489 xBin++;
490 }
491 xBin--;
492 }
493
494 if (yBin == 0) {
495
496 while (yValue >= _yBins[yBin]) {
497 yBin++;
498 }
499 yBin--;
500 }
501
502 if (!UOflow) {
503 // No under or over flow.
504 _histo[yBin][xBin] += (HistoType) weight;
505 return;
506 }
507 else {
508 // An under or over flow has happened.
509 if (xBin == -1) {
510 if (yBin == -1) {
511 _xUnderYUnderflow += (HistoType) weight;
512 return;
513 }
514 else {
515 if (yBin == (int) _yBins.size()) {
516 _xUnderYOverflow += (HistoType) weight;
517 return;
518 }
519 else {
520 _xUnderflow[yBin] += (HistoType) weight;
521 return;
522 }
523 }
524 }
525
526 if (xBin == (int) _xBins.size()) {
527 if (yBin == -1) {
528 _xOverYUnderflow += (HistoType) weight;
529 return;
530 }
531 else {
532 if (yBin == (int) _yBins.size()) {
533 _xOverYOverflow += (HistoType) weight;
534 return;
535 }
536 else {
537 _xOverflow[yBin] += (HistoType) weight;
538 return;
539 }
540 }
541 }
542
543 if (yBin == -1) {
544 _yUnderflow[xBin] += (HistoType) weight;
545 return;
546 }
547
548 if (yBin == (int) _yBins.size()) {
549 _yOverflow[xBin] += (HistoType) weight;
550 return;
551 }
552
553 }
554 }
555
556 template<class HistoType>
557 void Histo2DAction<HistoType>::Finalize() {
558
559 if (_outFileBase != "") {
560 // Write the ROOT file
561 if (_rootHisto && _outRoot) {
562 TFile outRootFile((_outFileBase + ".root"), _mode);
563 outRootFile.cd();
564 _rootHisto->Write();
565 outRootFile.Close();
566 }
567
568 //Write the text file
569 if (_outText) {
570 ofstream outTextFile((_outFileBase + ".txt").Data(), ios_base::out);
571 for (unsigned int i = 0; i < _histo.GetNRows(); i++) {
572 for (unsigned int j = 0; j < _histo.GetNCols(); j++) {
573 outTextFile << _histo[i][j] << " ";
574 }
575 outTextFile << "\n";
576 }
577 outTextFile << endl;
578 outTextFile.close();
579 }
580 }
581
582 }
583 #endif /* HISTO2DACTION_H_ */

  ViewVC Help
Powered by ViewVC 1.1.23