FairRoot/PandaRoot
Public Member Functions | Private Member Functions | Private Attributes | Friends | List of all members
GFDetPlane Class Reference

Detector plane genfit geometry class. More...

#include <GFDetPlane.h>

Inheritance diagram for GFDetPlane:

Public Member Functions

 GFDetPlane (GFAbsFinitePlane *finite=NULL)
 
 GFDetPlane (const TVector3 &o, const TVector3 &u, const TVector3 &v, GFAbsFinitePlane *finite=NULL)
 
 GFDetPlane (const TVector3 &o, const TVector3 &n, GFAbsFinitePlane *finite=NULL)
 
virtual ~GFDetPlane ()
 
 GFDetPlane (const GFDetPlane &)
 
GFDetPlaneoperator= (const GFDetPlane &)
 
TVector3 getO () const
 
TVector3 getU () const
 
TVector3 getV () const
 
void set (const TVector3 &o, const TVector3 &u, const TVector3 &v)
 
void setO (const TVector3 &o)
 
void setO (double, double, double)
 
void setU (const TVector3 &u)
 
void setU (double, double, double)
 
void setV (const TVector3 &v)
 
void setV (double, double, double)
 
void setUV (const TVector3 &u, const TVector3 &v)
 
void setON (const TVector3 &o, const TVector3 &n)
 
void setFinitePlane (GFAbsFinitePlane *finite)
 
TVector3 getNormal () const
 
void setNormal (TVector3 n)
 
void setNormal (double, double, double)
 
void setNormal (const double &theta, const double &phi)
 
TVector2 project (const TVector3 &x) const
 projecting a direction onto the plane: More...
 
TVector2 LabToPlane (const TVector3 &x) const
 transform from Lab system into plane More...
 
TVector3 toLab (const TVector2 &x) const
 transform from plane coordinates to lab system More...
 
TVector3 dist (const TVector3 &point) const
 
TVector2 straightLineToPlane (const TVector3 &point, const TVector3 &dir) const
 gives u,v coordinates of the intersection point of a straight line with plane More...
 
void Print () const
 
void getGraphics (double mesh, double length, TPolyMarker3D **pl, TPolyLine3D **plLine, TPolyLine3D **u, TPolyLine3D **v, TPolyLine3D **n=NULL)
 for poor attempts of making an event display. There is a lot of room for improvements. More...
 
double distance (TVector3 &) const
 
double distance (double, double, double) const
 
bool inActive (const TVector3 &point, const TVector3 &dir) const
 intersect in the active area? C.f. GFAbsFinitePlane More...
 
bool inActive (double u, double v) const
 inActive methods refer to finite plane. C.f. GFAbsFinitePlane More...
 
bool inActive (const TVector2 &v) const
 inActive methods refer to finite plane. C.f. GFAbsFinitePlane More...
 

Private Member Functions

void sane ()
 

Private Attributes

TVector3 fO
 
TVector3 fU
 
TVector3 fV
 
GFAbsFinitePlanefFinitePlane
 

Friends

bool operator== (const GFDetPlane &lhs, const GFDetPlane &rhs)
 
bool operator!= (const GFDetPlane &lhs, const GFDetPlane &rhs)
 returns NOT == More...
 

Detailed Description

Detector plane genfit geometry class.

A detector plane is the principle object to define coordinate systems for track fitting in genfit. Since a particle trajectory is a one-dimensional object (regardless of any specific parameterization) positions with repect to the track are always meassured in a plane.

Which plane is choosen depends on the type of detector. Fixed plane detectors have their detector plane defined by their mechanical setup. While wire chambers or time projection chambers might want to define a detector plane more flexibly.

This class parameterizes a plane in terms of an origin vector o and two plane-spanning directions u and v.

Definition at line 59 of file GFDetPlane.h.

Constructor & Destructor Documentation

GFDetPlane::GFDetPlane ( GFAbsFinitePlane finite = NULL)

Definition at line 37 of file GFDetPlane.cxx.

References fO, fU, fV, r, and sane().

38  :fFinitePlane(finite)
39 {
40  static TRandom3 r(0);
41  fO.SetXYZ(0.,0.,0.);
42  fU.SetXYZ(r.Uniform(),r.Uniform(),0.);
43  fV.SetXYZ(r.Uniform(),r.Uniform(),0.);
44  sane();
45 }
double r
Definition: RiemannTest.C:14
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
void sane()
Definition: GFDetPlane.cxx:217
TVector3 fU
Definition: GFDetPlane.h:153
GFDetPlane::GFDetPlane ( const TVector3 &  o,
const TVector3 &  u,
const TVector3 &  v,
GFAbsFinitePlane finite = NULL 
)
GFDetPlane::GFDetPlane ( const TVector3 &  o,
const TVector3 &  n,
GFAbsFinitePlane finite = NULL 
)

Definition at line 47 of file GFDetPlane.cxx.

References setNormal().

50  :fO(o),fFinitePlane(finite){
51  setNormal(n);
52 }
int n
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
TVector3 fO
Definition: GFDetPlane.h:151
void setNormal(TVector3 n)
Definition: GFDetPlane.cxx:156
GFDetPlane::~GFDetPlane ( )
virtual

Definition at line 54 of file GFDetPlane.cxx.

References fFinitePlane.

54  {
55  if(fFinitePlane!=NULL) delete fFinitePlane;
56 }
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
GFDetPlane::GFDetPlane ( const GFDetPlane rhs)

Definition at line 58 of file GFDetPlane.cxx.

References GFAbsFinitePlane::clone(), fFinitePlane, fO, fU, and fV.

58  : TObject(rhs) {
59  if(rhs.fFinitePlane != NULL) fFinitePlane = rhs.fFinitePlane->clone();
60  else fFinitePlane = NULL;
61  fO = rhs.fO;
62  fU = rhs.fU;
63  fV = rhs.fV;
64 }
virtual GFAbsFinitePlane * clone() const =0
Deep copy ctor for polymorphic class.
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
TVector3 fU
Definition: GFDetPlane.h:153

Member Function Documentation

TVector3 GFDetPlane::dist ( const TVector3 &  point) const

Definition at line 207 of file GFDetPlane.cxx.

References LabToPlane(), p, toLab(), and x.

Referenced by GFWirepointHitPolicy::detPlane(), GFWireHitPolicy::detPlane(), RKTrackRep::Extrap(), GeaneTrackRep::extrapolate(), and GeaneTrackRep::extrapolateToPoint().

208 {
209  TVector2 p=LabToPlane(x);
210  TVector3 xplane=toLab(p);
211  return xplane-x;
212 }
Double_t p
Definition: anasim.C:58
TVector3 toLab(const TVector2 &x) const
transform from plane coordinates to lab system
Definition: GFDetPlane.cxx:196
TVector2 LabToPlane(const TVector3 &x) const
transform from Lab system into plane
Definition: GFDetPlane.cxx:189
Double_t x
double GFDetPlane::distance ( TVector3 &  v) const

Definition at line 343 of file GFDetPlane.cxx.

References fO, fU, fV, s, and t.

Referenced by RKTrackRep::Extrap(), and RKTrackRep::RKutta().

343  {
344  double s = (v - fO)*fU;
345  double t = (v - fO)*fV;
346  TVector3 distanceVector = v - fO - (s*fU) - (t*fV);
347  return distanceVector.Mag();
348 }
TLorentzVector s
Definition: Pnd2DStar.C:50
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
__m128 v
Definition: P4_F32vec4.h:4
TVector3 fU
Definition: GFDetPlane.h:153
TTree * t
Definition: bump_analys.C:13
double GFDetPlane::distance ( double  x,
double  y,
double  z 
) const

Definition at line 349 of file GFDetPlane.cxx.

References fO, fU, fV, s, t, and v.

349  {
350  TVector3 v(x,y,z);
351  double s = (v - fO)*fU;
352  double t = (v - fO)*fV;
353  TVector3 distanceVector = v - fO - (s*fU) - (t*fV);
354  return distanceVector.Mag();
355 }
TLorentzVector s
Definition: Pnd2DStar.C:50
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
__m128 v
Definition: P4_F32vec4.h:4
TVector3 fU
Definition: GFDetPlane.h:153
Double_t z
Double_t x
TTree * t
Definition: bump_analys.C:13
Double_t y
void GFDetPlane::getGraphics ( double  mesh,
double  length,
TPolyMarker3D **  pl,
TPolyLine3D **  plLine,
TPolyLine3D **  u,
TPolyLine3D **  v,
TPolyLine3D **  n = NULL 
)

for poor attempts of making an event display. There is a lot of room for improvements.

Definition at line 283 of file GFDetPlane.cxx.

References fO, fU, fV, getNormal(), i, and vec.

283  {
284  *pl = new TPolyMarker3D(21*21,24);
285  (*pl)->SetMarkerSize(0.1);
286  (*pl)->SetMarkerColor(kBlue);
287  int nI=10;
288  int nJ=10;
289  *plLine = new TPolyLine3D(5);
290 
291  {
292  TVector3 linevec;
293  int i,j;
294  i=-1*nI;j=-1*nJ;
295  linevec=(fO+(mesh*i)*fU+(mesh*j)*fV);
296  (*plLine)->SetPoint(0,linevec.X(),linevec.Y(),linevec.Z());
297  i=-1*nI;j=1*nJ;
298  linevec=(fO+(mesh*i)*fU+(mesh*j)*fV);
299  (*plLine)->SetPoint(0,linevec.X(),linevec.Y(),linevec.Z());
300  i=1*nI;j=-1*nJ;
301  linevec=(fO+(mesh*i)*fU+(mesh*j)*fV);
302  (*plLine)->SetPoint(2,linevec.X(),linevec.Y(),linevec.Z());
303  i=1*nI;j=1*nJ;
304  linevec=(fO+(mesh*i)*fU+(mesh*j)*fV);
305  (*plLine)->SetPoint(1,linevec.X(),linevec.Y(),linevec.Z());
306  i=-1*nI;j=-1*nJ;
307  linevec=(fO+(mesh*i)*fU+(mesh*j)*fV);
308  (*plLine)->SetPoint(4,linevec.X(),linevec.Y(),linevec.Z());
309 
310  }
311  for (int i=-1*nI;i<=nI;++i){
312  for (int j=-1*nJ;j<=nJ;++j){
313  TVector3 vec(fO+(mesh*i)*fU+(mesh*j)*fV);
314  int id=(i+10)*21+j+10;
315  (*pl)->SetPoint(id,vec.X(),vec.Y(),vec.Z());
316  }
317  }
318 
319 
320  *u = new TPolyLine3D(2);
321  (*u)->SetPoint(0,fO.X(),fO.Y(),fO.Z());
322  (*u)->SetPoint(1,fO.X()+length*fU.X(),fO.Y()+length*fU.Y(),fO.Z()+length*fU.Z());
323  (*u)->SetLineWidth(2);
324  (*u)->SetLineColor(kGreen);
325 
326 
327  *v = new TPolyLine3D(2);
328  (*v)->SetPoint(0,fO.X(),fO.Y(),fO.Z());
329  (*v)->SetPoint(1,fO.X()+length*fV.X(),fO.Y()+length*fV.Y(),fO.Z()+length*fV.Z());
330  (*v)->SetLineWidth(2);
331  (*v)->SetLineColor(kRed);
332 
333  if(n!=NULL){
334  *n = new TPolyLine3D(2);
335  TVector3 _n=getNormal();
336  (*n)->SetPoint(0,fO.X(),fO.Y(),fO.Z());
337  (*n)->SetPoint(1,fO.X()+length*_n.X(),fO.Y()+length*_n.Y(),fO.Z()+length*_n.Z());
338  (*n)->SetLineWidth(2);
339  (*n)->SetLineColor(kBlue);
340  }
341 }
Int_t i
Definition: run_full.C:25
int n
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
__m128 v
Definition: P4_F32vec4.h:4
TVector3 fU
Definition: GFDetPlane.h:153
TVector3 getNormal() const
Definition: GFDetPlane.cxx:138
dble_vec_t vec[12]
Definition: ranlxd.cxx:380
TVector3 GFDetPlane::getNormal ( ) const

Definition at line 138 of file GFDetPlane.cxx.

References fU, and fV.

Referenced by RKTrackRep::extrapolate(), GenfitTrack2PndTrack(), getGraphics(), GeaneTrackRep::getMom(), RKTrackRep::getMom(), GeaneTrackRep::getPosMom(), RKTrackRep::getPosMom(), GeaneTrackRep::getPosMomCov(), Print(), RKTrackRep::RKutta(), sane(), and straightLineToPlane().

139 {
140  TVector3 result=fU.Cross(fV);
141  result.SetMag(1.);
142  return result;
143 }
TVector3 fV
Definition: GFDetPlane.h:154
TVector3 fU
Definition: GFDetPlane.h:153
TVector3 GFDetPlane::getO ( ) const
inline
TVector3 GFDetPlane::getU ( ) const
inline
TVector3 GFDetPlane::getV ( ) const
inline
bool GFDetPlane::inActive ( const TVector3 &  point,
const TVector3 &  dir 
) const
inline

intersect in the active area? C.f. GFAbsFinitePlane

Definition at line 132 of file GFDetPlane.h.

References straightLineToPlane().

Referenced by RKTrackRep::Extrap(), inActive(), and RKTrackRep::RKutta().

132  {
133  return this->inActive( this->straightLineToPlane(point,dir));
134  }
TClonesArray * point
Definition: anaLmdDigi.C:29
TVector2 straightLineToPlane(const TVector3 &point, const TVector3 &dir) const
gives u,v coordinates of the intersection point of a straight line with plane
Definition: GFDetPlane.cxx:357
bool inActive(const TVector3 &point, const TVector3 &dir) const
intersect in the active area? C.f. GFAbsFinitePlane
Definition: GFDetPlane.h:132
bool GFDetPlane::inActive ( double  u,
double  v 
) const
inline

inActive methods refer to finite plane. C.f. GFAbsFinitePlane

Definition at line 137 of file GFDetPlane.h.

References fFinitePlane, and GFAbsFinitePlane::inActive().

137  {
138  if(fFinitePlane==NULL) return true;
139  return fFinitePlane->inActive(u,v);
140  }
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
__m128 v
Definition: P4_F32vec4.h:4
virtual bool inActive(const double &u, const double &v) const =0
bool GFDetPlane::inActive ( const TVector2 &  v) const
inline

inActive methods refer to finite plane. C.f. GFAbsFinitePlane

Definition at line 143 of file GFDetPlane.h.

References inActive().

143  {
144  return inActive(v.X(),v.Y());
145  }
__m128 v
Definition: P4_F32vec4.h:4
bool inActive(const TVector3 &point, const TVector3 &dir) const
intersect in the active area? C.f. GFAbsFinitePlane
Definition: GFDetPlane.h:132
TVector2 GFDetPlane::LabToPlane ( const TVector3 &  x) const

transform from Lab system into plane

Definition at line 189 of file GFDetPlane.cxx.

References d, fO, and project().

Referenced by dist().

190 {
191  TVector3 d=x-fO;
192  return project(d);
193 }
TObjArray * d
TVector3 fO
Definition: GFDetPlane.h:151
TVector2 project(const TVector3 &x) const
projecting a direction onto the plane:
Definition: GFDetPlane.cxx:181
Double_t x
GFDetPlane & GFDetPlane::operator= ( const GFDetPlane rhs)

Definition at line 65 of file GFDetPlane.cxx.

References GFAbsFinitePlane::clone(), fFinitePlane, fO, fU, and fV.

65  {
66  assert(this!=&rhs);
67  if(fFinitePlane!=NULL) {
68  delete fFinitePlane;
69  }
70  if(rhs.fFinitePlane != NULL){
72  }
73  else{
74  fFinitePlane = NULL;
75  }
76  fO = rhs.fO;
77  fU = rhs.fU;
78  fV = rhs.fV;
79  return *this;
80 }
virtual GFAbsFinitePlane * clone() const =0
Deep copy ctor for polymorphic class.
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
TVector3 fU
Definition: GFDetPlane.h:153
void GFDetPlane::Print ( ) const

Definition at line 238 of file GFDetPlane.cxx.

References fFinitePlane, fO, fU, fV, getNormal(), and GFAbsFinitePlane::Print().

Referenced by GFBookkeeping::Print(), and GFAbsTrackRep::Print().

239 {
240  std::cout<<"GFDetPlane: "
241  <<"O("<<fO.X()<<","<<fO.Y()<<","<<fO.Z()<<") "
242  <<"u("<<fU.X()<<","<<fU.Y()<<","<<fU.Z()<<") "
243  <<"v("<<fV.X()<<","<<fV.Y()<<","<<fV.Z()<<") "
244  <<"n("<<getNormal().X()<<","<<getNormal().Y()<<","<<getNormal().Z()<<") "
245  <<std::endl;
246  std::cout << fFinitePlane << std::endl;
247  if(fFinitePlane!=NULL) fFinitePlane->Print();
248 
249 }
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
TVector3 fU
Definition: GFDetPlane.h:153
virtual void Print() const =0
TVector3 getNormal() const
Definition: GFDetPlane.cxx:138
TVector2 GFDetPlane::project ( const TVector3 &  x) const

projecting a direction onto the plane:

Definition at line 181 of file GFDetPlane.cxx.

References Double_t, fU, fV, and x.

Referenced by LabToPlane(), and straightLineToPlane().

182 {
183  Double_t xfU=fU*x;
184  Double_t xfV=fV*x;
185  return TVector2(xfU,xfV);
186 }
TVector3 fV
Definition: GFDetPlane.h:154
Double_t
TVector3 fU
Definition: GFDetPlane.h:153
Double_t x
void GFDetPlane::sane ( )
private

Definition at line 217 of file GFDetPlane.cxx.

References fU, fV, getNormal(), n, Pi, and v.

Referenced by GFDetPlane(), set(), setO(), setU(), setUV(), and setV().

217  {
218  assert(fU!=fV);
219  // ensure unit vectors
220  fU.SetMag(1.);
221  fV.SetMag(1.);
222  // ensure orthogonal system
223  // fV should be reached by
224  // rotating fU around _n in a counterclockwise direction
225 
226  TVector3 n=getNormal();
227  fV=n.Cross(fU);
228 
229  TVector3 v=fU;
230  v.Rotate(TMath::Pi()*0.5,n);
231  TVector3 null=v-fV;
232  assert(null.Mag()<1E-6);
233 
234 }
int n
TVector3 fV
Definition: GFDetPlane.h:154
__m128 v
Definition: P4_F32vec4.h:4
TVector3 fU
Definition: GFDetPlane.h:153
Double_t Pi
TVector3 getNormal() const
Definition: GFDetPlane.cxx:138
void GFDetPlane::set ( const TVector3 &  o,
const TVector3 &  u,
const TVector3 &  v 
)

Definition at line 83 of file GFDetPlane.cxx.

References fO, fU, fV, sane(), and v.

Referenced by GFAbsTrackRep::reset().

86 {
87  fO=o;fU=u;fV=v;
88  sane();
89 }
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
__m128 v
Definition: P4_F32vec4.h:4
void sane()
Definition: GFDetPlane.cxx:217
TVector3 fU
Definition: GFDetPlane.h:153
void GFDetPlane::setFinitePlane ( GFAbsFinitePlane finite)
inline

Optionally, set the finite plane definition. This is most important for avoiding fake intersection points in fitting of loopers. This should be implemented for silicon detectors most importantly.

Definition at line 96 of file GFDetPlane.h.

References fFinitePlane.

96 {fFinitePlane=finite;}
GFAbsFinitePlane * fFinitePlane
Definition: GFDetPlane.h:156
void GFDetPlane::setNormal ( TVector3  n)

Definition at line 156 of file GFDetPlane.cxx.

References fabs(), fU, and fV.

Referenced by GFSpacepointHitPolicy::detPlane(), GFDetPlane(), RKTrackRep::RKTrackRep(), setNormal(), setON(), and trackProximity().

156  {
157  n.SetMag(1.);
158  if( fabs(n.X()) > 0.1 ){
159  fU.SetXYZ(1./n.X()*(-1.*n.Y()-1.*n.Z()),1.,1.);
160  fU.SetMag(1.);
161  }
162  else {
163  if(fabs(n.Y()) > 0.1){
164  fU.SetXYZ(1.,1./n.Y()*(-1.*n.X()-1.*n.Z()),1.);
165  fU.SetMag(1.);
166  }
167  else{
168  fU.SetXYZ(1.,1.,1./n.Z()*(-1.*n.X()-1.*n.Y()));
169  fU.SetMag(1.);
170  }
171  }
172  fV=n.Cross(fU);
173 }
int n
TVector3 fV
Definition: GFDetPlane.h:154
TVector3 fU
Definition: GFDetPlane.h:153
friend F32vec4 fabs(const F32vec4 &a)
Definition: P4_F32vec4.h:47
void GFDetPlane::setNormal ( double  X,
double  Y,
double  Z 
)

Definition at line 151 of file GFDetPlane.cxx.

References setNormal().

151  {
152  TVector3 N(X,Y,Z);
153  setNormal(N);
154 }
double Y
Definition: anaLmdDigi.C:68
void setNormal(TVector3 n)
Definition: GFDetPlane.cxx:156
double X
Definition: anaLmdDigi.C:68
double Z
Definition: anaLmdDigi.C:68
void GFDetPlane::setNormal ( const double &  theta,
const double &  phi 
)

Definition at line 175 of file GFDetPlane.cxx.

References CAMath::Cos(), n, setNormal(), and CAMath::Sin().

175  {
177  setNormal(n);
178 }
static T Sin(const T &x)
Definition: PndCAMath.h:42
int n
static T Cos(const T &x)
Definition: PndCAMath.h:43
void setNormal(TVector3 n)
Definition: GFDetPlane.cxx:156
void GFDetPlane::setO ( const TVector3 &  o)

Definition at line 92 of file GFDetPlane.cxx.

References fO, and sane().

Referenced by GFSpacepointHitPolicy::detPlane(), RKTrackRep::extrapolateToLine(), RKTrackRep::RKTrackRep(), and trackProximity().

93 {
94  fO=o;
95  sane();
96 }
TVector3 fO
Definition: GFDetPlane.h:151
void sane()
Definition: GFDetPlane.cxx:217
void GFDetPlane::setO ( double  X,
double  Y,
double  Z 
)

Definition at line 98 of file GFDetPlane.cxx.

References fO, and sane().

99 {
100  fO.SetXYZ(X,Y,Z);
101  sane();
102 }
TVector3 fO
Definition: GFDetPlane.h:151
double Y
Definition: anaLmdDigi.C:68
void sane()
Definition: GFDetPlane.cxx:217
double X
Definition: anaLmdDigi.C:68
double Z
Definition: anaLmdDigi.C:68
void GFDetPlane::setON ( const TVector3 &  o,
const TVector3 &  n 
)

Definition at line 145 of file GFDetPlane.cxx.

References fO, and setNormal().

Referenced by RKTrackRep::extrapolateToPoint(), and RKTrackRep::stepalong().

145  {
146  fO = o;
147  setNormal(n);
148 }
int n
TVector3 fO
Definition: GFDetPlane.h:151
void setNormal(TVector3 n)
Definition: GFDetPlane.cxx:156
void GFDetPlane::setU ( const TVector3 &  u)

Definition at line 105 of file GFDetPlane.cxx.

References fU, and sane().

Referenced by RKTrackRep::extrapolateToLine().

106 {
107  fU=u;
108  sane();
109 }
void sane()
Definition: GFDetPlane.cxx:217
TVector3 fU
Definition: GFDetPlane.h:153
void GFDetPlane::setU ( double  X,
double  Y,
double  Z 
)

Definition at line 111 of file GFDetPlane.cxx.

References fU, and sane().

112 {
113  fU.SetXYZ(X,Y,Z);
114  sane();
115 }
double Y
Definition: anaLmdDigi.C:68
void sane()
Definition: GFDetPlane.cxx:217
TVector3 fU
Definition: GFDetPlane.h:153
double X
Definition: anaLmdDigi.C:68
double Z
Definition: anaLmdDigi.C:68
void GFDetPlane::setUV ( const TVector3 &  u,
const TVector3 &  v 
)

Definition at line 130 of file GFDetPlane.cxx.

References fU, fV, sane(), and v.

131 {
132  fU=u;
133  fV=v;
134  sane();
135 }
TVector3 fV
Definition: GFDetPlane.h:154
__m128 v
Definition: P4_F32vec4.h:4
void sane()
Definition: GFDetPlane.cxx:217
TVector3 fU
Definition: GFDetPlane.h:153
void GFDetPlane::setV ( const TVector3 &  v)

Definition at line 118 of file GFDetPlane.cxx.

References fV, sane(), and v.

Referenced by RKTrackRep::extrapolateToLine().

119 {
120  fV=v;
121  sane();
122 }
TVector3 fV
Definition: GFDetPlane.h:154
__m128 v
Definition: P4_F32vec4.h:4
void sane()
Definition: GFDetPlane.cxx:217
void GFDetPlane::setV ( double  X,
double  Y,
double  Z 
)

Definition at line 124 of file GFDetPlane.cxx.

References fV, and sane().

125 {
126  fV.SetXYZ(X,Y,Z);
127  sane();
128 }
double Y
Definition: anaLmdDigi.C:68
TVector3 fV
Definition: GFDetPlane.h:154
void sane()
Definition: GFDetPlane.cxx:217
double X
Definition: anaLmdDigi.C:68
double Z
Definition: anaLmdDigi.C:68
TVector2 GFDetPlane::straightLineToPlane ( const TVector3 &  point,
const TVector3 &  dir 
) const

gives u,v coordinates of the intersection point of a straight line with plane

Definition at line 357 of file GFDetPlane.cxx.

References fabs(), fO, getNormal(), point, project(), and t.

Referenced by inActive().

357  {
358  TVector3 dirNorm(dir);
359  dirNorm.SetMag(1.);
360  TVector3 normal = getNormal();
361  double dirTimesN = dirNorm*normal;
362  if(fabs(dirTimesN)<1.E-6){//straight line is parallel to plane, so return infinity
363  //doesnt parallel mean that they cross in infinity ;-)
364  return TVector2(1.E100,1.E100);
365  }
366  double t = 1/dirTimesN * ((fO-point)*normal);
367  return project(point - fO + t * dirNorm);
368 }
TVector3 fO
Definition: GFDetPlane.h:151
TClonesArray * point
Definition: anaLmdDigi.C:29
TVector2 project(const TVector3 &x) const
projecting a direction onto the plane:
Definition: GFDetPlane.cxx:181
friend F32vec4 fabs(const F32vec4 &a)
Definition: P4_F32vec4.h:47
TTree * t
Definition: bump_analys.C:13
TVector3 getNormal() const
Definition: GFDetPlane.cxx:138
TVector3 GFDetPlane::toLab ( const TVector2 &  x) const

transform from plane coordinates to lab system

Definition at line 196 of file GFDetPlane.cxx.

References d, fO, fU, and fV.

Referenced by dist().

197 {
198  TVector3 d(fO);
199  d+=x.X()*fU;
200  d+=x.Y()*fV;
201  return d;
202 }
TObjArray * d
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
TVector3 fU
Definition: GFDetPlane.h:153
Double_t x

Friends And Related Function Documentation

bool operator!= ( const GFDetPlane lhs,
const GFDetPlane rhs 
)
friend

returns NOT ==

Definition at line 278 of file GFDetPlane.cxx.

278  {
279  return !(lhs==rhs);
280 }
bool operator== ( const GFDetPlane lhs,
const GFDetPlane rhs 
)
friend

this operator is called very often in Kalman filtering. It checks equality of planes by comparing the 9 double values that define them.

Definition at line 259 of file GFDetPlane.cxx.

259  {
260  if(
261  fabs( (lhs.fO.X()-rhs.fO.X()) ) > DETPLANE_EPSILON ||
262  fabs( (lhs.fO.Y()-rhs.fO.Y()) ) > DETPLANE_EPSILON ||
263  fabs( (lhs.fO.Z()-rhs.fO.Z()) ) > DETPLANE_EPSILON
264  ) return false;
265  else if(
266  fabs( (lhs.fU.X()-rhs.fU.X()) ) > DETPLANE_EPSILON ||
267  fabs( (lhs.fU.Y()-rhs.fU.Y()) ) > DETPLANE_EPSILON ||
268  fabs( (lhs.fU.Z()-rhs.fU.Z()) ) > DETPLANE_EPSILON
269  ) return false;
270  else if(
271  fabs( (lhs.fV.X()-rhs.fV.X()) ) > DETPLANE_EPSILON ||
272  fabs( (lhs.fV.Y()-rhs.fV.Y()) ) > DETPLANE_EPSILON ||
273  fabs( (lhs.fV.Z()-rhs.fV.Z()) ) > DETPLANE_EPSILON
274  ) return false;
275  return true;
276 }
TVector3 fO
Definition: GFDetPlane.h:151
TVector3 fV
Definition: GFDetPlane.h:154
#define DETPLANE_EPSILON
Definition: GFDetPlane.cxx:258
TVector3 fU
Definition: GFDetPlane.h:153
friend F32vec4 fabs(const F32vec4 &a)
Definition: P4_F32vec4.h:47

Member Data Documentation

GFAbsFinitePlane* GFDetPlane::fFinitePlane
private

Definition at line 156 of file GFDetPlane.h.

Referenced by GFDetPlane(), inActive(), operator=(), Print(), setFinitePlane(), and ~GFDetPlane().

TVector3 GFDetPlane::fO
private
TVector3 GFDetPlane::fU
private
TVector3 GFDetPlane::fV
private

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