BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
MdcDedxCorrection.cxx
Go to the documentation of this file.
1// BesIII MDC dE/dx Reconstruction Module
2// Class: MdcDedxCorrection
3// Created by Dayong WANG 2003/11
4
5// #include "CLHEP/config/CLHEP.h"
6#include "MdcDedxCorrection.h"
7#include "MdcDedxParam.h"
8#include "MdcDedxTrk.h"
9
10#include <cmath>
11#include <iostream>
12#include <vector>
13
14using namespace std;
15
16extern "C" {
17float prob_( float*, int* );
18}
19
21#ifdef DEBUG
22 std::cout << "MdcDedxCorrection constructed!!!" << std::endl;
23#endif
24}
25
27#ifdef DEBUG
28 std::cout << "MdcDedxCorrection destructed!!!" << std::endl;
29#endif
30}
31
32/* ------------- calculate the expects of dE/dx -------------*/
33void MdcDedxCorrection::dedx_pid_exp_old( int landau, int runflag, float dedx, int Nohit,
34 float mom, float theta, float t0, float lsamp,
35 double dedx_exp[5], double ex_sigma[5],
36 double pid_prob[5], double chi_dedx[5] ) const {
37 double par[5], sigma_par[4], sigma_index_nhit, sigma_index_sin;
38
39 if ( runflag == 1 )
40 {
46 sigma_par[0] = MdcDedxParam::HV1_sigmap0;
47 sigma_par[1] = MdcDedxParam::HV1_sigmap1;
48 sigma_par[2] = MdcDedxParam::HV1_sigmap2;
49 sigma_par[3] = MdcDedxParam::HV1_sigmap3;
50 sigma_index_nhit = MdcDedxParam::HV1_index_nhit;
51 sigma_index_sin = MdcDedxParam::HV1_index_sin;
52 }
53 else // must be runflag == 2
54 {
60 sigma_par[0] = MdcDedxParam::HV2_sigmap0;
61 sigma_par[1] = MdcDedxParam::HV2_sigmap1;
62 sigma_par[2] = MdcDedxParam::HV2_sigmap2;
63 sigma_par[3] = MdcDedxParam::HV2_sigmap3;
64 sigma_index_nhit = MdcDedxParam::HV2_index_nhit;
65 sigma_index_sin = MdcDedxParam::HV2_index_sin;
66 }
67
68 const int par_cand( 5 );
69 const float Charge_Mass[par_cand] = { 0.00051100, 0.10566, 0.13957, 0.4937, 0.93827 };
70 double beta_G, beta, betterm, bethe_B, sig_param;
71
72 int Nmax_prob( 0 );
73 float max_prob( -0.01 );
74 int ndf;
75 float chi2;
76
77 for ( int it = 0; it < par_cand; it++ )
78 {
79 beta_G = mom / Charge_Mass[it];
80 beta = beta_G / sqrt( 1 + ( beta_G ) * ( beta_G ) );
81 betterm = par[1] - log( par[2] + pow( 1 / beta_G, par[4] ) );
82 bethe_B = par[0] / pow( beta, par[3] ) * betterm - par[0];
83
84 if ( Nohit > 0 )
85 {
86 dedx_exp[it] = bethe_B;
87 double sig_the = std::sin( (double)theta );
88
89 if ( runflag < 3 && runflag > 0 )
90 {
91 if ( landau == 0 )
92 {
93 sig_param = 1.6 * std::sin( (double)theta ) / ( lsamp * double( Nohit ) );
94 ex_sigma[it] = 0.05 * bethe_B * sqrt( 50.0 * sig_param );
95 }
96 else
97 {
98 // currently use one sigmap0
99 if ( beta_G < 4 )
100 { sig_param = sigma_par[1] + sigma_par[2] * std::pow( beta_G, sigma_par[3] ); }
101 else { sig_param = sigma_par[0]; }
102 // double sig_the=std::sin( (double)theta );
103 sig_the = std::pow( sig_the, sigma_index_sin );
104 double sig_n;
105 sig_n = 35.0 / double( Nohit );
106 sig_n = std::pow( sig_n, sigma_index_nhit );
107 ex_sigma[it] = sig_param * sig_the * sig_n;
108 }
109 }
110
111 MdcDedxTrk dedxtrk;
112 double dedx_correc;
113 if ( runflag == 2 ) dedx_correc = dedxtrk.SpaceChargeCorrec( theta, mom, it, dedx );
114 else dedx_correc = dedx;
115 chi_dedx[it] = ( dedx_correc - dedx_exp[it] ) / ex_sigma[it];
116 chi2 = chi_dedx[it] * chi_dedx[it];
117 ndf = 1;
118 pid_prob[it] = prob_( &chi2, &ndf );
119 // if(it ==0 ) cout<<"runflag: "<<runflag<<" dedx : "<<dedx<<" chi_dedx: "
120 //<<chi_dedx[it] <<" ptrk: "<<mom<<endl;
121 if ( it == -999 )
122 { // here a debug flag
123 std::cout << " mom = " << mom << "exp" << dedx_exp[it] << " sigma " << ex_sigma[it]
124 << " prob " << pid_prob[it] << std::endl;
125 }
126 if ( pid_prob[it] > max_prob )
127 {
128 max_prob = pid_prob[it];
129 Nmax_prob = it;
130 }
131 }
132 else
133 {
134 dedx_exp[it] = 0.0;
135 ex_sigma[it] = 1000.0;
136 pid_prob[it] = 0.0;
137 chi_dedx[it] = 999.0;
138 }
139 }
140 // std::cout<<"MdcDedxCorrection::dedx_pid_exp(blum)!!!"<<std::endl;
141}
142
143void MdcDedxCorrection::dedx_pid_exp( int vflag[3], float dedx, int trkalg, int Nohit,
144 float mom, float theta, float t0, float lsamp,
145 double dedx_exp[5], double ex_sigma[5],
146 double pid_prob[5], double chi_dedx[5],
147 std::vector<double>& par,
148 std::vector<double>& sig_par ) const {
149 const int par_cand( 5 );
150 const float Charge_Mass[par_cand] = { 0.00051100, 0.10566, 0.13957, 0.4937, 0.93827 };
151 double beta_G, beta, betterm, bethe_B = 0;
152
153 int dedxflag = vflag[0];
154 int sigmaflag = vflag[1];
155 bool ifMC = false;
156 if ( vflag[2] == 1 ) ifMC = true;
157
158 int Nmax_prob( 0 );
159 float max_prob( -0.01 );
160 int ndf;
161 float chi2;
162
163 for ( int it = 0; it < par_cand; it++ )
164 {
165 beta_G = mom / Charge_Mass[it];
166
167 if ( dedxflag == 1 )
168 {
169 beta = beta_G / sqrt( 1 + ( beta_G ) * ( beta_G ) );
170 betterm = par[1] - log( par[2] + pow( 1 / beta_G, par[4] ) );
171 bethe_B = par[0] / pow( beta, par[3] ) * betterm - par[0];
172 }
173 else if ( dedxflag == 2 )
174 {
175 double A = 0, B = 0, C = 0;
176 double x = beta_G;
177 if ( x < 4.5 ) A = 1;
178 else if ( x < 10 ) B = 1;
179 else C = 1;
180 double partA = par[0] * pow( sqrt( x * x + 1 ), par[2] ) / pow( x, par[2] ) *
181 ( par[1] - par[5] * log( pow( 1 / x, par[3] ) ) ) -
182 par[4] + exp( par[6] + par[7] * x );
183 double partB = par[8] * pow( x, 3 ) + par[9] * pow( x, 2 ) + par[10] * x + par[11];
184 double partC = -par[12] * log( par[15] + pow( 1 / x, par[13] ) ) + par[14];
185 bethe_B = 550 * ( A * partA + B * partB + C * partC );
186 // for fermi plateau ( the electron region) we just use 1.0
187 if ( beta_G > 100 ) bethe_B = 550 * 1.0;
188 }
189
190 if ( ifMC )
191 {
192 double A = 0, B = 0, C = 0;
193 double x = beta_G;
194 if ( x < 4.5 ) A = 1;
195 else if ( x < 10 ) B = 1;
196 else C = 1;
197 double partA = par[0] * pow( sqrt( x * x + 1 ), par[2] ) / pow( x, par[2] ) *
198 ( par[1] - par[5] * log( pow( 1 / x, par[3] ) ) ) -
199 par[4] + exp( par[6] + par[7] * x );
200 double partB = par[8] * pow( x, 3 ) + par[9] * pow( x, 2 ) + par[10] * x + par[11];
201 double partC = -par[12] * log( par[15] + pow( 1 / x, par[13] ) ) + par[14];
202 bethe_B = 550 * ( A * partA + B * partB + C * partC );
203 // for fermi plateau ( the electron region) we just use 1.0
204 if ( beta_G > 100 ) bethe_B = 550 * 1.0;
205 }
206
207 if ( Nohit > 0 )
208 {
209 dedx_exp[it] = bethe_B;
210 double sig_the = std::sin( (double)theta );
211 double f_betagamma, g_sinth, h_nhit, i_t0;
212
213 if ( ifMC )
214 {
215
216 if ( sigmaflag >= 6 )
217 {
218
219 double nhit = (double)Nohit;
220 double sigma_bg = 1.0;
221 double ndedx = bethe_B / 550;
222 sigma_bg = sig_par[0] + sig_par[1] * ndedx + sig_par[2] * ndedx * ndedx;
223
224 double cor_nhit = 1.0;
225 if ( nhit < 5 ) nhit = 5;
226 if ( nhit <= 35 )
227 cor_nhit = sig_par[3] * pow( nhit, 4 ) + sig_par[4] * pow( nhit, 3 ) +
228 sig_par[5] * pow( nhit, 2 ) + sig_par[6] * nhit + sig_par[7];
229
230 double cor_sin = 1.0;
231 if ( sig_the > 0.99 ) sig_the = 0.99;
232 cor_sin = sig_par[8] * pow( sig_the, 4 ) + sig_par[9] * pow( sig_the, 3 ) +
233 sig_par[10] * pow( sig_the, 2 ) + sig_par[11] * sig_the + sig_par[12];
234
235 // sigma vs t0
236 double cor_t0 = 1.0;
237
238 // calculate sigma
239 if ( trkalg == 1 ) ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0;
240 else ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0 * sig_par[13];
241 }
242
243 else
244 {
245 double x = beta_G;
246 double nhit = (double)Nohit;
247 double sigma_bg = 1.0;
248
249 if ( x > 0.3 )
250 sigma_bg = sig_par[0] * exp( sig_par[1] * x ) +
251 sig_par[2] * exp( sig_par[3] * pow( x, 0.25 ) ) + sig_par[4];
252 else sigma_bg = sig_par[5] * exp( sig_par[6] * x ) + sig_par[7];
253
254 double cor_nhit = 1.0;
255 if ( nhit < 5 ) nhit = 5;
256 if ( nhit <= 35 )
257 cor_nhit = sig_par[8] * pow( nhit, 4 ) + sig_par[9] * pow( nhit, 3 ) +
258 sig_par[10] * pow( nhit, 2 ) + sig_par[11] * nhit + sig_par[12];
259
260 double cor_sin = 1.0;
261 if ( sig_the > 0.99 ) sig_the = 0.99;
262 cor_sin = sig_par[13] * pow( sig_the, 4 ) + sig_par[14] * pow( sig_the, 3 ) +
263 sig_par[15] * pow( sig_the, 2 ) + sig_par[16] * sig_the + sig_par[17];
264
265 // sigma vs t0
266 double cor_t0 = 1.0;
267
268 // calculate sigma
269 if ( trkalg == 1 ) ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0;
270 else ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0 * sig_par[18];
271 }
272 } // end of MC
273 else
274 {
275 if ( sigmaflag == 1 )
276 {
277 f_betagamma = sig_par[0] * pow( beta_G, sig_par[1] ) + sig_par[2];
278 g_sinth = ( sig_par[3] * sig_the * sig_the + sig_par[4] ) /
279 ( sig_par[3] * sig_par[5] * sig_par[5] + sig_par[4] );
280 h_nhit =
281 ( sig_par[6] * Nohit * Nohit + sig_par[7] * Nohit + sig_par[8] ) /
282 ( sig_par[6] * sig_par[9] * sig_par[9] + sig_par[7] * sig_par[9] + sig_par[8] );
283 if ( sig_par[13] != 0 )
284 i_t0 = ( sig_par[10] * t0 * t0 + sig_par[11] * t0 + sig_par[12] ) /
285 ( sig_par[10] * sig_par[13] * sig_par[13] + sig_par[11] * sig_par[13] +
286 sig_par[12] );
287 else if ( sig_par[13] == 0 ) i_t0 = 1;
288 // cout<<"f_betagamma : "<<f_betagamma<<" g_sinth : "<<g_sinth<<" h_nhit
289 // : "
290 //<<h_nhit<<" i_t0 : "<<i_t0<<endl;
291 ex_sigma[it] = f_betagamma * g_sinth * h_nhit * i_t0;
292 }
293 else if ( sigmaflag == 2 )
294 {
295 double x = beta_G;
296 double nhit = (double)Nohit;
297 double sigma_bg = 1.0;
298 if ( x > 0.3 )
299 sigma_bg = sig_par[0] * exp( sig_par[1] * x ) +
300 sig_par[2] * exp( sig_par[3] * pow( x, 0.25 ) ) + sig_par[4];
301 else sigma_bg = sig_par[5] * exp( sig_par[6] * x ) + sig_par[7];
302
303 double cor_nhit = 1.0;
304 if ( nhit < 5 ) nhit = 5;
305 if ( nhit <= 35 )
306 cor_nhit = sig_par[8] * pow( nhit, 4 ) + sig_par[9] * pow( nhit, 3 ) +
307 sig_par[10] * pow( nhit, 2 ) + sig_par[11] * nhit + sig_par[12];
308
309 double cor_sin = 1.0;
310 if ( sig_the > 0.99 ) sig_the = 0.99;
311 cor_sin = sig_par[13] * pow( sig_the, 4 ) + sig_par[14] * pow( sig_the, 3 ) +
312 sig_par[15] * pow( sig_the, 2 ) + sig_par[16] * sig_the + sig_par[17];
313
314 double cor_t0 = 1;
315 if ( t0 > 1200 ) t0 = 1200;
316 if ( t0 > 800 ) cor_t0 = sig_par[18] * pow( t0, 2 ) + sig_par[19] * t0 + sig_par[20];
317
318 ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0;
319 }
320 else if ( sigmaflag == 3 )
321 {
322 double x = beta_G;
323 double nhit = (double)Nohit;
324 double sigma_bg = 1.0;
325 if ( x > 0.3 )
326 sigma_bg = sig_par[0] * exp( sig_par[1] * x ) +
327 sig_par[2] * exp( sig_par[3] * pow( x, 0.25 ) ) + sig_par[4];
328 else sigma_bg = sig_par[5] * exp( sig_par[6] * x ) + sig_par[7];
329
330 double cor_nhit = 1.0;
331 if ( nhit < 5 ) nhit = 5;
332 if ( nhit <= 35 )
333 cor_nhit = sig_par[8] * pow( nhit, 4 ) + sig_par[9] * pow( nhit, 3 ) +
334 sig_par[10] * pow( nhit, 2 ) + sig_par[11] * nhit + sig_par[12];
335
336 double cor_sin = 1.0;
337 if ( sig_the > 0.99 ) sig_the = 0.99;
338 cor_sin = sig_par[13] * pow( sig_the, 4 ) + sig_par[14] * pow( sig_the, 3 ) +
339 sig_par[15] * pow( sig_the, 2 ) + sig_par[16] * sig_the + sig_par[17];
340
341 ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin;
342 }
343 else if ( sigmaflag == 4 )
344 {
345 double x = beta_G;
346 double nhit = (double)Nohit;
347 double sigma_bg = 1.0;
348 if ( x > 0.3 )
349 sigma_bg = sig_par[0] * exp( sig_par[1] * x ) +
350 sig_par[2] * exp( sig_par[3] * pow( x, 0.25 ) ) + sig_par[4];
351 else sigma_bg = sig_par[5] * exp( sig_par[6] * x ) + sig_par[7];
352
353 double cor_nhit = 1.0;
354 if ( nhit < 5 ) nhit = 5;
355 if ( nhit <= 35 )
356 cor_nhit = sig_par[8] * pow( nhit, 4 ) + sig_par[9] * pow( nhit, 3 ) +
357 sig_par[10] * pow( nhit, 2 ) + sig_par[11] * nhit + sig_par[12];
358
359 double cor_sin = 1.0;
360 if ( sig_the > 0.99 ) sig_the = 0.99;
361 cor_sin = sig_par[13] * pow( sig_the, 4 ) + sig_par[14] * pow( sig_the, 3 ) +
362 sig_par[15] * pow( sig_the, 2 ) + sig_par[16] * sig_the + sig_par[17];
363
364 if ( trkalg == 1 ) ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin;
365 else ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * sig_par[18];
366 }
367 else if ( sigmaflag == 5 )
368 {
369 double x = beta_G;
370 double nhit = (double)Nohit;
371 double sigma_bg = 1.0;
372 if ( x > 0.3 )
373 sigma_bg = sig_par[0] * exp( sig_par[1] * x ) +
374 sig_par[2] * exp( sig_par[3] * pow( x, 0.25 ) ) + sig_par[4];
375 else sigma_bg = sig_par[5] * exp( sig_par[6] * x ) + sig_par[7];
376
377 double cor_nhit = 1.0;
378 if ( nhit < 5 ) nhit = 5;
379 if ( nhit <= 35 )
380 cor_nhit = sig_par[8] * pow( nhit, 4 ) + sig_par[9] * pow( nhit, 3 ) +
381 sig_par[10] * pow( nhit, 2 ) + sig_par[11] * nhit + sig_par[12];
382 double cor_sin = 1.0;
383 if ( sig_the > 0.99 ) sig_the = 0.99;
384 cor_sin = sig_par[13] * pow( sig_the, 4 ) + sig_par[14] * pow( sig_the, 3 ) +
385 sig_par[15] * pow( sig_the, 2 ) + sig_par[16] * sig_the + sig_par[17];
386
387 double cor_t0 = 1;
388 if ( t0 > 1200 ) t0 = 1200;
389 if ( t0 > 800 ) cor_t0 = sig_par[18] * pow( t0, 2 ) + sig_par[19] * t0 + sig_par[20];
390
391 if ( trkalg == 1 ) ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0;
392 else ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0 * sig_par[21];
393 }
394 else if ( sigmaflag == 6 )
395 {
396 double x = bethe_B / 550;
397 double nhit = (double)Nohit;
398 double sigma_bg = 1.0;
399
400 sigma_bg = sig_par[0] + sig_par[1] * x + sig_par[2] * x * x;
401
402 double cor_nhit = 1.0;
403 if ( nhit < 5 ) nhit = 5;
404 if ( nhit <= 35 )
405 cor_nhit = sig_par[3] * pow( nhit, 4 ) + sig_par[4] * pow( nhit, 3 ) +
406 sig_par[5] * pow( nhit, 2 ) + sig_par[6] * nhit + sig_par[7];
407
408 double cor_sin = 1.0;
409 if ( sig_the > 0.99 ) sig_the = 0.99;
410 cor_sin = sig_par[8] * pow( sig_the, 4 ) + sig_par[9] * pow( sig_the, 3 ) +
411 sig_par[10] * pow( sig_the, 2 ) + sig_par[11] * sig_the + sig_par[12];
412
413 if ( trkalg == 1 ) ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin;
414 else ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * sig_par[13];
415 }
416 else if ( sigmaflag == 7 )
417 {
418 // for(int i=0; i<17; i++) cout << i << " " << sig_par[i] << endl;
419 double x = bethe_B / 550;
420 double nhit = (double)Nohit;
421 double sigma_bg = 1.0;
422
423 sigma_bg = sig_par[0] + sig_par[1] * x + sig_par[2] * x * x;
424
425 double cor_nhit = 1.0;
426 if ( nhit < 5 ) nhit = 5;
427 if ( nhit <= 35 )
428 cor_nhit = sig_par[3] * pow( nhit, 4 ) + sig_par[4] * pow( nhit, 3 ) +
429 sig_par[5] * pow( nhit, 2 ) + sig_par[6] * nhit + sig_par[7];
430 double cor_sin = 1.0;
431 if ( sig_the > 0.99 ) sig_the = 0.99;
432 cor_sin = sig_par[8] * pow( sig_the, 4 ) + sig_par[9] * pow( sig_the, 3 ) +
433 sig_par[10] * pow( sig_the, 2 ) + sig_par[11] * sig_the + sig_par[12];
434
435 double cor_t0 = 1;
436 if ( t0 > 1200 ) t0 = 1200;
437 if ( t0 < 800 ) t0 = 800;
438 cor_t0 = sig_par[13] * pow( t0, 2 ) + sig_par[14] * t0 + sig_par[15];
439
440 if ( trkalg == 1 ) ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0;
441 else ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0 * sig_par[16];
442 // cout << "x " << x << " nhit " << nhit << " sigma_bg " << sigma_bg << "
443 // cor_nhit " << cor_nhit << " cor_sin " << cor_sin << " cor_t0 " << cor_t0 << "
444 // trkalg " << trkalg << endl; cout << "it " << it << " mom " << mom << " dedx "
445 // << dedx << " sinth " << sig_the << " t0 " << t0 << endl;
446 }
447 }
448
449 MdcDedxTrk dedxtrk;
450 double dedx_correc = dedx;
451 chi_dedx[it] = ( dedx_correc - dedx_exp[it] ) / ex_sigma[it];
452 chi2 = chi_dedx[it] * chi_dedx[it];
453 ndf = 1;
454 pid_prob[it] = prob_( &chi2, &ndf );
455 // if(it ==0 ) cout<<"runflag: "<<runflag<<" dedx : "<<dedx<<" chi_dedx: "
456 //<<chi_dedx[it] <<" ptrk: "<<mom<<endl;
457 if ( it == -999 )
458 { // here a debug flag
459 std::cout << " mom = " << mom << "exp" << dedx_exp[it] << " sigma " << ex_sigma[it]
460 << " prob " << pid_prob[it] << std::endl;
461 }
462 if ( pid_prob[it] > max_prob )
463 {
464 max_prob = pid_prob[it];
465 Nmax_prob = it;
466 }
467 }
468 else
469 {
470 dedx_exp[it] = 0.0;
471 ex_sigma[it] = 1000.0;
472 pid_prob[it] = 0.0;
473 chi_dedx[it] = 999.0;
474 } // if Nohit > 0
475 }
476 // std::cout<<"MdcDedxCorrection::dedx_pid_exp(blum)!!!"<<std::endl;
477}
478
479/*void
480 MdcDedxCorrection::getCalib(void) const {
481 Db2GeoDeDx obj1;
482 vector<float> VecA(43);
483 DbGeoDeDx obj2;
484 int i,j,m;
485
486
487 std::cout<<"success"<<std::endl;
488 obj1.Get_DbLayerGain(VecA);
489 std::cout<<"success!!"<<std::endl;
490
491 for(j=0;j<43;j++)
492 {
493 std::cout<<"VecA[" << j << "] is: "<<VecA[j]<<std::endl;
494 }
495
496 vector<float> VecB(6860);
497 std::cout<<"success"<<std::endl;
498 obj1.Get_DbWireGain(VecB);
499 for (i=0;i<6860;i++)
500 {
501 cout<<"VecB[" << i << "] is: "<<VecB[i]<<endl;
502 }
503 std::cout<<"success!!"<<std::endl;
504
505
506 vector<float> VecC(129);
507 std::cout<<"success"<<std::endl;
508 obj1.Get_DbDeDxAttach(VecC);
509 for(i=0;i<43;i++)
510 for(j=0;j<3;j++)
511 {
512 m=i*3+j;
513 cout<<"VecC[" << m << "] is "<<VecC[m]<<endl;
514 }
515 cout<<"success!!"<<endl;
516
517 vector<float> VecD(172);
518 std::cout<<"success"<<std::endl;
519 obj1.Get_DbDeDxSatur(VecD);
520 for(i=0;i<43;i++)
521 for(j=0;j<4;j++)
522 {
523 m=i*4+j;
524 cout<<"VecD[" << m << "] is "<<VecD[m]<<endl;
525 }
526 cout<<"success!!"<<endl;
527
528 }*/
float prob_(float *, int *)
EvtComplex exp(const EvtComplex &c)
float prob_(float *, int *)
***************************************************************************************Pseudo Class RRes *****************************************************************************************Parameters and physical constants **Maarten sept ************************************************************************DOUBLE PRECISION xsmu **************************************************************************PARTICLE DATA all others are from PDG *Only resonances with known widths into electron pairs are sept ************************************************************************C Declarations C
Definition RRes.h:29
void dedx_pid_exp_old(int landau, int runflag, float dedx, int Nhit, float mom, float theta, float t0, float lsamp, double dedx_ex[5], double ex_sigma[5], double pid_prob[5], double chi_dedx[5]) const
void dedx_pid_exp(int vflag[3], float dedx, int trkalg, int Nhit, float mom, float theta, float t0, float lsamp, double dedx_ex[5], double ex_sigma[5], double pid_prob[5], double chi_dedx[5], std::vector< double > &DedxCurve_Parameter, std::vector< double > &DedxSigma_Parameter) const
static double HV2_index_nhit
static double HV2_curvep3
static double HV1_curvep4
static double HV2_index_sin
static double HV1_sigmap2
static double HV1_curvep0
static double HV1_curvep1
static double HV2_curvep4
static double HV2_sigmap1
static double HV1_index_sin
static double HV2_sigmap0
static double HV2_curvep1
static double HV1_index_nhit
static double HV2_curvep2
static double HV2_curvep0
static double HV2_sigmap3
static double HV1_curvep3
static double HV1_sigmap1
static double HV1_sigmap0
static double HV1_sigmap3
static double HV2_sigmap2
static double HV1_curvep2
double SpaceChargeCorrec(double, double, int, double)