Getting a list of locations within a triangle in the form of x,y positions

Question

Let's say I have a triangle given by the three integer vertices (x1,y1), (x2,y2) and (x3,y3). What sort of algorithm can I use to return a comprehensive list of ALL (x,y) integer pairs that lie inside the triangle.

Answer 1

If you're dealing with rasterizing a fill, then we can skip more involved algorithms like Bresenham's line algorithm as we're not interested in getting "the best rasterization of a 1px thick line", we want all the coordinates located inside our triangular bounds, with some control of what to do with coordinates that are "on" the boundary itself (Do we exclude all of those? Or include all of them? Do we only include them if they're at least 50% "inside" the triangle?).

This simplifies things a great deal.

We can construct a sequence of x-intervals (rather than the y direction) so we can read consecutive stretches of values in a row in a 2D array, or simple stretches in a 1D flattened array, with the only challenge being "finding the right start and end points" for each stretch of coordinates to include:

1. Sort your three coordinates on y coordinate. Call the top point P1, the bottom point P3, and the middle point P2.
2. Compute the expressions for the three lines P1--P3, P1--P2, and P2--P3
   1. P1--P3 as a function of y is x = P1.x + (P3.x-P1.x)/(P3.y-P1.y) * (y - P1.y)
   2. P1--P2 as a function of y is x = P1.x + (P2.x-P1.x)/(P2.y-P1.y) * (y - P1.y)
   3. P2--P3 as a function of y is x = P2.x + (P3.x-P2.x)/(P3.y-P2.y) * (y - P2.y)
3. Add P1 to your set.
4. For i=1 to i=(P3.y - P1.y) do:
   1. set y = P1.y + i
   2. compute x1 on P1--P3 for the current y
   3. compute x2:
      1. if i < (P2.y - P1.y): compute x2 on P1--P2 for the current y
      2. else: compute x2 on P2--P3 for the current y instead
   4. add all coordinates from min(x1,x2) to max(x1,x2) to your set.

Of course, we need to turn those min and max values into integers rather than floats, and depending on our inclusion/exclusion policy for pixels on the bounary, we have different options:

• for points "inside the triangle to any degree", use floor(min(x1,x2)) and ceil(max(x1,x2)),
• for points "fully inside the triangle" use ceil(min(x1,x2)) and floor(max(x1,x2)), and
• for points "at least half inside the triangle", use round(min(x1,x2)) and round(max(x1,x2)).

const { random, PI, min, max, floor, ceil } = Math;
const cvs = document.querySelector(`canvas`);
const w = (cvs.width = 400);
const h = (cvs.height = 180);
const ctx = cvs.getContext(`2d`);

const pts = [
  { x: random() * w, y: random() * h },
  { x: random() * w, y: random() * h },
  { x: random() * w, y: random() * h },
];
const [P1, P2, P3] = pts.sort((a, b) => a.y - b.y);

const f12 = (y) => P1.x + ((P2.x - P1.x) / (P2.y - P1.y)) * (y - P1.y);
const f13 = (y) => P1.x + ((P3.x - P1.x) / (P3.y - P1.y)) * (y - P1.y);
const f23 = (y) => P2.x + ((P3.x - P2.x) / (P3.y - P2.y)) * (y - P2.y);

// Let's gather up our points!
const allPoints = [P1];

// We're skipping by 5, instead of 1, because otherwise we'll
// just end up coloring every pixel, which is not as useful
// as illustrative code example.
for (let i = 1; i < P3.y - P1.y; i += 4) {
  const y = P1.y + i;
  const x1 = f13(y);
  const f = i < P2.y - P1.y ? f12 : f23;
  const x2 = f(y);
  const x = floor(min(x1, x2));
  const X = ceil(max(x1, x2));
  // Same here, we're skipping by 5, rather than 1, purely
  // for illustrative purposes.
  for (let j = 0; j < X - x; j += 5) {
    allPoints.push({ x: x + j, y });
  }
}

// And then let's draw them out to see how we did.
ctx.strokeStyle = `#C00`;
ctx.beginPath();
ctx.moveTo(P1.x, P1.y);
allPoints.forEach((p) => {
  const { x, y } = p;
  ctx.moveTo(x, y);
  ctx.arc(x, y, 0.1, 0, 2 * PI);
});
ctx.stroke();

ctx.strokeStyle = `#F3F6`;
ctx.beginPath();
ctx.moveTo(P1.x, P1.y);
ctx.lineTo(P2.x, P2.y);
ctx.lineTo(P3.x, P3.y);
ctx.lineTo(P1.x, P1.y);
ctx.stroke();

canvas { 
  border: 1px solid black;
}

<canvas></canvas>

(note that the HTML canvas is a bit weird for illustrative purposes here, but it's the best we get: it does not have "pixels" but grid coordinates, so drawing something on (1,1) does not draw it on the first pixel, it draws it between the top-left four pixels. It is what it is)

And of course: remember to take the edge cases into account: if P1 and P2 are on the same row, or if P2 and P3 are on the same row, then the algorithm simplifies to one where we're not splitting "which function to use" based on which value i is: we simply always compute x2 using whichever of the two lines is not the flat one =)

Answer 2

The following algorithm should be appropriate:

1. Sort the triangle vertices by x coordinate in increasing order. Now we have two segments (1-2 and 2-3) on the one side (top or bottom), and one segment from the other one (1-3).

2. Compute coefficients of equations of lines (which contain the segments):

   A * x + B * y + C = 0
   A = y2 - y1
   B = x1 - x2
   C = x2 * y1 - x1 * y2
   

   There (x1, y1) and (x2, y2) are two points of the line.

3. For each of ranges [x1, x2), (x2, x3], and x2 (special case) iterate over integer points in ranges and do the following for every x:

   1. Find top bound as y_top = (- A1 * x - C1) div B1.
   2. Find bottom bound as y_bottom = (- A2 * y - C2 - 1) div B2 + 1.
   3. Add all points between (x, y_bottom) and (x, y_top) to the result.

This algorithm is for not strictly internal vertices. For strictly internal vertices items 3.1 and 3.2 slightly differ.

Answer 3

The proper name for this problem is triangle rasterization.

It's a well researched problem and there's variety of methods to do it. The two popular methods are:

1. Scan line by scan line.

   For each scan-line you require some basic geometry to recalculate the start and the end of the line. See Bresenham's Line drawing algorithm.

2. Test every pixel in the bounding box to see if it is in the triangle.

   This is usually done by using barycentric co-ordinates.

Most people assume method 1) is more efficient as you don't waste time testing pixels that can are outside the triangle, approximately half of all the pixels in the bounding box. However, 2) has a major advantage - it can be run in parallel far more easily and so for hardware is usually the much faster option. 2) is also simpler to code.

The original paper for describing exactly how to use method 2) is written by Juan Pineda in 1988 and is called "A Parallel Algorithm for Polygon Rasterization".

For triangles, it's conceptually very simple (if you learn barycentric co-ordindates). If you convert each pixel into triangle barycentric coordinates, alpha, beta and gamma - then the simple test is that alpha, beta and gamma must be between 0 and 1.

Answer 4

I like ray casting, nicely described in this Wikipedia article. Used it in my project for the same purpose. That method scales on other polygons too, including concave. Not sure about the performance, but it is easily coded, so you could try it yourself (I had no performance issues in my project)

Answer 5

I suppose you have a list of pairs you want to test (if this is not what your problem is about, please specify your question clearly). You should store the pairs into quad-tree or kd-tree structure first, in order to have a set of candidates which is small enough. If you have few points, this is probably not worth the hassle (but it won't scale well if you don't do it).

You can also narrow down candidates further by testing against a bounding box for your triangle.

Then, for each candidate pair (x, y), solve in a, b, c the system

a + b + c = 1
a x1 + b x2 + c x3 = x
a y2 + b y2 + c y3 = y

(I let you work this out), and the point is inside the triangle if a b and c are all positive.