The Wu line drawing algorithm for anti-aliased lines optimization

Question

I need to optimize the following function with low-level optimizations (No SIMD, Multithreading). I already applied a lot of optimizations and got it ~65% faster, but it is possible to get it 200% faster. I am practicing on code optimization, some explanation would be welcome.

The function describes the anti-aliased Wu line drawing algorithm. It draws colored lines. It takes into account the sensitivity of the human eye which can't be changed. Approximations are completely fine.

void DrawWuLine( Surface *screen, int X0, int Y0, int X1, int Y1, Color& clrLine )
{
    /* Make sure the line runs top to bottom */
    constexpr float pixelNormalization = 1.f / 255.f;

    uint primaryColor, secondaryColor;

    if (Y0 > Y1)
    {
        std::swap(Y0, Y1);
        std::swap(X0, X1);
    }

    int XDir = 1;
    int DeltaX = X1 - X0;
    int DeltaY = Y1 - Y0;

    if (DeltaX < 0)
    {
        XDir = -1;
        DeltaX = 0 - DeltaX;
    }
    
    int ErrorAdj, ErrorAcc = 0;
    int Weighting;

    int grayl = ((clrLine.r * 19595) + (clrLine.g * 38469) + (clrLine.b * 7471)) >> 16;

    screen->Plot(X0, Y0, clrLine.RGBA);
    
    /* Is this an X-major or Y-major line? */
    if (DeltaY > DeltaX)
    {
        ErrorAdj = (DeltaX << 16) / DeltaY;

        while (--DeltaY) {

            if (((ErrorAcc + ErrorAdj) & 0xFFFF) <= ErrorAcc) {
                X0 += XDir;
            }
            Y0++; 

            ErrorAcc = (ErrorAcc + ErrorAdj) & 0xFFFF;
            Weighting = (ErrorAcc >> 8) & 0xFF;
            
            uint clrBackGround = screen->pixels[X0 + Y0 * SCRWIDTH];

            int rb = clrBackGround & 0xFF;
            int gb = (clrBackGround >> 8) & 0xFF;
            int bb = (clrBackGround >> 16) & 0xFF;
            
            int grayb = ((rb * 19595) + (gb * 38469) + (bb * 7471)) >> 16;
            
            int weightXOR = (~Weighting) & 0xFF;

            bool isLighter = grayl < grayb;

            float horizontalWeighting = (isLighter ? Weighting : weightXOR) * pixelNormalization;

            uint rr = horizontalWeighting * abs(rb - clrLine.r) + std::min(rb, (int) clrLine.r);
            uint gr = horizontalWeighting * abs(gb - clrLine.g) + std::min(gb, (int) clrLine.g);
            uint br = horizontalWeighting * abs(bb - clrLine.b) + std::min(bb, (int) clrLine.b);
            
            primaryColor = (br << 16) | (gr << 8) | rr;

            clrBackGround = screen->pixels[X0 + XDir + Y0 * SCRWIDTH];
            
            rb = clrBackGround & 0xFF;
            gb = (clrBackGround >> 8) & 0xFF;
            bb = (clrBackGround >> 16) & 0xFF;
            
            grayb = ((rb * 19595) + (gb * 38469) + (bb * 7471)) >> 16;

            isLighter = grayl < grayb;

            float verticalWeighting = (isLighter ? weightXOR : Weighting) * pixelNormalization;
            
            rr = verticalWeighting * abs(rb - clrLine.r) + std::min(rb, (int)clrLine.r);
            gr = verticalWeighting * abs(gb - clrLine.g) + std::min(rb, (int)clrLine.g);
            br = verticalWeighting * abs(bb - clrLine.b) + std::min(rb, (int)clrLine.b);

            secondaryColor = (br << 16) | (gr << 8) | rr;
            
            screen->Plot(X0, Y0, primaryColor);
            screen->Plot( X0 + XDir, Y0, secondaryColor);
        }
      
        screen->Plot( X1, Y1, clrLine.RGBA );

        return;
    }
    
    ErrorAdj = (DeltaY << 16) / DeltaX;

    while (--DeltaX) {
       
        if (((ErrorAcc + ErrorAdj) & 0xFFFF) <= ErrorAcc) {
            Y0++;
        }
        
        X0 += XDir; 

        ErrorAcc = (ErrorAcc + ErrorAdj) & 0xFFFF;
        Weighting = (ErrorAcc >> 8) & 0xFF;
        
        uint clrBackGround = screen->pixels[X0 + Y0 * SCRWIDTH];
        int rb = clrBackGround & 0xFF;
        int gb = (clrBackGround >> 8) & 0xFF;
        int bb = (clrBackGround >> 16) & 0xFF;
        
        int grayb = ((rb * 19595) + (gb * 38469) + (bb * 7471)) >> 16;

        int weightXOR = (~Weighting) & 0xFF;

        bool isLighter = grayl < grayb;

        float horizontalWeighting = (isLighter ? Weighting : weightXOR) * pixelNormalization;

        uint rr = horizontalWeighting * abs(rb - clrLine.r) + std::min(rb, (int)clrLine.r);
        uint gr = horizontalWeighting * abs(gb - clrLine.g) + std::min(gb, (int)clrLine.g);
        uint br = horizontalWeighting * abs(bb - clrLine.b) + std::min(bb, (int)clrLine.b);

        primaryColor = (br << 16) | (gr << 8) | rr;
        
        clrBackGround = screen->pixels[X0 + (Y0 + 1) * SCRWIDTH];
        
        rb = clrBackGround & 0xFF;
        gb = (clrBackGround >> 8) & 0xFF;
        bb = (clrBackGround >> 16) & 0xFF;

        grayb = ((rb * 19595) + (gb * 38469) + (bb * 7471)) >> 16;

        isLighter = grayl < grayb;

        float verticalWeighting = (isLighter ? weightXOR : Weighting) * pixelNormalization;

        rr = verticalWeighting * abs(rb - clrLine.r) + std::min(rb, (int)clrLine.r);
        gr = verticalWeighting * abs(gb - clrLine.g) + std::min(rb, (int)clrLine.g);
        br = verticalWeighting * abs(bb - clrLine.b) + std::min(rb, (int)clrLine.b);

        secondaryColor = (br << 16) | (gr << 8) | rr;

        screen->Plot(X0, Y0, primaryColor);
        screen->Plot( X0, Y0 + 1, secondaryColor );
    }
    
    screen->Plot( X1, Y1, clrLine.RGBA );
}

Any tips or insights on how to get rid of branches inside of the while loops? Other optimizations are welcome.

Answer 1

This submission is about performance, yet it includes zero profiling, nor even throughput measurements (lines per second, for some uniform line length).

describe prior work

It would have been helpful Review Context to tell us what improvements are incorporated in the OP, and what initially promising changes didn't pan out.

For example, did manual loop unroll at the source code level have any useful effect? You would need to tell us about deployment environment and compiler flags when describing such attempts.

cite your reference

If the human physiology {19595, 38469, 7471} constants are a scaled version of one of these matrices, then tell us which one. As stated, it's unclear what spec. the OP code tries to conform to, and therefore it's not possible to say whether it is "correct". It can only be correct w.r.t. some spec.

algorithm

This is a color algorithm. It would be of interest to see its B&W analogue.

One such is presented on the corresponding wikipedia page. I am struck by the extreme simplicity of its 3-line inner loop, compared to the 31 SLOC in the OP loop. I am unwilling to believe that adding a "must plot in full color!" requirement would necessitate a 10x blow-up in number of statements. You may find that benching Wu B&W antialiased line drawing proves to be an educational experience.

In any event, do cite a reference when you tell us you're implementing Wu's algorithm, as it is written up in multiple places.

DRY

    int grayl = ((clrLine.r * 19595) + (clrLine.g * 38469) + (clrLine.b * 7471)) >> 16;
            ...           
            int grayb = ((rb * 19595) + (gb * 38469) + (bb * 7471)) >> 16;
            ...
            grayb = ((rb * 19595) + (gb * 38469) + (bb * 7471)) >> 16;

Looks like we hit the Rule of Three, there. Please clean this up with an inlined function or a macro.

branch misprediction

        ... = (isLighter ? weightXOR : Weighting) * ... ;

That looks like a bunch of pipeline stalls to me, since it is hard-to-predict. Consider using an unconditional expression. From arithmetic we have:
(isLighter * weightXOR + (1 - isLighter) * Weighting).
Or from bitwise operators we could assemble a similar expression.

more comments, please

The "8 bits per color" RGBA unpacks were perfectly clear. But it's not obvious why 16 bits matters in this expression:

            ErrorAcc = (ErrorAcc + ErrorAdj) & 0xFFFF;

A brief explanation, or citation of an author, would help future maintenance engineers.

framebuffer

All graphic updates are accomplished via screen->Plot(x, y, c). For "steep" lines we keep writing adjacent pixel pairs like this:

            screen->Plot(X0,        Y0, ... );
            screen->Plot(X0 + XDir, Y0, ... );

Suppose you had access to the underlying pixel storage. Would issuing a single store which updates both pixels be winning?