Python基于随机采样一至性实现拟合椭圆(优化版)

上次版本如果在没有找到轮廓或轮廓的点集数很小无法拟合椭圆或在RANSAC中寻找最优解时会死循环中,优化后的代码

import cv2
import os
import numpy as np
import matplotlib.pyplot as plt
import math
from Ransac_Process import RANSAC

def cul_area(x_mask, y_mask, r_circle, mask):
    mask_label = mask.copy()
    num_area = 0
    for xm in range(x_mask+r_circle-10, x_mask+r_circle+10):
        for ym in range(y_mask+r_circle-10, y_mask+r_circle+10):
            # print(mask[ym, xm])
            if (pow((xm-x_mask), 2) + pow((ym-y_mask), 2) - pow(r_circle,  2)) == 0 and mask[ym, xm][0] == 255:
                num_area += 1
                mask_label[ym, xm] = (0, 0, 255)
    cv2.imwrite('./test2/mask_label.png', mask_label)
    print(num_area)
    return num_area

def mainFigure(img, point0):

    point_center = []
    # cv2.imwrite('./test2/img_source.png', img)
    img_hsv = cv2.cvtColor(img, cv2.COLOR_BGR2HSV)
    # cv2.imwrite('./test2/img_hsv.png', img_hsv)
    w, h = img.shape[1], img.shape[0]
    w_hsv, h_hsv = img_hsv.shape[1], img_hsv.shape[0]
    for i_hsv in range(w_hsv):
        for j_hsv in range(h_hsv):
            if img_hsv[j_hsv, i_hsv][0] < 200 and img_hsv[j_hsv, i_hsv][1] < 130 and img_hsv[j_hsv, i_hsv][2] > 120:
                # if hsv[j_hsv, i_hsv][0] < 100 and hsv[j_hsv, i_hsv][1] < 200 and hsv[j_hsv, i_hsv][2] > 80:
                img_hsv[j_hsv, i_hsv] = 255, 255, 255
            else:
                img_hsv[j_hsv, i_hsv] = 0, 0, 0
    # cv2.imwrite('./test2/img_hsvhb.png', img_hsv)
    # cv2.imshow("hsv", img_hsv)
    # cv2.waitKey()

    # 灰度化处理图像
    grayImage = cv2.cvtColor(img_hsv, cv2.COLOR_BGR2GRAY)
    # mask = np.zeros((grayImage.shape[0], grayImage.shape[1]), np.uint8)
    # mask = cv2.cvtColor(mask, cv2.COLOR_GRAY2BGR)
    # cv2.imwrite('./mask.png', mask)

    # 尝试寻找轮廓
    contours, hierarchy = cv2.findContours(grayImage, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_NONE)

    # 合并轮廓
    if len(contours) > 1:
        contours_merge = np.vstack([contours[0], contours[1]])
        for i in range(2, len(contours)):
            contours_merge = np.vstack([contours_merge, contours[i]])
        # cv2.drawContours(img, contours_merge, -1, (0, 255, 255), 1)
        # # cv2.imwrite('./test2/img_res.png', img)
        # cv2.imshow("contours_merge", img)
        # cv2.waitKey()
    elif len(contours) == 1:
        contours_merge = contours[0]
    else:
        print("No contours!")
        return 0,0

    # RANSAC拟合
    points_data = np.reshape(contours_merge, (-1, 2))  # ellipse edge points set

    # print("points_data", len(points_data))
    # 2.Ransac fit ellipse param
    # Ransac = RANSAC(data=points_data, threshold=0.1, P=.99, S=.5, N=10)
    Ransac = RANSAC(data=points_data, threshold=0.5, P=.98, S=.6, N=10)

    ellipse_values = Ransac.execute_ransac()
    # 检测到轮廓里数量太少(<5)则无法拟合椭圆
    if ellipse_values is None:
        return 0,0

    (X, Y), (LAxis, SAxis), Angle = ellipse_values
    # print( (X, Y), (LAxis, SAxis))
    # 拟合圆
    cv2.ellipse(img, ((X, Y), (LAxis, SAxis), Angle), (0, 0, 255), 1, cv2.LINE_AA)  # 画圆
    cv2.circle(img, (int(X), int(Y)), 3, (0, 0, 255), -1)  # 画圆心
    point_center.append(int(X))
    point_center.append(int(Y))

    # 直接拟合
    # rrt = cv2.fitEllipse(contours_merge)  # x, y)代表椭圆中心点的位置(a, b)代表长短轴长度,应注意a、b为长短轴的直径,而非半径,angle 代表了中心旋转的角度
    # # print("rrt", rrt)
    # cv2.ellipse(img, rrt, (255, 0, 0), 1, cv2.LINE_AA)  # 画圆
    # x, y = rrt[0]
    # cv2.circle(img, (int(x), int(y)), 3, (255, 0, 0), -1)  # 画圆心
    # point_center.append(int(x))
    # point_center.append(int(y))
    # # print("no",(x,y))
    #
    # # 两种方法坐标的距离
    # dis_two_method = math.sqrt(math.pow(X - x, 2) + math.pow(Y - y, 2))
    # print("两种方法坐标的距离", dis_two_method)

    cv2.imshow("fit circle", img)
    cv2.waitKey(3)
    # cv2.imwrite("./test2/fitcircle.png", img)

    return point_center[0], point_center[1]

if __name__ == "__main__":

    # 测试所有图片
    mainFolder = "./Images/save_img"
    myFolders = os.listdir(mainFolder)
    print("myFolders", myFolders)
    myImageList = []
    path = ''
    for folder in myFolders:
        path = mainFolder + '/' + folder

        myImageList = os.listdir(path)
        # print(myImageList)
        # print(f'Tatal images deteted is  {len(myImageList)}')

        i = 0
        for imagN in myImageList:
            curImg = cv2.imread(f'{path}/{imagN}')
            # images.append(curImg)
            print(f'{path}/{imagN}')
            point0 = [0, 0]
            cir_x, cir_y = mainFigure(curImg, point0)
            print("This is ", i, "圆心为",(cir_x, cir_y))
            i += 1

    # # 测试2
    # for i in range(1,6):
    #     imageName = "s"
    #     imageName += str(i)
    #     path = './Images/danHoles/' + imageName + '.png'
    #     print(path)
    #     img = cv2.imread(path)
    #     point0 = [0, 0]
    #     cir_x, cir_y = mainFigure(img, point0)

    # # 测试1
    # img = cv2.imread('./Images/danHoles/s6.png')
    # point0 = [0, 0]
    # cir_x, cir_y = mainFigure(img, point0)

Ransac_Process.py

import cv2
import math
import random
import numpy as np
from numpy.linalg import inv, svd, det
import time

class RANSAC:
    def __init__(self, data, threshold, P, S, N):
        self.point_data = data  # 椭圆轮廓点集
        self.length = len(self.point_data)  # 椭圆轮廓点集长度
        self.error_threshold = threshold  # 模型评估误差容忍阀值

        self.N = N  # 随机采样数
        self.S = S  # 设定的内点比例
        self.P = P  # 采得N点去计算的正确模型概率
        self.max_inliers = self.length * self.S  # 设定最大内点阀值
        self.items = 8

        self.count = 0  # 内点计数器
        self.best_model = ((0, 0), (1e-6, 1e-6), 0)  # 椭圆模型存储器

    def random_sampling(self, n):
        # 这个部分有修改的空间,这样循环次数太多了,可以看看别人改进的ransac拟合椭圆的论文
        """随机取n个数据点"""
        all_point = self.point_data

        if len(all_point) >= n:
            select_point = np.asarray(random.sample(list(all_point), n))
            return select_point
        else:
            print("轮廓点数太少,数量为", len(all_point))
            return None

    def Geometric2Conic(self, ellipse):
        # 这个部分参考了GitHub中的一位大佬的,但是时间太久,忘记哪个人的了
        """计算椭圆方程系数"""
        # Ax ^ 2 + Bxy + Cy ^ 2 + Dx + Ey + F
        (x0, y0), (bb, aa), phi_b_deg = ellipse

        a, b = aa / 2, bb / 2  # Semimajor and semiminor axes
        phi_b_rad = phi_b_deg * np.pi / 180.0  # Convert phi_b from deg to rad
        ax, ay = -np.sin(phi_b_rad), np.cos(phi_b_rad)  # Major axis unit vector

        # Useful intermediates
        a2 = a * a
        b2 = b * b

        # Conic parameters
        if a2 > 0 and b2 > 0:
            A = ax * ax / a2 + ay * ay / b2
            B = 2 * ax * ay / a2 - 2 * ax * ay / b2
            C = ay * ay / a2 + ax * ax / b2
            D = (-2 * ax * ay * y0 - 2 * ax * ax * x0) / a2 + (2 * ax * ay * y0 - 2 * ay * ay * x0) / b2
            E = (-2 * ax * ay * x0 - 2 * ay * ay * y0) / a2 + (2 * ax * ay * x0 - 2 * ax * ax * y0) / b2
            F = (2 * ax * ay * x0 * y0 + ax * ax * x0 * x0 + ay * ay * y0 * y0) / a2 + \
                (-2 * ax * ay * x0 * y0 + ay * ay * x0 * x0 + ax * ax * y0 * y0) / b2 - 1
        else:
            # Tiny dummy circle - response to a2 or b2 == 0 overflow warnings
            A, B, C, D, E, F = (1, 0, 1, 0, 0, -1e-6)

        # Compose conic parameter array
        conic = np.array((A, B, C, D, E, F))
        return conic

    def eval_model(self, ellipse):
        # 这个地方也有很大修改空间,判断是否内点的条件在很多改进的ransac论文中有说明,可以多看点论文
        """评估椭圆模型,统计内点个数"""
        # this an ellipse ?
        a, b, c, d, e, f = self.Geometric2Conic(ellipse)
        E = 4 * a * c - b * b
        if E <= 0:
            # print('this is not an ellipse')
            return 0, 0

        #  which long axis ?
        (x, y), (LAxis, SAxis), Angle = ellipse
        LAxis, SAxis = LAxis / 2, SAxis / 2
        if SAxis > LAxis:
            temp = SAxis
            SAxis = LAxis
            LAxis = temp

        # calculate focus
        Axis = math.sqrt(LAxis * LAxis - SAxis * SAxis)
        f1_x = x - Axis * math.cos(Angle * math.pi / 180)
        f1_y = y - Axis * math.sin(Angle * math.pi / 180)
        f2_x = x + Axis * math.cos(Angle * math.pi / 180)
        f2_y = y + Axis * math.sin(Angle * math.pi / 180)

        # identify inliers points
        f1, f2 = np.array([f1_x, f1_y]), np.array([f2_x, f2_y])
        f1_distance = np.square(self.point_data - f1)
        f2_distance = np.square(self.point_data - f2)
        all_distance = np.sqrt(f1_distance[:, 0] + f1_distance[:, 1]) + np.sqrt(f2_distance[:, 0] + f2_distance[:, 1])

        Z = np.abs(2 * LAxis - all_distance)
        delta = math.sqrt(np.sum((Z - np.mean(Z)) ** 2) / len(Z))

        # Update inliers set
        inliers = np.nonzero(Z < 0.8 * delta)[0]
        inlier_pnts = self.point_data[inliers]

        return len(inlier_pnts), inlier_pnts

    def execute_ransac(self):
        Time_start = time.time()
        while math.ceil(self.items):
            # print(self.max_inliers)

            # 1.select N points at random
            select_points = self.random_sampling(self.N)
            # 当从轮廓中采集的点不够拟合椭圆时跳出循环
            if select_points is None or len(select_points) < 5:
                # print(select_points)
                return None
            else:
                # 2.fitting N ellipse points
                ellipse = cv2.fitEllipse(select_points)

                # 3.assess model and calculate inliers points
                inliers_count, inliers_set = self.eval_model(ellipse)

                # 4.number of new inliers points more than number of old inliers points ?
                if inliers_count > self.count:
                    if len(inliers_set) > 4:
                        ellipse_ = cv2.fitEllipse(inliers_set)  # fitting ellipse for inliers points
                        self.count = inliers_count  # Update inliers set
                        self.best_model = ellipse_  # Update best ellipse
                        # print("self.count", self.count)

                    # 5.number of inliers points reach the expected value
                    if self.count > self.max_inliers:
                        # print('the number of inliers: ', self.count)
                        break

                    # Update items
                    # print(math.log(1 - pow(inliers_count / self.length, self.N)))
                    if math.log(1 - pow(inliers_count / self.length, self.N)) != 0:
                        self.items = math.log(1 - self.P) / math.log(1 - pow(inliers_count / self.length, self.N))

            Time_end = time.time()
            # print(Time_end - Time_start )
            if Time_end - Time_start >= 4:
                # print("time is too long")
                break

        return self.best_model

if __name__ == '__main__':

    # 1.find ellipse edge line
    contours, hierarchy = cv2.findContours(grayImage, cv2.RETR_CCOMP, cv2.CHAIN_APPROX_NONE)

    # 2.Ransac fit ellipse param
    points_data = np.reshape(contours, (-1, 2))  # ellipse edge points set
    Ransac = RANSAC(data=points_data, threshold=0.5, P=.99, S=.618, N=10)
    (X, Y), (LAxis, SAxis), Angle = Ransac.execute_ransac()

以上就是Python基于随机采样一至性实现拟合椭圆(优化版)的详细内容,更多关于Python拟合椭圆的资料请关注我们其它相关文章!

(0)

相关推荐

  • python opencv常用图形绘制方法(线段、矩形、圆形、椭圆、文本)

    最近学了下 python opencv,分享下使用 opencv 在图片上绘制常用图形的方法. 案例中实现了在图片中添加线段.圆形.矩形.椭圆形以及添加文字的方法,使用 opencv2 实现的. 实现方法 1)画线段 cv.line 在图片中绘制一段直线 # 绘制线段 # 参数1:图片 # 参数2:起点 # 参数3:终点 # 参数4:BGR颜色 # 参数5:宽度 cv2.line(img, (60, 40), (90, 90), (255, 255, 255), 2); 参数说明 参数 值 说明

  • python+matplotlib绘制旋转椭圆实例代码

    旋转椭圆 实例代码: import matplotlib.pyplot as plt import numpy as np from matplotlib.patches import Ellipse delta = 45.0 # degrees angles = np.arange(0, 360 + delta, delta) ells = [Ellipse((1, 1), 4, 2, a) for a in angles] a = plt.subplot(111, aspect='equal

  • python编程实现随机生成多个椭圆实例代码

    椭圆演示: 代码示例: import matplotlib.pyplot as plt import numpy as np from matplotlib.patches import Ellipse NUM = 250 ells = [Ellipse(xy=np.random.rand(2) * 10, width=np.random.rand(), height=np.random.rand(), angle=np.random.rand() * 360) for i in range(N

  • python opencv圆、椭圆与任意多边形的绘制实例详解

    圆形的绘制 : OpenCV中使用circle(img,center,radius,color,thickness=None,lineType=None,shift=None)函数来绘制圆形 import cv2 import numpy as np image=np.zeros((400,400,3),np.uint8) cv2.circle(image,(200,200),50,(0,0,255),2) #画圆 ''' 参数2 center:必选参数.圆心坐标 参数3 radius:必选参数

  • Python基于随机采样一至性实现拟合椭圆

    检测这些圆,先找轮廓后通过轮廓点拟合椭圆 import cv2 import numpy as np import matplotlib.pyplot as plt import math from Ransac_Process import RANSAC def lj_img(img): wlj, hlj = img.shape[1], img.shape[0] lj_dis = 7 # 连接白色区域的判定距离 for ilj in range(wlj): for jlj in range(h

  • Python实现霍夫圆和椭圆变换代码详解

    在极坐标中,圆的表示方式为: x=x0+rcosθ y=y0+rsinθ 圆心为(x0,y0),r为半径,θ为旋转度数,值范围为0-359 如果给定圆心点和半径,则其它点是否在圆上,我们就能检测出来了.在图像中,我们将每个非0像素点作为圆心点,以一定的半径进行检测,如果有一个点在圆上,我们就对这个圆心累加一次.如果检测到一个圆,那么这个圆心点就累加到最大,成为峰值.因此,在检测结果中,一个峰值点,就对应一个圆心点. 霍夫圆检测的函数: skimage.transform.hough_circle

  • Python基于随机采样一至性实现拟合椭圆(优化版)

    上次版本如果在没有找到轮廓或轮廓的点集数很小无法拟合椭圆或在RANSAC中寻找最优解时会死循环中,优化后的代码 import cv2 import os import numpy as np import matplotlib.pyplot as plt import math from Ransac_Process import RANSAC def cul_area(x_mask, y_mask, r_circle, mask): mask_label = mask.copy() num_a

  • Python实现从概率分布中随机采样

    目录 1. 二项(binomial)/伯努利(Bernoulli)分布 2. 多项(multinomial)分布 3.均匀(uniform)分布 4. 狄利克雷(Dirichlet)分布 参考文献 在上一篇博文<Python中的随机采样和概率分布(一)>中,我们介绍了Python中最简单的随机采样函数.接下来我们更进一步,来看看如何从一个概率分布中采样,我们以几个机器学习中最常用的概率分布为例. 1. 二项(binomial)/伯努利(Bernoulli)分布 1.1 概率质量函数(pmf)

  • 浅析Python中的随机采样和概率分布

    目录 1. random.choice 2. random.choices(有放回) 3. numpy.sample(无放回) 4.rng.choices 和 rng.sample 5. numpy.random.choices 参考文献  Python(包括其包Numpy)中包含了了许多概率算法,包括基础的随机采样以及许多经典的概率分布生成.我们这个系列介绍几个在机器学习中常用的概率函数.先来看最基础的功能--随机采样. 1. random.choice 如果我们只需要从序列里采一个样本(所有

  • Python随机采样及概率分布(二)

    目录 1. 二项(binomial)/伯努利(Bernoulli)分布 1.1 概率质量函数(pmf) 1.2 函数原型 1.3 使用样例 2. 多项(multinomial)分布 2.1 概率质量函数(pmf) 2.2 函数原型 2.3 使用样例 3.均匀(uniform)分布 3.1 概率密度函数(pdf) 3.2 函数原型 3.3 使用样例 4. 狄利克雷(Dirichlet)分布 4.1 概率密度函数(pdf) 4.2 函数原型 4.3 使用样例 前言: 之前的<Python中的随机采样

  • 基于Python实现随机点名系统的示例代码

    目录 效果展示 代码展示 导入模块 子线程调用 应用初始化信息 姓名信息布局 开始信息布局 数据信息布局 整体布局 运行 大家好,我是了不起! 在某些难以抉择得时候,我们经常要用外力来帮助我们做出选择 比如,梁山出征方腊前沙场点兵,挑选先锋的场景 这个时候,有一个随机点名系统就非常好啦,毕竟我水泊梁山的名号~ 效果展示 创建一个这样的文件夹,然后把要随机点名的名字写在里面 导入后,这里就显示你导入了多少人员信息 点击开始点名后,会随机从导入名字里挑选一位幸运儿~ 效果大概就是这样,下面我们来看看

  • python基于机器学习预测股票交易信号

    引言 近年来,随着技术的发展,机器学习和深度学习在金融资产量化研究上的应用越来越广泛和深入.目前,大量数据科学家在Kaggle网站上发布了使用机器学习/深度学习模型对股票.期货.比特币等金融资产做预测和分析的文章.从金融投资的角度看,这些文章可能缺乏一定的理论基础支撑(或交易思维),大都是基于数据挖掘.但从量化的角度看,有很多值得我们学习参考的地方,尤其是Pyhton的深入应用.数据可视化和机器学习模型的评估与优化等.下面借鉴Kaggle上的一篇文章<Building an Asset Trad

  • Python基于纹理背景和聚类算法实现图像分割详解

    目录 一.基于纹理背景的图像分割 二.基于K-Means聚类算法的区域分割 三.总结 一.基于纹理背景的图像分割 该部分主要讲解基于图像纹理信息(颜色).边界信息(反差)和背景信息的图像分割算法.在OpenCV中,GrabCut算法能够有效地利用纹理信息和边界信息分割背景,提取图像目标物体.该算法是微软研究院基于图像分割和抠图的课题,它能有效地将目标图像分割提取,如图1所示[1]. GrabCut算法原型如下所示: mask, bgdModel, fgdModel = grabCut(img,

  • 用Python实现随机森林算法的示例

    拥有高方差使得决策树(secision tress)在处理特定训练数据集时其结果显得相对脆弱.bagging(bootstrap aggregating 的缩写)算法从训练数据的样本中建立复合模型,可以有效降低决策树的方差,但树与树之间有高度关联(并不是理想的树的状态). 随机森林算法(Random forest algorithm)是对 bagging 算法的扩展.除了仍然根据从训练数据样本建立复合模型之外,随机森林对用做构建树(tree)的数据特征做了一定限制,使得生成的决策树之间没有关联,

  • Python基于pygame实现的font游戏字体(附源码)

    本文实例讲述了Python基于pygame实现的font游戏字体.分享给大家供大家参考,具体如下: 在pygame游戏开发中,一个友好的UI中,漂亮的字体是少不了的 今天就给大伙带来有关pygame中字体的一些介绍说明 首先我们得判断一下我们的pygame中有没有font这个模块 复制代码 代码如下: if not pygame.font: print('Warning, fonts disabled') 如果有的话才可以进行接下来的操作:-) 我们可以这样使用pygame中的字体: 复制代码

随机推荐