Warning, /include/Geant4/tools/lina/geom2 is written in an unsupported language. File is not indexed.
0001 // Copyright (C) 2010, Guy Barrand. All rights reserved.
0002 // See the file tools.license for terms.
0003
0004 #ifndef tools_geom2
0005 #define tools_geom2
0006
0007 #include <vector>
0008 #include <cstddef>
0009
0010 namespace tools {
0011
0012 template <class VEC2>
0013 inline double is_left(const VEC2& P0,const VEC2& P1,const VEC2& P2){
0014 return ( (P1.v0() - P0.v0()) * (P2.v1() - P0.v1())
0015 - (P2.v0() - P0.v0()) * (P1.v1() - P0.v1()) );
0016 }
0017
0018 template <class VEC2>
0019 inline bool is_inside(const VEC2& a_P,const std::vector<VEC2>& a_V) {
0020 // V[] = vertex points of a polygon V[n+1] with V[n]=V[0]
0021
0022 // From :
0023 // http://softsurfer.com/Archive/algorithm_0103/algorithm_0103.htm
0024 // Copyright 2001, softSurfer (www.softsurfer.com)
0025 // This code may be freely used and modified for any purpose
0026 // providing that this copyright notice is included with it.
0027 // SoftSurfer makes no warranty for this code, and cannot be held
0028 // liable for any real or imagined damage resulting from its use.
0029 // Users of this code must verify correctness for their application.
0030
0031 size_t n = a_V.size()-1;
0032
0033 int wn = 0; // the winding number counter
0034
0035 // loop through all edges of the polygon
0036 for (size_t i=0; i<n; i++) { // edge from V[i] to V[i+1]
0037 if (a_V[i].v1() <= a_P.v1()) { // start y <= P[1]
0038 if (a_V[i+1].v1() > a_P.v1()) // an upward crossing
0039 if (is_left( a_V[i], a_V[i+1], a_P) > 0) // P left of edge
0040 ++wn; // have a valid up intersect
0041 } else { // start y > P[1] (no test needed)
0042 if (a_V[i+1].v1() <= a_P.v1()) // a downward crossing
0043 if (is_left( a_V[i], a_V[i+1], a_P) < 0) // P right of edge
0044 --wn; // have a valid down intersect
0045 }
0046 }
0047
0048 return ((wn!=0)?true:false);
0049 }
0050
0051 // the same done with std::pair.
0052
0053 template <class T>
0054 inline double is_left(const std::pair<T,T>& P0,const std::pair<T,T>& P1,const std::pair<T,T>& P2){
0055 return ( (P1.first - P0.first) * (P2.second - P0.second)
0056 - (P2.first - P0.first) * (P1.second - P0.second) );
0057 }
0058
0059 template <class T>
0060 inline bool inside(const std::pair<T,T>& a_P,const std::vector< std::pair<T,T> >& a_V) {
0061 // V[] = vertex points of a polygon V[n+1] with V[n]=V[0]
0062
0063 // From :
0064 // http://softsurfer.com/Archive/algorithm_0103/algorithm_0103.htm
0065 // Copyright 2001, softSurfer (www.softsurfer.com)
0066 // This code may be freely used and modified for any purpose
0067 // providing that this copyright notice is included with it.
0068 // SoftSurfer makes no warranty for this code, and cannot be held
0069 // liable for any real or imagined damage resulting from its use.
0070 // Users of this code must verify correctness for their application.
0071
0072 size_t n = a_V.size()-1;
0073
0074 int wn = 0; // the winding number counter
0075
0076 // loop through all edges of the polygon
0077 for (size_t i=0; i<n; i++) { // edge from V[i] to V[i+1]
0078 if (a_V[i].second <= a_P.second) { // start y <= P[1]
0079 if (a_V[i+1].second > a_P.second) // an upward crossing
0080 if (is_left( a_V[i], a_V[i+1], a_P) > 0) // P left of edge
0081 ++wn; // have a valid up intersect
0082 } else { // start y > P[1] (no test needed)
0083 if (a_V[i+1].second <= a_P.second) // a downward crossing
0084 if (is_left( a_V[i], a_V[i+1], a_P) < 0) // P right of edge
0085 --wn; // have a valid down intersect
0086 }
0087 }
0088
0089 return ((wn!=0)?true:false);
0090 }
0091
0092 template <class VEC2>
0093 inline bool intersect(const VEC2& P1,const VEC2& Q1,const VEC2& P2,const VEC2& Q2,VEC2& a_out) {
0094 // (P1,Q1) first line, (P2,Q2) second line and a_out intersection.
0095 // return false if no intersection. (For exa, parallel lines).
0096
0097 // P1+(Q1-P1)*r = P2+(Q2-P2)*s
0098 // (Q1-P1)*r - (Q2-P2)*s = P2-P1
0099 // r*(Q1-P1).v0 - s*(Q2-P2).v0 = (P2-P1).v0 //x
0100 // r*(Q1-P1).v1 - s*(Q2-P2).v1 = (P2-P1).v1 //y
0101
0102 // a*r + b*s = c
0103 // d*r + e*s = f
0104
0105 // r = (c*e-f*b)/(a*e-d*b)
0106 // s = (a*f-d*c)/(a*e-d*b)
0107
0108 typedef typename VEC2::elem_t T;
0109 VEC2 A = Q1-P1;
0110 VEC2 B = P2-Q2;
0111 VEC2 C = P2-P1;
0112
0113 T a = A.v0();
0114 T b = B.v0();
0115 T c = C.v0();
0116 T d = A.v1();
0117 T e = B.v1();
0118 T f = C.v1();
0119
0120 T det = a*e-d*b;
0121 if(det==T()) return false;
0122
0123 T r = (c*e-f*b)/det;
0124
0125 a_out = P1+A*r;
0126 return true;
0127 }
0128
0129 }
0130
0131 #endif
