BOSS 7.0.7
BESIII Offline Software System
Loading...
Searching...
No Matches
Gam4pikp Class Reference

#include <Gam4pikp.h>

+ Inheritance diagram for Gam4pikp:

Public Member Functions

 Gam4pikp (const std::string &name, ISvcLocator *pSvcLocator)
 
StatusCode initialize ()
 
StatusCode execute ()
 
StatusCode finalize ()
 

Detailed Description

Definition at line 9 of file Gam4pikp.h.

Constructor & Destructor Documentation

◆ Gam4pikp()

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

Definition at line 61 of file Gam4pikp.cxx.

61 :
62 Algorithm(name, pSvcLocator) {
63
64 //Declare the properties
65 declareProperty("skim4pi", m_skim4pi=false);
66 declareProperty("skim4k", m_skim4k=false);
67 declareProperty("skim2pi2pr", m_skim2pi2pr=false);
68 declareProperty("rootput", m_rootput=false);
69 declareProperty("Ecms", m_ecms=3.6862);
70 declareProperty("Vr0cut", m_vr0cut=1.0);
71 declareProperty("Vz0cut", m_vz0cut=5.0);
72 declareProperty("EnergyThreshold", m_energyThreshold=0.05);
73 declareProperty("GammaPhiCut", m_gammaPhiCut=20.0);
74 declareProperty("GammaThetaCut", m_gammaThetaCut=20.0);
75 declareProperty("GammaDangCut", m_gammadangcut=20.0);
76// declareProperty("Test4C", m_test4C = 1);
77// declareProperty("Test5C", m_test5C = 1);
78// declareProperty("CheckDedx", m_checkDedx = 1);
79// declareProperty("CheckTof", m_checkTof = 1);
80 }

Member Function Documentation

◆ execute()

StatusCode Gam4pikp::execute ( )

Definition at line 393 of file Gam4pikp.cxx.

393 {
394
395// std::cout << "execute()" << std::endl;
396 setFilterPassed(false);
397
398
399 MsgStream log(msgSvc(), name());
400 log << MSG::INFO << "in execute()" << endreq;
401
402 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
403 int run=eventHeader->runNumber();
404 int event=eventHeader->eventNumber();
405 log << MSG::DEBUG <<"run, evtnum = "
406 << run << " , "
407 << event <<endreq;
408// cout<<"run "<<run<<endl;
409// cout<<"event "<<event<<endl;
410if(m_rootput)
411{
412 m_run=run;
413 m_rec=event;
414}
415 Ncut[0]++;
416 SmartDataPtr<EvtRecEvent> evtRecEvent(eventSvc(), EventModel::EvtRec::EvtRecEvent);
417 log << MSG::DEBUG <<"ncharg, nneu, tottks = "
418 << evtRecEvent->totalCharged() << " , "
419 << evtRecEvent->totalNeutral() << " , "
420 << evtRecEvent->totalTracks() <<endreq;
421
422 SmartDataPtr<EvtRecTrackCol> evtRecTrkCol(eventSvc(), EventModel::EvtRec::EvtRecTrackCol);
423
424 if(m_rootput)
425{
426 Gam4pikp::InitVar();
427}
428
429 // for MC topo analysis
430if(m_rootput)
431{
432 if (eventHeader->runNumber()<0)
433 {
434 SmartDataPtr<Event::McParticleCol> mcParticleCol(eventSvc(), "/Event/MC/McParticleCol");
435 int m_numParticle = 0;
436 if (!mcParticleCol)
437 {
438// std::cout << "Could not retrieve McParticelCol" << std::endl;
439 return StatusCode::FAILURE;
440 }
441 else
442 {
443 bool psipDecay = false;
444 int rootIndex = -1;
445 Event::McParticleCol::iterator iter_mc = mcParticleCol->begin();
446 for (; iter_mc != mcParticleCol->end(); iter_mc++)
447 {
448 if ((*iter_mc)->primaryParticle()) continue;
449 if (!(*iter_mc)->decayFromGenerator()) continue;
450 if ((*iter_mc)->particleProperty()==100443)
451 {
452 psipDecay = true;
453 rootIndex = (*iter_mc)->trackIndex();
454 }
455 if (!psipDecay) continue;
456 int mcidx = ((*iter_mc)->mother()).trackIndex() - rootIndex;
457 int pdgid = (*iter_mc)->particleProperty();
458 m_pdgid[m_numParticle] = pdgid;
459 m_motheridx[m_numParticle] = mcidx;
460 m_numParticle += 1;
461 }
462 }
463 m_idxmc = m_numParticle;
464 }
465}
466
467 Vint iGood, ipip, ipim;
468 iGood.clear();
469 ipip.clear();
470 ipim.clear();
471 Vp4 ppip, ppim;
472 ppip.clear();
473 ppim.clear();
474 int nCharge = 0;
475 int ngdcgx=0;
476 int ncgx=0;
477 Vdouble pTrkCh;
478 pTrkCh.clear();
479 Hep3Vector v(0,0,0);
480 Hep3Vector vv(0,0,0);
481
482 IVertexDbSvc* vtxsvc;
483 Gaudi::svcLocator()->service("VertexDbSvc", vtxsvc);
484 if(vtxsvc->isVertexValid())
485 {
486 double* vertex = vtxsvc->PrimaryVertex();
487 double* vertexsigma = vtxsvc->SigmaPrimaryVertex();
488 v.setX(vertex[0]);
489 v.setY(vertex[1]);
490 v.setZ(vertex[2]);
491 vv.setX(vertexsigma[0]);
492 vv.setY(vertexsigma[1]);
493 vv.setZ(vertexsigma[2]);
494 }
495
496
497 if(run<0)
498 {
499 v[0]=0.;
500 v[1]=0.;
501 v[2]=0.;
502 vv[0]=0.;
503 vv[1]=0.;
504 vv[2]=0.;
505
506 }
507
508
509 for(int i = 0; i < evtRecEvent->totalCharged(); i++){
510 if(i>=evtRecTrkCol->size()) break;
511 ncgx++;
512 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + i;
513 if(!(*itTrk)->isMdcTrackValid()) continue;
514 if (!(*itTrk)->isMdcKalTrackValid()) continue;
515 ngdcgx++;
516 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
517
518 HepVector a = mdcTrk->helix();
519 HepSymMatrix Ea = mdcTrk->err();
520 HepPoint3D pivot(0.,0.,0.);
521 HepPoint3D IP(v[0],v[1],v[2]);
522 VFHelix helixp(pivot,a,Ea);
523 helixp.pivot(IP);
524 HepVector vec = helixp.a();
525 pTrkCh.push_back(mdcTrk->p()*mdcTrk->charge());
526 iGood.push_back((*itTrk)->trackId());
527 nCharge += mdcTrk->charge();
528 }
529
530
531 int nGood = iGood.size();
532 log << MSG::DEBUG << "ngood, totcharge = " << nGood << " , " << nCharge << endreq;
533 if((nGood < 4))
534 {
535 return StatusCode::SUCCESS;
536 }
537
538 Ncut[1]++;
539 Vint iGam;
540 int ngamx=0;
541 int ngdgamx=0;
542 iGam.clear();
543 Vint iGamnofit;
544 iGamnofit.clear();
545 for(int i = evtRecEvent->totalCharged(); i< evtRecEvent->totalTracks(); i++) {
546 ngamx++;
547 if(i>=evtRecTrkCol->size()) break;
548 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + i;
549 if(!(*itTrk)->isEmcShowerValid()) continue;
550 RecEmcShower *emcTrk = (*itTrk)->emcShower();
551 Hep3Vector emcpos(emcTrk->x(), emcTrk->y(), emcTrk->z());
552 // find the nearest charged track
553// double dthe = 200.;
554// double dphi = 200.;
555// double dang = 200.;
556 double dthe = 200.;
557 double dphi = 200.;
558 double dang = 200.;
559
560 ngdgamx++;
561 for(int j = 0; j < evtRecEvent->totalCharged(); j++) {
562 if(j>=evtRecTrkCol->size()) break;
563 EvtRecTrackIterator jtTrk = evtRecTrkCol->begin() + j;
564 if(!(*jtTrk)->isExtTrackValid()) continue;
565 RecExtTrack *extTrk = (*jtTrk)->extTrack();
566 if(extTrk->emcVolumeNumber() == -1) continue;
567 Hep3Vector extpos = extTrk->emcPosition();
568 // double ctht = extpos.cosTheta(emcpos);
569 double angd = extpos.angle(emcpos);
570 double thed = extpos.theta() - emcpos.theta();
571 double phid = extpos.deltaPhi(emcpos);
572 thed = fmod(thed+CLHEP::twopi+CLHEP::twopi+pi, CLHEP::twopi) - CLHEP::pi;
573 phid = fmod(phid+CLHEP::twopi+CLHEP::twopi+pi, CLHEP::twopi) - CLHEP::pi;
574
575// if(fabs(thed) < fabs(dthe)) dthe = thed;
576// if(fabs(phid) < fabs(dphi)) dphi = phid;
577// if(angd < dang) dang = angd;
578
579 if(angd < dang){ dang=angd;
580 dthe=thed;
581 dphi=phid;
582 }
583
584 }
585 // if(dang>=200) continue;
586 double eraw = emcTrk->energy();
587 dthe = dthe * 180 / (CLHEP::pi);
588 dphi = dphi * 180 / (CLHEP::pi);
589 dang = dang * 180 / (CLHEP::pi);
590// dthert = dthert * 180 / (CLHEP::pi);
591// dphirt = dphirt * 180 / (CLHEP::pi);
592// dangrt = dangrt * 180 / (CLHEP::pi);
593 // if(eraw < m_energyThreshold) continue;
594 iGam.push_back((*itTrk)->trackId());
595 if(eraw < m_energyThreshold) continue;
596 if(dang < 20.0) continue;
597 iGamnofit.push_back((*itTrk)->trackId());
598
599 }
600 // Finish Good Photon Selection
601 //
602 int nGam = iGam.size();
603 int nGamnofit = iGamnofit.size();
604
605 log << MSG::DEBUG << "num Good Photon " << nGam << " , " <<evtRecEvent->totalNeutral()<<endreq;
606 if(nGam<1){
607 return StatusCode::SUCCESS;
608 }
609 Ncut[2]++;
610
611 if(nGood>20||nGam>30) return StatusCode::SUCCESS;
612
613if(m_rootput)
614{
615m_idxmdc = nGood;
616m_idxemc=nGam;
617}
618//
619 Vp4 pGam;
620 pGam.clear();
621 for(int i = 0; i < nGam; i++) {
622 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + iGam[i];
623 RecEmcShower* emcTrk = (*itTrk)->emcShower();
624 double eraw = emcTrk->energy();
625 double phi = emcTrk->phi();
626 double the = emcTrk->theta();
627 HepLorentzVector ptrk;
628 ptrk.setPx(eraw*sin(the)*cos(phi));
629 ptrk.setPy(eraw*sin(the)*sin(phi));
630 ptrk.setPz(eraw*cos(the));
631 ptrk.setE(eraw);
632
633 // ptrk = ptrk.boost(-0.011,0,0);// boost to cms
634
635 pGam.push_back(ptrk);
636 }
637 //
638 // Assign 4-momentum to each charged track
639 //
641
642
643 for(int i = 0; i < nGood; i++) {
644 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + iGood[i];
645 RecMdcTrack* mdcTrk = (*itTrk)->mdcTrack();
646
647 RecMdcKalTrack* mdcKalTrk = (*itTrk)->mdcKalTrack();//After ParticleID, use RecMdcKalTrack substitute RecMdcTrack
648 RecMdcKalTrack::setPidType (RecMdcKalTrack::pion);//PID can set to electron, muon, pion, kaon and proton;The default setting is pion
649
650 if(mdcKalTrk->charge() >0) {
651 ipip.push_back(iGood[i]);
652 HepLorentzVector ptrk;
653 ptrk.setPx(mdcKalTrk->px());
654 ptrk.setPy(mdcKalTrk->py());
655 ptrk.setPz(mdcKalTrk->pz());
656 double p3 = ptrk.mag();
657 ptrk.setE(sqrt(p3*p3+mpi*mpi));
658
659
660 ppip.push_back(ptrk);
661 } else {
662 ipim.push_back(iGood[i]);
663 HepLorentzVector ptrk;
664 ptrk.setPx(mdcKalTrk->px());
665 ptrk.setPy(mdcKalTrk->py());
666 ptrk.setPz(mdcKalTrk->pz());
667 double p3 = ptrk.mag();
668 ptrk.setE(sqrt(p3*p3+mpi*mpi));
669
670 // ptrk = ptrk.boost(-0.011,0,0);//boost to cms
671
672 ppim.push_back(ptrk);
673 }
674 }
675 int npip = ipip.size();
676 int npim = ipim.size();
677 // if(npip*npim != 4) return SUCCESS;
678 if((npip < 2)||(npim < 2)) return SUCCESS;
679
680 Ncut[3]++;
681
682 //
683 // Loop each gamma pair, check ppi0 and pTot
684 //
685
686 int icgp1=-1;
687 int icgp2=-1;
688 int icgm1=-1;
689 int icgm2=-1;
690 int igam=-1;
691 double chisqtrk = 9999.;
692 WTrackParameter wcgp1;
693 WTrackParameter wcgp2;
694 WTrackParameter wcgm1;
695 WTrackParameter wcgm2;
696 int ipip1js=-1;
697 int ipip2js=-1;
698 int ipim1js=-1;
699 int ipim2js=-1;
700 double mcompall=9999;
701 double mppreclst=9999;
702 double meelst=9999;;
703 double mchic2kpilst=9999;
704 double chis4c2kpilst=9999;
705 int type2kpilst=0;
706 double dtpr2kpilst[4]={9999,9999,9999,9999};
707
708 double mc1=mpi,mc2=mpi,mc3=mpi,mc4=mpi;
709 double mchic01=0.0;
710 HepLorentzVector pgam=(0,0,0,0);
711 HepPoint3D xorigin(0.,0.,0.);
712 xorigin[0]=9999.;
713 xorigin[1]=9999.;
714 xorigin[2]=9999.;
715 HepSymMatrix xem(3,0);
716
717 int cl4pi=0;
718 int clajs=0;
719 HepLorentzVector p4psipjs(0.011*m_ecms, 0.0, 0.0, m_ecms);
720 double psipBetajs = (p4psipjs.vect()).mag()/(p4psipjs.e());
721 HepLorentzVector p4psip(0.011*m_ecms, 0.0, 0.0, m_ecms);
722 double psipBeta = (p4psip.vect()).mag()/(p4psip.e()); // beta = P/E
723
724 for(int ii = 0; ii < npip-1; ii++) {
725 RecMdcKalTrack *pip1Trk = (*(evtRecTrkCol->begin()+ipip[ii]))->mdcKalTrack();
726 for(int ij = ii+1; ij < npip; ij++) {
727 RecMdcKalTrack *pip2Trk = (*(evtRecTrkCol->begin()+ipip[ij]))->mdcKalTrack();
728 for(int ik = 0; ik < npim-1; ik++) {
729 RecMdcKalTrack *pim1Trk = (*(evtRecTrkCol->begin()+ipim[ik]))->mdcKalTrack();
730 for(int il = ik+1; il < npim; il++) {
731 RecMdcKalTrack *pim2Trk = (*(evtRecTrkCol->begin()+ipim[il]))->mdcKalTrack();
732 double squar[3]={9999.,9999.,9999.};
733 double squarkpi[6]={9999.,9999.,9999.,9999.,9999.,9999.};
734 WTrackParameter wvpip1Trk, wvpim1Trk, wvpip2Trk, wvpim2Trk;
735
736 Vint iGoodjs;
737 iGoodjs.clear();
738 Vdouble pTrkjs;
739 pTrkjs.clear();
740
741
742
743 pTrkjs.push_back(pip1Trk->p()*pip1Trk->charge());
744 pTrkjs.push_back(pip2Trk->p()*pip2Trk->charge());
745 pTrkjs.push_back(pim1Trk->p()*pim1Trk->charge());
746 pTrkjs.push_back(pim2Trk->p()*pim2Trk->charge());
747 iGoodjs.push_back(ipip[ii]);
748 iGoodjs.push_back(ipip[ij]);
749 iGoodjs.push_back(ipim[ik]);
750 iGoodjs.push_back(ipim[il]);
751
752 Gam4pikp::BubbleSort(pTrkjs, iGoodjs);
753 Vint jGoodjs;
754 jGoodjs.clear();
755 Vp4 p4chTrkjs;
756 p4chTrkjs.clear();
757 Vint i4cpip1js, i4cpip2js, i4cpim1js, i4cpim2js;
758 i4cpip1js.clear();
759 i4cpip2js.clear();
760 i4cpim1js.clear();
761 i4cpim2js.clear();
762 i4cpip1js.push_back(iGoodjs[2]);
763 i4cpip2js.push_back(iGoodjs[3]);
764 i4cpim1js.push_back(iGoodjs[1]);
765 i4cpim2js.push_back(iGoodjs[0]);
766 jGoodjs.push_back(i4cpip1js[0]);
767 jGoodjs.push_back(i4cpim1js[0]);
768 jGoodjs.push_back(i4cpip2js[0]);
769 jGoodjs.push_back(i4cpim2js[0]);
770
771 for (int i = 0; i < 4; i++)
772 {
773 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + jGoodjs[i];
774 if (!(*itTrk)->isMdcTrackValid()) continue;
775 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
776 if (!(*itTrk)->isMdcKalTrackValid()) continue;
777 RecMdcKalTrack *mdcKalTrk = (*itTrk)->mdcKalTrack();
778 double ptrk;
779 HepLorentzVector p4trk;
780 if (i<2)
781 {
783 ptrk = mdcKalTrk->p();
784 p4trk.setPx(mdcKalTrk->px());
785 p4trk.setPy(mdcKalTrk->py());
786 p4trk.setPz(mdcKalTrk->pz());
787 p4trk.setE(sqrt(ptrk*ptrk+mpi*mpi));
788 }
789 else{
791 ptrk = mdcKalTrk->p();
792 p4trk.setPx(mdcKalTrk->px());
793 p4trk.setPy(mdcKalTrk->py());
794 p4trk.setPz(mdcKalTrk->pz());
795 p4trk.setE(sqrt(ptrk*ptrk+xmass[0]*xmass[0]));
796 }
797 p4trk.boost(-1.0*psipBetajs, 0.0, 0.0);
798 p4chTrkjs.push_back(p4trk);
799 }
800 p4psipjs.boost(-1.0*psipBetajs, 0.0, 0.0);
801
802 HepLorentzVector p4pipijs = p4chTrkjs[0] + p4chTrkjs[1];
803 HepLorentzVector p4eejs = p4chTrkjs[2] + p4chTrkjs[3];
804 HepLorentzVector p4pipirecjs = p4psipjs - p4pipijs;
805 double mpprecjs=p4pipirecjs.m();
806 double mpipijs=p4pipijs.m();
807 double meejs=p4eejs.m();
808 double mcomp=sqrt((mpprecjs-3.097)*(mpprecjs-3.097)+(meejs-3.097)*(meejs-3.097));
809 if(mcomp<mcompall)
810 {
811 mcompall=mcomp;
812 ipip1js=i4cpip1js[0];
813 ipim1js=i4cpim1js[0];
814 ipip2js=i4cpip2js[0];
815 ipim2js=i4cpim2js[0];
816 mppreclst=mpprecjs;
817 meelst=meejs;
818 }
819
820if(m_rootput)
821{
822 m_mpprecall=mppreclst;
823 m_meeall=meelst;
824}
825 }
826
827 }
828 }
829 }
830
831 //{
832// HepLorentzVector p4psip(0.011*m_ecms, 0.0, 0.0, m_ecms);
833// double psipBeta = (p4psip.vect()).mag()/(p4psip.e()); // beta = P/E
834 Vint iGood4c;
835 iGood4c.clear();
836 Vdouble pTrk4c;
837 pTrk4c.clear();
838 Vint jGood;
839 jGood.clear();
840 Vint jsGood;
841 jsGood.clear();
842
843 if(mcompall>9997)
844 return StatusCode::SUCCESS;
845
846 jsGood.push_back(ipip1js);
847 jsGood.push_back(ipim1js);
848 jsGood.push_back(ipip2js);
849 jsGood.push_back(ipim2js);
850
851 for(int i = 0; i < evtRecEvent->totalCharged(); i++)
852 {
853 if(i>=evtRecTrkCol->size()) break;
854 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + i;
855 if(!(*itTrk)->isMdcTrackValid()) continue;
856 if (!(*itTrk)->isMdcKalTrackValid()) continue;
857 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
858 if((i!=ipip1js)&&(i!=ipim1js)&&(i!=ipip2js)&&(i!=ipim2js))
859 {
860 jsGood.push_back(i);
861 }
862 }
863
864 int njsGood=jsGood.size();
865 //ParticleID *pid = ParticleID::instance();
866if(m_rootput)
867{
868 for (int i = 0; i < njsGood; i++)
869 {
870 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + jsGood[i];
871 if (!(*itTrk)->isMdcTrackValid()) continue;
872 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
873 if (!(*itTrk)->isMdcKalTrackValid()) continue;
874 RecMdcKalTrack *mdcKalTrk = (*itTrk)->mdcKalTrack();
875 double ptrk;
876 ptrk = mdcKalTrk->p();
877 m_x0js[i] = mdcTrk->x();
878 m_y0js[i] = mdcTrk->y();
879 m_z0js[i] = mdcTrk->z();
880 m_r0js[i] = mdcTrk->r();
881 m_ppmdcjs[i] = mdcTrk->p();
882 m_pxmdcjs[i] = mdcTrk->px();
883 m_pymdcjs[i] = mdcTrk->py();
884 m_pzmdcjs[i] = mdcTrk->pz();
885 m_ppkaljs[i] = mdcKalTrk->p();
886 Hep3Vector p3jsi(mdcTrk->px(),mdcTrk->py(),mdcTrk->pz());
887 if(njsGood>4){
888 EvtRecTrackIterator itTrk5 = evtRecTrkCol->begin() + jsGood[4];
889 RecMdcTrack *mdcTrk5 = (*itTrk5)->mdcTrack();
890 Hep3Vector p3js5(mdcTrk5->px(),mdcTrk5->py(),mdcTrk5->pz());
891 m_angjs5[i]=p3jsi.angle(p3js5);
892 m_nearjs5[i]=p3jsi.howNear(p3js5);
893 }
894 if(njsGood>5){
895 EvtRecTrackIterator itTrk6 = evtRecTrkCol->begin() + jsGood[5];
896 RecMdcTrack *mdcTrk6 = (*itTrk6)->mdcTrack();
897 Hep3Vector p3js6(mdcTrk6->px(),mdcTrk6->py(),mdcTrk6->pz());
898 m_angjs6[i]=p3jsi.angle(p3js6);
899 m_nearjs6[i]=p3jsi.howNear(p3js6);
900 }
901
902
903 m_ptmdcjs[i] = mdcTrk->pxy();
904 m_ptkaljs[i] = mdcKalTrk->pxy();
905 double x0=mdcTrk->x();
906 double y0=mdcTrk->y();
907 double z0=mdcTrk->z();
908 double phi0=mdcTrk->helix(1);
909 double xv=v[0];
910 double yv=v[1];
911 double zv=v[2];
912 double Rxy=(x0-xv)*cos(phi0)+(y0-yv)*sin(phi0);
913 double Rz=z0-zv;
914 m_Rxyjs[i]=Rxy;
915 m_Rzjs[i]=Rz;
916 HepVector a = mdcTrk->helix();
917 HepSymMatrix Ea = mdcTrk->err();
918 HepPoint3D pivot(0.,0.,0.);
919 HepPoint3D IP(v[0],v[1],v[2]);
920 VFHelix helixp(pivot,a,Ea);
921 helixp.pivot(IP);
922 HepVector vec = helixp.a();
923 m_Rnxyjs[i]=vec[0];
924 m_Rnzjs[i]=vec[3];
925 m_phinjs[i]=vec[1];
926 //tof
927 if ((*itTrk)->isTofTrackValid())
928 { m_bjtofvaljs[i]=1;
929 m_cotof1js[i]=1;
930 SmartRefVector<RecTofTrack> tofTrkCol = (*itTrk)->tofTrack();
931 SmartRefVector<RecTofTrack>::iterator iter_tof=tofTrkCol.begin();
932 for (;iter_tof != tofTrkCol.end(); iter_tof++ ){
933 TofHitStatus *status = new TofHitStatus;
934 status->setStatus((*iter_tof)->status());
935 if(status->is_cluster()){
936 m_bjtofjs[i]=1;
937 m_counterjs[i] = status->is_counter();
938 m_barreljs[i] = status->is_barrel();
939 m_layertofjs[i] = status->layer();
940 m_readoutjs[i] = status->is_readout();
941 m_clusterjs[i] = status->is_cluster();
942 m_cotof2js[i]=2;
943 m_betajs[i]= (*iter_tof)->beta();
944 m_tofjs[i] = (*iter_tof)->tof();
945 m_tofpathjs[i] = (*iter_tof)->path();
946 m_zhitjs[i]= (*iter_tof)->zrhit();
947 m_texejs[i] = (*iter_tof)->texpElectron();
948 m_texmujs[i] = (*iter_tof)->texpMuon();
949 m_texpijs[i] = (*iter_tof)->texpPion();
950 m_texkjs[i] = (*iter_tof)->texpKaon();
951 m_texprjs[i] = (*iter_tof)->texpProton();
952 m_dtejs[i] = m_tofjs[i] - m_texejs[i];
953 m_dtmujs[i] = m_tofjs[i] - m_texmujs[i];
954 m_dtpijs[i] = m_tofjs[i] - m_texpijs[i];
955 m_dtkjs[i] = m_tofjs[i] - m_texkjs[i];
956 m_dtprjs[i] = m_tofjs[i] - m_texprjs[i];
957
958 m_sigmaetofjs[i]=(*iter_tof)->sigma(0);
959 m_sigmamutofjs[i]=(*iter_tof)->sigma(1);
960 m_sigmapitofjs[i]=(*iter_tof)->sigma(2);
961 m_sigmaktofjs[i]=(*iter_tof)->sigma(3);
962 m_sigmaprtofjs[i]=(*iter_tof)->sigma(4);
963 m_t0tofjs[i]=(*iter_tof)->t0();
964 m_errt0tofjs[i]=(*iter_tof)->errt0();
965
966 m_tofIDjs[i]=(*iter_tof)->tofID();
967 m_clusterIDjs[i]=(*iter_tof)->tofTrackID();
968 }
969 delete status;
970 }
971 }
972 //dedx
973 if ((*itTrk)->isMdcDedxValid())
974 {
975
976 RecMdcDedx* dedxTrk = (*itTrk)->mdcDedx();
977 m_chiejs[i] = dedxTrk->chiE();
978 m_chimujs[i] = dedxTrk->chiMu();
979 m_chipijs[i] = dedxTrk->chiPi();
980 m_chikjs[i] = dedxTrk->chiK();
981 m_chipjs[i] = dedxTrk->chiP();
982 m_ghitjs[i] = dedxTrk->numGoodHits();
983 m_thitjs[i] = dedxTrk->numTotalHits();
984 m_probphjs[i] = dedxTrk->probPH();
985 m_normphjs[i] = dedxTrk->normPH();
986 }
987
988
989
990 //emc
991 if ( (*itTrk)->isEmcShowerValid() )
992 {
993 RecEmcShower* emcTrk = (*itTrk)->emcShower();
994 m_bjemcjs[i]=1;
995 m_emcjs[i] = emcTrk->energy();
996 m_evpjs[i] = emcTrk->energy()/ptrk;
997 m_timecgjs[i] = emcTrk->time();
998 // totalEnergy += emcTrk->energy();
999 }
1000
1001 //muc
1002 if ( (*itTrk)->isMucTrackValid() )
1003 {
1004 RecMucTrack* mucTrk = (*itTrk)->mucTrack();
1005 double dpthp = mucTrk->depth()/25.0;//why?
1006 // m_depthmucjs[i] = dpthp<10.0 ? dpthp : 10.0;
1007 m_bjmucjs[i]=1;
1008 m_depthmucjs[i]=mucTrk->depth();
1009 m_layermucjs[i] = mucTrk->numLayers();
1010 }
1011
1012 m_cosmdcjs[i] = cos(mdcTrk->theta());
1013 m_phimdcjs[i] = mdcTrk->phi();
1014 //pid
1015 pid->init();
1016 pid->setMethod(pid->methodProbability());
1017 pid->setChiMinCut(4);
1018 pid->setRecTrack(*itTrk);
1019 pid->usePidSys(pid->useDedx() | pid->useTof1() | pid->useTof2());
1020 pid->identify(pid->onlyPion() | pid->onlyKaon() | pid->onlyProton());
1021 pid->identify(pid->onlyElectron() | pid->onlyMuon());
1022 pid->calculate();
1023 if(!(pid->IsPidInfoValid())) continue;
1024 // m_costpidjs[i] = cos(mdcTrk->theta());
1025 m_dedxpidjs[i] = pid->chiDedx(2);
1026 m_tof1pidjs[i] = pid->chiTof1(2);
1027 m_tof2pidjs[i] = pid->chiTof2(2);
1028 m_probejs[i] = pid->probElectron();
1029 m_probmujs[i] = pid->probMuon();
1030 m_probpijs[i] = pid->probPion();
1031 m_probkjs[i] = pid->probKaon();
1032 m_probprjs[i] = pid->probProton();
1033
1034
1035
1036 }
1037
1038}
1039
1040 Vint jGam2kpi;
1041 jGam2kpi.clear();
1042
1043 Vint iGood2kpi, ipip2kpi, ipim2kpi;
1044 iGood2kpi.clear();
1045 ipip2kpi.clear();
1046 ipim2kpi.clear();
1047 Vp4 ppip2kpi, ppim2kpi;
1048 ppip2kpi.clear();
1049 ppim2kpi.clear();
1050
1051 Vint ipipnofit, ipimnofit, ikpnofit, ikmnofit, ipropnofit, ipromnofit;
1052 ipipnofit.clear();
1053 ipimnofit.clear();
1054 ikpnofit.clear();
1055 ikmnofit.clear();
1056 ipropnofit.clear();
1057 ipromnofit.clear();
1058 Vp4 ppipnofit, ppimnofit, pkpnofit, pkmnofit, ppropnofit, ppromnofit;
1059 ppipnofit.clear();
1060 ppimnofit.clear();
1061 pkpnofit.clear();
1062 pkmnofit.clear();
1063 ppropnofit.clear();
1064 ppromnofit.clear();
1065
1066 Vdouble p3pip2kpi;
1067 p3pip2kpi.clear();
1068 Vdouble p3pim2kpi;
1069 p3pim2kpi.clear();
1070
1071 Vint itrak2kpi;
1072 itrak2kpi.clear();
1073 int cls2kpi;
1074 double chis2kpi=9999.;
1075 double m4piall=9999.;
1076 double m4kall=9999.;
1077 double mgam4piall=9999.;
1078 double mgam4kall=9999.;
1079 for(int i = 0; i < evtRecEvent->totalCharged(); i++){
1080 if(i>=evtRecTrkCol->size()) break;
1081 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + i;
1082 if(!(*itTrk)->isMdcTrackValid()) continue;
1083 if (!(*itTrk)->isMdcKalTrackValid()) continue;
1084 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
1085 double z02kpi = mdcTrk->z();
1086 double r02kpi = mdcTrk->r();
1087 HepVector a = mdcTrk->helix();
1088 HepSymMatrix Ea = mdcTrk->err();
1089 HepPoint3D pivot(0.,0.,0.);
1090 HepPoint3D IP(v[0],v[1],v[2]);
1091 VFHelix helixp(pivot,a,Ea);
1092 helixp.pivot(IP);
1093 HepVector vec = helixp.a();
1094 double Rnxy02kpi=vec[0];
1095 double Rnz02kpi=vec[3];
1096 // if(fabs(Rnxy02kpi>1.0)) continue;
1097 // if(fabs(Rnz02kpi>10.0)) continue;
1098 iGood2kpi.push_back((*itTrk)->trackId());
1099 }
1100 int nGood2kpi = iGood2kpi.size();
1101 if(nGood2kpi==4)
1102 {
1103 for(int i = 0; i < nGood2kpi; i++) {
1104 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + iGood2kpi[i];
1105 RecMdcTrack* mdcTrk = (*itTrk)->mdcTrack();
1106 RecMdcKalTrack* mdcKalTrk = (*itTrk)->mdcKalTrack();
1107 if(mdcKalTrk->charge() >0) {
1108 ipip2kpi.push_back(iGood2kpi[i]);
1109 p3pip2kpi.push_back(mdcKalTrk->p());
1110 HepLorentzVector ptrk;
1111 ptrk.setPx(mdcKalTrk->px());
1112 ptrk.setPy(mdcKalTrk->py());
1113 ptrk.setPz(mdcKalTrk->pz());
1114 double p3 = ptrk.mag();
1115 ptrk.setE(sqrt(p3*p3+mpi*mpi));
1116 ppip2kpi.push_back(ptrk);
1117 } else {
1118 ipim2kpi.push_back(iGood2kpi[i]);
1119 p3pim2kpi.push_back(mdcKalTrk->p());
1120 HepLorentzVector ptrk;
1121 ptrk.setPx(mdcKalTrk->px());
1122 ptrk.setPy(mdcKalTrk->py());
1123 ptrk.setPz(mdcKalTrk->pz());
1124 double p3 = ptrk.mag();
1125 ptrk.setE(sqrt(p3*p3+mpi*mpi));
1126 ppim2kpi.push_back(ptrk);
1127 }
1128 }
1129 int npip2kpi = ipip2kpi.size();
1130 int npim2kpi = ipim2kpi.size();
1131
1132
1133
1135
1136if(m_rootput)
1137{
1138 m_cy2kpi=6;
1139}
1140
1141 if((npip2kpi == 2)&&(npim2kpi == 2))
1142 {
1143 //if(nGamnofit>=1)
1144
1145
1146Ncut[4]++;
1147
1148 HepPoint3D xorigin2kpi(0.,0.,0.);
1149 xorigin2kpi[0]=9999.;
1150 xorigin2kpi[1]=9999.;
1151 xorigin2kpi[2]=9999.;
1152 HepSymMatrix xem2kpi(3,0);
1153
1154 Gam4pikp::BubbleSort(p3pip2kpi, ipip2kpi);
1155 Gam4pikp::BubbleSort(p3pim2kpi, ipim2kpi);
1156
1157 Mcut[1]++;
1158 RecMdcKalTrack *pip12kpiTrk = (*(evtRecTrkCol->begin()+ipip2kpi[0]))->mdcKalTrack();
1159 RecMdcKalTrack *pip22kpiTrk = (*(evtRecTrkCol->begin()+ipip2kpi[1]))->mdcKalTrack();
1160 RecMdcKalTrack *pim12kpiTrk = (*(evtRecTrkCol->begin()+ipim2kpi[0]))->mdcKalTrack();
1161 RecMdcKalTrack *pim22kpiTrk = (*(evtRecTrkCol->begin()+ipim2kpi[1]))->mdcKalTrack();
1162 double squar2kpi[10]={9999.,9999.,9999.,9999.,9999.,9999.,9999.,9999.,9999.,9999.};
1163 double mc12kpi,mc22kpi,mc32kpi,mc42kpi;
1164 // double chis2kpi=9999.;
1165 WTrackParameter wcgp12kpi;
1166 WTrackParameter wcgp22kpi;
1167 WTrackParameter wcgm12kpi;
1168 WTrackParameter wcgm22kpi;
1169 int icgp12kpi=-1;
1170 int icgp22kpi=-1;
1171 int icgm12kpi=-1;
1172 int icgm22kpi=-1;
1173 int igam2kpi=-1;
1174 double m2kpi=9999;
1175
1176 int n20=0;
1177 int n30=0;
1178 WTrackParameter wvpip12kpiTrk, wvpim12kpiTrk, wvpip22kpiTrk, wvpim22kpiTrk;
1179 for(int k=0;k<6;k++)
1180 {
1181 if(k==0){mc12kpi=mpi;mc22kpi=mpi;mc32kpi=mpi;mc42kpi=mpi;}
1182 if(k==1){mc12kpi=mk;mc22kpi=mk;mc32kpi=mk;mc42kpi=mk;}
1183 if(k==2){mc12kpi=mpi;mc22kpi=mpro;mc32kpi=mpi;mc42kpi=mpro;}
1184 if(k==3){mc12kpi=mpro;mc22kpi=mpi;mc32kpi=mpro;mc42kpi=mpi;}
1185 if(k==4){mc12kpi=mpi;mc22kpi=mpro;mc32kpi=mpro;mc42kpi=mpi; }
1186 if(k==5){mc12kpi=mpro;mc22kpi=mpi;mc32kpi=mpi;mc42kpi=mpro;}
1187
1188
1189
1190 wvpip12kpiTrk = WTrackParameter(mc12kpi, pip12kpiTrk->getZHelix(), pip12kpiTrk->getZError());
1191 wvpip22kpiTrk = WTrackParameter(mc22kpi, pip22kpiTrk->getZHelix(), pip22kpiTrk->getZError());
1192 wvpim12kpiTrk = WTrackParameter(mc32kpi, pim12kpiTrk->getZHelix(), pim12kpiTrk->getZError());
1193 wvpim22kpiTrk = WTrackParameter(mc42kpi, pim22kpiTrk->getZHelix(), pim22kpiTrk->getZError());
1194 HepPoint3D vx(0., 0., 0.);
1195 HepSymMatrix Evx(3, 0);
1196 double bx = 1E+6;
1197 double by = 1E+6;
1198 double bz = 1E+6;
1199 Evx[0][0] = bx*bx;
1200 Evx[1][1] = by*by;
1201 Evx[2][2] = bz*bz;
1202 VertexParameter vxpar;
1203 vxpar.setVx(vx);
1204 vxpar.setEvx(Evx);
1205 VertexFit* vtxfit = VertexFit::instance();
1206 vtxfit->init();
1207 vtxfit->AddTrack(0, wvpip12kpiTrk);
1208 vtxfit->AddTrack(1, wvpim12kpiTrk);
1209 vtxfit->AddTrack(2, wvpip22kpiTrk);
1210 vtxfit->AddTrack(3, wvpim22kpiTrk);
1211 vtxfit->AddVertex(0, vxpar,0, 1, 2, 3);
1212 if(!vtxfit->Fit(0)) continue;
1213 vtxfit->Swim(0);
1214 WTrackParameter wpip12kpi = vtxfit->wtrk(0);
1215 WTrackParameter wpim12kpi = vtxfit->wtrk(1);
1216 WTrackParameter wpip22kpi = vtxfit->wtrk(2);
1217 WTrackParameter wpim22kpi = vtxfit->wtrk(3);
1218 WTrackParameter wpip12kpiys = vtxfit->wtrk(0);
1219 WTrackParameter wpim12kpiys = vtxfit->wtrk(1);
1220 WTrackParameter wpip22kpiys = vtxfit->wtrk(2);
1221 WTrackParameter wpim22kpiys = vtxfit->wtrk(3);
1222 xorigin2kpi = vtxfit->vx(0);
1223 xem2kpi = vtxfit->Evx(0);
1225
1226 HepLorentzVector ecms(0.040547,0,0,3.68632);
1227
1228 double chisq = 9999.;
1229 int ig1 = -1;
1230 for(int i = 0; i < nGam; i++) {
1231 RecEmcShower *g1Trk = (*(evtRecTrkCol->begin()+iGam[i]))->emcShower();
1232 kmfit->init();
1233 kmfit->setBeamPosition(xorigin2kpi);
1234 kmfit->setVBeamPosition(xem2kpi);
1235 kmfit->AddTrack(0, wpip12kpi);
1236 kmfit->AddTrack(1, wpim12kpi);
1237 kmfit->AddTrack(2, wpip22kpi);
1238 kmfit->AddTrack(3, wpim22kpi);
1239 kmfit->AddTrack(4, 0.0, g1Trk);
1240 kmfit->AddFourMomentum(0, p4psip);
1241 // if(!kmfit->Fit(0)) continue;
1242 bool oksq = kmfit->Fit();
1243 if(oksq) {
1244 double chi2 = kmfit->chisq();
1245 if(chi2 < chisq) {
1246 chisq = chi2;
1247 squar2kpi[k]=chi2;
1248 ig1 = iGam[i];
1249
1250 }
1251 }
1252 }
1253 if(squar2kpi[k]<200&&squar2kpi[k]<chis2kpi)
1254 {
1255// m_comcs2kpi[k]=squar2kpi[k];
1256 chis2kpi=squar2kpi[k];
1257 if(squar2kpi[k]<20) n20=n20+1;
1258 if(squar2kpi[k]<30) n30=n30+1;
1259
1260 icgp12kpi=ipip2kpi[0];
1261 icgp22kpi=ipip2kpi[1];
1262 icgm12kpi=ipim2kpi[0];
1263 icgm22kpi=ipim2kpi[1];
1264 igam2kpi=ig1;
1265 wcgp12kpi=wpip12kpiys;
1266 wcgp22kpi=wpip22kpiys;
1267 wcgm12kpi=wpim12kpiys;
1268 wcgm22kpi=wpim22kpiys;
1269 cls2kpi=k;
1270
1271 if(k==0){
1272 itrak2kpi.push_back(ipip2kpi[0]); itrak2kpi.push_back(ipim2kpi[0]);
1273 itrak2kpi.push_back(ipip2kpi[1]); itrak2kpi.push_back(ipim2kpi[1]);
1274 }
1275
1276 if(k==1){
1277 itrak2kpi.push_back(ipip2kpi[0]); itrak2kpi.push_back(ipim2kpi[0]);
1278 itrak2kpi.push_back(ipip2kpi[1]); itrak2kpi.push_back(ipim2kpi[1]);
1279 }
1280
1281
1282 if(k==2){
1283 itrak2kpi.push_back(ipip2kpi[0]); itrak2kpi.push_back(ipim2kpi[0]);
1284 itrak2kpi.push_back(ipip2kpi[1]); itrak2kpi.push_back(ipim2kpi[1]);
1285 }
1286
1287 if(k==3){
1288 itrak2kpi.push_back(ipip2kpi[1]); itrak2kpi.push_back(ipim2kpi[1]);
1289 itrak2kpi.push_back(ipip2kpi[0]); itrak2kpi.push_back(ipim2kpi[0]);
1290 }
1291
1292 if(k==4){
1293 itrak2kpi.push_back(ipip2kpi[0]); itrak2kpi.push_back(ipim2kpi[1]);
1294 itrak2kpi.push_back(ipip2kpi[1]); itrak2kpi.push_back(ipim2kpi[0]);
1295 }
1296
1297 if(k==5){
1298 itrak2kpi.push_back(ipip2kpi[1]); itrak2kpi.push_back(ipim2kpi[0]);
1299 itrak2kpi.push_back(ipip2kpi[0]); itrak2kpi.push_back(ipim2kpi[1]);
1300 }
1301
1302
1303
1304
1305
1306 RecEmcShower *g1Trk = (*(evtRecTrkCol->begin()+igam2kpi))->emcShower();
1307 kmfit->init();
1308 kmfit->setBeamPosition(xorigin2kpi);
1309 kmfit->setVBeamPosition(xem2kpi);
1310 kmfit->AddTrack(0, wpip12kpi);
1311 kmfit->AddTrack(1, wpim12kpi);
1312 kmfit->AddTrack(2, wpip22kpi);
1313 kmfit->AddTrack(3, wpim22kpi);
1314 kmfit->AddTrack(4, 0.0, g1Trk);
1315 kmfit->AddFourMomentum(0, p4psip);
1316 bool oksq = kmfit->Fit();
1317 if(oksq){
1318 HepLorentzVector pchic2kpi = kmfit->pfit(0) + kmfit->pfit(1)+kmfit->pfit(2) + kmfit->pfit(3);
1319 HepLorentzVector ppsip2kpi = kmfit->pfit(0) + kmfit->pfit(1)+kmfit->pfit(2) + kmfit->pfit(3) + kmfit->pfit(4);
1320 mchic2kpilst=pchic2kpi.m();
1321 chis4c2kpilst=kmfit->chisq();
1322if(m_rootput)
1323{
1324 m_mchic2kpi = pchic2kpi.m();
1325 m_chis2kpi = kmfit->chisq();
1326 m_mpsip2kpi=ppsip2kpi.m();
1327}
1328
1329
1330
1331 }
1332
1333
1334 }
1335 }
1336
1337
1338 if(chis2kpi<200)
1339 {
1340Ncut[5]++;
1341 if(m_rootput)
1342 {
1343 m_ncy20=n20;
1344 m_ncy30=n30;
1345 m_cla2kpi=cls2kpi;
1346 }
1347 type2kpilst=cls2kpi;
1348
1349 Vp4 p2kpi;
1350 p2kpi.clear();
1351
1352 for (int i = 0; i < 4; i++)
1353 {
1354 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + itrak2kpi[i];
1355 if (!(*itTrk)->isMdcTrackValid()) continue;
1356 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
1357 if (!(*itTrk)->isMdcKalTrackValid()) continue;
1358 RecMdcKalTrack *mdcKalTrk = (*itTrk)->mdcKalTrack();
1359 double ptrk2kpi;
1360 HepLorentzVector p4trk2kpi;
1361 if(cls2kpi==1)
1362 {
1363 mdcKalTrk->setPidType(RecMdcKalTrack::kaon);
1364 ptrk2kpi = mdcKalTrk->p();
1365 p4trk2kpi.setPx(mdcKalTrk->px());
1366 p4trk2kpi.setPy(mdcKalTrk->py());
1367 p4trk2kpi.setPz(mdcKalTrk->pz());
1368 p4trk2kpi.setE(sqrt(ptrk2kpi*ptrk2kpi+mk*mk));
1369 p2kpi.push_back(p4trk2kpi);
1370 }
1371
1372 if(cls2kpi==2)
1373 {
1374 if (i<2)
1375 {
1376 mdcKalTrk->setPidType(RecMdcKalTrack::pion);
1377 ptrk2kpi = mdcKalTrk->p();
1378 p4trk2kpi.setPx(mdcKalTrk->px());
1379 p4trk2kpi.setPy(mdcKalTrk->py());
1380 p4trk2kpi.setPz(mdcKalTrk->pz());
1381 p4trk2kpi.setE(sqrt(ptrk2kpi*ptrk2kpi+mpi*mpi));
1382 p2kpi.push_back(p4trk2kpi);
1383
1384 }
1385 else{
1387 ptrk2kpi = mdcKalTrk->p();
1388 p4trk2kpi.setPx(mdcKalTrk->px());
1389 p4trk2kpi.setPy(mdcKalTrk->py());
1390 p4trk2kpi.setPz(mdcKalTrk->pz());
1391 p4trk2kpi.setE(sqrt(ptrk2kpi*ptrk2kpi+mpro*mpro));
1392 p2kpi.push_back(p4trk2kpi);
1393
1394 }
1395
1396 }
1397 if(cls2kpi!=1&&cls2kpi!=2)
1398 {
1399 mdcKalTrk->setPidType(RecMdcKalTrack::pion);
1400 ptrk2kpi = mdcKalTrk->p();
1401 p4trk2kpi.setPx(mdcKalTrk->px());
1402 p4trk2kpi.setPy(mdcKalTrk->py());
1403 p4trk2kpi.setPz(mdcKalTrk->pz());
1404 p4trk2kpi.setE(sqrt(ptrk2kpi*ptrk2kpi+mpi*mpi));
1405 p2kpi.push_back(p4trk2kpi);
1406 }
1407
1408
1409
1410
1411 }
1412
1413 for (int i = 0; i < 4; i++)
1414 {
1415 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + itrak2kpi[i];
1416 if (!(*itTrk)->isMdcTrackValid()) continue;
1417 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
1418 if (!(*itTrk)->isMdcKalTrackValid()) continue;
1419 RecMdcKalTrack *mdcKalTrk = (*itTrk)->mdcKalTrack();
1420 if ((*itTrk)->isTofTrackValid())
1421 {
1422 SmartRefVector<RecTofTrack> tofTrkCol = (*itTrk)->tofTrack();
1423 SmartRefVector<RecTofTrack>::iterator iter_tof=tofTrkCol.begin();
1424 for (;iter_tof != tofTrkCol.end(); iter_tof++ ){
1425 TofHitStatus *status = new TofHitStatus;
1426 status->setStatus((*iter_tof)->status());
1427 if(status->is_cluster()){ dtpr2kpilst[i]=((*iter_tof)->tof()-(*iter_tof)->texpProton());
1428 }
1429 delete status;
1430 }
1431 }
1432 }
1433
1434
1435 /*
1436 HepLorentzVector ecmsb(0.040547,0,0,3.68632);
1437 double chisq = 9999.;
1438 RecEmcShower *g1Trk = (*(evtRecTrkCol->begin()+igam2kpi))->emcShower();
1439 KalmanKinematicFit * kmfit = KalmanKinematicFit::instance();
1440 kmfit->init();
1441 kmfit->setBeamPosition(xorigin2kpi);
1442 kmfit->setVBeamPosition(xem2kpi);
1443 kmfit->AddTrack(0, wcgp12kpi);
1444 kmfit->AddTrack(1, wcgp22kpi);
1445 kmfit->AddTrack(2, wcgm12kpi);
1446 kmfit->AddTrack(3, wcgm22kpi);
1447 kmfit->AddTrack(4, 0.0, g1Trk);
1448 kmfit->AddFourMomentum(0, p4psip);
1449 bool oksq = kmfit->Fit();
1450 if(oksq) {
1451 Mcut[3]++;
1452 HepLorentzVector pchic4c2kpi = kmfit->pfit(0) + kmfit->pfit(1)+kmfit->pfit(2) + kmfit->pfit(3);
1453 HepLorentzVector ppsip4c2kpi = kmfit->pfit(0) + kmfit->pfit(1)+kmfit->pfit(2) + kmfit->pfit(3) + kmfit->pfit(4);
1454
1455 m_mchic4c2kpi = pchic4c2kpi.m();
1456// m2kpi=m_mchic4c2kpi;
1457 m_chis4c2kpi = kmfit->chisq();
1458 m_mpsip4c2kpi=ppsip4c2kpi.m();
1459
1460
1461
1462
1463 }
1464
1465*/
1466 Mcut[2]++;
1467if(m_rootput)
1468{
1469 for (int i = 0; i < 4; i++)
1470 {
1471 EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + itrak2kpi[i];
1472 if (!(*itTrk)->isMdcTrackValid()) continue;
1473 RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
1474 if (!(*itTrk)->isMdcKalTrackValid()) continue;
1475 RecMdcKalTrack *mdcKalTrk = (*itTrk)->mdcKalTrack();
1476 m_ppmdc2kpi[i]=mdcTrk->p();
1477 m_pxmdc2kpi[i]=mdcTrk->px();
1478 m_pymdc2kpi[i]=mdcTrk->py();
1479 m_pzmdc2kpi[i]=mdcTrk->pz();
1480 m_ppkal2kpi[i]=mdcKalTrk->p();
1481 m_charge2kpi[i]=mdcKalTrk->charge();
1482 double ptrk;
1483 ptrk=mdcKalTrk->p();
1484
1485 if (eventHeader->runNumber()<0)
1486 {double mcall=9999;
1487 SmartDataPtr<Event::McParticleCol> mcParticleCol(eventSvc(), "/Event/MC/McParticleCol");
1488 int m_numParticle = 0;
1489 if (!mcParticleCol)
1490 {
1491// std::cout << "Could not retrieve McParticelCol" << std::endl;
1492 return StatusCode::FAILURE;
1493 }
1494 else
1495 {
1496 Event::McParticleCol::iterator iter_mc = mcParticleCol->begin();
1497 for (; iter_mc != mcParticleCol->end(); iter_mc++)
1498 { double commc;
1499 int pdgcode = (*iter_mc)->particleProperty();
1500 float px,py,pz;
1501 px=(*iter_mc)->initialFourMomentum().x();
1502 py=(*iter_mc)->initialFourMomentum().y();
1503 pz=(*iter_mc)->initialFourMomentum().z();
1504
1505 commc=ptrk*ptrk-px*px-py*py-pz*pz;
1506 if(fabs(commc)<fabs(mcall))
1507 {
1508 mcall=commc;
1509 m_pdg[i]=pdgcode;
1510 m_cbmc[i]=commc;
1511 }
1512
1513 }
1514
1515 }
1516
1517 }
1518
1519 if ((*itTrk)->isTofTrackValid())
1520 { m_bjtofval2kpi[i]=1;
1521 m_cotof12kpi[i]=1;
1522 SmartRefVector<RecTofTrack> tofTrkCol = (*itTrk)->tofTrack();
1523 SmartRefVector<RecTofTrack>::iterator iter_tof=tofTrkCol.begin();
1524 for (;iter_tof != tofTrkCol.end(); iter_tof++ ){
1525 TofHitStatus *status = new TofHitStatus;
1526 status->setStatus((*iter_tof)->status());
1527
1528 if(status->is_cluster()){
1529
1530
1531 m_bjtof2kpi[i]=1;
1532 m_counter2kpi[i] = status->is_counter();
1533 m_barrel2kpi[i] = status->is_barrel();
1534 m_layertof2kpi[i] = status->layer();
1535 m_readout2kpi[i] = status->is_readout();
1536 m_cluster2kpi[i] = status->is_cluster();
1537 m_cotof22kpi[i]=2;
1538 m_beta2kpi[i]= (*iter_tof)->beta();
1539 m_tof2kpi[i] = (*iter_tof)->tof();
1540 m_tofpath2kpi[i] = (*iter_tof)->path();
1541 m_zhit2kpi[i]= (*iter_tof)->zrhit();
1542 m_texe2kpi[i] = (*iter_tof)->texpElectron();
1543 m_texmu2kpi[i] = (*iter_tof)->texpMuon();
1544 m_texpi2kpi[i] = (*iter_tof)->texpPion();
1545 m_texk2kpi[i] = (*iter_tof)->texpKaon();
1546 m_texpr2kpi[i] = (*iter_tof)->texpProton();
1547 m_dte2kpi[i] = m_tof2kpi[i] - m_texe2kpi[i];
1548 m_dtmu2kpi[i] = m_tof2kpi[i] - m_texmu2kpi[i];
1549 m_dtpi2kpi[i] = m_tof2kpi[i] - m_texpi2kpi[i];
1550 m_dtk2kpi[i] = m_tof2kpi[i] - m_texk2kpi[i];
1551 m_dtpr2kpi[i] = m_tof2kpi[i] - m_texpr2kpi[i];
1552 m_tofID2kpi[i]=(*iter_tof)->tofID();
1553 m_clusterID2kpi[i]=(*iter_tof)->tofTrackID();
1554 m_sigmaetof2kpi[i]=(*iter_tof)->sigma(0);
1555 m_sigmamutof2kpi[i]=(*iter_tof)->sigma(1);
1556 m_sigmapitof2kpi[i]=(*iter_tof)->sigma(2);
1557 m_sigmaktof2kpi[i]=(*iter_tof)->sigma(3);
1558 m_sigmaprtof2kpi[i]=(*iter_tof)->sigma(4);
1559 m_t0tof2kpi[i]=(*iter_tof)->t0();
1560 m_errt0tof2kpi[i]=(*iter_tof)->errt0();
1561
1562 }
1563 delete status;
1564 }
1565 }
1566
1567
1568
1569
1570
1571
1572 if ((*itTrk)->isMdcDedxValid())
1573 {
1574 RecMdcDedx* dedxTrk = (*itTrk)->mdcDedx();
1575 m_chie2kpi[i] = dedxTrk->chiE();
1576 m_chimu2kpi[i] = dedxTrk->chiMu();
1577 m_chipi2kpi[i] = dedxTrk->chiPi();
1578 m_chik2kpi[i] = dedxTrk->chiK();
1579 m_chip2kpi[i] = dedxTrk->chiP();
1580 m_ghit2kpi[i] = dedxTrk->numGoodHits();
1581 m_thit2kpi[i] = dedxTrk->numTotalHits();
1582 m_probph2kpi[i] = dedxTrk->probPH();
1583 m_normph2kpi[i] = dedxTrk->normPH();
1584 }
1585
1586 //emc
1587 if ( (*itTrk)->isEmcShowerValid() )
1588 {
1589 RecEmcShower* emcTrk = (*itTrk)->emcShower();
1590 m_bjemc2kpi[i]=1;
1591 m_emc2kpi[i] = emcTrk->energy();
1592 m_evp2kpi[i] = emcTrk->energy()/ptrk;
1593 m_timecg2kpi[i] = emcTrk->time();
1594 }
1595
1596
1597 //muc
1598 if ( (*itTrk)->isMucTrackValid() )
1599 {
1600 RecMucTrack* mucTrk = (*itTrk)->mucTrack();
1601 double dpthp = mucTrk->depth()/25.0;//why?
1602 m_bjmuc2kpi[i]=1;
1603 m_depthmuc2kpi[i]=mucTrk->depth();
1604 m_layermuc2kpi[i] = mucTrk->numLayers();
1605 }
1606
1607 m_cosmdc2kpi[i] = cos(mdcTrk->theta());
1608 m_phimdc2kpi[i] = mdcTrk->phi();
1609
1610 m_pidnum2kpi[i]=(*itTrk)->partId();
1611
1612 if(m_skim4pi)
1613 {
1614 if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==0)
1615 {
1616 if(i<4) (*itTrk)->setPartId(3);
1617 }
1618
1619if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==1&&dtpr2kpilst[0]<-0.4&&dtpr2kpilst[1]<-0.4&&dtpr2kpilst[2]<-0.4&&dtpr2kpilst[3]<-0.4)
1620 {
1621 if(i<4) (*itTrk)->setPartId(4);
1622 }
1623
1624
1625if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==2)
1626 {
1627 if(i==0||i==1) (*itTrk)->setPartId(3);
1628 if(i==2||i==3) (*itTrk)->setPartId(5);
1629 }
1630
1631 }
1632 //pid
1634 pid->init();
1635 pid->setMethod(pid->methodProbability());
1636 pid->setChiMinCut(4);
1637 pid->setRecTrack(*itTrk);
1638 pid->usePidSys(pid->useDedx() | pid->useTof1() | pid->useTof2());
1639 pid->identify(pid->onlyPion() | pid->onlyKaon() | pid->onlyProton());
1640 pid->identify(pid->onlyElectron() | pid->onlyMuon());
1641 pid->calculate();
1642 if(!(pid->IsPidInfoValid())) continue;
1643 m_costpid2kpi[i] = cos(mdcTrk->theta());
1644
1645
1646 m_probe2kpi[i] = pid->probElectron();
1647 m_probmu2kpi[i] = pid->probMuon();
1648 m_probpi2kpi[i] = pid->probPion();
1649 m_probk2kpi[i] = pid->probKaon();
1650 m_probpr2kpi[i] = pid->probProton();
1651
1652
1653
1654
1655
1656
1657 }
1658}
1659
1660 }
1661// Ncut[10]++;
1662
1663if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==0)
1664{
1665 jGam2kpi.push_back(igam2kpi);
1666}
1667
1668 for(int i = evtRecEvent->totalCharged(); i< evtRecEvent->totalTracks(); i++) {
1669 if(i>=evtRecTrkCol->size()) break;
1670 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + i;
1671 if(!(*itTrk)->isEmcShowerValid()) continue;
1672if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==0)
1673 {
1674 if(i!=igam2kpi) jGam2kpi.push_back((*itTrk)->trackId());
1675 }
1676 else{
1677 jGam2kpi.push_back((*itTrk)->trackId());
1678 }
1679
1680 }
1681
1682int ngam2kpi=jGam2kpi.size();
1683
1684if(m_rootput)
1685{
1686for(int i = 0; i< ngam2kpi; i++) {
1687 if(i>=evtRecTrkCol->size()) break;
1688 EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + jGam2kpi[i];
1689 if(!(*itTrk)->isEmcShowerValid()) continue;
1690 RecEmcShower *emcTrk = (*(evtRecTrkCol->begin()+jGam2kpi[i]))->emcShower();
1691 Hep3Vector emcpos(emcTrk->x(), emcTrk->y(), emcTrk->z());
1692 double dthe = 200.;
1693 double dphi = 200.;
1694 double dang = 200.;
1695// double dthert = 200.;
1696// double dphirt = 200.;
1697// double dangrt = 200.;
1698 for(int j = 0; j < evtRecEvent->totalCharged(); j++) {
1699 if(j>=evtRecTrkCol->size()) break;
1700 EvtRecTrackIterator jtTrk = evtRecTrkCol->begin() + j;
1701 if(!(*jtTrk)->isExtTrackValid()) continue;
1702 RecExtTrack *extTrk = (*jtTrk)->extTrack();
1703 if(extTrk->emcVolumeNumber() == -1) continue;
1704 Hep3Vector extpos = extTrk->emcPosition();
1705 double angd = extpos.angle(emcpos);
1706 double thed = extpos.theta() - emcpos.theta();
1707 double phid = extpos.deltaPhi(emcpos);
1708 thed = fmod(thed+CLHEP::twopi+CLHEP::twopi+pi, CLHEP::twopi) - CLHEP::pi;
1709 phid = fmod(phid+CLHEP::twopi+CLHEP::twopi+pi, CLHEP::twopi) - CLHEP::pi;
1710// if(fabs(thed) < fabs(dthe)) dthe = thed;
1711// if(fabs(phid) < fabs(dphi)) dphi = phid;
1712// if(angd < dang) dang = angd;
1713 if(angd < dang) { dang = angd;
1714 dthe = thed;
1715 dphi = phid;
1716 }
1717
1718 }
1719 if(dang>=200) continue;
1720 double eraw = emcTrk->energy();
1721 dthe = dthe * 180 / (CLHEP::pi);
1722 dphi = dphi * 180 / (CLHEP::pi);
1723 dang = dang * 180 / (CLHEP::pi);
1724// dthert = dthert * 180 / (CLHEP::pi);
1725// dphirt = dphirt * 180 / (CLHEP::pi);
1726// dangrt = dangrt * 180 / (CLHEP::pi);
1727 m_numHits[i]= emcTrk->numHits();
1728 m_secondmoment[i] = emcTrk->secondMoment();
1729 m_latmoment[i] = emcTrk->latMoment();
1730 m_timegm[i] = emcTrk->time();
1731 m_cellId[i]=emcTrk->cellId();
1732 m_module[i]=emcTrk->module();
1733 m_a20Moment[i]=emcTrk->a20Moment();
1734 m_a42Moment[i]=emcTrk->a42Moment();
1735 m_getShowerId[i]=emcTrk->getShowerId();
1736 m_getClusterId[i]=emcTrk->getClusterId();
1737 m_getEAll[i]=emcTrk->getEAll();
1738 m_x[i]= emcTrk->x();
1739 m_y[i]= emcTrk->y();
1740 m_z[i]= emcTrk->z();
1741 m_cosemc[i] = cos(emcTrk->theta());
1742 m_phiemc[i] = emcTrk->phi();
1743 m_energy[i] = emcTrk->energy();
1744 m_eSeed[i] = emcTrk->eSeed();
1745 m_e3x3[i] = emcTrk->e3x3();
1746 m_e5x5[i] = emcTrk->e5x5();
1747 m_dang4c[i] = dang;
1748 m_dthe4c[i] = dthe;
1749 m_dphi4c[i] = dphi;
1750// m_dang4crt[i] = dangrt;
1751// m_dthe4crt[i] = dthert;
1752// m_dphi4crt[i] = dphirt;
1753
1754}
1755}
1756}
1757
1758}
1759
1760
1761
1762
1763
1764 //4pi
1765 if(m_skim4pi)
1766 {
1767
1768if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==0)
1769 {m_subalg1->execute(); // save gam4pi data
1770 Ncut[6]++;
1771 }
1772 }
1773
1774 //4k
1775 if(m_skim4k)
1776 {
1777
1778if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==1&&dtpr2kpilst[0]<-0.4&&dtpr2kpilst[1]<-0.4&&dtpr2kpilst[2]<-0.4&&dtpr2kpilst[3]<-0.4)
1779 {m_subalg2->execute(); // save gam4k data
1780 Ncut[7]++;
1781 }
1782 }
1783
1784 if(m_skim2pi2pr)
1785 {
1786if(mppreclst<3.06&&chis4c2kpilst<20&&((mchic2kpilst>3.39&&mchic2kpilst<3.44)||(mchic2kpilst>3.5&&mchic2kpilst<3.57))&&type2kpilst==2)
1787 {m_subalg3->execute(); // save gam 2(pi p) data
1788 // pi+ pi- with low momentum and p pbar with high momentum.
1789 Ncut[8]++;
1790 }
1791 }
1792
1793
1794
1795//cout<<"chis4c2kpilst="<<chis4c2kpilst<<endl;
1796if(m_rootput)
1797{
1798
1799if((mppreclst<3.06&&chis4c2kpilst<40)||((meelst>3.06&&meelst<3.12)&&(fabs(mppreclst-3.097)<0.01)))
1800 { Ncut[9]++;
1801 m_tuple1->write();
1802
1803 }
1804}
1805
1806 return StatusCode::SUCCESS;
1807 }
double sin(const BesAngle a)
Definition: BesAngle.h:210
double cos(const BesAngle a)
Definition: BesAngle.h:213
EvtRecTrackCol::iterator EvtRecTrackIterator
Definition: EvtRecTrack.h:111
std::vector< HepLorentzVector > Vp4
Definition: Gam4pikp.cxx:53
const double xmass[5]
Definition: Gam4pikp.cxx:50
const double mpro
Definition: Gam4pikp.cxx:49
const double mk
Definition: Gam4pikp.cxx:48
std::vector< double > Vdouble
Definition: Gam4pikp.cxx:54
const double mpi
Definition: Gam4pikp.cxx:47
std::vector< int > Vint
Definition: Gam4pikp.cxx:52
**********Class see also m_nmax DOUBLE PRECISION m_amel DOUBLE PRECISION m_x2 DOUBLE PRECISION m_alfinv DOUBLE PRECISION m_Xenph INTEGER m_KeyWtm INTEGER m_idyfs DOUBLE PRECISION m_zini DOUBLE PRECISION m_q2 DOUBLE PRECISION m_Wt_KF DOUBLE PRECISION m_WtCut INTEGER m_KFfin *COMMON c_KarLud $ !Input CMS energy[GeV] $ !CMS energy after beam spread beam strahlung[GeV] $ !Beam energy spread[GeV] $ !z boost due to beam spread $ !electron beam mass *ff pair spectrum $ !minimum v
Definition: KarLud.h:35
IMessageSvc * msgSvc()
int cellId() const
Definition: DstEmcShower.h:32
double latMoment() const
Definition: DstEmcShower.h:52
double a42Moment() const
Definition: DstEmcShower.h:54
double eSeed() const
Definition: DstEmcShower.h:47
double theta() const
Definition: DstEmcShower.h:38
int module() const
Definition: DstEmcShower.h:33
double e3x3() const
Definition: DstEmcShower.h:48
double phi() const
Definition: DstEmcShower.h:39
double secondMoment() const
Definition: DstEmcShower.h:51
double x() const
Definition: DstEmcShower.h:35
double e5x5() const
Definition: DstEmcShower.h:49
double time() const
Definition: DstEmcShower.h:50
double z() const
Definition: DstEmcShower.h:37
int numHits() const
Definition: DstEmcShower.h:30
double a20Moment() const
Definition: DstEmcShower.h:53
double energy() const
Definition: DstEmcShower.h:45
double y() const
Definition: DstEmcShower.h:36
const Hep3Vector emcPosition() const
Definition: DstExtTrack.h:126
const int emcVolumeNumber() const
Definition: DstExtTrack.h:132
double probPH() const
Definition: DstMdcDedx.h:66
double chiE() const
Definition: DstMdcDedx.h:59
int numTotalHits() const
Definition: DstMdcDedx.h:65
int numGoodHits() const
Definition: DstMdcDedx.h:64
double normPH() const
Definition: DstMdcDedx.h:67
double chiPi() const
Definition: DstMdcDedx.h:61
double chiK() const
Definition: DstMdcDedx.h:62
double chiMu() const
Definition: DstMdcDedx.h:60
double chiP() const
Definition: DstMdcDedx.h:63
const double px() const
static void setPidType(PidType pidType)
const double pz() const
const double py() const
const double p() const
const int charge() const
const double pxy() const
const double theta() const
Definition: DstMdcTrack.h:59
const double r() const
Definition: DstMdcTrack.h:64
const double py() const
Definition: DstMdcTrack.h:56
const HepSymMatrix err() const
const int charge() const
Definition: DstMdcTrack.h:53
const double px() const
Definition: DstMdcTrack.h:55
const double phi() const
Definition: DstMdcTrack.h:60
const double pz() const
Definition: DstMdcTrack.h:57
const double pxy() const
Definition: DstMdcTrack.h:54
const HepVector helix() const
......
const double z() const
Definition: DstMdcTrack.h:63
const double p() const
Definition: DstMdcTrack.h:58
const double y() const
Definition: DstMdcTrack.h:62
const double x() const
Definition: DstMdcTrack.h:61
double depth() const
Definition: DstMucTrack.h:45
int numLayers() const
Definition: DstMucTrack.h:42
virtual bool isVertexValid()=0
virtual double * SigmaPrimaryVertex()=0
virtual double * PrimaryVertex()=0
void AddFourMomentum(int number, HepLorentzVector p4)
HepLorentzVector pfit(int n) const
static KalmanKinematicFit * instance()
int useTof2() const
int onlyProton() const
int methodProbability() const
int useDedx() const
int onlyMuon() const
int onlyKaon() const
int onlyElectron() const
int onlyPion() const
int useTof1() const
void setChiMinCut(const double chi=4)
void setRecTrack(EvtRecTrack *trk)
double probKaon() const
Definition: ParticleID.h:124
void setMethod(const int method)
Definition: ParticleID.h:94
double chiTof2(int n) const
void identify(const int pidcase)
Definition: ParticleID.h:103
double probMuon() const
Definition: ParticleID.h:122
double probElectron() const
Definition: ParticleID.h:121
void usePidSys(const int pidsys)
Definition: ParticleID.h:97
static ParticleID * instance()
Definition: ParticleID.cxx:22
bool IsPidInfoValid() const
double probPion() const
Definition: ParticleID.h:123
double chiTof1(int n) const
void calculate()
Definition: ParticleID.cxx:97
void init()
Definition: ParticleID.cxx:27
double probProton() const
Definition: ParticleID.h:125
double chiDedx(int n) const
RecEmcID getClusterId() const
Definition: RecEmcShower.h:59
RecEmcID getShowerId() const
Definition: RecEmcShower.h:55
RecEmcEnergy getEAll() const
Definition: RecEmcShower.h:92
const HepVector & getZHelix() const
const HepSymMatrix & getZError() const
bool is_barrel() const
Definition: TofHitStatus.h:26
unsigned int layer() const
Definition: TofHitStatus.h:28
bool is_cluster() const
Definition: TofHitStatus.h:25
void setStatus(unsigned int status)
bool is_counter() const
Definition: TofHitStatus.h:24
bool is_readout() const
Definition: TofHitStatus.h:23
void setVBeamPosition(const HepSymMatrix VBeamPosition)
Definition: TrackPool.h:144
void AddTrack(const int number, const double mass, const RecMdcTrack *trk)
Definition: TrackPool.cxx:22
void setBeamPosition(const HepPoint3D BeamPosition)
Definition: TrackPool.h:143
WTrackParameter wtrk(int n) const
Definition: VertexFit.h:78
void init()
Definition: VertexFit.cxx:29
void AddVertex(int number, VertexParameter vpar, std::vector< int > lis)
Definition: VertexFit.cxx:89
HepSymMatrix Evx(int n) const
Definition: VertexFit.h:86
HepPoint3D vx(int n) const
Definition: VertexFit.h:84
static VertexFit * instance()
Definition: VertexFit.cxx:15
void Swim(int n)
Definition: VertexFit.h:58
bool Fit()
Definition: VertexFit.cxx:301
void setEvx(const HepSymMatrix &eVx)
void setVx(const HepPoint3D &vx)
const double ecms
Definition: inclkstar.cxx:42
_EXTERN_ std::string EvtRecEvent
Definition: EventModel.h:116
_EXTERN_ std::string EvtRecTrackCol
Definition: EventModel.h:117
float ptrk
const float pi
Definition: vector3.h:133

◆ finalize()

StatusCode Gam4pikp::finalize ( )

Definition at line 1811 of file Gam4pikp.cxx.

1811 {
1812// cout<<"total number: "<<Ncut[0]<<endl;
1813// cout<<"nGood>=4, nCharge==0: "<<Ncut[1]<<endl;
1814// cout<<"nGam>=1: "<<Ncut[2]<<endl;
1815// cout<<"cgp>2 cgm>2: "<<Ncut[3]<<endl;
1816// cout<<"2+ 2-: "<<Ncut[4]<<endl;
1817// cout<<"all cg cycle, chisq<200: "<<Ncut[5]<<endl;
1818// cout<<"4pi ok: "<<Ncut[6]<<endl;
1819// cout<<"4k ok: "<<Ncut[7]<<endl;
1820// cout<<"2pi 2p ok: "<<Ncut[8]<<endl;
1821// cout<<"ntuple write: "<<Ncut[9]<<endl;
1822 MsgStream log(msgSvc(), name());
1823 log << MSG::INFO << "in finalize()" << endmsg;
1824 return StatusCode::SUCCESS;
1825 }

◆ initialize()

StatusCode Gam4pikp::initialize ( )

Definition at line 83 of file Gam4pikp.cxx.

83 {
84 MsgStream log(msgSvc(), name());
85
86 log << MSG::INFO << "in initialize()" << endmsg;
87// setFilterPassed(false);
88 StatusCode status;
89
90if(m_rootput)
91{
92 NTuplePtr nt1(ntupleSvc(), "FILE1/total4c");
93 if ( nt1 ) m_tuple1 = nt1;
94 else {
95 m_tuple1 = ntupleSvc()->book ("FILE1/total4c", CLID_ColumnWiseTuple, "ks N-Tuple example");
96 if ( m_tuple1 ) {
97
98 status = m_tuple1->addItem ("run", m_run );
99 status = m_tuple1->addItem ("rec", m_rec );
100 status = m_tuple1->addItem ("mpprecall", m_mpprecall );
101 status = m_tuple1->addItem ("meeall", m_meeall );
102 status = m_tuple1->addItem ("ncgjs", m_ncgjs );
103 status = m_tuple1->addItem ("cla2kpi", m_cla2kpi );
104 status = m_tuple1->addItem("indexmc", m_idxmc, 0, 100);
105 status = m_tuple1->addIndexedItem("pdgid", m_idxmc, m_pdgid);
106
107 status = m_tuple1->addIndexedItem("motheridx", m_idxmc, m_motheridx);
108 status = m_tuple1->addItem("indexmdc", m_idxmdc, 0, 5000);
109 status = m_tuple1->addIndexedItem ("x0js", m_idxmdc, m_x0js);
110 status = m_tuple1->addIndexedItem ("y0js", m_idxmdc, m_y0js);
111 status = m_tuple1->addIndexedItem ("z0js",m_idxmdc, m_z0js);
112 status = m_tuple1->addIndexedItem ("r0js",m_idxmdc, m_r0js);
113 status = m_tuple1->addIndexedItem ("Rxyjs",m_idxmdc, m_Rxyjs);
114 status = m_tuple1->addIndexedItem ("Rzjs",m_idxmdc, m_Rzjs);
115 status = m_tuple1->addIndexedItem ("Rnxyjs",m_idxmdc, m_Rnxyjs);
116 status = m_tuple1->addIndexedItem ("phinjs",m_idxmdc, m_phinjs);
117 status = m_tuple1->addIndexedItem ("Rnzjs",m_idxmdc, m_Rnzjs);
118 status = m_tuple1->addItem ("ncy20", m_ncy20);
119 status = m_tuple1->addItem ("ncy30", m_ncy30);
120 status = m_tuple1->addIndexedItem("angjs5", m_idxmdc, m_angjs5);
121 status = m_tuple1->addIndexedItem("nearjs5", m_idxmdc, m_nearjs5);
122 status = m_tuple1->addIndexedItem("angjs6", m_idxmdc, m_angjs6);
123 status = m_tuple1->addIndexedItem("nearjs6", m_idxmdc, m_nearjs6);
124 status = m_tuple1->addIndexedItem("ang4pi5", m_idxmdc, m_ang4pi5);
125 status = m_tuple1->addIndexedItem("near4pi5", m_idxmdc, m_near4pi5);
126 status = m_tuple1->addIndexedItem("ang4pi6", m_idxmdc, m_ang4pi6);
127 status = m_tuple1->addIndexedItem("near4pi6", m_idxmdc, m_near4pi6);
128 status = m_tuple1->addIndexedItem("ppmdcjs", m_idxmdc, m_ppmdcjs);
129 status = m_tuple1->addIndexedItem("pxmdcjs", m_idxmdc, m_pxmdcjs);
130 status = m_tuple1->addIndexedItem("pymdcjs", m_idxmdc, m_pymdcjs);
131 status = m_tuple1->addIndexedItem("pzmdcjs", m_idxmdc, m_pzmdcjs);
132 status = m_tuple1->addIndexedItem("ppkaljs", m_idxmdc, m_ppkaljs);
133 status = m_tuple1->addIndexedItem("ptmdcjs", m_idxmdc, m_ptmdcjs);
134 status = m_tuple1->addIndexedItem("ptkaljs", m_idxmdc, m_ptkaljs);
135 status = m_tuple1->addIndexedItem("ppmdc2kpi", m_idxmdc, m_ppmdc2kpi);
136 status = m_tuple1->addIndexedItem("pxmdc2kpi", m_idxmdc, m_pxmdc2kpi);
137 status = m_tuple1->addIndexedItem("pymdc2kpi", m_idxmdc, m_pymdc2kpi);
138 status = m_tuple1->addIndexedItem("pzmdc2kpi", m_idxmdc, m_pzmdc2kpi);
139 status = m_tuple1->addIndexedItem("ppkal2kpi", m_idxmdc, m_ppkal2kpi);
140 status = m_tuple1->addIndexedItem("ptmdc2kpi", m_idxmdc, m_ptmdc2kpi);
141 status = m_tuple1->addIndexedItem("charge2kpi", m_idxmdc, m_charge2kpi);
142 status = m_tuple1->addIndexedItem("ptkal2kpi", m_idxmdc, m_ptkal2kpi);
143 status = m_tuple1->addItem ("cy2pi", m_cy2kpi, 0, 100 );
144 status = m_tuple1->addIndexedItem("comcs2kpi", m_cy2kpi, m_comcs2kpi);
145 status = m_tuple1->addItem ("chiejs", m_idxmdc, m_chiejs);
146 status = m_tuple1->addItem ("chimujs", m_idxmdc, m_chimujs);
147 status = m_tuple1->addItem ("chipijs", m_idxmdc, m_chipijs);
148 status = m_tuple1->addItem ("chikjs", m_idxmdc, m_chikjs);
149 status = m_tuple1->addItem ("chipjs", m_idxmdc, m_chipjs);
150 status = m_tuple1->addItem ("ghitjs", m_idxmdc, m_ghitjs);
151 status = m_tuple1->addItem ("thitjs", m_idxmdc, m_thitjs);
152 status = m_tuple1->addIndexedItem("probphjs", m_idxmdc, m_probphjs);
153 status = m_tuple1->addIndexedItem("normphjs", m_idxmdc, m_normphjs);
154 status = m_tuple1->addItem ("pdg", m_idxmdc, m_pdg);
155 status = m_tuple1->addItem ("cbmc", m_idxmdc, m_cbmc);
156 status = m_tuple1->addIndexedItem("sigmaetof2kpi", m_idxmdc, m_sigmaetof2kpi);
157 status = m_tuple1->addIndexedItem("sigmamutof2kpi", m_idxmdc, m_sigmamutof2kpi);
158 status = m_tuple1->addIndexedItem("sigmapitof2kpi", m_idxmdc, m_sigmapitof2kpi);
159 status = m_tuple1->addIndexedItem("sigmaktof2kpi", m_idxmdc, m_sigmaktof2kpi);
160 status = m_tuple1->addIndexedItem("sigmaprtof2kpi", m_idxmdc, m_sigmaprtof2kpi);
161 status = m_tuple1->addIndexedItem("t0tof2kpi", m_idxmdc, m_t0tof2kpi);
162 status = m_tuple1->addIndexedItem("errt0tof2kpi", m_idxmdc, m_errt0tof2kpi);
163
164 status = m_tuple1->addItem ("chie2kpi", m_idxmdc, m_chie2kpi);
165 status = m_tuple1->addItem ("chimu2kpi", m_idxmdc, m_chimu2kpi);
166 status = m_tuple1->addItem ("chipi2kpi", m_idxmdc, m_chipi2kpi);
167 status = m_tuple1->addItem ("chik2kpi", m_idxmdc, m_chik2kpi);
168 status = m_tuple1->addItem ("chip2kpi", m_idxmdc, m_chip2kpi);
169 status = m_tuple1->addItem ("ghit2kpi", m_idxmdc, m_ghit2kpi);
170 status = m_tuple1->addItem ("thit2kpi", m_idxmdc, m_thit2kpi);
171 status = m_tuple1->addIndexedItem("probph2kpi", m_idxmdc, m_probph2kpi);
172 status = m_tuple1->addIndexedItem("normph2kpi", m_idxmdc, m_normph2kpi);
173 status = m_tuple1->addIndexedItem("pidnum2kpi", m_idxmdc, m_pidnum2kpi);
174 status = m_tuple1->addIndexedItem("bjmucjs", m_idxmdc, m_bjmucjs);
175 status = m_tuple1->addIndexedItem("bjmuc2kpi", m_idxmdc, m_bjmuc2kpi);
176 status = m_tuple1->addIndexedItem("bjemcjs", m_idxmdc, m_bjemcjs);
177 status = m_tuple1->addIndexedItem("bjemc2kpi", m_idxmdc, m_bjemc2kpi);
178 status = m_tuple1->addIndexedItem("bjtofjs", m_idxmdc, m_bjtofjs);
179 status = m_tuple1->addIndexedItem("bjtof2kpi", m_idxmdc, m_bjtof2kpi);
180 status = m_tuple1->addIndexedItem("bjtofvaljs", m_idxmdc, m_bjtofvaljs);
181 status = m_tuple1->addIndexedItem("bjtofval2kpi", m_idxmdc, m_bjtofval2kpi);
182
183 status = m_tuple1->addIndexedItem("emcjs", m_idxmdc, m_emcjs);
184 status = m_tuple1->addIndexedItem("evpjs", m_idxmdc, m_evpjs);
185 status = m_tuple1->addIndexedItem("timecgjs", m_idxmdc, m_timecgjs);
186 status = m_tuple1->addIndexedItem("depthjs", m_idxmdc, m_depthmucjs);
187 status = m_tuple1->addIndexedItem("layermucjs", m_idxmdc, m_layermucjs);
188
189 status = m_tuple1->addIndexedItem("emc2kpi", m_idxmdc, m_emc2kpi);
190 status = m_tuple1->addIndexedItem("evp2kpi", m_idxmdc, m_evp2kpi);
191 status = m_tuple1->addIndexedItem("timecg2kpi", m_idxmdc, m_timecg2kpi);
192 status = m_tuple1->addIndexedItem("depth2kpi", m_idxmdc, m_depthmuc2kpi);
193 status = m_tuple1->addIndexedItem("layermuc2kpi", m_idxmdc, m_layermuc2kpi);
194
195 status = m_tuple1->addIndexedItem("cotof1js", m_idxmdc, m_cotof1js);
196 status = m_tuple1->addIndexedItem("cotof2js", m_idxmdc, m_cotof2js);
197 status = m_tuple1->addIndexedItem("counterjs", m_idxmdc, m_counterjs);
198 status = m_tuple1->addIndexedItem("barreljs", m_idxmdc, m_barreljs);
199 status = m_tuple1->addIndexedItem("layertofjs", m_idxmdc, m_layertofjs);
200 status = m_tuple1->addIndexedItem("readoutjs", m_idxmdc, m_readoutjs);
201 status = m_tuple1->addIndexedItem("clusterjs", m_idxmdc, m_clusterjs);
202 status = m_tuple1->addIndexedItem("betajs", m_idxmdc, m_betajs);
203 status = m_tuple1->addIndexedItem("tofjs", m_idxmdc, m_tofjs);
204 status = m_tuple1->addIndexedItem("tofpathjs", m_idxmdc, m_tofpathjs);
205 status = m_tuple1->addIndexedItem("zhitjs", m_idxmdc, m_zhitjs);
206 status = m_tuple1->addIndexedItem("tofIDjs", m_idxmdc, m_tofIDjs);
207 status = m_tuple1->addIndexedItem("clusterIDjs", m_idxmdc, m_clusterIDjs);
208 status = m_tuple1->addIndexedItem("texejs", m_idxmdc, m_texejs);
209 status = m_tuple1->addIndexedItem("texmujs", m_idxmdc, m_texmujs);
210 status = m_tuple1->addIndexedItem("texpijs", m_idxmdc, m_texpijs);
211 status = m_tuple1->addIndexedItem("texkjs", m_idxmdc, m_texkjs);
212 status = m_tuple1->addIndexedItem("texprjs", m_idxmdc, m_texprjs);
213 status = m_tuple1->addIndexedItem("dtejs", m_idxmdc, m_dtejs);
214 status = m_tuple1->addIndexedItem("dtmujs", m_idxmdc, m_dtmujs);
215 status = m_tuple1->addIndexedItem("dtpijs", m_idxmdc, m_dtpijs);
216 status = m_tuple1->addIndexedItem("dtkjs", m_idxmdc, m_dtkjs);
217 status = m_tuple1->addIndexedItem("dtprjs", m_idxmdc, m_dtprjs);
218 status = m_tuple1->addIndexedItem("sigmaetofjs", m_idxmdc, m_sigmaetofjs);
219 status = m_tuple1->addIndexedItem("sigmamutofjs", m_idxmdc, m_sigmamutofjs);
220 status = m_tuple1->addIndexedItem("sigmapitofjs", m_idxmdc, m_sigmapitofjs);
221 status = m_tuple1->addIndexedItem("sigmaktofjs", m_idxmdc, m_sigmaktofjs);
222 status = m_tuple1->addIndexedItem("sigmaprtofjs", m_idxmdc, m_sigmaprtofjs);
223 status = m_tuple1->addIndexedItem("t0tofjs", m_idxmdc,m_t0tofjs);
224 status = m_tuple1->addIndexedItem("errt0tofjs", m_idxmdc,m_errt0tofjs);
225 status = m_tuple1->addIndexedItem("cotof12kpi", m_idxmdc, m_cotof12kpi);
226 status = m_tuple1->addIndexedItem("cotof22kpi", m_idxmdc, m_cotof22kpi);
227 status = m_tuple1->addIndexedItem("counter2kpi", m_idxmdc, m_counter2kpi);
228 status = m_tuple1->addIndexedItem("barrel2kpi", m_idxmdc, m_barrel2kpi);
229 status = m_tuple1->addIndexedItem("layertof2kpi", m_idxmdc, m_layertof2kpi);
230 status = m_tuple1->addIndexedItem("readout2kpi", m_idxmdc, m_readout2kpi);
231 status = m_tuple1->addIndexedItem("cluster2kpi", m_idxmdc, m_cluster2kpi);
232 status = m_tuple1->addIndexedItem("beta2kpi", m_idxmdc, m_beta2kpi);
233 status = m_tuple1->addIndexedItem("tof2kpi", m_idxmdc, m_tof2kpi);
234 status = m_tuple1->addIndexedItem("tofpath2kpi", m_idxmdc, m_tofpath2kpi);
235 status = m_tuple1->addIndexedItem("zhit2kpi", m_idxmdc, m_zhit2kpi);
236 status = m_tuple1->addIndexedItem("tofID2kpi", m_idxmdc, m_tofID2kpi);
237 status = m_tuple1->addIndexedItem("clusterID2kpi", m_idxmdc, m_clusterID2kpi);
238 status = m_tuple1->addIndexedItem("texe2kpi", m_idxmdc, m_texe2kpi);
239 status = m_tuple1->addIndexedItem("texmu2kpi", m_idxmdc, m_texmu2kpi);
240 status = m_tuple1->addIndexedItem("texpi2kpi", m_idxmdc, m_texpi2kpi);
241 status = m_tuple1->addIndexedItem("texk2kpi", m_idxmdc, m_texk2kpi);
242 status = m_tuple1->addIndexedItem("texpr2kpi", m_idxmdc, m_texpr2kpi);
243 status = m_tuple1->addIndexedItem("dte2kpi", m_idxmdc, m_dte2kpi);
244 status = m_tuple1->addIndexedItem("dtmu2kpi", m_idxmdc, m_dtmu2kpi);
245 status = m_tuple1->addIndexedItem("dtpi2kpi", m_idxmdc, m_dtpi2kpi);
246 status = m_tuple1->addIndexedItem("dtk2kpi", m_idxmdc, m_dtk2kpi);
247 status = m_tuple1->addIndexedItem("dtpr2kpi", m_idxmdc, m_dtpr2kpi);
248 status = m_tuple1->addIndexedItem("costpid2kpi", m_idxmdc, m_costpid2kpi);
249 status = m_tuple1->addIndexedItem("dedxpid2kpi", m_idxmdc, m_dedxpid2kpi);
250 status = m_tuple1->addIndexedItem("tof1pid2kpi", m_idxmdc, m_tof1pid2kpi);
251 status = m_tuple1->addIndexedItem("tof2pid2kpi", m_idxmdc, m_tof2pid2kpi);
252 status = m_tuple1->addIndexedItem("probe2kpi", m_idxmdc, m_probe2kpi);
253 status = m_tuple1->addIndexedItem("probmu2kpi", m_idxmdc, m_probmu2kpi);
254 status = m_tuple1->addIndexedItem("probpi2kpi", m_idxmdc, m_probpi2kpi);
255 status = m_tuple1->addIndexedItem("probk2kpi", m_idxmdc, m_probk2kpi);
256 status = m_tuple1->addIndexedItem("probpr2kpi", m_idxmdc, m_probpr2kpi);
257
258 status = m_tuple1->addIndexedItem("chipidxpid2kpi", m_idxmdc, m_chipidxpid2kpi);
259 status = m_tuple1->addIndexedItem("chipitof1pid2kpi", m_idxmdc, m_chipitof1pid2kpi);
260 status = m_tuple1->addIndexedItem("chipitof2pid2kpi", m_idxmdc, m_chipitof2pid2kpi);
261 status = m_tuple1->addIndexedItem("chipitofpid2kpi", m_idxmdc, m_chipitofpid2kpi);
262 status = m_tuple1->addIndexedItem("chipitofepid2kpi", m_idxmdc, m_chipitofepid2kpi);
263 status = m_tuple1->addIndexedItem("chipitofqpid2kpi", m_idxmdc, m_chipitofqpid2kpi);
264 status = m_tuple1->addIndexedItem("probpidxpid2kpi", m_idxmdc, m_probpidxpid2kpi);
265 status = m_tuple1->addIndexedItem("probpitofpid2kpi", m_idxmdc, m_probpitofpid2kpi);
266 status = m_tuple1->addIndexedItem("chikdxpid2kpi", m_idxmdc, m_chikdxpid2kpi);
267 status = m_tuple1->addIndexedItem("chiktof1pid2kpi", m_idxmdc, m_chiktof1pid2kpi);
268 status = m_tuple1->addIndexedItem("chiktof2pid2kpi", m_idxmdc, m_chiktof2pid2kpi);
269 status = m_tuple1->addIndexedItem("chiktofpid2kpi", m_idxmdc, m_chiktofpid2kpi);
270 status = m_tuple1->addIndexedItem("chiktofepid2kpi", m_idxmdc, m_chiktofepid2kpi);
271 status = m_tuple1->addIndexedItem("chiktofqpid2kpi", m_idxmdc, m_chiktofqpid2kpi);
272 status = m_tuple1->addIndexedItem("probkdxpid2kpi", m_idxmdc, m_probkdxpid2kpi);
273 status = m_tuple1->addIndexedItem("probktofpid2kpi", m_idxmdc, m_probktofpid2kpi);
274
275 status = m_tuple1->addIndexedItem("chiprdxpid2kpi", m_idxmdc, m_chiprdxpid2kpi);
276 status = m_tuple1->addIndexedItem("chiprtof1pid2kpi", m_idxmdc, m_chiprtof1pid2kpi);
277 status = m_tuple1->addIndexedItem("chiprtof2pid2kpi", m_idxmdc, m_chiprtof2pid2kpi);
278 status = m_tuple1->addIndexedItem("chiprtofpid2kpi", m_idxmdc, m_chiprtofpid2kpi);
279 status = m_tuple1->addIndexedItem("chiprtofepid2kpi", m_idxmdc, m_chiprtofepid2kpi);
280 status = m_tuple1->addIndexedItem("chiprtofqpid2kpi", m_idxmdc, m_chiprtofqpid2kpi);
281 status = m_tuple1->addIndexedItem("probprdxpid2kpi", m_idxmdc, m_probprdxpid2kpi);
282 status = m_tuple1->addIndexedItem("probprtofpid2kpi", m_idxmdc, m_probprtofpid2kpi);
283
284 status = m_tuple1->addIndexedItem("cosmdcjs", m_idxmdc, m_cosmdcjs);
285 status = m_tuple1->addIndexedItem("phimdcjs", m_idxmdc, m_phimdcjs);
286 status = m_tuple1->addIndexedItem("cosmdc2kpi", m_idxmdc, m_cosmdc2kpi);
287 status = m_tuple1->addIndexedItem("phimdc2kpi", m_idxmdc, m_phimdc2kpi);
288
289 status = m_tuple1->addIndexedItem("dedxpidjs", m_idxmdc, m_dedxpidjs);
290 status = m_tuple1->addIndexedItem("tof1pidjs", m_idxmdc, m_tof1pidjs);
291 status = m_tuple1->addIndexedItem("tof2pidjs", m_idxmdc, m_tof2pidjs);
292 status = m_tuple1->addIndexedItem("probejs", m_idxmdc, m_probejs);
293 status = m_tuple1->addIndexedItem("probmujs", m_idxmdc, m_probmujs);
294 status = m_tuple1->addIndexedItem("probpijs", m_idxmdc, m_probpijs);
295 status = m_tuple1->addIndexedItem("probkjs", m_idxmdc, m_probkjs);
296 status = m_tuple1->addIndexedItem("probprjs", m_idxmdc, m_probprjs);
297 status = m_tuple1->addItem ("mchic2kpi", m_mchic2kpi);
298 status = m_tuple1->addItem ("mpsip2kpi", m_mpsip2kpi);
299 status = m_tuple1->addItem ("chis2kpi", m_chis2kpi);
300 status = m_tuple1->addItem ("mchic4c2kpi", m_mchic4c2kpi);
301 status = m_tuple1->addItem ("mpsip4c2kpi", m_mpsip4c2kpi);
302 status = m_tuple1->addItem ("chis4c2kpi", m_chis4c2kpi);
303
304 status = m_tuple1->addItem("indexemc", m_idxemc, 0, 5000);
305 status = m_tuple1->addIndexedItem("numHits", m_idxemc, m_numHits);
306 status = m_tuple1->addIndexedItem("secmom", m_idxemc, m_secondmoment);
307 status = m_tuple1->addIndexedItem("latmom", m_idxemc, m_latmoment);
308 status = m_tuple1->addIndexedItem("timegm", m_idxemc, m_timegm);
309 status = m_tuple1->addIndexedItem("cellId", m_idxemc, m_cellId);
310 status = m_tuple1->addIndexedItem("module", m_idxemc, m_module);
311 status = m_tuple1->addIndexedItem("a20Moment", m_idxemc, m_a20Moment);
312 status = m_tuple1->addIndexedItem("a42Moment", m_idxemc, m_a42Moment);
313 status = m_tuple1->addIndexedItem("getEAll", m_idxemc, m_getEAll);
314 status = m_tuple1->addIndexedItem("getShowerId", m_idxemc, m_getShowerId);
315 status = m_tuple1->addIndexedItem("getClusterId", m_idxemc, m_getClusterId);
316 status = m_tuple1->addIndexedItem("x", m_idxemc, m_x);
317 status = m_tuple1->addIndexedItem("y", m_idxemc, m_y);
318 status = m_tuple1->addIndexedItem("z", m_idxemc, m_z);
319 status = m_tuple1->addIndexedItem("cosemc", m_idxemc, m_cosemc);
320 status = m_tuple1->addIndexedItem("phiemc", m_idxemc, m_phiemc);
321 status = m_tuple1->addIndexedItem("energy", m_idxemc, m_energy);
322 status = m_tuple1->addIndexedItem("e1", m_idxemc, m_eSeed);
323 status = m_tuple1->addIndexedItem("e9", m_idxemc, m_e3x3);
324 status = m_tuple1->addIndexedItem("e25", m_idxemc, m_e5x5);
325 status = m_tuple1->addIndexedItem("dang4c", m_idxemc, m_dang4c);
326 status = m_tuple1->addIndexedItem("dthe4c", m_idxemc, m_dthe4c);
327 status = m_tuple1->addIndexedItem("dphi4c", m_idxemc, m_dphi4c);
328 status = m_tuple1->addIndexedItem("dang4crt", m_idxemc, m_dang4crt);
329 status = m_tuple1->addIndexedItem("dthe4crt", m_idxemc, m_dthe4crt);
330 status = m_tuple1->addIndexedItem("dphi4crt", m_idxemc, m_dphi4crt);
331 status = m_tuple1->addIndexedItem("phtof", m_idxemc, 3, m_phgmtof,0.0,10000.0);
332 status = m_tuple1->addIndexedItem("phgmtof0", m_idxemc, m_phgmtof0);
333 status = m_tuple1->addIndexedItem("phgmtof1", m_idxemc, m_phgmtof1);
334 status = m_tuple1->addIndexedItem("phgmtof2", m_idxemc, m_phgmtof2);
335
336 }
337 else {
338 log << MSG::ERROR << " Cannot book N-tuple:" << long(m_tuple1) << endmsg;
339 return StatusCode::FAILURE;
340 }
341 }
342}
343
344StatusCode sc;
345
346 if(m_skim4pi)
347 {
348 sc = createSubAlgorithm( "EventWriter", "Selectgam4pi", m_subalg1);
349 if( sc.isFailure() ) {
350 log << MSG::ERROR << "Error creating Sub-Algorithm Selectgam4pi" <<endreq;
351 return sc;
352 } else {
353 log << MSG::INFO << "Success creating Sub-Algorithm Selectgam4pi" <<endreq;
354 }
355 }
356
357
358 if(m_skim4k)
359 {
360 sc = createSubAlgorithm( "EventWriter", "Selectgam4k", m_subalg2);
361 if( sc.isFailure() ) {
362 log << MSG::ERROR << "Error creating Sub-Algorithm Selectgam4k" <<endreq;
363 return sc;
364 } else {
365 log << MSG::INFO << "Success creating Sub-Algorithm Selectgam4k" <<endreq;
366 }
367 }
368
369 if(m_skim2pi2pr)
370 {
371 sc = createSubAlgorithm( "EventWriter", "Selectgam2pi2pr", m_subalg3);
372 if( sc.isFailure() ) {
373 log << MSG::ERROR << "Error creating Sub-Algorithm Selectgam2pi2pr" <<endreq;
374 return sc;
375 } else {
376 log << MSG::INFO << "Success creating Sub-Algorithm Selectgam2pi2pr" <<endreq;
377 }
378 }
379
380
381
382
383 //
384 //--------end of book--------
385 //
386
387 log << MSG::INFO << "successfully return from initialize()" <<endmsg;
388 return StatusCode::SUCCESS;
389
390}
INTupleSvc * ntupleSvc()

The documentation for this class was generated from the following files: