Editing Bresenham's line algorithm (section)

==Derivation==
To derive Bresenham's algorithm, two steps must be taken.  The first step is transforming the equation of a line from the typical slope-intercept form into something different; and then using this new equation to draw a line based on the idea of accumulation of error.

===Line equation===
[[File:Line 1.5x+1.svg|300px|thumb|y=f(x)=.5x+1 or f(x,y)=x-2y+2=0]]
[[File:Line 1.5x+1 -- planes.svg|300px|thumb|Positive and negative half-planes]]

The slope-intercept form of a line is written as

:<math>y = f(x) = mx + b</math>

where <math>m</math> is the slope and <math>b</math> is the y-intercept. Because this is a function of only <math>x</math>, it can't represent a vertical line. Therefore, it would be useful to make this equation written as a function of both <math>x</math> ''and'' <math>y</math>, to be able to draw lines at any angle. The angle (or slope) of a line can be stated as "rise over run", or <math>\Delta y/\Delta x</math>. Then, using algebraic manipulation,

:<math>
\begin{align}
y & = mx + b \\
y & = \frac{\Delta y}{\Delta x} x + b \\
(\Delta x) y & = (\Delta y) x + (\Delta x) b \\
0 & = (\Delta y) x - (\Delta x) y + (\Delta x) b
\end{align}
</math>

Letting this last equation be a function of <math>x</math> and <math>y</math>, it can be written as

:<math>f(x,y) := Ax + By + C = 0</math>

where the constants are
* <math>A = \Delta y = y_1 - y_0</math>
* <math>B = - \Delta x = - (x_1 - x_0)</math>
* <math>C = (\Delta x) b = (x_1 - x_0) b</math>

The line is then defined for some constants <math>A</math>, <math>B</math>, and <math>C</math> anywhere <math>f(x,y) = 0</math>.  That is, for any <math>(x,y)</math> not on the line, <math>f(x,y) \ne 0</math>. This form involves only integers if <math>x</math> and <math>y</math> are integers, since the constants <math>A</math>, <math>B</math>, and <math>C</math> are defined as integers.

As an example, the line <math display="inline">y=\frac{1}{2}x + 1</math> then this could be written as <math>f(x,y) = x - 2y + 2</math>.  The point (2,2) is on the line

:<math>f(2,2) = x - 2y + 2 = (2) - 2(2) + 2 = 2 - 4 + 2 = 0</math>

and the point (2,3) is not on the line

:<math>f(2,3) = (2) - 2(3) + 2 = 2 - 6 + 2 = -2</math>

and neither is the point (2,1)

:<math>f(2,1) = (2) - 2(1) + 2 = 2 - 2 + 2 = 2</math>

Notice that the points (2,1) and (2,3) are on opposite sides of the line and <math>f(x,y)</math> evaluates to positive or negative.  A line splits a plane into halves and the half-plane that has a negative <math>f(x,y)</math> can be called the negative half-plane, and the other half can be called the positive half-plane.  This observation is very important in the remainder of the derivation.

===Algorithm===
The starting point is on the line

:<math>f(x_0, y_0) = 0</math>

only because the line is defined to start and end on integer coordinates (though it is entirely reasonable to want to draw a line with non-integer end points).

[[File:Line 1.5x+1 -- candidates.svg|300px|thumb|Candidate point (2,2) in blue and two candidate points in green (3,2) and (3,3)]]

Keeping in mind that the slope is at most <math>1</math>, the problem now presents itself as to whether the next point should be at <math>(x_0 + 1, y_0)</math> or <math>(x_0 + 1, y_0 + 1)</math>. Perhaps intuitively, the point should be chosen based upon which is closer to the line at <math>x_0 + 1</math>.  If it is closer to the former then include the former point on the line, if the latter then the latter. To answer this, evaluate the line function at the midpoint between these two points:

:<math>f(x_0 + 1, y_0 + \tfrac 1 2)</math>

If the value of this is positive then the ideal line is below the midpoint and closer to the candidate point <math>(x_0+1,y_0+1)</math>; i.e. the y coordinate should increase.  Otherwise, the ideal line passes through or above the midpoint, and the y coordinate should stay the same; in which case the point <math>(x_0+1,y_0)</math> is chosen. The value of the line function at this midpoint is the sole determinant of which point should be chosen.

The adjacent image shows the blue point (2,2) chosen to be on the line with two candidate points in green (3,2) and (3,3).  The black point (3, 2.5) is the midpoint between the two candidate points.

=== Algorithm for integer arithmetic ===
Alternatively, the difference between points can be used instead of evaluating f(x,y) at midpoints. This alternative method allows for integer-only arithmetic, which is generally faster than using floating-point arithmetic. To derive the other method, define the difference to be as follows: 
:<math>
D_i = f(x_i+1,y_i+\tfrac 1 2) - f(x_0,y_0)
</math>
For the first decision, this formulation is equivalent to the midpoint method since <math>f(x_0,y_0)=0</math> at the starting point. Simplifying this expression yields:
:<math>\begin{array}{rclcl}
D_0 & = & \left[ A(x_0+1) + B \left(y_0+\frac{1}{2}\right) + C \right] & - & \left[ A x_0 + B y_0 + C \right] \\
& = & \left[ Ax_0 + B y_0+ C + A + \frac {1}{2} B\right]  & - & \left[ A x_0 + B y_0 + C \right] \\
& = & A + \frac{1}{2} B = \Delta y - \frac{1}{2} \Delta x
\end{array}</math>

Just as with the midpoint method, if <math>D_0</math> is positive, then choose <math>(x_0+1,y_0+1)</math>, otherwise choose <math>(x_0+1,y_0)</math>.

If <math>(x_0+1,y_0)</math> is chosen, the change in <math>D_i</math> will be: 
:<math>\begin{array}{lclcl}
\Delta D &=& f(x_0+2,y_0+\tfrac 1 2) - f(x_0+1,y_0+\tfrac 1 2) &=& A &=& \Delta y \\
\end{array}</math>
If <math>(x_0+1,y_0+1)</math> is chosen the change in <math>D_i</math> will be:
:<math>\begin{array}{lclcl}
\Delta D &=& f(x_0+2,y_0+\tfrac 3 2) - f(x_0+1,y_0+\tfrac 1 2) &=& A+B &=& \Delta y - \Delta x
\end{array}</math>

If the new D is positive then <math>(x_0+2,y_0+1)</math> is chosen, otherwise <math>(x_0+2,y_0)</math>. This decision can be generalized by accumulating the error on each subsequent point.

[[File:Line 1.5x+1 -- points.svg|300px|thumb|Plotting the line from (0,1) to (6,4) showing a plot of grid lines and pixels]]

All of the derivation for the algorithm is done.  One performance issue is the 1/2 factor in the initial value of D.  Since all of this is about the sign of the accumulated difference, then everything can be multiplied by 2 with no consequence.

This results in an algorithm that uses only integer arithmetic.

 plotLine(x0, y0, x1, y1)
     dx = x1 - x0
     dy = y1 - y0
     D = 2*dy - dx
     y = y0
 
     '''for''' x from x0 to x1
         plot(x, y)
         '''if''' D > 0
             y = y + 1
             D = D - 2*dx
         '''end if'''
         D = D + 2*dy

Running this algorithm for <math>f(x,y) = x-2y+2</math> from (0,1) to (6,4) yields the following differences with dx=6 and dy=3:

 D=2*3-6=0
 Loop from 0 to 6
  * x=0: '''plot(0, 1)''', D≤0: D=0+6=6
  * x=1: '''plot(1, 1)''', D>0: D=6-12=-6, y=1+1=2, D=-6+6=0
  * x=2: '''plot(2, 2)''', D≤0: D=0+6=6
  * x=3: '''plot(3, 2)''', D>0: D=6-12=-6, y=2+1=3, D=-6+6=0
  * x=4: '''plot(4, 3)''', D≤0: D=0+6=6
  * x=5: '''plot(5, 3)''', D>0: D=6-12=-6, y=3+1=4, D=-6+6=0
  * x=6: '''plot(6, 4)''', D≤0: D=0+6=6

The result of this plot is shown to the right.  The plotting can be viewed by plotting at the intersection of lines (blue circles) or filling in pixel boxes (yellow squares).  Regardless, the plotting is the same.

===All cases===
However, as mentioned above this only works for [[octant (plane geometry)|octant]] zero, that is lines starting at the origin with a slope between 0 and 1 where x increases by exactly 1 per iteration and y increases by 0 or 1.

The algorithm can be extended to cover slopes between 0 and -1 by checking whether y needs to increase or decrease (i.e. dy < 0)

 plotLineLow(x0, y0, x1, y1)
     dx = x1 - x0
     dy = y1 - y0
     yi = 1
     '''if''' dy < 0
         yi = -1
         dy = -dy
     '''end if'''
     D = (2 * dy) - dx
     y = y0
 
     '''for''' x from x0 to x1
         plot(x, y)
         '''if''' D > 0
             y = y + yi
             D = D + (2 * (dy - dx))
         '''else'''
             D = D + 2*dy
         '''end if'''

By switching the x and y axis an implementation for positive or negative steep slopes can be written as 

 plotLineHigh(x0, y0, x1, y1)
     dx = x1 - x0
     dy = y1 - y0
     xi = 1
     '''if''' dx < 0
         xi = -1
         dx = -dx
     '''end if'''
     D = (2 * dx) - dy
     x = x0
 
     '''for''' y from y0 to y1
         plot(x, y)
         '''if''' D > 0
             x = x + xi
             D = D + (2 * (dx - dy))
         '''else'''
             D = D + 2*dx
         '''end if'''

A complete solution would need to detect whether x1 > x0 or y1 > y0 and reverse the input coordinates before drawing, thus
 
 plotLine(x0, y0, x1, y1)
     '''if''' abs(y1 - y0) < abs(x1 - x0)
         '''if''' x0 > x1
             plotLineLow(x1, y1, x0, y0)
         '''else'''
             plotLineLow(x0, y0, x1, y1)
         '''end if'''
     '''else'''
         '''if''' y0 > y1
             plotLineHigh(x1, y1, x0, y0)
         '''else'''
             plotLineHigh(x0, y0, x1, y1)
         '''end if'''
     '''end if'''

In low level implementations which access the video memory directly, it would be typical for the special cases of vertical and horizontal lines to be handled separately as they can be highly optimized.

Some versions use Bresenham's principles of integer incremental error to perform all octant line draws, balancing the positive and negative error between the x and y coordinates.<ref name=Zingl/>

 plotLine(x0, y0, x1, y1)
     dx = abs(x1 - x0)
     sx = x0 < x1 ? 1 : -1
     dy = -abs(y1 - y0)
     sy = y0 < y1 ? 1 : -1
     error = dx + dy
     
     '''while''' true
         plot(x0, y0)
         e2 = 2 * error
         '''if''' e2 >= dy
             '''if''' x0 == x1 '''break'''
             error = error + dy
             x0 = x0 + sx
         '''end if'''
         '''if''' e2 <= dx
             '''if''' y0 == y1 '''break'''
             error = error + dx
             y0 = y0 + sy
         '''end if'''
     '''end while'''