Garfield++ v1r0
A toolkit for the detailed simulation of particle detectors based on ionisation measurement in gases and semiconductors
Loading...
Searching...
No Matches
ComponentAnalyticField.hh
Go to the documentation of this file.
1#ifndef G_COMPONENT_ANALYTIC_FIELD_H
2#define G_COMPONENT_ANALYTIC_FIELD_H
3
4#include <cmath>
5#include <complex>
6
7#include "ComponentBase.hh"
8#include "FundamentalConstants.hh"
9
10namespace Garfield {
11
12class ComponentAnalyticField : public ComponentBase {
13
14 public:
15 // Constructor
16 ComponentAnalyticField();
17 // Destructor
18 ~ComponentAnalyticField() {}
19
20 void ElectricField(const double x, const double y, const double z, double& ex,
21 double& ey, double& ez, Medium*& m, int& status);
22 void ElectricField(const double x, const double y, const double z, double& ex,
23 double& ey, double& ez, double& v, Medium*& m,
24 int& status);
25
26 bool GetVoltageRange(double& pmin, double& pmax);
27
28 void WeightingField(const double x, const double y, const double z,
29 double& wx, double& wy, double& wz,
30 const std::string label);
31 double WeightingPotential(const double x, const double y, const double z,
32 const std::string label);
33
34 bool GetBoundingBox(double& x0, double& y0, double& z0, double& x1,
35 double& y1, double& z1);
36
37 bool IsWireCrossed(double x0, double y0, double z0, double x1, double y1,
38 double z1, double& xc, double& yc, double& zc);
39
40 bool IsInTrapRadius(double x0, double y0, double z0, double& xw, double& yx,
41 double& rw);
42
43 // Add a wire at (x, y) .
44 void AddWire(const double x, const double y, const double diameter,
45 const double voltage, const std::string label,
46 const double length = 100., const double tension = 50.,
47 const double rho = 19.3, const int ntrap = 5);
48 // Add a tube.
49 void AddTube(const double radius, const double voltage, const int nEdges,
50 const std::string label);
51 // Add a plane at constant x
52 void AddPlaneX(const double x, const double voltage, const std::string label);
53 // Add a plane at constant y
54 void AddPlaneY(const double y, const double voltage, const std::string label);
55
56 void AddStripOnPlaneX(const char direction, const double x, const double smin,
57 const double smax, const std::string label,
58 const double gap = -1.);
59 void AddStripOnPlaneY(const char direction, const double y, const double smin,
60 const double smax, const std::string label,
61 const double gap = -1.);
62
63 // Set the periodic length [cm] in x/y direction
64 void SetPeriodicityX(const double s);
65 void SetPeriodicityY(const double s);
66 bool GetPeriodicityX(double& s);
67 bool GetPeriodicityY(double& s);
68
69 void AddCharge(const double x, const double y, const double z,
70 const double q);
71 void ClearCharges();
72 void PrintCharges();
73
74 std::string GetCellType() { return cellType; }
75 // Cells are classified according to the number
76 // and orientation of planes, the presence of
77 // periodicities and the location of the wires
78 // as one of the following types:
79 //
80 // A non-periodic cells with at most 1 x- and 1 y-plane
81 // B1X x-periodic cells without x-planes and at most 1 y-plane
82 // B1Y y-periodic cells without y-planes and at most 1 x-plane
83 // B2X cells with 2 x-planes and at most 1 y-plane
84 // B2Y cells with 2 y-planes and at most 1 x-plane
85 // C1 doubly periodic cells without planes
86 // C2X doubly periodic cells with x-planes
87 // C2Y doubly periodic cells with y-planes
88 // C3 double periodic cells with x- and y-planes
89 // D1 round tubes without axial periodicity
90 // D2 round tubes with axial periodicity
91 // D3 polygonal tubes without axial periodicity
92
93 void AddReadout(const std::string label);
94
95 void EnableChargeCheck() { chargeCheck = true; }
96 void DisableChargeCheck() { chargeCheck = false; }
97
98 int GetNumberOfWires() { return nWires; }
99 bool GetWire(const int i, double& x, double& y, double& diameter,
100 double& voltage, std::string& label, double& length,
101 double& charge, int& ntrap);
102
103 int GetNumberOfPlanesX();
104 int GetNumberOfPlanesY();
105 bool GetPlaneX(const int i, double& x, double& voltage, std::string& label);
106 bool GetPlaneY(const int i, double& y, double& voltage, std::string& label);
107
108 bool GetTube(double& r, double& voltage, int& nEdges, std::string& label);
109
110 private:
111 bool chargeCheck;
112
113 bool cellset;
114 bool sigset;
115
116 bool polar;
117
118 // Cell type (as string and number)
119 std::string cellType;
120 int iCellType;
121
122 // Bounding box
123 double xmin, xmax;
124 double ymin, ymax;
125 double zmin, zmax;
126
127 // Voltage range
128 double vmin, vmax;
129
130 // Periodicities
131 bool perx, pery;
132 double sx, sy;
133
134 // Signals
135 int nFourier;
136 std::string cellTypeFourier;
137 bool fperx, fpery;
138 int mxmin, mxmax, mymin, mymax;
139 int mfexp;
140
141 int nReadout;
142 std::vector<std::string> readout;
143
144 // Wires
145 int nWires;
146 struct wire {
147 // Location
148 double x;
149 double y;
150 // Diameter
151 double d;
152 // Potential
153 double v;
154 // Charge
155 double e;
156 // Label
157 std::string type;
158 // Length
159 double u;
160 // Readout group
161 int ind;
162 // Trap Radius - tracked particle "trapped" if within
163 // nTrap*radius of wire.
164 int nTrap;
165 };
166 std::vector<wire> w;
167
168 // Stretching weight
169 std::vector<double> weight;
170 // Density
171 std::vector<double> dens;
172 // Mirror charges for force calculations
173 std::vector<double> cnalso;
174
175 // Option for computation of dipole terms
176 bool dipole;
177 // Dipole angle and amplitude
178 std::vector<double> cosph2;
179 std::vector<double> sinph2;
180 std::vector<double> amp2;
181
182 // Parameters for B2 type cells
183 std::vector<double> b2sin;
184 // Parameters for C type cells
185 int mode;
186 std::complex<double> zmult;
187 double p1, p2, c1;
188 // Parameters for D3 type cells
189 // Conformal mapping in polygons
190 std::vector<std::complex<double> > wmap;
191 double kappa;
192 // Tables of coefficients
193 std::vector<std::vector<double> > cc1;
194 std::vector<std::vector<double> > cc2;
195
196 // Reference potential
197 double v0;
198 double corvta, corvtb, corvtc;
199
200 // Planes
201 // Existence
202 bool ynplan[4];
203 bool ynplax, ynplay;
204 // Coordinates
205 double coplan[4];
206 double coplax, coplay;
207 // Voltages
208 double vtplan[4];
209
210 struct strip {
211 // Label
212 std::string type;
213 // Readout group
214 int ind;
215 // Coordinates
216 double smin, smax;
217 double gap;
218 };
219
220 struct plane {
221 // Labels
222 std::string type;
223 // Readout group
224 int ind;
225 // Background weighting fields
226 double ewxcor;
227 double ewycor;
228 // x/y strips
229 int nStrips1;
230 std::vector<strip> strips1;
231 // z strips
232 int nStrips2;
233 std::vector<strip> strips2;
234 };
235
236 std::vector<plane> planes;
237
238 // Tube
239 bool tube;
240 int mtube, ntube;
241 double cotube;
242 double vttube;
243
244 // Capacitance matrix
245 std::vector<std::vector<double> > a;
246 // Signal matrix
247 std::vector<std::vector<std::complex<double> > > sigmat;
248 // Induced charges on planes
249 std::vector<std::vector<double> > qplane;
250
251 // Point charges
252 int n3d;
253 struct charge3d {
254 double x, y, z;
255 double e;
256 };
257 std::vector<charge3d> ch3d;
258 int nTermBessel;
259 int nTermPoly;
260
261 // Gravity
262 double down[3];
263
264 void UpdatePeriodicity();
265 void Reset() { CellInit(); }
266
267 void CellInit();
268 bool Prepare();
269 bool CellCheck();
270 bool CellType();
271 bool PrepareStrips();
272
273 bool PrepareSignals();
274 bool SetupWireSignals();
275 bool SetupPlaneSignals();
276
277 // Calculation of charges
278 bool Setup();
279 bool SetupA00();
280 bool SetupB1X();
281 bool SetupB1Y();
282 bool SetupB2X();
283 bool SetupB2Y();
284 bool SetupC10();
285 bool SetupC2X();
286 bool SetupC2Y();
287 bool SetupC30();
288 bool SetupD10();
289 bool SetupD20();
290 bool SetupD30();
291
292 bool IprA00(const int mx, const int my);
293 bool IprB2X(const int my);
294 bool IprB2Y(const int mx);
295 bool IprC2X();
296 bool IprC2Y();
297 bool IprC30();
298 bool IprD10();
299 bool IprD30();
300
301 bool SetupDipole() { return true; }
302
303 // Inversion of capacitance matrix
304 bool Charge();
305
306 // Evaluation of the electric field
307 int Field(const double xin, const double yin, const double zin, double& ex,
308 double& ey, double& ez, double& volt, const bool opt);
309 void FieldA00(const double xpos, const double ypos, double& ex, double& ey,
310 double& volt, const bool opt);
311 void FieldB1X(const double xpos, const double ypos, double& ex, double& ey,
312 double& volt, const bool opt);
313 void FieldB1Y(const double xpos, const double ypos, double& ex, double& ey,
314 double& volt, const bool opt);
315 void FieldB2X(const double xpos, const double ypos, double& ex, double& ey,
316 double& volt, const bool opt);
317 void FieldB2Y(const double xpos, const double ypos, double& ex, double& ey,
318 double& volt, const bool opt);
319 void FieldC10(const double xpos, const double ypos, double& ex, double& ey,
320 double& volt, const bool opt);
321 void FieldC2X(const double xpos, const double ypos, double& ex, double& ey,
322 double& volt, const bool opt);
323 void FieldC2Y(const double xpos, const double ypos, double& ex, double& ey,
324 double& volt, const bool opt);
325 void FieldC30(const double xpos, const double ypos, double& ex, double& ey,
326 double& volt, const bool opt);
327 void FieldD10(const double xpos, const double ypos, double& ex, double& ey,
328 double& volt, const bool opt);
329 void FieldD20(const double xpos, const double ypos, double& ex, double& ey,
330 double& volt, const bool opt);
331 void FieldD30(const double xpos, const double ypos, double& ex, double& ey,
332 double& volt, const bool opt);
333
334 // Field due to point charges
335 void Field3dA00(const double x, const double y, const double z, double& ex,
336 double& ey, double& ez, double& volt);
337 void Field3dB2X(const double x, const double y, const double z, double& ex,
338 double& ey, double& ez, double& volt);
339 void Field3dB2Y(const double x, const double y, const double z, double& ex,
340 double& ey, double& ez, double& volt);
341 void Field3dD10(const double x, const double y, const double z, double& ex,
342 double& ey, double& ez, double& volt);
343 // Evaluation of the weighting field
344 bool Wfield(const double xpos, const double ypos, const double zpos,
345 double& ex, double& ey, double& ez, double& volt, const int isw,
346 const bool opt);
347 void WfieldWireA00(const double xpos, const double ypos, double& ex,
348 double& ey, double& volt, const int mx, const int my,
349 const int sw, const bool opt);
350 void WfieldWireB2X(const double xpos, const double ypos, double& ex,
351 double& ey, double& volt, const int my, const int sw,
352 const bool opt);
353 void WfieldWireB2Y(const double xpos, const double ypos, double& ex,
354 double& ey, double& volt, const int mx, const int sw,
355 const bool opt);
356 void WfieldWireC2X(const double xpos, const double ypos, double& ex,
357 double& ey, double& volt, const int sw, const bool opt);
358 void WfieldWireC2Y(const double xpos, const double ypos, double& ex,
359 double& ey, double& volt, const int sw, const bool opt);
360 void WfieldWireC30(const double xpos, const double ypos, double& ex,
361 double& ey, double& volt, const int sw, const bool opt);
362 void WfieldWireD10(const double xpos, const double ypos, double& ex,
363 double& ey, double& volt, const int sw, const bool opt);
364 void WfieldWireD30(const double xpos, const double ypos, double& ex,
365 double& ey, double& volt, const int sw, const bool opt);
366 void WfieldPlaneA00(const double xpos, const double ypos, double& ex,
367 double& ey, double& volt, const int mx, const int my,
368 const int iplane, const bool opt);
369 void WfieldPlaneB2X(const double xpos, const double ypos, double& ex,
370 double& ey, double& volt, const int my, const int iplane,
371 const bool opt);
372 void WfieldPlaneB2Y(const double xpos, const double ypos, double& ex,
373 double& ey, double& volt, const int mx, const int iplane,
374 const bool opt);
375 void WfieldPlaneC2X(const double xpos, const double ypos, double& ex,
376 double& ey, double& volt, const int iplane,
377 const bool opt);
378 void WfieldPlaneC2Y(const double xpos, const double ypos, double& ex,
379 double& ey, double& volt, const int iplane,
380 const bool opt);
381 void WfieldPlaneC30(const double xpos, const double ypos, double& ex,
382 double& ey, double& volt, const int iplane,
383 const bool opt);
384 void WfieldPlaneD10(const double xpos, const double ypos, double& ex,
385 double& ey, double& volt, const int iplane,
386 const bool opt);
387 void WfieldPlaneD30(const double xpos, const double ypos, double& ex,
388 double& ey, double& volt, const int iplane,
389 const bool opt);
390 void WfieldStripZ(const double xpos, const double ypos, double& ex,
391 double& ey, double& volt, const int ip, const int is,
392 const bool opt);
393 void WfieldStripXy(const double xpos, const double ypos, const double zpos,
394 double& ex, double& ey, double& ez, double& volt,
395 const int ip, const int is, const bool opt);
396
397 // Auxiliary functions for C type cells
398 double Ph2(const double xpos, const double ypos);
399 double Ph2Lim(const double radius) {
400 return -log(abs(zmult) * radius * (1. - 3. * p1 + 5. * p2));
401 }
402 void E2Sum(const double xpos, const double ypos, double& ex, double& ey);
403
404 // Mapping function for D30 type cells
405 void ConformalMap(std::complex<double> z, std::complex<double>& ww,
406 std::complex<double>& wd);
407 void InitializeCoefficientTables();
408
409 bool InTube(const double x0, const double y0, const double a, const int n);
410
411 // Transformation between cartesian and polar coordinates
412 void Cartesian2Polar(const double x0, const double y0, double& r,
413 double& theta) {
414
415 if (x0 == 0. && y0 == 0.) {
416 r = theta = 0.;
417 return;
418 }
419 r = sqrt(x0 * x0 + y0 * y0);
420 theta = 180. * atan2(y0, x0) / Pi;
421 }
422
423 void Polar2Cartesian(const double r, const double theta, double& x0,
424 double& y0) {
425
426 x0 = r * cos(Pi * theta / 180.);
427 y0 = r * sin(Pi * theta / 180.);
428 }
429
430 // Transformation (r, theta) to (rho, phi) via the map
431 // (r, theta) = (exp(rho), 180 * phi / Pi).
432 void RTheta2RhoPhi(const double rho, const double phi, double& r,
433 double& theta) {
434
435 r = exp(rho);
436 theta = 180. * phi / Pi;
437 }
438};
439}
440
441#endif
cos
DoubleAc cos(const DoubleAc &f)
Definition: DoubleAc.cpp:431
sqrt
DoubleAc sqrt(const DoubleAc &f)
Definition: DoubleAc.cpp:313
sin
DoubleAc sin(const DoubleAc &f)
Definition: DoubleAc.cpp:383
exp
DoubleAc exp(const DoubleAc &f)
Definition: DoubleAc.cpp:376
Garfield::ComponentAnalyticField::GetVoltageRange
bool GetVoltageRange(double &pmin, double &pmax)
Definition: ComponentAnalyticField.cc:87
Garfield::ComponentAnalyticField::AddPlaneX
void AddPlaneX(const double x, const double voltage, const std::string label)
Definition: ComponentAnalyticField.cc:415
Garfield::ComponentAnalyticField::AddWire
void AddWire(const double x, const double y, const double diameter, const double voltage, const std::string label, const double length=100., const double tension=50., const double rho=19.3, const int ntrap=5)
Definition: ComponentAnalyticField.cc:317
Garfield::ComponentAnalyticField::IsInTrapRadius
bool IsInTrapRadius(double x0, double y0, double z0, double &xw, double &yx, double &rw)
Definition: ComponentAnalyticField.cc:256
Garfield::ComponentAnalyticField::AddCharge
void AddCharge(const double x, const double y, const double z, const double q)
Definition: ComponentAnalyticField.cc:682
Garfield::ComponentAnalyticField::AddReadout
void AddReadout(const std::string label)
Definition: ComponentAnalyticField.cc:1938
Garfield::ComponentAnalyticField::AddStripOnPlaneX
void AddStripOnPlaneX(const char direction, const double x, const double smin, const double smax, const std::string label, const double gap=-1.)
Definition: ComponentAnalyticField.cc:471
Garfield::ComponentAnalyticField::AddTube
void AddTube(const double radius, const double voltage, const int nEdges, const std::string label)
Definition: ComponentAnalyticField.cc:374
Garfield::ComponentAnalyticField::IsWireCrossed
bool IsWireCrossed(double x0, double y0, double z0, double x1, double y1, double z1, double &xc, double &yc, double &zc)
Definition: ComponentAnalyticField.cc:184
Garfield::ComponentAnalyticField::GetPlaneX
bool GetPlaneX(const int i, double &x, double &voltage, std::string &label)
Definition: ComponentAnalyticField.cc:760
Garfield::ComponentAnalyticField::WeightingField
void WeightingField(const double x, const double y, const double z, double &wx, double &wy, double &wz, const std::string label)
Definition: ComponentAnalyticField.cc:104
Garfield::ComponentAnalyticField::AddStripOnPlaneY
void AddStripOnPlaneY(const char direction, const double y, const double smin, const double smax, const std::string label, const double gap=-1.)
Definition: ComponentAnalyticField.cc:524
Garfield::ComponentAnalyticField::GetBoundingBox
bool GetBoundingBox(double &x0, double &y0, double &z0, double &x1, double &y1, double &z1)
Definition: ComponentAnalyticField.cc:165
Garfield::ComponentAnalyticField::AddPlaneY
void AddPlaneY(const double y, const double voltage, const std::string label)
Definition: ComponentAnalyticField.cc:443
Garfield::ComponentAnalyticField::WeightingPotential
double WeightingPotential(const double x, const double y, const double z, const std::string label)
Definition: ComponentAnalyticField.cc:134
Garfield::ComponentAnalyticField::GetWire
bool GetWire(const int i, double &x, double &y, double &diameter, double &voltage, std::string &label, double &length, double &charge, int &ntrap)
Definition: ComponentAnalyticField.cc:738
Garfield::ComponentAnalyticField::ElectricField
void ElectricField(const double x, const double y, const double z, double &ex, double &ey, double &ez, Medium *&m, int &status)
Definition: ComponentAnalyticField.cc:18
Garfield::ComponentAnalyticField::GetPlaneY
bool GetPlaneY(const int i, double &y, double &voltage, std::string &label)
Definition: ComponentAnalyticField.cc:775
Garfield::ComponentAnalyticField::GetTube
bool GetTube(double &r, double &voltage, int &nEdges, std::string &label)
Definition: ComponentAnalyticField.cc:790