Python function that shears image by n degrees while keeping the relationship between height and width

Question

I am trying to find a method to shear/skew an image horizontally or vertically by n degrees (n ranges from -90 to 90 excluding the terminals) so that the result would meet the following conditions:

1. The result would keep the size of the dimension that is not sheared, it means if you shear the image horizontally, the result would have the same height as the input and vice-versa.

2. The result would keep the relationship between the width and height of the original image.

The above means, if you name the width and height of the original image w and h, and shear the image horizontally by n degrees, the result would have be h tall, and the width w1 would be:

w1 = w / math.cos(math.radians(n)) + h * math.tan(math.radians(n))

So that the skewed width w2 and skewed height h2 would satisfy the following:

h2 = h / math.cos(math.radians(n))
w2 = w1 - h * math.tan(math.radians(n))
assert h2/w2 == h/w

I have written a custom function that does exactly this, except that it is inefficient and inelegant

---

Code

import math
import numpy as np
from itertools import product
from PIL import Image

def rhombize(arr, angle, dimension='vertical', reverse=False):
    assert 0 < angle < 90
    assert dimension in ('horizontal', 'vertical')
    height, width, dimensions = arr.shape
    ratio1 = math.cos(math.radians(angle))
    ratio2 = math.tan(math.radians(angle))
    if reverse:
        arr = np.flipud(arr) if dimension == 'vertical' else np.fliplr(arr)
    h, w = height, width
    
    if dimension == 'vertical':
        h = round(height / ratio1 + width * ratio2 )
    else:
        w = round(width / ratio1 + height * ratio2 )
    
    rhomb = np.zeros((h, w, dimensions), dtype=np.uint8)
    for y in range(height):
        for x in range(width):
            a, b = x, y
            a1, b1 = a+1, b+1
            if dimension == 'vertical':
                b = round(y / ratio1 + x * ratio2)
                b1 = round((y + 1) / ratio1 + (x + 1) * ratio2)
            else:
                a = round(x / ratio1 + y * ratio2)
                a1 = round((x + 1) / ratio1 + (y + 1) * ratio2)
            rhomb[b:b1, a:a1] = arr[y, x]
    if reverse:
        rhomb = np.flipud(rhomb) if dimension == 'vertical' else np.fliplr(rhomb)
    return rhomb

if __name__ == '__main__':
    arr = np.zeros((1024, 768, 3), dtype=np.uint8)
    for y, x in product(range(256), repeat=2):
        arr[y*4:y*4+4, x*3:x*3+3] = (x, 0, y)
    img1 = Image.fromarray(arr)
    img1.show()
    img1.save('D:/img1.jpg')
    rhombized = rhombize(arr, 30)
    img2 = Image.fromarray(rhombized)
    img2.show()
    img2.save('D:/img2.jpg')

First image:

Second image:

I am wondering what is the proper way to do this, I know about affine transformations and that seems indeed to be the right way to do this, but my math is not that good and none of the online tutorials I can find covers it.

So how can I properly implement the exactly same logic using affine transformations?

---

Edit

Maybe my description does not succinctly describe what I am really trying to achieve here, so here are some pictures that explains what I am really trying to do better:

Then I am trying to upscale that rectangle, so that the dimension that is not skewed equals to its original value.

---

My code did exactly that, shearing while keeping aspect ratio, I only wanted the shearing to happen in one dimension, either horizontal or vertical, but not both.

I have already discovered that for vertical shearing, the relationship between the new y coordinate y1 and the original y0 coordinate should satisfy the following in order to achieve intended result:

y1 = y0 / cos(a) + x0 * tan(a)
x1 = x0

For horizontal shearing:

x1 = x0 / cos(a) + y0 * tan(a)
y1 = y0

And the above is exactly what my code implements, and I have measured the results and they do indeed keep the aspect ratio.

I just can't write my idea into a transformation matrix. None of the online tutorials I can find covered how to convert math formulas into transformation matrixes, Google search is useless as usual.

Answer 1

It's actually more work to accept a dimension and a reverse. Just accept simultaneous horizontal and vertical angles, and if they're negative, assume reverse.

You don't need a cos; you only need tan in both dimensions. You also don't need flipud and fliplr.

Since you're interested in affine transformations, consider expressing your new h and w as a combined multiplication with a 2x2 transformation matrix.

In no case (sample image generation, transformation) should you be pixel-wise iterating.

You should not operate on the image as a Numpy array during transformation; you should use the pillow built-in affine transformation support.

Add PEP484 type hints.

Move your __main__ symbols into a function, since where they are now they still have global scope.

You've transposed your sample image from what you expected, in all likelihood. Numpy has height first, PIL has width first.

Suggested

This works for positive or negative horizontal or vertical shear. For simultaneous shear, the image is correct but the bounding box is not; I leave this edge case for you to figure out.

import numpy as np
from PIL import Image


def rhombize(
    img: Image,
    h_angle: float = 0,
    v_angle: float = 0,
):
    tans = np.tan(np.radians((v_angle, h_angle)))

    # Transform the old image size to the new image size
    size_transform = np.eye(2)
    size_transform[
        (0, 1),
        (1, 0),
    ] = np.abs(tans)
    new_size = img.size @ size_transform

    '''
    [ a b c ][ x ]   [ x']
    [ d e f ][ y ] = [ y']
    [ 0 0 1 ][ 1 ]   [ 1 ]
    '''
    shear = np.eye(2, 3)

    # Shift the x position based on the y position and vice versa
    shear[
        (1, 0),
        (0, 1),
    ] = tans

    # Shift the x and y positions everywhere; 
    # we need to be piecewise via clip() because in half of the cases
    # our origin needs to shift
    shear[(1, 0), 2] = np.clip(-tans * img.size, a_max=0, a_min=None)

    rhombised = img.transform(
        size=tuple(new_size.astype(int)),
        method=Image.AFFINE,
        data=shear.flatten(),
        # resample omitted - we don't really care about quality
    )

    return rhombised


def test() -> None:
    arr = np.zeros((768, 1024, 3), dtype=np.uint8)
    arr[:, :, 0] = np.linspace(0, 255, 1024, dtype=np.uint8)[np.newaxis, :]  # red ramp, horz
    arr[:, :, 2] = np.linspace(0, 255, 768, dtype=np.uint8)[:, np.newaxis]  # blue ramp, vert
    img1 = Image.fromarray(arr)
    # img1.show()

    img2 = rhombize(img1, h_angle=0, v_angle=-15)
    img2.show()


if __name__ == '__main__':
    test()

Strict equivalence

I still don't pretend to understand the geometric interpretation of what you're doing, but mathematically it's pretty easy to make an affine transformation to represent it. Note that due to a quirk in Pillow, the final affine transformation needs to be an inverse of the size transform. You still need to clip for the negative cases.

import numpy as np
from numpy.linalg import inv
from PIL import Image


def rhombize(
    img: Image,
    x_angle: float = 0,
    y_angle: float = 0,
):
    w0, h0 = img.size
    angles = np.radians((x_angle, y_angle))
    tanx, tany = np.tan(angles)
    cosx, cosy = np.cos(angles)

    '''
    Transform the old image coordinates to the new image coordinates
    [ a b c ][ x ]   [ x']
    [ d e f ][ y ] = [ y']
    [ 0 0 1 ][ 1 ]   [ 1 ]
    '''
    shear = np.array((
        # x col    y col      global col
        (1/cosx,   tanx, max(0, -h0*tanx)),  # -> x' row
        (  tany, 1/cosy, max(0, -w0*tany)),  # -> y' row
        (     0,      0,                1),  # -> 1  row
    ))

    size_transform = np.abs(shear[:2, :2])
    w1, h1 = (size_transform @ img.size).astype(int)

    print(shear)
    '''
    The original implementation was assigning old coordinates to new 
    coordinates on the left-hand side, so this needs to be inverted
    '''
    shear = inv(shear)
    print(shear)

    rhombised = img.transform(
        size=(w1, h1),
        method=Image.AFFINE,
        data=shear[:2, :].flatten(),
        # resample omitted - we don't really care about quality
    )

    return rhombised


def test() -> None:
    arr = np.zeros((768, 1024, 3), dtype=np.uint8)
    arr[:, :, 0] = np.linspace(0, 255, 1024, dtype=np.uint8)[np.newaxis, :]  # red ramp, horz
    arr[:, :, 2] = np.linspace(0, 255, 768, dtype=np.uint8)[:, np.newaxis]  # blue ramp, vert

    img1 = Image.fromarray(arr)
    img2 = rhombize(img1, y_angle=25, x_angle=-35)
    img2.show()


if __name__ == '__main__':
    test()

Simultaneous shear works:

though it appears that you don't want to allow that. The one-dimensional cases simplify your inverse and make it so that the special cases can be solved analytically, instead of numerically through linalg.inv().

import numpy as np
from PIL import Image


def rhombize(
    img: Image,
    x_angle: float = 0,
    y_angle: float = 0,
) -> Image:
    epsilon = 1e-12
    has_x, has_y = np.abs((x_angle, y_angle)) > epsilon
    if has_x and has_y:
        raise ValueError('Two-axis shear is not supported')
    if not (has_x or has_y):
        return img  # no-op

    w0, h0 = img.size
    angles = np.radians((x_angle, y_angle))
    sinx, siny = np.sin(angles)
    cosx, cosy = np.cos(angles)
    tanx, tany = np.tan(angles)

    '''
    Transform the old image coordinates to the new image coordinates
    [ a b c ][ x ]   [ x']
    [ d e f ][ y ] = [ y']
    [ 0 0 1 ][ 1 ]   [ 1 ]
    In the forward case (used for size transformation), the last row and last column are irrelevant.
    In the inverse case (used for image transformation), the last row is implied.
    '''
    shear_forward = np.array((
        (1/cosx,   tanx),
        (  tany, 1/cosy),
    ))
    shear_inverse = np.array((
        ( cosx, -sinx, min(0, h0*tanx)*cosx),
        (-siny,  cosy, min(0, w0*tany)*cosy),
    ))
    size_transform = np.abs(shear_forward)
    w1, h1 = (size_transform @ img.size).astype(int)

    rhombised = img.transform(
        size=(w1, h1),
        method=Image.AFFINE,
        data=shear_inverse.flatten(),
        # resample omitted - we don't really care about quality
    )

    return rhombised


def test() -> None:
    arr = np.zeros((768, 1024, 3), dtype=np.uint8)
    arr[:, :, 0] = np.linspace(0, 255, 1024, dtype=np.uint8)[np.newaxis, :]  # red ramp, horz
    arr[:, :, 2] = np.linspace(0, 255,  768, dtype=np.uint8)[:, np.newaxis]  # blue ramp, vert

    img1 = Image.fromarray(arr)
    img2 = rhombize(img1, y_angle=30)
    img2.show()


if __name__ == '__main__':
    test()