1---
2title: 'Two-point Conical Gradient'
3linkTitle: 'Two-point Conical Gradient'
4---
5
6<script type="text/x-mathjax-config">
7MathJax.Hub.Config({
8    tex2jax: {
9        inlineMath: [['$','$'], ['\\(','\\)']]
10    }
11});
12</script>
13
14<script src='https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.2/MathJax.js?config=TeX-MML-AM_CHTML'></script>
15
16(Please refresh the page if you see a lot of dollars instead of math symbols.)
17
18We present a fast shading algorithm (compared to bruteforcely solving the
19quadratic equation of gradient $t$) for computing the two-point conical gradient
20(i.e., `createRadialGradient` in
21[spec](https://html.spec.whatwg.org/multipage/canvas.html#dom-context-2d-createradialgradient)).
22It reduced the number of multiplications per pixel from ~10 down to 3, and
23brought a speedup of up to 26% in our nanobenches.
24
25This document has 3 parts:
26
271. [Problem Statement and Setup](#problem-statement)
282. [Algorithm](#algorithm)
293. [Appendix](#appendix)
30
31Part 1 and 2 are self-explanatory. Part 3 shows how to geometrically proves our
32Theorem 1 in part 2; it's more complicated but it gives us a nice picture about
33what's going on.
34
35## <span id="problem-statement">Problem Statement and Setup</span>
36
37Let two circles be $C_0, r_0$ and $C_1, r_1$ where $C$ is the center and $r$ is
38the radius. For any point $P = (x, y)$ we want the shader to quickly compute a
39gradient $t \in \mathbb R$ such that $p$ is on the linearly interpolated circle
40with center $C_t = (1-t) \cdot C_0 + t \cdot C_1$ and radius
41$r_t = (1-t) \cdot r_0 + t \cdot r_1 > 0$ (note that radius $r_t$ has to be
42_positive_). If there are multiple (at most 2) solutions of $t$, choose the
43bigger one.
44
45There are two degenerated cases:
46
471. $C_0 = C_1$ so the gradient is essentially a simple radial gradient.
482. $r_0 = r_1$ so the gradient is a single strip with bandwidth $2 r_0 = 2 r_1$.
49
50<!-- TODO maybe add some fiddle or images here to illustrate the two degenerated cases -->
51
52They are easy to handle so we won't cover them here. From now on, we assume
53$C_0 \neq C_1$ and $r_0
54\neq r_1$.
55
56As $r_0 \neq r_1$, we can find a focal point
57$C_f = (1-f) \cdot C_0 + f \cdot C_1$ where its corresponding linearly
58interpolated radius $r_f = (1-f) \cdot r_0 + f \cdot r_1 = 0$. Solving the
59latter equation gets us $f = r_0 / (r_0 - r_1)$.
60
61As $C_0 \neq C_1$, focal point $C_f$ is different from $C_1$ unless $r_1 = 0$.
62If $r_1 = 0$, we can swap $C_0, r_0$ with $C_1, r_1$, compute swapped gradient
63$t_s$ as if $r_1 \neq 0$, and finally set $t = 1 - t_s$. The only catch here is
64that with multiple solutions of $t_s$, we shall choose the smaller one (so $t$
65could be the bigger one).
66
67Assuming that we've done swapping if necessary so $C_1 \neq C_f$, we can then do
68a linear transformation to map $C_f, C_1$ to $(0, 0), (1, 0)$. After the
69transformation:
70
711. All centers $C_t = (x_t, 0)$ must be on the $x$ axis
722. The radius $r_t$ is $x_t r_1$.
733. Given $x_t$ , we can derive $t = f + (1 - f) x_t$
74
75From now on, we'll focus on how to quickly computes $x_t$. Note that $r_t > 0$
76so we're only interested positive solution $x_t$. Again, if there are multiple
77$x_t$ solutions, we may want to find the bigger one if $1 - f > 0$, and smaller
78one if $1 - f < 0$, so the corresponding $t$ is always the bigger one (note that
79$f \neq 1$, otherwise we'll swap $C_0, r_0$ with $C_1, r_1$).
80
81## <span id="algorithm">Algorithm</span>
82
83**Theorem 1.** The solution to $x_t$ is
84
851. $\frac{x^2 + y^2}{(1 + r_1) x} = \frac{x^2 + y^2}{2 x}$ if $r_1 = 1$
862. $\left(\sqrt{(r_1^2 - 1) y ^2 + r_1^2 x^2}  - x\right) / (r_1^2 - 1)$ if
87   $r_1 > 1$
883. $\left(\pm \sqrt{(r_1^2 - 1) y ^2 + r_1^2 x^2}  - x\right) / (r_1^2 - 1)$ if
89   $r_1 < 1$.
90
91Case 2 always produces a valid $x_t$. Case 1 and 3 requires $x > 0$ to produce
92valid $x_t > 0$. Case 3 may have no solution at all if
93$(r_1^2 - 1) y^2 + r_1^2 x^2 < 0$.
94
95_Proof._ Algebriacally, solving the quadratic equation
96$(x_t - x)^2 + y^2 = (x_t r_1)^2$ and eliminate negative $x_t$ solutions get us
97the theorem.
98
99Alternatively, we can also combine Corollary 2., 3., and Lemma 4. in the
100Appendix to geometrically prove the theorem. $\square$
101
102Theorem 1 by itself is not sufficient for our shader algorithm because:
103
1041. we still need to compute $t$ from $x_t$ (remember that $t = f + (1-f) x_t$);
1052. we still need to handle cases of choosing the bigger/smaller $x_t$;
1063. we still need to handle the swapped case (we swap $C_0, r_0$ with $C_1, r_1$
107   if $r_1 = 0$);
1084. there are way too many multiplications and divisions in Theorem 1 that would
109   slow our shader.
110
111Issue 2 and 3 are solved by generating different shader code based on different
112situations. So they are mainly correctness issues rather than performance
113issues. Issue 1 and 4 are performance critical, and they will affect how we
114handle issue 2 and 3.
115
116The key to handle 1 and 4 efficiently is to fold as many multiplications and
117divisions into the linear transformation matrix, which the shader has to do
118anyway (remember our linear transformation to map $C_f, C_1$ to
119$(0, 0), (1, 0)$).
120
121For example, let $\hat x, \hat y = |1-f|x, |1-f|y$. Computing $\hat x_t$ with
122respect to $\hat x,
123\hat y$ allow us to have
124$t = f + (1 - f)x_t = f + \text{sign}(1-f) \cdot \hat x_t$. That saves us one
125multiplication. Applying similar techniques to Theorem 1 gets us:
126
1271. If $r_1 = 1$, let $x' = x/2,~ y' = y/2$, then $x_t = (x'^2 + y'^2) / x'$.
1282. If $r_1 > 1$, let
129   $x' = r_1 / (r_1^2 - 1) x,~ y' = \frac{\sqrt{r_1^2 - 1}}{r_1^2 - 1} y$, then
130   $x_t = \sqrt{x'^2 + y'^2} - x' / r_1$
1313. If $r_1 < 1$, let
132   $x' = r_1 / (r_1^2 - 1) x,~ y' = \frac{\sqrt{1 - r_1^2}}{r_1^2 - 1} y$, then
133   $x_t = \pm\sqrt{x'^2 - y'^2} - x' / r_1$
134
135Combining it with the swapping, the equation $t = f + (1-f) x_t$, and the fact
136that we only want positive $x_t > 0$ and bigger $t$, we have our final
137algorithm:
138
139**Algorithm 1.**
140
1411. Let $C'_0, r'_0, C'_1, r'_1 = C_0, r_0, C_1, r_1$ if there is no swapping and
142   $C'_0,
143    r'_0, C'_1, r'_1 = C_1, r_1, C_0, r_0$ if there is swapping.
1442. Let $f = r'_0 / (r'_0 - r'_1)$ and $1 - f = r'_1 / (r'_1 - r'_0)$
1453. Let $x' = x/2,~ y' = y/2$ if $r_1 = 1$, and
146   $x' = r_1 / (r_1^2 - 1) x,~ y' = \sqrt{|r_1^2 - 1|} / (r_1^2 - 1) y$ if
147   $r_1 \neq 1$
1484. Let $\hat x = |1 - f|x', \hat y = |1 - f|y'$
1495. If $r_1 = 1$, let $\hat x_t = (\hat x^2 + \hat y^2) / \hat x$
1506. If $r_1 > 1$, let $\hat x_t = \sqrt{\hat x^2 + \hat y^2} - \hat x / r_1$
1517. If $r_1 < 1$
1528. return invalid if $\hat x^2 - \hat y^2 < 0$
1539. let $\hat x_t =  -\sqrt{\hat x^2 - \hat y^2} - \hat x / r_1$ if we've swapped
154   $r_0, r_1$, or if $1 - f < 0$
155
15610. let $\hat x_t =  \sqrt{\hat x^2 - \hat y^2} - \hat x / r_1$ otherwise
157
15811. $t$ is invalid if $\hat x_t < 0$ (this check is unnecessary if $r_1 > 1$)
15912. Let $t = f + \text{sign}(1 - f) \hat x_t$
16013. If swapped, let $t = 1 - t$
161
162In step 7, we try to select either the smaller or bigger $\hat x_t$ based on
163whether the final $t$ has a negative or positive relationship with $\hat x_t$.
164It's negative if we've swapped, or if $\text{sign}(1 - f)$ is negative (these
165two cannot both happen).
166
167Note that all the computations and if decisions not involving $\hat x, \hat y$
168can be precomputed before the shading stage. The two if decisions
169$\hat x^2 - \hat y^2 < 0$ and $\hat x^t < 0$ can also be omitted by precomputing
170the shading area that never violates those conditions.
171
172The number of operations per shading is thus:
173
174- 1 addition, 2 multiplications, and 1 division if $r_1 = 1$
175- 2 additions, 3 multiplications, and 1 sqrt for $r_1 \neq 1$ (count subtraction
176  as addition; dividing $r_1$ is multiplying $1/r_1$)
177- 1 more addition operation if $f \neq 0$
178- 1 more addition operation if swapped.
179
180In comparison, for $r_1 \neq 1$ case, our current raster pipeline shading
181algorithm (which shall hopefully soon be upgraded to the algorithm described
182here) mainly uses formula
183$$t = 0.5 \cdot
184(1/a) \cdot \left(-b \pm \sqrt{b^2 - 4ac}\right)$$ It precomputes
185$a = 1 - (r_1 - r_0)^2, 1/a, r1 -
186r0$. Number
187$b = -2 \cdot (x + (r1 - r0) \cdot r0)$ costs 2 multiplications and 1 addition.
188Number $c = x^2 + y^2 - r_0^2$ costs 3 multiplications and 2 additions. And the
189final $t$ costs 5 more multiplications, 1 more sqrt, and 2 more additions.
190That's a total of 5 additions, 10 multiplications, and 1 sqrt. (Our algorithm
191has 2-4 additions, 3 multiplications, and 1 sqrt.) Even if it saves the
192$0.5 \cdot (1/a), 4a, r_0^2$ and $(r_1 - r_0) r_0$ multiplications, there are
193still 6 multiplications. Moreover, it sends in 4 unitofmrs to the shader while
194our algorithm only needs 2 uniforms ($1/r_1$ and $f$).
195
196## <span id="appendix">Appendix</span>
197
198**Lemma 1.** Draw a ray from $C_f = (0, 0)$ to $P = (x, y)$. For every
199intersection points $P_1$ between that ray and circle $C_1 = (1, 0), r_1$, there
200exists an $x_t$ that equals to the length of segment $C_f P$ over length of
201segment $C_f P_1$. That is, $x_t = || C_f P || / ||C_f P_1||$
202
203_Proof._ Draw a line from $P$ that's parallel to $C_1 P_1$. Let it intersect
204with $x$-axis on point $C = (x', y')$.
205
206<img src="./lemma1.svg"/>
207
208Triangle $\triangle C_f C P$ is similar to triangle $\triangle C_f C_1 P_1$.
209Therefore $||P C|| = ||P_1 C_1|| \cdot (||C_f C|| / ||C_f C_1||) = r_1 x'$. Thus
210$x'$ is a solution to $x_t$. Because triangle $\triangle C_f C P$ and triangle
211$\triangle C_f C_1 P_1$ are similar,
212$x'
213= ||C_f C_1|| \cdot (||C_f P|| / ||C_f P_1||) = ||C_f P|| / ||C_f P_1||$.
214$\square$
215
216**Lemma 2.** For every solution $x_t$, if we extend/shrink segment $C_f P$ to
217$C_f P_1$ with ratio $1 / x_t$ (i.e., find $P_1$ on ray $C_f P$ such that
218$||C_f P_1|| / ||C_f P|| = 1 / x_t$), then $P_1$ must be on circle $C_1, r_1$.
219
220_Proof._ Let $C_t = (x_t, 0)$. Triangle $\triangle C_f C_t P$ is similar to
221$C_f C_1 P_1$. Therefore $||C_1 P_1|| = r_1$ and $P_1$ is on circle $C_1, r_1$.
222$\square$
223
224**Corollary 1.** By lemma 1. and 2., we conclude that the number of solutions
225$x_t$ is equal to the number of intersections between ray $C_f P$ and circle
226$C_1, r_1$. Therefore
227
228- when $r_1 > 1$, there's always one unique intersection/solution; we call this
229  "well-behaved"; this was previously known as the "inside" case;
230- when $r_1 = 1$, there's either one or zero intersection/solution (excluding
231  $C_f$ which is always on the circle); we call this "focal-on-circle"; this was
232  previously known as the "edge" case;
233
234<img src="./corollary2.2.1.svg"/>
235<img src="./corollary2.2.2.svg"/>
236
237- when $r_1 < 1$, there may be $0, 1$, or $2$ solutions; this was also
238  previously as the "outside" case.
239
240<img src="./corollary2.3.1.svg" width="30%"/>
241<img src="./corollary2.3.2.svg" width="30%"/>
242<img src="./corollary2.3.3.svg" width="30%"/>
243
244**Lemma 3.** When solution exists, one such solution is
245
246$$
247    x_t = {|| C_f P || \over ||C_f P_1||} = \frac{x^2 + y^2}{x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}}
248$$
249
250_Proof._ As $C_f = (0, 0), P = (x, y)$, we have $||C_f P|| = \sqrt(x^2 + y^2)$.
251So we'll mainly focus on how to compute $||C_f P_1||$.
252
253**When $x \geq 0$:**
254
255<img src="./lemma3.1.svg"/>
256
257Let $X_P = (x, 0)$ and $H$ be a point on $C_f P_1$ such that $C_1 H$ is
258perpendicular to $C_1
259P_1$. Triangle $\triangle C_1 H C_f$ is similar to triangle
260$\triangle P X_P C_f$. Thus
261$$||C_f H|| = ||C_f C_1|| \cdot (||C_f X_P|| / ||C_f P||) = x / \sqrt{x^2 + y^2}$$
262$$||C_1 H|| = ||C_f C_1|| \cdot (||P X_P|| / ||C_f P||) = y / \sqrt{x^2 + y^2}$$
263
264Triangle $\triangle C_1 H P_1$ is a right triangle with hypotenuse $r_1$. Hence
265$$ ||H P_1|| = \sqrt{r_1^2 - ||C_1 H||^2} = \sqrt{r_1^2 - y^2 / (x^2 + y^2)} $$
266
267We have \begin{align} ||C_f P_1|| &= ||C_f H|| + ||H P_1|| \\\\\\ &= x /
268\sqrt{x^2 + y^2} + \sqrt{r_1^2 - y^2 / (x^2 + y^2)} \\\\\\ &= \frac{x +
269\sqrt{r_1^2 (x^2 + y^2) - y^2}}{\sqrt{x^2 + y^2}} \\\\\\ &= \frac{x +
270\sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}}{\sqrt{x^2 + y^2}} \end{align}
271
272**When $x < 0$:**
273
274Define $X_P$ and $H$ similarly as before except that now $H$ is on ray $P_1 C_f$
275instead of $C_f P_1$.
276
277<img src="./lemma3.2.svg"/>
278
279As before, triangle $\triangle C_1 H C_f$ is similar to triangle
280$\triangle P X_P C_f$, and triangle $\triangle C_1 H P_1$ is a right triangle,
281so we have
282$$||C_f H|| = ||C_f C_1|| \cdot (||C_f X_P|| / ||C_f P||) = -x / \sqrt{x^2 + y^2}$$
283$$||C_1 H|| = ||C_f C_1|| \cdot (||P X_P|| / ||C_f P||) = y / \sqrt{x^2 + y^2}$$
284$$ ||H P_1|| = \sqrt{r_1^2 - ||C_1 H||^2} = \sqrt{r_1^2 - y^2 / (x^2 + y^2)} $$
285
286Note that the only difference is changing $x$ to $-x$ because $x$ is negative.
287
288Also note that now $||C_f P_1|| = -||C_f H|| + ||H P_1||$ and we have
289$-||C_f H||$ instead of $||C_f H||$. That negation cancels out the negation of
290$-x$ so we get the same equation of $||C_f P_1||$ for both $x \geq 0$ and
291$x < 0$ cases:
292
293$$
294    ||C_f P_1|| = \frac{x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}}{\sqrt{x^2 + y^2}}
295$$
296
297Finally
298
299$$
300    x_t = \frac{||C_f P||}{||C_f P_1||} = \frac{\sqrt{x^2 + y^2}}{||C_f P_1||}
301        = \frac{x^2 + y^2}{x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}}
302$$ $\square$
303
304**Corollary 2.** If $r_1 = 1$, then the solution $x_t = \frac{x^2 + y^2}{(1 + r_1) x}$, and
305it's valid (i.e., $x_t > 0$) iff $x > 0$.
306
307*Proof.* Simply plug $r_1 = 1$ into the formula of Lemma 3. $\square$
308
309**Corollary 3.** If $r_1 > 1$, then the unique solution is
310$x_t = \left(\sqrt{(r_1^2 - 1) y ^2 + r_1^2 x^2}  - x\right) / (r_1^2 - 1)$.
311
312*Proof.* From Lemma 3., we have
313
314\begin{align}
315    x_t &= \frac{x^2 + y^2}{x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}} \\\\\\
316        &=  {
317                (x^2 + y^2) \left ( -x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2} \right )
318            \over
319                \left (x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2} \right )
320                \left (-x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2} \right )
321            } \\\\\\
322        &=  {
323                (x^2 + y^2) \left ( -x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2} \right )
324            \over
325                -x^2 + (r_1^2 - 1) y^2 + r_1^2 x^2
326            } \\\\\\
327        &=  {
328                (x^2 + y^2) \left ( -x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2} \right )
329            \over
330                (r_1^2 - 1) (x^2 + y^2)
331            } \\\\\\
332        &=  \left(\sqrt{(r_1^2 - 1) y ^2 + r_1^2 x^2}  - x\right) / (r_1^2 - 1)
333\end{align}
334
335The transformation above (multiplying $-x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}$ to enumerator and
336denomenator) is always valid because $r_1 > 1$ and it's the unique solution due to Corollary 1.
337$\square$
338
339**Lemma 4.** If $r_1 < 1$, then
340
3411. there's no solution to $x_t$ if $(r_1^2 - 1) y^2 + r_1^2 x^2 < 0$
3422. otherwise, the solutions are
343    $x_t = \left(\sqrt{(r_1^2 - 1) y ^2 + r_1^2 x^2}  - x\right) / (r_1^2 - 1)$,
344    or
345    $x_t = \left(-\sqrt{(r_1^2 - 1) y ^2 + r_1^2 x^2}  - x\right) / (r_1^2 - 1)$.
346
347(Note that solution $x_t$ still has to be nonnegative to be valid; also note that
348$x_t > 0 \Leftrightarrow x > 0$ if the solution exists.)
349
350*Proof.* Case 1 follows naturally from Lemma 3. and Corollary 1.
351
352<img src="./lemma4.svg"/>
353
354For case 2, we notice that $||C_f P_1||$ could be
355
3561. either $||C_f H|| + ||H P_1||$ or $||C_f H|| - ||H P_1||$ if $x \geq 0$,
3572. either $-||C_f H|| + ||H P_1||$ or $-||C_f H|| - ||H P_1||$ if $x < 0$.
358
359By analysis similar to Lemma 3., the solution to $x_t$ does not depend on the sign of $x$ and
360they are either $\frac{x^2 + y^2}{x + \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}}$
361or $\frac{x^2 + y^2}{x - \sqrt{(r_1^2 - 1) y^2 + r_1^2 x^2}}$.
362
363As $r_1 \neq 1$, we can apply the similar transformation in Corollary 3. to get the two
364formula in the lemma.
365$\square$
366
367
368$$
369