BOSS 7.0.5
BESIII Offline Software System
Loading...
Searching...
No Matches
Public Member Functions | List of all members
MdcCalib Class Referenceabstract

#include <MdcCalib.h>

+ Inheritance diagram for MdcCalib:

Public Member Functions

 MdcCalib ()
 
virtual ~MdcCalib ()
 
virtual void initialize (TObjArray *hlist, IMdcGeomSvc *mdcGeomSvc, IMdcCalibFunSvc *mdcFunSvc, IMdcUtilitySvc *mdcUtilitySvc)=0
 
virtual void setParam (MdcCalParams &param)=0
 
virtual int fillHist (MdcCalEvent *event)=0
 
virtual int updateConst (MdcCalibConst *calconst)=0
 
virtual void clear ()=0
 
 MdcCalib ()
 
virtual ~MdcCalib ()
 
virtual void initialize (TObjArray *hlist, IMdcGeomSvc *mdcGeomSvc, IMdcCalibFunSvc *mdcFunSvc, IMdcUtilitySvc *mdcUtilitySvc)=0
 
virtual void setParam (MdcCalParams &param)=0
 
virtual int fillHist (MdcCalEvent *event)=0
 
virtual int updateConst (MdcCalibConst *calconst)=0
 
virtual void clear ()=0
 

Detailed Description

Definition at line 36 of file InstallArea/include/MdcCalibAlg/MdcCalibAlg/MdcCalib.h.

Constructor & Destructor Documentation

◆ MdcCalib() [1/2]

MdcCalib::MdcCalib ( )

Definition at line 36 of file MdcCalib.cxx.

36 {
37 m_nEvt=0;
38 m_cut1=0;
39 m_cut2=0;
40 m_cut3=0;
41 m_cut4=0;
42 m_cut5=0;
43 m_cut6=0;
44 m_nGrPoint = 0;
45 fgReadWireEff = false;
46
47 int lay;
48 for(lay=0; lay<43; lay++){
49 if(lay < 15) m_nEntr[lay] = 1;
50 else m_nEntr[lay] = 2;
51 }
52 m_dwid = 0.5; // mm
53 m_fgIni = false;
54
55 m_phiWid = PI2 / (double)NPhiBin;
56 m_theWid = 2.0 / (double)NThetaBin;
57
58 m_nEvtNtuple = 0;
59
60 for(lay=0; lay<MdcCalNLayer; lay++){
61 if(lay < 8) m_nBin[lay] = 12;
62 else m_nBin[lay] = 16;
63 }
64
65 // setting boundary layer flags
66 for(lay=0; lay<MdcCalNLayer; lay++){
67 if((0==lay) || (7==lay) || (8==lay) || (19==lay) || (20==lay) ||
68 (35==lay) || (36==lay) || (42==lay) ) m_layBound[lay] = true;
69 else m_layBound[lay] = false;
70 }
71}

◆ ~MdcCalib() [1/2]

MdcCalib::~MdcCalib ( )
virtual

Definition at line 73 of file MdcCalib.cxx.

73 {
74}

◆ MdcCalib() [2/2]

MdcCalib::MdcCalib ( )

◆ ~MdcCalib() [2/2]

virtual MdcCalib::~MdcCalib ( )
virtual

Member Function Documentation

◆ clear() [1/2]

void MdcCalib::clear ( )
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

Definition at line 76 of file MdcCalib.cxx.

76 {
77 int lay;
78 for(lay=0; lay<m_nlayer; lay++){
79 delete m_htraw[lay];
80 delete m_htdr[lay];
81 delete m_hresInc[lay];
82 delete m_hresExc[lay];
83 delete m_hresAve[lay];
84 delete m_hadc[lay];
85 for (int lr=0; lr<2; lr++){
86 delete m_htdrlr[lay][lr];
87 delete m_hreslrInc[lay][lr];
88 delete m_hreslrExc[lay][lr];
89 delete m_hreslrAve[lay][lr];
90 }
91 }
92
93 delete m_effNtrk;
94 delete m_effNtrkRecHit;
95 delete m_hresAllInc;
96 delete m_hresAllExc;
97 delete m_hresAllAve;
98 for(int i=0; i<14; i++){
99 delete m_hresAveAllQ[i];
100 delete m_hresAveOutQ[i];
101 }
102 for(lay=0; lay<43; lay++){
103 for(int i=0; i<14; i++) delete m_hresAveLayQ[lay][i];
104 }
105 delete m_hresInnInc;
106 delete m_hresInnExc;
107 delete m_hresStpInc;
108 delete m_hresStpExc;
109 delete m_hresOutInc;
110 delete m_hresOutExc;
111 for(int iEs=0; iEs<m_param.nEsFlag; iEs++) delete m_hTes[iEs];
112 delete m_hbbTrkFlg;
113 delete m_hTesAll;
114 delete m_hTesGood;
115 delete m_hTesAllFlag;
116 delete m_hTesRec;
117 delete m_hTesCalFlag;
118 delete m_hTesCalUse;
119 delete m_hnRawHit;
120 delete m_hpt;
121 delete m_hpMax;
122 delete m_hpMaxCms;
123 delete m_hptPos;
124 delete m_hptNeg;
125 delete m_hp;
126 delete m_hp_cms;
127 delete m_hpPos;
128 delete m_hpNeg;
129 delete m_hpPoscms;
130 delete m_hpNegcms;
131 delete m_hp_cut;
132 delete m_hchisq;
133 delete m_hnTrk;
134 delete m_hnTrkCal;
135 delete m_hnhitsRec;
136 delete m_hnhitsRecInn;
137 delete m_hnhitsRecStp;
138 delete m_hnhitsRecOut;
139 delete m_hnhitsCal;
140 delete m_hnhitsCalInn;
141 delete m_hnhitsCalStp;
142 delete m_hnhitsCalOut;
143 delete m_wirehitmap;
144 delete m_layerhitmap;
145 delete m_hnoisephi;
146 delete m_hnoiselay;
147 delete m_hnoisenhits;
148 delete m_hratio;
149 delete m_hdr;
150 delete m_hphi0;
151 delete m_hkap;
152 delete m_hdz;
153 delete m_htanl;
154 delete m_hcosthe;
155 delete m_hcostheNeg;
156 delete m_hcosthePos;
157 delete m_hx0;
158 delete m_hy0;
159 delete m_hdelZ0;
160 delete m_grX0Y0;
161 delete m_hitEffAll;
162 delete m_hitEffRaw;
163 delete m_hitEffRec;
164 int bin;
165 int thbin;
166 for(bin=0; bin<NPhiBin; bin++){
167 delete m_ppPhi[bin];
168 delete m_pnPhi[bin];
169 delete m_ppPhiCms[bin];
170 delete m_pnPhiCms[bin];
171 for(thbin=0; thbin<NThetaBin; thbin++){
172 delete m_ppThePhi[thbin][bin];
173 delete m_pnThePhi[thbin][bin];
174 delete m_ppThePhiCms[thbin][bin];
175 delete m_pnThePhiCms[thbin][bin];
176 }
177 }
178 for(thbin=0; thbin<NThetaBin; thbin++){
179 delete m_ppThe[thbin];
180 delete m_pnThe[thbin];
181 delete m_ppTheCms[thbin];
182 delete m_pnTheCms[thbin];
183 }
184
185 for(unsigned i=0; i<m_hr2dInc.size(); i++){
186 delete m_hr2dInc[i];
187 delete m_hr2dExc[i];
188 }
189 m_hr2dInc.clear();
190 m_hr2dExc.clear();
191 m_mapr2d.clear();
192
193 delete m_fdTime;
194 delete m_fdAdc;
195 delete m_fdres;
196 delete m_fdresAve;
197 delete m_fdres2d;
198 delete m_fdcom;
199 delete m_fdResQ;
200}
bin
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per bin
Definition: FoamA.h:85

Referenced by GrXtMdcCalib::clear(), PreT0MdcCalib::clear(), QtMdcCalib::clear(), T0MdcCalib::clear(), Wr2dMdcCalib::clear(), WrMdcCalib::clear(), XtInteMdcCalib::clear(), XtMdcCalib::clear(), and MdcCalibAlg::finalize().

◆ clear() [2/2]

virtual void MdcCalib::clear ( )
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

◆ fillHist() [1/2]

int MdcCalib::fillHist ( MdcCalEvent *  event)
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

Definition at line 692 of file MdcCalib.cxx.

692 {
693 IMessageSvc* msgSvc;
694 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
695 MsgStream log(msgSvc, "MdcCalib");
696 log << MSG::DEBUG << "MdcCalib::fillHist()" << endreq;
697
698 int i;
699 int k;
700 int lay;
701 int cel;
702 int wir;
703 int lr;
704 int stat;
705
706 int hid;
707 int key;
708 int iEntr;
709 int bin;
710
711 int phiBin;
712 int phiBinCms;
713 int theBin;
714 int theBinCms;
715 double lamda;
716 double theta;
717
718 double qhit;
719 double traw;
720 double tdr;
721 double doca;
722 double resiInc;
723 double resiExc;
724 double entr;
725 double pt;
726 double p;
727 double p_cms;
728 double chisq;
729 double ecm = m_param.ecm;
730 double xboost = m_param.boostPar[0] * ecm;
731 double yboost = m_param.boostPar[1];
732 double zboost = m_param.boostPar[2];
733 HepLorentzVector p4;
734
735 double dr;
736 double phi0;
737 double dz;
738 double kap;
739 double tanl;
740
741 double x0;
742 double y0;
743 double zminus = 9999.0;
744 double zplus = -9999.0;
745
746 double hitphi;
747 double philab;
748 double phicms;
749 double thetacms;
750 double costheCMS;
751
752 double dphi;
753 double wphi;
754 double xx;
755 double yy;
756 double rr;
757
758 int nhitlay;
759 bool fgHitLay[MdcCalNLayer];
760 bool fgTrk;
761
762 int ntrk = event -> getNTrk();
763 int nhitRec;
764 int nhitRecInn;
765 int nhitRecStp;
766 int nhitRecOut;
767 int nhitCal;
768 int nhitCalInn;
769 int nhitCalStp;
770 int nhitCalOut;
771 MdcCalRecTrk* rectrk;
772 MdcCalRecHit* rechit;
773
774 int fgAdd[43]; // for calculating layer efficiency
775
776 // read dead wire
777 if(!fgReadWireEff){
778 for(lay=0; lay<MdcCalNLayer; lay++){
779 int ncel = m_mdcGeomSvc->Layer(lay)->NCell();
780 for(cel=0; cel<ncel; cel++){
781 double eff = m_mdcFunSvc->getWireEff(lay, cel);
782 if(eff > 0.5) m_fgGoodWire[lay][cel] = true;
783 else m_fgGoodWire[lay][cel] = false;
784 if(eff<0.9) cout << "dead channel: " << setw(5) << lay << setw(5) << cel << endl;
785 }
786 }
787 fgReadWireEff = true;
788 }
789
790 int nRawHit = event->getNRawHitTQ();
791 m_hnRawHit->Fill(nRawHit);
792
793 IDataProviderSvc* eventSvc = NULL;
794 Gaudi::svcLocator()->service("EventDataSvc", eventSvc);
795
796 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc,"/Event/EventHeader");
797 if (!eventHeader) {
798 log << MSG::FATAL << "Could not find Event Header" << endreq;
799 return( StatusCode::FAILURE);
800 }
801 int iRun = eventHeader->runNumber();
802 int iEvt = eventHeader->eventNumber();
803
804 int esTimeflag = event->getEsFlag();
805 double tes = event->getTes();
806 bool esCutFg = event->getEsCutFlag();
807 int iEs = event->getNesCutFlag();
808 //calculate the efficiency of Bhabha event
809 if (-1 != esTimeflag) {
810 SmartDataPtr<RecMdcTrackCol> newtrkCol(eventSvc, "/Event/Recon/RecMdcTrackCol");
811 if(!newtrkCol){
812 log << MSG::ERROR << "Could not find RecMdcTrackCol" << endreq;
813 return ( StatusCode::FAILURE );
814 }
815 int nGoodTrk = 0;
816 int icharge = 0;
817 Vp4 p4p;
818 Vp4 p4m;
819 RecMdcTrackCol::iterator it_trk = newtrkCol->begin();
820 for(; it_trk != newtrkCol->end(); it_trk++){
821 double mass = 0.000511; // only for eletron
822 HepLorentzVector p4 = (*it_trk)->p4(mass);
823 icharge = (*it_trk)->charge();
824 if (icharge > 0) p4p.push_back(p4);
825 else p4m.push_back(p4);
826 }
827 if (1 == p4p.size() * p4m.size()){
828 double dang = p4p[0].vect().angle(p4m[0].vect());
829 m_hbbTrkFlg->Fill(1);
830 if (dang > 2.94) {
831 m_hbbTrkFlg->Fill(2);
832 }
833 }
834
835 }
836 m_hTesAll->Fill(tes);
837// cout << "tes " << tes << endl;
838 if (-1 != esTimeflag) m_hTesGood->Fill(tes);
839 m_hTesAllFlag->Fill(esTimeflag);
840 if(ntrk > 0) m_hTesRec->Fill(tes);
841 if( (iEs >= 0) && (iEs < m_param.nEsFlag) ) m_hTes[iEs]->Fill(tes);
842 if( esCutFg ) m_hTesCalFlag->Fill(tes);
843 else return -1;
844
845 if(! m_fgIni){
846 for(lay=0; lay<MdcCalNLayer; lay++){
847 if(lay < 8) m_docaMax[lay] = m_param.maxDocaInner;
848 else m_docaMax[lay] = m_param.maxDocaOuter;
849 }
850 m_fgIni = true;
851 }
852
853 bool trkFlag[200]; // for calculating hit efficiency without impact of track fitting
854 for(i=0; i<200; i++) trkFlag[i] = false;
855
856 int ntrkCal = 0;
857 double pTrk[2];
858 double pTrkcms[2];
859 double hitphiplus = 9999.0;
860 double hitthetaplus = 9999.0;
861 double hitphiminus = -9999.0;
862 double hitthetaminus = -9999.0;
863 Vp4 pp;
864 Vp4 pm;
865 pp.clear();
866 pm.clear();
867 int phibinp;
868 int phibinm;
869
870 m_hnTrk->Fill(ntrk);
871 if((ntrk < m_param.nTrkCut[0]) || (ntrk > m_param.nTrkCut[1])){
872 m_cut1++;
873 return -1;
874 }
875// if(ntrk > 2) return -1;
876 for(i=0; i<ntrk; i++){
877 fgTrk = true;
878 rectrk = event -> getRecTrk(i);
879 nhitRec = rectrk -> getNHits();
880 m_hnhitsRec -> Fill( nhitRec );
881
882 for(lay=0; lay<MdcCalNLayer; lay++){
883 fgHitLay[lay] = false;
884 }
885
886 nhitRecInn = 0;
887 nhitRecStp = 0;
888 nhitRecOut = 0;
889 for(k=0; k<nhitRec; k++){
890 rechit = rectrk -> getRecHit(k);
891 lay = rechit -> getLayid();
892 doca = rechit -> getDocaExc();
893 resiExc = rechit -> getResiExc();
894 fgHitLay[lay] = true;
895
896 if(lay < 8) nhitRecInn++;
897 else if(lay < 20) nhitRecStp++;
898 else nhitRecOut++;
899 }
900 m_hnhitsRecInn->Fill(nhitRecInn);
901 m_hnhitsRecStp->Fill(nhitRecStp);
902 m_hnhitsRecOut->Fill(nhitRecOut);
903
904 // get momentum
905 pt = rectrk -> getPt();
906 p = rectrk -> getP();
907
908 // boost P to the cms
909 p4 = rectrk->getP4();
910 HepLorentzVector psip(xboost, yboost, zboost, ecm);
911 Hep3Vector boostv = psip.boostVector();
912 p4.boost(- boostv);
913 p_cms = p4.rho();
914 phicms = p4.phi();
915 if(phicms < 0) phicms += PI2;
916 thetacms = p4.theta();
917 costheCMS = cos(thetacms);
918 if (pt < 0) p_cms *= -1.0;
919 p4.boost(boostv);
920// cout << setw(15) << p << setw(15) << p_cms << setw(15) << costheCMS << endl;
921
922 // cos(theta) cut
923 if( (costheCMS < m_param.costheCut[0]) || (costheCMS > m_param.costheCut[1]) ){
924 m_cut2++;
925 continue;
926 }
927
928 // dr cut
929 dr = rectrk->getDr();
930 if(fabs(dr) > m_param.drCut){
931 m_cut3++;
932 continue;
933 }
934
935 // dz cut
936 dz = rectrk->getDz();
937 if(fabs(dz) > m_param.dzCut){
938 m_cut4++;
939 continue;
940 }
941
942// if(! fgTrk) continue;
943
944 // hit layer cut
945 nhitlay = 0;
946 for(lay=0; lay<MdcCalNLayer; lay++){
947 if(fgHitLay[lay]) nhitlay++;
948 }
949 if(nhitlay < m_param.nHitLayCut){
950 m_cut5++;
951 continue;
952 }
953
954 // nhit cut
955 if(nhitRec < m_param.nHitCut){
956 m_cut6++;
957 continue;
958 }
959
960// bool fgNoise = rectrk->getFgNoiseRatio();
961// if(m_param.noiseCut && (!fgNoise)) continue;
962// cout << setw(10) << iRun << setw(15) << iEvt << setw(5) << fgNoise << endl;
963
964// if(! ((fgHitLay[0]||fgHitLay[1]) && (fgHitLay[41]||fgHitLay[42])) ){
965// continue;
966// }
967
968 // calculate cell on the track
969 int cellTrkPass[43];
970 bool fgPass = getCellTrkPass(event, i, cellTrkPass);
971 for(lay=0; lay<m_nlayer; lay++){
972 fgAdd[lay] = 0;
973// if((16==lay) || (18==lay) || (19==lay) || (41==lay)){ // hv2200 2009-3
974 if((15==lay) || (16==lay) || (18==lay) || (19==lay) || (40==lay) || (41==lay) || (42==lay)){
975 int iCell = cellTrkPass[lay];
976 if(fgPass && (iCell >= 0) && m_fgGoodWire[lay][iCell]) m_effNtrk->Fill(lay);
977 else fgAdd[lay] = 1;
978 } else{
979 m_effNtrk->Fill(lay);
980 }
981 }
982
983 chisq = rectrk -> getChisq();
984 m_hchisq -> Fill( chisq );
985
986 if(pt > 0){
987 m_hpt -> Fill(pt);
988 m_hptPos -> Fill(pt);
989 m_hp -> Fill(p);
990 m_hp_cms -> Fill(p_cms);
991 m_hpPos -> Fill(p);
992 m_hpPoscms -> Fill(p_cms);
993 } else{
994 m_hpt -> Fill(-1.0*pt);
995 m_hptNeg -> Fill(-1.0*pt);
996 m_hp -> Fill(-1.0*p);
997 m_hp_cms -> Fill(-1.0*p_cms);
998 m_hpNeg -> Fill(-1.0*p);
999 m_hpNegcms -> Fill(-1.0*p_cms);
1000 }
1001 if(2 == ntrk){
1002 pTrk[i] = fabs(p);
1003 pTrkcms[i] = fabs(p_cms);
1004 }
1005
1006 dr = rectrk -> getDr();
1007 phi0 = rectrk -> getPhi0();
1008 kap = rectrk -> getKappa();
1009 dz = rectrk -> getDz();
1010 tanl = rectrk -> getTanLamda();
1011 lamda = atan(tanl);
1012 theta = HFPI - lamda;
1013
1014 m_hdr -> Fill(dr);
1015 m_hphi0 -> Fill(phi0);
1016 m_hkap -> Fill(kap);
1017 m_hdz -> Fill(dz);
1018 m_htanl -> Fill(tanl);
1019 m_hcosthe -> Fill(cos(theta));
1020 if(pt > 0) m_hcosthePos->Fill(cos(theta));
1021 else m_hcostheNeg->Fill(cos(theta));
1022
1023 philab = phi0 + HFPI;
1024 if(philab > PI2) philab -= PI2;
1025// cout << setw(15) << phi0 << setw(15) << philab << setw(15) << phicms << endl;
1026
1027 phiBin = (int)(philab / m_phiWid);
1028 phiBinCms = (int)(phicms / m_phiWid);
1029 theBin = (int)((cos(theta) + 1.0) / m_theWid);
1030 theBinCms = (int)((cos(thetacms) + 1.0) / m_theWid);
1031 if(phiBin < 0) phiBin = 0;
1032 if(phiBin >= NPhiBin) phiBin = NPhiBin-1;
1033 if(phiBinCms < 0) phiBinCms = 0;
1034 if(phiBinCms >= NPhiBin) phiBinCms = NPhiBin-1;
1035 if(theBin < 0) theBin = 0;
1036 if(theBin >= NThetaBin) theBin = NThetaBin-1;
1037 if(theBinCms < 0) theBinCms = 0;
1038 if(theBinCms >= NThetaBin) theBinCms = NThetaBin-1;
1039
1040 if(pt > 0){
1041 m_ppPhi[phiBin]->Fill(p);
1042 m_ppPhiCms[phiBinCms]->Fill(p_cms);
1043 m_ppThe[theBin]->Fill(p);
1044 m_ppTheCms[theBinCms]->Fill(p_cms);
1045 m_ppThePhi[theBin][phiBin]->Fill(p);
1046 m_ppThePhiCms[theBinCms][phiBinCms]->Fill(p_cms);
1047 } else{
1048 m_pnPhi[phiBin]->Fill(-1.0*p);
1049 m_pnPhiCms[phiBinCms]->Fill(-1.0*p_cms);
1050 m_pnThe[theBin]->Fill(-1.0*p);
1051 m_pnTheCms[theBinCms]->Fill(-1.0*p_cms);
1052 m_pnThePhi[theBin][phiBin]->Fill(-1.0*p);
1053 m_pnThePhiCms[theBinCms][phiBinCms]->Fill(-1.0*p_cms);
1054 }
1055
1056 x0 = dr * cos(phi0);
1057 y0 = dr * sin(phi0);
1058 m_hx0 -> Fill(x0);
1059 m_hy0 -> Fill(y0);
1060 if(m_nGrPoint < 10000){
1061 m_grX0Y0->SetPoint(m_nGrPoint, x0, y0);
1062 m_nGrPoint++;
1063 }
1064
1065 if(kap < 0) {
1066 zminus = dz;
1067 pm.push_back(p4);
1068 phibinm = phiBinCms;
1069 } else {
1070 zplus = dz;
1071 pp.push_back(p4);
1072 phibinp = phiBinCms;
1073 }
1074
1075// cout << "phi = " << setw(15) << philab << setw(15) << philab*180./3.14159 << setw(15) << p << endl;
1076 ntrkCal++;
1077 trkFlag[i] = true;
1078 nhitCal = 0;
1079 nhitCalInn = 0;
1080 nhitCalStp = 0;
1081 nhitCalOut = 0;
1082 for(k=0; k<nhitRec; k++){
1083 rechit = rectrk -> getRecHit(k);
1084
1085 lay = rechit -> getLayid();
1086 cel = rechit -> getCellid();
1087 lr = rechit -> getLR();
1088 stat = rechit -> getStat();
1089 doca = rechit -> getDocaExc();
1090 resiInc = rechit -> getResiIncLR();
1091 resiExc = rechit -> getResiExcLR();
1092 entr = rechit -> getEntra();
1093 tdr = rechit -> getTdrift();
1094 traw = (rechit -> getTdc()) * MdcCalTdcCnv;
1095 wir = m_mdcGeomSvc -> Wire(lay, cel) -> Id();
1096
1097 m_cel[lay] = (long)cel;
1098 m_lr[lay] = (long)lr;
1099 m_run[lay] = iRun;
1100 m_evt[lay] = iEvt;
1101 m_doca[lay] = doca;
1102 m_dm[lay] = rechit -> getDmeas();
1103 m_tdr[lay] = tdr;
1104 m_tdc[lay] = traw;
1105 m_entr[lay] = entr*180.0/3.14;
1106 m_zhit[lay] = rechit -> getZhit();
1107 m_qhit[lay] = rechit -> getQhit();
1108 m_p[lay] = p;
1109 m_pt[lay] = pt;
1110 m_phi0[lay] = phi0;
1111 m_tanl[lay] = tanl;
1112 qhit = rechit -> getQhit();
1113
1114 // calculating hitphi
1115 xx = (m_zhit[lay] - m_zw[wir]) * (m_xe[wir] - m_xw[wir]) /
1116 (m_ze[wir] - m_zw[wir]) + m_xw[wir];
1117 yy = (m_zhit[lay] - m_zw[wir]) * (m_ye[wir] - m_yw[wir]) /
1118 (m_ze[wir] - m_zw[wir]) + m_yw[wir];
1119 rr = sqrt( (xx * xx) + (yy * yy) );
1120 dphi = fabs(doca) / m_radii[lay];
1121
1122 if( yy >= 0 ) wphi = acos(xx / rr);
1123 else wphi = PI2 - acos(xx / rr);
1124 if(1 == lr) hitphi = wphi + dphi; // mention
1125 else hitphi = wphi - dphi;
1126 if(hitphi < 0) hitphi += PI2;
1127 else if(hitphi > PI2) hitphi -= PI2;
1128
1129 m_hitphi[lay] = hitphi;
1130
1131 if( (fabs(doca) > m_docaMax[lay]) ||
1132 (fabs(resiExc) > m_param.resiCut[lay]) ){
1133 continue;
1134 }
1135
1136 if(m_param.fgAdjacLayerCut){
1137 if(0 == lay){
1138 if( ! fgHitLay[1] ) continue;
1139 } else if(42 == lay){
1140 if( ! fgHitLay[41] ) continue;
1141 } else{
1142 if( (!fgHitLay[lay-1]) && (!fgHitLay[lay+1]) ) continue;
1143
1144 // for boundary layers
1145 if( m_param.fgBoundLayerCut && m_layBound[lay] &&
1146 ((!fgHitLay[lay-1]) || (!fgHitLay[lay+1])) ) continue;
1147 }
1148 }
1149
1150 if((1 == m_param.hitStatCut) && (1 != stat)) continue;
1151
1152 // fill xtplot tree
1153 if((1 == m_param.fillNtuple) && (m_nEvtNtuple < m_param.nEvtNtuple)){
1154 m_xtTuple[lay] -> write();
1155 }
1156
1157 if(1 == m_param.hitStatCut){
1158 if( (0 == fgAdd[lay]) && (1 == stat) ){
1159 m_effNtrkRecHit->Fill(lay);
1160 fgAdd[lay] = 1;
1161 }
1162 } else{
1163 if(0 == fgAdd[lay]){
1164 m_effNtrkRecHit->Fill(lay);
1165 fgAdd[lay] = 1;
1166 }
1167 }
1168
1169 nhitCal++;
1170 if(lay < 8) nhitCalInn++;
1171 else if(lay < 20) nhitCalStp++;
1172 else nhitCalOut++;
1173
1174 m_wirehitmap -> Fill(wir);
1175 m_layerhitmap -> Fill(lay);
1176
1177 m_htraw[lay] -> Fill(traw);
1178 m_htdr[lay] -> Fill(tdr);
1179 m_htdrlr[lay][lr]->Fill(tdr);
1180 m_hadc[lay] -> Fill(qhit);
1181
1182 m_hresAllInc -> Fill(resiInc);
1183 m_hresAllExc -> Fill(resiExc);
1184 double resiAve = 0.5 * (resiInc + resiExc);
1185 m_hresAllAve -> Fill(resiAve);
1186
1187 if(lay < 8){
1188 m_hresInnInc -> Fill(resiInc);
1189 m_hresInnExc -> Fill(resiExc);
1190 } else if(lay < 20){
1191 m_hresStpInc -> Fill(resiInc);
1192 m_hresStpExc -> Fill(resiExc);
1193 } else{
1194 m_hresOutInc -> Fill(resiInc);
1195 m_hresOutExc -> Fill(resiExc);
1196 }
1197
1198 int qbin = (int)((qhit-100.0)/100.0);
1199 if(qbin>=0 && qbin<14){
1200 m_hresAveAllQ[qbin]->Fill(resiAve);
1201 m_hresAveLayQ[lay][qbin]->Fill(resiAve);
1202 if(lay > 7) m_hresAveOutQ[qbin]->Fill(resiAve);
1203 }
1204
1205 m_hresInc[lay] -> Fill(resiInc);
1206 m_hreslrInc[lay][lr]->Fill(resiInc);
1207 m_hresExc[lay] -> Fill(resiExc);
1208 m_hreslrExc[lay][lr]->Fill(resiExc);
1209 m_hresAve[lay] -> Fill(resiAve);
1210 m_hreslrAve[lay][lr]->Fill(resiAve);
1211
1212 int iPhi = (int)(hitphi*20.0/PI2);
1213 if(iPhi>=20) iPhi = 19;
1214 m_hresphi[lay][iPhi]->Fill((resiInc+resiExc)*0.5);
1215
1216// bin = (int)(fabs(doca) / m_dwid);
1217 bin = (int)(fabs(rechit->getDmeas()) / m_dwid);
1218 iEntr = m_mdcFunSvc -> getSdEntrIndex(entr);
1219 if(1 == m_nEntr[lay]){
1220 iEntr = 0;
1221 } else if(2 == m_nEntr[lay]){
1222 if(entr > 0.0) iEntr = 1;
1223 else iEntr = 0;
1224 }
1225 if((iEntr < MdcCalNENTRSD) && (bin < MdcCalSdNBIN)){
1226 key = getHresId(lay, iEntr, lr, bin);
1227 if( 1 == m_mapr2d.count(key) ){
1228 hid = m_mapr2d[key];
1229 m_hr2dInc[hid] -> Fill(resiInc);
1230 m_hr2dExc[hid] -> Fill(resiExc);
1231 }
1232 }
1233
1234 if((tdr>0) && (tdr<750)){
1235 if(tdr<300) bin = (int)(tdr/10.0);
1236 else bin = (int)((tdr-300.0)/30.0) + 29;
1237 m_hr2t[lay][iEntr][lr][bin]->Fill(resiExc);
1238 }
1239 } // loop of nhits
1240 m_nEvtNtuple++;
1241 m_hnhitsCal->Fill(nhitCal);
1242 m_hnhitsCalInn->Fill(nhitCalInn);
1243 m_hnhitsCalStp->Fill(nhitCalStp);
1244 m_hnhitsCalOut->Fill(nhitCalOut);
1245 } // end of track loop
1246 m_hnTrkCal->Fill(ntrkCal);
1247 if(2 == ntrkCal){
1248 if(pTrk[0] > pTrk[1]) m_hpMax->Fill(pTrk[0]);
1249 else m_hpMax->Fill(pTrk[1]);
1250
1251 if(pTrkcms[0] > pTrkcms[1]) m_hpMaxCms->Fill(pTrkcms[0]);
1252 else m_hpMaxCms->Fill(pTrkcms[1]);
1253 }
1254 if(ntrkCal > 0) m_hTesCalUse->Fill(tes);
1255
1256 double delZ0;
1257 if((fabs(zminus) < 9000.0) && (fabs(zplus) < 9000.0)) delZ0 = zplus - zminus;
1258 m_hdelZ0 -> Fill(delZ0);
1259
1260 if (1 == pp.size() * pm.size()){
1261 HepLorentzVector ptot = pp[0] + pm[0];
1262 bool fourmomcut = false;
1263// fourmomcut = (ptot.x()>0.02 && ptot.x()<0.06) && (fabs(ptot.y()) < 0.02)
1264// && (ptot.z()>-0.01 && ptot.z()<0.03) && (ptot.e()>3.4 && ptot.e()<4.0);
1265 fourmomcut = (fabs(ptot.x()-0.04)<0.026) && (fabs(ptot.y()) < 0.026)
1266 && (fabs(ptot.z()-0.005)<0.036) && (fabs(ptot.e()-ecm)<0.058);
1267 //cout << "x = " << ptot.x() << ", y = " << ptot.y() << ", z = " << ptot.z() << ", e = " << ptot.e() << endl;
1268 if (fourmomcut) {
1269 HepLorentzVector psip(xboost, yboost, zboost, ecm);
1270 Hep3Vector boostv = psip.boostVector();
1271 pp[0].boost(- boostv);
1272 pm[0].boost(- boostv);
1273 m_hp_cut->Fill(pp[0].rho());
1274 m_hp_cut->Fill(pm[0].rho());
1275 }
1276 }
1277
1278 if(2==ntrk) for(i=0; i<ntrk; i++) pTrk[i] = (event -> getRecTrk(i)) -> getP();
1279 if((5==m_param.particle) && (2==ntrk) && (fabs(pTrk[0])<5) && (fabs(pTrk[1])<5)){
1280// if(1==ntrk) p = (event->getRecTrk(0)) -> getP();
1281// if((5==m_param.particle) && (1==ntrk) && (fabs(p)<5)){
1282 m_tescos = tes;
1283 m_tesFlagcos = esTimeflag;
1284 for(i=0; i<ntrk; i++){
1285 rectrk = event -> getRecTrk(i);
1286 phi0 = rectrk -> getPhi0();
1287 phi0 = ((phi0+HFPI) > PI2) ? (phi0+HFPI-PI2) : (phi0+HFPI);
1288
1289 tanl = rectrk -> getTanLamda();
1290 lamda = atan(tanl);
1291 theta = HFPI - lamda;
1292
1293 if(phi0 < (2.0*HFPI)){
1294 m_nhitUpcos = rectrk -> getNHits();
1295 m_pUpcos = rectrk -> getP();
1296 m_ptUpcos = rectrk -> getPt();
1297 m_phiUpcos = phi0;
1298 m_drUpcos = rectrk->getDr();
1299 m_dzUpcos = rectrk->getDz();
1300 m_ctheUpcos = cos(theta);
1301 } else{
1302 m_nhitDwcos = rectrk -> getNHits();
1303 m_pDwcos = rectrk -> getP();
1304 m_ptDwcos = rectrk -> getPt();
1305 m_phiDwcos = phi0;
1306 m_drDwcos = rectrk->getDr();
1307 m_dzDwcos = rectrk->getDz();
1308 m_ctheDwcos = cos(theta);
1309
1310 if(m_pDwcos > 0) m_chargecos = 1;
1311 else m_chargecos = 0;
1312 }
1313 }
1314 m_cosmic->write();
1315 }
1316
1317 if(1 == m_param.fgCalDetEffi) calDetEffi();
1318
1319 return 1;
1320}
mass
double mass
Definition: CosmicGenerator.cxx:138
Fill
curve Fill()
p_cms
const HepLorentzVector p_cms(0.034067, 0.0, 0.0, 3.097)
Vp4
std::vector< HepLorentzVector > Vp4
Definition: Gam4pikp.cxx:53
sin
double sin(const BesAngle a)
Definition: InstallArea/include/MdcGeom/MdcGeom/BesAngle.h:210
cos
double cos(const BesAngle a)
Definition: InstallArea/include/MdcGeom/MdcGeom/BesAngle.h:213
key
*************DOUBLE PRECISION m_pi *DOUBLE PRECISION m_HvecTau2 DOUBLE PRECISION m_HvClone2 DOUBLE PRECISION m_gamma1 DOUBLE PRECISION m_gamma2 DOUBLE PRECISION m_thet1 DOUBLE PRECISION m_thet2 INTEGER m_IFPHOT *COMMON c_Taupair $ !Spin Polarimeter vector first Tau $ !Spin Polarimeter vector second Tau $ !Clone Spin Polarimeter vector first Tau $ !Clone Spin Polarimeter vector second Tau $ !Random Euler angle for cloning st tau $ !Random Euler angle for cloning st tau $ !Random Euler angle for cloning st tau $ !Random Euler angle for cloning nd tau $ !Random Euler angle for cloning nd tau $ !Random Euler angle for cloning nd tau $ !phi of HvecTau1 $ !theta of HvecTau1 $ !phi of HvecTau2 $ !theta of HvecTau2 $ !super key
Definition: Taupair.h:42
IMdcCalibFunSvc::getWireEff
virtual double getWireEff(int layid, int cellid) const =0
IMdcGeomSvc::Layer
virtual const MdcGeoLayer *const Layer(unsigned id)=0

Referenced by GrXtMdcCalib::fillHist(), QtMdcCalib::fillHist(), T0MdcCalib::fillHist(), Wr2dMdcCalib::fillHist(), WrMdcCalib::fillHist(), XtInteMdcCalib::fillHist(), and XtMdcCalib::fillHist().

◆ fillHist() [2/2]

virtual int MdcCalib::fillHist ( MdcCalEvent *  event)
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

◆ initialize() [1/2]

void MdcCalib::initialize ( TObjArray *  hlist,
IMdcGeomSvc *  mdcGeomSvc,
IMdcCalibFunSvc *  mdcFunSvc,
IMdcUtilitySvc *  mdcUtilitySvc 
)
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

Definition at line 202 of file MdcCalib.cxx.

203 {
204 IMessageSvc* msgSvc;
205 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
206 MsgStream log(msgSvc, "MdcCalib");
207 log << MSG::INFO << "MdcCalib::initialize()" << endreq;
208
209 m_hlist = hlist;
210 m_mdcGeomSvc = mdcGeomSvc;
211 m_mdcFunSvc = mdcFunSvc;
212 m_mdcUtilitySvc = mdcUtilitySvc;
213
214 int lay;
215 int iEntr;
216 int lr;
217 int bin;
218 char hname[200];
219
220 m_nlayer = m_mdcGeomSvc -> getLayerSize();
221
222 for(lay=0; lay<m_nlayer; lay++){
223 m_radii[lay] = m_mdcGeomSvc->Layer(lay)->Radius();
224 }
225 ofstream fwpc("wirelog.txt");
226 for(int wir=0; wir<MdcCalTotCell; wir++){
227 m_xe[wir] = m_mdcGeomSvc->Wire(wir)->Backward().x();
228 m_ye[wir] = m_mdcGeomSvc->Wire(wir)->Backward().y();
229 m_ze[wir] = m_mdcGeomSvc->Wire(wir)->Backward().z();
230 m_xw[wir] = m_mdcGeomSvc->Wire(wir)->Forward().x();
231 m_yw[wir] = m_mdcGeomSvc->Wire(wir)->Forward().y();
232 m_zw[wir] = m_mdcGeomSvc->Wire(wir)->Forward().z();
233 fwpc << setw(6) << wir << setw(15) << m_xe[wir] << setw(15) << m_ye[wir]
234 << setw(15) << m_ze[wir] << setw(15) << m_xw[wir]
235 << setw(15) << m_yw[wir] << setw(15) << m_zw[wir] << endl;
236 }
237 fwpc.close();
238
239 m_fdcom = new TFolder("common", "common");
240 m_hlist -> Add(m_fdcom);
241
242 m_effNtrk = new TH1F("effNtrk", "", 43, -0.5, 42.5);
243 m_fdcom->Add(m_effNtrk);
244
245 m_effNtrkRecHit = new TH1F("effNtrkRecHit", "", 43, -0.5, 42.5);
246 m_fdcom->Add(m_effNtrkRecHit);
247
248 m_hresAllInc = new TH1F("HResAllInc", "", 200, -1.0, 1.0);
249 m_fdcom -> Add(m_hresAllInc);
250
251 m_hresAllExc = new TH1F("HResAllExc", "", 200, -1.0, 1.0);
252 m_fdcom -> Add(m_hresAllExc);
253
254 m_hresAllAve = new TH1F("HResAllAve", "", 200, -1.0, 1.0);
255 m_fdcom -> Add(m_hresAllAve);
256
257 m_hresInnInc = new TH1F("HResInnInc", "", 200, -1.0, 1.0);
258 m_fdcom -> Add(m_hresInnInc);
259
260 m_hresInnExc = new TH1F("HResInnExc", "", 200, -1.0, 1.0);
261 m_fdcom -> Add(m_hresInnExc);
262
263 m_hresStpInc = new TH1F("HResStpInc", "", 200, -1.0, 1.0);
264 m_fdcom -> Add(m_hresStpInc);
265
266 m_hresStpExc = new TH1F("HResStpExc", "", 200, -1.0, 1.0);
267 m_fdcom -> Add(m_hresStpExc);
268
269 m_hresOutInc = new TH1F("HResOutInc", "", 200, -1.0, 1.0);
270 m_fdcom -> Add(m_hresOutInc);
271
272 m_hresOutExc = new TH1F("HResOutExc", "", 200, -1.0, 1.0);
273 m_fdcom -> Add(m_hresOutExc);
274
275 m_fdResQ = new TFolder("ResQ", "ResQ");
276 m_hlist->Add(m_fdResQ);
277 for(int i=0; i<14; i++){
278 sprintf(hname, "resoAll_qbin%02d", i);
279 m_hresAveAllQ[i] = new TH1F(hname, "", 200, -1, 1);
280 m_fdResQ->Add(m_hresAveAllQ[i]);
281
282 sprintf(hname, "resoOut_qbin%02d", i);
283 m_hresAveOutQ[i] = new TH1F(hname, "", 200, -1, 1);
284 m_fdResQ->Add(m_hresAveOutQ[i]);
285 }
286 for(lay=0; lay<43; lay++){
287 for(int i=0; i<14; i++){
288 sprintf(hname, "resoLay%02d_qbin%02d", lay, i);
289 m_hresAveLayQ[lay][i] = new TH1F(hname, "", 200, -1, 1);
290 m_fdResQ->Add(m_hresAveLayQ[lay][i]);
291 }
292 }
293
294 for(int iEs=0; iEs<m_param.nEsFlag; iEs++){
295 sprintf(hname, "Tes_%d", m_param.esFlag[iEs]);
296 m_hTes[iEs] = new TH1F(hname, "", 750, 0, 1500);
297 m_fdcom->Add(m_hTes[iEs]);
298 }
299
300 m_hbbTrkFlg = new TH1F("BbTrkFlg", "", 100, 0, 6);
301 m_fdcom -> Add(m_hbbTrkFlg);
302
303 m_hTesAll = new TH1F("TesAll", "", 1000, 0, 2000);
304 m_fdcom -> Add(m_hTesAll);
305
306 m_hTesGood = new TH1F("TesGood", "", 1000, 0, 2000);
307 m_fdcom -> Add(m_hTesGood);
308
309 m_hTesAllFlag = new TH1F("TesAllFlag", "", 300, -0.5, 299.5);
310 m_fdcom -> Add(m_hTesAllFlag);
311
312 m_hTesRec = new TH1F("TesRec", "", 1000, 0, 2000);
313 m_fdcom -> Add(m_hTesRec);
314
315 m_hTesCalFlag = new TH1F("TesCalFlag", "", 1000, 0, 2000);
316 m_fdcom -> Add(m_hTesCalFlag);
317
318 m_hTesCalUse = new TH1F("TesCalUse", "", 1000, 0, 2000);
319 m_fdcom -> Add(m_hTesCalUse);
320
321 m_hnRawHit = new TH1F("nRawHit", "", 6797, -0.5, 6796.5);
322 m_fdcom -> Add(m_hnRawHit);
323
324 m_hpt = new TH1F("HPt", "", 800, 0, 3);
325 m_fdcom -> Add(m_hpt);
326
327 m_hptPos = new TH1F("HPtPos", "", 800, 0, 3);
328 m_fdcom -> Add(m_hptPos);
329
330 m_hptNeg = new TH1F("HPtNeg", "", 800, 0, 3);
331 m_fdcom -> Add(m_hptNeg);
332
333 m_hp = new TH1F("HP", "", 800, 0, 3);
334 m_fdcom -> Add(m_hp);
335
336 m_hp_cms = new TH1F("HPCMS", "", 800, 0, 3);
337 m_fdcom -> Add(m_hp_cms);
338
339 m_hpMax = new TH1F("HPMax", "", 800, 0, 3);
340 m_fdcom -> Add(m_hpMax);
341
342 m_hpMaxCms = new TH1F("HPMax_Cms", "", 800, 0, 3);
343 m_fdcom -> Add(m_hpMaxCms);
344
345 m_hpPos = new TH1F("HP_Pos", "", 800, 0, 3);
346 m_fdcom -> Add(m_hpPos);
347
348 m_hpNeg = new TH1F("HP_Neg", "", 800, 0, 3);
349 m_fdcom -> Add(m_hpNeg);
350
351 m_hpPoscms = new TH1F("HP_Pos_cms", "", 800, 0, 3);
352 m_fdcom -> Add(m_hpPoscms);
353
354 m_hpNegcms = new TH1F("HP_Neg_cms", "", 800, 0, 3);
355 m_fdcom -> Add(m_hpNegcms);
356
357 m_hp_cut = new TH1F("HPCut", "", 800, 0, 3);
358 m_fdcom -> Add(m_hp_cut);
359
360 m_hchisq = new TH1F("Chisq", "", 10, 0, 100);
361 m_fdcom -> Add(m_hchisq);
362
363 m_hnTrk = new TH1F("HNtrack", "HNtrack", 10, -0.5, 9.5);
364 m_fdcom -> Add(m_hnTrk);
365
366 m_hnTrkCal = new TH1F("HNtrackCal", "HNtrackCal", 10, -0.5, 9.5);
367 m_fdcom -> Add(m_hnTrkCal);
368
369 m_hnhitsRec = new TH1F("HNhitsRec", "", 100, -0.5, 99.5);
370 m_fdcom -> Add(m_hnhitsRec);
371
372 m_hnhitsRecInn = new TH1F("HNhitsInnRec", "", 60, 0.5, 60.5);
373 m_fdcom -> Add(m_hnhitsRecInn);
374
375 m_hnhitsRecStp = new TH1F("HNhitsStpRec", "", 60, 0.5, 60.5);
376 m_fdcom -> Add(m_hnhitsRecStp);
377
378 m_hnhitsRecOut = new TH1F("HNhitsOutRec", "", 60, 0.5, 60.5);
379 m_fdcom -> Add(m_hnhitsRecOut);
380
381 m_hnhitsCal = new TH1F("HNhitsCal", "", 100, -0.5, 99.5);
382 m_fdcom -> Add(m_hnhitsCal);
383
384 m_hnhitsCalInn = new TH1F("HNhitsCalInn", "", 60, 0.5, 60.5);
385 m_fdcom -> Add(m_hnhitsCalInn);
386
387 m_hnhitsCalStp = new TH1F("HNhitsCalStp", "", 60, 0.5, 60.5);
388 m_fdcom -> Add(m_hnhitsCalStp);
389
390 m_hnhitsCalOut = new TH1F("HNhitsCalOut", "", 60, 0.5, 60.5);
391 m_fdcom -> Add(m_hnhitsCalOut);
392
393 m_wirehitmap = new TH1F("Wire_HitMap", "Wire_HitMap", 6796, -0.5, 6795.5);
394 m_fdcom -> Add(m_wirehitmap);
395
396 m_layerhitmap = new TH1F("Layer_HitMap", "Layer_HitMap", 43, -0.5, 42.5);
397 m_fdcom -> Add(m_layerhitmap);
398
399 m_hnoisephi = new TH1F("phi_noise", "", 100, 0, 6.284);
400 m_fdcom -> Add(m_hnoisephi);
401
402 m_hnoiselay = new TH1F("Layer_noise", "Layer_noise", 43, -0.5, 42.5);
403 m_fdcom -> Add(m_hnoiselay);
404
405 m_hnoisenhits = new TH1F("nhits_noise", "nhits_noise", 6796, -0.5, 6795.5);
406 m_fdcom -> Add(m_hnoisenhits);
407
408 m_hratio = new TH1F("ratio", "", 100, 0, 1);
409 m_fdcom -> Add(m_hratio);
410
411 m_hdr = new TH1F("dr", "", 500, -500, 500);
412 m_fdcom -> Add(m_hdr);
413
414 m_hphi0 = new TH1F("phi0", "", 100, 0, 6.284);
415 m_fdcom -> Add(m_hphi0);
416
417 m_hkap = new TH1F("kappa", "", 400, -50, 50);
418 m_fdcom -> Add(m_hkap);
419
420 m_hdz = new TH1F("dz", "", 500, -1000, 1000);
421 m_fdcom -> Add(m_hdz);
422
423 m_htanl = new TH1F("tanl", "", 200, -5, 5);
424 m_fdcom -> Add(m_htanl);
425
426 m_hcosthe = new TH1F("costheta", "", 200, -1, 1);
427 m_fdcom -> Add(m_hcosthe);
428
429 m_hcostheNeg = new TH1F("costhetaNeg", "", 200, -1, 1);
430 m_fdcom -> Add(m_hcostheNeg);
431
432 m_hcosthePos = new TH1F("costhetaPos", "", 200, -1, 1);
433 m_fdcom -> Add(m_hcosthePos);
434
435 m_hx0 = new TH1F("x0", "", 100, -10, 10);
436 m_fdcom -> Add(m_hx0);
437
438 m_hy0 = new TH1F("y0", "", 100, -10, 10);
439 m_fdcom -> Add(m_hy0);
440
441 m_hdelZ0 = new TH1F("delta_z0", "", 100, -50, 50);
442 m_fdcom -> Add(m_hdelZ0);
443
444 m_grX0Y0 = new TGraph();
445 m_grX0Y0->SetName("x0y0");
446 m_fdcom -> Add(m_grX0Y0);
447
448 m_hitEffAll = new TH1F("hitEffAll", "", 6800, -0.5, 6799.5);
449 m_fdcom -> Add(m_hitEffAll);
450
451 m_hitEffRaw = new TH1F("hitEffRaw", "", 6800, -0.5, 6799.5);
452 m_fdcom -> Add(m_hitEffRaw);
453
454 m_hitEffRec = new TH1F("hitEffRec", "", 6800, -0.5, 6799.5);
455 m_fdcom -> Add(m_hitEffRec);
456
457 // histograms for drift time
458 m_fdTime = new TFolder("time", "time");
459 m_hlist -> Add(m_fdTime);
460
461 for(lay=0; lay<m_nlayer; lay++){
462 sprintf(hname, "Traw%02d", lay);
463 m_htraw[lay] = new TH1F(hname, "", 1000, 0, 2000);
464 m_fdTime -> Add(m_htraw[lay]);
465
466 sprintf(hname, "Tdr%02d", lay);
467 m_htdr[lay] = new TH1F(hname, "", 510, -10, 500);
468 m_fdTime -> Add(m_htdr[lay]);
469
470 for (lr=0; lr<2; lr++){
471 sprintf(hname, "Tdr%02d_lr%01d", lay, lr);
472 m_htdrlr[lay][lr] = new TH1F(hname, "", 510, -10, 500);
473 m_fdTime -> Add(m_htdrlr[lay][lr]);
474 }
475 }
476
477 // histograms of adc
478 m_fdAdc = new TFolder("adc", "adc");
479 m_hlist -> Add(m_fdAdc);
480
481 for(lay=0; lay<m_nlayer; lay++){
482 sprintf(hname, "adc%02d", lay);
483 m_hadc[lay] = new TH1F(hname, "", 1500, 0, 3000);
484 m_fdAdc -> Add(m_hadc[lay]);
485 }
486 // histograms for resolution
487 m_fdres = new TFolder("resolution", "resolution");
488 m_hlist -> Add(m_fdres);
489
490 m_fdresAve = new TFolder("resAve", "resAve");
491 m_hlist -> Add(m_fdresAve);
492
493 for(lay=0; lay<m_nlayer; lay++){
494 sprintf(hname, "Reso%02dInc", lay);
495 m_hresInc[lay] = new TH1F(hname, "", 1000, -5, 5);
496 m_fdres -> Add(m_hresInc[lay]);
497
498 sprintf(hname, "Reso%02dExc", lay);
499 m_hresExc[lay] = new TH1F(hname, "", 1000, -5, 5);
500 m_fdres -> Add(m_hresExc[lay]);
501
502 sprintf(hname, "Reso%02d", lay);
503 m_hresAve[lay] = new TH1F(hname, "", 1000, -5, 5);
504 m_fdresAve -> Add(m_hresAve[lay]);
505
506 for (lr=0; lr<2; lr++){
507 sprintf(hname, "Reso%02dInc_lr%01d", lay, lr);
508// m_hreslrInc[lay][lr] = new TH1F(hname, "", 200, -1, 1);
509 m_hreslrInc[lay][lr] = new TH1F(hname, "", 1000, -5, 5);
510 m_fdres->Add(m_hreslrInc[lay][lr]);
511
512 sprintf(hname, "Reso%02dExc_lr%01d", lay, lr);
513// m_hreslrExc[lay][lr] = new TH1F(hname, "", 200, -1, 1);
514 m_hreslrExc[lay][lr] = new TH1F(hname, "", 1000, -5, 5);
515 m_fdres->Add(m_hreslrExc[lay][lr]);
516
517 sprintf(hname, "Reso%02d_lr%01d", lay, lr);
518// m_hreslrAve[lay][lr] = new TH1F(hname, "", 200, -1, 1);
519 m_hreslrAve[lay][lr] = new TH1F(hname, "", 1000, -5, 5);
520 m_fdresAve->Add(m_hreslrAve[lay][lr]);
521 }
522 for(int phi=0; phi<20; phi++){
523 sprintf(hname, "ResoPhi%02d_phi%02d", lay, phi);
524 m_hresphi[lay][phi] = new TH1F(hname, "", 200, -1, 1);
525 m_fdres->Add(m_hresphi[lay][phi]);
526 }
527 }
528
529 /* histograms for momentum vs phi */
530 m_fdmomPhi = new TFolder("momPhi", "momPhi");
531 m_hlist -> Add(m_fdmomPhi);
532
533 int thbin;
534 for(bin=0; bin<NPhiBin; bin++){
535 sprintf(hname, "hPpos_phi%02d", bin);
536 m_ppPhi[bin] = new TH1F(hname, "", 400, 1.0, 2.5);
537 m_fdmomPhi->Add(m_ppPhi[bin]);
538
539 sprintf(hname, "hPneg_phi%02d", bin);
540 m_pnPhi[bin] = new TH1F(hname, "", 400, 1.0, 2.5);
541 m_fdmomPhi->Add(m_pnPhi[bin]);
542
543 sprintf(hname, "hPpos_phi_cms%02d", bin);
544 m_ppPhiCms[bin] = new TH1F(hname, "", 400, 1.0, 2.5);
545 m_fdmomPhi->Add(m_ppPhiCms[bin]);
546
547 sprintf(hname, "hPneg_phi_cms%02d", bin);
548 m_pnPhiCms[bin] = new TH1F(hname, "", 400, 1.0, 2.5);
549 m_fdmomPhi->Add(m_pnPhiCms[bin]);
550
551 for(thbin=0; thbin<NThetaBin; thbin++){
552 sprintf(hname, "hPpos_theta%02d_phi%02d", thbin, bin);
553 m_ppThePhi[thbin][bin] = new TH1F(hname, "", 400, 1.0, 2.5);
554 m_fdmomPhi->Add(m_ppThePhi[thbin][bin]);
555
556 sprintf(hname, "hPneg_theta%02d_phi%02d", thbin, bin);
557 m_pnThePhi[thbin][bin] = new TH1F(hname, "", 400, 1.0, 2.5);
558 m_fdmomPhi->Add(m_pnThePhi[thbin][bin]);
559
560 sprintf(hname, "hPposCms_theta%02d_phi%02d", thbin, bin);
561 m_ppThePhiCms[thbin][bin] = new TH1F(hname, "", 400, 1.0, 2.5);
562 m_fdmomPhi->Add(m_ppThePhiCms[thbin][bin]);
563
564 sprintf(hname, "hPnegCms_theta%02d_phi%02d", thbin, bin);
565 m_pnThePhiCms[thbin][bin] = new TH1F(hname, "", 400, 1.0, 2.5);
566 m_fdmomPhi->Add(m_pnThePhiCms[thbin][bin]);
567 }
568 }
569 for(thbin=0; thbin<NThetaBin; thbin++){
570 sprintf(hname, "hPpos_the%02d", thbin);
571 m_ppThe[thbin] = new TH1F(hname, "", 400, 1.0, 2.5);
572 m_fdmomPhi->Add(m_ppThe[thbin]);
573
574 sprintf(hname, "hPneg_the%02d", thbin);
575 m_pnThe[thbin] = new TH1F(hname, "", 400, 1.0, 2.5);
576 m_fdmomPhi->Add(m_pnThe[thbin]);
577
578 sprintf(hname, "hPposCms_the%02d", thbin);
579 m_ppTheCms[thbin] = new TH1F(hname, "", 400, 1.0, 2.5);
580 m_fdmomPhi->Add(m_ppTheCms[thbin]);
581
582 sprintf(hname, "hPnegCms_the%02d", thbin);
583 m_pnTheCms[thbin] = new TH1F(hname, "", 400, 1.0, 2.5);
584 m_fdmomPhi->Add(m_pnTheCms[thbin]);
585 }
586
587 // histograms for resolution vs distance
588 m_fdres2d = new TFolder("res2d", "res2d");
589 m_hlist -> Add(m_fdres2d);
590
591 int hid = 0;
592 int key;
593 TH1F* hist;
594 for(lay=0; lay<m_nlayer; lay++){
595 for(iEntr=0; iEntr<MdcCalNENTRSD; iEntr++){
596 for(lr=0; lr<2; lr++){
597 for(bin=0; bin<MdcCalSdNBIN; bin++){
598 sprintf(hname, "r2d%02d_%02d_%01d_%02dInc", lay, iEntr, lr, bin);
599 hist = new TH1F(hname, "", 200, -1, 1);
600 m_hr2dInc.push_back(hist);
601 m_fdres2d -> Add(hist);
602
603 sprintf(hname, "r2d%02d_%02d_%01d_%02dExc", lay, iEntr, lr, bin);
604 hist = new TH1F(hname, "", 200, -1, 1);
605 m_hr2dExc.push_back(hist);
606 m_fdres2d -> Add(hist);
607
608 key = getHresId(lay, iEntr, lr, bin);
609 m_mapr2d.insert( valType3(key, hid) );
610 hid++;
611 }
612 }
613 }
614 } // end of layer loop
615
616 m_fdres2t = new TFolder("res2t", "res2t");
617// m_hlist -> Add(m_fdres2t);
618
619 for(lay=0; lay<m_nlayer; lay++){
620 for(iEntr=0; iEntr<MdcCalNENTRXT; iEntr++){
621 for(lr=0; lr<2; lr++){
622 for(bin=0; bin<45; bin++){
623 sprintf(hname, "r2t%02d_%02d_%01d_%02d", lay, iEntr, lr, bin);
624 m_hr2t[lay][iEntr][lr][bin] = new TH1F(hname, "", 600, -3, 3);
625 m_fdres2t -> Add(m_hr2t[lay][iEntr][lr][bin]);
626 }
627 }
628 }
629 }
630
631 INTupleSvc* ntupleSvc;
632 Gaudi::svcLocator() -> service("NTupleSvc", ntupleSvc);
633 for(lay=0; lay<m_nlayer; lay++){
634 sprintf(hname, "FILE136/xt%02d", lay);
635 NTuplePtr nt(ntupleSvc, hname);
636 if ( nt ) m_xtTuple[lay] = nt;
637 else{
638 m_xtTuple[lay] = ntupleSvc->book(hname, CLID_ColumnWiseTuple, "MdcXtNtuple");
639 if( m_xtTuple[lay] ){
640 m_xtTuple[lay]->addItem("cel", m_cel[lay]);
641 m_xtTuple[lay]->addItem("lr", m_lr[lay]);
642 m_xtTuple[lay]->addItem("run", m_run[lay]);
643 m_xtTuple[lay]->addItem("evt", m_evt[lay]);
644 m_xtTuple[lay]->addItem("doca", m_doca[lay]);
645 m_xtTuple[lay]->addItem("dm", m_dm[lay]);
646 m_xtTuple[lay]->addItem("tdr", m_tdr[lay]);
647 m_xtTuple[lay]->addItem("tdc", m_tdc[lay]);
648 m_xtTuple[lay]->addItem("entr", m_entr[lay]);
649 m_xtTuple[lay]->addItem("zhit", m_zhit[lay]);
650 m_xtTuple[lay]->addItem("qhit", m_qhit[lay]);
651 m_xtTuple[lay]->addItem("p", m_p[lay]);
652 m_xtTuple[lay]->addItem("pt", m_pt[lay]);
653 m_xtTuple[lay]->addItem("phi0", m_phi0[lay]);
654 m_xtTuple[lay]->addItem("tanl", m_tanl[lay]);
655 m_xtTuple[lay]->addItem("hitphi", m_hitphi[lay]);
656 } else{
657 log << MSG::FATAL << "Cannot book Xt N-tuple:"
658 << long(m_xtTuple[lay]) << endreq;
659 }
660 }
661 }
662
663 if(5==m_param.particle){
664 sprintf(hname, "FILE136/cosmic");
665 NTuplePtr nt(ntupleSvc, hname);
666 if ( nt ) m_cosmic = nt;
667 else{
668 m_cosmic = ntupleSvc->book(hname, CLID_ColumnWiseTuple, "MdcXtNtuple");
669 if( m_cosmic ){
670 m_cosmic->addItem("pUp", m_pUpcos);
671 m_cosmic->addItem("pDw", m_pDwcos);
672 m_cosmic->addItem("ptUp", m_ptUpcos);
673 m_cosmic->addItem("ptDw", m_ptDwcos);
674 m_cosmic->addItem("phiUp", m_phiUpcos);
675 m_cosmic->addItem("phiDw", m_phiDwcos);
676 m_cosmic->addItem("drUp", m_drUpcos);
677 m_cosmic->addItem("drDw", m_drDwcos);
678 m_cosmic->addItem("dzUp", m_dzUpcos);
679 m_cosmic->addItem("dzDw", m_dzDwcos);
680 m_cosmic->addItem("ctheUp", m_ctheUpcos);
681 m_cosmic->addItem("ctheDw", m_ctheDwcos);
682 m_cosmic->addItem("nhitUp", m_nhitUpcos);
683 m_cosmic->addItem("nhitDw", m_nhitDwcos);
684 m_cosmic->addItem("char", m_chargecos);
685 m_cosmic->addItem("tesfg", m_tesFlagcos);
686 m_cosmic->addItem("tes", m_tescos);
687 }
688 }
689 }
690}
valType3
map< int, int >::value_type valType3
Definition: MdcCalib.cxx:33
IMdcGeomSvc::Wire
virtual const MdcGeoWire *const Wire(unsigned id)=0
sprintf
sprintf(cut,"kal_costheta0_em>-0.93&&kal_costheta0_em<0.93&&kal_pxy0_em>=0.05+%d*0.1&&kal_pxy0_em<0.15+%d*0.1&&NGch>=2", j, j)
Add
mg Add(gr3)

Referenced by GrXtMdcCalib::initialize(), QtMdcCalib::initialize(), T0MdcCalib::initialize(), Wr2dMdcCalib::initialize(), WrMdcCalib::initialize(), XtInteMdcCalib::initialize(), and XtMdcCalib::initialize().

◆ initialize() [2/2]

virtual void MdcCalib::initialize ( TObjArray *  hlist,
IMdcGeomSvc *  mdcGeomSvc,
IMdcCalibFunSvc *  mdcFunSvc,
IMdcUtilitySvc *  mdcUtilitySvc 
)
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

◆ setParam() [1/2]

void MdcCalib::setParam ( MdcCalParams &  param)
inlinepure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

Definition at line 293 of file InstallArea/include/MdcCalibAlg/MdcCalibAlg/MdcCalib.h.

293 {
294 m_param = param;
295}

Referenced by GrXtMdcCalib::setParam(), IniMdcCalib::setParam(), PreT0MdcCalib::setParam(), PreXtMdcCalib::setParam(), QtMdcCalib::setParam(), T0MdcCalib::setParam(), Wr2dMdcCalib::setParam(), WrMdcCalib::setParam(), XtInteMdcCalib::setParam(), and XtMdcCalib::setParam().

◆ setParam() [2/2]

virtual void MdcCalib::setParam ( MdcCalParams &  param)
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

◆ updateConst() [1/2]

int MdcCalib::updateConst ( MdcCalibConst *  calconst)
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

Definition at line 1322 of file MdcCalib.cxx.

1322 {
1323 IMessageSvc* msgSvc;
1324 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
1325 MsgStream log(msgSvc, "MdcCalib");
1326 log << MSG::DEBUG << "MdcCalib::updateConst()" << endreq;
1327
1328 cout << "Tot " << m_hnTrk->GetEntries()
1329 << ", nTrkCut " << m_cut1 << ", cos(theta)_cut " << m_cut2 << ", drCut " << m_cut3
1330 << ", dzCut " << m_cut4 << ", nHitLayer_cut " << m_cut5 << ", nHit_cut " << m_cut6 << endl;
1331
1332 int lay;
1333 double effi;
1334 double effErr;
1335
1336 int nGoodAll = 0;
1337 int nGoodInn = 0;
1338 int nGoodStp = 0;
1339 int nGoodOut = 0;
1340 int nTotAll = 0;
1341 int nTotInn = 0;
1342 int nTotStp = 0;
1343 int nTotOut = 0;
1344 ofstream feffi("MdcLayerEffi.dat");
1345 for(lay=0; lay<m_nlayer; lay++){
1346 double effNtrk = m_effNtrk->GetBinContent(lay+1);
1347 double effGoodHit = m_effNtrkRecHit->GetBinContent(lay+1);
1348 nGoodAll += effGoodHit;
1349 if(lay < 8) nGoodInn += effGoodHit;
1350 else if (lay < 20) nGoodStp += effGoodHit;
1351 else nGoodOut += effGoodHit;
1352
1353 nTotAll += effNtrk;
1354 if(lay < 8) nTotInn += effNtrk;
1355 else if (lay < 20) nTotStp += effNtrk;
1356 else nTotOut += effNtrk;
1357
1358 effi = (double)effGoodHit / (double)effNtrk;
1359 effErr = sqrt(effi * (1-effi) / (double)effNtrk);
1360 feffi << setw(5) << lay << setw(15) << effi << setw(15) << effErr
1361 << setw(15) << effGoodHit << setw(15) << effNtrk << endl;
1362 }
1363 double effiAll = (double)nGoodAll / (double)(nTotAll);
1364 double errAll = sqrt(effiAll * (1-effiAll) / (double)(nTotAll));
1365 double effiInn = (double)nGoodInn / (double)(nTotInn);
1366 double errInn = sqrt(effiInn * (1-effiInn) / (double)(nTotInn));
1367 double effiStp = (double)nGoodStp / (double)(nTotStp);
1368 double errStp = sqrt(effiStp * (1-effiStp) / (double)(nTotStp));
1369 double effiOut = (double)nGoodOut / (double)(nTotOut);
1370 double errOut = sqrt(effiOut * (1-effiOut) / (double)(nTotOut));
1371 feffi << endl << "EffiAll: " << setw(15) << effiAll << setw(15) << errAll
1372 << setw(15) << nGoodAll << setw(15) << nTotAll << endl;
1373 feffi << endl << "EffiInn: " << setw(15) << effiInn << setw(15) << errInn
1374 << setw(15) << nGoodInn << setw(15) << nTotInn << endl;
1375 feffi << endl << "EffiStp: " << setw(15) << effiStp << setw(15) << errStp
1376 << setw(15) << nGoodStp << setw(15) << nTotStp << endl;
1377 feffi << endl << "EffiOut: " << setw(15) << effiOut << setw(15) << errOut
1378 << setw(15) << nGoodOut << setw(15) << nTotOut << endl;
1379 feffi.close();
1380
1381 // calculate efficiency without the impact of track fitting
1382 if(0 != m_param.fgCalDetEffi){
1383 int nHitAll[] = {0, 0};
1384 int nHitInn[] = {0, 0};
1385 int nHitStp[] = {0, 0};
1386 int nHitOut[] = {0, 0};
1387 ofstream feff2("MdcHitEffi.dat");
1388 for(lay=0; lay<m_nlayer; lay++){
1389 nHitAll[0] += m_hitNum[lay][0];
1390 nHitAll[1] += m_hitNum[lay][1];
1391 if(lay < 8){
1392 nHitInn[0] += m_hitNum[lay][0];
1393 nHitInn[1] += m_hitNum[lay][1];
1394 } else if (lay < 20){
1395 nHitStp[0] += m_hitNum[lay][0];
1396 nHitStp[1] += m_hitNum[lay][1];
1397 } else{
1398 nHitOut[0] += m_hitNum[lay][0];
1399 nHitOut[1] += m_hitNum[lay][1];
1400 }
1401
1402 effi = (double)m_hitNum[lay][1] / (double)m_hitNum[lay][0];
1403 effErr = sqrt(effi * (1-effi) / (double)m_hitNum[lay][0]);
1404 feff2 << setw(5) << lay << setw(15) << effi << setw(15) << effErr
1405 << setw(15) << m_hitNum[lay][1] << setw(15) << m_hitNum[lay][0] << endl;
1406 }
1407 effiAll = (double)nHitAll[1] / (double)nHitAll[0];
1408 errAll = sqrt(effiAll * (1-effiAll)) / (double)nHitAll[0];
1409 effiInn = (double)nHitInn[1] / (double)nHitInn[0];
1410 errInn = sqrt(effiInn * (1-effiInn)) / (double)nHitInn[0];
1411 effiStp = (double)nHitStp[1] / (double)nHitStp[0];
1412 errStp = sqrt(effiStp * (1-effiStp)) / (double)nHitStp[0];
1413 effiOut = (double)nHitOut[1] / (double)nHitOut[0];
1414 errOut = sqrt(effiOut * (1-effiOut)) / (double)nHitOut[0];
1415 feff2 << endl << "EffiAll: " << setw(15) << effiAll << setw(15) << errAll
1416 << setw(15) << nHitAll[1] << setw(15) << nHitAll[0] << endl;
1417 feff2 << endl << "EffiInn: " << setw(15) << effiInn << setw(15) << errInn
1418 << setw(15) << nHitInn[1] << setw(15) << nHitInn[0] << endl;
1419 feff2 << endl << "EffiStp: " << setw(15) << effiStp << setw(15) << errStp
1420 << setw(15) << nHitStp[1] << setw(15) << nHitStp[0] << endl;
1421 feff2 << endl << "EffiOut: " << setw(15) << effiOut << setw(15) << errOut
1422 << setw(15) << nHitOut[1] << setw(15) << nHitOut[0] << endl;
1423 feff2.close();
1424 }
1425
1426 // get resolution
1427 int i;
1428 int iEntr;
1429 int lr;
1430 int bin;
1431 int key;
1432 int hid;
1433
1434 Stat_t entry;
1435 double sigm[MdcCalSdNBIN];
1436 if(m_param.calSigma){
1437 ofstream fr2d("logr2d.dat");
1438 for(lay=0; lay<m_nlayer; lay++){
1439 for(iEntr=0; iEntr<m_nEntr[lay]; iEntr++){
1440 for(lr=0; lr<2; lr++){
1441 fr2d << setw(3) << lay << setw(3) << iEntr << setw(3) << lr << endl;
1442 for(bin=0; bin<m_nBin[lay]; bin++){
1443 key = getHresId(lay, iEntr, lr, bin);
1444 hid = m_mapr2d[key];
1445
1446 if(1 == m_param.resiType){
1447 entry = m_hr2dExc[hid] -> GetEntries();
1448 if(entry > 500){
1449 m_hr2dExc[hid] -> Fit("gaus", "Q");
1450 sigm[bin] = m_hr2dExc[hid]->GetFunction("gaus")->GetParameter(2);
1451 } else if(entry > 100){
1452 sigm[bin] = m_hr2dExc[hid]->GetRMS();
1453 } else{
1454 sigm[bin] = 0.2;
1455 }
1456 } else{
1457 entry = m_hr2dInc[hid] -> GetEntries();
1458 if(entry > 500){
1459 m_hr2dInc[hid] -> Fit("gaus", "Q");
1460 sigm[bin] = m_hr2dInc[hid]->GetFunction("gaus")->GetParameter(2);
1461 } else if(entry > 100){
1462 sigm[bin] = m_hr2dInc[hid]->GetRMS();
1463 } else{
1464 sigm[bin] = 0.2;
1465 }
1466 }
1467 if(sigm[bin] < 0.05) sigm[bin] = 0.05; // for boundary layers
1468 } // end of bin loop
1469
1470 for(bin=m_nBin[lay]; bin<MdcCalSdNBIN; bin++){
1471 sigm[bin] = sigm[m_nBin[lay]-1];
1472 }
1473
1474 for(bin=0; bin<MdcCalSdNBIN; bin++){
1475 if(1 == m_param.fgCalib[lay]){
1476// calconst -> resetSdpar(lay, iEntr, lr, bin, sigm[bin]);
1477 if(1 == m_nEntr[lay]){
1478 for(i=0; i<6; i++) calconst -> resetSdpar(lay, i, lr, bin, sigm[bin]);
1479 } else if(2 == m_nEntr[lay]){
1480 if(0 == iEntr){
1481 for(i=0; i<3; i++){ // entr<0
1482 calconst -> resetSdpar(lay, i, lr, bin, sigm[bin]);
1483 }
1484 } else{
1485 for(i=3; i<6; i++){ // entr>0
1486 calconst -> resetSdpar(lay, i, lr, bin, sigm[bin]);
1487 }
1488 }
1489 }
1490 } else{
1491 sigm[bin] = calconst->getSdpar(lay, iEntr, lr, bin);
1492 }
1493 fr2d << setw(5) << bin << setw(15) << sigm[bin] << endl;
1494 } // end of bin loop
1495 }
1496 } // end of entr loop
1497 }
1498 fr2d.close();
1499 }
1500
1501 return 1;
1502}
Fit
void Fit()
Definition: Eangle1D/Fit.cxx:3
MdcCalibConst::getSdpar
double getSdpar(int lay, int entr, int lr, int bin)
Definition: MdcCalibConst.cpp:105

Referenced by GrXtMdcCalib::updateConst(), QtMdcCalib::updateConst(), T0MdcCalib::updateConst(), WrMdcCalib::updateConst(), XtInteMdcCalib::updateConst(), and XtMdcCalib::updateConst().

◆ updateConst() [2/2]

virtual int MdcCalib::updateConst ( MdcCalibConst *  calconst)
pure virtual

Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, XtMdcCalib, GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.

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