[Math] Ray-Sphere
Author
I wish to determine the parameters tin and tout where a ray enters and leaves a sphere. Is there a formula to obtain these values without computing a square root?
no. If there would be, then we would be able to use it to compute square root without computing square root. [lol]
btw if you need code,
(and of course you can compute that not using sqrt and using something else instead, like det=exp(log(det)/2.0)) but it's square root anyway)
As you can see from formula , the "det" is kind of squared distance from middle between intersection points, so maybe you can still avoid square root somehow.
edit: By the way, it doesn't compute square root for rays that do not intersect sphere, so if most rays do not intersect your spheres it is almost as good as if it wouldn't compute square root at all.
btw if you need code,
//// Line-sphere intersection. // p=(line origin position - sphere position),// d=line direction,// r=sphere radius,// i1=first intersection distance,// i2=second intersection distance// i1<=i2//// returns true if intersection found,false otherwise.// bool LineSphereIntersect(const Vec3d &p, const Vec3d &d,real64 r, real64 &i1, real64 &i2){ real64 det,b; b = -DotProd(p,d); det=(b*b) - (p.x*p.x + p.y*p.y + p.z*p.z) + r*r; if (det<0){ return false; }//else det= sqrt(det); i1= b - det; i2= b + det; return true;};(and of course you can compute that not using sqrt and using something else instead, like det=exp(log(det)/2.0)) but it's square root anyway)
As you can see from formula , the "det" is kind of squared distance from middle between intersection points, so maybe you can still avoid square root somehow.
edit: By the way, it doesn't compute square root for rays that do not intersect sphere, so if most rays do not intersect your spheres it is almost as good as if it wouldn't compute square root at all.
Author
Quote:Original post by Dmytry
no. If there would be, then we would be able to use it to compute square root without computing square root.[lol]
Excellent answer [lol]
Quote:btw if you need code
*** Source Snippet Removed ***
Thank you, but your code does not apply for three reasons: I'm not using a C family language, my ray direction vector is not normalized, and I'm actually solving an equivalent problem (moving sphere vs moving sphere collision).
Snippet of my (untested) code:
(** Compute the two tangential contact times for a two-sphere system. @param pA The initial (x,y) position of the first sphere. @param vA The velocity of the first sphere. @param rA The radius of the first sphere. @param pB The position of the second sphere. @param vB The velocity of the second sphere. @param rB The radius of the second sphere. @return The times (t0,t1) of tangential contact. @raise NoCollision If there is no contact. *)let sphereContact (pA) (vA) (rA) (pB) (vB) (rB) = let px = (fst pA) -. (fst pB) and py = (snd pA) -. (snd pB) in let vx = (fst vA) -. (fst vB) and vy = (snd vA) -. (snd vB) in let r = rA +. rB in let sq x = x *. x in let delta = let sub = px *. vy -. py *. vx in (sq r) *. ((sq px) +. (sq py)) -. (sq sub) in if delta < 0. then raise NoCollision; let d = ((sq vx) +. (sq vy)) in let delta = (sqrt delta) /. d in let x = -. (px *. vx +. py *. vy) /. d in (x -. delta),(x +. delta);;
Ray-sphere collision detection is pretty easy and can be found here. Moving sphere-sphere collision isn't that bad either and can be found here. The latter link only checks to see if two spheres will intersect, but it provides you with the second-order equation you would solve to get intersection times so you can modify the algorithm.
Quote:Original post by ToohrVykQuote:Original post by Dmytry
no. If there would be, then we would be able to use it to compute square root without computing square root.[lol]
Excellent answer [lol]Quote:btw if you need code
*** Source Snippet Removed ***
Thank you, but your code does not apply for three reasons: I'm not using a C family language, my ray direction vector is not normalized,
modded for non-normalized direction vector:
//// Line-sphere intersection. // p=(line origin position - sphere position),// d=line direction,// (line equation is p+d*t)// r=sphere radius,// i1=first intersection,// i2=second intersection// i1<=i2//// returns true if intersection found,false otherwise.// bool LineSphereIntersect(const Vec3d &p, const Vec3d &d, real64 r, real64 &i1, real64 &i2){ real64 det,b; b = -DotProd(p,d); real64 scale = 1.0/DotProd(d,d); b*=scale; det=(b*b) - ((p.x*p.x + p.y*p.y + p.z*p.z) + r*r)*scale; if (det<0){ return false; }//else det= sqrt(det); i1= b - det; i2= b + det; return true;};Not tested, but I think should work.
Quote:
and I'm actually solving an equivalent problem (moving sphere vs moving sphere collision).
Yes, that does reduce to ray-sphere intersection for sphere with radius equal to sum of radiuses.
Quote:
Snippet of my (untested) code:
*** Source Snippet Removed ***
OCaml? something else ml-based?
no offense, but it's approximately as readable as assembler (of some procesor unknown to me). Well, i can kind of figure what is going on, but, why "+." ? Is it to show that it's floating point addition?
edit: fixed typo in code
[Edited by - Dmytry on November 27, 2005 3:48:29 PM]
Quote:Original post by Dmytry
OCaml? something else ml-based?
no offense, but it's approximately as readable as assembler (of some procesor unknown to me). Well, i can kind of figure what is going on, but, why "+." ? Is it to show that it's floating point addition?
edit: fixed typo in code
You're right on both accounts. It's OCaml, and +. is the floating-point variant for +
[sorry for offtopic] Actually I just seen some sample code with Ocaml a while ago that had these strange "let" "in" "." things, but had to guess what dots everywhere is for. I have few questions: what happens if . is missed? "omg_wtf_cast" to integer? (i suppose idea is that + can't or should not be overloaded to work with integers, reals, vectors, etc.) I googled for Ocamn now, must admit some language features does look nice, but syntax looks bit like "parallel evolution" to c++ , in sense of having strange strings of symbols that is hard to decipher. (and i don't like c/c++ sintax either)
[Edited by - Dmytry on November 28, 2005 2:40:05 PM]
[Edited by - Dmytry on November 28, 2005 2:40:05 PM]
Author
OCaml is a lot more readable with decent syntax coloring and enough coding experience. And my code isn't exactly adhering to common beauty guidelines.
The difference between (+) and (+.) is their type:
val (+) : int -> int -> int
val (+.) : float -> float -> float
And since all casts in O'Caml must be explicit, the compiler will throw a fit if you use one instead of the other. There is no operator overloading in O'Caml, because the type of the operator is used to infer the type of the objects it is applied to (and not the reverse).
Here, perhaps a more readable version:
The difference between (+) and (+.) is their type:
val (+) : int -> int -> int
val (+.) : float -> float -> float
And since all casts in O'Caml must be explicit, the compiler will throw a fit if you use one instead of the other. There is no operator overloading in O'Caml, because the type of the operator is used to infer the type of the objects it is applied to (and not the reverse).
Here, perhaps a more readable version:
let vdot (ax,ay) (bx,by) = (ax*.bx+.ay*.by);;let vmin (ax,ay) (bx,by) = (ax-.bx),(ay-.by);;(** Compute the parameter of the collision between a circle and a ray. @param p Circle center. @param r Circle radius. @param a Ray origin. @param v Ray direction (not normalized) @return The time of first contact. @raise NoCollision if there is no contact in the future. *)let rayCircleContact (p) (r) (a) (v) = let p = vmin a p in let sq x = x *. x in (* a = vx² + vy² *) (* b = px.vx + py.vy --- b_h = b/2 *) (* c = px² + py² - r² *) let a = vdot v v in let b_h = vdot p v in let c = vdot p p -. (sq r) in (* det = b² - 4ac --- det_f = det/4 *) let det_f = (sq b_h) -. a *. c in if det_f < 0. then raise NoCollision; (* d = sqrt( det ) --- d_ha = d/2a *) (* m = - b/2a *) let d_ha = (sqrt det_f) /. a in let m = -. b_h /. a in (* t = - (b + d) / 2a *) let t = (m -. d_ha) in if t > 0. then t else raise NoCollision;; (** Compute the first tangential contact time for a two-circle system. @param pA The initial (x,y) position of the first circle. @param vA The velocity of the first circle. @param rA The radius of the first circle. @param pB The position of the second circle. @param vB The velocity of the second circle. @param rB The radius of the second circle. @return The first time in the future when collision occurs. @raise NoCollision If there is no contact in the future. *)let circleCircleContact (pA) (vA) (rA) (pB) (vB) (rB) = let p = vmin pA pB in let v = vmin vA vB in let r = rA +. rB in rayCircleContact (0.0,0.0) (r) (p) (v);;This topic is closed to new replies.
Advertisement
Popular Topics
Advertisement
