/* ********************************************************************** sphere.cpp atelier iebele abel - 2001 atelier@iebele.nl http://www.iebele.nl sphere_line_intersection function adapted from: http://astronomy.swin.edu.au/pbourke/geometry/sphereline Paul Bourke pbourke@swin.edu.au ********************************************************************** */ #include #include #include "sphere.h" // "sphere.h" // float * sphere_line_intersection ( // float x1, float y1 , float z1, // float x2, float y2 , float z2, // float x3, float y3 , float z3, float r ); float square( float f ) { return (f*f) ;}; float * sphere_line_intersection ( float x1, float y1 , float z1, float x2, float y2 , float z2, float x3, float y3 , float z3, float r ) { // x1,y1,z1 P1 coordinates (point of line) // x2,y2,z2 P2 coordinates (point of line) // x3,y3,z3, r P3 coordinates and radius (sphere) // x,y,z intersection coordinates // // This function returns a pointer array which first index indicates // the number of intersection point, followed by coordinate pairs. float x , y , z; float a, b, c, mu, i ; float * p; // array containing #points and its coordinates p = (float*) malloc(7*sizeof(float)); a = square(x2 - x1) + square(y2 - y1) + square(z2 - z1); b = 2* ( (x2 - x1)*(x1 - x3) + (y2 - y1)*(y1 - y3) + (z2 - z1)*(z1 - z3) ) ; c = square(x3) + square(y3) + square(z3) + square(x1) + square(y1) + square(z1) - 2* ( x3*x1 + y3*y1 + z3*z1 ) - square(r) ; i = b * b - 4 * a * c ; if ( i < 0.0 ) { // no intersection p[0] = 0.0; return(p); } if ( i == 0.0 ) { // one intersection p[0] = 1.0; mu = -b/(2*a) ; p[1] = x1 + mu*(x2-x1); p[2] = y1 + mu*(y2-y1); p[3] = z1 + mu*(z2-z1); return(p); } if ( i > 0.0 ) { // two intersections p[0] = 2.0; // first intersection mu = (-b + sqrt( square(b) - 4*a*c )) / (2*a); p[1] = x1 + mu*(x2-x1); p[2] = y1 + mu*(y2-y1); p[3] = z1 + mu*(z2-z1); // second intersection mu = (-b - sqrt(square(b) - 4*a*c )) / (2*a); p[4] = x1 + mu*(x2-x1); p[5] = y1 + mu*(y2-y1); p[6] = z1 + mu*(z2-z1); return(p); } } // ************************************************************************ // Example use of function sphere_line_intersection (...) // main() { int i, num ; float * p; p=sphere_line_intersection ( 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 , 3.0 ); num = (int) *p++; cout << " Number of intersection points is " << num << endl; for ( i=0 ; i < num ; i++ ) { cout << " X-coordinate of point " << i << " : " << *p++ << endl ; cout << " Y-coordinate of point " << i << " : " << *p++ << endl ; cout << " Z-coordinate of point " << i << " : " << *p++ << endl ; } free(&p); } // ********************************************************************