BOSS 7.0.8
BESIII Offline Software System
Loading...
Searching...
No Matches
MucRec2DRoad Class Reference

#include <MucRec2DRoad.h>

+ Inheritance diagram for MucRec2DRoad:

Public Member Functions

 MucRec2DRoad (const int &part=0, const int &seg=0, const int &orient=0, const float &xVertex=0.0, const float &yVertex=0.0, const float &zVertex=0.0, const int &fittingMethod=0)
 Constructor.
 
MucRec2DRoadoperator= (const MucRec2DRoad &orig)
 Assignment constructor.
 
 MucRec2DRoad (const MucRec2DRoad &source)
 Copy constructor.
 
 ~MucRec2DRoad ()
 Destructor.
 
void SetIndex (const int &index)
 Set the index for this road.
 
void AttachHit (MucRecHit *hit)
 Attach the given hit to this road.
 
void AttachHitNoFit (MucRecHit *hit)
 Attach the given hit to this road, but not fit this road.
 
void SetMaxNSkippedGaps (const int &nGaps)
 
int SimpleFit (float &slope, float &intercept, float &sigmaSlope, float &sigmaIntercept, float &chisq, int &ndof)
 Calculate the best-fit straight line with "simple" weights.
 
int SimpleFitNoCurrentGap (int currentgap, float &slope, float &intercept, float &sigmaSlope, float &sigmaIntercept, float &chisq, int &ndof)
 Calculate the best-fit straight line with "simple" weights. not use current gap!!!
 
int Fit (const float &x, const float &y, const float &z, float &slope, float &intercept, float &sigmaSlope, float &sigmaIntercept, float &chisq, int &ndof)
 Fit (with weights) the hit center to a straight line.
 
void Project (const int &gap, float &x, float &y, float &z, float &x2, float &y2, float &z2)
 Where does the trajectory of this road intersect a specific gap?
 
float WeightFunc (const float &chisq, const float &distance) const
 
float WindowFunc (const float &chisq, const float &distance) const
 
int GetIndex () const
 A unique identifier for this road in the current event.
 
int GetPart () const
 In which part was this road found?
 
int GetSeg () const
 In which segment was this road found?
 
int GetOrient () const
 In which orientation was this road found?
 
void GetVertexPos (float &x, float &y, float &z) const
 Position of the vertex point.
 
int GetLastGap () const
 Which gap is the last one with hits attached to this road?
 
int GetNGapsWithHits () const
 How many gaps provide hits attached to this road?
 
int GetTotalHits () const
 How many hits in all does this road contain?
 
int GetHitsPerGap (const int &gap) const
 How many hits per gap does this road contain?
 
int GetMaxHitsPerGap () const
 How many hits were attached in the gap with the most attached hits?
 
bool HasHitInGap (const int &gap) const
 Does this road contain any hits in the given gap?
 
int GetNSharedHits (const MucRec2DRoad *road) const
 How many hits do two roads share?
 
float GetSlope () const
 Slope of trajectory.
 
float GetIntercept () const
 Intercept of trajectory.
 
int GetDegreesOfFreedom () const
 How many degrees of freedom in the trajectory fit?
 
float GetReducedChiSquare () const
 Chi-square parameter (per degree of freedom) of the trajectory fit.
 
void GetSimpleFitParams (float &slope, float &intercept, float &sigmaSlope, float &sigmaIntercept, float &chisq, int &ndof) const
 Get the parameters from the simple fit.
 
void GetSimpleFitParams (float &a, float &b, float &c, int &whichhalf, float &sigmaa, float &sigmab, float &sigmac, float &chisq, int &ndof) const
 Get the parameters from the simple quad fit.
 
bool GetQuadFitOk ()
 
void GetPosSigma (float &possigma) const
 
MucRecHitGetHit (const int &gap) const
 Get a pointer to the first found hit attached in a particular gap.
 
float GetHitDistance (const MucRecHit *hit)
 Distance to a hit.
 
float GetSearchWindowSize (const int &gap) const
 Determine the size of the search window in the given gap.
 
bool HasHit (MucRecHit *hit) const
 Does the hit exist in the road .
 
vector< IdentifierGetHitsID () const
 Get indices of all hits in the road.
 
vector< MucRecHit * > GetHits () const
 
void PrintHitsInfo () const
 Print Hits Infomation.
 

Detailed Description

Describes a "two-dimensional" road found in the muon chamber. part and orient decides dimension of the road. part orient dimension 1 0 Z-R 1 1 Phi-R 0,2 0 Y-Z 0,2 1 X-Z A full "three-dimensional" road is composed of two 2D roads with different orients on the same part.

Author
Zhengyun You \URL{youzy.nosp@m.@hep.nosp@m..pku..nosp@m.cn}

Definition at line 38 of file MucRec2DRoad.h.

Constructor & Destructor Documentation

◆ MucRec2DRoad() [1/2]

MucRec2DRoad::MucRec2DRoad ( const int &  part = 0,
const int &  seg = 0,
const int &  orient = 0,
const float &  xVertex = 0.0,
const float &  yVertex = 0.0,
const float &  zVertex = 0.0,
const int &  fittingMethod = 0 
)

Constructor.

Definition at line 28 of file MucRec2DRoad.cxx.

35 : m_VertexPos(xVertex, yVertex, zVertex),
36 m_VertexSigma(0.0, 0.0, 0.0),
37 m_Part(part), m_Seg(seg), m_Orient(orient),
38 m_Chi2(0.0), m_DOF(0),
39 m_MaxHitsPerGap(0), m_LastGap(0), m_TotalHits(0),
40 m_FitOK(false), m_QuadFitOK(false),
41 m_HitDistance(5, float(muckInfinity)),
42 m_pHits(0), m_fittingMethod(fittingMethod)
43{ }
const double muckInfinity

◆ MucRec2DRoad() [2/2]

MucRec2DRoad::MucRec2DRoad ( const MucRec2DRoad source)

Copy constructor.

Definition at line 72 of file MucRec2DRoad.cxx.

73 : m_VertexPos(source.m_VertexPos),
74 m_VertexSigma(source.m_VertexSigma),
75 m_Index(source.m_Index),
76 m_Part(source.m_Part), m_Seg(source.m_Seg), m_Orient(source.m_Orient),
77 m_Chi2(source.m_Chi2), m_DOF(source.m_DOF),
78 m_MaxHitsPerGap(source.m_MaxHitsPerGap),
79 m_LastGap(source.m_LastGap), m_TotalHits(source.m_TotalHits),
80 m_FitOK(source.m_FitOK),
81 m_HitDistance(source.m_HitDistance),
82 m_pHits(source.m_pHits),
83 m_fittingMethod(source.m_fittingMethod)
84{ }

◆ ~MucRec2DRoad()

MucRec2DRoad::~MucRec2DRoad ( )

Destructor.

Definition at line 87 of file MucRec2DRoad.cxx.

88{ }

Member Function Documentation

◆ AttachHit()

void MucRec2DRoad::AttachHit ( MucRecHit hit)

Attach the given hit to this road.

Definition at line 100 of file MucRec2DRoad.cxx.

101{
102 // cout << "MucRec2DRoad::AttachHit-I0 hit = " << hit << endl;
103 if (!hit) {
104 cout << "MucRec2DRoad::AttachHit-E1 null hit pointer!" << endl;
105 return ;
106 }
107
108 int gap = hit->Gap();
109 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
110 // The gap number of the hit is out of range.
111 cout << "MucRec2DRoad::AttachHit-W3 bad gap number = " << gap
112 << endl;
113 return;
114 }
115
116 // Attach the hit to the road.
117 //bool hitExist = false;
118 for (int i = 0; i < (int)m_pHits.size(); i++) {
119 if (m_pHits[i] == hit) return;
120 }
121 m_pHits.push_back(hit);
122 //cout << "before Count " << m_pHits.size() << endl;
123 // m_HitDistance[gap] = dX;
124
125 // Now recalculate the total number of hits and the max. number of
126 // hits per gap.
127 CountHits();
128 //cout << "after Count " << m_pHits.size() << endl;
129
130 // Redo the "simple" least-squares fit to the positions ...
131 float a, sa, b, sb, chi;
132 int n;
133 int status = SimpleFit(a, b, sa, sb, chi, n);
134 if (status < 0) {
135 //cout << "MucRec2DRoad::AttachHit-E5 simple fit fail status = " << status << endl;
136 }
137
138 //cout << "Hit center = " << hit->GetCenterPos() << endl;
139 //cout << "After attach hit, slope = " << a << " intercept = " << b << endl;
140}
const Int_t n
static value_type getGapMax()
Definition: MucID.cxx:195
int SimpleFit(float &slope, float &intercept, float &sigmaSlope, float &sigmaIntercept, float &chisq, int &ndof)
Calculate the best-fit straight line with "simple" weights.
int Gap() const
Get Gap.
Definition: MucRecHit.h:77
const double b
Definition: slope.cxx:9

Referenced by MucRecRoadFinder::execute(), and RecMucTrack::LineFit().

◆ AttachHitNoFit()

void MucRec2DRoad::AttachHitNoFit ( MucRecHit hit)

Attach the given hit to this road, but not fit this road.

Definition at line 145 of file MucRec2DRoad.cxx.

146{
147 // cout << "MucRec2DRoad::AttachHit-I0 hit = " << hit << endl;
148 if (!hit) {
149 cout << "MucRec2DRoad::AttachHit-E1 null hit pointer!" << endl;
150 return ;
151 }
152
153 int gap = hit->Gap();
154 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
155 // The gap number of the hit is out of range.
156 cout << "MucRec2DRoad::AttachHit-W3 bad gap number = " << gap
157 << endl;
158 return;
159 }
160
161 // Attach the hit to the road.
162 //bool hitExist = false;
163 for (int i = 0; i < (int)m_pHits.size(); i++) {
164 if (m_pHits[i] == hit) return;
165 }
166 m_pHits.push_back(hit);
167 //cout << "before Count " << m_pHits.size() << endl;
168 // m_HitDistance[gap] = dX;
169
170 // Now recalculate the total number of hits and the max. number of
171 // hits per gap.
172 CountHits();
173 //cout << "after Count " << m_pHits.size() << endl;
174
175}

Referenced by MucRecRoadFinder::execute().

◆ Fit()

int MucRec2DRoad::Fit ( const float &  x,
const float &  y,
const float &  z,
float &  slope,
float &  intercept,
float &  sigmaSlope,
float &  sigmaIntercept,
float &  chisq,
int &  ndof 
)

Fit (with weights) the hit center to a straight line.

Definition at line 575 of file MucRec2DRoad.cxx.

581{
582 int status;
583
584 // If the "simple" fit has not been done yet, do it now.
585 if (!m_SimpleFitOK) {
586 // Least-squares fit to the positions ...
587 float a, sa, b, sb, chi;
588 int n;
589 status = SimpleFit(a, b, sa, sb, chi, n);
590 if (status < 0) {
591 //cout << "MucRec2DRoad::Fit-E2 simple fit fail status = "
592 // << status << endl;
593 return status;
594 }
595 }
596
597 // Assign to temporary arrays to be used in fit.
598 float px[100];
599 float py[100];
600 float pw[100];
601 int npoints = 0;
602 float dx, dy, dr;
603
604
605 // include vertex point when do the fancy fitting
606 //if (m_Orient == kHORIZ) {
607 // px[npoints] = m_VertexPos.y();
608 // dx = px[npoints] - y;
609 //} else {
610 // px[npoints] = m_VertexPos.x();
611 // dx = px[npoints] - x;
612 //}
613 //pz[npoints] = m_VertexPos.z();
614 //dz = pz[npoints] - z;
615 //dr = sqrt(dx*dx + dz*dz);
616 //pw[npoints] = WeightFunc(m_SimpleChi2,dr);
617 //npoints++;
618
619 vector<MucRecHit*>::const_iterator iHit;
620
621 // Determine the weights based on the chi-square of the fit
622 // (the scatter of the points should be roughly related to the
623 // momentum) and the distance from the center to the estimated
624 // location.
625
626 //cout << m_pHits.size() << endl;
627 for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
628 if (*iHit) { // Check for a null pointer.
629
630 Hep3Vector center = (*iHit)->GetCenterPos();
631 Hep3Vector sigma = (*iHit)->GetCenterSigma();
632
633 if (m_Part == 1) { // change from 0 to 1 at 2006.11.30
634 if ( m_Orient == 0) {
635 px[npoints] = center.z();
636 dx = px[npoints] - z;
637 py[npoints] = sqrt(center.x()*center.x()
638 + center.y()*center.y());
639 if(m_Seg==2) py[npoints] = center.y(); //deal with seg2 seperately! because there is a hole here. 2006.11.9
640 dy = py[npoints] - sqrt(x*x + y*y);
641 }
642 else {
643 px[npoints] = center.x();
644 dx = px[npoints] - x;
645 py[npoints] = center.y();
646 dy = py[npoints] - y;
647 }
648 }
649 else {
650 if ( m_Orient == 0) {
651 px[npoints] = center.z();
652 dx = px[npoints] - z;
653 py[npoints] = center.y();
654 dy = py[npoints] - y;
655 }
656 else {
657 px[npoints] = center.z();
658 dx = px[npoints] - z;
659 py[npoints] = center.x();
660 dy = py[npoints] - x;
661 }
662 }
663
664 dr = sqrt(dx*dx + dy*dy);
665 pw[npoints] = WeightFunc(m_SimpleChi2, dr);
666 //pw[npoints] = WeightFunc(m_SimpleChi2,dr);
667
668 // cout << " " << npoints << " "
669 // << px[npoints] << " " << py[npoints] << " " << pw[npoints]
670 // << " " << dx << " " << dy << " " << dr
671 // << endl;
672
673 npoints++;
674
675 if(npoints > 99)
676 { cout<<"MucRec2DRoad: to many hits in this track, ignore it!"<<endl;
677 return -1;
678 }
679
680 }
681 }
682
683 // Refit ...
684 ndof = npoints - 2;
685 if (npoints == 1) {
686 if (m_Part == 1) { // change from 0 to 1 at 2006.11.30
687 if ( m_Orient == 0) {
688 px[npoints] = m_VertexPos.z();
689 py[npoints] = sqrt(m_VertexPos.x()*m_VertexPos.x()
690 + m_VertexPos.y()*m_VertexPos.y());
691 if(m_Seg==2) py[npoints] = m_VertexPos.y();
692 }
693 else {
694 px[npoints] = m_VertexPos.x();
695 py[npoints] = m_VertexPos.y();
696 }
697 }
698 else {
699 if ( m_Orient == 0) {
700 px[npoints] = m_VertexPos.z();
701 py[npoints] = m_VertexPos.y();
702 }
703 else {
704 px[npoints] = m_VertexPos.z();
705 py[npoints] = m_VertexPos.x();
706 }
707 }
708 pw[npoints] = 1.0;
709 npoints++;
710 }
711 else {
712 if (npoints == 0) {
713 return -1;
714 }
715 }
716
717 MucRecLineFit fit;
718 //cout << "npoints " << npoints << endl;
719 status = fit.LineFit(px, py, pw, npoints,
720 &slope, &intercept, &chisq,
721 &sigmaSlope, &sigmaIntercept);
722 m_DOF = ndof;
723 m_Chi2 = chisq;
724
725 if (status < 0) {
726 //cout << "MucRec2DRoad::Fit-E3 fit fail status = " << status << endl;
727 }
728
729 return status;
730}
TTree * sigma
Double_t x[10]
float WeightFunc(const float &chisq, const float &distance) const
int LineFit(float x[], float y[], float w[], int n, float *a, float *b, float *chisq, float *siga, float *sigb)
double y[1000]
void slope()
Definition: slope.cxx:12

◆ GetDegreesOfFreedom()

int MucRec2DRoad::GetDegreesOfFreedom ( ) const

How many degrees of freedom in the trajectory fit?

Definition at line 1064 of file MucRec2DRoad.cxx.

1064{ return m_SimpleDOF; }

Referenced by MucRecRoadFinder::execute(), and MucRec3DRoad::MucRec3DRoad().

◆ GetHit()

MucRecHit * MucRec2DRoad::GetHit ( const int &  gap) const

Get a pointer to the first found hit attached in a particular gap.

Definition at line 1131 of file MucRec2DRoad.cxx.

1132{
1133 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
1134 cout << "MucRec2DRoad::Hit-E1 invalid gap = " << gap << endl;
1135 return 0;
1136 }
1137
1138 vector<MucRecHit*>::const_iterator iHit;
1139
1140 for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
1141 if (*iHit) { // Check for a null pointer.
1142 if ( (*iHit)->Gap() == gap) {
1143 return (*iHit);
1144 }
1145 }
1146 }
1147
1148 return 0L;
1149}

Referenced by MucRec3DRoad::GetHit(), and MucRec3DRoad::MucRec3DRoad().

◆ GetHitDistance()

float MucRec2DRoad::GetHitDistance ( const MucRecHit hit)

Distance to a hit.

Definition at line 1154 of file MucRec2DRoad.cxx.

1155{
1156 // Use straight-line projection?
1157 if (!hit) {
1158 cout << "MucRec2DRoad::GetHitDistance-E1 null hit pointer!" << endl;
1159 return muckInfinity;
1160 }
1161
1162 int gap = hit->Gap();
1163 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
1164 cout << "MucRec2DRoad::HitDistance-W2 bad gap number = " << gap << endl;
1165 return muckInfinity;
1166 }
1167
1168 if ( hit->GetGap()->Orient() != m_Orient) {
1169 // The orientation of the hit is different from that of the road.
1170 cout << "MucRec2DRoad::GetHitDistance-W3 "
1171 << " hit has wrong orientation = "
1172 << hit->GetGap()->Orient()
1173 << endl;
1174 return muckInfinity;
1175 }
1176
1177 HepPoint3D r = hit->GetCenterPos();
1178 HepPoint3D rl = hit->GetGap()->TransformToGap(r);
1179 // cout << "hit center " << r << endl;
1180 // cout << "hit center local " << rl << endl;
1181
1182 // Determine the projection point of the "simple" fit to the desired gap.
1183 // FIXME(?): need to use fit with fancy weights instead?
1184 float delta, delta1,delta2;
1185 float x0, y0, z0;
1186 float x2, y2, z2;
1187 //float sigmaX0, sigmaY0, sigmaZ0;
1188
1189 //cout << "y:x = " << m_SimpleSlope << endl;
1190
1191 Project(gap, x0, y0, z0, x2, y2, z2);
1192 // cout << " gap = " << gap << " x0 = " << x0
1193 // << " y0 = " << y0 << " z0 = " << z0
1194 // << endl;
1195
1196 if(x0==y0&&x0==z0&&x0==0) {x0+=-9999;y0+=-9999;z0+=-9999;}//cout<<"wrong intersection"<<endl;}
1197 if(x2==y2&&x2==z2&&x2==0) {x2+=-9999;y2+=-9999;z2+=-9999;}//cout<<"wrong intersection"<<endl;} //wrong intersection!!!
1198 HepPoint3D s = HepPoint3D(x0, y0, z0);
1199 HepPoint3D s_2 = HepPoint3D(x2, y2, z2);
1200 HepPoint3D sl = hit->GetGap()->TransformToGap(s);
1201 HepPoint3D s_2l = hit->GetGap()->TransformToGap(s_2);
1202 //cout << "project center " << s << endl;
1203
1204// cout << "project center local sl= " << sl << endl;
1205// cout << "project center local sl= " << s_2l << endl;
1206// cout << "project center local rl= " << rl << endl;
1207 if ( m_Part == 0 ) {
1208 if ( m_Orient == 0 ) {
1209 delta1 = fabs( sl.y() - rl.y() );
1210 delta2= fabs( s_2l.y() - rl.y() );
1211 }
1212 else {
1213 delta1 = fabs( sl.x() - rl.x() );
1214 delta2= fabs( s_2l.x() - rl.x() );
1215 }
1216 }
1217 else {
1218 if ( m_Orient == 0 ) {
1219 delta1 = fabs( sl.y() - rl.y() );
1220 delta2= fabs( s_2l.y() - rl.y() );
1221 }
1222 else {
1223 delta1 = fabs( sl.x() - rl.x() );
1224 delta2= fabs( s_2l.x() - rl.x() );
1225 }
1226 }
1227// if(which==0) {
1228// if(delta1<delta2)delta = delta1;
1229// else delta = delta2;
1230// }
1231 //cout<<"which = "<<which<<" distance = "<<delta1<<" "<<delta2<<endl;
1232
1233 if(delta1 < delta2) return delta1;
1234 else return delta2;
1235}
const double delta
HepGeom::Point3D< double > HepPoint3D
Definition: Gam4pikp.cxx:37
XmlRpcServer s
Definition: HelloServer.cpp:11
HepPoint3D TransformToGap(const HepPoint3D gPoint) const
Transform a point from global coordinate to gap coordinate.
Definition: MucGeoGap.cxx:627
int Orient() const
Definition: MucGeoGap.h:108
void Project(const int &gap, float &x, float &y, float &z, float &x2, float &y2, float &z2)
Where does the trajectory of this road intersect a specific gap?
MucGeoGap * GetGap() const
Get geometry data for the gap containing this hit.
Definition: MucRecHit.h:83
Hep3Vector GetCenterPos() const
Get Center position of the strip (global coords).
Definition: MucRecHit.cxx:104

Referenced by MucRecRoadFinder::execute(), and MucRec3DRoad::MucRec3DRoad().

◆ GetHits()

vector< MucRecHit * > MucRec2DRoad::GetHits ( ) const

Definition at line 1439 of file MucRec2DRoad.cxx.

1440{
1441return m_pHits;
1442
1443}

Referenced by MucRecRoadFinder::execute().

◆ GetHitsID()

vector< Identifier > MucRec2DRoad::GetHitsID ( ) const

Get indices of all hits in the road.

Definition at line 1418 of file MucRec2DRoad.cxx.

1419{
1420 vector<Identifier> idCon;
1421
1422 vector<MucRecHit*>::const_iterator iHit;
1423 for ( iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
1424 if (*iHit) { // Check for a null pointer.
1425 Identifier id = (*iHit)->GetID();
1426 idCon.push_back(id);
1427 /*
1428 cout << " MucRec2DRoad::HitIndices gap orientation twopack= "
1429 << (*iHit)->ChannelID().Plane() << " "
1430 << (*iHit)->ChannelID().Orient() << " "
1431 << (*iHit)->ChannelID().TwoPack() << endl;
1432 */
1433 }
1434 }
1435 return idCon;
1436}

Referenced by MucRec3DRoad::GetHitsID().

◆ GetHitsPerGap()

int MucRec2DRoad::GetHitsPerGap ( const int &  gap) const

How many hits per gap does this road contain?

Definition at line 977 of file MucRec2DRoad.cxx.

978{
979 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
980 cout << "MucRec2DRoad::GetHitsPerGap-E1 invalid gap = " << gap << endl;
981 return 0;
982 }
983
984 vector<MucRecHit*>::const_iterator iHit;
985 int hitsInGap = 0;
986
987 for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
988
989 if ( !(*iHit) ) {
990 cout << "MucRec2DRoad::GetHitsPerGap()-E2 null hit pointer !" << endl;
991 return 0;
992 }
993 else {
994 if( gap == (*iHit)->Gap() ) hitsInGap += 1;
995 }
996 }
997
998 return hitsInGap;
999}

Referenced by MucRec3DRoad::GetHitsPerGap(), and MucRec3DRoad::GetNGapsWithHits().

◆ GetIndex()

int MucRec2DRoad::GetIndex ( ) const

A unique identifier for this road in the current event.

Definition at line 909 of file MucRec2DRoad.cxx.

910{
911 return m_Index;
912}

◆ GetIntercept()

float MucRec2DRoad::GetIntercept ( ) const

Intercept of trajectory.

Definition at line 1060 of file MucRec2DRoad.cxx.

1060{ return m_SimpleIntercept; }

Referenced by MucRecRoadFinder::execute(), and RecMucTrack::LineFit().

◆ GetLastGap()

int MucRec2DRoad::GetLastGap ( ) const

Which gap is the last one with hits attached to this road?

Definition at line 939 of file MucRec2DRoad.cxx.

939{ return m_LastGap; }

Referenced by MucRecRoadFinder::execute(), MucRec3DRoad::GetLastGapDelta(), and MucRec3DRoad::MucRec3DRoad().

◆ GetMaxHitsPerGap()

int MucRec2DRoad::GetMaxHitsPerGap ( ) const

How many hits were attached in the gap with the most attached hits?

Definition at line 943 of file MucRec2DRoad.cxx.

943{ return m_MaxHitsPerGap; }

Referenced by MucRec3DRoad::GetMaxHitsPerGap().

◆ GetNGapsWithHits()

int MucRec2DRoad::GetNGapsWithHits ( ) const

How many gaps provide hits attached to this road?

Definition at line 951 of file MucRec2DRoad.cxx.

952{
953 const int ngap = (int)MucID::getGapMax();
954 int gap, count = 0;
955 vector<int> firedGap;
956 for ( gap = 0; gap < ngap; gap++) {
957 firedGap.push_back(0);
958 }
959
960 vector<MucRecHit*>::const_iterator iHit;
961 for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
962 if (*iHit) { // Check for a null pointer.
963 gap = (*iHit)->Gap();
964 firedGap[gap] = 1;
965 }
966 }
967
968 for ( gap = 0; gap < ngap; gap++) {
969 count += firedGap[gap];
970 }
971
972 return count;
973}
DOUBLE_PRECISION count[3]

Referenced by MucRecRoadFinder::execute().

◆ GetNSharedHits()

int MucRec2DRoad::GetNSharedHits ( const MucRec2DRoad road) const

How many hits do two roads share?

Definition at line 1026 of file MucRec2DRoad.cxx.

1027{
1028 if (!road2) {
1029 return 0;
1030 }
1031
1032 int count = 0;
1033 vector<MucRecHit*>::const_iterator iHit1;
1034 vector<MucRecHit*>::const_iterator iHit2;
1035 MucRecHit *hit1, *hit2;
1036
1037 for( iHit1 = m_pHits.begin(); iHit1 != m_pHits.end(); iHit1++){
1038 for( iHit2 = road2->m_pHits.begin();
1039 iHit2 != road2->m_pHits.end(); iHit2++){
1040 hit1 = (*iHit1);
1041 hit2 = (*iHit2);
1042
1043 if ( (hit1 != 0) && (hit2 != 0) ) {
1044 if (hit1->GetID() == hit2->GetID()) {
1045 count++;
1046 }
1047 }
1048 }
1049 }
1050
1051 return count;
1052}
Identifier GetID() const
Get soft identifier of this hit.
Definition: MucRecHit.h:68

Referenced by MucRec3DRoad::GetNSharedHits().

◆ GetOrient()

int MucRec2DRoad::GetOrient ( ) const

In which orientation was this road found?

Definition at line 924 of file MucRec2DRoad.cxx.

924{ return m_Orient; }

Referenced by MucRecRoadFinder::execute().

◆ GetPart()

int MucRec2DRoad::GetPart ( ) const

In which part was this road found?

Definition at line 916 of file MucRec2DRoad.cxx.

916{ return m_Part; }

Referenced by MucRecRoadFinder::execute(), and MucRec3DRoad::MucRec3DRoad().

◆ GetPosSigma()

void MucRec2DRoad::GetPosSigma ( float &  possigma) const

Definition at line 1099 of file MucRec2DRoad.cxx.

1100{
1101 possigma = m_SimplePosSigma;
1102
1103}

Referenced by MucRec3DRoad::ProjectNoCurrentGap(), MucRec3DRoad::ProjectWithSigma(), and MucRec3DRoad::RefitNoCurrentGap().

◆ GetQuadFitOk()

bool MucRec2DRoad::GetQuadFitOk ( )
inline

Definition at line 160 of file MucRec2DRoad.h.

160{return m_QuadFitOK;}

Referenced by MucRec3DRoad::Project(), and MucRec3DRoad::ProjectNoCurrentGap().

◆ GetReducedChiSquare()

float MucRec2DRoad::GetReducedChiSquare ( ) const

Chi-square parameter (per degree of freedom) of the trajectory fit.

Definition at line 1068 of file MucRec2DRoad.cxx.

1069{
1070 if ( (!m_SimpleFitOK) || (m_SimpleDOF < 0) ) {
1071 return -1.0;
1072 }
1073 else if (m_SimpleDOF == 0) {
1074 return 0.0;
1075 }
1076 else {
1077 return (m_SimpleChi2 / m_SimpleDOF);
1078 }
1079}

Referenced by MucRecRoadFinder::execute().

◆ GetSearchWindowSize()

float MucRec2DRoad::GetSearchWindowSize ( const int &  gap) const

Determine the size of the search window in the given gap.

Definition at line 1239 of file MucRec2DRoad.cxx.

1240{
1241 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
1242 cout << "MucRec2DRoad::GetSearchWindowSize-E1 invalid gap = " << gap << endl;
1243 return 0.0;
1244 }
1245
1246 // Determine the projection point of the "simple" fit to the last
1247 // gap and the desired gap.
1248 // FIXME? the "delta-x" variable is calculated as the scalar
1249 // difference between the positions obtained by projecting to the
1250 // last gap and to the desired gap.
1251
1253 float x1, y1, z1, x2, y2, z2, dx, dy, dr;
1254 float sigmaX, sigmaY, sigmaZ;
1255
1256 if (m_Part == 0) {
1257 if (m_Orient == 0) {
1258 geom->FindIntersection(m_Part, 0, m_LastGap,
1259 1.0, m_SimpleSlope, 0.0,
1260 0.0, m_SimpleIntercept, 0.0,
1261 0.0, m_SimpleSlopeSigma, 0.0,
1262 0.0, m_SimpleInterceptSigma, 0.0,
1263 x1, y1, z1,
1264 sigmaX, sigmaY, sigmaZ);
1265 geom->FindIntersection(m_Part, 0, gap,
1266 1.0, m_SimpleSlope, 0.0,
1267 0.0, m_SimpleIntercept, 0.0,
1268 0.0, m_SimpleSlopeSigma, 0.0,
1269 0.0, m_SimpleInterceptSigma, 0.0,
1270 x2, y2, z2,
1271 sigmaX, sigmaY, sigmaZ);
1272 dx = z2 - z1;
1273 dy = sqrt(x2*x2 + y2*y2) - sqrt(x1*x1 + y1*y1);
1274 }
1275 else {
1276 geom->FindIntersection(m_Part, m_Seg, m_LastGap,
1277 m_SimpleSlope, 0.0, 1.0,
1278 m_SimpleIntercept, 0.0, 0.0,
1279 m_SimpleSlopeSigma, 0.0, 0.0,
1280 m_SimpleInterceptSigma, 0.0, 0.0,
1281 x1, y1, z1,
1282 sigmaX, sigmaY, sigmaZ);
1283 geom->FindIntersection(m_Part, m_Seg, gap,
1284 m_SimpleSlope, 0.0, 1.0,
1285 m_SimpleIntercept, 0.0, 0.0,
1286 m_SimpleSlopeSigma, 0.0, 0.0,
1287 m_SimpleInterceptSigma, 0.0, 0.0,
1288 x2, y2, z2,
1289 sigmaX, sigmaY, sigmaZ);
1290 dx = x2 - x1;
1291 dy = y2 - y1;
1292 }
1293 }
1294 else {
1295 if (m_Orient == 0) {
1296 geom->FindIntersection(m_Part, m_Seg, m_LastGap,
1297 0.0, m_SimpleSlope, 1.0,
1298 0.0, m_SimpleIntercept, 0.0,
1299 0.0, m_SimpleSlopeSigma, 0.0,
1300 0.0, m_SimpleInterceptSigma, 0.0,
1301 x1, y1, z1,
1302 sigmaX, sigmaY, sigmaZ);
1303 geom->FindIntersection(m_Part, m_Seg, gap,
1304 0.0, m_SimpleSlope, 1.0,
1305 0.0, m_SimpleIntercept, 0.0,
1306 0.0, m_SimpleSlopeSigma, 0.0,
1307 0.0, m_SimpleInterceptSigma, 0.0,
1308 x2, y2, z2,
1309 sigmaX, sigmaY, sigmaZ);
1310 dx = z2 - z1;
1311 dy = y2 - y1;
1312 }
1313 else {
1314 geom->FindIntersection(m_Part, m_Seg, m_LastGap,
1315 m_SimpleSlope, 0.0, 1.0,
1316 m_SimpleIntercept, 0.0, 0.0,
1317 m_SimpleSlopeSigma, 0.0, 0.0,
1318 m_SimpleInterceptSigma, 0.0, 0.0,
1319 x1, y1, z1,
1320 sigmaX, sigmaY, sigmaZ);
1321 geom->FindIntersection(m_Part, m_Seg, gap,
1322 m_SimpleSlope, 0.0, 1.0,
1323 m_SimpleIntercept, 0.0, 0.0,
1324 m_SimpleSlopeSigma, 0.0, 0.0,
1325 m_SimpleInterceptSigma, 0.0, 0.0,
1326 x2, y2, z2,
1327 sigmaX, sigmaY, sigmaZ);
1328 dx = z2 - z1;
1329 dy = x2 - x1;
1330 }
1331 }
1332
1333 dr = sqrt(dx*dx + dy*dy);
1334
1335 return WindowFunc(m_SimpleChi2,dr);
1336}
static MucGeoGeneral * Instance()
Get a pointer to the single instance of MucGeoGeneral.
void FindIntersection(const int part, const int seg, const int gap, const float vx, const float vy, const float vz, const float x0, const float y0, const float z0, const float sigmaVx, const float sigmaVy, const float sigmaVz, const float sigmaX0, const float sigmaY0, const float sigmaZ0, float &x, float &y, float &z, float &sigmaX, float &sigmaY, float &sigmaZ)
float WindowFunc(const float &chisq, const float &distance) const

◆ GetSeg()

int MucRec2DRoad::GetSeg ( ) const

In which segment was this road found?

Definition at line 920 of file MucRec2DRoad.cxx.

920{ return m_Seg; }

Referenced by MucRecRoadFinder::execute(), and MucRec3DRoad::MucRec3DRoad().

◆ GetSimpleFitParams() [1/2]

void MucRec2DRoad::GetSimpleFitParams ( float &  a,
float &  b,
float &  c,
int &  whichhalf,
float &  sigmaa,
float &  sigmab,
float &  sigmac,
float &  chisq,
int &  ndof 
) const

Get the parameters from the simple quad fit.

Definition at line 1108 of file MucRec2DRoad.cxx.

1111{
1112 a = m_SimpleQuad_a;
1113 b = m_SimpleQuad_b;
1114 c = m_SimpleQuad_c;
1115 whichhalf = m_SimpleQuad_whichhalf;
1116
1117 sigmaa = m_SimpleQuad_aSigma;
1118 sigmab = m_SimpleQuad_bSigma;
1119 sigmac = m_SimpleQuad_cSigma;
1120
1121 chisq = m_SimpleChi2;
1122 ndof = m_SimpleDOF;
1123
1124 return;
1125}

◆ GetSimpleFitParams() [2/2]

void MucRec2DRoad::GetSimpleFitParams ( float &  slope,
float &  intercept,
float &  sigmaSlope,
float &  sigmaIntercept,
float &  chisq,
int &  ndof 
) const

Get the parameters from the simple fit.

Definition at line 1083 of file MucRec2DRoad.cxx.

1086{
1087 slope = m_SimpleSlope;
1088 intercept = m_SimpleIntercept;
1089 sigmaSlope = m_SimpleSlopeSigma;
1090 sigmaIntercept = m_SimpleInterceptSigma;
1091 chisq = m_SimpleChi2;
1092 ndof = m_SimpleDOF;
1093
1094 return;
1095}

Referenced by Project(), MucRec3DRoad::Project(), MucRec3DRoad::ProjectNoCurrentGap(), and MucRec3DRoad::ProjectWithSigma().

◆ GetSlope()

float MucRec2DRoad::GetSlope ( ) const

Slope of trajectory.

Definition at line 1056 of file MucRec2DRoad.cxx.

1056{ return m_SimpleSlope; }

Referenced by MucRecRoadFinder::execute(), and RecMucTrack::LineFit().

◆ GetTotalHits()

int MucRec2DRoad::GetTotalHits ( ) const

How many hits in all does this road contain?

Definition at line 947 of file MucRec2DRoad.cxx.

947{ return m_TotalHits; }

Referenced by MucRecRoadFinder::execute(), MucRec3DRoad::GetTotalHits(), and MucRec3DRoad::GetTotalHitsDelta().

◆ GetVertexPos()

void MucRec2DRoad::GetVertexPos ( float &  x,
float &  y,
float &  z 
) const

Position of the vertex point.

Definition at line 928 of file MucRec2DRoad.cxx.

929{
930 x = m_VertexPos.x();
931 y = m_VertexPos.y();
932 z = m_VertexPos.z();
933
934 return;
935}

Referenced by MucRec3DRoad::MucRec3DRoad().

◆ HasHit()

bool MucRec2DRoad::HasHit ( MucRecHit hit) const

Does the hit exist in the road .

Definition at line 1395 of file MucRec2DRoad.cxx.

1396{
1397
1398 vector<MucRecHit*>::const_iterator iHit;
1399 bool HitExist = false;
1400
1401 // Also, if a track hits overlap region of two panels, we avoid
1402 // to double count hits in two panels
1403
1404 Identifier id = hit->GetID();
1405
1406 for ( iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++ ) {
1407 if ( *iHit ) { // Check for a null pointer.
1408 if ( (*iHit)->GetID() == id ) HitExist = true;
1409 }
1410 if (HitExist) break;
1411 }
1412
1413 return HitExist;
1414}

◆ HasHitInGap()

bool MucRec2DRoad::HasHitInGap ( const int &  gap) const

Does this road contain any hits in the given gap?

Definition at line 1003 of file MucRec2DRoad.cxx.

1004{
1005 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
1006 cout << "MucRec2DRoad::HasHitInGap-E2 invalid gap = " << gap << endl;
1007 return false;
1008 }
1009
1010 bool found = false;
1011 vector<MucRecHit*>::const_iterator iHit;
1012
1013 for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
1014 if (*iHit) { // Check for a null pointer.
1015 if ( (*iHit)->Gap() == gap ) {
1016 found = true;
1017 }
1018 }
1019 }
1020
1021 return found;
1022}

Referenced by MucRec3DRoad::HasHitInGap().

◆ operator=()

MucRec2DRoad & MucRec2DRoad::operator= ( const MucRec2DRoad orig)

Assignment constructor.

Definition at line 47 of file MucRec2DRoad.cxx.

48{
49 // Assignment operator.
50 if ( this != &orig ) { // Watch out for self-assignment!
51 m_VertexPos = orig.m_VertexPos;
52 m_VertexSigma = orig.m_VertexSigma;
53 m_Index = orig.m_Index;
54 m_Part = orig.m_Part;
55 m_Seg = orig.m_Seg;
56 m_Orient = orig.m_Orient;
57 m_Chi2 = orig.m_Chi2;
58 m_DOF = orig.m_DOF;
59 m_FitOK = orig.m_FitOK;
60 m_MaxHitsPerGap = orig.m_MaxHitsPerGap;
61 m_LastGap = orig.m_LastGap;
62 m_TotalHits = orig.m_TotalHits;
63 m_HitDistance = orig.m_HitDistance;
64 m_pHits = orig.m_pHits;
65 m_fittingMethod = orig.m_fittingMethod;
66 }
67
68 return *this;
69}

◆ PrintHitsInfo()

void MucRec2DRoad::PrintHitsInfo ( ) const

Print Hits Infomation.

Definition at line 1449 of file MucRec2DRoad.cxx.

1450{
1451 vector<MucRecHit*>::const_iterator iHit;
1452 for ( iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
1453 if (*iHit) { // Check for a null pointer.
1454 float xl, yl, zl;
1455 (*iHit)->GetStrip()->GetCenterPos(xl, yl, zl);
1456 HepPoint3D vl(xl, yl, zl);
1457 HepPoint3D vg = (*iHit)->GetGap()->TransformToGlobal(vl);
1458
1459 cout << " orient " << m_Orient
1460 << " part " << (*iHit)->Part()
1461 << " seg " << (*iHit)->Seg()
1462 << " gap " << (*iHit)->Gap()
1463 << " strip " << (*iHit)->Strip()
1464 << " pos (" << vg.x() << ", " << vg.y() << ", " << vg.z() << ")"
1465 << endl;
1466 }
1467 }
1468
1469}

Referenced by MucRec3DRoad::PrintHitsInfo().

◆ Project()

void MucRec2DRoad::Project ( const int &  gap,
float &  x,
float &  y,
float &  z,
float &  x2,
float &  y2,
float &  z2 
)

Where does the trajectory of this road intersect a specific gap?

Definition at line 734 of file MucRec2DRoad.cxx.

736{
737 float sigmaX, sigmaY, sigmaZ;
738
739 x = 0.0; sigmaX = 0.0; x2 = 0.0;
740 y = 0.0; sigmaY = 0.0; y2 = 0.0;
741 z = 0.0; sigmaZ = 0.0; z2 = 0.0;
742
743 if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
744 cout << "MucRec2DRoad::Project-E1 invalid gap = " << gap << endl;
745 return;
746 }
747
748 // Determine the projection point of the "simple" fit to the desired gap.
749 float x0, y0, z0, sigmaX0, sigmaY0, sigmaZ0;
750 float phi = m_Seg*0.25*kPi;
752
753 if (m_Part == 1) {
754 if (m_Orient == 0) {
755 geom->FindIntersection(m_Part, m_Seg, gap,
756 m_SimpleSlope*cos(phi), m_SimpleSlope*sin(phi), 1.0,
757 m_SimpleIntercept*cos(phi),m_SimpleIntercept*sin(phi), 0.0,
758 m_SimpleSlopeSigma, 0.0, 0.0,
759 m_SimpleInterceptSigma, 0.0, 0.0,
760 x0, y0, z0,
761 sigmaX0, sigmaY0, sigmaZ0);
762
763 if(m_SimpleSlope>10000){ //the line is in right center of this segment, we can not get intersection in y coordinate, in this case, m_SimpleIntercept is in z coordinate.
764 geom->FindIntersection(m_Part, m_Seg, gap,
765 m_SimpleSlope*cos(phi), m_SimpleSlope*sin(phi), 1.0,
766 0.0, 0.0, m_SimpleIntercept,
767 m_SimpleSlopeSigma, 0.0, 0.0,
768 m_SimpleInterceptSigma, 0.0, 0.0,
769 x0, y0, z0,
770 sigmaX0, sigmaY0, sigmaZ0);
771
772 }
773
774 }
775 else {
776 geom->FindIntersection(m_Part, m_Seg, gap,
777 1.0, m_SimpleSlope, 0.0,
778 0.0, m_SimpleIntercept, 0.0,
779 0.0, m_SimpleSlopeSigma, 0.0,
780 0.0, m_SimpleInterceptSigma, 0.0,
781 x0, y0, z0,
782 sigmaX0, sigmaY0, sigmaZ0);
783 //cout<<"in MucRec2DRoad line Project xyz0 = "<<x0<<" "<<y0<<" "<<z0<<endl;
784
785 }
786 }
787 else {
788 if (m_Orient == 0) {
789 geom->FindIntersection(m_Part, m_Seg, gap,
790 m_SimpleSlope, m_SimpleSlope, 1.0,
791 0.0, m_SimpleIntercept, 0.0,
792 0.0, m_SimpleSlopeSigma, 0.0,
793 0.0, m_SimpleInterceptSigma, 0.0,
794 x0, y0, z0,
795 sigmaX0, sigmaY0, sigmaZ0);
796 }
797 else {
798 geom->FindIntersection(m_Part, m_Seg, gap,
799 m_SimpleSlope, m_SimpleSlope, 1.0,
800 m_SimpleIntercept, 0.0, 0.0,
801 m_SimpleSlopeSigma, 0.0, 0.0,
802 m_SimpleInterceptSigma, 0.0, 0.0,
803 x0, y0, z0,
804 sigmaX0, sigmaY0, sigmaZ0);
805 }
806 //cout << " part " << m_Part
807 // << " seg " << m_Seg
808 // << " gap " << gap
809 // << " orient " << m_Orient
810 // << " slope = " << m_SimpleSlope
811 // << endl;
812 }
813
814 //cout << "In find intersection x0 = " << x0 << " y0 = " << y0 << " z0 = " << z0 << endl;
815
816 float a,b,sa,sb,chi2;
817 int ndof;
818
819 int status = Fit(x0, y0, z0, a, b, sa, sb, chi2, ndof);
820
821// m_FitOK = (status == 0) && (chi2<1000.0);
822// if (!fFitOK) {
823// cout << "MucRec2DRoad::Project-E2 fit fail status = "
824// << status << " npoints = " << npoints << " chi-sq = "
825// << chi2 << endl;
826// }
827
828// // Assign to fields of TMuiRoad object.
829// m_DOF = npoints - 2;
830// m_Chi2 = chi2;
831
832 if (m_Part == 1) { //change from 0 to 1 at 2006.11.30
833 if (m_Orient == 0) {
834 geom->FindIntersection(m_Part, m_Seg, gap,
835 a*cos(phi), a*sin(phi), 1.0,
836 //a, 0.0, 1.0,
837 b*cos(phi), b*sin(phi), 0.0,
838 sa, 0.0, 0.0,
839 sb, 0.0, 0.0,
840 x, y, z,
841 sigmaX, sigmaY, sigmaZ);
842
843 if(fabs(a)>10000){ ///
844 geom->FindIntersection(m_Part, m_Seg, gap,
845 a*cos(phi), a*sin(phi), 1.0,
846 0.0, 0.0, b,
847 //a, 0.0, 1.0,
848 //b, 0.0, 0.0,
849 sa, 0.0, 0.0,
850 sb, 0.0, 0.0,
851 x, y, z,
852 sigmaX, sigmaY, sigmaZ);
853
854 }
855 }
856 else {
857 geom->FindIntersection(m_Part, m_Seg, gap,
858 1.0, a, 0.0,
859 0.0, b, 0.0,
860 0.0, sa, 0.0,
861 0.0, sb, 0.0,
862 x, y, z,
863 sigmaX, sigmaY, sigmaZ);
864
865 if(m_fittingMethod == 2 && m_QuadFitOK ){
866 float a, b, c;
867 float sa, sb, sc, chi2x; int ydof; int whichhalf;
868
869 GetSimpleFitParams(a, b, c, whichhalf, sa, sb, sc, chi2x, ydof);
870 geom->FindIntersection(m_Part, m_Seg, gap,
871 10E30, 0.0, m_SimpleIntercept, 0.0, //vy = infinite
872 a, b, c, //y = a*x*x + b*x +c
873 whichhalf,
874 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
875 x, y, z, x2, y2, z2,
876 sigmaX, sigmaY, sigmaZ);
877
878
879 }
880
881 }
882 }
883 else {
884 if (m_Orient == 0) {
885 geom->FindIntersection(m_Part, m_Seg, gap,
886 a, a, 1.0,
887 0.0, b, 0.0,
888 0.0, sa, 0.0,
889 0.0, sb, 0.0,
890 x, y, z,
891 sigmaX, sigmaY, sigmaZ);
892 }
893 else {
894 geom->FindIntersection(m_Part, m_Seg, gap,
895 a, a, 1.0,
896 b, 0.0, 0.0,
897 sa, 0.0, 0.0,
898 sb, 0.0, 0.0,
899 x, y, z,
900 sigmaX, sigmaY, sigmaZ);
901 }
902 }
903
904 return;
905}
double sin(const BesAngle a)
Definition: BesAngle.h:210
double cos(const BesAngle a)
Definition: BesAngle.h:213
void Fit()
Definition: Eangle1D/Fit.cxx:3
const double kPi
Definition: MucConstant.h:6
void GetSimpleFitParams(float &slope, float &intercept, float &sigmaSlope, float &sigmaIntercept, float &chisq, int &ndof) const
Get the parameters from the simple fit.

Referenced by GetHitDistance().

◆ SetIndex()

void MucRec2DRoad::SetIndex ( const int &  index)

Set the index for this road.

Definition at line 92 of file MucRec2DRoad.cxx.

93{
94 if (index >= 0) m_Index = index;
95}

◆ SetMaxNSkippedGaps()

void MucRec2DRoad::SetMaxNSkippedGaps ( const int &  nGaps)

Max number of consecutive gaps allowed with no hits attached. This parameter affects the calculation of the last gap.

Definition at line 180 of file MucRec2DRoad.cxx.

181{
182 m_MaxNSkippedGaps = numGaps;
183 CountHits(); // recalculate the last gap and the hit counts.
184}

Referenced by MucRecRoadFinder::execute().

◆ SimpleFit()

int MucRec2DRoad::SimpleFit ( float &  slope,
float &  intercept,
float &  sigmaSlope,
float &  sigmaIntercept,
float &  chisq,
int &  ndof 
)

Calculate the best-fit straight line with "simple" weights.

Definition at line 188 of file MucRec2DRoad.cxx.

194{
195//hits in one track can not more than 100.
196 if(m_pHits.size()>100){
197 cout<<"MucRec2DRoad: too many hits in this track!"<<endl;
198 return -1;
199 }
200 // Assign to temporary arrays to be used in fit.
201 float px[100];
202 float py[100];
203 float pw[100];
204 int npoints = 0;
205
206 vector<MucRecHit*>::const_iterator iHit;
207
208 float weight[100];
209
210// for (int i = 0; i < m_pHits.size(); i++) {
211// cout<<"info: "<<m_pHits[i]->Seg()<<" "<<m_pHits[i]->Gap()<<" "<< m_pHits[i]->Strip()<<endl;
212// }
213 for (int i = 0; i < m_pHits.size(); i++) {
214 weight[i] = 1;
215
216 for(int j = 0; j < m_pHits.size(); j++){
217
218 if(j == i) continue;
219 if(m_pHits[i]->Part() == m_pHits[j]->Part() &&
220 m_pHits[i]->Seg() == m_pHits[j]->Seg() &&
221 m_pHits[i]->Gap() == m_pHits[j]->Gap() )
222 {
223 int deltaStrip = fabs(m_pHits[i]->Strip()- m_pHits[j]->Strip());
224
225 //cout<<i<<" "<<m_pHits[i]->Seg()<<" "<<m_pHits[i]->Gap()<<" "<<m_pHits[i]->Strip()<<" - "<<m_pHits[j]->Strip()<<endl;
226 if(deltaStrip == 0){
227 cout<<"deltaStrip == 0 ? check it"<<endl;
228 }
229 else{
230 weight[i] *= (deltaStrip+1) * (deltaStrip+1);
231 }
232
233 }
234 }
235 }
236
237 for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
238 if (*iHit) { // Check for a null pointer.
239
240 /*
241 float localx, localy, localz;
242 (*iHit)->GetStrip()->GetCenterPos(localx, localy, localz);
243 if ( m_Orient == 0) {
244 px[npoints] = localy;
245 py[npoints] = localz;
246 }
247 else {
248 px[npoints] = localx;
249 py[npoints] = localz;
250 }
251 npoints++;
252 }
253 }
254 */
255
256
257 Hep3Vector center = (*iHit)->GetCenterPos();
258 Hep3Vector sigma = (*iHit)->GetCenterSigma();
259 //Hep3Vector sigma(1.0, 1.0, 1.0);
260
261 if (m_Part == 1) {
262 if ( m_Orient == 0) {
263 px[npoints] = center.z();
264 py[npoints] = sqrt(center.x()*center.x()
265 + center.y()*center.y());
266 if(m_Seg==2) py[npoints] = center.y(); //deal with seg2 seperately! because there is a hole here. 2006.11.9
267
268 pw[npoints] = 4.0 * 1.0/(sigma.y()*sigma.y()) / weight[npoints];
269 }
270 else {
271 px[npoints] = center.x();
272 py[npoints] = center.y();
273 pw[npoints] = 4.0 * 1.0/(sigma.x()*sigma.x()) / weight[npoints];
274 }
275 }
276 else {
277 if ( m_Orient == 0) {
278 px[npoints] = center.z();
279 py[npoints] = center.y();
280 pw[npoints] = 4.0 * 1.0/(sigma.y()*sigma.y()) / weight[npoints];
281 }
282 else {
283 px[npoints] = center.z();
284 py[npoints] = center.x();
285 pw[npoints] = 4.0 * 1.0/(sigma.x()*sigma.x()) / weight[npoints];
286 }
287 }
288
289 //cout << " in MucRec2DRoad ID: " <<(*iHit)->Part()<<" "<<(*iHit)->Seg()<<" "<<(*iHit)->Gap()<<" "<< (*iHit)->Strip()<<" "<<px[npoints] << " " << py[npoints] << " " << pw[npoints] << endl;
290 npoints++;
291
292 if(npoints > 99)
293 { cout<<"MucRec2DRoad: to many hits in this track, ignore it!"<<endl;
294 return -1;
295 }
296 }
297}
298/*
299float px2[100];
300float py2[100];
301for (int i = 0; i < m_pHits.size(); i++) {
302 int hitsInGap = 1;
303 px2[i] = -999; py2[i] = -999;
304 for(int j = 0; j < m_pHits.size(); j++){
305
306 if(j == i) continue;
307 if(m_pHits[i]->Part() == m_pHits[j]->Part() &&
308 m_pHits[i]->Seg() == m_pHits[j]->Seg() &&
309 m_pHits[i]->Gap() == m_pHits[j]->Gap() )
310 {
311 hitsInGap++;
312 px2[i] = (px[i]*(hitsInGap-1) + px[j])/hitsInGap;
313 py2[i] = (py[i]*(hitsInGap-1) + py[j])/hitsInGap;
314 cout<<hitsInGap<<" "<<px[i]<<" "<<px[j]<<" "<<px2[i]<<endl;
315 }
316 }
317}
318
319for(int i = 0; i < m_pHits.size(); i++){
320 if(px2[i] != -999&&py2[i] != -999){
321 px[i] = px2[i]; py[i] = py2[i];
322 }
323}
324*/
325
326ndof = npoints - 2;
327if (ndof < 0) return -1;
328
329if (npoints == 1) {
330 if (m_Part == 1) {
331 if ( m_Orient == 0) {
332 px[npoints] = m_VertexPos.z();
333 py[npoints] = sqrt(m_VertexPos.x()*m_VertexPos.x()
334 + m_VertexPos.y()*m_VertexPos.y());
335 if(m_Seg==2) py[npoints] = m_VertexPos.y();
336 }
337 else {
338 px[npoints] = m_VertexPos.x();
339 py[npoints] = m_VertexPos.y();
340 }
341 }
342 else {
343 if ( m_Orient == 0) {
344 px[npoints] = m_VertexPos.z();
345 py[npoints] = m_VertexPos.y();
346 }
347 else {
348 px[npoints] = m_VertexPos.z();
349 py[npoints] = m_VertexPos.x();
350 }
351 }
352 pw[npoints] = 1.0;
353 npoints++;
354}
355else {
356 if (npoints == 0 ) {
357 return -1;
358 }
359}
360
361// Do the fits here.
362MucRecLineFit fit;
363 int status = fit.LineFit(px, py, pw, npoints,
364 &slope, &intercept, &chisq,
365 &sigmaSlope, &sigmaIntercept);
366
367
368 float tempslope, tempintercept,tempchisq, tempsigmaslope, sigmaPos;
369 int status4 = fit.LineFit(px, py, pw,m_Part,m_Seg,m_Orient, npoints,
370 &tempslope, &tempintercept, &tempchisq,
371 &tempsigmaslope, &sigmaPos);
372
373 MucRecQuadFit quadfit;
374 float quad_a, quad_b, quad_c, sigmaquad_a, sigmaquad_b, sigmaquad_c, chisq_quad;
375 int whichhalf, status2;
376
377 if(m_fittingMethod == 2){
378 status2 = quadfit.QuadFit(px, py, pw, npoints,
379 &quad_a, &quad_b, &quad_c, &whichhalf, &chisq_quad,
380 &sigmaquad_a, &sigmaquad_b, &sigmaquad_c);
381
382 }
383 //cout << " in MucRec2DRoad slope " << slope << " "<<intercept<<endl;
384
385 if (slope > 1.0e10 || slope < -1.0e10) {
386 if (m_Seg > 4) slope *= -1.0; // to avoid wrong direction
387 }
388
389 m_SimpleSlope = slope;
390 m_SimpleSlopeSigma = sigmaSlope;
391 m_SimpleIntercept = intercept;
392 m_SimpleInterceptSigma = sigmaIntercept;
393 m_SimplePosSigma = sigmaPos; //new 20071227
394 m_SimpleChi2 = chisq;
395 m_SimpleDOF = ndof;
396 m_SimpleFitOK = (status == 0) && (chisq < 1000.0);
397 m_QuadFitOK = (status2 == 1);
398
399 m_SimpleQuad_a = quad_a;
400 m_SimpleQuad_b = quad_b;
401 m_SimpleQuad_c = quad_c;
402 m_SimpleQuad_whichhalf = whichhalf;
403 m_SimpleQuad_aSigma = sigmaquad_a;
404 m_SimpleQuad_bSigma = sigmaquad_b;
405 m_SimpleQuad_cSigma = sigmaquad_c;
406
407 return status;
408}
*********Class see also m_nmax DOUBLE PRECISION m_MasPhot DOUBLE PRECISION m_phsu DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_r2 DOUBLE PRECISION m_WtMass INTEGER m_nmax INTEGER m_Nevgen INTEGER m_IsFSR INTEGER m_MarTot *COMMON c_KarFin $ !Output file $ !Event serial number $ !alpha QED at Thomson limit $ !minimum energy at CMS for remooval $ !infrared dimensionless $ !dummy photon IR regulator $ !crude photon multiplicity enhancement factor *EVENT $ !MC crude volume of PhhSpace *Sfactors $ !YFS formfactor IR part only $ !YFS formfactor non IR finite part $ !mass weight
Definition: KarFin.h:34
int QuadFit(float x[], float y[], float w[], int n, float *a, float *b, float *c, int *half, float *chisq, float *siga, float *sigb, float *sigc)

Referenced by AttachHit(), and Fit().

◆ SimpleFitNoCurrentGap()

int MucRec2DRoad::SimpleFitNoCurrentGap ( int  currentgap,
float &  slope,
float &  intercept,
float &  sigmaSlope,
float &  sigmaIntercept,
float &  chisq,
int &  ndof 
)

Calculate the best-fit straight line with "simple" weights. not use current gap!!!

Definition at line 414 of file MucRec2DRoad.cxx.

422{
423 // Assign to temporary arrays to be used in fit.
424 float px[100];
425 float py[100];
426 float pw[100];
427 int npoints = 0;
428
429 int notused = 0;
430
431 vector<MucRecHit*>::const_iterator iHit;
432
433 float pw2[100];
434 for(int i = 0; i< 100; i++) {
435 pw2[i] = 1;
436 //----if( m_pHits[i]->Gap()==currentgap ) pw2[i] = 9999;
437 }
438
439 for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
440 if (*iHit) { // Check for a null pointer.
441 if( (*iHit)->Gap()==currentgap ) continue;
442
443 Hep3Vector center = (*iHit)->GetCenterPos();
444 Hep3Vector sigma = (*iHit)->GetCenterSigma();
445 //Hep3Vector sigma(1.0, 1.0, 1.0);
446
447 if (m_Part == 1) {
448 if ( m_Orient == 0) {
449 px[npoints] = center.z();
450 py[npoints] = sqrt(center.x()*center.x()
451 + center.y()*center.y());
452 if(m_Seg==2) py[npoints] = center.y(); //deal with seg2 seperately! because there is a hole here. 2006.11.9
453
454 pw[npoints] = 1.0/(sigma.y()*sigma.y())/pw2[npoints]/pw2[npoints];
455 }
456 else {
457 px[npoints] = center.x();
458 py[npoints] = center.y();
459 pw[npoints] = 1.0/(sigma.x()*sigma.x())/pw2[npoints]/pw2[npoints];
460 }
461 }
462 else {
463 if ( m_Orient == 0) {
464 px[npoints] = center.z();
465 py[npoints] = center.y();
466 pw[npoints] = 1.0/(sigma.y()*sigma.y())/pw2[npoints]/pw2[npoints];
467 }
468 else {
469 px[npoints] = center.z();
470 py[npoints] = center.x();
471 pw[npoints] = 1.0/(sigma.x()*sigma.x())/pw2[npoints]/pw2[npoints];
472 }
473 }
474
475 //cout << " in MucRec2DRoad ID: " <<(*iHit)->Part()<<" "<<(*iHit)->Seg()<<" "<<(*iHit)->Gap()<<" "<< (*iHit)->Strip()<<" "<<px[npoints] << " " << py[npoints] << " " << pw[npoints] << endl;
476 //cout<<"process ngap: "<<currentgap<<" current: "<<(*iHit)->Gap()<<" x: "<<px[npoints]<<" y: "<<py[npoints]<<" weight: "<<pw[npoints]<<endl;
477 npoints++;
478
479 if(npoints > 99)
480 { cout<<"MucRec2DRoad: to many hits in this track, ignore it!"<<endl;
481 return -1;
482 }
483 }
484 }
485
486 ndof = npoints - 2;
487 if (ndof < 0) return -1;
488
489 if (npoints == 1) {
490 if (m_Part == 1) {
491 if ( m_Orient == 0) {
492 px[npoints] = m_VertexPos.z();
493 py[npoints] = sqrt(m_VertexPos.x()*m_VertexPos.x()
494 + m_VertexPos.y()*m_VertexPos.y());
495 if(m_Seg==2) py[npoints] = m_VertexPos.y();
496 }
497 else {
498 px[npoints] = m_VertexPos.x();
499 py[npoints] = m_VertexPos.y();
500 }
501 }
502 else {
503 if ( m_Orient == 0) {
504 px[npoints] = m_VertexPos.z();
505 py[npoints] = m_VertexPos.y();
506 }
507 else {
508 px[npoints] = m_VertexPos.z();
509 py[npoints] = m_VertexPos.x();
510 }
511 }
512 pw[npoints] = 1.0;
513 npoints++;
514 }
515else {
516 if (npoints == 0 ) {
517 return -1;
518 }
519
520}
521
522// Do the fits here.
523MucRecLineFit fit;
524int status = fit.LineFit(px, py, pw, npoints,
525 &slope, &intercept, &chisq,
526 &sigmaSlope, &sigmaIntercept);
527
528MucRecQuadFit quadfit;
529float quad_a, quad_b, quad_c, sigmaquad_a, sigmaquad_b, sigmaquad_c, chisq_quad;
530int whichhalf, status2;
531
532if(m_fittingMethod == 2){
533 status2 = quadfit.QuadFit(px, py, pw, npoints,
534 &quad_a, &quad_b, &quad_c, &whichhalf, &chisq_quad,
535 &sigmaquad_a, &sigmaquad_b, &sigmaquad_c);
536
537}
538//cout << " in MucRec2DRoad slope " << slope << " "<<intercept<<endl;
539
540if (slope > 1.0e10 || slope < -1.0e10) {
541 if (m_Seg > 4) slope *= -1.0; // to avoid wrong direction
542}
543
544////////2009-03-12
545 m_SimpleSlope = slope;
546 m_SimpleSlopeSigma = sigmaSlope;
547 m_SimpleIntercept = intercept;
548 m_SimpleInterceptSigma = sigmaIntercept;
549 //m_SimplePosSigma = sigmaPos; //new 20071227
550 m_SimpleChi2 = chisq;
551 m_SimpleDOF = ndof;
552 m_SimpleFitOK = (status == 0) && (chisq < 1000.0);
553 m_QuadFitOK = (status2 == 1);
554
555 m_SimpleQuad_a = quad_a;
556 m_SimpleQuad_b = quad_b;
557 m_SimpleQuad_c = quad_c;
558 m_SimpleQuad_whichhalf = whichhalf;
559 m_SimpleQuad_aSigma = sigmaquad_a;
560 m_SimpleQuad_bSigma = sigmaquad_b;
561 m_SimpleQuad_cSigma = sigmaquad_c;
562
563
564return status;
565}

Referenced by MucRec3DRoad::ProjectNoCurrentGap(), and MucRec3DRoad::RefitNoCurrentGap().

◆ WeightFunc()

float MucRec2DRoad::WeightFunc ( const float &  chisq,
const float &  distance 
) const

Definition at line 1473 of file MucRec2DRoad.cxx.

1474{
1475 return 1.0;
1476}

Referenced by Fit().

◆ WindowFunc()

float MucRec2DRoad::WindowFunc ( const float &  chisq,
const float &  distance 
) const

Definition at line 1480 of file MucRec2DRoad.cxx.

1481{
1482 return 1.0;
1483}

Referenced by GetSearchWindowSize().


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