CbmRoot
CbmBlastWave.cxx
Go to the documentation of this file.
1 /*
2  *====================================================================
3  *
4  * CBM Blast Wave model analysis
5  *
6  * Authors: V.Vovchenko
7  *
8  * e-mail : new
9  *
10  *====================================================================
11  *
12  * Blast Wave parameters extraction
13  *
14  *====================================================================
15  */
16 
17 #include "CbmBlastWave.h"
18 #include "CbmL1Def.h"
19 
20 
21 //#include "KFParticleFinder.h"
22 //#include "KFParticleSIMD.h"
23 #include "CbmKFTrack.h"
24 #include "CbmKFVertex.h"
25 #include "CbmStsTrack.h"
26 
27 #include "TClonesArray.h"
28 #include "TDirectory.h"
29 #include "TFile.h"
30 #include "TMath.h"
31 
32 #include "TCanvas.h"
33 #include "TGraphErrors.h"
34 #include "TH1F.h"
35 #include "TH2F.h"
36 #include "TSpline.h"
37 
38 #include "L1Field.h"
39 
40 #include "CbmVertex.h"
41 
42 #include "TStopwatch.h"
43 #include <cmath>
44 #include <cstdio>
45 #include <iostream>
46 #include <map>
47 
48 //#include "CbmTrackMatch.h"
49 #include "CbmMCTrack.h"
50 #include "CbmTrackMatchNew.h"
51 //for particle ID from global track
52 #include "CbmGlobalTrack.h"
53 #include "CbmRichRing.h"
54 #include "CbmTofHit.h"
55 #include "CbmTrdTrack.h"
56 #include "TDatabasePDG.h"
57 //for RICH identification
58 #include "TSystem.h"
59 // #include "CbmRichElectronIdAnn.h"
60 
61 #include "CbmL1PFFitter.h"
62 
63 #include "KFPTrackVector.h"
64 #include "KFParticleTopoReconstructor.h"
65 
66 #include "BlastWave.h"
67 #include "BlastWaveLongitudinal.h"
68 #include "RapidityFitBlastWave.h"
69 
70 
71 using std::ios;
72 using std::vector;
73 
74 using namespace std;
75 
77 
78  CbmBlastWave::CbmBlastWave(Int_t recoLevel,
79  Int_t /*iVerbose*/,
80  TString Mode,
81  Int_t PDG,
82  TString pname,
83  Int_t EventStats,
84  KFParticleTopoReconstructor* tr,
85  Float_t ekin_,
86  Double_t Tlong)
87  : CbmModelBase(tr)
88  , name(pname)
89  ,
90  //ekin(ekin_),
91  //fusePID(usePID),
92  ekin(ekin_)
93  , p0cm(5.)
94  , ycm(2.)
95  , fUpdate(true)
96  , fusePID(true)
97  , fRecoLevel(recoLevel)
98  , fTrackNumber(1)
99  , fEventStats(EventStats)
100  , events(0)
101  , fModeName(Mode)
102  , outfileName("")
103  , histodir(0)
104  , flistMCTracks(0)
105  , IndexT(0)
106  , IndexMt(0)
107  , IndexModelMt(0)
108  , IndexMt2(0)
109  , IndexModelMt2(0)
110  , IndexModelMt4Pi(0)
111  , histodndy(0)
112  , histodndymodel(0)
113  , histoeta(0)
114  , histodndymodel2(0)
115  , histomult(0)
116  , histomultmodel(0)
117  , histo1DIntervals(0)
118  , grTy(0)
119  , grdndyReco(0)
120  , pullT(0)
121  , Ts()
122  , kProtonMass(0.938271998)
123  ,
124  //fTrackNumber(trackNumber),
125  //flistStsTracks(0),
126  //flistStsTracksMatch(0),
127  //fPrimVtx(0),
128  //flsitGlobalTracks(0),
129  //flistTofHits(0),
130  fPDGID(PDG)
131  , fMass(TDatabasePDG::Instance()->GetParticle(fPDGID)->Mass())
132  , fYminv()
133  , fYmaxv()
134  , paramGlobal(0)
135  , paramGlobalInterval(0)
136  , param2GlobalInterval(0)
137  , paramLocal(0)
138  , paramLocalInterval(0)
139  , totalLocal(0)
140  , totalGlobal(0)
141  , totalGlobalInterval(0)
142  , totalLocalInterval(0)
143  , y2Local(0.)
144  , y2Global(0.)
145  , y4Local(0.)
146  , y4Global(0.)
147  , totalEvents(0)
148  , model(0)
149  , modelmc(0)
150  , modelL(0)
151  , modelmcL(0)
152  , modelsY(0)
153 // flistRichRings(0),
154 // flistTrdTracks(0),
155 
156 {
157  // fModeName = Mode;
158  // fEventStats = EventStats;
159  // fMass = TDatabasePDG::Instance()->GetParticle(fPDGID)->Mass();
160 
161  // events = 0;
162  Ts.resize(0);
163 
164  //PPDG = 2212;
165  // kProtonMass = 0.938271998;
166 
167  double pbeam = sqrt((kProtonMass + ekin) * (kProtonMass + ekin)
169  double betacm = pbeam / (2. * kProtonMass + ekin);
170  ycm = 0.5 * log((1. + betacm) / (1. - betacm));
171 
172  //PDGtoIndex.clear();
173 
174  TDirectory* currentDir = gDirectory;
175 
176  gDirectory->cd("Models");
177 
178  histodir = gDirectory;
179 
180  char ccc[200];
181  sprintf(ccc, "BlastWave %s", name.Data());
182  gDirectory->mkdir(ccc);
183  gDirectory->cd(ccc);
184  gDirectory->mkdir(fModeName);
185  gDirectory->cd(fModeName);
186  TString tname = "PerEvent";
187  if (fEventStats != 1)
188  tname =
189  TString("Each ") + TString::Itoa(fEventStats, 10) + TString(" events");
190  gDirectory->mkdir(tname);
191  gDirectory->cd(tname);
192  int CurrentIndex = 0;
193 
194  IndexT = CurrentIndex;
195  histo1D[CurrentIndex] = new TH1F("T", "Event-by-event T", 100, 0., 0.4);
196  histo1D[CurrentIndex]->SetXTitle("T (GeV)");
197  histo1D[CurrentIndex]->SetYTitle("Entries");
198  CurrentIndex++;
199 
200  IndexMt = CurrentIndex;
201  histo1D[CurrentIndex] = new TH1F("f(m_{T})", "mt distribution", 200, 0., 2.5);
202  histo1D[CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
203  histo1D[CurrentIndex]->SetYTitle("Entries");
204  CurrentIndex++;
205 
206  IndexModelMt = CurrentIndex;
207  histo1D[CurrentIndex] =
208  new TH1F("Model f(m_{T})", "Model mt distribution", 200, 0., 2.5);
209  histo1D[CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
210  histo1D[CurrentIndex]->SetYTitle("Entries");
211  CurrentIndex++;
212 
213  IndexMt2 = CurrentIndex;
214  histo1D[CurrentIndex] =
215  new TH1F("f2(m_{T})", "mt2 distribution", 40, 0., 2.5);
216  histo1D[CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
217  histo1D[CurrentIndex]->SetYTitle(
218  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
219  CurrentIndex++;
220 
221  IndexModelMt2 = CurrentIndex;
222  histo1D[CurrentIndex] =
223  new TH1F("Model f2(m_{T})", "Model mt2 distribution", 200, 0., 2.5);
224  histo1D[CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
225  histo1D[CurrentIndex]->SetYTitle(
226  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
227  CurrentIndex++;
228 
229  IndexModelMt4Pi = CurrentIndex;
230  histo1D[CurrentIndex] =
231  new TH1F("Model f2(m_{T}) 4Pi", "Model mt2 distribution", 200, 0., 2.5);
232  histo1D[CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
233  histo1D[CurrentIndex]->SetYTitle(
234  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
235  CurrentIndex++;
236 
237  histodndy =
238  new TH1F("dN/dy", "Rapidity distribution", 40, -2. * ycm, 2. * ycm);
239  histodndy->SetXTitle("y");
240  histodndy->SetYTitle("dN/dy");
241 
242  histodndymodel = new TH1F(
243  "dN/dy model", "Rapidity distribution from model", 40, -2. * ycm, 2. * ycm);
244  histodndymodel->SetXTitle("y");
245  histodndymodel->SetYTitle("dN/dy");
246 
247  histoeta = new TH1F("etamax", "Event-by-event etamax", 50, 0., ycm);
248  histoeta->SetXTitle("#eta_{max}");
249  histoeta->SetYTitle("Entries");
250 
251  histodndymodel2 = new TH1F("dN/dy model 2",
252  "Rapidity distribution from model",
253  40,
254  -2. * ycm,
255  2. * ycm);
256  histodndymodel2->SetXTitle("y");
257  histodndymodel2->SetYTitle("dN/dy");
258 
259  histomult =
260  new TH1F("Multiplicity", "Event-by-event multiplicity", 1000, 0.5, 1000.5);
261  histomult->SetXTitle("Event");
262  histomult->SetYTitle("N");
263 
264  histomultmodel = new TH1F("Multiplicity from model",
265  "Event-by-event multiplicity",
266  1000,
267  0.5,
268  1000.5);
269  histomultmodel->SetXTitle("Event");
270  histomultmodel->SetYTitle("N");
271 
272  grdndyReco = new TGraphErrors();
273  grdndyReco->SetTitle("Rapidity distribution");
274  grdndyReco->GetXaxis()->SetTitle("y");
275  grdndyReco->GetYaxis()->SetTitle("dN/dy");
276  grdndyReco->SetName(TString("dNdy-") + fModeName);
277  gDirectory->Add(grdndyReco);
278 
279 
280  if (fYminv.size() > 0) {
281  histo1DIntervals = new TH1F**[fYminv.size()];
282  for (unsigned int ind = 0; ind < fYminv.size(); ++ind) {
283  histo1DIntervals[ind] = new TH1F*[nHisto1D];
284  char cc[200], cc2[200];
285  CurrentIndex = 0;
286 
287  sprintf(cc, "Event-by-event T, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
288  sprintf(cc2, "T, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
289  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 100, 0., 1.);
290  histo1DIntervals[ind][CurrentIndex]->SetXTitle("T (GeV)");
291  histo1DIntervals[ind][CurrentIndex]->SetYTitle("Entries");
292  CurrentIndex++;
293 
294  sprintf(cc, "mt distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
295  sprintf(cc2, "f(m_{T}, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
296  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
297  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
298  histo1DIntervals[ind][CurrentIndex]->SetYTitle("Entries");
299  CurrentIndex++;
300 
301  sprintf(
302  cc, "Model mt distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
303  sprintf(cc2, "Model f(m_{T}), %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
304  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
305  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
306  histo1DIntervals[ind][CurrentIndex]->SetYTitle("Entries");
307  CurrentIndex++;
308 
309  sprintf(cc, "mt2 distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
310  sprintf(cc2, "f2(m_{T}, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
311  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 40, 0., 2.5);
312  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
313  histo1DIntervals[ind][CurrentIndex]->SetYTitle(
314  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
315  CurrentIndex++;
316 
317  sprintf(
318  cc, "Model mt2 distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
319  sprintf(cc2, "Model f2(m_{T}), %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
320  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
321  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
322  histo1DIntervals[ind][CurrentIndex]->SetYTitle(
323  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
324  CurrentIndex++;
325 
326  sprintf(cc,
327  "Model mt2 distribution 4Pi, %.2lf<y<%.2lf",
328  fYminv[ind],
329  fYmaxv[ind]);
330  sprintf(
331  cc2, "Model f2(m_{T}) 4Pi, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
332  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
333  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
334  histo1DIntervals[ind][CurrentIndex]->SetYTitle(
335  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
336  CurrentIndex++;
337  }
338  }
339 
340 
341  gDirectory->cd("..");
342  gDirectory->cd("..");
343  gDirectory->cd("..");
344 
345  gDirectory = currentDir;
346 
347 
348  double ssqrt = sqrt(2. * kProtonMass * (ekin + 2. * kProtonMass));
349  p0cm = sqrt(0.25 * ssqrt * ssqrt - kProtonMass * kProtonMass);
350 
351  std::cout << "ekin = " << ekin << "\n";
352  std::cout << "ycm = " << ycm << "\n";
353 
354  events = 0;
355 
356  paramGlobal = paramLocal = 0.;
357  y2Global = y2Local = 0.;
358  y4Global = y4Local = 0.;
359  totalGlobal = totalLocal = 0;
360  paramGlobalInterval.resize(0);
361  param2GlobalInterval.resize(0);
362  paramLocalInterval.resize(0);
363  totalGlobalInterval.resize(0);
364  totalLocalInterval.resize(0);
365  modelsY.resize(0);
366 
367 
368  if (fRecoLevel == -1 || fRecoLevel > 10)
369  model = new BlastWave(fMass, fPDGID, false, -2. * ycm, 2. * ycm, ycm, 0.);
370  else
371  model = new BlastWave(fMass, fPDGID, true, -2. * ycm, 2. * ycm, ycm, 0.);
372  modelmc = new BlastWave(fMass, fPDGID, true, -2. * ycm, 2. * ycm, ycm, 0.);
373 
374  if (fRecoLevel == -1 || fRecoLevel > 10)
375  modelL = new BlastWaveLongitudinal(
376  fMass, fPDGID, false, -2. * ycm, 2. * ycm, ycm, Tlong);
377  else
378  modelL = new BlastWaveLongitudinal(
379  fMass, fPDGID, true, -2. * ycm, 2. * ycm, ycm, Tlong);
380  modelmcL = new BlastWaveLongitudinal(
381  fMass, fPDGID, false, -2. * ycm, 2. * ycm, ycm, Tlong);
382 
383  //cout << "Level=" << recoLevel << "\n";
384 }
385 
387  if (model != NULL) delete model;
388  for (unsigned int ind = 0; ind < fYminv.size(); ++ind)
389  if (modelsY[ind] != NULL) delete modelsY[ind];
390 }
391 
392 void CbmBlastWave::AddRapidityInterval(double ymin, double ymax) {
393  fYminv.push_back(ymin);
394  fYmaxv.push_back(ymax);
395  paramGlobalInterval.push_back(0.);
396  param2GlobalInterval.push_back(0.);
397  paramLocalInterval.push_back(0.);
398  totalGlobalInterval.push_back(0);
399  totalLocalInterval.push_back(0);
400  if (fRecoLevel == -1 || fRecoLevel > 10)
401  modelsY.push_back(new BlastWave(fMass, fPDGID, false, ymin, ymax, ycm, 0.));
402  else
403  modelsY.push_back(new BlastWave(fMass, fPDGID, true, ymin, ymax, ycm, 0.));
404 }
405 
407  TDirectory* currentDir = gDirectory;
408 
409  gDirectory->cd("Models");
410 
411  histodir = gDirectory;
412 
413  char ccc[200];
414  sprintf(ccc, "BlastWave %s", name.Data());
415  //gDirectory->mkdir(ccc);
416  gDirectory->cd(ccc);
417  //gDirectory->mkdir(fModeName);
418  gDirectory->cd(fModeName);
419  TString tname = "PerEvent";
420  if (fEventStats != 1)
421  tname =
422  TString("Each ") + TString::Itoa(fEventStats, 10) + TString(" events");
423  //gDirectory->mkdir(tname);
424  gDirectory->cd(tname);
425  int CurrentIndex = 0;
426 
427  if (fYminv.size() > 0) {
428  histo1DIntervals = new TH1F**[fYminv.size()];
429  for (unsigned int ind = 0; ind < fYminv.size(); ++ind) {
430  char cc3[200];
431  sprintf(cc3, "%.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
432  gDirectory->mkdir(cc3);
433  gDirectory->cd(cc3);
434 
435  histo1DIntervals[ind] = new TH1F*[nHisto1D];
436  char cc[200], cc2[200];
437  CurrentIndex = 0;
438 
439  sprintf(cc, "Event-by-event T, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
440  sprintf(cc2, "T, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
441  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 100, 0., 1.);
442  histo1DIntervals[ind][CurrentIndex]->SetXTitle("T (GeV)");
443  histo1DIntervals[ind][CurrentIndex]->SetYTitle("Entries");
444  CurrentIndex++;
445 
446  sprintf(cc, "mt distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
447  sprintf(cc2, "f(m_{T}, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
448  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
449  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
450  histo1DIntervals[ind][CurrentIndex]->SetYTitle("Entries");
451  CurrentIndex++;
452 
453  sprintf(
454  cc, "Model mt distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
455  sprintf(cc2, "Model f(m_{T}), %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
456  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
457  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
458  histo1DIntervals[ind][CurrentIndex]->SetYTitle("Entries");
459  CurrentIndex++;
460 
461  sprintf(cc, "mt2 distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
462  sprintf(cc2, "f2(m_{T}, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
463  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 40, 0., 2.5);
464  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
465  histo1DIntervals[ind][CurrentIndex]->SetYTitle(
466  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
467  CurrentIndex++;
468 
469  sprintf(
470  cc, "Model mt2 distribution, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
471  sprintf(cc2, "Model f2(m_{T}), %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
472  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
473  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
474  histo1DIntervals[ind][CurrentIndex]->SetYTitle(
475  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
476  CurrentIndex++;
477 
478  sprintf(cc,
479  "Model mt2 distribution 4Pi, %.2lf<y<%.2lf",
480  fYminv[ind],
481  fYmaxv[ind]);
482  sprintf(
483  cc2, "Model f2(m_{T}) 4Pi, %.2lf<y<%.2lf", fYminv[ind], fYmaxv[ind]);
484  histo1DIntervals[ind][CurrentIndex] = new TH1F(cc2, cc, 200, 0., 2.5);
485  histo1DIntervals[ind][CurrentIndex]->SetXTitle("m_{T} - m_{0} (GeV)");
486  histo1DIntervals[ind][CurrentIndex]->SetYTitle(
487  "#frac{1}{N} #frac{dN}{m_{T} dm_{T}} [(GeV)^{-2}]");
488  CurrentIndex++;
489 
490  gDirectory->cd("..");
491  }
492  grTy = new TGraphErrors();
493  grTy->SetTitle("T(y)");
494  grTy->GetXaxis()->SetTitle("y");
495  grTy->GetYaxis()->SetTitle("T(y) (MeV)");
496  grTy->SetName(TString("T(y)-") + fModeName);
497  gDirectory->Add(grTy);
498  }
499 
500 
501  gDirectory->cd("..");
502  gDirectory->cd("..");
503  gDirectory->cd("..");
504 
505  gDirectory = currentDir;
506 }
507 
508 void CbmBlastWave::ReInit(FairRootManager* fManger) {
509  flistMCTracks = dynamic_cast<TClonesArray*>(fManger->GetObject("MCTrack"));
510 }
511 
513  //ReInit();
514 }
515 
517  if (fRecoLevel == -1 && !flistMCTracks) return;
518  if (fRecoLevel != -1 && !fTopoReconstructor) return;
519 
520  if (fRecoLevel != -1 && fTopoReconstructor != NULL
521  && fTopoReconstructor->GetTracks()[2].PDG().size() == 0)
522  return;
523 
525 
526  events++;
527  totalEvents++;
528 
529  if (events % fEventStats == 0) {
530 
531  //std::cout << paramLocal/ totalLocal << "\n";
532  double T = model->GetT(paramLocal / totalLocal);
533 
534  histo1D[IndexT]->Fill(T);
535  //histo1D[IndexSigz]->Fill(sigz);
536  //histo1D[IndexQz]->Fill(p0cm-qz);
537  double etamax = modelL->Geteta(y2Local / totalLocal);
538  histoeta->Fill(etamax);
539 
540  //std::cout << "T = " << T << "\n";
541  histomult->SetBinContent(totalEvents,
542  totalLocal / static_cast<double>(fEventStats));
543  histomult->SetBinError(
544  totalEvents, 0. * sqrt(totalLocal) / static_cast<double>(fEventStats));
545  double s1 = y2Local;
546  double s2 = y4Local;
547  int Numb = totalLocal;
548  double erretamax =
550  * TMath::Sqrt((s2 - s1 * s1 / Numb) / (Numb - 1.) / Numb);
551  histomultmodel->SetBinContent(
552  totalEvents,
553  modelL->GetA(totalLocal / static_cast<double>(fEventStats), etamax));
554  histomultmodel->SetBinError(
555  totalEvents,
556  modelL->GetAerror(totalLocal / static_cast<double>(fEventStats),
557  etamax,
558  0. * TMath::Sqrt(totalLocal)
559  / static_cast<double>(fEventStats),
560  erretamax));
561 
562  events = 0;
563  paramLocal = 0.;
564  totalLocal = 0;
565  y2Local = y4Local = 0.;
566 
567  for (unsigned int i = 0; i < fYminv.size(); ++i) {
568  if (totalLocalInterval[i] > 0)
569  histo1DIntervals[i][IndexT]->Fill(
571  //std::cout << events << " " << paramLocalInterval[i] / totalLocalInterval[i] << "\n";
572  paramLocalInterval[i] = 0.;
573  totalLocalInterval[i] = 0;
574  }
575  }
576 }
577 
579 
580  double T = model->GetT(paramGlobal / totalGlobal);
581  std::cout << "T = " << T << "\t<m_T> = " << paramGlobal / totalGlobal << "\n";
582 
583  histo1D[IndexMt]->Sumw2();
584  histo1D[IndexMt2]->Sumw2();
585  histo1D[IndexMt]->Scale(1. / totalGlobal
586  / histo1D[IndexMt]->GetXaxis()->GetBinWidth(1));
587  histo1D[IndexMt2]->Scale(1. / totalGlobal
588  / histo1D[IndexMt2]->GetXaxis()->GetBinWidth(1));
589  for (int n = 1; n < histo1D[IndexMt2]->GetNbinsX(); n++) {
590  histo1D[IndexMt2]->SetBinContent(
591  n,
592  histo1D[IndexMt2]->GetBinContent(n)
593  / (histo1D[IndexMt2]->GetXaxis()->GetBinCenter(n) + fMass));
594  histo1D[IndexMt2]->SetBinError(
595  n,
596  histo1D[IndexMt2]->GetBinError(n)
597  / (histo1D[IndexMt2]->GetXaxis()->GetBinCenter(n) + fMass));
598  }
599  histodndy->Sumw2();
600  histodndy->Scale(1. / totalEvents / histodndy->GetXaxis()->GetBinWidth(1));
601 
602  for (int n = 1; n < histo1D[IndexModelMt]->GetNbinsX(); n++) {
603  histo1D[IndexModelMt]->SetBinContent(
604  n,
605  model->fmt(histo1D[IndexModelMt]->GetXaxis()->GetBinCenter(n) + fMass,
606  T));
607  histo1D[IndexModelMt2]->SetBinContent(
608  n,
609  model->fmt(histo1D[IndexModelMt2]->GetXaxis()->GetBinCenter(n) + fMass, T)
610  / (histo1D[IndexModelMt2]->GetXaxis()->GetBinCenter(n) + fMass));
611  histo1D[IndexModelMt4Pi]->SetBinContent(
612  n,
613  modelmc->fmt(
614  histo1D[IndexModelMt4Pi]->GetXaxis()->GetBinCenter(n) + fMass, T)
615  / (histo1D[IndexModelMt4Pi]->GetXaxis()->GetBinCenter(n) + fMass));
616  }
617 
618  std::vector<double> ys(0), dndys(0), dndyerrs(0);
619  int grindex = 0;
620  for (unsigned int ind = 0; ind < fYminv.size(); ++ind) {
621  if (totalGlobalInterval[ind] > 0)
622  T =
623  modelsY[ind]->GetT(paramGlobalInterval[ind] / totalGlobalInterval[ind]);
624  else
625  continue;
626  histo1DIntervals[ind][IndexMt]->Sumw2();
627  histo1DIntervals[ind][IndexMt2]->Sumw2();
628  histo1DIntervals[ind][IndexMt]->Scale(
629  1. / totalGlobalInterval[ind]
630  / histo1DIntervals[ind][IndexMt]->GetXaxis()->GetBinWidth(1));
631  histo1DIntervals[ind][IndexMt2]->Scale(
632  1. / totalGlobalInterval[ind]
633  / histo1DIntervals[ind][IndexMt2]->GetXaxis()->GetBinWidth(1));
634  for (int n = 1; n < histo1DIntervals[ind][IndexMt2]->GetNbinsX(); n++) {
635  histo1DIntervals[ind][IndexMt2]->SetBinContent(
636  n,
637  histo1DIntervals[ind][IndexMt2]->GetBinContent(n)
638  / (histo1DIntervals[ind][IndexMt2]->GetXaxis()->GetBinCenter(n)
639  + fMass));
640  histo1DIntervals[ind][IndexMt2]->SetBinError(
641  n,
642  histo1DIntervals[ind][IndexMt2]->GetBinError(n)
643  / (histo1DIntervals[ind][IndexMt2]->GetXaxis()->GetBinCenter(n)
644  + fMass));
645  }
646  double avmt = 0.;
647  if (T == T && totalGlobalInterval[ind] > 0) {
648  for (int n = 1; n < histo1DIntervals[ind][IndexModelMt]->GetNbinsX();
649  n++) {
650  histo1DIntervals[ind][IndexModelMt]->SetBinContent(
651  n,
652  modelsY[ind]->fmt(
653  histo1DIntervals[ind][IndexModelMt]->GetXaxis()->GetBinCenter(n)
654  + fMass,
655  T));
656  histo1DIntervals[ind][IndexModelMt2]->SetBinContent(
657  n,
658  modelsY[ind]->fmt(
659  histo1DIntervals[ind][IndexModelMt2]->GetXaxis()->GetBinCenter(n)
660  + fMass,
661  T)
662  / (histo1DIntervals[ind][IndexModelMt2]->GetXaxis()->GetBinCenter(n)
663  + fMass));
664  histo1DIntervals[ind][IndexModelMt4Pi]->SetBinContent(
665  n,
666  modelmc->fmt(
667  histo1DIntervals[ind][IndexModelMt4Pi]->GetXaxis()->GetBinCenter(n)
668  + fMass,
669  T)
670  / (histo1DIntervals[ind][IndexModelMt4Pi]->GetXaxis()->GetBinCenter(
671  n)
672  + fMass));
673  avmt +=
674  (histo1DIntervals[ind][IndexModelMt]->GetXaxis()->GetBinCenter(n)
675  + fMass)
676  * modelsY[ind]->fmt(
677  histo1DIntervals[ind][IndexModelMt]->GetXaxis()->GetBinCenter(n)
678  + fMass,
679  T);
680  //std::cout << histo1DIntervals[ind][IndexModelMt]->GetXaxis()->GetBinCenter(n) << " " << modelsY[ind]->fmt(histo1DIntervals[ind][IndexModelMt]->GetXaxis()->GetBinCenter(n), T) <<
681  //" " << modelmc->fmt(histo1DIntervals[ind][IndexModelMt]->GetXaxis()->GetBinCenter(n), T) << "\n";
682  }
683  }
684  avmt *= histo1DIntervals[ind][IndexModelMt]->GetXaxis()->GetBinWidth(1);
685  std::cout << fYminv[ind] << "<y<" << fYmaxv[ind] << "\tT = " << T
686  << "\t<m_T> = "
688  << "\t<m_T>_2 = " << avmt << "\n";
689  double s1 = paramGlobalInterval[ind];
690  double s2 = param2GlobalInterval[ind];
691  int Numb = totalGlobalInterval[ind];
692  double errT =
693  modelsY[ind]->dTdmt(paramGlobalInterval[ind] / totalGlobalInterval[ind])
694  * TMath::Sqrt((s2 - s1 * s1 / Numb) / (Numb - 1.) / Numb);
695  //std::cout << T << " " << errT << " " << (s2 - s1*s1 / Numb) << "\n";
696  if (T != T || errT != errT || errT > T / 3.) continue;
697  grTy->SetPoint(grindex, 0.5 * (fYminv[ind] + fYmaxv[ind]), T * 1.e3);
698  grTy->SetPointError(
699  grindex, 0. * 0.5 * (fYmaxv[ind] - fYminv[ind]), errT * 1.e3);
700  //std::cout << T << " " << errT << "\n";
701 
702  double A = modelsY[ind]->GetA(
703  totalGlobalInterval[ind] / static_cast<double>(totalEvents), T);
704  double errA = modelsY[ind]->GetAerror(
705  totalGlobalInterval[ind] / static_cast<double>(totalEvents),
706  T,
707  sqrt(totalGlobalInterval[ind]) / static_cast<double>(totalEvents),
708  errT);
709  std::cout << "A = " << A << " error = " << errA << "\n";
710  grdndyReco->SetPoint(grindex,
711  0.5 * (fYminv[ind] + fYmaxv[ind]),
712  A / (fYmaxv[ind] - fYminv[ind]));
713  grdndyReco->SetPointError(grindex,
714  0. * 0.5 * (fYmaxv[ind] - fYminv[ind]),
715  errA / (fYmaxv[ind] - fYminv[ind]));
716 
717  ys.push_back(0.5 * (fYminv[ind] + fYmaxv[ind]));
718  dndys.push_back(A / (fYmaxv[ind] - fYminv[ind]));
719  dndyerrs.push_back(errA / (fYmaxv[ind] - fYminv[ind]));
720 
721  grindex++;
722  }
723  if (fYminv.size() > 0) {
724  grTy->GetXaxis()->SetLimits(-3., 3.);
725  grTy->GetYaxis()->SetLimits(0., 300.);
726  //grTy->SetMinimum(0.);
727  //grTy->SetMaximum(300.);
728  grTy->GetXaxis()->SetTitle("y");
729  grTy->GetYaxis()->SetTitle("T (MeV)");
730  grdndyReco->GetXaxis()->SetLimits(-3., 3.);
731  //grdndyReco->GetYaxis()->SetLimits(0., 300.);
732  //grdndyReco->SetMinimum(0.);
733  //grdndyReco->SetMaximum(300.);
734  grdndyReco->GetXaxis()->SetTitle("y");
735  grdndyReco->GetYaxis()->SetTitle("dN/dy");
736  }
737 
739  RapidityFitBlastWave RFit(ys, dndys, dndyerrs, T, fMass);
741 
742  for (int n = 1; n < histodndymodel->GetNbinsX(); n++) {
743  histodndymodel->SetBinContent(n,
745  histodndymodel->GetXaxis()->GetBinCenter(n),
746  params.A.value,
747  params.eta.value,
748  T,
749  fMass));
750  }
751 
752  std::cout << "eta = " << params.eta.value << " "
753  << "Error = " << params.eta.error << "\n";
754 
755  double etamax = modelL->Geteta(y2Global / totalGlobal);
756  double s1 = y2Global;
757  double s2 = y4Global;
758  int Numb = totalGlobal;
759  double erretamax = modelL->detady2(y2Global / totalGlobal)
760  * TMath::Sqrt((s2 - s1 * s1 / Numb) / (Numb - 1.) / Numb);
761  //std::cout << "<y2> = " << y2Global / totalGlobal << " <y2>2 = " << modelL->y2Av(etamax) << " <y2>3 = " << modelmcL->y2Av2(etamax) << "\n";
762  std::cout << "etamax = " << etamax << " "
763  << "Error = " << erretamax << "\n";
764 
765  //double y2av = 0., y2norm = 0.;
766  for (int n = 1; n < histodndymodel2->GetNbinsX(); n++) {
767  histodndymodel2->SetBinContent(
768  n,
769  modelL->GetA(totalGlobal / static_cast<double>(totalEvents), etamax)
770  * modelmcL->fy(histodndymodel2->GetXaxis()->GetBinCenter(n), etamax));
771  }
772 
773  std::cout << "<" << name << "> = "
774  << modelL->GetA(totalGlobal / static_cast<double>(totalEvents),
775  etamax)
776  << " error = "
777  << modelL->GetAerror(totalGlobal / static_cast<double>(totalEvents),
778  etamax,
779  TMath::Sqrt(totalGlobal)
780  / static_cast<double>(totalEvents),
781  erretamax)
782  << "\n";
783 
784  //std::cout << "<y2>4 = " << y2av / y2norm << "\n";
785 }
786 
787 void CbmBlastWave::CalculateAveragesInEvent(int RecoLevel, bool UpdateGlobal) {
788  if (RecoLevel == -1) {
789  vector<CbmMCTrack> vRTracksMC;
790  int nTracksMC = flistMCTracks->GetEntries();
791  std::cout << "MC tracks: " << nTracksMC << "\n";
792  vRTracksMC.resize(nTracksMC);
793  for (int iTr = 0; iTr < nTracksMC; iTr++)
794  // vRTracksMC[iTr] = *( (CbmMCTrack*) flistMCTracks->At(iTr));
795  vRTracksMC[iTr] = *(static_cast<CbmMCTrack*>(flistMCTracks->At(iTr)));
796 
797  for (int iTr = 0; iTr < nTracksMC; iTr++) {
798  if (vRTracksMC[iTr].GetPdgCode() == fPDGID
799  && vRTracksMC[iTr].GetMotherId() == -1) {
800  totalLocal++;
801  double pt = vRTracksMC[iTr].GetPt();
802  double ty = vRTracksMC[iTr].GetRapidity() - ycm;
803  double mt = TMath::Sqrt(
804  vRTracksMC[iTr].GetMass() * vRTracksMC[iTr].GetMass() + pt * pt);
805  //double pz = mt * TMath::SinH(ty);
806  paramLocal += mt;
807  y2Local += ty * ty;
808  y4Local += ty * ty * ty * ty;
809  if (UpdateGlobal) {
810  totalGlobal++;
811  paramGlobal += mt;
812  y2Global += ty * ty;
813  y4Global += ty * ty * ty * ty;
814  }
815  histo1D[IndexMt]->Fill(mt - fMass);
816  histo1D[IndexMt2]->Fill(mt - fMass);
817  histodndy->Fill(ty);
818  for (unsigned int ind = 0; ind < fYminv.size(); ++ind) {
819  if (ty >= fYminv[ind] && ty <= fYmaxv[ind]) {
820  if (UpdateGlobal) totalGlobalInterval[ind]++;
821  totalLocalInterval[ind]++;
822  if (UpdateGlobal) {
823  paramGlobalInterval[ind] += mt;
824  param2GlobalInterval[ind] += mt * mt;
825  }
826  paramLocalInterval[ind] += mt;
827  histo1DIntervals[ind][IndexMt]->Fill(mt - fMass);
828  histo1DIntervals[ind][IndexMt2]->Fill(mt - fMass);
829  }
830  }
831  }
832  }
833  } else {
834  for (int itype = 2; itype <= 3; ++itype) {
835  const KFPTrackVector& tr = fTopoReconstructor->GetTracks()[itype];
836  const kfvector_int& pdgs = tr.PDG();
837  for (unsigned int ind = 0; ind < pdgs.size(); ++ind) {
838  int iPDG = pdgs[ind];
839  //for(unsigned int part=0;part<ParticlePDGsTrack.size();++part) {
840  if (iPDG == fPDGID) {
841  totalLocal++;
842  double pt = tr.Pt(ind);
843  double p = tr.P(ind);
844  double m = TDatabasePDG::Instance()->GetParticle(iPDG)->Mass();
845  double p0 = TMath::Sqrt(m * m + p * p);
846  double pz = TMath::Sqrt(p * p - pt * pt);
847  double ty = 0.5 * log((p0 + pz) / (p0 - pz)) - ycm;
848  pz = TMath::Sqrt(m * m + pt * pt) * TMath::SinH(ty);
849  double mt = TMath::Sqrt(m * m + pt * pt);
850 
851  paramLocal += mt;
852  y2Local += ty * ty;
853  y4Local += ty * ty * ty * ty;
854  if (UpdateGlobal) {
855  totalGlobal++;
856  paramGlobal += mt;
857  y2Global += ty * ty;
858  y4Global += ty * ty * ty * ty;
859  }
860  histo1D[IndexMt]->Fill(mt - fMass);
861  histo1D[IndexMt2]->Fill(mt - fMass);
862  histodndy->Fill(ty);
863  for (unsigned int ind2 = 0; ind2 < fYminv.size(); ++ind2) {
864  if (ty >= fYminv[ind2] && ty <= fYmaxv[ind2]) {
865  if (UpdateGlobal) totalGlobalInterval[ind2]++;
866  totalLocalInterval[ind2]++;
867  if (UpdateGlobal) {
868  paramGlobalInterval[ind2] += mt;
869  param2GlobalInterval[ind2] += mt * mt;
870  }
871  paramLocalInterval[ind2] += mt;
872  histo1DIntervals[ind2][IndexMt]->Fill(mt - fMass);
873  histo1DIntervals[ind2][IndexMt2]->Fill(mt - fMass);
874  }
875  }
876  }
877  //}
878  }
879  }
880  }
881 }
CbmBlastWave::fRecoLevel
Int_t fRecoLevel
Definition: CbmBlastWave.h:77
CbmBlastWave::modelsY
std::vector< BlastWave * > modelsY
Definition: CbmBlastWave.h:129
CbmBlastWave::model
BlastWave * model
Definition: CbmBlastWave.h:127
CbmBlastWave::paramGlobalInterval
std::vector< double > paramGlobalInterval
Definition: CbmBlastWave.h:114
CbmBlastWave::IndexModelMt2
int IndexModelMt2
Definition: CbmBlastWave.h:93
CbmModelBase::fTopoReconstructor
KFParticleTopoReconstructor * fTopoReconstructor
Definition: CbmModelBase.h:42
CbmBlastWave::y2Global
double y2Global
Definition: CbmBlastWave.h:120
RapidityFitBlastWaveParameters::A
RapidityFitBlastWaveParameter A
Definition: RapidityFitBlastWave.h:29
CbmVertex.h
CbmBlastWave::y2Local
double y2Local
Definition: CbmBlastWave.h:120
CbmBlastWave::IndexMt2
int IndexMt2
Definition: CbmBlastWave.h:93
CbmBlastWave
Definition: CbmBlastWave.h:39
RapidityFitBlastWave.h
CbmBlastWave::name
TString name
Definition: CbmBlastWave.h:42
BlastWaveLongitudinal::Geteta
double Geteta(double ay2)
Definition: BlastWaveLongitudinal.h:41
CbmBlastWave::Init
virtual void Init()
Definition: CbmBlastWave.cxx:512
CbmBlastWave::histodndymodel
TH1F * histodndymodel
Definition: CbmBlastWave.h:96
L1Field.h
CbmBlastWave::param2GlobalInterval
std::vector< double > param2GlobalInterval
Definition: CbmBlastWave.h:114
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
BlastWave.h
CbmBlastWave::events
Int_t events
Definition: CbmBlastWave.h:80
ClassImp
ClassImp(CbmBlastWave) CbmBlastWave
Definition: CbmBlastWave.cxx:76
CbmBlastWave::histo1D
TH1F * histo1D[nHisto1D]
Definition: CbmBlastWave.h:95
BlastWaveLongitudinal::detady2
double detady2(double ay2)
Definition: BlastWaveLongitudinal.h:76
CbmBlastWave::histo1DIntervals
TH1F *** histo1DIntervals
Definition: CbmBlastWave.h:99
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
BlastWaveLongitudinal::GetAerror
double GetAerror(double multiplicity, double eta, double multerr, double etaerr)
Definition: BlastWaveLongitudinal.h:57
CbmBlastWave::histodir
TDirectory * histodir
Definition: CbmBlastWave.h:86
CbmBlastWave::paramGlobal
double paramGlobal
Definition: CbmBlastWave.h:113
CbmBlastWave::totalGlobal
int totalGlobal
Definition: CbmBlastWave.h:117
CbmBlastWave::flistMCTracks
TClonesArray * flistMCTracks
Definition: CbmBlastWave.h:88
CbmBlastWave::histoeta
TH1F * histoeta
Definition: CbmBlastWave.h:97
CbmGlobalTrack.h
CbmBlastWave::fYmaxv
std::vector< double > fYmaxv
Definition: CbmBlastWave.h:112
BlastWave::fmt
double fmt(double amt, double T)
Definition: BlastWave.cxx:106
BlastWaveLongitudinal::GetA
double GetA(double multiplicity, double eta)
Definition: BlastWaveLongitudinal.h:49
CbmBlastWave::IndexModelMt4Pi
int IndexModelMt4Pi
Definition: CbmBlastWave.h:93
CbmBlastWave::fMass
double fMass
Definition: CbmBlastWave.h:111
CbmRichRing.h
CbmBlastWave::CalculateAveragesInEvent
void CalculateAveragesInEvent(int RecoLevel, bool UpdateGlobal=0)
Definition: CbmBlastWave.cxx:787
CbmBlastWave::histomult
TH1F * histomult
Definition: CbmBlastWave.h:98
CbmBlastWave::y4Global
double y4Global
Definition: CbmBlastWave.h:120
RapidityFitBlastWaveParameters
Definition: RapidityFitBlastWave.h:28
CbmBlastWave::totalEvents
int totalEvents
Definition: CbmBlastWave.h:122
CbmBlastWave::totalLocalInterval
std::vector< int > totalLocalInterval
Definition: CbmBlastWave.h:118
CbmL1Def.h
CbmBlastWave.h
RapidityFitBlastWaveParameter::value
double value
Definition: RapidityFitBlastWave.h:10
BlastWaveLongitudinal
Definition: BlastWaveLongitudinal.h:12
CbmStsTrack.h
Data class for STS tracks.
CbmBlastWave::Exec
virtual void Exec()
Definition: CbmBlastWave.cxx:516
log
friend F32vec4 log(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:135
RapidityFitBlastWaveNamespace::dndy
double dndy(double y, double A, double eta, double T, double mass)
Definition: RapidityFitBlastWave.cxx:20
CbmTrackMatchNew.h
RapidityFitBlastWaveParameter::error
double error
Definition: RapidityFitBlastWave.h:11
CbmL1PFFitter.h
RapidityFitBlastWave
Definition: RapidityFitBlastWave.h:68
CbmBlastWave::fYminv
std::vector< double > fYminv
Definition: CbmBlastWave.h:112
RapidityFitBlastWaveParameters::eta
RapidityFitBlastWaveParameter eta
Definition: RapidityFitBlastWave.h:29
CbmKFTrack.h
CbmBlastWave::modelmcL
BlastWaveLongitudinal * modelmcL
Definition: CbmBlastWave.h:128
CbmBlastWave::CbmBlastWave
CbmBlastWave(Int_t recoLevel=-1, Int_t iVerbose=1, TString Mode="MC", Int_t PDG=-211, TString pname="pi-", Int_t EventStats=1, KFParticleTopoReconstructor *tr=0, Float_t ekin_=25., Double_t Tlong=0.132)
CbmBlastWave::totalLocal
int totalLocal
Definition: CbmBlastWave.h:117
CbmBlastWave::fModeName
TString fModeName
Definition: CbmBlastWave.h:83
CbmModelBase
Definition: CbmModelBase.h:25
CbmBlastWave::ycm
Float_t ycm
Definition: CbmBlastWave.h:74
CbmBlastWave::histodndymodel2
TH1F * histodndymodel2
Definition: CbmBlastWave.h:97
CbmBlastWave::paramLocal
double paramLocal
Definition: CbmBlastWave.h:115
CbmBlastWave::paramLocalInterval
std::vector< double > paramLocalInterval
Definition: CbmBlastWave.h:116
CbmBlastWave::totalGlobalInterval
std::vector< int > totalGlobalInterval
Definition: CbmBlastWave.h:118
CbmBlastWave::AddHistos
void AddHistos()
Definition: CbmBlastWave.cxx:406
CbmBlastWave::fPDGID
int fPDGID
Definition: CbmBlastWave.h:110
CbmBlastWave::~CbmBlastWave
~CbmBlastWave()
Definition: CbmBlastWave.cxx:386
CbmBlastWave::Finish
virtual void Finish()
Definition: CbmBlastWave.cxx:578
CbmBlastWave::fEventStats
Int_t fEventStats
Definition: CbmBlastWave.h:79
CbmBlastWave::ReInit
virtual void ReInit(FairRootManager *fManger)
Definition: CbmBlastWave.cxx:508
CbmMCTrack.h
CbmBlastWave::grTy
TGraphErrors * grTy
Definition: CbmBlastWave.h:101
CbmMCTrack
Definition: CbmMCTrack.h:34
m
__m128 m
Definition: L1/vectors/P4_F32vec4.h:26
CbmBlastWave::modelmc
BlastWave * modelmc
Definition: CbmBlastWave.h:127
CbmBlastWave::nHisto1D
static const int nHisto1D
Definition: CbmBlastWave.h:92
BlastWaveLongitudinal.h
RapidityFitBlastWave::PerformFit
RapidityFitBlastWaveParameters PerformFit()
Definition: RapidityFitBlastWave.cxx:76
CbmBlastWave::histodndy
TH1F * histodndy
Definition: CbmBlastWave.h:96
CbmBlastWave::IndexT
int IndexT
Definition: CbmBlastWave.h:93
CbmBlastWave::IndexModelMt
int IndexModelMt
Definition: CbmBlastWave.h:93
CbmBlastWave::IndexMt
int IndexMt
Definition: CbmBlastWave.h:93
CbmBlastWave::AddRapidityInterval
void AddRapidityInterval(double ymin, double ymax)
Definition: CbmBlastWave.cxx:392
CbmBlastWave::grdndyReco
TGraphErrors * grdndyReco
Definition: CbmBlastWave.h:102
CbmTrdTrack.h
kProtonMass
const Double_t kProtonMass
Definition: CbmPhsdGenerator.cxx:31
CbmBlastWave::histomultmodel
TH1F * histomultmodel
Definition: CbmBlastWave.h:98
BlastWave
Definition: BlastWave.h:12
CbmKFVertex.h
CbmBlastWave::y4Local
double y4Local
Definition: CbmBlastWave.h:120
BlastWave::GetT
double GetT(double amt)
Definition: BlastWave.h:38
BlastWaveLongitudinal::fy
double fy(double y, double eta)
Definition: BlastWaveLongitudinal.cxx:122
CbmBlastWave::modelL
BlastWaveLongitudinal * modelL
Definition: CbmBlastWave.h:128