StraightLine.cxx

Go to the documentation of this file.

#include <iostream>
#include <iomanip>
#include <cmath>
#include <cfloat>
 
#include "CgemGeomSvc/StraightLine.h"
 
StraightLine::StraightLine(const HepPoint3D & pivot, const HepVector & a)
 :m_pivot(pivot),m_a(a)
{
}
 
StraightLine::StraightLine(double dr, double phi0, double dz, double tanl)
 :m_a(HepVector(4,0))
{
     m_a[0]=dr;
     m_a[1]=phi0;
     m_a[2]=dz;
     m_a[3]=tanl;
 
     HepPoint3D pivot(0,0,0);
     m_pivot=pivot;
}
 
StraightLine::StraightLine(HepPoint3D pos1, HepPoint3D pos2)
 :m_a(HepVector(4,0))
{
     // cout << "debug int StraightLine" << endl;
     // cout << "2points: " << setw(15) << pos1.x() << setw(15) << pos1.y() << setw(15) << pos1.z()
     //       << setw(15) << pos2.x() << setw(15) << pos2.y() << setw(15) << pos2.z() << endl;
     HepPoint3D pivot(0,0,0);
     m_pivot=pivot;
 
     if(pos1.y()<pos2.y()) {
      HepPoint3D pos_tmp=pos1;
      pos1=pos2;
      pos2=pos_tmp;
     }
     else if(pos1.y()==pos2.y()&&pos1.x()>pos2.x()) {
      HepPoint3D pos_tmp = pos1;
      pos1=pos2;
      pos2=pos_tmp;
     }
     //cout<<__FILE__<<": pos1, pos2 = "<<pos1<<", "<<pos2<<endl;
     HepPoint3D vec = pos2-pos1;
     double phi_vec = vec.phi();// (pi,2pi]
     //cout<<__FILE__<<": phi_vec="<<phi_vec<<endl;
     double phi0=phi_vec+M_PI/2.;// automatically in (-pi/2,pi/2)
     //cout<<__FILE__<<": phi0   ="<<phi0   <<endl;
     double dphi1=pos1.phi()-phi0;
     double dr = pos1.perp()*cos(dphi1);
 
     double s1(0.),s2(0.);
     if(pos1.y()!=pos2.y()) {
      s1=(pos1.y()-dr*sin(phi0))/sin(phi_vec);
      s2=(pos2.y()-dr*sin(phi0))/sin(phi_vec);
     }
     else {
      s1=(pos1.x()-dr*cos(phi0))/cos(phi_vec);
      s2=(pos2.x()-dr*cos(phi0))/cos(phi_vec);
     }
    
     double dz(0.),tanl(0.);
     if(s1>=s2) {
      cout<<"StraightLine: s1>=s2!!!"<<endl;
     }
     else {
      tanl=(pos2.z()-pos1.z())/(s2-s1);
      dz=-s1*tanl+pos1.z();
     }
     m_a[0]=dr;
     m_a[1]=phi0;
     m_a[2]=dz;
     m_a[3]=tanl;
     // cout << "line: " << setw(15) << dr << setw(15) << phi0 << setw(15) << dz << setw(15) << tanl << endl;
}
 
StraightLine::~StraightLine()
{}
 
HepPoint3D StraightLine::x(double s) const
{
     double x=m_a[0]*cos(m_a[1])+s*cos(m_a[1]-M_PI/2.);
     double y=m_a[0]*sin(m_a[1])+s*sin(m_a[1]-M_PI/2.);
     double z=m_a[2]+s*m_a[3];
     return HepPoint3D(x,y,z);
}
 
double StraightLine::sAtR(double R, int direction) const
{
     // if pivot at (0,0,0)
     double ss=999999.;
     if(fabs(m_a[0])<=R) {
      ss=sqrt(R*R-m_a[0]*m_a[0]);
      ss*=direction;
     }
     return ss;
}
 
HepPoint3D StraightLine::xAtR(double R, int direction) const
{
     HepPoint3D pos(999999,999999,999999);
     double ss=sAtR(R, direction);
     if(ss!=999999.) pos=x(ss);
     return pos;
}
 
double StraightLine::sAtPhi(double phi) const
{
     // if pivot at (0,0,0)
     double ss=999999.;
     //double       phi0_global  = m_a[1];//(-pi/2,pi/2]
     if(m_a[0]<0) {
      //phi0_global+ = M_PI;  //(pi/2,3/2*pi]
      phi+=M_PI;
     }
     double dphi=phi0()-phi;
     while(dphi> M_PI) dphi-=2*M_PI;
     while(dphi<-M_PI) dphi+=2*M_PI;
     if(fabs(dphi)<M_PI/2)
     {
      ss=m_a[0]*tan(dphi);
     }
     return ss;
}
 
HepPoint3D StraightLine::xAtPhi(double phi) const
{
     HepPoint3D pos(999999,999999,999999);
     double ss=sAtPhi(phi);
     if(ss!=999999.) pos=x(ss);
     return pos;
}
 
void StraightLine::set(const HepVector & a)
{
     m_a=a;
}
 
void StraightLine::set(const HepPoint3D & pivot)
{
     m_pivot = pivot;
     // m_a to be updated // FIXME
}
 
void StraightLine::set(const HepPoint3D & pivot, const HepVector & a)
{
     m_pivot = pivot;
     m_a=a;
}
 
void StraightLine::print()
{
     cout<<"line pivot: "<<m_pivot<<endl;
     cout<<"     par: "<<m_a<<endl;
}