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.
MucRec2DRoad &	operator= (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
MucRecHit *	GetHit (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< Identifier >	GetHitsID () 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.

  : m_VertexPos(xVertex, yVertex, zVertex),
    m_VertexSigma(0.0, 0.0, 0.0),
    m_Part(part), m_Seg(seg), m_Orient(orient),
    m_Chi2(0.0), m_DOF(0),
    m_MaxHitsPerGap(0), m_LastGap(0), m_TotalHits(0),
    m_FitOK(false), m_QuadFitOK(false),
    m_HitDistance(5, float(muckInfinity)),
    m_pHits(0), m_fittingMethod(fittingMethod)
{ }

MucRec2DRoad() [2/2]

MucRec2DRoad::MucRec2DRoad

(

const MucRec2DRoad &

source

)

Copy constructor.

Definition at line 72 of file MucRec2DRoad.cxx.

  : m_VertexPos(source.m_VertexPos),
    m_VertexSigma(source.m_VertexSigma),
    m_Index(source.m_Index),
    m_Part(source.m_Part), m_Seg(source.m_Seg), m_Orient(source.m_Orient),
    m_Chi2(source.m_Chi2), m_DOF(source.m_DOF),
    m_MaxHitsPerGap(source.m_MaxHitsPerGap),
    m_LastGap(source.m_LastGap), m_TotalHits(source.m_TotalHits),
    m_FitOK(source.m_FitOK), 
    m_HitDistance(source.m_HitDistance),
    m_pHits(source.m_pHits),
    m_fittingMethod(source.m_fittingMethod)
{ }

~MucRec2DRoad()

MucRec2DRoad::~MucRec2DRoad

(

)

Destructor.

Definition at line 87 of file MucRec2DRoad.cxx.

{ }

Member Function Documentation

AttachHit()

void MucRec2DRoad::AttachHit

(

MucRecHit *

hit

)

Attach the given hit to this road.

Definition at line 100 of file MucRec2DRoad.cxx.

{
  //  cout << "MucRec2DRoad::AttachHit-I0  hit = " << hit << endl;
  if (!hit) {
    cout << "MucRec2DRoad::AttachHit-E1  null hit pointer!" << endl;
    return ;
  }
  
  int gap = hit->Gap();
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    // The gap number of the hit is out of range.
    cout << "MucRec2DRoad::AttachHit-W3  bad gap number = " << gap
     << endl;
    return;
  }
 
  // Attach the hit to the road.
  //bool hitExist = false;
  for (int i = 0; i < (int)m_pHits.size(); i++) {
    if (m_pHits[i] == hit) return;
  }
  m_pHits.push_back(hit);
  //cout << "before Count " << m_pHits.size() << endl;
   //  m_HitDistance[gap] = dX;
 
  // Now recalculate the total number of hits and the max. number of
  // hits per gap.
  CountHits();
  //cout << "after Count " << m_pHits.size() << endl;
 
  // Redo the "simple" least-squares fit to the positions ...
  float a, sa, b, sb, chi;
  int n;
  int status = SimpleFit(a, b, sa, sb, chi, n);
  if (status < 0) {
    //cout << "MucRec2DRoad::AttachHit-E5  simple fit fail status = " << status << endl;
  }
 
  //cout << "Hit center = " << hit->GetCenterPos() << endl;
  //cout << "After attach hit, slope = " << a << "  intercept = " << b << endl;
}

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.

{
  //  cout << "MucRec2DRoad::AttachHit-I0  hit = " << hit << endl;
  if (!hit) {
    cout << "MucRec2DRoad::AttachHit-E1  null hit pointer!" << endl;
    return ;
  }
  
  int gap = hit->Gap();
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    // The gap number of the hit is out of range.
    cout << "MucRec2DRoad::AttachHit-W3  bad gap number = " << gap
     << endl;
    return;
  }
 
  // Attach the hit to the road.
  //bool hitExist = false;
  for (int i = 0; i < (int)m_pHits.size(); i++) {
    if (m_pHits[i] == hit) return;
  }
  m_pHits.push_back(hit);
  //cout << "before Count " << m_pHits.size() << endl;
   //  m_HitDistance[gap] = dX;
 
  // Now recalculate the total number of hits and the max. number of
  // hits per gap.
  CountHits();
  //cout << "after Count " << m_pHits.size() << endl;
 
}

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.

{
  int status;
 
  // If the "simple" fit has not been done yet, do it now.
  if (!m_SimpleFitOK) {
    // Least-squares fit to the positions ...
    float a, sa, b, sb, chi;
    int n;
    status = SimpleFit(a, b, sa, sb, chi, n);
    if (status < 0) {
      //cout << "MucRec2DRoad::Fit-E2  simple fit fail status = "
      //  << status << endl;
      return status;
    }
  }
 
  // Assign to temporary arrays to be used in fit.
  float px[100];
  float py[100];
  float pw[100];
  int npoints = 0;
  float dx, dy, dr;
 
 
  // include vertex point when do the fancy fitting
  //if (m_Orient == kHORIZ) {
  //  px[npoints] = m_VertexPos.y();
  //  dx = px[npoints] - y;
  //} else {
  //  px[npoints] = m_VertexPos.x();
  //  dx = px[npoints] - x;
  //}
  //pz[npoints] = m_VertexPos.z();
  //dz = pz[npoints] - z;
  //dr = sqrt(dx*dx + dz*dz);
  //pw[npoints] = WeightFunc(m_SimpleChi2,dr);  
  //npoints++;
 
  vector<MucRecHit*>::const_iterator iHit;
 
  // Determine the weights based on the chi-square of the fit
  // (the scatter of the points should be roughly related to the
  // momentum) and the distance from the center to the estimated
  // location.
 
  //cout << m_pHits.size() << endl;
  for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
 
      Hep3Vector center = (*iHit)->GetCenterPos();
      Hep3Vector sigma  = (*iHit)->GetCenterSigma();
 
      if (m_Part == 1) {    // change from 0 to 1 at 2006.11.30
        if ( m_Orient == 0) {
          px[npoints] = center.z();
          dx = px[npoints] - z;
          py[npoints] = sqrt(center.x()*center.x() 
              + center.y()*center.y());
          if(m_Seg==2) py[npoints] = center.y();  //deal with seg2 seperately! because there is a hole here. 2006.11.9
          dy = py[npoints] - sqrt(x*x + y*y);
        } 
        else {
          px[npoints] = center.x();
          dx = px[npoints] - x;
          py[npoints] = center.y();
          dy = py[npoints] - y;
        }
      }
      else {
        if ( m_Orient == 0) {
          px[npoints] = center.z();
          dx = px[npoints] - z;
          py[npoints] = center.y();
          dy = py[npoints] - y;
        }
        else {
          px[npoints] = center.z();
          dx = px[npoints] - z;
          py[npoints] = center.x();
          dy = py[npoints] - x;
        }
      }  
 
      dr = sqrt(dx*dx + dy*dy);
      pw[npoints] = WeightFunc(m_SimpleChi2, dr);
      //pw[npoints] = WeightFunc(m_SimpleChi2,dr);
 
      //    cout << "       " << npoints  << "     "
      // << px[npoints] << "  " << py[npoints] << "  " << pw[npoints]
      // << "         " << dx << "  " << dy << "  " << dr
      // << endl;
 
      npoints++;
 
      if(npoints > 99)
      {  cout<<"MucRec2DRoad: to many hits in this track, ignore it!"<<endl;
        return -1;
      }   
 
    }
  }
 
  // Refit ...
  ndof = npoints - 2;
  if (npoints == 1) {
    if (m_Part == 1) {     // change from 0 to 1 at 2006.11.30
      if ( m_Orient == 0) {
        px[npoints] = m_VertexPos.z();
        py[npoints] = sqrt(m_VertexPos.x()*m_VertexPos.x() 
            + m_VertexPos.y()*m_VertexPos.y());
        if(m_Seg==2) py[npoints] = m_VertexPos.y();
      } 
      else {
        px[npoints] = m_VertexPos.x();
        py[npoints] = m_VertexPos.y();
      }
    }
    else {
      if ( m_Orient == 0) {
        px[npoints] = m_VertexPos.z();
        py[npoints] = m_VertexPos.y();
      }
      else {
        px[npoints] = m_VertexPos.z();
        py[npoints] = m_VertexPos.x();
      }
    }  
    pw[npoints] = 1.0;
    npoints++;
  } 
  else {
    if (npoints == 0) {
      return -1;
    }
  }
 
  MucRecLineFit fit;
  //cout << "npoints " << npoints << endl;
  status = fit.LineFit(px, py, pw, npoints,
      &slope, &intercept, &chisq,
      &sigmaSlope, &sigmaIntercept);
  m_DOF  = ndof;
  m_Chi2 = chisq;
 
  if (status < 0) {
    //cout << "MucRec2DRoad::Fit-E3  fit fail status = " << status << endl;
  }
 
  return status;
}

GetDegreesOfFreedom()

int MucRec2DRoad::GetDegreesOfFreedom

(

)

const

How many degrees of freedom in the trajectory fit?

Definition at line 1064 of file MucRec2DRoad.cxx.

{ 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.

{
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    cout << "MucRec2DRoad::Hit-E1  invalid gap = " << gap << endl;
    return 0;
  }
  
  vector<MucRecHit*>::const_iterator iHit;
  
  for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
      if ( (*iHit)->Gap() == gap) { 
    return (*iHit);
      }
    }
  }
  
  return 0L;
}

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.

{
  // Use straight-line projection?
  if (!hit) {
    cout << "MucRec2DRoad::GetHitDistance-E1  null hit pointer!" << endl;
    return muckInfinity;
  }
 
  int gap = hit->Gap();
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    cout << "MucRec2DRoad::HitDistance-W2  bad gap number = " << gap << endl;
    return muckInfinity;
  }
 
  if ( hit->GetGap()->Orient() != m_Orient) {
    // The orientation of the hit is different from that of the road.
    cout << "MucRec2DRoad::GetHitDistance-W3 "
     << " hit has wrong orientation = "
     << hit->GetGap()->Orient()
     << endl;
    return muckInfinity;
  }
 
  HepPoint3D r  = hit->GetCenterPos();
  HepPoint3D rl = hit->GetGap()->TransformToGap(r);
  // cout << "hit center " << r << endl;
  // cout << "hit center local " << rl << endl;
 
  // Determine the projection point of the "simple" fit to the desired gap.
  // FIXME(?):  need to use fit with fancy weights instead?
  float delta, delta1,delta2;
  float x0, y0, z0;
  float x2, y2, z2;
  //float sigmaX0, sigmaY0, sigmaZ0;
 
  //cout << "y:x = " << m_SimpleSlope << endl;
 
  Project(gap, x0, y0, z0, x2, y2, z2);
  // cout << " gap = " << gap << " x0 = " << x0
  //      << " y0 = "  << y0  << " z0 = " << z0
  //      << endl;
 
  if(x0==y0&&x0==z0&&x0==0) {x0+=-9999;y0+=-9999;z0+=-9999;}//cout<<"wrong intersection"<<endl;}
  if(x2==y2&&x2==z2&&x2==0) {x2+=-9999;y2+=-9999;z2+=-9999;}//cout<<"wrong intersection"<<endl;}   //wrong intersection!!!
  HepPoint3D s  = HepPoint3D(x0, y0, z0);
  HepPoint3D s_2  = HepPoint3D(x2, y2, z2);
  HepPoint3D sl = hit->GetGap()->TransformToGap(s);
  HepPoint3D s_2l = hit->GetGap()->TransformToGap(s_2);
  //cout << "project center       " << s  << endl;
  
//   cout << "project center local sl= " << sl << endl;
//   cout << "project center local sl= " << s_2l << endl;
//   cout << "project center local rl= " << rl << endl;
  if ( m_Part == 0 ) {
    if ( m_Orient == 0 ) {
      delta1 = fabs( sl.y() - rl.y() );
      delta2= fabs( s_2l.y() - rl.y() );
    } 
    else {
      delta1 = fabs( sl.x() - rl.x() );
      delta2= fabs( s_2l.x() - rl.x() );
    }
  }
  else {
    if ( m_Orient == 0 ) {
      delta1 = fabs( sl.y() - rl.y() );
      delta2= fabs( s_2l.y() - rl.y() );
    } 
    else {
      delta1 = fabs( sl.x() - rl.x() );
      delta2= fabs( s_2l.x() - rl.x() );
    }
  }
//   if(which==0) {
//     if(delta1<delta2)delta = delta1;
//     else delta = delta2;
//   }
  //cout<<"which = "<<which<<"  distance = "<<delta1<<" "<<delta2<<endl;
  
  if(delta1 < delta2)  return delta1;
  else return delta2;
}

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

GetHits()

vector< MucRecHit * > MucRec2DRoad::GetHits

(

)

const

Definition at line 1439 of file MucRec2DRoad.cxx.

{
return m_pHits;
 
}

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.

{
  vector<Identifier> idCon;
  
  vector<MucRecHit*>::const_iterator iHit;
  for ( iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
      Identifier id = (*iHit)->GetID();
      idCon.push_back(id);
      /*
      cout << " MucRec2DRoad::HitIndices  gap  orientation  twopack= "
       <<  (*iHit)->ChannelID().Plane() << "   "
       <<  (*iHit)->ChannelID().Orient()  << "   "
       <<  (*iHit)->ChannelID().TwoPack() << endl; 
      */
    }
  }
  return idCon;
}

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.

{
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    cout << "MucRec2DRoad::GetHitsPerGap-E1  invalid gap = " << gap << endl;
    return 0;
  } 
 
  vector<MucRecHit*>::const_iterator iHit;
  int hitsInGap = 0;
 
  for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    
    if ( !(*iHit) ) {
      cout << "MucRec2DRoad::GetHitsPerGap()-E2 null hit pointer !"  << endl;
      return 0;
    }
    else {
      if( gap == (*iHit)->Gap() ) hitsInGap += 1;
    }
  }
 
  return hitsInGap;
}

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.

{
  return m_Index;
}

GetIntercept()

float MucRec2DRoad::GetIntercept

(

)

const

Intercept of trajectory.

Definition at line 1060 of file MucRec2DRoad.cxx.

{ 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.

{ 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.

{ 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.

{
  const int ngap = (int)MucID::getGapMax();
  int gap, count = 0;
  vector<int> firedGap;
  for ( gap = 0; gap < ngap; gap++) {
    firedGap.push_back(0);
  }
 
  vector<MucRecHit*>::const_iterator iHit;
  for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
      gap = (*iHit)->Gap();
      firedGap[gap] = 1;
    }
  }
  
  for ( gap = 0; gap < ngap; gap++) {
    count += firedGap[gap];
  }
 
  return count;
}

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.

{
  if (!road2) {
    return 0;
  }
 
  int count = 0;
  vector<MucRecHit*>::const_iterator iHit1;
  vector<MucRecHit*>::const_iterator iHit2;
  MucRecHit *hit1, *hit2;
 
  for( iHit1 = m_pHits.begin(); iHit1 != m_pHits.end(); iHit1++){
    for( iHit2 = road2->m_pHits.begin(); 
     iHit2 != road2->m_pHits.end(); iHit2++){
      hit1 = (*iHit1);
      hit2 = (*iHit2);
 
      if ( (hit1 != 0) && (hit2 != 0) ) {
         if (hit1->GetID() == hit2->GetID()) {
       count++;
     }
      }
    }
  }
  
  return count;
}

Referenced by MucRec3DRoad::GetNSharedHits().

GetOrient()

int MucRec2DRoad::GetOrient

(

)

const

In which orientation was this road found?

Definition at line 924 of file MucRec2DRoad.cxx.

{ 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.

{ 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.

{
  possigma = m_SimplePosSigma;
 
}

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

GetQuadFitOk()

bool MucRec2DRoad::GetQuadFitOk ( )

inline

Definition at line 160 of file MucRec2DRoad.h.

{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.

{
  if ( (!m_SimpleFitOK) || (m_SimpleDOF < 0) ) {
    return -1.0;
  }
  else if (m_SimpleDOF == 0) {
    return 0.0;
  }
  else {
    return (m_SimpleChi2 / m_SimpleDOF);
  }
}

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.

{
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    cout << "MucRec2DRoad::GetSearchWindowSize-E1  invalid gap = " << gap << endl;
    return 0.0;
  }
 
  // Determine the projection point of the "simple" fit to the last
  // gap and the desired gap.
  // FIXME?  the "delta-x" variable is calculated as the scalar
  // difference between the positions obtained by projecting to the
  // last gap and to the desired gap.
 
  MucGeoGeneral *geom = MucGeoGeneral::Instance();
  float x1, y1, z1, x2, y2, z2, dx, dy, dr;
  float sigmaX, sigmaY, sigmaZ;
 
  if (m_Part == 0) {
    if (m_Orient == 0) {
      geom->FindIntersection(m_Part, 0, m_LastGap,
                 1.0, m_SimpleSlope,          0.0,
                 0.0, m_SimpleIntercept,      0.0,
                 0.0, m_SimpleSlopeSigma,     0.0,
                 0.0, m_SimpleInterceptSigma, 0.0, 
                 x1, y1, z1,
                 sigmaX, sigmaY, sigmaZ);
      geom->FindIntersection(m_Part, 0, gap,
                 1.0, m_SimpleSlope,          0.0,
                 0.0, m_SimpleIntercept,      0.0,
                 0.0, m_SimpleSlopeSigma,     0.0,
                 0.0, m_SimpleInterceptSigma, 0.0, 
                 x2, y2, z2,
                 sigmaX, sigmaY, sigmaZ);
      dx = z2 - z1;
      dy = sqrt(x2*x2 + y2*y2) - sqrt(x1*x1 + y1*y1);
    }
    else {
      geom->FindIntersection(m_Part, m_Seg, m_LastGap,
                 m_SimpleSlope,          0.0, 1.0,
                 m_SimpleIntercept,      0.0, 0.0,
                 m_SimpleSlopeSigma,     0.0, 0.0,
                 m_SimpleInterceptSigma, 0.0, 0.0, 
                 x1, y1, z1,
                 sigmaX, sigmaY, sigmaZ);
      geom->FindIntersection(m_Part, m_Seg, gap,
                 m_SimpleSlope,          0.0, 1.0,
                 m_SimpleIntercept,      0.0, 0.0,
                 m_SimpleSlopeSigma,     0.0, 0.0,
                 m_SimpleInterceptSigma, 0.0, 0.0, 
                 x2, y2, z2,
                 sigmaX, sigmaY, sigmaZ);
      dx = x2 - x1;
      dy = y2 - y1;
    }
  }
  else {
    if (m_Orient == 0) {
      geom->FindIntersection(m_Part, m_Seg, m_LastGap,
                 0.0, m_SimpleSlope,          1.0,
                 0.0, m_SimpleIntercept,      0.0,
                 0.0, m_SimpleSlopeSigma,     0.0,
                 0.0, m_SimpleInterceptSigma, 0.0, 
                 x1, y1, z1,
                 sigmaX, sigmaY, sigmaZ);
      geom->FindIntersection(m_Part, m_Seg, gap,
                 0.0, m_SimpleSlope,          1.0,
                 0.0, m_SimpleIntercept,      0.0,
                 0.0, m_SimpleSlopeSigma,     0.0,
                 0.0, m_SimpleInterceptSigma, 0.0, 
                 x2, y2, z2,
                 sigmaX, sigmaY, sigmaZ);
      dx = z2 - z1;
      dy = y2 - y1;
    }
    else {
      geom->FindIntersection(m_Part, m_Seg, m_LastGap,
                 m_SimpleSlope,          0.0, 1.0,
                 m_SimpleIntercept,      0.0, 0.0,
                 m_SimpleSlopeSigma,     0.0, 0.0,
                 m_SimpleInterceptSigma, 0.0, 0.0,
                 x1, y1, z1,
                 sigmaX, sigmaY, sigmaZ);
      geom->FindIntersection(m_Part, m_Seg, gap,
                 m_SimpleSlope,          0.0, 1.0,
                 m_SimpleIntercept,      0.0, 0.0,
                 m_SimpleSlopeSigma,     0.0, 0.0,
                 m_SimpleInterceptSigma, 0.0, 0.0,
                 x2, y2, z2,
                 sigmaX, sigmaY, sigmaZ);
      dx = z2 - z1;
      dy = x2 - x1;
    }
  }
  
  dr = sqrt(dx*dx + dy*dy);
 
  return WindowFunc(m_SimpleChi2,dr);
}

GetSeg()

int MucRec2DRoad::GetSeg

(

)

const

In which segment was this road found?

Definition at line 920 of file MucRec2DRoad.cxx.

{ 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.

{
  a           = m_SimpleQuad_a;
  b           = m_SimpleQuad_b;
  c           = m_SimpleQuad_c;
  whichhalf   = m_SimpleQuad_whichhalf;
 
  sigmaa      = m_SimpleQuad_aSigma; 
  sigmab      = m_SimpleQuad_bSigma;
  sigmac      = m_SimpleQuad_cSigma;
 
  chisq           = m_SimpleChi2;
  ndof            = m_SimpleDOF;
 
  return;
}

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.

{
  slope           = m_SimpleSlope;
  intercept       = m_SimpleIntercept;
  sigmaSlope      = m_SimpleSlopeSigma;
  sigmaIntercept  = m_SimpleInterceptSigma;
  chisq           = m_SimpleChi2;
  ndof            = m_SimpleDOF;
 
  return;
}

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.

{ 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.

{ 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.

{
  x = m_VertexPos.x();
  y = m_VertexPos.y();
  z = m_VertexPos.z();
 
  return;
}

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.

{
 
  vector<MucRecHit*>::const_iterator iHit;
  bool HitExist = false;
 
  // Also, if a track hits overlap region of two panels, we avoid
  // to double count hits in two panels
  
  Identifier id = hit->GetID();
 
  for ( iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++ ) {
    if ( *iHit ) {  // Check for a null pointer.
      if ( (*iHit)->GetID() == id ) HitExist = true;    
    }
    if (HitExist) break;
  }
  
  return HitExist;
}

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.

{
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    cout << "MucRec2DRoad::HasHitInGap-E2  invalid gap = " << gap << endl;
    return false;
  }
 
  bool found = false;
  vector<MucRecHit*>::const_iterator iHit;
 
  for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
      if ( (*iHit)->Gap() == gap ) {
    found = true;
      }
    }
  }
 
  return found;
}

Referenced by MucRec3DRoad::HasHitInGap().

operator=()

MucRec2DRoad & MucRec2DRoad::operator=

(

const MucRec2DRoad &

orig

)

Assignment constructor.

Definition at line 47 of file MucRec2DRoad.cxx.

{
  // Assignment operator.
  if ( this != &orig ) {             // Watch out for self-assignment!
    m_VertexPos     = orig.m_VertexPos;
    m_VertexSigma   = orig.m_VertexSigma;
    m_Index         = orig.m_Index;
    m_Part          = orig.m_Part;
    m_Seg           = orig.m_Seg;
    m_Orient        = orig.m_Orient;
    m_Chi2          = orig.m_Chi2; 
    m_DOF           = orig.m_DOF;
    m_FitOK         = orig.m_FitOK;
    m_MaxHitsPerGap = orig.m_MaxHitsPerGap;
    m_LastGap       = orig.m_LastGap;
    m_TotalHits     = orig.m_TotalHits;
    m_HitDistance   = orig.m_HitDistance;
    m_pHits         = orig.m_pHits;
    m_fittingMethod = orig.m_fittingMethod;
  }
  
  return *this;
}

PrintHitsInfo()

void MucRec2DRoad::PrintHitsInfo

(

)

const

Print Hits Infomation.

Definition at line 1449 of file MucRec2DRoad.cxx.

{
  vector<MucRecHit*>::const_iterator iHit;
  for ( iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
      float xl, yl, zl;
      (*iHit)->GetStrip()->GetCenterPos(xl, yl, zl);
      HepPoint3D vl(xl, yl, zl);
      HepPoint3D vg = (*iHit)->GetGap()->TransformToGlobal(vl);
 
      cout << " orient " << m_Orient
       << " part "   << (*iHit)->Part() 
       << " seg  "   << (*iHit)->Seg()
       << " gap  "   << (*iHit)->Gap()
       << " strip "  << (*iHit)->Strip()
       << " pos ("   << vg.x() << ", " << vg.y() << ", " << vg.z() << ")"  
       << endl;
    }
  }
 
}

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.

{
  float sigmaX, sigmaY, sigmaZ;
 
  x = 0.0; sigmaX = 0.0; x2 = 0.0;
  y = 0.0; sigmaY = 0.0; y2 = 0.0;
  z = 0.0; sigmaZ = 0.0; z2 = 0.0;
 
  if ( (gap < 0) || (gap >= (int)MucID::getGapMax()) ) {
    cout << "MucRec2DRoad::Project-E1  invalid gap = " << gap << endl;
    return;
  }
  
  // Determine the projection point of the "simple" fit to the desired gap.
  float x0, y0, z0, sigmaX0, sigmaY0, sigmaZ0;
  float phi = m_Seg*0.25*kPi;
  MucGeoGeneral *geom = MucGeoGeneral::Instance();
 
  if (m_Part == 1) {
    if (m_Orient == 0) {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 m_SimpleSlope*cos(phi), m_SimpleSlope*sin(phi), 1.0,
                 m_SimpleIntercept*cos(phi),m_SimpleIntercept*sin(phi), 0.0,
                 m_SimpleSlopeSigma,     0.0, 0.0,
                 m_SimpleInterceptSigma, 0.0, 0.0, 
                 x0, y0, z0,
           sigmaX0, sigmaY0, sigmaZ0);
 
      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.                                               
        geom->FindIntersection(m_Part, m_Seg, gap,
        m_SimpleSlope*cos(phi), m_SimpleSlope*sin(phi), 1.0,
        0.0, 0.0, m_SimpleIntercept,
        m_SimpleSlopeSigma,     0.0, 0.0,
        m_SimpleInterceptSigma, 0.0, 0.0,
        x0, y0, z0,
        sigmaX0, sigmaY0, sigmaZ0);
 
      }
 
    }
    else {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 1.0, m_SimpleSlope,          0.0,
                 0.0, m_SimpleIntercept,      0.0,
                 0.0, m_SimpleSlopeSigma,     0.0,
                 0.0, m_SimpleInterceptSigma, 0.0, 
                 x0, y0, z0,
                 sigmaX0, sigmaY0, sigmaZ0);
      //cout<<"in MucRec2DRoad line Project xyz0 = "<<x0<<" "<<y0<<" "<<z0<<endl;
 
    }
  }
  else {
    if (m_Orient == 0) {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 m_SimpleSlope, m_SimpleSlope, 1.0,
                 0.0, m_SimpleIntercept,      0.0,
                 0.0, m_SimpleSlopeSigma,     0.0,
                 0.0, m_SimpleInterceptSigma, 0.0, 
                 x0, y0, z0,
                 sigmaX0, sigmaY0, sigmaZ0);
    }
    else {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 m_SimpleSlope, m_SimpleSlope, 1.0,
                 m_SimpleIntercept,      0.0, 0.0,
                 m_SimpleSlopeSigma,     0.0, 0.0,
                 m_SimpleInterceptSigma, 0.0, 0.0,
                 x0, y0, z0,
                 sigmaX0, sigmaY0, sigmaZ0);
    }
    //cout << " part " << m_Part 
    // << " seg "  << m_Seg 
    // << " gap "  << gap
    // << " orient "  << m_Orient
    // << " slope = " << m_SimpleSlope
    // << endl;
  }
    
  //cout << "In find intersection x0 = " << x0 << " y0 = " << y0 << " z0 = " << z0 << endl; 
 
  float a,b,sa,sb,chi2;
  int ndof;
 
  int status = Fit(x0, y0, z0, a, b, sa, sb, chi2, ndof);
 
//   m_FitOK = (status == 0) && (chi2<1000.0);
//   if (!fFitOK) {
//     cout << "MucRec2DRoad::Project-E2  fit fail status = "
//   << status << "  npoints = " << npoints << "  chi-sq = "
//   << chi2 << endl;
//   }
 
//   // Assign to fields of TMuiRoad object.
//   m_DOF  = npoints - 2;
//   m_Chi2 = chi2;
 
  if (m_Part == 1) {    //change from 0 to 1 at 2006.11.30
    if (m_Orient == 0) {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 a*cos(phi), a*sin(phi), 1.0,
                 //a,  0.0, 1.0,
                 b*cos(phi), b*sin(phi), 0.0,
                 sa, 0.0, 0.0,
                 sb, 0.0, 0.0, 
                 x,  y,   z,
                 sigmaX, sigmaY, sigmaZ);
 
      if(fabs(a)>10000){ /// 
        geom->FindIntersection(m_Part, m_Seg, gap,
            a*cos(phi), a*sin(phi), 1.0,
            0.0, 0.0, b,
            //a,  0.0, 1.0,
            //b,  0.0, 0.0,
            sa, 0.0, 0.0,
            sb, 0.0, 0.0,
            x,  y,   z,
            sigmaX, sigmaY, sigmaZ);
      
    }
    }
    else {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 1.0, a,  0.0,
                 0.0, b,  0.0,
                 0.0, sa, 0.0,
                 0.0, sb, 0.0, 
                 x,   y,  z,
                 sigmaX, sigmaY, sigmaZ);
 
      if(m_fittingMethod == 2 && m_QuadFitOK ){
    float a, b, c;
    float sa, sb, sc, chi2x; int ydof; int whichhalf;
    
    GetSimpleFitParams(a, b, c, whichhalf, sa, sb, sc, chi2x, ydof);
    geom->FindIntersection(m_Part, m_Seg, gap,
                   10E30, 0.0, m_SimpleIntercept,  0.0,   //vy = infinite
                   a, b, c,  //y = a*x*x + b*x +c
                   whichhalf,
                   0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
                   x, y, z, x2, y2, z2,
                   sigmaX, sigmaY, sigmaZ);
 
 
      }
 
    }
  }
  else {
    if (m_Orient == 0) {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 a,   a,  1.0,
                 0.0, b,  0.0,
                 0.0, sa, 0.0,
                 0.0, sb, 0.0, 
                 x,  y,   z,
                 sigmaX, sigmaY, sigmaZ);
    }
    else {
      geom->FindIntersection(m_Part, m_Seg, gap,
                 a,  a, 1.0,
                 b,  0.0, 0.0,
                 sa, 0.0, 0.0,
                 sb, 0.0, 0.0,
                 x,  y,   z,
                 sigmaX, sigmaY, sigmaZ);
    }
  }
 
  return;
}

Referenced by GetHitDistance().

SetIndex()

void MucRec2DRoad::SetIndex

(

const int &

index

)

Set the index for this road.

Definition at line 92 of file MucRec2DRoad.cxx.

{
  if (index >= 0) m_Index = index;
}

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.

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

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.

{
//hits in one track can not more than 100.
  if(m_pHits.size()>100){
    cout<<"MucRec2DRoad: too many hits in this track!"<<endl;
    return -1;
    }
  // Assign to temporary arrays to be used in fit.
  float px[100];
  float py[100];
  float pw[100];
  int npoints = 0;
 
  vector<MucRecHit*>::const_iterator iHit;
 
  float weight[100];
 
//    for (int i = 0; i < m_pHits.size(); i++) {
//  cout<<"info: "<<m_pHits[i]->Seg()<<" "<<m_pHits[i]->Gap()<<" "<< m_pHits[i]->Strip()<<endl;
//    }
  for (int i = 0; i < m_pHits.size(); i++) {
     weight[i] = 1;
      
     for(int j = 0; j < m_pHits.size(); j++){
 
    if(j == i) continue;
        if(m_pHits[i]->Part() == m_pHits[j]->Part() &&
       m_pHits[i]->Seg()  == m_pHits[j]->Seg() &&
       m_pHits[i]->Gap()  == m_pHits[j]->Gap() )
     {
       int deltaStrip = fabs(m_pHits[i]->Strip()- m_pHits[j]->Strip());
 
       //cout<<i<<"  "<<m_pHits[i]->Seg()<<"  "<<m_pHits[i]->Gap()<<"  "<<m_pHits[i]->Strip()<<" - "<<m_pHits[j]->Strip()<<endl;
       if(deltaStrip == 0){
             cout<<"deltaStrip == 0 ? check it"<<endl;
       } 
       else{
         weight[i] *=  (deltaStrip+1) * (deltaStrip+1);
       }
 
     }       
     }
  }
  
  for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
      
      /*
    float localx, localy, localz;
      (*iHit)->GetStrip()->GetCenterPos(localx, localy, localz);
      if ( m_Orient == 0) {
    px[npoints] = localy;
    py[npoints] = localz;
      }
      else {
    px[npoints] = localx;
    py[npoints] = localz;
      }
      npoints++;
      }
      }
      */
      
      
      Hep3Vector center = (*iHit)->GetCenterPos();
      Hep3Vector sigma  = (*iHit)->GetCenterSigma();
      //Hep3Vector sigma(1.0, 1.0, 1.0);
 
      if (m_Part == 1) {
        if ( m_Orient == 0) {
          px[npoints] = center.z();
          py[npoints] = sqrt(center.x()*center.x() 
              + center.y()*center.y());
          if(m_Seg==2) py[npoints] = center.y();  //deal with seg2 seperately! because there is a hole here. 2006.11.9
 
          pw[npoints] = 4.0 * 1.0/(sigma.y()*sigma.y()) / weight[npoints];
        } 
        else {
          px[npoints] = center.x();
          py[npoints] = center.y();
          pw[npoints] = 4.0 * 1.0/(sigma.x()*sigma.x()) / weight[npoints];
        }
      }
      else {
        if ( m_Orient == 0) {
          px[npoints] = center.z();
          py[npoints] = center.y();
          pw[npoints] = 4.0 * 1.0/(sigma.y()*sigma.y()) / weight[npoints];
        }
        else {
          px[npoints] = center.z();
          py[npoints] = center.x();
          pw[npoints] = 4.0 * 1.0/(sigma.x()*sigma.x()) / weight[npoints];
        }
      }  
 
      //cout << " in MucRec2DRoad  ID: " <<(*iHit)->Part()<<" "<<(*iHit)->Seg()<<" "<<(*iHit)->Gap()<<" "<< (*iHit)->Strip()<<" "<<px[npoints] << "  " << py[npoints] << "  " << pw[npoints] << endl;
      npoints++;
 
      if(npoints > 99) 
      {  cout<<"MucRec2DRoad: to many hits in this track, ignore it!"<<endl;
        return -1;
      }   
  }
}
/*
float px2[100];
float py2[100];
for (int i = 0; i < m_pHits.size(); i++) {
    int hitsInGap = 1;
    px2[i] = -999; py2[i] = -999;
    for(int j = 0; j < m_pHits.size(); j++){
 
        if(j == i) continue;
        if(m_pHits[i]->Part() == m_pHits[j]->Part() &&
                m_pHits[i]->Seg()  == m_pHits[j]->Seg() &&
                m_pHits[i]->Gap()  == m_pHits[j]->Gap() )
        {
            hitsInGap++;
            px2[i] = (px[i]*(hitsInGap-1) + px[j])/hitsInGap;
            py2[i] = (py[i]*(hitsInGap-1) + py[j])/hitsInGap;
            cout<<hitsInGap<<"  "<<px[i]<<"  "<<px[j]<<"  "<<px2[i]<<endl;
        }
    }
}
 
for(int i = 0; i < m_pHits.size(); i++){
    if(px2[i] != -999&&py2[i] != -999){
            px[i] = px2[i]; py[i] = py2[i];
    }
}
*/
 
ndof = npoints - 2;
if (ndof < 0) return -1;
 
if (npoints == 1) {
  if (m_Part == 1) {
    if ( m_Orient == 0) {
      px[npoints] = m_VertexPos.z();
      py[npoints] = sqrt(m_VertexPos.x()*m_VertexPos.x() 
          + m_VertexPos.y()*m_VertexPos.y());
      if(m_Seg==2) py[npoints] = m_VertexPos.y();
    } 
    else {
      px[npoints] = m_VertexPos.x();
      py[npoints] = m_VertexPos.y();
    }
  }
  else {
    if ( m_Orient == 0) {
      px[npoints] = m_VertexPos.z();
      py[npoints] = m_VertexPos.y();
    }
    else {
      px[npoints] = m_VertexPos.z();
      py[npoints] = m_VertexPos.x();
    }
  }  
  pw[npoints] = 1.0;
  npoints++;
} 
else {
  if (npoints == 0 ) {
    return -1;
  }
}
 
// Do the fits here.
MucRecLineFit fit;
  int status = fit.LineFit(px, py, pw, npoints,
               &slope, &intercept, &chisq,
               &sigmaSlope, &sigmaIntercept);
 
 
  float tempslope, tempintercept,tempchisq, tempsigmaslope, sigmaPos;
  int status4 = fit.LineFit(px, py, pw,m_Part,m_Seg,m_Orient, npoints,
      &tempslope, &tempintercept, &tempchisq,
      &tempsigmaslope, &sigmaPos);
  
  MucRecQuadFit quadfit;
  float quad_a, quad_b, quad_c, sigmaquad_a, sigmaquad_b, sigmaquad_c, chisq_quad;
  int whichhalf, status2;
 
  if(m_fittingMethod == 2){
  status2 = quadfit.QuadFit(px, py, pw, npoints,
                &quad_a, &quad_b, &quad_c, &whichhalf, &chisq_quad,
                &sigmaquad_a, &sigmaquad_b, &sigmaquad_c);
 
  }
  //cout << " in MucRec2DRoad slope " << slope << " "<<intercept<<endl;
  
  if (slope > 1.0e10 || slope < -1.0e10) {
    if (m_Seg > 4) slope *= -1.0;  // to avoid wrong direction
  }
 
  m_SimpleSlope          = slope;
  m_SimpleSlopeSigma     = sigmaSlope;
  m_SimpleIntercept      = intercept;
  m_SimpleInterceptSigma = sigmaIntercept;
  m_SimplePosSigma       = sigmaPos;   //new 20071227
  m_SimpleChi2           = chisq;
  m_SimpleDOF            = ndof;
  m_SimpleFitOK          = (status == 0) && (chisq < 1000.0);
  m_QuadFitOK            = (status2 == 1);
 
  m_SimpleQuad_a         = quad_a;
  m_SimpleQuad_b         = quad_b;
  m_SimpleQuad_c         = quad_c;
  m_SimpleQuad_whichhalf = whichhalf;
  m_SimpleQuad_aSigma    = sigmaquad_a;
  m_SimpleQuad_bSigma    = sigmaquad_b;
  m_SimpleQuad_cSigma    = sigmaquad_c;
 
  return status;
}

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.

{
  // Assign to temporary arrays to be used in fit.
  float px[100];
  float py[100];
  float pw[100];
  int npoints = 0;
 
  int notused = 0;
 
  vector<MucRecHit*>::const_iterator iHit;
 
  float pw2[100];
  for(int i = 0; i< 100; i++) {
    pw2[i] = 1;
    //----if( m_pHits[i]->Gap()==currentgap ) pw2[i] = 9999;
  }
 
  for (iHit = m_pHits.begin(); iHit != m_pHits.end(); iHit++) {
    if (*iHit) {  // Check for a null pointer.
      if( (*iHit)->Gap()==currentgap ) continue;
 
      Hep3Vector center = (*iHit)->GetCenterPos();
      Hep3Vector sigma  = (*iHit)->GetCenterSigma();
      //Hep3Vector sigma(1.0, 1.0, 1.0);
 
      if (m_Part == 1) {
        if ( m_Orient == 0) {
          px[npoints] = center.z();
          py[npoints] = sqrt(center.x()*center.x()
              + center.y()*center.y());
          if(m_Seg==2) py[npoints] = center.y();  //deal with seg2 seperately! because there is a hole here. 2006.11.9
 
          pw[npoints] = 1.0/(sigma.y()*sigma.y())/pw2[npoints]/pw2[npoints];
        }
        else {
          px[npoints] = center.x();
          py[npoints] = center.y();
          pw[npoints] = 1.0/(sigma.x()*sigma.x())/pw2[npoints]/pw2[npoints];
        }
      }
      else {
        if ( m_Orient == 0) {
          px[npoints] = center.z();
          py[npoints] = center.y();
          pw[npoints] = 1.0/(sigma.y()*sigma.y())/pw2[npoints]/pw2[npoints];
        }
        else {
          px[npoints] = center.z();
          py[npoints] = center.x();
          pw[npoints] = 1.0/(sigma.x()*sigma.x())/pw2[npoints]/pw2[npoints];
        }
      }
 
      //cout << " in MucRec2DRoad  ID: " <<(*iHit)->Part()<<" "<<(*iHit)->Seg()<<" "<<(*iHit)->Gap()<<" "<< (*iHit)->Strip()<<" "<<px[npoints] << "  " << py[npoints] << "  " << pw[npoints] << endl;
      //cout<<"process ngap: "<<currentgap<<" current: "<<(*iHit)->Gap()<<"  x: "<<px[npoints]<<"  y: "<<py[npoints]<<"  weight: "<<pw[npoints]<<endl;
      npoints++;
 
      if(npoints > 99)
      {  cout<<"MucRec2DRoad: to many hits in this track, ignore it!"<<endl;
        return -1;
      }
    }
  }
 
  ndof = npoints - 2;
  if (ndof < 0) return -1;
 
  if (npoints == 1) {
    if (m_Part == 1) {
      if ( m_Orient == 0) {
        px[npoints] = m_VertexPos.z();
        py[npoints] = sqrt(m_VertexPos.x()*m_VertexPos.x()
            + m_VertexPos.y()*m_VertexPos.y());
        if(m_Seg==2) py[npoints] = m_VertexPos.y();
      }
      else {
        px[npoints] = m_VertexPos.x();
        py[npoints] = m_VertexPos.y();
      }
    }
    else {
      if ( m_Orient == 0) {
        px[npoints] = m_VertexPos.z();
        py[npoints] = m_VertexPos.y();
      }
      else {
        px[npoints] = m_VertexPos.z();
        py[npoints] = m_VertexPos.x();
      }
    }
    pw[npoints] = 1.0;
    npoints++;
  }
else {
  if (npoints == 0 ) {
    return -1;
  }
 
}
 
// Do the fits here.
MucRecLineFit fit;
int status = fit.LineFit(px, py, pw, npoints,
    &slope, &intercept, &chisq,
    &sigmaSlope, &sigmaIntercept);
 
MucRecQuadFit quadfit;
float quad_a, quad_b, quad_c, sigmaquad_a, sigmaquad_b, sigmaquad_c, chisq_quad;
int whichhalf, status2;
 
if(m_fittingMethod == 2){
  status2 = quadfit.QuadFit(px, py, pw, npoints,
      &quad_a, &quad_b, &quad_c, &whichhalf, &chisq_quad,
      &sigmaquad_a, &sigmaquad_b, &sigmaquad_c);
 
}
//cout << " in MucRec2DRoad slope " << slope << " "<<intercept<<endl;
 
if (slope > 1.0e10 || slope < -1.0e10) {
  if (m_Seg > 4) slope *= -1.0;  // to avoid wrong direction
}
 
////////2009-03-12  
  m_SimpleSlope          = slope;
  m_SimpleSlopeSigma     = sigmaSlope;
  m_SimpleIntercept      = intercept;
  m_SimpleInterceptSigma = sigmaIntercept;
  //m_SimplePosSigma       = sigmaPos;   //new 20071227
  m_SimpleChi2           = chisq;
  m_SimpleDOF            = ndof;
  m_SimpleFitOK          = (status == 0) && (chisq < 1000.0);
  m_QuadFitOK            = (status2 == 1);
 
  m_SimpleQuad_a         = quad_a;
  m_SimpleQuad_b         = quad_b;
  m_SimpleQuad_c         = quad_c;
  m_SimpleQuad_whichhalf = whichhalf;
  m_SimpleQuad_aSigma    = sigmaquad_a;
  m_SimpleQuad_bSigma    = sigmaquad_b;
  m_SimpleQuad_cSigma    = sigmaquad_c;
 
 
return status;
}

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.

{
  return 1.0;
}

Referenced by Fit().

WindowFunc()

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

Definition at line 1480 of file MucRec2DRoad.cxx.

{
  return 1.0;
}

Referenced by GetSearchWindowSize().

---

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

• MucRec2DRoad.h
• MucRec2DRoad.cxx