Motivation
EE5806: Topics in Digital Image Processing Lecture Notes: Landmark registration
• Wish to compare source image with target image (Fig. 1).
• Need to align (i.e., register) source with target before comparison.
Approach
• We will deal with cases where the source image is a rotated, translated and scaled version of the target. Therefore, we need to find a geometric transformation consisting of rotation, translation and scaling to align the source with the target.
• In landmark registration, we identify landmarks in the source image and corresponding landmarks in the target image.
o Landmarks:Distinctivefeaturesthatcanbeaccuratelyandreproduciblyidentified in both the source and target images.
• Landmarks are used to find an “optimal” geometric transformation to register the two images. In the following, we will define “optimality” mathematically.
Fig. 1
Page 1 of 5
Mathematical Derivation
Given: N landmarks (n = 1, 2, …., N) in source image and corresponding landmarks in
the target image (n = 1, 2, …., N). In 2D images:
and
Problem: Find a geometric transformation (specifically, one consisting of rotation, translation and scaling) that aligns the two images based on the two sets of landmarks. We will assume the scaling in i and j directions are the same (i.e., )
Solution:
Fig. 2. Source image before and after transformation and the target image. denotes the landmarks chosen in the source image (blue crosses in left panel) and denotes
the corresponding landmarks chosen in the target image (red crosses in the right panel). After the transformation T has been applied, becomes (blue crosses in the right panel).
So, there are infinitely many transformation, how do we choose a “best” transformation to perform? Before answering this question, you would ask: “After applying a certain transformation T, how close are the transformed source landmarks to the target landmarks .” The “closeness” is quantified by the following cost function:
Define a cost function, C, that is the mean-squared distance between the target landmarks, , and transformed source landmarks, . Our goal is to minimize this cost function:
where denotes the magnitude of the vector v, and
Page 2 of 5
where and .
Now, we express C in terms of :
To minimize this cost function, compute its derivative with respect to a, b, ti and tj and set to zero. Example for a:
Performing differentiation with respect to b, ti and tj as well gives the following set of linear equations:
where
This system of linear equations can be expressed in matrix notation where
,,
Solving for v will define our optimal geometric transformation T. Page 3 of 5
Python implementation
import numpy as np
import cv2
import matplotlib.pyplot as plt
def landmark_register(fnameS=’imS.jpg’, fnameT=’imT.jpg’): “””
landmark_register: Landmark registration of two images. Takes into account translation, rotation and scaling only.
:param fnameS: a string containing the filename of the SOURCE image :param fnameT: a string containing the filename of the TARGET image Note: If you use Pycharm and cannot get the landmarks successfully,
please uncheck ‘Show Plots in Tool window’
Windows: Settings \ Tools \ Python Scientific \ Show Plots in Tool
window
MacOS: Preferences \ Tools \ Python Scientific \ Show Plots in Tool
window “””
# Read in images and display side-by-side; get landmarks source = cv2.imread(fnameS, cv2.IMREAD_GRAYSCALE)
target = cv2.imread(fnameT, cv2.IMREAD_GRAYSCALE) height, width = target.shape
plt.ion()
plt.subplot(131)
plt.imshow(source, cmap=’gray’, vmin = 0, vmax = 255) plt.axis(‘off’)
plt.title(‘Source image’, {‘fontsize’: 15, ‘fontweight’: print(‘Select landmarks in source image. Press ENTER key ps = plt.ginput(n=-1, timeout=0, show_clicks=True)
ps = np.asarray(ps)
plt.plot(ps[:, 0], ps[:, 1], ‘r+’)
plt.subplot(132)
plt.imshow(target, cmap=’gray’, vmin = 0, vmax = 255) plt.axis(‘off’)
plt.title(‘Target image’, {‘fontsize’: 15, ‘fontweight’: print(‘Select landmarks in target image. Press ENTER key pt = plt.ginput(n=-1, timeout=0, show_clicks=True)
pt = np.asarray(pt)
plt.plot(pt[:, 0], pt[:, 1], ‘r+’)
plt.ioff()
# Set up the matrix B and the vector k
N = len(ps)
Sj = np.sum(ps[:, 0])
Si = np.sum(ps[:, 1])
Sjj = np.sum(ps[:, 0] ** 2)
Sii = np.sum(ps[:, 1] ** 2)
Siip = np.sum(ps[:, 1] * pt[:, 1])
Sjjp = np.sum(ps[:, 0] * pt[:, 0])
Sjip = np.sum(ps[:, 0] * pt[:, 1])
Sjpi = np.sum(pt[:, 0] * ps[:, 1])
Sjp = np.sum(pt[:, 0])
Sip = np.sum(pt[:, 1])
B = np.array([[Sii + Sjj, 0, Si, Sj], [0, Sii + Sjj, Sj,
Si, 0, N], [Si, Sj, N, 0]])
k = np.array([Siip + Sjjp, Sjip – Sjpi, Sjp, Sip]).T
# Solve least-squares problem: Bv = k using the inv(B)*k v = np.dot(np.linalg.inv(B), k)
a = v[0]
b = v[1]
ti = v[2]
‘bold’})
when done.’)
‘bold’})
when done.’)
Page 4 of 5
-Si], [Sj, –
tj = v[3]
# Perform the transformation using built-in function. # Step 1: Define affine matrix, T
T = np.array([[a, -b, tj], [b, a, ti]])
# Step 2: Perform the transformation.
im2 = cv2.warpAffine(source, T, (width, height))
plt.subplot(133)
plt.imshow(im2, cmap=’gray’, vmin = 0, vmax = 255)
plt.axis(‘off’)
plt.title(‘Registered image’, {‘fontsize’: 15, ‘fontweight’: ‘bold’}) plt.show()
Page 5 of 5