From 39e77678eedd335aacef173f52aebf522a9ea7af Mon Sep 17 00:00:00 2001
From: Josh Holtrop <jholtrop+git@gmail.com>
Date: Thu, 19 Feb 2009 18:03:23 +0000
Subject: [PATCH] updated Cyl intersection routine

git-svn-id: svn://anubis/fart/trunk@132 7f9b0f55-74a9-4bce-be96-3c2cd072584d
---
 shapes/Cyl.cc | 57 ++++++++++++++++++++++++++++++++++++++++++++++-----
 1 file changed, 52 insertions(+), 5 deletions(-)

diff --git a/shapes/Cyl.cc b/shapes/Cyl.cc
index 66f23e5..cca5ca3 100644
--- a/shapes/Cyl.cc
+++ b/shapes/Cyl.cc
@@ -22,9 +22,44 @@ Cyl::Cyl(double bottom_radius, double top_radius, double height)
 Shape::IntersectList Cyl::intersect(const Ray & ray)
 {
     Ray ray_inv = m_inverse.transform_ray(ray);
-
     IntersectList res;
+
+    /* First intersect with the bottom plane, if it has positive area */
+    if (m_bottom_radius > 0.0)
+    {
+        LinearSolver solver(-ray_inv.getDirection()[2],
+                            -ray_inv.getOrigin()[2]);
+        Solver::Result solutions = solver.solve();
+        if (solutions.numResults > 0)
+        {
+            Vector isect_point = ray_inv[solutions.results[0]];
+            if (isect_point[0]*isect_point[0] + isect_point[1]*isect_point[1]
+                    < m_bottom_radius_2)
+            {
+                res.push_back(m_transform.transform_point(isect_point));
+            }
+        }
+    }
+
+    /* Same for the top plane */
+    if (m_top_radius > 0.0)
+    {
+        LinearSolver solver(ray_inv.getDirection()[2],
+                            ray_inv.getOrigin()[2] + m_height);
+        Solver::Result solutions = solver.solve();
+        if (solutions.numResults > 0)
+        {
+            Vector isect_point = ray_inv[solutions.results[0]];
+            if (isect_point[0]*isect_point[0] + isect_point[1]*isect_point[1]
+                    < m_top_radius_2)
+            {
+                res.push_back(m_transform.transform_point(isect_point));
+            }
+        }
+    }
+
     /*
+     * Now see if the ray hit the side of the cylinder/cone thingy
      * Ray equation:   R = R0 + tRd
      *                 x = R0x + tRdx
      *                 y = R0y + tRdy
@@ -47,15 +82,27 @@ Shape::IntersectList Cyl::intersect(const Ray & ray)
     double R0x = ray_inv.getOrigin()[0];
     double R0y = ray_inv.getOrigin()[1];
     double R0z = ray_inv.getOrigin()[2];
-    double m = m_m;
+    double m = m_slope;
     double m2 = m*m;
+    double b = m_bottom_radius;
     QuadraticSolver solver(Rdx * Rdx + Rdy * Rdy - m2 * Rdz * Rdz,
-                           2.0 * (R0x*Rdx + R0y*Rdy - m2*R0z*Rdz - m*Rdz*m_b),
+                           2.0 * (R0x*Rdx + R0y*Rdy - m2*R0z*Rdz - m*Rdz*b),
                            R0x*R0x + R0y*R0y
-                             - m2*R0z*R0z - m_b*m_b - 2*m*R0z*m_b
+                             - m2*R0z*R0z - b*b - 2*m*R0z*b
                           );
 
-    Solver::Results results = solver.solve();
+    Solver::Result solutions = solver.solve();
+    for (int i = 0; i < solutions.numResults; i++)
+    {
+        if (solutions.results[i] >= 0.0)
+        {
+            Vector isect_point = ray_inv[solutions.results[i]];
+            if (isect_point[2] >= 0.0 && isect_point[2] <= m_height)
+            {
+                res.push_back(m_transform.transform_point(isect_point));
+            }
+        }
+    }
 
     return res;
 }