BoundingSphere ,RayIntersection implementation ...
RuiTheAlmeida wrote at 2012-11-29 11:48:
Hi Objo and everybody,
First of all I want to congratulate you for this great Helix 3D tookit
I had to use one of the methods of class BoundingSphere and as it was not implemented so I did it. Feel free to test and include the code I'm sending here.
Best Regards,
Rui
public bool RayIntersection(Ray3D ray, out Point3D[] result) { result = null; bool inter = false; double h = this.center.X ; double j = this.center.Y; double k = this.center.Z; double ra = this.radius; double x1 = ray.Origin.X ; double y1 = ray.Origin.Y; double z1 = ray.Origin.Z; double a = x1 - ray.Origin.X * ray.Direction.X; double b = y1 - ray.Origin.Y * ray.Direction.Y; double c = z1 - ray.Origin.Z * ray.Direction.Z; //Quadratic solving double aQ = a * a + b * b + c * c; double bQ = 2 * a * (x1 - h) + 2 * b * (y1 - j) + 2 * c * (z1 - k); double cQ = x1 * x1 + y1 * y1 + z1 * z1 + h * h + k * k + j * j - 2 * (j * y1 + k * z1 + h * x1) - ra * ra; //Dscriminant double Q = bQ * bQ - 4 * aQ * cQ; if (Q >= 0) { double rx1 = (-bQ + Math.Sqrt(bQ * bQ - 4 * aQ * cQ)) / (2 * aQ); //First Root double rx2 = (-bQ - Math.Sqrt(bQ * bQ - 4 * aQ * cQ)) / (2 * aQ); //Second Root //first Intersection Point3D I1 = new Point3D(x1 + a * rx1,y1 + b * rx1 , z1 + c * rx1); //Second Intersection Point3D I2 = new Point3D(x1 + a * rx2, y1 + b * rx2, z1 + c * rx2); //ray inside sphere if ((ray.Origin - this.center).Length <= this.radius) { result = new Point3D[1]; //Check if the vectors have the same direction if (Vector3D.DotProduct(ray.Direction, I1 - ray.Origin) >= 0) { result[0] = I1; } else { result[0] = I2; } } else { result = new Point3D[2]; result[0] = I1; result[1] = I2; } inter = true; } else { inter = false; } return inter; }
objo wrote at 2012-12-04 19:45:
thanks for the code! Did you make this yourself, or do you have a reference?
I created some unit tests for this code, but could not get my last case (no intersection) to work:
[Test, Ignore] public void RayInterSection_NoIntersections() { var ray = new Ray3D(new Point3D(0.2, 0.3, 1), new Vector3D(1, 0.1, 0.1)); var sphere = new BoundingSphere(new Point3D(0.2, 0.3, 0), 1.0); Point3D[] result; Assert.IsFalse(sphere.RayIntersection(ray, out result)); Assert.AreEqual(sphere.Radius, sphere.Center.DistanceTo(result[0]), 1e-6, "Point is not on sphere."); Assert.AreEqual(0, ray.GetNearest(result[0]).DistanceTo(result[0]), 1e-6, "Point is not on ray."); }
RuiTheAlmeida wrote at 2012-12-07 05:12:
Hi objo,
I did the all code myself, and yes I have some references from Wikipedia and other common links that I will give.
I'm also finishing a detailed presentation with the all theoretical stuff, and I will post it all here.
Regards,
Rui
objo wrote at 2012-12-07 08:01:
Good! Can you also look at the "no intersections" case? Is there an error in my test?
I found some efficient implementations at
- http://wiki.cgsociety.org/index.php/Ray_Sphere_Intersection
- http://www.siggraph.org/education/materials/HyperGraph/raytrace/rtinter1.htm
- http://paulbourke.net/geometry/circlesphere/
- http://www.lighthouse3d.com/tutorials/maths/ray-sphere-intersection/
- http://www.scratchapixel.com/lessons/3d-basic-lessons/lesson-7-intersecting-simple-shapes/ray-sphere-intersection/
- http://www.ccs.neu.edu/home/fell/CSU540/programs/RayTracingFormulas.htm
RuiTheAlmeida wrote at 2012-12-09 19:11:
Hi Objo,
I've made some improvements, and you can find the new code below. I'm also posting some test cases. (I've changed your "no intersections" case )
You can also find my detailed explanation and references here
thanks,
Rui
public bool RayIntersection(Ray3D ray, out Point3D[] result) { result = new Point3D[2]; result[0] = new Point3D(double.PositiveInfinity, double.PositiveInfinity, double.PositiveInfinity); result[1] = new Point3D(double.PositiveInfinity, double.PositiveInfinity, double.PositiveInfinity); bool inter = false; double h = this.center.X; double j = this.center.Y; double k = this.center.Z; double ra = this.radius; double x1 = ray.Origin.X; double y1 = ray.Origin.Y; double z1 = ray.Origin.Z; //[a,b,c] = Ray.direction double a = ray.Direction.X; double b = ray.Direction.Y; double c = ray.Direction.Z; //Quadratic solving double aQ = a * a + b * b + c * c; double bQ = 2 * a * (x1 - h) + 2 * b * (y1 - j) + 2 * c * (z1 - k); double cQ = x1 * x1 + y1 * y1 + z1 * z1 + h * h + k * k + j * j - 2 * (j * y1 + k * z1 + h * x1) - ra * ra; //Dscriminant double Q = bQ * bQ - 4 * aQ * cQ; if (Q >= 0) //We have at least one possible Intersection { double rx1 = (-bQ + Math.Sqrt(bQ * bQ - 4 * aQ * cQ)) / (2 * aQ); //First Root double rx2 = (-bQ - Math.Sqrt(bQ * bQ - 4 * aQ * cQ)) / (2 * aQ); //Second Root //first Possible Intersection Point3D I1 = new Point3D(x1 + a * rx1, y1 + b * rx1, z1 + c * rx1); //Second Possible Intersection Point3D I2 = new Point3D(x1 + a * rx2, y1 + b * rx2, z1 + c * rx2); //Check if I1 is intersection by laying in the ray if (Vector3D.DotProduct(ray.Direction, I1 - ray.Origin) >= 0) { result[0] = I1; inter = true; } //Check if I2 is intersection by laying in the ray if (Vector3D.DotProduct(ray.Direction, I2 - ray.Origin) >= 0) { result[1] = I2; inter = true; } } return inter; }
TESTS:
[Test] public void RayInterSection_HasIntersections() { var ray = new Ray3D(new Point3D(0.2, 0.3, 0), new Vector3D(1, 0, 0)); var sphere = new BoundingSphere(new Point3D(0.2, 0.3, 0), 1.0); Point3D[] result; Assert.IsTrue(sphere.RayIntersection(ray, out result)); Assert.AreEqual(sphere.Radius, sphere.Center.DistanceTo(result[0]), 1e-6, "Point is not on sphere."); Assert.AreEqual(0, ray.GetNearest(result[0]).DistanceTo(result[0]), 1e-6, "Point "+ result[0] +" is not on ray."); } [Test] public void RayInterSection_NoIntersections() { var ray = new Ray3D(new Point3D(0.2, 0.3, 1), new Vector3D(1, 0.1, 0.1)); var sphere = new BoundingSphere(new Point3D(0.2, 0.3, 0), 1.0); Point3D[] result; Assert.IsFalse(sphere.RayIntersection(ray, out result)); Assert.IsTrue( double.IsInfinity(sphere.Center.DistanceTo(result[0])), "There is a finite solution 1."); Assert.AreNotEqual(0, ray.GetNearest(result[0]).DistanceTo(result[0]), "Point 1 is on ray."); Assert.IsTrue(double.IsInfinity(sphere.Center.DistanceTo(result[1])), "There is a finite solution 2."); Assert.AreNotEqual(0, ray.GetNearest(result[0]).DistanceTo(result[1]), "Point 2 is on ray."); } [Test] public void RayInterSection_TwoIntersections() { var sphere = new BoundingSphere(new Point3D(0.2, 0.3, 0), 1.0); var ray = new Ray3D(new Point3D(12, 23, 32), sphere.Center - new Point3D(12, 23, 32)); // vector away and pointing to the center of sphere Point3D[] result; Assert.IsTrue(sphere.RayIntersection(ray, out result)); Assert.AreEqual(sphere.Radius, sphere.Center.DistanceTo(result[0]), 1e-6, "Point 1 is not on sphere."); Assert.AreEqual(0, ray.GetNearest(result[0]).DistanceTo(result[0]), 1e-6, "Point 1" + result[0] + " is not on ray."); Assert.AreEqual(sphere.Radius, sphere.Center.DistanceTo(result[1]), 1e-6, "Point 2 is not on sphere."); Assert.AreEqual(0, ray.GetNearest(result[1]).DistanceTo(result[1]), 1e-6, "Point 2 " + result[0] + " is not on ray."); } [Test] public void RayInterSection_OnlyOneIntersection() { var sphere = new BoundingSphere(new Point3D(0.2, 0.3, 0), 1.0); var ray = new Ray3D(sphere.Center, new Vector3D(1, 0.1, 0.1)); // vector origin in the center of sphere pointing away Point3D[] result; Assert.IsTrue(sphere.RayIntersection(ray, out result), "No intersection"); Assert.IsTrue(sphere.Center.DistanceTo(ray.Origin ) <= sphere.Radius,"Ray origin outside sphere"); Assert.IsFalse( double.IsInfinity(sphere.Center.DistanceTo(result[0])) && double.IsInfinity(sphere.Center.DistanceTo(result[1])), "No Intersections"); Assert.IsFalse(!double.IsInfinity(sphere.Center.DistanceTo(result[0])) && !double.IsInfinity(sphere.Center.DistanceTo(result[1])), "Two Intersections"); } [Test] public void RayInterSection_TwoTangentIntersections() { var sphere = new BoundingSphere(new Point3D(0.2, 0.3, 0), 1.0); var ray = new Ray3D(new Point3D(sphere.Center.X + sphere.Radius, sphere.Center.Y + sphere.Radius, sphere.Center.Z + 50), new Point3D(sphere.Center.X + sphere.Radius, sphere.Center.Y, sphere.Center.Z) - new Point3D(sphere.Center.X + sphere.Radius, sphere.Center.Y + sphere.Radius, sphere.Center.Z + 50)); // ray tangent to sphere Point3D[] result; Assert.IsTrue(sphere.RayIntersection(ray, out result), "No intersection"); Assert.IsTrue(!double.IsInfinity(sphere.Center.DistanceTo(result[0])) && !double.IsInfinity(sphere.Center.DistanceTo(result[1])), "Less then two Intersections"); Assert.IsTrue(result[0]==result[1], "Solutions are not equal"); }
objo wrote at 2012-12-13 23:08:
thanks for the unit tests. I have done some modifications and submitted the changes!
Customer support service by UserEcho