out of core rendering (not working)

Rendering Engine (Not Working)

    A. First Edition

    This is just a temporary version of rendering framework which includes preprocessing and basic view frustum

    culling.

    ¡¡

    B.The problem

1. Convert original "big-endian" format to "little-endian" with supplied software, merge total 1100 model files into one big file. Also during this pass calculate the bounding of whole data for octree construction.
¡¡
2. Using two-pass to sort data model into octree by loading segments by segments of whole data into memory. (I think it is impossible to do it in one pass because the sequence of vertex is changed and the triangle's index array can only be modified until vertex is sorted.) So, the first pass to sort vertex and second pass to sort triangle.
¡¡
3. Precalculate all neighbours of each octree node and store this infomation in octree nodes. I think this is the best solution for "neighbour-finding" in octree.
¡¡
4. All data structure will be saved as memory image in file so that at run time just simple loading and update of starting address of array is needed. To do that, all "pointers" are just index numbers which can used when reloading.
¡¡
5. The triangles at boundary which run across multiple octree nodes will be calculated multiple times in each node. Therefore the number of triangle will inflate after sorting about 30%.
¡¡
6. The preprocess of 11million vertex and 12million triangles take a little more than one hour. (After preprocessing, the number of triangle becomes 16million.) There are total 44414 octree nodes and total memory usage will be nearly 9M and will be saved as memory image in file. The maximum depth of octree is 12 and total pre-calculated neighbour nodes for all octree nodes are 265k.
¡¡
7. I use 1000 vertex as threshold for octree leaf node. Maybe it is a little bit too big? What do you think? The model is "power plant" of which the length of bounding box is more than 80000. Should I decide this "leaf node" size according to each model?

****************************************************************************************************************************************************

1. I changed my strategy because the triangles which run across octree node are impossible for rendering because we load data page by page, or node by node.

2. The index of triangles are now changed to local index. i.e. index for local octree node because absolute index has no meaning for individual octree nodes.

3. Data of both vertex and triangles are now NOT consecutive. Instead they are aligned with page size or 64k.

4. After first preprocessing, the octree can be saved directly in disk and loaded from disk by first reading a small "info.data" file.

5. The preprocessing of 12M triangles and 11M vertices takes about 40Mins in machine in lab, but should be longer in home.

6. File mapping limit is only about size which maybe up to nearly 2G, depending on whole image size.

¡¡

¡¡

¡¡

C.The idea of program
¡¡

There is one subtle problem which is worth mentioning. Originally the octree is designed to be an "self-incremental" array.

And I am also very careful about the invalid addressing after reallocating array. i.e.

void Octree::createNewSonNode() {

Octee*ptr;

ptr=root+getNewNodeIndex(); //====>this kind of thing I am very careful! you see the "root" might be changed in

//getNewNodeIndex which "re-allocate" the array. So, the "ptr" will be invalid.

ptr->level=level+1; //====>>Big Problem!! Because "member data" "level" becomes invalid after "octree array" is moved!

¡¡

¡¡

However, this is NOT enough! Even those "member data" is also becoming invalid. This is really subtle and tricky!

    ¡¡

    D.The major functions

    ¡¡

    E.Further improvement

    ¡¡

    F.File listing

    1. plyModel.h

    2. Octree.h

    3. Octree.cpp

    4. rply.h   (downloaded and modified a bit)

    5. rply.c    (downloaded and modified a bit)

    6. frustum.h  (downloaded)

    7. frustum.cpp  (downloaded)

    8. plyconverter.cpp  (downloaded and modified a bit)

    9. octreeLoader.h

    10. octreeLoader.cpp

    11. renderEngine.cpp (main)

    ¡¡

    file name: plyModel.h

    #ifndef PLYMODEL_H
    #define PLYMODEL_H

    typedef float FVector[3];
    typedef unsigned char UCVector[3];
    typedef unsigned int UIVector[3];
    typedef double DVector[3];

    struct PlyVertex
    {
    	FVector coord;
    	UCVector color;
    	FVector normal;
    };

    struct PlyFace
    {
    	unsigned char number;
    	int* indexList;
    	UCVector color;
    };

    struct PlyModel
    {
    	int vNumber;//internal usage only
    	int fNumber;//internal usage only
    	PlyVertex* vertex;
    	UIVector* triangle;
    	int vertexCount;
    	int faceCount;
    	FVector maxCoord;
    	FVector minCoord;
    	float maxEdge;

    	//PlyFace* faceList;
    	void calcLongestEdge(UIVector tri);
    	bool analysis(char* buf);
    	int readFile(const char* fileName);
    	int convertFile(const char* fileName);
    	void readDirectory(char* dir);

    };

    void sortTriangle(UIVector triangle);

    #endif

    ¡¡

    ¡¡

    file name: octree.h

    /**********************************************************************
    E:\working>octree
    finished vertex processing and now start face processing...
    finished face processing, now start octree calculation...
    total vertex is 11070509 face is 16364639
    total process time 4623141
    octree info: nodeCount=44414, nodeSize=204
    totalNeighbourCount=265602
    max depth of octree is 12

    E:\working>

    *********************************************************************/
    #ifndef OCTREE_H
    #define OCTREE_H

    #include "plyModel.h"
    #include <windows.h>
    #include <cstdio>
    #include <cstdlib>
    #include <stack>

    using namespace std;

    const int PageSize=64*1024;
    //const int MaxNodeArrayLength=1024*256;
    const int MaxNodeArrayLength=64*1024;
    const int MemoryThreshold=1000;
    //const int TempFileStartOffset=sizeof(int);
    //const int NodeArrayLengthIncrement=10;
    const int MaxNeighbourCount=26;

    struct InfoStruct
    {
    	int nodeCount;
    	int nodeSize;
    	int memoryThreshold;
    	float maxEdge;
    	int totalNeighbourCount;
    	int maxDepth;
    };

    struct Octree
    {
    	enum OctreeNodeEnum
    	{
    		TOP_LEFT_FRONT =2, TOP_LEFT_BACK =3,
    		TOP_RIGHT_BACK =7, TOP_RIGHT_FRONT =6,
    		BOTTOM_LEFT_FRONT =0, BOTTOM_LEFT_BACK =1,
    		BOTTOM_RIGHT_BACK =5, BOTTOM_RIGHT_FRONT =4
    	};
    	static int nodeCount;
    	static int maxDepth;
    	static Octree* root;
    	static int nodeArrayLength;
    	static int totalNeighbourCount;
    	static float maxEdge;
    	static int currentVPageIndex;
    	static int currentFPageIndex;

    	////////////////////////////////////////////////////
    	//here goes the preprocessing methods and they will be used only once
    	static void createOctreeFromFile(char* vertexFileName, char* faceFileName,
    		char* infoFileName, char* outVertexFileName, char* outFaceFileName,
    		char* octreeFileName);

    	static void calculatePageIndex();
    	static void calculateNeighbour();
    	void doCalculateNeighbour(Octree* theNode);

    	void doCreateOctree(HANDLE outVertexHandle);

    	int doTraverseOctree(HANDLE outFaceHandle);

    	/////////////////////////////////////////////
    	//utility methods
    	int relativePosition(FVector vPoint, FVector vCtr);
    	int relativePosition(UIVector triangle);
    	void setSonNodeCentre(int nodeID, FVector vNodeCenter);
    	void setupNewSonNode(int nodeID, int sonNodePos);

    	/////////////////////////////////////////////
    	//constructor
    	int createNewOctreeNode(int nodeID);

    	////////////////////////////////////////////////////
    	//the following is all data and try to make sure they are smaller than 64bytes
    	FVector center;
    	float width;//from center to edge, the half width of edge
    	int vNumber;
    	int fNumber;
    	int vStartIndex;
    	int fStartIndex;
    	int myNodeIndex;
    	int parentIndex;
    	int sons[8];
    	int vStartPageIndex;//vertex
    	int vPageCount;
    	int fStartPageIndex;//face
    	int fPageCount;
    	void* reserved1;
    	void* reserved2;
    	LPVOID vAddress;
    	LPVOID fAddress;
    	unsigned int vOffset;
    	unsigned int fOffset;
    	unsigned char sonCount;
    	unsigned char level;
    	unsigned char neighbourCount;
    	int neighbours[MaxNeighbourCount];//maximumly 26 neighbours

    };

    #endif

    file name: octree.cpp

    ///////////////////////////////////////////////////////////////
    //The molton code scheme is as simple as following:
    //1. using an unsigned int (32) for code.
    //2. The first 4 bits represents level, , root node starting from 0x0 (this limits our deepest level from 16
    //3. each level using 3 bits representing son index.
    // I forget this after only two weeks
    //
    ////////////////////////////////////////////////////////////////
    ////we use the largest radius for radius of octree node
    //
    //////////////////////////////////////////////////////////////

    #include <stdio.h>
    #include <math.h>
    #include "Octree.h"
    #include <windows.h>
    #include <winbase.h>

    #include "PlyModel.h"

    //////////////////////////////////////////////

    #define CHECK_NEW(x) if (x==NULL){printf("memory alloc failed\n"); exit(1);}

    const double Epsilon=0.001;

    int Octree::maxDepth=0;
    int Octree::nodeCount=0;
    Octree* Octree::root=NULL;
    int Octree::nodeArrayLength=0;
    int* hashTable;
    int vertexIndexCounter=0;
    int Octree::totalNeighbourCount=0;
    float Octree::maxEdge=0.0;
    int Octree::currentVPageIndex=0;
    int Octree::currentFPageIndex=0;

    HANDLE dumpHandle;

    void dumpEveryThing()
    {
    	unsigned long length;
    	WriteFile(dumpHandle, &Octree::nodeCount, sizeof(int), &length, NULL);
    	WriteFile(dumpHandle, &vertexIndexCounter, sizeof(int), &length, NULL);
    	CloseHandle(dumpHandle);
    }

    //return the new node index
    int Octree::createNewOctreeNode(int nodeID)
    {
    	Octree *ptr;
    	//lazy
    	if (nodeCount==MaxNodeArrayLength)
    	{
    		printf("max node array length reached!\n");
    		exit(78);

    		/*
    		newRoot=new Octree[nodeArrayLength+NodeArrayLengthIncrement];
    		if (newRoot==NULL)
    		{
    			printf("allocate new octree node array failed\n");
    			exit(1);
    		}
    		memcpy(newRoot, root, sizeof(Octree)*nodeCount);
    		delete[] root;
    		root=newRoot;
    		nodeArrayLength+=NodeArrayLengthIncrement;
    		*/
    	}

    	ptr=root+nodeCount;
    	///////////////////
    	//parent
    	sonCount++;
    	sons[nodeID]=nodeCount;

    	setSonNodeCentre(nodeID, ptr->center);
    	////////////////////////////////
    	//son
    	ptr->width=width/2.0;
    	ptr->myNodeIndex=nodeCount;//this will not change
    	ptr->level=level+1;
    	ptr->parentIndex=myNodeIndex;

    	//these are initialization to make things neat
    	ptr->fNumber=ptr->vNumber=ptr->fStartIndex=ptr->vStartIndex=-1;
    	ptr->vPageCount=ptr->fPageCount=ptr->vStartPageIndex=ptr->fStartPageIndex=-1;
    	ptr->vAddress=ptr->fAddress=NULL;
    	ptr->vOffset=ptr->fOffset=0;////////////////////now it is 0
    	ptr->neighbourCount=0;

    	if (ptr->level>maxDepth)
    	{
    		maxDepth=ptr->level;
    	}
    	for (int i=0; i<8; i++)
    	{
    		ptr->sons[i]=-1;
    	}
    	ptr->sonCount=0;

    	return 	nodeCount++;;
    }

    /////////////////////////////////////////

    void Octree::setSonNodeCentre(int nodeID, FVector vNodeCenter)
    {
    	switch(nodeID)
        {
    	case TOP_LEFT_FRONT:
    		//vNodeCenter = TPoint(vCtr[0] - width/4, vCtr[1] + width/4, vCtr[2] + width/4);
    		vNodeCenter[0] = center[0] - width/4;
    		vNodeCenter[1] = center[1] + width/4;
    		vNodeCenter[2] = center[2] + width/4;
    		break;

    	case TOP_LEFT_BACK:
    		vNodeCenter[0] = center[0] - width/4;
    		vNodeCenter[1] = center[1] + width/4;
    		vNodeCenter[2] = center[2] - width/4;
    		break;

    	case TOP_RIGHT_BACK:
    		//vNodeCenter = TPoint(vCtr[0] + width/4, vCtr[1] + width/4, vCtr[2] - width/4);
    		vNodeCenter[0] = center[0] + width/4;
    		vNodeCenter[1] = center[1] + width/4;
    		vNodeCenter[2] = center[2] - width/4;
    		break;

    	case TOP_RIGHT_FRONT:
    		//vNodeCenter = TPoint(vCtr[0] + width/4, vCtr[1] + width/4, vCtr[2] + width/4);
    		vNodeCenter[0] = center[0] + width/4;
    		vNodeCenter[1] = center[1] + width/4;
    		vNodeCenter[2] = center[2] + width/4;
    		break;

    	case BOTTOM_LEFT_FRONT:
    		//vNodeCenter = TPoint(vCtr[0] - width/4, vCtr[1] - width/4, vCtr[2] + width/4);
    		vNodeCenter[0] = center[0] - width/4;
    		vNodeCenter[1] = center[1] - width/4;
    		vNodeCenter[2] = center[2] + width/4;
    		break;

    	case BOTTOM_LEFT_BACK:
    		//vNodeCenter = TPoint(vCtr[0] - width/4, vCtr[1] - width/4, vCtr[2] - width/4);
    		vNodeCenter[0] = center[0] - width/4;
    		vNodeCenter[1] = center[1] - width/4;
    		vNodeCenter[2] = center[2] - width/4;
    		break;

    	case BOTTOM_RIGHT_BACK:
    		//vNodeCenter = TPoint(vCtr[0] + width/4, vCtr[1] - width/4, vCtr[2] - width/4);
    		vNodeCenter[0] = center[0] + width/4;
    		vNodeCenter[1] = center[1] - width/4;
    		vNodeCenter[2] = center[2] - width/4;
    		break;

    	case BOTTOM_RIGHT_FRONT:
    		//vNodeCenter = TPoint(vCtr[0] + width/4, vCtr[1] - width/4, vCtr[2] + width/4);
    		vNodeCenter[0] = center[0] + width/4;
    		vNodeCenter[1] = center[1] - width/4;
    		vNodeCenter[2] = center[2] + width/4;
    		break;
        }
    }

    int Octree::relativePosition(FVector vPoint, FVector vCtr)
    {

    	if( (vPoint[0] <= vCtr[0]) && (vPoint[1] >= vCtr[1]) && (vPoint[2] >= vCtr[2]) )
    	{
    		return TOP_LEFT_FRONT;
    	}

    	if( (vPoint[0] <= vCtr[0]) && (vPoint[1] >= vCtr[1]) && (vPoint[2] <= vCtr[2]) )
    	{
    		return TOP_LEFT_BACK;
    	}
    	if( (vPoint[0] >= vCtr[0]) && (vPoint[1] >= vCtr[1]) && (vPoint[2] <= vCtr[2]) )
    	{
    		return TOP_RIGHT_BACK;
    	}

    	if( (vPoint[0] >= vCtr[0]) && (vPoint[1] >= vCtr[1]) && (vPoint[2] >= vCtr[2]) )
    	{
    		return TOP_RIGHT_FRONT;
    	}
    	if( (vPoint[0] <= vCtr[0]) && (vPoint[1] <= vCtr[1]) && (vPoint[2] >= vCtr[2]) )
    	{
    		return BOTTOM_LEFT_FRONT;
    	}

    	if( (vPoint[0] <= vCtr[0]) && (vPoint[1] <= vCtr[1]) && (vPoint[2] <= vCtr[2]) )
    	{
    		return BOTTOM_LEFT_BACK;
    	}

    	if( (vPoint[0] >= vCtr[0]) && (vPoint[1] <= vCtr[1]) && (vPoint[2] <= vCtr[2]) )
    	{
    		return BOTTOM_RIGHT_BACK;
    	}

    	if( (vPoint[0] >= vCtr[0]) && (vPoint[1] <= vCtr[1]) && (vPoint[2] >= vCtr[2]) )
    	{
    		return BOTTOM_RIGHT_FRONT;
    	}
    	//we should never reach this
    	return -1;
    }

    //hard coding page size=64k
    //this is deprecated!!! because we allocate page index during construction
    //such that data is alignment with page index
    void Octree::calculatePageIndex()
    {
    	int vertexSize=sizeof(PlyVertex);
    	int faceSize=sizeof(UIVector);
    	int startOffset, endOffset;
    	Octree*ptr;

    	for (int i=0; i<nodeCount; i++)
    	{
    		ptr=root+i;

    		startOffset=ptr->vStartIndex*vertexSize;
    		endOffset=(ptr->vStartIndex+ptr->vNumber)*vertexSize;
    		ptr->vStartPageIndex=startOffset/PageSize;

    		ptr->vOffset=startOffset-ptr->vStartPageIndex*PageSize;

    		ptr->vPageCount=endOffset/PageSize-startOffset/PageSize +1;
    		if (ptr->vPageCount==0)
    		{
    			ptr->vPageCount=1;//at least one
    		}

    		//repeat

    		startOffset=ptr->fStartIndex*faceSize;
    		endOffset=(ptr->fStartIndex+ptr->fNumber)*faceSize;
    		ptr->fStartPageIndex=startOffset/PageSize;

    		ptr->fOffset=startOffset-ptr->fStartPageIndex*PageSize;

    		ptr->fPageCount=endOffset/PageSize-startOffset/PageSize +1;

    	}
    }

    void Octree::doCalculateNeighbour(Octree* theNode)
    {
    	Octree* ptr;
    	int i;
    	float distance;
    	if (theNode==this)
    	{
    		return;
    	}
    	if (level<theNode->level)
    	{
    		for (i=0; i<8; i++)
    		{
    			if (sons[i]!=-1)
    			{
    				ptr=root+sons[i];
    				ptr->doCalculateNeighbour(theNode);
    			}
    		}
    	}
    	else
    	{
    		if (level==theNode->level)
    		{
    			//check all sides to see if they colinear
    			for (i=0; i<3; i++)
    			{
    				//all coord must be exact one width away or colinear
    				distance=fabs(center[i]-theNode->center[i]);
    				if (distance>Epsilon)
    				{
    					distance=fabs(distance-width);
    					if (distance>Epsilon)
    					{
    						//no! this is far away
    						return;
    					}
    				}
    			}
    			//got it
    			theNode->neighbours[theNode->neighbourCount]=myNodeIndex;
    			theNode->neighbourCount++;
    			totalNeighbourCount++;
    			if (theNode->neighbourCount>MaxNeighbourCount)
    			{
    				printf("corrupted data or neighbour calculation mistake\n");
    				exit(77);
    			}
    		}
    	}
    	//in all case we end
    }

    void Octree::calculateNeighbour()
    {
    	Octree* ptr;
    	for (int i=0; i<nodeCount; i++)
    	{
    		ptr=root+i;
    		root->doCalculateNeighbour(ptr);
    	}
    }

    void Octree::createOctreeFromFile(char* vertexFileName, char* faceFileName, char* infoFileName,
    		char* outVertexFileName, char* outFaceFileName, char* octreeFileName)
    {
    	HANDLE infoHandle, inVertexHandle, outVertexHandle, inFaceHandle, outFaceHandle, octreeHandle;
    	Octree* ptr;
    	PlyModel ply;
    	int i;
    	float sideWidth;
    	unsigned long length;

    	infoHandle=CreateFile(infoFileName, GENERIC_READ|GENERIC_WRITE, 0, NULL, OPEN_EXISTING,
    		FILE_ATTRIBUTE_NORMAL , NULL);

    	if (infoHandle==INVALID_HANDLE_VALUE)
    	{
    		printf("open file failed\n");
    		exit(3);
    	}
    	dumpHandle=CreateFile("dumpfile.data", GENERIC_WRITE, 0, NULL, CREATE_ALWAYS,
    		FILE_ATTRIBUTE_NORMAL , NULL);

    	//MUST OVER WRITE
    	octreeHandle=CreateFile(octreeFileName, GENERIC_WRITE, 0, NULL, CREATE_ALWAYS,
    		FILE_ATTRIBUTE_NORMAL , NULL);

    	if (ReadFile(infoHandle, &ply, sizeof(PlyModel), &length, NULL)==0)
    	{
    		printf("read info failed\n");
    		exit(4);
    	}
    	inVertexHandle=CreateFile(vertexFileName, GENERIC_READ, 0, NULL, OPEN_EXISTING,
    		FILE_ATTRIBUTE_NORMAL , NULL);

    	outVertexHandle=CreateFile(outVertexFileName, GENERIC_WRITE, 0, NULL, CREATE_ALWAYS,
    		FILE_ATTRIBUTE_NORMAL , NULL);

    	inFaceHandle=CreateFile(faceFileName, GENERIC_READ, 0, NULL, OPEN_EXISTING,
    		FILE_ATTRIBUTE_NORMAL , NULL);

    	outFaceHandle=CreateFile(outFaceFileName, GENERIC_WRITE, 0, NULL, CREATE_ALWAYS,
    		FILE_ATTRIBUTE_NORMAL , NULL);
    	if (outVertexHandle==INVALID_HANDLE_VALUE||inVertexHandle==INVALID_HANDLE_VALUE
    		||outFaceHandle==INVALID_HANDLE_VALUE||inFaceHandle==INVALID_HANDLE_VALUE
    		||octreeHandle==INVALID_HANDLE_VALUE||dumpHandle==INVALID_HANDLE_VALUE)
    	{
    		printf("cannot create output vertex file\n");
    		exit(9);
    	}

    	//nodeArrayLength+=NodeArrayLengthIncrement;
    	//root=new Octree[nodeArrayLength];
    	//CHECK_NEW(root);
    	hashTable=new int[ply.vertexCount];
    	CHECK_NEW(hashTable);
    	root=new Octree[MaxNodeArrayLength];
    	CHECK_NEW(root);
    	//calculate the center and width of root
    	ptr=root;
    	nodeCount=1;
    	maxEdge=ply.maxEdge;
    	ptr->width=0.0;
    	for (i=0; i<3; i++)
    	{
    		ptr->center[i]=(ply.minCoord[i]+ply.maxCoord[i])/2.0;
    		sideWidth=ply.maxCoord[i]-ply.minCoord[i];
    		if (ptr->width<sideWidth)
    		{
    			ptr->width=sideWidth;
    		}
    	}
    	ptr->width;
    	ptr->level=0;
    	ptr->parentIndex=-1;
    	ptr->fStartIndex=0;
    	ptr->fNumber=ply.faceCount;
    	ptr->vStartIndex=0;
    	ptr->vNumber=ply.vertexCount;
    	ptr->neighbourCount=0;
    	for(i=0; i<8; i++)
    	{
    		ptr->sons[i]=-1;
    	}
    	ptr->sonCount=0;
    	//hard coding this;
    	ptr->myNodeIndex=0;
    	ptr->vAddress=ptr->fAddress=NULL;
    	ptr->fOffset=ptr->vOffset=-1;
    	ptr->vPageCount=ptr->fPageCount=-1;
    	ptr->vStartPageIndex=ptr->fStartPageIndex=-1;

    	if (ply.vertexCount>MemoryThreshold)
    	{
    		ptr->reserved1=inVertexHandle;
    		ptr->reserved2=INVALID_HANDLE_VALUE;
    		ptr->doCreateOctree(outVertexHandle);
    	}

    	printf("finished vertex processing and now start face processing...\n");
    	if (ptr->sonCount>0)
    	{
    		ptr->reserved1=inFaceHandle;
    		ptr->fNumber=ptr->doTraverseOctree(outFaceHandle);
    	}

    	ply.faceCount=ptr->fNumber;

    	delete[] hashTable;
    	printf("finished face processing, now start octree calculation...\n");
    	//calculatePageIndex();//deprecated

    	calculateNeighbour();

    	InfoStruct infoStruct;
    	infoStruct.nodeCount=Octree::nodeCount;
    	infoStruct.nodeSize=sizeof(Octree);
    	infoStruct.totalNeighbourCount=Octree::totalNeighbourCount;
    	infoStruct.maxDepth=Octree::maxDepth;
    	infoStruct.memoryThreshold=MemoryThreshold;
    	//important!!!
    	SetFilePointer(infoHandle, 0, NULL, FILE_BEGIN);
    	if (WriteFile(infoHandle, &infoStruct, sizeof(InfoStruct), &length, NULL)==0)
    	{
    		printf("write back info file failed\n");
    		exit(82);
    	}

    	printf("total vertex is %d face is %d\n", ply.vertexCount, ply.faceCount);

    	if (WriteFile(octreeHandle, Octree::root, sizeof(Octree)*Octree::nodeCount, &length, NULL)==0)
    	{
    		printf("write octree file failed\n");
    		exit(82);
    	}
    	delete[] root;
    	CloseHandle(inVertexHandle);
    	CloseHandle(octreeHandle);
    	CloseHandle(outVertexHandle);
    	CloseHandle(inFaceHandle);
    	CloseHandle(outFaceHandle);
    	CloseHandle(infoHandle);
    }

    /*
    void Octree::setupNewSonNode(int nodeID, int sonNodePos)
    {
    	Octree* ptr;
    	ptr=root+ sonNodePos;
    	sons[nodeID]=sonNodePos;
    	ptr->parentIndex=myNodeIndex;
    	setSonNodeCentre(nodeID, ptr->center);
    	ptr->width=width/2.0;
    }
    */

    //this will return the actual face after sorting
    int Octree::doTraverseOctree(HANDLE outFaceHandle)
    {
    	HANDLE tempFileHandle[8];
    	char tempFileName[8][MAX_PATH];
    	int i, j, k;
    	unsigned long length;
    	int tempFaceCounter[8];
    	UIVector triangle;
    	HANDLE faceHandle;

    	faceHandle=(HANDLE)reserved1;

    	//remember to read from beginning
    	SetFilePointer(faceHandle, 0, NULL, FILE_BEGIN);

    	if (sonCount==0)
    	{
    		///////////////////////////////////////
    		//we also want to allocate face data alignment with page index
    		//and we want to leave 40% space for extra triangle in future
    		fStartPageIndex=currentFPageIndex;
    		fPageCount=fNumber*sizeof(UIVector)*7/5/PageSize+1;
    		currentFPageIndex+=fPageCount;
    		SetFilePointer(outFaceHandle, fStartPageIndex*PageSize, NULL, FILE_BEGIN);

    		/////////////////////////////////////////////////////////////////////////

    		for (i=0; i<fNumber; i++)
    		{
    			if (ReadFile(faceHandle, triangle, sizeof(UIVector), &length, NULL)==0)
    			{
    				printf("read temp vertex file failed\n");
    				exit(11);
    			}
    			for (j=0; j<3; j++)
    			{
    				triangle[j]-=vStartIndex;
    			}

    			if (WriteFile(outFaceHandle, triangle, sizeof(UIVector), &length, NULL)==0)
    			{
    				printf("write output vertex file failed\n");
    				exit(11);
    			}
    		}
    		//////////////////////////////////////
    		return fNumber;
    	}

    	//else
    	for (i=0; i<8; i++)
    	{
    		tempFaceCounter[i]=0;
    		if (GetTempFileName(".", "face", 0, tempFileName[i])==0)
    		{
    			printf("get temp file name failed\n");
    			exit(44);
    		}
    		tempFileHandle[i]=CreateFile(tempFileName[i], GENERIC_READ|GENERIC_WRITE, 0, NULL,
    			OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL);

    		if (tempFileHandle[i]==INVALID_HANDLE_VALUE)
    		{
    			printf("open temp file failed %d\n", GetLastError());
    			exit(55);
    		}
    	}

    	for (i=0; i<fNumber; i++)
    	{
    		if (ReadFile(faceHandle, triangle, sizeof(UIVector), &length, NULL)==0)
    		{
    			printf("read file fails\n");
    			exit(2);
    		}
    		//update vertex index
    		if (level==0)
    		{
    			for (j=0; j<3; j++)
    			{
    				triangle[j]=hashTable[triangle[j]];
    			}
    			//in this case, we don't have to sort triangle because it is still
    			//possible in different nodes, but it makes neat
    			sortTriangle(triangle);
    		}
    		for (j=0; j<8; j++)
    		{
    			Octree* ptr;
    			bool isOK;
    			if (sons[j]!=-1)
    			{
    				ptr=root+sons[j];
    				isOK=true;
    				for (k=0; k<3; k++)
    				{
    					//this will count the face possibly in multiple nodes
    					if (triangle[k]<ptr->vStartIndex||triangle[k]>=ptr->vStartIndex+ptr->vNumber)
    					{
    						isOK=false;
    						break;
    					}
    				}
    				//nodeID=j;
    				if (isOK)
    				{
    					if (WriteFile(tempFileHandle[j], triangle, sizeof(UIVector),
    							&length, NULL)==0)
    					{
    						printf("write temp file failed with error code %d\n", GetLastError());
    						exit(4);
    					}
    					tempFaceCounter[j]++;
    				}
    						//don't count twice
    			}
    		}
    	}
    	//so, we have two kinds of face number "inflation". the following is our sibling at same
    	//level. The below is our children nodes. So, here we don't have to update face number
    	//because it is still not right
    	//first we need to update the "face number" because some face is overlapping across nodes

    	//we calculate from parent starting index
    	int currentFaceCount=0;

    	fPageCount=0;
    	fStartPageIndex=-1;
    	for (i=0; i<8; i++)
    	{
    		if (tempFaceCounter[i]>0)
    		{
    			Octree* ptr;
    			ptr=root+sons[i];
    			ptr->reserved1=tempFileHandle[i];
    			ptr->fStartIndex=fStartIndex+currentFaceCount;
    			////////////////////////////////////
    			//calculate the start offset of face

    			///////////////////////////////////////
    			ptr->fNumber=tempFaceCounter[i];
    			//the son node will update inflated face number itself
    			currentFaceCount+=ptr->doTraverseOctree(outFaceHandle);
    			if (fStartPageIndex==-1)
    			{
    				fStartPageIndex=root[sons[i]].fStartPageIndex;
    			}
    			fPageCount+=root[sons[i]].fPageCount;
    		}
    		//this is of course the safe point to delete and close file
    		//and in both case, we need to close file.
    		if (CloseHandle(tempFileHandle[i])==0)
    		{
    			printf("close handle failed\n");
    			exit(9);
    		}
    		if (DeleteFile(tempFileName[i])==0)
    		{
    			printf("delete file failed\n");
    			exit(83);
    		}
    	}
    	//here update the total face number
    	fNumber=currentFaceCount;
    	//reporting inflated face number to parent for updating
    	return fNumber;
    }

    void Octree::doCreateOctree(HANDLE outVertexHandle)
    {
    	HANDLE tempFileHandle[8], tempIndexFileHandle[8];
    	char tempFileName[8][MAX_PATH], tempIndexFileName[8][MAX_PATH];
    	int i, j;
    	unsigned long length;
    	int currentIndex;
    	int tempVertexCounter[8];
    	PlyVertex plyVertex;
    	int nodeID;
    	HANDLE vertexHandle, indexHandle;

    	vertexHandle=(HANDLE)reserved1;
    	indexHandle=(HANDLE)reserved2;
    	//remember to read from beginning
    	SetFilePointer(vertexHandle, 0, NULL, FILE_BEGIN);
    	if (indexHandle!=INVALID_HANDLE_VALUE)
    	{
    		SetFilePointer(indexHandle, 0, NULL, FILE_BEGIN);
    	}

    	if (vNumber<=MemoryThreshold)
    	{
    		//////////////////////////////////////////
    		//here we want to make sure each node starts alignment with page index
    		//so, we move file pointer so that it is

    		vStartPageIndex=currentVPageIndex;
    		vPageCount=vNumber*sizeof(PlyVertex)*7/5/PageSize+1;
    		//still we want to allocate with 40% of space for extra vertex
    		currentVPageIndex+=vPageCount;

    		//because of writing, we move file pointer
    		SetFilePointer(outVertexHandle, vStartPageIndex*PageSize, NULL, FILE_BEGIN);
    		///////////////////////////////////////////////////////////////////////////

    		for (j=0; j<vNumber; j++)
    		{
    			if (ReadFile(vertexHandle, &plyVertex, sizeof(PlyVertex), &length, NULL)==0)
    			{
    				printf("read temp vertex file failed\n");
    				exit(11);
    			}
    			if (ReadFile(indexHandle, &currentIndex, sizeof(int), &length, NULL)==0)
    			{
    				printf("read temp index file failed\n");
    				exit(2);
    			}
    			hashTable[currentIndex]=vertexIndexCounter;

    			//////////////////////////////////////////
    			//here we want to make sure each node starts alignment with page index
    			if (WriteFile(outVertexHandle, &plyVertex, sizeof(PlyVertex), &length, NULL)==0)
    			{
    				printf("write output vertex file failed\n");
    				exit(11);
    			}
    			vertexIndexCounter++;
    		}
    		return;

    		//////////////////////////////////////////
    		//here we want to make sure each node starts alignment with page index
    		//so, we move file pointer so that it is
    	}

    	for (i=0; i<8; i++)
    	{
    		tempVertexCounter[i]=0;
    		if (GetTempFileName(".", "vertex", 0, tempFileName[i])==0)
    		{
    			printf("get temp file name failed\n");
    			exit(44);
    		}
    		if (GetTempFileName(".", "index", 0, tempIndexFileName[i])==0)
    		{
    			printf("get temp file name failed\n");
    			exit(44);
    		}
    		tempFileHandle[i]=CreateFile(tempFileName[i], GENERIC_READ|GENERIC_WRITE, 0, NULL,
    			OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL);
    		tempIndexFileHandle[i]=CreateFile(tempIndexFileName[i], GENERIC_READ|GENERIC_WRITE,
    			0, NULL, OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL);
    		if (tempFileHandle[i]==INVALID_HANDLE_VALUE
    			||tempIndexFileHandle[i]==INVALID_HANDLE_VALUE)
    		{
    			printf("open temp file failed %d\n", GetLastError());
    			exit(55);
    		}
    	}

    	for (i=0; i<vNumber; i++)
    	{
    		if (ReadFile(vertexHandle, &plyVertex, sizeof(PlyVertex), &length, NULL)==0)
    		{
    			printf("read file fails\n");
    			exit(2);
    		}
    		if (indexHandle!=INVALID_HANDLE_VALUE)
    		{
    			if (ReadFile(indexHandle, &currentIndex, sizeof(int), &length, NULL)==0)
    			{
    				printf("read index file failed\n");
    				exit(9);
    			}
    		}
    		else
    		{
    			currentIndex=i;
    		}

    		nodeID=relativePosition(plyVertex.coord, center);
    		if (WriteFile(tempFileHandle[nodeID], &plyVertex, sizeof(PlyVertex),
    			&length, NULL)==0)
    		{
    			printf("write temp file failed with error code %d\n", GetLastError());
    			exit(4);
    		}
    		if (WriteFile(tempIndexFileHandle[nodeID], &currentIndex, sizeof(int), &length,
    			NULL)==0)
    		{
    			printf("write index file failed code %d\n", GetLastError());
    			exit(3);
    		}
    		tempVertexCounter[nodeID]++;
    	}
    	//we calculate from parent starting index
    	int currentVertexIndex=vStartIndex;
    	vStartPageIndex=-1;
    	vPageCount=0;
    	for (i=0; i<8; i++)
    	{
    		int sonNodePos;
    		if (tempVertexCounter[i]>0)
    		{
    			sonNodePos=createNewOctreeNode(i);
    			root[sonNodePos].vNumber=tempVertexCounter[i];
    			root[sonNodePos].vStartIndex=currentVertexIndex;
    			currentVertexIndex+=tempVertexCounter[i];

    			root[sonNodePos].reserved1=tempFileHandle[i];
    			root[sonNodePos].reserved2=tempIndexFileHandle[i];
    			root[sonNodePos].doCreateOctree(outVertexHandle);
    			if (vStartPageIndex==-1)
    			{
    				vStartPageIndex=root[sonNodePos].vStartPageIndex;
    			}
    			vPageCount+=root[sonNodePos].vPageCount;
    			//here is just counting, we don't move file pointers
    			//after this recursion, the temp file will not be used any more
    		}

    		//this is of course the safe point to delete and close file
    		if (CloseHandle(tempFileHandle[i])==0||CloseHandle(tempIndexFileHandle[i])==0)
    		{
    			printf("close handle failed\n");
    			exit(9);
    		}
    		if (DeleteFile(tempFileName[i])==0||DeleteFile(tempIndexFileName[i])==0)
    		{
    			printf("delete file failed\n");
    			exit(81);
    		}
    	}

    	//for all interanl nodes, we have to calculate sum of page of sons
    }

    file name: rply.h  (downloaded and modified a bit)

    #ifndef PLY_H
    #define PLY_H
    /* ----------------------------------------------------------------------
     * RPly library, read/write PLY files
     * Diego Nehab, Princeton University
     * http://www.cs.princeton.edu/~diego/professional/rply
     *
     * This library is distributed under the MIT License. See notice
     * at the end of this file.
     * ---------------------------------------------------------------------- */

    #ifdef __cplusplus
    extern "C" {
    #endif

    #define RPLY_VERSION   "RPly 1.01"
    #define RPLY_COPYRIGHT "Copyright (C) 2003-2005 Diego Nehab"
    #define RPLY_AUTHORS   "Diego Nehab"

    /* ----------------------------------------------------------------------
     * Types
     * ---------------------------------------------------------------------- */
    /* structures are opaque */
    typedef struct t_ply_ *p_ply;
    typedef struct t_ply_element_ *p_ply_element;
    typedef struct t_ply_property_ *p_ply_property;
    typedef struct t_ply_argument_ *p_ply_argument;

    /* ply format mode type */
    typedef enum e_ply_storage_mode_ {
        PLY_BIG_ENDIAN,
        PLY_LITTLE_ENDIAN,
        PLY_ASCII,
        PLY_DEFAULT      /* has to be the last in enum */
    } e_ply_storage_mode; /* order matches ply_storage_mode_list */

    /* ply data type */
    typedef enum e_ply_type {
        PLY_INT8, PLY_UINT8, PLY_INT16, PLY_UINT16,
        PLY_INT32, PLY_UIN32, PLY_FLOAT32, PLY_FLOAT64,
        PLY_CHAR, PLY_UCHAR, PLY_SHORT, PLY_USHORT,
        PLY_INT, PLY_UINT, PLY_FLOAT, PLY_DOUBLE,
        PLY_LIST    /* has to be the last in enum */
    } e_ply_type;   /* order matches ply_type_list */

    /* ----------------------------------------------------------------------
     * Property reading callback prototype
     *
     * message: error message
     * ---------------------------------------------------------------------- */
    typedef void (*p_ply_error_cb)(const char *message);

    /* ----------------------------------------------------------------------
     * Opens a ply file for reading (fails if file is not a ply file)
     *
     * error_cb: error callback function
     * name: file name
     *
     * Returns 1 if successful, 0 otherwise
     * ---------------------------------------------------------------------- */
    p_ply ply_open(const char *name, p_ply_error_cb error_cb);

    /* ----------------------------------------------------------------------
     * Reads and parses the header of a ply file returned by ply_open
     *
     * ply: handle returned by ply_open
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_read_header(p_ply ply);

    /* ----------------------------------------------------------------------
     * Property reading callback prototype
     *
     * argument: parameters for property being processed when callback is called
     *
     * Returns 1 if should continue processing file, 0 if should abort.
     * ---------------------------------------------------------------------- */
    typedef int (*p_ply_read_cb)(p_ply_argument argument);

    /* ----------------------------------------------------------------------
     * Sets up callbacks for property reading after header was parsed
     *
     * ply: handle returned by ply_open
     * element_name: element where property is
     * property_name: property to associate element with
     * read_cb: function to be called for each property value
     * pdata/idata: user data that will be passed to callback
     *
     * Returns 0 if no element or no property in element, returns the
     * number of element instances otherwise.
     * ---------------------------------------------------------------------- */
    long ply_set_read_cb(p_ply ply, const char *element_name,
            const char *property_name, p_ply_read_cb read_cb,
            void *pdata, long idata);

    /* ----------------------------------------------------------------------
     * Returns information about the element originating a callback
     *
     * argument: handle to argument
     * element: receives a the element handle (if non-null)
     * instance_index: receives the index of the current element instance
     *     (if non-null)
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_get_argument_element(p_ply_argument argument,
            p_ply_element *element, long *instance_index);

    /* ----------------------------------------------------------------------
     * Returns information about the property originating a callback
     *
     * argument: handle to argument
     * property: receives the property handle (if non-null)
     * length: receives the number of values in this property (if non-null)
     * value_index: receives the index of current property value (if non-null)
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_get_argument_property(p_ply_argument argument,
            p_ply_property *property, long *length, long *value_index);

    /* ----------------------------------------------------------------------
     * Returns user data associated with callback
     *
     * pdata: receives a copy of user custom data pointer (if non-null)
     * idata: receives a copy of user custom data integer (if non-null)
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_get_argument_user_data(p_ply_argument argument, void **pdata,
            long *idata);

    /* ----------------------------------------------------------------------
     * Returns the value associated with a callback
     *
     * argument: handle to argument
     *
     * Returns the current data item
     * ---------------------------------------------------------------------- */
    double ply_get_argument_value(p_ply_argument argument);

    /* ----------------------------------------------------------------------
     * Reads all elements and properties calling the callbacks defined with
     * calls to ply_set_read_cb
     *
     * ply: handle returned by ply_open
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_read(p_ply ply);

    /* ----------------------------------------------------------------------
     * Iterates over all elements by returning the next element.
     * Call with NULL to return handle to first element.
     *
     * ply: handle returned by ply_open
     * last: handle of last element returned (NULL for first element)
     *
     * Returns element if successfull or NULL if no more elements
     * ---------------------------------------------------------------------- */
    p_ply_element ply_get_next_element(p_ply ply, p_ply_element last);

    /* ----------------------------------------------------------------------
     * Iterates over all comments by returning the next comment.
     * Call with NULL to return pointer to first comment.
     *
     * ply: handle returned by ply_open
     * last: pointer to last comment returned (NULL for first comment)
     *
     * Returns comment if successfull or NULL if no more comments
     * ---------------------------------------------------------------------- */
    const char *ply_get_next_comment(p_ply ply, const char *last);

    /* ----------------------------------------------------------------------
     * Iterates over all obj_infos by returning the next obj_info.
     * Call with NULL to return pointer to first obj_info.
     *
     * ply: handle returned by ply_open
     * last: pointer to last obj_info returned (NULL for first obj_info)
     *
     * Returns obj_info if successfull or NULL if no more obj_infos
     * ---------------------------------------------------------------------- */
    const char *ply_get_next_obj_info(p_ply ply, const char *last);

    /* ----------------------------------------------------------------------
     * Returns information about an element
     *
     * element: element of interest
     * name: receives a pointer to internal copy of element name (if non-null)
     * ninstances: receives the number of instances of this element (if non-null)
     *
     * Returns 1 if successfull or 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_get_element_info(p_ply_element element, const char** name,
            long *ninstances);

    /* ----------------------------------------------------------------------
     * Iterates over all properties by returning the next property.
     * Call with NULL to return handle to first property.
     *
     * element: handle of element with the properties of interest
     * last: handle of last property returned (NULL for first property)
     *
     * Returns element if successfull or NULL if no more properties
     * ---------------------------------------------------------------------- */
    p_ply_property ply_get_next_property(p_ply_element element,
            p_ply_property last);

    /* ----------------------------------------------------------------------
     * Returns information about a property
     *
     * property: handle to property of interest
     * name: receives a pointer to internal copy of property name (if non-null)
     * type: receives the property type (if non-null)
     * length_type: for list properties, receives the scalar type of
     *     the length field (if non-null)
     * value_type: for list properties, receives the scalar type of the value
     *     fields  (if non-null)
     *
     * Returns 1 if successfull or 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_get_property_info(p_ply_property property, const char** name,
            e_ply_type *type, e_ply_type *length_type, e_ply_type *value_type);

    /* ----------------------------------------------------------------------
     * Creates new ply file
     *
     * name: file name
     * storage_mode: file format mode
     *
     * Returns handle to ply file if successfull, NULL otherwise
     * ---------------------------------------------------------------------- */
    p_ply ply_create(const char *name, e_ply_storage_mode storage_mode,
            p_ply_error_cb error_cb);

    /* ----------------------------------------------------------------------
     * Adds a new element to the ply file created by ply_create
     *
     * ply: handle returned by ply_create
     * name: name of new element
     * ninstances: number of element of this time in file
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_add_element(p_ply ply, const char *name, long ninstances);

    /* ----------------------------------------------------------------------
     * Adds a new property to the last element added by ply_add_element
     *
     * ply: handle returned by ply_create
     * name: name of new property
     * type: property type
     * length_type: scalar type of length field of a list property
     * value_type: scalar type of value fields of a list property
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_add_property(p_ply ply, const char *name, e_ply_type type,
            e_ply_type length_type, e_ply_type value_type);

    /* ----------------------------------------------------------------------
     * Adds a new list property to the last element added by ply_add_element
     *
     * ply: handle returned by ply_create
     * name: name of new property
     * length_type: scalar type of length field of a list property
     * value_type: scalar type of value fields of a list property
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_add_list_property(p_ply ply, const char *name,
            e_ply_type length_type, e_ply_type value_type);

    /* ----------------------------------------------------------------------
     * Adds a new property to the last element added by ply_add_element
     *
     * ply: handle returned by ply_create
     * name: name of new property
     * type: property type
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_add_scalar_property(p_ply ply, const char *name, e_ply_type type);

    /* ----------------------------------------------------------------------
     * Adds a new comment item
     *
     * ply: handle returned by ply_create
     * comment: pointer to string with comment text
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_add_comment(p_ply ply, const char *comment);

    /* ----------------------------------------------------------------------
     * Adds a new obj_info item
     *
     * ply: handle returned by ply_create
     * comment: pointer to string with obj_info data
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_add_obj_info(p_ply ply, const char *obj_info);

    /* ----------------------------------------------------------------------
     * Writes the ply file header after all element and properties have been
     * defined by calls to ply_add_element and ply_add_property
     *
     * ply: handle returned by ply_create
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_write_header(p_ply ply);

    /* ----------------------------------------------------------------------
     * Writes one property value, in the order they should be written to the
     * file. For each element type, write all elements of that type in order.
     * For each element, write all its properties in order. For scalar
     * properties, just write the value. For list properties, write the length
     * and then each of the values.
     *
     * ply: handle returned by ply_create
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_write(p_ply ply, double value);

    /* ----------------------------------------------------------------------
     * Closes a ply file handle. Releases all memory used by handle
     *
     * ply: handle to be closed.
     *
     * Returns 1 if successfull, 0 otherwise
     * ---------------------------------------------------------------------- */
    int ply_close(p_ply ply);

    #ifdef __cplusplus
    }
    #endif

    #endif /* RPLY_H */

    /* ----------------------------------------------------------------------
     * Copyright (C) 2003-2005 Diego Nehab. All rights reserved.
     *
     * Permission is hereby granted, free of charge, to any person obtaining
     * a copy of this software and associated documentation files (the
     * "Software"), to deal in the Software without restriction, including
     * without limitation the rights to use, copy, modify, merge, publish,
     * distribute, sublicense, and/or sell copies of the Software, and to
     * permit persons to whom the Software is furnished to do so, subject to
     * the following conditions:
     *
     * The above copyright notice and this permission notice shall be
     * included in all copies or substantial portions of the Software.
     *
     * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
     * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
     * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
     * IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
     * CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
     * TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
     * SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
     * ---------------------------------------------------------------------- */

    ¡¡

    ¡¡

    file name: rply.c (downloaded and modified a bit)

    /* ----------------------------------------------------------------------
     * RPly library, read/write PLY files
     * Diego Nehab, Princeton University
     * http://www.cs.princeton.edu/~diego/professional/rply
     *
     * This library is distributed under the MIT License. See notice
     * at the end of this file.
     * ---------------------------------------------------------------------- */
    #include <stdio.h>
    #include <ctype.h>
    #include <assert.h>
    #include <string.h>
    #include <limits.h>
    #include <float.h>
    #include <stdarg.h>
    #include <stdlib.h>
    #include <stddef.h>

    #include "rply.h"

    /* ----------------------------------------------------------------------
     * Constants
     * ---------------------------------------------------------------------- */
    #define WORDSIZE 256
    #define LINESIZE 1024
    #define BUFFERSIZE (8*1024)

    typedef enum e_ply_io_mode_ {
        PLY_READ,
        PLY_WRITE
    } e_ply_io_mode;

    static const char *const ply_storage_mode_list[] = {
        "binary_big_endian", "binary_little_endian", "ascii", NULL
    };     /* order matches e_ply_storage_mode enum */

    static const char *const ply_type_list[] = {
        "int8", "uint8", "int16", "uint16",
        "int32", "uint32", "float32", "float64",
        "char", "uchar", "short", "ushort",
        "int", "uint", "float", "double",
        "list", NULL
    };     /* order matches e_ply_type enum */

    /* ----------------------------------------------------------------------
     * Property reading callback argument
     *
     * element: name of element being processed
     * property: name of property being processed
     * nelements: number of elements of this kind in file
     * instance_index: index current element of this kind being processed
     * length: number of values in current list (or 1 for scalars)
     * value_index: index of current value int this list (or 0 for scalars)
     * value: value of property
     * pdata/idata: user data defined with ply_set_cb
     *
     * Returns handle to ply file if succesful, NULL otherwise.
     * ---------------------------------------------------------------------- */
    typedef struct t_ply_argument_ {
        p_ply_element element;
        long instance_index;
        p_ply_property property;
        long length, value_index;
        double value;
        void *pdata;
        long idata;
    } t_ply_argument;

    /* ----------------------------------------------------------------------
     * Property information
     *
     * name: name of this property
     * type: type of this property (list or type of scalar value)
     * length_type, value_type: type of list property count and values
     * read_cb: function to be called when this property is called
     *
     * Returns 1 if should continue processing file, 0 if should abort.
     * ---------------------------------------------------------------------- */
    typedef struct t_ply_property_ {
        char name[WORDSIZE];
        e_ply_type type, value_type, length_type;
        p_ply_read_cb read_cb;
        void *pdata;
        long idata;
    } t_ply_property;

    /* ----------------------------------------------------------------------
     * Element information
     *
     * name: name of this property
     * ninstances: number of elements of this type in file
     * property: property descriptions for this element
     * nproperty: number of properties in this element
     *
     * Returns 1 if should continue processing file, 0 if should abort.
     * ---------------------------------------------------------------------- */
    typedef struct t_ply_element_ {
        char name[WORDSIZE];
        long ninstances;
        p_ply_property property;
        long nproperties;
    } t_ply_element;

    /* ----------------------------------------------------------------------
     * Input/output driver
     *
     * Depending on file mode, different functions are used to read/write
     * property fields. The drivers make it transparent to read/write in ascii,
     * big endian or little endian cases.
     * ---------------------------------------------------------------------- */
    typedef int (*p_ply_ihandler)(p_ply ply, double *value);
    typedef int (*p_ply_ichunk)(p_ply ply, void *anydata, size_t size);
    typedef struct t_ply_idriver_ {
        p_ply_ihandler ihandler[16];
        p_ply_ichunk ichunk;
        const char *name;
    } t_ply_idriver;
    typedef t_ply_idriver *p_ply_idriver;

    typedef int (*p_ply_ohandler)(p_ply ply, double value);
    typedef int (*p_ply_ochunk)(p_ply ply, void *anydata, size_t size);
    typedef struct t_ply_odriver_ {
        p_ply_ohandler ohandler[16];
        p_ply_ochunk ochunk;
        const char *name;
    } t_ply_odriver;
    typedef t_ply_odriver *p_ply_odriver;

    /* ----------------------------------------------------------------------
     * Ply file handle.
     *
     * io_mode: read or write (from e_ply_io_mode)
     * storage_mode: mode of file associated with handle (from e_ply_storage_mode)
     * element: elements description for this file
     * nelement: number of different elements in file
     * comment: comments for this file
     * ncomments: number of comments in file
     * obj_info: obj_info items for this file
     * nobj_infos: number of obj_info items in file
     * fp: file pointer associated with ply file
     * c: last character read from ply file
     * buffer: last word/chunck of data read from ply file
     * buffer_first, buffer_last: interval of untouched good data in buffer
     * buffer_token: start of parsed token (line or word) in buffer
     * idriver, odriver: input driver used to get property fields from file
     * argument: storage space for callback arguments
     * welement, wproperty: element/property type being written
     * winstance_index: index of instance of current element being written
     * wvalue_index: index of list property value being written
     * wlength: number of values in list property being written
     * error_cb: callback for error messages
     * ---------------------------------------------------------------------- */
    typedef struct t_ply_ {
        e_ply_io_mode io_mode;
        e_ply_storage_mode storage_mode;
        p_ply_element element;
        long nelements;
        char *comment;
        long ncomments;
        char *obj_info;
        long nobj_infos;
        FILE *fp;
        int c;
        char buffer[BUFFERSIZE];
        size_t buffer_first, buffer_token, buffer_last;
        p_ply_idriver idriver;
        p_ply_odriver odriver;
        t_ply_argument argument;
        long welement, wproperty;
        long winstance_index, wvalue_index, wlength;
        p_ply_error_cb error_cb;
    } t_ply;

    /* ----------------------------------------------------------------------
     * I/O functions and drivers
     * ---------------------------------------------------------------------- */
    static t_ply_idriver ply_idriver_ascii;
    static t_ply_idriver ply_idriver_binary;
    static t_ply_idriver ply_idriver_binary_reverse;
    static t_ply_odriver ply_odriver_ascii;
    static t_ply_odriver ply_odriver_binary;
    static t_ply_odriver ply_odriver_binary_reverse;

    static int ply_read_word(p_ply ply);
    static int ply_check_word(p_ply ply);
    static int ply_read_line(p_ply ply);
    static int ply_check_line(p_ply ply);
    static int ply_read_chunk(p_ply ply, void *anybuffer, size_t size);
    static int ply_read_chunk_reverse(p_ply ply, void *anybuffer, size_t size);
    static int ply_write_chunk(p_ply ply, void *anybuffer, size_t size);
    static int ply_write_chunk_reverse(p_ply ply, void *anybuffer, size_t size);
    static void ply_reverse(void *anydata, size_t size);

    /* ----------------------------------------------------------------------
     * String functions
     * ---------------------------------------------------------------------- */
    static int ply_find_string(const char *item, const char* const list[]);
    static p_ply_element ply_find_element(p_ply ply, const char *name);
    static p_ply_property ply_find_property(p_ply_element element,
            const char *name);

    /* ----------------------------------------------------------------------
     * Header parsing
     * ---------------------------------------------------------------------- */
    static int ply_read_header_format(p_ply ply);
    static int ply_read_header_comment(p_ply ply);
    static int ply_read_header_obj_info(p_ply ply);
    static int ply_read_header_property(p_ply ply);
    static int ply_read_header_element(p_ply ply);

    /* ----------------------------------------------------------------------
     * Error handling
     * ---------------------------------------------------------------------- */
    static void ply_error_cb(const char *message);
    static void ply_error(p_ply ply, const char *fmt, ...);

    /* ----------------------------------------------------------------------
     * Memory allocation and initialization
     * ---------------------------------------------------------------------- */
    static void ply_init(p_ply ply);
    static void ply_element_init(p_ply_element element);
    static void ply_property_init(p_ply_property property);
    static p_ply ply_alloc(void);
    static p_ply_element ply_grow_element(p_ply ply);
    static p_ply_property ply_grow_property(p_ply ply, p_ply_element element);
    static void *ply_grow_array(p_ply ply, void **pointer, long *nmemb, long size);

    /* ----------------------------------------------------------------------
     * Special functions
     * ---------------------------------------------------------------------- */
    static e_ply_storage_mode ply_arch_endian(void);
    static int ply_type_check(void);

    /* ----------------------------------------------------------------------
     * Auxiliary read functions
     * ---------------------------------------------------------------------- */
    static int ply_read_element(p_ply ply, p_ply_element element,
            p_ply_argument argument);
    static int ply_read_property(p_ply ply, p_ply_element element,
            p_ply_property property, p_ply_argument argument);
    static int ply_read_list_property(p_ply ply, p_ply_element element,
            p_ply_property property, p_ply_argument argument);
    static int ply_read_scalar_property(p_ply ply, p_ply_element element,
            p_ply_property property, p_ply_argument argument);

    /* ----------------------------------------------------------------------
     * Buffer support functions
     * ---------------------------------------------------------------------- */
    /* pointers to tokenized word and line in buffer */
    #define BWORD(p) (p->buffer + p->buffer_token)
    #define BLINE(p) (p->buffer + p->buffer_token)

    /* pointer to start of untouched bytes in buffer */
    #define BFIRST(p) (p->buffer + p->buffer_first)

    /* number of bytes untouched in buffer */
    #define BSIZE(p) (p->buffer_last - p->buffer_first)

    /* consumes data from buffer */
    #define BSKIP(p, s) (p->buffer_first += s)

    /* refills the buffer */
    static int BREFILL(p_ply ply) {
        /* move untouched data to beginning of buffer */
        size_t size = BSIZE(ply);
        memmove(ply->buffer, BFIRST(ply), size);
        ply->buffer_last = size;
        ply->buffer_first = ply->buffer_token = 0;
        /* fill remaining with new data */
        size = fread(ply->buffer+size, 1, BUFFERSIZE-size-1, ply->fp);
        /* place sentinel so we can use str* functions with buffer */
        ply->buffer[BUFFERSIZE-1] = '\0';
        /* check if read failed */
        if (size <= 0) return 0;
        /* increase size to account for new data */
        ply->buffer_last += size;
        return 1;
    }

    /* ----------------------------------------------------------------------
     * Exported functions
     * ---------------------------------------------------------------------- */
    /* ----------------------------------------------------------------------
     * Read support functions
     * ---------------------------------------------------------------------- */
    p_ply ply_open(const char *name, p_ply_error_cb error_cb) {
        char magic[5] = "    ";
        FILE *fp = NULL;
        p_ply ply = NULL;
        if (error_cb == NULL) error_cb = ply_error_cb;
        if (!ply_type_check()) {
            error_cb("Incompatible type system");
            return NULL;
        }
        assert(name);
        fp = fopen(name, "rb");
        if (!fp) {
            error_cb("Unable to open file");
            return NULL;
        }
        if (fread(magic, 1, 4, fp) < 4) {
            error_cb("Error reading from file");
            fclose(fp);
            return NULL;
        }
        if (strcmp(magic, "ply\n")) {
            fclose(fp);
            error_cb("Not a PLY file. Expected magic number 'ply\\n'");
            return NULL;
        }
        ply = ply_alloc();
        if (!ply) {
            error_cb("Out of memory");
            fclose(fp);
            return NULL;
        }
        ply->fp = fp;
        ply->io_mode = PLY_READ;
        ply->error_cb = error_cb;
        return ply;
    }

    int ply_read_header(p_ply ply) {
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        if (!ply_read_word(ply)) return 0;
        /* parse file format */
        if (!ply_read_header_format(ply)) {
            ply_error(ply, "Invalid file format");
            return 0;
        }
        /* parse elements, comments or obj_infos until the end of header */
        while (strcmp(BWORD(ply), "end_header")) {
            if (!ply_read_header_comment(ply) &&
                    !ply_read_header_element(ply) &&
                    !ply_read_header_obj_info(ply)) {
                ply_error(ply, "Unexpected token '%s'", BWORD(ply));
                return 0;
            }
        }
        return 1;
    }

    long ply_set_read_cb(p_ply ply, const char *element_name,
            const char* property_name, p_ply_read_cb read_cb,
            void *pdata, long idata) {
        p_ply_element element = NULL;
        p_ply_property property = NULL;
        assert(ply && element_name && property_name);
        element = ply_find_element(ply, element_name);
        if (!element) return 0;
        property = ply_find_property(element, property_name);
        if (!property) return 0;
        property->read_cb = read_cb;
        property->pdata = pdata;
        property->idata = idata;
        return (int) element->ninstances;
    }

    int ply_read(p_ply ply) {
        long i;
        p_ply_argument argument;
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        argument = &ply->argument;
        /* for each element type */
        for (i = 0; i < ply->nelements; i++) {
            p_ply_element element = &ply->element[i];
            argument->element = element;
            if (!ply_read_element(ply, element, argument))
                return 0;
        }
        return 1;
    }

    /* ----------------------------------------------------------------------
     * Write support functions
     * ---------------------------------------------------------------------- */
    p_ply ply_create(const char *name, e_ply_storage_mode storage_mode,
            p_ply_error_cb error_cb) {
        FILE *fp = NULL;
        p_ply ply = NULL;
        if (error_cb == NULL) error_cb = ply_error_cb;
        if (!ply_type_check()) {
            error_cb("Incompatible type system");
            return NULL;
        }
        assert(name && storage_mode <= PLY_DEFAULT);
        fp = fopen(name, "wb");
        if (!fp) {
            error_cb("Unable to create file");
            return NULL;
        }
        ply = ply_alloc();
        if (!ply) {
            fclose(fp);
            error_cb("Out of memory");
            return NULL;
        }
        ply->io_mode = PLY_WRITE;
        if (storage_mode == PLY_DEFAULT) storage_mode = ply_arch_endian();
        if (storage_mode == PLY_ASCII) ply->odriver = &ply_odriver_ascii;
        else if (storage_mode == ply_arch_endian())
            ply->odriver = &ply_odriver_binary;
        else ply->odriver = &ply_odriver_binary_reverse;
        ply->storage_mode = storage_mode;
        ply->fp = fp;
        ply->error_cb = error_cb;
        return ply;
    }

    int ply_add_element(p_ply ply, const char *name, long ninstances) {
        p_ply_element element = NULL;
        assert(ply && ply->fp && ply->io_mode == PLY_WRITE);
        assert(name && strlen(name) < WORDSIZE && ninstances >= 0);
        if (strlen(name) >= WORDSIZE || ninstances < 0) {
            ply_error(ply, "Invalid arguments");
            return 0;
        }
        element = ply_grow_element(ply);
        if (!element) return 0;
        strcpy(element->name, name);
        element->ninstances = ninstances;
        return 1;
    }

    int ply_add_scalar_property(p_ply ply, const char *name, e_ply_type type) {
        p_ply_element element = NULL;
        p_ply_property property = NULL;
        assert(ply && ply->fp && ply->io_mode == PLY_WRITE);
        assert(name && strlen(name) < WORDSIZE);
        assert(type < PLY_LIST);
        if (strlen(name) >= WORDSIZE || type >= PLY_LIST) {
            ply_error(ply, "Invalid arguments");
            return 0;
        }
        element = &ply->element[ply->nelements-1];
        property = ply_grow_property(ply, element);
        if (!property) return 0;
        strcpy(property->name, name);
        property->type = type;
        return 1;
    }

    int ply_add_list_property(p_ply ply, const char *name,
            e_ply_type length_type, e_ply_type value_type) {
        p_ply_element element = NULL;
        p_ply_property property = NULL;
        assert(ply && ply->fp && ply->io_mode == PLY_WRITE);
        assert(name && strlen(name) < WORDSIZE);
        if (strlen(name) >= WORDSIZE) {
            ply_error(ply, "Invalid arguments");
            return 0;
        }
        assert(length_type < PLY_LIST);
        assert(value_type < PLY_LIST);
        if (length_type >= PLY_LIST || value_type >= PLY_LIST) {
            ply_error(ply, "Invalid arguments");
            return 0;
        }
        element = &ply->element[ply->nelements-1];
        property = ply_grow_property(ply, element);
        if (!property) return 0;
        strcpy(property->name, name);
        property->type = PLY_LIST;
        property->length_type = length_type;
        property->value_type = value_type;
        return 1;
    }

    int ply_add_property(p_ply ply, const char *name, e_ply_type type,
            e_ply_type length_type, e_ply_type value_type) {
        if (type == PLY_LIST)
            return ply_add_list_property(ply, name, length_type, value_type);
        else
            return ply_add_scalar_property(ply, name, type);
    }

    int ply_add_comment(p_ply ply, const char *comment) {
        char *new_comment = NULL;
        assert(ply && comment && strlen(comment) < LINESIZE);
        if (!comment || strlen(comment) >= LINESIZE) {
            ply_error(ply, "Invalid arguments");
            return 0;
        }
        new_comment = (char *) ply_grow_array(ply, (void **) &ply->comment,
                &ply->ncomments, LINESIZE);
        if (!new_comment) return 0;
        strcpy(new_comment, comment);
        return 1;
    }

    int ply_add_obj_info(p_ply ply, const char *obj_info) {
        char *new_obj_info = NULL;
        assert(ply && obj_info && strlen(obj_info) < LINESIZE);
        if (!obj_info || strlen(obj_info) >= LINESIZE) {
            ply_error(ply, "Invalid arguments");
            return 0;
        }
        new_obj_info = (char *) ply_grow_array(ply, (void **) &ply->obj_info,
                &ply->nobj_infos, LINESIZE);
        if (!new_obj_info) return 0;
        strcpy(new_obj_info, obj_info);
        return 1;
    }

    int ply_write_header(p_ply ply) {
        long i, j;
        assert(ply && ply->fp && ply->io_mode == PLY_WRITE);
        assert(ply->element || ply->nelements == 0);
        assert(!ply->element || ply->nelements > 0);
        if (fprintf(ply->fp, "ply\nformat %s 1.0\n",
                    ply_storage_mode_list[ply->storage_mode]) <= 0) goto error;
        for (i = 0; i < ply->ncomments; i++)
            if (fprintf(ply->fp, "comment %s\n", ply->comment + LINESIZE*i) <= 0)
                goto error;
        for (i = 0; i < ply->nobj_infos; i++)
            if (fprintf(ply->fp, "obj_info %s\n", ply->obj_info + LINESIZE*i) <= 0)
                goto error;
        for (i = 0; i < ply->nelements; i++) {
            p_ply_element element = &ply->element[i];
            assert(element->property || element->nproperties == 0);
            assert(!element->property || element->nproperties > 0);
            if (fprintf(ply->fp, "element %s %ld\n", element->name,
                        element->ninstances) <= 0) goto error;
            for (j = 0; j < element->nproperties; j++) {
                p_ply_property property = &element->property[j];
                if (property->type == PLY_LIST) {
                    if (fprintf(ply->fp, "property list %s %s %s\n",
                                ply_type_list[property->length_type],
                                ply_type_list[property->value_type],
                                property->name) <= 0) goto error;
                } else {
                    if (fprintf(ply->fp, "property %s %s\n",
                                ply_type_list[property->type],
                                property->name) <= 0) goto error;
                }
            }
        }
        return fprintf(ply->fp, "end_header\n") > 0;
    error:
        ply_error(ply, "Error writing to file");
        return 0;
    }

    int ply_write(p_ply ply, double value) {
        p_ply_element element = NULL;
        p_ply_property property = NULL;
        int type = -1;
        int breakafter = 0;
        if (ply->welement > ply->nelements) return 0;
        element = &ply->element[ply->welement];
        if (ply->wproperty > element->nproperties) return 0;
        property = &element->property[ply->wproperty];
        if (property->type == PLY_LIST) {
            if (ply->wvalue_index == 0) {
                type = property->length_type;
                ply->wlength = (long) value;
            } else type = property->value_type;
        } else {
            type = property->type;
            ply->wlength = 0;
        }
        if (!ply->odriver->ohandler[type](ply, value)) {
            ply_error(ply, "Failed writing %s of %s %d (%s: %s)",
                        property->name, element->name,
                        ply->winstance_index,
                        ply->odriver->name, ply_type_list[type]);
            return 0;
        }
        ply->wvalue_index++;
        if (ply->wvalue_index > ply->wlength) {
            ply->wvalue_index = 0;
            ply->wproperty++;
        }
        if (ply->wproperty >= element->nproperties) {
            ply->wproperty = 0;
            ply->winstance_index++;
            if (ply->storage_mode == PLY_ASCII) breakafter = 1;
        }
        if (ply->winstance_index >= element->ninstances) {
            ply->winstance_index = 0;
            ply->welement++;
        }
        return !breakafter || putc('\n', ply->fp) > 0;
    }

    int ply_close(p_ply ply) {
        long i;
        assert(ply && ply->fp);
        assert(ply->element || ply->nelements == 0);
        assert(!ply->element || ply->nelements > 0);
        /* write last chunk to file */
        if (ply->io_mode == PLY_WRITE &&
          fwrite(ply->buffer, 1, ply->buffer_last, ply->fp) < ply->buffer_last) {
            ply_error(ply, "Error closing up");
            return 0;
        }
        fclose(ply->fp);
        /* free all memory used by handle */
        if (ply->element) {
            for (i = 0; i < ply->nelements; i++) {
                p_ply_element element = &ply->element[i];
                if (element->property) free(element->property);
            }
            free(ply->element);
        }
        if (ply->obj_info) free(ply->obj_info);
        if (ply->comment) free(ply->comment);
        free(ply);
        return 1;
    }

    /* ----------------------------------------------------------------------
     * Query support functions
     * ---------------------------------------------------------------------- */
    p_ply_element ply_get_next_element(p_ply ply,
            p_ply_element last) {
        assert(ply);
        if (!last) return ply->element;
        last++;
        if (last < ply->element + ply->nelements) return last;
        else return NULL;
    }

    int ply_get_element_info(p_ply_element element, const char** name,
            long *ninstances) {
        assert(element);
        if (name) *name = element->name;
        if (ninstances) *ninstances = (long) element->ninstances;
        return 1;
    }

    p_ply_property ply_get_next_property(p_ply_element element,
            p_ply_property last) {
        assert(element);
        if (!last) return element->property;
        last++;
        if (last < element->property + element->nproperties) return last;
        else return NULL;
    }

    int ply_get_property_info(p_ply_property property, const char** name,
            e_ply_type *type, e_ply_type *length_type, e_ply_type *value_type) {
        assert(property);
        if (name) *name = property->name;
        if (type) *type = property->type;
        if (length_type) *length_type = property->length_type;
        if (value_type) *value_type = property->value_type;
        return 1;

    }

    const char *ply_get_next_comment(p_ply ply, const char *last) {
        assert(ply);
        if (!last) return ply->comment;
        last += LINESIZE;
        if (last < ply->comment + LINESIZE*ply->ncomments) return last;
        else return NULL;
    }

    const char *ply_get_next_obj_info(p_ply ply, const char *last) {
        assert(ply);
        if (!last) return ply->obj_info;
        last += LINESIZE;
        if (last < ply->obj_info + LINESIZE*ply->nobj_infos) return last;
        else return NULL;
    }

    /* ----------------------------------------------------------------------
     * Callback argument support functions
     * ---------------------------------------------------------------------- */
    int ply_get_argument_element(p_ply_argument argument,
            p_ply_element *element, long *instance_index) {
        assert(argument);
        if (!argument) return 0;
        if (element) *element = argument->element;
        if (instance_index) *instance_index = argument->instance_index;
        return 1;
    }

    int ply_get_argument_property(p_ply_argument argument,
            p_ply_property *property, long *length, long *value_index) {
        assert(argument);
        if (!argument) return 0;
        if (property) *property = argument->property;
        if (length) *length = argument->length;
        if (value_index) *value_index = argument->value_index;
        return 1;
    }

    int ply_get_argument_user_data(p_ply_argument argument, void **pdata,
            long *idata) {
        assert(argument);
        if (!argument) return 0;
        if (pdata) *pdata = argument->pdata;
        if (idata) *idata = argument->idata;
        return 1;
    }

    double ply_get_argument_value(p_ply_argument argument) {
        assert(argument);
        if (!argument) return 0.0;
        return argument->value;
    }

    /* ----------------------------------------------------------------------
     * Internal functions
     * ---------------------------------------------------------------------- */
    static int ply_read_list_property(p_ply ply, p_ply_element element,
            p_ply_property property, p_ply_argument argument) {
        int l;
        p_ply_read_cb read_cb = property->read_cb;
        p_ply_ihandler *driver = ply->idriver->ihandler;
        /* get list length */
        p_ply_ihandler handler = driver[property->length_type];
        double length;
        if (!handler(ply, &length)) {
            ply_error(ply, "Error reading '%s' of '%s' number %d",
                    property->name, element->name, argument->instance_index);
            return 0;
        }
        /* invoke callback to pass length in value field */
        argument->length = (long) length;
        argument->value_index = -1;
        argument->value = length;
        if (read_cb && !read_cb(argument)) {
            ply_error(ply, "Aborted by user");
            return 0;
        }
        /* read list values */
        handler = driver[property->value_type];
        /* for each value in list */
        for (l = 0; l < (long) length; l++) {
            /* read value from file */
            argument->value_index = l;
            if (!handler(ply, &argument->value)) {
                ply_error(ply, "Error reading value number %d of '%s' of "
                        "'%s' number %d", l+1, property->name,
                        element->name, argument->instance_index);
                return 0;
            }
            /* invoke callback to pass value */
            if (read_cb && !read_cb(argument)) {
                ply_error(ply, "Aborted by user");
                return 0;
            }
        }
        return 1;
    }

    static int ply_read_scalar_property(p_ply ply, p_ply_element element,
            p_ply_property property, p_ply_argument argument) {
        p_ply_read_cb read_cb = property->read_cb;
        p_ply_ihandler *driver = ply->idriver->ihandler;
        p_ply_ihandler handler = driver[property->type];
        argument->length = 1;
        argument->value_index = 0;
        if (!handler(ply, &argument->value)) {
            ply_error(ply, "Error reading '%s' of '%s' number %d",
                    property->name, element->name, argument->instance_index);
            return 0;
        }
        if (read_cb && !read_cb(argument)) {
            ply_error(ply, "Aborted by user");
            return 0;
        }
        return 1;
    }

    static int ply_read_property(p_ply ply, p_ply_element element,
            p_ply_property property, p_ply_argument argument) {
        if (property->type == PLY_LIST)
            return ply_read_list_property(ply, element, property, argument);
        else
            return ply_read_scalar_property(ply, element, property, argument);
    }

    static int ply_read_element(p_ply ply, p_ply_element element,
            p_ply_argument argument) {
        long j, k;
        /* for each element of this type */
        for (j = 0; j < element->ninstances; j++) {
            argument->instance_index = j;
            /* for each property */
            for (k = 0; k < element->nproperties; k++) {
                p_ply_property property = &element->property[k];
                argument->property = property;
                argument->pdata = property->pdata;
                argument->idata = property->idata;
                if (!ply_read_property(ply, element, property, argument))
                    return 0;
            }
        }
        return 1;
    }

    static int ply_find_string(const char *item, const char* const list[]) {
        int i;
        assert(item && list);
        for (i = 0; list[i]; i++)
            if (!strcmp(list[i], item)) return i;
        return -1;
    }

    static p_ply_element ply_find_element(p_ply ply, const char *name) {
        p_ply_element element;
        int i, nelements;
        assert(ply && name);
        element = ply->element;
        nelements = ply->nelements;
        assert(element || nelements == 0);
        assert(!element || nelements > 0);
        for (i = 0; i < nelements; i++)
            if (!strcmp(element[i].name, name)) return &element[i];
        return NULL;
    }

    static p_ply_property ply_find_property(p_ply_element element,
            const char *name) {
        p_ply_property property;
        int i, nproperties;
        assert(element && name);
        property = element->property;
        nproperties = element->nproperties;
        assert(property || nproperties == 0);
        assert(!property || nproperties > 0);
        for (i = 0; i < nproperties; i++)
            if (!strcmp(property[i].name, name)) return &property[i];
        return NULL;
    }

    static int ply_check_word(p_ply ply) {
        if (strlen(BLINE(ply)) >= WORDSIZE) {
            ply_error(ply, "Word too long");
            return 0;
        }
        return 1;
    }

    static int ply_read_word(p_ply ply) {
        size_t t = 0;
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        /* skip leading blanks */
        while (1) {
            t = strspn(BFIRST(ply), " \n\r\t");
            /* check if all buffer was made of blanks */
            if (t >= BSIZE(ply)) {
                if (!BREFILL(ply)) {
                    ply_error(ply, "Unexpected end of file");
                    return 0;
                }
            } else break;
        }
        BSKIP(ply, t);
        /* look for a space after the current word */
        t = strcspn(BFIRST(ply), " \n\r\t");
        /* if we didn't reach the end of the buffer, we are done */
        if (t < BSIZE(ply)) {
            ply->buffer_token = ply->buffer_first;
            BSKIP(ply, t);
            *BFIRST(ply) = '\0';
            BSKIP(ply, 1);
            return ply_check_word(ply);
        }
        /* otherwise, try to refill buffer */
        if (!BREFILL(ply)) {
            ply_error(ply, "Unexpected end of file");
            return 0;
        }
        /* keep looking from where we left */
        t += strcspn(BFIRST(ply) + t, " \n\r\t");
        /* check if the token is too large for our buffer */
        if (t >= BSIZE(ply)) {
            ply_error(ply, "Token too large");
            return 0;
        }
        /* we are done */
        ply->buffer_token = ply->buffer_first;
        BSKIP(ply, t);
        *BFIRST(ply) = '\0';
        BSKIP(ply, 1);
        return ply_check_word(ply);
    }

    static int ply_check_line(p_ply ply) {
        if (strlen(BLINE(ply)) >= LINESIZE) {
            ply_error(ply, "Line too long");
            return 0;
        }
        return 1;
    }

    static int ply_read_line(p_ply ply) {
        const char *end = NULL;
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        /* look for a end of line */
        end = strchr(BFIRST(ply), '\n');
        /* if we didn't reach the end of the buffer, we are done */
        if (end) {
            ply->buffer_token = ply->buffer_first;
            BSKIP(ply, end - BFIRST(ply));
            *BFIRST(ply) = '\0';
            BSKIP(ply, 1);
            return ply_check_line(ply);
        } else {
            end = ply->buffer + BSIZE(ply);
            /* otherwise, try to refill buffer */
            if (!BREFILL(ply)) {
                ply_error(ply, "Unexpected end of file");
                return 0;
            }
        }
        /* keep looking from where we left */
        end = strchr(end, '\n');
        /* check if the token is too large for our buffer */
        if (!end) {
            ply_error(ply, "Token too large");
            return 0;
        }
        /* we are done */
        ply->buffer_token = ply->buffer_first;
        BSKIP(ply, end - BFIRST(ply));
        *BFIRST(ply) = '\0';
        BSKIP(ply, 1);
        return ply_check_line(ply);
    }

    static int ply_read_chunk(p_ply ply, void *anybuffer, size_t size) {
        char *buffer = (char *) anybuffer;
        size_t i = 0;
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        assert(ply->buffer_first <= ply->buffer_last);
        while (i < size) {
            if (ply->buffer_first < ply->buffer_last) {
                buffer[i] = ply->buffer[ply->buffer_first];
                ply->buffer_first++;
                i++;
            } else {
                ply->buffer_first = 0;
                ply->buffer_last = fread(ply->buffer, 1, BUFFERSIZE, ply->fp);
                if (ply->buffer_last <= 0) return 0;
            }
        }
        return 1;
    }

    static int ply_write_chunk(p_ply ply, void *anybuffer, size_t size) {
        char *buffer = (char *) anybuffer;
        size_t i = 0;
        assert(ply && ply->fp && ply->io_mode == PLY_WRITE);
        assert(ply->buffer_last <= BUFFERSIZE);
        while (i < size) {
            if (ply->buffer_last < BUFFERSIZE) {
                ply->buffer[ply->buffer_last] = buffer[i];
                ply->buffer_last++;
                i++;
            } else {
                ply->buffer_last = 0;
                if (fwrite(ply->buffer, 1, BUFFERSIZE, ply->fp) < BUFFERSIZE)
                    return 0;
            }
        }
        return 1;
    }

    static int ply_write_chunk_reverse(p_ply ply, void *anybuffer, size_t size) {
        int ret = 0;
        ply_reverse(anybuffer, size);
        ret = ply_write_chunk(ply, anybuffer, size);
        ply_reverse(anybuffer, size);
        return ret;
    }

    static int ply_read_chunk_reverse(p_ply ply, void *anybuffer, size_t size) {
        if (!ply_read_chunk(ply, anybuffer, size)) return 0;
        ply_reverse(anybuffer, size);
        return 1;
    }

    static void ply_reverse(void *anydata, size_t size) {
        char *data = (char *) anydata;
        char temp;
        size_t i;
        for (i = 0; i < size/2; i++) {
            temp = data[i];
            data[i] = data[size-i-1];
            data[size-i-1] = temp;
        }
    }

    static void ply_init(p_ply ply) {
        ply->c = ' ';
        ply->element = NULL;
        ply->nelements = 0;
        ply->comment = NULL;
        ply->ncomments = 0;
        ply->obj_info = NULL;
        ply->nobj_infos = 0;
        ply->idriver = NULL;
        ply->odriver = NULL;
        ply->buffer[0] = '\0';
        ply->buffer_first = ply->buffer_last = ply->buffer_token = 0;
        ply->welement = 0;
        ply->wproperty = 0;
        ply->winstance_index = 0;
        ply->wlength = 0;
        ply->wvalue_index = 0;
    }

    static void ply_element_init(p_ply_element element) {
        element->name[0] = '\0';
        element->ninstances = 0;
        element->property = NULL;
        element->nproperties = 0;
    }

    static void ply_property_init(p_ply_property property) {
        property->name[0] = '\0';
        property->type = -1;
        property->length_type = -1;
        property->value_type = -1;
        property->read_cb = (p_ply_read_cb) NULL;
        property->pdata = NULL;
        property->idata = 0;
    }

    static p_ply ply_alloc(void) {
        p_ply ply = (p_ply) malloc(sizeof(t_ply));
        if (!ply) return NULL;
        ply_init(ply);
        return ply;
    }

    static void *ply_grow_array(p_ply ply, void **pointer,
            long *nmemb, long size) {
        void *temp = *pointer;
        long count = *nmemb + 1;
        if (!temp) temp = malloc(count*size);
        else temp = realloc(temp, count*size);
        if (!temp) {
            ply_error(ply, "Out of memory");
            return NULL;
        }
        *pointer = temp;
        *nmemb = count;
        return (char *) temp + (count-1) * size;
    }

    static p_ply_element ply_grow_element(p_ply ply) {
        p_ply_element element = NULL;
        assert(ply);
        assert(ply->element || ply->nelements == 0);
        assert(!ply->element || ply->nelements > 0);
        element = (p_ply_element) ply_grow_array(ply, (void **) &ply->element,
                &ply->nelements, sizeof(t_ply_element));
        if (!element) return NULL;
        ply_element_init(element);
        return element;
    }

    static p_ply_property ply_grow_property(p_ply ply, p_ply_element element) {
        p_ply_property property = NULL;
        assert(ply);
        assert(element);
        assert(element->property || element->nproperties == 0);
        assert(!element->property || element->nproperties > 0);
        property = (p_ply_property) ply_grow_array(ply,
                (void **) &element->property,
                &element->nproperties, sizeof(t_ply_property));
        if (!property) return NULL;
        ply_property_init(property);
        return property;
    }

    static int ply_read_header_format(p_ply ply) {
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        if (strcmp(BWORD(ply), "format")) return 0;
        if (!ply_read_word(ply)) return 0;
        ply->storage_mode = ply_find_string(BWORD(ply), ply_storage_mode_list);
        if (ply->storage_mode == (e_ply_storage_mode) (-1)) return 0;
        if (ply->storage_mode == PLY_ASCII) ply->idriver = &ply_idriver_ascii;
        else if (ply->storage_mode == ply_arch_endian())
            ply->idriver = &ply_idriver_binary;
        else ply->idriver = &ply_idriver_binary_reverse;
        if (!ply_read_word(ply)) return 0;
        if (strcmp(BWORD(ply), "1.0")) return 0;
        if (!ply_read_word(ply)) return 0;
        return 1;
    }

    static int ply_read_header_comment(p_ply ply) {
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        if (strcmp(BWORD(ply), "comment")) return 0;
        if (!ply_read_line(ply)) return 0;
        if (!ply_add_comment(ply, BLINE(ply))) return 0;
        if (!ply_read_word(ply)) return 0;
        return 1;
    }

    static int ply_read_header_obj_info(p_ply ply) {
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        if (strcmp(BWORD(ply), "obj_info")) return 0;
        if (!ply_read_line(ply)) return 0;
        if (!ply_add_obj_info(ply, BLINE(ply))) return 0;
        if (!ply_read_word(ply)) return 0;
        return 1;
    }

    static int ply_read_header_property(p_ply ply) {
        p_ply_element element = NULL;
        p_ply_property property = NULL;
        /* make sure it is a property */
        if (strcmp(BWORD(ply), "property")) return 0;
        element = &ply->element[ply->nelements-1];
        property = ply_grow_property(ply, element);
        if (!property) return 0;
        /* get property type */
        if (!ply_read_word(ply)) return 0;
        property->type = ply_find_string(BWORD(ply), ply_type_list);
        if (property->type == (e_ply_type) (-1)) return 0;
        if (property->type == PLY_LIST) {
            /* if it's a list, we need the base types */
            if (!ply_read_word(ply)) return 0;
            property->length_type = ply_find_string(BWORD(ply), ply_type_list);
            if (property->length_type == (e_ply_type) (-1)) return 0;
            if (!ply_read_word(ply)) return 0;
            property->value_type = ply_find_string(BWORD(ply), ply_type_list);
            if (property->value_type == (e_ply_type) (-1)) return 0;
        }
        /* get property name */
        if (!ply_read_word(ply)) return 0;
        strcpy(property->name, BWORD(ply));
        if (!ply_read_word(ply)) return 0;
        return 1;
    }

    static int ply_read_header_element(p_ply ply) {
        p_ply_element element = NULL;
        long dummy;
        assert(ply && ply->fp && ply->io_mode == PLY_READ);
        if (strcmp(BWORD(ply), "element")) return 0;
        /* allocate room for new element */
        element = ply_grow_element(ply);
        if (!element) return 0;
        /* get element name */
        if (!ply_read_word(ply)) return 0;
        strcpy(element->name, BWORD(ply));
        /* get number of elements of this type */
        if (!ply_read_word(ply)) return 0;
        if (sscanf(BWORD(ply), "%ld", &dummy) != 1) {
            ply_error(ply, "Expected number got '%s'", BWORD(ply));
            return 0;
        }
        element->ninstances = dummy;
        /* get all properties for this element */
        if (!ply_read_word(ply)) return 0;
        while (ply_read_header_property(ply) ||
            ply_read_header_comment(ply) || ply_read_header_obj_info(ply))
            /* do nothing */;
        return 1;
    }

    static void ply_error_cb(const char *message) {
        fprintf(stderr, "RPly: %s\n", message);
    }

    static void ply_error(p_ply ply, const char *fmt, ...) {
        char buffer[1024];
        va_list ap;
        va_start(ap, fmt);
        vsprintf(buffer, fmt, ap);
        va_end(ap);
        ply->error_cb(buffer);
    }

    static e_ply_storage_mode ply_arch_endian(void) {
        unsigned long i = 1;
        unsigned char *s = (unsigned char *) &i;
        if (*s == 1) return PLY_LITTLE_ENDIAN;
        else return PLY_BIG_ENDIAN;
    }

    static int ply_type_check(void) {
        assert(sizeof(char) == 1);
        assert(sizeof(unsigned char) == 1);
        assert(sizeof(short) == 2);
        assert(sizeof(unsigned short) == 2);
        assert(sizeof(long) == 4);
        assert(sizeof(unsigned long) == 4);
        assert(sizeof(float) == 4);
        assert(sizeof(double) == 8);
        if (sizeof(char) != 1) return 0;
        if (sizeof(unsigned char) != 1) return 0;
        if (sizeof(short) != 2) return 0;
        if (sizeof(unsigned short) != 2) return 0;
        if (sizeof(long) != 4) return 0;
        if (sizeof(unsigned long) != 4) return 0;
        if (sizeof(float) != 4) return 0;
        if (sizeof(double) != 8) return 0;
        return 1;
    }

    /* ----------------------------------------------------------------------
     * Output handlers
     * ---------------------------------------------------------------------- */
    static int oascii_int8(p_ply ply, double value) {
        if (value > CHAR_MAX || value < CHAR_MIN) return 0;
        return fprintf(ply->fp, "%d ", (char) value) > 0;
    }

    static int oascii_uint8(p_ply ply, double value) {
        if (value > UCHAR_MAX || value < 0) return 0;
        return fprintf(ply->fp, "%d ", (unsigned char) value) > 0;
    }

    static int oascii_int16(p_ply ply, double value) {
        if (value > SHRT_MAX || value < SHRT_MIN) return 0;
        return fprintf(ply->fp, "%d ", (short) value) > 0;
    }

    static int oascii_uint16(p_ply ply, double value) {
        if (value > USHRT_MAX || value < 0) return 0;
        return fprintf(ply->fp, "%d ", (unsigned short) value) > 0;
    }

    static int oascii_int32(p_ply ply, double value) {
        if (value > LONG_MAX || value < LONG_MIN) return 0;
        return fprintf(ply->fp, "%d ", (int) value) > 0;
    }

    static int oascii_uint32(p_ply ply, double value) {
        if (value > ULONG_MAX || value < 0) return 0;
        return fprintf(ply->fp, "%d ", (unsigned int) value) > 0;
    }

    static int oascii_float32(p_ply ply, double value) {
        if (value < -FLT_MAX || value > FLT_MAX) return 0;
        return fprintf(ply->fp, "%g ", (float) value) > 0;
    }

    static int oascii_float64(p_ply ply, double value) {
        if (value < -DBL_MAX || value > DBL_MAX) return 0;
        return fprintf(ply->fp, "%g ", value) > 0;
    }

    static int obinary_int8(p_ply ply, double value) {
        char int8 = (char) value;
        if (value > CHAR_MAX || value < CHAR_MIN) return 0;
        return ply->odriver->ochunk(ply, &int8, sizeof(int8));
    }

    static int obinary_uint8(p_ply ply, double value) {
        unsigned char uint8 = (unsigned char) value;
        if (value > UCHAR_MAX || value < 0) return 0;
        return ply->odriver->ochunk(ply, &uint8, sizeof(uint8));
    }

    static int obinary_int16(p_ply ply, double value) {
        short int16 = (short) value;
        if (value > SHRT_MAX || value < SHRT_MIN) return 0;
        return ply->odriver->ochunk(ply, &int16, sizeof(int16));
    }

    static int obinary_uint16(p_ply ply, double value) {
        unsigned short uint16 = (unsigned short) value;
        if (value > USHRT_MAX || value < 0) return 0;
        return ply->odriver->ochunk(ply, &uint16, sizeof(uint16));
    }

    static int obinary_int32(p_ply ply, double value) {
        long int32 = (long) value;
        if (value > LONG_MAX || value < LONG_MIN) return 0;
        return ply->odriver->ochunk(ply, &int32, sizeof(int32));
    }

    static int obinary_uint32(p_ply ply, double value) {
        unsigned long uint32 = (unsigned long) value;
        if (value > ULONG_MAX || value < 0) return 0;
        return ply->odriver->ochunk(ply, &uint32, sizeof(uint32));
    }

    static int obinary_float32(p_ply ply, double value) {
        float float32 = (float) value;
        if (value > FLT_MAX || value < -FLT_MAX) return 0;
        return ply->odriver->ochunk(ply, &float32, sizeof(float32));
    }

    static int obinary_float64(p_ply ply, double value) {
        return ply->odriver->ochunk(ply, &value, sizeof(value));
    }

    /* ----------------------------------------------------------------------
     * Input  handlers
     * ---------------------------------------------------------------------- */
    static int iascii_int8(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtol(BWORD(ply), &end, 10);
        if (*end || *value > CHAR_MAX || *value < CHAR_MIN) return 0;
        return 1;
    }

    static int iascii_uint8(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtol(BWORD(ply), &end, 10);
        if (*end || *value > UCHAR_MAX || *value < 0) return 0;
        return 1;
    }

    static int iascii_int16(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtol(BWORD(ply), &end, 10);
        if (*end || *value > SHRT_MAX || *value < SHRT_MIN) return 0;
        return 1;
    }

    static int iascii_uint16(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtol(BWORD(ply), &end, 10);
        if (*end || *value > USHRT_MAX || *value < 0) return 0;
        return 1;
    }

    static int iascii_int32(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtol(BWORD(ply), &end, 10);
        if (*end || *value > LONG_MAX || *value < LONG_MIN) return 0;
        return 1;
    }

    static int iascii_uint32(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtol(BWORD(ply), &end, 10);
        if (*end || *value < 0) return 0;
        return 1;
    }

    static int iascii_float32(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtod(BWORD(ply), &end);
        if (*end || *value < -FLT_MAX || *value > FLT_MAX) return 0;
        return 1;
    }

    static int iascii_float64(p_ply ply, double *value) {
        char *end;
        if (!ply_read_word(ply)) return 0;
        *value = strtod(BWORD(ply), &end);
        if (*end || *value < -DBL_MAX || *value > DBL_MAX) return 0;
        return 1;
    }

    static int ibinary_int8(p_ply ply, double *value) {
        char int8;
        if (!ply->idriver->ichunk(ply, &int8, 1)) return 0;
        *value = int8;
        return 1;
    }

    static int ibinary_uint8(p_ply ply, double *value) {
        unsigned char uint8;
        if (!ply->idriver->ichunk(ply, &uint8, 1)) return 0;
        *value = uint8;
        return 1;
    }

    static int ibinary_int16(p_ply ply, double *value) {
        short int16;
        if (!ply->idriver->ichunk(ply, &int16, sizeof(int16))) return 0;
        *value = int16;
        return 1;
    }

    static int ibinary_uint16(p_ply ply, double *value) {
        unsigned short uint16;
        if (!ply->idriver->ichunk(ply, &uint16, sizeof(uint16))) return 0;
        *value = uint16;
        return 1;
    }

    static int ibinary_int32(p_ply ply, double *value) {
        long int32;
        if (!ply->idriver->ichunk(ply, &int32, sizeof(int32))) return 0;
        *value = int32;
        return 1;
    }

    static int ibinary_uint32(p_ply ply, double *value) {
        unsigned long uint32;
        if (!ply->idriver->ichunk(ply, &uint32, sizeof(uint32))) return 0;
        *value = uint32;
        return 1;
    }

    static int ibinary_float32(p_ply ply, double *value) {
        float float32;
        if (!ply->idriver->ichunk(ply, &float32, sizeof(float32))) return 0;
        *value = float32;
        ply_reverse(&float32, sizeof(float32));
        return 1;
    }

    static int ibinary_float64(p_ply ply, double *value) {
        return ply->idriver->ichunk(ply, value, sizeof(double));
    }

    /* ----------------------------------------------------------------------
     * Constants
     * ---------------------------------------------------------------------- */
    static t_ply_idriver ply_idriver_ascii = {
        {   iascii_int8, iascii_uint8, iascii_int16, iascii_uint16,
            iascii_int32, iascii_uint32, iascii_float32, iascii_float64,
            iascii_int8, iascii_uint8, iascii_int16, iascii_uint16,
            iascii_int32, iascii_uint32, iascii_float32, iascii_float64
        }, /* order matches e_ply_type enum */
        NULL,
        "ascii input"
    };

    static t_ply_idriver ply_idriver_binary = {
        {   ibinary_int8, ibinary_uint8, ibinary_int16, ibinary_uint16,
            ibinary_int32, ibinary_uint32, ibinary_float32, ibinary_float64,
            ibinary_int8, ibinary_uint8, ibinary_int16, ibinary_uint16,
            ibinary_int32, ibinary_uint32, ibinary_float32, ibinary_float64
        }, /* order matches e_ply_type enum */
        ply_read_chunk,
        "binary input"
    };

    static t_ply_idriver ply_idriver_binary_reverse = {
        {   ibinary_int8, ibinary_uint8, ibinary_int16, ibinary_uint16,
            ibinary_int32, ibinary_uint32, ibinary_float32, ibinary_float64,
            ibinary_int8, ibinary_uint8, ibinary_int16, ibinary_uint16,
            ibinary_int32, ibinary_uint32, ibinary_float32, ibinary_float64
        }, /* order matches e_ply_type enum */
        ply_read_chunk_reverse,
        "reverse binary input"
    };

    static t_ply_odriver ply_odriver_ascii = {
        {   oascii_int8, oascii_uint8, oascii_int16, oascii_uint16,
            oascii_int32, oascii_uint32, oascii_float32, oascii_float64,
            oascii_int8, oascii_uint8, oascii_int16, oascii_uint16,
            oascii_int32, oascii_uint32, oascii_float32, oascii_float64
        }, /* order matches e_ply_type enum */
        NULL,
        "ascii output"
    };

    static t_ply_odriver ply_odriver_binary = {
        {   obinary_int8, obinary_uint8, obinary_int16, obinary_uint16,
            obinary_int32, obinary_uint32, obinary_float32, obinary_float64,
            obinary_int8, obinary_uint8, obinary_int16, obinary_uint16,
            obinary_int32, obinary_uint32, obinary_float32, obinary_float64
        }, /* order matches e_ply_type enum */
        ply_write_chunk,
        "binary output"
    };

    static t_ply_odriver ply_odriver_binary_reverse = {
        {   obinary_int8, obinary_uint8, obinary_int16, obinary_uint16,
            obinary_int32, obinary_uint32, obinary_float32, obinary_float64,
            obinary_int8, obinary_uint8, obinary_int16, obinary_uint16,
            obinary_int32, obinary_uint32, obinary_float32, obinary_float64
        }, /* order matches e_ply_type enum */
        ply_write_chunk_reverse,
        "reverse binary output"
    };

    /* ----------------------------------------------------------------------
     * Copyright (C) 2003 Diego Nehab.  All rights reserved.
     *
     * Permission is hereby granted, free of charge, to any person obtaining
     * a copy of this software and associated documentation files (the
     * "Software"), to deal in the Software without restriction, including
     * without limitation the rights to use, copy, modify, merge, publish,
     * distribute, sublicense, and/or sell copies of the Software, and to
     * permit persons to whom the Software is furnished to do so, subject to
     * the following conditions:
     *
     * The above copyright notice and this permission notice shall be
     * included in all copies or substantial portions of the Software.
     *
     * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
     * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
     * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
     * IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
     * CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
     * TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
     * SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
     * ---------------------------------------------------------------------- */

    ¡¡

    ¡¡

    file name: frustum.h (downloaded and modified a bit)

    #ifndef _FRUSTUM_H
    #define _FRUSTUM_H

    // This will allow us to create an object to keep track of our frustum
    class CFrustum {

    public:

    	// Call this every time the camera moves to update the frustum
    	void CalculateFrustum();

    	// This takes a 3D point and returns TRUE if it's inside of the frustum
    	bool PointInFrustum(float x, float y, float z);

    	// This takes a 3D point and a radius and returns TRUE if the sphere is inside of the frustum
    	bool SphereInFrustum(float x, float y, float z, float radius);

    	// This takes the center and half the length of the cube.
    	int CubeInFrustum( float x, float y, float z, float size );

    private:

    	// This holds the A B C and D values for each side of our frustum.
    	float m_Frustum[6][4];
    };

    #endif //_FRUSTUM_H

    //
    // Ben Humphrey (DigiBen)
    // Game Programmer
    // DigiBen@GameTutorials.com
    // Co-Web Host of www.GameTutorials.com
    //
    //

    ¡¡

    file name: frustum.c (downloaded and modified a bit)

    //***********************************************************************//
    //																		 //
    //		- "Talk to me like I'm a 3 year old!" Programming Lessons -		 //
    //                                                                       //
    //		$Author:		DigiBen		digiben@gametutorials.com			 //
    //																		 //
    //		$Program:		Frustum Culling									 //
    //																		 //
    //		$Description:	Demonstrates checking if shapes are in view		 //
    //																		 //
    //		$Date:			8/28/01											 //
    //																		 //
    //***********************************************************************//

    										// Header File For The OpenGL32 Library

    #include <cstdio>
    #include <cstdlib>
    #include <cmath>
    #include <GL/glut.h>
    #include "Frustum.h"

    // We create an enum of the sides so we don't have to call each side 0 or 1.
    // This way it makes it more understandable and readable when dealing with frustum sides.
    enum FrustumSide
    {
    	RIGHT	= 0,		// The RIGHT side of the frustum
    	LEFT	= 1,		// The LEFT	 side of the frustum
    	BOTTOM	= 2,		// The BOTTOM side of the frustum
    	TOP		= 3,		// The TOP side of the frustum
    	BACK	= 4,		// The BACK	side of the frustum
    	FRONT	= 5			// The FRONT side of the frustum
    };

    // Like above, instead of saying a number for the ABC and D of the plane, we
    // want to be more descriptive.
    enum PlaneData
    {
    	A = 0,				// The X value of the plane's normal
    	B = 1,				// The Y value of the plane's normal
    	C = 2,				// The Z value of the plane's normal
    	D = 3				// The distance the plane is from the origin
    };

    ///////////////////////////////// NORMALIZE PLANE \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*
    /////
    /////	This normalizes a plane (A side) from a given frustum.
    /////
    ///////////////////////////////// NORMALIZE PLANE \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*

    void NormalizePlane(float frustum[6][4], int side)
    {
    	// Here we calculate the magnitude of the normal to the plane (point A B C)
    	// Remember that (A, B, C) is that same thing as the normal's (X, Y, Z).
    	// To calculate magnitude you use the equation:  magnitude = sqrt( x^2 + y^2 + z^2)
    	float magnitude = (float)sqrt( frustum[side][A] * frustum[side][A] +
    								   frustum[side][B] * frustum[side][B] +
    								   frustum[side][C] * frustum[side][C] );

    	// Then we divide the plane's values by it's magnitude.
    	// This makes it easier to work with.
    	frustum[side][A] /= magnitude;
    	frustum[side][B] /= magnitude;
    	frustum[side][C] /= magnitude;
    	frustum[side][D] /= magnitude;
    }

    ///////////////////////////////// CALCULATE FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*
    /////
    /////	This extracts our frustum from the projection and modelview matrix.
    /////
    ///////////////////////////////// CALCULATE FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*

    void CFrustum::CalculateFrustum()
    {
    	float   proj[16];								// This will hold our projection matrix
    	float   modl[16];								// This will hold our modelview matrix
    	float   clip[16];								// This will hold the clipping planes

    	// glGetFloatv() is used to extract information about our OpenGL world.
    	// Below, we pass in GL_PROJECTION_MATRIX to abstract our projection matrix.
    	// It then stores the matrix into an array of [16].
    	glGetFloatv( GL_PROJECTION_MATRIX, proj );

    	// By passing in GL_MODELVIEW_MATRIX, we can abstract our model view matrix.
    	// This also stores it in an array of [16].
    	glGetFloatv( GL_MODELVIEW_MATRIX, modl );

    	// Now that we have our modelview and projection matrix, if we combine these 2 matrices,
    	// it will give us our clipping planes.  To combine 2 matrices, we multiply them.

    	clip[ 0] = modl[ 0] * proj[ 0] + modl[ 1] * proj[ 4] + modl[ 2] * proj[ 8] + modl[ 3] * proj[12];
    	clip[ 1] = modl[ 0] * proj[ 1] + modl[ 1] * proj[ 5] + modl[ 2] * proj[ 9] + modl[ 3] * proj[13];
    	clip[ 2] = modl[ 0] * proj[ 2] + modl[ 1] * proj[ 6] + modl[ 2] * proj[10] + modl[ 3] * proj[14];
    	clip[ 3] = modl[ 0] * proj[ 3] + modl[ 1] * proj[ 7] + modl[ 2] * proj[11] + modl[ 3] * proj[15];

    	clip[ 4] = modl[ 4] * proj[ 0] + modl[ 5] * proj[ 4] + modl[ 6] * proj[ 8] + modl[ 7] * proj[12];
    	clip[ 5] = modl[ 4] * proj[ 1] + modl[ 5] * proj[ 5] + modl[ 6] * proj[ 9] + modl[ 7] * proj[13];
    	clip[ 6] = modl[ 4] * proj[ 2] + modl[ 5] * proj[ 6] + modl[ 6] * proj[10] + modl[ 7] * proj[14];
    	clip[ 7] = modl[ 4] * proj[ 3] + modl[ 5] * proj[ 7] + modl[ 6] * proj[11] + modl[ 7] * proj[15];

    	clip[ 8] = modl[ 8] * proj[ 0] + modl[ 9] * proj[ 4] + modl[10] * proj[ 8] + modl[11] * proj[12];
    	clip[ 9] = modl[ 8] * proj[ 1] + modl[ 9] * proj[ 5] + modl[10] * proj[ 9] + modl[11] * proj[13];
    	clip[10] = modl[ 8] * proj[ 2] + modl[ 9] * proj[ 6] + modl[10] * proj[10] + modl[11] * proj[14];
    	clip[11] = modl[ 8] * proj[ 3] + modl[ 9] * proj[ 7] + modl[10] * proj[11] + modl[11] * proj[15];

    	clip[12] = modl[12] * proj[ 0] + modl[13] * proj[ 4] + modl[14] * proj[ 8] + modl[15] * proj[12];
    	clip[13] = modl[12] * proj[ 1] + modl[13] * proj[ 5] + modl[14] * proj[ 9] + modl[15] * proj[13];
    	clip[14] = modl[12] * proj[ 2] + modl[13] * proj[ 6] + modl[14] * proj[10] + modl[15] * proj[14];
    	clip[15] = modl[12] * proj[ 3] + modl[13] * proj[ 7] + modl[14] * proj[11] + modl[15] * proj[15];

    	// Now we actually want to get the sides of the frustum.  To do this we take
    	// the clipping planes we received above and extract the sides from them.

    	// This will extract the RIGHT side of the frustum
    	m_Frustum[RIGHT][A] = clip[ 3] - clip[ 0];
    	m_Frustum[RIGHT][B] = clip[ 7] - clip[ 4];
    	m_Frustum[RIGHT][C] = clip[11] - clip[ 8];
    	m_Frustum[RIGHT][D] = clip[15] - clip[12];

    	// Now that we have a normal (A,B,C) and a distance (D) to the plane,
    	// we want to normalize that normal and distance.

    	// Normalize the RIGHT side
    	NormalizePlane(m_Frustum, RIGHT);

    	// This will extract the LEFT side of the frustum
    	m_Frustum[LEFT][A] = clip[ 3] + clip[ 0];
    	m_Frustum[LEFT][B] = clip[ 7] + clip[ 4];
    	m_Frustum[LEFT][C] = clip[11] + clip[ 8];
    	m_Frustum[LEFT][D] = clip[15] + clip[12];

    	// Normalize the LEFT side
    	NormalizePlane(m_Frustum, LEFT);

    	// This will extract the BOTTOM side of the frustum
    	m_Frustum[BOTTOM][A] = clip[ 3] + clip[ 1];
    	m_Frustum[BOTTOM][B] = clip[ 7] + clip[ 5];
    	m_Frustum[BOTTOM][C] = clip[11] + clip[ 9];
    	m_Frustum[BOTTOM][D] = clip[15] + clip[13];

    	// Normalize the BOTTOM side
    	NormalizePlane(m_Frustum, BOTTOM);

    	// This will extract the TOP side of the frustum
    	m_Frustum[TOP][A] = clip[ 3] - clip[ 1];
    	m_Frustum[TOP][B] = clip[ 7] - clip[ 5];
    	m_Frustum[TOP][C] = clip[11] - clip[ 9];
    	m_Frustum[TOP][D] = clip[15] - clip[13];

    	// Normalize the TOP side
    	NormalizePlane(m_Frustum, TOP);

    	// This will extract the BACK side of the frustum
    	m_Frustum[BACK][A] = clip[ 3] - clip[ 2];
    	m_Frustum[BACK][B] = clip[ 7] - clip[ 6];
    	m_Frustum[BACK][C] = clip[11] - clip[10];
    	m_Frustum[BACK][D] = clip[15] - clip[14];

    	// Normalize the BACK side
    	NormalizePlane(m_Frustum, BACK);

    	// This will extract the FRONT side of the frustum
    	m_Frustum[FRONT][A] = clip[ 3] + clip[ 2];
    	m_Frustum[FRONT][B] = clip[ 7] + clip[ 6];
    	m_Frustum[FRONT][C] = clip[11] + clip[10];
    	m_Frustum[FRONT][D] = clip[15] + clip[14];

    	// Normalize the FRONT side
    	NormalizePlane(m_Frustum, FRONT);
    }

    // The code below will allow us to make checks within the frustum.  For example,
    // if we want to see if a point, a sphere, or a cube lies inside of the frustum.
    // Because all of our planes point INWARDS (The normals are all pointing inside the frustum)
    // we then can assume that if a point is in FRONT of all of the planes, it's inside.

    ///////////////////////////////// POINT IN FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*
    /////
    /////	This determines if a point is inside of the frustum
    /////
    ///////////////////////////////// POINT IN FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*

    bool CFrustum::PointInFrustum( float x, float y, float z )
    {
    	// If you remember the plane equation (A*x + B*y + C*z + D = 0), then the rest
    	// of this code should be quite obvious and easy to figure out yourself.
    	// In case don't know the plane equation, it might be a good idea to look
    	// at our Plane Collision tutorial at www.GameTutorials.com in OpenGL Tutorials.
    	// I will briefly go over it here.  (A,B,C) is the (X,Y,Z) of the normal to the plane.
    	// They are the same thing... but just called ABC because you don't want to say:
    	// (x*x + y*y + z*z + d = 0).  That would be wrong, so they substitute them.
    	// the (x, y, z) in the equation is the point that you are testing.  The D is
    	// The distance the plane is from the origin.  The equation ends with "= 0" because
    	// that is true when the point (x, y, z) is ON the plane.  When the point is NOT on
    	// the plane, it is either a negative number (the point is behind the plane) or a
    	// positive number (the point is in front of the plane).  We want to check if the point
    	// is in front of the plane, so all we have to do is go through each point and make
    	// sure the plane equation goes out to a positive number on each side of the frustum.
    	// The result (be it positive or negative) is the distance the point is front the plane.

    	// Go through all the sides of the frustum
    	for(int i = 0; i < 6; i++ )
    	{
    		// Calculate the plane equation and check if the point is behind a side of the frustum
    		if(m_Frustum[i][A] * x + m_Frustum[i][B] * y + m_Frustum[i][C] * z + m_Frustum[i][D] <= 0)
    		{
    			// The point was behind a side, so it ISN'T in the frustum
    			return false;
    		}
    	}

    	// The point was inside of the frustum (In front of ALL the sides of the frustum)
    	return true;
    }

    ///////////////////////////////// SPHERE IN FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*
    /////
    /////	This determines if a sphere is inside of our frustum by it's center and radius.
    /////
    ///////////////////////////////// SPHERE IN FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*

    bool CFrustum::SphereInFrustum( float x, float y, float z, float radius )
    {
    	// Now this function is almost identical to the PointInFrustum(), except we
    	// now have to deal with a radius around the point.  The point is the center of
    	// the radius.  So, the point might be outside of the frustum, but it doesn't
    	// mean that the rest of the sphere is.  It could be half and half.  So instead of
    	// checking if it's less than 0, we need to add on the radius to that.  Say the
    	// equation produced -2, which means the center of the sphere is the distance of
    	// 2 behind the plane.  Well, what if the radius was 5?  The sphere is still inside,
    	// so we would say, if(-2 < -5) then we are outside.  In that case it's false,
    	// so we are inside of the frustum, but a distance of 3.  This is reflected below.

    	// Go through all the sides of the frustum
    	for(int i = 0; i < 6; i++ )
    	{
    		// If the center of the sphere is farther away from the plane than the radius
    		if( m_Frustum[i][A] * x + m_Frustum[i][B] * y + m_Frustum[i][C] * z + m_Frustum[i][D] <= -radius )
    		{
    			// The distance was greater than the radius so the sphere is outside of the frustum
    			return false;
    		}
    	}

    	// The sphere was inside of the frustum!
    	return true;
    }

    ///////////////////////////////// CUBE IN FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*
    /////
    /////	This determines if a cube is in or around our frustum by it's center and 1/2 it's length
    /////
    ///////////////////////////////// CUBE IN FRUSTUM \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*
    ///////////////////////////////////////////////////////////////////////////////////
    //Let's define the return value like this:
    // 1 means the cube is completely inside frustum
    // 0 means the cube is intersects with frustum
    // -1 means it is completely outside frustum
    //the code below seems a bit clumnsy, but I think it should be the better way of coding.

    int CFrustum::CubeInFrustum( float x, float y, float z, float size )
    {
    	// This test is a bit more work, but not too much more complicated.
    	// Basically, what is going on is, that we are given the center of the cube,
    	// and half the length.  Think of it like a radius.  Then we checking each point
    	// in the cube and seeing if it is inside the frustum.  If a point is found in front
    	// of a side, then we skip to the next side.  If we get to a plane that does NOT have
    	// a point in front of it, then it will return false.

    	// *Note* - This will sometimes say that a cube is inside the frustum when it isn't.
    	// This happens when all the corners of the bounding box are not behind any one plane.
    	// This is rare and shouldn't effect the overall rendering speed.
    	bool isInside=false, isOutside=false;
    	for(int i = 0; i < 6; i++ )
    		// A violation of all the six conditions means its completely inside the frustum. @RG (I hope so :))
    	{
    		if(m_Frustum[i][A] * (x - size) + m_Frustum[i][B] * (y - size) + m_Frustum[i][C] * (z - size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}
    		//////////////////////////////////////////////////////

    		if(m_Frustum[i][A] * (x + size) + m_Frustum[i][B] * (y - size) + m_Frustum[i][C] * (z - size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}

    		if(m_Frustum[i][A] * (x - size) + m_Frustum[i][B] * (y + size) + m_Frustum[i][C] * (z - size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}

    		if(m_Frustum[i][A] * (x + size) + m_Frustum[i][B] * (y + size) + m_Frustum[i][C] * (z - size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}

    		if(m_Frustum[i][A] * (x - size) + m_Frustum[i][B] * (y - size) + m_Frustum[i][C] * (z + size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}

    		if(m_Frustum[i][A] * (x + size) + m_Frustum[i][B] * (y - size) + m_Frustum[i][C] * (z + size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}

    		if(m_Frustum[i][A] * (x - size) + m_Frustum[i][B] * (y + size) + m_Frustum[i][C] * (z + size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}

    		if(m_Frustum[i][A] * (x + size) + m_Frustum[i][B] * (y + size) + m_Frustum[i][C] * (z + size) + m_Frustum[i][D] > 0)
    		{
    			isInside=true;
    		}
    		else
    		{
    			isOutside=true;
    		}

    		if (isOutside&&isInside)
    		{
    			return 0; //intersection
    		}
    	}

    	//if we reach here it means they are not same
    	if (isInside)
    	{
    		return 1;
    	}

    	if (isOutside)
    	{
    		return -1;
    	}

    }

    /////////////////////////////////////////////////////////////////////////////////
    //
    // * QUICK NOTES *
    //
    // WOZZERS!  That seemed like an incredible amount to look at, but if you break it
    // down, it's not.  Frustum culling is a VERY useful thing when it comes to 3D.
    // If you want a large world, there is no way you are going to send it down the
    // 3D pipeline every frame and let OpenGL take care of it for you.  That would
    // give you a 0.001 frame rate.  If you hit '+' and bring the sphere count up to
    // 1000, then take off culling, you will see the HUGE difference it makes.
    // Also, you wouldn't really be rendering 1000 spheres.  You would most likely
    // use the sphere code for larger objects.  Let me explain.  Say you have a bunch
    // of objects, well... all you need to do is give the objects a radius, and then
    // test that radius against the frustum.  If that sphere is in the frustum, then you
    // render that object.  Also, you won't be rendering a high poly sphere so it won't
    // be so slow.  This goes for bounding box's too (Cubes).  If you don't want to
    // do a cube, it is really easy to convert the code for rectangles.  Just pass in
    // a width and height, instead of just a length.  Remember, it's HALF the length of
    // the cube, not the full length.  So it would be half the width and height for a rect.
    //
    // This is a perfect starter for an octree tutorial.  Wrap you head around the concepts
    // here and then see if you can apply this to making an octree.  Hopefully we will have
    // a tutorial up and running for this subject soon.  Once you have frustum culling,
    // the next step is getting space partitioning.  Either it being a BSP tree of an Octree.
    //
    // Let's go over a brief overview of the things we learned here:
    //
    // 1) First we need to abstract the frustum from OpenGL.  To do that we need the
    //    projection and modelview matrix.  To get the projection matrix we use:
    //
    //			glGetFloatv( GL_PROJECTION_MATRIX, /* An Array of 16 floats */ );
    //    Then, to get the modelview matrix we use:
    //
    //			glGetFloatv( GL_MODELVIEW_MATRIX, /* An Array of 16 floats */ );
    //
    //	  These 2 functions gives us an array of 16 floats (The matrix).
    //
    // 2) Next, we need to combine these 2 matrices.  We do that by matrix multiplication.
    //
    // 3) Now that we have the 2 matrixes combined, we can abstract the sides of the frustum.
    //    This will give us the normal and the distance from the plane to the origin (ABC and D).
    //
    // 4) After abstracting a side, we want to normalize the plane data.  (A B C and D).
    //
    // 5) Now we have our frustum, and we can check points against it using the plane equation.
    //    Once again, the plane equation (A*x + B*y + C*z + D = 0) says that if, point (X,Y,Z)
    //    times the normal of the plane (A,B,C), plus the distance of the plane from origin,
    //    will equal 0 if the point (X, Y, Z) lies on that plane.  If it is behind the plane
    //    it will be a negative distance, if it's in front of the plane (the way the normal is facing)
    //    it will be a positive number.
    //
    //
    // If you need more help on the plane equation and why this works, download our
    // Ray Plane Intersection Tutorial at www.GameTutorials.com.
    //
    // That's pretty much it with frustums.  There is a lot more we could talk about, but
    // I don't want to complicate this tutorial more than I already have.
    //
    // I want to thank Mark Morley for his tutorial on frustum culling.  Most of everything I got
    // here comes from his teaching.  If you want more in-depth, visit his tutorial at:
    //
    // http://www.markmorley.com/opengl/frustumculling.html
    //
    // Good luck!
    //
    //
    // Ben Humphrey (DigiBen)
    // Game Programmer
    // DigiBen@GameTutorials.com
    // Co-Web Host of www.GameTutorials.com
    //
    //

    ¡¡

    file name: plyConverter.cpp (downloaded and modified a bit)

    #include <stdio.h>
    #include <vector>
    #include "rply.h"
    #include <cmath>
    #include <windows.h>
    #include <winbase.h>
    #include "PlyModel.h"

    using namespace std;

    PlyModel ply;

    static int callback(p_ply_argument argument)
    {
        void *pdata;
        /* just pass the value from the input file to the output file */
        ply_get_argument_user_data(argument, &pdata, NULL);
        ply_write((p_ply) pdata, ply_get_argument_value(argument));
        return 1;
    }

    static int setup_callbacks(p_ply iply, p_ply oply)
    {
        p_ply_element element = NULL;
        /* iterate over all elements in input file */
        while ((element = ply_get_next_element(iply, element))) {
            p_ply_property property = NULL;
            long nelems = 0;
            const char *element_name;
            ply_get_element_info(element, &element_name, &nelems);
            /* add this element to output file */
            if (!ply_add_element(oply, element_name, nelems)) return 0;
            /* iterate over all properties of current element */
            while ((property = ply_get_next_property(element, property))) {
                const char *property_name;
                e_ply_type type, length_type, value_type;
                ply_get_property_info(property, &property_name, &type,
                        &length_type, &value_type);
                /* setup input callback for this property */
                if (!ply_set_read_cb(iply, element_name, property_name, callback,
                        oply, 0)) return 0;
                /* add this property to output file */
                if (!ply_add_property(oply, property_name, type, length_type,
                        value_type)) return 0;
            }
        }
        return 1;
    }

    int convert(char* fileName)
    {
        const char *value;
    	char buf[32];
    	p_ply iply, oply;

    	GetTempFileName(".", "whatever", 0, buf);
    	if (MoveFileEx(fileName, buf, MOVEFILE_REPLACE_EXISTING)==0)
    	{
    		printf("rename fails with error code %d\n", GetLastError());
    		return 2;
    	}

        iply = ply_open(buf, NULL);
        if (!iply) return 1;
        if (!ply_read_header(iply)) return 1;

        oply = ply_create(fileName, PLY_LITTLE_ENDIAN, NULL);
        if (!oply) return 1;
        if (!setup_callbacks(iply, oply)) return 1;
        //pass comments and obj_infos from input to output
        value = NULL;
        while ((value = ply_get_next_comment(iply, value)))
            if (!ply_add_comment(oply, value)) return 1;
        value = NULL;
        while ((value = ply_get_next_obj_info(iply, value)))
            if (!ply_add_obj_info(oply, value)) return 1;;
        // write output header
        if (!ply_write_header(oply)) return 1;
        // read input file generating callbacks that pass data to output file /
        if (!ply_read(iply)) return 1;
        // close up, we are done
        if (!ply_close(iply)) return 1;
        if (!ply_close(oply)) return 1;

    	if (DeleteFile(buf)==0)
    	{
    		printf("delete file failed with code %d\n", GetLastError());
    		exit(22);
    	}
        return 0;
    }

    void sortTriangle(unsigned int triangle[3])
    {
    	int swap;
    	if (triangle[0]>triangle[1])
    	{
    		swap=triangle[0];
    		triangle[0]=triangle[1];
    		triangle[1]=swap;
    	}
    	if (triangle[0]>triangle[2])
    	{
    		swap=triangle[0];
    		triangle[0]=triangle[2];
    		triangle[2]=swap;
    	}
    	if (triangle[1]>triangle[2])
    	{
    		swap=triangle[1];
    		triangle[1]=triangle[2];
    		triangle[2]=swap;
    	}
    }

    float calcEdge(FVector v1, FVector v2)
    {
    	float result=0;
    	for (int i=0; i<3; i++)
    	{
    		result+=(v1[i]-v2[i])*(v1[i]-v2[i]);
    	}
    	return sqrt(result);
    }

    void PlyModel::calcLongestEdge(UIVector tri)
    {
    	float result=0, temp;
    	float edge[3];
    	int first, second;
    	for (int i=0; i<3; i++)
    	{
    		temp=calcEdge(vertex[tri[i]].coord, vertex[tri[(i+1)%3]].coord);
    		if (temp>maxEdge)
    		{
    			maxEdge=temp;
    		}
    	}
    }

    void retrieveAllFile(char* dir, vector<string>& fileNameVector)
    {
    	HANDLE handle;
    	char curFileName[256];
    	char wildFileName[256];

    	WIN32_FIND_DATA ffd;

    	sprintf(wildFileName, "%s\\*.*", dir);
    	handle=FindFirstFile(wildFileName, &ffd);

    	if (handle==INVALID_HANDLE_VALUE)
    	{
    		printf("findfirst failed of error code =%d\n", GetLastError());
    		exit(19);
    	}

    	while (FindNextFile(handle, &ffd))
    	{
    		if (strcmp(ffd.cFileName, "..")!=0)
    		{
    			sprintf(curFileName, "%s\\%s", dir, ffd.cFileName);
    			if  (GetFileAttributes(curFileName)==FILE_ATTRIBUTE_DIRECTORY)
    			{
    				retrieveAllFile(curFileName, fileNameVector);
    			}
    			else
    			{
    				fileNameVector.push_back(string(curFileName));
    			}
    		}
    	}
    	FindClose(handle);
    }

    void PlyModel::readDirectory(char* dir)
    {
    	int i, j, k;
    	bool firstTime=true;
    	unsigned long length;
    	vector<string> fileNameVector;
    	HANDLE vertexHandle, faceHandle;
    	string fileName;

    	vertexCount=faceCount=0;

    	if ((vertexHandle=CreateFile("vertex.data", GENERIC_WRITE, 0, NULL, CREATE_ALWAYS,
    		FILE_ATTRIBUTE_NORMAL, NULL))==NULL)
    	{
    		printf("create file failed\n");
    		exit(15);
    	}
    	if ((faceHandle=CreateFile("face.data", GENERIC_WRITE, 0, NULL, CREATE_ALWAYS,
    		FILE_ATTRIBUTE_NORMAL, NULL))==NULL)
    	{
    		printf("create file failed\n");
    		exit(15);
    	}

    	retrieveAllFile(dir, fileNameVector);

    	for (i=0; i<fileNameVector.size(); i++)
    	{
    		fileName=fileNameVector[i];
    		if (convertFile(fileName.c_str())!=0)
    		{
    			printf("convert file %s failed\n", fileName.c_str());
    			exit(12);
    		}

    		readFile(fileName.c_str());

    		if (firstTime)
    		{
    			//initialize the min and max coord
    			memcpy(minCoord, vertex[0].coord, sizeof(FVector));
    			memcpy(maxCoord, vertex[0].coord, sizeof(FVector));
    			firstTime=false;
    		}

    		for (j=0; j<vNumber; j++)
    		{
    			for (k=0; k<3; k++)
    			{
    				if (vertex[j].coord[k]<minCoord[k])
    				{
    					minCoord[k]=vertex[j].coord[k];
    				}
    				if (vertex[j].coord[k]>maxCoord[k])
    				{
    					maxCoord[k]=vertex[j].coord[k];
    				}
    			}
    		}

    		if (WriteFile(vertexHandle, vertex, sizeof(PlyVertex)*vNumber,
    			&length, NULL)==0)
    		{
    			printf("write vertex file failed\n");
    			exit(17);
    		}
    		if (length!=vNumber*sizeof(PlyVertex))
    		{
    			printf("write vertex file failed\n");
    			exit(17);
    		}

    		for (j=0; j<fNumber; j++)
    		{
    			sortTriangle(triangle[j]);
    			triangle[j][0]+=vertexCount;
    			triangle[j][1]+=vertexCount;
    			triangle[j][2]+=vertexCount;
    		}
    		WriteFile(faceHandle, triangle, sizeof(FVector)*fNumber, &length, NULL);
    		vertexCount+=vNumber;
    		faceCount+=fNumber;
    		delete[]vertex;
    		delete[]triangle;
    	}

    	CloseHandle(vertexHandle);
    	CloseHandle(faceHandle);
    }

    bool PlyModel::analysis(char* buf)
    {
    	char* token;
    	token=strtok(buf, " \n");
    	while (token!=NULL)
    	{
    		if (strcmp(token, "vertex")==0)
    		{
    			token=strtok(NULL, " \n");
    			sscanf(token, "%d", &vNumber);
    			break;
    		}
    		if (strcmp(token, "face")==0)
    		{
    			token=strtok(NULL, " \n");
    			sscanf(token, "%d", &fNumber);
    			break;
    		}
    		if (strcmp(token, "end_header" )==0)
    		{
    			return true;
    		}
    		token=strtok(NULL, " \n");
    	}
    	return false;
    }

    int PlyModel::convertFile(const char* fileName)
    {
    	return convert(const_cast<char*>(fileName));
    }

    int PlyModel::readFile(const char* fileName)
    {
    	FILE* file;
    	char buffer[256];
    	int i;
    	unsigned char cNumber;
    	UCVector dummy;

    	file=fopen(fileName, "r+b");

    	while (true)
    	{
    		fgets(buffer, 256, file);
    		//i=fread(buffer, 1, 256, file);

    		if (analysis(buffer))
    		{
    			break;
    		}
    	}
    	vertex=new PlyVertex[vNumber];
    	triangle=new UIVector[fNumber];
    	if (vertex==NULL||triangle==NULL)
    	{
    		printf("allocate memory fails\n");
    		exit(13);
    	}

    	for (i=0; i<vNumber; i++)
    	{
    		fread(vertex[i].coord, sizeof(float), 3, file);
    		fread(vertex[i].color, sizeof(unsigned char), 3, file);
    		fread(vertex[i].normal, sizeof(float), 3, file);
    	}

    	for (i=0; i<fNumber; i++)
    	{
    		fread(&cNumber, sizeof(unsigned char), 1, file);
    		if (cNumber!=3)
    		{
    			printf("found trouble, some polygon is not triangle\n");
    			exit(14);
    		}
    		fread(triangle[i], sizeof(int), cNumber, file);
    		calcLongestEdge(triangle[i]);
    		//discard face color
    		fread(dummy, sizeof(UCVector), 1, file);
    	}
    	fclose(file);
    	return 0;
    }

    void convertPlyFile(char* dirName)
    {
    	DWORD timeCount;
    	HANDLE handle;
    	unsigned long length;
    	ply.maxEdge=0.0;

    	timeCount=GetTickCount();
    	ply.readDirectory(dirName);

    	timeCount=GetTickCount()- timeCount;
    	printf("convert file process takes time %d seconds\n", timeCount/1000);

    	handle=CreateFile("info.data", GENERIC_READ|GENERIC_WRITE, 0, NULL,
    		CREATE_ALWAYS, FILE_ATTRIBUTE_NORMAL, NULL);
    	if (handle==INVALID_HANDLE_VALUE)
    	{
    		printf("open info file failed\n");
    		exit(89);
    	}

    	if (WriteFile(handle, &ply, sizeof(PlyModel),&length, NULL)==0)
    	{
    		printf("error of writing info and only write %d bytes %d\n", length,  GetLastError());
    		exit(88);
    	}

    	printf("vcount=%d, fcount=%d\n", ply.vertexCount, ply.faceCount);
    	for(int i=0; i<3; i++)
    	{
    		printf("max[%c]=%f, min[%c]=%f\n", 'x'+i, ply.maxCoord[i], 'x'+i, ply.minCoord[i]);
    	}
    	printf("max edge is %f\n", ply.maxEdge);
    	CloseHandle(handle);

    }

    ¡¡

    file name: octreeLoader.h

    #ifndef OCTREELOADER_H
    #define OCTREELOADER_H

    #include "octree.h"
    //#include "range.h"
    #include "frustum.h"
    #include <cstdio>
    #include <set>

    using namespace std;

    typedef set<int> OctreeNodeSet;

    struct OctreeLoader
    {
    	HANDLE vertexHandle, faceHandle, vertexMapHandle, faceMapHandle;
    	//PageManager vertexManager, faceManager;
    	int currentIndex, nextIndex;
    	OctreeNodeSet nodeSet[2];
    	void loadOctree(char* infoFileName, char* octreeFileName);
    	void initData(char* infoFileName, char* octreeFileName,
    		char* vertexFileName, char* faceFileName);
    	Octree* root;
    	InfoStruct infoStruct;
    	CFrustum frustum;
    	void calculateOctree(int nodeIndex);

    	void displayByOctreeNode(int nodeIndex);

    	void swapBuffer();

    };

    #endif

    ¡¡

    ¡¡

    file name: octreeLoader.cpp

    #include<GL/glut.h>
    #include "octree.h"
    #include "octreeLoader.h"
    #include "frustum.h"
    #include <cstdio>
    #include <vector>

    using namespace std;

    extern void getLastError();

    void OctreeLoader::loadOctree(char* infoFileName, char* octreeFileName)
    {
    	HANDLE handle;
    	unsigned long length;
    	handle=CreateFile(infoFileName, GENERIC_READ, FILE_SHARE_READ, NULL,
    		OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL);
    	if (handle==INVALID_HANDLE_VALUE)
    	{
    		printf("read octree file failed");
    		exit(6);
    	}
    	if (ReadFile(handle, &infoStruct, sizeof(InfoStruct), &length, NULL)==0)
    	{
    		printf("read info struct failed %d\n", GetLastError());
    		exit(3);
    	}
    	CloseHandle(handle);
    	root=new Octree[infoStruct.nodeCount];

    	if (root==NULL)
    	{
    		printf("allocate octree memory failed");
    		exit(5);
    	}
    	handle=CreateFile(octreeFileName, GENERIC_READ, FILE_SHARE_READ, NULL,
    		OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL);
    	if (handle==INVALID_HANDLE_VALUE)
    	{
    		printf("read octree file failed");
    		exit(6);
    	}
    	if (ReadFile(handle, root, sizeof(Octree)*infoStruct.nodeCount, &length, NULL)==0)
    	{
    		printf("read info struct failed %d\n", GetLastError());
    		exit(3);
    	}
    	CloseHandle(handle);

    	Octree::root=root;
    	Octree::nodeCount=infoStruct.nodeCount;
    	Octree::maxEdge=infoStruct.maxEdge;

    }

    void OctreeLoader::initData(char* infoFileName, char* octreeFileName,
    		char* vertexFileName, char* faceFileName)
    {
    	loadOctree(infoFileName, octreeFileName);
    	vertexHandle=CreateFile(vertexFileName, GENERIC_READ, FILE_SHARE_READ, NULL, OPEN_EXISTING,
    		FILE_ATTRIBUTE_NORMAL, NULL);
    	if (vertexHandle==INVALID_HANDLE_VALUE)
    	{
    		printf("create failed\n");
    		exit(1);
    	}
    	vertexMapHandle=CreateFileMapping(vertexHandle, NULL, PAGE_READONLY, 0, 0, NULL);
    	if (vertexMapHandle==NULL)
    	{
    		printf("error code=%d\n", GetLastError());
    		exit(2);
    	}

    	faceHandle=CreateFile(faceFileName, GENERIC_READ, FILE_SHARE_READ, NULL, OPEN_EXISTING,
    		FILE_ATTRIBUTE_NORMAL, NULL);
    	if (faceHandle==INVALID_HANDLE_VALUE)
    	{
    		printf("create failed\n");
    		exit(1);
    	}
    	faceMapHandle=CreateFileMapping(faceHandle, NULL, PAGE_READONLY, 0, 0, NULL);
    	if (faceMapHandle==NULL)
    	{
    		printf("error code=%d\n", GetLastError());
    		exit(2);
    	}

    	currentIndex=0;
    	nextIndex=1;

    }

    void OctreeLoader::calculateOctree(int nodeIndex)
    {
    	Octree* ptr=root+nodeIndex;
    	int i;
    	int result=frustum.CubeInFrustum(ptr->center[0], ptr->center[1],
    		ptr->center[2], ptr->width);

    	switch (result)
    	{
    	case 1:
    		//inside
    		nodeSet[nextIndex].insert(nodeIndex);
    		break;
    	case 0:
    		//intersect
    		if (ptr->sonCount==0)
    		{
    			nodeSet[nextIndex].insert(nodeIndex);
    		}
    		else
    		{
    			for (i=0; i<8; i++)
    			{
    				if (ptr->sons[i]!=-1)
    				{
    					calculateOctree(ptr->sons[i]);
    				}
    			}
    		}
    		break;
    	case -1:
    		//outside
    		break;
    	}
    }

    void OctreeLoader::displayByOctreeNode(int index)
    {
    	int mapSize;
    	if (root[index].sonCount>0)
    	{
    		for (int i=0; i<8; i++)
    		{
    			if (root[index].sons[i]!=-1)
    			{
    				displayByOctreeNode(root[index].sons[i]);
    			}
    		}

    	}
    	else
    	{
    		if (root[index].vAddress==NULL)
    		{
    			mapSize=root[index].vOffset+root[index].vNumber*sizeof(PlyVertex);
    			if ((root[index].vAddress=MapViewOfFile(vertexMapHandle, FILE_MAP_READ, 0,
    				root[index].vStartPageIndex*PageSize, mapSize))==0)
    			{
    				printf("map vertex failed %d\n", GetLastError());
    				exit(4);
    			}

    		}
    		/////////////////////////////

    		//////////////////////////////

    		if (root[index].fAddress==NULL&&root[index].fNumber>0)
    		{
    			mapSize=root[index].fOffset+root[index].fNumber*sizeof(UIVector);
    			if((root[index].fAddress=MapViewOfFile(faceMapHandle, FILE_MAP_READ, 0,
    				root[index].fStartPageIndex*PageSize, mapSize))==0)
    			{
    				printf("map face failed\n");
    				exit(4);
    			}

    		}
    		/////////////////////////////

    		/////////////////////////////

    		PlyVertex* vPtr;
    		UIVector* fPtr;
    		vPtr=(PlyVertex*)((char*)(root[index].vAddress) +root[index].vOffset);
    		fPtr=(UIVector*)((char*)(root[index].fAddress)+root[index].fOffset);

    		glEnableClientState (GL_VERTEX_ARRAY);
    		glEnableClientState (GL_NORMAL_ARRAY);
    		glEnableClientState (GL_COLOR_ARRAY);
    		getLastError();
    		glVertexPointer(3, GL_FLOAT, sizeof(PlyVertex), &(vPtr->coord));
    		glNormalPointer(GL_FLOAT, sizeof(PlyVertex), &(vPtr->normal));
    		glColorPointer(3, GL_UNSIGNED_BYTE, sizeof(PlyVertex), &(vPtr->color));

    		getLastError();

    		if (root[index].fNumber>0)
    		{
    			glDrawElements(GL_TRIANGLES, root[index].fNumber, GL_UNSIGNED_INT, fPtr);
    		}

    		getLastError();

    		glDisableClientState (GL_VERTEX_ARRAY);
    		glDisableClientState (GL_NORMAL_ARRAY);
    		glDisableClientState (GL_COLOR_ARRAY);
    		glFinish();

    		/////////////////////////////////////////
    		if (root[index].vAddress!=NULL)
    		{
    			if (UnmapViewOfFile(root[index].vAddress)==0)
    			{
    				printf("unmap failed\n");
    				exit(3);
    			}
    			root[index].vAddress=NULL;
    		}
    		if (root[index].fAddress!=NULL)
    		{
    			if (UnmapViewOfFile(root[index].fAddress)==0)
    			{
    				printf("unmap failed\n");
    				exit(3);
    			}
    			root[index].fAddress=NULL;
    		}
    		/////////////////////////////////////
    	}

    }

    void OctreeLoader::swapBuffer()
    {
    	int i;
    	OctreeNodeSet::const_iterator it;

    	frustum.CalculateFrustum();
    	calculateOctree(0);

    /*
    	for (it=nodeSet[currentIndex].begin(); it!=nodeSet[currentIndex].end(); ++it)
    	{
    		if (nodeSet[nextIndex].find(*it)==nodeSet[nextIndex].end())
    		{
    			//unmap
    			if (UnmapViewOfFile(root[*it].vAddress)==0)
    			{
    				printf("unmap failed\n");
    				exit(3);
    			}
    			root[*it].vAddress=NULL;
    			if (UnmapViewOfFile(root[*it].fAddress)==0)
    			{
    				printf("unmap failed\n");
    				exit(3);
    			}
    			root[*it].fAddress=NULL;
    		}
    	}
    */

    	for (it=nodeSet[nextIndex].begin(); it!=nodeSet[nextIndex].end(); it++)
    	{
    		int index=*it;
    		int mapSize;

    		displayByOctreeNode(index);

    	}

    	//int swap=nextIndex;
    	//nextIndex=currentIndex;
    	//currentIndex=swap;
    }

    ¡¡

    ¡¡

    file name: renderEngine.cpp

    /***********************************************************************************
    C:\Program Files\Microsoft Visual Studio\MyProjects\octree\Debug>octree
    total process time 3109
    vcount=216734, fcount=222728
    max[x]=35146.398438, min[x]=-30148.000000
    max[y]=85618.000000, min[y]=6990.450195
    max[z]=92353.500000, min[z]=-18330.500000
    finished vertex processing and now start face processing...
    finished face processing, now start octree calculation...
    total vertex is 216734 face is 256275
    total process time 144062
    octree info: nodeCount=948, nodeSize=220
    totalNeighbourCount=6818
    max depth of octree is 7

    ************************************************************************************/

    #include <GL/glut.h>
    #include "PlyModel.h"
    #include "octree.h"
    #include "octreeloader.h"
    #include <sys/stat.h>
    #include <cstdio>
    #include <cmath>
    #include <windows.h>

    //#define GL_ARB_shader_objects 1

    int scrW=800;
    int scrH=800;

    char* vertexFileName="vertex.data";
    char* faceFileName="face.data";
    char* infoFileName="info.data";
    char* outVertexFileName="outVertex.data";
    char* outFaceFileName="outFace.data";
    char* octreeFileName="octree.data";
    char* plyDirName="small";

    OctreeLoader octreeLoader;

    struct RenderContext
    {
    	DVector eye;
    	DVector center;
    	DVector up;
    	double	fovy;
    	double	aspect;
    	double	zNear;
    	double	zFar;
    };

    RenderContext renderContext;

    GLfloat preMatrix[16]=
    {
    	1,0,0,0,
    	0,1,0,0,
    	0,0,1,0,
    	0,0,0,1
    };

    float* LightPosition;
    float* LightColor;
    float* TeapotColor;

    const float Offset=0.05;

    float pos_x=0.0, pos_y=0.0, pos_z=0.0;
    float angle=0.0;

    float dist=10.0;
    float zoom=1.0;
    float nearPlane=dist/100.0f;
    float farPlane=dist*10.0f;

    char buffer[256];
    char* sourceCode;

    GLuint myTexture;

    typedef GLfloat Color[4];

    enum ColorName
    {
    	BLACK, RED, GREEN, BLUE, YELLOW, MAGENTA, PINK, WHITE
    };
    int length;
    const int MaxLightNumber=8;
    Color color[MaxLightNumber] =
    {
    	{ 0.0, 0.0, 0.0 , 1.0},	// BLACK
    	{ 1.0, 0.0, 0.0 , 1.0},	// RED
    	{ 0.0, 1.0, 0.0 , 1.0},	// GREEN
    	{ 0.0, 0.0, 1.0 , 1.0},	// BLUE
    	{ 1.0, 1.0, 0.0 , 1.0},	// YEWLLOW
    	{ 0.0, 1.0, 1.0 , 1.0},	// MAGENTA
    	{ 1.0, 0.0, 1.0 , 1.0},	// PINK
    	{ 1.0, 1.0, 1.0 , 1.0},	// WHITE

    };						// Definition of the COLORS for the unit cube model.

    GLfloat no_mat[] = { 0.0, 0.0, 0.0, 1.0 };
    GLfloat half_mat[] = { 1.0, 1.0, 1.0, 1.0 };
    GLfloat light_diffuse[] = {0.9, 0.9, 0.9, 1.0};	//{ 1.0, 1.0, 1.0, 1.0 };
    GLfloat light_specular[]={1.0, 1.0, 1.0, 1.0};

    GLint location;

    GLfloat light_position[MaxLightNumber][4] =
    {
    	{1.0, 3.0, -16.0, 1.0 },
    	{0.0, 15.0, 0.0, 1.0 },
    	{0.0, 16.0, 0.0, 1.0 },
    	{0.0, 21.0, 0.0, 1.0 },
    	{0.0, -10.0, 0.0, 1.0 },
    	{0.0, -14.0, 0.0, 1.0 },
    	{0.0, -16.0, 0.0, 1.0 },
    	{0.0, -11.0, 0.0, 1.0 }
    };

    GLfloat light_direction[MaxLightNumber][3] =
    {
    	{0.0, -1.0, 0.0 },
    	{0.0, -1.0, 0.0 },
    	{0.0, -1.0, 1.0 },
    	{0.0, -1.0, 0.0 },
    	{0.0, 1.0, 0.0 },
    	{0.0, 1.0, 0.0 },
    	{0.0, 1.0, 0.0 },
    	{0.0, 1.0, 0.0 }
    };

    GLfloat fire[] = { 1.0f, 0.6f, 0.2f, 0.8f };

    void getLastError()
    {
    	GLenum err;
    	char msg[256];
    	err=glGetError();
    	switch(err)
    	{
    	case GL_NO_ERROR:
    		sprintf(msg, "GL_NO_ERROR");
    		break;
    	case GL_INVALID_ENUM:
    		sprintf(msg, "GL_INVALID_ENUM");
    		break;
    	case GL_INVALID_VALUE:
    		sprintf(msg, "GL_INVALID_VALUE");
    		break;
    	case GL_INVALID_OPERATION:
    		sprintf(msg, "GL_INVALID_OPERATION");
    		break;
    	case GL_STACK_OVERFLOW:
    		sprintf(msg, "GL_STACK_OVERFLOW");
    		break;
    	case GL_STACK_UNDERFLOW:
    		sprintf(msg, "GL_STACK_UNDERFLOW");
    		break;
    	case GL_OUT_OF_MEMORY:
    		sprintf(msg, "GL_OUT_OF_MEMORY");
    		break;
    	}
    	if (err!=GL_NO_ERROR)
    	{
    		printf("error message %s\n", msg);
    		exit(err);
    	}
    }

    void updateRenderContext()
    {
    	renderContext.center[0]=pos_x;
    	renderContext.center[1]=1.8;
    	renderContext.center[2]=pos_z;
    	renderContext.eye[0]=renderContext.center[0] + 10*sin(angle/3.14);
    	renderContext.eye[1]=renderContext.center[1];
    	renderContext.eye[2]=renderContext.center[2]+10.0*cos(angle/3.14);
    	renderContext.up[0]=renderContext.center[0];
    	renderContext.up[1]=renderContext.center[1]+10;
    	renderContext.up[2]=renderContext.center[2];
    	gluLookAt(renderContext.eye[0], renderContext.eye[1], renderContext.eye[2],
    		renderContext.center[0], renderContext.center[1], renderContext.center[2],
    		renderContext.up[0], renderContext.up[1] , renderContext.up[2]);

    }

    void init()
    {
    	renderContext.fovy=1;
    	renderContext.aspect=1.0;
    	renderContext.zFar=90000;
    	renderContext.zNear=0.1;

    	glShadeModel(GL_FLAT);
    	glMatrixMode(GL_PROJECTION);
    	glLoadIdentity();
    	gluPerspective(renderContext.fovy*40.0f, (float)scrW/(float)scrH,
    		renderContext.zNear, renderContext.zFar);
    	glMatrixMode(GL_MODELVIEW);
    	glLoadIdentity();
    	updateRenderContext();
    	octreeLoader.initData(infoFileName, octreeFileName, outVertexFileName, outFaceFileName);

    }

    void initLights()
    					// Initialize te common lights related data.
    {
    	GLenum lightEnum[MaxLightNumber]=
    	{
    		GL_LIGHT0, GL_LIGHT1, GL_LIGHT2,  GL_LIGHT3, GL_LIGHT4, GL_LIGHT5,
    			GL_LIGHT6,GL_LIGHT7
    	};

    	int i;
    	//glLightfv(GL_LIGHT0, GL_AMBIENT, no_mat);
    	for (i=0; i<MaxLightNumber; i++)
    	{
    		glLightfv(lightEnum[i], GL_AMBIENT, no_mat);
    		//glLightfv(lightEnum[i], GL_DIFFUSE, light_diffuse);
    		glLightfv(lightEnum[i], GL_DIFFUSE, color[i]);
    		//glLightfv(lightEnum[i], GL_SPECULAR, light_specular);
    		glLightfv(lightEnum[i], GL_SPECULAR, color[i]);

    		glLightf(lightEnum[i], GL_CONSTANT_ATTENUATION, 1.0);
    		glLightf(lightEnum[i], GL_LINEAR_ATTENUATION, 0.0);
    		glLightf(lightEnum[i], GL_QUADRATIC_ATTENUATION, 0.0);

    		glLightfv(lightEnum[i], GL_POSITION, light_position[i]);
    		glLightfv(lightEnum[i], GL_SPOT_DIRECTION, light_direction[i]);
    		glLightf(lightEnum[i], GL_SPOT_CUTOFF, 85.0);
    	}

    	glMaterialfv(GL_FRONT, GL_AMBIENT_AND_DIFFUSE, color[WHITE]);
    	glMaterialfv(GL_FRONT, GL_SPECULAR, color[WHITE]);
    	glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, 50.0);

    	//glEnable(GL_LIGHT0);
    	glColorMaterial(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE);

    	glEnable(GL_COLOR_MATERIAL);

    	glEnable(GL_LIGHTING);
    	for (i=0; i<MaxLightNumber; i++)
    	{
    		glEnable(lightEnum[i]);
    	}

    }

    void displayCallbackProc()
    {
    	static bool firstTime=true;
    	glClearColor(1.0, 1.0, 1.0, 1.0);
    	glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);

    	octreeLoader.swapBuffer();

    	glutSwapBuffers();
    }

    void keyboardCallbackProc(unsigned char key, int x, int y)
    {
    	switch (key)
    	{
    	case 27:
    		exit(0);
    		break;
    	case 'l':
    		pos_x-=Offset;
    		break;
    	case 'r':
    		pos_x+=Offset;
    		break;
    	case 'u':
    		pos_y+=Offset;
    		break;
    	case 'd':
    		pos_y-=Offset;
    		break;
    	case '1':
    		angle+=5.0;
    		break;
    	case '2':
    		angle-=5.0;
    		break;

    	}
    	updateRenderContext();
    	glutPostRedisplay();
    }

    extern void dumpEveryThing();
    extern void convertPlyFile(char* dirName);

    void createFiles(char* dirName=plyDirName);

    void createFiles(char* dirName)
    {
    	printf("now starting merging ply files...\n");
    	convertPlyFile(dirName);

    	printf("now starting creating octree files...\n");
    	DWORD start=GetTickCount();
    	Octree::createOctreeFromFile(vertexFileName, faceFileName, infoFileName,
    		outVertexFileName, outFaceFileName, octreeFileName);

    	start=GetTickCount()-start;
    	printf("create octree process takes time %d seconds\n", start/1000);
    	printf("octree info: nodeCount=%d, nodeSize=%d\n", Octree::nodeCount, sizeof(Octree));
    	printf("totalNeighbourCount=%d\n", Octree::totalNeighbourCount);
    	printf("max depth of octree is %d and max edge is%f\n", Octree::maxDepth, Octree::maxEdge);
    }

    int main(int argc, char** argv)
    {
    	if (argc==1)
    	{
    		glutInit(&argc, argv);

    		glutInitDisplayMode(GLUT_DOUBLE | GLUT_RGBA | GLUT_DEPTH | GLUT_ALPHA);
    		glutInitWindowSize(800, 800);
    		glutInitWindowPosition(10, 10);
    		glutCreateWindow("out of core rendering");
    		glutDisplayFunc(displayCallbackProc);

    		glutKeyboardFunc(keyboardCallbackProc);
    		//glutReshapeFunc(reShapeCallbackProc);
    		//glutIdleFunc(idleCallbackProc);
    		//glutMouseFunc(mouseCallbackProc);
    		//glutMotionFunc(trackMotionCallbackProc);

    		//createFiles();
    		init();
    		glutMainLoop();
    	}
    	else
    	{
    		createFiles(argv[1]);
    	}
    	return 0;
    }

    /*********************************************************
    middle files.

    convert file process takes time 6 seconds
    vcount=319923, fcount=382062
    max[x]=27049.550781, min[x]=-26959.589844
    max[y]=88348.546875, min[y]=8529.599609
    max[z]=111630.500000, min[z]=-14434.150391
    max edge is 47640.019531
    finished vertex processing and now start face processing...
    finished face processing, now start octree calculation...
    total vertex is 319923 face is 350371
    create octree process takes time 199875
    octree info: nodeCount=1186, nodeSize=220
    totalNeighbourCount=5952
    max depth of octree is 8 and max edge is47640.019531
    Press any key to continue

    **************************************************************/

    ¡¡

    ¡¡

    ¡¡

    ///////////////////////////////

    The input files are downloaded from here:

    C:\data\running\second>octree powerplant
    now starting merging ply files...
    convert file process takes time 308 seconds
    vcount=11070509, fcount=12748510
    max[x]=406335.000000, min[x]=-205000.000000
    max[y]=249000.000000, min[y]=0.000000
    max[z]=160098.453125, min[z]=-25646.591797
    max edge is 214669.968750
    now starting creating octree files...
    finished vertex processing and now start face processing...
    finished face processing, now start octree calculation...
    total vertex is 11070509 face is 9294341
    create octree process takes time 2233 seconds
    octree info: nodeCount=44414, nodeSize=220
    totalNeighbourCount=265602
    max depth of octree is 12 and max edge is214669.968750

    C:\data\running\second>

    ¡¡