Line | Branch | Exec | Source |
---|---|---|---|
1 | #include "Algorithm.h" | ||
2 | #include <cmath> | ||
3 | |||
4 | namespace iguana::clas12 { | ||
5 | |||
6 | REGISTER_IGUANA_ALGORITHM(MomentumCorrection); | ||
7 | |||
8 | 4 | void MomentumCorrection::Start(hipo::banklist& banks) | |
9 | { | ||
10 |
1/2✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
|
4 | b_particle = GetBankIndex(banks, "REC::Particle"); |
11 |
1/2✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
|
4 | b_sector = GetBankIndex(banks, "REC::Particle::Sector"); |
12 |
1/2✓ Branch 0 taken 4 times.
✗ Branch 1 not taken.
|
4 | b_config = GetBankIndex(banks, "RUN::config"); |
13 | 4 | } | |
14 | |||
15 | |||
16 | 2100 | void MomentumCorrection::Run(hipo::banklist& banks) const | |
17 | { | ||
18 |
1/2✓ Branch 0 taken 2100 times.
✗ Branch 1 not taken.
|
2100 | auto& particleBank = GetBank(banks, b_particle, "REC::Particle"); |
19 |
1/2✓ Branch 0 taken 2100 times.
✗ Branch 1 not taken.
|
2100 | auto& sectorBank = GetBank(banks, b_sector, "REC::Particle::Sector"); |
20 |
1/2✓ Branch 0 taken 2100 times.
✗ Branch 1 not taken.
|
4200 | auto& configBank = GetBank(banks, b_config, "RUN::config"); |
21 |
1/2✓ Branch 0 taken 2100 times.
✗ Branch 1 not taken.
|
2100 | ShowBank(particleBank, Logger::Header("INPUT PARTICLES")); |
22 | |||
23 | 2100 | auto torus = configBank.getFloat("torus", 0); | |
24 | |||
25 |
2/2✓ Branch 0 taken 9684 times.
✓ Branch 1 taken 2100 times.
|
11784 | for(auto const& row : particleBank.getRowList()) { |
26 | |||
27 |
1/2✓ Branch 0 taken 9684 times.
✗ Branch 1 not taken.
|
19368 | auto [px, py, pz] = Transform( |
28 | 9684 | particleBank.getFloat("px", row), | |
29 | 9684 | particleBank.getFloat("py", row), | |
30 | 9684 | particleBank.getFloat("pz", row), | |
31 | sectorBank.getInt("sector", row), | ||
32 | particleBank.getInt("pid", row), | ||
33 | torus); | ||
34 |
1/2✓ Branch 0 taken 9684 times.
✗ Branch 1 not taken.
|
9684 | particleBank.putFloat("px", row, px); |
35 |
1/2✓ Branch 0 taken 9684 times.
✗ Branch 1 not taken.
|
9684 | particleBank.putFloat("py", row, py); |
36 |
1/2✓ Branch 0 taken 9684 times.
✗ Branch 1 not taken.
|
9684 | particleBank.putFloat("pz", row, pz); |
37 | } | ||
38 | |||
39 |
1/2✓ Branch 0 taken 2100 times.
✗ Branch 1 not taken.
|
2100 | ShowBank(particleBank, Logger::Header("OUTPUT PARTICLES")); |
40 | 2100 | } | |
41 | |||
42 | |||
43 | 9877 | Momentum3 MomentumCorrection::Transform(vector_element_t const px, vector_element_t const py, vector_element_t const pz, int const sec, int const pid, float const torus) const | |
44 | { | ||
45 | // energy loss correction | ||
46 | auto e_cor = torus < 0 | ||
47 |
1/2✓ Branch 0 taken 9877 times.
✗ Branch 1 not taken.
|
9877 | ? EnergyLossInbending(px, py, pz, pid) |
48 | ✗ | : EnergyLossOutbending(px, py, pz, pid); | |
49 | // momentum correction | ||
50 | auto p_cor = torus < 0 | ||
51 |
1/2✓ Branch 0 taken 9877 times.
✗ Branch 1 not taken.
|
9877 | ? CorrectionInbending(e_cor * px, e_cor * py, e_cor * pz, sec, pid) |
52 | ✗ | : CorrectionOutbending(e_cor * px, e_cor * py, e_cor * pz, sec, pid); | |
53 | // return the corrected momentum | ||
54 | return { | ||
55 | 9877 | e_cor * p_cor * px, | |
56 | 9877 | e_cor * p_cor * py, | |
57 | 9877 | e_cor * p_cor * pz}; | |
58 | } | ||
59 | |||
60 | |||
61 | 9877 | double MomentumCorrection::CorrectionInbending(vector_element_t const px, vector_element_t const py, vector_element_t const pz, int const sec, int const pid) const | |
62 | { | ||
63 | |||
64 | // skip the correction if it's not defined | ||
65 |
4/4✓ Branch 0 taken 6295 times.
✓ Branch 1 taken 3582 times.
✓ Branch 2 taken 1500 times.
✓ Branch 3 taken 4795 times.
|
9877 | if(!(pid == particle::electron || pid == particle::pi_plus || pid == particle::pi_minus || pid == particle::proton)) |
66 | return 1.0; | ||
67 | |||
68 | // Momentum Magnitude | ||
69 | 5082 | double pp = sqrt(px * px + py * py + pz * pz); | |
70 | |||
71 | // Initializing the correction factor | ||
72 | double dp = 0; | ||
73 | |||
74 | // Defining Phi Angle | ||
75 | 5082 | double Phi = (180 / M_PI) * atan2(py, px); | |
76 | |||
77 | // (Initial) Shift of the Phi Angle (done to realign sectors whose data is separated when plotted from ±180˚) | ||
78 |
8/8✓ Branch 0 taken 1006 times.
✓ Branch 1 taken 4076 times.
✓ Branch 2 taken 844 times.
✓ Branch 3 taken 162 times.
✓ Branch 4 taken 1023 times.
✓ Branch 5 taken 3897 times.
✓ Branch 6 taken 1011 times.
✓ Branch 7 taken 12 times.
|
5082 | if(((sec == 4 || sec == 3) && Phi < 0) || (sec > 4 && Phi < 90)) { |
79 | 1173 | Phi += 360; | |
80 | } | ||
81 | |||
82 | // Getting Local Phi Angle | ||
83 | 5082 | double PhiLocal = Phi - (sec - 1) * 60; | |
84 | |||
85 | // Applying Shift Functions to Phi Angles (local shifted phi = phi) | ||
86 | double phi = PhiLocal; | ||
87 | |||
88 | // For Electron Shift | ||
89 |
2/2✓ Branch 0 taken 1444 times.
✓ Branch 1 taken 3638 times.
|
5082 | if(pid == particle::electron) { |
90 | 1444 | phi = PhiLocal - 30 / pp; | |
91 | } | ||
92 | |||
93 | // For π+ Pion/Proton Shift | ||
94 |
2/2✓ Branch 0 taken 2652 times.
✓ Branch 1 taken 2430 times.
|
5082 | if(pid == particle::pi_plus || pid == particle::proton) { |
95 | 2652 | phi = PhiLocal + (32 / (pp - 0.05)); | |
96 | } | ||
97 | |||
98 | // For π- Pion Shift | ||
99 |
2/2✓ Branch 0 taken 986 times.
✓ Branch 1 taken 4096 times.
|
5082 | if(pid == particle::pi_minus) { |
100 | 986 | phi = PhiLocal - (32 / (pp - 0.05)); | |
101 | } | ||
102 | |||
103 | ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
104 | //==================================================================================================================================// | ||
105 | //=======================//=======================// Electron Corrections //=======================//=======================// | ||
106 | //==================================================================================================================================// | ||
107 | ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
108 |
2/2✓ Branch 0 taken 1444 times.
✓ Branch 1 taken 3638 times.
|
5082 | if(pid == particle::electron) { |
109 |
2/2✓ Branch 0 taken 121 times.
✓ Branch 1 taken 1323 times.
|
1444 | if(sec == 1) { |
110 | // The CONTINUOUS QUADRATIC function predicted for ∆p_{El} for [Cor = Uncorrected][Sector 1] is: | ||
111 | 121 | dp = ((-4.3303e-06) * phi * phi + (1.1006e-04) * phi + (-5.7235e-04)) * pp * pp + ((3.2555e-05) * phi * phi + (-0.0014559) * phi + (0.0014878)) * pp + ((-1.9577e-05) * phi * phi + (0.0017996) * phi + (0.025963)); | |
112 | } | ||
113 |
2/2✓ Branch 0 taken 57 times.
✓ Branch 1 taken 1387 times.
|
1444 | if(sec == 2) { |
114 | // The CONTINUOUS QUADRATIC function predicted for ∆p_{El} for [Cor = Uncorrected][Sector 2] is: | ||
115 | 57 | dp = ((-9.8045e-07) * phi * phi + (6.7395e-05) * phi + (-4.6757e-05)) * pp * pp + ((-1.4958e-05) * phi * phi + (-0.0011191) * phi + (-0.0025143)) * pp + ((1.2699e-04) * phi * phi + (0.0033121) * phi + (0.020819)); | |
116 | } | ||
117 |
2/2✓ Branch 0 taken 53 times.
✓ Branch 1 taken 1391 times.
|
1444 | if(sec == 3) { |
118 | // The CONTINUOUS QUADRATIC function predicted for ∆p_{El} for [Cor = Uncorrected][Sector 3] is: | ||
119 | 53 | dp = ((-5.9459e-07) * phi * phi + (-2.8289e-05) * phi + (-4.3541e-04)) * pp * pp + ((-1.5025e-05) * phi * phi + (5.7730e-04) * phi + (-0.0077582)) * pp + ((7.3348e-05) * phi * phi + (-0.001102) * phi + (0.057052)); | |
120 | } | ||
121 |
2/2✓ Branch 0 taken 23 times.
✓ Branch 1 taken 1421 times.
|
1444 | if(sec == 4) { |
122 | // The CONTINUOUS QUADRATIC function predicted for ∆p_{El} for [Cor = Uncorrected][Sector 4] is: | ||
123 | 23 | dp = ((-2.2714e-06) * phi * phi + (-3.0360e-05) * phi + (-8.9322e-04)) * pp * pp + ((2.9737e-05) * phi * phi + (5.1142e-04) * phi + (0.0045641)) * pp + ((-1.0582e-04) * phi * phi + (-5.6852e-04) * phi + (0.027506)); | |
124 | } | ||
125 |
2/2✓ Branch 0 taken 35 times.
✓ Branch 1 taken 1409 times.
|
1444 | if(sec == 5) { |
126 | // The CONTINUOUS QUADRATIC function predicted for ∆p_{El} for [Cor = Uncorrected][Sector 5] is: | ||
127 | 35 | dp = ((-1.1490e-06) * phi * phi + (-6.2147e-06) * phi + (-4.7235e-04)) * pp * pp + ((3.7039e-06) * phi * phi + (-1.5943e-04) * phi + (-8.5238e-04)) * pp + ((4.4069e-05) * phi * phi + (0.0014152) * phi + (0.031933)); | |
128 | } | ||
129 |
2/2✓ Branch 0 taken 40 times.
✓ Branch 1 taken 1404 times.
|
1444 | if(sec == 6) { |
130 | // The CONTINUOUS QUADRATIC function predicted for ∆p_{El} for [Cor = Uncorrected][Sector 6] is: | ||
131 | 40 | dp = ((1.1076e-06) * phi * phi + (4.0156e-05) * phi + (-1.6341e-04)) * pp * pp + ((-2.8613e-05) * phi * phi + (-5.1861e-04) * phi + (-0.0056437)) * pp + ((1.2419e-04) * phi * phi + (4.9084e-04) * phi + (0.049976)); | |
132 | } | ||
133 | } | ||
134 | |||
135 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
136 | //====================================================================================================================================// | ||
137 | //=========================//=========================// π+ Corrections //=========================//=========================// | ||
138 | //====================================================================================================================================// | ||
139 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
140 |
2/2✓ Branch 0 taken 2138 times.
✓ Branch 1 taken 2944 times.
|
5082 | if(pid == particle::pi_plus) { |
141 |
2/2✓ Branch 0 taken 386 times.
✓ Branch 1 taken 1752 times.
|
2138 | if(sec == 1) { |
142 | 386 | dp = ((-5.4904e-07) * phi * phi + (-1.4436e-05) * phi + (3.1534e-04)) * pp * pp + ((3.8231e-06) * phi * phi + (3.6582e-04) * phi + (-0.0046759)) * pp + ((-5.4913e-06) * phi * phi + (-4.0157e-04) * phi + (0.010767)); | |
143 | 386 | dp = dp + ((6.1103e-07) * phi * phi + (5.5291e-06) * phi + (-1.9120e-04)) * pp * pp + ((-3.2300e-06) * phi * phi + (1.5377e-05) * phi + (7.5279e-04)) * pp + ((2.1434e-06) * phi * phi + (-6.9572e-06) * phi + (-7.9333e-05)); | |
144 | 386 | dp = dp + ((-1.3049e-06) * phi * phi + (1.1295e-05) * phi + (4.5797e-04)) * pp * pp + ((9.3122e-06) * phi * phi + (-5.1074e-05) * phi + (-0.0030757)) * pp + ((-1.3102e-05) * phi * phi + (2.2153e-05) * phi + (0.0040938)); | |
145 | } | ||
146 |
2/2✓ Branch 0 taken 278 times.
✓ Branch 1 taken 1860 times.
|
2138 | if(sec == 2) { |
147 | 278 | dp = ((-1.0087e-06) * phi * phi + (2.1319e-05) * phi + (7.8641e-04)) * pp * pp + ((6.7485e-06) * phi * phi + (7.3716e-05) * phi + (-0.0094591)) * pp + ((-1.1820e-05) * phi * phi + (-3.8103e-04) * phi + (0.018936)); | |
148 | 278 | dp = dp + ((8.8155e-07) * phi * phi + (-2.8257e-06) * phi + (-2.6729e-04)) * pp * pp + ((-5.4499e-06) * phi * phi + (3.8397e-05) * phi + (0.0015914)) * pp + ((6.8926e-06) * phi * phi + (-5.9386e-05) * phi + (-0.0021749)); | |
149 | 278 | dp = dp + ((-2.0147e-07) * phi * phi + (1.1061e-05) * phi + (3.8827e-04)) * pp * pp + ((4.9294e-07) * phi * phi + (-6.0257e-05) * phi + (-0.0022087)) * pp + ((9.8548e-07) * phi * phi + (5.9047e-05) * phi + (0.0022905)); | |
150 | } | ||
151 |
2/2✓ Branch 0 taken 337 times.
✓ Branch 1 taken 1801 times.
|
2138 | if(sec == 3) { |
152 | 337 | dp = ((8.6722e-08) * phi * phi + (-1.7975e-05) * phi + (4.8118e-05)) * pp * pp + ((2.6273e-06) * phi * phi + (3.1453e-05) * phi + (-0.0015943)) * pp + ((-6.4463e-06) * phi * phi + (-5.8990e-05) * phi + (0.0041703)); | |
153 | 337 | dp = dp + ((9.6317e-07) * phi * phi + (-1.7659e-06) * phi + (-8.8318e-05)) * pp * pp + ((-5.1346e-06) * phi * phi + (8.3318e-06) * phi + (3.7723e-04)) * pp + ((3.9548e-06) * phi * phi + (-6.9614e-05) * phi + (2.1393e-04)); | |
154 | 337 | dp = dp + ((5.6438e-07) * phi * phi + (8.1678e-06) * phi + (-9.4406e-05)) * pp * pp + ((-3.9074e-06) * phi * phi + (-6.5174e-05) * phi + (5.4218e-04)) * pp + ((6.3198e-06) * phi * phi + (1.0611e-04) * phi + (-4.5749e-04)); | |
155 | } | ||
156 |
2/2✓ Branch 0 taken 263 times.
✓ Branch 1 taken 1875 times.
|
2138 | if(sec == 4) { |
157 | 263 | dp = ((4.3406e-07) * phi * phi + (-4.9036e-06) * phi + (2.3064e-04)) * pp * pp + ((1.3624e-06) * phi * phi + (3.2907e-05) * phi + (-0.0034872)) * pp + ((-5.1017e-06) * phi * phi + (2.4593e-05) * phi + (0.0092479)); | |
158 | 263 | dp = dp + ((6.0218e-07) * phi * phi + (-1.4383e-05) * phi + (-3.1999e-05)) * pp * pp + ((-1.1243e-06) * phi * phi + (9.3884e-05) * phi + (-4.1985e-04)) * pp + ((-1.8808e-06) * phi * phi + (-1.2222e-04) * phi + (0.0014037)); | |
159 | 263 | dp = dp + ((-2.5490e-07) * phi * phi + (-8.5120e-07) * phi + (7.9109e-05)) * pp * pp + ((2.5879e-06) * phi * phi + (8.6108e-06) * phi + (-5.1533e-04)) * pp + ((-4.4521e-06) * phi * phi + (-1.7012e-05) * phi + (7.4848e-04)); | |
160 | } | ||
161 |
2/2✓ Branch 0 taken 335 times.
✓ Branch 1 taken 1803 times.
|
2138 | if(sec == 5) { |
162 | 335 | dp = ((2.4292e-07) * phi * phi + (8.8741e-06) * phi + (2.9482e-04)) * pp * pp + ((3.7229e-06) * phi * phi + (7.3215e-06) * phi + (-0.0050685)) * pp + ((-1.1974e-05) * phi * phi + (-1.3043e-04) * phi + (0.0078836)); | |
163 | 335 | dp = dp + ((1.0867e-06) * phi * phi + (-7.7630e-07) * phi + (-4.4930e-05)) * pp * pp + ((-5.6564e-06) * phi * phi + (-1.3417e-05) * phi + (2.5224e-04)) * pp + ((6.8460e-06) * phi * phi + (9.0495e-05) * phi + (-4.6587e-04)); | |
164 | 335 | dp = dp + ((8.5720e-07) * phi * phi + (-6.7464e-06) * phi + (-4.0944e-05)) * pp * pp + ((-4.7370e-06) * phi * phi + (5.8808e-05) * phi + (1.9047e-04)) * pp + ((5.7404e-06) * phi * phi + (-1.1105e-04) * phi + (-1.9392e-04)); | |
165 | } | ||
166 |
2/2✓ Branch 0 taken 272 times.
✓ Branch 1 taken 1866 times.
|
2138 | if(sec == 6) { |
167 | 272 | dp = ((2.1191e-06) * phi * phi + (-3.3710e-05) * phi + (2.5741e-04)) * pp * pp + ((-1.2915e-05) * phi * phi + (2.3753e-04) * phi + (-2.6882e-04)) * pp + ((2.2676e-05) * phi * phi + (-2.3115e-04) * phi + (-0.001283)); | |
168 | 272 | dp = dp + ((6.0270e-07) * phi * phi + (-6.8200e-06) * phi + (1.3103e-04)) * pp * pp + ((-1.8745e-06) * phi * phi + (3.8646e-05) * phi + (-8.8056e-04)) * pp + ((2.0885e-06) * phi * phi + (-3.4932e-05) * phi + (4.5895e-04)); | |
169 | 272 | dp = dp + ((4.7349e-08) * phi * phi + (-5.7528e-06) * phi + (-3.4097e-06)) * pp * pp + ((1.7731e-06) * phi * phi + (3.5865e-05) * phi + (-5.7881e-04)) * pp + ((-9.7008e-06) * phi * phi + (-4.1836e-05) * phi + (0.0035403)); | |
170 | } | ||
171 | } | ||
172 | |||
173 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
174 | //====================================================================================================================================// | ||
175 | //==================//==================// π- Corrections (Updated as of 01-13-2023) //==================//==================// | ||
176 | //====================================================================================================================================// | ||
177 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
178 |
2/2✓ Branch 0 taken 986 times.
✓ Branch 1 taken 4096 times.
|
5082 | if(pid == particle::pi_minus) { |
179 |
2/2✓ Branch 0 taken 145 times.
✓ Branch 1 taken 841 times.
|
986 | if(sec == 1) { |
180 | |||
181 | 145 | dp = (-9.2163E-07 * phi * phi + 3.1862E-06 * phi + 2.9805E-03) * pp * pp + (1.0435E-05 * phi * phi + -8.7298E-05 * phi + -1.7730E-02) * pp + -1.5154E-05 * phi * phi + -1.3716E-04 * phi + 2.2410E-02; | |
182 | } | ||
183 |
2/2✓ Branch 0 taken 124 times.
✓ Branch 1 taken 862 times.
|
986 | if(sec == 2) { |
184 | |||
185 | 124 | dp = (-1.9656E-06 * phi * phi + 9.7389E-05 * phi + 4.1250E-03) * pp * pp + (1.6439E-06 * phi * phi + -4.6007E-04 * phi + -1.9809E-02) * pp + 3.5794E-07 * phi * phi + 4.8250E-04 * phi + 1.7333E-02; | |
186 | } | ||
187 | |||
188 |
2/2✓ Branch 0 taken 116 times.
✓ Branch 1 taken 870 times.
|
986 | if(sec == 3) { |
189 | |||
190 | 116 | dp = (2.5351E-06 * phi * phi + 4.1043E-05 * phi + 3.1157E-03) * pp * pp + (-1.3573E-05 * phi * phi + -1.7609E-04 * phi + -1.6759E-02) * pp + 1.4647E-05 * phi * phi + 1.7484E-04 * phi + 1.3805E-02; | |
191 | } | ||
192 | |||
193 |
2/2✓ Branch 0 taken 116 times.
✓ Branch 1 taken 870 times.
|
986 | if(sec == 4) { |
194 | |||
195 | 116 | dp = (2.3500E-06 * phi * phi + -7.7894E-05 * phi + 4.4837E-03) * pp * pp + (-9.7915E-06 * phi * phi + 4.6576E-04 * phi + -2.6809E-02) * pp + 1.3819E-05 * phi * phi + -5.6017E-04 * phi + 3.0320E-02; | |
196 | } | ||
197 | |||
198 |
2/2✓ Branch 0 taken 92 times.
✓ Branch 1 taken 894 times.
|
986 | if(sec == 5) { |
199 | |||
200 | 92 | dp = (-2.1809E-06 * phi * phi + 2.4948E-05 * phi + 2.7995E-03) * pp * pp + (6.3908E-06 * phi * phi + -6.5122E-05 * phi + -1.7571E-02) * pp + -1.9146E-06 * phi * phi + -6.3799E-05 * phi + 2.0877E-02; | |
201 | } | ||
202 | |||
203 |
2/2✓ Branch 0 taken 151 times.
✓ Branch 1 taken 835 times.
|
986 | if(sec == 6) { |
204 | |||
205 | 151 | dp = (-9.3043E-06 * phi * phi + 6.2678E-05 * phi + 5.9660E-03) * pp * pp + (4.0581E-05 * phi * phi + -3.0537E-04 * phi + -3.1485E-02) * pp + -3.8345E-05 * phi * phi + 2.0267E-04 * phi + 3.3363E-02; | |
206 | } | ||
207 | } | ||
208 | |||
209 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
210 | //====================================================================================================================================// | ||
211 | //=======================//=======================// All Proton Corrections //=======================//=======================// | ||
212 | //====================================================================================================================================// | ||
213 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
214 |
2/2✓ Branch 0 taken 514 times.
✓ Branch 1 taken 4568 times.
|
5082 | if(pid == particle::proton) { |
215 |
2/2✓ Branch 0 taken 56 times.
✓ Branch 1 taken 458 times.
|
514 | if(sec == 1) { |
216 | 56 | dp = ((1 + std::copysign(1, (pp - 1.4))) / 2) * ((4.4034e-03) * pp + (-0.01703)) + ((1 + std::copysign(1, -(pp - 1.4))) / 2) * ((-0.10898) * (pp - 1.4) * (pp - 1.4) + (-0.09574) * (pp - 1.4) + ((4.4034e-03) * 1.4 + (-0.01703))); | |
217 | } | ||
218 |
2/2✓ Branch 0 taken 50 times.
✓ Branch 1 taken 464 times.
|
514 | if(sec == 2) { |
219 | 50 | dp = ((1 + std::copysign(1, (pp - 1.5))) / 2) * ((0.01318) * pp + (-0.03403)) + ((1 + std::copysign(1, -(pp - 1.5))) / 2) * ((-0.09829) * (pp - 1.5) * (pp - 1.5) + (-0.0986) * (pp - 1.5) + ((0.01318) * 1.5 + (-0.03403))); | |
220 | } | ||
221 |
2/2✓ Branch 0 taken 48 times.
✓ Branch 1 taken 466 times.
|
514 | if(sec == 3) { |
222 | 48 | dp = ((1 + std::copysign(1, (pp - 1.05))) / 2) * ((-4.7052e-03) * pp + (1.2410e-03)) + ((1 + std::copysign(1, -(pp - 1.05))) / 2) * ((-0.22721) * (pp - 1.05) * (pp - 1.05) + (-0.09702) * (pp - 1.05) + ((-4.7052e-03) * 1.05 + (1.2410e-03))); | |
223 | } | ||
224 |
2/2✓ Branch 0 taken 50 times.
✓ Branch 1 taken 464 times.
|
514 | if(sec == 4) { |
225 | 50 | dp = ((1 + std::copysign(1, (pp - 1.4))) / 2) * ((-1.0900e-03) * pp + (-4.0573e-03)) + ((1 + std::copysign(1, -(pp - 1.4))) / 2) * ((-0.09236) * (pp - 1.4) * (pp - 1.4) + (-0.073) * (pp - 1.4) + ((-1.0900e-03) * 1.4 + (-4.0573e-03))); | |
226 | } | ||
227 |
2/2✓ Branch 0 taken 52 times.
✓ Branch 1 taken 462 times.
|
514 | if(sec == 5) { |
228 | 52 | dp = ((1 + std::copysign(1, (pp - 1.5))) / 2) * ((7.3965e-03) * pp + (-0.02428)) + ((1 + std::copysign(1, -(pp - 1.5))) / 2) * ((-0.09539) * (pp - 1.5) * (pp - 1.5) + (-0.09263) * (pp - 1.5) + ((7.3965e-03) * 1.5 + (-0.02428))); | |
229 | } | ||
230 |
2/2✓ Branch 0 taken 46 times.
✓ Branch 1 taken 468 times.
|
514 | if(sec == 6) { |
231 | 46 | dp = ((1 + std::copysign(1, (pp - 1.15))) / 2) * ((-7.6214e-03) * pp + (8.1014e-03)) + ((1 + std::copysign(1, -(pp - 1.15))) / 2) * ((-0.12718) * (pp - 1.15) * (pp - 1.15) + (-0.06626) * (pp - 1.15) + ((-7.6214e-03) * 1.15 + (8.1014e-03))); | |
232 | } | ||
233 | } | ||
234 | |||
235 | 5082 | return dp / pp + 1; | |
236 | } | ||
237 | |||
238 | |||
239 | ✗ | double MomentumCorrection::CorrectionOutbending(vector_element_t const px, vector_element_t const py, vector_element_t const pz, int const sec, int const pid) const | |
240 | { | ||
241 | |||
242 | // skip the correction if it's not defined | ||
243 | ✗ | if(!(pid == particle::electron || pid == particle::pi_plus || pid == particle::pi_minus)) | |
244 | return 1.0; | ||
245 | |||
246 | // Momentum Magnitude | ||
247 | ✗ | double pp = sqrt(px * px + py * py + pz * pz); | |
248 | |||
249 | // Initializing the correction factor | ||
250 | double dp = 0; | ||
251 | |||
252 | // Defining Phi Angle | ||
253 | ✗ | double Phi = (180 / M_PI) * atan2(py, px); | |
254 | |||
255 | // (Initial) Shift of the Phi Angle (done to realign sectors whose data is separated when plotted from ±180˚) | ||
256 | ✗ | if(((sec == 4 || sec == 3) && Phi < 0) || (sec > 4 && Phi < 90)) { | |
257 | ✗ | Phi += 360; | |
258 | } | ||
259 | |||
260 | // Getting Local Phi Angle | ||
261 | ✗ | double PhiLocal = Phi - (sec - 1) * 60; | |
262 | |||
263 | // Applying Shift Functions to Phi Angles (local shifted phi = phi) | ||
264 | double phi = PhiLocal; | ||
265 | // For Electron Shift | ||
266 | ✗ | if(pid == particle::electron) { | |
267 | ✗ | phi = PhiLocal - 30 / pp; | |
268 | } | ||
269 | // For π+ Pion/Proton Shift | ||
270 | ✗ | if(pid == particle::pi_plus || pid == particle::proton) { | |
271 | ✗ | phi = PhiLocal + (32 / (pp - 0.05)); | |
272 | } | ||
273 | // For π- Pion Shift | ||
274 | ✗ | if(pid == particle::pi_minus) { | |
275 | ✗ | phi = PhiLocal - (32 / (pp - 0.05)); | |
276 | } | ||
277 | |||
278 | ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
279 | //=======================//=======================// Electron Corrections //=======================//=======================// | ||
280 | ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
281 | ✗ | if(pid == particle::electron) { | |
282 | ✗ | if(sec == 1) { | |
283 | ✗ | dp = ((1.3189e-06) * phi * phi + (4.26057e-05) * phi + (-0.002322628)) * pp * pp + ((-1.1409e-05) * phi * phi + (2.2188e-05) * phi + (0.02878927)) * pp + ((2.4950e-05) * phi * phi + (1.6170e-06) * phi + (-0.061816275)); | |
284 | } | ||
285 | ✗ | if(sec == 2) { | |
286 | ✗ | dp = ((-2.9240e-07) * phi * phi + (3.2448e-07) * phi + (-0.001848308)) * pp * pp + ((4.4500e-07) * phi * phi + (4.76324e-04) * phi + (0.02219469)) * pp + ((6.9220e-06) * phi * phi + (-0.00153517) * phi + (-0.0479058)); | |
287 | } | ||
288 | ✗ | if(sec == 3) { | |
289 | ✗ | dp = ((2.71911e-06) * phi * phi + (1.657148e-05) * phi + (-0.001822211)) * pp * pp + ((-4.96814e-05) * phi * phi + (-3.761117e-04) * phi + (0.02564148)) * pp + ((1.97748e-04) * phi * phi + (9.58259e-04) * phi + (-0.05818292)); | |
290 | } | ||
291 | ✗ | if(sec == 4) { | |
292 | ✗ | dp = ((1.90966e-06) * phi * phi + (-2.4761e-05) * phi + (-0.00231562)) * pp * pp + ((-2.3927e-05) * phi * phi + (2.25262e-04) * phi + (0.0291831)) * pp + ((8.0515e-05) * phi * phi + (-6.42098e-04) * phi + (-0.06159197)); | |
293 | } | ||
294 | ✗ | if(sec == 5) { | |
295 | ✗ | dp = ((-3.6760323e-06) * phi * phi + (4.04398e-05) * phi + (-0.0021967515)) * pp * pp + ((4.90857e-05) * phi * phi + (-4.37437e-04) * phi + (0.02494339)) * pp + ((-1.08257e-04) * phi * phi + (0.00146111) * phi + (-0.0648485)); | |
296 | } | ||
297 | ✗ | if(sec == 6) { | |
298 | ✗ | dp = ((-6.2488e-08) * phi * phi + (2.23173e-05) * phi + (-0.00227522)) * pp * pp + ((1.8372e-05) * phi * phi + (-7.5227e-05) * phi + (0.032636)) * pp + ((-6.6566e-05) * phi * phi + (-2.4450e-04) * phi + (-0.072293)); | |
299 | } | ||
300 | } | ||
301 | |||
302 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
303 | //=========================//=========================// π+ Corrections //=========================//=========================// | ||
304 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
305 | ✗ | if(pid == particle::pi_plus) { | |
306 | ✗ | if(sec == 1) { | |
307 | ✗ | dp = ((-1.7334e-06) * phi * phi + (1.45112e-05) * phi + (0.00150721)) * pp * pp + ((6.6234e-06) * phi * phi + (-4.81191e-04) * phi + (-0.0138695)) * pp + ((-3.23625e-06) * phi * phi + (2.79751e-04) * phi + (0.027726)); | |
308 | } | ||
309 | ✗ | if(sec == 2) { | |
310 | ✗ | dp = ((-4.475464e-06) * phi * phi + (-4.11573e-05) * phi + (0.00204557)) * pp * pp + ((2.468278e-05) * phi * phi + (9.3590e-05) * phi + (-0.015399)) * pp + ((-1.61547e-05) * phi * phi + (-2.4206e-04) * phi + (0.0231743)); | |
311 | } | ||
312 | ✗ | if(sec == 3) { | |
313 | ✗ | dp = ((-8.0374e-07) * phi * phi + (2.8728e-06) * phi + (0.00152163)) * pp * pp + ((5.1347e-06) * phi * phi + (3.71709e-04) * phi + (-0.0165735)) * pp + ((4.0105e-06) * phi * phi + (-5.289869e-04) * phi + (0.02175395)); | |
314 | } | ||
315 | ✗ | if(sec == 4) { | |
316 | ✗ | dp = ((-3.8790e-07) * phi * phi + (-4.78445e-05) * phi + (0.002324725)) * pp * pp + ((6.80543e-06) * phi * phi + (5.69358e-04) * phi + (-0.0199162)) * pp + ((-1.30264e-05) * phi * phi + (-5.91606e-04) * phi + (0.03202088)); | |
317 | } | ||
318 | ✗ | if(sec == 5) { | |
319 | ✗ | dp = ((2.198518e-06) * phi * phi + (-1.52535e-05) * phi + (0.001187761)) * pp * pp + ((-1.000264e-05) * phi * phi + (1.63976e-04) * phi + (-0.01429673)) * pp + ((9.4962e-06) * phi * phi + (-3.86691e-04) * phi + (0.0303695)); | |
320 | } | ||
321 | ✗ | if(sec == 6) { | |
322 | ✗ | dp = ((-3.92944e-07) * phi * phi + (1.45848e-05) * phi + (0.00120668)) * pp * pp + ((3.7899e-06) * phi * phi + (-1.98219e-04) * phi + (-0.0131312)) * pp + ((-3.9961e-06) * phi * phi + (-1.32883e-04) * phi + (0.0294497)); | |
323 | } | ||
324 | } | ||
325 | |||
326 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
327 | //=======================//=======================// π- Corrections //=======================//=======================// | ||
328 | //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// | ||
329 | ✗ | if(pid == particle::pi_minus) { | |
330 | ✗ | if(sec == 1) { | |
331 | ✗ | dp = (7.8044E-06 * phi * phi + -9.4703E-05 * phi + 4.6696E-03) * pp * pp + (-3.4668E-05 * phi * phi + 6.2280E-04 * phi + -2.4273E-02) * pp + 2.3566E-05 * phi * phi + -5.8519E-04 * phi + 3.9226E-02; | |
332 | } | ||
333 | ✗ | if(sec == 2) { | |
334 | ✗ | dp = (-4.6611E-06 * phi * phi + -8.1637E-05 * phi + 7.5013E-03) * pp * pp + (1.7616E-05 * phi * phi + 3.5439E-04 * phi + -3.7122E-02) * pp + -1.6286E-05 * phi * phi + -2.6545E-04 * phi + 4.5659E-02; | |
335 | } | ||
336 | ✗ | if(sec == 3) { | |
337 | ✗ | dp = (4.5270E-06 * phi * phi + 2.2578E-04 * phi + 5.9214E-03) * pp * pp + (-1.6419E-05 * phi * phi + -8.1776E-04 * phi + -3.2776E-02) * pp + 1.3734E-05 * phi * phi + 6.6125E-04 * phi + 4.5784E-02; | |
338 | } | ||
339 | ✗ | if(sec == 4) { | |
340 | ✗ | dp = (-1.3141E-06 * phi * phi + 1.9648E-04 * phi + 7.6109E-03 - 0.006) * pp * pp + (8.0912E-06 * phi * phi + -8.2672E-04 * phi + -4.0495E-02 + 0.03) * pp + -3.1380E-06 * phi * phi + 6.2211E-04 * phi + 5.3361E-02 - 0.04; | |
341 | } | ||
342 | ✗ | if(sec == 5) { | |
343 | ✗ | dp = (-5.4065E-06 * phi * phi + -1.6325E-05 * phi + 1.2269E-02 - 0.002) * pp * pp + (1.9512E-05 * phi * phi + 1.0228E-04 * phi + -6.2351E-02 + 0.01) * pp + -9.5023E-06 * phi * phi + -3.7997E-05 * phi + 7.1061E-02 - 0.02; | |
344 | } | ||
345 | ✗ | if(sec == 6) { | |
346 | ✗ | dp = (-1.1882E-05 * phi * phi + 2.0101E-04 * phi + 1.1635E-02 - 0.01) * pp * pp + (5.8488E-05 * phi * phi + -6.4709E-04 * phi + -5.3833E-02 + 0.05) * pp + -4.4462E-05 * phi * phi + 3.7529E-04 * phi + 6.2130E-02 - 0.06; | |
347 | } | ||
348 | } | ||
349 | |||
350 | ✗ | return dp / pp + 1; | |
351 | } | ||
352 | |||
353 | |||
354 | 9877 | double MomentumCorrection::EnergyLossInbending(vector_element_t const px, vector_element_t const py, vector_element_t const pz, int const pid) const | |
355 | { | ||
356 | |||
357 | // The following code is for the Energy Loss Corrections for the proton | ||
358 |
2/2✓ Branch 0 taken 514 times.
✓ Branch 1 taken 9363 times.
|
9877 | if(pid != particle::proton) |
359 | return 1.0; | ||
360 | |||
361 | double dE_loss = 0; | ||
362 | 514 | auto pro = sqrt(px * px + py * py + pz * pz); | |
363 | 514 | auto proth = atan2(sqrt(px * px + py * py), pz) * (180 / M_PI); | |
364 | |||
365 | // Inbending Energy Loss Correction // | ||
366 |
2/2✓ Branch 0 taken 222 times.
✓ Branch 1 taken 292 times.
|
514 | if(proth < 27) { |
367 | 222 | dE_loss = exp(-2.739 - 3.932 * pro) + 0.002907; | |
368 | } | ||
369 | else { | ||
370 | 292 | dE_loss = exp(-1.2 - 4.228 * pro) + 0.007502; | |
371 | } | ||
372 | 514 | return (pro + dE_loss) / pro; | |
373 | } | ||
374 | |||
375 | |||
376 | ✗ | double MomentumCorrection::EnergyLossOutbending(vector_element_t const px, vector_element_t const py, vector_element_t const pz, int const pid) const | |
377 | { | ||
378 | |||
379 | // The following code is for the Energy Loss Corrections for the proton | ||
380 | ✗ | if(pid != particle::proton) | |
381 | return 1.0; | ||
382 | |||
383 | double dE_loss = 0; | ||
384 | ✗ | auto pro = sqrt(px * px + py * py + pz * pz); | |
385 | ✗ | auto proth = atan2(sqrt(px * px + py * py), pz) * (180 / M_PI); | |
386 | |||
387 | // Outbending Energy Loss Correction // | ||
388 | ✗ | if(proth > 27) { | |
389 | ✗ | dE_loss = exp(-1.871 - 3.063 * pro) + 0.007517; | |
390 | } | ||
391 | ✗ | return (pro + dE_loss) / pro; | |
392 | } | ||
393 | |||
394 | |||
395 | 3 | void MomentumCorrection::Stop() | |
396 | { | ||
397 | 3 | } | |
398 | |||
399 | } | ||
400 |