从紧耦合线和噪声曲线中寻找直线
Finding straight lines from tightly coupled lines and noise curvy lines
我有这张树线裁剪图片。我需要找到裁剪对齐的大体方向。我正在尝试获取图像的霍夫线，然后找到角度分布模式。
我一直在关注 this tutorial 作物线，但是在那一个上，作物线很稀疏。在这里，它们密集地堆积在一起，经过灰度化、模糊处理和使用精明的边缘检测后，这就是我得到的
import cv2
import numpy as np
import matplotlib.pyplot as plt
img = cv2.imread('drive/MyDrive/tree/sample.jpg')
gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY)
gauss = cv2.GaussianBlur(gray, (3,3), 3)
plt.figure(figsize=(15,15))
plt.subplot(1,2,1)
plt.imshow(gauss)
gscale = cv2.Canny(gauss, 80, 140)
plt.subplot(1,2,2)
plt.imshow(gscale)
plt.show()
（左边一张没有canny的模糊图片，左边一张经过canny预处理）
之后，我按照教程对预处理后的图像进行了“骨架化”
size = np.size(gscale)
skel = np.zeros(gscale.shape, np.uint8)
ret, gscale = cv2.threshold(gscale, 128, 255,0)
element = cv2.getStructuringElement(cv2.MORPH_CROSS, (3,3))
done = False
while not done:
eroded = cv2.erode(gscale, element)
temp = cv2.dilate(eroded, element)
temp = cv2.subtract(gscale, temp)
skel = cv2.bitwise_or(skel, temp)
gscale = eroded.copy()
zeros = size - cv2.countNonZero(gscale)
if zeros==size:
done = True
给我
如你所见，仍然有一堆弯曲的线条。在上面使用HoughLines算法时，到处散落着11k条线
lines = cv2.HoughLinesP(skel,1,np.pi/180,130)
a,b,c = lines.shape
for i in range(a):
rho = lines[i][0][0]
theta = lines[i][0][1]
a = np.cos(theta)
b = np.sin(theta)
x0 = a*rho
y0 = b*rho
x1 = int(x0 + 1000*(-b))
y1 = int(y0 + 1000*(a))
x2 = int(x0 - 1000*(-b))
y2 = int(y0 - 1000*(a))
cv2.line(img,(x1,y1),(x2,y2),(0,0,255),2, cv2.LINE_AA)#showing the results:
plt.figure(figsize=(15,15))
plt.subplot(121)#OpenCV reads images as BGR, this corrects so it is displayed as RGB
plt.plot()
plt.imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
plt.title('Row Detection')
plt.xticks([])
plt.yticks([])
plt.subplot(122)
plt.plot()
plt.imshow(skel,cmap='gray')
plt.title('Skeletal Image')
plt.xticks([])
plt.yticks([])
plt.show()
我是 cv2 的新手，所以我不知道该怎么做。搜索并尝试了很多东西，但 none 有效。我怎样才能去除稍大的点，去除波浪线？
您可以使用 2D FFT 来找到裁剪对齐的大致方向（正如 mozway 在评论中所建议的那样）。这个想法是，当输入包含许多相同方向的线时，可以很容易地从出现在幅度谱 中的 中心光束射线中提取大体方向。您可以在 this previous post 中找到有关其工作原理的更多信息。它直接与输入图像一起工作，但最好应用高斯 + Canny 过滤器。
这里是过滤后的灰度图像的幅度谱中有趣的部分：
可以很容易地看到主光束。您可以通过以递增的角度迭代多条线来提取它的角度，并对每条线上的幅度值求和，如下图所示：
这里是根据线的角度（以弧度为单位）绘制的每条线的幅度总和：
在此基础上，你只需要找到使计算和最大化的角度。
这是结果代码：
def computeAngle(arr):
# Naive inefficient algorithm
n, m = arr.shape
yCenter, xCenter = (n-1, m//2-1)
lineLen = m//2-2
sMax = 0.0
bestAngle = np.nan
for angle in np.arange(0, math.pi, math.pi/300):
i = np.arange(lineLen)
y, x = (np.sin(angle) * i + 0.5).astype(np.int_), (np.cos(angle) * i + 0.5).astype(np.int_)
s = np.sum(arr[yCenter-y, xCenter+x])
if s > sMax:
bestAngle = angle
sMax = s
return bestAngle
# Load the image in gray
img = cv2.imread('lines.jpg')
gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY)
# Apply some filters
gauss = cv2.GaussianBlur(gray, (3,3), 3)
gscale = cv2.Canny(gauss, 80, 140)
# Compute the 2D FFT of real values
freqs = np.fft.rfft2(gscale)
# Shift the frequencies (centering) and select the low frequencies
upperPart = freqs[:freqs.shape[0]//4,:freqs.shape[1]//2]
lowerPart = freqs[-freqs.shape[0]//4:,:freqs.shape[1]//2]
filteredFreqs = np.vstack((lowerPart, upperPart))
# Compute the magnitude spectrum
magnitude = np.log(np.abs(filteredFreqs))
# Correct the angle
magnitude = np.rot90(magnitude).copy()
# Find the major angle
bestAngle = computeAngle(magnitude)
为了完整起见，我也想 post Sobel 梯度角方式。
总体思路：
- 对于每个像素，计算 X 和 Y 梯度（例如使用 Sobel）
- 根据 X 和 Y 梯度信息及其分布计算角度
- 找到模式，例如带有直方图和 select 最高的
用 C++ 编写，但可能很容易转换为 python:
int main()
{
try
{
cv::Mat img = cv::imread("C:/data/Whosebug/cropLines/lines.jpg", cv::IMREAD_GRAYSCALE);
// tests with artificial lines:
//img = cv::Mat::zeros(img.size(), CV_8UC1);
//img = cv::Mat(img.size(), CV_8UC1, cv::Scalar::all(255));
//cv::line(img, cv::Point(0, img.rows), cv::Point(img.cols, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols, img.rows), cv::Point(0, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols, img.rows/2), cv::Point(0, img.rows/2), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols/2, img.rows), cv::Point(img.cols/2, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols / 2, img.rows), cv::Point(img.cols / 2, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols / 2, img.rows), cv::Point(img.cols / 2, 0), cv::Scalar::all(0), 10); // sample to check angles
cv::imshow("img", img);
cv::Mat gradX, gradY, mag, angle;
cv::Sobel(img, gradX, CV_32F, 1, 0, 3);
cv::Sobel(img, gradY, CV_32F, 0, 1, 3);
cv::cartToPolar(gradX, gradY, mag, angle, true);
cv::Mat magMask = mag > 0; // dont use pixels where angle is 0 just because there is no gradient.
float scaleX = 3;
float scaleY = 0.15;
float maxValueY = 3000;
cv::Mat chart = cv::Mat(maxValueY * scaleY, 360 * scaleX, CV_8UC3);
cv::Point prev(-1, -1);
double window = 5.0; // window size 1 is much more noisy but still works
// this loop can probably be optimized with an optimized histogram compuation from any library
for (int i = 0; i <= 360; i = i + 1)
{
double amount = cv::countNonZero((angle >= i-window/2 & angle < i + window/2) & (magMask));
std::cout << i << "°: " << amount << std::endl;
cv::Point current(i*scaleX, chart.rows - amount*scaleY/window);
if (prev.x >= 0) cv::line(chart, prev, current, cv::Scalar(0, 0, 255), 1);
prev = current;
}
cv::line(chart, cv::Point(45 * scaleX, 0), cv::Point(45 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(90 * scaleX, 0), cv::Point(90 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(135 * scaleX, 0), cv::Point(135 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(180 * scaleX, 0), cv::Point(180 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(225 * scaleX, 0), cv::Point(225 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(270 * scaleX, 0), cv::Point(270 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(315 * scaleX, 0), cv::Point(315 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(360 * scaleX, 0), cv::Point(360 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::imshow("chart", chart);
cv::imwrite("C:/data/Whosebug/cropLines/chart.png", chart);
cv::imwrite("C:/data/Whosebug/cropLines/input.png", img);
cv::waitKey(0);
}
catch (std::exception& e)
{
std::cout << e.what() << std::endl;
}
}
给出这个结果，图像顶部的每条蓝线都是 45°。最大值为 52°（及其 180° 的倍数），其中梯度相对于直线方向旋转 90°。所以该图像中的线方向是从 x 轴顺时针旋转 142° 或逆时针旋转 38°。
我有这张树线裁剪图片。我需要找到裁剪对齐的大体方向。我正在尝试获取图像的霍夫线，然后找到角度分布模式。
我一直在关注 this tutorial 作物线，但是在那一个上，作物线很稀疏。在这里，它们密集地堆积在一起，经过灰度化、模糊处理和使用精明的边缘检测后，这就是我得到的
import cv2
import numpy as np
import matplotlib.pyplot as plt
img = cv2.imread('drive/MyDrive/tree/sample.jpg')
gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY)
gauss = cv2.GaussianBlur(gray, (3,3), 3)
plt.figure(figsize=(15,15))
plt.subplot(1,2,1)
plt.imshow(gauss)
gscale = cv2.Canny(gauss, 80, 140)
plt.subplot(1,2,2)
plt.imshow(gscale)
plt.show()
（左边一张没有canny的模糊图片，左边一张经过canny预处理）
之后，我按照教程对预处理后的图像进行了“骨架化”
size = np.size(gscale)
skel = np.zeros(gscale.shape, np.uint8)
ret, gscale = cv2.threshold(gscale, 128, 255,0)
element = cv2.getStructuringElement(cv2.MORPH_CROSS, (3,3))
done = False
while not done:
eroded = cv2.erode(gscale, element)
temp = cv2.dilate(eroded, element)
temp = cv2.subtract(gscale, temp)
skel = cv2.bitwise_or(skel, temp)
gscale = eroded.copy()
zeros = size - cv2.countNonZero(gscale)
if zeros==size:
done = True
给我
如你所见，仍然有一堆弯曲的线条。在上面使用HoughLines算法时，到处散落着11k条线
lines = cv2.HoughLinesP(skel,1,np.pi/180,130)
a,b,c = lines.shape
for i in range(a):
rho = lines[i][0][0]
theta = lines[i][0][1]
a = np.cos(theta)
b = np.sin(theta)
x0 = a*rho
y0 = b*rho
x1 = int(x0 + 1000*(-b))
y1 = int(y0 + 1000*(a))
x2 = int(x0 - 1000*(-b))
y2 = int(y0 - 1000*(a))
cv2.line(img,(x1,y1),(x2,y2),(0,0,255),2, cv2.LINE_AA)#showing the results:
plt.figure(figsize=(15,15))
plt.subplot(121)#OpenCV reads images as BGR, this corrects so it is displayed as RGB
plt.plot()
plt.imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
plt.title('Row Detection')
plt.xticks([])
plt.yticks([])
plt.subplot(122)
plt.plot()
plt.imshow(skel,cmap='gray')
plt.title('Skeletal Image')
plt.xticks([])
plt.yticks([])
plt.show()
我是 cv2 的新手，所以我不知道该怎么做。搜索并尝试了很多东西，但 none 有效。我怎样才能去除稍大的点，去除波浪线？
您可以使用 2D FFT 来找到裁剪对齐的大致方向（正如 mozway 在评论中所建议的那样）。这个想法是，当输入包含许多相同方向的线时，可以很容易地从出现在幅度谱 中的 中心光束射线中提取大体方向。您可以在 this previous post 中找到有关其工作原理的更多信息。它直接与输入图像一起工作，但最好应用高斯 + Canny 过滤器。
这里是过滤后的灰度图像的幅度谱中有趣的部分：
可以很容易地看到主光束。您可以通过以递增的角度迭代多条线来提取它的角度，并对每条线上的幅度值求和，如下图所示：
这里是根据线的角度（以弧度为单位）绘制的每条线的幅度总和：
在此基础上，你只需要找到使计算和最大化的角度。
这是结果代码：
def computeAngle(arr):
# Naive inefficient algorithm
n, m = arr.shape
yCenter, xCenter = (n-1, m//2-1)
lineLen = m//2-2
sMax = 0.0
bestAngle = np.nan
for angle in np.arange(0, math.pi, math.pi/300):
i = np.arange(lineLen)
y, x = (np.sin(angle) * i + 0.5).astype(np.int_), (np.cos(angle) * i + 0.5).astype(np.int_)
s = np.sum(arr[yCenter-y, xCenter+x])
if s > sMax:
bestAngle = angle
sMax = s
return bestAngle
# Load the image in gray
img = cv2.imread('lines.jpg')
gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY)
# Apply some filters
gauss = cv2.GaussianBlur(gray, (3,3), 3)
gscale = cv2.Canny(gauss, 80, 140)
# Compute the 2D FFT of real values
freqs = np.fft.rfft2(gscale)
# Shift the frequencies (centering) and select the low frequencies
upperPart = freqs[:freqs.shape[0]//4,:freqs.shape[1]//2]
lowerPart = freqs[-freqs.shape[0]//4:,:freqs.shape[1]//2]
filteredFreqs = np.vstack((lowerPart, upperPart))
# Compute the magnitude spectrum
magnitude = np.log(np.abs(filteredFreqs))
# Correct the angle
magnitude = np.rot90(magnitude).copy()
# Find the major angle
bestAngle = computeAngle(magnitude)
为了完整起见，我也想 post Sobel 梯度角方式。
总体思路：
- 对于每个像素，计算 X 和 Y 梯度（例如使用 Sobel）
- 根据 X 和 Y 梯度信息及其分布计算角度
- 找到模式，例如带有直方图和 select 最高的
用 C++ 编写，但可能很容易转换为 python:
int main()
{
try
{
cv::Mat img = cv::imread("C:/data/Whosebug/cropLines/lines.jpg", cv::IMREAD_GRAYSCALE);
// tests with artificial lines:
//img = cv::Mat::zeros(img.size(), CV_8UC1);
//img = cv::Mat(img.size(), CV_8UC1, cv::Scalar::all(255));
//cv::line(img, cv::Point(0, img.rows), cv::Point(img.cols, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols, img.rows), cv::Point(0, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols, img.rows/2), cv::Point(0, img.rows/2), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols/2, img.rows), cv::Point(img.cols/2, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols / 2, img.rows), cv::Point(img.cols / 2, 0), cv::Scalar::all(255), 10); // sample to check angles
//cv::line(img, cv::Point(img.cols / 2, img.rows), cv::Point(img.cols / 2, 0), cv::Scalar::all(0), 10); // sample to check angles
cv::imshow("img", img);
cv::Mat gradX, gradY, mag, angle;
cv::Sobel(img, gradX, CV_32F, 1, 0, 3);
cv::Sobel(img, gradY, CV_32F, 0, 1, 3);
cv::cartToPolar(gradX, gradY, mag, angle, true);
cv::Mat magMask = mag > 0; // dont use pixels where angle is 0 just because there is no gradient.
float scaleX = 3;
float scaleY = 0.15;
float maxValueY = 3000;
cv::Mat chart = cv::Mat(maxValueY * scaleY, 360 * scaleX, CV_8UC3);
cv::Point prev(-1, -1);
double window = 5.0; // window size 1 is much more noisy but still works
// this loop can probably be optimized with an optimized histogram compuation from any library
for (int i = 0; i <= 360; i = i + 1)
{
double amount = cv::countNonZero((angle >= i-window/2 & angle < i + window/2) & (magMask));
std::cout << i << "°: " << amount << std::endl;
cv::Point current(i*scaleX, chart.rows - amount*scaleY/window);
if (prev.x >= 0) cv::line(chart, prev, current, cv::Scalar(0, 0, 255), 1);
prev = current;
}
cv::line(chart, cv::Point(45 * scaleX, 0), cv::Point(45 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(90 * scaleX, 0), cv::Point(90 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(135 * scaleX, 0), cv::Point(135 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(180 * scaleX, 0), cv::Point(180 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(225 * scaleX, 0), cv::Point(225 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(270 * scaleX, 0), cv::Point(270 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(315 * scaleX, 0), cv::Point(315 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::line(chart, cv::Point(360 * scaleX, 0), cv::Point(360 * scaleX, 100 * scaleY), cv::Scalar(255, 0, 0), 1);
cv::imshow("chart", chart);
cv::imwrite("C:/data/Whosebug/cropLines/chart.png", chart);
cv::imwrite("C:/data/Whosebug/cropLines/input.png", img);
cv::waitKey(0);
}
catch (std::exception& e)
{
std::cout << e.what() << std::endl;
}
}
给出这个结果，图像顶部的每条蓝线都是 45°。最大值为 52°（及其 180° 的倍数），其中梯度相对于直线方向旋转 90°。所以该图像中的线方向是从 x 轴顺时针旋转 142° 或逆时针旋转 38°。