Annotation of researchv10dc/cmd/sml/doc/examples/ray.sml, revision 1.1.1.1

1.1       root        1: To: research!dbm
                      2: Status: R
                      3: 
                      4: (* Find the intersection of a ray with a cylinder in three-space.  If
                      5:    there is more than one intersection, return the parameter of
                      6:    intersection closest to the start of the ray.  If there is no
                      7:    intersection, return NONE.
                      8: 
                      9:    The ray is specified by its starting point (rayStart) and direction
                     10:    vector (rayVec).  The cylinder is specified by a starting point (cStart)
                     11:    a length vector (cVec), and a radius (cRad).  The cylinder is finite;
                     12:    it does not extend beyond the direction vector at either end.
                     13: 
                     14: *)
                     15: 
                     16: local open Real
                     17:       val EP = 0.0000001
                     18: 
                     19:   fun newton f x =
                     20:     let val y = f x 
                     21:      in if abs(y)>EP then newton f (x+EP*y/(y-f(x+EP))) else x 
                     22:     end
                     23: 
                     24:   infix 7 @   fun (a,b,c) @ (d,e,f) = a*d+b*e+c*f     (* dot product *)
                     25:   infix 7 *!  fun (a,b,c) *! t = (a*t,b*t,c*t)     (* vector * scalar *)
                     26:   infix 6 +!  fun (a,b,c) +! (d,e,f) = (a+d,b+e,c+f)  (* vector addition *)
                     27:   infix 6 -!  fun (a,b,c) -! (d,e,f) = (a-d,b-e,c-f)  (* vector subtraction *)
                     28: 
                     29:   fun along (C,D) t = C +! (D *! t)
                     30: 
                     31: in 
                     32: 
                     33: fun cylhit (ray as (rayStart, rayVec), cyl as (cStart, cVec), cRad) =
                     34: let val rayStart = rayStart -! cStart
                     35:     fun cylpos p = (cVec@p) / (cVec@cVec)
                     36:         (* parameter of closest point on cVec to the point p *)
                     37: 
                     38:     fun cyldist p = let val v = p -! (along cyl (cylpos p)) in v@v end
                     39:         (* distance of p from closest point on cVec *)
                     40: 
                     41:     fun cyldist' p = p*!(cVec@cVec) -! cVec*!(p@cVec)
                     42:        (* gradient of the cyldist function *)
                     43: 
                     44:     val tf = newton(fn t => cyldist' (along ray t) @ rayVec) 0.0
                     45:        (* parameter along the ray of the closest point to cVec *)
                     46: 
                     47:  in if cyldist (along ray tf) > cRad*cRad then NONE else
                     48:    let val t0 = newton(fn t => cyldist(along ray t)-(cRad*cRad)) 0.0
                     49:        (* one of the two intersection points *)
                     50:        val t1 = newton(fn t => cyldist(along ray t)-(cRad*cRad)) (tf+tf-t0)
                     51:        (* the other intersection point *)
                     52:        fun inRange v = v >= 0.0 andalso v <= 1.0
                     53:     in case (t0 > 0.0 andalso inRange (cylpos (along ray t0)),
                     54:             t1 > 0.0 andalso inRange (cylpos (along ray t1)))
                     55:        of (false,false) => NONE
                     56:         | (true,false) => SOME t0
                     57:         | (false,true) => SOME t1
                     58:         | (true,true) => SOME (if (t0<t1) then t0 else t1)
                     59:    end
                     60: end
                     61: 
                     62: end

unix.superglobalmegacorp.com

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.