The following steps assume you start with two points – A & B – and are trying to determine the point C to be used to form a triangle:
a. Create a member function that will determine if a given point C is left or right of the line formed by two points A and B. Hint: to do this, take the cross product of the vector between the two points A and B and the vector between A and C. Since the cross product is proportional to the sine of the angle between two vectors, it will be a positive value for an angle between 0 and 180 degrees (i.e., if point C is to the left of the line from A to B).
b. Create a member function that will determine if a given point is inside the circle formed by three other points (in this case, the three points will be the points of a triangle and the resulting circle will be the circumcircle). Hint: a very elegant implementation of this function can be found on the Wikipedia entry for Delaunay Triangulation. If you borrow this implementation, make sure to test it thoroughly, reference it, and explain how it works in your lab report.
c. Create a member function that, given two points, finds the point for the next Delaunay triangle. This function can search the whole list (not the most efficient implementation, but sufficient for this lab) and most likely will call the functions defined in 4a and 4b above.
d. Create a recursive member function, Delaunay(Point A, Point B), that calls the function from 4c above the find the next point C. When point C is found, this function should insert a new triangle into the triangle vector and also update the Boolean variable in the point list. This function should then recursively call itself using points A and C, and again using points C and B. Make sure the function also has two base cases; one for when the point C found has already been used, in which case the triangle should be added, but the recursive call should not occur, and; another case where no point C is found (because A and B form an outside edge of the dataset)
I have completed the code but seem to keep coming up with error in the debugging process and I believe I've narrowed it down to the read function so here is my code:
PointList::PointList()
{
totPt = 0;
totTri = 0;
}
void PointList::read(const char* filename)
{
FILE* file;
file = fopen ( filename, "r" );
fscanf ( file, "%d \n", &totPt);
datapoint.resize( totPt );
for ( unsigned int i=0; i < datapoint.size(); i++ )
{
fscanf ( file, " %d %lf %lf \n", &datapoint.at(i).ptnum, &datapoint.at(i).xCoord, &datapoint.at(i).yCoord );
}
}
void PointList::TriPrint()
{
cout << "# of triangles created: " << triPoint.size() << endl;
cout << "Triangle # : Vertice1 Vertice2 Vertice3" << endl;
for ( unsigned int i = 0; i < triPoint.size() ; i++)
{
cout << triPoint.at(i).triNum << " : " << triPoint.at(i).vert1.ptnum << " " << triPoint.at(i).vert2.ptnum << " " << triPoint.at(i).vert3.ptnum << endl;
}
return;
}
//TASK 4a: will determine if a given point C is left or right of the line forme by the two points a and b
bool PointList::isRight( double a, double b, double c )
{
Point A = datapoint.at(a-1);
Point B = datapoint.at(b-1);
Point C = datapoint.at(c-1);
//Cross product of AB and AC
Point AB;
Point AC;
AB.xCoord = B.xCoord - A.xCoord;
AB.yCoord = B.yCoord - A.yCoord;
AC.xCoord = C.xCoord - A.xCoord;
AC.yCoord = C.yCoord - A.yCoord;
double det = (AB.xCoord*AC.yCoord) - (AC.xCoord*AB.yCoord);
if ( det >= 0 )
{
return true;
}
else
{
return false;
}
}
//Task 4b: determine if a given point is inside circle formed by three other points
bool PointList::inCircle ( double a, double b, double c, double d )
{
bool inCirc = false;
Point A = datapoint.at(a-1);
Point B = datapoint.at(b-1);
Point C = datapoint.at(c-1);
Point D = datapoint.at(d-1);
vector <vector<double>> Matrix;
Matrix.resize(3);
for ( int i = 0; i < 3 ; i++)
{
Matrix.at(i).resize(3);
}
Matrix.at(0).at(0) = A.xCoord - D.xCoord;
Matrix.at(1).at(0) = B.xCoord - D.xCoord;
Matrix.at(2).at(0) = C.xCoord - D.xCoord;
Matrix.at(0).at(1) = A.yCoord - D.yCoord;
Matrix.at(1).at(1) = B.yCoord - D.yCoord;
Matrix.at(2).at(1) = C.yCoord - D.yCoord;
Matrix.at(0).at(2) = ((A.xCoord*A.xCoord) - (D.xCoord*D.xCoord)) + ((A.yCoord*A.yCoord) - (D.yCoord*D.yCoord));
Matrix.at(1).at(2) = ((B.xCoord*B.xCoord) - (D.xCoord*D.xCoord)) + ((B.yCoord*B.yCoord) - (D.yCoord*D.yCoord));
Matrix.at(2).at(2) = ((C.xCoord*C.xCoord) - (D.xCoord*D.xCoord)) + ((C.yCoord*C.yCoord) - (D.yCoord*D.yCoord));
double det = 0;
det = (Matrix.at(0).at(0) * Matrix.at(1).at(1) * Matrix.at(2).at(2)) + (Matrix.at(0).at(1) * Matrix.at(1).at(2) * Matrix.at(2).at(0)) + (Matrix.at(0).at(2) * Matrix.at(1).at(0) * Matrix.at(2).at(1));
det = det - (Matrix.at(2).at(0) * Matrix.at(1).at(1) * Matrix.at(0).at(2)) - (Matrix.at(2).at(1) * Matrix.at(1).at(2) * Matrix.at(0).at(0)) - (Matrix.at(2).at(2) * Matrix.at(1).at(0) * Matrix.at(0).at(1));
if ( det >= 0 ) // determinant is positive if and only if D lies inside the circumcircle
{
inCirc = false;
}
else
{
inCirc = true;
}
return inCirc;
}
//Task 4c: Given two points, finds the point for the next Delaunay triangle
double PointList::nextDelaunay ( double a, double b )
{
bool incircle = false;
bool isright = false;
bool noRight = false;
int f = -1;
//check all points in file
for ( int i = 0; i < datapoint.size() ; i++)
{
if ( i == a || i == b)
{
continue;
}
else
{
isright = isRight( a, b, i); //checks to see if its to the right
if(isright)
{
noRight = false;
for ( int j = 1; j < datapoint.size(); j++ )
{
if ( j == a || j == b || j == i )
{
continue;
}
else
{
incircle = inCircle ( a, b, i, j );
//checks to see if point in circle
if(incircle)
{
break;
}
}
}
if ( !incircle)
{
return i;
}
}
else
{
continue;
}
}
}
if (noRight)
{
return f;
}
}
//Task 4d: a recursive member function, Delaunay, the calls from 4c above to find the next point C
void PointList::Delaunay( int a, int b )
{
Triangle x;
Point A = datapoint.at(a - 1);
Point B = datapoint.at(b - 1);
int c = nextDelaunay(a,b);
if ( c == -1)
{
return;
}
else
{
Point C = datapoint.at(c-1);
if ( C.usedPt == false)
{
x.vert1.ptnum = a;
x.vert1.xCoord = A.xCoord;
x.vert1.yCoord = A.yCoord;
x.vert2.ptnum = b;
x.vert2.xCoord = B.xCoord;
x.vert2.yCoord = B.yCoord;
x.vert3.ptnum = c;
x.vert3.xCoord = C.xCoord;
x.vert3.yCoord = C.yCoord;
x.triNum = triPoint.size()+1;
triPoint.push_back(x);
datapoint.at(c-1).usedPt = true;
Delaunay( a, c );
Delaunay( c, b );
return;
}
if ( C.usedPt == true )
{
x.vert1.ptnum = a;
x.vert1.xCoord = A.xCoord;
x.vert1.yCoord = A.yCoord;
x.vert2.ptnum = b;
x.vert2.xCoord = B.xCoord;
x.vert2.yCoord = B.yCoord;
x.vert3.ptnum = c;
x.vert3.xCoord = C.xCoord;
x.vert3.yCoord = C.yCoord;
x.triNum = triPoint.size()+1;
triPoint.push_back(x);
return;
}
}
}
read
, it's not immediately clear what you think the problem might be with that function. Have you checked the return value fromfscanf
to ensure it's doing as many conversions as you expect? – Joe Z