How do I generalise Bresenham's line algorithm to floating-point endpoints?

I'm trying to combine two things. I'm writing a game and I need to determine the grid squares lying on a line with the floating-point endpoints.

Moreover I need it to include all the grid squares it touches (i.e. not just Bresenham's line but the blue one):

Can someone offer me any insight on how to do it? The obvious solution is to use naive line algorithm, but is there something more optimized (faster)?

• In case the link goes offline, just google for "A faster voxel traversal algorithm for raytracing" – Gustavo Maciel Aug 3 '14 at 16:49

You are looking for a grid traversal algorithm. This paper gives a good implementation;

Here's the basic implementation in 2D found on the paper:

loop {
if(tMaxX < tMaxY) {
tMaxX= tMaxX + tDeltaX;
X= X + stepX;
} else {
tMaxY= tMaxY + tDeltaY;
Y= Y + stepY;
}
NextVoxel(X,Y);
}


There's also a 3D ray-casting version on the paper.

In case the link rots, you can find many mirrors with its name: A faster voxel traversal algorithm for raytracing.

• Well, awkward. I guess, I'll switch answer to you and up-vote ltjax. Because I solved based on your link to that paper. – SmartK8 Aug 4 '14 at 21:09

Blue's idea is good, but the implementation is a bit clumsy. In fact, you can easily do it without sqrt. Let's assume for the moment that you exclude degenerate cases (BeginX==EndX || BeginY==EndY) and focus only on line directions in the first quadrant, so BeginX < EndX && BeginY < EndY. You'll have to implement a version for at least one other quadrant too, but that's very similar to the version for the first quadrant - you only check other edges. In C'ish pseudo code:

int cx = floor(BeginX); // Begin/current cell coords
int cy = floor(BeginY);
int ex = floor(EndX); // End cell coords
int ey = floor(EndY);

// Delta or direction
double dx = EndX-BeginX;
double dy = EndY-BeginY;

while (cx < ex && cy < ey)
{
// find intersection "time" in x dir
float t0 = (ceil(BeginX)-BeginX)/dx;
float t1 = (ceil(BeginY)-BeginY)/dy;

visit_cell(cx, cy);

if (t0 < t1) // cross x boundary first=?
{
++cx;
BeginX += t0*dx;
BeginY += t0*dy;
}
else
{
++cy;
BeginX += t1*dx;
BeginY += t1*dy;
}
}


Now for other quadrants, you just change the ++cx or ++cy and the loop condition. If you use this for collision, you probably have to implement all 4 versions, otherwise you can get away with two by appropriately swapping begin and end points.

• The algorithm Gustavo Maciel provided is a bit more efficient. It only determines first Ts and then just adds 1 to vertical or horizontal and shifts Ts by a size of cell. But as he didn't converted it to an answer I'll accept this one as the nearest answer. – SmartK8 Aug 4 '14 at 15:29

Your assumption isn't necessarily to find the cells but the lines it cross on this grid.

For example taking your image we can highlight not the cells, but lines of the grid it crosses:

This then shows that if it crosses a grid line that the cells either side of this line are those that are filled.

You can use an intersection algorithm to find if your floating point line will cross these by scaling your points to pixels. If you have a 1.0:1 ratio of floating coordinates:pixels then you're sorted and you can just translate it directly. Using the Line segment intersection algorithm you can check if your lower left line (1,7)(2,7) intersects with your line (1.3,6.2)(6.51,2.9). http://alienryderflex.com/intersect/

Some translation from c to C# will be needed but you can get the idea from that paper. I'll put the code below in-case the link breaks.

//  public domain function by Darel Rex Finley, 2006

//  Determines the intersection point of the line defined by points A and B with the
//  line defined by points C and D.
//
//  Returns YES if the intersection point was found, and stores that point in X,Y.
//  Returns NO if there is no determinable intersection point, in which case X,Y will
//  be unmodified.

bool lineIntersection(
double Ax, double Ay,
double Bx, double By,
double Cx, double Cy,
double Dx, double Dy,
double *X, double *Y) {

double  distAB, theCos, theSin, newX, ABpos ;

//  Fail if either line is undefined.
if (Ax==Bx && Ay==By || Cx==Dx && Cy==Dy) return NO;

//  (1) Translate the system so that point A is on the origin.
Bx-=Ax; By-=Ay;
Cx-=Ax; Cy-=Ay;
Dx-=Ax; Dy-=Ay;

//  Discover the length of segment A-B.
distAB=sqrt(Bx*Bx+By*By);

//  (2) Rotate the system so that point B is on the positive X axis.
theCos=Bx/distAB;
theSin=By/distAB;
newX=Cx*theCos+Cy*theSin;
Cy  =Cy*theCos-Cx*theSin; Cx=newX;
newX=Dx*theCos+Dy*theSin;
Dy  =Dy*theCos-Dx*theSin; Dx=newX;

//  Fail if the lines are parallel.
if (Cy==Dy) return NO;

//  (3) Discover the position of the intersection point along line A-B.
ABpos=Dx+(Cx-Dx)*Dy/(Dy-Cy);

//  (4) Apply the discovered position to line A-B in the original coordinate system.
*X=Ax+ABpos*theCos;
*Y=Ay+ABpos*theSin;

//  Success.
return YES; }


If you need to find out only when (and where) the line segments intersect, you can modify the function as follows:

//  public domain function by Darel Rex Finley, 2006

//  Determines the intersection point of the line segment defined by points A and B
//  with the line segment defined by points C and D.
//
//  Returns YES if the intersection point was found, and stores that point in X,Y.
//  Returns NO if there is no determinable intersection point, in which case X,Y will
//  be unmodified.

bool lineSegmentIntersection(
double Ax, double Ay,
double Bx, double By,
double Cx, double Cy,
double Dx, double Dy,
double *X, double *Y) {

double  distAB, theCos, theSin, newX, ABpos ;

//  Fail if either line segment is zero-length.
if (Ax==Bx && Ay==By || Cx==Dx && Cy==Dy) return NO;

//  Fail if the segments share an end-point.
if (Ax==Cx && Ay==Cy || Bx==Cx && By==Cy
||  Ax==Dx && Ay==Dy || Bx==Dx && By==Dy) {
return NO; }

//  (1) Translate the system so that point A is on the origin.
Bx-=Ax; By-=Ay;
Cx-=Ax; Cy-=Ay;
Dx-=Ax; Dy-=Ay;

//  Discover the length of segment A-B.
distAB=sqrt(Bx*Bx+By*By);

//  (2) Rotate the system so that point B is on the positive X axis.
theCos=Bx/distAB;
theSin=By/distAB;
newX=Cx*theCos+Cy*theSin;
Cy  =Cy*theCos-Cx*theSin; Cx=newX;
newX=Dx*theCos+Dy*theSin;
Dy  =Dy*theCos-Dx*theSin; Dx=newX;

//  Fail if segment C-D doesn't cross line A-B.
if (Cy<0. && Dy<0. || Cy>=0. && Dy>=0.) return NO;

//  (3) Discover the position of the intersection point along line A-B.
ABpos=Dx+(Cx-Dx)*Dy/(Dy-Cy);

//  Fail if segment C-D crosses line A-B outside of segment A-B.
if (ABpos<0. || ABpos>distAB) return NO;

//  (4) Apply the discovered position to line A-B in the original coordinate system.
*X=Ax+ABpos*theCos;
*Y=Ay+ABpos*theSin;

//  Success.
return YES; }

• Hi, the grid traversal is exactly for the purpose of optimizing thousands of line intersections all over the grid. This cannot be solved by thousands of line intersections. I have a map in a game with ground lines that player cannot cross. There can be thousands of these. I need to determine which one to calculate expensive intersection for. To determine these I only want calculate the intersections of those in line of player movement (or light from light source). In your case I would need to determine intersections with ~256x256x2 line segments each round. That would not be optimized at all. – SmartK8 Aug 3 '14 at 19:46
• But still thank you for you answer. Technically it works and is correct. But just not feasible for me. – SmartK8 Aug 3 '14 at 20:07
float difX = end.x - start.x;
float difY = end.y - start.y;
float dist = abs(difX) + abs(difY);

float dx = difX / dist;
float dy = difY / dist;

for (int i = 0, int x, int y; i <= ceil(dist); i++) {
x = floor(start.x + dx * i);
y = floor(start.y + dy * i);
draw(x,y);
}
return true;


JS Demo :

//C# stackoverflow

function verifyLoS(start, end) {
var difX = end.x - start.x;
var difY = end.y - start.y;
var dist = Math.abs(difX) + Math.abs(difY);

var dx = difX / dist;
var dy = difY / dist;

for (var i = 0, x, y; i <= Math.ceil(dist); i++) {
x = Math.floor(start.x + dx * i);
y = Math.floor(start.y + dy * i);
draw(x,y);
}
return true;
}

var HEIGHT = 14,
WIDTH = 14,
PX = 32;

var canvas, ctx;

function main () {
canvas = document.getElementById("canvas");
ctx = canvas.getContext("2d");

canvas.height = HEIGHT * PX;
canvas.width = WIDTH * PX;

loop();
}

function loop() {
for (var x=0;x<HEIGHT; x++) {
for (var y=0; y<WIDTH; y++) {
ctx.fillStyle = ((x+y)%2===0)?"#AFA":"#AEA";
ctx.fillRect(x*PX,y*PX,PX,PX);
}
}
var a = {x:1.5, y:9.5},
b = {x:12.5, y:1.5};
verifyLoS(a, b);
ctx.beginPath();
ctx.moveTo(a.x*PX,a.y*PX);
ctx.lineTo(b.x*PX,b.y*PX);
ctx.stroke();
//window.requestAnimationFrame(loop);
}

function draw(x, y) {
ctx.fillStyle = "rgba(255,0,0,0.2)";
ctx.fillRect(x*PX,y*PX,PX,PX);
}

main();
<!DOCTYPE html>
<html>
<meta charset="utf-8">
<title>JS Bin</title>
<body>
<canvas id="canvas"></canvas>
</body>
</html>

• This failed for me due to floating point numerical errors (the loop will do an extra iteration for the most tiny fraction over the next integer which will push the line end-point beyond the 'end' location). The simple fix is to calculate dist as a ceil in the first place so dx,dy are divided by the integer number of iterations of the loop (this means you can lose the ceil(dist) in the for loop). – PeteB Apr 6 '16 at 10:54

I ran into the same problem today and made a pretty big mountain of spaghetti out of a mole hill but ended up with something that works: https://github.com/SnpM/Pan-Line-Algorithm.