Search code examples
pythonopencvcomputer-visionhomography

Projecting ellipses onto circles using homographies


I am relatively new to computer vision, and I am experimenting with extracting a clock/dial face from an image.

I managed to fit an ellipse onto the dial, and I wanted to fix the perspective such that the dial would face the camera directly. Essentially, I wanted to project an ellipse onto a circle.

I found online that homographies are normally used to accomplish this task.

I tried following the OpenCV tutorial on homographies, but I am having problems with the source poorly matching the destination.

Reading similar questions, a perfect projection seems impossible since perspective-projected circles are not really ellipses. However, I am not looking for a mathematically exact projection, but even the simplest cases I have seem to produce bad results (see example 4 in the images below).

Here are four example input images on which I annotated the ellipses (red), unit circles (blue), source points (yellow), and destination points (cyan) (note: the centers are annotated but are not used to compute the homography). Annotated input image

However, after applying the homography computed using OpenCV (I got the same exact results using Scikit), the ellipses are not projected well onto the unit circles. Homography results

I also tried using more than 4 points to compute the homography, but again, the results were identical.

Here is the code I used:

def persp_transform(orig):
    img = orig.copy()

    ellipses = find_ellipses(img)
    ellipse = ellipses[0] # Use the largest one

    center_x = ellipse[0]
    center_y = ellipse[1]
    center = np.array([center_x, center_y])
    width = ellipse[2]
    height = ellipse[3]
    angle = ellipse[4]
    minor_semiaxis = min(width, height)/2
    major_semiaxis = max(width, height)/2

    # Translate the image to center the ellipse
    img_center = np.array([img.shape[0]//2, img.shape[1]//2])
    translation = center - img_center
    M = np.float32([[1, 0, -translation[1]], [0, 1, -translation[0]]])
    img = cv.warpAffine(img, M, (img.shape[1], img.shape[0]))

    # Draw ellipse before projection
    rr, cc = ellipse_perimeter(img_center[0], img_center[1], int(major_semiaxis), int(minor_semiaxis), orientation=angle, shape=img.shape)
    draw_ellipse(img, rr, cc, color=0, thickness=10)
    
    def sin(a): return np.sin(np.deg2rad(a))
    def cos(a): return np.cos(np.deg2rad(a))
    
    # Source points around the ellipse
    num_points = 5
    rotation = np.array([[np.cos(angle), -np.sin(angle)], [np.sin(angle), np.cos(angle)]])
    from_points = np.array([[major_semiaxis*sin(a), minor_semiaxis*cos(a)] for a in np.linspace(0, 360, num_points)]) @ rotation + img_center
    from_points = from_points.astype(np.float32)

    # Destination points around a centered circle
    to_radius = int((min(img.shape[:2]) / 2) * 0.8)
    to_points = np.array([[to_radius*sin(a), to_radius*cos(a)] for a in np.linspace(0, 360, num_points)]) @ rotation + img_center
    to_points = to_points.astype(np.float32)

    # Draw ellipse center and source points before projection
    cv.circle(img, (img_center[1], img_center[0]), 20, (255, 255, 0), -1)
    for fp in from_points:
        cv.circle(img, (int(fp[1]), int(fp[0])), 20, (255, 255, 0), -1)
    
    # Compute homography and project
    M, _ = cv.findHomography(from_points, to_points, cv.RANSAC, 5.0)
    img = cv.warpPerspective(img, M, (img.shape[1], img.shape[0]))

    # Draw target unit circle and destination points after projection
    cv.circle(img, (img_center[1], img_center[0]), to_radius, (0, 0, 255), 20)
    cv.circle(img, (img_center[1], img_center[0]), 20, (0, 255, 255), -1)
    for tp in to_points:
        cv.circle(img, (int(tp[1]), int(tp[0])), 20, (0, 255, 255), -1)

    return img

EDIT 1

Here are the 4 raw example images I used as input: GDrive Link

To fit the ellipses, I use the method proposed in AAMED, which can be found on GH here.

I manually played with the input parameters until I found satisfactory ones since I didn't need to find all the ellipses in the image, only the main ones.

My ellipse finding code looks like this:

def find_ellipses(orig, scale=0.3, theta_fsa=20, length_fsa=5.4, T_val=0.9):
    img = cv.resize(orig.copy(), None, fx=scale, fy=scale)
    gray_image = cv.cvtColor(img, cv.COLOR_BGR2GRAY)

    aamed = AAMED(img.shape[0]+1, img.shape[1]+1)
    aamed.setParameters(np.deg2rad(theta_fsa), length_fsa, T_val)
    ellipses = aamed.run_AAMED(gray_image)
    aamed.release()

    ellipses = sorted(ellipses, key=lambda e: ellipse_area(e), reverse=True)
    ellipses = [np.array(e) / scale for e in ellipses]

    return ellipses

Then I use skimage.draw.ellipse_perimeter(...) to impose the largest ellipse on the image.

Edit 2

I tried to better visualize the homography and narrowed down the issue possibly to cv.warpPerspective(...).

Firstly, I color-coded each point pair to check that I was mapping the right pairs: Color-coded points

Then, I computed the homography and individually projected each source point using cv.perspectiveTransform(...) without warping the image to see where they would end up, and they seemed to be projected properly: Projected points

However, when I then use cv.warpPerspective(...) to project the image, the projection is inconsistent with the result from cv.perspectiveTransform(...): Wrong projection


Solution

  • Here is one attempt in Python/OpenCV. The idea is the get the ellipse. Then get the 4 points on the ellipse that corresponds to the min and max semi-major radii. Then get the 4 corresponding points on a circle or radius equal to the max semi-major radii. Then compute the homography for the two sets of 4 corresponding points. Finally warp the image.

    • Read the input
    • Read the ellipse mask
    • Get the contour of the mask
    • Fit and ellipse to the contour
    • Compute the 4 points on the ellipse from the ellipse parameters
    • Compute the 4 points on the circle using the max semi-major radius
    • Get the homography between the two sets of points
    • Warp the input image
    • Save the results

    Input:

    enter image description here

    Ellipse Mask (fit manually in Photoshop):

    enter image description here

    import cv2
    import numpy as np
    import math
    
    # read the input
    img = cv2.imread('ex_4.jpg')
    
    # read ellipse mask as grayscale
    mask = cv2.imread('ex_4_mask.png', cv2.IMREAD_GRAYSCALE)
    
    # get ellipse contour
    contours = cv2.findContours(mask, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    contours = contours[0] if len(contours) == 2 else contours[1]
    big_contour = max(contours, key=cv2.contourArea)
    
    # draw contour on black background
    mask2 = np.zeros_like(mask, dtype=np.uint8)
    cv2.drawContours(mask2, [big_contour], 0, 255, 1)
    
    # fit ellipse to mask contour
    points = np.argwhere(mask2.transpose()>0)
    ellipse = cv2.fitEllipse(points)
    (cx,cy), (d1,d2), angle = ellipse
    
    # draw ellipse
    ellipse = img.copy()
    cv2.ellipse(ellipse, (int(cx),int(cy)), (int(d1/2),int(d2/2)), angle, 0, 360, (0,0,255), 4)
    
    # correct angle
    if angle > 90:
        angle = angle - 90
    else:
        angle = angle + 90
    print('center: (', cx,cy, ')', 'diameters: (', d1, d2, ')', 'angle:', angle)
    
    # get 4 points from min and max semi-major radii
    r1 = min(d1,d2)/2
    r2 = max(d1,d2)/2
    
    x1 = np.intp( cx + math.cos(math.radians(angle))*r2 )
    y1 = np.intp( cy + math.sin(math.radians(angle))*r2 )
    x2 = np.intp( cx + math.cos(math.radians(angle+180))*r2 )
    y2 = np.intp( cy + math.sin(math.radians(angle+180))*r2 )
    x3 = np.intp( cx + math.cos(math.radians(angle+90))*r1 )
    y3 = np.intp( cy + math.sin(math.radians(angle+90))*r1 )
    x4 = np.intp( cx + math.cos(math.radians(angle+270))*r1 )
    y4 = np.intp( cy + math.sin(math.radians(angle+270))*r1 )
    
    print('x1,y1:', x1, y1)
    print('x2,y2:', x1, y1)
    print('x3,y3:', x1, y1)
    print('x4,y4:', x1, y1)
    
    # input points from ellipse
    input_pts = np.float32([[x1,y1], [x2,y2], [x3,y3], [x4,y4]])
    
    # output points from circle of radius = max semi-major radii
    ox1 = np.intp( cx + math.cos(math.radians(angle))*r2 )
    oy1 = np.intp( cy + math.sin(math.radians(angle))*r2 )
    ox2 = np.intp( cx + math.cos(math.radians(angle+180))*r2 )
    oy2 = np.intp( cy + math.sin(math.radians(angle+180))*r2 )
    ox3 = np.intp( cx + math.cos(math.radians(angle+90))*r2 )
    oy3 = np.intp( cy + math.sin(math.radians(angle+90))*r2 )
    ox4 = np.intp( cx + math.cos(math.radians(angle+270))*r2 )
    oy4 = np.intp( cy + math.sin(math.radians(angle+270))*r2 )
    
    cx = np.intp(cx)
    cy = np.intp(cy)
    r2 = np.intp(r2)
    ellipse_pts = ellipse.copy()
    
    # draw white circle on copy of ellipse
    cv2.circle(ellipse_pts, (cx,cy), r2, (255,255,255), 4)
    
    output_pts = np.float32([[ox1,oy1], [ox2,oy2], [ox3,oy3], [ox4,oy4]])
    
    # draw output points on copy of ellipse image
    cv2.circle(ellipse_pts, (ox1,oy1), 16, (0,255,0), 4)
    cv2.circle(ellipse_pts, (ox2,oy2), 16, (0,255,0), 4)
    cv2.circle(ellipse_pts, (ox3,oy3), 16, (0,255,0), 4)
    cv2.circle(ellipse_pts, (ox4,oy4), 16, (0,255,0), 4)
    
    
    # draw input points on copy of ellipse image
    cv2.circle(ellipse_pts, (x1,y1), 12, (255,0,0), -1)
    cv2.circle(ellipse_pts, (x2,y2), 12, (255,0,0), -1)
    cv2.circle(ellipse_pts, (x3,y3), 12, (255,0,0), -1)
    cv2.circle(ellipse_pts, (x4,y4), 12, (255,0,0), -1)
    
    # get homography when only 4 pts
    h = cv2.getPerspectiveTransform(input_pts, output_pts)
    
    # warp image
    result = cv2.warpPerspective(img, h, (img.shape[1],img.shape[0]))
    
    
    # save results (compress to fit 2Mb limit on this forum)
    cv2.imwrite('ex_4_ellipse.jpg', ellipse, [int(cv2.IMWRITE_JPEG_QUALITY), 85])
    cv2.imwrite('ex_4_ellipse_pts.jpg', ellipse_pts, [int(cv2.IMWRITE_JPEG_QUALITY), 85])
    cv2.imwrite('ex_4_ellipse2circle.jpg', result, [int(cv2.IMWRITE_JPEG_QUALITY), 85])
    
    # show results
    cv2.imshow('ellipse', ellipse)
    cv2.imshow('ellipse_pts', ellipse_pts)
    cv2.imshow('result', result)
    cv2.waitKey(0)
    

    Ellipse on input:

    enter image description here

    Ellipse with points and circle:

    enter image description here

    Warped Result:

    enter image description here