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