13 #define N(n) index_->nodes_[(n)]
14 #define NC(n,m) index_->nodes_[(n)].childID_[(m)] // the children n->m
15 #define NV(m) index_->nodes_[id].v_[(m)] // the vertices of n
16 #define V(m) index_->vertices_[(m)] // the vertex vector m
18 #define SGN(x) ( (x)<0? -1: ( (x)>0? 1:0 ) ) // signum
49 if(s1 < 0.0L) a1 = (-1) * a1 ;
50 if(s2 < 0.0L) a2 = (-1) * a2 ;
51 if(s3 < 0.0L) a3 = (-1) * a3 ;
76 for (i = 0, k = 0; i < 4 ; i++)
77 for (j = i+1; j < 4; j++, k++) {
78 d[k] = (*v[i]) ^ (*v[j]);
80 for (l = 0, m = 0; l < 4; l++)
81 if(l != i && l != j)s[k][m++] = d[k] * (*v[l]);
90 for(i = 0; i < 6; i++)
91 if(s[i][0] * s[i][1] > 0.0)
93 d[i] : (-1 * d[i])), 0.0));
101 for(i = 0; i < 6; i++)
102 if(s[i][0] == 0.0 || s[i][1] == 0.0) {
104 d[i] : (-1 * d[i])), 0.0));
121 for (
size_t i =
constraints_.size() - 1; i > 0; i-- ) {
185 std::vector<size_t> cornerConstr1, cornerConstr2, removeConstr;
186 std::vector<SpatialVector> corner;
212 bool ruledout =
true;
225 bool vi1ok =
true, vi2ok =
true;
229 if(k == i || k == j)
continue;
230 if(vi1ok && vi1 *
constraints_[k].a_ <= 0.0) vi1ok =
false;
231 if(vi2ok && vi2 *
constraints_[k].a_ <= 0.0) vi2ok =
false;
232 if(!vi1ok && !vi2ok)
break;
235 corner.push_back(vi1);
236 cornerConstr1.push_back(i);
237 cornerConstr2.push_back(j);
241 corner.push_back(vi2);
242 cornerConstr1.push_back(i);
243 cornerConstr2.push_back(j);
249 if(ruledout) removeConstr.push_back(i);
263 i = cornerConstr1[0];
264 j = cornerConstr2[0];
265 size_t c1(0),c2(0),k1(0),k2(0);
269 for( k = 1; k < cornerConstr1.size(); k ++) {
270 if(cornerConstr1[k] == i) {
272 c1 = cornerConstr2[k];
275 if(cornerConstr2[k] == i) {
277 c1 = cornerConstr1[k];
280 if(cornerConstr1[k] == j) {
282 c2 = cornerConstr2[k];
285 if(cornerConstr2[k] == j) {
287 c2 = cornerConstr1[k];
298 size_t c,currentCorner;
299 if( ((vi1 - corner[0]) ^
constraints_[i].a_) * corner[0] > 0 ) {
318 while( currentCorner ) {
319 for (k = 0; k < cornerConstr1.size(); k++) {
320 if(k == currentCorner)
continue;
321 if(cornerConstr1[k] == c) {
322 if( (currentCorner = k) == 0)
break;
324 c = cornerConstr2[k];
327 if(cornerConstr2[k] == c) {
328 if( (currentCorner = k) == 0)
break;
330 c = cornerConstr1[k];
336 for ( i = 0; i < removeConstr.size(); i++)
345 for(i = 0; i <
corners_.size(); i++)
346 for(j = i+1; j <
corners_.size(); j++)
347 for(k = j+1; k <
corners_.size(); k++) {
393 bool redundancy =
true;
399 for(i = 0; i < clen; i++) {
400 for(j = 0; j < i; j++) {
450 if(redundancy)
break;
497 phi = (phi <= -1.0L +
gEpsilon ?
gPi : acos(phi)) ;
503 if ( phi > a1 + a2 )
return -1;
504 if ( a1 > phi + a2 )
return 1;
505 if ( a2 > phi + a1 )
return 2;
517 addlevel_ = idx->maxlevel_ - idx->buildlevel_;
524 for(
uint32 i = 1; i <= 8; i++){
536 int level, i, shifts;
541 IDHIGHBIT = IDHIGHBIT << 63;
542 IDHIGHBIT2 = IDHIGHBIT2 << 62;
545 for(i = 0; i <
IDSIZE; i+=2) {
546 if ( (htmid << i) & IDHIGHBIT )
break;
549 level = (IDSIZE-i) >> 1;
555 shifts = (
olevel - level) << 1;
556 lo = htmid << shifts;
557 hi = lo + ((
uint64) 1 << shifts) -1;
578 const struct SpatialIndex::QuadNode *indexNode = &
index_->nodes_[id];
618 childID = indexNode->childID_[0];
622 childID = indexNode->childID_[0];
testTrixel(childID);
623 childID = indexNode->childID_[1];
testTrixel(childID);
624 childID = indexNode->childID_[2];
testTrixel(childID);
625 childID = indexNode->childID_[3];
testTrixel(childID);
669 ids[0] = id0 =
id << 2;
690 if ((level-- <= 0) || ((P == 4) || (F >= 2) || (P == 3 && F == 1) || (P > 1 && PPrev == 3))){
695 for(
int i=0; i<4; i++){
719 const struct SpatialIndex::QuadNode *indexNode = &
index_->nodes_[id];
721 m = indexNode->v_[0];
722 v0 = &
index_->vertices_[m];
724 m = indexNode->v_[1];
725 v1 = &
index_->vertices_[m];
727 m = indexNode->v_[2];
728 v2 = &
index_->vertices_[m];
772 if(vsum == 1 || vsum == 2)
return pARTIAL;
954 edge[0].e = v0 ^ v1; edge[0].e1 = &v0; edge[0].e2 = &v1;
955 edge[1].e = v1 ^ v2; edge[1].e1 = &v1; edge[1].e2 = &v2;
956 edge[2].e = v2 ^ v0; edge[2].e1 = &v2; edge[2].e2 = &v0;
957 edge[0].l = acos(v0 * v1);
958 edge[1].l = acos(v1 * v2);
959 edge[2].l = acos(v2 * v0);
961 for(
size_t i = 0; i <
corners_.size(); i++) {
963 if(i <
corners_.size() - 1) j = i+1;
969 for (
size_t iedge = 0; iedge < 3; iedge++) {
976 for(
size_t k = 0; k < 2; k++) {
980 l1 = acos( *(edge[iedge].e1) * a1 );
981 l2 = acos( *(edge[iedge].e2) * a1 );
982 if( l1 - edge[iedge].l <=
gEpsilon &&
1005 if (
eSolve(v0, v1, i) )
return true;
1006 if (
eSolve(v1, v2, i) )
return true;
1007 if (
eSolve(v2, v0, i) )
return true;
1025 float64 u2 = (1 - mu) / (1 + mu);
1028 float64 b = gamma1 * ( u2 - 1 ) + gamma2 * ( u2 + 1 );
1031 float64 D = b * b - 4 * a * c;
1033 if( D < 0.0L )
return false;
1037 float64 q = -0.5L * ( b + (
SGN(b) * sqrt(D) ) );
1048 if (i == 0)
return false;
1049 if ( root1 >= 0.0L && root1 <= 1.0L )
return true;
1050 if ( i == 2 && ( (root1 >= 0.0L && root1 <= 1.0L ) ||
1051 (root2 >= 0.0L && root2 <= 1.0L ) ) )
return true;
1105 if (
eSolve(v0, v1, cIndex) )
return true;
1106 if (
eSolve(v1, v2, cIndex) )
return true;
1107 if (
eSolve(v2, v0, cIndex) )
return true;
1151 if( ( (( v0 ^ v1 ) * v) < 0 ) ||
1152 ( (( v1 ^ v2 ) * v) < 0 ) ||
1153 ( (( v2 ^ v0 ) * v) < 0 ) )
All constraints negative or zero.
The Spatial Index is a quad tree of spherical triangles.
const SpatialIndex * index_
int testConstraints(size_t i, size_t j)
SpatialMarkup testNode(uint64 id)
void intersect(const SpatialIndex *index, HtmRange *hr)
Intersect with index.
SpatialConstraint boundingCircle_
Triangle partially intersected.
void add(SpatialConstraint &)
Add a constraint.
RangeConvex()
Default Constructor.
Triangle is outside the queried area.
bool eSolve(const SpatialVector &v1, const SpatialVector &v2, size_t cIndex)
bool testVectorInside(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2, SpatialVector &v)
bool testBoundingCircle(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2)
The Constraint is really a cone on the sky-sphere.
bool testEdge0(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2)
std::vector< SpatialVector > corners_
bool testEdgeConstraint(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2, size_t cIndex)
size_t testOtherPosNone(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2)
The SpatialVector is a 3D vector usually living on the surface of the sphere.
bool testHole(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2)
void mergeRange(const Key lo, const Key hi)
bool testEdge(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2)
void testPartial(size_t level, uint64 id, const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2, int previousPartials)
All of the triangle is inside queried area.
std::vector< SpatialConstraint > constraints_
All constraints positive or zero.
SpatialMarkup testTrixel(uint64 nodeIndex)
void saveTrixel(uint64 htmid)
unsigned long long uint64
void simplify()
Simplify the convex, remove redundancies.
bool testConstraintInside(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2, size_t cIndex)
SpatialMarkup testTriangle(const SpatialVector &v0, const SpatialVector &v1, const SpatialVector &v2, int vsum)
int testVertex(const SpatialVector &v)
void normalize()
Normalize vector length to 1.