Heavy ion Analysis Libriares
Loading...
Searching...
No Matches
FemtoSHCF.cxx
1/*
2 * FemtoSHCF.cxx
3 *
4 * Created on: 28 kwi 2016
5 * Author: Daniel Wielanek
6 * E-mail: daniel.wielanek@gmail.com
7 * Warsaw University of Technology, Faculty of Physics
8 */
9
10#include "FemtoSHCF.h"
11
12#include <TAttLine.h>
13#include <TAxis.h>
14#include <TBrowser.h>
15#include <TCanvas.h>
16#include <TCollection.h>
17#include <TH1.h>
18#include <TH3.h>
19#include <TList.h>
20#include <TMathBase.h>
21#include <TNamed.h>
22#include <TObjString.h>
23#include <TROOT.h>
24#include <TRandom.h>
25#include <TSystem.h>
26#ifndef DISABLE_GSL
27#include <gsl/gsl_blas.h>
28#include <gsl/gsl_cblas.h>
29#include <gsl/gsl_matrix_double.h>
30#include <gsl/gsl_vector_double.h>
31#endif
32#include <iostream>
33#ifdef __CIA__
34#include "CorrFitSHCF.h"
35#endif
36#include "Cout.h"
37#include "FemtoPair.h"
38#include "FemtoSHCF.h"
39#include "FemtoSHCFPainter.h"
40#include "FemtoSHSlice.h"
41#include "FemtoSerializationInterfaceSH.h"
42#include "FemtoYlmSolver.h"
43#include "HistoStyle.h"
44#include "HtmlCore.h"
45#include "HtmlFile.h"
46#include "HtmlObject.h"
47#include "HtmlTable.h"
48#include "LineStyle.h"
49#include "MarkerStyle.h"
50#include "PadStyle.h"
51#include "Std.h"
52#include "StdHist.h"
53#include "StdString.h"
54#include "Style.h"
55
56//#define FULL_CALC
57//#define _FINISH_DEBUG_
58namespace Hal {
61 fMaxJM(0),
62 fFrame(Femto::EKinematics::kLCMS),
63 fNormPurity(0),
64 fNormRadius(0),
65 fNormBohr(0),
66 fLmVals(FemtoYlmIndexes(1)) {
67 gSystem->Load("libgsl.so");
68 gSystem->Load("libgslcblas.so");
69 }
70
71 FemtoSHCF::FemtoSHCF(TString name, Int_t maxL, Int_t bins, Double_t min, Double_t max, Femto::EKinematics kinematics) :
72 DividedHisto1D(name, bins, min, max),
73 fMaxJM((TMath::Min(maxL, 5) + 1) * (TMath::Min(maxL, 5) + 1)),
74 fFrame(kinematics),
75 fNormPurity(0),
76 fNormRadius(0),
77 fNormBohr(0),
78 fLmVals(FemtoYlmIndexes(TMath::Min(maxL, 5))),
79 fLmMath() {
80 SetNorm(0, 0.5, 0);
81 if (maxL > 5) { Hal::Cout::PrintInfo("MaxL > 5 not supported", EInfo::kError); }
82 fNumReal = new TH1D*[fMaxJM];
83 fNumImag = new TH1D*[fMaxJM];
84 fDenReal = new TH1D*[fMaxJM];
85 fDenImag = new TH1D*[fMaxJM];
86
87 TString x_name = Hal::Femto::KinematicsToAxisLabel(fFrame, 0, 1);
88 TString y_name = Hal::Femto::KinematicsToAxisLabel(fFrame, 1, 1);
89
90 for (int ihist = 0; ihist < fMaxJM; ihist++) {
91 TString hname = Form("NumReYlm%i%i", fLmVals.GetElsi(ihist), fLmVals.GetEmsi(ihist));
92 fNumReal[ihist] = new TH1D(hname, hname, bins, min, max);
93 hname = Form("NumImYlm%i%i", fLmVals.GetElsi(ihist), fLmVals.GetEmsi(ihist));
94 fNumImag[ihist] = new TH1D(hname, hname, bins, min, max);
95 hname = Form("DenReYlm%i%i", fLmVals.GetElsi(ihist), fLmVals.GetEmsi(ihist));
96 fDenReal[ihist] = new TH1D(hname, hname, bins, min, max);
97 hname = Form("DenImYlm%i%i", (int) fLmVals.GetElsi(ihist), fLmVals.GetEmsi(ihist));
98 fDenImag[ihist] = new TH1D(hname, hname, bins, min, max);
99 fNumReal[ihist]->GetXaxis()->SetTitle(x_name);
100 fNumReal[ihist]->GetYaxis()->SetTitle(y_name);
101 fNumImag[ihist]->GetXaxis()->SetTitle(x_name);
102 fNumImag[ihist]->GetYaxis()->SetTitle(y_name);
103 fDenReal[ihist]->GetXaxis()->SetTitle(x_name);
104 fDenReal[ihist]->GetYaxis()->SetTitle(y_name);
105 fDenImag[ihist]->GetXaxis()->SetTitle(x_name);
106 fDenImag[ihist]->GetYaxis()->SetTitle(y_name);
107 fNumReal[ihist]->Sumw2();
108 fNumImag[ihist]->Sumw2();
109 fDenReal[ihist]->Sumw2();
110 fDenImag[ihist]->Sumw2();
111 }
112
113 fCovNum.MakeBigger(bins, fMaxJM * 2, fMaxJM * 2);
114 fCovDen.MakeBigger(bins, fMaxJM * 2, fMaxJM * 2);
115 fCovCf.MakeBigger(bins, fMaxJM * 2, fMaxJM * 2);
116
117
118 // gSystem->Load("libgsl.so");
119 // gSystem->Load("libgslcblas.so");
120 AddLabel(Hal::Femto::KinematicsToLabel(fFrame));
122 }
123
125 DividedHisto1D(other),
126 fMaxJM(other.fMaxJM),
127 fFrame(other.fFrame),
128 fNormPurity(other.fNormPurity),
129 fNormRadius(other.fNormRadius),
130 fNormBohr(other.fNormBohr),
131 fLmVals(other.fLmVals) {
132 if (other.fNumReal) {
133 fNumReal = new TH1D*[fMaxJM];
134 fNumImag = new TH1D*[fMaxJM];
135 fDenReal = new TH1D*[fMaxJM];
136 fDenImag = new TH1D*[fMaxJM];
137 for (int i = 0; i < fMaxJM; i++) {
138 fNumReal[i] = (TH1D*) other.fNumReal[i]->Clone();
139 fNumImag[i] = (TH1D*) other.fNumImag[i]->Clone();
140 fDenReal[i] = (TH1D*) other.fDenReal[i]->Clone();
141 fDenImag[i] = (TH1D*) other.fDenImag[i]->Clone();
142 }
143 }
144
145 if (other.fCFReal) {
146 fCFReal = new TH1D*[fMaxJM];
147 fCFImag = new TH1D*[fMaxJM];
148 for (int i = 0; i < fMaxJM; i++) {
149 fCFReal[i] = (TH1D*) other.fCFReal[i]->Clone();
150 fCFImag[i] = (TH1D*) other.fCFImag[i]->Clone();
151 }
152 }
153 fCovNum = other.fCovNum;
154 fCovDen = other.fCovDen;
155 fCovCf = other.fCovCf;
156 if (other.fCfcov) { fCfcov = (TH3D*) other.fCfcov->Clone(); }
157 UnpackCovariances();
158 RecalculateCF();
159 }
160
161 TH1D* FemtoSHCF::GetCFRe(int el, int em) const {
162 if (fLmVals.GetIndex(el, em) >= 0 && fCFReal != nullptr) {
163 return fCFReal[fLmVals.GetIndex(el, em)];
164 } else {
165 return nullptr;
166 }
167 }
168
169 TH1D* FemtoSHCF::GetCFIm(int el, int em) const {
170 if (fLmVals.GetIndex(el, em) >= 0 && fCFImag != nullptr) {
171 return fCFImag[fLmVals.GetIndex(el, em)];
172 } else {
173 return nullptr;
174 }
175 }
176
177 TH1D* FemtoSHCF::GetNumRe(int el, int em) const {
178 if (fLmVals.GetIndex(el, em) >= 0)
179 return fNumReal[fLmVals.GetIndex(el, em)];
180 else
181 return nullptr;
182 }
183
184 TH1D* FemtoSHCF::GetNumIm(int el, int em) const {
185 if (fLmVals.GetIndex(el, em) >= 0)
186 return fNumImag[fLmVals.GetIndex(el, em)];
187 else
188 return nullptr;
189 }
190
191 TH1D* FemtoSHCF::GetDenRe(int el, int em) const {
192 if (fLmVals.GetIndex(el, em) >= 0)
193 return fDenReal[fLmVals.GetIndex(el, em)];
194 else
195 return nullptr;
196 }
197
198 TH1D* FemtoSHCF::GetDenIm(int el, int em) const {
199 if (fLmVals.GetIndex(el, em) >= 0)
200 return fDenImag[fLmVals.GetIndex(el, em)];
201 else
202 return nullptr;
203 }
204
205 void FemtoSHCF::SetNumRe(TH1D** histograms, Bool_t clone) {
206 if (clone) {
207 if (fNumReal == nullptr) fNumReal = new TH1D*[fMaxJM];
208 for (int i = 0; i < fMaxJM; i++) {
209 fNumReal[i] = (TH1D*) histograms[i]->Clone();
210 }
211 } else {
212 fNumReal = histograms;
213 }
214 }
215
216 void FemtoSHCF::SetNumIm(TH1D** histograms, Bool_t clone) {
217 if (clone) {
218 if (fNumImag == nullptr) fNumImag = new TH1D*[fMaxJM];
219 for (int i = 0; i < fMaxJM; i++) {
220 fNumImag[i] = (TH1D*) histograms[i]->Clone();
221 }
222 } else {
223 fNumImag = histograms;
224 }
225 }
226
227 void FemtoSHCF::SetDenRe(TH1D** histograms, Bool_t clone) {
228 if (clone) {
229 if (fDenReal == nullptr) fDenReal = new TH1D*[fMaxJM];
230 for (int i = 0; i < fMaxJM; i++) {
231 fDenReal[i] = (TH1D*) histograms[i]->Clone();
232 }
233 } else {
234 fDenReal = histograms;
235 }
236 }
237
238 void FemtoSHCF::SetDenIm(TH1D** histograms, Bool_t clone) {
239 if (clone) {
240 if (fDenImag == nullptr) fDenImag = new TH1D*[fMaxJM];
241 for (int i = 0; i < fMaxJM; i++) {
242 fDenImag[i] = (TH1D*) histograms[i]->Clone();
243 }
244 } else {
245 fDenImag = histograms;
246 }
247 }
248
249 void FemtoSHCF::Draw(Option_t* opt) {
250 TString option = opt;
251 if (fPainter) {
252 fPainter->SetOption(option);
253 fPainter->Paint();
254 } else {
255 fPainter = new Hal::FemtoSHCFPainter(this);
256 fPainter->SetOption(option);
257 fPainter->Paint();
258 }
259 }
260
261 void FemtoSHCF::PackCfcCovariance() {
262 TString bufName;
263 if (fCfcov) delete fCfcov;
264 bufName = Form("CovCfc%s", fNumReal[0]->GetName());
265 fCfcov = new TH3D(bufName,
266 bufName,
267 fCFReal[0]->GetNbinsX(),
268 fCFReal[0]->GetXaxis()->GetXmin(),
269 fCFReal[0]->GetXaxis()->GetXmax(),
270 GetMaxJM() * 2,
271 -0.5,
272 GetMaxJM() * 2 - 0.5,
273 GetMaxJM() * 2,
274 -0.5,
275 GetMaxJM() * 2 - 0.5);
276
277 double tK;
278 for (int ibin = 1; ibin <= fCfcov->GetNbinsX(); ibin++) {
279 const int binCov = ibin - 1;
280 for (int ilmz = 0; ilmz < GetMaxJM() * 2; ilmz++)
281 for (int ilmp = 0; ilmp < GetMaxJM() * 2; ilmp++) {
282 tK = fCovCf[binCov][ilmz][ilmp];
283 fCfcov->SetBinContent(ibin, ilmz + 1, ilmp + 1, tK);
284 }
285 }
286 fCfcov->SetBinContent(0, 0, 0, 1.0);
287 }
288
289 void FemtoSHCF::Browse(TBrowser* b) {
290 if (fCFReal == nullptr) {
291 Cout::PrintInfo("No correlation functions, call RecalculateCF", EInfo::kError);
292 return;
293 }
294 for (int i = 0; i < fMaxJM; i++) {
295 b->Add(fNumReal[i]);
296 b->Add(fNumImag[i]);
297 b->Add(fDenReal[i]);
298 b->Add(fDenImag[i]);
299 }
300 for (int i = 0; i < fMaxJM; i++) {
301 b->Add(fCFReal[i]);
302 b->Add(fCFImag[i]);
303 }
304 b->Add(fCfcov);
305 b->Add(fNum);
306 b->Add(fDen);
307 NormalizeBy(0, 0);
308 Draw();
309 }
310
311 void FemtoSHCF::NormalizeBy(int el, int em, Option_t* opt) {
312 TString option = opt;
313 TH1D* num = nullptr;
314 TH1D* den = nullptr;
315 if (option.EqualTo("re")) {
316 num = GetNumRe(el, em);
317 den = GetDenRe(el, em);
318 } else {
319 num = GetNumIm(el, em);
320 den = GetDenIm(el, em);
321 if (num) num->SetLineColor(kRed);
322 }
323 fScale = 1.0;
324 return;
325 Double_t n_integral, d_integral;
326 if (fNormMax == fNormMin) {
327 n_integral = num->Integral(1, num->GetNbinsX());
328 d_integral = den->Integral(1, den->GetNbinsX());
329 } else {
330 n_integral = num->Integral(num->GetXaxis()->FindBin(fNormMin[0]), num->GetXaxis()->FindBin(fNormMax[0]));
331 d_integral = den->Integral(den->GetXaxis()->FindBin(fNormMin[0]), den->GetXaxis()->FindBin(fNormMax[0]));
332 }
333 fScale = d_integral / n_integral;
334 }
335
336 void FemtoSHCF::AddLabel(TString label) {
337 if (fLabels == nullptr) {
338 fLabels = new TList();
339 fLabels->SetName("Labels");
340 }
341 fLabels->AddLast(new TObjString(label));
342 }
343
344 void FemtoSHCF::Add(const Object* pack) {
345 FemtoSHCF* cf = (FemtoSHCF*) pack;
346 if (cf->fMaxJM != fMaxJM) {
347 Cout::PrintInfo("Not compatible maxL numbers during merging FemtoSHCF", EInfo::kError);
348 return;
349 }
350 if (cf->fFrame != this->fFrame) {
351 Cout::PrintInfo("Incompatible frames !", EInfo::kError);
352 return;
353 }
354 for (int i = 0; i < fMaxJM; i++) {
355 fNumReal[i]->Add(cf->fNumReal[i]);
356 fNumImag[i]->Add(cf->fNumImag[i]);
357 fDenReal[i]->Add(cf->fDenReal[i]);
358 fDenImag[i]->Add(cf->fDenImag[i]);
359 }
360 fNum->Add(cf->fNum);
361 fDen->Add(cf->fDen);
362 for (int i = 0; i < fCovNum.GetSize(); i++) {
363 for (int j = 0; j < fCovNum[0].GetSize(); j++) {
364 for (int k = 0; k < fCovNum[0][0].GetSize(); k++) {
365 fCovNum[i][j][k] = fCovNum[i][j][k] + cf->fCovNum[i][j][k];
366 fCovDen[i][j][k] = fCovDen[i][j][k] + cf->fCovDen[i][j][k];
367 }
368 }
369 }
370 // RecalculateCF(); // DO NOT MERGE !!!
371 }
372
373 TH1D* FemtoSHCF::GetHisto(int el, int em, Bool_t norm, Option_t* opt) const {
374 TString option = opt;
375 TH1D* cf = nullptr;
376 if (option.EqualTo("re")) {
377 cf = GetCFRe(el, em);
378 } else {
379 cf = GetCFIm(el, em);
380 }
381 cf = (TH1D*) cf->Clone();
382 cf->GetYaxis()->SetTitle(Form("C^{%i}_{%i}", el, em));
383 cf->SetTitle(Form("CF^{%i}_{%i}", el, em));
384 cf->SetMinimum(0);
385 cf->SetMaximum(cf->GetBinContent(cf->GetMaximumBin()) * 1.1);
386 if (norm && fScale != 0) cf->Scale(fScale);
387 return cf;
388 }
389
390 TH1D* FemtoSHCF::GetCubic(Option_t* opt, Bool_t R) const {
391 if (fCFReal == nullptr) {
392 Cout::PrintInfo("No correlation functions, call RecalculateCF", EInfo::kError);
393 return nullptr;
394 }
395 TString option = opt;
396 if (option.Length() > GetLMax()) return nullptr;
397 TH1D* cf = nullptr;
398 TH1D* cf2 = nullptr;
399 TH1D* cf3 = nullptr;
400 TH1D* cf4 = nullptr;
401 TH1D* cf5 = nullptr;
402 if (option.Length() == 1) {
403 Double_t scale1 = 2.0 * Sqr(1, 3);
404 Double_t scale2 = Sqr(2, 3);
405 // Double_t scale4 = Sqr(1,4);
406 if (option == "z") {
407 cf = Histo(1, 0, "re", scale1);
408 } else if (option == "y") {
409 cf = Histo(1, 1, "im", scale2);
410 cf2 = Histo(1, -1, "im", scale2);
411 } else if (option == "x") {
412 cf = Histo(1, 1, "re", -scale2);
413 cf2 = Histo(1, -1, "re", scale2);
414 } else if (option == "0") {
415 cf = Histo(0, 0, "re", 1.0);
416 }
417 }
418 if (option.Length() == 2) {
419 Double_t scale = Sqr(2, 15);
420 if (option == "zz") {
421 cf = Histo(2, 0, "re", 4.0 / 3.0 * Sqr(1, 5));
422 } else if (option == "yz") {
423 cf = Histo(2, 1, "im", scale);
424 cf2 = Histo(2, -1, "im", scale);
425 } else if (option == "yy") {
426 cf = Histo(2, 2, "re", -scale);
427 cf2 = Histo(2, 0, "re", -2.0 / 3.0 * Sqr(1, 5));
428 cf3 = Histo(2, -2, "re", -scale);
429 } else if (option == "xz") {
430 cf = Histo(2, 1, "re", -scale);
431 cf2 = Histo(2, -1, "re", scale);
432 } else if (option == "xy") {
433 cf = Histo(2, 2, "im", -scale);
434 cf2 = Histo(2, -2, "im", scale);
435 } else if (option == "xx") {
436 cf = Histo(2, 2, "re", scale);
437 cf2 = Histo(2, 0, "re", -2.0 / 3.0 * Sqr(1, 5));
438 cf3 = Histo(2, -2, "re", scale);
439 }
440 }
441 if (option.Length() == 3) {
442 Double_t scale7 = 1.0 / 5.0 * Sqr(1, 7);
443 Double_t scale21 = 1.0 / 5.0 * Sqr(1, 21);
444 Double_t scale105 = Sqr(2, 105);
445 Double_t scale35 = Sqr(1, 35);
446 if (option == "zzz") {
447 cf = Histo(3, 0, "re", 4.0 * scale7);
448 } else if (option == "yzz") {
449 cf = Histo(3, 1, "im", 4.0 * scale21);
450 cf2 = Histo(3, -1, "im", 4.0 * scale21);
451 } else if (option == "yyz") {
452 cf = Histo(3, 2, "re", -scale105);
453 cf2 = Histo(3, 0, "re", -2.0 * scale7);
454 cf3 = Histo(3, -2, "re", -scale105);
455 } else if (option == "yyy") {
456 cf = Histo(3, 3, "im", -scale35);
457 cf2 = Histo(3, 1, "im", -Sqr(3, 7) / 5.0);
458 cf3 = Histo(3, -1, "im", -Sqr(3, 7) / 5.0);
459 cf4 = Histo(3, -3, "im", -scale35);
460 } else if (option == "xzz") {
461 cf = Histo(3, 1, "re", -4.0 * scale21);
462 cf2 = Histo(3, -1, "re", 4.0 * scale21);
463 } else if (option == "xyz") {
464 cf = Histo(3, 2, "im", -scale105);
465 cf2 = Histo(3, -2, "im", scale105);
466 } else if (option == "xyy") {
467 cf = Histo(3, 3, "re", scale35);
468 cf2 = Histo(3, 1, "re", scale21);
469 cf3 = Histo(3, -1, "re", -scale21);
470 cf4 = Histo(3, -3, "re", -scale35);
471 } else if (option == "xxz") {
472 cf = Histo(3, 2, "re", scale105);
473 cf2 = Histo(3, 0, "re", -2.0 * scale7);
474 cf3 = Histo(3, -2, "re", scale105);
475 } else if (option == "xxy") {
476 cf = Histo(3, 3, "im", scale35);
477 cf2 = Histo(3, 1, "im", -scale21);
478 cf3 = Histo(3, -1, "im", -scale21);
479 cf4 = Histo(3, -3, "im", scale35);
480 } else if (option == "xxx") {
481 cf = Histo(3, 3, "re", -scale35);
482 cf2 = Histo(3, 1, "re", scale7);
483 cf3 = Histo(3, -1, "re", -scale7);
484 cf4 = Histo(3, -3, "re", scale35);
485 }
486 }
487 if (option.Length() == 4) {
488 Double_t scale5 = Sqr(1, 5) / 21.0;
489 Double_t scale25 = Sqr(2, 5) / 21.0;
490 Double_t scale35 = Sqr(1, 35) / 3.0;
491 Double_t scale235 = Sqr(2, 35) / 3.0;
492 Double_t scale57 = Sqr(1, 5) / 7.0;
493 Double_t pis = TMath::Sqrt(TMath::Pi());
494 if (option == "zzzz") {
495 cf = Histo(4, 0, "re", TMath::Sqrt(TMath::Pi()) * 16.0 / 105.0);
496 } else if (option == "yzzz") {
497 cf = Histo(4, 1, "im", scale5 * 4);
498 cf2 = Histo(4, -1, "im", scale5 * 4.0);
499 } else if (option == "yyzz") {
500 cf = Histo(4, 2, "re", -2.0 * scale25);
501 cf2 = Histo(4, 0, "re", -8.0 / 105.0 * pis);
502 cf3 = Histo(4, -2, "re", -2.0 * scale25);
503 } else if (option == "yyyz") {
504 cf = Histo(4, 3, "im", -scale35);
505 cf2 = Histo(4, 1, "im", -scale57);
506 cf3 = Histo(4, -1, "im", -scale57);
507 cf4 = Histo(4, -3, "im", -scale35);
508 } else if (option == "yyyy") {
509 cf = Histo(4, 4, "re", scale235);
510 cf2 = Histo(4, 2, "re", scale25 * 2);
511 cf3 = Histo(4, 0, "re", 2.0 / 35.0 * pis);
512 cf4 = Histo(4, 2, "re", scale25 * 2);
513 cf5 = Histo(4, -4, "re", scale235);
514 } else if (option == "xzzz") {
515 cf = Histo(4, 1, "re", -4.0 * scale5);
516 cf2 = Histo(4, -1, "re", 4.0 * scale5);
517 } else if (option == "xyzz") {
518 cf = Histo(4, 2, "im", -2.0 * scale25);
519 cf2 = Histo(4, -2, "im", 2.0 * scale25);
520 } else if (option == "xyzz") {
521 cf = Histo(4, 3, "re", scale35);
522 cf2 = Histo(4, 1, "re", scale5);
523 cf3 = Histo(4, -1, "re", -scale5);
524 cf4 = Histo(4, -3, "re", -scale35);
525 } else if (option == "xyyz") {
526 cf = Histo(4, 3, "re", scale35);
527 cf2 = Histo(4, 1, "re", scale5);
528 cf3 = Histo(4, -1, "re", -scale5);
529 cf4 = Histo(4, -3, "re", -scale35);
530 } else if (option == "xyyy") {
531 cf = Histo(4, 4, "im", scale235);
532 cf2 = Histo(4, 2, "im", scale25);
533 cf3 = Histo(4, -2, "im", -scale25);
534 cf4 = Histo(4, -4, "im", -scale235);
535 } else if (option == "xxzz") {
536 cf = Histo(4, 2, "re", 2.0 * scale25);
537 cf2 = Histo(4, 0, "re", -8.0 / 105.0 * pis);
538 cf3 = Histo(4, -2, "re", 2.0 * scale25);
539 } else if (option == "xxyz") {
540 cf = Histo(4, 3, "im", scale35);
541 cf2 = Histo(4, 1, "im", -scale5);
542 cf3 = Histo(4, -1, "im", -scale5);
543 cf4 = Histo(4, -3, "im", scale35);
544 } else if (option == "xxyy") {
545 cf = Histo(4, 4, "re", -scale235);
546 cf2 = Histo(4, 0, "re", 2.0 / 105. * pis);
547 cf3 = Histo(4, 4, "re", -scale235);
548 } else if (option == "xxxz") {
549 cf = Histo(4, 3, "re", -scale35);
550 cf2 = Histo(4, 1, "re", scale57);
551 cf3 = Histo(4, -1, "re", -scale57);
552 cf4 = Histo(4, -3, "re", scale35);
553 } else if (option == "xxxy") {
554 cf = Histo(4, 4, "im", -scale235);
555 cf2 = Histo(4, 2, "im", scale25);
556 cf3 = Histo(4, -2, "im", -scale25);
557 cf4 = Histo(4, -4, "im", scale235);
558 } else if (option == "xxxx") {
559 cf = Histo(4, 4, "re", scale235);
560 cf2 = Histo(4, 2, "re", -2.0 * scale25);
561 cf3 = Histo(4, 0, "re", 2.0 / 35.0 * pis);
562 cf4 = Histo(4, -2, "re", -2.0 * scale25);
563 cf5 = Histo(4, -4, "re", scale235);
564 }
565 }
566 if (cf == nullptr) return nullptr;
567 if (cf2) {
568 cf->Add(cf2);
569 delete cf2;
570 if (cf3) {
571 cf->Add(cf3);
572 delete cf3;
573 if (cf4) {
574 cf->Add(cf4);
575 delete cf4;
576 if (cf5) {
577 cf->Add(cf5);
578 delete cf5;
579 }
580 }
581 }
582 }
583 cf->SetTitle(option);
584 cf->SetName(option);
585 if (R) {
586 cf->SetTitle(option + " R(q)");
587 for (int i = 1; i <= cf->GetNbinsX(); i++) {
588 cf->SetBinContent(i, cf->GetBinContent(i) - fScale);
589 }
590 }
591 return cf;
592 }
593
594 Double_t FemtoSHCF::Sqr(Double_t val1, Double_t val2) const { return TMath::Sqrt(TMath::Pi() * val1 / val2); }
595
596 TH1D* FemtoSHCF::Histo(int lm, int em, Option_t* opt, Double_t scale) const {
597 TH1D* cf;
598 TString option = opt;
599 if (option.EqualTo("im")) {
600 cf = (TH1D*) GetCFIm(lm, TMath::Abs(em))->Clone();
601 } else {
602 cf = (TH1D*) GetCFRe(lm, TMath::Abs(em))->Clone();
603 }
604 if (option == "im") {
605 cf->Scale(-scale); // multiply by i
606 } else {
607 cf->Scale(scale);
608 }
609 return cf;
610 }
611
612 TH3D* FemtoSHCF::GetCovCF() const { return fCfcov; }
613
614 void FemtoSHCF::RecalculateCF(Int_t debug, Bool_t suwm) {
615#ifndef DISABLE_GSL
616 if (fDenImag == nullptr) {
617 Cout::PrintInfo("No imaginary denominators!", EInfo::kError);
618 return;
619 }
620 if (fDenReal == nullptr) {
621 Cout::PrintInfo("No real denominators", EInfo::kError);
622 return;
623 }
624 if (fNumReal == nullptr) {
625 Cout::PrintInfo("No real numeraotrs", EInfo::kError);
626 return;
627 }
628 if (fNumImag == nullptr) {
629 Cout::PrintInfo("No imaginary numerators", EInfo::kError);
630 return;
631 }
632
633 Cout::PrintInfo("Computing SF correlation functions", EInfo::kDebugInfo);
634 if (fCFReal) {
635 for (int i = 0; i < fMaxJM; i++)
636 delete fCFReal[i];
637 delete[] fCFReal;
638 }
639 if (fCFImag) {
640 for (int i = 0; i < fMaxJM; i++)
641 delete fCFImag[i];
642 delete[] fCFImag;
643 }
644 fCFReal = new TH1D*[fMaxJM];
645 fCFImag = new TH1D*[fMaxJM];
646 // UnpackCovariances(); do nothing
647 for (int i = 0; i < fMaxJM; i++) {
648 TString name = fNumReal[i]->GetName();
649 name.ReplaceAll("Num", "CF");
650 Int_t nbins;
651 Double_t min, max;
652 Std::GetAxisPar(*fNumReal[i], nbins, min, max, "x");
653 fCFReal[i] = new TH1D(name, name, nbins, min, max);
654 name = fNumImag[i]->GetName();
655 name.ReplaceAll("Num", "CF");
656 fCFImag[i] = new TH1D(name, name, nbins, min, max);
657 fCFImag[i]->SetLineColor(kRed);
658 fCFImag[i]->SetMarkerColor(kRed);
659 int el, em;
660 el = fLmVals.GetElsi(i);
661 em = fLmVals.GetEmsi(i);
662 fCFReal[i]->GetYaxis()->SetTitle(Form("C^{%i}_{%i}", el, em));
663 fCFReal[i]->GetXaxis()->SetTitle(Hal::Femto::KinematicsToAxisLabel(fFrame, 0, 1));
664 fCFReal[i]->SetTitle(Form("CF^{%i}_{%i} Re", el, em));
665 fCFImag[i]->GetYaxis()->SetTitle(Form("C^{%i}_{%i}", el, em));
666 fCFImag[i]->SetTitle(Form("CF^{%i}_{%i} Im", el, em));
667 fCFImag[i]->GetXaxis()->SetTitle(Hal::Femto::KinematicsToAxisLabel(fFrame, 0, 1));
668 }
669 if ((fCovNum[0][0][0] > 0.0)) {
670 std::cout << "Detected calculated covariance matrix. Do not recalculate !!!" << std::endl;
671 // recalccov = 0;
672 }
673 if (fNum->GetEntries() == 0) {
674 std::cout << "Emtpy numerator, stop calculating cF" << std::endl;
675 return; // no data, CF is empdy
676 }
677
678 //=====================================================
679 FemtoYlmSolver solver(GetLMax(), this, suwm);
680 solver.SetDebugBin(debug);
681 solver.SetNormalizationArea(GetNormMin(0), GetNormMax(0));
682 solver.Solve(kTRUE);
683 PackCfcCovariance();
684#else
685 std::cout << "GLS NOT ENABLED!" << std::endl;
686#endif
687 }
688
689 Long64_t FemtoSHCF::Merge(TCollection* collection) {
690 if (collection) {
691 FemtoSHCF* pack = nullptr;
692 TIter iterator(collection);
693 while ((pack = (FemtoSHCF*) iterator())) {
694 FastAdd(pack);
695 }
696 }
697 RecalculateCF();
698 return 1;
699 }
700
701 void FemtoSHCF::FastAdd(const FemtoSHCF* obj) {
702 FemtoSHCF* cf = (FemtoSHCF*) obj;
703 if (cf->fMaxJM != fMaxJM) {
704 Cout::PrintInfo("Not compatible maxL numbers during merging FemtoSHCF", EInfo::kError);
705 return;
706 }
707 if (cf->fFrame != this->fFrame) {
708 Cout::PrintInfo("Incompatible frames !", EInfo::kError);
709 return;
710 }
711 for (int i = 0; i < fMaxJM; i++) {
712 fNumReal[i]->Add(cf->fNumReal[i]);
713 fNumImag[i]->Add(cf->fNumImag[i]);
714 fDenReal[i]->Add(cf->fDenReal[i]);
715 fDenImag[i]->Add(cf->fDenImag[i]);
716 }
717 fCovCf += obj->fCovCf;
718 fCovNum += obj->fCovNum;
719 fCovDen += obj->fCovDen;
720 fNum->Add(cf->fNum);
721 fDen->Add(cf->fDen);
722 }
723
724 FemtoSHCF::~FemtoSHCF() {
725 auto cleanArray = [&](TH1D** x) {
726 if (x == nullptr) return;
727 for (int i = 0; i < fMaxJM; i++)
728 delete x[i];
729 delete[] x;
730 };
731 cleanArray(fNumImag);
732 cleanArray(fNumReal);
733 cleanArray(fDenReal);
734 cleanArray(fDenImag);
735 cleanArray(fCFReal);
736 cleanArray(fCFImag);
737 if (fCfcov) delete fCfcov;
738 if (fPainter) delete fPainter; // KURWA
739 }
740
741 TString FemtoSHCF::HTMLExtract(Int_t counter, TString dir) const {
742 gSystem->MakeDirectory(Form("%s/divided_%i", dir.Data(), counter));
743 TString filename = Form("%s/divided_%i/divided.html", dir.Data(), counter);
744 HtmlFile file(filename, kFALSE);
745 HtmlTable table("", "haltable", "");
746 HtmlRow row1("", "dark_blue", "");
747 row1.AddContent(HtmlCell("label"));
748 row1.AddContent(HtmlCellCol("value", 2));
749 table.AddContent(row1);
750 HtmlRow row2("", "dark_blue", "");
751 row2.AddContent(HtmlCell("Name"));
752 row2.AddContent(HtmlCellCol(GetName(), 2));
753 table.AddContent(row2);
754 HtmlRow row3("", "dark_blue", "");
755 row3.AddContent(HtmlCell("Comment"));
756 row3.AddContent(HtmlCellCol(fComment, 2));
757 table.AddContent(row3);
758 if (fScale == 0) {
759 HtmlRow row4("", "dark_blue", "");
760 row4.AddContent(HtmlCell("Scale"));
761 row4.AddContent(HtmlCellCol("0(disabled)", 2));
762 table.AddContent(row4);
763 } else {
764 HtmlRow row4("", "dark_blue", "");
765 row4.AddContent(HtmlCell("Norm on xaxis"));
766 row4.AddContent(HtmlCell(Form("%4.2f", fNormMin[0])));
767 row4.AddContent(HtmlCell(Form("%4.2f", fNormMax[0])));
768 table.AddContent(row4);
769 }
770 file.AddContent(table);
771 HtmlTable table2("", "haltable", "");
772 HtmlRow row5("", "dark_blue", "");
773 row5.AddContent(HtmlCell("No."));
774 row5.AddContent(HtmlCell("Name/Value"));
775 table2.AddContent(row5);
776
777 if (fLabels) // for backward compatilibiyt
778 for (int i = 0; i < fLabels->GetEntries(); i++) {
779 TObjString* obj = (TObjString*) fLabels->At(i);
780 HtmlRow row6("", "light_blue", "");
781 row6.AddContent(HtmlCell(Form("%i", i)));
782 row6.AddContent(HtmlCell(obj->GetString()));
783 table2.AddContent(row6);
784 }
785 file.AddContent(table2);
786 // file<<"<img border=\"0\" src=\""<<Form("divided.png")<<"\" width=\"996\"
787 // height=\"1472\">"<<std::endl;
788 HtmlTable table3;
789 Bool_t batch = gROOT->IsBatch();
790 gROOT->SetBatch(kTRUE);
791 Int_t jmax = GetLMax() * 2 + 1;
792 Int_t width = 966 / jmax;
793 Int_t height = 966 / GetLMax();
794 if (height > width) height = width;
795 for (int i = 0; i <= GetLMax(); i++) {
796 HtmlRow row7("", "dark_blue", "");
797
798 int I = i;
799 for (int j = 0; j <= GetLMax(); j++) {
800 int J = j - i;
801 if (j <= i) I = i;
802 TH1D* cfr = GetHisto(I, J, kFALSE, "re");
803 TH1D* cfi = GetHisto(I, J, kFALSE, "im");
804 if (cfi) cfi->SetLineColor(kRed);
805 HtmlCell cell;
806 if (cfi && cfr) {
807 TCanvas* c1 = new TCanvas();
808 cfr->Draw();
809 cfi->Draw("SAME");
810 c1->SaveAs(Form("%s/divided_%i/C_%i%i.png", dir.Data(), counter, I, J));
811 TString cell_content = Form("<a href=\"C_%i%i.png\">", I, J);
812 cell_content = cell_content
813 + Form("<img border=\"0\" src=\"C_%i%i.png\" "
814 "width=\"%i\" height=\"%i\"></a>",
815 I,
816 J,
817 width,
818 height);
819 delete c1;
820 cell.SetStringContent(cell_content);
821 }
822 row7.AddContent(cell);
823 I++;
824 }
825 file.AddContent(row7);
826 }
827
828 file.AddContent(table3);
829 file.Save();
830 gROOT->SetBatch(batch);
831 return HtmlCore::GetUrl(Form("divided_%i/divided.html", counter), this->ClassName());
832 }
833
834 void FemtoSHCF::FillNumObj(TObject* obj) {
835 FemtoPair* pair = (FemtoPair*) obj;
836 Double_t kv = TMath::Sqrt(pair->GetX() * pair->GetX() + pair->GetY() * pair->GetY() + pair->GetZ() * pair->GetZ());
837 std::complex<double>* YlmBuffer = fLmMath.YlmUpToL(fLmVals.GetMaxL(), pair->GetX(), pair->GetY(), pair->GetZ());
838 for (int ilm = 0; ilm < GetMaxJM(); ilm++) {
839 fNumReal[ilm]->Fill(kv, real(YlmBuffer[ilm]) * pair->GetWeight());
840 fNumImag[ilm]->Fill(kv, -imag(YlmBuffer[ilm]) * pair->GetWeight());
841 fNum->Fill(kv, 1.0); //?
842 }
843
844 // Fill in the error matrix
845 int nqbin = fNum->GetXaxis()->FindFixBin(kv) - 1;
846 if (nqbin < fNum->GetNbinsX()) {
847 Double_t weight2 = pair->GetWeight() * pair->GetWeight();
848 for (int ilmzero = 0; ilmzero < GetMaxJM(); ilmzero++) {
849 const int twoilmzero = ilmzero * 2;
850 for (int ilmprim = 0; ilmprim < GetMaxJM(); ilmprim++) {
851 const int twoilmprim = ilmprim * 2;
852 fCovNum[nqbin][twoilmzero][twoilmprim] += real(YlmBuffer[ilmzero]) * real(YlmBuffer[ilmprim]) * weight2;
853 fCovNum[nqbin][twoilmzero][twoilmprim + 1] += real(YlmBuffer[ilmzero]) * -imag(YlmBuffer[ilmprim]) * weight2;
854 fCovNum[nqbin][twoilmzero + 1][twoilmprim] -= imag(YlmBuffer[ilmzero]) * real(YlmBuffer[ilmprim]) * weight2;
855 fCovNum[nqbin][twoilmzero + 1][twoilmprim + 1] -= imag(YlmBuffer[ilmzero]) * -imag(YlmBuffer[ilmprim]) * weight2;
856 }
857 }
858 }
859 }
860
861 void FemtoSHCF::FillDenObj(TObject* obj) {
862 FemtoPair* pair = (FemtoPair*) obj;
863 Double_t kv = TMath::Sqrt(pair->GetX() * pair->GetX() + pair->GetY() * pair->GetY() + pair->GetZ() * pair->GetZ());
864 std::complex<double>* YlmBuffer = fLmMath.YlmUpToL(fLmVals.GetMaxL(), pair->GetX(), pair->GetY(), pair->GetZ());
865 for (int ilm = 0; ilm < GetMaxJM(); ilm++) {
866 fDenReal[ilm]->Fill(kv, real(YlmBuffer[ilm]) * pair->GetWeight());
867 fDenImag[ilm]->Fill(kv, -imag(YlmBuffer[ilm]) * pair->GetWeight());
868 fDen->Fill(kv, 1.0);
869 }
870
871 // Fill in the error matrix
872 int nqbin = fDen->GetXaxis()->FindFixBin(kv) - 1;
873 if (nqbin < fDen->GetNbinsX()) {
874 Double_t weight2 = pair->GetWeight() * pair->GetWeight();
875 for (int ilmzero = 0; ilmzero < GetMaxJM(); ilmzero++) {
876 const int twoilmzero = ilmzero * 2;
877 for (int ilmprim = 0; ilmprim < GetMaxJM(); ilmprim++) {
878 const int twoilmprim = ilmprim * 2;
879 fCovDen[nqbin][twoilmzero][twoilmprim] += real(YlmBuffer[ilmzero]) * real(YlmBuffer[ilmprim]) * weight2;
880 fCovDen[nqbin][twoilmzero][twoilmprim + 1] += real(YlmBuffer[ilmzero]) * -imag(YlmBuffer[ilmprim]) * weight2;
881 fCovDen[nqbin][twoilmzero + 1][twoilmprim] -= imag(YlmBuffer[ilmzero]) * real(YlmBuffer[ilmprim]) * weight2;
882 fCovDen[nqbin][twoilmzero + 1][twoilmprim + 1] -= imag(YlmBuffer[ilmzero]) * -imag(YlmBuffer[ilmprim]) * weight2;
883 }
884 }
885 }
886 }
887
888 Array_1<Float_t>* FemtoSHCF::ExportToFlatNum() const {
889 Array_1<Float_t>* data = new Array_1<Float_t>(fNum->GetNbinsX() * fMaxJM * 2);
890 ExportIntoToFlatNum(data);
891 return data;
892 }
893
894 void FemtoSHCF::ExportIntoToFlatNum(Array_1<Float_t>* output) const {
895 Int_t bin = 0;
896 Int_t nbins = GetNum()->GetNbinsX();
897 output->MakeBigger(fNum->GetNbinsX() * fMaxJM * 2);
898 for (int ibin = 1; ibin <= nbins; ibin++) {
899 for (int l = 0; l <= GetLMax(); l++) {
900 for (int m = -l; m <= l; m++) {
901 TH1D* h = GetCFRe(l, m);
902 (*output)[bin++] = h->GetBinContent(ibin);
903 h = GetCFIm(l, m);
904 (*output)[bin++] = h->GetBinContent(ibin);
905 }
906 }
907 }
908 }
909
910 void FemtoSHCF::ExportIntoToFlatNumValkyria(Array_1<Float_t>* output) const {
911 Int_t bin = 0;
912 Int_t nbins = GetNum()->GetNbinsX();
913 output->MakeBigger(fNum->GetNbinsX() * fMaxJM * 2);
914 for (int ibin = 1; ibin <= nbins; ibin++) {
915 for (int l = 0; l <= GetLMax(); l++) {
916 for (int m = -l; m <= l; m++) {
917 TH1D* h = GetNumRe(l, m);
918 (*output)[bin++] = h->GetBinContent(ibin);
919 h = GetNumIm(l, m);
920 (*output)[bin++] = h->GetBinContent(ibin);
921 }
922 }
923 }
924 }
925
926 void FemtoSHCF::MakeDummyCov() {
927 Cout::Text("FemtoSHCF::MakeDummyCov", "M", kRed);
928 double nent = fNum->GetEntries();
929 double nentd = fDen->GetEntries();
930 for (int iBin = 1; iBin <= fNum->GetNbinsX(); iBin++) {
931 const int bin = iBin - 1;
932 for (int ilmx = 0; ilmx < GetMaxJM(); ilmx++) {
933 const int ilmx2 = ilmx * 2;
934 for (int ilmy = 0; ilmy < GetMaxJM(); ilmy++) {
935 const int ilmy2 = ilmy * 2;
936 double t1t2rr = fNumReal[ilmx]->GetBinContent(iBin) * fNumReal[ilmy]->GetBinContent(iBin) / nent / nent;
937 double t1t2ri = fNumReal[ilmx]->GetBinContent(iBin) * fNumImag[ilmy]->GetBinContent(iBin) / nent / nent;
938 double t1t2ir = fNumImag[ilmx]->GetBinContent(iBin) * fNumReal[ilmy]->GetBinContent(iBin) / nent / nent;
939 double t1t2ii = fNumImag[ilmx]->GetBinContent(iBin) * fNumImag[ilmy]->GetBinContent(iBin) / nent / nent;
940 if (ilmx == ilmy) {
941 fCovNum[bin][ilmx2][ilmy2] = nent * (TMath::Power(fNumReal[ilmx]->GetBinError(iBin) / nent, 2) * (nent - 1) + t1t2rr);
942 fCovNum[bin][ilmx2][ilmy2 + 1] = nent * t1t2ri;
943 fCovNum[bin][ilmx2 + 1][ilmy2] = nent * t1t2ir;
944 fCovNum[bin][ilmx2 + 1][ilmy2 + 1] =
945 nent * (TMath::Power(fNumImag[ilmx]->GetBinError(iBin) / nent, 2) * (nent - 1) + t1t2rr);
946 } else {
947 fCovNum[bin][ilmx2][ilmy2] = nent * t1t2rr;
948 fCovNum[bin][ilmx2][ilmy2 + 1] = nent * t1t2ri;
949 fCovNum[bin][ilmx2 + 1][ilmy2] = nent * t1t2ir;
950 fCovNum[bin][ilmx2 + 1][ilmy2 + 1] = nent * t1t2ii;
951 }
952 t1t2rr = fDenReal[ilmx]->GetBinContent(iBin) * fDenReal[ilmy]->GetBinContent(iBin) / nentd / nentd;
953 t1t2ri = fDenReal[ilmx]->GetBinContent(iBin) * fDenImag[ilmy]->GetBinContent(iBin) / nentd / nentd;
954 t1t2ir = fDenImag[ilmx]->GetBinContent(iBin) * fDenReal[ilmy]->GetBinContent(iBin) / nentd / nentd;
955 t1t2ii = fDenImag[ilmx]->GetBinContent(iBin) * fDenImag[ilmy]->GetBinContent(iBin) / nentd / nentd;
956
957
958 fCovDen[bin][ilmx2 + 0][ilmy2 + 0] = nentd * t1t2rr;
959 fCovDen[bin][ilmx2 + 0][ilmy2 + 1] = nentd * t1t2ri;
960 fCovDen[bin][ilmx2 + 1][ilmy2 + 0] = nentd * t1t2ir;
961 fCovDen[bin][ilmx2 + 1][ilmy2 + 1] = nentd * t1t2ii;
962 }
963 }
964 }
965
966
967 return;
968 for (int i = 0; i < fNum->GetNbinsX(); i++) {
969 for (int ilmzero = 0; ilmzero < GetMaxJM(); ilmzero++) {
970 for (int ilmzero2 = 0; ilmzero2 < GetMaxJM(); ilmzero2++) {
971 // Int_t gbin = GetBin(i, ilmzero, 0, ilmzero2, 1);
972 // (*covmnum)[gbin] = 0;
973 }
974 }
975 for (int ilmzero = 0; ilmzero < GetMaxJM(); ilmzero++) {
976 // Int_t gbin = GetBin(i, ilmzero, 0, ilmzero, 1);
977 // (*covmnum)[gbin] = 1;
978 }
979 }
980 }
981
982 void FemtoSHCF::Rebin(Int_t ngroup, Option_t* opt) {
983
984 int oldBins = fNumReal[0]->GetNbinsX();
985 for (int i = 0; i < fMaxJM; i++) {
986 fNumReal[i]->Rebin(ngroup);
987 fNumImag[i]->Rebin(ngroup);
988 fDenReal[i]->Rebin(ngroup);
989 fDenImag[i]->Rebin(ngroup);
990 }
991 Int_t newBins = oldBins / ngroup;
992 Array_3<Double_t> CovNum, CovDen, CovCf;
993 CovNum.MakeBigger(oldBins, fMaxJM * 2, fMaxJM * 2);
994 CovDen.MakeBigger(oldBins, fMaxJM * 2, fMaxJM * 2);
995 CovCf.MakeBigger(oldBins, fMaxJM * 2, fMaxJM * 2);
996 for (int i = 0; i < newBins; i++) {
997 for (int a = 0; a < fMaxJM * 2; a++) {
998 for (int b = 0; b < fMaxJM * 2; b++) {
999 for (int j = 0; j < ngroup; j++) {
1000 CovNum[i][a][b] = CovNum[i][a][b] + fCovNum[i * ngroup + j][a][b];
1001 CovDen[i][a][b] = CovDen[i][a][b] + fCovDen[i * ngroup + j][a][b];
1002 CovCf[i][a][b] = CovCf[i][a][b] + fCovCf[i * ngroup + j][a][b];
1003 }
1004 }
1005 }
1006 }
1007 fCovNum = CovNum;
1008 fCovDen = CovDen;
1009 fCovCf = CovCf;
1010 if (fCFReal[0]) { RecalculateCF(); }
1011 }
1012#ifdef __CIA__
1013 void FemtoSHCF::Fit(CorrFitSHCF* fit) { fit->Fit(this); }
1014
1015 void FemtoSHCF::FitDummy(CorrFitSHCF* fit) { fit->FitDummy(this); }
1016
1017 CorrFitMaskSH FemtoSHCF::MakeEmptyMask() const {
1018 CorrFitMaskSH mask;
1019 mask.Build(*this);
1020 return mask;
1021 }
1022#endif
1023
1024 void FemtoSHCF::ImportSlice(Array_1<Float_t>* array, Int_t bin) {
1025 if (bin < 0) {
1026 return;
1027 } else {
1028 int count = 0;
1029 fNum->SetBinContent(bin, array->Get(count++));
1030 fDen->SetBinContent(bin, array->Get(count++));
1031 for (int i = 0; i < fMaxJM; i++) {
1032 fNumReal[i]->SetBinContent(bin, array->Get(count++));
1033 fDenReal[i]->SetBinContent(bin, array->Get(count++));
1034 fNumImag[i]->SetBinContent(bin, array->Get(count++));
1035 fDenImag[i]->SetBinContent(bin, array->Get(count++));
1036 }
1037 for (int i = 0; i < fMaxJM * 2; i++) {
1038 for (int j = 0; j < fMaxJM * 2; j++) {
1039 fCovNum[bin - 1][i][j] = array->Get(count++);
1040 }
1041 }
1042 }
1043 }
1044
1045 void FemtoSHCF::ApplyStyle(const Hal::HistoStyle& h) {
1046 Hal::HistoStyle anti = h;
1047 anti.SetAntiColor();
1049 Hal::HistoStyle copy = h; // because of const
1050 for (int i = 0; i < fMaxJM; i++) {
1051 h.Apply(*fNumReal[i]);
1052 h.Apply(*fDenReal[i]);
1053 anti.Apply(*fNumImag[i]);
1054 anti.Apply(*fDenImag[i]);
1055 if (fCFReal[i]) h.Apply(*fCFReal[i]);
1056 if (fCFImag[i]) anti.Apply(*fCFImag[i]);
1057 if (i == 0
1058 && (copy.GetMarkerStyle().Find(Hal::MarkerStyle::kColor)
1059 || copy.GetLineStyle().Find(Hal::LineStyle::kColor))) { // cols set by style, overwrite colz that are use by draw
1060 fColzSet = kTRUE;
1061 fColRe = fCFReal[i]->GetMarkerColor();
1062 fColIm = fCFImag[i]->GetMarkerColor();
1063 }
1064 }
1065 }
1066
1067 TObject* FemtoSHCF::GetSpecial(TString opt) const {
1068 if (opt == "serialization") return new FemtoSerializationInterfaceSH();
1069 if (opt == "painter") return fPainter;
1070 return nullptr;
1071 }
1072
1073} // namespace Hal
void MakeBigger(Int_t new_dim)
Definition Array.cxx:5
void MakeBigger(Int_t sizeA, Int_t sizeB, Int_t sizeC)
Definition Array.cxx:263
static void PrintInfo(TString text, Hal::EInfo status)
Definition Cout.cxx:370
virtual void ApplyStyle(const HistoStyle &h)
void SetNorm(Double_t min, Double_t max, Int_t axis=0)
Double_t GetY() const
Definition FemtoPair.h:312
Double_t GetWeight() const
Definition FemtoPair.h:282
Double_t GetZ() const
Definition FemtoPair.h:317
Double_t GetX() const
Definition FemtoPair.h:307
void RecalculateCF(Int_t debugBin=-1, Bool_t suwm=kFALSE)
void AddLabel(TString label)
Int_t GetEmsi(Int_t i) const
Int_t GetElsi(Int_t i) const
void SetAntiColor(Bool_t safe=kFALSE)
virtual void AddContent(const HtmlObject &obj)
Definition HtmlTable.cxx:29
virtual void AddContent(const HtmlObject &obj)
Definition HtmlTable.cxx:17
virtual void Add(const Object *pack)
Definition Object.cxx:24
Bool_t Find(Int_t bit) const
Definition Style.h:55