GCC Code Coverage Report


Directory: ./
File: src/iguana/algorithms/clas12/MomentumCorrection/Algorithm.cc
Date: 2025-01-05 09:03:17
Exec Total Coverage
Lines: 115 177 65.0%
Functions: 6 8 75.0%
Branches: 94 172 54.7%

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