Arbitrary angle rotation through shear (Paeth algorithm)

Question

I'm trying to write an Java implementation of the 3 shear rotation algorithm described by Alan Paeth. The problem is not the the calculation of the values but to fit the rotated points onto the image grid. In the paper, the rotation is performed by 3 consecutive shears given by the following calculation:

1. x = x + alpha * y
2. y= y + beta * x
3. x = x + alpha * y

Alpha and beta are calculated from the given angle (theta;in radians) by the following formulas:

• beta = sin(theta)
• alpha = - tan(theta/2)

Using these formulas, the points are rotated around the center of the coordinate system.

To correct the negative values, i add the minimal calculated coordinate for the respective axis to each point so that the minimal value will always be 0.

My Java implementation so far:

ShiftPoint[] val = new ShiftPoint[m*n];
double minX = 0,minY = 0, maxX = 0, maxY = 0;
double alpha = -1d*  Math.tan(Math.toRadians(theta)/2d);
double beta = Math.sin(Math.toRadians(theta));
for(int a = 0; a < m; a++) {
    for(int b = 0; b < n; b++) {
        ShiftPoint temp = new ShiftPoint(a, b, values[a][b]);
        double newX = b + alpha * a;    //first shear
        double newY = a + beta * newX;  //second shear
        newX += alpha * newY;           //third shear
        temp.setX(newX);
        temp.setY(newY);
        val[m * b + b] = temp;
    }
}

Note: ShiftPoint is a simple self-written class to hold the specific coordinates and the value inside the matrix (in case of image processing: the rgb value of the pixel). Here's a graphic representation of the calculations:

The problem: Whilst the calculated values seem to be correct and the graphic representation shows that the rotation actually works, i am not sure how to fit the calculated values on the fixed grid an image (or a 2d array) has without distorting it. Also i don't fully understand the implementation (for the x-axis shear) given in Paeths paper:

I get that skewi is the integer part of the calculated value and skewf is the fractional part, but what are width, height, oleft and left supposed to be? Also: Why does he add 0.5 to the y value and does not take the x value into account in his first calculation?

Note: I'm aware of the fact that Java offers simple ways to rotate images, but i'm trying to implement this specific algorithm just for the fun of it. I'm also aware of the 3 - 5 websites which can be found via websearch (such as #1 and #2) and try to explain that algorithm, but firstly they do not use java and secondly they mostly reference the example implementation by Paeth, so they are not terribly useful.

Answer 1

The essential principle behind this approach to image rotation is twofold:

• Define a coordinate transformation which rotates the x & y axes in the rotated image onto the x & y axes of the original image
• Compute each pixel in the rotated image by interpolating between values near the corresponding point in the original image

The first step is generally the easier one, and will involve simple linear combinations of the x & y coordinates. Shear transformations (which aren't strictly rotations, because they don't conserve the area of each pixel) just involve something like x -> x + alpha * y.

The algorithm given in Paeth's paper (which dates from 1986) seems to be a carefully optimized approach to doing the second (interpolation) step for a shear transformation. I think this boils-down to a piecewise-linear interpolation along the x-axis, but written in a form that doesn't require more than one array lookup for each pixel in the output image. A clearer (and slightly less efficient) approach might involve something like skewf * pixel(x-skewi-1, y) + (1-skewf) * pixel(x-skewi, y).

This particular algorithm is obviously highly specialized for a single-axis skew. For a general rotation, you might need something more like a bilinear interpolation over each 2x2 square of pixels that surround the unrotated location corresponding to the centre of each pixel in the rotated image. (Computing this centre-pixel value may be the origin of the y+0.5 in Paeth's code.)

Given how much compute power we now have relative 1986, I suspect your code would be much easier to understand if you included an explicit formula for this bilinear interpolation, even if it does use more array lookups than Paeth's approach, unless you really need maximum performance.