GeodesicGrid.cpp   StelGeodesicGrid.cpp 
/* /*
GeodesicGrid: a library for dividing the sphere into triangle zones StelGeodesicGrid: a library for dividing the sphere into triangle zones
by subdividing the icosahedron by subdividing the icosahedron
Author and Copyright: Johannes Gajdosik, 2006 Author and Copyright: Johannes Gajdosik, 2006
This library requires a simple Vector library, This library requires a simple Vector library,
which may have different copyright and license, which may have different copyright and license,
for example Vec3d from vecmath.h. for example Vec3d from VecMath.hpp.
In the moment I choose to distribute the library under the GPL, In the moment I choose to distribute the library under the GPL,
later I may choose to additionally distribute it under a more later I may choose to additionally distribute it under a more
relaxed license like the LGPL. If you want to have the library relaxed license like the LGPL. If you want to have the library
under another license, please ask me. under another license, please ask me.
This library is free software; you can redistribute it and/or This library is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License modify it under the terms of the GNU General Public License
as published by the Free Software Foundation; either version 2 as published by the Free Software Foundation; either version 2
of the License, or (at your option) any later version. of the License, or (at your option) any later version.
skipping to change at line 33 skipping to change at line 33
but WITHOUT ANY WARRANTY; without even the implied warranty of but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details. GNU General Public License for more details.
You should have received a copy of the GNU General Public License You should have received a copy of the GNU General Public License
along with this library; if not, write to the Free Software along with this library; if not, write to the Free Software
Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
*/ */
#include "GeodesicGrid.hpp" #include "StelGeodesicGrid.hpp"
#include <QDebug> #include <QDebug>
#include <cmath> #include <cmath>
#include <cstdlib> #include <cstdlib>
static const double icosahedron_G = 0.5*(1.0+sqrt(5.0)); static const double icosahedron_G = 0.5*(1.0+sqrt(5.0));
static const double icosahedron_b = 1.0/sqrt(1.0+icosahedron_G*icosahedron_ G); static const double icosahedron_b = 1.0/sqrt(1.0+icosahedron_G*icosahedron_ G);
static const double icosahedron_a = icosahedron_b*icosahedron_G; static const double icosahedron_a = icosahedron_b*icosahedron_G;
static const Vec3d icosahedron_corners[12] = static const Vec3d icosahedron_corners[12] =
skipping to change at line 88 skipping to change at line 88
{{10,11, 4}}, // 10 {{10,11, 4}}, // 10
{{ 2, 4,11}}, // 19 {{ 2, 4,11}}, // 19
{{ 5, 4, 2}}, // 5 {{ 5, 4, 2}}, // 5
{{ 2, 8, 5}}, // 15 {{ 2, 8, 5}}, // 15
{{ 4, 1,10}}, // 17 {{ 4, 1,10}}, // 17
{{ 4, 5, 1}}, // 4 {{ 4, 5, 1}}, // 4
{{ 5, 9, 1}}, // 13 {{ 5, 9, 1}}, // 13
{{ 8, 9, 5}} // 8 {{ 8, 9, 5}} // 8
}; };
GeodesicGrid::GeodesicGrid(const int lev) : maxLevel(lev<0?0:lev), lastMaxS earchlevel(-1) StelGeodesicGrid::StelGeodesicGrid(const int lev) : maxLevel(lev<0?0:lev), lastMaxSearchlevel(-1)
{ {
if (maxLevel > 0) if (maxLevel > 0)
{ {
triangles = new Triangle*[maxLevel+1]; triangles = new Triangle*[maxLevel+1];
int nr_of_triangles = 20; int nr_of_triangles = 20;
for (int i=0;i<maxLevel;i++) for (int i=0;i<maxLevel;i++)
{ {
triangles[i] = new Triangle[nr_of_triangles]; triangles[i] = new Triangle[nr_of_triangles];
nr_of_triangles *= 4; nr_of_triangles *= 4;
} }
skipping to change at line 115 skipping to change at line 115
icosahedron_corners[corners[2]]); icosahedron_corners[corners[2]]);
} }
} }
else else
{ {
triangles = 0; triangles = 0;
} }
cacheSearchResult = new GeodesicSearchResult(*this); cacheSearchResult = new GeodesicSearchResult(*this);
} }
GeodesicGrid::~GeodesicGrid(void) StelGeodesicGrid::~StelGeodesicGrid(void)
{ {
if (maxLevel > 0) if (maxLevel > 0)
{ {
for (int i=maxLevel-1;i>=0;i--) delete[] triangles[i]; for (int i=maxLevel-1;i>=0;i--) delete[] triangles[i];
delete[] triangles; delete[] triangles;
} }
delete cacheSearchResult; delete cacheSearchResult;
cacheSearchResult = NULL; cacheSearchResult = NULL;
} }
void GeodesicGrid::getTriangleCorners(int lev,int index, void StelGeodesicGrid::getTriangleCorners(int lev,int index,
Vec3d &h0, Vec3d &h0,
Vec3d &h1, Vec3d &h1,
Vec3d &h2) const Vec3d &h2) const
{ {
if (lev <= 0) if (lev <= 0)
{ {
const int *const corners = icosahedron_triangles[index].corn ers; const int *const corners = icosahedron_triangles[index].corn ers;
h0 = icosahedron_corners[corners[0]]; h0 = icosahedron_corners[corners[0]];
h1 = icosahedron_corners[corners[1]]; h1 = icosahedron_corners[corners[1]];
h2 = icosahedron_corners[corners[2]]; h2 = icosahedron_corners[corners[2]];
skipping to change at line 181 skipping to change at line 181
break; break;
case 3: case 3:
h0 = t.e0; h0 = t.e0;
h1 = t.e1; h1 = t.e1;
h2 = t.e2; h2 = t.e2;
break; break;
} }
} }
} }
int GeodesicGrid::getPartnerTriangle(int lev, int index) const int StelGeodesicGrid::getPartnerTriangle(int lev, int index) const
{ {
if (lev==0) if (lev==0)
{ {
Q_ASSERT(index<20); Q_ASSERT(index<20);
return (index&1) ? index+1 : index-1; return (index&1) ? index+1 : index-1;
} }
switch(index&7) switch(index&7)
{ {
case 2: case 2:
case 6: case 6:
skipping to change at line 210 skipping to change at line 210
case 4: case 4:
return (lev==1) ? index-3 : (getPartnerTriangle(lev-1, index >>2)<<2)+1; return (lev==1) ? index-3 : (getPartnerTriangle(lev-1, index >>2)<<2)+1;
case 5: case 5:
return (lev==1) ? index-5 : (getPartnerTriangle(lev-1, index >>2)<<2)+0; return (lev==1) ? index-5 : (getPartnerTriangle(lev-1, index >>2)<<2)+0;
default: default:
Q_ASSERT(0); Q_ASSERT(0);
} }
return 0; return 0;
} }
void GeodesicGrid::initTriangle(int lev,int index, void StelGeodesicGrid::initTriangle(int lev,int index,
const Vec3d &c0, const Vec3d &c0,
const Vec3d &c1, const Vec3d &c1,
const Vec3d &c2) const Vec3d &c2)
{ {
Q_ASSERT((c0^c1)*c2 >= 0.0); Q_ASSERT((c0^c1)*c2 >= 0.0);
Triangle &t(triangles[lev][index]); Triangle &t(triangles[lev][index]);
t.e0 = c1+c2; t.e0 = c1+c2;
t.e0.normalize(); t.e0.normalize();
t.e1 = c2+c0; t.e1 = c2+c0;
t.e1.normalize(); t.e1.normalize();
skipping to change at line 234 skipping to change at line 234
if (lev < maxLevel) if (lev < maxLevel)
{ {
index *= 4; index *= 4;
initTriangle(lev,index+0,c0,t.e2,t.e1); initTriangle(lev,index+0,c0,t.e2,t.e1);
initTriangle(lev,index+1,t.e2,c1,t.e0); initTriangle(lev,index+1,t.e2,c1,t.e0);
initTriangle(lev,index+2,t.e1,t.e0,c2); initTriangle(lev,index+2,t.e1,t.e0,c2);
initTriangle(lev,index+3,t.e0,t.e1,t.e2); initTriangle(lev,index+3,t.e0,t.e1,t.e2);
} }
} }
void GeodesicGrid::visitTriangles(int maxVisitLevel, void StelGeodesicGrid::visitTriangles(int maxVisitLevel,
VisitFunc *func, VisitFunc *func,
void *context) const void *context) const
{ {
if (func && maxVisitLevel >= 0) if (func && maxVisitLevel >= 0)
{ {
if (maxVisitLevel > maxLevel) maxVisitLevel = maxLevel; if (maxVisitLevel > maxLevel) maxVisitLevel = maxLevel;
for (int i=0;i<20;i++) for (int i=0;i<20;i++)
{ {
const int *const corners = icosahedron_triangles[i]. corners; const int *const corners = icosahedron_triangles[i]. corners;
visitTriangles(0,i, visitTriangles(0,i,
icosahedron_corners[corners[0]], icosahedron_corners[corners[0]],
icosahedron_corners[corners[1]], icosahedron_corners[corners[1]],
icosahedron_corners[corners[2]], icosahedron_corners[corners[2]],
maxVisitLevel,func,context); maxVisitLevel,func,context);
} }
} }
} }
void GeodesicGrid::visitTriangles(int lev,int index, void StelGeodesicGrid::visitTriangles(int lev,int index,
const Vec3d &c0, const Vec3d &c0,
const Vec3d &c1, const Vec3d &c1,
const Vec3d &c2, const Vec3d &c2,
int maxVisitLevel, int maxVisitLevel,
VisitFunc *func, VisitFunc *func,
void *context) const void *context) const
{ {
(*func)(lev,index,c0,c1,c2,context); (*func)(lev,index,c0,c1,c2,context);
Triangle &t(triangles[lev][index]); Triangle &t(triangles[lev][index]);
lev++; lev++;
if (lev <= maxVisitLevel) if (lev <= maxVisitLevel)
{ {
index *= 4; index *= 4;
visitTriangles(lev,index+0,c0,t.e2,t.e1,maxVisitLevel,func,c ontext); visitTriangles(lev,index+0,c0,t.e2,t.e1,maxVisitLevel,func,c ontext);
visitTriangles(lev,index+1,t.e2,c1,t.e0,maxVisitLevel,func,c ontext); visitTriangles(lev,index+1,t.e2,c1,t.e0,maxVisitLevel,func,c ontext);
visitTriangles(lev,index+2,t.e1,t.e0,c2,maxVisitLevel,func,c ontext); visitTriangles(lev,index+2,t.e1,t.e0,c2,maxVisitLevel,func,c ontext);
visitTriangles(lev,index+3,t.e0,t.e1,t.e2,maxVisitLevel,func ,context); visitTriangles(lev,index+3,t.e0,t.e1,t.e2,maxVisitLevel,func ,context);
} }
} }
int GeodesicGrid::searchZone(const Vec3d &v,int searchLevel) const int StelGeodesicGrid::searchZone(const Vec3d &v,int searchLevel) const
{ {
for (int i=0;i<20;i++) for (int i=0;i<20;i++)
{ {
const int *const corners = icosahedron_triangles[i].corners; const int *const corners = icosahedron_triangles[i].corners;
const Vec3d &c0(icosahedron_corners[corners[0]]); const Vec3d &c0(icosahedron_corners[corners[0]]);
const Vec3d &c1(icosahedron_corners[corners[1]]); const Vec3d &c1(icosahedron_corners[corners[1]]);
const Vec3d &c2(icosahedron_corners[corners[2]]); const Vec3d &c2(icosahedron_corners[corners[2]]);
if (((c0^c1)*v >= 0.0) && ((c1^c2)*v >= 0.0) && ((c2^c0)*v > = 0.0)) if (((c0^c1)*v >= 0.0) && ((c1^c2)*v >= 0.0) && ((c2^c0)*v > = 0.0))
{ {
// v lies inside this icosahedron triangle // v lies inside this icosahedron triangle
skipping to change at line 309 skipping to change at line 309
i += 2; i += 2;
} }
else else
{ {
i += 3; i += 3;
} }
} }
return i; return i;
} }
} }
qWarning() << "ERROR GeodesicGrid::searchZone - not found"; qWarning() << "ERROR StelGeodesicGrid::searchZone - not found";
exit(-1); exit(-1);
// shut up the compiler // shut up the compiler
return -1; return -1;
} }
void GeodesicGrid::searchZones(const StelGeom::ConvexS& convex, void StelGeodesicGrid::searchZones(const StelGeom::ConvexS& convex,
int **inside_list,int **border_list, int **inside_list,int **border_list,
int maxSearchLevel) const int maxSearchLevel) const
{ {
if (maxSearchLevel < 0) maxSearchLevel = 0; if (maxSearchLevel < 0) maxSearchLevel = 0;
else if (maxSearchLevel > maxLevel) maxSearchLevel = maxLevel; else if (maxSearchLevel > maxLevel) maxSearchLevel = maxLevel;
#if defined __STRICT_ANSI__ || !defined __GNUC__ #if defined __STRICT_ANSI__ || !defined __GNUC__
int *halfs_used = new int[convex.size()]; int *halfs_used = new int[convex.size()];
#else #else
int halfs_used[convex.size()]; int halfs_used[convex.size()];
#endif #endif
skipping to change at line 356 skipping to change at line 356
corner_inside[icosahedron_triangles[i].corners[1 ]], corner_inside[icosahedron_triangles[i].corners[1 ]],
corner_inside[icosahedron_triangles[i].corners[2 ]], corner_inside[icosahedron_triangles[i].corners[2 ]],
inside_list,border_list,maxSearchLevel); inside_list,border_list,maxSearchLevel);
} }
#if defined __STRICT_ANSI__ || !defined __GNUC__ #if defined __STRICT_ANSI__ || !defined __GNUC__
delete[] halfs_used; delete[] halfs_used;
for(int ci=0; ci < 12; ci++) delete[] corner_inside[ci]; for(int ci=0; ci < 12; ci++) delete[] corner_inside[ci];
#endif #endif
} }
void GeodesicGrid::searchZones(int lev,int index, void StelGeodesicGrid::searchZones(int lev,int index,
const StelGeom::ConvexS& convex, const StelGeom::ConvexS& convex,
const int *indexOfUsedHalfSpaces, const int *indexOfUsedHalfSpaces,
const int halfSpacesUsed, const int halfSpacesUsed,
const bool *corner0_inside, const bool *corner0_inside,
const bool *corner1_inside, const bool *corner1_inside,
const bool *corner2_inside, const bool *corner2_inside,
int **inside_list,int **border_list, int **inside_list,int **border_list,
const int maxSearchLevel) const const int maxSearchLevel) const
{ {
#if defined __STRICT_ANSI__ || !defined __GNUC__ #if defined __STRICT_ANSI__ || !defined __GNUC__
skipping to change at line 455 skipping to change at line 455
} }
} }
#if defined __STRICT_ANSI__ || !defined __GNUC__ #if defined __STRICT_ANSI__ || !defined __GNUC__
delete[] halfs_used; delete[] halfs_used;
#endif #endif
} }
/************************************************************************* /*************************************************************************
Return a search result matching the given spatial region Return a search result matching the given spatial region
*************************************************************************/ *************************************************************************/
const GeodesicSearchResult* GeodesicGrid::search(const StelGeom::ConvexS& c onvex, int maxSearchLevel) const const GeodesicSearchResult* StelGeodesicGrid::search(const StelGeom::Convex S& convex, int maxSearchLevel) const
{ {
// Try to use the cached version // Try to use the cached version
if (maxSearchLevel==lastMaxSearchlevel && convex==lastSearchRegion) if (maxSearchLevel==lastMaxSearchlevel && convex==lastSearchRegion)
{ {
return cacheSearchResult; return cacheSearchResult;
} }
// Else recompute it and update cache parameters // Else recompute it and update cache parameters
lastMaxSearchlevel = maxSearchLevel; lastMaxSearchlevel = maxSearchLevel;
lastSearchRegion = convex; lastSearchRegion = convex;
cacheSearchResult->search(convex, maxSearchLevel); cacheSearchResult->search(convex, maxSearchLevel);
return cacheSearchResult; return cacheSearchResult;
} }
/************************************************************************* /*************************************************************************
Return a search result matching the given spatial region Return a search result matching the given spatial region
*************************************************************************/ *************************************************************************/
const GeodesicSearchResult* GeodesicGrid::search(const Vec3d &e0,const Vec3 d &e1,const Vec3d &e2,const Vec3d &e3,int maxSearchLevel) const const GeodesicSearchResult* StelGeodesicGrid::search(const Vec3d &e0,const Vec3d &e1,const Vec3d &e2,const Vec3d &e3,int maxSearchLevel) const
{ {
StelGeom::ConvexS c(e0, e1, e2, e3); StelGeom::ConvexS c(e0, e1, e2, e3);
return search(c,maxSearchLevel); return search(c,maxSearchLevel);
} }
GeodesicSearchResult::GeodesicSearchResult(const GeodesicGrid &grid) GeodesicSearchResult::GeodesicSearchResult(const StelGeodesicGrid &grid)
:grid(grid), :grid(grid),
zones(new int*[grid.getMaxLevel()+1]), zones(new int*[grid.getMaxLevel()+1]),
inside(new int*[grid.getMaxLevel()+1]), inside(new int*[grid.getMaxLevel()+1]),
border(new int*[grid.getMaxLevel()+1]) border(new int*[grid.getMaxLevel()+1])
{ {
for (int i=0;i<=grid.getMaxLevel();i++) for (int i=0;i<=grid.getMaxLevel();i++)
{ {
zones[i] = new int[GeodesicGrid::nrOfZones(i)]; zones[i] = new int[StelGeodesicGrid::nrOfZones(i)];
} }
} }
GeodesicSearchResult::~GeodesicSearchResult(void) GeodesicSearchResult::~GeodesicSearchResult(void)
{ {
for (int i=grid.getMaxLevel();i>=0;i--) for (int i=grid.getMaxLevel();i>=0;i--)
{ {
delete[] zones[i]; delete[] zones[i];
} }
delete[] border; delete[] border;
delete[] inside; delete[] inside;
delete[] zones; delete[] zones;
} }
void GeodesicSearchResult::search(const StelGeom::ConvexS& convex, void GeodesicSearchResult::search(const StelGeom::ConvexS& convex,
int maxSearchLevel) int maxSearchLevel)
{ {
for (int i=grid.getMaxLevel();i>=0;i--) for (int i=grid.getMaxLevel();i>=0;i--)
{ {
inside[i] = zones[i]; inside[i] = zones[i];
border[i] = zones[i]+GeodesicGrid::nrOfZones(i); border[i] = zones[i]+StelGeodesicGrid::nrOfZones(i);
} }
grid.searchZones(convex,inside,border,maxSearchLevel); grid.searchZones(convex,inside,border,maxSearchLevel);
} }
void GeodesicSearchInsideIterator::reset(void) void GeodesicSearchInsideIterator::reset(void)
{ {
level = 0; level = 0;
maxCount = 1<<(maxLevel<<1); // 4^maxLevel maxCount = 1<<(maxLevel<<1); // 4^maxLevel
indexP = r.zones[0]; indexP = r.zones[0];
endP = r.inside[0]; endP = r.inside[0];
 End of changes. 19 change blocks. 
19 lines changed or deleted 19 lines changed or added

This html diff was produced by rfcdiff 1.41. The latest version is available from http://tools.ietf.org/tools/rfcdiff/