超像素聚类

Definition: 

 定义： 

A superclass is a class that has been extended by another class. It allows the extending class to inherit its state and behaviors.

 超类是被另一个类扩展的类。 它允许扩展类继承其状态和行为。 

Also Known As: base class, parent class

 也称为：基类，父类 

Examples: 

 例子： 

Imagine you define a Person class:

 假设您定义了一个Person类： 
public class Person
{
}

A new class can be created by extending this class:

 可以通过扩展此类来创建一个新类： 
public class Employee extends Person
{
}

The Person class is said to be the superclass of the Employee class.

 据说Person类是Employee类的超类。 
翻译自: https://www.thoughtco.com/superclass-2034298
超像素聚类

使用K-means进行像素聚类



python通过使用K-means对像素聚类以此进行图像分割

K-means聚类简介：

K-means算法是最为经典的基于划分的聚类方法，是十大经典数据挖掘算法之一。K-means算法的基本思想是：以空间中k个点为中心进行聚类，对最靠近他们的对象归类。通过迭代的方法，逐次更新各聚类中心的值，直至得到最好的聚类结果。

因此，简单来说Kmeans是一种将输入数据划分为k个族的简单的聚类算法。

算法步骤（基本K-means）：



repeat  
    将每个点指派到最近的质心，形成K个簇  
    重新计算每个簇的质心  
until 簇不发生变化或达到最大迭代次数 

步骤

（1）适当选择k个类的初始中心；一般是采用随机或猜测的方式初始化类中心，(注：采用随机的方法实现简单，但是族的质量往往比较差，所以有好几种关于中心选取的解决方案，比如先使用层次聚类进行聚类，从层次聚类中提取K个簇，并用这些簇的质心作为初始质心。也有通过使类内总方差最小的方式，选择方差最小的类中心。这篇博客不是讨论算法本身，所以这里就不过多叙述，感兴趣的同学可以去查一些资料。)
（2）对任意一个样本，求其到k个中心的距离，将该样本归到距离最短的中心所在的类；(注：常用的距离度量方法包括：欧几里得距离和余弦相似度。两者都是评定个体间差异的大小的。欧几里得距离度量会受指标不同单位刻度的影响，所以一般需要先进行标准化。余弦相似度倾向给出更优解)
（3）对所有属于该类的数据点求平均，将平均值作为新的类中心；（取平均即向量各维取平均）
（4）重复步骤（2）（3）直到收敛。
该算法的最大优势在于简洁和快速，最重要的是容易实现与可并行。最大的缺陷在于要预先设定聚类数k，选择不当则会导致聚类效果很差。算法的关键在于初始中心的选择和距离公式。

算法复杂度：

时间复杂度：O(tKmn)，其中，t为迭代次数，K为簇的数目，m为记录数，n为维数

空间复杂度：O((m+K)n)，其中，K为簇的数目，m为记录数，n为维数

使用scipy中的K-means：

尽 管 K-means 算 法 很 容 易 实 现, 但 我 们 没 有 必 要 自 己 实 现 它。 SciPy 矢 量 量 化 包 

scipy.cluster.vq 中有 K-means 的实现,下面是使用方法。

scipy K-means 的简单实例:



# -*- coding: utf-8 -*-

#导入scipy中K-means的相关工具
from scipy.cluster.vq import *

#randn是NumPy中的一个函数
from numpy import *

from pylab import *

#生成简单的二维数据:生成两类二维正态分布数据。
class1 = 1.5 * randn(100,2)
class2 = randn(100,2) + array([5,5])
features = vstack((class1,class2))

#用 k=2 对这些数据进行聚类:
centroids,variance = kmeans(features,2)

"""

由于 SciPy 中实现的 K-means 会计算若干次(默认为 20 次),并为我们选择方差最
小的结果,所以这里返回的方差并不是我们真正需要的。

"""

#用 SciPy 包中的矢量量化函数对每个数据点进行归类:通过得到的 code ,我们可以检查是否有归类错误

code,distance = vq(features,centroids)

#可视化结果：画出这些数据点及最终的聚类中心:函数 where() 给出每个类的索引

figure()

ndx = where(code==0)[0]
plot(features[ndx,0],features[ndx,1],'*')
ndx = where(code==1)[0]
plot(features[ndx,0],features[ndx,1],'r.')
plot(centroids[:,0],centroids[:,1],'go')
axis('off')
show()


绘制结果：



K-means像素聚类

将图像区域或像素合并成有意义的部分称为图像分割。单纯在像素水平上应用 K-means可以用于一些简单图像的图像分割，但是对于复杂图像得出的结果往往是毫无意义的。对于复杂图像的图像分割往往需要更复杂的类模型而非平均像素色彩或空间一致性。

下面在 RGB 三通道的像素值上运用 K-means 进行聚类：

 K-means RGB三通道像素聚类

对像素进行聚类首先要降分辨率处理，否则像素的值太多，运行效率会相当差，噪点也一定很多，一般情况下传统的方式会以steps×steps个方格为单位对每个方格中的像素值进行均化以达到降低分辨率的目的，但是如果图像不是正方形或者step值不合适的话，结果就会出现严重的变形，因此我将原本step×step的方式改为了stepX*stepY并且根据输入的step自动计算合理step值，成功解决了这个问题，以下代码中clusterpixels_square为传统的像素聚类，clusterpixels_rectangular是我改进之后的聚类方法，希望对大家有所帮助：



# -*- coding: utf-8 -*-

"""
对像素进行聚类。
在像素级水平进行聚类可以用在一些很简单的图像

载入图像，并将其下采样到一个较低的分辨率，然后对这些区域用k-means进行聚类

K-means 的输入是一个有 stepsX × stepsY 行的数组,数组的每一行有 3 列,各列分别为区域块 R、G、B 三个通道的像素平均值。

为可视化最后的结果 , 我们用 SciPy 的imresize() 函数在原图像坐标中显示这幅图像。

参数 interp 指定插值方法;我们在这里采用最近邻插值法,以便在类间进行变换时不需要引入新的像素值。

"""

from scipy.cluster.vq import *
from scipy.misc import imresize

from pylab import *

from PIL import Image

#steps*steps像素聚类
def clusterpixels_square(infile, k, steps):

    im = array(Image.open(infile))    

    #im.shape[0] 高 im.shape[1] 宽
    dx = im.shape[0] / steps
    dy = im.shape[1] / steps
    # 计算每个区域的颜色特征
    features = []
    for x in range(steps):
        for y in range(steps):
            R = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 0])
            G = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 1])
            B = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 2])
            features.append([R, G, B])
    features = array(features, 'f')     # 变为数组
    # 聚类， k是聚类数目
    centroids, variance = kmeans(features, k)
    code, distance = vq(features, centroids)

    # 用聚类标记创建图像
    codeim = code.reshape(steps, steps)
    codeim = imresize(codeim, im.shape[:2], 'nearest')
    return codeim

#stepsX*stepsY像素聚类
def clusterpixels_rectangular(infile, k, stepsX):

    im = array(Image.open(infile))

    stepsY = stepsX * im.shape[1] / im.shape[0]

    #im.shape[0] 高 im.shape[1] 宽
    dx = im.shape[0] / stepsX
    dy = im.shape[1] / stepsY
    # 计算每个区域的颜色特征
    features = []
    for x in range(stepsX):
        for y in range(stepsY):
            R = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 0])
            G = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 1])
            B = mean(im[x * dx:(x + 1) * dx, y * dy:(y + 1) * dy, 2])
            features.append([R, G, B])
    features = array(features, 'f')     # 变为数组
    # 聚类， k是聚类数目
    centroids, variance = kmeans(features, k)
    code, distance = vq(features, centroids)
      # 用聚类标记创建图像
    codeim = code.reshape(stepsX, stepsY)
    codeim = imresize(codeim, im.shape[:2], 'nearest')
    return codeim



#计算最优steps 为保证速度以及减少噪点 最大值为maxsteps 其值为最接近且小于maxsteps 的x边长的约数
def getfirststeps(img,maxsteps):

    msteps = img.shape[0]

    n = 2

    while(msteps>maxsteps):

        msteps =   img.shape[0]/n
        n = n + 1    

    return msteps



#Test


#图像文件 路径
infile = './data/10.jpg'

im = array(Image.open(infile))

#参数
m_k = 10

m_maxsteps = 128

#显示原图empire.jpg
figure()

subplot(131)

title('source')

imshow(im)

# 用改良矩形块对图片的像素进行聚类
codeim= clusterpixels_rectangular(infile, m_k,getfirststeps(im,m_maxsteps))

subplot(132)

title('New steps = '+str(getfirststeps(im,m_maxsteps))+' K = '+str(m_k));

imshow(codeim)

#方形块对图片的像素进行聚类
codeim= clusterpixels_square(infile, 15, 200)

subplot(133)

title('Old steps = 200 K = '+str(m_k));

imshow(codeim)

show()



聚类结果：



我们看到，使用xy方向步数相等的方式出现了严重的图像变形缺失，使用我改良的方法则解决了这个问题





聚类效果方面，对于人物这种简单图片效果还是不错的。





但是对于风景这种复杂图片效果就比较一般了，可能采取分水岭这种专门的分割算法能达到一个比较好的效果。

0.摘要
在双色反射模型的基础上，最近的镜面反射高光去除技术通常会估计和群集照明色度值，​​以从单个图像中分离出漫反射和镜面反射分量。尽管这些技术能够获得视觉上令人愉悦的结果，但它们的聚类算法存在初始化错误或成本太高而无法实时计算的问题。在本文中，我们提出了一种高质量的像素聚类方法，该方法可以实时从单个图像中去除镜面高光。我们遵循先前的工作，并估计每个像素的最小和最大色度值。然后，我们分析了这些值在最小最大色度空间中的分布模式，以提出一种有效的像素聚类方法。然后，我们估计每个群集的强度比，以分离漫反射和镜面反射分量。最后，我们提出了优化策略，以针对CPU和GPU架构有效地实施我们的方法。在可用数据集中评估的实验结果表明，所提出的方法不仅更准确，而且比仅在CPU上运行的最新技术快两倍。在GPU上运行，我们证明了我们的方法需要约24毫秒的时间才能去除分辨率为3840×2160（4k）的图像中的镜面高光，这使得我们的GPU实现比当前状态快了一个数量级（20×）。先进的4k分辨率图像，同时准确提供所需的效果
1.说明
计算机视觉应用程序（例如内在图像分解）通常会假定场景中存在的材质纯粹是朗伯（即没有镜面反射），并将镜面反射区域视为噪声或异常值。在这种情况下，如果图像中存在镜面高光的较大区域，则应用程序的精度可能会严重下降。
为了提高这些应用程序的精度，可以将镜面高光去除用作输入图像中的预处理步骤。虽然大多数现有技术（例如[3] – [5]）仍以非交互式帧速率执行此任务，但是诸如对象跟踪[6]之类的应用需要实时解决方案以通过去除镜面高光来提高跟踪精度。
为了减轻镜面高光去除的任务，二向色反射模型[7]是最常用的假设，用于指导从单个图像中分离漫反射和镜面反射分量。在这样的模型中，颜色可以通过漫反射和镜面反射分量之和或漫反射和镜面色度的线性组合来描述在像素处观察到的光。通常，大多数现有的精确技术（例如，[5]，[8]，[9]）将输入像素的颜色转换为基于色度的空间，该空间允许通过像素聚类的方式来分离漫反射分量和镜面反射分量。尽管这些技术能够最小化图像中镜面反射高光的存在，但它们中的大多数都使用了聚类策略，这些策略对簇种子的初始化很差，或者在处理时间方面太昂贵。
在本文中，我们提出了一种算法从静止图像中实时，高质量地去除镜面高光，比相关工作更准确，更快捷，同时即使对于高分辨率图像也能提供实时性能。受沈和郑[8]的启发，我们提出了一种有效且准确的算法来对最小-最大色度空间中的像素进行聚类。然后，我们估计每个聚类的强度比，以分离输入图像的漫反射和镜面反射分量。我们进一步讨论了减少提议解决方案的处理时间，优化CPU和GPU架构获得的性能的策略。

在这项工作中，我们的主要贡献有以下三个方面：

1）一种有效的像素聚类算法，可改善聚类种子的初始化并比相关结果更准确工作；

2）一种优化策略，使所提出的解决方案的CPU实现速度比相关工作快两倍；

3）在GPU上实现所提出的解决方案的效果，与高工作量相比，其结果要比相关工作快一个数量级。分辨率图像，同时仍获得较高的准确率；
2.相关工作
在本节中，我们将介绍从单幅图像中去除镜面高光的领域中建议的相关工作。 对于现有的镜面高光去除技术的深入研究，我们请读者阅读Artusiet 等人作品。



图（1）    （a）的输入图像，我们为每个像素计算最小值（b），最大值（c）和范围（d）值。然后，我们从最小的图像中减去输入图像，以估算出无伪镜面反射的图像（e）。之后，我们估计无伪镜面反射图像的每个像素的色度值，并将相关像素的最小（f）和最大（g）色度通道存储到单独的图像中。建立了最小最大色度空间，并且仅在三个簇（h）中有效地对其进行了簇聚。接下来，对于每个聚类，我们使用最大值和范围值来估计单个强度比。最后，先前估计的三个强度比用于通过镜面高光（j）的分离来计算真实的无镜面图像（i）。
tan等人[11]和沉等。[12]率先使用无伪镜面反射图像（仅扩散色度图像的近似值）来指导优化框架并协助镜面反射高光去除的任务。由于这种框架通常需要5秒钟以上的时间才能去除低分辨率（例如480p）图像的镜面高光，因此其使用不适用于实时应用。
Shen和Cai [13]和Yanget al。[14]分别使用了连接的组件标签和迭代双边过滤来去除镜面高光。两种方法都比以前的方法快了两个数量级，可以实现交互式帧速率，但是仍然需要超过1秒的时间才能去除高分辨率（例如2160p）图像的镜面高光。
Shen和Zheng [8]将像素聚类到最小最大色度空间。然后，对于每个聚类，估计每个像素的最大和范围色度之间的强度比，并为整个聚类分配单个强度比，以准确地分离镜面反射像素和漫射像素。
而不是依赖于最大色度-强度空间[11]或最小-最大色度空间[8]，不同的颜色空间（例如Ch-Cv和HSI空间[15]，a * -b *空间[16]，色相空间[4]）和像素聚类策略（例如区域增长，均值位移，k均值，亲和力传播）已用于镜面高光去除。使用颜色分割算法（例如，区域增长，均值偏移）的技术[15]，[17]容易产生噪声伪像，并且不适用于纹理化或多色表面。使用k均值或亲和力传播作为聚类策略的技术[3] – [5]通常需要大量迭代才能收敛到准确的解决方案，从而产生远非实时的帧速率。同时，在交互处理时能够达到较高准确率的技术[8]，[9]可能会遭受初始簇种子选择错误的困扰。
在这项工作中，我们利用了Shen and Zheng [8]提出了一种改进，高效的像素聚类方案，该方案比原始方法更准确，更快。我们还提出了一种优化算法来估计每个群集的单个强度比，该算法用于分离漫反射和镜面反射强度。最后，我们展示了如何在GPU上有效地实现所提出的方法，据我们所知，在文献中对于镜面高光去除单个图像的任务，可以获得最佳的精度和处理时间。
3.实时特殊高光去除
在本节中，我们介绍实时镜面高光去除的方法。拟议方法的概述如图1所示。我们首先计算最小值（图1-（b）），最大值（图1-（c））和范围（图1-（d））的每个像素的值。然后，我们计算出无伪镜面反射的图像（图1-（e）），该图像用于估计每个像素的最小（图1-（f））和最大（图1-（g））色度值。我们分析了最小-最大色度空间中这些值的分布模式，并提出了仅将色度空间划分为三个簇的阳极聚类方案（图1-（h））。接下来，对于每个聚类，我们提出一个新算法来估计每个聚类的最大值和范围值的单一强度比，以实时生成无反射镜图像（图1-（i））。
让我们假设每个像素x在I（x）= [ Ir（x），Ig（x），Ib（x）] T分别存储红色，绿色和蓝色通道（图1-（a））。根据双色反射模型[7]，像素  I（x）可以由其漫反射D（x）和镜面反射S（x）反射分量之和组成:

I（x）= D（x）+ S（x）

或者，通过漫反射Λ和镜面反射Γ色度值:

I（x）= wd（x）Λ（x）+ ws（x）Γ

的线性组合，其中wd和ws是与表面几何上的漫反射和镜面反射有关的权重:

Λ（x） = [Λr（x），Λg（x），Λb（x）] T，

Γ   通常假定是[5]，[8]，[9]，并且对于输入图像Γ= [1/3,1/3,1/3] T，如本文其余部分所述。
对于相同的颜色表面，具有镜面反射高光的图像通常包含仅具有漫反射分量的像素，以及具有漫反射和镜面反射分量的像素。因此，要去除图像的镜面反射高光，我们需要定义一个色彩空间，在该空间中我们可以聚类并进一步将纯漫反射与镜面像素分开。接下来的相关工作[8]，让我们使用
强度比
的概念来进行漫反射和镜面反射分量分离。Imin和Imax是存储最小值的单通道图像（图1-（b））

Imin（x）= min（Ir（x），Ig （x），Ib（x））

= wd（x）Λmin（x）+ ws（x）Γ， ------------------------------------                           （1）

和最大值（图1-（c））

Imax（x）= max（Ir（x） ，Ig（x），Ib（x））

= wd（x）Λmax（x）+ ws（x）Γ，                      ------------------------------------（2）

每个像素的红色，绿色和蓝色强度，其中

Λmin（x）= min（Λr（ x），Λg（x），Λb（x））           --------------       （3）

和

Λmax（x）= max（Λr（x），Λg（x），Λb（x））。          -----------（4）

让Irange为存储Imax和Imin（图1-（d））

Irange（x）= Imax（x）-Imin（x） = wd（x）Λrange（x）     ---- （5）

之间的减法的单通道图像。

防止Irange（x）为0 ，可以将此项加一个小值≥0。
从（5）可以很容易地看出I range是无镜面的图像，因为该图像不再包含双色反射模型（2）的镜面项ws（x）Γ。在图1-（d）中也可以看到这种效果。现在，让我们将强度比Iratio（x）定义为最大值与范围值之间的比率

Iratio（x）= Imax（x）/ Irange（x）。（6）

对于具有纯漫反射，强度比由最大和范围漫反射色度之间的比率明确定义

Iratio（x）= wd（x）Λmax  /  wd（x）（Λmax-Λmin）

=Λmax / (Λmax-Λmin)。（7）

但是，对于具有漫反射和镜面反射，强度比为

Iratio（x）= wd（x）Λmax+ ws（x） Γ  /   wd（x）（Λmax-Λmin）（8）

从（7）和（8），我们可以说，

对于具有在相同颜色的表面或相同的漫反射色度上，镜面反射像素的强度比高于纯漫射像素的强度比。从这个意义上讲，强度比可以用作有效地分离漫反射和镜面像素的度量。

为了使用比率来分离漫反射和镜面像素，我们需要首先对具有几乎相同的漫射色度的像素进行聚类。

为了估算像素的漫反射色度，让我们使用无伪镜面反射的image Ipsf（图1-（e）），可以通过从Imin减去输入图像I轻松获得:

Ipsf（x）= I（x）-Imin（ x）= wd（x）Λpsf（x）。（9）

类似地，Ipsf不包含双色反射模型的镜面项（2）。

与Irange不同，Ipsf是类似于原始输入图像的三通道图像，允许估计红色，绿色和蓝色强度的漫反射色度。

为了进一步防止Ipsf比I暗，并提高图像噪声的鲁棒性，我们添加了均值



到Ipsf [13]的每个像素的Imin。

一旦有了Ipsf，每个像素的漫反射色度值Λ psf（x）估计为

Λ psf（x）= Ipsf（x） / (Ipsf r （x）+ Ipsf g （x）+ Ipsf b（x）)。（10） （7）和（8），基于最小和最大漫射色度值来计算像素的强度比。然后，我们利用Ipsf（9）和Λpsf（10）的优势来计算最小值:

Λpsfmin（x）= min（Λpsfr（x），Λpsfg（x），Λpsfb（x））

（图1-（f））和最大值

Λpsfmax（ x）= max（Λpsfr（x），Λpsfg（x），Λpsfb（x））

（图1-（g））



图 2.我们针对最小（x轴）最大（y轴）色度空间的聚类方法。 在算法（a）的第一次迭代中，我们选择位于最高最小值（蓝色圆圈），最高最大值（绿色圆圈）和最低最小色度值（红色圆圈）处的三个簇种子（箭头指向的大圆圈）。 然后，我们运行两次k均值迭代，以将每个像素与欧几里得空间中最近的簇种子相关联，根据其相关像素的质心更新簇种子的位置，并根据更新后的簇种子将每个像素相关联（ b）。
伪漫射色度值仅适用于相关可能位于具有镜面反射高光的区域中的像素。 如下确定像素是否与镜面高光去除有关


从这个意义上讲，背景像素和强度值低于平均最小强度Imin的像素由于可能不包含镜面像素而被进一步计算，因此我们将其剔除，以将具有相同扩散色度的像素聚类，将Ipsf的相关像素投影到最小-最大色度中空间，如图2所示。该空间中的任何像素都以三角形的形式投影（此陈述的证明可以在补充文档的第一节中看到），如图2所示。而不是选择数字和位置在初始簇种子中，我们利用最小最大色度空间样本位于其中的三角形形状，为在该空间中进行像素聚类提出了一种有效的方法。我们已经看到，而不是将空间分成一个或一个两个集群，分离三个群集中的最小最大色度空间为相关像素的群集提供了准确的解决方案。为此，我们首先选择其像素在最小-最大色度空间中存储最高最小，最高最大和最低色度值的像素作为聚类种子，如图2-（a）中的大圆圈所示。这些像素在我们的算法中被选择为簇种子，因为它们可以很容易地从Λpsfmin和Λpsfmax中检测出来，而它们的投影是在最小最大色度空间中形成的三角形形状的顶点。然后，我们仅运行两次k均值迭代[18]，以将输入图像的每个相关像素与最小-最大色度空间中最近的群集种子相关联，以将每个群集种子的位置更新为其关联像素的中心。


图3.给定每个群集的强度比的无序数组（a），我们可以对数组（b）进行排序，并选择中值（（b）中的蓝色矩形）作为所选强度比，因为此值将数组相等地扩散（绿色矩形）和镜面反射（红色矩形）强度比率（c）。另一种方法是从无序数组中迭代选择一个值（（d）中的蓝色矩形），以更好地分离群集中的漫反射强度和镜面反射强度比（e） 。
在最小最大色度空间中（请参见图2-（b）中的簇种子的新位置），并将每个像素与最近的更新簇种子相关联（图2-（b））。使用此算法，我们能够为了消除随机性以选择初始聚类种子，提供了一种改进的初始化方案，该方案对于图像的方向是不变的，可以轻松地从Ipsf计算得出，并且需要对k均值进行两次迭代才能收敛到一个精确的解决方案。在相同的漫反射色度下，我们需要为每个聚类选择基于比例（7）和（8）将漫反射像素与镜面像素分开的强度比率（图3）。每个群集计算机的一个聪明的想法是，以升序对与相同群集种子相关联的每个像素的强度比进行排序（图3-（b）），然后将位于已排序阵列中心位置的强度比作为选定对象valuer（图3-（b）中的1.6），因为它几乎将簇分为漫射像素和镜面像素（图3-（c））。请注意，根据输入数组的大小，排序步骤可能也是如此由于处理时间效率低下，我们提出了另一种估算方法，不需要每个聚类进行任何排序。算法中的伪代码在算法1中列出。在我们的建议中，我们首先以计算机为每个群集的平均强度比（算法1的第6行，图3-（d）），然后确定群集中其强度的像素数量比率低于或高于r（算法1的第7、9和10行）。如果不能将像素数量几乎均等地划分为漫反射和镜面反射组，或者，按顺序，如果其中一个组包含的像素多于该群集中像素总数的β分数（算法1的第11行和第13行），则按预先定义的方式进行更改步长γ（算法3的第12和14行），并重复算法直到达到最大迭代次数α（算法1的第8行）或更新后的r正确地划分了漫反射像素和镜面像素（算法3的第15行，图3-（e））。最后，我们为与选择强度比为r的聚类相关联的每个像素x分配Iratio（x）= r（算法3的第18-22行）。

一旦确定了每个聚类的强度比，就可以如下计算真实的镜面反射图像S（1）（图1-（j））。让我们从（1）和（2）中回想

S（x）= ws（x）Γ。

现在，从（4）和（6）可以看出，对于具有纯漫反射的像素

Imax（x）= wd（x）Λmax（x）

和

Imax（x）= Iratio（x）Irange（x） 。

因此，

wd（x）^ max（x）= Iratio（x）Irange（x）。

然后，以最大强度图像（4）计算S



（12）从（6），在（12）中S（x）将为零，因为

Imax（x）= Iratio（x）Irange（x）。

但是，如算法1中所列，在计算了（6）之后，Iratio的值发生了更改，使得每个群集使用一个强度比，使（6）成为近似方程式。最后，仅扩散图像D（图1-（i）） 可以通过输入图像和S之间的减法轻松计算。
4.基于GPU的DSP高亮度删除

在本节中，我们描述图1所示管道的GPU实现。默认情况下，我们假设内核中的每个线程都基于图像中相应像素执行其计算。此外，我们的GPU算法处理的所有图像都会创建并存储在内存中，并通过合并的内存访问来优化读写操作。首先，我们计算Imin（图1-（b）），Imax（图1-（c））和Irange （图1-（d））在GPU上的单个内核中。每个线程负责并行计算最小（3），最大（4）和范围（5）强度。接下来，我们估计Imin（Imin的平均值），以便对相关像素和无关像素进行分类（11）。为此，我们在Imin上运行并行前缀和[19]，以便求和Imin的最小强度值，并估计其平均最小强度。然后，在另一个内核中，我们计算Ipsf（9）（图1-（e））， pspsf（10）并在单独的图像中存储最小值Λpsfmin（图1-（f））和最大值Λpsfmax（图1-（g））伪漫射色度值。对于此任务，每个线程并行处理其对应的像素。在最小-最大色度空间中对像素进行聚类之前，我们在GPU上从“ psfmin”和“ psfmax”中定位三个初始簇种子。我们使用并行前缀和的优化扫描结构[19]来恢复具有最低和最高最小值，以及最高最大伪漫射色度值的像素。如图2-（a）所示，它们是我们聚类算法的初始聚类种子。将它们上传到恒定的设备内存后，我们将每个相关像素与其最近的群集种子相关联，该操作与像素无关。因此，在单独的内核中，每个线程在最小-最大色度空间中投影其对应的像素，测量投影像素与三个现有群集种子的距离，并将最近的群集种子的标识符存储到像素。对于三个群集中的每个群集，我们利用扫描结构并行计算与每个群集关联的像素数量。接下来，我们运行并行前缀总和以求和这些像素的最小和最大弥散色度值，并更新簇种子的质心。将这些新的簇种子位置复制到常量存储器后，我们并行执行另一个分配步骤，以将每个相关像素与其更新的最近簇种子相关联（图2-（b））。如图3所示，有两种计算强度比率的方法，可将每个群集分为漫射像素和镜面像素。在这两种方法中，我们首先使用线程并行计算Iratio（6）。接下来，我们运行扫描运算符以并行计算与每个群集关联的像素数。对于基于排序的方法，我们将与同一簇种子关联的像素复制到一个单独的紧凑数组中。然后，我们运行并行基数排序[19]来对强度比进行排序（图3-（b）），并为Iratio（x）的每个像素x复制所选择的强度比。非分类（a）原始（b）地面真相（c）OursPNSR37.5（d）Shen和Zheng [8] PNSR37.3（e）Yanget等人[14] PNSR37.2（f）Shen和Cai [13]算法1中列出的基于PNSR34.8的方法也可以很容易地并行化。对于每个聚类，我们运行并行前缀总和以估计与相同聚类种子关联的像素的平均强度比（算法1的第6-7行）。使用scan运算符，我们计算强度比低于r的像素（算法1的第9行）。在CPU上的单线程中，我们控制迭代和更新overr（算法1的8,11-16行）。接下来，我们并行地为Iratio（x）的每个像素x复制选择的强度比（算法1第18-19行）。最后，每个GPU线程利用每个簇重新估计的强度比来并行计算D和S的每个像素（图根据等式（12）的1-（i）），（图1-（j）），即像素独立且易于并行化。
5.结果和讨论
在本节中，我们通过文献中提出的独特的镜面高光去除技术来评估视觉质量和处理时间。我们将评估限于技术[5]，[9]，这些技术在Shen和Zheng [8]提出的标准数据集上报告了他们的结果，或者其源代码是开放且免费的[8]，[11] – [14]，[ 20]，或由各自的作者共享[4]。

A.实验设置

为了衡量我们的方法以及相关工作提出的方法的处理时间，我们使用了NVIDIA GeForce GTX Titan X图形卡和Intel CoreTMi7-3770K CPU（3.50 GHz），8GB RAM。 OpenCV 2.3.1 [21]和CUDA 8.0 [22]分别用于图像和并行处理任务。此外，我们在Thrust库中使用了并行前缀和和基数排序的优化实现

1.我们的最终源代码是开放的，可以免费获得

2，在准确性方面，我们将基于排序的（SB）（图3-（b，c））和基于非排序的（NSB）（图3-（d，e））的变体进行了比较，关于相关工作的测量，以测量标准镜面高光去除数据集中可用的四个图像上的峰值信噪比（PSNR）。
在性能方面，我们针对算法480p，720p，1080p和2160p的典型图像分辨率，测量了我们的方法的CPU和GPU版本以及相关工作的CPU实现所获得的平均处理时间。 ，我们凭经验选择了α= 3以保持算法的实时性能，β= 0.55和γ= 0.025，以使算法能够更快地收敛到准确的解决方案。如本节其余部分所示，这些参数足以为我们的测试提供准确的结果。
B.准确性评估

在本节中，我们定量和定性地评估了我们的方法和相关工作获得的准确性。在表I中，我们比较了我们的方法和相关工作对于标准数据集中使用的四个图像的准确性[8]。由于我们无法获得J.Yanget等人[23]，Suoet等人[5]和Renet等人[9]的技术源代码，因此，我们仅显示了其原始论文中报告的准确性结果。对于标准数据集中可用的四个图像中的三个，建议的方法比最先进的技术更准确，这表明我们的聚类方案能够准确去除镜面高光。从表I中还可以看出，我们方法的排序和非排序变体都提供了几乎相同的精度，能够获得比相关工作更好的精度。


在图4、5和6中，我们提供了一种方法与最准确的三项相关工作之间的直观比较，但其源代码尚未公开的技术除外[5]，[9]，[15]。我们请读者参阅补充文件的第二部分，以对结果进行更完整的视觉分析。 Tanet等人[11]，Shenet等人[12]提出的技术。 Shen和Cai [13]减少了原始图像中出现的镜面反射高光，但是它们也改变了图像仅扩散区域的颜色强度。这种伪影主要在图4-（f）的绿色特写中可见，鸭子的红色变暗。这是因为Tan等人的迭代框架[11]。这项工作是基于在色彩和几何形状方面都不一致的无镜面图像，而其他工作[12]，[13]易受噪声伪影的影响。Yanget等人[14] Shen and Zheng [8]倾向于降低位于物体边界附近的明亮扩散区域的颜色强度。这个事实可以在图4-（d，e）的红色特写中看到。还有，诸如Shen andCai [13]，Yanget等人[14]的技术。和Liuet等人。[4]不能准确地减少在纹理对象中发现的镜面反射高光，如图5-（e，f）和图6-（d）所示。从这个意义上讲，我们的方法能够更好地处理纹理对象的镜面高光去除（图5-（c）和6-（c）），并保留在图像的仅扩散区域中发现的明亮强度（图4-（c ））。
C.处理时间评估
表II中，我们提供了我们的方法与相关工作之间针对不同图像分辨率的处理时间比较。 Liuet等人的工作[4]是最慢的一种，因为他们的算法使用了昂贵的优化步骤来去除镜面高光。 Tanet等人[11]和Shenet等人[12]的方法。使用大图像分辨率无法缩放的迭代算法，从而提供远非实时的性能。 Q. Yanget等人[14]提出的算法。为镜面反射高光去除提供了可伸缩性，但是需要两次到四个迭代的双边滤波才能收敛到一个好的解决方案，需要2秒钟以上的时间才能去除具有2160分辨率的图像的镜面反射高光。 Shen和Cai [13]利用连接的组件标记算法来检测输入图像中高光的最大区域，扩展该区域以检测周围的扩散像素，并解决最小二乘问题以去除镜面高光。这种简单的方法可实现低分辨率图像的实时性能，但不能扩展为高分辨率图像。如表II所示，用于镜面高光去除的最快技术是基于最小最大色度空间聚类的技术，即Shen和Zheng [8]和我们的方法。但是，我们表明，通过使用基于非排序的方法来估计将每个聚类的漫反射和镜面像素分开的单个强度比，我们可以将最快的相关工作获得的处理时间减少一半[8]，获得最快，最准确的基于CPU的镜面高光去除算法。但是，如表II所示，对于1080p和2160p的图像分辨率，我们基于CPU的算法不能实时运行。从这个意义上讲，我们表明基于GPU的实现可提供图像分辨率方面的可伸缩性，比Shen和Zheng [8]的工作至少快20倍，并且仅需24毫秒即可去除2160的镜面高光，从而实现实时低分辨率和高分辨率图像的性能。
6.结论与未来工作
在本文中，我们提出了一种实时的高品质算法，用于去除单个图像的镜面高光。使用简单而高效的像素聚类算法，并结合优化方法将漫反射和镜面像素分开，我们证明了我们的方法能够与现有技术相比，处理时间要快两倍，同时提供现有数据集中可用于镜面高光去除的四个图像中的三个图像的最准确结果。此外，借助基于GPU的实现，我们能够提供镜面高光去除的实时性能，即使对于2160p图像分辨率也是如此。从这个意义上讲，我们可以为实时应用（例如增强现实和实时图像处理）使用精确的镜面反射高光去除。我们认为镜面反射高光去除领域仍需要开发更大的镜面反射高光去除数据集，以实现更好的镜面高光去除数据集。普遍性评估

OpenCV and Python K-Means Color Clustering（图像像素聚类）
英文原文
translator：aaron-clark-aic
简单讲，聚类就是将一组数据按照相识度分割成几类，也就是归类算法
k-means算法是一种简单的迭代型聚类算法
Python opencv中的k-mean聚类算法有现成的实现
需要注意的是，对于图像聚类我们需要将图像矩阵转换为像素的列表
# load the image and show it
image_base = cv.imread(".././img/base.jpg")

# reshape the image to be a list of pixels
image = image_base.reshape((image_base.shape[0] * image_base.shape[1], 3))

调用聚类算法的实现，实现聚类
k = 5 #聚类的类别个数
iterations = 4 #并发数4
iteration = 200 #聚类最大循环次数

clt = KMeans(n_clusters = k, n_jobs = iterations, max_iter = iteration) 
clt.fit(image)


调用直方图算法的实现，统计出图像的聚类分布
def centroid_histogram(clt):
    # grab the number of different clusters and create a histogram
    # based on the number of pixels assigned to each cluster
    numLabels = np.arange(0, len(np.unique(clt.labels_)) + 1)
    (hist, _) = np.histogram(clt.labels_, bins=numLabels)

    # normalize the histogram, such that it sums to one
    hist = hist.astype("float")
    hist /= hist.sum()

    # return the histogram
    return hist


def plot_colors(hist, centroids):
    # initialize the bar chart representing the relative frequency
    # of each of the colors
    bar = np.zeros((50, 300, 3), dtype="uint8")
    startX = 0

    # loop over the percentage of each cluster and the color of
    # each cluster
    for (percent, color) in zip(hist, centroids):
        # plot the relative percentage of each cluster
        endX = startX + (percent * 300)
        cv.rectangle(bar, (int(startX), 0), (int(endX), 50),
                      color.astype("uint8").tolist(), -1)
        startX = endX

    # return the bar chart
    return bar

hist = centroid_histogram(clt)
bar = plot_colors(hist, clt.cluster_centers_)

完整代码

# import the necessary packages
from matplotlib import pyplot as plt
from sklearn.cluster import KMeans
import cv2 as cv
import numpy as np


def centroid_histogram(clt):
    # grab the number of different clusters and create a histogram
    # based on the number of pixels assigned to each cluster
    numLabels = np.arange(0, len(np.unique(clt.labels_)) + 1)
    (hist, _) = np.histogram(clt.labels_, bins=numLabels)

    # normalize the histogram, such that it sums to one
    hist = hist.astype("float")
    hist /= hist.sum()

    # return the histogram
    return hist


def plot_colors(hist, centroids):
    # initialize the bar chart representing the relative frequency
    # of each of the colors
    bar = np.zeros((50, 300, 3), dtype="uint8")
    startX = 0

    # loop over the percentage of each cluster and the color of
    # each cluster
    for (percent, color) in zip(hist, centroids):
        # plot the relative percentage of each cluster
        endX = startX + (percent * 300)
        cv.rectangle(bar, (int(startX), 0), (int(endX), 50),
                      color.astype("uint8").tolist(), -1)
        startX = endX

    # return the bar chart
    return bar

# load the image and show it
image_base = cv.imread(".././img/base.jpg")

# reshape the image to be a list of pixels
image = image_base.reshape((image_base.shape[0] * image_base.shape[1], 3))

k = 5 #聚类的类别个数
iterations = 4 #并发数4
iteration = 200 #聚类最大循环次数

clt = KMeans(n_clusters = k, n_jobs = iterations, max_iter = iteration)
clt.fit(image)

hist = centroid_histogram(clt)
bar = plot_colors(hist, clt.cluster_centers_)

# show our color bart
fig = plt.figure()
ax = fig.add_subplot(211)
ax.imshow(image_base)
ax = fig.add_subplot(212)
ax.imshow(bar)
plt.show()

cv.waitKey(0)