Heavy ion Analysis Libriares
Loading...
Searching...
No Matches
StdHist.cxx
1/*
2 * Hal::StdHist.cxx
3 *
4 * Created on: 24 lut 2019
5 * Author: Daniel Wielanek
6 * E-mail: daniel.wielanek@gmail.com
7 * Warsaw University of Technology, Faculty of Physics
8 */
9
10#include "StdHist.h"
11#include "StdString.h"
12
13
14#include <TArray.h>
15#include <TAxis.h>
16#include <TCollection.h>
17#include <TColor.h>
18#include <TDirectory.h>
19#include <TGaxis.h>
20#include <TH1.h>
21#include <TH2.h>
22#include <TH3.h>
23#include <TList.h>
24#include <TMath.h>
25#include <TMatrixDfwd.h>
26#include <TMatrixT.h>
27#include <TNamed.h>
28#include <TObjArray.h>
29#include <TPaletteAxis.h>
30#include <TROOT.h>
31#include <TRandom.h>
32#include <TStyle.h>
33#include <TVirtualPad.h>
34#include <iostream>
35#include <utility>
36
37#include "Cout.h"
38#include "Splines.h"
39#include "Std.h"
40
41NamespaceImp(Hal::Std)
42
43 namespace Hal {
44 namespace Std {
45 void RemoveNan(TH1* h, Double_t fill, Double_t fill_e) {
46 if (h->InheritsFrom("TH3")) {
47 for (int i = 0; i <= h->GetNbinsX() + 1; i++) {
48 for (int j = 0; j <= h->GetNbinsY() + 1; j++) {
49 for (int k = 0; k <= h->GetNbinsZ() + 1; k++) {
50 if (TMath::IsNaN(h->GetBinContent(i, j, k))) {
51 h->SetBinContent(i, j, k, fill);
52 h->SetBinError(i, j, k, fill_e);
53 }
54 }
55 }
56 }
57 } else if (h->InheritsFrom("TH2")) {
58 for (int i = 0; i <= h->GetNbinsX() + 1; i++) {
59 for (int j = 0; j <= h->GetNbinsY() + 1; j++) {
60 if (TMath::IsNaN(h->GetBinContent(i, j))) {
61 h->SetBinContent(i, j, fill);
62 h->SetBinError(i, j, fill_e);
63 }
64 }
65 }
66 } else {
67 for (int i = 0; i <= h->GetNbinsX() + 1; i++) {
68 if (TMath::IsNaN(h->GetBinContent(i))) {
69 h->SetBinContent(i, fill);
70 h->SetBinError(i, fill_e);
71 }
72 }
73 }
74 }
75
76 TH1D* GetProjection1D(const TH3* histo, Double_t min1, Double_t max1, Double_t min2, Double_t max2, Option_t* opt) {
77 TString option = opt;
78 Bool_t sumw = Hal::Std::FindParam(option, "sumw", kTRUE);
79 const TAxis *axisA, *axisB, *axisC;
80 TH1D* projection;
81 Int_t nbins = 1;
82 if (option.Contains("z")) {
83 axisA = histo->GetXaxis();
84 axisB = histo->GetYaxis();
85 axisC = histo->GetZaxis();
86 TString name = "_pz";
87 if (!option.Contains("noautoname")) { name = Form("%i_pz", anonymCounter++); }
88 if (option.Contains("bins")) {
89 projection = histo->ProjectionZ(name, min1, max1, min2, max2, "e");
90 } else {
91 projection = histo->ProjectionZ(
92 name, axisA->FindFixBin(min1), axisA->FindFixBin(max1), axisB->FindFixBin(min2), axisB->FindFixBin(max2), "o");
93 }
94 } else if (option.Contains("y")) {
95 axisA = histo->GetXaxis();
96 axisB = histo->GetZaxis();
97 axisC = histo->GetYaxis();
98 TString name = "_py";
99 if (!option.Contains("noautoname")) { name = Form("%i_py", anonymCounter++); }
100 if (option.Contains("bins")) {
101 projection = histo->ProjectionY(name, min1, max1, min2, max2, "e");
102 } else {
103 projection = histo->ProjectionY(
104 name, axisA->FindFixBin(min1), axisA->FindFixBin(max1), axisB->FindFixBin(min2), axisB->FindFixBin(max2), "o");
105 }
106 } else {
107 axisA = histo->GetYaxis();
108 axisB = histo->GetZaxis();
109 axisC = histo->GetXaxis();
110 TString name = "_px";
111 if (!option.Contains("noautoname")) { name = Form("%i_px", anonymCounter++); }
112 if (option.Contains("bins")) {
113 projection = histo->ProjectionX(name, min1, max1, min2, max2, "e");
114 } else {
115 projection = histo->ProjectionX(
116 name, axisA->FindFixBin(min1), axisA->FindFixBin(max1), axisB->FindFixBin(min2), axisB->FindFixBin(max2), "o");
117 }
118 }
119 if (option.Contains("bins")) {
120 nbins = (max1 - min1 + 1) * (max2 - min2 + 1);
121 } else {
122 nbins = axisA->FindFixBin(max1) - axisA->FindFixBin(min1) + 1;
123 nbins = nbins * (axisB->FindFixBin(max2) - axisB->FindFixBin(min2) + 1);
124 }
125 if (sumw) projection->Sumw2();
126 if (option.Contains("scale")) {
127 if (nbins != 0) projection->Scale(1.0 / nbins);
128 }
129 projection->GetXaxis()->SetTitle(axisC->GetTitle());
130 // projection->SetDirectory(0);
131 return projection;
132 }
133
134 TH1D* GetProjection1D(const TH2* histo, Double_t min, Double_t max, Option_t* opt) {
135 TString option = opt;
136 TH1D* projection;
137 Bool_t sumw = Hal::Std::FindParam(option, "sumw", kTRUE);
138 Double_t nbins = 0;
139 if (option.Contains("y")) {
140 TString name = "_py";
141 if (!option.Contains("noautoname")) { name = Form("%i_py", anonymCounter++); }
142 if (option.Contains("bins")) {
143 projection = histo->ProjectionY(name, min, max);
144 nbins = max - min;
145 } else {
146 projection = histo->ProjectionY(name, histo->GetXaxis()->FindBin(min), histo->GetXaxis()->FindBin(max));
147 nbins = histo->GetXaxis()->FindBin(max) - histo->GetXaxis()->FindBin(min);
148 }
149 projection->GetXaxis()->SetTitle(histo->GetYaxis()->GetTitle());
150 } else {
151 TString name = "_px";
152 if (!option.Contains("noautoname")) { name = Form("%i_px", anonymCounter++); }
153 if (option.Contains("bins")) {
154 projection = histo->ProjectionX(name, min, max);
155 nbins = max - min;
156 } else {
157 projection = histo->ProjectionX(name, histo->GetYaxis()->FindBin(min), histo->GetYaxis()->FindBin(max));
158 nbins = histo->GetYaxis()->FindBin(max) - histo->GetYaxis()->FindBin(min);
159 }
160 projection->GetXaxis()->SetTitle(histo->GetXaxis()->GetTitle());
161 }
162 if (sumw) projection->Sumw2();
163 if (option.Contains("scale")) {
164 if (nbins != 0) projection->Scale(1.0 / (nbins + 1));
165 }
166 projection->SetDirectory(0);
167 return projection;
168 }
169
170 TH2D* GetProjection2D(const TH3* histo, Double_t min, Double_t max, Option_t* opt) {
171 TString option = opt;
172 TString title = histo->GetTitle();
173 TH3* histo_cloned = (TH3*) histo->Clone();
174 TH2D* projection;
175 TAxis* axis;
176 Bool_t sumw = Hal::Std::FindParam(option, "sumw", kTRUE);
177 // TAxis *axisA, *axisB;
178 Int_t nbins = 0;
179 TString projection_option;
180 if (Hal::Std::FindParam(option, "yz")) {
181 axis = histo_cloned->GetXaxis();
182 projection_option = "yze";
183 } else if (Hal::Std::FindParam(option, "zy")) {
184 axis = histo_cloned->GetXaxis();
185 projection_option = "zye";
186 } else if (Hal::Std::FindParam(option, "xz")) {
187 axis = histo_cloned->GetYaxis();
188 projection_option = "xze";
189 // axisA = histo->GetXaxis();
190 // axisB = histo->GetZaxis();
191 } else if (Hal::Std::FindParam(option, "zx")) {
192 axis = histo_cloned->GetYaxis();
193 projection_option = "zxe";
194 } else if (Hal::Std::FindParam(option, "yx")) {
195 axis = histo_cloned->GetZaxis();
196 projection_option = "yxe";
197 } else {
198 axis = histo_cloned->GetZaxis();
199 projection_option = "xye";
200 // axisA = histo->GetXaxis();
201 // axisB = histo->GetYaxis();
202 }
203 if (Hal::Std::FindParam(option, "bins")) {
204 axis->SetRange(min, max);
205 nbins = 1 + max - min;
206 } else {
207 if (min == max) { // same values, this makes problem
208 max = min + axis->GetBinWidth(1) * 0.1;
209 }
210 axis->SetRangeUser(min, max);
211 nbins = 1 + axis->FindBin(max) - axis->FindBin(min);
212 }
213 projection = (TH2D*) histo_cloned->Project3D(projection_option);
214 delete histo_cloned;
215 // axis->SetRange(0, 0);
216 if (sumw) projection->Sumw2();
217 if (option.Contains("scale")) { projection->Scale(1.0 / nbins); }
218 if (!option.Contains("noautoname")) {
219 TString name = Form("%i_2dproj", anonymCounter++);
220 projection->SetName(name);
221 }
222 projection->SetDirectory(0);
223 return projection;
224 }
225
226 TH1D* SmearHistogram(TH1D* input_histogram, TH2D* smear_matrix, Option_t* opt) {
227 TString option = opt;
228 Bool_t bad_map = kFALSE;
229 if (input_histogram->GetNbinsX() != smear_matrix->GetNbinsX()) {
230 if (!option.Contains("forced")) Hal::Cout::PrintInfo("Incompatible histograms for smearing", Hal::EInfo::kWarning);
231 bad_map = kTRUE;
232 }
233 if (input_histogram->GetNbinsX() != smear_matrix->GetNbinsY()) {
234 if (!option.Contains("forced")) Hal::Cout::PrintInfo("Incompatible histograms for smearing", Hal::EInfo::kWarning);
235 bad_map = kTRUE;
236 }
237 if (option.Contains("forced") && bad_map) {
238 Hal::Spline2D* m_map = new Hal::Spline2D(smear_matrix);
239 m_map->Refit();
240 TH2D* smear_new = new TH2D("smear_temp",
241 "smear_temp",
242 input_histogram->GetNbinsX(),
243 input_histogram->GetXaxis()->GetXmin(),
244 input_histogram->GetXaxis()->GetXmax(),
245 input_histogram->GetNbinsX(),
246 input_histogram->GetXaxis()->GetXmin(),
247 input_histogram->GetXaxis()->GetXmax());
248 for (int i = 0; i <= input_histogram->GetNbinsX() + 1; i++) {
249 Double_t x = smear_new->GetXaxis()->GetBinCenter(i);
250 for (int j = 0; j <= input_histogram->GetNbinsX() + 1; j++) {
251 Double_t y = smear_new->GetYaxis()->GetBinCenter(j);
252 if (m_map->Eval(x, y) >= 0) {
253 smear_new->SetBinContent(i, j, m_map->Eval(x, y));
254 smear_new->SetBinError(i, j, m_map->GetError(x, y));
255 }
256 }
257 }
258 smear_matrix = smear_new;
259 delete m_map;
260 } else {
261 if (bad_map) return NULL;
262 }
263
264 TH1D* cloned = (TH1D*) input_histogram->Clone(Form("%s_smeared", input_histogram->GetName()));
265 cloned->Reset();
266 Int_t N = smear_matrix->GetNbinsX() + 2;
267
268 TH2D* norm_matrix = smear_matrix; //(TH2D*)smear_matrix->Clone("temp");
269 for (int i = 0; i <= N - 1; i++) {
270 Double_t val = 0;
271 for (int j = 0; j < N; j++) {
272 val += smear_matrix->GetBinContent(i, j);
273 }
274 for (int j = 0; j < N; j++) {
275 if (val != 0) {
276 norm_matrix->SetBinContent(i, j, smear_matrix->GetBinContent(i, j) / val);
277 norm_matrix->SetBinError(i, j, smear_matrix->GetBinError(i, j) / val);
278 } else {
279 Double_t diag = 0.0;
280 if (j == i) diag = 1.0;
281 norm_matrix->SetBinContent(i, j, diag);
282 norm_matrix->SetBinError(i, j, 0);
283 }
284 }
285 }
286
287 TMatrixD AM(N, N);
288 TMatrixD AME(N, N);
289 TMatrixD XM(N, 1);
290 TMatrixD XME(N, 1);
291 TMatrixD YM(N, 1), YME(N, 1);
292
293 for (int i = 0; i < N; i++) {
294 XM[i][0] = input_histogram->GetBinContent(i);
295 XME[i][0] = input_histogram->GetBinError(i);
296 for (int j = 0; j < N; j++) {
297 AM[i][j] = norm_matrix->GetBinContent(j, i);
298 AME[i][j] = norm_matrix->GetBinError(j, i);
299 }
300 }
301
302 if (option.Contains("REV") || option.Contains("rev")) {
303 YME = AM * XME; // TODO approximation !
304 AM.Invert();
305 YM = AM * XM;
306 for (int i = 0; i < N; i++) {
307 cloned->SetBinContent(i, YM[i][0]);
308 cloned->SetBinError(i, YME[i][0]);
309 }
310 } else {
311 YM = AM * XM;
312 YME = AM * XME;
313 TMatrixD YME2 = AME * XM;
314 for (int i = 0; i < N; i++) {
315 cloned->SetBinContent(i, YM[i][0]);
316 cloned->SetBinError(i, YME[i][0]);
317 }
318 }
319 if (bad_map) delete smear_matrix;
320 // delete norm_matrix;
321 return cloned;
322 }
323
324 TH1* MakeHisto1D(TString name, TString title, TVector3 Xaxis, Char_t type) {
325 switch (type) {
326 case 'D': return new TH1D(name, title, (int) Xaxis.X(), Xaxis.Y(), Xaxis.Z()); break;
327 case 'F': return new TH1F(name, title, (int) Xaxis.X(), Xaxis.Y(), Xaxis.Z()); break;
328 case 'I': return new TH1I(name, title, (int) Xaxis.X(), Xaxis.Y(), Xaxis.Z()); break;
329 }
330 return nullptr;
331 }
332
333 TH2* MakeHisto2D(TString name, TString title, TVector3 Xaxis, TVector3 Yaxis, Char_t type) {
334 switch (type) {
335 case 'D':
336 return new TH2D(name, title, (int) Xaxis.X(), Xaxis.Y(), Xaxis.Z(), (int) Yaxis.X(), Yaxis.Y(), Yaxis.Z());
337 break;
338 case 'F':
339 return new TH2F(name, title, (int) Xaxis.X(), Xaxis.Y(), Xaxis.Z(), (int) Yaxis.X(), Yaxis.Y(), Yaxis.Z());
340 break;
341 case 'I':
342 return new TH2I(name, title, (int) Xaxis.X(), Xaxis.Y(), Xaxis.Z(), (int) Yaxis.X(), Yaxis.Y(), Yaxis.Z());
343 break;
344 }
345 return nullptr;
346 }
347
348 TH3* MakeHisto3D(TString name, TString title, TVector3 Xaxis, TVector3 Yaxis, TVector3 Zaxis, Char_t type) {
349 switch (type) {
350 case 'D':
351 return new TH3D(name,
352 title,
353 (int) Xaxis.X(),
354 Xaxis.Y(),
355 Xaxis.Z(),
356 (int) Yaxis.X(),
357 Yaxis.Y(),
358 Yaxis.Z(),
359 (int) Zaxis.X(),
360 Zaxis.Y(),
361 Zaxis.Z());
362 break;
363 case 'F':
364 return new TH3F(name,
365 title,
366 (int) Xaxis.X(),
367 Xaxis.Y(),
368 Xaxis.Z(),
369 (int) Yaxis.X(),
370 Yaxis.Y(),
371 Yaxis.Z(),
372 (int) Zaxis.X(),
373 Zaxis.Y(),
374 Zaxis.Z());
375 break;
376 case 'I':
377 return new TH3I(name,
378 title,
379 (int) Xaxis.X(),
380 Xaxis.Y(),
381 Xaxis.Z(),
382 (int) Yaxis.X(),
383 Yaxis.Y(),
384 Yaxis.Z(),
385 (int) Zaxis.X(),
386 Zaxis.Y(),
387 Zaxis.Z());
388 break;
389 }
390 return nullptr;
391 }
392
393 void HistogramEdges(TH1* h, TString option, Double_t value) {
394 Char_t axis = 'x';
395 enum eFill { kByVal, kByNeighbour };
396 enum eBin { kFirst, kLast };
397 if (option.Contains("y")) { axis = 'y'; }
398 if (option.Contains("z")) { axis = 'z'; }
399 eBin bin = kFirst;
400 if (option.Contains("ov")) { bin = kLast; }
401 eFill fill = kByNeighbour;
402 if (option.Contains("val")) {
403 fill = kByVal;
404 } else if (h->InheritsFrom("TH3")) {
405 if (axis == 'z') {
406 for (int i = 0; i <= h->GetNbinsX(); i++) {
407 for (int j = 0; j < h->GetNbinsY(); j++) {
408 switch (bin) {
409 case kFirst: {
410 switch (fill) {
411 case kByNeighbour: {
412 h->SetBinContent(i, j, 0, h->GetBinContent(i, j, 1));
413 } break;
414 case kByVal: {
415 h->SetBinContent(i, j, 0, value);
416 } break;
417 }
418 } break;
419 case kLast: {
420 switch (fill) {
421 case kByNeighbour: {
422 h->SetBinContent(i, j, h->GetNbinsX() + 1, h->GetBinContent(i, j, 1));
423 } break;
424 case kByVal: {
425 h->SetBinContent(i, j, h->GetNbinsX() + 1, value);
426 } break;
427 }
428 } break;
429 }
430 }
431 }
432 } else if (axis == 'y') {
433 for (int i = 0; i <= h->GetNbinsX(); i++) {
434 for (int j = 0; j < h->GetNbinsZ(); j++) {
435 switch (bin) {
436 case kFirst: {
437 switch (fill) {
438 case kByNeighbour: {
439 h->SetBinContent(i, 0, j, h->GetBinContent(i, 1, j));
440 } break;
441 case kByVal: {
442 h->SetBinContent(i, 0, j, value);
443 } break;
444 }
445 } break;
446 case kLast: {
447 switch (fill) {
448 case kByNeighbour: {
449 h->SetBinContent(i, h->GetNbinsX() + 1, j, h->GetBinContent(i, 1, j));
450 } break;
451 case kByVal: {
452 h->SetBinContent(i, h->GetNbinsX() + 1, j, value);
453 } break;
454 }
455 } break;
456 }
457 }
458 }
459 } else { // x
460 for (int i = 0; i <= h->GetNbinsY(); i++) {
461 for (int j = 0; j < h->GetNbinsZ(); j++) {
462 switch (bin) {
463 case kFirst: {
464 switch (fill) {
465 case kByNeighbour: {
466 h->SetBinContent(0, i, j, h->GetBinContent(1, i, j));
467 } break;
468 case kByVal: {
469 h->SetBinContent(0, i, j, value);
470 } break;
471 }
472 } break;
473 case kLast: {
474 switch (fill) {
475 case kByNeighbour: {
476 h->SetBinContent(h->GetNbinsX() + 1, i, j, h->GetBinContent(1, i, j));
477 } break;
478 case kByVal: {
479 h->SetBinContent(h->GetNbinsX() + 1, i, j, value);
480 } break;
481 }
482 } break;
483 }
484 }
485 }
486 }
487 } else if (h->InheritsFrom("TH2")) {
488 if (axis == 'y') {
489 for (int i = 0; i <= h->GetNbinsX() + 1; i++) {
490 switch (bin) {
491 case kFirst: {
492 switch (fill) {
493 case kByNeighbour: {
494 h->SetBinContent(i, 0, h->GetBinContent(i, 1));
495 } break;
496 case kByVal: {
497 h->SetBinContent(i, 0, value);
498 } break;
499 }
500 } break;
501 case kLast: {
502 switch (fill) {
503 case kByNeighbour: {
504 h->SetBinContent(i, h->GetNbinsY() + 1, h->GetBinContent(i, h->GetNbinsY()));
505 } break;
506 case kByVal: {
507 h->SetBinContent(i, h->GetNbinsY() + 1, value);
508 } break;
509 }
510 } break;
511 }
512 }
513 } else { // axis x
514 for (int i = 0; i <= h->GetNbinsY() + 1; i++) {
515 switch (bin) {
516 case kFirst: {
517 switch (fill) {
518 case kByNeighbour: {
519 h->SetBinContent(0, i, h->GetBinContent(1, i));
520 } break;
521 case kByVal: {
522 h->SetBinContent(0, i, value);
523 } break;
524 }
525 } break;
526 case kLast: {
527 switch (fill) {
528 case kByNeighbour: {
529 h->SetBinContent(h->GetNbinsX() + 1, i, h->GetBinContent(h->GetNbinsX(), i));
530 } break;
531 case kByVal: {
532 h->SetBinContent(h->GetNbinsX() + 1, i, value);
533 } break;
534 }
535 } break;
536 }
537 }
538 }
539
540 } else { // TH1
541 if (bin == kFirst) {
542 if (fill == kByVal) { h->SetBinContent(0, value); }
543 if (fill == kByNeighbour) { h->SetBinContent(0, h->GetBinContent(1)); }
544 }
545 if (bin == kLast) {
546 if (fill == kByVal) { h->SetBinContent(h->GetNbinsX() + 1, value); }
547 if (fill == kByNeighbour) { h->SetBinContent(h->GetNbinsX() + 1, h->GetBinContent(h->GetNbinsX())); }
548 }
549 }
550 }
551
552 void HistogramExtend(TH1* h, Char_t axis, Double_t factor) {
553 TAxis* x = NULL;
554 switch (axis) {
555 case 'y': {
556 x = h->GetYaxis();
557 } break;
558 case 'z': {
559 x = h->GetZaxis();
560 } break;
561 default: {
562 x = h->GetXaxis();
563 } break;
564 }
565 if (x->IsVariableBinSize()) {
566 const TArrayD* ar = x->GetXbins();
567 Int_t size = ar->GetSize();
568 Double_t* edges = new Double_t[size];
569 for (int i = 0; i < ar->GetSize(); i++) {
570 edges[i] = ar->At(i) * factor;
571 }
572 x->Set(size - 1, edges);
573 delete[] edges;
574 } else {
575 x->Set(x->GetNbins(), x->GetXmin() * factor, x->GetXmax() * factor);
576 }
577 h->ResetStats();
578 }
579
580 void CopyAxisProp(const TAxis* from, TAxis* to, TString option) {
581 TAxis default_axis = TAxis();
582 if (default_axis.GetCenterTitle() != from->GetCenterTitle()) to->CenterTitle(from->GetCenterTitle());
583 if (default_axis.GetNdivisions() != from->GetNdivisions()) to->SetNdivisions(from->GetNdivisions());
584 if (default_axis.GetTitleColor() != from->GetTitleColor()) to->SetTitleColor(from->GetTitleColor());
585 if (default_axis.GetTitleFont() != from->GetTitleFont()) to->SetTitleFont(from->GetTitleFont());
586 if (default_axis.GetTitleOffset() != from->GetTitleOffset()) to->SetTitleOffset(from->GetTitleOffset());
587 if (0.035 != from->GetTitleSize()) to->SetTitleSize(from->GetTitleSize()); // why 0.35 is not default? don't know
588 if (default_axis.GetLabelColor() != from->GetLabelColor()) to->SetLabelColor(from->GetLabelColor());
589 if (default_axis.GetLabelFont() != from->GetLabelFont()) to->SetLabelFont(from->GetLabelFont());
590 if (default_axis.GetLabelOffset() != from->GetLabelOffset()) to->SetLabelOffset(from->GetLabelOffset());
591 if (default_axis.GetLabelSize() != from->GetLabelSize()) to->SetLabelSize(from->GetLabelSize());
592 if (!Hal::Std::FindParam(option, "!tit", kTRUE)) { to->SetTitle(from->GetTitle()); }
593 }
594
595 TH1D* GetDiagonalProjection1D(TH3* h, TString dir, Double_t start, Double_t start2) {
596 if (h->GetNbinsX() != h->GetNbinsY() || h->GetNbinsX() != h->GetNbinsZ()) {
597 Hal::Cout::PrintInfo("Cannot make projection in nonsymetric histo", Hal::EInfo::kLowWarning);
598 return nullptr;
599 }
600 TString temp = dir;
601 Double_t xmin = h->GetXaxis()->GetBinLowEdge(1);
602 Double_t xmax = h->GetXaxis()->GetBinUpEdge(h->GetNbinsX());
603 Int_t binsTot = h->GetXaxis()->GetNbins();
604 TString titleX = h->GetXaxis()->GetTitle();
605 TString titleY = h->GetYaxis()->GetTitle();
606 TString titleZ = h->GetZaxis()->GetTitle();
607 TString randname = Form("%i", (int) gRandom->Uniform(1E+9));
608 TH1D* res = new TH1D(randname, "res", binsTot, xmin, xmax);
609 res->SetMarkerStyle(h->GetMarkerStyle());
610 res->SetMarkerColor(h->GetMarkerColor());
611 res->SetMarkerSize(h->GetMarkerSize());
612 res->SetLineColor(h->GetLineColor());
613 res->SetLineStyle(h->GetLineStyle());
614 res->SetLineWidth(h->GetLineWidth());
615 TString newTitle = "";
616 std::vector<int> binsVecX(binsTot), binsVecY(binsTot), binsVecZ(binsTot);
617 if (temp.Length() == 6 && temp.BeginsWith("xyz")) {
618 for (int i = 1; i <= binsTot; i++) {
619 binsVecX[i - 1] = i;
620 binsVecY[i - 1] = i;
621 binsVecZ[i - 1] = i;
622 }
623 if (temp[3] == '-') std::sort(binsVecX.begin(), binsVecX.end(), std::greater<int>());
624 if (temp[4] == '-') std::sort(binsVecY.begin(), binsVecY.end(), std::greater<int>());
625 if (temp[5] == '-') std::sort(binsVecZ.begin(), binsVecZ.end(), std::greater<int>());
626 newTitle = titleX + temp[3] + titleY + temp[4] + titleZ + temp[5];
627 } else if (temp.Length() == 4) {
628 if (temp.BeginsWith("yz")) {
629 Int_t binMid = h->GetXaxis()->FindBin(start);
630 for (int i = 1; i <= binsTot; i++) {
631 binsVecX[i - 1] = binMid;
632 binsVecY[i - 1] = i;
633 binsVecZ[i - 1] = i;
634 }
635 if (temp[2] == '-') std::sort(binsVecY.begin(), binsVecY.end(), std::greater<int>());
636 if (temp[3] == '-') std::sort(binsVecZ.begin(), binsVecZ.end(), std::greater<int>());
637 newTitle = titleY + temp[2] + titleZ + temp[3];
638 } else if (temp.BeginsWith("xz")) {
639 Int_t binMid = h->GetYaxis()->FindBin(start);
640 for (int i = 1; i <= binsTot; i++) {
641 binsVecX[i - 1] = i;
642 binsVecY[i - 1] = binMid;
643 binsVecZ[i - 1] = i;
644 }
645 if (temp[2] == '-') std::sort(binsVecX.begin(), binsVecX.end(), std::greater<int>());
646 if (temp[3] == '-') std::sort(binsVecZ.begin(), binsVecZ.end(), std::greater<int>());
647 newTitle = titleX + temp[2] + titleZ + temp[3];
648 } else { // xy
649 Int_t binMid = h->GetZaxis()->FindBin(start);
650 for (int i = 1; i <= binsTot; i++) {
651 binsVecX[i - 1] = i;
652 binsVecY[i - 1] = i;
653 binsVecZ[i - 1] = binMid;
654 }
655 if (temp[2] == '-') std::sort(binsVecX.begin(), binsVecX.end(), std::greater<int>());
656 if (temp[3] == '-') std::sort(binsVecY.begin(), binsVecY.end(), std::greater<int>());
657 newTitle = titleX + temp[2] + titleY + temp[3];
658 }
659 } else if (temp.Length() == 2 || temp.Length() == 1) {
660 TString signFlag = "";
661 if (temp.Length() == 2) {
662 if (temp[1] == '-') signFlag = "-";
663 if (temp[1] == '+') signFlag = "+";
664 }
665 Bool_t rev = kFALSE;
666 if (signFlag.EqualTo("-")) rev = kTRUE;
667 if (temp[0] == 'y') {
668 Int_t bin1 = h->GetXaxis()->FindBin(start);
669 Int_t bin2 = h->GetZaxis()->FindBin(start2);
670 for (int i = 1; i <= binsTot; i++) {
671 binsVecX[i - 1] = bin1;
672 binsVecY[i - 1] = i;
673 binsVecZ[i - 1] = bin2;
674 }
675 if (rev) std::sort(binsVecY.begin(), binsVecY.end(), std::greater<int>());
676 newTitle = titleY + signFlag;
677 } else if (temp[0] == 'z') {
678 Int_t bin1 = h->GetXaxis()->FindBin(start);
679 Int_t bin2 = h->GetYaxis()->FindBin(start2);
680 for (int i = 1; i <= binsTot; i++) {
681 binsVecX[i - 1] = bin1;
682 binsVecY[i - 1] = bin2;
683 binsVecZ[i - 1] = i;
684 }
685 if (rev) std::sort(binsVecZ.begin(), binsVecZ.end(), std::greater<int>());
686 newTitle = titleZ + signFlag;
687 } else {
688 Int_t bin1 = h->GetYaxis()->FindBin(start);
689 Int_t bin2 = h->GetZaxis()->FindBin(start2);
690 for (int i = 1; i <= binsTot; i++) {
691 binsVecX[i - 1] = i;
692 binsVecY[i - 1] = bin1;
693 binsVecZ[i - 1] = bin2;
694 }
695 if (rev) std::sort(binsVecX.begin(), binsVecX.end(), std::greater<int>());
696 newTitle = titleX + signFlag;
697 }
698 } else {
699 Hal::Cout::PrintInfo("Std::GetDiagonalProjection1D unknown option", Hal::EInfo::kLowWarning);
700 delete res;
701 return nullptr;
702 }
703
704 for (int i = 1; i <= binsTot; i++) {
705 res->SetBinContent(i, h->GetBinContent(binsVecX[i - 1], binsVecY[i - 1], binsVecZ[i - 1]));
706 res->SetBinError(i, h->GetBinError(binsVecX[i - 1], binsVecY[i - 1], binsVecZ[i - 1]));
707 }
708 res->GetXaxis()->SetTitle(newTitle);
709 return res;
710 } // namespace Hal::Std
711
712 Bool_t AreSimilar(const TH1* x, const TH1* y, Bool_t classes) {
713 if (classes) {
714 TString classname1 = x->ClassName();
715 TString classname2 = y->ClassName();
716 if (!classname1.EqualTo(classname2)) return kFALSE;
717 }
718 if (x->InheritsFrom("TH3") && y->InheritsFrom("TH3")) {
719 const TH3* h1 = static_cast<const TH3*>(x);
720 const TH3* h2 = static_cast<const TH3*>(y);
721 if (h1->GetNbinsX() != h2->GetNbinsX()) return kFALSE;
722 if (h1->GetNbinsY() != h2->GetNbinsY()) return kFALSE;
723 if (h1->GetNbinsZ() != h2->GetNbinsZ()) return kFALSE;
724 if (h1->GetXaxis()->GetXmin() != h2->GetXaxis()->GetXmin()) return kFALSE;
725 if (h1->GetXaxis()->GetXmax() != h2->GetXaxis()->GetXmax()) return kFALSE;
726 if (h1->GetYaxis()->GetXmin() != h2->GetYaxis()->GetXmin()) return kFALSE;
727 if (h1->GetYaxis()->GetXmax() != h2->GetYaxis()->GetXmax()) return kFALSE;
728 if (h1->GetZaxis()->GetXmin() != h2->GetZaxis()->GetXmin()) return kFALSE;
729 if (h1->GetZaxis()->GetXmax() != h2->GetZaxis()->GetXmax()) return kFALSE;
730 } else if (x->InheritsFrom("TH2") && y->InheritsFrom("TH2")) {
731 const TH2* h1 = static_cast<const TH2*>(x);
732 const TH2* h2 = static_cast<const TH2*>(y);
733 if (h1->GetNbinsX() != h2->GetNbinsX()) return kFALSE;
734 if (h1->GetNbinsY() != h2->GetNbinsY()) return kFALSE;
735 if (h1->GetXaxis()->GetXmin() != h2->GetXaxis()->GetXmin()) return kFALSE;
736 if (h1->GetXaxis()->GetXmax() != h2->GetXaxis()->GetXmax()) return kFALSE;
737 if (h1->GetYaxis()->GetXmin() != h2->GetYaxis()->GetXmin()) return kFALSE;
738 if (h1->GetYaxis()->GetXmax() != h2->GetYaxis()->GetXmax()) return kFALSE;
739 } else {
740 const TH1* h1 = static_cast<const TH1*>(x);
741 const TH1* h2 = static_cast<const TH1*>(y);
742 if (h1->GetNbinsX() != h2->GetNbinsX()) return kFALSE;
743 if (h1->GetXaxis()->GetXmin() != h2->GetXaxis()->GetXmin()) return kFALSE;
744 if (h1->GetXaxis()->GetXmax() != h2->GetXaxis()->GetXmax()) return kFALSE;
745 }
746 return kTRUE;
747 }
748
749 void GetAxisPar(const TH1& obj, Int_t& nbins, Double_t& min, Double_t& max, Option_t* opt) {
750 TString option = opt;
751 if (option == "x" || option == "X") {
752 const TAxis* x = obj.GetXaxis();
753 nbins = x->GetNbins();
754 min = x->GetBinLowEdge(1);
755 max = x->GetBinUpEdge(nbins);
756 } else if (option == "y" || option == "Y") {
757 const TAxis* x = obj.GetYaxis();
758 nbins = x->GetNbins();
759 min = x->GetBinLowEdge(1);
760 max = x->GetBinUpEdge(nbins);
761 } else if (option == "z" || option == "Z") {
762 const TAxis* x = obj.GetZaxis();
763 nbins = x->GetNbins();
764 min = x->GetBinLowEdge(1);
765 max = x->GetBinUpEdge(nbins);
766 }
767 }
768
769 std::vector<double> GetAxisCenters(const TH1& obj, Option_t* opt) {
770 TString option = opt;
771 auto getAxis = [](TString& opti, const TH1& hst) {
772 if (Hal::Std::FindParam(opti, "x", kTRUE)) { return hst.GetXaxis(); }
773 if (Hal::Std::FindParam(opti, "y", kTRUE)) { return hst.GetYaxis(); }
774 if (Hal::Std::FindParam(opti, "z", kTRUE)) { return hst.GetZaxis(); }
775 return (const TAxis*) nullptr;
776 };
777 auto axis = (const TAxis*) getAxis(option, obj);
778 if (axis == nullptr) return std::vector<double>();
779
780 int bins = axis->GetNbins();
781 std::vector<double> vec;
782 int startbin = 1;
783 if (Hal::Std::FindParam(option, "o", kTRUE)) {
784 startbin = 0;
785 bins++;
786 }
787 for (int i = startbin; i <= bins; i++) {
788 vec.push_back(axis->GetBinCenter(i));
789 }
790 return vec;
791 }
792
793 namespace { // anonymous namespace this is for computing axes of histogram
794
795 std::vector<int> FoldAxis(const TAxis& x, Double_t val) {
796 const Int_t nBins = x.GetNbins();
797 Int_t binFold = x.FindBin(val);
798 Double_t center = x.GetBinCenter(binFold);
799 // Double_t low = x.GetBinLowEdge(binFold);
800 std::vector<int> foldBinMap(nBins + 1);
801 Int_t shift = 0;
802 if (val != center) shift = 1;
803 for (int i = 1; i <= nBins; i++) {
804 Int_t delta = binFold - i + shift;
805 Int_t bin = binFold + delta;
806 if (bin < 1) bin = bin + nBins;
807 if (bin > nBins) bin = bin - nBins;
808 foldBinMap[i] = bin;
809 }
810 return foldBinMap;
811 }
812
813 void CropAxis(const TAxis* x,
814 Int_t& nbins,
815 std::pair<int, int>& binId,
816 std::pair<double, double>& val,
817 Double_t lmin,
818 Double_t lmax,
819 TString opt) {
820 if (opt == "vals") {
821 binId.first = x->FindBin(lmin);
822 binId.second = x->FindBin(lmax);
823 } else {
824 binId.first = lmin;
825 binId.second = lmax;
826 }
827 binId.first = TMath::Max(1, binId.first);
828 binId.second = TMath::Min(binId.second, x->GetNbins());
829 val.first = x->GetBinLowEdge(binId.first);
830 val.second = x->GetBinUpEdge(binId.second);
831 nbins = binId.second - binId.first + 1;
832 };
833
834 } // namespace
835
836 void Fold1D(Double_t val, TH1& h) {
837 TH1* tempCopy = (TH1*) h.Clone();
838 const Int_t nBins = h.GetXaxis()->GetNbins();
839 std::vector<int> foldBinMap = FoldAxis(*h.GetXaxis(), val);
840 h.Reset();
841 for (int i = 1; i <= nBins; i++) {
842 Double_t old = h.GetBinContent(i);
843 h.SetBinContent(i, old + tempCopy->GetBinContent(foldBinMap[i]));
844 Double_t sumv = tempCopy->GetSumw2()->GetAt(foldBinMap[i]);
845 h.GetSumw2()->SetAt(h.GetSumw2()->At(i) + sumv, i);
846 }
847 delete tempCopy;
848 }
849
850 void Fold2D(Double_t val, TH2& h, TString opt) {
851 TH2* tempCopy = (TH2*) h.Clone();
852 const Int_t nBinsX = h.GetXaxis()->GetNbins();
853 const Int_t nBinsY = h.GetYaxis()->GetNbins();
854 if (Hal::Std::FindParam(opt, "x", kTRUE)) {
855 std::vector<int> foldBinMap = FoldAxis(*h.GetXaxis(), val);
856 for (int i = 1; i <= nBinsX; i++) {
857 for (int j = 1; j <= nBinsY; j++) {
858 Double_t old = h.GetBinContent(i, j);
859 Double_t add = tempCopy->GetBinContent(foldBinMap[i], j);
860 h.SetBinContent(i, j, old + add);
861 }
862 }
863 } else if (Hal::Std::FindParam(opt, "y", kTRUE)) {
864 std::vector<int> foldBinMap = FoldAxis(*h.GetYaxis(), val);
865 for (int i = 1; i <= nBinsX; i++) {
866 for (int j = 1; j <= nBinsY; j++) {
867 Double_t old = h.GetBinContent(i, j);
868 Double_t add = tempCopy->GetBinContent(i, foldBinMap[j]);
869 h.SetBinContent(i, j, old + add);
870 }
871 }
872 }
873 h.ResetStats();
874 delete tempCopy;
875 }
876
877 void Fold3D(Double_t val, TH2& h, TString opt) {
878 TH3* tempCopy = (TH3*) h.Clone();
879 const Int_t nBinsX = h.GetXaxis()->GetNbins();
880 const Int_t nBinsY = h.GetYaxis()->GetNbins();
881 const Int_t nBinsZ = h.GetZaxis()->GetNbins();
882 if (Hal::Std::FindParam(opt, "x", kTRUE)) {
883 std::vector<int> foldBinMap = FoldAxis(*h.GetXaxis(), val);
884 for (int i = 1; i <= nBinsX; i++) {
885 for (int j = 1; j <= nBinsY; j++) {
886 for (int k = 1; k <= nBinsZ; k++) {
887 Double_t old = h.GetBinContent(i, j, k);
888 Double_t add = tempCopy->GetBinContent(foldBinMap[i], j, k);
889 h.SetBinContent(i, j, k, old + add);
890 }
891 }
892 }
893 } else if (Hal::Std::FindParam(opt, "y", kTRUE)) {
894 std::vector<int> foldBinMap = FoldAxis(*h.GetYaxis(), val);
895 for (int i = 1; i <= nBinsX; i++) {
896 for (int j = 1; j <= nBinsY; j++) {
897 for (int k = 1; k <= nBinsZ; k++) {
898 Double_t old = h.GetBinContent(i, j, k);
899 Double_t add = tempCopy->GetBinContent(i, foldBinMap[j], k);
900 h.SetBinContent(i, j, k, old + add);
901 }
902 }
903 }
904 } else if (Hal::Std::FindParam(opt, "z", kTRUE)) {
905 std::vector<int> foldBinMap = FoldAxis(*h.GetZaxis(), val);
906 for (int i = 1; i <= nBinsX; i++) {
907 for (int j = 1; j <= nBinsY; j++) {
908 for (int k = 1; k <= nBinsZ; k++) {
909 Double_t old = h.GetBinContent(i, j, k);
910 Double_t add = tempCopy->GetBinContent(i, j, foldBinMap[k]);
911 h.SetBinContent(i, j, k, old + add);
912 }
913 }
914 }
915 }
916 h.ResetStats();
917 delete tempCopy;
918 }
919
920 void SetColor(TH1& h, Color_t color) {
921 h.SetLineColor(color);
922 h.SetMarkerColor(color);
923 }
924
925 void SetColorAndMarker(TH1& h, Color_t color, Marker_t m, Size_t s) {
926 h.SetLineColor(color);
927 h.SetMarkerColor(color);
928 h.SetMarkerStyle(m);
929 h.SetMarkerSize(s);
930 }
931
932 void MakeBeautiful() {
933 gStyle->SetPalette(kRainBow);
934 TGaxis::SetMaxDigits(3);
935 gStyle->SetCanvasPreferGL(kTRUE);
936 }
937
938 void SetRainbow(TH2& h, Double_t x1, Double_t y1, Double_t x2, Double_t y2) {
939 gPad->Update();
940 TPaletteAxis* palette = (TPaletteAxis*) h.GetListOfFunctions()->FindObject("palette");
941 palette->SetX1NDC(x1);
942 palette->SetX2NDC(x2);
943 palette->SetY1NDC(y1);
944 palette->SetY2NDC(y2);
945 gPad->Modified();
946 gPad->Update();
947 }
948
949 Int_t GetAntiColor(Int_t col) {
950 if (col < 0) return -1;
951
952 // Get list of all defined colors
953 TObjArray* colors = (TObjArray*) ::ROOT::GetROOT()->GetListOfColors();
954 Int_t ncolors = colors->GetSize();
955 // Get existing color at index n
956 TColor* color = nullptr;
957 if (col < ncolors) color = (TColor*) colors->At(col);
958 if (!color) return -1;
959
960 // Get the rgb of the new bright color corresponding to color n
961 Float_t r, g, b;
962
963 color->GetRGB(r, g, b);
964 r = 1.0 - r;
965 b = 1.0 - b;
966 g = 1.0 - g;
967
968 // Build the bright color (unless the slot nb is already used)
969 Int_t nb = col + 150;
970 TColor* colorb = nullptr;
971 if (nb < ncolors) colorb = (TColor*) colors->At(nb);
972 if (colorb) return nb;
973 colorb = new TColor(nb, r, g, b);
974 colorb->SetName(Form("%s_bright", color->GetName()));
975 colors->AddAtAndExpand(colorb, nb);
976 return nb;
977 }
978
979 Int_t GetListOfSubPads(TVirtualPad* pad) {
980 TList* l = pad->GetListOfPrimitives();
981 Int_t subpads = 0;
982 for (int i = 0; i < l->GetEntries(); i++) {
983 if (dynamic_cast<TVirtualPad*>(l->At(i))) subpads++;
984 }
985 return subpads;
986 }
987
988 TH1D* Crop1D(const TH1& h, Double_t min, Double_t max, TString option) {
989 std::pair<int, int> binsX;
990 std::pair<double, double> valsX;
991 Int_t nBins = 0;
992 CropAxis(h.GetXaxis(), nBins, binsX, valsX, min, max, option);
993
994 TH1D* nh = new TH1D(h.GetName(), h.GetTitle(), nBins, valsX.first, valsX.second);
995 int count = 0;
996 for (int i = binsX.first; i <= binsX.second; i++) {
997 nh->SetBinContent(++count, h.GetBinContent(i));
998 nh->SetBinError(count, h.GetBinError(i));
999 }
1000 return nh;
1001 }
1002
1003 TH2D* Crop2D(const TH2& h, Double_t minX, Double_t maxX, Double_t minY, Double_t maxY, TString option) {
1004 Int_t nBinsX, nBinsY;
1005 std::pair<int, int> binsX, binsY;
1006 std::pair<double, double> valsX, valsY;
1007 CropAxis(h.GetXaxis(), nBinsX, binsX, valsX, minX, maxX, option);
1008 CropAxis(h.GetYaxis(), nBinsY, binsY, valsY, minY, maxY, option);
1009 TH2D* nh = new TH2D(h.GetName(), h.GetTitle(), nBinsX, valsX.first, valsX.second, nBinsY, valsY.first, valsY.second);
1010 int countX = 0;
1011 int countY = 0;
1012 for (int i = binsX.first; i <= binsX.second; i++) {
1013 ++countX;
1014 for (int j = binsY.first; j <= binsY.second; j++) {
1015 ++countY;
1016 nh->SetBinContent(countX, countY, h.GetBinContent(i, j));
1017 nh->SetBinError(countX, countY, h.GetBinError(i, j));
1018 }
1019 }
1020 return nh;
1021 }
1022
1023 TH3D* Crop3D(const TH3& h,
1024 Double_t minX,
1025 Double_t maxX,
1026 Double_t minY,
1027 Double_t maxY,
1028 Double_t minZ,
1029 Double_t maxZ,
1030 TString option) {
1031 std::pair<int, int> X, Y, Z;
1032 std::pair<double, double> valsX, valsY, valsZ;
1033
1034 Int_t nBinsX, nBinsY, nBinsZ;
1035 CropAxis(h.GetXaxis(), nBinsX, X, valsX, minX, maxX, option);
1036 CropAxis(h.GetYaxis(), nBinsY, Y, valsY, minY, maxY, option);
1037 CropAxis(h.GetZaxis(), nBinsZ, Z, valsZ, minZ, maxZ, option);
1038 TH3D* nh = new TH3D(h.GetName(),
1039 h.GetTitle(),
1040 nBinsX,
1041 valsX.first,
1042 valsX.second,
1043 nBinsY,
1044 valsY.first,
1045 valsY.second,
1046 nBinsZ,
1047 valsZ.first,
1048 valsZ.second);
1049 int countX = 0;
1050 int countY = 0;
1051 int countZ = 0;
1052 for (int i = X.first; i <= X.second; i++) {
1053 ++countX;
1054 countY = 0;
1055 for (int j = Y.first; j <= Y.second; j++) {
1056 ++countY;
1057 countZ = 0;
1058 for (int k = Z.first; k <= Z.second; k++) {
1059 ++countZ;
1060 nh->SetBinContent(countX, countY, countZ, h.GetBinContent(i, j, k));
1061 nh->SetBinError(countX, countY, countZ, h.GetBinError(i, j, k));
1062 }
1063 }
1064 }
1065 return nh;
1066 }
1067 std::vector<TObject*> GetPadChildren(TString objName, TString className, TVirtualPad* pad) {
1068 if (!pad) pad = gPad;
1069 std::vector<TObject*> list;
1070 if (!pad) return list;
1071 TList* l = gPad->GetListOfPrimitives();
1072 for (int i = 0; i < l->GetEntries(); i++) {
1073 auto obj = l->At(i);
1074 TString classNameObj = obj->ClassName();
1075 TString objectName = obj->GetName();
1076 Bool_t objNameSame = ((objName == objectName) || objName.Length() == 0);
1077 Bool_t classNameSame = (className == classNameObj || className.Length() == 0);
1078 if (objNameSame && classNameSame) { list.push_back(obj); }
1079 }
1080 return list;
1081 }
1082 Bool_t CheckHistogramData(const TH1& h1, const TH1& h2, Double_t thres, Option_t* opt) {
1083 TString option = opt;
1084 Bool_t debug = Hal::Std::FindParam(option, "print", kTRUE);
1085 Bool_t skiperr = Hal::Std::FindParam(option, "skiperr", kTRUE);
1086 Bool_t relative = Hal::Std::FindParam(option, "rel", kTRUE);
1087 if (!AreSimilar(&h1, &h2, kFALSE)) {
1088 if (debug) std::cout << __FILE__ << __LINE__ << " histograms are not similiar" << std::endl;
1089 return kFALSE;
1090 }
1091 auto compare = [&](double x1, double e1, double x2, double e2) {
1092 if (relative) {
1093 double relval = TMath::Abs((x1 - x2) / x1);
1094 if (relval > thres) return kFALSE;
1095 } else {
1096 if (TMath::Abs(x1 - x2) > thres) return kFALSE;
1097 }
1098
1099 if (skiperr) return kTRUE;
1100 if (relative) {
1101 double relval = TMath::Abs((e1 - e2) / e1);
1102 if (relval > thres) return kFALSE;
1103 } else {
1104 if (TMath::Abs(e1 - e2) > thres) return kFALSE;
1105 }
1106 return kTRUE;
1107 };
1108 if (dynamic_cast<const TH3*>(&h1)) {
1109 auto H1 = static_cast<const TH3*>(&h1);
1110 auto H2 = static_cast<const TH3*>(&h2);
1111 for (int i = 1; i <= H1->GetNbinsX(); i++) {
1112 for (int j = 1; j <= H1->GetNbinsY(); j++) {
1113 for (int k = 1; k <= H1->GetNbinsZ(); k++) {
1114 bool compared = compare(
1115 H1->GetBinContent(i, j, k), H1->GetBinError(i, j, k), H2->GetBinContent(i, j, k), H2->GetBinError(i, j, k));
1116 if (!compared) {
1117 if (debug) { std::cout << "Detected difference at bin " << i << " " << j << " " << k << std::endl; }
1118 return kFALSE;
1119 }
1120 }
1121 }
1122 }
1123 } else if (dynamic_cast<const TH2*>(&h1)) {
1124 auto H1 = static_cast<const TH2*>(&h1);
1125 auto H2 = static_cast<const TH2*>(&h2);
1126 for (int i = 1; i <= H1->GetNbinsX(); i++) {
1127 for (int j = 1; j <= H1->GetNbinsY(); j++) {
1128 bool compared =
1129 compare(H1->GetBinContent(i, j), H1->GetBinError(i, j), H2->GetBinContent(i, j), H2->GetBinError(i, j));
1130 if (!compared) {
1131 if (debug) { std::cout << "Detected difference at bin " << i << " " << j << " " << std::endl; }
1132 return kFALSE;
1133 }
1134 }
1135 }
1136 } else {
1137 auto H1 = static_cast<const TH1*>(&h1);
1138 auto H2 = static_cast<const TH1*>(&h2);
1139 for (int i = 1; i <= H1->GetNbinsX(); i++) {
1140 bool compared = compare(H1->GetBinContent(i), H1->GetBinError(i), H2->GetBinContent(i), H2->GetBinError(i));
1141 if (!compared) {
1142 if (debug) { std::cout << "Detected difference at bin " << i << " " << std::endl; }
1143 return kFALSE;
1144 }
1145 }
1146 }
1147 return kTRUE;
1148 }
1149 void CopyHistProp(const TH1& from, TH1& to, TString opt) {
1150 auto d3 = static_cast<const TH3*>(&from);
1151
1152 CopyAxisProp(from.GetXaxis(), to.GetXaxis(), opt);
1153 CopyAxisProp(from.GetYaxis(), to.GetYaxis(), opt);
1154 if (d3) CopyAxisProp(from.GetZaxis(), to.GetZaxis(), opt);
1155 to.SetMarkerSize(from.GetMarkerSize());
1156 to.SetMarkerStyle(from.GetMarkerStyle());
1157 to.SetMarkerColor(from.GetMarkerColor());
1158 to.SetLineStyle(from.GetLineStyle());
1159 to.SetLineWidth(from.GetLineWidth());
1160 to.SetLineColor(from.GetLineColor());
1161 to.SetFillColor(from.GetFillColor());
1162 to.SetFillStyle(from.GetFillStyle());
1163 }
1164 } // namespace Std
1165} // namespace Hal
Hal::Cout::PrintInfo
static void PrintInfo(TString text, Hal::EInfo status)
Definition Cout.cxx:370
Hal::Spline2D::Eval
Double_t Eval(Double_t x, Double_t y) const
Definition Splines.cxx:327
Hal::Spline2D::GetError
Double_t GetError(Double_t x, Double_t y) const
Definition Splines.cxx:481