BOSS 7.0.5
BESIII Offline Software System
Loading...
Searching...
No Matches
HoughPeak.cxx
Go to the documentation of this file.
1#include "MdcHoughFinder/HoughPeak.h"
2#include<fstream>
3#include<sstream>
4
5HoughPeak::HoughPeak(){
6}
7HoughPeak::~HoughPeak(){
8// cout<<" HoughPeak destruct" <<endl;
9 for(int i =0;i<_houghPeakHitList.size();i++){
10 // delete _houghPeakHitList[i];
11 _houghPeakHitList[i]=NULL;
12 }
13 _houghPeakHitList.clear();
14}
15HoughPeak& HoughPeak::operator=(const HoughPeak& other)
16{
17 if( this == &other ) return *this;
18 _height=other._height;
19 _theta=other._theta;
20 _rho=other._rho;
21 _thetaBin=other._thetaBin;
22 _rhoBin=other._rhoBin;
23 _isCandiTrack=other._isCandiTrack;
24 _peakNum=other._peakNum;
25 _charge=other._charge;
26 if(_houghPeakHitList.size() != 0 ){
27 for(int i =0;i<_houghPeakHitList.size();i++){
28 _houghPeakHitList[i]=NULL;
29 }
30 _houghPeakHitList.clear();
31 }
32 _houghPeakHitList = other._houghPeakHitList;
33 //for(int i =0;i<other._houghPeakHitList.size();i++){
34 // const HoughHit* p_hit= new HoughHit( *(other._houghPeakHitList[i]));
35 // _houghPeakHitList.push_back(p_hit);
36 //}
37 return *this;
38}
39
40HoughPeak::HoughPeak(const HoughPeak& other)
41{
42 if(_houghPeakHitList.size() != 0 ){
43 for(int i =0;i<_houghPeakHitList.size();i++){
44 _houghPeakHitList[i]=NULL;
45 }
46 _houghPeakHitList.clear();
47 }
48 _height=other._height;
49 _theta=other._theta;
50 _rho=other._rho;
51 _thetaBin=other._thetaBin;
52 _rhoBin=other._rhoBin;
53 _houghPeakHitList = other._houghPeakHitList;
54 _isCandiTrack=other._isCandiTrack;
55 _peakNum=other._peakNum;
56 _charge=other._charge;
57 //for(int i =0;i<other._houghPeakHitList.size();i++){
58 // const HoughHit* p_hit= new HoughHit( *(other._houghPeakHitList[i]));
59 // _houghPeakHitList.push_back(p_hit);
60 //}
61}
62
63HoughPeak::HoughPeak(int itheta ,int irho, double theta,double rho,vector< const HoughHit* >** mapHitList,bool isCandiTrack,int peakNum){
64 _thetaBin=itheta;
65 _rhoBin=irho;
66 _theta=theta;
67 _rho=rho;
68 _houghPeakHitList = mapHitList[itheta][irho];
69 _isCandiTrack =isCandiTrack;
70 _peakNum=peakNum;
71 // int t_size=mapHitList[itheta][irho].size();
72 // for(int i =0;i<t_size;i++){
73 // const HoughHit* p_hit= new HoughHit( *(mapHitList[itheta][irho][i]));
74 // _houghPeakHitList.push_back(p_hit);
75 // }
76}
77
78HoughPeak::HoughPeak(int height,int itheta ,int irho, double theta,double rho,bool isCandiTrack,int peakNum,int charge){
79 _height=height;
80 _thetaBin=itheta;
81 _rhoBin=irho;
82 _theta=theta;
83 _rho=rho;
84// _houghPeakHitList = mapHitList[itheta][irho];
85 _isCandiTrack =isCandiTrack;
86 _peakNum=peakNum;
87 _charge=charge;
88// collectHits(hitlist);
89}
90
91void HoughPeak::printAllHit() const{
92 int size=_houghPeakHitList.size();
93 if( size ==0 ) return ;
94 cout<<"Print hitlist in HoughPeak " <<endl;
95 cout<<"at center ("<<_theta<<" ,"<<_rho<<"): "<<size<<endl;
96 for(int i=0;i<size;i++){
97 int layer= _houghPeakHitList[i]->getLayerId();
98 int wire = _houghPeakHitList[i]->getWireId();
99 std::cout<<"Peak ("<<layer<<","<<wire<<") "<< _houghPeakHitList[i]<<std::endl;
100 }
101}
102
103int HoughPeak::getHitNum(int select ) const{
104 int size=_houghPeakHitList.size();
105 int n0,n1,n2,n3,n4,n5,n6;
106 n0=n1=n2=n3=n4=n5=n6=0;
107 for(int i=0;i<size;i++){
108 int cir= _houghPeakHitList[i]->getCirList();
109 int index = _houghPeakHitList[i]->digi()->getTrackIndex();
110 int style= _houghPeakHitList[i]->getStyle();
111 // cout<<"peak compare ("<<_houghPeakHitList[i]->getLayerId()<<" ," <<_houghPeakHitList[i]->getWireId()<<" )"<<endl;
112 n0++;
113 if( style == 1 ) n1++;
114 if( style == 2 ) n2++;
115 if( index < 0 ) n3++;
116 if( index >=0 && cir<=1 && style==0 ) n4++;
117 if( index <0 || cir>1 ) n5++;
118 if( index >=0 && cir==0 && style==0 ) n6++;
119 }
120 if( select == 0 ) return n0;
121 if( select == 1 ) return n1;
122 if( select == 2 ) return n2;
123 if( select == 3 ) return n3;
124 if( select == 4 ) return n4;
125 if( select == 5 ) return n5;
126 if( select == 6 ) return n6;
127}
128int HoughPeak::getHitNumA(int select ) const{
129 int size=_houghPeakHitList.size();
130 int n0,n1,n2,n3,n4,n5,n6;
131 n0=n1=n2=n3=n4=n5=n6=0;
132 for(int i=0;i<size;i++){
133 int cir= _houghPeakHitList[i]->getCirList();
134 int index = _houghPeakHitList[i]->digi()->getTrackIndex();
135 int type = _houghPeakHitList[i]->getSlayerType();
136 int style= _houghPeakHitList[i]->getStyle();
137 if( type ==0 ) n0++;
138 if( type==0 && style == 1 ) n1++;
139 if( type==0 && style == 2 ) n2++;
140 if( type==0 && index < 0 ) n3++;
141 if( type==0 && index >=0 && cir<=1 && style==0) n4++;
142 if( type==0 && (index <0 || cir>1) ) n5++;
143 if( type==0 && index >=0 && cir==0 && style==0) n6++;
144 }
145 if( select == 0 ) return n0;
146 if( select == 1 ) return n1;
147 if( select == 2 ) return n2;
148 if( select == 3 ) return n3;
149 if( select == 4 ) return n4;
150 if( select == 5 ) return n5;
151 if( select == 6 ) return n6;
152}
153int HoughPeak::getHitNumS(int select ) const{
154 int size=_houghPeakHitList.size();
155 int n0,n1,n2,n3,n4,n5,n6;
156 n0=n1=n2=n3=n4=n5=n6=0;
157 for(int i=0;i<size;i++){
158 int cir= _houghPeakHitList[i]->getCirList();
159 int index = _houghPeakHitList[i]->digi()->getTrackIndex();
160 int type = _houghPeakHitList[i]->getSlayerType();
161 int style= _houghPeakHitList[i]->getStyle();
162 if( type !=0 ) n0++;
163 if( type!=0 && style == 1 ) n1++;
164 if( type!=0 && style == 2 ) n2++;
165 if( type!=0 && index < 0 ) n3++;
166 if( type!=0 && index >=0 && cir<=1 && style==0) n4++;
167 if( type!=0 && (index <0 || cir>1) ) n5++;
168 if( type!=0 && index >=0 && cir==0 && style==0) n6++;
169 }
170 if( select == 0 ) return n0;
171 if( select == 1 ) return n1;
172 if( select == 2 ) return n2;
173 if( select == 3 ) return n3;
174 if( select == 4 ) return n4;
175 if( select == 5 ) return n5;
176 if( select == 6 ) return n6;
177}
178
179int HoughPeak::collectHits(const HoughHitList& hitlist){
180
181 double cirr_track = fabs(1./(_rho));
182 double cirx_track = (1./_rho)*cos(_theta);
183 double ciry_track = (1./_rho)*sin(_theta);
184 std::vector<HoughHit>::const_iterator iter= hitlist.getList().begin();
185 for(; iter!= hitlist.getList().end(); iter++){
186 const HoughHit *hit = &(*iter);
187 int layer = hit->getLayerId();
188 if( hit->getSlayerType()!=0 ) continue;
189 double cirx_hit = hit->getMidX();
190 double ciry_hit = hit->getMidY();
191 double cirr_hit = hit->getDriftDist();
192 if( _charge == 1 && (cirx_track*ciry_hit - ciry_track*cirx_hit>=0) ) continue; // suppose charge -1 FIXME
193 if( _charge ==-1 && (cirx_track*ciry_hit - ciry_track*cirx_hit<=0) ) continue; // suppose charge -1 FIXME
194 double l1l2 = sqrt( (cirx_hit-cirx_track)*(cirx_hit-cirx_track)+(ciry_hit-ciry_track)*(ciry_hit-ciry_track) );
195 double deltaD ;
196 if( l1l2>cirr_track ) deltaD = l1l2-cirr_track-cirr_hit;
197 if( l1l2<=cirr_track ) deltaD = l1l2-cirr_track+cirr_hit;
198
199 //cal flight length
200
201 double theta_temp;
202 double l_temp;
203 if(cirx_track==0||cirx_hit-cirx_track==0){
204 theta_temp=0;
205 }
206 else{
207 theta_temp=M_PI-atan2(ciry_hit-ciry_track,cirx_hit-cirx_track)+atan2(ciry_track,cirx_track);
208 if(theta_temp>2*M_PI){
209 theta_temp=theta_temp-2*M_PI;
210 }
211 if(theta_temp<0){
212 theta_temp=theta_temp+2*M_PI;
213 }
214 }
215 //cout<<" charge "<<_charge <<" "<<theta_temp<<endl;
216 if(_charge==-1) {
217 l_temp=cirr_track*theta_temp;
218 }
219 if(_charge==1) {
220 theta_temp=2*M_PI-theta_temp;
221 l_temp=cirr_track*(theta_temp);
222 }
223 double pt = fabs( (1./_rho)/333.567 );
224 //if( fabs(deltaD)<0.2 && l_temp<=35 ) _houghPeakHitList.push_back(hit); //suppose 60MeV FIXME
225 //if( fabs(deltaD)<0.2 && l_temp>35 && l_temp<=45) _houghPeakHitList.push_back(hit);
226
227 if( pt<0.06 && fabs(deltaD)<0.1&& l_temp<=35 ) _houghPeakHitList.push_back(hit);
228 if( 0.06<pt&&pt<0.07 && fabs(deltaD)<0.1 && l_temp<=35 ) _houghPeakHitList.push_back(hit);
229 if( 0.07<pt&&pt<0.08 && fabs(deltaD)<0.1 && l_temp<=43 ) _houghPeakHitList.push_back(hit);
230 if( 0.08<pt&&pt<0.09 && fabs(deltaD)<0.1 && l_temp<=43 ) _houghPeakHitList.push_back(hit);
231 if( 0.09<pt&&pt<0.10 && fabs(deltaD)<0.1 && l_temp<=41 ) _houghPeakHitList.push_back(hit);
232 if( 0.10<pt&&pt<0.11 && fabs(deltaD)<0.1 && l_temp<=41 ) _houghPeakHitList.push_back(hit);
233 if( 0.11<pt&&pt<0.12 && fabs(deltaD)<0.1 && l_temp<=41 ) _houghPeakHitList.push_back(hit);
234
235 if( pt<0.06 && fabs(deltaD)<0.1&& l_temp>35&&l_temp<45 ) _houghPeakHitList.push_back(hit);
236 if( 0.06<pt&&pt<0.07 && fabs(deltaD)<0.1 && l_temp>35&& l_temp<=45) _houghPeakHitList.push_back(hit);
237 if( 0.07<pt&&pt<0.08 && fabs(deltaD)<0.1 && l_temp>43&& l_temp<=50) _houghPeakHitList.push_back(hit);
238 if( 0.08<pt&&pt<0.09 && fabs(deltaD)<0.1 && l_temp>43&& l_temp<=50) _houghPeakHitList.push_back(hit);
239 // if( 0.09<pt&&pt<0.10 && fabs(deltaD)<0. && l_temp>41&& l_temp<=50) _houghPeakHitList.push_back(hit);
240 // if( 0.10<pt&&pt<0.11 && fabs(deltaD)<0. && l_temp>80&& l_temp<=105) _houghPeakHitList.push_back(hit);
241 // if( 0.11<pt&&pt<0.12 && fabs(deltaD)<0. && l_temp>80&& l_temp<=105) _houghPeakHitList.push_back(hit);
242 if( pt>0.12 && fabs(deltaD)<0.1) _houghPeakHitList.push_back(hit);
243
244
245 }
246 return _houghPeakHitList.size();
247
248}
iter
EvtStreamInputIterator< typename Generator::result_type > iter(Generator gen, int N=0)
Definition: EvtStreamInputIterator.hh:100
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
n2
int n2
Definition: SD0Tag.cxx:55
n1
int n1
Definition: SD0Tag.cxx:54
M_PI
#define M_PI
Definition: TConstant.h:4
HoughHitList::getList
const std::vector< HoughHit > & getList() const
Definition: HoughHitList.cxx:344
HoughPeak::getHitNumS
int getHitNumS(int) const
Definition: HoughPeak.cxx:153
HoughPeak::getHitNum
int getHitNum(int) const
Definition: HoughPeak.cxx:103
HoughPeak::~HoughPeak
~HoughPeak()
Definition: HoughPeak.cxx:7
HoughPeak::operator=
HoughPeak & operator=(const HoughPeak &other)
Definition: HoughPeak.cxx:15
HoughPeak::printAllHit
void printAllHit() const
Definition: HoughPeak.cxx:91
HoughPeak::collectHits
int collectHits(const HoughHitList &)
Definition: HoughPeak.cxx:179
HoughPeak::getHitNumA
int getHitNumA(int) const
Definition: HoughPeak.cxx:128
charge
float charge
Definition: prepare_data.cxx:28