Belle II Software development
Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
EDepInGas Class Reference

The Class for Energy deposit in the gas. More...

#include <EDepInGas.h>

Public Member Functions

virtual ~EDepInGas ()
 Destructor.
 
void initialize ()
 Initialize theclass.
 
double getEDepInGas (int mode, int pdg, double p, double dx, double e3) const
 Return the energy deosite in the gas.
 
double getMostProbabEDep (double p, double mass, double zi, double dx, double z, double a, double i, double rho) const
 Return the energy deosite in the material.
 

Static Public Member Functions

static EDepInGas & getInstance ()
 Static method to get a reference to the EDepInGas instance.
 

Private Member Functions

 EDepInGas ()
 Singleton class.
 
 EDepInGas (const EDepInGas &)
 Singleton class.
 
EDepInGas & operator= (const EDepInGas &)
 Singleton class.
 

Private Attributes

double m_z1 = 0
 Z of He-C2H6 gas.
 
double m_a1 = 0
 A of He-C2H6 gas.
 
double m_i1 = 0
 I of He-C2H6 gas.
 
double m_rho1 = 0
 rho of He-C2H6 gas
 
double m_z2 = 0
 Z of "wire gas".
 
double m_a2 = 0
 A of "wire gas".
 
double m_i2 = 0
 I of "wire gas".
 
double m_rho2 = 0
 rho of "wire gas"
 
double m_z0 = 0
 Z of wire+gas.
 
double m_a0 = 0
 A of wire+gas.
 
double m_i0 = 0
 I of wire+gas.
 
double m_rho0 = 0
 rho of wire+gas
 
const double m_massE = 0.51099895e6
 electron mass
 
double m_ak1 [6][2] = {{ }}
 coeffs a for density effect
 
double m_bk1 [6][2] = {{ }}
 coeffs b for density effect
 
double m_ck1 [6][2] = {{ }}
 coeffs c for density effect
 

Static Private Attributes

static EDepInGas * m_pntr = 0
 Pointer that saves the instance of this class.
 

Detailed Description

The Class for Energy deposit in the gas.

This class provides the energy deposit in the gas.

Definition at line 20 of file EDepInGas.h.

Constructor & Destructor Documentation

◆ ~EDepInGas()

~EDepInGas ( )
virtual

Destructor.

Definition at line 34 of file EDepInGas.cc.

34{}

◆ EDepInGas()

EDepInGas ( )
private

Singleton class.

Definition at line 32 of file EDepInGas.cc.

32{}

Member Function Documentation

◆ getEDepInGas()

double getEDepInGas ( int  mode,
int  pdg,
double  p,
double  dx,
double  e3 
) const

Return the energy deosite in the gas.

Parameters
mode0: simple scaling; 1: scaling based on mpe; 2: based on probab.
pdgpdg code of incoming particle
pabsolute momentum (GeV)
dxthickness of material (cm)
e3energy deposit in gas+wire

Definition at line 137 of file EDepInGas.cc.

138{
139 double xi1 = m_rho1 * (m_z1 / m_a1);
140 double xi2 = m_rho2 * (m_z2 / m_a2);
141 // B2DEBUG(m_dbg, (xi1+xi2)/xi1 <<" "<< e3);
142 if (mode == 0) return e3 * xi1 / (xi1 + xi2);
143
144 static TDatabasePDG* ptrTDatabasePDG = TDatabasePDG::Instance();
145 const double mass = ptrTDatabasePDG->GetParticle(pdg)->Mass();
146 const double chrg = fabs(ptrTDatabasePDG->GetParticle(pdg)->Charge() / 3);
147
148 const double me1 = getMostProbabEDep(mom, mass, chrg, dx, m_z1, m_a1, m_i1, m_rho1);
149 if (mode == 2) {
150 const double me0 = getMostProbabEDep(mom, mass, chrg, dx, m_z0, m_a0, m_i0, m_rho0);
151 return e3 * me1 / me0;
152 }
153 const double me2 = getMostProbabEDep(mom, mass, chrg, dx, m_z2, m_a2, m_i2, m_rho2);
154 if (mode == 1) return e3 * me1 / (me1 + me2);
155
156 // extract energy deposit (e1) based on probab., L(e1)*L(e3-e1),
157 // where L is the analytic approx. of Landau dist.
158 // find lm1 which maximizes probability
159 const double mom2 = mom * mom;
160 const double beta2 = mom2 / (mass * mass + mom2);
161 const double buf = 0.1535e-3 * dx * chrg * chrg / beta2;
162 xi1 *= buf;
163 xi2 *= buf;
164 const double xi1i = 1 / xi1;
165 const double xi2i = 1 / xi2;
166 const double lm1min = -me1 * xi1i;
167 const double lm1max = (e3 - me1) * xi1i;
168 double lm1 = (e3 * me1 / (me1 + me2) - me1) * xi1i; // initial guess
169 lm1 = min(lm1, lm1max);
170 lm1 = max(lm1, lm1min);
171 // double edm = 1;
172 int itry = 0;
173 // B2DEBUG(m_dbg, "e3,me1,me2,xi1i=" << e3 <<" "<< me1 <<" "<< me2 <<" "<< xi1i);
174 // B2DEBUG(m_dbg, "lm1min,lm1max,lm1=" << lm1min <<" "<< lm1max <<" "<< lm1);
175 const double lm2b = (e3 - me1 - me2) * xi2i;
176 const double r = xi1 * xi2i;
177 double deltal = 0;
178 double lm1old = lm1;
179 while (itry <= 2 || fabs(lm1 - lm1old) > 0.001) {
180 ++itry;
181 lm1old = lm1;
182 if (itry > 50) {
183 // B2DEBUG(m_dbg, "itry > 50" << "itry,lm1min,lm1max,lm1,edm=" << itry <<" "<< lm1min <<" "<< lm1max <<" "<< lm1);
184 break;
185 }
186 lm1 -= deltal;
187 lm1 = max(lm1min, lm1);
188 lm1 = min(lm1max, lm1);
189 double lm2 = lm2b - r * lm1;
190 double emlm1 = exp(-lm1);
191 double emlm2 = exp(-lm2);
192 double dldlm1 = 1 - emlm1 - r + r * emlm2;
193 double dl2dlm12 = emlm1 + r * r * emlm2;
194 deltal = dldlm1 / dl2dlm12;
195 // edm = dldlm1 * deltal; //drop 0.5* for speed up
196 }
197
198 if (mode == 3) return xi1 * lm1 + me1;
199
200 double lm2 = lm2b - r * lm1;
201 double lmin = 0.5 * (lm1 + lm2 + exp(-lm1) + exp(-lm2));
202 const double ymax = exp(-lmin);
203 // if (ymax < 1.e-10) return e3*me1/(me1+me2);
204
205 // generate e1 according to probab.
206 const double lm1w = lm1max - lm1min;
207 double y = 10;
208 double p = 1;
209 int jtry = 0;
210 double urndm[2];
211 while (y > p) {
212 ++jtry;
213 if (jtry > 500) {
214 // B2DEBUG(m_dbg, "jtry > 500" <<"ymax" << ymax);
215 lm1 = (e3 * me1 / (me1 + me2) - me1) * xi1i;
216 break;
217 }
218 gRandom->RndmArray(2, urndm);
219 lm1 = lm1min + lm1w * urndm[0];
220 lm2 = lm2b - r * lm1;
221 double l = 0.5 * (lm1 + lm2 + exp(-lm1) + exp(-lm2));
222 p = exp(-l);
223 y = ymax * urndm[1];
224 }
225
226 return xi1 * lm1 + me1;
227}
Belle2::CDC::EDepInGas::m_i1
double m_i1
I of He-C2H6 gas.
Definition: EDepInGas.h:71
Belle2::CDC::EDepInGas::m_i0
double m_i0
I of wire+gas.
Definition: EDepInGas.h:81
Belle2::CDC::EDepInGas::m_z1
double m_z1
Z of He-C2H6 gas.
Definition: EDepInGas.h:69
Belle2::CDC::EDepInGas::m_i2
double m_i2
I of "wire gas".
Definition: EDepInGas.h:76
Belle2::CDC::EDepInGas::m_rho0
double m_rho0
rho of wire+gas
Definition: EDepInGas.h:82
Belle2::CDC::EDepInGas::m_rho1
double m_rho1
rho of He-C2H6 gas
Definition: EDepInGas.h:72
Belle2::CDC::EDepInGas::m_rho2
double m_rho2
rho of "wire gas"
Definition: EDepInGas.h:77
Belle2::CDC::EDepInGas::m_a1
double m_a1
A of He-C2H6 gas.
Definition: EDepInGas.h:70
Belle2::CDC::EDepInGas::m_z2
double m_z2
Z of "wire gas".
Definition: EDepInGas.h:74
Belle2::CDC::EDepInGas::m_z0
double m_z0
Z of wire+gas.
Definition: EDepInGas.h:79
Belle2::CDC::EDepInGas::m_a0
double m_a0
A of wire+gas.
Definition: EDepInGas.h:80
Belle2::CDC::EDepInGas::getMostProbabEDep
double getMostProbabEDep(double p, double mass, double zi, double dx, double z, double a, double i, double rho) const
Return the energy deosite in the material.
Definition: EDepInGas.cc:229
Belle2::CDC::EDepInGas::m_a2
double m_a2
A of "wire gas".
Definition: EDepInGas.h:75

◆ getInstance()

EDepInGas & getInstance ( )
static

Static method to get a reference to the EDepInGas instance.

Returns
A reference to an instance of this class.

Definition at line 23 of file EDepInGas.cc.

24{
25 if (!m_pntr) {
26 m_pntr = new EDepInGas();
27 m_pntr->initialize();
28 }
29 return *m_pntr;
30}
Belle2::CDC::EDepInGas::EDepInGas
EDepInGas()
Singleton class.
Definition: EDepInGas.cc:32
Belle2::CDC::EDepInGas::initialize
void initialize()
Initialize theclass.
Definition: EDepInGas.cc:36
Belle2::CDC::EDepInGas::m_pntr
static EDepInGas * m_pntr
Pointer that saves the instance of this class.
Definition: EDepInGas.h:89

◆ getMostProbabEDep()

double getMostProbabEDep ( double  p,
double  mass,
double  zi,
double  dx,
double  z,
double  a,
double  i,
double  rho 
) const

Return the energy deosite in the material.

Parameters
pabsolute momentum (GeV)
masslmass of incoming particle (GeV)
ziz of incoming particle
dxthickness of material (cm)
zatomic number of material
aatomic weight of material
iionization potential of material (eV)
rhodensity of material

Definition at line 229 of file EDepInGas.cc.

231{
232// Calculate most probable energy deposit in GeV.
233// Ref. J.Va'vra et al., NIM 203 (1982) 109.
234// ----------------------------------------------------------------------
235// <Input>
236// p : total momentum (GeV)
237// mass : mass of particle (GeV)
238// zi : charge of particle (e)
239// dx : thickness of gas (cm)
240// Z : atomic no. of gas
241// A : atomic weight of gas
242// I : mean ionization potential of gas (eV)
243// rho : density of gas (g/cm**3)
244// ----------------------------------------------------------------------
245 const double psq = p * p;
246 const double masq = mass * mass;
247 const double tote = psq + masq;
248 const double beta2 = psq / tote;
249 const double bgam2 = psq / masq;
250 const double gamma = sqrt(tote) / mass;
251 const double alft = 0.1535e6 * rho * dx * Z / A;
252 const double k = log(m_massE * alft / I / I) + 0.891;
253 const double k1 = k - log(dx);
254 const double alngam = log(gamma);
255 const double alggam = log10(gamma);
256 // calculate density correction term
257 int irgn = 0;
258 double aa, bb, cc;
259 if (k1 < 6.) {
260 if (alggam > 3.6) irgn = 1;
261 aa = (m_ak1[1][irgn] - m_ak1[0][irgn]) * (6. - k1) + m_ak1[0][irgn];
262 bb = (m_bk1[1][irgn] - m_bk1[0][irgn]) * (6. - k1) + m_bk1[0][irgn];
263 cc = (m_ck1[1][irgn] - m_ck1[0][irgn]) * (6. - k1) + m_ck1[0][irgn];
264 } else if (k1 < 7.) {
265 if (alggam > 3.6) irgn = 1;
266 aa = (m_ak1[1][irgn] - m_ak1[0][irgn]) * (k1 - 6.) + m_ak1[0][irgn];
267 bb = (m_bk1[1][irgn] - m_bk1[0][irgn]) * (k1 - 6.) + m_bk1[0][irgn];
268 cc = (m_ck1[1][irgn] - m_ck1[0][irgn]) * (k1 - 6.) + m_ck1[0][irgn];
269 } else if (k1 < 8.) {
270 if (alggam > 3.6) irgn = 1;
271 aa = (m_ak1[2][irgn] - m_ak1[1][irgn]) * (k1 - 7.) + m_ak1[1][irgn];
272 bb = (m_bk1[2][irgn] - m_bk1[1][irgn]) * (k1 - 7.) + m_bk1[1][irgn];
273 cc = (m_ck1[2][irgn] - m_ck1[1][irgn]) * (k1 - 7.) + m_ck1[1][irgn];
274 } else if (k1 < 9.) {
275 if (alggam > 3.6) irgn = 1;
276 aa = (m_ak1[3][irgn] - m_ak1[2][irgn]) * (k1 - 8.) + m_ak1[2][irgn];
277 bb = (m_bk1[3][irgn] - m_bk1[2][irgn]) * (k1 - 8.) + m_bk1[2][irgn];
278 cc = (m_ck1[3][irgn] - m_ck1[2][irgn]) * (k1 - 8.) + m_ck1[2][irgn];
279 } else if (k1 < 10.) {
280 if (alggam > 3.6) irgn = 1;
281 aa = (m_ak1[4][irgn] - m_ak1[3][irgn]) * (k1 - 9.) + m_ak1[3][irgn];
282 bb = (m_bk1[4][irgn] - m_bk1[3][irgn]) * (k1 - 9.) + m_bk1[3][irgn];
283 cc = (m_ck1[4][irgn] - m_ck1[3][irgn]) * (k1 - 9.) + m_ck1[3][irgn];
284 } else if (k1 < 11.) {
285 if (alggam > 4.0) irgn = 1;
286 aa = (m_ak1[5][irgn] - m_ak1[4][irgn]) * (k1 - 10.) + m_ak1[4][irgn];
287 bb = (m_bk1[5][irgn] - m_bk1[4][irgn]) * (k1 - 10.) + m_bk1[4][irgn];
288 cc = (m_ck1[5][irgn] - m_ck1[4][irgn]) * (k1 - 10.) + m_ck1[4][irgn];
289 } else {
290 if (alggam > 4.0) irgn = 1;
291 aa = (m_ak1[5][irgn] - m_ak1[4][irgn]) * (k1 - 11.) + m_ak1[5][irgn];
292 bb = (m_bk1[5][irgn] - m_bk1[4][irgn]) * (k1 - 11.) + m_bk1[5][irgn];
293 cc = (m_ck1[5][irgn] - m_ck1[4][irgn]) * (k1 - 11.) + m_ck1[5][irgn];
294 }
295
296 double delta = aa + alggam * (bb + cc * alggam);
297 delta = max(delta, 0.);
298 delta = min(delta, fabs(log(bgam2)));
299
300 double val = zi * zi * alft * (k + 2 * alngam - beta2 - delta) / beta2;
301 val = max(val, 0.);
302 val *= Unit::eV;
303 return val;
304}
Belle2::CDC::EDepInGas::m_bk1
double m_bk1[6][2]
coeffs b for density effect
Definition: EDepInGas.h:86
Belle2::CDC::EDepInGas::m_ak1
double m_ak1[6][2]
coeffs a for density effect
Definition: EDepInGas.h:85
Belle2::CDC::EDepInGas::m_massE
const double m_massE
electron mass
Definition: EDepInGas.h:84
Belle2::CDC::EDepInGas::m_ck1
double m_ck1[6][2]
coeffs c for density effect
Definition: EDepInGas.h:87
Belle2::Unit::eV
static const double eV
[electronvolt]
Definition: Unit.h:112
Belle2::sqrt
double sqrt(double a)
sqrt for double
Definition: beamHelpers.h:28

◆ initialize()

void initialize ( )

Initialize theclass.

Definition at line 36 of file EDepInGas.cc.

37{
38// gas
39 const double zHe = 2;
40 const double aHe = 4.003;
41 double fHe = 0.1108;
42
43 const double zC = 6;
44 const double aC = 12.011;
45 double fC = 0.2223;
46
47 const double zH = 1;
48 const double aH = 1.008;
49 double fH = 0.6669;
50
51 double f = fHe + fC + fH;
52 fHe /= f;
53 fC /= f;
54 fH /= f;
55 m_z1 = fHe * zHe + fC * zC + fH * zH;
56 m_a1 = fHe * aHe + fC * aC + fH * aH;
57 m_i1 = 49.229;
58 m_rho1 = 0.703e-3;
59 // B2DEBUG(m_dbg, m_z1 <<" "<< m_a1 <<" "<< m_i1 <<" "<< m_rho1*m_z1/m_a1);
60
61 // wire
62 const double zW = 74;
63 const double aW = 183.842;
64 double fW = 0.0198;
65
66 const double zAl = 13;
67 const double aAl = 26.982;
68 double fAl = 0.9802;
69
70 f = fW + fAl;
71 fW /= f;
72 fAl /= f;
73 m_z2 = fW * zW + fAl * zAl;
74 m_a2 = fW * aW + fAl * aAl;
75 m_i2 = 193.301;
76 m_rho2 = 0.419e-3;
77 // B2DEBUG(m_dbg, m_z2 <<" "<< m_a2 <<" "<< m_i2 <<" "<< m_rho2*m_z2/m_a2);
78
79 fHe = 0.01031;
80 fC = 0.02068;
81 fH = 0.06204;
82 fW = 0.0014;
83 fAl = 0.0684;
84 f = fHe + fC + fH + fW + fAl;
85 fHe /= f;
86 fC /= f;
87 fH /= f;
88 fW /= f;
89 fAl /= f;
90 m_z0 = fHe * zHe + fC * zC + fH * zH + fW * zW + fAl * zAl;
91 m_a0 = fHe * aHe + fC * aC + fH * aH + fW * aW + fAl * aAl;
92 m_i0 = 76.693;
93 // m_rho0 = m_rho1 + m_rho2;
94 m_rho0 = 1.121e-3;
95
96 //Define coefficients of density correction term
97 m_ak1[0][0] = 10.065;
98 m_ak1[0][1] = -13.63;
99 m_ak1[1][0] = 4.1334;
100 m_ak1[1][1] = -12.4;
101 m_ak1[2][0] = 0.8135;
102 m_ak1[2][1] = -11.33;
103 m_ak1[3][0] = 0.00402;
104 m_ak1[3][1] = -10.36;
105 m_ak1[4][0] = - 1.7564;
106 m_ak1[4][1] = - 9.611;
107 m_ak1[5][0] = - 1.5065;
108 m_ak1[5][1] = - 8.6753;
109
110 m_bk1[0][0] = - 8.5619;
111 m_bk1[0][1] = 4.6;
112 m_bk1[1][0] = - 4.5857;
113 m_bk1[1][1] = 4.6;
114 m_bk1[2][0] = - 2.1447;
115 m_bk1[2][1] = 4.6;
116 m_bk1[3][0] = - 1.1582;
117 m_bk1[3][1] = 4.6;
118 m_bk1[4][0] = 0.67256;
119 m_bk1[4][1] = 4.6;
120 m_bk1[5][0] = 1.0156;
121 m_bk1[5][1] = 4.6;
122
123 m_ck1[0][0] = 1.82804;
124 m_ck1[0][1] = 0.0;
125 m_ck1[1][0] = 1.27579;
126 m_ck1[1][1] = 0.0;
127 m_ck1[2][0] = 0.93676;
128 m_ck1[2][1] = 0.0;
129 m_ck1[3][0] = 0.79975;
130 m_ck1[3][1] = 0.0;
131 m_ck1[4][0] = 0.49093;
132 m_ck1[4][1] = 0.0;
133 m_ck1[5][0] = 0.44805;
134 m_ck1[5][1] = 0.0;
135}

Member Data Documentation

◆ m_a0

double m_a0 = 0
private

A of wire+gas.

Definition at line 80 of file EDepInGas.h.

◆ m_a1

double m_a1 = 0
private

A of He-C2H6 gas.

Definition at line 70 of file EDepInGas.h.

◆ m_a2

double m_a2 = 0
private

A of "wire gas".

Definition at line 75 of file EDepInGas.h.

◆ m_ak1

double m_ak1[6][2] = {{ }}
private

coeffs a for density effect

Definition at line 85 of file EDepInGas.h.

◆ m_bk1

double m_bk1[6][2] = {{ }}
private

coeffs b for density effect

Definition at line 86 of file EDepInGas.h.

◆ m_ck1

double m_ck1[6][2] = {{ }}
private

coeffs c for density effect

Definition at line 87 of file EDepInGas.h.

◆ m_i0

double m_i0 = 0
private

I of wire+gas.

Definition at line 81 of file EDepInGas.h.

◆ m_i1

double m_i1 = 0
private

I of He-C2H6 gas.

Definition at line 71 of file EDepInGas.h.

◆ m_i2

double m_i2 = 0
private

I of "wire gas".

Definition at line 76 of file EDepInGas.h.

◆ m_massE

const double m_massE = 0.51099895e6
private

electron mass

Definition at line 84 of file EDepInGas.h.

◆ m_pntr

EDepInGas * m_pntr = 0
staticprivate

Pointer that saves the instance of this class.

Definition at line 89 of file EDepInGas.h.

◆ m_rho0

double m_rho0 = 0
private

rho of wire+gas

Definition at line 82 of file EDepInGas.h.

◆ m_rho1

double m_rho1 = 0
private

rho of He-C2H6 gas

Definition at line 72 of file EDepInGas.h.

◆ m_rho2

double m_rho2 = 0
private

rho of "wire gas"

Definition at line 77 of file EDepInGas.h.

◆ m_z0

double m_z0 = 0
private

Z of wire+gas.

Definition at line 79 of file EDepInGas.h.

◆ m_z1

double m_z1 = 0
private

Z of He-C2H6 gas.

Definition at line 69 of file EDepInGas.h.

◆ m_z2

double m_z2 = 0
private

Z of "wire gas".

Definition at line 74 of file EDepInGas.h.

The documentation for this class was generated from the following files: