Loading...
Searching...
No Matches
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 ¶m)=0 |
| virtual int | fillHist (MdcCalEvent *event)=0 |
| virtual int | updateConst (MdcCalibConst *calconst)=0 |
| virtual void | printCut () const =0 |
| virtual void | clear ()=0 |
Detailed Description
Definition at line 36 of file MdcCalib.h.
Constructor & Destructor Documentation
◆ MdcCalib()
| MdcCalib::MdcCalib | ( | ) |
Definition at line 36 of file MdcCalib.cxx.
36 {
37 m_nEvt=0;
38 m_cut1=0;
39 m_nTrkAfTrkCut=0;
40 m_cut2=0;
41 m_cut3=0;
42 m_cut4=0;
43 m_cut5=0;
44 m_cut6=0;
45 m_cut7=0;
46 m_nTrkCal=0;
47 m_nGrPoint = 0;
48 fgReadWireEff = false;
49
50 int lay;
51 for(lay=0; lay<43; lay++){
52 if(lay < 15) m_nEntr[lay] = 1;
53 else m_nEntr[lay] = 2;
54 // m_nEntr[lay] = 3;
55 }
56 m_dwid = 0.5; // mm
57 m_fgIni = false;
58
59 m_phiWid = PI2 / (double)NPhiBin;
60 m_theWid = 2.0 / (double)NThetaBin;
61
62 m_nEvtNtuple = 0;
63
65 if(lay < 8) m_nBin[lay] = 12;
66 else m_nBin[lay] = 16;
67 }
68
69 // setting boundary layer flags
71 if((0==lay) || (7==lay) || (8==lay) || (19==lay) || (20==lay) ||
72 (35==lay) || (36==lay) || (42==lay) ) m_layBound[lay] = true;
73 else m_layBound[lay] = false;
74 }
75}
◆ ~MdcCalib()
|
virtual |
Definition at line 77 of file MdcCalib.cxx.
77 {
78}
Member Function Documentation
◆ clear()
|
pure virtual |
Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.
Definition at line 80 of file MdcCalib.cxx.
80 {
81 int lay;
82 for(lay=0; lay<m_nlayer; lay++){
83 delete m_htraw[lay];
84 delete m_htdr[lay];
85 delete m_hresInc[lay];
86 delete m_hresExc[lay];
87 delete m_hresAve[lay];
88 delete m_hadc[lay];
89 for (int lr=0; lr<2; lr++){
90 delete m_htdrlr[lay][lr];
91 delete m_hreslrInc[lay][lr];
92 delete m_hreslrExc[lay][lr];
93 delete m_hreslrAve[lay][lr];
94 }
95 }
96
97 delete m_effNtrk;
98 delete m_effNtrkRecHit;
99 delete m_hresAllInc;
100 delete m_hresAllExc;
101 delete m_hresAllAve;
102 for(int i=0; i<14; i++){
103 delete m_hresAveAllQ[i];
104 delete m_hresAveOutQ[i];
105 }
106 for(lay=0; lay<43; lay++){
107 for(int i=0; i<14; i++) delete m_hresAveLayQ[lay][i];
108 }
109 delete m_hresInnInc;
110 delete m_hresInnExc;
111 delete m_hresStpInc;
112 delete m_hresStpExc;
113 delete m_hresOutInc;
114 delete m_hresOutExc;
115 delete m_hresAllExcUp;
116 delete m_hresAllExcDw;
117 delete m_hresInnExcUp;
118 delete m_hresInnExcDw;
119 delete m_hresStpExcUp;
120 delete m_hresStpExcDw;
121 delete m_hresOutExcUp;
122 delete m_hresOutExcDw;
124 delete m_hbbTrkFlg;
125 delete m_hTesAll;
126 delete m_hTesGood;
127 delete m_hTesAllFlag;
128 delete m_hTesRec;
129 delete m_hTesCalFlag;
130 delete m_hTesCalUse;
131 delete m_hnRawHit;
132 delete m_hpt;
133 delete m_hpMax;
134 delete m_hpMaxCms;
135 delete m_hptPos;
136 delete m_hptNeg;
137 delete m_hp;
138 delete m_hp_cms;
139 delete m_hpPos;
140 delete m_hpNeg;
141 delete m_hpPoscms;
142 delete m_hpNegcms;
143 delete m_hp_cut;
144 delete m_hchisq;
145 delete m_hnTrk;
146 delete m_hnTrkCal;
147 delete m_hnhitsRec;
148 delete m_hnhitsRecInn;
149 delete m_hnhitsRecStp;
150 delete m_hnhitsRecOut;
151 delete m_hnhitsCal;
152 delete m_hnhitsCalInn;
153 delete m_hnhitsCalStp;
154 delete m_hnhitsCalOut;
155 delete m_wirehitmap;
156 delete m_layerhitmap;
157 delete m_hnoisephi;
158 delete m_hnoiselay;
159 delete m_hnoisenhits;
160 delete m_hratio;
161 delete m_hdr;
162 delete m_hphi0;
163 delete m_hkap;
164 delete m_hdz;
165 delete m_htanl;
166 delete m_hcosthe;
167 delete m_hcostheNeg;
168 delete m_hcosthePos;
169 delete m_hx0;
170 delete m_hy0;
171 delete m_hdelZ0;
172 delete m_grX0Y0;
173 delete m_hitEffAll;
174 delete m_hitEffRaw;
175 delete m_hitEffRec;
177 int thbin;
183 for(thbin=0; thbin<NThetaBin; thbin++){
188 }
189 }
190 for(thbin=0; thbin<NThetaBin; thbin++){
191 delete m_ppThe[thbin];
192 delete m_pnThe[thbin];
193 delete m_ppTheCms[thbin];
194 delete m_pnTheCms[thbin];
195 }
196
197 for(unsigned i=0; i<m_hr2dInc.size(); i++){
198 delete m_hr2dInc[i];
199 delete m_hr2dExc[i];
200 }
201 m_hr2dInc.clear();
202 m_hr2dExc.clear();
203 m_mapr2d.clear();
204
205 delete m_fdTime;
206 delete m_fdAdc;
207 delete m_fdres;
208 delete m_fdresAve;
209 delete m_fdres2d;
210 delete m_fdcom;
211 delete m_fdResQ;
212}
*******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(), and XtMdcCalib::clear().
◆ fillHist()
|
pure virtual |
Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.
Definition at line 730 of file MdcCalib.cxx.
730 {
731 IMessageSvc* msgSvc;
734 log << MSG::DEBUG << "MdcCalib::fillHist()" << endreq;
735
736 int i;
737 int k;
738 int lay;
739 int cel;
740 int wir;
741 int lr;
742 int stat;
743
744 int hid;
746 int iEntr;
748
749 int phiBin;
750 int phiBinCms;
751 int theBin;
752 int theBinCms;
753 double lamda;
754 double theta;
755
756 double qhit;
757 double traw;
758 double tdr;
759 double doca;
760 double resiInc;
761 double resiExc;
762 double entr;
763 double pt;
764 double p;
766 double chisq;
771 HepLorentzVector p4;
772
773 double dr;
774 double phi0;
775 double dz;
776 double kap;
777 double tanl;
778
779 double x0;
780 double y0;
781 double zminus = 9999.0;
782 double zplus = -9999.0;
783
784 double hitphi;
785 double philab;
786 double phicms;
787 double thetacms;
788 double costheCMS;
789
790 double dphi;
791 double wphi;
792 double xx;
793 double yy;
794 double rr;
795
796 int nhitlay;
798 bool fgTrk;
799
800 int ntrk = event -> getNTrk();
801 int nhitRec;
802 int nhitRecInn;
803 int nhitRecStp;
804 int nhitRecOut;
805 int nhitCal;
806 int nhitCalInn;
807 int nhitCalStp;
808 int nhitCalOut;
809 MdcCalRecTrk* rectrk;
810 MdcCalRecHit* rechit;
811
812 int fgAdd[43]; // for calculating layer efficiency
813
814 // read dead wire
815 if(!fgReadWireEff){
818 for(cel=0; cel<ncel; cel++){
820 if(eff > 0.5) m_fgGoodWire[lay][cel] = true;
821 else m_fgGoodWire[lay][cel] = false;
822 int wireid = m_mdcGeomSvc -> Wire(lay, cel) -> Id();
823 if(eff<0.9) cout << "dead channel: " << setw(5) << lay << setw(5) << cel << setw(8) << wireid << endl;
824 }
825 }
826 fgReadWireEff = true;
827
836 fwpc << setw(6) << wir << setw(15) << m_xe[wir] << setw(15) << m_ye[wir]
837 << setw(15) << m_ze[wir] << setw(15) << m_xw[wir]
838 << setw(15) << m_yw[wir] << setw(15) << m_zw[wir] << endl;
839 }
840 fwpc.close();
841 }
842
843 int nRawHit = event->getNRawHitTQ();
844 m_hnRawHit->Fill(nRawHit);
845
846 IDataProviderSvc* eventSvc = NULL;
847 Gaudi::svcLocator()->service("EventDataSvc", eventSvc);
848
849 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc,"/Event/EventHeader");
850 if (!eventHeader) {
851 log << MSG::FATAL << "Could not find Event Header" << endreq;
852 return( StatusCode::FAILURE);
853 }
854 int iRun = eventHeader->runNumber();
855 int iEvt = eventHeader->eventNumber();
856
857 int esTimeflag = event->getEsFlag();
858 double tes = event->getTes();
859 bool esCutFg = event->getEsCutFlag();
860 int iEs = event->getNesCutFlag();
861 //calculate the efficiency of Bhabha event
862 if (-1 != esTimeflag) {
863 SmartDataPtr<RecMdcTrackCol> newtrkCol(eventSvc, "/Event/Recon/RecMdcTrackCol");
864 if(!newtrkCol){
865 log << MSG::ERROR << "Could not find RecMdcTrackCol" << endreq;
866 return ( StatusCode::FAILURE );
867 }
868 int nGoodTrk = 0;
869 int icharge = 0;
870 Vp4 p4p;
871 Vp4 p4m;
872 RecMdcTrackCol::iterator it_trk = newtrkCol->begin();
873 for(; it_trk != newtrkCol->end(); it_trk++){
876 icharge = (*it_trk)->charge();
879 }
880 if (1 == p4p.size() * p4m.size()){
881 double dang = p4p[0].vect().angle(p4m[0].vect());
882 m_hbbTrkFlg->Fill(1);
883 if (dang > 2.94) {
884 m_hbbTrkFlg->Fill(2);
885 }
886 }
887
888 }
889 m_hTesAll->Fill(tes);
890// cout << "tes " << tes << endl;
891 // if (-1 != esTimeflag) m_hTesGood->Fill(tes);
892 if (-1 == esTimeflag) return -1;
893 m_hTesGood->Fill(tes);
894 m_hTesAllFlag->Fill(esTimeflag);
895 if(ntrk > 0) m_hTesRec->Fill(tes);
897
898 if( esCutFg ) m_hTesCalFlag->Fill(tes);
899 else return -1;
900
901 if(! m_fgIni){
905 }
906 m_fgIni = true;
907 }
908
909 bool trkFlag[200]; // for calculating hit efficiency without impact of track fitting
910 for(i=0; i<200; i++) trkFlag[i] = false;
911
912 int ntrkCal = 0;
913 double pTrk[2];
914 double pTrkcms[2];
915 double hitphiplus = 9999.0;
916 double hitthetaplus = 9999.0;
917 double hitphiminus = -9999.0;
918 double hitthetaminus = -9999.0;
919 Vp4 pp;
920 Vp4 pm;
921 pp.clear();
922 pm.clear();
923 int phibinp;
924 int phibinm;
925
926 m_hnTrk->Fill(ntrk);
928 m_cut1++;
929 return -1;
930 }
931// if(ntrk > 2) return -1;
932 for(i=0; i<ntrk; i++){
933 m_nTrkAfTrkCut++;
934 fgTrk = true;
935 rectrk = event -> getRecTrk(i);
936 nhitRec = rectrk -> getNHits();
937 m_hnhitsRec -> Fill( nhitRec );
938
940 fgHitLay[lay] = false;
941 }
942
943 nhitRecInn = 0;
944 nhitRecStp = 0;
945 nhitRecOut = 0;
946 for(k=0; k<nhitRec; k++){
947 rechit = rectrk -> getRecHit(k);
948 lay = rechit -> getLayid();
949 doca = rechit -> getDocaExc();
950 resiExc = rechit -> getResiExc();
951 fgHitLay[lay] = true;
952
953 if(lay < 8) nhitRecInn++;
954 else if(lay < 20) nhitRecStp++;
955 else nhitRecOut++;
956 }
957 m_hnhitsRecInn->Fill(nhitRecInn);
958 m_hnhitsRecStp->Fill(nhitRecStp);
959 m_hnhitsRecOut->Fill(nhitRecOut);
960
961 // get momentum
962 pt = rectrk -> getPt();
963 p = rectrk -> getP();
964
965 // boost P to the cms
967 HepLorentzVector psip(xboost, yboost, zboost, ecm);
968 Hep3Vector boostv = psip.boostVector();
969 p4.boost(- boostv);
971 phicms = p4.phi();
973 thetacms = p4.theta();
974 costheCMS = cos(thetacms);
976 p4.boost(boostv);
977// cout << setw(15) << p << setw(15) << p_cms << setw(15) << costheCMS << endl;
978
979 // cos(theta) cut
981 m_cut2++;
982 continue;
983 }
984
985 // dr cut
986 dr = rectrk->getDr();
988 m_cut3++;
989 continue;
990 }
991
992 // dz cut
993 dz = rectrk->getDz();
995 m_cut4++;
996 continue;
997 }
998
999 // momentum cut
1001 m_cut5++;
1002 continue;
1003 }
1004
1005 // charge cut, to check +/- asymmetry
1009 continue;
1010 }
1011
1012// if(! fgTrk) continue;
1013
1014 // hit layer cut
1015 nhitlay = 0;
1017 if(fgHitLay[lay]) nhitlay++;
1018 }
1020 m_cut6++;
1021 continue;
1022 }
1023
1024 // nhit cut
1026 m_cut7++;
1027 continue;
1028 }
1029
1030// bool fgNoise = rectrk->getFgNoiseRatio();
1031// if(m_param.noiseCut && (!fgNoise)) continue;
1032// cout << setw(10) << iRun << setw(15) << iEvt << setw(5) << fgNoise << endl;
1033
1034// if(! ((fgHitLay[0]||fgHitLay[1]) && (fgHitLay[41]||fgHitLay[42])) ){
1035// continue;
1036// }
1037
1038 // calculate cell on the track
1039 int cellTrkPass[43];
1040 bool fgPass = getCellTrkPass(event, i, cellTrkPass);
1041 for(lay=0; lay<m_nlayer; lay++){
1042 fgAdd[lay] = 0;
1043// if((16==lay) || (18==lay) || (19==lay) || (41==lay)){ // hv2200 2009-3
1044 if((15==lay) || (16==lay) || (18==lay) || (19==lay) || (40==lay) || (41==lay) || (42==lay)){
1045 int iCell = cellTrkPass[lay];
1046 if(fgPass && (iCell >= 0) && m_fgGoodWire[lay][iCell]) m_effNtrk->Fill(lay);
1047 else fgAdd[lay] = 1;
1048 } else{
1049 m_effNtrk->Fill(lay);
1050 }
1051 }
1052
1053 m_nTrkCal++;
1054 chisq = rectrk -> getChisq();
1055 m_hchisq -> Fill( chisq );
1056
1057 if(pt > 0){
1058 m_hpt -> Fill(pt);
1059 m_hptPos -> Fill(pt);
1060 m_hp -> Fill(p);
1062 m_hpPos -> Fill(p);
1064 } else{
1065 m_hpt -> Fill(-1.0*pt);
1066 m_hptNeg -> Fill(-1.0*pt);
1067 m_hp -> Fill(-1.0*p);
1069 m_hpNeg -> Fill(-1.0*p);
1071 }
1072 if(2 == ntrk){
1073 pTrk[i] = fabs(p);
1074 pTrkcms[i] = fabs(p_cms);
1075 }
1076
1077 dr = rectrk -> getDr();
1078 phi0 = rectrk -> getPhi0();
1079 kap = rectrk -> getKappa();
1080 dz = rectrk -> getDz();
1081 tanl = rectrk -> getTanLamda();
1082 lamda = atan(tanl);
1083 theta = HFPI - lamda;
1084
1085 m_hdr -> Fill(dr);
1086 m_hphi0 -> Fill(phi0);
1087 m_hkap -> Fill(kap);
1088 m_hdz -> Fill(dz);
1089 m_htanl -> Fill(tanl);
1093
1094 philab = phi0 + HFPI;
1096// cout << setw(15) << phi0 << setw(15) << philab << setw(15) << phicms << endl;
1097
1098 phiBin = (int)(philab / m_phiWid);
1099 phiBinCms = (int)(phicms / m_phiWid);
1100 theBin = (int)((cos(theta) + 1.0) / m_theWid);
1101 theBinCms = (int)((cos(thetacms) + 1.0) / m_theWid);
1102 if(phiBin < 0) phiBin = 0;
1103 if(phiBin >= NPhiBin) phiBin = NPhiBin-1;
1104 if(phiBinCms < 0) phiBinCms = 0;
1105 if(phiBinCms >= NPhiBin) phiBinCms = NPhiBin-1;
1106 if(theBin < 0) theBin = 0;
1107 if(theBin >= NThetaBin) theBin = NThetaBin-1;
1108 if(theBinCms < 0) theBinCms = 0;
1109 if(theBinCms >= NThetaBin) theBinCms = NThetaBin-1;
1110
1111 if(pt > 0){
1112 m_ppPhi[phiBin]->Fill(p);
1113 m_ppPhiCms[phiBinCms]->Fill(p_cms);
1114 m_ppThe[theBin]->Fill(p);
1115 m_ppTheCms[theBinCms]->Fill(p_cms);
1116 m_ppThePhi[theBin][phiBin]->Fill(p);
1117 m_ppThePhiCms[theBinCms][phiBinCms]->Fill(p_cms);
1118 } else{
1119 m_pnPhi[phiBin]->Fill(-1.0*p);
1120 m_pnPhiCms[phiBinCms]->Fill(-1.0*p_cms);
1121 m_pnThe[theBin]->Fill(-1.0*p);
1122 m_pnTheCms[theBinCms]->Fill(-1.0*p_cms);
1123 m_pnThePhi[theBin][phiBin]->Fill(-1.0*p);
1124 m_pnThePhiCms[theBinCms][phiBinCms]->Fill(-1.0*p_cms);
1125 }
1126
1127 x0 = dr * cos(phi0);
1128 y0 = dr * sin(phi0);
1129 m_hx0 -> Fill(x0);
1130 m_hy0 -> Fill(y0);
1131 if(m_nGrPoint < 10000){
1132 m_grX0Y0->SetPoint(m_nGrPoint, x0, y0);
1133 m_nGrPoint++;
1134 }
1135
1136 if(kap < 0) {
1137 zminus = dz;
1138 pm.push_back(p4);
1139 phibinm = phiBinCms;
1140 } else {
1141 zplus = dz;
1142 pp.push_back(p4);
1143 phibinp = phiBinCms;
1144 }
1145
1146// cout << "phi = " << setw(15) << philab << setw(15) << philab*180./3.14159 << setw(15) << p << endl;
1147 ntrkCal++;
1148 trkFlag[i] = true;
1149 nhitCal = 0;
1150 nhitCalInn = 0;
1151 nhitCalStp = 0;
1152 nhitCalOut = 0;
1153 for(k=0; k<nhitRec; k++){
1154 rechit = rectrk -> getRecHit(k);
1155
1156 lay = rechit -> getLayid();
1157 cel = rechit -> getCellid();
1158 lr = rechit -> getLR();
1159 stat = rechit -> getStat();
1160 doca = rechit -> getDocaExc();
1161 resiInc = rechit -> getResiIncLR();
1162 resiExc = rechit -> getResiExcLR();
1163 entr = rechit -> getEntra();
1164 tdr = rechit -> getTdrift();
1165 traw = (rechit -> getTdc()) * MdcCalTdcCnv;
1166 wir = m_mdcGeomSvc -> Wire(lay, cel) -> Id();
1167
1168 m_cel[lay] = (long)cel;
1169 m_lr[lay] = (long)lr;
1170 m_run[lay] = iRun;
1171 m_evt[lay] = iEvt;
1172 m_doca[lay] = doca;
1173 m_dm[lay] = rechit -> getDmeas();
1174 m_tdr[lay] = tdr;
1175 m_tdc[lay] = traw;
1176 m_entr[lay] = entr*180.0/3.14;
1177 m_zhit[lay] = rechit -> getZhit();
1178 m_qhit[lay] = rechit -> getQhit();
1179 m_p[lay] = p_cms;
1180 m_pt[lay] = pt;
1181 m_phi0[lay] = phi0;
1182 m_tanl[lay] = tanl;
1183 qhit = rechit -> getQhit();
1184// cout << "layer " << setw(5) << lay << setw(5) << cel << setw(5) << lr << endl;
1185 // calculating hitphi
1186 xx = (m_zhit[lay] - m_zw[wir]) * (m_xe[wir] - m_xw[wir]) /
1187 (m_ze[wir] - m_zw[wir]) + m_xw[wir];
1188 yy = (m_zhit[lay] - m_zw[wir]) * (m_ye[wir] - m_yw[wir]) /
1189 (m_ze[wir] - m_zw[wir]) + m_yw[wir];
1190 rr = sqrt( (xx * xx) + (yy * yy) );
1191 dphi = fabs(doca) / m_radii[lay];
1192
1193 if( yy >= 0 ) wphi = acos(xx / rr);
1195 if(1 == lr) hitphi = wphi + dphi; // mention
1196 else hitphi = wphi - dphi;
1199
1200 m_hitphi[lay] = hitphi;
1201
1202 if( (fabs(doca) > m_docaMax[lay]) ||
1203 (fabs(resiExc) > m_param.resiCut[lay]) ){
1204 continue;
1205 }
1206
1208 if(0 == lay){
1209 if( ! fgHitLay[1] ) continue;
1210 } else if(42 == lay){
1211 if( ! fgHitLay[41] ) continue;
1212 } else{
1213 if( (!fgHitLay[lay-1]) && (!fgHitLay[lay+1]) ) continue;
1214
1215 // for boundary layers
1217 ((!fgHitLay[lay-1]) || (!fgHitLay[lay+1])) ) continue;
1218 }
1219 }
1220
1222
1223 // fill xtplot tree
1225 m_xtTuple[lay] -> write();
1226 }
1227
1229 if( (0 == fgAdd[lay]) && (1 == stat) ){
1230 m_effNtrkRecHit->Fill(lay);
1231 fgAdd[lay] = 1;
1232 }
1233 } else{
1234 if(0 == fgAdd[lay]){
1235 m_effNtrkRecHit->Fill(lay);
1236 fgAdd[lay] = 1;
1237 }
1238 }
1239
1240 nhitCal++;
1241 if(lay < 8) nhitCalInn++;
1242 else if(lay < 20) nhitCalStp++;
1243 else nhitCalOut++;
1244
1245 m_wirehitmap -> Fill(wir);
1246 m_layerhitmap -> Fill(lay);
1247
1248 m_htraw[lay] -> Fill(traw);
1249 m_htdr[lay] -> Fill(tdr);
1250 m_htdrlr[lay][lr]->Fill(tdr);
1251 m_hadc[lay] -> Fill(qhit);
1252
1253 m_hresAllInc -> Fill(resiInc);
1254 m_hresAllExc -> Fill(resiExc);
1255 double resiAve = 0.5 * (resiInc + resiExc);
1256 m_hresAllAve -> Fill(resiAve);
1257 if(yy > 0) m_hresAllExcUp->Fill(resiExc);
1258 else m_hresAllExcDw->Fill(resiExc);
1259
1260 if(lay < 8){
1261 m_hresInnInc -> Fill(resiInc);
1262 m_hresInnExc -> Fill(resiExc);
1263 if(yy > 0) m_hresInnExcUp->Fill(resiExc);
1264 else m_hresInnExcDw->Fill(resiExc);
1265 } else if(lay < 20){
1266 m_hresStpInc -> Fill(resiInc);
1267 m_hresStpExc -> Fill(resiExc);
1268 if(yy > 0) m_hresStpExcUp->Fill(resiExc);
1269 else m_hresStpExcDw->Fill(resiExc);
1270 } else{
1271 m_hresOutInc -> Fill(resiInc);
1272 m_hresOutExc -> Fill(resiExc);
1273 if(yy > 0) m_hresOutExcUp->Fill(resiExc);
1274 else m_hresOutExcDw->Fill(resiExc);
1275 }
1276
1277 int qbin = (int)((qhit-100.0)/100.0);
1278 if(qbin>=0 && qbin<14){
1279 m_hresAveAllQ[qbin]->Fill(resiAve);
1280 m_hresAveLayQ[lay][qbin]->Fill(resiAve);
1281 if(lay > 7) m_hresAveOutQ[qbin]->Fill(resiAve);
1282 }
1283
1285 m_hresInc[lay] -> Fill(resiInc);
1286 m_hreslrInc[lay][lr]->Fill(resiInc);
1287 m_hresExc[lay] -> Fill(resiExc);
1288 m_hreslrExc[lay][lr]->Fill(resiExc);
1289 m_hresAve[lay] -> Fill(resiAve);
1290 m_hreslrAve[lay][lr]->Fill(resiAve);
1291 }
1292
1294 if(iPhi>=20) iPhi = 19;
1295 m_hresphi[lay][iPhi]->Fill((resiInc+resiExc)*0.5);
1296
1297// bin = (int)(fabs(doca) / m_dwid);
1299 iEntr = m_mdcFunSvc -> getSdEntrIndex(entr);
1300 if(1 == m_nEntr[lay]){
1301 iEntr = 0;
1302 } else if(2 == m_nEntr[lay]){
1303 if(entr > 0.0) iEntr = 1;
1304 else iEntr = 0;
1305 } else if(3 == m_nEntr[lay]){
1306 double entrBinWid = 0.174532925; // 10 degree
1307 iEntr = (int)(fabs(entr)/entrBinWid);
1308 if(iEntr > 2) iEntr = 2;
1309 }
1313 hid = m_mapr2d[key];
1314 m_hr2dInc[hid] -> Fill(resiInc);
1315 m_hr2dExc[hid] -> Fill(resiExc);
1316 }
1317 }
1318
1319 if((tdr>0) && (tdr<750)){
1322 m_hr2t[lay][iEntr][lr][bin]->Fill(resiExc);
1323 }
1324 } // loop of nhits
1325 m_nEvtNtuple++;
1326 m_hnhitsCal->Fill(nhitCal);
1327 m_hnhitsCalInn->Fill(nhitCalInn);
1328 m_hnhitsCalStp->Fill(nhitCalStp);
1329 m_hnhitsCalOut->Fill(nhitCalOut);
1330 } // end of track loop
1331 m_hnTrkCal->Fill(ntrkCal);
1332 if(2 == ntrkCal){
1333 if(pTrk[0] > pTrk[1]) m_hpMax->Fill(pTrk[0]);
1334 else m_hpMax->Fill(pTrk[1]);
1335
1336 if(pTrkcms[0] > pTrkcms[1]) m_hpMaxCms->Fill(pTrkcms[0]);
1337 else m_hpMaxCms->Fill(pTrkcms[1]);
1338 }
1339 if(ntrkCal > 0) m_hTesCalUse->Fill(tes);
1340
1341 double delZ0;
1342 if((fabs(zminus) < 9000.0) && (fabs(zplus) < 9000.0)) delZ0 = zplus - zminus;
1343 m_hdelZ0 -> Fill(delZ0);
1344
1345 if (1 == pp.size() * pm.size()){
1346 HepLorentzVector ptot = pp[0] + pm[0];
1347 bool fourmomcut = false;
1348// fourmomcut = (ptot.x()>0.02 && ptot.x()<0.06) && (fabs(ptot.y()) < 0.02)
1349// && (ptot.z()>-0.01 && ptot.z()<0.03) && (ptot.e()>3.4 && ptot.e()<4.0);
1350 fourmomcut = (fabs(ptot.x()-0.04)<0.026) && (fabs(ptot.y()) < 0.026)
1351 && (fabs(ptot.z()-0.005)<0.036) && (fabs(ptot.e()-ecm)<0.058);
1352 //cout << "x = " << ptot.x() << ", y = " << ptot.y() << ", z = " << ptot.z() << ", e = " << ptot.e() << endl;
1353 if (fourmomcut) {
1354 HepLorentzVector psip(xboost, yboost, zboost, ecm);
1355 Hep3Vector boostv = psip.boostVector();
1356 pp[0].boost(- boostv);
1357 pm[0].boost(- boostv);
1358 m_hp_cut->Fill(pp[0].rho());
1359 m_hp_cut->Fill(pm[0].rho());
1360 }
1361 }
1362
1363 if(2==ntrk) for(i=0; i<ntrk; i++) pTrk[i] = (event -> getRecTrk(i)) -> getP();
1365// if(1==ntrk) p = (event->getRecTrk(0)) -> getP();
1366// if((5==m_param.particle) && (1==ntrk) && (fabs(p)<5)){
1367 m_tescos = tes;
1368 m_tesFlagcos = esTimeflag;
1369 for(i=0; i<ntrk; i++){
1370 rectrk = event -> getRecTrk(i);
1371 phi0 = rectrk -> getPhi0();
1373
1374 tanl = rectrk -> getTanLamda();
1375 lamda = atan(tanl);
1376 theta = HFPI - lamda;
1377
1379 m_nhitUpcos = rectrk -> getNHits();
1380 m_pUpcos = rectrk -> getP();
1381 m_ptUpcos = rectrk -> getPt();
1382 m_phiUpcos = phi0;
1383 m_drUpcos = rectrk->getDr();
1384 m_dzUpcos = rectrk->getDz();
1385 m_ctheUpcos = cos(theta);
1386 } else{
1387 m_nhitDwcos = rectrk -> getNHits();
1388 m_pDwcos = rectrk -> getP();
1389 m_ptDwcos = rectrk -> getPt();
1390 m_phiDwcos = phi0;
1391 m_drDwcos = rectrk->getDr();
1392 m_dzDwcos = rectrk->getDz();
1393 m_ctheDwcos = cos(theta);
1394
1395 if(m_pDwcos > 0) m_chargecos = 1;
1396 else m_chargecos = 0;
1397 }
1398 }
1399 m_cosmic->write();
1400 }
1401
1403
1404 return 1;
1405}
curve Fill()
const HepLorentzVector p_cms(0.034067, 0.0, 0.0, 3.097)
*************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
virtual double getWireEff(int layid, int cellid) const =0
virtual const MdcGeoLayer *const Layer(unsigned id)=0
virtual const MdcGeoWire *const Wire(unsigned id)=0
Definition: MdcCalRecHit.h:18
Definition: MdcCalRecTrk.h:19
Referenced by GrXtMdcCalib::fillHist(), QtMdcCalib::fillHist(), T0MdcCalib::fillHist(), Wr2dMdcCalib::fillHist(), WrMdcCalib::fillHist(), XtInteMdcCalib::fillHist(), and XtMdcCalib::fillHist().
◆ initialize()
|
pure virtual |
Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.
Definition at line 214 of file MdcCalib.cxx.
215 {
216 IMessageSvc* msgSvc;
219 log << MSG::INFO << "MdcCalib::initialize()" << endreq;
220
221 m_hlist = hlist;
222 m_mdcGeomSvc = mdcGeomSvc;
223 m_mdcFunSvc = mdcFunSvc;
224 m_mdcUtilitySvc = mdcUtilitySvc;
225
226 int lay;
227 int iEntr;
228 int lr;
230 char hname[200];
231
232 int nbinMom = 2400;
233
234 m_nlayer = m_mdcGeomSvc -> getLayerSize();
235
236 for(lay=0; lay<m_nlayer; lay++){
238 }
247 fwpc << setw(6) << wir << setw(15) << m_xe[wir] << setw(15) << m_ye[wir]
248 << setw(15) << m_ze[wir] << setw(15) << m_xw[wir]
249 << setw(15) << m_yw[wir] << setw(15) << m_zw[wir] << endl;
250 }
251 fwpc.close();
252
253 m_fdcom = new TFolder("common", "common");
254 m_hlist -> Add(m_fdcom);
255
256 m_effNtrk = new TH1F("effNtrk", "", 43, -0.5, 42.5);
257 m_fdcom->Add(m_effNtrk);
258
259 m_effNtrkRecHit = new TH1F("effNtrkRecHit", "", 43, -0.5, 42.5);
260 m_fdcom->Add(m_effNtrkRecHit);
261
262 m_hresAllInc = new TH1F("HResAllInc", "", 200, -1.0, 1.0);
263 m_fdcom -> Add(m_hresAllInc);
264
265 m_hresAllExc = new TH1F("HResAllExc", "", 200, -1.0, 1.0);
266 m_fdcom -> Add(m_hresAllExc);
267
268 m_hresAllAve = new TH1F("HResAllAve", "", 200, -1.0, 1.0);
269 m_fdcom -> Add(m_hresAllAve);
270
271 m_hresInnInc = new TH1F("HResInnInc", "", 200, -1.0, 1.0);
272 m_fdcom -> Add(m_hresInnInc);
273
274 m_hresInnExc = new TH1F("HResInnExc", "", 200, -1.0, 1.0);
275 m_fdcom -> Add(m_hresInnExc);
276
277 m_hresStpInc = new TH1F("HResStpInc", "", 200, -1.0, 1.0);
278 m_fdcom -> Add(m_hresStpInc);
279
280 m_hresStpExc = new TH1F("HResStpExc", "", 200, -1.0, 1.0);
281 m_fdcom -> Add(m_hresStpExc);
282
283 m_hresOutInc = new TH1F("HResOutInc", "", 200, -1.0, 1.0);
284 m_fdcom -> Add(m_hresOutInc);
285
286 m_hresOutExc = new TH1F("HResOutExc", "", 200, -1.0, 1.0);
287 m_fdcom -> Add(m_hresOutExc);
288
289 m_hresAllExcUp = new TH1F("HResAllExcUp", "", 200, -1.0, 1.0);
290 m_fdcom -> Add(m_hresAllExcUp);
291
292 m_hresAllExcDw = new TH1F("HResAllExcDw", "", 200, -1.0, 1.0);
293 m_fdcom -> Add(m_hresAllExcDw);
294
295 m_hresInnExcUp = new TH1F("HResInnExcUp", "", 200, -1.0, 1.0);
296 m_fdcom -> Add(m_hresInnExcUp);
297
298 m_hresInnExcDw = new TH1F("HResInnExcDw", "", 200, -1.0, 1.0);
299 m_fdcom -> Add(m_hresInnExcDw);
300
301 m_hresStpExcUp = new TH1F("HResStpExcUp", "", 200, -1.0, 1.0);
302 m_fdcom -> Add(m_hresStpExcUp);
303
304 m_hresStpExcDw = new TH1F("HResStpExcDw", "", 200, -1.0, 1.0);
305 m_fdcom -> Add(m_hresStpExcDw);
306
307 m_hresOutExcUp = new TH1F("HResOutExcUp", "", 200, -1.0, 1.0);
308 m_fdcom -> Add(m_hresOutExcUp);
309
310 m_hresOutExcDw = new TH1F("HResOutExcDw", "", 200, -1.0, 1.0);
311 m_fdcom -> Add(m_hresOutExcDw);
312
313 m_fdResQ = new TFolder("ResQ", "ResQ");
314 m_hlist->Add(m_fdResQ);
315 for(int i=0; i<14; i++){
317 m_hresAveAllQ[i] = new TH1F(hname, "", 200, -1, 1);
318 m_fdResQ->Add(m_hresAveAllQ[i]);
319
321 m_hresAveOutQ[i] = new TH1F(hname, "", 200, -1, 1);
322 m_fdResQ->Add(m_hresAveOutQ[i]);
323 }
324 for(lay=0; lay<43; lay++){
325 for(int i=0; i<14; i++){
327 m_hresAveLayQ[lay][i] = new TH1F(hname, "", 200, -1, 1);
328 m_fdResQ->Add(m_hresAveLayQ[lay][i]);
329 }
330 }
331
334 m_hTes[iEs] = new TH1F(hname, "", 750, 0, 1500);
335 m_fdcom->Add(m_hTes[iEs]);
336 }
337
338 m_hbbTrkFlg = new TH1F("BbTrkFlg", "", 100, 0, 6);
339 m_fdcom -> Add(m_hbbTrkFlg);
340
341 m_hTesAll = new TH1F("TesAll", "", 2000, 0, 2000);
342 m_fdcom -> Add(m_hTesAll);
343
344 m_hTesGood = new TH1F("TesGood", "", 2000, 0, 2000);
345 m_fdcom -> Add(m_hTesGood);
346
347 m_hTesAllFlag = new TH1F("TesAllFlag", "", 300, -0.5, 299.5);
348 m_fdcom -> Add(m_hTesAllFlag);
349
350 m_hTesRec = new TH1F("TesRec", "", 2000, 0, 2000);
351 m_fdcom -> Add(m_hTesRec);
352
353 m_hTesCalFlag = new TH1F("TesCalFlag", "", 2000, 0, 2000);
354 m_fdcom -> Add(m_hTesCalFlag);
355
356 m_hTesCalUse = new TH1F("TesCalUse", "", 2000, 0, 2000);
357 m_fdcom -> Add(m_hTesCalUse);
358
359 m_hnRawHit = new TH1F("nRawHit", "", 6797, -0.5, 6796.5);
360 m_fdcom -> Add(m_hnRawHit);
361
362 m_hpt = new TH1F("HPt", "", nbinMom, 0, 3);
363 m_fdcom -> Add(m_hpt);
364
365 m_hptPos = new TH1F("HPtPos", "", nbinMom, 0, 3);
366 m_fdcom -> Add(m_hptPos);
367
368 m_hptNeg = new TH1F("HPtNeg", "", nbinMom, 0, 3);
369 m_fdcom -> Add(m_hptNeg);
370
371 m_hp = new TH1F("HP", "", nbinMom, 0, 3);
372 m_fdcom -> Add(m_hp);
373
374 m_hp_cms = new TH1F("HPCMS", "", nbinMom, 0, 3);
375 m_fdcom -> Add(m_hp_cms);
376
377 m_hpMax = new TH1F("HPMax", "", nbinMom, 0, 3);
378 m_fdcom -> Add(m_hpMax);
379
380 m_hpMaxCms = new TH1F("HPMax_Cms", "", nbinMom, 0, 3);
381 m_fdcom -> Add(m_hpMaxCms);
382
383 m_hpPos = new TH1F("HP_Pos", "", nbinMom, 0, 3);
384 m_fdcom -> Add(m_hpPos);
385
386 m_hpNeg = new TH1F("HP_Neg", "", nbinMom, 0, 3);
387 m_fdcom -> Add(m_hpNeg);
388
389 m_hpPoscms = new TH1F("HP_Pos_cms", "", nbinMom, 0, 3);
390 m_fdcom -> Add(m_hpPoscms);
391
392 m_hpNegcms = new TH1F("HP_Neg_cms", "", nbinMom, 0, 3);
393 m_fdcom -> Add(m_hpNegcms);
394
395 m_hp_cut = new TH1F("HPCut", "", nbinMom, 0, 3);
396 m_fdcom -> Add(m_hp_cut);
397
398 m_hchisq = new TH1F("Chisq", "", 10, 0, 100);
399 m_fdcom -> Add(m_hchisq);
400
401 m_hnTrk = new TH1F("HNtrack", "HNtrack", 10, -0.5, 9.5);
402 m_fdcom -> Add(m_hnTrk);
403
404 m_hnTrkCal = new TH1F("HNtrackCal", "HNtrackCal", 10, -0.5, 9.5);
405 m_fdcom -> Add(m_hnTrkCal);
406
407 m_hnhitsRec = new TH1F("HNhitsRec", "", 100, -0.5, 99.5);
408 m_fdcom -> Add(m_hnhitsRec);
409
410 m_hnhitsRecInn = new TH1F("HNhitsInnRec", "", 60, 0.5, 60.5);
411 m_fdcom -> Add(m_hnhitsRecInn);
412
413 m_hnhitsRecStp = new TH1F("HNhitsStpRec", "", 60, 0.5, 60.5);
414 m_fdcom -> Add(m_hnhitsRecStp);
415
416 m_hnhitsRecOut = new TH1F("HNhitsOutRec", "", 60, 0.5, 60.5);
417 m_fdcom -> Add(m_hnhitsRecOut);
418
419 m_hnhitsCal = new TH1F("HNhitsCal", "", 100, -0.5, 99.5);
420 m_fdcom -> Add(m_hnhitsCal);
421
422 m_hnhitsCalInn = new TH1F("HNhitsCalInn", "", 60, 0.5, 60.5);
423 m_fdcom -> Add(m_hnhitsCalInn);
424
425 m_hnhitsCalStp = new TH1F("HNhitsCalStp", "", 60, 0.5, 60.5);
426 m_fdcom -> Add(m_hnhitsCalStp);
427
428 m_hnhitsCalOut = new TH1F("HNhitsCalOut", "", 60, 0.5, 60.5);
429 m_fdcom -> Add(m_hnhitsCalOut);
430
431 m_wirehitmap = new TH1F("Wire_HitMap", "Wire_HitMap", 6796, -0.5, 6795.5);
432 m_fdcom -> Add(m_wirehitmap);
433
434 m_layerhitmap = new TH1F("Layer_HitMap", "Layer_HitMap", 43, -0.5, 42.5);
435 m_fdcom -> Add(m_layerhitmap);
436
437 m_hnoisephi = new TH1F("phi_noise", "", 100, 0, 6.284);
438 m_fdcom -> Add(m_hnoisephi);
439
440 m_hnoiselay = new TH1F("Layer_noise", "Layer_noise", 43, -0.5, 42.5);
441 m_fdcom -> Add(m_hnoiselay);
442
443 m_hnoisenhits = new TH1F("nhits_noise", "nhits_noise", 6796, -0.5, 6795.5);
444 m_fdcom -> Add(m_hnoisenhits);
445
446 m_hratio = new TH1F("ratio", "", 100, 0, 1);
447 m_fdcom -> Add(m_hratio);
448
449 m_hdr = new TH1F("dr", "", 2000, -100, 100);
450 m_fdcom -> Add(m_hdr);
451
452 m_hphi0 = new TH1F("phi0", "", 100, 0, 6.284);
453 m_fdcom -> Add(m_hphi0);
454
455 m_hkap = new TH1F("kappa", "", 400, -50, 50);
456 m_fdcom -> Add(m_hkap);
457
458 m_hdz = new TH1F("dz", "", 2000, -500, 500);
459 m_fdcom -> Add(m_hdz);
460
461 m_htanl = new TH1F("tanl", "", 200, -5, 5);
462 m_fdcom -> Add(m_htanl);
463
464 m_hcosthe = new TH1F("costheta", "", 200, -1, 1);
465 m_fdcom -> Add(m_hcosthe);
466
467 m_hcostheNeg = new TH1F("costhetaNeg", "", 200, -1, 1);
468 m_fdcom -> Add(m_hcostheNeg);
469
470 m_hcosthePos = new TH1F("costhetaPos", "", 200, -1, 1);
471 m_fdcom -> Add(m_hcosthePos);
472
473 m_hx0 = new TH1F("x0", "", 100, -10, 10);
474 m_fdcom -> Add(m_hx0);
475
476 m_hy0 = new TH1F("y0", "", 100, -10, 10);
477 m_fdcom -> Add(m_hy0);
478
479 m_hdelZ0 = new TH1F("delta_z0", "", 100, -50, 50);
480 m_fdcom -> Add(m_hdelZ0);
481
482 m_grX0Y0 = new TGraph();
483 m_grX0Y0->SetName("x0y0");
484 m_fdcom -> Add(m_grX0Y0);
485
486 m_hitEffAll = new TH1F("hitEffAll", "", 6800, -0.5, 6799.5);
487 m_fdcom -> Add(m_hitEffAll);
488
489 m_hitEffRaw = new TH1F("hitEffRaw", "", 6800, -0.5, 6799.5);
490 m_fdcom -> Add(m_hitEffRaw);
491
492 m_hitEffRec = new TH1F("hitEffRec", "", 6800, -0.5, 6799.5);
493 m_fdcom -> Add(m_hitEffRec);
494
495 // histograms for drift time
496 m_fdTime = new TFolder("time", "time");
497 m_hlist -> Add(m_fdTime);
498
499 for(lay=0; lay<m_nlayer; lay++){
501 m_htraw[lay] = new TH1F(hname, "", 1000, 0, 2000);
502 m_fdTime -> Add(m_htraw[lay]);
503
505 m_htdr[lay] = new TH1F(hname, "", 510, -10, 500);
506 m_fdTime -> Add(m_htdr[lay]);
507
508 for (lr=0; lr<2; lr++){
510 m_htdrlr[lay][lr] = new TH1F(hname, "", 510, -10, 500);
511 m_fdTime -> Add(m_htdrlr[lay][lr]);
512 }
513 }
514
515 // histograms of adc
516 m_fdAdc = new TFolder("adc", "adc");
517 m_hlist -> Add(m_fdAdc);
518
519 for(lay=0; lay<m_nlayer; lay++){
521 m_hadc[lay] = new TH1F(hname, "", 1500, 0, 3000);
522 m_fdAdc -> Add(m_hadc[lay]);
523 }
524 // histograms for resolution
525 m_fdres = new TFolder("resolution", "resolution");
526 m_hlist -> Add(m_fdres);
527
528 m_fdresAve = new TFolder("resAve", "resAve");
529 m_hlist -> Add(m_fdresAve);
530
531 for(lay=0; lay<m_nlayer; lay++){
533 m_hresInc[lay] = new TH1F(hname, "", 1000, -5, 5);
534 m_fdres -> Add(m_hresInc[lay]);
535
537 m_hresExc[lay] = new TH1F(hname, "", 1000, -5, 5);
538 m_fdres -> Add(m_hresExc[lay]);
539
541 m_hresAve[lay] = new TH1F(hname, "", 1000, -5, 5);
542 m_fdresAve -> Add(m_hresAve[lay]);
543
544 for (lr=0; lr<2; lr++){
546// m_hreslrInc[lay][lr] = new TH1F(hname, "", 200, -1, 1);
547 m_hreslrInc[lay][lr] = new TH1F(hname, "", 1000, -5, 5);
548 m_fdres->Add(m_hreslrInc[lay][lr]);
549
551// m_hreslrExc[lay][lr] = new TH1F(hname, "", 200, -1, 1);
552 m_hreslrExc[lay][lr] = new TH1F(hname, "", 1000, -5, 5);
553 m_fdres->Add(m_hreslrExc[lay][lr]);
554
556// m_hreslrAve[lay][lr] = new TH1F(hname, "", 200, -1, 1);
557 m_hreslrAve[lay][lr] = new TH1F(hname, "", 1000, -5, 5);
558 m_fdresAve->Add(m_hreslrAve[lay][lr]);
559 }
560 for(int phi=0; phi<20; phi++){
562 m_hresphi[lay][phi] = new TH1F(hname, "", 200, -1, 1);
563 m_fdres->Add(m_hresphi[lay][phi]);
564 }
565 }
566
567 /* histograms for momentum vs phi */
568 m_fdmomPhi = new TFolder("momPhi", "momPhi");
569 m_hlist -> Add(m_fdmomPhi);
570
571 int thbin;
575 m_fdmomPhi->Add(m_ppPhi[bin]);
576
579 m_fdmomPhi->Add(m_pnPhi[bin]);
580
583 m_fdmomPhi->Add(m_ppPhiCms[bin]);
584
587 m_fdmomPhi->Add(m_pnPhiCms[bin]);
588
589 for(thbin=0; thbin<NThetaBin; thbin++){
592 m_fdmomPhi->Add(m_ppThePhi[thbin][bin]);
593
596 m_fdmomPhi->Add(m_pnThePhi[thbin][bin]);
597
600 m_fdmomPhi->Add(m_ppThePhiCms[thbin][bin]);
601
604 m_fdmomPhi->Add(m_pnThePhiCms[thbin][bin]);
605 }
606 }
607 for(thbin=0; thbin<NThetaBin; thbin++){
609 m_ppThe[thbin] = new TH1F(hname, "", nbinMom, 0.0, 3.0);
610 m_fdmomPhi->Add(m_ppThe[thbin]);
611
613 m_pnThe[thbin] = new TH1F(hname, "", nbinMom, 0.0, 3.0);
614 m_fdmomPhi->Add(m_pnThe[thbin]);
615
617 m_ppTheCms[thbin] = new TH1F(hname, "", nbinMom, 0.0, 3.0);
618 m_fdmomPhi->Add(m_ppTheCms[thbin]);
619
621 m_pnTheCms[thbin] = new TH1F(hname, "", nbinMom, 0.0, 3.0);
622 m_fdmomPhi->Add(m_pnTheCms[thbin]);
623 }
624
625 // histograms for resolution vs distance
626 m_fdres2d = new TFolder("res2d", "res2d");
627 m_hlist -> Add(m_fdres2d);
628
629 int hid = 0;
631 TH1F* hist;
632 for(lay=0; lay<m_nlayer; lay++){
634 for(lr=0; lr<2; lr++){
637 hist = new TH1F(hname, "", 200, -1, 1);
638 m_hr2dInc.push_back(hist);
639 m_fdres2d -> Add(hist);
640
642 hist = new TH1F(hname, "", 200, -1, 1);
643 m_hr2dExc.push_back(hist);
644 m_fdres2d -> Add(hist);
645
648 hid++;
649 }
650 }
651 }
652 } // end of layer loop
653
654 m_fdres2t = new TFolder("res2t", "res2t");
655// m_hlist -> Add(m_fdres2t);
656
657 for(lay=0; lay<m_nlayer; lay++){
659 for(lr=0; lr<2; lr++){
664 }
665 }
666 }
667 }
668
669 INTupleSvc* ntupleSvc;
671 for(lay=0; lay<m_nlayer; lay++){
673 NTuplePtr nt(ntupleSvc, hname);
674 if ( nt ) m_xtTuple[lay] = nt;
675 else{
677 if( m_xtTuple[lay] ){
678 m_xtTuple[lay]->addItem("cel", m_cel[lay]);
679 m_xtTuple[lay]->addItem("lr", m_lr[lay]);
680 m_xtTuple[lay]->addItem("run", m_run[lay]);
681 m_xtTuple[lay]->addItem("evt", m_evt[lay]);
682 m_xtTuple[lay]->addItem("doca", m_doca[lay]);
683 m_xtTuple[lay]->addItem("dm", m_dm[lay]);
684 m_xtTuple[lay]->addItem("tdr", m_tdr[lay]);
685 m_xtTuple[lay]->addItem("tdc", m_tdc[lay]);
686 m_xtTuple[lay]->addItem("entr", m_entr[lay]);
687 m_xtTuple[lay]->addItem("zhit", m_zhit[lay]);
688 m_xtTuple[lay]->addItem("qhit", m_qhit[lay]);
689 m_xtTuple[lay]->addItem("p", m_p[lay]);
690 m_xtTuple[lay]->addItem("pt", m_pt[lay]);
691 m_xtTuple[lay]->addItem("phi0", m_phi0[lay]);
692 m_xtTuple[lay]->addItem("tanl", m_tanl[lay]);
693 m_xtTuple[lay]->addItem("hitphi", m_hitphi[lay]);
694 } else{
695 log << MSG::FATAL << "Cannot book Xt N-tuple:"
696 << long(m_xtTuple[lay]) << endreq;
697 }
698 }
699 }
700
703 NTuplePtr nt(ntupleSvc, hname);
704 if ( nt ) m_cosmic = nt;
705 else{
707 if( m_cosmic ){
708 m_cosmic->addItem("pUp", m_pUpcos);
709 m_cosmic->addItem("pDw", m_pDwcos);
710 m_cosmic->addItem("ptUp", m_ptUpcos);
711 m_cosmic->addItem("ptDw", m_ptDwcos);
712 m_cosmic->addItem("phiUp", m_phiUpcos);
713 m_cosmic->addItem("phiDw", m_phiDwcos);
714 m_cosmic->addItem("drUp", m_drUpcos);
715 m_cosmic->addItem("drDw", m_drDwcos);
716 m_cosmic->addItem("dzUp", m_dzUpcos);
717 m_cosmic->addItem("dzDw", m_dzDwcos);
718 m_cosmic->addItem("ctheUp", m_ctheUpcos);
719 m_cosmic->addItem("ctheDw", m_ctheDwcos);
720 m_cosmic->addItem("nhitUp", m_nhitUpcos);
721 m_cosmic->addItem("nhitDw", m_nhitDwcos);
722 m_cosmic->addItem("char", m_chargecos);
723 m_cosmic->addItem("tesfg", m_tesFlagcos);
724 m_cosmic->addItem("tes", m_tescos);
725 }
726 }
727 }
728}
Referenced by GrXtMdcCalib::initialize(), QtMdcCalib::initialize(), T0MdcCalib::initialize(), Wr2dMdcCalib::initialize(), WrMdcCalib::initialize(), XtInteMdcCalib::initialize(), and XtMdcCalib::initialize().
◆ printCut()
|
pure virtual |
Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.
Definition at line 1407 of file MdcCalib.cxx.
1407 {
1408 cout << "Events " << m_hTesAll->GetEntries() << ", Tes cut " << m_hTesCalFlag->GetEntries()
1409 << ", nTrkCut " << m_cut1 << endl;
1410 cout << "Tot Tracks " << m_nTrkAfTrkCut
1411 << ", cos(theta)_cut " << m_cut2 << ", drCut " << m_cut3
1412 << ", dzCut " << m_cut4 << ", momCut" << m_cut5
1413 << ", nHitLayer_cut " << m_cut6 << ", nHit_cut " << m_cut7 << endl;
1414 cout << "nTrack after all cuts: " << m_nTrkCal << endl;
1415}
Referenced by GrXtMdcCalib::printCut(), QtMdcCalib::printCut(), T0MdcCalib::printCut(), Wr2dMdcCalib::printCut(), WrMdcCalib::printCut(), XtInteMdcCalib::printCut(), and XtMdcCalib::printCut().
◆ setParam()
|
inlinepure virtual |
Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.
Definition at line 306 of file MdcCalib.h.
306 {
307 m_param = param;
308}
Referenced by GrXtMdcCalib::setParam(), IniMdcCalib::setParam(), PreT0MdcCalib::setParam(), PreXtMdcCalib::setParam(), QtMdcCalib::setParam(), T0MdcCalib::setParam(), Wr2dMdcCalib::setParam(), WrMdcCalib::setParam(), XtInteMdcCalib::setParam(), and XtMdcCalib::setParam().
◆ updateConst()
|
pure virtual |
Implemented in GrXtMdcCalib, IniMdcCalib, PreT0MdcCalib, PreXtMdcCalib, QtMdcCalib, T0MdcCalib, Wr2dMdcCalib, WrMdcCalib, XtInteMdcCalib, and XtMdcCalib.
Definition at line 1417 of file MdcCalib.cxx.
1417 {
1418 IMessageSvc* msgSvc;
1421 log << MSG::DEBUG << "MdcCalib::updateConst()" << endreq;
1422
1423 int lay;
1424 double effi;
1425 double effErr;
1426
1427 int nGoodAll = 0;
1428 int nGoodInn = 0;
1429 int nGoodStp = 0;
1430 int nGoodOut = 0;
1431 int nTotAll = 0;
1432 int nTotInn = 0;
1433 int nTotStp = 0;
1434 int nTotOut = 0;
1436 for(lay=0; lay<m_nlayer; lay++){
1437 double effNtrk = m_effNtrk->GetBinContent(lay+1);
1438 double effGoodHit = m_effNtrkRecHit->GetBinContent(lay+1);
1439 nGoodAll += effGoodHit;
1440 if(lay < 8) nGoodInn += effGoodHit;
1441 else if (lay < 20) nGoodStp += effGoodHit;
1442 else nGoodOut += effGoodHit;
1443
1444 nTotAll += effNtrk;
1445 if(lay < 8) nTotInn += effNtrk;
1446 else if (lay < 20) nTotStp += effNtrk;
1447 else nTotOut += effNtrk;
1448
1449 effi = (double)effGoodHit / (double)effNtrk;
1450 effErr = sqrt(effi * (1-effi) / (double)effNtrk);
1451 feffi << setw(5) << lay << setw(15) << effi << setw(15) << effErr
1452 << setw(15) << effGoodHit << setw(15) << effNtrk << endl;
1453 }
1454 double effiAll = (double)nGoodAll / (double)(nTotAll);
1455 double errAll = sqrt(effiAll * (1-effiAll) / (double)(nTotAll));
1456 double effiInn = (double)nGoodInn / (double)(nTotInn);
1457 double errInn = sqrt(effiInn * (1-effiInn) / (double)(nTotInn));
1458 double effiStp = (double)nGoodStp / (double)(nTotStp);
1459 double errStp = sqrt(effiStp * (1-effiStp) / (double)(nTotStp));
1460 double effiOut = (double)nGoodOut / (double)(nTotOut);
1461 double errOut = sqrt(effiOut * (1-effiOut) / (double)(nTotOut));
1462 feffi << endl << "EffiAll: " << setw(15) << effiAll << setw(15) << errAll
1463 << setw(15) << nGoodAll << setw(15) << nTotAll << endl;
1464 feffi << endl << "EffiInn: " << setw(15) << effiInn << setw(15) << errInn
1465 << setw(15) << nGoodInn << setw(15) << nTotInn << endl;
1466 feffi << endl << "EffiStp: " << setw(15) << effiStp << setw(15) << errStp
1467 << setw(15) << nGoodStp << setw(15) << nTotStp << endl;
1468 feffi << endl << "EffiOut: " << setw(15) << effiOut << setw(15) << errOut
1469 << setw(15) << nGoodOut << setw(15) << nTotOut << endl;
1470 feffi.close();
1471
1472 // calculate efficiency without the impact of track fitting
1474 int nHitAll[] = {0, 0};
1475 int nHitInn[] = {0, 0};
1476 int nHitStp[] = {0, 0};
1477 int nHitOut[] = {0, 0};
1479 for(lay=0; lay<m_nlayer; lay++){
1480 nHitAll[0] += m_hitNum[lay][0];
1481 nHitAll[1] += m_hitNum[lay][1];
1482 if(lay < 8){
1483 nHitInn[0] += m_hitNum[lay][0];
1484 nHitInn[1] += m_hitNum[lay][1];
1485 } else if (lay < 20){
1486 nHitStp[0] += m_hitNum[lay][0];
1487 nHitStp[1] += m_hitNum[lay][1];
1488 } else{
1489 nHitOut[0] += m_hitNum[lay][0];
1490 nHitOut[1] += m_hitNum[lay][1];
1491 }
1492
1493 effi = (double)m_hitNum[lay][1] / (double)m_hitNum[lay][0];
1494 effErr = sqrt(effi * (1-effi) / (double)m_hitNum[lay][0]);
1495 feff2 << setw(5) << lay << setw(15) << effi << setw(15) << effErr
1496 << setw(15) << m_hitNum[lay][1] << setw(15) << m_hitNum[lay][0] << endl;
1497 }
1498 effiAll = (double)nHitAll[1] / (double)nHitAll[0];
1499 errAll = sqrt(effiAll * (1-effiAll)) / (double)nHitAll[0];
1500 effiInn = (double)nHitInn[1] / (double)nHitInn[0];
1501 errInn = sqrt(effiInn * (1-effiInn)) / (double)nHitInn[0];
1502 effiStp = (double)nHitStp[1] / (double)nHitStp[0];
1503 errStp = sqrt(effiStp * (1-effiStp)) / (double)nHitStp[0];
1504 effiOut = (double)nHitOut[1] / (double)nHitOut[0];
1505 errOut = sqrt(effiOut * (1-effiOut)) / (double)nHitOut[0];
1506 feff2 << endl << "EffiAll: " << setw(15) << effiAll << setw(15) << errAll
1507 << setw(15) << nHitAll[1] << setw(15) << nHitAll[0] << endl;
1508 feff2 << endl << "EffiInn: " << setw(15) << effiInn << setw(15) << errInn
1509 << setw(15) << nHitInn[1] << setw(15) << nHitInn[0] << endl;
1510 feff2 << endl << "EffiStp: " << setw(15) << effiStp << setw(15) << errStp
1511 << setw(15) << nHitStp[1] << setw(15) << nHitStp[0] << endl;
1512 feff2 << endl << "EffiOut: " << setw(15) << effiOut << setw(15) << errOut
1513 << setw(15) << nHitOut[1] << setw(15) << nHitOut[0] << endl;
1514 feff2.close();
1515 }
1516
1517 // get resolution
1518 int i;
1519 int iEntr;
1520 int lr;
1523 int hid;
1524
1525 Stat_t entry;
1529 for(lay=0; lay<m_nlayer; lay++){
1530 for(iEntr=0; iEntr<m_nEntr[lay]; iEntr++){
1531 for(lr=0; lr<2; lr++){
1532 fr2d << setw(3) << lay << setw(3) << iEntr << setw(3) << lr << endl;
1535 hid = m_mapr2d[key];
1536
1538 entry = m_hr2dExc[hid] -> GetEntries();
1539 if(entry > 500){
1542 } else if(entry > 100){
1543 sigm[bin] = m_hr2dExc[hid]->GetRMS();
1544 } else{
1545 sigm[bin] = 0.2;
1546 }
1547 } else{
1548 entry = m_hr2dInc[hid] -> GetEntries();
1549 if(entry > 500){
1552 } else if(entry > 100){
1553 sigm[bin] = m_hr2dInc[hid]->GetRMS();
1554 } else{
1555 sigm[bin] = 0.2;
1556 }
1557 }
1559 } // end of bin loop
1560
1562 sigm[bin] = sigm[m_nBin[lay]-1];
1563 }
1564
1567// calconst -> resetSdpar(lay, iEntr, lr, bin, sigm[bin]);
1568 if(1 == m_nEntr[lay]){
1570 } else if(2 == m_nEntr[lay]){
1571 if(0 == iEntr){
1572 for(i=0; i<3; i++){ // entr<0
1574 }
1575 } else{
1576 for(i=3; i<6; i++){ // entr>0
1578 }
1579 }
1580 }
1581 } else{
1583 }
1585 } // end of bin loop
1586 }
1587 } // end of entr loop
1588 }
1589 fr2d.close();
1590 }
1591
1592 return 1;
1593}
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().
The documentation for this class was generated from the following files:
