/*
 
StelGeodesicGrid: a library for dividing the sphere into triangle zones
by subdividing the icosahedron
 
Author and Copyright: Johannes Gajdosik, 2006
 
This library requires a simple Vector library,
which may have different copyright and license,
for example Vec3f from VecMath.hpp.
 
In the moment I choose to distribute the library under the GPL,
later I may choose to additionally distribute it under a more
relaxed license like the LGPL. If you want to have the library
under another license, please ask me.
 
 
 
This library is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License
as published by the Free Software Foundation; either version 2
of the License, or (at your option) any later version.
 
This library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.
 
You should have received a copy of the GNU General Public License
along with this library; if not, write to the Free Software
Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
 
*/

#ifndef _STELGEODESICGRID_HPP_
#define _STELGEODESICGRID_HPP_

#include "StelSphereGeometry.hpp"

class GeodesicSearchResult;

//! @class StelGeodesicGrid
//! Grid of triangles (zones) on the sphere with radius 1, generated by subdividing the icosahedron.
//! level 0: just the icosahedron, 20 zones
//! level 1: 80 zones, level 2: 320 zones, ...
//! Each zone has a unique integer number in the range [0,getNrOfZones()-1].
class StelGeodesicGrid
{
public:
	StelGeodesicGrid(int maxLevel);
	~StelGeodesicGrid(void);
	
	int getMaxLevel(void) const {return maxLevel;}
	
	static int nrOfZones(int level) {return (20<<(level<<1));} // 20*4^level
	
	int getNrOfZones(void) const {return nrOfZones(maxLevel);}
	
	typedef void (VisitFunc)(int lev,int index,
							 const Vec3f &c0,
							 const Vec3f &c1,
							 const Vec3f &c2,
	                         void *context);
	void visitTriangles(int maxVisitLevel, VisitFunc *func,void *context) const;

	//! Find the zone number in which a given point lies.
	//! prerequisite: v*v==1
	//! When the point lies on the border of two or more zones,
	//! one such zone is returned (always the same one,
	//! because the algorithm is deterministic).
	int getZoneNumberForPoint(const Vec3f &v,int searchLevel) const;

	//! Return the position of the 3 corners for the triangle at the given level and index
	void getTriangleCorners(int lev, int index, Vec3f& c0, Vec3f& c1, Vec3f& c2) const;
	
	//! Return the index of the partner triangle with which to form a parallelogram
	int getPartnerTriangle(int lev, int index) const;
	
	//! Return a search result matching the given spatial region
	//! The result is cached, meaning that it is very fast to search the same region consecutively
	//! @return a GeodesicSearchResult instance which must be used with GeodesicSearchBorderIterator and GeodesicSearchInsideIterator
	const GeodesicSearchResult* search(const QVector<SphericalCap>& convex, int maxSearchLevel) const;

private:
	friend class GeodesicSearchResult;
	
	//! Find all zones that lie fully(inside) or partly(border)
	//! in the intersection of the given half spaces.
	//! The result is accurate when (0,0,0) lies on the border of
	//! each half space. If this is not the case,
	//! the result may be inaccurate, because it is assumed, that
	//! a zone lies in a half space when its 3 corners lie in this half space.
	//! inside[l] points to the begin of an integer array of size nrOfZones(l),
	//! border[l] points to one after the end of the same integer array.
	//! The array will be filled from the beginning with the inside zone numbers
	//! and from the end with the border zone numbers of the given level l
	//! for 0<=l<=getMaxLevel().
	//! inside[l] will not contain zones that are already contained
	//! in inside[l1] for some l1 < l.
	//! In order to restrict search depth set maxSearchLevel < maxLevel,
	//! for full search depth set maxSearchLevel = maxLevel,
	void searchZones(const QVector<SphericalCap>& convex,
					 int **inside,int **border,int maxSearchLevel) const;
	
	const Vec3f& getTriangleCorner(int lev, int index, int cornerNumber) const;
	void initTriangle(int lev,int index,
					  const Vec3f &c0,
					  const Vec3f &c1,
					  const Vec3f &c2);
	void visitTriangles(int lev,int index,
						const Vec3f &c0,
						const Vec3f &c1,
						const Vec3f &c2,
	                    int maxVisitLevel,
	                    VisitFunc *func,
	                    void *context) const;
	void searchZones(int lev,int index,
	                 const QVector<SphericalCap>& convex,
	                 const int *indexOfUsedSphericalCaps,
	                 const int halfSpacesUsed,
	                 const bool *corner0_inside,
	                 const bool *corner1_inside,
	                 const bool *corner2_inside,
	                 int **inside,int **border,int maxSearchLevel) const;

	const int maxLevel;
	struct Triangle
	{
		Vec3f e0,e1,e2;   // Seitenmittelpunkte
	};
	Triangle **triangles;
	// 20*(4^0+4^1+...+4^n)=20*(4*(4^n)-1)/3 triangles total
	// 2+10*4^n corners
	
	//! A cached search result used to avoid doing twice the same search
	mutable GeodesicSearchResult* cacheSearchResult;
	mutable int lastMaxSearchlevel;
	mutable QVector<SphericalCap> lastSearchRegion;
};

class GeodesicSearchResult
{
public:
	GeodesicSearchResult(const StelGeodesicGrid &grid);
	~GeodesicSearchResult(void);
	void print(void) const;
private:
	friend class GeodesicSearchInsideIterator;
	friend class GeodesicSearchBorderIterator;
	friend class StelGeodesicGrid;
	
	void search(const QVector<SphericalCap>& convex, int maxSearchLevel);
	
	const StelGeodesicGrid &grid;
	int **const zones;
	int **const inside;
	int **const border;
};

class GeodesicSearchBorderIterator
{
public:
	GeodesicSearchBorderIterator(const GeodesicSearchResult &ar,int alevel)
		: r(ar),level((alevel<0)?0:(alevel>ar.grid.getMaxLevel())
			             ?ar.grid.getMaxLevel():alevel),
			end(ar.zones[GeodesicSearchBorderIterator::level]+
			    StelGeodesicGrid::nrOfZones(GeodesicSearchBorderIterator::level))
	{reset();}
	void reset(void) {index = r.border[level];}
	int next(void) // returns -1 when finished
	{if (index < end) {return *index++;} return -1;}
private:
	const GeodesicSearchResult &r;
	const int level;
	const int *const end;
	const int *index;
};


class GeodesicSearchInsideIterator
{
public:
	GeodesicSearchInsideIterator(const GeodesicSearchResult &ar,int alevel)
		: 	r(ar), 
			maxLevel((alevel<0)?0:(alevel>ar.grid.getMaxLevel())?ar.grid.getMaxLevel():alevel)
	{reset();}
	void reset(void);
	int next(void); // returns -1 when finished
private:
	const GeodesicSearchResult &r;
	const int maxLevel;
	int level;
	int maxCount;
	int *indexP;
	int *endP;
	int index;
	int count;
};

#endif // _STELGEODESICGRID_HPP_