BOSS 7.1.1
BESIII Offline Software System
Loading...
Searching...
No Matches
MdcTunningSvc.cc
Go to the documentation of this file.
1#include <math.h>
2#include <iostream>
3#include <iomanip>
4#include <fstream>
5#include <iostream>
6
7#include "MdcTunningSvc.h"
8#include "GaudiKernel/Kernel.h"
9#include "GaudiKernel/IInterface.h"
10#include "GaudiKernel/StatusCode.h"
11#include "GaudiKernel/SvcFactory.h"
12#include "GaudiKernel/MsgStream.h"
13
14#include "GaudiKernel/IIncidentSvc.h"
15#include "GaudiKernel/Incident.h"
16#include "GaudiKernel/IIncidentListener.h"
17
18#include "GaudiKernel/ISvcLocator.h"
19#include "GaudiKernel/Bootstrap.h"
20
21
22#include "EventModel/EventModel.h"
23#include "EventModel/Event.h"
24#include "EventModel/EventHeader.h"
25#include "GaudiKernel/SmartDataPtr.h"
26
27#include "DatabaseSvc/DatabaseSvc.h"
28
29using namespace std;
30
31DECLARE_COMPONENT(MdcTunningSvc)
32
33MdcTunningSvc::MdcTunningSvc( const string& name, ISvcLocator* svcloc) :
34 base_class (name, svcloc){
35 m_BesMdcRes=0;
36 // declare properties
37 declareProperty("UseDatabase",m_dbFlag = false);
38 declareProperty("UseEndcapTuning",m_EndcapTuning = 1); // 0:no endcap Tuning , 1: using endcap Tuning
39 declareProperty("EffFile", m_effFile = std::string("no path"));
40 declareProperty("ResFile", m_resFile = std::string("no path"));
41 declareProperty("EffFile_endcap", m_effFile_endcap = std::string("no path"));
42 declareProperty("ResFile_endcap", m_resFile_endcap = std::string("no path"));
43 declareProperty("path_mdc", m_path = std::string("no path"));
44 declareProperty("Host" , host = std::string("bes3db2.ihep.ac.cn"));
45 declareProperty("DbName" , dbName = std::string("offlinedb"));
46 declareProperty("UserName" , userName = std::string("guest"));
47 declareProperty("Password" , password = std::string("guestpass"));
48 declareProperty("SerialNo" , serialNo = 0);
49 declareProperty("fromDB", m_fromDB = true);
50 declareProperty("ParBossVer", m_ParBossVer = std::string("unknown"));
51
52 int no[43]={40,44,48,56,64,72,80,80,76,76,88,88,100,100,112,112,128,128,140,140,160,160,160,160,176,176,176,176,208,208,208,208,240,240,240,240,256,256,256,256,288,288,288};
53
54 for(int i=0;i<43;i++){
55 cellNo[i]=no[i];
56 }
57
58}
59
60MdcTunningSvc::~MdcTunningSvc(){
61 if(m_BesMdcRes) delete m_BesMdcRes;
62}
63
64/*StatusCode MdcTunningSvc::queryInterface(const InterfaceID& riid, void** ppvInterface){
65 if( IID_IMdcTunningSvc.versionMatch(riid) ){
66 *ppvInterface = static_cast<IMdcTunningSvc*> (this);
67 } else{
68 return Service::queryInterface(riid, ppvInterface);
69 }
70 return StatusCode::SUCCESS;
71}*/
72
73
74StatusCode MdcTunningSvc::initialize(){
75 MsgStream log(messageService(), name());
76 log << MSG::INFO << "========== MdcTunningSvc::initialize() ==========" << endreq;
77
78 m_ParBossVer=getenv("BES_RELEASE");
79
80 StatusCode sc = Service::initialize();
81 if( sc.isFailure() ) return sc;
82
83 if(m_fromDB)
84 cout << " MdcTunningSvc read from database. " << endl;
85 else if(!m_fromDB)
86 cout << " MdcTunningSvc read from localfile. " << endl;
87
88 if(m_fromDB==true){
89
90 sc = serviceLocator()->service("DatabaseSvc", m_dbsvc, true);
91 if (sc .isFailure() ) {
92 log << MSG::ERROR << " ERROR: unable to find DatabaseSvc " << endreq;
93 return sc;
94 }
95
96 // MYSQL *conn;
97 // char *opt_host_name = "202.122.33.53";
98 // char *opt_user_name = "maqm";
99 // char *opt_password = "maqm_offline";
100 //// unsigned int opt_port_num = 3306;
101 //// char *opt_socket_name = NULL;
102 // char *opt_db_name = "offlinedb";
103 ////conn = mysql_init(NULL);
104 ////unsigned int opt_flags = 0;
105 ////mysql_real_connect(conn, host.c_str(), userName.c_str(), password.c_str(),
106 //// dbName.c_str(), opt_port_num, opt_socket_name, opt_flags);
107 // mysql_real_connect(conn, opt_host_name, opt_user_name, opt_password,
108 // opt_db_name, opt_port_num, opt_socket_name, opt_flags);
109
110
111 IIncidentSvc* incsvc;
112 sc = service("IncidentSvc", incsvc);
113 int priority = 100;
114 if( sc.isSuccess() ){
115 incsvc -> addListener(this, "NewRun", priority);
116 }
117 sc = serviceLocator()->service("EventDataSvc", m_eventSvc, true);
118 if (sc .isFailure() ) {
119 log << MSG::ERROR << " ERROR: unable to find EventDataSvc " << endreq;
120 return sc;
121 }
122 }
123
124 if(m_fromDB!=true){
125 bool initStat = initTuningConst();
126 // if(m_path!=std::string("no path")) setMdcRes(m_path);
127 if(!initStat){
128 cout << "========== MdcTunningSvc::initialize() failure ! ==========" << endl;
129 return StatusCode::FAILURE;
130 }
131 }
132
133 return StatusCode::SUCCESS;
134}
135
136StatusCode MdcTunningSvc::finalize(){
137 MsgStream log(messageService(), name());
138 log << MSG::INFO << "========== MdcTunningSvc::finalize() ==========" << endreq;
139 //// mysql_close(conn);
140
141 return StatusCode::SUCCESS;
142}
143
144void MdcTunningSvc::handle(const Incident& inc){
145 cout << "========== MdcTunningSvc::handle() ==========" << endl;
146
147 MsgStream log( messageService(), name() );
148 log << MSG::DEBUG << "handle: " << inc.type() << endreq;
149
150 if ( inc.type() == "NewRun" ){
151 log << MSG::DEBUG << "NewRun" << endreq;
152 if(m_fromDB==true) {
153 cout << " Start getMdcTuningTableInfo. " << endl;
154 StatusCode sc = getMdcTuningTableInfo();
155 if(sc.isFailure()) {
156 cout << " ERROR: can not get MdcTuning data from the database. " << endl;
157 log << MSG::ERROR << " ERROR: can not get MdcTuning data from the database. " << endreq;
158 exit(1);
159 }
160 }
161 }
162}
163
164bool MdcTunningSvc::initTuningConst(){
165 std::string FilePath = getenv("MDCTUNNINGSVCROOT");
166
167 if(m_effFile==std::string("no path")){
168 if(!m_fromDB) cout << " ERROR: no mdc tuning eff file, please check the input! " << endl;
169 //m_effFile=FilePath+"/dat/mc_eff.dat";
170 //m_effFile = FilePath+"/dat/mc_eff_psipp.dat";
171 return false;
172 }
173 bool setMcEffStat = setMcEff(m_effFile);
174
175 if(m_resFile==std::string("no path")){
176 if(!m_fromDB) cout << " ERROR: no mdc tuning res file, please check the input! " << endl;
177 //m_resFile=FilePath+"/dat/mc_res.dat";
178 //m_resFile = FilePath+"/dat/mc_res_psipp.dat";
179 return false;
180 }
181 //setMcRes2(m_resFile);
182 bool setMcRes3Stat = setMcRes3(m_resFile);
183
184 if(!(setMcRes3Stat&&setMcEffStat))
185 return false;
186
187 return true;
188}
189
190bool MdcTunningSvc::setMcEff(std::string eff_con){
191 int i,j;
192 string line;
193 double lay,bin,expect,hit;
194 //ifstream readMCEff(eff_con.c_str());
195 std::istringstream readMCEff;
196 // (eff_con);
197 if(m_fromDB){
198 readMCEff.str(eff_con);
199 }
200 if(!m_fromDB){
201 ifstream in(eff_con.c_str());
202 //char ch;
203 //string hhh;
204 //while(ift.get(ch)) hhh.putback(ch);
205 //std::cout<<"hhh:"<<hhh<<std::endl;
206 // stringstream strBuf ;
207 istreambuf_iterator<char> iter(in) ;
208 string strCache = string( iter, (istreambuf_iterator<char>()) );
209 readMCEff.str(strCache);
210 //std::cout<<"strCache:"<<strCache<<std::endl;
211 }
212
213
214 ifstream fin(eff_con.c_str());
215 if(!m_fromDB)
216 if(!fin){
217 cout << " ERROR: the mdc tunning eff file " << m_effFile << " not exist, please check the input! " << endl;
218 return false;
219 }
220 // ifstream readMCEff(eff_con);
221 if(!readMCEff.good()){
222 cout << " ERROR: mdc tuning eff file " << m_effFile << " not exist. " << endl;
223 return false;
224 }else{
225 if(!m_fromDB)
226 if(fin)
227 cout << " Open mdc tuning eff file: " << m_effFile << endl;
228 for(i=0;i<43;i++){
229 readMCEff>>lay;
230 getline(readMCEff,line);
231 for(j=0;j<docaNo;j++){
232 readMCEff>>bin>>docaEff[i][j]>>expect>>hit;
233 }
234 readMCEff>>lay;
235 getline(readMCEff,line);
236 for(j=0;j<thetaNo;j++){
237 readMCEff>>bin>>thetaEff[i][j]>>expect>>hit;
238 }
239 readMCEff>>lay;
240 getline(readMCEff,line);
241 for(j=0;j<cellNo[i];j++){
242 readMCEff>>bin>>cellEff[i][j]>>expect>>hit;
243
244 }
245 }
246 for(i=0;i<43;i++){
247 readMCEff>>lay;
248 getline(readMCEff,line);
249 for(j=0;j<docaNo;j++){
250 readMCEff>>bin>>docaEff_2[i][j]>>expect>>hit;
251 }
252 readMCEff>>lay;
253 getline(readMCEff,line);
254 for(j=0;j<thetaNo;j++){
255 readMCEff>>bin>>thetaEff_2[i][j]>>expect>>hit;
256
257 }
258 readMCEff>>lay;
259 getline(readMCEff,line);
260 for(j=0;j<cellNo[i];j++){
261 readMCEff>>bin>>cellEff_2[i][j]>>expect>>hit;
262
263 }
264 }
265 }
266 return true;
267}
268
269bool MdcTunningSvc::setMcRes(){
270
271
272 int i,j;
273 string line;
274 double lay,bin;
275 ifstream readMCRes(m_resFile.c_str());
276 if(!readMCRes.good()){
277 cout << " ERROR: mdc tuning file: " << m_resFile << " not exist. " << endl;
278 return false;
279 }else{
280 cout << " Open mdc tuning file: " << m_resFile << endl;
281 for(i=0;i<43;i++){
282 readMCRes>>lay;
283 getline(readMCRes,line);
284 getline(readMCRes,line);
285 for(j=0;j<docaNo;j++){
286 readMCRes>>bin>>docaRes[i][j][0][0]>>docaRes[i][j][0][1]; //entranceAngle<0
287 }
288 readMCRes>>lay;
289 getline(readMCRes,line);
290 getline(readMCRes,line);
291 for(j=0;j<docaNo;j++){
292 readMCRes>>bin>>docaRes[i][j][1][0]>>docaRes[i][j][1][1]; //entranceAngle>0
293 /*
294 if(i==0||i==42){
295 cout<<"lay "<<i<<" docaNo "<<j<<" <0 mean "<<docaRes[i][j][0][0]<<" sigma "<<docaRes[i][j][0][1]<<" >0 mean "<<docaRes[i][j][1][0]<<" sigma "<<docaRes[i][j][1][1]<<endl;
296 }
297 */
298 }
299 }
300 }
301 return true;
302}
303
304bool MdcTunningSvc::setMcRes2(std::string res_con){
305
306 int i,j;
307 string line;
308 double lay,bin;
309 //ifstream readMCRes(m_resFile.c_str());
310 std::istringstream readMCRes;
311 if(m_fromDB)
312 {
313 readMCRes.str(res_con);
314 }
315 if(!m_fromDB)
316 {
317 ifstream in(res_con.c_str());
318 istreambuf_iterator<char> iter(in) ;
319 string strCache = string( iter, (istreambuf_iterator<char>()) );
320 readMCRes.str(strCache);
321 }
322 if(!readMCRes.good()){
323 cout << " ERROR: mdc tuning file: " << m_resFile << " not exist. " << endl;
324 return false;
325 }else{
326 cout << " MdcTunningSvc::setMcRes2() Open mdc tuning resfile " << m_resFile << endl;
327 for(i=0;i<43;i++){
328 readMCRes>>lay;
329 getline(readMCRes,line);
330 getline(readMCRes,line);
331 for(j=0;j<docaNo;j++){
332 readMCRes>>bin>>docaF[i][j][0]>>docaMean1[i][j][0]>>docaSigma1[i][j][0]>>docaMean2[i][j][0]>>docaSigma2[i][j][0]; //entranceAngle<0
333 }
334 readMCRes>>lay;
335 getline(readMCRes,line);
336 getline(readMCRes,line);
337 for(j=0;j<docaNo;j++){
338 readMCRes>>bin>>docaF[i][j][1]>>docaMean1[i][j][1]>>docaSigma1[i][j][1]>>docaMean2[i][j][1]>>docaSigma2[i][j][1]; //entranceAngle>0
339 }
340 }
341 for(i=0;i<43;i++){
342 readMCRes>>lay;
343 getline(readMCRes,line);
344 getline(readMCRes,line);
345 for(j=0;j<docaNo;j++){
346 readMCRes>>bin>>docaF_2[i][j][0]>>docaMean1_2[i][j][0]>>docaSigma1_2[i][j][0]>>docaMean2_2[i][j][0]>>docaSigma2_2[i][j][0]; //entranceAngle<0
347 }
348 readMCRes>>lay;
349 getline(readMCRes,line);
350 getline(readMCRes,line);
351 for(j=0;j<docaNo;j++){
352 readMCRes>>bin>>docaF_2[i][j][1]>>docaMean1_2[i][j][1]>>docaSigma1_2[i][j][1]>>docaMean2_2[i][j][1]>>docaSigma2_2[i][j][1]; //entranceAngle>0
353 }
354 }
355 }
356 return true;
357}
358
359bool MdcTunningSvc::setMcRes3(std::string res_con){
360
361 int i,j;
362 string line;
363 double lay,bin;
364 //ifstream readMCRes(m_resFile.c_str());
365 std::istringstream readMCRes;
366 if(m_fromDB){
367 readMCRes.str(res_con);
368 }
369 if(!m_fromDB){
370 ifstream in(res_con.c_str());
371 istreambuf_iterator<char> iter(in) ;
372 string strCache = string( iter, (istreambuf_iterator<char>()) );
373 readMCRes.str(strCache);
374 }
375
376 ifstream fin(res_con.c_str());
377
378 if(!m_fromDB)
379 if(!fin){
380 cout << " ERROR: the mdc tunning res file " << m_resFile << " not exist, please check the input!" << endl;
381 return false;
382 }
383 if(!readMCRes.good()){
384 cout << " ERROR: the mdc tuning res file: " << m_resFile << " not exist, please check the input! " << endl;
385 return false;
386 }else{
387 if(!m_fromDB)
388 if(fin)
389 cout << " MdcTunningSvc::setMcRes3() Open mdc tuning resfile: " << m_resFile << endl;
390 for(i=0;i<43;i++){
391 readMCRes>>lay;
392 getline(readMCRes,line);
393 getline(readMCRes,line);
394 for(j=0;j<docaNo;j++){
395 readMCRes>>bin>>docaF[i][j][0]>>docaMean1[i][j][0]>>docaSigma1[i][j][0]>>docaMean2[i][j][0]>>docaSigma2[i][j][0]>>resLargest[i][j][0]>>resSmallest[i][j][0]>>resRatio[i][j][0]; //entranceAngle<0
396 }
397 readMCRes>>lay;
398 getline(readMCRes,line);
399 getline(readMCRes,line);
400 for(j=0;j<docaNo;j++){
401 readMCRes>>bin>>docaF[i][j][1]>>docaMean1[i][j][1]>>docaSigma1[i][j][1]>>docaMean2[i][j][1]>>docaSigma2[i][j][1]>>resLargest[i][j][1]>>resSmallest[i][j][1]>>resRatio[i][j][1]; //entranceAngle>0
402
403 }
404 }
405 for(i=0;i<43;i++){
406 readMCRes>>lay;
407 getline(readMCRes,line);
408 getline(readMCRes,line);
409 for(j=0;j<docaNo;j++){
410 readMCRes>>bin>>docaF_2[i][j][0]>>docaMean1_2[i][j][0]>>docaSigma1_2[i][j][0]>>docaMean2_2[i][j][0]>>docaSigma2_2[i][j][0]>>resLargest_2[i][j][0]>>resSmallest_2[i][j][0]>>resRatio_2[i][j][0]; //entranceAngle<0
411 }
412 readMCRes>>lay;
413 getline(readMCRes,line);
414 getline(readMCRes,line);
415 for(j=0;j<docaNo;j++){
416
417 readMCRes>>bin>>docaF_2[i][j][1]>>docaMean1_2[i][j][1]>>docaSigma1_2[i][j][1]>>docaMean2_2[i][j][1]>>docaSigma2_2[i][j][1]>>resLargest_2[i][j][1]>>resSmallest_2[i][j][1]>>resRatio_2[i][j][1]; //entranceAngle>0
418
419
420 }
421 }
422 }
423
424 return true;
425}
426
427
428
429BesMdcRes* MdcTunningSvc::getMdcRes(){
430 return m_BesMdcRes;
431}
432
433void MdcTunningSvc::setMdcRes(std::string path){
434 if(m_BesMdcRes) delete m_BesMdcRes;
435 m_BesMdcRes = new BesMdcRes(path);
436}
437
438double MdcTunningSvc::NewSig(int layerId,double driftD){
439 int bindD = 0;
440 double mindD = 0.0 ;
441 double maxdD = 9.0 ;
442 int maxbin =8.0 ;
443
444 if((driftD<mindD)||(driftD>maxdD)){
445 bindD = maxbin ;
446 }else {
447 for(int kk = 0; kk < 9; kk++){
448 if((driftD>=(double)kk)&&(driftD<(double)(kk+1))){
449 bindD = kk ; }
450 }
451 }
452
453 double sigma1 = 0 ;
454
455 sigma1 = (m_BesMdcRes -> getD_dD(layerId ,bindD)) ;
456
457 return sigma1;
458}
459
460
461double MdcTunningSvc::DeldriftD(int layerId,double driftD){
462 int bindD = 0;
463 int maxbin =8 ;
464 double mindD = 0.0 ;
465 double maxdD = 9.0 ;
466
467 for(int jj =0;jj<9;jj++){
468 if((driftD<mindD)||(driftD>maxdD)){
469 bindD = maxbin;
470 }else if((driftD>=jj)&&(driftD<(jj+1))){
471 bindD = jj ;
472 }
473 }
474 double y0D = (m_BesMdcRes -> getD_dD(layerId ,bindD)) ;
475 double y1D = (m_BesMdcRes -> getD_dD(layerId, bindD+1)) ;
476 double yD = y0D + (y1D-y0D)*(driftD - bindD); // calculate the data
477 double y0M = (m_BesMdcRes -> getM_dD(layerId ,bindD)) ;
478 double y1M = (m_BesMdcRes -> getM_dD(layerId ,bindD+1));
479 double yM = y0M + (y1M-y0M)*(driftD - bindD); // calculate the mc data
480 double dely = yD - yM ;
481
482
483 /*if((bindD>=1)&&(bindD<=4)){
484 // cout<<"layerId :"<<layerId<<" dritfD :"<<driftD<<" dely*0.618:"<<dely*0.618<<endl;
485 return dely*0.618 ;
486 }else if((bindD>=5)&&(bindD<7)||(bindD ==8)){
487 // cout<<"layerId :"<<layerId<<" dritfD :"<<driftD<<" dely*0.418:"<<dely*0.418<<endl;
488 return dely*0.418 ;
489 }else {
490 // cout<<"layerId :"<<layerId<<" dritfD :"<<driftD<<" dely:"<<dely<<endl;
491 return dely;
492 }*/
493 return dely;
494
495}
496
497double MdcTunningSvc::Delcostta(int layerId,double costta){
498 int binTa = 0;
499 int maxTa = 15;
500 double minCos = -0.8 ;
501 double minCos2 = -0.7 ;
502 double maxCos = 0.8 ;
503
504 for(int ii = 0; ii <16; ii++){
505 if((costta<minCos)||(costta>maxCos)){
506 binTa = maxTa;
507 }else if((costta>=(minCos + ii*0.1))&&(costta<(minCos2 + ii*0.1))){
508 binTa = ii;
509 }
510 }
511
512 double y0D = (m_BesMdcRes -> getD_theta(layerId ,binTa));
513 double y1D = (m_BesMdcRes -> getD_theta(layerId ,binTa+1));
514 double y0M = (m_BesMdcRes -> getM_theta(layerId,binTa));
515 double y1M = (m_BesMdcRes -> getM_theta(layerId,binTa+1));
516
517 double yD[16],yM[16],Del[16];
518 for(int ll =0;ll<16;ll++){
519 yD[ll] = y0D + (y1D - y0D)/0.1*(costta - (minCos + ll*0.1));
520 yM[ll] = y0M + (y1M - y0D)/0.1*(costta - (minCos + ll*0.1));
521 Del[ll] = yD[ll] - yM[ll] ;
522 }
523
524 double delTha = 0;
525
526 if((binTa>=0)&&(binTa<=5)){
527 delTha = Del[binTa] * 0.118 ;
528 }else if((binTa>5)&&(binTa<=10)){
529 delTha = Del[binTa] * 0.518 ;
530 }else if((binTa>10)&&(binTa<=15)){
531 delTha = Del[binTa] * 0.218 ;
532 }
533
534 return delTha ;
535
536}
537
538
539double MdcTunningSvc::GetEff(int layerId,int cellId,double driftD,double cosTheta,int posFlag){
540 driftD=fabs(driftD);
541 if(driftD>12){
542 //std::cout<<"MdcTuningSvc: driftD overflow "<<driftD<<std::endl;
543 driftD=12;
544 }
545 if(posFlag==0)driftD *= -1;
546
547 if(layerId<0 || layerId>42)std::cout<<" MdcTuningSvc:wrong LayerId "<<layerId<<std::endl;
548 if(cellId<0 || cellId>=cellNo[layerId])std::cout<<"MdcTuningSvc:wrong cellId "<<cellId<<std::endl;
549
550 if(fabs(cosTheta)>1){
551 std::cout<<"MdcTuningSvc:wrong coseTheta "<<cosTheta<<std::endl;
552 cosTheta=1;
553 }
554 double eff;
555 int thetaBin=(int)floor((cosTheta+1)*thetaNo/2.);
556 //debug
557 // std::cout<<"cosTheta "<<cosTheta<<" caled "<<(cosTheta+1)*thetaNo/2.<<" floor "<<thetaBin<<std::endl;
558 int docaBin=(int)floor((driftD+12)*docaNo/24.);
559 if(m_EndcapTuning==0)
560 eff = docaEff[layerId][docaBin] * thetaEff[layerId][thetaBin] * cellEff[layerId][cellId];
561 else {
562 if(fabs(cosTheta)<=0.83)
563 eff = docaEff[layerId][docaBin] * thetaEff[layerId][thetaBin] * cellEff[layerId][cellId];
564 else
565 eff = docaEff_2[layerId][docaBin] * thetaEff_2[layerId][thetaBin] * cellEff_2[layerId][cellId];
566
567 }
568 //debug
569 //std::cout<<"tuning svc layer "<<layerId<<"doca"<<docaBin<<" theta "<<thetaBin<<"cell"<<cellId<<" eff "<<eff<<std::endl;
570
571 return eff;
572}
573
574
575
576double MdcTunningSvc::GetRes(int layerId,int cellId,double driftD,double cosTheta,int posFlag,double entranceAngle,double& mean,double& sigma){
577
578 driftD=fabs(driftD);
579 if(driftD>12){
580 //std::cout<<"MdcTuningSvc: driftD overflow "<<driftD<<std::endl;
581 driftD=12;
582 }
583 if(posFlag==0)driftD *= -1;
584
585 if(layerId<0 || layerId>42)std::cout<<" MdcTuningSvc:wrong LayerId "<<layerId<<std::endl;
586 if(cellId<0 || cellId>=cellNo[layerId])std::cout<<"MdcTuningSvc:wrong cellId "<<cellId<<std::endl;
587
588 if(fabs(cosTheta)>1){
589 std::cout<<"MdcTuningSvc:wrong cosTheta "<<cosTheta<<std::endl;
590 cosTheta=1;
591 }
592
593 //// int thetaBin=(int)floor((cosTheta+1)*thetaNo/2.);
594 //debug
595 // std::cout<<"cosTheta "<<cosTheta<<" caled "<<(cosTheta+1)*thetaNo/2.<<" floor "<<thetaBin<<std::endl;
596 int docaBin=(int)floor((driftD+12)*docaNo/24.);
597 if(entranceAngle<0){
598 mean=docaRes[layerId][docaBin][0][0];
599 sigma=docaRes[layerId][docaBin][0][1];
600 }else{
601 mean=docaRes[layerId][docaBin][1][0];
602 sigma=docaRes[layerId][docaBin][1][1];
603 }
604
605 //debug
606 //std::cout<<"tuning svc layer "<<layerId<<" driftD "<<driftD<<" posFlag "<<posFlag<<" doca "<<docaBin<<" theta "<<thetaBin<<" angle "<<entranceAngle<<" mean "<<mean<<" sigma "<<sigma<<std::endl;
607
608 return 1;
609}
610
611double MdcTunningSvc::GetRes2(int layerId,int cellId,double driftD,double cosTheta,int posFlag,double entranceAngle,double& f,double& mean1,double& sigma1,double& mean2,double& sigma2){
612
613 driftD=fabs(driftD);
614 if(driftD>12){
615 //std::cout<<"MdcTuningSvc: driftD overflow "<<driftD<<std::endl;
616 driftD=12;
617 }
618 if(posFlag==0)driftD *= -1;
619
620 if(layerId<0 || layerId>42)std::cout<<" MdcTuningSvc:wrong LayerId "<<layerId<<std::endl;
621 if(cellId<0 || cellId>=cellNo[layerId])std::cout<<"MdcTuningSvc:wrong cellId "<<cellId<<std::endl;
622
623 if(fabs(cosTheta)>1){
624 std::cout<<"MdcTuningSvc:wrong cosTheta "<<cosTheta<<std::endl;
625 cosTheta=1;
626 }
627
628 //// int thetaBin=(int)floor((cosTheta+1)*thetaNo/2.);
629 //debug
630 // std::cout<<"cosTheta "<<cosTheta<<" caled "<<(cosTheta+1)*thetaNo/2.<<" floor "<<thetaBin<<std::endl;
631 int docaBin=(int)floor((driftD+12)*docaNo/24.);
632 if(m_EndcapTuning==0) {
633 if(entranceAngle<0){
634 f=docaF[layerId][docaBin][0];
635 mean1=docaMean1[layerId][docaBin][0];
636 sigma1=docaSigma1[layerId][docaBin][0];
637 mean2=docaMean2[layerId][docaBin][0];
638 sigma2=docaSigma2[layerId][docaBin][0];
639 }else{
640 f=docaF[layerId][docaBin][1];
641 mean1=docaMean1[layerId][docaBin][1];
642 sigma1=docaSigma1[layerId][docaBin][1];
643 mean2=docaMean2[layerId][docaBin][1];
644 sigma2=docaSigma2[layerId][docaBin][1];
645 }
646 }else {
647 if(fabs(cosTheta)<=0.83) {
648 if(entranceAngle<0){
649 f=docaF[layerId][docaBin][0];
650 mean1=docaMean1[layerId][docaBin][0];
651 sigma1=docaSigma1[layerId][docaBin][0];
652 mean2=docaMean2[layerId][docaBin][0];
653 sigma2=docaSigma2[layerId][docaBin][0];
654 }else{
655 f=docaF[layerId][docaBin][1];
656 mean1=docaMean1[layerId][docaBin][1];
657 sigma1=docaSigma1[layerId][docaBin][1];
658 mean2=docaMean2[layerId][docaBin][1];
659 sigma2=docaSigma2[layerId][docaBin][1];
660 }
661 } else {
662 if(entranceAngle<0){
663 f=docaF_2[layerId][docaBin][0];
664 mean1=docaMean1_2[layerId][docaBin][0];
665 sigma1=docaSigma1_2[layerId][docaBin][0];
666 mean2=docaMean2_2[layerId][docaBin][0];
667 sigma2=docaSigma2_2[layerId][docaBin][0];
668 }else{
669 f=docaF_2[layerId][docaBin][1];
670 mean1=docaMean1_2[layerId][docaBin][1];
671 sigma1=docaSigma1_2[layerId][docaBin][1];
672 mean2=docaMean2_2[layerId][docaBin][1];
673 sigma2=docaSigma2_2[layerId][docaBin][1];
674 }
675 }
676 }
677
678 //debug
679 //std::cout<<"tuning svc layer "<<layerId<<" driftD "<<driftD<<" posFlag "<<posFlag<<" doca "<<docaBin<<" theta "<<thetaBin<<" angle "<<entranceAngle<<" f "<<f<<" mean1 "<<mean1<<" sigma1 "<<sigma1<<" mean2 "<<mean2<<" sigma2 "<<sigma2<<std::endl;
680
681 return 1;
682}
683
684double MdcTunningSvc::GetRes3(int layerId,int cellId,double driftD,double cosTheta,int posFlag,double entranceAngle,double& f,double& mean1,double& sigma1,double& mean2,double& sigma2,double& ResLargest,double& ResSmallest,double& ResRatio){
685
686 driftD=fabs(driftD);
687 if(driftD>12){
688 //std::cout<<"MdcTuningSvc: driftD overflow "<<driftD<<std::endl;
689 driftD=12;
690 }
691 if(posFlag==0)driftD *= -1;
692
693 if(layerId<0 || layerId>42)std::cout<<" MdcTuningSvc:wrong LayerId "<<layerId<<std::endl;
694 if(cellId<0 || cellId>=cellNo[layerId])std::cout<<"MdcTuningSvc:wrong cellId "<<cellId<<std::endl;
695
696 if(fabs(cosTheta)>1){
697 std::cout<<"MdcTuningSvc:wrong cosTheta "<<cosTheta<<std::endl;
698 cosTheta=1;
699 }
700
701 ////int thetaBin=(int)floor((cosTheta+1)*thetaNo/2.);
702 //debug
703 // std::cout<<"cosTheta "<<cosTheta<<" caled "<<(cosTheta+1)*thetaNo/2.<<" floor "<<thetaBin<<std::endl;
704 int docaBin=(int)floor((driftD+12)*docaNo/24.);
705 if(m_EndcapTuning==0) {
706 if(entranceAngle<0){
707 f=docaF[layerId][docaBin][0];
708 mean1=docaMean1[layerId][docaBin][0];
709 sigma1=docaSigma1[layerId][docaBin][0];
710 mean2=docaMean2[layerId][docaBin][0];
711 sigma2=docaSigma2[layerId][docaBin][0];
712 ResLargest=resLargest[layerId][docaBin][0];
713 ResSmallest=resSmallest[layerId][docaBin][0];
714 ResRatio=resRatio[layerId][docaBin][0];
715
716 }else{
717 f=docaF[layerId][docaBin][1];
718 mean1=docaMean1[layerId][docaBin][1];
719 sigma1=docaSigma1[layerId][docaBin][1];
720 mean2=docaMean2[layerId][docaBin][1];
721 sigma2=docaSigma2[layerId][docaBin][1];
722 ResLargest=resLargest[layerId][docaBin][1];
723 ResSmallest=resSmallest[layerId][docaBin][1];
724 ResRatio=resRatio[layerId][docaBin][1];
725 }
726 }else {
727 if(fabs(cosTheta)<=0.83) {
728 if(entranceAngle<0){
729 f=docaF[layerId][docaBin][0];
730 mean1=docaMean1[layerId][docaBin][0];
731 sigma1=docaSigma1[layerId][docaBin][0];
732 mean2=docaMean2[layerId][docaBin][0];
733 sigma2=docaSigma2[layerId][docaBin][0];
734 ResLargest=resLargest[layerId][docaBin][0];
735 ResSmallest=resSmallest[layerId][docaBin][0];
736 ResRatio=resRatio[layerId][docaBin][0];
737 }else{
738 f=docaF[layerId][docaBin][1];
739 mean1=docaMean1[layerId][docaBin][1];
740 sigma1=docaSigma1[layerId][docaBin][1];
741 mean2=docaMean2[layerId][docaBin][1];
742 sigma2=docaSigma2[layerId][docaBin][1];
743 ResLargest=resLargest[layerId][docaBin][1];
744 ResSmallest=resSmallest[layerId][docaBin][1];
745 ResRatio=resRatio[layerId][docaBin][1];
746 }
747 } else {
748 if(entranceAngle<0){
749 f=docaF_2[layerId][docaBin][0];
750 mean1=docaMean1_2[layerId][docaBin][0];
751 sigma1=docaSigma1_2[layerId][docaBin][0];
752 mean2=docaMean2_2[layerId][docaBin][0];
753 sigma2=docaSigma2_2[layerId][docaBin][0];
754 ResLargest=resLargest_2[layerId][docaBin][0];
755 ResSmallest=resSmallest_2[layerId][docaBin][0];
756 ResRatio=resRatio_2[layerId][docaBin][0];
757 }else{
758 f=docaF_2[layerId][docaBin][1];
759 mean1=docaMean1_2[layerId][docaBin][1];
760 sigma1=docaSigma1_2[layerId][docaBin][1];
761 mean2=docaMean2_2[layerId][docaBin][1];
762 sigma2=docaSigma2_2[layerId][docaBin][1];
763 ResLargest=resLargest_2[layerId][docaBin][1];
764 ResSmallest=resSmallest_2[layerId][docaBin][1];
765 ResRatio=resRatio_2[layerId][docaBin][1];
766
767 }
768 }
769 }
770
771 //debug
772 //std::cout<<"tuning svc layer "<<layerId<<" driftD "<<driftD<<" posFlag "<<posFlag<<" doca "<<docaBin<<" theta "<<thetaBin<<" angle "<<entranceAngle<<" f "<<f<<" mean1 "<<mean1<<" sigma1 "<<sigma1<<" mean2 "<<mean2<<" sigma2 "<<sigma2<<std::endl;
773 //debug information
774 //std::cout<<"MdcTunningSvc::GetRes3() debug Info."<<endl
775 // <<"layer docaBin thetaBin entranceAngle f mean1 sigma1 mean2 sigma2 largest smallest ratio"<<endl
776 // <<setw(5)<<layerId<<setw(5)<<docaBin<<setw(5)<<thetaBin<<setw(18)<<entranceAngle<<setw(15)<<f<<setw(15)<<mean1<<setw(15)<<sigma1<<setw(15)<<mean2<<setw(15)<<sigma2<<setw(15)<<ResLargest<<setw(15)<<ResSmallest<<ResRatio<<endl;
777
778 return 1;
779}
780
781double MdcTunningSvc::ResvEntr(int layerId,double enterA,int ilr ,double driftD){
782 int bindD =0;
783 int maxbin = 17;
784 int iEntr = 0;
785 int ll = -1;
786 double mindD = -9.;
787 double maxdD = 9.;
788 double sigmaE = 0.0;
789
790 if(enterA < 0){
791 iEntr = 0 ;
792 }else{ iEntr = 1 ;}
793
794 if(ilr == 0 ){
795 driftD = -1.*driftD;
796 }else{ driftD = driftD ;}
797
798 if( (driftD < mindD) || (driftD > maxdD) ){
799 bindD = maxbin ;
800 }else{
801 for(double dd=-9.;dd<9.;dd++){
802 ll++;
803 if( (driftD>= dd ) && (driftD < (dd+1.)) ){
804 bindD = ll ;
805 }
806 }
807 }
808
809 sigmaE = (m_BesMdcRes -> getD_iEntr(layerId,iEntr,bindD) );
810 //cout<<"Svc lay : "<<layerId<<" iEntr : "<<iEntr<<" bindD : "<<bindD<<" sigmaE :"<<sigmaE<<endl;
811 return sigmaE;
812}
813
814double MdcTunningSvc::DelEtr_Sig(int lay,double enterA,int ilr,double driftD){
815 int bindD = 0;
816 int maxbin =17;
817 int iEntr =0;
818 int ll = -1;
819 double mindD = -9.;
820 double maxdD =9.0;
821 //// double delsigE =0;
822
823 if(enterA < 0){
824 iEntr = 0;
825 }else {iEntr = 1;}
826
827 if(ilr == 0 ){
828 driftD = -1.*driftD;
829 }else{ driftD = driftD ;}
830
831 if( (driftD < mindD) || (driftD > maxdD) ){
832 bindD = maxbin;
833 }else {
834 for(double dd =-9.; dd<9.;dd++){
835 ll++;
836 if( (driftD>=dd ) && (driftD < (dd+1.)) ){
837 bindD = ll;
838 dD[bindD] = dd;
839 }
840 }
841 }
842
843 double y0D = (m_BesMdcRes -> getD_iEntr(lay,iEntr,bindD) );
844 double y1D = (m_BesMdcRes -> getD_iEntr(lay,iEntr,bindD+1) );
845 double yD = y0D + (y1D-y0D)*(driftD - dD[bindD]); // calculate the data
846 double y0M = (m_BesMdcRes -> getM_iEntr(lay,iEntr,bindD) );
847 double y1M = (m_BesMdcRes -> getM_iEntr(lay,iEntr,bindD+1));
848 double yM = y0M + (y1M-y0M)*(driftD - dD[bindD]); // calculate the mc data
849 double dely = yD - yM ;
850
851 //cout<<"Svc lay:"<<lay<<" iEntr :"<<iEntr<<" bindD :"<<bindD
852 // <<" dD["<<bindD<<"] : "<<dD[bindD]
853 // <<" y0D : "<<y0D<<" y1D : "<<y1D<<" yD :"<<yD
854 // <<" y0M : "<<y0M<<" y1M : "<<y1M<<" yM :"<<yM<<" dely : "<<dely<<endl;
855
856 return dely;
857
858}
859
860StatusCode MdcTunningSvc::getMdcTuningTableInfo(){
861 MsgStream log(messageService(), name());
862 SmartDataPtr<Event::EventHeader> eventHeader(m_eventSvc,"/Event/EventHeader");
863 int run=eventHeader->runNumber();
864
865 log << MSG::INFO << "MdcTuningSvc::getMdcTuningTableInfo() run =" << run << endreq;
866 if(m_ParBossVer==std::string("unknown"))
867 cout << "MdcTuningSvc::getMdcTuningTableInfo() : ERROR: there is no ParBossVer " << endl;
868 else log << MSG::INFO << "MdcTuningSvc::getMdcTuningTableInfo() : ParBossVer = " << m_ParBossVer << endl;
869
870 run = fabs(run);
871 ////unsigned long *lengths;
872 //// MYSQL_RES *res_set;
873 //// MYSQL_ROW row;
874 int j=0;
875 int cnt=0;
876 DatabaseRecordVector result;
877 for(int i=0;i<1000;i++){
878 char stmt1[200];
879 cnt = i;
880 int run1 = run+i;
881 cout << " ==================== " << endl;
882 //sprintf(stmt1,"select MdcRes,MdcEff,dEdxTuning from MdcTuning where RunFrom <= %d and RunTo >= %d ",run1,run1);
883 sprintf(stmt1,"select MdcRes,MdcEff from MdcTuning where RunFrom <= %d and RunTo >= %d and SftVer = \"%s\"",run1,run1,m_ParBossVer.c_str());
884 cout << " stmt1: " << stmt1 << endl;
885
886 result.clear();
887 int status = m_dbsvc->query("offlinedb",stmt1,result);
888 if(status<0){
889 cout << " ERROR: can not read MdcRes, MdcEff from the MdcTuning table " << endl;
890 log << MSG::ERROR << " ERROR Read MdcRes, MdcEff from the MdcTuning table " << endreq;
891 return StatusCode::FAILURE;
892 }
893
894 if(result.size()>=1){
895
896 //// int aaa = -1;
897 ////aaa = mysql_real_query(conn, stmt1, strlen(stmt1));
898
899 ////std::cout<<" mysql_real_query: "<<aaa<<std::endl;
900 //// if(aaa){
901 // if(mysql_real_query(conn, stmt1, strlen(stmt1)))
902 ////fprintf(stderr, "query error\n");
903 ////return StatusCode::FAILURE;
904 ////}
905 ////res_set = mysql_store_result (conn);
906 ////row = mysql_fetch_row(res_set);
907 ////int row_no = mysql_num_rows(res_set);
908 ////std::cout<<"row_no="<<row_no<<std::endl;
909 ////if (row_no) {}
910 j=1;
911 break;
912 }
913
914 int run2 = run-i;
915 char stmt2[200];
916 //sprintf(stmt2,"select MdcRes,MdcEff.dEdxTuning from MdcTuning where RunFrom <= %d and RunTo >= %d ",run2,run2);
917 sprintf(stmt2,"select MdcRes,MdcEff from MdcTuning where RunFrom <= %d and RunTo >= %d and SftVer = \"%s\"",run2,run2,m_ParBossVer.c_str());
918
919 result.clear();
920 status = m_dbsvc->query("offlinedb",stmt2,result);
921 if(status<0){
922 log << MSG::ERROR << " ERROR Read MdcRes, MdcEff.dEdxTuning from the MdcTuning table " << endreq;
923 return StatusCode::FAILURE;
924 }
925
926
927 if(result.size()>=1){
928 //// mysql_real_query(conn, stmt2, strlen(stmt2));
929 ////res_set = mysql_store_result (conn);
930 ////row = mysql_fetch_row(res_set);
931 ////row_no = mysql_num_rows(res_set);
932 ////std::cout<<"row_no="<<row_no<<std::endl;
933 ////if (row_no) {}
934 j=-1;
935 break;
936 }
937 }
938
939 if(cnt!=0&&cnt!=1000) {
940 log << MSG::INFO << " get MDC tuning data from run " << run + cnt*j << " instead of run " << run<< endreq;
941 //// std::cout<<"get MDC tuning data form run " <<run+cnt*j<<" instread of run"<< run<<std::endl;
942 }
943 cout << " cnt = " << cnt << endl;
944 if(cnt==1000){
945 log << MSG::ERROR << " can not read Data from DB" << endreq;
946 //// mysql_close(conn);
947 return StatusCode::FAILURE;
948 }
949
950 //// mysql_field_seek (res_set, 0);
951
952 //// lengths = mysql_fetch_lengths(res_set);
953 // string fff = std::string(row[0]);
954 // m_BesMdcRes->setMdcRes(fff);
955 ////string ggg = std::string(row[1]);
956 ////string fff = std::string(row[0]);
957
958 //get last row
959 int row = result.size()-1;
960 cout << " row = " << row << endl;
961
962 if(row<0){
963 cout << " ERROR: can not read Data from DB, please check MdcTunningSvc Version. " << endl;
964 return StatusCode::FAILURE;
965 }
966
967 string ggg = result[row]->GetString("MdcEff");
968 string fff = result[row]->GetString("MdcRes");
969
970 log << MSG::DEBUG << " MdcTunning Data: MdcEff: " << ggg << " MdcRes: " << fff << endreq;
971 //cout << " MdcTunning Data: MdcEff: " << ggg << " MdcRes: " << fff << endl;
972
973 bool stRes = setMcRes3(fff);
974 //bool stRes = setMcRes2(fff);
975 bool stEff = setMcEff(ggg);
976
977 if(!stEff) return StatusCode::FAILURE;
978 if(!stRes) return StatusCode::FAILURE;
979 return StatusCode::SUCCESS;
980}
no
const int no
Definition CheckRung/check/dEdxVsCos.cxx:17
sigma
TTree * sigma
Definition CurConst_dEdXCal.cxx:6
sprintf
sprintf(cut,"kal_costheta0_em>-0.93&&kal_costheta0_em<0.93&&kal_pxy0_em>=0.05+%d*0.1&&kal_pxy0_em<0.15+%d*0.1&&NGch>=2", j, j)
f
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")
iter
EvtStreamInputIterator< typename Generator::result_type > iter(Generator gen, int N=0)
Definition EvtStreamInputIterator.hh:100
bin
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per bin
Definition FoamA.h:85
IDatabaseSvc::query
virtual int query(const std::string &dbName, const std::string &sql, DatabaseRecordVector &res)=0
MdcTunningSvc::DelEtr_Sig
double DelEtr_Sig(int lay, double enterA, int ilr, double driftD)
Definition MdcTunningSvc.cc:814
MdcTunningSvc::setMdcRes
void setMdcRes(std::string path)
Definition MdcTunningSvc.cc:433
MdcTunningSvc::setMcRes2
bool setMcRes2(std::string res_con)
Definition MdcTunningSvc.cc:304
MdcTunningSvc::setMcRes3
bool setMcRes3(std::string res_con)
Definition MdcTunningSvc.cc:359
MdcTunningSvc::getMdcTuningTableInfo
StatusCode getMdcTuningTableInfo()
Definition MdcTunningSvc.cc:860
MdcTunningSvc::Delcostta
double Delcostta(int layerId, double costta)
Definition MdcTunningSvc.cc:497
MdcTunningSvc::handle
void handle(const Incident &)
Definition MdcTunningSvc.cc:144
MdcTunningSvc::getMdcRes
BesMdcRes * getMdcRes()
Definition MdcTunningSvc.cc:429
MdcTunningSvc::ResvEntr
double ResvEntr(int layerId, double enterA, int ilr, double driftD)
Definition MdcTunningSvc.cc:781
MdcTunningSvc::setMcEff
bool setMcEff(std::string eff_con)
Definition MdcTunningSvc.cc:190
MdcTunningSvc::NewSig
double NewSig(int layerId, double driftD)
Definition MdcTunningSvc.cc:438
MdcTunningSvc::DeldriftD
double DeldriftD(int layerId, double driftD)
Definition MdcTunningSvc.cc:461
MdcTunningSvc::GetRes
double GetRes(int layerId, int cellId, double driftD, double cosTheta, int posFlag, double entranceAngle, double &mean, double &sigma)
Definition MdcTunningSvc.cc:576
MdcTunningSvc::GetEff
double GetEff(int layerId, int cellId, double driftD, double cosTheta, int posFlag)
Definition MdcTunningSvc.cc:539
MdcTunningSvc::finalize
virtual StatusCode finalize()
Definition MdcTunningSvc.cc:136
MdcTunningSvc::GetRes2
double GetRes2(int layerId, int cellId, double driftD, double cosTheta, int posFlag, double entranceAngle, double &f, double &mean1, double &sigma1, double &mean2, double &sigma2)
Definition MdcTunningSvc.cc:611
MdcTunningSvc::GetRes3
double GetRes3(int layerId, int cellId, double driftD, double cosTheta, int posFlag, double entranceAngle, double &f, double &mean1, double &sigma1, double &mean2, double &sigma2, double &ResLargest, double &ResSmallest, double &ResRatio)
Definition MdcTunningSvc.cc:684
MdcTunningSvc::initialize
virtual StatusCode initialize()
Definition MdcTunningSvc.cc:74