CbmRoot
CbmLitCheckBrem.cxx
Go to the documentation of this file.
2 
6 
7 #include "TAxis.h"
8 #include "TCanvas.h"
9 #include "TGraph.h"
10 #include "TLegend.h"
11 #include "TMultiGraph.h"
12 #include "TPad.h"
13 #include "TStyle.h"
14 
15 #include <cmath>
16 #include <cstdlib>
17 
19  fNofMom = 10;
20  fMom.resize(fNofMom);
21  fMom[0] = 0.05; // GeV
22  fMom[1] = 0.06; // GeV
23  fMom[2] = 0.08; // GeV
24  fMom[3] = 0.1; // GeV
25  fMom[4] = 0.2; // GeV
26  fMom[5] = 0.5; // GeV
27  fMom[6] = 1.; // GeV
28  fMom[7] = 2.; // GeV
29  fMom[8] = 5.; // GeV
30  fMom[9] = 10.; // GeV
31 
32  fNofMaterials = 3;
33  fTable.resize(fNofMaterials);
34  fCalc.resize(fNofMaterials);
35  fMaterials.resize(fNofMaterials);
36 
37  //N2 gas
39  m0.fA = 14.00674;
40  m0.fZ = 7;
41  m0.fRho = 1.250; //* 1e-3; // [g/cm3]
42  m0.fX0 = 47.1; // [cm]
43  m0.fPHIRAD.resize(fNofMom);
44  m0.fPHIRAD[0] = 17.63;
45  m0.fPHIRAD[1] = 17.99;
46  m0.fPHIRAD[2] = 18.51;
47  m0.fPHIRAD[3] = 18.87;
48  m0.fPHIRAD[4] = 19.81;
49  m0.fPHIRAD[5] = 20.60;
50  m0.fPHIRAD[6] = 20.96;
51  m0.fPHIRAD[7] = 21.18;
52  m0.fPHIRAD[8] = 21.36;
53  m0.fPHIRAD[9] = 21.43;
54  fMaterials[0] = m0;
55 
56  // Cu
58  m1.fA = 63.546;
59  m1.fZ = 29;
60  m1.fRho = 8.96; // [g/cm3]
61  m1.fX0 = 1.43; // [cm]
62  m1.fPHIRAD.resize(fNofMom);
63  m1.fPHIRAD[0] = 14.36;
64  m1.fPHIRAD[1] = 14.60;
65  m1.fPHIRAD[2] = 14.95;
66  m1.fPHIRAD[3] = 15.19;
67  m1.fPHIRAD[4] = 15.81;
68  m1.fPHIRAD[5] = 16.30;
69  m1.fPHIRAD[6] = 16.52;
70  m1.fPHIRAD[7] = 16.66;
71  m1.fPHIRAD[8] = 16.76;
72  m1.fPHIRAD[9] = 16.80;
73  fMaterials[1] = m1;
74 
75  // Sn
77  m2.fA = 118.710;
78  m2.fZ = 50;
79  m2.fRho = 7.31; // [g/cm3]
80  m2.fX0 = 1.21; // [cm]
81  m2.fPHIRAD.resize(fNofMom);
82  m2.fPHIRAD[0] = 13.31;
83  m2.fPHIRAD[1] = 13.53;
84  m2.fPHIRAD[2] = 13.83;
85  m2.fPHIRAD[3] = 14.05;
86  m2.fPHIRAD[4] = 14.58;
87  m2.fPHIRAD[5] = 15.01;
88  m2.fPHIRAD[6] = 15.20;
89  m2.fPHIRAD[7] = 15.32;
90  m2.fPHIRAD[8] = 15.41;
91  m2.fPHIRAD[9] = 15.44;
92  fMaterials[2] = m2;
93 }
94 
96 
98  gStyle->SetCanvasColor(kWhite);
99  gStyle->SetFrameFillColor(kWhite);
100  gStyle->SetPadColor(kWhite);
101  gStyle->SetStatColor(kWhite);
102  gStyle->SetTitleFillColor(kWhite);
103  gStyle->SetPalette(1);
104 
105  CreateGraphs();
106  FillGraphs();
107  DrawGraphs();
108 }
109 
111  TCanvas* c1 = new TCanvas("brem_loss", "c1", 800, 800);
112  c1->SetGrid();
113 
114  for (int i = 0; i < fNofMaterials; i++) {
115  fTable[i]->SetLineStyle(3);
116  fTable[i]->SetLineColor(2);
117  fTable[i]->SetMarkerColor(2);
118  fTable[i]->SetLineWidth(3);
119  fTable[i]->SetMarkerSize(2);
120 
121  fCalc[i]->SetLineStyle(1);
122  fCalc[i]->SetLineColor(4);
123  fCalc[i]->SetMarkerColor(4);
124  fCalc[i]->SetLineWidth(3);
125  fCalc[i]->SetMarkerSize(2);
126  }
127  fTable[0]->SetMarkerStyle(20);
128  fTable[1]->SetMarkerStyle(26);
129  fTable[2]->SetMarkerStyle(27);
130  // fTable[3]->SetMarkerStyle(25);
131 
132  fCalc[0]->SetMarkerStyle(20);
133  fCalc[1]->SetMarkerStyle(26);
134  fCalc[2]->SetMarkerStyle(27);
135  // fCalc[3]->SetMarkerStyle(25);
136 
137  TMultiGraph* mg = new TMultiGraph();
138  for (int i = 0; i < fNofMaterials; i++) {
139  mg->Add(fTable[i]);
140  mg->Add(fCalc[i]);
141  }
142  mg->SetMinimum(0.001);
143  mg->SetMaximum(10.);
144  gPad->SetLogx();
145  gPad->SetLogy();
146 
147  mg->Draw("ALP");
148 
149  mg->GetXaxis()->SetTitle("Electron momentum [GeV/c]");
150  mg->GetYaxis()->SetTitle("Brehmstrahlung energy loss [GeV/cm]");
151  mg->GetXaxis()->SetLimits(0.04, 11);
152 
153  TLegend* l1 = new TLegend(0.20, 0.97, 0.9, 0.7);
154  l1->SetFillColor(kWhite);
155  l1->SetHeader("Energy losses for muons in iron vs. momentum");
156  l1->AddEntry(fTable[0], "N2 (table)", "lp");
157  l1->AddEntry(fCalc[0], "N2 (calculation)", "lp");
158  l1->AddEntry(fTable[1], "Cu (table)", "lp");
159  l1->AddEntry(fCalc[1], "Cu (calculation)", "lp");
160  l1->AddEntry(fTable[2], "Sn (table)", "lp");
161  l1->AddEntry(fCalc[2], "Sn (calculation)", "lp");
162  l1->Draw();
163 
164  c1->SaveAs("brem_loss.gif");
165  c1->SaveAs("brem_loss.eps");
166  c1->SaveAs("brem_loss.svg");
167 }
168 
170  for (Int_t i = 0; i < fNofMom; i++) {
171  Double_t p = fMom[i];
172 
173  for (Int_t j = 0; j < fNofMaterials; j++) {
175  // Double_t brem_calc = p / (m.fX0);
176  // Double_t Z2 = (m.fZ*m.fZ);
177  // Double_t X0 = 1./(m.fA * 716.408) * (Z2 * (std::log(184.15/std::pow(Z, 1./3.)) ) + Z*)
178  Double_t X0 =
179  (716.4 * m.fA) / (m.fZ * (m.fZ + 1) * std::log(287. / std::sqrt(m.fZ)));
180  std::cout << "X0=" << m.fX0 * m.fRho << " X0calc=" << X0 << std::endl;
181  Double_t brem_calc = p / X0;
182  fCalc[j]->SetPoint(i, p, brem_calc);
183 
184  Double_t brem_table = (6.022045e23 / m.fA) * 5.794661e-28 * (m.fZ * m.fZ)
185  * (p + 0.000511) * m.fPHIRAD[i]; // * m.fRho;
186  fTable[j]->SetPoint(i, p, brem_table);
187  }
188  }
189 }
190 
192  for (int i = 0; i < fNofMaterials; i++) {
193  fTable[i] = new TGraph();
194  fCalc[i] = new TGraph();
195  }
196 }
197 
CbmLitTrackParam.h
Data class for track parameters.
CbmLitCheckBrem
Definition: CbmLitCheckBrem.h:19
CbmLitCheckBrem::fNofMaterials
Int_t fNofMaterials
Definition: CbmLitCheckBrem.h:35
CbmLitCheckBrem::fMaterials
std::vector< CbmLitSimpleMaterial > fMaterials
Definition: CbmLitCheckBrem.h:38
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
CbmLitSimpleMaterial::fZ
Double_t fZ
Definition: CbmLitCheckBrem.h:12
CbmLitCheckBrem::fCalc
std::vector< TGraph * > fCalc
Definition: CbmLitCheckBrem.h:37
CbmLitCheckBrem::CbmLitCheckBrem
CbmLitCheckBrem()
Definition: CbmLitCheckBrem.cxx:18
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmLitCheckBrem::CreateGraphs
void CreateGraphs()
Definition: CbmLitCheckBrem.cxx:191
CbmLitSimpleMaterial::fX0
Double_t fX0
Definition: CbmLitCheckBrem.h:15
CbmLitCheckBrem::fTable
std::vector< TGraph * > fTable
Definition: CbmLitCheckBrem.h:36
CbmLitCheckBrem::~CbmLitCheckBrem
virtual ~CbmLitCheckBrem()
Definition: CbmLitCheckBrem.cxx:95
log
friend F32vec4 log(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:135
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmLitSimpleMaterial::fRho
Double_t fRho
Definition: CbmLitCheckBrem.h:14
CbmLitCheckBrem::fNofMom
Int_t fNofMom
Definition: CbmLitCheckBrem.h:33
m
__m128 m
Definition: L1/vectors/P4_F32vec4.h:26
CbmLitSimpleMaterial::fPHIRAD
std::vector< Double_t > fPHIRAD
Definition: CbmLitCheckBrem.h:16
CbmLitCheckBrem::FillGraphs
void FillGraphs()
Definition: CbmLitCheckBrem.cxx:169
CbmLitCheckBrem.h
CbmLitSimpleMaterial::fA
Double_t fA
Definition: CbmLitCheckBrem.h:13
CbmLitMaterialInfo.h
CbmLitCheckBrem::Check
virtual void Check()
Definition: CbmLitCheckBrem.cxx:97
CbmLitSimpleMaterial
Definition: CbmLitCheckBrem.h:9
CbmLitCheckBrem::fMom
std::vector< double > fMom
Definition: CbmLitCheckBrem.h:34
CbmLitCheckBrem::DrawGraphs
void DrawGraphs()
Definition: CbmLitCheckBrem.cxx:110
CbmLitMaterialEffectsImp.h
Calculation of multiple scattering and energy loss.