5. Barycentric Coordinates
5.1. How to determine if a point is inside or outside a triangle?
Given a triangle on a 2d-plane with vertices [x1,y1], [x2,y2] and [x3,y3], and given a point [x,y], how can you find out whether the point is inside or outside the triangle? Before we address how to solve this problem, think about how you would approach the it for a minute or two.
This problem seemingly has no connection with linear algebra but, as we will see below, can be solved elegantly with the tools of linear algebra that we have learnt so far.
5.1.1. Let's start with an example
Let's say the vertices of the triangle are [3,0], [2,3] and [-4,-2]. And we want to know whether or not the point [3,0] (shown in amber color) is inside the triangle.
We will also keep implementing a function to do the same as we make progress:
def is_point_in_triangle(vertices, point):
'''
vertices: [[x1,y1],[x2,y2],[x3,y3]]
point: [p,q]
'''
...
5.2. Move (or translate) the triangle and the point to a better location
Consider the two cases below: the only difference between them is that in the second case, one of the vertices of the triangle ([3,0] in the left figure) lies on the origin. Everything else is the same - the position of the point relative to the vertices of the triangle is still the same.
We want to find whether or not the point [2,1] lies in the triangle with vertices [3,0], [2,3] and [-4,-2]
We want to find whether or not the point [-1,1] lies in the triangle with vertices [0,0], [-1,3] and [-7,-2]
These two problems are actually the same problem. But it is easier to find the answer of the second problem, which is also the answer of the first problem. So we will now focus on how to solve the second problem.
5.2.1. How to go from problem on the left to problem on the right?
This is easy: if we want to move the vertex [x,y] to the origin, we subtract [x,y] from all three vertices and the point.
def is_point_in_triangle(vertices, point):
'''
vertices: [[x1,y1],[x2,y2],[x3,y3]]
point: [p,q]
'''
x1,y1 = vertices[0]
x2,y2 = vertices[1]
x3,y3 = vertices[2]
shifted_vertices = [[x1-x1,y1-y1],[x2-x1,y2-y1],[x3-x1,y3-y1]]
shifted_point = [point[0]-x1, point[1]-y1]5.3. First change in viewpoint: look at sides of the triangle as vectors
Let's name the vertices and the point first: A = [0,0], B = [-1,3], C = [-7,-2] and P = [-1,1].
We want to look at sides of the triangle AB, AC and the point P as vectors. Next we would like to get to point P by taking weighted sums of vectors AB and AC.
The side AB can be seen as a vector \(
\begin{bmatrix}
-1 \\ 3
\end{bmatrix}
\) in the 2-dimensional space. Similarly AC can be seen as a vector \(
\begin{bmatrix}
-7 \\ -2
\end{bmatrix}
\). The point P is the vector \(
\begin{bmatrix}
-1 \\ 1
\end{bmatrix}
\)
This shift in perspective can be shown visually as well:
Vertices of the triangle ABC and the input point P
Sides of a triangle AB, AC and the point Pwhen seen as vectors in a 2d space
5.3.1. Let's walk in the 2d space by taking weighted sums of these two vectors
Remember that this walk can be represented algebriacally as a matrix-vector multiplication where the matrix is \(
\begin{bmatrix}
-1 & -7 \\ 3 & -2
\end{bmatrix}
\) and the vector is the weights for each vector \(
\begin{bmatrix}
w_{1} \\ w_{2}
\end{bmatrix}
\).
We want to observe for what values of \(w_{1}, w_{2}\) does the output point lie inside the triangle. It is this insight that leads to the solution. Before we reveal the pattern for the weights, let's play with it a bit and get a feel for it:
We take weighted sum of the vectors AB and AC.
Change \(w_{1}\): 0.2
Change \(w_{2}\): 0.2
\(
\begin{bmatrix}
-1 & -7 \\ 3 & -2
\end{bmatrix}
\)\(
\begin{bmatrix}
0.2 \\ 0.2
\end{bmatrix}
\) = \(
\begin{bmatrix}
-1.6 \\ 0.2
\end{bmatrix}
\)
Output point is inside the triangle!
Let's also write a function for matrix multiplication:
def scale_vector(scale, vector):
return [i*scale for i in vector]
def add_vectors(array_of_vectors):
dimensionality = len(array_of_vectors[0])
output = [0]*dimensionality
for vector in array_of_vectors:
output = [i+j for (i,j) in zip(output, vector)]
return output
def matMul(columns_of_matrix, weights):
scaled_columns = [scale_vector(scale, vector)
for (scale, vector) in zip(weights, columns_of_matrix)]
return add_vectors(scaled_columns)5.3.2. A glimpse of the solution
Here are some observations about the output point of the matrix multiplication \(
\begin{bmatrix}
-1 & -7 \\ 3 & -2
\end{bmatrix}
\)\(
\begin{bmatrix}
w_{1} \\ w_{2}
\end{bmatrix}
\).
- It lies inside the triangle only when \(w_{1} > 0, w_{2}> 0\) and \(w_{1} + w_{2} <= 1\)
- When \(w_{1}= 0, w_{2} = 0\), the input vector \(
\begin{bmatrix}
0 \\ 0
\end{bmatrix}
\) gives as output the point A (which is at the origin).
- When \(w_{1}= 1, w_{2} = 0\), the input vector \(
\begin{bmatrix}
1 \\ 0
\end{bmatrix}
\) gives as output the point B.
- When \(w_{1}= 0, w_{2} = 1\), the input vector \(
\begin{bmatrix}
0 \\ 1
\end{bmatrix}
\) gives as output the point C.
Now the question is, if we want to find out whether the point \(
\begin{bmatrix}
-1 \\ 1
\end{bmatrix}
\) lies inside the triangle or not, we need to find out which values of \(w_{1}, w_{2}\) does \(
\begin{bmatrix}
-1 \\ 1
\end{bmatrix}
\) = \(
\begin{bmatrix}
-1 & -7 \\ 3 & -2
\end{bmatrix}
\)\(
\begin{bmatrix}
w_{1} \\ w_{2}
\end{bmatrix}
\).
5.4. The final solution
This brings us to the final solution. We know that \(
\begin{bmatrix}
-1 & -7 \\ 3 & -2
\end{bmatrix}
\) maps \(
\begin{bmatrix}
w_{1} \\ w_{2}
\end{bmatrix}
\) to \(
\begin{bmatrix}
-1 \\ 1
\end{bmatrix}
\). Then the matrix that maps \(
\begin{bmatrix}
-1 \\ 1
\end{bmatrix}
\) back to \(
\begin{bmatrix}
w_{1} \\ w_{2}
\end{bmatrix}
\) is the inverse of \(
\begin{bmatrix}
-1 & -7 \\ 3 & -2
\end{bmatrix}
\).
We saw how to invert a 2-by-2 matrix in the last chapter. This is how a function to do so looks like:
def get_inverse_of(matrix2d):
x1,y1 = matrix2d[0]
x2,y2 = matrix2d[1]
scale = 1/(x1*y2 - x2*y1)
return [[y2*scale, -y1*scale],[-x2*scale, x1*scale]]Let's see what the inverse of this matrix \(M^{-1}\) is: Inverse of \(
\begin{bmatrix}
-1 & -7 \\ 3 & -2
\end{bmatrix}
\) = \(
\begin{bmatrix}
-0.09 & 0.3 \\ -0.13 & -0.04
\end{bmatrix}
\)
It is important to understand this part: Multiplying this inverse matrix \(M^{-1}\) with a vector \(
\begin{bmatrix}
x \\ y
\end{bmatrix}
\) gives us values (or weights) \(w_{1}, w_{2}\) for which the weighted sum of vectors AB and AC = \(
\begin{bmatrix}
x \\ y
\end{bmatrix}
\).
Now if these weights follow these three rules, it means that the point \(
\begin{bmatrix}
x \\ y
\end{bmatrix}
\) lies inside the triangle:
1. \(w_{1} >= 0\)
2. \(w_{2} >= 0\)
3. \(w_{1} + w_{2} <= 1\)
Let's write a function to check for these conditions:
def weights_are_valid(w1,w2):
return w1 >= 0 and w2 >= 0 and w1+w2 <= 1Note that in this example we have \(
\begin{bmatrix}
-0.09 & 0.3 \\ -0.13 & -0.04
\end{bmatrix}
\)\(
\begin{bmatrix}
-1 \\ 1
\end{bmatrix}
\) = \(
\begin{bmatrix}
0.39 \\ 0.09
\end{bmatrix}
\). Thus \(w_{1} = 0.39\) and \(w_{2} = 0.09\). Since these values satisfy the three conditions above, the point [-1,1] indeed lies inside the triangle.
Below is a full implementation of a function that takes as input vertices of a triangle and a point. It returns True if the point is inside the triangle and False otherwise.
def is_point_in_triangle(vertices, point):
'''
vertices: [[x1,y1],[x2,y2],[x3,y3]]
point: [p,q]
'''
x1,y1 = vertices[0]
x2,y2 = vertices[1]
x3,y3 = vertices[2]
shifted_vertices = [[x1-x1,y1-y1],[x2-x1,y2-y1],[x3-x1,y3-y1]]
shifted_point = [point[0]-x1, point[1]-y1]
# create a matrix with sides AB and AC as vectors
matrix = [shifted_vertices[1], shifted_vertices[2]]
inverse_matrix = get_inverse_of(matrix)
# get the weights needed to scale the vectors AB and AC
w1, w2 = matMul(inverse_matrix, shifted_point)
return weights_are_valid(w1,w2)5.4.1. Extending the idea to polygons
Once you figure out how to determine whether or not a point is inside a given triangle, you can easily extend this idea to polygons. All you need to do is to break down the polygon into a set of triangles. Then a given point is inside the polygon if it is inside any of these triangles. An example is shown below:
Problem: figure out whether or not a point [x,y] is inside a pentagon.
This problem can be broken down into whether or not the point [x,y] is inside a triangle. You need to solve this problem for three triangles.
This is what a function to do so would look like:
def is_point_in_polygon(vertices, point):
triangles = decompose_polygon_into_triangles(vertices)
for triangle in triangles:
if is_point_in_triangle(triangle, point):
return True
return FalseSome fun questions: how will you implement the function decompose_polygon_into_triangles? Does the function work for all polygons? What if a polygon is not convex?