MarlinUtil  1.12.1
B_Util.h
1 #ifndef __B_UTIL_H
2 #define __B_UTIL_H
3 
4 #include <math.h>
5 
6 //-----------------------------------------------------------------------
7 //-----------------------------------------------------------------------
8 // Common independent classes
9 //-----------------------------------------------------------------------
10 //-----------------------------------------------------------------------
11 
12 class Point2D {
13 //-----------------------------------------------------------------------
14  public:
15  double x; // Cartesian coordinate
16  double y; // Cartesian coordinate
17 
18  double rz;
19  double rz2;
20 
21  void calc() { rz2=x*x+y*y; rz=sqrt(rz2); }
22  void calc(double _x,double _y) { x=_x; y=_y; rz2=x*x+y*y; rz=sqrt(rz2); }
23  void calc(Point2D const &a) { x=a.x; y=a.y; rz2=a.rz2; rz=a.rz; }
24  double dist2(Point2D &a) { return (x-a.x)*(x-a.x)+(y-a.y)*(y-a.y); }
25  double dist(Point2D &a) { return sqrt(dist2(a)); }
26 
27  Point2D() : x(0.),y(0.),rz(0.),rz2(0.) { }
28  Point2D(double _x,double _y) { calc(_x,_y); }
29  Point2D(Point2D const & a) { calc(a); }
30  Point2D& operator=(const Point2D &a) { calc(a); return *this; }
31 };
32 
33 class Vector3D;
34 
35 //-----------------------------------------------------------------------
36 class Point3D : public Point2D {
37 //-----------------------------------------------------------------------
38  public:
39  double z; // Cartesian coordinate
40  double r; // Radius-vector length
41 
42  void calc() { Point2D::calc(); r=sqrt(rz2+z*z); }
43  void calc(double _x,double _y,double _z) { Point2D::calc(_x,_y); z=_z; r=sqrt(rz2+z*z); }
44  void calc(Point3D const &_a) { Point2D::calc(_a); z=_a.z; r=_a.r; }
45 
46  double dist2(Point3D &a) { return Point2D::dist2(a)+(z-a.z)*(z-a.z); }
47  double dist(Point3D &a) { return sqrt(dist2(a)); }
48 
49  // Point3D operator=(Point3D &a) { return a; }
50 
51  Point3D operator+(Point3D &a){
52  Point3D t(x+a.x,y+a.y,z+a.z);
53  return t;
54  }
55 
56  Point3D() : Point2D(),z(0.),r(0.) {}
57  Point3D(double _x,double _y,double _z) : Point2D(_x,_y) { z=_z; r=sqrt(rz2+z*z); }
58  Point3D(Point3D const &a) : Point2D(a) { z=a.z; r=a.r; }
59 
60  double rdist(Point3D &a) { return fabs(a.r-r); }
61 
62  bool IsBetween(Point3D &p0,Point3D &p1);
63  bool IsNear(float px,float py,float pz,float tolerance=5){
64  if(fabs(px-x)+fabs(py-y)+fabs(pz-z)<tolerance)
65  return true;
66  return false;
67  }
68  Point3D& operator=(const Point3D &_a) { calc(_a); return *this; }
69 };
70 
71 //-----------------------------------------------------------------------
72 class Vector3D : public Point3D {
73 //-----------------------------------------------------------------------
74  public:
75 
76  Vector3D(Point3D &a,Point3D &b) : Point3D(b.x-a.x,b.y-a.y,b.z-a.z) { }
77  Vector3D(Point3D const &a) : Point3D(a) { }
78  Vector3D() : Point3D() { }
79  Vector3D(double _x,double _y,double _z) : Point3D(_x,_y,_z) { }
80  double dot(Vector3D &a) { return x*a.x+y*a.y+z*a.z; }
81  double module() { return sqrt(x*x+y*y+z*z); }
82  void divide(double a) { x/=a; y/=a; z/=a; }
83  void mult(double a) { x*=a; y*=a; z*=a; }
84 
85  double length2() const { return x*x+y*y+z*z; }
86 
87  Point3D linear(double a,Point3D &b) {
88  Point3D p(a*x+b.x,a*y+b.y,a*z+b.z);
89  return p;
90  }
91 
92  Point3D add(const Point3D &_a){
93  Point3D p(x+_a.x,y+_a.y,z+_a.z);
94  return p;
95  }
96 
97  Point3D rotate(const Point3D &v1,const Point3D &v2,const Point3D &v3){
98  Point3D p(x*v1.x+y*v2.x+z*v3.x,
99  x*v1.y+y*v2.y+z*v3.y,
100  x*v1.z+y*v2.z+z*v3.z);
101  return p;
102  }
103 };
104 
105 Point3D Projection(Point3D &p,Point3D &p0,Point3D &p1);
106 
107 //-----------------------------------------------------------------------
108 class Sphere3D {
109 //-----------------------------------------------------------------------
110  public:
111  static const double RADDEG;// = 57.2957795130823209;
112 
113  double r; // Radius
114  double t; // Theta
115  double p; // Phi
116  Sphere3D() : r(0.),t(0.),p(0.) {}
117  void calc(Point3D &a) {
118 // Fill Spherical coordinate system
119  double rsq = hypot(a.x,a.y);
120  r = hypot(rsq,a.z);
121  t = atan2(rsq,a.z);
122  p = atan2(a.y,a.x);
123  if (p < 0.0)
124  p = 2*M_PI + p;
125  }
126 };
127 //-----------------------------------------------------------------------
129 //-----------------------------------------------------------------------
130 public:
131  double sum_x,sum_y,sum_x2,sum_y2,sum_xy;
132  unsigned n;
133 
134  double a;
135  double b;
136  double md; // mean displacement
137 
138  void init(){
139  sum_x=sum_y=sum_x2=sum_y2=sum_xy=0.;
140  a=b=md=0.;
141  n=0;
142  }
143  void add(double x,double y){
144  sum_x+=x;
145  sum_x2+=x*x;
146  sum_y+=y;
147  sum_xy+=x*y;
148  n++;
149  }
150  void calc(){
151  if(n>1){
152  a=(n*sum_xy-sum_x*sum_y)/(n*sum_x2-sum_x*sum_x);
153  b=(sum_y-a*sum_x)/n;
154  md=(sum_y2+a*a*sum_x2+n*b*b+2*a*b*sum_x-2*a*sum_xy-2*b*sum_y)/n;
155  }
156  }
157 };
158 
159 //-----------------------------------------------------------------------
160 template <class _Tp>
161 //-----------------------------------------------------------------------
162 class MinMax {
163 //-----------------------------------------------------------------------
164 public:
165  _Tp min;
166  _Tp max;
167  MinMax() : min(100000), max(-100000) { }
168  void init() { min=100000; max=-100000; }
169  void set(_Tp v){
170  if(v<min)
171  min=v;
172  if(v>max)
173  max=v;
174  }
175 };
176 
177 #define ARSIZE(x) (sizeof(x)/sizeof(x[0]))
178 
179 //-----------------------------------------------------------------------
180 //-----------------------------------------------------------------------
181 // Some useful functions
182 //-----------------------------------------------------------------------
183 //-----------------------------------------------------------------------
184 double angle_dist(double t1,double p1,double t2,double p2);
185 
186 #endif /* __B_UTIL_H */
Definition: B_Util.h:36
Definition: B_Util.h:108
Definition: B_Util.h:72
Definition: B_Util.h:162
Definition: B_Util.h:12
Definition: B_Util.h:128