Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Rossegger.h
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file Rossegger.h
1 #ifndef __ROSSEGGER_H__
2 #define __ROSSEGGER_H__
3 
4 //
5 // Hello Space Charge Fans: (although you should hate space charge)
6 //
7 // This is a code that implements the calculations of space charge contributions
8 // In a cylindrical TPC. It uses the solutions discovered by Stefan Rossegger for
9 // the ALICE experiment. These solutions use various Bessel functions as a series
10 // solution to the "point charge in a conducting cylinder" problem imposed by
11 // The typical configuration of a TPC found at a collider.
12 //
13 // These calculations have only a single dimension [length] and we shall choose
14 // the units of cm throughout the calculations.
15 //
16 // TKH
17 // 12-3-2015
18 //
19 
20 #include <cmath>
21 #include <cstdio>
22 #include <map>
23 #include <string>
24 
25 class TH2;
26 class TH3;
27 
28 #define NumberOfOrders 15 // Convergence problems after 15; Rossegger used 30
29 
30 class Rossegger
31 {
32  public:
33  explicit Rossegger(std::string filename);
34  Rossegger(double a = 30, double b = 80, double L = 80, double epsilon = 1E-4);
35  virtual ~Rossegger() {}
36 
37  void Verbosity(int v)
38  {
39  printf("verbosity set to %d. was %d\n", v, verbosity);
40  verbosity = v;
41  return;
42  };
43  double Rmn(int m, int n, double r); // Rmn function from Rossegger
44  double Rmn_for_zeroes(int m, double x); // Rmn function from Rossegger, as used to find Betamn zeroes.
45  double Rmn1(int m, int n, double r); // Rmn1 function from Rossegger
46  double Rmn2(int m, int n, double r); // Rmn2 function from Rossegger
47  double RPrime(int m, int n, double a, double r); // RPrime function from Rossegger
48 
49  double Rnk(int n, int k, double r); // Rnk function from Rossegger
50  double Rnk_for_zeroes(int n, double mu); // Rnk function from Rossegger, as used to find munk zeroes.
51 
52  double Limu(double mu, double x); // Bessel functions of purely imaginary order
53  double Kimu(double mu, double x); // Bessel functions of purely imaginary order
54 
55  double Ez(double r, double phi, double z, double r1, double phi1, double z1);
56  double Er(double r, double phi, double z, double r1, double phi1, double z1);
57  double Ephi(double r, double phi, double z, double r1, double phi1, double z1);
58 
59  // alternate versions that don't use precalc constants.
60  double Rmn_(int m, int n, double r); // Rmn function from Rossegger
61  // Rmn_for_zeroes doesn't have a way to speed it up with precalcs.
62  double Rmn1_(int m, int n, double r); // Rmn1 function from Rossegger
63  double Rmn2_(int m, int n, double r); // Rmn2 function from Rossegger
64  double RPrime_(int m, int n, double a, double r); // RPrime function from Rossegger
65  double Rnk_(int n, int k, double r); // Rnk function from Rossegger
66  double Rnk_for_zeroes_(int n, double mu); // Rnk function from Rossegger, as used to find munk zeroes.
67  double Ez_(double r, double phi, double z, double r1, double phi1, double z1);
68  double Er_(double r, double phi, double z, double r1, double phi1, double z1);
69  double Ephi_(double r, double phi, double z, double r1, double phi1, double z1);
70 
71  protected:
72  bool fByFile = false;
73  double a = NAN;
74  double b = NAN;
75  double L = NAN; // InnerRadius, OuterRadius, Length of 1/2 the TPC.
76  int verbosity = 0;
77  double pi = M_PI;
78  double epsilon = NAN; // precision.
79 
80  bool tweak = false;
81 
82  double MinimumDR = NAN;
83  double MinimumDPHI = NAN;
84  double MinimumDZ = NAN;
85 
86  double FindNextZero(double xstart, double epsilon, int order, double (Rossegger::*func)(int, double)); // Routine to find zeroes of func.
87  void FindBetamn(double epsilon); // Routine used to fill the Betamn array with resolution epsilon...
88  void FindMunk(double epsilon); // Routine used to fill the Munk array with resolution epsilon...
89  bool CheckZeroes(double epsilon); // confirm that the zeroes match to the desired precision.
90 
91  void LoadZeroes(const char* destfile);
92  void SaveZeroes(const char* destfile);
93 
94  double Betamn[NumberOfOrders][NumberOfOrders]{}; // Betamn array from Rossegger
95  double N2mn[NumberOfOrders][NumberOfOrders]{}; // N2mn array from Rossegger
96  double Munk[NumberOfOrders][NumberOfOrders]{}; // Munk array from Rossegger
97  double N2nk[NumberOfOrders][NumberOfOrders]{}; // N2nk array from Rossegger
98 
99  void PrecalcFreeConstants(); // Routine used to fill the arrays of other values used repeatedly:
100  // helpful values to precompute for frequent use:
101  double BetaN[NumberOfOrders]{}; // BetaN=(n+1)*pi/L as used in eg 5.32, .46
102  double BetaN_a[NumberOfOrders]{}; // BetaN*a as used in eg 5.32, .46
103  double BetaN_b[NumberOfOrders]{}; // BetaN*b as used in eg 5.32, .46
104  double km_BetaN_a[NumberOfOrders][NumberOfOrders]{}; // kn(m,BetaN*a) as used in Rossegger 5.32
105  double im_BetaN_a[NumberOfOrders][NumberOfOrders]{}; // in(m,BetaN*a) as used in Rossegger 5.32
106  double km_BetaN_b[NumberOfOrders][NumberOfOrders]{}; // kn(m,BetaN*b) as used in Rossegger 5.33
107  double im_BetaN_b[NumberOfOrders][NumberOfOrders]{}; // in(m,BetaN*b) as used in Rossegger 5.33
108  double bessel_denominator[NumberOfOrders][NumberOfOrders]{}; // BesselI(m,BetaN*a)*BesselK(m,BetaN*b)-BesselI(m,BetaN*b)*BesselK(m,BetaN*a) as in Rossegger 5.65
109 
110  void PrecalcDerivedConstants(); // Routine used to fill repeated values that depend on the Betamn and Munk zeroes:
111  double ym_Betamn_a[NumberOfOrders][NumberOfOrders]{}; // yn(m,Betamn[m][n]*a) as used in Rossegger 5.11
112  double jm_Betamn_a[NumberOfOrders][NumberOfOrders]{}; // jn(m,Betamn[m][n]*a) as used in Rossegger 5.11
113  double liMunk_BetaN_a[NumberOfOrders][NumberOfOrders]{}; // limu(Munk[n][k],BetaN*a) as used in Rossegger 5.45
114  double kiMunk_BetaN_a[NumberOfOrders][NumberOfOrders]{}; // kimu(Munk[n][k],BetaN*a) as used in Rossegger 5.45
115  double sinh_Betamn_L[NumberOfOrders][NumberOfOrders]{}; // sinh(Betamn[m][n]*L) as in Rossegger 5.64
116  double sinh_pi_Munk[NumberOfOrders][NumberOfOrders]{}; // sinh(pi*Munk[n][k]) as in Rossegger 5.66
117 
118  TH2* Tags = nullptr;
119  std::map<std::string, TH3*> Grid;
120 };
121 
122 #endif /* __SPACECHARGE_H__ */