在skimage中绘制渐变椭圆

Question

我想在渐变颜色的skimage中绘制一个椭圆形蒙版.颜色从椭圆内部开始变化,并在外部椭圆处结束.如何用skimage或open-cv绘制它？

如下图所示:

Answer 1

介绍

让我们首先详细描述样本图像.

• 它是一个4通道图像(RGB + alpha透明度),但它只使用灰色阴影.
• 图像非常适合绘图,形状周围只有极小的余量.
• 有一个填充的,抗锯齿的,旋转的外椭圆,周围是透明的背景.
• 所述外椭圆内填充有旋转椭圆渐变(同心,并用相同的旋转为椭圆),这是黑色的边缘,在中心直线前进到白色.
• 所述外椭圆覆盖有一个同心的,填充的,反走样,旋转内椭圆(再次相同的旋转,这两个轴都以相同的比例缩放).填充颜​​色为白色.

此外,让:

• a并且b是椭圆的半长轴和半短轴
• theta是椭圆周围椭圆的旋转角度(a沿x轴为0 点)
• inner_scale是内部与外部椭圆的相应轴之间的比率(即0.5表示内部按比例缩小50%)
• h并且k是椭圆中心的(x,y)坐标

在演示的代码中,我们将使用以下导入:

import cv2
import numpy as np
import math

并设置定义我们的绘图的参数(我们将计算h和k)为:

a, b = (360.0, 200.0) # Semi-major and semi-minor axis
theta = math.radians(40.0) # Ellipse rotation (radians)
inner_scale = 0.6 # Scale of the inner full-white ellipse

步骤1

为了生成这样的图像,我们需要采取的第一步是计算我们需要的"画布"(我们将要绘制的图像)的大小.为此,我们可以计算旋转外椭圆的边界框,并在其周围添加一些小边距.

我不知道现有的OpenCV函数可以有效地执行此操作,但StackOverflow可以节省时间 - 已经有一个相关问题,答案链接到讨论问题的有用文章.我们可以使用这些资源来实现以下Python实现:

def ellipse_bbox(h, k, a, b, theta):
    ux = a * math.cos(theta)
    uy = a * math.sin(theta)
    vx = b * math.cos(theta + math.pi / 2)
    vy = b * math.sin(theta + math.pi / 2)
    box_halfwidth = np.ceil(math.sqrt(ux**2 + vx**2))
    box_halfheight = np.ceil(math.sqrt(uy**2 + vy**2))
    return ((int(h - box_halfwidth), int(k - box_halfheight))
        , (int(h + box_halfwidth), int(k + box_halfheight)))

注意:我将浮点大小向上舍入,因为我们必须覆盖整个像素,并将左上角和右下角作为整数(x,y)对返回.

然后我们可以按以下方式使用该函数:

# Calculate the image size needed to draw this and center the ellipse
_, (h, k) = ellipse_bbox(0, 0, a, b, theta) # Ellipse center
h += 2 # Add small margin
k += 2 # Add small margin
width, height = (h*2+1, k*2+1) # Canvas size

第2步

第二步是生成透明层.这是单通道8位图像,其中黑色(0)表示完全透明,白色(255)表示完全不透明的像素.这个任务很简单,因为我们可以使用cv2.ellipse.

我们可以将外椭圆定义为一个RotatedRect结构(一个与椭圆紧密配合的旋转矩形).在Python中,这表示为包含以下内容的元组:

• 表示旋转矩形中心的元组(x和y坐标)
• 表示旋转矩形大小的元组(宽度和高度)
• 旋转角度

这是代码:

ellipse_outer = ((h,k), (a*2, b*2), math.degrees(theta))

transparency = np.zeros((height, width), np.uint8)
cv2.ellipse(transparency, ellipse_outer, 255, -1, cv2.LINE_AA)

......以及它产生的图像:

第3步

作为第三步,我们创建一个包含我们所需的旋转椭圆梯度的单通道(灰度或强度)图像.但首先,我们怎么连数学在我们的形象的直角来定义一个旋转的椭圆(x, y)坐标,用我们的(a, b),theta(θ)和(h, k)参数？

这次数学StackExchange是一个拯救一天的人:有一个问题与我们的问题完全匹配,并提供了这个有用的等式的答案:

观察到,对于我们从椭圆中心取的任何方向,左侧在椭圆的周长处评估为1.它在中心处为0,并且在周边处朝向1线性增加,然后进一步超过它.

weight由于缺乏更好的术语,让我们称之为右手边.由于它从中心向外扩展得非常好,我们可以用它来计算我们想要的梯度.我们的公式在外部给出了白色(浮点图像为1.0),中心为黑色(0.0).我们想要逆,所以我们简单地weight从1.0结果中减去并剪切结果[0.0, 1.0].

让我们从一个简单的Python实现开始(如手动迭代numpy.array代表我们图像的各个元素)来计算权重.然而,由于我们是懒惰的程序员,我们将使用Numpy将计算的weights转换为分级图像,使用矢量化减法,以及numpy.clip.

这是代码:

def make_gradient_v1(width, height, h, k, a, b, theta):
    # Precalculate constants
    st, ct =  math.sin(theta), math.cos(theta)
    aa, bb = a**2, b**2

    weights = np.zeros((height, width), np.float64)    
    for y in range(height):
        for x in range(width):
            weights[y,x] = ((((x-h) * ct + (y-k) * st) ** 2) / aa
                + (((x-h) * st - (y-k) * ct) ** 2) / bb)

    return np.clip(1.0 - weights, 0, 1)

......以及它产生的图像:

这很好,花花公子,但是因为我们迭代每个像素并在Python解释器中进行计算,所以它也是aaawwwfffuuullllllyyy ssslllooowww ....可能需要一秒钟,但我们正在使用Numpy,所以如果我们可以做得更好利用它.这意味着我们可以做任何事情.

首先,让我们注意到,唯一不同的输入是每个给定像素的坐标.这意味着要对我们的算法进行矢量化,我们需要输入两个数组(与图像大小相同),保持每个像素的坐标x和y坐标.幸运的是,Numpy为我们提供了生成此类数组的工具 - numpy.mgrid.我们可以写

y,x = np.mgrid[:height,:width]

生成我们需要的输入数组.然而,让我们看到,我们从来不使用x和y直接-而我们总是通过不断弥补他们.让我们通过生成x-h和y-k...来避免偏移操作

y,x = np.mgrid[-k:height-k,-h:width-h]

我们可以再次预先计算4个常数,除此之外,其余的只是矢量化加法,减法,乘法,除法和幂,它们都是由Numpy提供的矢量化操作(即更快).

def make_gradient_v2(width, height, h, k, a, b, theta):
    # Precalculate constants
    st, ct =  math.sin(theta), math.cos(theta)
    aa, bb = a**2, b**2

    # Generate (x,y) coordinate arrays
    y,x = np.mgrid[-k:height-k,-h:width-h]
    # Calculate the weight for each pixel
    weights = (((x * ct + y * st) ** 2) / aa) + (((x * st - y * ct) ** 2) / bb)

    return np.clip(1.0 - weights, 0, 1)

与仅使用Python的脚本相比,使用此版本的脚本大约需要30%的时间.没有什么令人惊叹的,但它产生了相同的结果,这项任务似乎是你不必经常做的事情,所以这对我来说已经足够了.

如果您[读者]知道更快的方式,请将其作为答案发布.

现在我们有一个浮点图像,强度范围在0.0到1.0之间.为了生成我们的结果,我们希望得到一个8位图像,其值介于0到255之间.

intensity = np.uint8(make_gradient_v2(width, height, h, k, a, b, theta) * 255)

第4步

第四步 - 绘制内椭圆.这很简单,我们之前已经完成了.我们只需要适当地缩放轴.

ellipse_inner = ((h,k), (a*2*inner_scale, b*2*inner_scale), math.degrees(theta))

cv2.ellipse(intensity, ellipse_inner, 255, -1, cv2.LINE_AA)

这给了我们以下强度图像:

第五步

第五步 - 我们快到了.我们所要做的就是将强度和透明度层组合成BGRA图像,然后将其保存为PNG.

result = cv2.merge([intensity, intensity, intensity, transparency])

注意:对红色,绿色和蓝色使用相同的强度只给我们灰色阴影.

保存结果后,我们得到以下图像:

结论

鉴于我猜测你用来生成样本图像的参数,我会说我的脚本的结果非常接近.它运行得相当快 - 如果你想要更好的东西,你可能无法避免接近裸机(C,C++等).更智能的方法,或者GPU也许可以做得更好.值得尝试......

总而言之,这里有一点证明这个代码也适用于其他轮换:

我用来写这个完整的脚本:

import cv2
import numpy as np
import math

# ============================================================================

def ellipse_bbox(h, k, a, b, theta):
    ux = a * math.cos(theta)
    uy = a * math.sin(theta)
    vx = b * math.cos(theta + math.pi / 2)
    vy = b * math.sin(theta + math.pi / 2)
    box_halfwidth = np.ceil(math.sqrt(ux**2 + vx**2))
    box_halfheight = np.ceil(math.sqrt(uy**2 + vy**2))
    return ((int(h - box_halfwidth), int(k - box_halfheight))
        , (int(h + box_halfwidth), int(k + box_halfheight)))

# ----------------------------------------------------------------------------

# Rotated elliptical gradient - slow, Python-only approach
def make_gradient_v1(width, height, h, k, a, b, theta):
    # Precalculate constants
    st, ct =  math.sin(theta), math.cos(theta)
    aa, bb = a**2, b**2

    weights = np.zeros((height, width), np.float64)    
    for y in range(height):
        for x in range(width):
            weights[y,x] = ((((x-h) * ct + (y-k) * st) ** 2) / aa
                + (((x-h) * st - (y-k) * ct) ** 2) / bb)

    return np.clip(1.0 - weights, 0, 1)

# ----------------------------------------------------------------------------

# Rotated elliptical gradient - faster, vectorized numpy approach
def make_gradient_v2(width, height, h, k, a, b, theta):
    # Precalculate constants
    st, ct =  math.sin(theta), math.cos(theta)
    aa, bb = a**2, b**2

    # Generate (x,y) coordinate arrays
    y,x = np.mgrid[-k:height-k,-h:width-h]
    # Calculate the weight for each pixel
    weights = (((x * ct + y * st) ** 2) / aa) + (((x * st - y * ct) ** 2) / bb)

    return np.clip(1.0 - weights, 0, 1)

# ============================================================================ 

def draw_image(a, b, theta, inner_scale, save_intermediate=False):
    # Calculate the image size needed to draw this and center the ellipse
    _, (h, k) = ellipse_bbox(0,0,a,b,theta) # Ellipse center
    h += 2 # Add small margin
    k += 2 # Add small margin
    width, height = (h*2+1, k*2+1) # Canvas size

    # Parameters defining the two ellipses for OpenCV (a RotatedRect structure)
    ellipse_outer = ((h,k), (a*2, b*2), math.degrees(theta))
    ellipse_inner = ((h,k), (a*2*inner_scale, b*2*inner_scale), math.degrees(theta))

    # Generate the transparency layer -- the outer ellipse filled and anti-aliased
    transparency = np.zeros((height, width), np.uint8)
    cv2.ellipse(transparency, ellipse_outer, 255, -1, cv2.LINE_AA)
    if save_intermediate:
        cv2.imwrite("eligrad-t.png", transparency) # Save intermediate for demo

    # Generate the gradient and scale it to 8bit grayscale range
    intensity = np.uint8(make_gradient_v1(width, height, h, k, a, b, theta) * 255)
    if save_intermediate:
        cv2.imwrite("eligrad-i1.png", intensity) # Save intermediate for demo

    # Draw the inter ellipse filled and anti-aliased
    cv2.ellipse(intensity, ellipse_inner, 255, -1, cv2.LINE_AA)
    if save_intermediate:
        cv2.imwrite("eligrad-i2.png", intensity) # Save intermediate for demo

    # Turn it into a BGRA image
    result = cv2.merge([intensity, intensity, intensity, transparency])
    return result

# ============================================================================ 

a, b = (360.0, 200.0) # Semi-major and semi-minor axis
theta = math.radians(40.0) # Ellipse rotation (radians)
inner_scale = 0.6 # Scale of the inner full-white ellipse

cv2.imwrite("eligrad.png", draw_image(a, b, theta, inner_scale, True))

# ============================================================================ 

rows = []
for j in range(0, 4, 1):
    cols = []
    for i in range(0, 90, 10):
        tile = np.zeros((170, 170, 4), np.uint8)
        image = draw_image(80.0, 50.0, math.radians(i + j * 90), 0.6)
        tile[:image.shape[0],:image.shape[1]] = image
        cols.append(tile)
    rows.append(np.hstack(cols))

cv2.imwrite("eligrad-m.png", np.vstack(rows))

注意:如果您发现任何愚蠢的错误,令人困惑的术语或此帖的任何其他问题,请随意留下建设性的评论,或者直接编辑答案以使其更好.我知道有很多方法可以进一步优化这一点 - 让读者做这个练习(也许提供和补充答案).