BOSS_7.0.2: MdcDedxRecon Class Reference

MdcDedxRecon::MdcDedxRecon	(	const std::string &	name,
		ISvcLocator *	pSvcLocator
	)

Definition at line 52 of file MdcDedxRecon.cxx.

References calib_flag, doNewGlobal, ex_calib, landau, m_alg, m_rootfile, ntpFlag, ParticleFlag, recalg, deljobs::string, and truncate.

00052                                                                           :  Algorithm(name, pSvcLocator)
00053 {
00054     ex_calib=0;  // get MdcDedxCorrection
00055     calib_flag = 0x7F; //all calib on
00056     landau = 1;  //0: gauss distribution;  1:landau distribution;
00057     ntpFlag = 1;
00058     doNewGlobal = 0;
00059     recalg = 0; //reconstruction method: 0:RecMdcTrack; 1:RecMdcKalTrack; 
00060     //2:use RecMdcTrack when no RecMdcKalTrack
00061     ParticleFlag = -1;  //e:0,  mu:1,  pi:2,  K:3,  p:4,  cosmic:5
00062     m_alg = 1;  //calculation method of dE/dx of one track; 1: truncation and get average;
00063     truncate = 0.70; // truncation rate
00064 
00065     // Declare the properties
00066     declareProperty("CalibFlag",calib_flag);
00067     declareProperty("NTupleFlag",ntpFlag);
00068     declareProperty("RecAlg",recalg);
00069     declareProperty("ParticleType",ParticleFlag);
00070     declareProperty("dEdxMethod",m_alg);
00071     declareProperty("TruncRate",truncate);
00072     declareProperty("RootFile",m_rootfile = std::string("no rootfile")); 
00073 }

MdcDedxRecon::~MdcDedxRecon ( ) [inline]

Definition at line 33 of file MdcDedxRecon.h.

00033 {};

StatusCode MdcDedxRecon::beginRun ( )

Definition at line 237 of file MdcDedxRecon.cxx.

References Bes_Common::DEBUG, calibUtil::ERROR, geosvc, Bes_Common::INFO, and msgSvc().

00238 {
00239     MsgStream log(msgSvc(), name());
00240     log << MSG::DEBUG <<"in MdcDedxRecon::beginrun()!!!"<< endreq;
00241 
00242     StatusCode gesc = service("MdcGeomSvc", geosvc);
00243     if (gesc ==  StatusCode::SUCCESS) 
00244     {
00245         log << MSG::INFO <<"MdcGeomSvc is running"<<endreq;
00246     } 
00247     else 
00248     {
00249         log << MSG::ERROR <<"MdcGeomSvc is failed"<<endreq;
00250         return StatusCode::SUCCESS;
00251     }
00252 
00253     return StatusCode::SUCCESS;
00254 }

void MdcDedxRecon::clearTables ( )

Definition at line 1115 of file MdcDedxRecon.cxx.

01115 {}

StatusCode MdcDedxRecon::execute ( )

Definition at line 286 of file MdcDedxRecon.cxx.

References abs, cos(), Bes_Common::DEBUG, dedxcursvc, dedxhit, dedxhitList, dedxhitrefvec, dedxList, doNewGlobal, electron, calibUtil::ERROR, eventNo, ex_calib, ex_trks, exsvc, Bes_Common::FATAL, geosvc, IDedxCurSvc::getCurve(), IDedxCurSvc::getCurveSize(), DstMdcKalTrack::getFError(), DstMdcKalTrack::getFHelix(), IDedxCurSvc::getSigma(), IDedxCurSvc::getSigmaSize(), genRecEmupikp::i, Iner_DriftDistCut, Bes_Common::INFO, kaltrackrec(), landau, IMdcGeomSvc::Layer(), MdcID::layer(), m_alg, m_charge, m_charge1, m_chidedx, m_chidedxe, m_chidedxk, m_chidedxmu, m_chidedxp, m_chidedxpi, m_costheta, m_costheta1, m_dEdx_hit, m_dEdx_meas, m_doca_ex, m_doca_in, m_driftdist, m_driftT, m_eangle, m_evtNO, m_expect, m_expectid, m_exraw, m_hitNo_wire, m_layer, m_localwid, m_nhit, m_nhit_hit, m_nphlisthit, m_nt1, m_nt2, m_parttype, m_pathL, m_phraw, m_phtm, M_PI, m_poca_x, m_poca_y, m_poca_z, m_prob, m_probpid, m_ptrk, m_ptrk1, m_ptrk_hit, m_recalg, m_rootfile, m_runNO, m_runNO1, m_sigma, m_sigmaid, m_sintheta, m_t0, m_t01, m_tdc_raw, m_trackId, m_trackId1, m_trkalg, m_trkalgs, m_usedhit, m_wire, MaxHistValue, msgSvc(), ntpFlag, Out_DriftDistCut, ParticleFlag, IDedxCorrecSvc::PathL(), phi0, recalg, EventModel::Recon::RecMdcDedxCol, EventModel::Recon::RecMdcDedxHitCol, resdedx, RPhi_PathMaxCut, MdcDedxParam::RUN0, MdcDedxParam::RUN1, MdcDedxParam::RUN2, MdcDedxParam::RUN3, MdcDedxParam::RUN4, RunNO1, RunNO2, MdcDedxTrk::set_phlist(), MdcDedxTrk::set_phlist_hit(), DstMdcDedx::setChi(), RecMdcDedxHit::setDedx(), RecMdcDedx::setDedxEsat(), RecMdcDedx::setDedxExpect(), RecMdcDedx::setDedxHit(), RecMdcDedx::setDedxMoment(), RecMdcDedx::setDedxNoRun(), RecMdcDedxHit::setFlagLR(), RecMdcDedxHit::setMdcHit(), RecMdcDedxHit::setMdcId(), RecMdcDedxHit::setMdcKalHelixSeg(), RecMdcDedx::setMdcKalTrack(), RecMdcDedx::setMdcTrack(), DstMdcDedx::setNormPH(), DstMdcDedx::setNumGoodHits(), DstMdcDedx::setNumTotalHits(), DstMdcDedx::setParticleId(), RecMdcDedxHit::setPathLength(), RecMdcDedx::setPidProb(), DstMdcDedx::setProbPH(), RecMdcDedx::setSigmaDedx(), DstMdcDedx::setStatus(), DstMdcDedx::setTrackId(), RecMdcDedxHit::setTrkId(), DstMdcDedx::setTruncAlg(), RecMdcDedx::setVecDedxHits(), MdcDedxTrk::setVecDedxHits(), sin(), IDedxCorrecSvc::StandardHitCorrec(), IDedxCorrecSvc::StandardTrackCorrec(), switchtomdctrack(), trackNO1, trackNO2, trackNO3, MdcDedxTrk::trk_id(), truncate, Bes_Common::WARNING, MdcID::wire(), and MdcGeoLayer::Wirst().

00287 {
00288     MsgStream log(msgSvc(), name());
00289     log << MSG::INFO << "in execute()" << endreq;
00290 
00291     vector<double> Curve_Para, Sigma_Para;
00292     int vFlag[3];//vFlag[0]:dedx curve version;  vFlag[1]:dedx sigma version;  vFlag[2]: 0:data; 1:MC
00293 
00294     for(int i=0; i< dedxcursvc->getCurveSize(); i++) // get the dedx curve parameters
00295     {
00296         if(i==0)  vFlag[0] = (int) dedxcursvc->getCurve(i);  //first element is dedx curve version
00297         else  Curve_Para.push_back( dedxcursvc->getCurve(i) );  //dedx curve parameters
00298     }
00299     for(int i=0; i< dedxcursvc->getSigmaSize(); i++)
00300     {
00301         if(i==0) vFlag[1] = (int) dedxcursvc->getSigma(i);  //dedx sigma version: 2: psip; 3:jpsi
00302         else  Sigma_Para.push_back( dedxcursvc->getSigma(i) );  //dedx sigma parameters
00303     }
00304 
00305     //check if size of parameters is right
00306     if(vFlag[0] ==1 && Curve_Para.size() != 5)    //version 1: 5 parameters 652a psip data 
00307       cout<<" size of dedx curve paramters for version 1 is not right!"<<endl;
00308     if(vFlag[0] ==2 && Curve_Para.size() != 16)   //version 2: increase from 5 parameters of 652 to 16
00309       cout<<" size of dedx curve paramters for version 2 is not right!"<<endl;
00310 
00311     if(vFlag[1] ==1 && Sigma_Para.size() != 14)   //vesion 1: 14 parameters 652a psip data (old way)
00312       cout<<" size of dedx sigma paramters for version 1 is not right!"<<endl;
00313     if(vFlag[1] ==2 && Sigma_Para.size() != 21)   //version 2: include t0 correction (for psip09 data)
00314       cout<<" size of dedx sigma paramters for version 2 is not right!"<<endl;
00315     if(vFlag[1] ==3 && Sigma_Para.size() != 18)   //version 3: no t0 correction (for jpsi09 data and beyond)
00316       cout<<" size of dedx sigma paramters for version 3 is not right!"<<endl;
00317     if(vFlag[1] ==4 && Sigma_Para.size() != 19)  //version 4: data with mdc track defaulted when kal track not ok(no t0 correction) 
00318       cout<<" size of dedx sigma paramters for version 4 is not right!"<<endl;
00319     if(vFlag[1] ==5 && Sigma_Para.size() != 22)  //version 5: data with mdc track defaulted when kal track not ok(with t0 correction) 
00320       cout<<" size of dedx sigma paramters for version 5 is not right!"<<endl;
00321 
00322 
00323     //---------- register RecMdcDedxCol and RecMdcDedxHitCol to tds---------//
00324     DataObject *aReconEvent;
00325     eventSvc()->findObject("/Event/Recon",aReconEvent);
00326     if(aReconEvent==NULL) 
00327     {        
00328         aReconEvent = new ReconEvent();
00329         StatusCode sc = eventSvc()->registerObject("/Event/Recon",aReconEvent);
00330         if(sc!=StatusCode::SUCCESS) 
00331         {
00332             log << MSG::FATAL << "Could not register ReconEvent" <<endreq;            
00333             return( StatusCode::FAILURE);
00334         } 
00335     }
00336 
00337     SmartIF<IDataManagerSvc> dataManSvc(eventSvc());
00338 
00339     DataObject *aDedxcol;
00340     eventSvc()->findObject("/Event/Recon/RecMdcDedxCol",aDedxcol);    
00341     if(aDedxcol != NULL) 
00342     {
00343         dataManSvc->clearSubTree("/Event/Recon/RecMdcDedxCol");
00344         eventSvc()->unregisterObject("/Event/Recon/RecMdcDedxCol");
00345     }
00346     dedxList = new RecMdcDedxCol;
00347     StatusCode dedxsc;
00348     dedxsc =  eventSvc()->registerObject(EventModel::Recon::RecMdcDedxCol, dedxList);
00349     if(!dedxsc.isSuccess()) 
00350     {
00351         log << MSG::FATAL << " Could not register Mdc dedx collection" <<endreq;
00352         return ( StatusCode::FAILURE);
00353     }
00354 
00355     DataObject *aDedxhitcol;
00356     eventSvc()->findObject("/Event/Recon/RecMdcDedxHitCol",aDedxhitcol);
00357     if(aDedxhitcol != NULL) 
00358     {
00359         dataManSvc->clearSubTree("/Event/Recon/RecMdcDedxHitCol");
00360         eventSvc()->unregisterObject("/Event/Recon/RecMdcDedxHitCol");
00361     }
00362     dedxhitList = new RecMdcDedxHitCol;
00363     StatusCode dedxhitsc;
00364     dedxhitsc =  eventSvc()->registerObject(EventModel::Recon::RecMdcDedxHitCol, dedxhitList);
00365     if(!dedxhitsc.isSuccess()) 
00366     {
00367         log << MSG::FATAL << " Could not register Mdc dedx collection" <<endreq;
00368         return ( StatusCode::FAILURE);
00369     }    
00370 
00371     SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
00372     if (!eventHeader) 
00373     {
00374         log << MSG::INFO << "Could not find Event Header" << endreq;
00375         return( StatusCode::FAILURE);
00376     }
00377     int eventNO = eventHeader->eventNumber();
00378     int runNO = eventHeader->runNumber();
00379     log << MSG::INFO << "now begin the event:  runNO= "<<runNO<<" eventNO= "<< eventNO<< endreq;
00380     RunNO2= runNO; 
00381     if(RunNO1==0) RunNO1=runNO;
00382     if(RunNO2 != RunNO1)
00383     {
00384         cout<<"RunNO2 = "<<RunNO2 <<"      RunNO1 = "<<RunNO1<<endl;
00385     }
00386     RunNO1 = runNO;
00387     int runFlag; //data type flag
00388     if(runNO<MdcDedxParam::RUN0) runFlag =0;   //MC: less than 0
00389     else if(runNO>=MdcDedxParam::RUN1 && runNO<MdcDedxParam::RUN2) runFlag =1;
00390     else if(runNO>=MdcDedxParam::RUN2 && runNO<MdcDedxParam::RUN3) runFlag =2;
00391     else if(runNO>=MdcDedxParam::RUN3 && runNO<MdcDedxParam::RUN4) runFlag =3;
00392     else  runFlag =4;
00393 
00394     //vFlag[2] identify MC or data: 1:Mc;  0:data
00395     if(runNO < 0)   vFlag[2]=1;
00396     else            vFlag[2]=0;
00397 
00398     double tes = -999.0;
00399     int esTimeflag = -1;
00400     SmartDataPtr<RecEsTimeCol> estimeCol(eventSvc(),"/Event/Recon/RecEsTimeCol");
00401     if( ! estimeCol) 
00402     {
00403         log << MSG::WARNING << "Could not find EvTimeCol" << endreq;
00404     }
00405     else
00406     {
00407         RecEsTimeCol::iterator iter_evt = estimeCol->begin();
00408         for(; iter_evt!=estimeCol->end(); iter_evt++)
00409         {
00410             tes =  (*iter_evt)->getTest();
00411             esTimeflag = (*iter_evt)->getStat();
00412         }
00413     }   
00414     //cout<<"estime : "<<tes<<endl;
00415 
00416 
00417     Identifier mdcid; 
00418     int ntrk;
00419     ex_trks.clear();  // clear the vector of MdcDedxTrk,when read a new event  
00420     m_trkalgs.clear();
00421     if( !doNewGlobal )
00422     {
00423         log << MSG::DEBUG << "recalg:  "<<recalg<<endreq;
00424 
00425         //---------using kal algorithm by default, switch to seek mdc track if no kaltrack hits //
00426         if(recalg ==2 )
00427         {
00428             //retrieve  MdcKalTrackCol from TDS
00429             SmartDataPtr<RecMdcKalTrackCol> newtrkCol(eventSvc(),"/Event/Recon/RecMdcKalTrackCol");
00430             if (!newtrkCol) 
00431             {
00432                 log << MSG::WARNING << "Could not find RecMdcKalTrackCol" << endreq;
00433                 return StatusCode::SUCCESS;
00434             }
00435             if(ntpFlag>0) eventNo++;
00436             ntrk = newtrkCol->size();
00437             log << MSG::DEBUG << "track collection size:  " << newtrkCol->size() <<endreq;
00438 
00439             RecMdcKalTrackCol::iterator trk = newtrkCol->begin();
00440             //cout<<"MdcDedxRecon  newtrkCol.size()  "<<newtrkCol->size()<<endl;
00441             for( ; trk != newtrkCol->end(); trk++) 
00442             {
00443                 if(ntpFlag>0) trackNO1++;
00444 
00445                 HelixSegRefVec gothits= (*trk)->getVecHelixSegs(ParticleFlag);  
00446                 //if not set ParticleFlag, we will get the last succefully hypothesis;
00447                 if(gothits.size()==0)
00448                 {
00449                     m_trkalg=0;
00450                     int trkid=(*trk)->trackId();
00451                     switchtomdctrack(trkid, mdcid, tes, runNO, runFlag, log);
00452                 }
00453                 else
00454                 {
00455                     m_trkalg=1;
00456                     if(gothits.size()<2) continue;
00457                     kaltrackrec(trk, mdcid, tes, runNO, runFlag, log );
00458                 }
00459             }//end of track loop
00460         }//end of recalg==2
00461 
00462         //------------------------ Use information of MDC track rec --------------------------//
00463         else if(recalg ==0 )
00464         {
00465             m_trkalg=0; 
00466 
00467             //retrieve  MdcTrackCol from TDS     
00468             SmartDataPtr<RecMdcTrackCol> newtrkCol(eventSvc(),"/Event/Recon/RecMdcTrackCol");
00469             if (!newtrkCol) 
00470             {
00471                 log << MSG::WARNING << "Could not find RecMdcTrackCol" << endreq;
00472                 return StatusCode::SUCCESS;
00473             }
00474             if(ntpFlag>0) eventNo++;
00475             ntrk = newtrkCol->size();
00476             log << MSG::DEBUG << "track collection size:  " << newtrkCol->size() <<endreq;
00477 
00478             vector<double> phlist;  //dE/dx only after hit level correction
00479             vector<double> phlist_hit;  //dE/dx after hit and track level correction
00480             double m_d0E=0,m_phi0E=0,m_cpaE=0,m_z0E=0,phi0=0,costheta=0,sintheta=0,m_Pt=0,m_P=0;
00481             int charge=0,layid=0,localwid=0,w0id=0,wid=0,lr=-2;
00482             double adc_raw=0,tdc_raw=0,zhit=0,driftd=0,driftD=0,driftT=0,dd=0,eangle=0;
00483             int Nhits=0;
00484             double ph=0,ph_hit=0,pathlength=0,rphi_path=0;
00485 
00486             RecMdcTrackCol::iterator trk = newtrkCol->begin(); 
00487             for( ; trk != newtrkCol->end(); trk++)
00488             {
00489                 if(ntpFlag>0) trackNO1++;
00490 
00491                 MdcDedxTrk trk_ex( **trk);
00492                 log << MSG::DEBUG <<" %%%%% trk id = "<<trk_ex.trk_id() <<endreq;
00493                 log << MSG::DEBUG <<" %%%%% initial id = "<<(*trk)->trackId() <<endreq;
00494 
00495                 HepVector a = (*trk)->helix();
00496                 HepSymMatrix tkErrM = (*trk)->err();
00497                 m_d0E   = tkErrM[0][0];
00498                 m_phi0E = tkErrM[1][1];
00499                 m_cpaE  = tkErrM[2][2];
00500                 m_z0E   = tkErrM[3][3];
00501 
00502                 HepPoint3D IP(0,0,0);
00503                 Dedx_Helix exhel(IP, a);  //initialize class Dedx_Helix for DedxCorrecSvc
00504                 log << MSG::DEBUG <<"MDC Helix 5 parameters: "<<a[0]<<"     "<<a[1]<<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endreq;
00505                 //cout<<"MDC Helix 5 parameters: "<<a[0]<<"     "<<a[1]<<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endl;
00506                 phi0= a[1];
00507                 costheta = cos(M_PI_2-atan(a[4]));
00508                 //cout<<"track theta: "<<M_PI_2-atan(a[4])<<" extrack costheta: "<<trk_ex.theta()<<endl;
00509                 sintheta = sin(M_PI_2-atan(a[4]));
00510                 m_Pt = 1.0/fabs( a[2] );
00511                 m_P = m_Pt*sqrt(1 + a[4]*a[4]);
00512                 charge = ( a[2] > 0 )? 1 : -1;
00513                 //cout<<"track charge: "<<charge<<" extrack charge: "<<trk_ex.charge()<<endl;
00514                 dedxhitrefvec = new DedxHitRefVec;
00515 
00516 
00517                 HitRefVec gothits= (*trk)->getVecHits();
00518                 Nhits = gothits.size();
00519                 log << MSG::DEBUG << "hits size  = " <<gothits.size()<<endreq;
00520                 if(gothits.size()<2) 
00521                 {
00522                     delete dedxhitrefvec;
00523                     continue;
00524                 }
00525                 //if(gothits.size()<15) cout<<"eventNO : "<<eventNO<<"  hits : "<<gothits.size()<<endl;
00526 
00527                 RecMdcKalHelixSeg* mdcKalHelixSeg = 0;
00528                 HitRefVec::iterator it_gothit = gothits.begin();
00529                 for( ; it_gothit != gothits.end(); it_gothit++) 
00530                 {
00531                     mdcid = (*it_gothit)->getMdcId();
00532                     layid =  MdcID::layer(mdcid);
00533                     localwid = MdcID::wire(mdcid); 
00534                     w0id = geosvc->Layer(layid)->Wirst();
00535                     wid = w0id + localwid;
00536                     adc_raw = (*it_gothit)->getAdc();
00537                     tdc_raw = (*it_gothit)->getTdc();
00538                     log << MSG::INFO <<"hit layer id  " <<layid<< " wire id: " <<localwid<< " adc_raw: " <<adc_raw<<" tdc_raw: "<<tdc_raw<<endreq;
00539                     zhit = (*it_gothit)->getZhit();
00540                     lr = (*it_gothit)->getFlagLR();
00541                     if(lr == 2) cout<<"lr = "<<lr<<endl;
00542                     if(lr == 0 || lr == 2) driftD = (*it_gothit)->getDriftDistLeft();
00543                     else  driftD = (*it_gothit)->getDriftDistRight();   
00544                     //cout<<"lr: "<<lr<<"  driftD: "<<driftD<<endl;
00545                     driftd = abs(driftD);
00546                     dd = (*it_gothit)->getDoca();
00547                     if(lr == 0 || lr == 2 ) dd = -abs(dd);
00548                     if(lr == 1 ) dd = abs(dd);
00549                     driftT = (*it_gothit)->getDriftT();
00550 
00551                     eangle = (*it_gothit)->getEntra();
00552                     eangle = eangle/M_PI;
00553                     pathlength=exsvc->PathL( ntpFlag, exhel, layid, localwid, zhit);
00554                     rphi_path=pathlength * sintheta;
00555 
00556                     //cout<<"ph para check: "<<"adc_raw: "<<adc_raw<<"  dd: "<<dd<<" eangle: "<<eangle<<" zhit: "<<zhit<<" costheta: "<<costheta<<endl;
00557                     ph=exsvc->StandardHitCorrec(0,runFlag,ntpFlag,runNO,pathlength,wid,layid,adc_raw,dd,eangle,costheta);
00558                     ph_hit = exsvc->StandardTrackCorrec(0,runFlag, ntpFlag, runNO, ph, costheta,tes );
00559                     //if(pathlength == -1)
00560                     //cout<<"parameter0: "<<"eventNO: "<<eventNO<<" layid: "<<layid<<" localwid: "<<localwid<<" driftd: "<<driftd<<" rphi_path: "<<rphi_path<<" pathlength: "<<pathlength<<" ph: "<<ph<<" ph_hit: "<<ph_hit<<endl;
00561 
00562                     if(runNO<=5911 && runNO>=5854 && layid ==14) continue;
00563                     if(runNO<0)
00564                     {
00565                         if (layid<8)
00566                         {  
00567                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;
00568                         }
00569                         else if(layid >= 8)
00570                         {
00571                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue ||rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue;                     
00572                         }
00573                     }
00574                     else if(runNO>=0)
00575                     {
00576                         if(layid <8)
00577                         {  
00578                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;
00579                         }
00580                         else if(layid >= 8)
00581                         {
00582                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue;
00583                         }
00584                     }
00585                     //cout<<"recAlg=0 para mdc:  m_phraw: "<<adc_raw<<" m_exraw: "<<ph_hit<<" m_doca: "<<dd<<" m_eangle: "<<eangle<<" m_costheta: "<<costheta<<endl;
00586 
00587                     if (ph<0.0||ph==0) continue;
00588                     else
00589                     {
00590                         //-----------------------put data into TDS-----------------------------//
00591                         dedxhit = new RecMdcDedxHit;
00592                         dedxhit->setMdcHit(*it_gothit);
00593                         dedxhit->setMdcKalHelixSeg(mdcKalHelixSeg );
00594                         dedxhit->setMdcId(mdcid);
00595                         dedxhit->setFlagLR(lr);
00596                         dedxhit->setTrkId(trk_ex.trk_id());
00597                         dedxhit->setDedx(ph_hit);
00598                         dedxhit->setPathLength(pathlength);
00599 
00600                         //---------------------------Fill root file----------------------------//
00601                         if(m_rootfile!= "no rootfile")
00602                         {
00603                             m_hitNo_wire->Fill(wid);
00604                         }
00605 
00606                         //-------------------------Fill ntuple n102---------------------------//
00607                         if ( ntpFlag ==2 )
00608                         {
00609                             m_charge1 = (*trk)->charge();
00610                             m_t01 = tes;
00611                             m_driftT = driftT;
00612                             m_tdc_raw = tdc_raw;
00613                             m_phraw = adc_raw;
00614                             m_exraw = ph_hit;
00615                             m_localwid = localwid;
00616                             m_wire = wid;
00617                             m_runNO1 = runNO;
00618                             m_nhit_hit = Nhits;
00619                             m_doca_in = dd;
00620                             m_doca_ex = dd;
00621                             m_driftdist = driftD;
00622                             m_eangle = eangle;
00623                             m_costheta1 = costheta;
00624                             m_pathL = pathlength;
00625                             m_layer = layid;
00626                             m_ptrk1 = m_P;
00627                             //cout<<"adc_raw : "<<adc_raw<<"     "<<ph_hit<<"     "<<dd_in<<"  "<<dd_ex<<"  "<<eangle<<endl;
00628                             m_trackId1 = trk_ex.trk_id();
00629                             StatusCode status2 = m_nt2->write();
00630                             if ( status2.isFailure() )
00631                               log << MSG::ERROR << "Cannot fill Ntuple n102" << endreq;
00632                         }
00633                         if(layid>3)   
00634                         {
00635                             phlist.push_back(ph);
00636                             phlist_hit.push_back(ph_hit);
00637                         }
00638                     }
00639                     dedxhitList->push_back( dedxhit );
00640                     dedxhitrefvec->push_back( dedxhit );
00641                 }//end of hit loop
00642                 trk_ex.set_phlist( phlist );
00643                 trk_ex.set_phlist_hit( phlist_hit );
00644                 trk_ex.setVecDedxHits( *dedxhitrefvec );
00645                 ex_trks.push_back(trk_ex );
00646                 m_trkalgs.push_back(m_trkalg);
00647 
00648                 delete dedxhitrefvec;
00649                 phlist.clear();
00650                 phlist_hit.clear();
00651                 if(ntpFlag>0) trackNO2++;
00652             }//end of track loop
00653             log << MSG::DEBUG << "size in processing:  " << ex_trks.size() <<endreq; 
00654         }//end of recalg==0
00655 
00656         //------------------------ Use information of MDC KAL track rec -----------------------//
00657         else if(recalg ==1 )
00658         {
00659             m_trkalg=1;
00660 
00661             //retrieve  MdcKalTrackCol from TDS
00662             SmartDataPtr<RecMdcKalTrackCol> newtrkCol(eventSvc(),"/Event/Recon/RecMdcKalTrackCol");
00663             if (!newtrkCol) 
00664             {
00665                 log << MSG::WARNING << "Could not find RecMdcKalTrackCol" << endreq;
00666                 return  StatusCode::SUCCESS;
00667             }
00668             if(ntpFlag>0) eventNo++; 
00669             ntrk = newtrkCol->size();
00670             log << MSG::DEBUG << "track collection size:  " << newtrkCol->size() <<endreq;
00671 
00672             vector<double> phlist; //dE/dx only after hit level correction
00673             vector<double> phlist_hit;  //dE/dx after hit and track level correction
00674             double m_d0E=0,m_phi0E=0,m_cpaE=0,m_z0E=0,phi0=0,costheta=0,sintheta=0,m_Pt=0,m_P=0;
00675             int charge=0,layid=0,localwid=0,w0id=0,wid=0,lr=-2;
00676             double p_hit=0,adc_raw=0,tdc_raw=0,zhit=0,driftd=0,driftD=0,driftT=0,dd_in=0,dd_ex=0,eangle=0;
00677             int Nhits=0;
00678             double ph=0,ph_hit=0,pathlength=0,rphi_path=0;
00679 
00680 
00681             RecMdcKalTrackCol::iterator trk = newtrkCol->begin();
00682             for( ; trk != newtrkCol->end(); trk++)
00683             {
00684                 if(ntpFlag>0) trackNO1++;
00685 
00686                 MdcDedxTrk trk_ex( **trk, ParticleFlag);
00687                 log << MSG::DEBUG <<" %%%%% trk id = "<<trk_ex.trk_id() <<endreq;
00688                 log << MSG::DEBUG <<" %%%%% initial id = "<<(*trk)->trackId() <<endreq;
00689 
00690                 HepVector a;
00691                 HepSymMatrix tkErrM;
00692                 if(ParticleFlag>-1&&ParticleFlag<5) 
00693                 {
00694                     DstMdcKalTrack* dstTrk = *trk;
00695                     a = dstTrk->getFHelix(ParticleFlag);
00696                     tkErrM = dstTrk->getFError(ParticleFlag);
00697                 }
00698                 else 
00699                 {
00700                     a = (*trk)->getFHelix();
00701                     tkErrM = (*trk)->getFError();
00702                 }
00703                 //cout<<"FHelix 5 parameters: "<<a[0]<<"     "<<a[1] <<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endl;  //getFHelix is first layer of MdcKalTrack;
00704                 log << MSG::DEBUG <<"FHelix 5 parameters: "<<a[0]<<"     "<<a[1] <<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endreq;
00705 
00706                 m_d0E   = tkErrM[0][0];
00707                 m_phi0E = tkErrM[1][1];
00708                 m_cpaE  = tkErrM[2][2];
00709                 m_z0E   = tkErrM[3][3];
00710 
00711                 phi0= a[1];
00712                 costheta = cos(M_PI_2-atan(a[4]));
00713                 //cout<<"track theta: "<<M_PI_2-atan(a[4])<<" extrack theta: "<<trk_ex.theta()<<endl;  //theta in trk_ex is got by getFHelix();
00714                 sintheta = sin(M_PI_2-atan(a[4]));
00715                 m_Pt = 1.0/fabs( a[2] );
00716                 m_P = m_Pt*sqrt(1 + a[4]*a[4]);
00717                 //cout<<"track ptrk: "<<m_P<<" extrack ptrk: "<<trk_ex.P()<<endl;
00718                 charge = ( a[2] > 0 )? 1 : -1;
00719                 //cout<<"track charge: "<<charge<<" extrack charge: "<<trk_ex.charge()<<endl;
00720                 dedxhitrefvec = new DedxHitRefVec;
00721 
00722 
00723                 HelixSegRefVec gothits= (*trk)->getVecHelixSegs(ParticleFlag);
00724                 //HelixSegRefVec gothits= (*trk)->getVecHelixSegs();
00725                 //if not set ParticleFlag, we will get the last succefully hypothesis;
00726                 Nhits = gothits.size();
00727                 log << MSG::DEBUG << "hits size  = " <<gothits.size()<<endreq;
00728                 if(gothits.size()<2) 
00729                 {
00730                     delete dedxhitrefvec;
00731                     continue;
00732                 }
00733                 //if(gothits.size()<15) cout<<"eventNO : "<<eventNO<<"  hits : "<<gothits.size()<<endl;
00734 
00735                 RecMdcHit* mdcHit = 0;
00736                 HelixSegRefVec::iterator it_gothit = gothits.begin();
00737                 for( ; it_gothit != gothits.end(); it_gothit++)
00738                 {
00739                     HepVector a_hit_in = (*it_gothit)->getHelixIncl();
00740                     //HepVector a_hit_ex = (*it_gothit)->getHelixExcl(); //same with getHelixIncl()
00741                     HepPoint3D IP(0,0,0);
00742                     Dedx_Helix exhel_hit_in(IP, a_hit_in);
00743                     //Dedx_Helix exhel_hit_ex(IP, a_hit_ex);
00744                     p_hit = 1.0/fabs(a_hit_in(3))*sqrt(1+a_hit_in(5)*a_hit_in(5));   
00745                     //cout<<"getHelixIncl 5 para: "<<a_hit_in[0]<<"     "<<a_hit_in[1] <<"      "<<a_hit_in[2]<<"      "<<a_hit_in[3]<<"     "<<a_hit_in[4]<<endl;
00746 
00747                     mdcid = (*it_gothit)->getMdcId();
00748                     layid =  MdcID::layer(mdcid);
00749                     localwid = MdcID::wire(mdcid);
00750                     w0id = geosvc->Layer(layid)->Wirst();
00751                     wid = w0id + localwid;
00752                     adc_raw = (*it_gothit)->getAdc();
00753                     tdc_raw = (*it_gothit)->getTdc();
00754                     log << MSG::INFO <<"hit layer id  " <<layid<< " wire id: " <<localwid<< " adc_raw: " <<adc_raw<<" tdc_raw: "<<tdc_raw<<endreq;
00755                     zhit = (*it_gothit)->getZhit();
00756                     lr = (*it_gothit)->getFlagLR();
00757                     if(lr == 2) cout<<"lr = "<<lr<<endl;
00758                     driftD = (*it_gothit)->getDD();
00759                     driftd = abs(driftD);
00760                     driftT = (*it_gothit)->getDT();
00761                     dd_in = (*it_gothit)->getDocaIncl(); //getDocaIncl() include fit unused hit
00762                     dd_ex = (*it_gothit)->getDocaExcl(); //getDocaExcl() exclude fit unused hit
00763 
00764                     if(lr==-1 || lr == 2) {dd_in = dd_in; dd_ex = dd_ex;}
00765                     else if(lr ==1) {dd_in = -dd_in; dd_ex = -dd_ex;}
00766                     //dd = dd/doca_norm[layid];
00767                     //cout<<"lr "<<lr<<" dd_in: "<<dd_in<<" dd_ex: "<<dd_ex<<endl;
00768 
00769                     eangle = (*it_gothit)->getEntra();
00770                     eangle = eangle/M_PI;
00771                     pathlength=exsvc->PathL( ntpFlag, exhel_hit_in, layid, localwid, zhit);
00772                     rphi_path=pathlength * sintheta;
00773                     //cout<<"ph para check: "<<"adc_raw: "<<adc_raw<<"  dd_in: "<<dd_in<<" eangle: "<<eangle<<" zhit: "<<zhit<<" costheta: "<<costheta<<endl;
00774                     ph=exsvc->StandardHitCorrec(0,runFlag,ntpFlag,runNO,pathlength,wid,layid,adc_raw,dd_in,eangle,costheta); 
00775                     ph_hit = exsvc->StandardTrackCorrec(0,runFlag, ntpFlag, runNO, ph, costheta,tes );
00776                     //if(pathlength == -1)  
00777                     //cout<<"parameter1: "<<" layid: "<<layid<<" localwid: "<<localwid<<" driftd: "<<driftd<<" rphi_path: "<<rphi_path<<" pathlength: "<<pathlength<<" ph: "<<ph<<"  ph_hit: "<<ph_hit<<endl;
00778 
00779                     if(runNO<=5911 && runNO>=5854 && layid ==14) continue;
00780                     if(runNO<0)
00781                     {
00782                         if (layid<8)
00783                         {  
00784                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;                         
00785                         }
00786                         else if(layid >= 8)
00787                         {  
00788                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue ||rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue;  
00789                         }
00790                     }    
00791                     else if(runNO>=0)
00792                     {
00793                         if(layid <8)
00794                         {
00795                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;
00796                         }
00797                         else if(layid >= 8)
00798                         {  
00799                             if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue;
00800                         }
00801                     }
00802                     //cout<<"recAlg=1 para kal:  m_phraw: "<<adc_raw<<" m_exraw: "<<ph_hit<<" m_doca: "<<dd_in<<" m_eangle: "<<eangle<<" m_costheta: "<<costheta<<endl;
00803 
00804                     if (ph<0.0||ph==0) continue;
00805                     else
00806                     {
00807                         //-----------------------put data into TDS-----------------------------//
00808                         dedxhit = new RecMdcDedxHit;
00809                         dedxhit->setMdcHit(mdcHit);
00810                         dedxhit->setMdcKalHelixSeg(*it_gothit);
00811                         dedxhit->setMdcId(mdcid);
00812                         dedxhit->setFlagLR(lr);
00813                         dedxhit->setTrkId(trk_ex.trk_id());
00814                         dedxhit->setDedx(ph_hit);
00815                         dedxhit->setPathLength(pathlength); 
00816 
00817                         //---------------------------Fill root file----------------------------//
00818                         if(m_rootfile!= "no rootfile")                    
00819                         {
00820                             m_hitNo_wire->Fill(wid);
00821                         }
00822 
00823                         //-------------------------Fill ntuple n102---------------------------//
00824                         if ( ntpFlag ==2 )
00825                         {
00826                             m_charge1 = (*trk)->charge();
00827                             m_t01 = tes;
00828                             m_driftT = driftT;
00829                             m_tdc_raw = tdc_raw;
00830                             m_phraw = adc_raw;
00831                             m_exraw = ph_hit;
00832                             m_localwid = localwid;
00833                             m_wire = wid;
00834                             m_runNO1 = runNO;
00835                             m_nhit_hit = Nhits;
00836                             m_doca_in = dd_in;
00837                             m_doca_ex = dd_ex;
00838                             m_driftdist = driftD;
00839                             m_eangle = eangle;
00840                             m_costheta1 = costheta;
00841                             m_pathL = pathlength;
00842                             m_layer = layid;
00843                             m_ptrk1 = m_P;
00844                             m_ptrk_hit = p_hit;
00845                             //cout<<"adc_raw : "<<adc_raw<<"     "<<ph_hit<<"     "<<dd_in<<"  "<<dd_ex<<"  "<<eangle<<endl;
00846                             m_trackId1 = trk_ex.trk_id();
00847                             StatusCode status2 = m_nt2->write();
00848                             if ( status2.isFailure() )  
00849                               log << MSG::ERROR << "Cannot fill Ntuple n102" << endreq;
00850                         }
00851                         if(layid>3)   
00852                         {
00853                             phlist.push_back(ph);
00854                             phlist_hit.push_back(ph_hit);
00855                         }
00856                     }
00857                     dedxhitList->push_back( dedxhit );
00858                     dedxhitrefvec->push_back( dedxhit );
00859                 }//end of hit loop
00860                 trk_ex.set_phlist( phlist );
00861                 trk_ex.set_phlist_hit( phlist_hit );
00862                 trk_ex.setVecDedxHits( *dedxhitrefvec );
00863                 ex_trks.push_back(trk_ex );
00864                 m_trkalgs.push_back(m_trkalg); 
00865 
00866                 delete dedxhitrefvec;
00867                 phlist.clear();
00868                 phlist_hit.clear();
00869                 if(ntpFlag>0) trackNO3++;
00870             }//end of track loop
00871             log << MSG::DEBUG << "size in processing:  " << ex_trks.size() <<endreq;
00872         }//end of recalg==1
00873         //--------------------------------- Hit level finished ---------------------------------//
00874 
00875         //-------------------------------- Track level begin --------------------------------//
00876         if( ntrk != ex_trks.size())
00877           log << MSG::DEBUG <<"ntrkCol size: "<<ntrk<<" dedxtrk size: "<<ex_trks.size()<< endreq;
00878 
00879         double poca_x=0,poca_y=0,poca_z=0;
00880         float sintheta=0,costheta=0,ptrk=0,charge=0;
00881         int Nhit=0,Nphlisthit=0;
00882         int usedhit = 0;
00883         double phtm=0,median=0,geometric=0,geometric_trunc=0,harmonic=0,harmonic_trunc=0,transform=0,logtrunc=0;
00884 
00885         enum pid_dedx parType_temp(electron);
00886         float probpid_temp=-0.01,expect_temp=-0.01,sigma_temp=-0.01,chidedx_temp=-0.01;
00887 
00888         double dEdx_meas_hit=0;
00889         double dEdx_meas=0,dEdx_meas_esat=0,dEdx_meas_norun=0; 
00890         int trk_recalg=-1;
00891 
00892         int idedxid = 0;
00893         for(std::vector<MdcDedxTrk>::iterator it = ex_trks.begin(); it != ex_trks.end(); it++, idedxid++ )
00894         {
00895             log << MSG::DEBUG <<"-------------------------------------------------------"<< endreq;
00896             log << MSG::DEBUG <<" trk id ="<< it->trk_id()<<" : P ="<<it->P() <<"  Pt ="<<it->Pt()<<" : theta ="<<it->theta() <<" : phi ="<<it->phi()<< " : charge = "<<it->charge()<<endreq;
00897             log << MSG::DEBUG <<"all hits on track: "<<it->nsample()<<"  phlist size: "<<it->get_phlist().size()<<endreq;
00898 
00899             if(it->trk_ptr_kal()!=0)
00900             {
00901                 poca_x = it->trk_ptr_kal()->x();  //get poca, default pid is pion; change pid using setPidType();
00902                 poca_y = it->trk_ptr_kal()->y();
00903                 poca_z = it->trk_ptr_kal()->z();
00904             }
00905             else if(it->trk_ptr()!=0)
00906             {
00907                 poca_x = it->trk_ptr()->x();
00908                 poca_y = it->trk_ptr()->y();
00909                 poca_z = it->trk_ptr()->z();
00910             }
00911             //cout<<"poca_x: "<<poca_x<<" poca_y: "<<poca_y<<" poca_z: "<<poca_z<<endl;
00912 
00913             sintheta = sin(it->theta());
00914             costheta = cos(it->theta());
00915             ptrk = it->P();
00916             charge = it->charge();
00917             Nhit = it->nsample();   //total hits on track used as sample; 
00918             Nphlisthit = it->get_phlist_hit().size();  //hits in phlist_hit, exclude first 4 layers;
00919             //usedhit: hits after truncting phlist and used in cal dE/dx value;
00920 
00921             //--------------------------extrk truncation--------------------------------// 
00922             phtm = it->cal_dedx( truncate );
00923             //cout<<"phtm:  "<<phtm<<endl; 
00924             //median = it->cal_dedx_median( truncate );
00925             //geometric = it->cal_dedx_geometric( truncate );
00926             //geometric_trunc = it->cal_dedx_geometric_trunc( truncate );
00927             //harmonic = it->cal_dedx_harmonic( truncate );
00928             //harmonic_trunc = it->cal_dedx_harmonic_trunc( truncate );
00929             //transform = it->cal_dedx_transform( 0 );
00930             //logtrunc = it->cal_dedx_log( 1.0, 0);
00931 
00932             if(m_alg==1) dEdx_meas_hit = it->cal_dedx_bitrunc(truncate, 0, usedhit);
00933             else if(m_alg==2) dEdx_meas_hit = it-> cal_dedx_transform(1);
00934             else if(m_alg==3) dEdx_meas_hit = it->cal_dedx_log(1.0, 1);
00935             else cout<<"-------------Truncate Algorithm Error!!!------------------------"<<endl;
00936             if(m_alg==1 && usedhit==0)
00937             {
00938                 cout<<"***************bad extrk with no hits!!!!!******************"<<endl;
00939                 continue;
00940             }
00941             // force to use the same definition of usedhit in TDS and what used in setDedx() function
00942             usedhit = (int)(Nphlisthit*truncate);  
00943 
00944             //--------------------track level correction of extrk dE/dx Value---------------// 
00945             dEdx_meas = exsvc->StandardTrackCorrec(0,runFlag, ntpFlag, runNO, dEdx_meas_hit,  cos(it->theta()),  tes );
00946             dEdx_meas_esat = exsvc->StandardTrackCorrec(1,runFlag, ntpFlag, runNO, dEdx_meas_hit,  cos(it->theta()),  tes );
00947             dEdx_meas_norun = exsvc->StandardTrackCorrec(2,runFlag, ntpFlag, runNO, dEdx_meas_hit,  cos(it->theta()),  tes );
00948             log << MSG::INFO << "===================== " << endreq << endreq;
00949             log << MSG::DEBUG <<"dEdx_meas_hit: "<<dEdx_meas_hit<<" dEdx_meas: "<<dEdx_meas<<" dEdx_meas_esat: "<<dEdx_meas_esat<<" dEdx_meas_norun: "<<dEdx_meas_norun<<" ptrk: "<<it->P()<<endreq;
00950             log << MSG::INFO  << "ptrk " << ptrk << "  costhe " << costheta << "  runno " << runNO << "  evtno " << eventNO << endreq;
00951             //cout<<"dEdx_meas_hit: "<<dEdx_meas_hit<<" dEdx_meas: "<<dEdx_meas<<" dEdx_meas_esat: "<<dEdx_meas_esat<<" dEdx_meas_norun: "<<dEdx_meas_norun<<" ptrk: "<<it->P()<<endl;
00952 
00953             trk_recalg = m_trkalgs[idedxid];
00954 
00955             if(dEdx_meas<0 || dEdx_meas==0)  continue;
00956             it->set_dEdx( landau , dEdx_meas, trk_recalg, runFlag, vFlag , tes, Curve_Para, Sigma_Para, ex_calib); // calculate expect
00957             parType_temp = electron;
00958             probpid_temp=-0.01; 
00959             expect_temp=-0.01;
00960             sigma_temp=-0.01;
00961             chidedx_temp=-0.01;
00962             for(int i=0;i<5;i++)
00963             {   
00964                 if(it->pprob()[i]>probpid_temp)
00965                 {                        
00966                     parType_temp = (enum pid_dedx) i;  //e:0,  mu:1,  pi:2,  K:3,  p:4
00967                     probpid_temp = it->pprob()[i];
00968                     expect_temp = it->pexpect()[i];
00969                     sigma_temp = it->pexp_sigma()[i];
00970                     chidedx_temp = it->pchi_dedx()[i];
00971                 }
00972             }
00973             log<< MSG::INFO<<"expect dE/dx: type: "<<parType_temp<<" probpid: "<<probpid_temp<<" expect: "<<expect_temp<<" sigma: "<<sigma_temp<<" chidedx: "<<chidedx_temp<<endreq;
00974             //cout<<"dEdx_meas: "<<dEdx_meas<<endl;
00975 
00976             //-----------------------put data into TDS-----------------------------//
00977             resdedx = new RecMdcDedx;
00978             resdedx->setDedxHit(dEdx_meas_hit);
00979             resdedx->setDedxEsat(dEdx_meas_esat);
00980             resdedx->setDedxNoRun(dEdx_meas_norun);
00981             resdedx->setDedxMoment(it->P());
00982             resdedx->setTrackId( it->trk_id() );
00983             resdedx->setProbPH( dEdx_meas );
00984             resdedx->setNormPH( dEdx_meas/550.0 );
00985             resdedx->setDedxExpect( it->pexpect() );
00986             resdedx->setSigmaDedx( it->pexp_sigma() );
00987             resdedx->setPidProb( it->pprob() );
00988             resdedx->setChi( it->pchi_dedx() );
00989             //cout<<"recdedx: "<<"  "<<resdedx->getPidProb(parType_temp)<<"  "<<resdedx->getDedxExpect(parType_temp)<<"  "<<resdedx->getSigmaDedx(parType_temp)<<"  "<<resdedx->chi(parType_temp)<<endl;
00990             resdedx->setNumTotalHits(it->nsample() );  //all hits on track;
00991             resdedx->setNumGoodHits(usedhit);   //hits after truncing the phlist
00992             //phlist are all hits on track excluding first 4 layers; 
00993             //resdedx->setStatus( it->quality() );
00994             resdedx->setStatus(trk_recalg );
00995             //cout<<"trk_recalg: "<<trk_recalg<<"  trk stat: "<<it->quality()<<endl;
00996             resdedx->setTruncAlg( m_alg );
00997             resdedx->setParticleId(parType_temp);
00998             //cout<<"ParticleType: "<<parType_temp<<"  "<<resdedx->particleType()<<endl;
00999             resdedx->setVecDedxHits(it->getVecDedxHits());
01000             resdedx->setMdcTrack(it->trk_ptr());
01001             resdedx->setMdcKalTrack(it->trk_ptr_kal());
01002 
01003             //-------------------------Fill ntuple n103---------------------------//
01004             if(ntpFlag ==2)
01005             {
01006                 m_phtm=phtm;
01007                 //m_median=median;
01008                 //m_geometric=geometric;
01009                 //m_geometric_trunc=geometric_trunc;
01010                 //m_harmonic=harmonic;
01011                 //m_harmonic_trunc=harmonic_trunc;
01012                 //m_transform=transform;
01013                 //m_logtrunc=logtrunc;
01014                 m_dEdx_meas = dEdx_meas;
01015 
01016                 m_poca_x = poca_x;
01017                 m_poca_y = poca_y;
01018                 m_poca_z = poca_z;
01019                 m_ptrk=ptrk;
01020                 m_sintheta=sintheta;
01021                 m_costheta=costheta;
01022                 m_charge=charge;
01023                 m_runNO = runNO;
01024                 m_evtNO = eventNO;
01025 
01026                 m_t0 = tes;
01027                 m_trackId = it->trk_id();
01028                 m_recalg = trk_recalg;
01029 
01030                 m_nhit=Nhit;
01031                 m_nphlisthit=Nphlisthit; 
01032                 m_usedhit=usedhit;
01033                 for(int i=0; i<Nphlisthit; i++)  m_dEdx_hit[i] = it->get_phlist_hit()[i]; 
01034 
01035                 m_parttype = parType_temp;
01036                 m_prob = probpid_temp;
01037                 m_expect = expect_temp;
01038                 m_sigma = sigma_temp;
01039                 m_chidedx = chidedx_temp;
01040                 m_chidedxe=it->pchi_dedx()[0];
01041                 m_chidedxmu=it->pchi_dedx()[1];
01042                 m_chidedxpi=it->pchi_dedx()[2];
01043                 m_chidedxk=it->pchi_dedx()[3];
01044                 m_chidedxp=it->pchi_dedx()[4];
01045                 for(int i=0;i<5;i++)
01046                 {
01047                     m_probpid[i]= it->pprob()[i];
01048                     m_expectid[i]= it->pexpect()[i];
01049                     m_sigmaid[i]= it->pexp_sigma()[i];
01050                 }
01051 
01052                 StatusCode status = m_nt1->write();
01053                 if ( status.isFailure() )  
01054                 {
01055                     log << MSG::ERROR << "Cannot fill Ntuple n103" << endreq;
01056                 } 
01057             }
01058 
01059             log<< MSG::INFO<<"-----------------Summary of this ExTrack----------------"<<endreq
01060                 <<"dEdx_mean: "<<dEdx_meas<<"  type: "<<parType_temp<<" probpid: "<<probpid_temp
01061                 <<" expect: "<<expect_temp<<" sigma: "<<sigma_temp<<" chidedx: "<<chidedx_temp<<endreq<<endreq;
01062 
01063             dedxList->push_back( resdedx );            
01064         }//ExTrk loop end
01065     } //doNewGlobal==0 END . . .
01066 
01067 
01068     // check MdcDedxCol registered
01069     SmartDataPtr<RecMdcDedxCol> newexCol(eventSvc(),"/Event/Recon/RecMdcDedxCol");
01070     if (!newexCol) 
01071     { 
01072         log << MSG::FATAL << "Could not find RecMdcDedxCol" << endreq;
01073         return( StatusCode::SUCCESS);
01074     }       
01075     log << MSG::DEBUG << "----------------Begin to check RecMdcDedxCol-----------------"<<endreq; 
01076     RecMdcDedxCol::iterator it_dedx = newexCol->begin();
01077     for( ; it_dedx != newexCol->end(); it_dedx++) 
01078     {
01079         log << MSG::INFO << "retrieved MDC dE/dx:" << endreq
01080             << "dEdx Id: " << (*it_dedx)->trackId()
01081             << "   part Id: " << (*it_dedx)->particleType()
01082             << "   measured dEdx:  " << (*it_dedx)->probPH()
01083             << "   dEdx std:  " << (*it_dedx)->normPH() << endreq
01084             << "hits on track:   "<<(*it_dedx)->numTotalHits()
01085             << "   usedhits:  " << (*it_dedx)->numGoodHits() << endreq
01086             << "Electron: expect: " << (*it_dedx)->getDedxExpect(0)
01087             << "   sigma: " << (*it_dedx)->getSigmaDedx(0)
01088             << "   chi: " << (*it_dedx)->chi(0)
01089             << "   prob: " << (*it_dedx)->getPidProb(0) << endreq
01090             << "Muon: expect: " << (*it_dedx)->getDedxExpect(1)
01091             << "   sigma: " << (*it_dedx)->getSigmaDedx(1)  
01092             << "   chi: " << (*it_dedx)->chi(1)
01093             << "   prob: " << (*it_dedx)->getPidProb(1) << endreq
01094             << "Pion: expect: " << (*it_dedx)->getDedxExpect(2)
01095             << "   sigma: " << (*it_dedx)->getSigmaDedx(2)
01096             << "   chi: " << (*it_dedx)->chi(2)
01097             << "   prob: " << (*it_dedx)->getPidProb(2) << endreq
01098             << "Kaon: expect: " << (*it_dedx)->getDedxExpect(3)
01099             << "   sigma: " << (*it_dedx)->getSigmaDedx(3)
01100             << "   chi: " << (*it_dedx)->chi(3)
01101             << "   prob: " << (*it_dedx)->getPidProb(3) << endreq
01102             << "Proton: expect: " << (*it_dedx)->getDedxExpect(4)
01103             << "   sigma: " << (*it_dedx)->getSigmaDedx(4)
01104             << "   chi: " << (*it_dedx)->chi(4)
01105             << "   prob: " << (*it_dedx)->getPidProb(4) << endreq;
01106     }
01107     return StatusCode::SUCCESS;
01108 }

StatusCode MdcDedxRecon::finalize ( )

Definition at line 257 of file MdcDedxRecon.cxx.

References Bes_Common::DEBUG, eventNo, ex_calib, ex_trks, f1, Bes_Common::INFO, m_hitlevel, m_hitNo_wire, m_hlist, m_rootfile, m_trkalgs, msgSvc(), trackNO1, trackNO2, and trackNO3.

00258 {
00259     MsgStream log(msgSvc(), name());
00260     log << MSG::INFO << "in finalize()" << endreq;
00261 
00262     ex_trks.clear();
00263     m_trkalgs.clear();
00264 
00265     if(m_rootfile != "no rootfile")
00266     {
00267         TFile *f1 = new TFile(m_rootfile.c_str(),"recreate");
00268         m_hlist->Write();
00269         f1->Close();
00270         delete f1;
00271         delete m_hitNo_wire;
00272         delete m_hitlevel;
00273         delete m_hlist;
00274     }
00275     delete ex_calib;
00276 
00277     cout<<"total event number is : "<<eventNo<<endl;
00278     cout<<"total track number is : "<<trackNO1
00279         <<"     RecMdcTrack number is : "<<trackNO2
00280         <<"     RecMdcKalTrack number is :"<<trackNO3<<endl;
00281     log << MSG::DEBUG <<"MdcDedxRecon terminated!!!"<< endreq;
00282     return StatusCode::SUCCESS;
00283 }

StatusCode MdcDedxRecon::initialize ( )

Definition at line 76 of file MdcDedxRecon.cxx.

References calib_flag, Bes_Common::DEBUG, dedxcursvc, doNewGlobal, calibUtil::ERROR, ex_calib, exsvc, Bes_Common::FATAL, Bes_Common::INFO, landau, m_alg, m_charge, m_charge1, m_chidedx, m_chidedxe, m_chidedxk, m_chidedxmu, m_chidedxp, m_chidedxpi, m_costheta, m_costheta1, m_dEdx_hit, m_dEdx_meas, m_doca_ex, m_doca_in, m_driftdist, m_driftT, m_eangle, m_evtNO, m_expect, m_expectid, m_exraw, m_hitlevel, m_hitNo_wire, m_hlist, m_layer, m_localwid, m_nhit, m_nhit_hit, m_nphlisthit, m_nt1, m_nt2, m_parttype, m_pathL, m_phraw, m_phtm, m_poca_x, m_poca_y, m_poca_z, m_prob, m_probpid, m_ptrk, m_ptrk1, m_ptrk_hit, m_recalg, m_rootfile, m_runNO, m_runNO1, m_sigma, m_sigmaid, m_sintheta, m_t0, m_t01, m_tdc_raw, m_trackId, m_trackId1, m_usedhit, m_wire, msgSvc(), ntpFlag, ntupleSvc(), recalg, IDedxCorrecSvc::set_flag(), and truncate.

00077 {
00078     MsgStream log(msgSvc(), name());
00079     log << MSG::INFO << "in initialize()" << endreq;
00080 
00081     log << MSG::INFO <<"--------------------<  MdcDedxRecon  V2.1   >---------------------"<< endreq;
00082     log << MSG::INFO <<"#######  correction for the wire current dependence   #######"<< endreq;
00083     log << MSG::INFO <<"#######      correction for the new z dependence      #######"<< endreq;
00084     log << MSG::INFO <<"-------------------------------------------------------------"<< endreq;
00085     log << MSG::INFO <<"++++++++++++++++++++[  initialization  ]+++++++++++++++++++++"<< endreq;
00086     log << MSG::INFO << "MdcDedxRecon has been initialized with calib_flag: "<< calib_flag <<endreq;
00087     log << MSG::INFO << "MdcDedxRecon has been initialized with landau: "<< landau <<  endreq;
00088     if( landau == 0 ) {truncate = 0.7;}
00089     log << MSG::INFO << "MdcDedxRecon has been initialized with ntpFlag: "<< ntpFlag <<  endreq;
00090     log << MSG::INFO << "MdcDedxRecon has been initialized with recalg: "<< recalg <<  endreq;
00091     log << MSG::INFO << "MdcDedxRecon has been initialized with dE/dx cal method m_alg: "<< m_alg << endreq;
00092     log << MSG::INFO << "MdcDedxRecon has been initialized with truncate: "<< truncate <<endreq;
00093     log << MSG::INFO << "MdcDedxRecon has been initialized with doNewGlobal: "<<doNewGlobal<<endreq;
00094     log << MSG::DEBUG  << "+++++++++++ MdcDedxRecon initialization end ++++++++++++ "<< endreq;
00095     if( truncate <= 0.0 || 1.0 < truncate )
00096     {
00097         log << MSG::FATAL <<" Initialization ERROR of truncate  = "<<truncate<< endreq;
00098         log << MSG::FATAL <<" truncate must be within 0.0 to 1.0 ! "<< endreq;
00099         log << MSG::FATAL <<" Please stop exec. "<< endreq;
00100     }
00101     log << MSG::DEBUG <<"-------------------------------------------------------------"<< endreq;
00102     log << MSG::DEBUG <<"MdcDedxRecon init 2nd part!!!"<< endreq;
00103 
00104     StatusCode scex = service("DedxCorrecSvc", exsvc);
00105     if (scex ==  StatusCode::SUCCESS) 
00106     {
00107         log << MSG::INFO <<"Hi, DedxCorrectSvc is running"<<endreq;
00108     } 
00109     else 
00110     {
00111         log << MSG::ERROR <<"DedxCorrectSvc is failed"<<endreq;
00112         return StatusCode::SUCCESS;
00113     }
00114     exsvc->set_flag( calib_flag );
00115 
00116     StatusCode cursvc = service("DedxCurSvc", dedxcursvc);
00117     if (cursvc ==  StatusCode::SUCCESS) 
00118     {
00119         log << MSG::INFO <<"DedxCurSvc is running"<<endreq;
00120     } 
00121     else 
00122     {
00123         log << MSG::ERROR <<"DedxCurSvc is failed"<<endreq;
00124         return StatusCode::SUCCESS;
00125     }
00126 
00127 
00128     if( !ex_calib ) ex_calib = new MdcDedxCorrection;
00129     //#ifdef DBCALIB
00130     //    ex_calib->getCalib();  //cleate MdcDedxWire and MdcDedxLayer.
00131     //#endif
00132 
00133     //------------------------------produce histogram root files-----------------------------//
00134     if(m_rootfile!="no rootfile")
00135     {
00136         const char* preDir = gDirectory->GetPath();
00137         m_hlist  = new TObjArray(0); 
00138         m_hitlevel = new TFolder("dedx_hitlevel","hitlevel");
00139         m_hlist -> Add(m_hitlevel);
00140         m_hitNo_wire = new TH1F("dedx_HitNo_wire","dedx hitNo vs wire",6797, -0.5, 6796.5);
00141         m_hitlevel -> Add(m_hitNo_wire);
00142         gDirectory->cd(preDir);
00143     }
00144 
00145     //------------------------------produce ntuple files-------------------------------------//
00146     if( ntpFlag ==2 )
00147     {
00148         StatusCode status;
00149         NTuplePtr  nt1(ntupleSvc(),"FILE103/n103");
00150         if ( nt1 )
00151           m_nt1 = nt1;
00152         else 
00153         {
00154             m_nt1= ntupleSvc()->book("FILE103/n103",CLID_ColumnWiseTuple,"dEdx vs momentum");
00155             if ( m_nt1 )
00156             {
00157                 status = m_nt1->addItem("phtm",m_phtm);
00158                 //status = m_nt1->addItem("median",m_median);
00159                 //status = m_nt1->addItem("geom",m_geometric);
00160                 //status = m_nt1->addItem("geom_tr",m_geometric_trunc);
00161                 //status = m_nt1->addItem("harm",m_harmonic);
00162                 //status = m_nt1->addItem("harm_tr",m_harmonic_trunc);
00163                 //status = m_nt1->addItem("transf",m_transform);
00164                 //status = m_nt1->addItem("logex",m_logtrunc);
00165                 status = m_nt1->addItem("dEdx_meas", m_dEdx_meas);
00166 
00167                 status = m_nt1->addItem("ptrk",m_ptrk);
00168                 status = m_nt1->addItem("sintheta",m_sintheta);
00169                 status = m_nt1->addItem("costheta",m_costheta);
00170                 status = m_nt1->addItem("charge",m_charge);
00171                 status = m_nt1->addItem("runNO",m_runNO);
00172                 status = m_nt1->addItem("evtNO",m_evtNO);
00173                 status = m_nt1->addItem("t0",m_t0);
00174                 status = m_nt1->addItem("trackId",m_trackId); 
00175                 status = m_nt1->addItem("poca_x",m_poca_x);
00176                 status = m_nt1->addItem("poca_y",m_poca_y);
00177                 status = m_nt1->addItem("poca_z",m_poca_z);
00178                 status = m_nt1->addItem("recalg",m_recalg);
00179 
00180                 status = m_nt1->addItem("nhit",m_nhit);
00181                 status = m_nt1->addItem("usedhit",m_usedhit);
00182                 status = m_nt1->addItem("ndedxhit",m_nphlisthit,0,100);
00183                 status = m_nt1->addIndexedItem("dEdx_hit",m_nphlisthit,m_dEdx_hit);
00184 
00185                 status = m_nt1->addItem("type",m_parttype);
00186                 status = m_nt1->addItem("prob",m_prob);
00187                 status = m_nt1->addItem("expect",m_expect);
00188                 status = m_nt1->addItem("sigma",m_sigma);
00189                 status = m_nt1->addItem("chidedx",m_chidedx);
00190                 status = m_nt1->addItem("chidedx_e",m_chidedxe);
00191                 status = m_nt1->addItem("chidedx_mu",m_chidedxmu);
00192                 status = m_nt1->addItem("chidedx_pi",m_chidedxpi);
00193                 status = m_nt1->addItem("chidedx_k",m_chidedxk);
00194                 status = m_nt1->addItem("chidedx_p",m_chidedxp);
00195                 status = m_nt1->addItem("partid",5,m_probpid);
00196                 status = m_nt1->addItem("expectid",5,m_expectid);
00197                 status = m_nt1->addItem("sigmaid",5,m_sigmaid);
00198             }
00199         }
00200 
00201         NTuplePtr  nt2(ntupleSvc(),"FILE103/n102");
00202         if ( nt2 )  m_nt2 = nt2;
00203         else 
00204         {
00205             m_nt2= ntupleSvc()->book("FILE103/n102",CLID_RowWiseTuple,"pulae height raw");
00206             if ( m_nt2 )
00207             {
00208                 status = m_nt2->addItem("charge",m_charge1);
00209                 status = m_nt2->addItem("adc_raw",m_phraw);
00210                 status = m_nt2->addItem("exraw",m_exraw);
00211                 status = m_nt2->addItem("runNO1",m_runNO1);
00212                 status = m_nt2->addItem("nhit_hit", m_nhit_hit);
00213                 status = m_nt2->addItem("wire",m_wire);
00214                 //status = m_nt2->addItem("doca",m_doca);
00215                 status = m_nt2->addItem("doca_in",m_doca_in);
00216                 status = m_nt2->addItem("doca_ex",m_doca_ex);
00217                 status = m_nt2->addItem("driftdist",m_driftdist);
00218                 status = m_nt2->addItem("eangle",m_eangle);
00219                 status = m_nt2->addItem("costheta1",m_costheta1);
00220                 status = m_nt2->addItem("path_rphi",m_pathL);
00221                 status = m_nt2->addItem("layer",m_layer);
00222                 status = m_nt2->addItem("ptrk1",m_ptrk1);
00223                 status = m_nt2->addItem("ptrk_hit",m_ptrk_hit);
00224                 status = m_nt2->addItem("t01",m_t01); 
00225                 status = m_nt2->addItem("tdc_raw",m_tdc_raw);  
00226                 status = m_nt2->addItem("driftT",m_driftT);
00227                 status = m_nt2->addItem("localwid",m_localwid);
00228                 status = m_nt2->addItem("trackId1",m_trackId1);
00229             }
00230         }
00231     }
00232 
00233     return StatusCode::SUCCESS;
00234 }

void MdcDedxRecon::kaltrackrec	(	RecMdcKalTrackCol::iterator	trk,
		Identifier	mdcid,
		double	tes,
		int	RunNO,
		int	runFlag,
		MsgStream	log
	)

Definition at line 1286 of file MdcDedxRecon.cxx.

References abs, cos(), Bes_Common::DEBUG, dedxhit, dedxhitList, dedxhitrefvec, calibUtil::ERROR, ex_trks, exsvc, geosvc, DstMdcKalTrack::getFError(), DstMdcKalTrack::getFHelix(), Iner_DriftDistCut, Bes_Common::INFO, IMdcGeomSvc::Layer(), MdcID::layer(), m_charge1, m_costheta1, m_doca_ex, m_doca_in, m_driftdist, m_driftT, m_eangle, m_exraw, m_hitNo_wire, m_layer, m_localwid, m_nhit_hit, m_nt2, m_pathL, m_phraw, M_PI, m_ptrk1, m_ptrk_hit, m_rootfile, m_runNO1, m_t01, m_tdc_raw, m_trackId1, m_trkalg, m_trkalgs, m_wire, MaxHistValue, ntpFlag, Out_DriftDistCut, ParticleFlag, IDedxCorrecSvc::PathL(), phi0, RPhi_PathMaxCut, MdcDedxTrk::set_phlist(), MdcDedxTrk::set_phlist_hit(), RecMdcDedxHit::setDedx(), RecMdcDedxHit::setFlagLR(), RecMdcDedxHit::setMdcHit(), RecMdcDedxHit::setMdcId(), RecMdcDedxHit::setMdcKalHelixSeg(), RecMdcDedxHit::setPathLength(), RecMdcDedxHit::setTrkId(), MdcDedxTrk::setVecDedxHits(), sin(), IDedxCorrecSvc::StandardHitCorrec(), IDedxCorrecSvc::StandardTrackCorrec(), trackNO3, MdcDedxTrk::trk_id(), MdcID::wire(), and MdcGeoLayer::Wirst().

Referenced by execute().

01287 {
01288     vector<double> phlist;  //dE/dx only after hit level correction
01289     vector<double> phlist_hit; //dE/dx after hit and track level correction
01290     double m_d0E=0,m_phi0E=0,m_cpaE=0,m_z0E=0,phi0=0,costheta=0,sintheta=0,m_Pt=0,m_P=0;
01291     int charge=0,layid=0,localwid=0,w0id=0,wid=0,lr=-2;
01292     double p_hit=0,adc_raw=0,tdc_raw=0,zhit=0,driftd=0,driftD=0,driftT=0,dd_in=0,dd_ex=0,eangle=0;
01293     int Nhits=0;
01294     double ph=0,ph_hit=0,pathlength=0,rphi_path=0;
01295 
01296     MdcDedxTrk trk_ex( **trk, ParticleFlag);
01297     log << MSG::DEBUG <<" %%%%% trk id = "<<trk_ex.trk_id() <<endreq;
01298     log << MSG::DEBUG <<" %%%%% initial id = "<<(*trk)->trackId() <<endreq;
01299 
01300     HepVector a;
01301     HepSymMatrix tkErrM;
01302     if(ParticleFlag>-1&&ParticleFlag<5)
01303     {
01304         DstMdcKalTrack* dstTrk = *trk;
01305         a = dstTrk->getFHelix(ParticleFlag);
01306         tkErrM = dstTrk->getFError(ParticleFlag);
01307     }
01308     else
01309     {
01310         a = (*trk)->getFHelix();
01311         tkErrM = (*trk)->getFError();
01312     }
01313     log << MSG::DEBUG <<"FHelix 5 parameters: "<<a[0]<<"     "<<a[1] <<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endreq;
01314     //cout<<"FHelix 5 parameters: "<<a[0]<<"     "<<a[1] <<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endl;  //getFHelix is first layer of MdcKalTrack;
01315 
01316     m_d0E   = tkErrM[0][0];
01317     m_phi0E = tkErrM[1][1];
01318     m_cpaE  = tkErrM[2][2];
01319     m_z0E   = tkErrM[3][3];
01320 
01321     phi0= a[1];
01322     costheta = cos(M_PI_2-atan(a[4]));
01323     //cout<<"track theta: "<<M_PI_2-atan(a[4])<<" extrack theta: "<<trk_ex.theta()<<endl;  //theta in trk_ex is got by getFHelix();
01324     sintheta = sin(M_PI_2-atan(a[4]));
01325     m_Pt = 1.0/fabs( a[2] );
01326     m_P = m_Pt*sqrt(1 + a[4]*a[4]);
01327     //cout<<"track ptrk: "<<m_P<<" extrack ptrk: "<<trk_ex.P()<<endl;
01328     charge = ( a[2] > 0 )? 1 : -1;
01329     //cout<<"track charge: "<<charge<<" extrack charge: "<<(*trk)->charge()<<endl;
01330     dedxhitrefvec = new DedxHitRefVec;
01331 
01332 
01333     HelixSegRefVec gothits= (*trk)->getVecHelixSegs(ParticleFlag);
01334     //HelixSegRefVec gothits= (*trk)->getVecHelixSegs();
01335     //if not set ParticleFlag, we will get the last succefully hypothesis;
01336     Nhits = gothits.size();
01337     log << MSG::DEBUG << "hits size  = " <<gothits.size()<<endreq;
01338     //if(gothits.size()<15) cout<<"eventNO : "<<eventNO<<"  hits : "<<gothits.size()<<endl;
01339 
01340     RecMdcHit* mdcHit = 0;
01341     HelixSegRefVec::iterator it_gothit = gothits.begin();
01342     for( ; it_gothit != gothits.end(); it_gothit++)
01343     {
01344         HepVector a_hit_in = (*it_gothit)->getHelixIncl();
01345         //HepVector a_hit_ex = (*it_gothit)->getHelixExcl(); //same with getHelixIncl()
01346         HepPoint3D IP(0,0,0);
01347         Dedx_Helix exhel_hit_in(IP, a_hit_in);
01348         //Dedx_Helix exhel_hit_ex(IP, a_hit_ex);
01349         p_hit = 1.0/fabs(a_hit_in(3))*sqrt(1+a_hit_in(5)*a_hit_in(5));
01350         //cout<<"getHelixIncl 5 para: "<<a_hit_in[0]<<"     "<<a_hit_in[1] <<"      "<<a_hit_in[2]<<"      "<<a_hit_in[3]<<"     "<<a_hit_in[4]<<endl;
01351 
01352         mdcid = (*it_gothit)->getMdcId();
01353         layid =  MdcID::layer(mdcid);
01354         localwid = MdcID::wire(mdcid);                    
01355         w0id = geosvc->Layer(layid)->Wirst();
01356         wid = w0id + localwid;
01357         adc_raw = (*it_gothit)->getAdc();
01358         tdc_raw = (*it_gothit)->getTdc();
01359         log << MSG::INFO <<"hit layer id  " <<layid<< " wire id: " <<localwid<< " adc_raw: " <<adc_raw<<" tdc_raw: "<<tdc_raw<<endreq;
01360         zhit = (*it_gothit)->getZhit();
01361         lr = (*it_gothit)->getFlagLR();
01362         if(lr == 2) cout<<"lr = "<<lr<<endl;
01363         driftD = (*it_gothit)->getDD();
01364         driftd = abs(driftD);
01365         driftT = (*it_gothit)->getDT();
01366         dd_in = (*it_gothit)->getDocaIncl(); //getDocaIncl() include fit unused hit
01367         dd_ex = (*it_gothit)->getDocaExcl(); //getDocaExcl() exclude fit unused hit
01368 
01369         if(lr==-1 || lr == 2) {dd_in = dd_in; dd_ex = dd_ex;}
01370         else if(lr ==1) {dd_in = -dd_in; dd_ex = -dd_ex;}
01371         //dd = dd/doca_norm[layid];
01372         //cout<<"lr: "<<lr<<" dd_in: "<<dd_in<<" dd_ex: "<<dd_ex<<endl;
01373 
01374         eangle = (*it_gothit)->getEntra();
01375         eangle = eangle/M_PI;
01376         pathlength=exsvc->PathL( ntpFlag, exhel_hit_in, layid, localwid, zhit);
01377         rphi_path=pathlength * sintheta;
01378 
01379         //cout<<"ph para check: "<<"runFlag: "<<runFlag<<" runNO: "<<runNO<<" pathlength: "<<pathlength<<" wid: "<<wid<<" layid: "<<layid<<" adc_raw: "<<adc_raw<<" dd_in: "<<dd_in<<" eangle: "<<eangle<<" zhit: "<<zhit<<" costheta: "<<costheta<<endl;
01380         ph=exsvc->StandardHitCorrec(0,runFlag,ntpFlag,runNO,pathlength,wid,layid,adc_raw,dd_in,eangle,costheta);
01381         //cout<<"tes= "<<tes<<endl;
01382         ph_hit = exsvc->StandardTrackCorrec(0,runFlag, ntpFlag, runNO, ph, costheta,tes );
01383         //cout<<"adc_raw=  "<<adc_raw<<" ph= "<<ph<<endl;
01384         //cout<<"adc_raw = "<<adc_raw<<" ph_hit= "<<ph_hit<<endl;
01385         //if(pathlength == -1)  
01386         //cout<<"parameter1: "<<" layid: "<<layid<<" localwid: "<<localwid<<" driftd: "<<driftd<<" rphi_path: "<<rphi_path<<" pathlength: "<<pathlength<<" ph: "<<ph<<"  ph_hit: "<<ph_hit<<endl;
01387 
01388         if(runNO<=5911 && runNO>=5854 && layid ==14) continue;
01389         if(runNO<0)
01390         {
01391             if (layid<8)
01392             {
01393                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;
01394             }
01395             else if(layid >= 8)
01396             {
01397                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue ||rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue;
01398             }
01399         }
01400         else if(runNO>=0)
01401         {
01402             if(layid <8)
01403             {
01404                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;
01405             }
01406             else if(layid >= 8)
01407             {
01408                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue;
01409             }
01410         }
01411         //cout<<"recAlg=2 para kal:  m_phraw: "<<adc_raw<<" m_exraw: "<<ph_hit<<" m_doca: "<<dd_in<<" m_eangle: "<<eangle<<" m_costheta: "<<costheta<<endl;   
01412 
01413         if (ph<0.0||ph==0) continue;
01414         else
01415         {
01416             //-----------------------put data into TDS-----------------------------//
01417             dedxhit = new RecMdcDedxHit;
01418             dedxhit->setMdcHit(mdcHit);
01419             dedxhit->setMdcKalHelixSeg(*it_gothit);
01420             dedxhit->setMdcId(mdcid);
01421             dedxhit->setFlagLR(lr);
01422             dedxhit->setTrkId(trk_ex.trk_id());
01423             dedxhit->setDedx(ph_hit);
01424             dedxhit->setPathLength(pathlength);
01425 
01426             //---------------------------Fill root file----------------------------//
01427             if(m_rootfile!= "no rootfile")
01428             {
01429                 m_hitNo_wire->Fill(wid);
01430             }
01431 
01432             //-------------------------Fill ntuple n102---------------------------//
01433             if ( ntpFlag ==2 )
01434             {
01435                 m_charge1 = (*trk)->charge();
01436                 m_t01 = tes;
01437                 m_driftT = driftT;
01438                 m_tdc_raw = tdc_raw;
01439                 m_phraw = adc_raw;
01440                 m_exraw = ph_hit;
01441                 m_localwid = localwid;
01442                 m_wire = wid;
01443                 m_runNO1 = runNO;
01444                 m_nhit_hit = Nhits;
01445                 m_doca_in = dd_in;
01446                 m_doca_ex = dd_ex;
01447                 m_driftdist = driftD;
01448                 m_eangle = eangle;
01449                 m_costheta1 = costheta;
01450                 m_pathL = pathlength;
01451                 m_layer = layid;
01452                 m_ptrk1 = m_P;
01453                 m_ptrk_hit = p_hit;
01454                 //cout<<"adc_raw : "<<adc_raw<<"     "<<ph_hit<<"     "<<dd_in<<"  "<<dd_ex<<"  "<<eangle<<endl;
01455                 m_trackId1 = trk_ex.trk_id();
01456                 StatusCode status2 = m_nt2->write();
01457                 if ( status2.isFailure() )
01458                   log << MSG::ERROR << "Cannot fill Ntuple n102" << endreq;
01459             }
01460             if(layid>3)   
01461             {
01462                 phlist.push_back(ph);
01463                 phlist_hit.push_back(ph_hit);
01464             }
01465         }
01466         dedxhitList->push_back( dedxhit );
01467         dedxhitrefvec->push_back( dedxhit );
01468     }//end of hit loop
01469     trk_ex.set_phlist( phlist );
01470     trk_ex.set_phlist_hit( phlist_hit );
01471     trk_ex.setVecDedxHits( *dedxhitrefvec );
01472     ex_trks.push_back(trk_ex );
01473     m_trkalgs.push_back(m_trkalg);
01474 
01475     delete dedxhitrefvec;
01476     phlist.clear();
01477     phlist_hit.clear();
01478     if(ntpFlag>0) trackNO3++;
01479 }

void MdcDedxRecon::mdctrackrec	(	RecMdcTrackCol::iterator	trk,
		Identifier	mdcid,
		double	tes,
		int	RunNO,
		int	runFlag,
		MsgStream	log
	)

Definition at line 1117 of file MdcDedxRecon.cxx.

References abs, cos(), Bes_Common::DEBUG, dedxhit, dedxhitList, dedxhitrefvec, calibUtil::ERROR, ex_trks, exsvc, geosvc, Iner_DriftDistCut, Bes_Common::INFO, IMdcGeomSvc::Layer(), MdcID::layer(), m_charge1, m_costheta1, m_doca_ex, m_doca_in, m_driftdist, m_driftT, m_eangle, m_exraw, m_hitNo_wire, m_layer, m_localwid, m_nhit_hit, m_nt2, m_pathL, m_phraw, M_PI, m_ptrk1, m_rootfile, m_runNO1, m_t01, m_tdc_raw, m_trackId1, m_trkalg, m_trkalgs, m_wire, MaxHistValue, ntpFlag, Out_DriftDistCut, IDedxCorrecSvc::PathL(), phi0, RPhi_PathMaxCut, MdcDedxTrk::set_phlist(), MdcDedxTrk::set_phlist_hit(), RecMdcDedxHit::setDedx(), RecMdcDedxHit::setFlagLR(), RecMdcDedxHit::setMdcHit(), RecMdcDedxHit::setMdcId(), RecMdcDedxHit::setMdcKalHelixSeg(), RecMdcDedxHit::setPathLength(), RecMdcDedxHit::setTrkId(), MdcDedxTrk::setVecDedxHits(), sin(), IDedxCorrecSvc::StandardHitCorrec(), IDedxCorrecSvc::StandardTrackCorrec(), trackNO2, MdcDedxTrk::trk_id(), MdcID::wire(), and MdcGeoLayer::Wirst().

Referenced by switchtomdctrack().

01118 {
01119     vector<double> phlist;  //dE/dx only after hit level correction
01120     vector<double> phlist_hit; //dE/dx after hit and track level correction
01121     double m_d0E=0,m_phi0E=0,m_cpaE=0,m_z0E=0,phi0=0,costheta=0,sintheta=0,m_Pt=0,m_P=0;
01122     int charge=0,layid=0,localwid=0,w0id=0,wid=0,lr=-2;
01123     double adc_raw=0,tdc_raw=0,zhit=0,driftd=0,driftD=0,driftT=0,dd=0,eangle=0;
01124     int Nhits=0;
01125     double ph=0,ph_hit=0,pathlength=0,rphi_path=0;
01126 
01127     MdcDedxTrk trk_ex( **trk);
01128     log << MSG::DEBUG <<" %%%%% trk id = "<<trk_ex.trk_id() <<endreq;
01129     log << MSG::DEBUG <<" %%%%% initial id = "<<(*trk)->trackId() <<endreq;
01130 
01131     HepVector a = (*trk)->helix();
01132     HepSymMatrix tkErrM = (*trk)->err();
01133     m_d0E   = tkErrM[0][0];
01134     m_phi0E = tkErrM[1][1];
01135     m_cpaE  = tkErrM[2][2];
01136     m_z0E   = tkErrM[3][3];
01137 
01138     HepPoint3D IP(0,0,0);
01139     Dedx_Helix exhel(IP, a);  //initialize class Dedx_Helix for DedxCorrecSvc
01140     log << MSG::DEBUG <<"MDC Helix 5 parameters: "<<a[0]<<"     "<<a[1]<<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endreq;
01141     //cout<<"MDC Helix 5 parameters: "<<a[0]<<"     "<<a[1]<<"      "<<a[2]<<"      "<<a[3]<<"     "<<a[4]<<endl;
01142     phi0= a[1];
01143     costheta = cos(M_PI_2-atan(a[4]));
01144     //cout<<"track theta: "<<M_PI_2-atan(a[4])<<" extrack theta: "<<trk_ex.theta()<<endl;
01145     sintheta = sin(M_PI_2-atan(a[4]));
01146     m_Pt = 1.0/fabs( a[2] );
01147     m_P = m_Pt*sqrt(1 + a[4]*a[4]);
01148     charge = ( a[2] > 0 )? 1 : -1;
01149     //cout<<"track charge: "<<charge<<" extrack charge: "<<trk_ex.charge()<<endl;
01150     dedxhitrefvec = new DedxHitRefVec;
01151 
01152 
01153     HitRefVec gothits= (*trk)->getVecHits();
01154     Nhits = gothits.size();
01155     log << MSG::DEBUG << "hits size  = " <<gothits.size()<<endreq;
01156     //if(gothits.size()<15) cout<<"eventNO : "<<eventNO<<"  hits : "<<gothits.size()<<endl;
01157 
01158     RecMdcKalHelixSeg* mdcKalHelixSeg = 0;
01159     HitRefVec::iterator it_gothit = gothits.begin();
01160     for( ; it_gothit != gothits.end(); it_gothit++)
01161     {
01162         mdcid = (*it_gothit)->getMdcId();
01163         layid =  MdcID::layer(mdcid);
01164         localwid = MdcID::wire(mdcid);
01165         w0id = geosvc->Layer(layid)->Wirst();
01166         wid = w0id + localwid;
01167         adc_raw = (*it_gothit)->getAdc();
01168         tdc_raw = (*it_gothit)->getTdc();
01169         log << MSG::INFO <<"hit layer id  " <<layid<< " wire id: " <<localwid<< " adc_raw: " <<adc_raw<<" tdc_raw: "<<tdc_raw<<endreq;
01170         zhit = (*it_gothit)->getZhit();
01171         lr = (*it_gothit)->getFlagLR();
01172         if(lr == 2) cout<<"lr = "<<lr<<endl;
01173         if(lr == 0 || lr == 2) driftD = (*it_gothit)->getDriftDistLeft();
01174         else  driftD = (*it_gothit)->getDriftDistRight();
01175         //cout<<"lr: "<<lr<<"  driftD: "<<driftD<<endl;
01176         driftd = abs(driftD);
01177         dd = (*it_gothit)->getDoca();
01178         if(lr == 0 || lr == 2 ) dd = -abs(dd);
01179         if(lr == 1 ) dd = abs(dd);
01180         driftT = (*it_gothit)->getDriftT();
01181 
01182         eangle = (*it_gothit)->getEntra();
01183         eangle = eangle/M_PI;
01184         pathlength=exsvc->PathL( ntpFlag, exhel, layid, localwid, zhit);
01185         rphi_path=pathlength * sintheta;
01186 
01187         //cout<<"ph para check: "<<"adc_raw: "<<adc_raw<<"  dd: "<<dd<<" eangle: "<<eangle<<" zhit: "<<zhit<<" costheta: "<<costheta<<endl;
01188         ph=exsvc->StandardHitCorrec(0,runFlag,ntpFlag,runNO,pathlength,wid,layid,adc_raw,dd,eangle,costheta);
01189         ph_hit = exsvc->StandardTrackCorrec(0,runFlag, ntpFlag, runNO, ph, costheta,tes );
01190         //if(pathlength == -1)
01191         //cout<<"parameter0: "<<"eventNO: "<<eventNO<<" layid: "<<layid<<" localwid: "<<localwid<<" driftd: "<<driftd<<" rphi_path: "<<rphi_path<<" pathlength: "<<pathlength<<" ph: "<<ph<<" ph_hit: "<<ph_hit<<endl;
01192 
01193         if(runNO<=5911 && runNO>=5854 && layid ==14) continue;
01194         if(runNO<0)
01195         {
01196             if (layid<8)
01197             {
01198                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;
01199             }
01200             else if(layid >= 8)
01201             {
01202                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue ||rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue; 
01203             }
01204         }
01205         else if(runNO>=0)
01206         {
01207             if(layid <8)
01208             {
01209                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Iner_RPhi_PathMinCut || driftd>Iner_DriftDistCut) continue;
01210             }
01211             else if(layid >= 8)
01212             {
01213                 if(adc_raw<MinHistValue || adc_raw>MaxHistValue || rphi_path>RPhi_PathMaxCut || rphi_path<Out_RPhi_PathMinCut || driftd>Out_DriftDistCut ) continue;
01214             }
01215         }
01216         //cout<<"recAlg=2 para mdc:  m_phraw: "<<adc_raw<<" m_exraw: "<<ph_hit<<" m_doca: "<<dd<<" m_eangle: "<<eangle<<" m_costheta: "<<costheta<<endl;
01217 
01218         if (ph<0.0||ph==0) continue;
01219         else
01220         {
01221             //-----------------------put data into TDS-----------------------------//
01222             dedxhit = new RecMdcDedxHit;
01223             dedxhit->setMdcHit(*it_gothit);
01224             dedxhit->setMdcKalHelixSeg(mdcKalHelixSeg);
01225             dedxhit->setMdcId(mdcid);
01226             dedxhit->setFlagLR(lr);
01227             dedxhit->setTrkId(trk_ex.trk_id());
01228             dedxhit->setDedx(ph_hit);
01229             dedxhit->setPathLength(pathlength);
01230 
01231             //---------------------------Fill root file----------------------------//
01232             if(m_rootfile!= "no rootfile")
01233             {
01234                 m_hitNo_wire->Fill(wid);
01235             }
01236 
01237             //-------------------------Fill ntuple n102---------------------------//
01238             if ( ntpFlag ==2 )
01239             {
01240                 m_charge1 = (*trk)->charge();
01241                 m_t01 = tes;
01242                 m_driftT = driftT;
01243                 m_tdc_raw = tdc_raw;
01244                 m_phraw = adc_raw;
01245                 m_exraw = ph_hit;
01246                 m_localwid = localwid;
01247                 m_wire = wid;
01248                 m_runNO1 = runNO;
01249                 m_nhit_hit = Nhits;
01250                 m_doca_in = dd;
01251                 m_doca_ex = dd;
01252                 m_driftdist = driftD;
01253                 m_eangle = eangle;
01254                 m_costheta1 = costheta;
01255                 m_pathL = pathlength;
01256                 m_layer = layid;
01257                 m_ptrk1 = m_P;
01258                 //cout<<"adc_raw : "<<adc_raw<<"     "<<ph_hit<<"     "<<dd_in<<"  "<<dd_ex<<"  "<<eangle<<endl;
01259                 m_trackId1 = trk_ex.trk_id();
01260                 StatusCode status2 = m_nt2->write();
01261                 if ( status2.isFailure() )
01262                   log << MSG::ERROR << "Cannot fill Ntuple n102" << endreq;
01263             }
01264             if(layid>3)   
01265             {
01266                 phlist.push_back(ph);
01267                 phlist_hit.push_back(ph_hit);
01268             }
01269         }
01270         dedxhitList->push_back( dedxhit );
01271         dedxhitrefvec->push_back( dedxhit );
01272     }//end of hit loop
01273     trk_ex.set_phlist( phlist );
01274     trk_ex.set_phlist_hit( phlist_hit );
01275     trk_ex.setVecDedxHits( *dedxhitrefvec );
01276     ex_trks.push_back(trk_ex );
01277     m_trkalgs.push_back(m_trkalg);
01278 
01279     delete dedxhitrefvec;
01280     phlist.clear();
01281     phlist_hit.clear();
01282     if(ntpFlag>0) trackNO2++;
01283 }

void MdcDedxRecon::switchtomdctrack	(	int	trkid,
		Identifier	mdcid,
		double	tes,
		int	RunNO,
		int	runFlag,
		MsgStream	log
	)

Definition at line 1481 of file MdcDedxRecon.cxx.

References mdctrackrec(), and Bes_Common::WARNING.

Referenced by execute().

01482 {
01483     //retrieve  MdcTrackCol from TDS
01484     SmartDataPtr<RecMdcTrackCol> newtrkCol(eventSvc(),"/Event/Recon/RecMdcTrackCol");
01485     if (!newtrkCol) 
01486     {
01487         log << MSG::WARNING << "Could not find RecMdcTrackCol in switchtomdctrack" << endreq;
01488         return ;
01489     }
01490 
01491     RecMdcTrackCol::iterator trk = newtrkCol->begin();
01492     for( ; trk != newtrkCol->end(); trk++) 
01493     {
01494         if( (*trk)->trackId()==trkid)
01495         {
01496             HitRefVec gothits= (*trk)->getVecHits();
01497             if(gothits.size()>=2)
01498               mdctrackrec(trk,mdcid,tes,runNO,runFlag,log);
01499         }
01500     }
01501 }

const std::vector< MdcDedxTrk > & MdcDedxRecon::tracks ( void ) const

Definition at line 1110 of file MdcDedxRecon.cxx.

References ex_trks.

01111 {
01112     return ex_trks;
01113 }

int MdcDedxRecon::calib_flag [private]

Definition at line 49 of file MdcDedxRecon.h.

Referenced by initialize(), and MdcDedxRecon().

IDedxCurSvc* MdcDedxRecon::dedxcursvc [private]

Definition at line 60 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

RecMdcDedxHit* MdcDedxRecon::dedxhit [private]

Definition at line 71 of file MdcDedxRecon.h.

Referenced by execute(), kaltrackrec(), and mdctrackrec().

RecMdcDedxHitCol* MdcDedxRecon::dedxhitList [private]

Definition at line 70 of file MdcDedxRecon.h.

Referenced by execute(), kaltrackrec(), and mdctrackrec().

DedxHitRefVec* MdcDedxRecon::dedxhitrefvec [private]

Definition at line 72 of file MdcDedxRecon.h.

Referenced by execute(), kaltrackrec(), and mdctrackrec().

RecMdcDedxCol* MdcDedxRecon::dedxList [private]

Definition at line 68 of file MdcDedxRecon.h.

Referenced by execute().

int MdcDedxRecon::doNewGlobal [private]

Definition at line 52 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), and MdcDedxRecon().

MdcDedxCorrection* MdcDedxRecon::ex_calib [private]

Definition at line 63 of file MdcDedxRecon.h.

Referenced by execute(), finalize(), initialize(), and MdcDedxRecon().

std::vector<MdcDedxTrk> MdcDedxRecon::ex_trks [private]

Definition at line 64 of file MdcDedxRecon.h.

Referenced by execute(), finalize(), kaltrackrec(), mdctrackrec(), and tracks().

IDedxCorrecSvc* MdcDedxRecon::exsvc [private]

Definition at line 62 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

IMdcGeomSvc* MdcDedxRecon::geosvc [private]

Definition at line 61 of file MdcDedxRecon.h.

Referenced by beginRun(), execute(), kaltrackrec(), and mdctrackrec().

int MdcDedxRecon::landau [private]

Definition at line 50 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), and MdcDedxRecon().

int MdcDedxRecon::m_alg [private]

Definition at line 56 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), and MdcDedxRecon().

NTuple::Item<float> MdcDedxRecon::m_charge [private]

Definition at line 94 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<double> MdcDedxRecon::m_charge1 [private]

Definition at line 125 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_chidedx [private]

Definition at line 116 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_chidedxe [private]

Definition at line 117 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_chidedxk [private]

Definition at line 120 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_chidedxmu [private]

Definition at line 118 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_chidedxp [private]

Definition at line 121 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_chidedxpi [private]

Definition at line 119 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_costheta [private]

Definition at line 92 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_costheta1 [private]

Definition at line 136 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Array<double> MdcDedxRecon::m_dEdx_hit [private]

Definition at line 107 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_dEdx_meas [private]

Definition at line 89 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_doca_ex [private]

Definition at line 133 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_doca_in [private]

Definition at line 132 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_driftdist [private]

Definition at line 134 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_driftT [private]

Definition at line 143 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_eangle [private]

Definition at line 135 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_evtNO [private]

Definition at line 96 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_expect [private]

Definition at line 114 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Array<float> MdcDedxRecon::m_expectid [private]

Definition at line 110 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<double> MdcDedxRecon::m_exraw [private]

Definition at line 127 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

TFolder* MdcDedxRecon::m_hitlevel [private]

Definition at line 76 of file MdcDedxRecon.h.

Referenced by finalize(), and initialize().

TH1F* MdcDedxRecon::m_hitNo_wire [private]

Definition at line 77 of file MdcDedxRecon.h.

Referenced by execute(), finalize(), initialize(), kaltrackrec(), and mdctrackrec().

TObjArray* MdcDedxRecon::m_hlist [private]

Definition at line 75 of file MdcDedxRecon.h.

Referenced by finalize(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_layer [private]

Definition at line 138 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_localwid [private]

Definition at line 129 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_nhit [private]

Definition at line 104 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_nhit_hit [private]

Definition at line 131 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<long> MdcDedxRecon::m_nphlisthit [private]

Definition at line 106 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Tuple* MdcDedxRecon::m_nt1 [private]

Definition at line 80 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Tuple* MdcDedxRecon::m_nt2 [private]

Definition at line 124 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_parttype [private]

Definition at line 112 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_pathL [private]

Definition at line 137 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<double> MdcDedxRecon::m_phraw [private]

Definition at line 126 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_phtm [private]

Definition at line 81 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_poca_x [private]

Definition at line 99 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_poca_y [private]

Definition at line 100 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_poca_z [private]

Definition at line 101 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_prob [private]

Definition at line 113 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Array<float> MdcDedxRecon::m_probpid [private]

Definition at line 109 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_ptrk [private]

Definition at line 91 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_ptrk1 [private]

Definition at line 139 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_ptrk_hit [private]

Definition at line 140 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), and kaltrackrec().

NTuple::Item<float> MdcDedxRecon::m_recalg [private]

Definition at line 102 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

std::string MdcDedxRecon::m_rootfile [private]

Definition at line 58 of file MdcDedxRecon.h.

Referenced by execute(), finalize(), initialize(), kaltrackrec(), MdcDedxRecon(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_runNO [private]

Definition at line 95 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<double> MdcDedxRecon::m_runNO1 [private]

Definition at line 130 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_sigma [private]

Definition at line 115 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Array<float> MdcDedxRecon::m_sigmaid [private]

Definition at line 111 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_sintheta [private]

Definition at line 93 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_t0 [private]

Definition at line 97 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_t01 [private]

Definition at line 141 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_tdc_raw [private]

Definition at line 142 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_trackId [private]

Definition at line 98 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_trackId1 [private]

Definition at line 144 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

int MdcDedxRecon::m_trkalg [private]

Definition at line 66 of file MdcDedxRecon.h.

Referenced by execute(), kaltrackrec(), and mdctrackrec().

std::vector<int> MdcDedxRecon::m_trkalgs [private]

Definition at line 65 of file MdcDedxRecon.h.

Referenced by execute(), finalize(), kaltrackrec(), and mdctrackrec().

NTuple::Item<float> MdcDedxRecon::m_usedhit [private]

Definition at line 105 of file MdcDedxRecon.h.

Referenced by execute(), and initialize().

NTuple::Item<float> MdcDedxRecon::m_wire [private]

Definition at line 128 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), and mdctrackrec().

int MdcDedxRecon::ntpFlag [private]

Definition at line 51 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), kaltrackrec(), MdcDedxRecon(), and mdctrackrec().

int MdcDedxRecon::ParticleFlag [private]

Definition at line 55 of file MdcDedxRecon.h.

Referenced by execute(), kaltrackrec(), and MdcDedxRecon().

int MdcDedxRecon::recalg [private]

Definition at line 53 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), and MdcDedxRecon().

RecMdcDedx* MdcDedxRecon::resdedx [private]

Definition at line 69 of file MdcDedxRecon.h.

Referenced by execute().

float MdcDedxRecon::truncate [private]

Definition at line 57 of file MdcDedxRecon.h.

Referenced by execute(), initialize(), and MdcDedxRecon().

MdcDedxRecon Class Reference

Public Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation


Public Member Functions
	MdcDedxRecon (const std::string &name, ISvcLocator *pSvcLocator)
	~MdcDedxRecon ()
StatusCode	initialize ()
StatusCode	execute ()
StatusCode	finalize ()
StatusCode	beginRun ()
const std::vector< MdcDedxTrk > &	tracks (void) const
void	clearTables ()
void	mdctrackrec (RecMdcTrackCol::iterator trk, Identifier mdcid, double tes, int RunNO, int runFlag, MsgStream log)
void	kaltrackrec (RecMdcKalTrackCol::iterator trk, Identifier mdcid, double tes, int RunNO, int runFlag, MsgStream log)
void	switchtomdctrack (int trkid, Identifier mdcid, double tes, int RunNO, int runFlag, MsgStream log)
Private Attributes
int	calib_flag
int	landau
int	ntpFlag
int	doNewGlobal
int	recalg
int	ParticleFlag
int	m_alg
float	truncate
std::string	m_rootfile
IDedxCurSvc *	dedxcursvc
IMdcGeomSvc *	geosvc
IDedxCorrecSvc *	exsvc
MdcDedxCorrection *	ex_calib
std::vector< MdcDedxTrk >	ex_trks
std::vector< int >	m_trkalgs
int	m_trkalg
RecMdcDedxCol *	dedxList
RecMdcDedx *	resdedx
RecMdcDedxHitCol *	dedxhitList
RecMdcDedxHit *	dedxhit
DedxHitRefVec *	dedxhitrefvec
TObjArray *	m_hlist
TFolder *	m_hitlevel
TH1F *	m_hitNo_wire
NTuple::Tuple *	m_nt1
NTuple::Item< float >	m_phtm
NTuple::Item< float >	m_dEdx_meas
NTuple::Item< float >	m_ptrk
NTuple::Item< float >	m_costheta
NTuple::Item< float >	m_sintheta
NTuple::Item< float >	m_charge
NTuple::Item< float >	m_runNO
NTuple::Item< float >	m_evtNO
NTuple::Item< float >	m_t0
NTuple::Item< float >	m_trackId
NTuple::Item< float >	m_poca_x
NTuple::Item< float >	m_poca_y
NTuple::Item< float >	m_poca_z
NTuple::Item< float >	m_recalg
NTuple::Item< float >	m_nhit
NTuple::Item< float >	m_usedhit
NTuple::Item< long >	m_nphlisthit
NTuple::Array< double >	m_dEdx_hit
NTuple::Array< float >	m_probpid
NTuple::Array< float >	m_expectid
NTuple::Array< float >	m_sigmaid
NTuple::Item< float >	m_parttype
NTuple::Item< float >	m_prob
NTuple::Item< float >	m_expect
NTuple::Item< float >	m_sigma
NTuple::Item< float >	m_chidedx
NTuple::Item< float >	m_chidedxe
NTuple::Item< float >	m_chidedxmu
NTuple::Item< float >	m_chidedxpi
NTuple::Item< float >	m_chidedxk
NTuple::Item< float >	m_chidedxp
NTuple::Tuple *	m_nt2
NTuple::Item< double >	m_charge1
NTuple::Item< double >	m_phraw
NTuple::Item< double >	m_exraw
NTuple::Item< float >	m_wire
NTuple::Item< float >	m_localwid
NTuple::Item< double >	m_runNO1
NTuple::Item< float >	m_nhit_hit
NTuple::Item< float >	m_doca_in
NTuple::Item< float >	m_doca_ex
NTuple::Item< float >	m_driftdist
NTuple::Item< float >	m_eangle
NTuple::Item< float >	m_costheta1
NTuple::Item< float >	m_pathL
NTuple::Item< float >	m_layer
NTuple::Item< float >	m_ptrk1
NTuple::Item< float >	m_ptrk_hit
NTuple::Item< float >	m_t01
NTuple::Item< float >	m_tdc_raw
NTuple::Item< float >	m_driftT
NTuple::Item< float >	m_trackId1