# Reflect a ray off a circle so it hits another point

my problem is the following:

I have two points ($e$ and $p$) in a 2D space and I am trying to figure out where on the circle is the reflection of $p$ as seen from $e$.  So the way I approached this is by looking for the vector from the center of the circle to the point, where the ray is reflected.

[For simplicity the circle's center is located at $\color{darkturquoise}{(0,0)}$ and it has a radius with unit length. The two points are represented by the vectors $\color{deeppink}{e}$ and $\color{gold}{p}$] In the picture were looking for the vector $\color{red}{r}$, which coincides with the normal vector of the surface of the circle. So this vector has to fulfill 2 properties:

1. $||\color{red}{r}|| = 1$
2. The angle between $\color{red}{r}$ and $\color{green}{p-r}$ should be equal to the angle between $\color{blue}{e-r}$ and $\color{red}{r}$.

The equation for the fist property is pretty straight forward. For the second property I came up with the following equation: $$\Bigg|\Bigg|\frac{\color{blue}{e-r}}{||\color{blue}{e-r}||}\times\color{red}{r}\Bigg|\Bigg| = \Bigg|\Bigg|\color{red}{r}\times\frac{\color{green}{p-r}}{||\color{green}{p-r}||}\Bigg|\Bigg|\\ \frac{||(\color{blue}{e-r})\times\color{red}{r}||}{||\color{blue}{e-r}||} =\frac{||\color{red}{r}\times(\color{green}{p-r})||}{||\color{green}{p-r}||}\\ \frac{1}{||\color{blue}{e-r}||}\Big|\Big|(\color{deeppink}{e}\times\color{red}{r})-\big((-\color{red}{r})\times\color{red}{r}\big)\Big|\Big|= \frac{1}{||\color{green}{p-r}||}\Big|\Big|(\color{red}{r}\times\color{gold}{p})-\big(\color{red}{r}\times(-\color{red}{r})\big)\Big|\Big|\\ \frac{||\color{deeppink}{e}\times\color{red}{r}||}{||\color{blue}{e-r}||}= \frac{||\color{red}{r}\times\color{gold}{p}||}{||\color{green}{p-r}||}\\ \frac{||\color{green}{p-r}||}{||\color{blue}{e-r}||}= \frac{||\color{red}{r}\times\color{gold}{p}||}{||\color{deeppink}{e}\times\color{red}{r}||}\\ \frac{\sqrt{1 + ||\color{deeppink}{e}||^{2} - 2(\color{red}{r}\cdot\color{pink}{e})}}{\sqrt{1 + ||\color{gold}{p}||^{2} - 2(\color{red}{r}\cdot\color{gold}{p})}}= \frac{||\color{red}{r}\times\color{gold}{p}||}{||\color{deeppink}{e}\times\color{red}{r}||}$$ After that I use the fact that the cross product of 2 2D-Vectors has a single non-zero value and it determines the length of that cross product. $$\frac{\sqrt{1 + ||\color{deeppink}{e}||^{2} - 2(\color{red}{r}_{1}\color{pink}{e}_{1} + \color{red}{r}_{2}\color{pink}{e}_{2})}}{\sqrt{1 + ||\color{gold}{p}||^{2} - 2(\color{red}{r}_{1}\color{gold}{p}_{1} + \color{red}{r}_{2}\color{gold}{p}_{2})}}= \frac{\color{red}{r}_{1}\color{gold}{p}_{2}-\color{red}{r}_{2}\color{gold}{p}_{1}}{\color{deeppink}{e}_{1}\color{red}{r}_{2}-\color{deeppink}{e}_{2}\color{red}{r}_{1}}$$ The left hand side is positive. I want both angles to have equal signs, because they are both in the same general direction as $\color{red}{r}$. This means that the right hand side is also positive. I square both sides to get rid of the roots on the left. $$\frac{1 + ||\color{deeppink}{e}||^{2} - 2(\color{red}{r}_{1}\color{deeppink}{e}_{1} + \color{red}{r}_{2}\color{pink}{e}_{2})}{1 + ||\color{gold}{p}||^{2} - 2(\color{red}{r}_{1}\color{gold}{p}_{1} + \color{red}{r}_{2}\color{gold}{p}_{2})}= \frac{(\color{red}{r}_{1}\color{gold}{p}_{2}-\color{red}{r}_{2}\color{gold}{p}_{1})^{2}}{(\color{deeppink}{e}_{1}\color{red}{r}_{2}-\color{deeppink}{e}_{2}\color{red}{r}_{1})^{2}}$$ And this is as far as I can get. I don't think substituting $p_{2}$ with $\sqrt{1-p_{1}^2}$ will get me anywhere and also avoid it, since I am going to use that in a graphic application and I am trying to make as few calculations as possible. Hope I explained the problem well enough and followed all conventions on this forum :) This is my first post here, so I accept any critique on it.

• Welcome to MSE, what about an approach by analytical geometry? – user Feb 12 '18 at 22:39
• I am not a matematitian, so I'm not sure what you mean. Can you give me an idea what to read? I opened up the Wiki page for analytical geometry, but I'm not sure I'm ready to start learning an entirely new to me field of math for this one problem :D – Орлин Митев Feb 12 '18 at 22:47
• I did some investigation myself, using a computer algebra system. It appears that the problem involves solving quartics, and trying this even for some integer points produces expressions that are not at all simple. While an analytic solution appears to be possible, it looks likely that such a solution will be nasty enough not to be worth investigating or implementing. :-( – Theo Bendit Feb 12 '18 at 23:55
• Yeah it seemed to me like something, that could be elegantly computed, but maybe that is not the case :/ Thanks a lot though! At least now I know I'm not just lacking the knowledge to solve it :) – Орлин Митев Feb 13 '18 at 0:57
• BTW, even though this was your first post you were very clear, LaTeX'ed your equations, and even provided diagrams. Nicely done! – JonathanZ supports MonicaC Feb 16 '18 at 17:37

## 4 Answers

Made a program which seems to work very well, i.e. it has passed all the examples I have tried. It uses no iteration but does require solving a cubic polynomial. I'll provide the code if you are interested.

My idea is the following: If I rotate the points E and P around the circle center, then at some rotation angle $\alpha$ the point of reflection must be at the top of the circle, i.e. at the point $(0,1)$. And a way to know I have found the right rotation angle is if the slope of a line from E' to $(0,1)$ is the negative equal of the slope from P' to $(0,1)$. If I could find the angle in this way, then the point of reflection, in the original setup, would just be located at an angle of $\frac{\pi}{2} - \alpha$. So, if we let $E=(a,b)$ and $P=(c,d)$, then the rotated position of $E$, call it $E'=(a_1, b_1)$, would be:$$a_1 = a \cos(\alpha)-b \sin(\alpha)$$ $$b_1 = a \sin(\alpha)+b \cos(\alpha)$$ Same goes for $P'=(c_1,d_1)$.

If the slopes from $E'$ and $P'$ to (0,1) are negative equals, then $$\frac{b_1-1}{a_1} = -\frac{d_1-1}{c_1}$$

Inserting the equations for $a_1, b_1,$ etc into the equation above and simplifying, we get $$(ac-bd)\sin(2\alpha) + (bc+ad) \cos(2\alpha)+(b+d) \sin(\alpha) -(a+c)\cos(\alpha) = 0$$ Not so easy to solve. But, I decided to try an approximation, simply substituting the trig functions with the first two terms of their Taylor series. These are quite accurate up to about $40^\circ$, which is enough as I'll later talk about.

Inserting the substitutions and simplifying we get $$-\frac{8e+g}{6}\alpha^3 + (\frac{h}{2}-2f)\alpha^2+(2e+g)\alpha+(f-h) = 0$$ where $e=(ac-bd), f=(bc+ad), g=(b+d)$ and $h=(a+c)$.

The roots of this cubic equation can now be found in the usual manner.

This was my first version of the program and it worked somewhat, but could be off by $5-6 ^\circ$ which is too much. The problem is that the approximations get inaccurate for large angles. The solution I found was to do an initial rotation of E and P such that the angle bisector of $\angle E'CP'$ is at $90^\circ$. This puts E and P close to where they should be and means the final adjustment angle to be found via the cubic equation, is small.

Based on the examples I've tried (and checked via GeoGebra drawings), the accuracy of the found angle is now close to 2 decimal places. An example below: The input was $E=(-1,-1.5), P=(-2,2)$ and the ouput was the angle $-168.05^\circ$.

• What I find interesting about your work, @Jens, is that your final equation (if you substitute $\sin 2x = 2\sin x \cos x$ and $\cos 2x = \cos^2 x - \sin ^2 x) ends up with the same degree of complexity as mine, a "trig quadratic". I can't tell from what you posted, but can you tell me if the double angles arose naturally from your initial slopes equation, or did you get it into that form by using the reverse of the double angle formulas? – JonathanZ supports MonicaC Feb 19 '18 at 0:40 • I used the reverse of the double angle formulas. Have you managed to make a working program? – Jens Feb 19 '18 at 15:26 • No, I've stopped at the point of getting that equation. I'm a little more interested in looking for general solutions, as opposed to your practical approximations which are probably more useful to the original poster. :--} – JonathanZ supports MonicaC Feb 19 '18 at 16:08 • @JonathanZ - Sorry, I simply assumed you were the OP without checking. :-) – Jens Feb 19 '18 at 16:25 I'm going to report the result I got from some quick calculations, in case someone can take them further. The result is kind of neat, but it doesn't seem to lead to any nice closed answer so the OP, who probably values speed of calculation over absolute accuracy, should probably go with other methods. Call the points$p$and$e$as$\vec a$and$\vec b$. (It doesn't matter which is which, as the configuration is symmetric in$\vec a$and$\vec b$.) Let$\vec u$be the unit vector to the point on the circle and$r$be the radius of the circle. I calculated the direction of the line you get by bouncing$\vec a - r\vec u$off the point$r\vec u$, and then wrote the equation requiring$\vec b$to lie on the line in that direction passing through$r \vec u$. If we write$\vec v ^\perp$for the vector$\vec v $rotated 90 degrees (counter-clockwise, to pick a definite direction) the equation simplifies to $$2(\vec a\cdot \vec u)(\vec b \cdot \vec u ^\perp) - r(\vec b \cdot \vec u ^\perp) + r(\vec a ^ \perp \cdot \vec u) = \vec b \cdot \vec a ^ \perp.$$ If you were to write$\vec u$as$(\cos\theta, \sin \theta)$and plug it in to the equation you'd get a "trigonometric quadratic", i.e. it would have$\cos$,$\sin$,$\cos^2$,$\cos\cdot\sin$, and$\sin^2$terms. I'm sure you could find an algorithm for solving these types of equations, but I'm also guessing it'd be to messy to be useful in your situation. I have made a Geogebra model (https://www.geogebra.org/m/r75cgAfw ) that works according to your original equations and seems to give the right answers as far as tested. Starting with your originals: The angle between$r$and$p-r$should be equal to the angle between$r$and$e-r$$$acos\left(\frac{\mathbf{r-(p-r)}}{\mathbf{|r|\cdot|p-r|}}\right)=acos\left(\frac{\mathbf{r-(e-r)}}{\mathbf{|r|\cdot|e-r|}}\right)$$ $$|r|=1$$ Next I need to simplify some, so 1) expand the numerator and 2)replace both |r| and r\cdot r with 1, and 3)drop the acos from both sides. to get: $$\mathbf{\frac{r\cdot p-1}{|p-r|}}=\mathbf{\frac{r\cdot e-1}{|e-r|}}$$ Next, I substitute the vector (x,y) for$\mathbf r$remembering that both$\mathbf e$and$\mathbf p$are known vectors. $$\frac{p_{x}x+p_{y}y-1}{\sqrt{(p_{x}-x)^{2}+(p_{y}-y)^{2}}}=\frac{e_{x}x+e_{y}y-1}{\sqrt{(e_{x}-x)+(e_{y}-y)^{2}}}$$ Finally, substitue in$x=\sqrt{-y^2+1}\$ for x. That creates a minor issue with where on the circle the reflection is, since it causes 2 answers, one on each side of the circle. Still, the CAS equation to Geogebra is $$\frac{p_{y} y + \sqrt{-y^{2} + 1} p_{x} - 1}{\sqrt{{\left(p_{x} -\sqrt{-y^{2} + 1}\right)}^{2} + {\left(p_{y} - y\right)}^{2}}} =\frac{e_{y} y + \sqrt{-y^{2} + 1} e_{x} - 1}{\sqrt{{\left(e_{x} -\sqrt{-y^{2} + 1}\right)}^{2} + {\left(e_{y} - y\right)}^{2}}}$$ which is 1 equation in 1 unknown, thus allowing a relatively fast numeric solution. Fast enough that one can move around points E and P and get various solutions.

• The vector that I'm looking for will coincide with the bisector of ECP only if E and P are equally distanced from the center. And in the common case they're not :( – Орлин Митев Feb 14 '18 at 23:16
• Take a look at this Geogebra sheet. Grab P with you mouse and move it around. geogebra.org/worksheet/edit/id/E5gCDg48 – Narlin Feb 15 '18 at 2:49
• Cannot open the link :( Can you tell me the name of the material, so I can look it up myself? – Орлин Митев Feb 15 '18 at 7:52
• Try the link again. I changed it from shared to public. The material is named Mild Vectors, but it may not be google indexed. – Narlin Feb 15 '18 at 12:53
• Yeah i found it :) But sadly that's really not what I'm looking for. I need the vector b to be the bisector of EBP, which it is not. Try putting P at e.g. (1, 1) and E at (-2, 2) and you can clearly see, that w is not the direction in which a ray from E would be reflected off the circle. But thanks a lot! It's really cool that you made the effort of making this graphic. – Орлин Митев Feb 15 '18 at 22:17

Fermat’s principle suggests an approach to this problem. The idea is to find the point $$\mathbf r$$ on the circle that minimizes the distance $$\|\mathbf r-\mathbf p\|+\|\mathbf r-\mathbf e\|.$$ Using the parameterization $$\mathbf r=(r\cos t,r\sin t)$$, this distance is $$\sqrt{(\cos t-p_x)^2+(\sin t-p_y)^2}+\sqrt{(\cos t-e_x)^2+(\sin t-e_y)^2}$$ and the problem becomes that of solving $${p_x\sin t-p_y\cos t \over \sqrt{(\cos t-p_x)^2+(\sin t-p_y)^2}}+{e_x\sin t-e_y\cos t \over \sqrt{(\cos t-e_x)^2+(\sin t-e_y)^2}} = 0.$$ Unfortunately, this doesn’t appear to lend itself to a nice (or even not so nice) closed-form solution. Rahul’s suggestion of binary search or perhaps a Newton’s method approximation to the root of the equation seems best for a practical solution that is to be coded. The angle bisector of $$\mathbf p$$ and $$\mathbf e$$ is a reasonable first guess for $$t$$ (equiv. $$\mathbf r$$) for any such iterative approximation.

A perhaps more promising approach that comes to mind is that the point $$\mathbf r$$ is such that $$\mathbf p$$, $$\mathbf r$$ and the reflection of $$\mathbf e$$ in the tangent at $$\mathbf r$$ are colinear. (By symmetry, you could reflect $$\mathbf p$$ instead, of course.) A normal vector to the tangent at $$\mathbf r$$ is $$\mathbf r$$ itself, so the reflection of $$\mathbf e$$ is $$\mathbf e' = \mathbf e + 2(1-\mathbf e\cdot\mathbf r)\,\mathbf r.$$ The colinearity condition can be expressed via the vanishing of the scalar triple product (in homogeneous coordinates) $$[\mathbf r,\mathbf p,\mathbf e']=0$$. Expanding by coordinates, this translates into the equation $$(p_x e_y+p_y e_x)(x^2-y^2) - 2(p_x e_x-p_y e_y)xy - (p_y+e_y)x + (p_x+e_x)y = 0. \tag{*}$$ Together with the constraint $$x^2+y^2=1$$, this makes the problem that of computing the intersection(s) of this hyperbola with the unit circle. There’s an unavoidable cubic equation that needs to be solved, but the entire computation can be mechanized via an algorithm described in Richter-Gebert’s Perspectives on Projective Geometry and sketched out in this answer. Since one of the conics is a unit circle, there are likely some simplifications that can be made to the general computation. There can be up to four intersections, so you’ll still need to select the correct one, but equation (*) is sufficient to mock up a solution in a graphical tool such as GeoGebra.