Garfield++ v1r0
A toolkit for the detailed simulation of particle detectors based on ionisation measurement in gases and semiconductors
Loading...
Searching...
No Matches
Public Member Functions | List of all members
Parabol Class Reference

#include <parabol.h>

Public Member Functions

double a (void) const
 
double b (void) const
 
double c (void) const
 
void put_a (double fa)
 
void put_b (double fb)
 
void put_c (double fc)
 
 Parabol (void)
 
 Parabol (const Parabol &f)
 
 Parabol (double fa, double fb, double fc)
 
 Parabol (double x[3], double y[3])
 
 Parabol (double x[3], double y[3], int)
 
 Parabol (double x1, double x2, double x3, double y1, double y2, double y3)
 
double y (double x)
 
int find_zero (double xzero[2]) const
 
double find_maxmin (void)
 
double determinant (void) const
 

Detailed Description

Definition at line 32 of file parabol.h.

Constructor & Destructor Documentation

◆ Parabol() [1/6]

Parabol::Parabol ( void  )
inline

Definition at line 54 of file parabol.h.

54: da(0.0), db(0.0), dc(0.0), s_det(0), s_dxzero(0) {}

◆ Parabol() [2/6]

Parabol::Parabol ( const Parabol &  f)

Definition at line 19 of file parabol.cpp.

19{ *this = f; }

◆ Parabol() [3/6]

Parabol::Parabol ( double  fa,
double  fb,
double  fc 
)
inline

Definition at line 56 of file parabol.h.

57 : da(fa), db(fb), dc(fc), s_det(0), s_dxzero(0) {}

◆ Parabol() [4/6]

Parabol::Parabol ( double  x[3],
double  y[3] 
)

Definition at line 21 of file parabol.cpp.

21 : s_det(0), s_dxzero(0) {
22 mfunname("Parabol::Parabol(double x[3], double y[3])");
23
24 check_econd12a(x[0], ==, x[1],
25 "x[2]=" << x[2] << " y[0]=" << y[0] << " y[1]=" << y[1]
26 << " y[2]=" << y[2] << '\n',
27 mcerr);
28 check_econd12a(x[0], ==, x[2],
29 "x[1]=" << x[1] << " y[0]=" << y[0] << " y[1]=" << y[1]
30 << " y[2]=" << y[2] << '\n',
31 mcerr);
32 check_econd12a(x[1], ==, x[2],
33 "x[0]=" << x[0] << " y[0]=" << y[0] << " y[1]=" << y[1]
34 << " y[2]=" << y[2] << '\n',
35 mcerr);
36 DynArr<DoubleAc> mat(3, 3);
37 DynLinArr<DoubleAc> par(3);
38 DynLinArr<DoubleAc> f(3);
39 for (int i = 0; i < 3; ++i) {
40 f[i] = y[i];
41 mat.ac(i, 2) = 1.0;
42 mat.ac(i, 1) = x[i];
43 mat.ac(i, 0) = x[i] * x[i];
44 }
45 int ierr;
46 int szero;
47 DynArr<DoubleAc> mat_inv;
48 inverse_DynArr_prot(mat, mat_inv, szero, ierr);
49 //check_econd11a( ierr, != 0 , "should never happen\n", mcerr );
50 if (ierr == 0) {
51 par = mat_inv * f;
52 //Iprintdla_DoubleAc(mcout, par, 3);
53 //if(fabs(par[0]) == 0.0)
54 // da=0.0;
55 //else
56 // da=par[0];
57 da = par[0];
58 db = par[1];
59 dc = par[2];
60 } else {
61 da = 0.0;
62 DynLinArr<int> s_var(3);
63 s_var[0] = 0;
64 s_var[1] = 1;
65 s_var[2] = 1;
66 DynArr<DoubleAc> mat_inv1(3, 3);
67 //int ierr1;
68 //inverse_DynArr(mat, s_var, mat_inv, ierr, mat_inv1, ierr1);
69 inverse_DynArr_prot(mat, s_var, mat_inv, szero, ierr);
70 if (ierr != 0) {
71 // what if x1 and x2 are the same but the both differ from x0
72 // Then the following should help:
73 mat.ac(1, 1) = mat.ac(0, 1);
74 mat.ac(1, 2) = mat.ac(0, 2);
75 f[1] = f[0];
76 s_var[0] = 0;
77 s_var[1] = 1;
78 s_var[2] = 1;
79 inverse_DynArr_prot(mat, s_var, mat_inv, szero, ierr);
80 check_econd11a(ierr, != 0,
81 "should never happen\nmat=" << mat << "\ns_var=" << s_var
82 << "\nmat_inv=" << mat_inv,
83 mcerr);
84 }
85 par = mat_inv * f;
86 db = par[1];
87 dc = par[2];
88 }
89
90}
check_econd11a
#define check_econd11a(a, signb, add, stream)
Definition: FunNameStack.h:395
check_econd12a
#define check_econd12a(a, sign, b, add, stream)
Definition: FunNameStack.h:411
mfunname
#define mfunname(string)
Definition: FunNameStack.h:67
Parabol::y
double y(double x)
Definition: parabol.h:66
inverse_DynArr_prot
void inverse_DynArr_prot(const DynArr< DoubleAc > &mi, DynArr< DoubleAc > &mr, int &szero, int &serr, int s_stop)
Definition: inverse.cpp:15
mcerr
#define mcerr
Definition: prstream.h:135

◆ Parabol() [5/6]

Parabol::Parabol ( double  x[3],
double  y[3],
int   
)

Definition at line 169 of file parabol.cpp.

169 : s_det(0), s_dxzero(0) {
170 mfunname("Parabol::Parabol(double x[3], double y[3], int)");
171
172 check_econd12(x[0], ==, x[1], mcerr);
173 //check_econd12( x[0] , == , x[2] , mcerr);
174 //check_econd12( x[1] , == , x[2] , mcerr);
175
176 DynArr<DoubleAc> mat(3, 3);
177 DynLinArr<DoubleAc> par(3);
178 DynLinArr<DoubleAc> f(3);
179 for (int i = 0; i < 3; ++i)
180 f[i] = y[i];
181 for (int i = 0; i < 2; ++i) {
182 mat.ac(i, 2) = 1.0;
183 mat.ac(i, 1) = x[i];
184 //Iprintdan(mcout, mat.ac(i,1));
185 mat.ac(i, 0) = x[i] * x[i];
186 }
187 mat.ac(2, 2) = 0.0;
188 mat.ac(2, 1) = 1.0;
189 mat.ac(2, 0) = 2.0 * x[2];
190 int ierr;
191 int szero;
192 DynArr<DoubleAc> mat_inv;
193 inverse_DynArr_prot(mat, mat_inv, szero, ierr);
194 check_econd11a(ierr, != 0, "should never happen\n", mcerr);
195 /*
196 if (ierr != 0) {
197 da = 0.0;
198 DynLinArr<int> s_var(3);
199 s_var[0] = 0;
200 s_var[1] = 1;
201 s_var[2] = 1;
202 inverse_DynArr(mat, mat_inv, ierr);
203 check_econd11a( ierr, != 0 , "should never happen\n", mcerr );
204 par = mat_inv * f;
205 db=par[1]; dc=par[2];
206 } else {
207 */
208 par = mat_inv * f;
209 if (fabs(par[0]) == 0.0) {
210 da = 0.0;
211 } else {
212 da = par[0];
213 db = par[1];
214 dc = par[2];
215 }
216 /*
217 HepMatrix mat(3,3,0);
218 HepVector par(3,0);
219 HepVector f(3,0);
220 for (int i = 0; i < 3; i++) f[i]=y[i];
221 for (int i = 0; i < 2; i++) {
222 mat[i][2] = 1.0;
223 mat[i][1] = x[i];
224 mat[i][0] = x[i] * x[i];
225 }
226 mat[2][2] = 0.0;
227 mat[2][1] = 1.0;
228 mat[2][0] = 2.0 * x[2];
229 par = solve(mat, f);
230 da = par[0]; db = par[1]; dc = par[2];
231 */
232}
fabs
DoubleAc fabs(const DoubleAc &f)
Definition: DoubleAc.h:616
check_econd12
#define check_econd12(a, sign, b, stream)
Definition: FunNameStack.h:380

◆ Parabol() [6/6]

Parabol::Parabol ( double  x1,
double  x2,
double  x3,
double  y1,
double  y2,
double  y3 
)

Definition at line 92 of file parabol.cpp.

94 : s_det(0), s_dxzero(0) {
95 mfunname("Parabol::Parabol(double x[3], double y[3])");
96
97 check_econd12a(x1, ==, x2, "x3=" << x3 << " y1=" << y1 << " y2=" << y2
98 << " y3=" << y3 << '\n',
99 mcerr);
100 check_econd12a(x1, ==, x3, "x2=" << x2 << " y1=" << y1 << " y2=" << y2
101 << " y3=" << y3 << '\n',
102 mcerr);
103 check_econd12a(x2, ==, x3, "x1=" << x1 << " y1=" << y1 << " y2=" << y2
104 << " y3=" << y3 << '\n',
105 mcerr);
106 DynArr<DoubleAc> mat(3, 3);
107 DynLinArr<DoubleAc> par(3);
108 DynLinArr<DoubleAc> f(3);
109 f[0] = y1;
110 mat.ac(0, 2) = 1.0;
111 mat.ac(0, 1) = x1;
112 mat.ac(0, 0) = x1 * x1;
113 f[1] = y2;
114 mat.ac(1, 2) = 1.0;
115 mat.ac(1, 1) = x2;
116 mat.ac(1, 0) = x2 * x2;
117 f[2] = y3;
118 mat.ac(2, 2) = 1.0;
119 mat.ac(2, 1) = x3;
120 mat.ac(2, 0) = x3 * x3;
121
122 int ierr;
123 int szero;
124 DynArr<DoubleAc> mat_inv;
125 inverse_DynArr_prot(mat, mat_inv, szero, ierr);
126 //check_econd11a( ierr, != 0 , "should never happen\n", mcerr );
127 if (ierr == 0) {
128 par = mat_inv * f;
129 //Iprintdla_DoubleAc(mcout, par, 3);
130 //if(fabs(par[0]) == 0.0)
131 // da=0.0;
132 //else
133 // da=par[0];
134 da = par[0];
135 db = par[1];
136 dc = par[2];
137 } else {
138 da = 0.0;
139 DynLinArr<int> s_var(3);
140 s_var[0] = 0;
141 s_var[1] = 1;
142 s_var[2] = 1;
143 DynArr<DoubleAc> mat_inv1(3, 3);
144 //int ierr1;
145 //inverse_DynArr(mat, s_var, mat_inv, ierr, mat_inv1, ierr1);
146 inverse_DynArr_prot(mat, s_var, mat_inv, szero, ierr);
147 if (ierr != 0) {
148 // what if x1 and x2 are the same but the both differ from x0
149 // Then the following should help:
150 mat.ac(1, 1) = mat.ac(0, 1);
151 mat.ac(1, 2) = mat.ac(0, 2);
152 f[1] = f[0];
153 s_var[0] = 0;
154 s_var[1] = 1;
155 s_var[2] = 1;
156 inverse_DynArr_prot(mat, s_var, mat_inv, szero, ierr);
157 check_econd11a(ierr, != 0,
158 "should never happen\nmat=" << mat << "\ns_var=" << s_var
159 << "\nmat_inv=" << mat_inv,
160 mcerr);
161 }
162 par = mat_inv * f;
163 db = par[1];
164 dc = par[2];
165 }
166
167}

Member Function Documentation

◆ a()

double Parabol::a ( void  ) const
inline

Definition at line 35 of file parabol.h.

35{ return da; }

Referenced by operator<<().

◆ b()

double Parabol::b ( void  ) const
inline

Definition at line 36 of file parabol.h.

36{ return db; }

Referenced by operator<<().

◆ c()

double Parabol::c ( void  ) const
inline

Definition at line 37 of file parabol.h.

37{ return dc; }

Referenced by operator<<().

◆ determinant()

double Parabol::determinant ( void  ) const
inline

Definition at line 72 of file parabol.h.

72 {
73 convmut(Parabol);
74 if (s_det == 0) {
75 t.s_det = 1;
76 t.det = db * db - 4 * da * dc;
77 }
78 return det;
79 }
convmut
#define convmut(type)
Definition: emul_new_stand.h:47

Referenced by find_zero().

◆ find_maxmin()

double Parabol::find_maxmin ( void  )

Definition at line 279 of file parabol.cpp.

279 {
280 mfunname("double Parabol::find_maxmin(void)");
281 check_econd11(da, == 0, mcerr);
282 return -db / (2.0 * da);
283}
check_econd11
#define check_econd11(a, signb, stream)
Definition: FunNameStack.h:366

◆ find_zero()

int Parabol::find_zero ( double  xzero[2]) const

Definition at line 234 of file parabol.cpp.

234 {
235 mfunnamep("int Parabol::find_zero(double xzero[2]) const");
236 convmut(Parabol);
237 if (s_dxzero == 0) {
238 //mcout<<"Parabol::find_zero: s_dxzero == 0\n";
239 t.s_dxzero = 1;
240 if (da == 0.0) {
241 if (db == 0.0) {
242 funnw.ehdr(mcerr);
243 mcerr << "can not find zero\n";
244 spexit(mcerr);
245 } else {
246 t.qdxzero = 1;
247 t.dxzero[0] = -dc / db;
248 }
249 } else {
250 if (determinant() < 0.0) {
251 t.qdxzero = 0;
252 t.dxzero[0] = 0;
253 t.dxzero[1] = 0;
254 } else if (determinant() == 0.0) {
255 t.qdxzero = 1;
256 t.dxzero[0] = -db / (2.0 * da);
257 } else {
258 const double sq = sqrt(determinant());
259 t.qdxzero = 2;
260 if (da > 0.0) {
261 t.dxzero[0] = (-db - sq) / (2.0 * da);
262 t.dxzero[1] = (-db + sq) / (2.0 * da);
263 } else {
264 t.dxzero[1] = (-db - sq) / (2.0 * da);
265 t.dxzero[0] = (-db + sq) / (2.0 * da);
266 }
267 //mcout<<"Parabol::find_zero: t.dxzero[0]="<<t.dxzero[0]
268 // <<" dxzero[0]="<<dxzero[0]<<'\n';
269 //mcout<<"Parabol::find_zero: t.dxzero[1]="<<t.dxzero[1]
270 // <<" dxzero[1]="<<dxzero[1]<<'\n';
271 }
272 }
273 }
274 xzero[0] = dxzero[0];
275 xzero[1] = dxzero[1];
276 return qdxzero;
277}
sqrt
DoubleAc sqrt(const DoubleAc &f)
Definition: DoubleAc.cpp:313
mfunnamep
#define mfunnamep(string)
Definition: FunNameStack.h:77
spexit
#define spexit(stream)
Definition: FunNameStack.h:536
Parabol::determinant
double determinant(void) const
Definition: parabol.h:72

Referenced by Cubic::find_maxmin(), and operator<<().

◆ put_a()

void Parabol::put_a ( double  fa)
inline

Definition at line 38 of file parabol.h.

38 {
39 da = fa;
40 s_det = 0;
41 s_dxzero = 0;
42 }

◆ put_b()

void Parabol::put_b ( double  fb)
inline

Definition at line 43 of file parabol.h.

43 {
44 db = fb;
45 s_det = 0;
46 s_dxzero = 0;
47 }

◆ put_c()

void Parabol::put_c ( double  fc)
inline

Definition at line 48 of file parabol.h.

48 {
49 dc = fc;
50 s_det = 0;
51 s_dxzero = 0;
52 }

◆ y()

double Parabol::y ( double  x)
inline

Definition at line 66 of file parabol.h.

66{ return da * x * x + db * x + dc; }

Referenced by Parabol().

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