Draw an arc by using end-points and bulge distance. In OpenCV or PIL
Working on a script to convert dxf to png, I need to draw arc which are having only three parameters, i.e. Start point of arc, end point of arc and bulge distance.
I have checked OpenCV and PIL both, and they require start and end angle to draw this arc. I can find out those angles using some geometry, but would like to know if there is any other solution out there which I am missing out.
You have three pieces of information defining your circular arc: two points on the circle (defining a chord of that circle) and the bulge distance (called the sagitta of the circular arc).
See the following graphic:
Here s is the sagitta, l is half the chord-length, and r is of course the radius. Other important non-marked positions are the points which the chord intersects the circle, the point where the sagitta intersects the circle, and the center of the circle, from which the radius is extending.
For OpenCV's ellipse() function, we would use the following prototype:
cv2.ellipse(img, center, axes, angle, startAngle, endAngle, color[, thickness[, lineType[, shift]]]) → img
where most of the parameters are described by the following graphic:
Since we're drawing a circular, not elliptical arc, the major/minor axes will have the same size and there's no difference rotating it, so the axes will just be (radius, radius) and the angle should be zero to simplify. Then the only parameters we need are the center of the circle, the radius, and the start angle and ending angle of drawing, corresponding to the points of the chord. The angles are easy to calculate (they're just some angles on the circle). So ultimately we need to find the radius and center of the circle.
Finding the radius and the center is the same as finding the equation of the circle, so there are a ton of ways to do it. But since we're programming here, the easiest way IMO is to define a third point on the circle by where the sagitta touches the circle, and then solve for the circle from those three points.
So first we'll need to get the midpoint of the chord, get a perpendicular line to that midpoint, and extend it the length of the sagitta to get to that third point, but that's easy enough. I'll start given pt1 = (x1, y1) and pt2 = (x2, y2) as my two points on the circle and sagitta is the 'bulge depth' (i.e. the parameters you have):
# extract point coordinates
x1, y1 = pt1
x2, y2 = pt2
# find normal from midpoint, follow by length sagitta
n = np.array([y2 - y1, x1 - x2])
n_dist = np.sqrt(np.sum(n**2))
if np.isclose(n_dist, 0):
# catch error here, d(pt1, pt2) ~ 0
print('Error: The distance between pt1 and pt2 is too small.')
n = n/n_dist
x3, y3 = (np.array(pt1) + np.array(pt2))/2 + sagitta * n
Now we've got the third point on the circle. Note that the sagitta is just some length, so it could go either direction---if sagitta were negative, it'd go one direction from the chord, and if it were positive, it'd go the other direction. Not sure if this is how that distance is given to you or not.
Then we can simply use determinants to solve for the radius and center.
# calculate the circle from three points
# see https://math.stackexchange.com/a/1460096/246399
A = np.array([
[x1**2 + y1**2, x1, y1, 1],
[x2**2 + y2**2, x2, y2, 1],
[x3**2 + y3**2, x3, y3, 1]])
M11 = np.linalg.det(A[:, (1, 2, 3)])
M12 = np.linalg.det(A[:, (0, 2, 3)])
M13 = np.linalg.det(A[:, (0, 1, 3)])
M14 = np.linalg.det(A[:, (0, 1, 2)])
if np.isclose(M11, 0):
# catch error here, the points are collinear (sagitta ~ 0)
print('Error: The third point is collinear.')
cx = 0.5 * M12/M11
cy = -0.5 * M13/M11
radius = np.sqrt(cx**2 + cy**2 + M14/M11)
Then lastly, as we need the starting and ending angles to draw the ellipse with OpenCV, we can just use atan2() to get the angles from the center to the initial points:
# calculate angles of pt1 and pt2 from center of circle
pt1_angle = 180*np.arctan2(y1 - cy, x1 - cx)/np.pi
pt2_angle = 180*np.arctan2(y2 - cy, x2 - cx)/np.pi
So I packaged this all up into one function:
def convert_arc(pt1, pt2, sagitta):
# extract point coordinates
x1, y1 = pt1
x2, y2 = pt2
# find normal from midpoint, follow by length sagitta
n = np.array([y2 - y1, x1 - x2])
n_dist = np.sqrt(np.sum(n**2))
if np.isclose(n_dist, 0):
# catch error here, d(pt1, pt2) ~ 0
print('Error: The distance between pt1 and pt2 is too small.')
n = n/n_dist
x3, y3 = (np.array(pt1) + np.array(pt2))/2 + sagitta * n
# calculate the circle from three points
# see https://math.stackexchange.com/a/1460096/246399
A = np.array([
[x1**2 + y1**2, x1, y1, 1],
[x2**2 + y2**2, x2, y2, 1],
[x3**2 + y3**2, x3, y3, 1]])
M11 = np.linalg.det(A[:, (1, 2, 3)])
M12 = np.linalg.det(A[:, (0, 2, 3)])
M13 = np.linalg.det(A[:, (0, 1, 3)])
M14 = np.linalg.det(A[:, (0, 1, 2)])
if np.isclose(M11, 0):
# catch error here, the points are collinear (sagitta ~ 0)
print('Error: The third point is collinear.')
cx = 0.5 * M12/M11
cy = -0.5 * M13/M11
radius = np.sqrt(cx**2 + cy**2 + M14/M11)
# calculate angles of pt1 and pt2 from center of circle
pt1_angle = 180*np.arctan2(y1 - cy, x1 - cx)/np.pi
pt2_angle = 180*np.arctan2(y2 - cy, x2 - cx)/np.pi
return (cx, cy), radius, pt1_angle, pt2_angle
With these values you can then the arc with OpenCV's ellipse() function. However, these are all floating point values. ellipse() does let you plot floating point values with the shift argument, but if you're not familiar with it it's a little weird, so instead we can borrow the solution from this answer to define a function
def draw_ellipse(
img, center, axes, angle,
startAngle, endAngle, color,
thickness=1, lineType=cv2.LINE_AA, shift=10):
# uses the shift to accurately get sub-pixel resolution for arc
# taken from https://stackoverflow.com/a/44892317/5087436
center = (
int(round(center[0] * 2**shift)),
int(round(center[1] * 2**shift))
)
axes = (
int(round(axes[0] * 2**shift)),
int(round(axes[1] * 2**shift))
)
return cv2.ellipse(
img, center, axes, angle,
startAngle, endAngle, color,
thickness, lineType, shift)
Then to use these functions it's as simple as:
img = np.zeros((500, 500), dtype=np.uint8)
pt1 = (50, 50)
pt2 = (350, 250)
sagitta = 50
center, radius, start_angle, end_angle = convert_arc(pt1, pt2, sagitta)
axes = (radius, radius)
draw_ellipse(img, center, axes, 0, start_angle, end_angle, 255)
cv2.imshow('', img)
cv2.waitKey()
And again note that a negative sagitta gives an arc the other direction:
center, radius, start_angle, end_angle = convert_arc(pt1, pt2, sagitta)
axes = (radius, radius)
draw_ellipse(img, center, axes, 0, start_angle, end_angle, 255)
center, radius, start_angle, end_angle = convert_arc(pt1, pt2, -sagitta)
axes = (radius, radius)
draw_ellipse(img, center, axes, 0, start_angle, end_angle, 127)
cv2.imshow('', img)
cv2.waitKey()
Lastly just to expand, I broke out two error cases in the convert_arc() function. First:
if np.isclose(n_dist, 0):
# catch error here, d(pt1, pt2) ~ 0
print('Error: The distance between pt1 and pt2 is too small.')
The error catch here is because we need to get a unit vector, so we need to divide by the length which can't be zero. Of course, this will only happen if pt1 and pt2 are the same point, so you can just check that they're unique at the top of the function instead of checking here.
Second:
if np.isclose(M11, 0):
# catch error here, the points are collinear (sagitta ~ 0)
print('Error: The third point is collinear.')
Here only happens if the three points are collinear, which only happens if the sagitta is 0. So again, you can check this at the top of your function (and maybe say, OK, if it is 0, then just draw a line from pt1 to pt2 or whatever you want to do).