← All notebooks

June 16, 2022 16 min read

Image Warping

Implementations of various image warping methods including fish eye, swirl, waves, cylinder anamorphosis, radial blur, bilinear warping, and perspective mapping effects.

Image Warping¶

In [1]:

import numpy as npimport imageioimport matplotlib.pyplot as pltimport scipy.ndimage as img

In [2]:

def imageRead(imgname, pilmode='L', arrtype=np.float):    """    pilmode: str        for luminance / intesity images use 'L'        for RGB color images use 'RGB'    arrtype: numpy dtype        use np.float, np.uint8, ...    """    return imageio.imread(imgname, pilmode=pilmode).astype(arrtype)def imageWrite(arrF, imgname, arrtype=np.uint8):    imageio.imwrite(imgname, arrF.astype(arrtype))

/tmp/ipykernel_2184682/3362905083.py:1: DeprecationWarning: `np.float` is a deprecated alias for the builtin `float`. To silence this warning, use `float` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.float64` here.
Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations
  def imageRead(imgname, pilmode='L', arrtype=np.float):

In [3]:

arrF = imageRead("Data/asterix.png")_ = plt.imshow(arrF, cmap="gray")

A transformation T alters the geometry of an image. If function f[x,y] denotes the original of an image and g[u,v] denotes the warped image and if $ \big[u,v \big] = \big[T_x(x,y), T_y(x,y) \big] $, then $g\big[u\big] = g\big[T(x)\big] = f\big[x\big]$

Forward (PUSH) Warps¶

In [4]:

def push_warp(arrF):    M, N = arrF.shape    arrG = np.zeros_like(arrF)    for y in range(M):        for x in range(N):            u = int(4.0 * x)            v = y            if u<N:                arrG[v, u] = arrF[y, x]    return arrGarrF = imageRead("Data/asterix.png")arrG = push_warp(arrF)_ = plt.imshow(arrG, cmap="gray")

Push warping caused leapfrogs over the coordinates in the destination image. Some coordinates did not get values from the source image.

To solve this problem, the corresponding value for each coordinate in the destination image is calculated from the source image.

$ u = T(x) \Leftrightarrow x = T^{-1}(u) $

Backward (PULL) Warps¶

$ x = T^{-1}(u) $

Leapfrog problem does not occur with this procedure. Corresponding value from the original function is obtained for each pixel in the destination image. However, some values might have to be read from a location situated in between the discrete pixel coordinates of the original image. Therefore, pull warping entails the use of interpolation methods.

In [5]:

def pull_warp(arrF, interpolation="nearest"):    M, N = arrF.shape    arrG = np.zeros_like(arrF)    for v in range(M):        for u in range(N):            x = u / 4.0            y = v                        if interpolation=="nearest":                nearest_x = np.round(x)                  arrG[v,u] = arrF[y, int(nearest_x)]                        elif interpolation=="linear":                x_up = np.ceil(x)                x_down = np.floor(x)                if x_up != x_down:                    delta = (x - x_down) * (arrF[y, int(x_up)] - arrF[y, int(x_down)]) / (x_up - x_down)                    arrG[v,u] = arrF[y, int(x_down)] + int(delta)                else:                    arrG[v,u] = arrF[y, int(x)]                return arrGarrF = imageRead("Data/asterix.png")arrG = pull_warp(arrF, interpolation="linear")_ = plt.imshow(arrG, cmap="gray")

In [6]:

# Without For Loopdef pull_warp_2(arrF):    M, N = arrF.shape    u, v = np.meshgrid(M, N)        arrG = np.zeros_like(arrF)         x = u / 4.0    y = v    # Nearest-neighbor interpolation    nearest_x = np.round(x)      arrG[v,u] = arrF[y, int(nearest_x)]                        return arrGarrF = imageRead("Data/asterix.png")arrG = pull_warp(arrF, interpolation="linear")_ = plt.imshow(arrG, cmap="gray")

Fisheye Effect¶

Notations:
c -> center of force
$ u = T(x)$
$ x = T^{-1}(u) = W(u)$

Approach (1D):

• Understanding $u \ge 0$ as the distance to the center of the lens, we require a warp function where
  $ W(u) \simeq u $ , if u is large
  $ W(u) < u $ , if u is small

• We may either add something to u which is negative for small u and equals 0 for large u: $ W(u) = u + \delta^+(u) $

• Or, we may multiply u by something smaller than 1 for small u and equal to 1 for large u: $ W(u) = u . \delta^.(u) $

The naive approach (2D):

1. Chose the center c of the virtual lens
2. Compute polar coordinates of u w.r.t. c
   $ [u,v] \rightarrow [r, \varphi]$
3. Warp the radial coordinate
   $ r' = W(r) = r + \delta^+ (r) $
4. Compute the warped Euclidean coordinates x
   $ [r', \varphi] \rightarrow [x,y] $

In the following code, instead of following the naive approach (2D), we follow the much better approach. We think like a physicist and understand W to behave like a central force.

Much better approach (2D):

• We consider warps of the following form: $ W(u,c,\sigma) = u + \delta^+(r, \sigma) $
  where $r = c-u$ denotes the vector pointing from u towards the center of the force c
• If we follow this idea, then u will be displaced into the direction r
• The magnitude of this displacement is proportional to $ r = ||r||$
• In the following, we consider the effect of different functions $\delta^+ (r, \sigma) = r . \delta(r,\sigma)$

---

In the following fisheye effect function, one of the following delta ($\delta$) functions is used. The formulas of the delta functions are provided below.

\begin{equation*} \delta_1 = \begin{cases} 1 - \frac{r}{\sigma} & \text{ if } r \le \sigma \\ 0 & \text{ otherwise} \end{cases} \end{equation*}

\begin{equation*} \delta_2 = \begin{cases} 1 - \frac{r^2}{\sigma^2} & \text{ if } r \le \sigma \\ 0 & \text{ otherwise} \end{cases} \end{equation*}

\begin{equation*} \delta_3 = \begin{cases} \frac{\sqrt{\sigma^2-r^2}}{\sigma} & \text{ if } r \le \sigma \\ 0 & \text{ otherwise} \end{cases} \end{equation*}

$$ \delta_4 = 1 -tanh\bigg(\frac{r}{\sigma}\bigg) $$

$$ \delta_5 = exp \bigg(-\frac{r^2}{2 \sigma^2} \bigg) $$

In [7]:

def delta1(r, sigma):    return np.where(r < sigma, 1 - r/sigma, 0)def delta2(r, sigma):    return np.where(r < sigma, 1 - r**2/sigma**2, 0)def delta3(r, sigma):    return np.where(r < sigma, np.sqrt(np.abs(sigma**2 - r**2)) / sigma, 0)def delta4(r, sigma):    return 1 - np.tanh(r/sigma)def delta5(r, sigma):    return np.exp(-0.5 * (r/sigma)**2)

In [8]:

def fisheye(arrF, vecC, sigma=100., dfct=delta1):    vecC = np.array(vecC).reshape(-1,1)    M, N = arrF.shape    u, v = np.meshgrid(np.arange(N), np.arange(M))    matX = np.vstack((v.flatten(), u.flatten())).astype(float) # u    matR = vecC - matX                      # vectors pointing to center # r = c-u    dist = np.sqrt(np.sum(matR**2, axis=0)) # distances to center # ||r||    matX = matX + matR * dfct(dist, sigma)  # W(u,c,sigma)= u + r*delta        arrG = img.map_coordinates(arrF, matX)  # matX MUST be float    arrG = arrG.reshape(M,N)    return arrG

In [9]:

white_cols = np.linspace(0, arrF.shape[1]-16, 16).astype(int)white_rows = np.linspace(0, arrF.shape[0]-16, 16).astype(int)black_img = np.zeros_like(arrF)for i in range(4):    black_img[white_rows+i]=255    black_img[:, white_cols+i]=255#_ = plt.imshow(black_img / 255, cmap='gray')fig, axs = plt.subplots(2,3, figsize=(20,10))axs = axs.flatten()axs[0].imshow(black_img / 255, cmap="gray", vmin=0, vmax=1)axs[0].title.set_text('Original Image')out = fisheye(black_img,(black_img.shape[0]//2,black_img.shape[1]//2), sigma=300, dfct=delta1)axs[1].title.set_text('Delta 1'); axs[1].imshow(out / 255, cmap="gray", vmin=0, vmax=1)out = fisheye(black_img,(black_img.shape[0]//2,black_img.shape[1]//2), sigma=300, dfct=delta2)axs[2].title.set_text('Delta 2'); axs[2].imshow(out / 255, cmap="gray", vmin=0, vmax=1)out = fisheye(black_img,(black_img.shape[0]//2,black_img.shape[1]//2), sigma=300, dfct=delta3)axs[3].title.set_text('Delta 2'); axs[3].imshow(out / 255, cmap="gray", vmin=0, vmax=1)out = fisheye(black_img,(black_img.shape[0]//2,black_img.shape[1]//2), sigma=300, dfct=delta4)axs[4].title.set_text('Delta 4'); axs[4].imshow(out / 255, cmap="gray", vmin=0, vmax=1)out = fisheye(black_img,(black_img.shape[0]//2,black_img.shape[1]//2), sigma=300, dfct=delta5)axs[5].title.set_text('Delta 5'); _ = axs[5].imshow(out / 255, cmap="gray", vmin=0, vmax=1)fig.tight_layout()

In [10]:

arrF = imageRead('Data/asterix.png')fig, axs = plt.subplots(1, 2, figsize=(15,15))out = fisheye(arrF, (400, 650))axs[0].imshow(out / 255, cmap="gray", vmin=0, vmax=1)axs[0].title.set_text('Center of force: (400,650)')out = fisheye(arrF, (350, 200), sigma=300)axs[1].imshow(out / 255, cmap="gray", vmin=0, vmax=1)_ = axs[1].title.set_text('Center of force: (350,200)')

The $L_p$-norm $\lVert \vec{x} \rVert_p$ of a vector $\vec{x} \in \mathbb{R}^m$ is defined as $$ \begin{equation} \label{eq:lpnorm} \lVert \vec{x} \rVert_p = \left( \sum_{i=1}^m \,\lvert x_i \rvert^p \right)^{\frac{1}{p}} \end{equation} $$

In the following, we will use it to compute PULL warps of an image where we consider $\vec{x} = W \bigl(\vec{u}, \vec{c}, \sigma, p \bigr)$ with

\begin{align} W \bigl(\vec{u}, \vec{c}, \sigma, p \bigr) & = \vec{u} + \bigl[ \vec{c} - \vec{u} \bigr] \cdot \exp \left( - \frac{\lVert \vec{c} - \vec{u} \rVert_p^2}{2 \, \sigma^2} \right) \label{eq:warp} \end{align}

In [11]:

def lpNorm(matX, p):    return np.power(np.sum(np.power(np.abs(matX), p), axis=0), 1/p)def imageWarp(arrF, vecC, sigma, p):    vecC = np.array(vecC).reshape(-1,1)    M, N = arrF.shape    u, v = np.meshgrid(np.arange(N), np.arange(M))        matX = np.vstack((v.flatten(), u.flatten())).astype(float)     matR = vecC - matX     # vectors pointing to center    matX = matX + matR * np.exp(-lpNorm(matR, p)**2 / (2*sigma**2))        arrG = img.map_coordinates(arrF, matX)    arrG = arrG.reshape(M,N)    return np.clip(arrG,0,255)arrF = imageRead("Data/lena.png")vecC = np.array([350,300])  #np.array([150,128])sigma = [50, 50, 100]p_values = [3, 12, 12]fig, axs = plt.subplots(1, 3, figsize=(15,15))for i, (s,p) in enumerate(zip(sigma, p_values)):    out = imageWarp(arrF, vecC=vecC, sigma=s, p=p)    axs[i].imshow(out/255, cmap="gray")

In [12]:

arrF = imageRead("Data/cat1.png")print("image.shape: ", arrF.shape)vecC = np.array([128,116])sigma = [50, 75, 100]p_values = [3, 12, 24]fig, axs = plt.subplots(1, 3, figsize=(15,15))for i, (s,p) in enumerate(zip(sigma, p_values)):    out = imageWarp(arrF, vecC=vecC, sigma=s, p=p)    axs[i].imshow(out/255, cmap="gray")

image.shape:  (256, 232)

Swirl Effect¶

• Compute the polar coordinates of u w.r.t. c
  $$ [u,v] \rightarrow [r, \varphi]$$

• Warp the angular coordinate
  $$ \varphi' = W(\varphi, r) = \varphi + \delta^+(\varphi,r)$$

• Compute warped Euclidean coordinates $$ [r, \varphi'] \rightarrow [x,y] $$

In [13]:

def custom_swirl_effect(arrF):    M, N = arrF.shape    u, v = np.meshgrid(np.arange(N), np.arange(M))    matX = np.stack((v, u)).astype(float)        # compute polar coordinates with respect to vecC    vecC = np.array([350, 200]).reshape(2,1,1)    diff = matX - vecC     r = np.linalg.norm(diff, axis=0)    angle = np.arctan2(diff[1], diff[0])        dist = r/r.max()    sigma = 0.1    gaussian = np.exp(- dist**2 / (2*(sigma**2)))    plt.imshow(gaussian)    plt.title("gaussian")    plt.show()        magnitude = 6    angle += magnitude*gaussian        # compute euclidian coordinates with respect to image zero    matX = np.stack([r * np.cos(angle) , r * np.sin(angle)])    matX += vecC    arrG = img.map_coordinates(arrF, matX) # matX MUST be float    arrG = arrG.reshape(M,N)    return arrGarrF = imageRead('Data/asterix.png')out = custom_swirl_effect(arrF)_ = plt.imshow(out / 255, cmap="gray")

Waves Effect¶

Transform the image such that it appears as if on the bottom of a swimming pool whose water is in motion

Approach (1D):

$ W(u) = u + \alpha . sin(v u -\phi) $ where the parameters $ \alpha, v$ and $ \phi$ denote amplitude, frequency, and phase of the water wave question

Appoach (2D): We follow the idea of a central force and choose:

$ W(u,c,\alpha, v, \phi) = u + \alpha . sin(v ||r|| - \phi) . \frac{r}{||r||} $
$ W(u,c,\alpha, v, \phi) = u + \delta^+(r, \alpha, v, \phi) $

In [14]:

def waves_effect(arrF, amplitude, frequency, phase, debug=False):    M, N = arrF.shape    u, v = np.meshgrid(np.arange(N+amplitude[1]*2), np.arange(M+amplitude[0]*2))    matX = np.stack((v, u)).astype(float)          # for j axis    b = amplitude[1] * np.sin(matX[0]/frequency[1] + phase[1])    matX[1] += b - amplitude[1]        # for i axis    a = amplitude[0] * np.sin(matX[1]/frequency[0] + phase[0])     matX[0] += a - amplitude[0]        if debug:        plt.imshow(a); plt.show()        plt.imshow(b); plt.show()        arrG = img.map_coordinates(arrF, matX) # matX MUST be float    arrG = arrG.reshape(matX.shape[1],matX.shape[2])    return arrG

In [15]:

arrF = imageRead('Data/lena.png')out = waves_effect(arrF, amplitude=[10,7], frequency=[10.0,6.5], phase=[0,2], debug=True)plt.imshow(out / 255, cmap="gray")plt.axis("off")plt.show()

In [16]:

arrF = imageRead('Data/lena.png')out1 = waves_effect(arrF, amplitude=[10,10], frequency=[20,20], phase=[5,0])out2 = waves_effect(arrF, amplitude=[10,7], frequency=[10.0,6.5], phase=[0,3.5])out3 = waves_effect(arrF, amplitude=[10,10], frequency=[5.0,40], phase=[0,np.pi])fig, axs = plt.subplots(1, 3, figsize=(15,15))axs[0].imshow(out1 / 255, cmap="gray"); axs[0].axis("off")axs[1].imshow(out2 / 255, cmap="gray"); axs[1].axis("off")axs[2].imshow(out3 / 255, cmap="gray"); _ = axs[2].axis("off")

In [17]:

arrF = imageRead('Data/clock.jpg')out1 = waves_effect(arrF, amplitude=[10,10], frequency=[20,20], phase=[5,0])out2 = waves_effect(arrF, amplitude=[10,7], frequency=[10.0,6.5], phase=[0,3.5])out3 = waves_effect(arrF, amplitude=[10,10], frequency=[5.0,40], phase=[0,np.pi])fig, axs = plt.subplots(1, 4, figsize=(15,10))axs[0].imshow(arrF / 255, cmap="gray"); axs[0].axis("off")axs[1].imshow(out1 / 255, cmap="gray"); axs[1].axis("off")axs[2].imshow(out2 / 255, cmap="gray"); axs[2].axis("off")axs[3].imshow(out3 / 255, cmap="gray"); _ = axs[3].axis("off")#To save the images#imageWrite(out1, save_dir + "Task6-2.png")

Cylinder Anamorphosis¶

In [18]:

def cylinder(arrF, debug=False):    M, N = arrF.shape    u, v = np.meshgrid(np.arange(N), np.arange(M))    matX = np.stack((v, u)).astype(float)        # compute polar coordinates with respect to vecC    vecC = np.array([arrF.shape[0]//2, arrF.shape[1]//2]).reshape(2,1,1)    diff = matX - vecC         r = np.linalg.norm(diff, axis=0)    #y = (1 - r/(arrF.shape[0]//2)) * (arrF.shape[0]-1)    y = (r/(arrF.shape[0]//2)) * (arrF.shape[0]-1)         angle = np.arctan2(diff[0], diff[1])    angle = angle-angle.min() # min angle is 0 with this line    angle = angle/angle.max() # angle is normalized to 0-1    x = angle * (arrF.shape[1]-1)        if debug:        plt.imshow(y); plt.show()        plt.imshow(x); plt.show()    matX = np.stack([y, x])        #arrG = img.map_coordinates(arrF, matX) # matX MUST be float    arrG = img.map_coordinates(np.flipud(arrF), matX) # matX MUST be float    arrG = arrG.reshape(M,N)    return arrGarrF = imageRead("Data/flower.png")out = cylinder(arrF, debug=True)plt.imshow(out / 255, cmap="gray")plt.show()

In [19]:

def cylinder_with_hole(arrF, debug=False):    M, N = arrF.shape    u, v = np.meshgrid(np.arange(min(M, N)), np.arange(min(M, N)))    matX = np.stack((v, u)).astype(float)        # compute polar coordinates with respect to vecC    vecC = np.array([min(M, N)//2, min(M, N)//2]).reshape(2,1,1)    diff = matX - vecC         r = np.linalg.norm(diff, axis=0)    #r = (1 - r/(arrF.shape[0]//2)) * (arrF.shape[0]-1)+60    #r = (1 - r/(arrF.shape[0]//2-30)) * (arrF.shape[0]-1) + 60    r = (1 - r/(min(M, N)//2-30)) * (min(M, N)-1) + 60        angle = np.arctan2(diff[0], diff[1])    angle = angle-angle.min()    angle = angle/angle.max() * (arrF.shape[1]-1)        if debug:        fig, axs = plt.subplots(1, 3, figsize=(15,15))        axs[0].imshow(r); axs[0].title.set_text("r")        axs[1].imshow(angle); axs[1].title.set_text("angle")            matX = np.stack([r, angle])        arrG = img.map_coordinates(arrF, matX) # matX MUST be float    arrG = arrG.reshape(min(M, N),min(M, N))    return arrGarrF_ = imageRead('Data/asterix.png')out = cylinder_with_hole(arrF_, debug=True)plt.imshow(out / 255, cmap="gray"); plt.title("warped image")plt.show()

In [20]:

arrF = imageRead('Data/flower.png')out = cylinder_with_hole(arrF)fig, axs = plt.subplots(1, 2, figsize=(10,10))axs[0].imshow(arrF / 255, cmap="gray"); plt.axis("off"); plt.axis("off"); axs[0].title.set_text("original image")axs[1].imshow(out / 255, cmap="gray"); plt.axis("off"); plt.axis("off"); axs[1].title.set_text("warped image")

Radial Blur Effect¶

In [21]:

def xy2rphi(x, y):    r = np.sqrt(x**2 + y**2)    phi = np.arctan2(y, x)    return r, phidef rphi2xy(r, phi):    x = r * np.cos(phi)    y = r * np.sin(phi)    return x, ydef to_r_phi_plane_(f, m, n, rmax, phimax):        rs, phis = np.meshgrid(np.linspace(0, rmax, n), np.linspace(0, phimax, m), sparse=True)        xs, ys = rphi2xy(rs, phis)    xs, ys = xs.reshape(-1), ys.reshape(-1)        coords = np.vstack((ys, xs))    #print(coords.shape)        vecC = np.array([f.shape[0]//2, f.shape[1]//2]).reshape(2,1)    coords += vecC        g = img.map_coordinates(f, coords, order=3)    g = g.reshape(m, n)    return np.flipud(g)def from_r_phi_plane_V2_(g, m, n, rmax, phimax):    xs, ys = np.meshgrid(np.arange(n), np.arange(m), sparse=True)        xs -= n//2    ys -= m//2        rs, phis = xy2rphi(xs, ys)    #print(rs, phis)    phis += np.pi    rs, phis = rs.reshape(-1), phis.reshape(-1)         iis = phis / phimax * (m-1)    jjs = rs / rmax * (n-1)    coords = np.vstack((iis, jjs))     #print(iis, jjs)        h = img.map_coordinates(g, coords, order=3)    h = h.reshape(m, n)    return np.fliplr(np.flipud(h))

In [22]:

fig, axs = plt.subplots(1, 4, figsize=(15,15))f = imageRead('Data/clock.jpg')axs[0].imshow(f / 255, cmap="gray"); axs[0].axis("off")f = np.flipud(f)m, n = f.shapermax = np.sqrt((m/2)**2 + (n/2)**2)phimax = 2 * np.pig = to_r_phi_plane_(f, m, n, rmax, phimax)axs[1].imshow(g / 255, cmap="gray"); axs[1].axis("off")blurred_g = img.gaussian_filter1d(g, sigma=6, axis=0, mode="wrap")axs[2].imshow(blurred_g / 255, cmap="gray"); axs[2].axis("off")#g = np.flipud(g)h = from_r_phi_plane_V2_(blurred_g, m, n, rmax, phimax)axs[3].imshow(h / 255, cmap="gray"); _ = axs[3].axis("off")

In [23]:

f = imageRead('Data/clock.jpg')f = np.flipud(f)m, n = f.shapermax = np.sqrt((m/2)**2 + (n/2)**2)phimax = 2 * np.pig = to_r_phi_plane_(f, m, n, rmax, phimax)fig, axs = plt.subplots(1, 4, figsize=(15,10))blurred_g1 = img.gaussian_filter1d(g, sigma=2, axis=0, mode="wrap")h1 = from_r_phi_plane_V2_(blurred_g1, m, n, rmax, phimax)_ = axs[0].imshow(h1 / 255, cmap="gray"); _ = axs[0].axis("off")blurred_g2 = img.gaussian_filter1d(g, sigma=6, axis=0, mode="wrap")h2 = from_r_phi_plane_V2_(blurred_g2, m, n, rmax, phimax)_ = axs[1].imshow(h2 / 255, cmap="gray"); _ = axs[1].axis("off")blurred_g3 = img.gaussian_filter1d(g, sigma=12, axis=0, mode="wrap")h3 = from_r_phi_plane_V2_(blurred_g3, m, n, rmax, phimax)_ = axs[2].imshow(h3 / 255, cmap="gray"); _ = axs[2].axis("off")blurred_g4 = img.gaussian_filter1d(g, sigma=24, axis=0, mode="wrap")h4 = from_r_phi_plane_V2_(blurred_g4, m, n, rmax, phimax)_ = axs[3].imshow(h4 / 255, cmap="gray"); _ = axs[3].axis("off")#To save the images#imageWrite(h2, save_dir + "Task6-4.png")

Bilinear Warping¶

In [24]:

def bilinear_warp(arrF, arrH, u_ul, u_ur, u_ll, u_lr):    M, N = arrH.shape    u, v = np.meshgrid(np.arange(N), np.arange(M))    matX = np.stack((u, v)).astype(float)        u1, v1 = u_ul    u2, v2 = u_ur    u3, v3 = u_ll    u4, v4 = u_lr        x1, y1 = 0, 0    x2, y2 = arrF.shape[1], 0    x3, y3 = 0, arrF.shape[0]    x4, y4 = arrF.shape[1], arrF.shape[0]        A = [[u1*v1, u1, v1, 1, 0, 0, 0, 0],         [u2*v2, u2, v2, 1, 0, 0, 0, 0],         [u3*v3, u3, v3, 1, 0, 0, 0, 0],         [u4*v4, u4, v4, 1, 0, 0, 0, 0],         [0, 0, 0, 0, u1*v1, u1, v1, 1],         [0, 0, 0, 0, u2*v2, u2, v2, 1],         [0, 0, 0, 0, u3*v3, u3, v3, 1],         [0, 0, 0, 0, u4*v4, u4, v4, 1]]    A = np.array(A)    b = np.array([x1, x2, x3, x4, y1, y2, y3, y4]).T    X, _, _, _ = np.linalg.lstsq(A, b, rcond=None)    a,b,c,d,e,f,g,h = X        X_ = a*matX[0]*matX[1] + b*matX[0] + c*matX[1] + d    Y_ = e*matX[0]*matX[1] + f*matX[0] + g*matX[1] + h        matX = np.stack((Y_, X_))        arrG = img.map_coordinates(arrF, matX, cval=-1) # matX MUST be float    arrG = arrG.reshape(matX.shape[1],matX.shape[2])        mask = arrG != -1    newArr = arrH.copy()    newArr[mask] = arrG[mask]    return newArrarrF = imageRead('Data/asterix.png')arrH = imageRead('Data/isle.jpg')out = bilinear_warp(arrF, arrH, (215, 56), (365, 10), (218, 258), (364, 296))plt.subplots(figsize=(10,10))plt.imshow(out / 255, cmap="gray"); plt.axis("off"); plt.show()

Perspective Mapping¶

In [25]:

def perspective_mapping(arrF, arrH, u_ul, u_ur, u_ll, u_lr, debug=False):    M, N = arrH.shape    u, v = np.meshgrid(np.arange(N), np.arange(M))    matX = np.stack((u, v)).astype(float)        u1, v1 = u_ul    u2, v2 = u_ur    u3, v3 = u_ll    u4, v4 = u_lr        x1, y1 = 0, 0    x2, y2 = arrF.shape[1], 0    x3, y3 = 0, arrF.shape[0]    x4, y4 = arrF.shape[1], arrF.shape[0]        A = [[u1, v1, 1, 0, 0, 0, -u1*x1, -v1*x1],         [u2, v2, 1, 0, 0, 0, -u2*x2, -v2*x2],         [u3, v3, 1, 0, 0, 0, -u3*x3, -v3*x3],         [u4, v4, 1, 0, 0, 0, -u4*x4, -v4*x4],         [0, 0, 0, u1, v1, 1, -u1*y1, -v1*y1],         [0, 0, 0, u2, v2, 1, -u2*y2, -v2*y2],         [0, 0, 0, u3, v3, 1, -u3*y3, -v3*y3],         [0, 0, 0, u4, v4, 1, -u4*y4, -v4*y4]]        A = np.array(A)    b = np.array([x1, x2, x3, x4, y1, y2, y3, y4]).T    X, _, _, _ = np.linalg.lstsq(A, b, rcond=None)    a,b,c,d,e,f,g,h = X        X_ = (a*matX[0] + b*matX[1] + c) / (g*matX[0] + h*matX[1] + 1)    Y_ = (d*matX[0] + e*matX[1] + f) / (g*matX[0] + h*matX[1] + 1)        matX = np.stack((Y_, X_))        arrG = img.map_coordinates(arrF, matX, cval=-1) # matX MUST be float    arrG = arrG.reshape(matX.shape[1],matX.shape[2])        mask = arrG != -1    #mask = np.bitwise_and(arrG != -1, arrG<230)        if debug:        plt.imshow(arrG, cmap="gray"); plt.title("Transformed Image"); plt.show()        plt.imshow(mask, cmap="gray"); plt.title("Mask")    newArr = arrH.copy()    newArr[mask] = arrG[mask]    return newArrdef perspective_mapping_transparent(arrF, arrH, u_ul, u_ur, u_ll, u_lr, debug=False):    M, N = arrH.shape    u, v = np.meshgrid(np.arange(N), np.arange(M))    matX = np.stack((u, v)).astype(float)        u1, v1 = u_ul    u2, v2 = u_ur    u3, v3 = u_ll    u4, v4 = u_lr        x1, y1 = 0, 0    x2, y2 = arrF.shape[1], 0    x3, y3 = 0, arrF.shape[0]    x4, y4 = arrF.shape[1], arrF.shape[0]        A = [[u1, v1, 1, 0, 0, 0, -u1*x1, -v1*x1],         [u2, v2, 1, 0, 0, 0, -u2*x2, -v2*x2],         [u3, v3, 1, 0, 0, 0, -u3*x3, -v3*x3],         [u4, v4, 1, 0, 0, 0, -u4*x4, -v4*x4],         [0, 0, 0, u1, v1, 1, -u1*y1, -v1*y1],         [0, 0, 0, u2, v2, 1, -u2*y2, -v2*y2],         [0, 0, 0, u3, v3, 1, -u3*y3, -v3*y3],         [0, 0, 0, u4, v4, 1, -u4*y4, -v4*y4]]        A = np.array(A)    b = np.array([x1, x2, x3, x4, y1, y2, y3, y4]).T    X, _, _, _ = np.linalg.lstsq(A, b, rcond=None)    a,b,c,d,e,f,g,h = X         X_ = (a*matX[0] + b*matX[1] + c) / (g*matX[0] + h*matX[1] + 1)    Y_ = (d*matX[0] + e*matX[1] + f) / (g*matX[0] + h*matX[1] + 1)        matX = np.stack((Y_, X_))        arrG = img.map_coordinates(arrF, matX, cval=-1) # matX MUST be float    arrG = arrG.reshape(matX.shape[1],matX.shape[2])        #mask = arrG != -1    mask = np.bitwise_and(arrG != -1, arrG<230)        if debug:        plt.imshow(arrG, cmap="gray"); plt.title("Transformed Image"); plt.show()        plt.imshow(mask, cmap="gray"); plt.title("Mask")    newArr = arrH.copy()    newArr[mask] = arrG[mask]    return newArrarrF_ = imageRead('Data/asterix.png')arrH_ = imageRead('Data/isle.jpg')out = perspective_mapping_transparent(arrF_, arrH_, (215, 56), (365, 10), (218, 258), (364, 296))out2 = perspective_mapping(arrF_, arrH_, (215, 56), (365, 10), (218, 258), (364, 296))fig, axs = plt.subplots(1, 2, figsize=(10,10))axs[0].imshow(out / 255, cmap="gray")axs[1].imshow(out2 / 255, cmap="gray")plt.show()

In [26]:

arrF_ = imageRead('Data/clock.jpg')arrH_ = imageRead('Data/isle.jpg')out = perspective_mapping_transparent(arrF_, arrH_, (215, 56), (365, 10), (218, 258), (364, 296))out2 = perspective_mapping(arrF_, arrH_, (215, 56), (365, 10), (218, 258), (364, 296))fig, axs = plt.subplots(1, 2, figsize=(10,10))axs[0].imshow(out / 255, cmap="gray")axs[1].imshow(out2 / 255, cmap="gray")plt.show()#imageWrite(out, save_dir + "Task6-5.png")