#include <MdcDedxCorrection.h>
Public Member Functions | |
| MdcDedxCorrection () | |
| ~MdcDedxCorrection () | |
| 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 |
| void | getCalib (void) const |
Definition at line 12 of file MdcDedxCorrection.h.
| MdcDedxCorrection::MdcDedxCorrection | ( | ) |
Definition at line 18 of file MdcDedxCorrection.cxx.
00019 { 00020 #ifdef DEBUG 00021 std::cout<<"MdcDedxCorrection constructed!!!"<<std::endl; 00022 #endif 00023 }
| MdcDedxCorrection::~MdcDedxCorrection | ( | ) |
Definition at line 25 of file MdcDedxCorrection.cxx.
00026 { 00027 #ifdef DEBUG 00028 std::cout<<"MdcDedxCorrection destructed!!!"<<std::endl; 00029 #endif 00030 }
| void MdcDedxCorrection::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 |
Definition at line 138 of file MdcDedxCorrection.cxx.
References EvtCyclic3::A, EvtCyclic3::B, EvtCyclic3::C, exp(), prob_(), sin(), and x.
Referenced by MdcDedxTrk::set_dEdx().
00143 { 00144 const int par_cand( 5 ); 00145 const float Charge_Mass[par_cand] = {0.00051100, 0.10566, 0.13957, 0.4937, 0.93827 }; 00146 double beta_G, beta, betterm, bethe_B; 00147 00148 int dedxflag = vflag[0]; 00149 int sigmaflag = vflag[1]; 00150 bool ifMC = false; 00151 if(vflag[2] == 1) ifMC = true; 00152 00153 int Nmax_prob(0); 00154 float max_prob(-0.01); 00155 int ndf; 00156 float chi2; 00157 00158 for( int it = 0; it < par_cand; it++ ) { 00159 beta_G = mom/Charge_Mass[it]; 00160 00161 if(dedxflag == 1){ 00162 beta = beta_G/sqrt(1+(beta_G)*(beta_G)); 00163 betterm = par[1]-log(par[2]+pow(1/beta_G,par[4])); 00164 bethe_B = par[0]/pow(beta,par[3])*betterm-par[0]; 00165 } 00166 else if(dedxflag == 2) { 00167 double A=0, B=0,C=0; 00168 double x = beta_G; 00169 if(x<4.5) 00170 A=1; 00171 else if(x<10) 00172 B=1; 00173 else 00174 C=1; 00175 double partA = par[0]*pow(sqrt(x*x+1),par[2])/pow(x,par[2])*(par[1]-par[5]*log(pow(1/x,par[3])) ) - par[4]+exp(par[6]+par[7]*x); 00176 double partB = par[8]*pow(x,3)+par[9]*pow(x,2)+par[10]*x+par[11]; 00177 double partC = -par[12]*log(par[15]+pow(1/x,par[13]))+par[14]; 00178 bethe_B = 550*(A*partA+B*partB+C*partC); 00179 //for fermi plateau ( the electron region) we just use 1.0 00180 if(beta_G> 100) bethe_B=550*1.0; 00181 } 00182 00183 if (ifMC) { 00184 double A=0, B=0,C=0; 00185 double x = beta_G; 00186 if(x<4.5) 00187 A=1; 00188 else if(x<10) 00189 B=1; 00190 else 00191 C=1; 00192 double partA = par[0]*pow(sqrt(x*x+1),par[2])/pow(x,par[2])*(par[1]-par[5]*log(pow(1/x,par[3])) ) - par[4]+exp(par[6]+par[7]*x); 00193 double partB = par[8]*pow(x,3)+par[9]*pow(x,2)+par[10]*x+par[11]; 00194 double partC = -par[12]*log(par[15]+pow(1/x,par[13]))+par[14]; 00195 bethe_B = 550*(A*partA+B*partB+C*partC); 00196 //for fermi plateau ( the electron region) we just use 1.0 00197 if(beta_G> 100) bethe_B=550*1.0; 00198 } 00199 00200 00201 if (Nohit > 0) { 00202 dedx_exp[it] = bethe_B; 00203 double sig_the = std::sin((double)theta); 00204 double f_betagamma, g_sinth, h_nhit, i_t0; 00205 00206 if (ifMC) { 00207 00208 if(sigmaflag >= 6) { 00209 00210 double nhit = (double)Nohit; 00211 double sigma_bg=1.0; 00212 double ndedx=bethe_B/550; 00213 sigma_bg= sig_par[0]+sig_par[1]*ndedx+sig_par[2]*ndedx*ndedx; 00214 00215 double cor_nhit = 1.0; 00216 if (nhit < 5) nhit = 5; 00217 if (nhit <= 35) 00218 cor_nhit = sig_par[3]*pow(nhit,4)+sig_par[4]*pow(nhit,3)+sig_par[5]*pow(nhit,2) + 00219 sig_par[6]*nhit+sig_par[7]; 00220 00221 double cor_sin= 1.0; 00222 if(sig_the>0.99) sig_the=0.99; 00223 cor_sin = sig_par[8]*pow(sig_the,4)+sig_par[9]*pow(sig_the,3) + 00224 sig_par[10]*pow(sig_the,2)+sig_par[11]*sig_the+sig_par[12]; 00225 00226 //sigma vs t0 00227 double cor_t0 = 1.0; 00228 00229 //calculate sigma 00230 if (trkalg == 1) 00231 ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0; 00232 else 00233 ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0 * sig_par[13]; 00234 00235 } 00236 00237 else{ 00238 double x= beta_G; 00239 double nhit = (double)Nohit; 00240 double sigma_bg=1.0; 00241 00242 if (x > 0.3) 00243 sigma_bg = sig_par[0]*exp(sig_par[1]*x)+sig_par[2]*exp(sig_par[3]*pow(x,0.25))+sig_par[4]; 00244 else 00245 sigma_bg = sig_par[5]*exp(sig_par[6]*x)+ sig_par[7]; 00246 00247 00248 double cor_nhit = 1.0; 00249 if (nhit < 5) nhit = 5; 00250 if (nhit <= 35) 00251 cor_nhit = sig_par[8]*pow(nhit,4)+sig_par[9]*pow(nhit,3)+sig_par[10]*pow(nhit,2) + 00252 sig_par[11]*nhit+sig_par[12]; 00253 00254 double cor_sin= 1.0; 00255 if(sig_the>0.99) sig_the=0.99; 00256 cor_sin = sig_par[13]*pow(sig_the,4)+sig_par[14]*pow(sig_the,3) + 00257 sig_par[15]*pow(sig_the,2)+sig_par[16]*sig_the+sig_par[17]; 00258 00259 //sigma vs t0 00260 double cor_t0 = 1.0; 00261 00262 //calculate sigma 00263 if (trkalg == 1) 00264 ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0; 00265 else 00266 ex_sigma[it] = 550 * sigma_bg * cor_nhit * cor_sin * cor_t0 * sig_par[18]; 00267 00268 } 00269 }//end of MC 00270 else { 00271 if(sigmaflag == 1) { 00272 f_betagamma=sig_par[0]*pow(beta_G,sig_par[1])+sig_par[2]; 00273 g_sinth=(sig_par[3]*sig_the*sig_the+sig_par[4])/(sig_par[3]*sig_par[5]*sig_par[5]+sig_par[4]); 00274 h_nhit=(sig_par[6]*Nohit*Nohit+sig_par[7]*Nohit+sig_par[8]) / 00275 (sig_par[6]*sig_par[9]*sig_par[9]+sig_par[7]*sig_par[9]+sig_par[8]); 00276 if(sig_par[13] != 0) 00277 i_t0 = (sig_par[10]*t0*t0+sig_par[11]*t0+sig_par[12]) / 00278 (sig_par[10]*sig_par[13]*sig_par[13]+sig_par[11]*sig_par[13]+sig_par[12]); 00279 else if(sig_par[13] == 0) 00280 i_t0 =1; 00281 //cout<<"f_betagamma : "<<f_betagamma<<" g_sinth : "<<g_sinth<<" h_nhit : " 00282 //<<h_nhit<<" i_t0 : "<<i_t0<<endl; 00283 ex_sigma[it]= f_betagamma* g_sinth * h_nhit * i_t0; 00284 } 00285 else if(sigmaflag == 2) { 00286 double x = beta_G; 00287 double nhit = (double)Nohit; 00288 double sigma_bg=1.0; 00289 if (x > 0.3) 00290 sigma_bg = sig_par[0]*exp(sig_par[1]*x)+sig_par[2]*exp(sig_par[3]*pow(x,0.25))+sig_par[4]; 00291 else 00292 sigma_bg = sig_par[5]*exp(sig_par[6]*x)+ sig_par[7]; 00293 00294 double cor_nhit=1.0; 00295 if(nhit<5) nhit=5; 00296 if (nhit <= 35) 00297 cor_nhit = sig_par[8]*pow(nhit,4)+sig_par[9]*pow(nhit,3)+sig_par[10]*pow(nhit,2) + sig_par[11]*nhit+sig_par[12]; 00298 00299 double cor_sin= 1.0; 00300 if(sig_the>0.99) sig_the = 0.99; 00301 cor_sin = sig_par[13]*pow(sig_the,4)+sig_par[14]*pow(sig_the,3) + 00302 sig_par[15]*pow(sig_the,2)+sig_par[16]*sig_the+sig_par[17]; 00303 00304 double cor_t0 = 1; 00305 if (t0 > 1200) t0 = 1200; 00306 if (t0 > 800) 00307 cor_t0 = sig_par[18]*pow(t0,2)+sig_par[19]*t0+sig_par[20]; 00308 00309 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin*cor_t0; 00310 } 00311 else if(sigmaflag == 3) { 00312 double x= beta_G; 00313 double nhit = (double)Nohit; 00314 double sigma_bg=1.0; 00315 if (x > 0.3) 00316 sigma_bg = sig_par[0]*exp(sig_par[1]*x)+sig_par[2]*exp(sig_par[3]*pow(x,0.25))+sig_par[4]; 00317 else 00318 sigma_bg = sig_par[5]*exp(sig_par[6]*x)+ sig_par[7]; 00319 00320 double cor_nhit = 1.0; 00321 if (nhit < 5) nhit = 5; 00322 if (nhit <= 35) 00323 cor_nhit = sig_par[8]*pow(nhit,4)+sig_par[9]*pow(nhit,3)+sig_par[10]*pow(nhit,2) + sig_par[11]*nhit+sig_par[12]; 00324 00325 double cor_sin= 1.0; 00326 if(sig_the>0.99) sig_the=0.99; 00327 cor_sin = sig_par[13]*pow(sig_the,4)+sig_par[14]*pow(sig_the,3) + 00328 sig_par[15]*pow(sig_the,2)+sig_par[16]*sig_the+sig_par[17]; 00329 00330 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin; 00331 } 00332 else if(sigmaflag == 4) { 00333 double x= beta_G; 00334 double nhit = (double)Nohit; 00335 double sigma_bg=1.0; 00336 if (x > 0.3) 00337 sigma_bg = sig_par[0]*exp(sig_par[1]*x)+sig_par[2]*exp(sig_par[3]*pow(x,0.25))+sig_par[4]; 00338 else 00339 sigma_bg = sig_par[5]*exp(sig_par[6]*x)+ sig_par[7]; 00340 00341 double cor_nhit = 1.0; 00342 if (nhit < 5) nhit = 5; 00343 if (nhit <= 35) 00344 cor_nhit = sig_par[8]*pow(nhit,4)+sig_par[9]*pow(nhit,3)+sig_par[10]*pow(nhit,2) + sig_par[11]*nhit+sig_par[12]; 00345 00346 double cor_sin= 1.0; 00347 if(sig_the>0.99) sig_the=0.99; 00348 cor_sin = sig_par[13]*pow(sig_the,4)+sig_par[14]*pow(sig_the,3) + 00349 sig_par[15]*pow(sig_the,2)+sig_par[16]*sig_the+sig_par[17]; 00350 00351 if(trkalg==1) 00352 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin; 00353 else 00354 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin*sig_par[18]; 00355 } 00356 else if(sigmaflag == 5) { 00357 double x = beta_G; 00358 double nhit = (double)Nohit; 00359 double sigma_bg=1.0; 00360 if (x > 0.3) 00361 sigma_bg = sig_par[0]*exp(sig_par[1]*x)+sig_par[2]*exp(sig_par[3]*pow(x,0.25))+sig_par[4]; 00362 else 00363 sigma_bg = sig_par[5]*exp(sig_par[6]*x)+ sig_par[7]; 00364 00365 double cor_nhit=1.0; 00366 if(nhit<5) nhit=5; 00367 if (nhit <= 35) 00368 cor_nhit = sig_par[8]*pow(nhit,4)+sig_par[9]*pow(nhit,3)+sig_par[10]*pow(nhit,2) + 00369 sig_par[11]*nhit+sig_par[12]; 00370 double cor_sin= 1.0; 00371 if(sig_the>0.99) sig_the = 0.99; 00372 cor_sin = sig_par[13]*pow(sig_the,4)+sig_par[14]*pow(sig_the,3) + 00373 sig_par[15]*pow(sig_the,2)+sig_par[16]*sig_the+sig_par[17]; 00374 00375 double cor_t0 = 1; 00376 if (t0 > 1200) t0 = 1200; 00377 if (t0 > 800) 00378 cor_t0 = sig_par[18]*pow(t0,2)+sig_par[19]*t0+sig_par[20]; 00379 00380 if(trkalg==1) 00381 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin*cor_t0; 00382 else 00383 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin*cor_t0*sig_par[21]; 00384 } 00385 else if(sigmaflag == 6) { 00386 double x= bethe_B/550; 00387 double nhit = (double)Nohit; 00388 double sigma_bg=1.0; 00389 00390 sigma_bg= sig_par[0]+sig_par[1]*x+sig_par[2]*x*x; 00391 00392 double cor_nhit = 1.0; 00393 if (nhit < 5) nhit = 5; 00394 if (nhit <= 35) 00395 cor_nhit = sig_par[3]*pow(nhit,4)+sig_par[4]*pow(nhit,3)+sig_par[5]*pow(nhit,2) + sig_par[6]*nhit+sig_par[7]; 00396 00397 double cor_sin= 1.0; 00398 if(sig_the>0.99) sig_the=0.99; 00399 cor_sin = sig_par[8]*pow(sig_the,4)+sig_par[9]*pow(sig_the,3) + 00400 sig_par[10]*pow(sig_the,2)+sig_par[11]*sig_the+sig_par[12]; 00401 00402 if(trkalg==1) 00403 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin; 00404 else 00405 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin*sig_par[13]; 00406 } 00407 else if(sigmaflag == 7) { 00408 // for(int i=0; i<17; i++) cout << i << " " << sig_par[i] << endl; 00409 double x= bethe_B/550; 00410 double nhit = (double)Nohit; 00411 double sigma_bg=1.0; 00412 00413 sigma_bg= sig_par[0]+sig_par[1]*x+sig_par[2]*x*x; 00414 00415 00416 double cor_nhit=1.0; 00417 if(nhit<5) nhit=5; 00418 if (nhit <= 35) 00419 cor_nhit = sig_par[3]*pow(nhit,4)+sig_par[4]*pow(nhit,3)+sig_par[5]*pow(nhit,2) + 00420 sig_par[6]*nhit+sig_par[7]; 00421 double cor_sin= 1.0; 00422 if(sig_the>0.99) sig_the = 0.99; 00423 cor_sin = sig_par[8]*pow(sig_the,4)+sig_par[9]*pow(sig_the,3) + 00424 sig_par[10]*pow(sig_the,2)+sig_par[11]*sig_the+sig_par[12]; 00425 00426 double cor_t0 = 1; 00427 if (t0 > 1200) t0 = 1200; 00428 if (t0 < 800) t0 = 800; 00429 cor_t0 = sig_par[13]*pow(t0,2)+sig_par[14]*t0+sig_par[15]; 00430 00431 if(trkalg==1) 00432 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin*cor_t0; 00433 else 00434 ex_sigma[it]=550*sigma_bg*cor_nhit*cor_sin*cor_t0*sig_par[16]; 00435 // cout << "x " << x << " nhit " << nhit << " sigma_bg " << sigma_bg << " cor_nhit " << cor_nhit << " cor_sin " << cor_sin << " cor_t0 " << cor_t0 << " trkalg " << trkalg << endl; 00436 // cout << "it " << it << " mom " << mom << " dedx " << dedx << " sinth " << sig_the << " t0 " << t0 << endl; 00437 } 00438 } 00439 00440 MdcDedxTrk dedxtrk; 00441 double dedx_correc = dedx; 00442 chi_dedx[it] = (dedx_correc - dedx_exp[it])/ex_sigma[it]; 00443 chi2 = chi_dedx[it]*chi_dedx[it]; 00444 ndf=1; 00445 pid_prob[it] = prob_(&chi2,&ndf); 00446 //if(it ==0 ) cout<<"runflag: "<<runflag<<" dedx : "<<dedx<<" chi_dedx: " 00447 //<<chi_dedx[it] <<" ptrk: "<<mom<<endl; 00448 if (it == -999) { // here a debug flag 00449 std::cout << " mom = " << mom <<"exp"<< dedx_exp[it] 00450 << " sigma "<<ex_sigma[it] <<" prob "<<pid_prob[it] 00451 << std::endl; 00452 } 00453 if( pid_prob[it] > max_prob ){ 00454 max_prob = pid_prob[it]; 00455 Nmax_prob = it; 00456 } 00457 } 00458 else{ 00459 dedx_exp[it] = 0.0; 00460 ex_sigma[it] = 1000.0; 00461 pid_prob[it] = 0.0; 00462 chi_dedx[it] = 999.0; 00463 } //if Nohit > 0 00464 } 00465 // std::cout<<"MdcDedxCorrection::dedx_pid_exp(blum)!!!"<<std::endl; 00466 }
| void MdcDedxCorrection::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 |
Definition at line 34 of file MdcDedxCorrection.cxx.
References MdcDedxParam::HV1_curvep0, MdcDedxParam::HV1_curvep1, MdcDedxParam::HV1_curvep2, MdcDedxParam::HV1_curvep3, MdcDedxParam::HV1_curvep4, MdcDedxParam::HV1_index_nhit, MdcDedxParam::HV1_index_sin, MdcDedxParam::HV1_sigmap0, MdcDedxParam::HV1_sigmap1, MdcDedxParam::HV1_sigmap2, MdcDedxParam::HV1_sigmap3, MdcDedxParam::HV2_curvep0, MdcDedxParam::HV2_curvep1, MdcDedxParam::HV2_curvep2, MdcDedxParam::HV2_curvep3, MdcDedxParam::HV2_curvep4, MdcDedxParam::HV2_index_nhit, MdcDedxParam::HV2_index_sin, MdcDedxParam::HV2_sigmap0, MdcDedxParam::HV2_sigmap1, MdcDedxParam::HV2_sigmap2, MdcDedxParam::HV2_sigmap3, prob_(), sin(), and MdcDedxTrk::SpaceChargeCorrec().
Referenced by MdcDedxTrk::set_dEdx().
00038 { 00039 double par[5], sigma_par[4], sigma_index_nhit, sigma_index_sin; 00040 00041 if(runflag==1) { 00042 par[0]= MdcDedxParam::HV1_curvep0; 00043 par[1]= MdcDedxParam::HV1_curvep1; 00044 par[2]= MdcDedxParam::HV1_curvep2; 00045 par[3]= MdcDedxParam::HV1_curvep3; 00046 par[4]= MdcDedxParam::HV1_curvep4; 00047 sigma_par[0] = MdcDedxParam::HV1_sigmap0; 00048 sigma_par[1] = MdcDedxParam::HV1_sigmap1; 00049 sigma_par[2] = MdcDedxParam::HV1_sigmap2; 00050 sigma_par[3] = MdcDedxParam::HV1_sigmap3; 00051 sigma_index_nhit = MdcDedxParam::HV1_index_nhit; 00052 sigma_index_sin = MdcDedxParam::HV1_index_sin; 00053 } 00054 else if(runflag==2) { 00055 par[0]= MdcDedxParam::HV2_curvep0; 00056 par[1]= MdcDedxParam::HV2_curvep1; 00057 par[2]= MdcDedxParam::HV2_curvep2; 00058 par[3]= MdcDedxParam::HV2_curvep3; 00059 par[4]= MdcDedxParam::HV2_curvep4; 00060 sigma_par[0] = MdcDedxParam::HV2_sigmap0; 00061 sigma_par[1] = MdcDedxParam::HV2_sigmap1; 00062 sigma_par[2] = MdcDedxParam::HV2_sigmap2; 00063 sigma_par[3] = MdcDedxParam::HV2_sigmap3; 00064 sigma_index_nhit = MdcDedxParam::HV2_index_nhit; 00065 sigma_index_sin = MdcDedxParam::HV2_index_sin; 00066 } 00067 00068 const int par_cand( 5 ); 00069 const float Charge_Mass[par_cand] = {0.00051100, 0.10566, 0.13957, 0.4937, 0.93827 }; 00070 double beta_G, beta, betterm, bethe_B, sig_param; 00071 00072 int Nmax_prob( 0 ); 00073 float max_prob( -0.01 ); 00074 int ndf; 00075 float chi2; 00076 00077 for( int it = 0; it < par_cand; it++ ) { 00078 beta_G = mom/Charge_Mass[it]; 00079 beta = beta_G/sqrt(1+(beta_G)*(beta_G)); 00080 betterm = par[1]-log(par[2]+pow(1/beta_G,par[4])); 00081 bethe_B = par[0]/pow(beta,par[3])*betterm-par[0]; 00082 00083 if( Nohit >0 ) { 00084 dedx_exp[it] = bethe_B; 00085 double sig_the=std::sin( (double)theta ); 00086 00087 if(runflag <3 && runflag>0){ 00088 if(landau == 0) { 00089 sig_param = 1.6*std::sin( (double) theta )/(lsamp*double(Nohit)); 00090 ex_sigma[it] = 0.05*bethe_B*sqrt( 50.0*sig_param ); 00091 } 00092 else { 00093 //currently use one sigmap0 00094 if(beta_G < 4) { 00095 sig_param=sigma_par[1]+sigma_par[2]*std::pow(beta_G,sigma_par[3]); 00096 } else { 00097 sig_param= sigma_par[0]; 00098 } 00099 //double sig_the=std::sin( (double)theta ); 00100 sig_the=std::pow(sig_the,sigma_index_sin); 00101 double sig_n; 00102 sig_n=35.0/double(Nohit); 00103 sig_n=std::pow(sig_n,sigma_index_nhit); 00104 ex_sigma[it]=sig_param*sig_the*sig_n; 00105 } 00106 } 00107 00108 MdcDedxTrk dedxtrk; 00109 double dedx_correc; 00110 if(runflag == 2 ) dedx_correc = dedxtrk.SpaceChargeCorrec(theta, mom, it, dedx); 00111 else dedx_correc = dedx; 00112 chi_dedx[it] = (dedx_correc - dedx_exp[it])/ex_sigma[it]; 00113 chi2 = chi_dedx[it]*chi_dedx[it]; 00114 ndf = 1; 00115 pid_prob[it] = prob_(&chi2,&ndf); 00116 //if(it ==0 ) cout<<"runflag: "<<runflag<<" dedx : "<<dedx<<" chi_dedx: " 00117 //<<chi_dedx[it] <<" ptrk: "<<mom<<endl; 00118 if( it == -999 ){ // here a debug flag 00119 std::cout << " mom = " << mom <<"exp"<< dedx_exp[it] 00120 << " sigma "<<ex_sigma[it] <<" prob "<<pid_prob[it] 00121 << std::endl; 00122 } 00123 if( pid_prob[it] > max_prob ) { 00124 max_prob = pid_prob[it]; 00125 Nmax_prob = it; 00126 } 00127 } 00128 else { 00129 dedx_exp[it] = 0.0; 00130 ex_sigma[it] = 1000.0; 00131 pid_prob[it] = 0.0; 00132 chi_dedx[it] = 999.0; 00133 } 00134 } 00135 // std::cout<<"MdcDedxCorrection::dedx_pid_exp(blum)!!!"<<std::endl; 00136 }
| void MdcDedxCorrection::getCalib | ( | void | ) | const |
1.4.7