本发明涉及计算机图像处理技术领域，具体地说是一种用于人体实时非刚性三维重建的模型顶点的并行处理方法、系统、求解器、芯片和介质。
背景技术：
人体的实时非刚性三维重建是计算机视觉、计算机图形学的重点问题。此技术广泛用于视觉特效、虚拟现实、远程教育和全息通信等多个领域，具有很高的商用价值。
非刚性实时三维重建一般采取彩色图和深度图信息作为综合输入，通过优化的方法同时重建模型并追踪运动情况。已有的工作常见通过深度传感器如kinect直接获取彩色和深度信息，或者用双目彩色相机在采集彩色信息的基础上计算得到深度信息。不同的深度信息获取途径使得深度信息的质量有差异，这要求求解器能够容忍输入深度信息的瑕疵。
由于非刚性实时三位重建的复杂性，目前还没有通用的商用级别的产品，系统还没有一个完全统一的框架，但是已有的非刚性实时三维重建系统框架的一个共同特点是包含追踪和融合两部分。融合部分维护并更新待建三维模型的几何，追踪部分对模型进行变形从而与实时输入的深度信息进行关联和对齐，这种对应关系继而又被融合模块利用来更新模型的几何。由此可见求解器作为求解追踪部分的计算引擎，直接决定了追踪部分对齐的效率和效果，进而决定了融合模块的效果，求解器是非刚性实时三维重建系统的重要核心。
然而目前从计算硬件到软件都没有一个通用框架来解决实时非刚性的三维重建，同时由于计算的复杂程度与实时的苛刻要求，整个重建的框架需要结合具体的计算硬件环境进行具体的设计，而作为重要核心的求解器占了绝大部分的计算量。仅仅调用通用的函数以及接口，例如英伟达显卡配套的cuda框架提供的计算函数，远远无法满足实时性的要求。
技术实现要素：
本发明为解决现有的问题，旨在提供一种模型顶点的并行处理方法、系统、求解器、芯片和介质。
为了达到上述目的，本发明采用的技术方案提供一种模型顶点的并行处理方法，包括：
s1，将输入矩阵j的分块矩阵设为母矩阵，表示为行数组、列数组和值数组；s2，获得列数组的非零列，对其元素排列组合后得出相乘后位于结果矩阵jtj中的行列位置；
s3，将元素与其行列位置所对应的两个分块矩阵进行并行计算，并排序去重；即得结果矩阵jtj的行数组和列数组的信息。
进一步地，s1中，值数组按照行优先的顺序依次放入母矩阵非零块；行数组的第j个数母矩阵第j行的首个非零块在值数组中的索引；列数组中存储母矩阵中每一行中非零块的列索引。
进一步地，s2中，对于输入矩阵j的第i行，首先从列数组得到所有非零列ci＝{c0,c1,...,ck}，对ci中的元素做排列组合，计算出相乘后位于结果矩阵jtj的行列位置：di＝{(c0,c0),(c0,c1),...,(ck,ck)}。
进一步地，s3中，di中每一个元素对应两个在输入矩阵j中用来计算的两个分块矩阵的位置，将所有行计算得到d＝{d0,d1,...dm}。
进一步地，s3中，排序方法为并行基数排序。
进一步地，s3还包括对若干结果矩阵jtj进行加法计算：输入矩阵j的非零块直接对结果矩阵jtj的非零块进行贡献，并键值排序去重。
进一步地，s3还包括转置：对d排列后，并行地将输入矩阵j中非零块复制到结果矩阵jtj的值数组。
进一步地，还包括s4，通过预条件共轭梯度法求解jtjδω＝-jtf(ω)，并得到变形参数的更新。
本发明还提供一种并行处理系统，包括存储模块，用于将输入矩阵j的分块矩阵设为母矩阵，存储为行数组、列数组和值数组；即针对非刚性实时三维重建设计的快速的块稀疏矩阵存储。
运算模块，用于计算块稀疏矩阵的行数组和列数组，以及计算每一块的具体数值并填充值数组，即执行块系数矩阵的快速代数运算。
求解模块，用于求解半正定稀疏线性方程以得到变形参数的更新，可以针对英伟达gpu和cuda框架的快速稀疏线性系统求解。
本发明还提供一种求解器，包括处理器、以及用于存储处理器的可执行指令的存储器，所述处理器运行时执行上述的并行处理方法。
本发明还提供一种芯片，包括处理器，用于从存储器中调用并运行计算机程序，使得安装有所述芯片的设备执行上述的并行处理方法。
本发明还提供一种计算机可读介质，其上存储有计算机程序指令，所述计算机程序指令被处理执行时，实现上述的并行处理方法。
和现有技术相比，本发明实现了一种用于人体实时非刚性三维重建的快速的并行处理方法，其中针对块结构的稀疏矩阵数据结构利用三维重建中数据的局部性，提高了数据的访问效率，在此基础上开发的矩阵运算将传统矩阵计算转化为排序加局部运算，提高了并行度，适合gpu快速计算。本发明所提出的线性系统减少调用gpu核函数的次数，提高了运行效率。
附图说明
为了更清楚地说明本发明实施例的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1是本发明实施例采用的实时三维重建装置的示意图；
图2是本发明在实时三维重建系统整体架构中的逻辑结构图；
图3是块稀疏矩阵存储的示意图；
图4是块稀疏矩阵的jtj的乘法示意图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图对本发明实施例中的技术方案进行清楚、完整地描述；显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。
参见图1和图2，本实施例用于非刚性实时三维重建的追踪模块，给定当前的模型几何以及非刚性变形参数ω，在输入新的深度和颜色信息时，通过给定的代价函数e来根据当前误差值求解模型的非刚性变形参数。代价函数定义在模型的顶点v上,表达为fv(ω)则可以得到代价函数的表达：
将所有顶点的误差排成列向量f(ω)我们可以得到向量表达：
e(ω)＝1/2f(ω)tf(ω)，通过对f求关于ω的雅可比矩阵
jtjδω＝-jtf(ω)
求解这个方程就能得到变形参数的更新δω。
上述理论是目前所有实时非刚性三维重建求解变形参数和追踪模块的理论基础。其中，雅可比矩阵j和残差f(ω)是输入，而变形参数的更新δω为输出。本发明主要针对非刚性实时三维重建设计j的数据存储，稀疏矩阵的代数运算以及稀疏线性系统的求解。
本发明提出的求解器用于非刚性实时三维重建，由于实时性的苛刻要求，求解器的设计与计算硬件的配置本身关系密切。由于模型各个顶点之间的独立性，很多运算可以并行处理，因此求解器的并行程度是提高计算速度的重要因素。求解器绝大部分的运算都在英伟达显卡(gpu)上完成，对gpu的操作通过英伟达官方的cuda库完成。本实施例的所有测试和部署在1080ti型号显卡上完成，因此对显存大小和计算兼容性(computecompatibility)的要求需要不低于1080ti。
在实时非刚性三维重建中，雅可比矩阵j和jtj都具有稀疏度高的性质。如果使用普通的密集存储形式，假设模型的顶点数为10万且每个顶点的误差函数为3维，控制点数量为1000，由于参数化形式的区别，每个控制点对应的变形参数为6至12个，那么一个雅可比矩阵的大小将大于1.8gb。这对于显存来讲是不可承受的存储压力，何况实际上顶点数目多于10万是常见的情况，因此引入稀疏矩阵来存储j和jtj是必然的。而结合实时三维重建的实际情况，j和jtj具有明显的分块结构：每个控制点的变形参数有6个，在计算雅克比矩阵时，从一个顶点p利用链式法则求导数求到对应的控制点c，最终雅克比矩阵j对应该顶点p和该控制点c确定的雅克比数值组成了一个3x6的块，计算得到jtj的对应位置则是6x6的块，并且在计算的时候往往每块的数据同时使用。我们提出保持数据局部性的分块稀疏矩阵格式，首先将每个分块矩阵当作一个整体元素，写成一个母矩阵，对这个母矩阵采取压缩稀疏行格式(compressedsparserowformat，csr)表示，即使用三个数组来表示母矩阵，值数组按照行优先的顺序依次放入母矩阵非零块，行数组的第j个数母矩阵第j行的首个非零块在值数组中的索引，而列数组中存储母矩阵中每一行中非零块的列索引。注意这里值数组存储的母矩阵每一个非零块实际上是将块矩阵继续按照行优先顺序展开成数组得到的，但是行数组和列数组中对于每个块矩阵是按照一个整体去计数及考量。
由于上述稀疏矩阵存储格式为自行定义，cuda无法支持所设定的格式的矩阵运算，则实现的块稀疏矩阵的每一种运算可以分为两个部分：第一步是计算块稀疏矩阵的行数组和列数组，第二步是计算每一块的具体数值并填充值数组。与传统的矩阵计算定义，本实施例的块稀疏矩阵的快速计算主要借助排序完成，而gpu具有非常快的并行排序算法，这保证大型稀疏矩阵的代数运算得以快速完成。本实施例中由输入矩阵j计算结果矩阵jtj是一个特殊的矩阵乘法，理论基础是使用jt的每一列乘j的每一行然后相加。参见图3和图4，具体实现的做法是对于j的第i行，首先从列数组得到所有非零列ci＝{c0,c1,...,ck}，由此对ci中的元素做排列组合能计算出相乘后位于jtj中的行列位置di＝{(c0,c0),(c0,c1),...,(ck,ck)}，di中每一个元素对应着两个j中用来计算的两个块矩阵的位置，将所有行计算得到d＝{d0,d1,...dm}作为键值进行排序并去重，就得到jtj中所有非零块的位置以及每个非零块由j中哪些块计算得来，也就得到了jtj的稀疏结构即行数组和列数组的信息，并能够对于所有非零块并行计算最终的乘积结果。其中计算d＝{d0,d1,...dm}可以对矩阵不同行进行并行计算，排序由并行基数排序(radixsort)完成，即所有步骤均有很大的并行度，保证了算法的高效。
作为优选，本实施例的方法对于不同误差来源的jtj之间可能需要加法计算。与乘法相比加法更加简单，输入矩阵的非零块直接对结果矩阵的非零块进行贡献，因此加法运算中针对结果矩阵的行列位置d＝{d0，d1，...dm}是恒等映射，同样通过键值排序去重，就能计算得到加法结果。
作为优选，不同条件下计算得到的雅可比矩阵会出现需要转置的情形。类似地，对于输入矩阵的每个非零块，经过计算到结果矩阵对应的行列位置d＝{d0，d1，...dm}，不同于乘法和加法，转置操作只需要将单个输入矩阵对应的非零块复制到相应位置即可，因此对d排列后无需去重操作，直接并行地将输入矩阵中非零块复制到结果矩阵的值数组即可完成转置操作。
为了求解jtjδω＝-jtf(ω)得到变形参数的更新，本实施例需要求解一系列大型的半正定稀疏线性方程，其理论框架是预条件共轭梯度法(preconditionedconjugategradientmethod，pcg)。按照其标准计算框架，涉及的代数运算主要包含稀疏矩阵与向量的乘法以及向量之间的点乘，其算法流程如下：
r0:＝b-ax0
p0：＝z0
k:＝0
repeat
xk 1：＝xk αkpk
rk 1:＝rk-αkapk
ifrk 1issufficientlysmallthenexitloopendif
pk 1:＝zk 1 βkpk
k:＝k 1
endrepeat
theresultisxk 1
按照标准框架用cuda框架实现需要调用9次核函数，其中向量乘法需要一次kernel调用，例如
例如，我们根据方法中涉及的中间变量之间的迭代更新依赖关系，将调用核函数的次数减少到4次，即将apk和
除此之外，本发明还去除了cpu与gpu之间的数据传输成本，并按照稀疏矩阵的定义按照矩阵和向量分块乘法的法则实现了矩阵向量的乘法，提升了pcg方法求解线性系统的效率。
本实施例的方法的特点在于：
能够实时求解基于颜色和深度输入的实时非刚性三维重建中的追踪问题；
适用于多种深度信息采集方式作为输入；
与重建系统其他软件模块解耦，在计算硬件不变的前提下，能够单独析出或者接入到其他任务中提供求解功能。
本发明实施例还提供一种并行处理系统，包括存储模块，用于将输入矩阵j的分块矩阵设为母矩阵，存储为行数组、列数组和值数组；即针对非刚性实时三维重建设计的快速的块稀疏矩阵存储。
运算模块，用于计算块稀疏矩阵的行数组和列数组，以及计算每一块的具体数值并填充值数组；即执行块系数矩阵的快速代数运算。
求解模块，用于求解半正定稀疏线性方程以得到变形参数的更新，可以针对英伟达gpu和cuda框架的快速稀疏线性系统求解。
本发明还提供一种求解器，包括处理器、以及用于存储处理器的可执行指令的存储器，所述处理器运行时执行上述的并行处理方法。
例如，存储器可以包括随机存储器、闪存、只读存储器、可编程只读存储器、非易失性存储器或寄存器等。处理器可以是中央处理器(centralprocessingunit，cpu)等。或者是图像处理器(graphicprocessingunit,gpu)存储器可以存储可执行指令。处理器可以执行在存储器中存储的可执行指令，从而实现本文描述的各个过程。
可以理解，本实施例中的存储器可以是易失性存储器或非易失性存储器，或可包括易失性和非易失性存储器两者。其中，非易失性存储器可以是rom(read-onlymemory，只读存储器)、prom(programmablerom，可编程只读存储器)、eprom(erasableprom，可擦除可编程只读存储器)、eeprom(electricallyeprom，电可擦除可编程只读存储器)或闪存。易失性存储器可以是ram(randomaccessmemory，随机存取存储器)，其用作外部高速缓存。通过示例性但不是限制性说明，许多形式的ram可用，例如sram(staticram，静态随机存取存储器)、dram(dynamicram，动态随机存取存储器)、sdram(synchronousdram，同步动态随机存取存储器)、ddrsdram(doubledataratesdram，双倍数据速率同步动态随机存取存储器)、esdram(enhancedsdram，增强型同步动态随机存取存储器)、sldram(synchlinkdram，同步连接动态随机存取存储器)和drram(directrambusram，直接内存总线随机存取存储器)。本文描述的存储器旨在包括但不限于这些和任意其它适合类型的存储器。
在一些实施方式中，存储器存储了如下的元素，升级包、可执行单元或者数据结构，或者他们的子集，或者他们的扩展集：操作系统和应用程序。
其中，操作系统，包含各种系统程序，例如框架层、核心库层、驱动层等，用于实现各种基础业务以及处理基于硬件的任务。应用程序，包含各种应用程序，用于实现各种应用业务。实现本发明实施例方法的程序可以包含在应用程序中。
在本发明实施例中，处理器通过调用存储器存储的程序或指令，具体的，可以是应用程序中存储的程序或指令，处理器用于执行上述方法步骤。
本发明实施例还提供一种芯片，用于执行上述的方法。具体地，该芯片包括：处理器，用于从存储器中调用并运行计算机程序，使得安装有该芯片的设备用于执行上述方法。
本发明还提供一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序，计算机程序被处理器执行时实现本发明上述的方法步骤。
例如，机器可读存储介质可以包括但不限于各种已知和未知类型的非易失性存储器。
本发明实施例还提供一种计算机程序产品，包括计算机程序指令，该计算机程序指令使得计算机执行上述方法。
本领域技术人员可以明白的是，结合本文中所公开的实施例描述的各示例的单元及算法步骤能够以电子硬件、或者软件和电子硬件的结合来实现。这些功能是以硬件还是软件方式来实现，取决于技术方案的特定应用和设计约束条件。本领域技术人员可以针对每个特定的应用，使用不同的方式来实现所描述的功能，但是这种实现并不应认为超出本申请的范围。
在本申请实施例中，所公开的系统、电子设备和方法可以通过其它方式来实现。例如，单元的划分仅仅为一种逻辑功能划分，在实际实现时还可以有另外的划分方式。例如，多个单元或组件可以进行组合或者可以集成到另一个系统中。
另外，各个单元之间的耦合可以是直接耦合或间接耦合。另外，在本申请实施例中的各功能单元可以集成在一个处理单元中，也可以是单独的物理存在等等。
应理解，在本申请的各种实施例中，各过程的序号的大小并不意味着执行顺序的先后，各过程的执行顺序应以其功能和内在逻辑确定，而不应对本申请的实施例的实施过程构成任何限定。
所述功能如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在机器可读存储介质中。因此，本申请的技术方案可以以软件产品的形式来体现，该软件产品可以存储在机器可读存储介质中，其可以包括若干指令用以使得电子设备执行本申请实施例所描述的技术方案的全部或部分过程。上述存储介质可以包括rom、ram、可移动盘、硬盘、磁盘或者光盘等各种可以存储程序代码的介质。
上面结合附图及实施例描述了本发明的实施方式，实施例给出的并不构成对本发明的限制，本领域内熟练的技术人员可依据需要做出调整，在所附权利要求的范围内做出各种变形或修改均在保护范围内。
1.一种模型顶点的并行处理方法，其特征在于包括：
s1，将输入矩阵j的分块矩阵设为母矩阵，表示为行数组、列数组和值数组；
s2，获得列数组的非零列，对其元素排列组合后得出相乘后位于结果矩阵jtj中的行列位置；
s3，将元素与其行列位置所对应的两个分块矩阵进行并行计算，并排序去重；
即得结果矩阵jtj的行数组和列数组的信息。
2.根据权利要求1所述的模型顶点的并行处理方法，其特征在于：s1中，值数组按照行优先的顺序依次放入母矩阵非零块；行数组的第j个数母矩阵第j行的首个非零块在值数组中的索引；列数组中存储母矩阵中每一行中非零块的列索引。
3.根据权利要求2所述的模型顶点的并行处理方法，其特征在于：s2中，对于输入矩阵j的第i行，首先从列数组得到所有非零列ci＝{c0,c1,...,ck}，对ci中的元素做排列组合，计算出相乘后位于结果矩阵jtj的行列位置：
di＝{(c0,c0),(c0,c1),...,(ck,ck)}。
4.根据权利要求3所述的模型顶点的并行处理方法，其特征在于：s3中，di中每一个元素对应两个在输入矩阵j中用来计算的两个分块矩阵的位置，将所有行计算得到d＝{d0,d1,...dm}。
5.根据权利要求1所述的模型顶点的并行处理方法，其特征在于：s3中，排序方法为并行基数排序。
6.根据权利要求1-5任一所述的模型顶点的并行处理方法，其特征在于：s3还包括对若干结果矩阵jtj进行加法计算：输入矩阵j的非零块直接对结果矩阵jtj的非零块进行贡献，并键值排序去重。
7.根据权利要求1-5任一所述的模型顶点的并行处理方法，其特征在于：s3还包括转置：对d排列后，并行地将输入矩阵j中非零块复制到结果矩阵jtj的值数组。
8.根据权利要求1-5任一所述的模型顶点的并行处理方法，其特征在于：还包括s4，通过预条件共轭梯度法求解jtjδω＝-jtf(ω)，并得到变形参数的更新。
9.一种并行处理系统，其特征在于：
包括存储模块，用于将输入矩阵j的分块矩阵设为母矩阵，存储为行数组、列数组和值数组；
运算模块，用于计算块稀疏矩阵的行数组和列数组，以及计算每一块的具体数值并填充值数组；
求解模块，用于求解半正定稀疏线性方程以得到变形参数的更新。
10.一种计算机可读介质，其特征在于：其上存储有计算机程序指令，所述计算机程序指令被处理执行时，实现权利要求1-8中任一所述的并行处理方法。
技术总结