Geometric warp of image in python

Question

I would like to use python to perform a geometric transform over an image, to 'straighten' or rectify an image along a given curve. It seems that scikit-image ProjectiveTransform() and warp() are very good for this, but the documentation is sparse. I followed the documentation here, but I couldn't get it to work properly for a sample case.

Here's an example: I'll create an image with two concentric circles, and the goal is to rectify one quarter of these circles, so that the resulting image are two parallel lines. Here is the sample data:

import numpy as np
a = np.zeros((500, 500))

# create two concentric circles with a thickness of a few pixels:
for i in range(500):
    for j in range(500):
        r = np.sqrt((i - 250)**2 + (j - 250)**2) 
        if r > 50 and r < 52:
            a[i, j] = 10
        if r > 100 and r < 102:
            a[i, j] = 10
# now create the coordinates of the control points in the original image:
(x0, y0) = (250, 250)
r = 30   # inner circle
x = np.linspace(250 - r, 250, 50)
y = np.sqrt(r ** 2 - (x - x0) ** 2) + x0
r2 = 120   # outer circle
x2 = np.linspace(250 - r2, 250, 50)
y2 = np.sqrt(r2 ** 2 - (x2 - x0) ** 2) + x0
dst = np.concatenate((np.array([x, y]).T, np.array([x2, y2]).T))

And this can be plotted, e.g.:

imshow(a, cmap='gist_gray_r')
plot(x, y, 'r.')
plot(x2, y2, 'r.')

So my goal is to rectify the image in the quadrant given by the red control points. (In this case, this is the same as a Cartesian to polar transformation.) Using scikit image from the documentation example, I've done:

# create corresponding coordinates for control points in final image:
xi = np.linspace(0, 100, 50)
yi = np.zeros(50)
xi2 = xi
yi2 = yi + (r2 - r)
src = np.concatenate((np.array([xi, yi]).T, np.array([xi2, yi2]).T))

# transform image
from skimage import transform, data
tform3 = transform.ProjectiveTransform()
tform3.estimate(src, dst)
warped = transform.warp(a, tform3)

I was expecting this warped image to show two parallel lines, but instead I get:

What am I doing wrong here?

Note that while in this case it is a Cartesian to polar transform, in the most general case I'm looking for a transformation from some arbitrary curve. If someone knows of a better way using some other package, please let me know. I can solve this problem by using ndimage.map_coordinates for a bunch of radial lines, but was looking for something more elegant.

Answer 1

A ProjectiveTransform is a linear transformation, and cannot match your deformation scheme. There may be better options, but for arbitrary curves you can make it work with a PiecewiseAffineTransform, which will match anything you throw at it by tessellating linear transformations. If you simply change the name of the transform in your code, this is the output I get:

So you'll probably need to tweak it a little bit to get what you are after, but at least it produces the two parallel lines you were expecting in the area where your transformation is well defined.