1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
|
# Copyright (c) 2016-2017 Weitian LI <weitian@aaronly.me>
# MIT license
"""
Image (only gray-scale image, i.e., matrix) transformation utilities.
References
----------
- Leptonica: Rotation
http://www.leptonica.com/rotation.html
- Image rotation by MATLAB without using imrotate
https://stackoverflow.com/a/19687481/4856091
https://stackoverflow.com/a/19689081/4856091
- Stackoverflow: Python: Rotating greyscale images
https://codereview.stackexchange.com/a/41903
"""
import numpy as np
import numba as nb
from scipy import ndimage
@nb.jit([nb.float64[:, :](nb.int64[:, :], nb.float64, nb.boolean,
nb.boolean, nb.float64),
nb.float64[:, :](nb.float64[:, :], nb.float64, nb.boolean,
nb.boolean, nb.float64)],
nopython=True)
def rotate_center(imgin, angle, interp=True, reshape=True, fill_value=0.0):
"""
Rotate the input image (only gray-scale image currently supported)
by a given angle about its center.
Parameters
----------
imgin : 2D `~numpy.ndarray`
Input gray-scale image to be rotated
angle : float
Rotation angle (unit: [ degree ])
interp : bool, optional
Use the area mapping of the 4 closest input pixels (``interp=True``),
which is also the same as "bilinear interpolation",
or use the nearest neighbor (``interp=False``) for rotated pixels.
reshape : bool, optional
Whether adapt the output shape so that the input image is contained
completely in the output?
fill_value : float, optional
Value used to fill pixels in the rotated image that corresponding
pixels outside the boundaries of the input image.
"""
nrow, ncol = imgin.shape
# Rotation transformation image
angle = np.deg2rad(angle)
mrotate = np.zeros((2, 2), dtype=np.float64)
mrotate[0, 0] = np.cos(angle)
mrotate[0, 1] = np.sin(angle)
mrotate[1, 0] = -np.sin(angle)
mrotate[1, 1] = np.cos(angle)
# Determine the shape of rotated image
corner00 = np.array((0, 0))
corner01 = np.array((0, ncol-1))
corner10 = np.array((nrow-1, 0))
corner11 = np.array((nrow-1, ncol-1))
corners = np.vstack((corner00, corner01, corner10, corner11))
if reshape:
dest = np.dot(corners.astype(np.float64), mrotate)
# XXX: ``numba`` does not support ``np.max()`` with arguments
minr = np.min(dest[:, 0])
minc = np.min(dest[:, 1])
maxr = np.max(dest[:, 0])
maxc = np.max(dest[:, 1])
nr = int(maxr - minr + 0.5)
nc = int(maxc - minc + 0.5)
else:
# Constraint to be same shape
nr = nrow
nc = ncol
imgout = np.ones((nr, nc)) * fill_value
#
# Calculate the offset, for easier transformation of rotated pixels
# back to input image.
#
# NOTE:
# Notations:
# P_out : (r_out, c_out) a pixel in the output rotated image
# C_out : center position of the output rotated image
# P_in : (r_in, c_in) a pixel in the input image
# C_in : center position of the input image
# R : rotation matrix
# R_T : transposed rotation matrix
# The rotation relation between pixel pair is (about the center):
# (P_in - C_in) * R = P_out - C_out
# Then:
# (P_in - C_in) = (P_out - C_out) * R_T
# And then:
# P_in = C_in + (P_out-C_out) * R_T = P_out*R_T + (C_in - C_out*R_T)
# Thus can define the "offset" as:
# offset = C_in - C_out * R_T
# Then the transformation back to input image is simply given by:
# P_in = P_out * R_T + offset
#
center_in = np.array((nrow/2.0 - 0.5, ncol/2.0 - 0.5))
center_out = np.array((nr/2.0 - 0.5, nc/2.0 - 0.5))
mrotate_T = mrotate.transpose()
offset = center_in - np.dot(center_out, mrotate_T)
# Map pixels of the rotated image to the input image
for rr in range(nr):
for cc in range(nc):
p_out = np.array((rr, cc))
p_in = np.dot(p_out.astype(np.float64), mrotate_T) + offset
if np.all((p_in > corner00) & (p_in < corner11)):
# Calculate the pixel value for the rotated image
if interp:
# Use area mapping of the 4 closest input pixels
idx_rf, idx_cf = np.floor(p_in).astype(np.int64)
idx_rc, idx_cc = np.ceil(p_in).astype(np.int64)
# NOTE:
# It is possible that ``p_in[0]`` and/or ``p_in[1]``
# are just integers, which cause ``idx_rf == idx_rc``
# and/or ``idx_cf == idx_cc``, which further lead to
# the calculated pixel value ``p_val = 0``.
if idx_rf == idx_rc:
idx_rc += 1
if idx_cf == idx_cc:
idx_cc += 1
# Calculate the overlapping areas
p_r, p_c = p_in
p4_area = np.array([(idx_rc - p_r) * (idx_cc - p_c),
(idx_rc - p_r) * (p_c - idx_cf),
(p_r - idx_rf) * (idx_cc - p_c),
(p_r - idx_rf) * (p_c - idx_cf)])
p4_val = np.array((imgin[idx_rf, idx_cf],
imgin[idx_rf, idx_cc],
imgin[idx_rc, idx_cf],
imgin[idx_rc, idx_cc]))
p_val = np.sum(p4_area * p4_val)
else:
# Use the nearest neighbor as the rotated value
idx_r = round(p_in[0])
idx_c = round(p_in[1])
p_val = imgin[idx_r, idx_c]
#
imgout[rr, cc] = p_val
return imgout
def circle2ellipse(imgcirc, bfraction, rotation=None):
"""
Shrink the input circle image with respect to the center along the
column (axis) to transform the circle to an ellipse, and then rotate
around the image center.
Parameters
----------
imgcirc : 2D `~numpy.ndarray`
Input image grid containing a circle at the center
bfraction : float
The fraction of the semi-minor axis w.r.t. the semi-major axis
(i.e., the half width of the input image), to determine the
shrunk size (height) of the output image.
Should be a fraction within [0, 1]
rotation : float, optional
Rotation angle (unit: [deg])
Default: ``None`` (i.e., no rotation)
Returns
-------
imgout : 2D `~numpy.ndarray`
Image of the same size as the input circle image.
"""
nrow, ncol = imgcirc.shape
# Shrink the circle to be elliptical
nrow2 = nrow * bfraction
nrow2 = int(nrow2 / 2) * 2 + 1 # be odd
# NOTE: zoom() calculate the output shape with round() instead of int();
# fix the warning about they may be different.
zoom = ((nrow2+0.1)/nrow, 1)
img2 = ndimage.zoom(imgcirc, zoom=zoom, order=1)
# Pad the shrunk image to have the same size as input
imgout = np.zeros(shape=(nrow, ncol))
r1 = int((nrow - nrow2) / 2)
imgout[r1:(r1+nrow2), :] = img2
if rotation:
imgout = ndimage.rotate(imgout, angle=rotation, reshape=False, order=1)
return imgout
|
