计算二值图像的圆度

Question

我将为您翻译一段关于IT技术的内容。该内容涉及计算给定二进制图像的圆度。经过一些研究，我明确了圆度公式如下：

4π*area/perimeter^2

值应该在0到1之间，其中1表示最圆。

给定二进制矩阵 im

计算面积很简单

面积 = sum(im)

我按照以下规则计算周长: 如果一个像素不为零且连接至少一个零值像素，则它是周长的一部分

per = matrix(0, nrow(im), ncol(im))
for(i in 2:(nrow(im)-1)){
  for(j in 2:(ncol(im)-1)){
    if(im[i,j] != 0){
      x=c(im[i-1,j],im[i+1,j],im[i,j-1], im[i,j+1])
      if(0 %in% x) per[i,j] = 1
    }
  }
}
perimeter = sum(per)

然后我这样计算循环性：

circ = (4*pi*area)/(perimiter^2)

然而，有时我会得到大于1的值，而且结果不正确。例如：

这张图片给出了 circ=1.155119 而这张图片给出了 circ=1.148728 你知道发生了什么吗？ 这些值不应该更接近于0.95和0.7吗？

- Omar Wagih

Answer 1

您对“二进制周长”的定义并不能很好地近似平滑周长。

# Sample data
n <- 100
im <- matrix(0, 3*n, 3*n+1)
x <- ( col(im) - 1.5*n ) / n
y <- ( row(im) - 1.5*n ) / n
im[ x^2 + y^2 <= 1 ] <- 1
image(im)

# Shift the image in one direction
s1 <- function(z) cbind(rep(0,nrow(z)), z[,-ncol(z)] )
s2 <- function(z) cbind(z[,-1], rep(0,nrow(z)) )
s3 <- function(z) rbind(rep(0,ncol(z)), z[-nrow(z),] )
s4 <- function(z) rbind(z[-1,], rep(0,ncol(z)) )

# Area, perimeter and circularity
area <- function(z) sum(z)
perimeter <- function(z) sum( z != 0 & s1(z)*s2(z)*s3(z)*s4(z) == 0)
circularity <- function(z) 4*pi*area(z) / perimeter(z)^2

circularity(im)
# [1] 1.241127

area(im)
# [1] 31417
n^2*pi
# [1] 31415.93

perimeter(im)
# [1] 564
2*pi*n
# [1] 628.3185

一个令人担忧的特性是其周长不具备旋转不变性：当你将一个边长为1（与坐标轴平行）的正方形旋转45度时，它的面积保持不变，但其周长被除以sqrt(2)...

square1 <- -1 <= x & x <= 1 & -1 <= y & y <= 1
c( perimeter(square1), area(square1) )
# [1]   800 40401

square2 <- abs(x) + abs(y) <= sqrt(2)
c( perimeter(square2), area(square2) )
# [1]   564 40045

这里是更好的周长近似值。 对于周长上的每个点， 查看其8邻域中哪些点也在周长上; 如果它们形成一个垂直或水平线段， 则该对点对周长的贡献为1， 如果它们在对角线上，则贡献为sqrt(2)。

edge <- function(z) z & !(s1(z)&s2(z)&s3(z)&s4(z))
perimeter <- function(z) {
  e <- edge(z)
  ( 
    # horizontal and vertical segments
    sum( e & s1(e) ) + sum( e & s2(e) ) + sum( e & s3(e) ) + sum( e & s4(e) ) + 
    # diagonal segments
    sqrt(2)*( sum(e & s1(s3(e))) + sum(e & s1(s4(e))) + sum(e & s2(s3(e))) + sum(e & s2(s4(e))) )
  ) / 2  # Each segment was counted twice, once for each end
}

perimeter(im)
# [1] 661.7544
c( perimeter(square1), area(square1) )
# [1]   805.6569 40401.0000
c( perimeter(square2), area(square2) )
# [1]   797.6164 40045.0000

circularity(im)
# [1] 0.9015315
circularity(square1)
# [1] 0.7821711
circularity(square2)
# [1] 0.7909881

- Vincent Zoonekynd

Answer 2

让我建议一个不同的算法。

仅仅通过计算像素数量，您应该可以非常准确地获取 blob 的面积。您还可以通过对每个内部像素坐标取平均值来找到中心位置。现在您可以使用 sqrt(area/pi) 找到半径。有了半径和中心位置，您就可以画出一个完美的圆形，其面积几乎相同 - 计算同时属于 blob 和完美圆形的像素数量，并除以先前计算出的面积。

- Mark Ransom

Answer 3

如果我没记错的话，对于给定的周长，一个完美的圆形可以围住最大的面积。因此，如果一个圆的比率是1.0，那么根据你的公式，任何其他形状的比率都会大于1.0。

- Nishanth

---

这是反过来的：对于一个非圆形，面积会更小，因此面积/周长^2的比率将小于1。例如，对于边长为1的正方形，面积=1，周长=4，圆度=pi/4 < 1。 - Vincent Zoonekynd