分类目录归档:技术

全景视频技术的产品化之路

版权声明:原创作品,欢迎转载,但转载请以超链接形式注明文章来源(planckscale.info)、作者信息和本声明,否则将追究法律责任。

甚嚣尘上的VR炒作终于在今年平静了,这大概意味着VR技术开始进入技术成熟度曲线的第三个时期:行业的公众关注进入低谷,人们开始冷静客观评估技术的适用范围和潜力,并逐步发现有效的经营模式。

技术成熟度曲线

VR时代的到来是不可避免的,或者说它已经到来,只是还没有推到大众面前。另外,真正具有想象力和冲击力的新技术乃是紧随其后的AR,这一点可能并不像公众预期的那样。这个时代需要由一系列扎实漂亮的产品撑起(不是概念,不是Demo,不是DIY,是产品),我们这次来谈谈全景摄像机的产品化之路上有哪些曲折和挑战。当然,全景摄像机本身并非仅限于VR应用,我们也要包括安防应用。

安防监控领域

泛泛来说有两种全景视频实时拼接方案,即前端(机内)拼接后端(PC/手机)拼接。在安防领域也是如此。前端拼接直接由全景摄像机输出拼接完成的全景帧,具有很好的兼容性,可以直接像一台普通IPC一样接入旧有系统;而后端拼接是将全景摄像机看做独立的多路IPC,同时接入监控PC服务器,由PC完成实时拼接和监看。后端拼接的优势在于可以完成极高分辨率(目前我们的后端方案全景监控分辨率最高已经有9600万像素)的全景监控,但兼容性不好,需要将全景拼接SDK嵌入平台软件,不能做到“即插即用”。

从实现上来说,大概有如下几种:FPGA/DSP/CUDA/OpenGL/CPU. 前两种用于前端拼接,FPGA的开发和维护都有较高代价,CUDA和OpenGL方案具有最高的处理能力,CPU方案除非无法选择否则是应该排除的。在前端拼接方案里,还要考虑编码问题,全景帧动辄数千万的分辨率编码并不是一个简单问题。这里我们主要谈我们自己比较熟悉的CUDA/OpenGL方案。

安防监控领域对于全景摄像机有一些特殊需求。对于后端拼接全景,其拼接参数应当保存在设备之中,由设备传给平台软件完成实时监看的初始化流程,而平台软件上则对实时拼接的效率,全景模块与其他设备如球机的互动都有颇多要求,我们简单罗列如下。

  1. 拼接参数应该是一个很小(几k到几十k)的文件,方便写入设备及在网络上传输;
  2. 灵活的裁剪/融合算法。安防全景细分需求繁多,催生大量不同类型的设备,不同目数,不同镜头,不同摆位都会导致非常不一样的应用场景,算法需要在任何场景下都能够完美融合,不产生图像瑕疵。
  3. 全景拼接模块应当支持实时的子码流/主码流切换。平台软件实时监看几十上百路网络摄像机,区分子码流/主码流非常重要,这样在小视图模式下采用子码流,而大视图下自动切换到主码流,既保证了性能又保证了操作体验。
  4. 全景拼接模块在子码流输入下能够同时完成多达几十路的实时拼接播放。
  5. 应该有多种全景投影模式。除了常见的球面/柱面展开,碗形交互式展开在安防监控领域颇受欢迎,如图所示。各种投影模式之间应该能够实时切换。
  6. ROI局域放大。
  7. 与球机联动,全景纵览全局,球机实现局域放大。
  8. 全景像素坐标到输入图像像素坐标的正反向投影。

碗形交互式全景

解决以上需求的任务并非平凡。简单一例,子/主码流实时切换中,除非子码流与主码流具有同样的视野,否则无法在不重新初始化算法的前提下完成切换,这要求算法具备瞬间初始化完毕的性能。同样的性能要求也出现在实时全景投影模式切换中。

从一个算法到一个成熟产品的道路是长远的。行业里很多学术型团队最终败在懂算法不懂软件工程,无法将一个Demo级的算法提升成一款结构良好,功能灵活,充分解决行业内需求的算法产品,令人扼腕。

对于前端拼接来说,要支持交互式全景类型如碗形、柱面等,同样也需要将一个全景播放模块嵌入平台软件,此种情况下,上述需求中除1、2、8外都仍然需要满足。

前端拼接技术的一个需要解决的问题是,以条带展开型全景作为全景帧类型做编码传输,如何节省带宽?将一个全景球展开为一个平面图像就如同将柚子皮拍扁在桌面上,总会像全球地图那样产生一个畸变,这是无法避免的。投影类型选择不好可能会导致相当大的畸变,比如球面两极地区一个像素被拉伸成一行像素。更好的投影类型应该是立方体展开

立方体展开

这一展开方式可以将畸变控制到很小的程度,但它一定程度上损失了条带全景图那种一览无余的直观性,需要特殊的全景播放器将它重新贴图到全景球或全景展开平面上才能够还原全景。当然也有其他采用更高级的数学方案设计的展开方式,这里不再提。

全景摄像机与全景直播

这里特意避免了提“VR全景”这一概念,因为严格来说VR全景和普通的全景摄像机并非同一概念,前者要求具有视觉深度感,后者只是个普通的2D曲面,沉浸感不强。但由于普通全景摄像机技术较前者简单,所以目前市面上大都为此类产品。

我的个人观点是,目前全景摄像机难以普及的一个关键是没有标准格式。并不像传统数码相机,全景输出格式杂乱无标准,全景视频播放器无法自动化决定采用何种投影类型播放,使得全景视频成了少数geek一族的玩具。但在真正的行业标准出现之前,让自己的产品对各种不同的输出格式都做好准备不失一个办法,而且不难。

这类产品中低端以前端拼的双鱼眼为主,高端以后端拼的多目摄像机为主,但迄今几乎没有很让人满意的产品出现。

双鱼眼方案的优势在于廉价且可以极小化拼缝。在所有可能的基于拼接算法的方案里,双鱼眼的拼缝是最小的。拼缝大小取决于多个摄像机投影中心的距离,摄像机的投影中心位置大致在sensor中心向后一个焦距远的地方,通常这是个很短的距离。理论上只有各个摄像机的投影中心重合于一点才能够产生出无缝的全景图,但这种情况下相机的体积需要压缩到极限,几乎不可达到,通常只能将尺寸压缩到极小以期更好的拼缝效果。除了双鱼眼方案,它是可以真的做到投影中心重合的。

所以,对于做基于拼接算法的全景摄像机的厂商,一个忠告是,将相机尺寸做小

全景直播机似乎有很长一段时间卡在很高的软件授权费和拼接服务器价格上,但这是比较奇怪的,因为这一技术并不困难——至少在安防领域,四年前就已经有公司做到了上千万像素的全景监控。像安防领域一样,最高性能且具有很好平台兼容性的方案就是OpenGL方案,现在的显卡处理几千万像素的全景拼接融合如同砍瓜切菜,顺便搞个硬编码做推流是不难的——我们自己的技术在这方面早已验证过。

全景直播机通常并不是多个摄像机拼一块儿这么简单粗暴,它需要解决两个基本问题,一是摄像机之间的帧同步,一是摄像机之间的成像参数同步。前者保证人通过拼缝时不会出现消失又出现这种诡异效果,后者使得全景画面亮度、色彩具有一致性,不出现尖锐的过渡。

但实际上,我们并不真的需要成像参数同步。理论上,多个摄像机各自自动曝光,可以实现HDR(高动态范围)全景,因而目前在硬件上做成像参数同步只是一个过渡方案,将来为了生成HDR效果全景,这一机制是必然要废弃的。

要有更好的拼接质量,可以选择CUDA或OpenCL,它比OpenGL提供更多控制力,使得开发者可以采用更复杂的图像处理算法。我们目前就在基于CUDA尝试HDR全景算法的开发。

VR全景

VR全景是万众期待众望所归,出于不可描述之原因,这一技术被视作新时代的宅男福利。但一定要冷静!因为我们真还有很多技术问题要解决。

实现3D效果,目前主要有基于传统拼接算法拼左右眼全景图(参见我们的文章《DIY 3D全景摄像机》)和光流算法(Google Jump/Facebook Surround360等)两种。

基于拼接算法基本是没有前途的(所以我们直接做成了DIY教程-_-!)。这一方案的死结在于,3D全景中深度感最强的近景,正好是拼缝最大的,而且你不能够通过缩小设备尺寸来解决,因为它至少应该有人的瞳距(~6.2cm)那么大,否则你戴上眼镜后,会发现自己缩小了——周边的一切都大了一遭。

光流算法是目前给出效果最好的,光流刻画了两个图像的像素是如何对应的,算法利用光流来插值计算没有被相机所采集的光线之颜色,从而产生出完全无缝的全景效果。但目前效率不高,关键是光流本身的计算是相当繁重的,而且算法对于每队图像还需要计算正反向两个光流,再考虑上光流在时间轴上的一致性,带来了非常大的计算开销。

实现VR视频采集,本质上是通过有限个相机采集几个点上的物理光线,然后用这些光线来猜测、插值出其他空间位置上任意光线。这在计算机视觉领域早已研究多年(想想黑客帝国里的子弹时间镜头是怎么来的),这个方向叫做”Image-Based Rendering”.

理想自然是通过采集有限个点上的光线就能够计算出一个邻域上的光场。这一定程度上做得到,而且有很好的工作,但付诸应用仍然有距离。

所以,仅就目前的情况来说,基于光流算法来做后期,做高质量近景VR视频是没问题的,但想要直播,还得等等。

DIY 3D全景摄像机

版权声明:原创作品,欢迎转载,但转载请以超链接形式注明文章来源(planckscale.info)、作者信息和本声明,否则将追究法律责任。

目前市面上的VR全景多是二维,没有深度感,若真想拥有身临其境般的体验,深度信息必不可少。诺基亚的OZO,Google的Jump,Facebook开源的Surround 360,都是为3D全景而设计。OZO设备8个鱼眼售价高达三十多万,Surround 360搭载的是Point Grey的相机,硬件成本二十多万,Jump也要搭载十几个GoPro,硬件成本少说也要几万,普通玩家真心想玩也要思量的。如果不追求那么高大上,其实自己就可以DIY出一台3D全景相机。

组件:

千兆交换机+网络摄像头模组+180°鱼眼镜头+线材

结构:

要实现深度感,结构是关键。

相机可以理解为光线采集设备,采集到的光线与成像平面的交点即像点。

通常的二维全景要求采集到的所有光线汇于中心点,即视点,以视点为中心的球面或圆柱面为成像面,所有光线交于成像面形成全景图,如下图(a)。二维全景相机要求所有相机共中心摆位,即所有相机的光轴相交于视点。

3D全景即左右眼各对应一个全景图。

两只眼睛分别对应两个不同的视点位置,当转头360度时,两只眼睛转过的轨迹即一个以瞳距为直径的圆,称之为Viewing Circle,3D全景要求采集到的所有光线相切于此圆。左右眼采集到的光线分别与成像面相交形成左眼全景图和右眼全景图, 如下图(b)(c)。

两只眼睛所在的视点位置投影出的图像称为一个立体对,如图(b)左眼光线1和图(c)右眼光线1即可看成一个立体对,同理左眼中的光线2,3,4,5等与右眼中的2,3,4,5等分别构成不同的立体对。

该图引自文献:Stereo Panorama with a Single Camera.

至于相机摆位一般有两种方案,如下图:

该图引自文献:Jump: Virtual Reality Video.

切向摆位如上图(a):每个相机的光轴相切于Viewing Circle,此种方案一半的相机用于左眼全景(图(a)中绿色相机),而另一半的相机用于右眼全景(图(a)中红色相机)。

径向摆位如上图(b): 每个相机的光轴沿Viewing Circle半径方向,此种方案不区分左右眼相机,每个相机都对两眼的全景图有贡献,因此对每个相机水平视场角有更高要求:R越小,要求每个相机的水平视场角越大。

一般R设计比较大时采用径向摆位,R较小时采用切向摆位。

切向摆位最简单的结构设计即正多边形,每条边上放置两个camera,其sensor中心的距离设为瞳距。如果镜头的视场角足够大,可以设计一个正三角形,用六个camera来实现3D全景。本文介绍的是正四边形八个camera的方案,用Solidworks设计一个简单的支架,预留出上camera的安装孔位。

结构设计及3D打印: 

组装:

效果:

原始视频截图

上下3D格式

VR眼镜观看3D效果

目前能实现3D全景的技术无外乎几种:

  • 拼接方案。左右眼的视野分别做拼接融合以达到3D全景的效果,这种方案最简单,可以实时化,其缺点是拼缝难消除。
  • 光流方案。Google Jump以及Facebook开源的算法即此方案,能很好的消除拼缝,但实时化比较困难,适合做后期处理。
  • 光场重建。用有限个相机重建光场,给出真正的3D效果,使用户拥有更多活动自由度,这是终极的VR视频采集方案,但即便在理论上也有很多困难之处。

本文中所用的实时3D全景拼接软件是奇点视觉实时全景拼接方案,演示中为2700万像素实时拼接融合效果。我们目前可以轻易的实现超高分辨率实时拼接和直播,但仍有以下几个问题:

  1. 硬件使用了安防用网络摄像头模组,不具备同步曝光功能,因此全景图明暗不均较严重;
  2. 近景拼缝明显。3D全景要求设备直径至少等于人眼瞳距,但过大的直径容易导致更严重的拼缝。这是基于拼接方案做3D全景的一个本质困难。

但至少到目前,拼接方案仍然是唯一能做到低开销高分辨率直播的。我们正在研发新一代的实时光流/光场重建算法,希望能解决高质量近景VR直播的问题。

奇点视觉是一个致力于计算机视觉技术的研发和产品化的团队,专注于算法,为了发挥我们的优势,没有去做产品,而是给有能力做产品的公司提供技术解决方案。过去两年,我们专注于安防全景,毫不夸张的说,安防全景技术我们做到了世界顶级水平,已经由客户厂商实现产品化并销往海内外。目前我们致力于把全景技术迁移到VR应用中,并作技术升级,希望能够为更高质量的VR内容贡献一份力量,敬请关注奇点视觉。

玩玩三维重建

版权声明:本文系本站作者自己翻译整理,欢迎转载,但转载请以超链接形式注明文章来源(planckscale.info)、作者信息和本声明,否则将追究法律责任。

我们在实时三维重建方面的工作今年已经密集展开。或许不久后某一天,你会在本站看到带有SLAM(即时定位与地图构建)功能的四轴飞行器,或者让你在书桌上打一场现代战争的增强现实应用。在敲锣打鼓欢天喜地亮出我们自己的三维重建实现前,先拿别人的东西给大家打打牙祭。

中科大刘利刚教授的3D建模软件与处理软件简介介绍了N多实用的3D相关软件。而基于照片的快速建模软件并不多,之前玩过123D Catch,很赞。围着你要建模的物体拍摄一圈,用123D Catch加载拍摄的图像,经过其强大的处理能力,生成具有纹理的3D模型。下图是我重建的我的蒙奇奇。你要做的只是拍照、上传、等待而已,相当简单。

蒙奇奇

但是123D Catch也存在一些局限,完全的黑盒子,对重建过程没有任何操控力。

要想了解从照片如何一步步重建出三维模型,并能操控某些过程,可用的免费开源软件也不少,较常用的是VisualSFM和Meshlab:

第一步:VisualSFM

VisualSFM软件允许我们上传一系列图像,它从这些图像中找到每一个图像的特定特征,利用这些特征信息重建出3D模型的稀疏点云,而后还可进行稠密点云重建。

输入:围着要重建对象拍摄的一系列照片;

输出:一个 .out文件,存储着每个相机的位置及重建出的稀疏点云;

一个.ply文件,存储着由稀疏点云重建出的稠密点云。

第二步:Meshlab

可用Meshlab对3D网格/点云做各种操作。输入VisualSFM的生成文件,Meshlab通过一系列操作可创建出包含纹理的、干净的、高分辨率的网格,并自动计算UV映射及创建纹理图像。

输入:VisualSFM的生成文件,.out文件和list.txt文件(存储照片序列); 以及.ply文件;

输出:一个.obj文件,3D模型的网格;

一个.png文件,任意大小的纹理图;

完整的流程见下图:

liucheng

 

第一步:运行VisualSFM

1

1. 输入一系列图片

拍照注意事项:切忌不要站在原地,仅转动身体去拍:相机共中心能拼接全景,但是给不出三维重建的深度信息。要以待重建的对象为中心,围着它每转10-20度拍一张,这样转一圈,有不同高度信息更好。VisualSFM没有照片数量限制,照片越多,重建出的细节越丰富,但重建过程花费时间越长。QQ图片20150314232349  

2.  特征检测及匹配

因照片可能存在旋转、缩放或亮度变化,此过程利用SIFT算法提取、描述特征,用 RANSAC算法过滤掉误匹配。此过程亦可利用GPU加速。工作状态实时显示在侧边的log窗口。

QQ图片20150314232955QQ图片20150314233141

3. 利用SFM进行稀疏3D重建

利用 SFM 方法,通过迭代求解出相机参数和三维点坐标。即重建出3D模型的稀疏点云。若有“bad”相机(位置错误或朝向错误),结合工具栏上的“3+”按钮和手型按钮即可删除之,使结果更准确。

QQ图片20150314233508

4. 利用  CMVS/PMVS 进行稠密3D重建

CMVS/PMVS需自己下载,编译,也可直接下载exe文件。而后把pmvs2.exe/cmvs.exe/genOption.exe文件放到VisualSFM.exe的同目录下。

通过 CMVS 对照片进行聚类,以减少稠密重建数据量,而后利用PMVS从3D模型的稀疏点云开始,在局部光度一致性和全局可见性约束下,经过匹配、扩散、过滤 生成带真实颜色的稠密点云。(下图为用Meshlab查看效果图)

6

 

至此,VisualSFM的工作告一段落,结果都已存盘。若因图片匹配失败或图片较少导致某区域重建失败或重建出的某区域细节不足,可以返回添加一些这个区域的照片,重新来过(本人较懒,未作补充,谅解)。但因特征检测和匹配的结果已存盘( 每张图像对应的.sift 和 .mat文件),所以已经匹配好的图像不必再次匹配,会更快完成。

第二步:运行Meshlab

11

1. 打开bundle.rd.out 文件

a. 按钮1,打开由 VisualSFM生成的存储在xx.nvm.cmvs文件夹下的 bundle.rd.out 文件。随后会询问是否选择照片列表文件,选择同文件夹下的 “list.txt”即可。这一步会把相机及对应的照片导入进来,对后续的纹理处理至关重要。

3

b. 按钮2,打开显示层目录,检测相机载入是否正确, Render –> Show Camera,因可视化相机的尺寸比网格尺寸大得多,所以需调整相机的缩放因子,scale factor可以从0.001开始调小,直到相机位置清晰可见。

45

 2. 稠密点云代替稀疏点云

a.  按钮3,隐藏可视的稀疏点云;

b. File –> Import Mesh加载稠密点云(xx/00/models/option-0000.ply);VisualSFM生成多个.ply文件时,需合并成一个mesh。在载入的任何一个.ply上右键选“Flatter Visible Layers”。

6

3. 清除杂点

按钮4选中杂点区,按钮5删除之。大致清了桌前的一些杂点。

QQ图片20150314234454

4. 网格化

Filter –> Point Set–> Surface Reconstruction: Poisson.

利用Poisson Surface Reconstruction算法由稠密点云生成多边形网格表面。

参数可调, Octree Depth:控制着网格的细节,此值越大细节越丰富但占内存越大运行起来慢,一般设10,可慢慢调大。

7

Poisson表面重建算法会生成一个“不漏水”气泡,把所有场景对象包裹在其中。即模型是封闭的。可以移除多余的面Filters –> Selection –> Select faces with edges longer than,而后删除。

QQ图片20150314223022 QQ图片20150314223134

保存(整个project和mesh)。

5. 修复流形边缘

后续的纹理处理要求网格化的模型必须是流形(MANIFOLD)的,因此需删除非流形边(简单讲就是任何由多面共享的边)。Filters –> Selection –> Select Non-Manifold edges,而后删除之。

QQ图片20150314234928

6. 参数化(Parameterization)

Filter –> Texture –> Parameterization from registered rasters。

根据相机投影关系创建UV映射。

QQ图片20150314235004

保存 (整个project和mesh)。

7.  投影纹理

Filter –> Texture –> Project active rasters color to current mesh, filling the texture。

可设置任意分辨率(512的2的二次方倍:512 / 1024 / 2048 / 4096 / 8192…)的纹理图。

QQ图片20150314235126

6和7可以合为一步: Filter –> Texturing –> Parameterization + texturing from registered rasters.

 QQ图片20150314235200

8. 完成、导出

当你调整满意了之后,File –> Save mesh as… a .obj文件。即可便有了一个包含你选定分辨率纹理的obj文件。

QQ图片20150314225720 QQ图片20150314225733

收官啦。而后关乎应用,就是拼想象的时候了!

更多细节参见:We Did Stuff 

CUDA, 软件抽象的幻影背后 之三

版权声明:原创作品,欢迎转载,但转载请以超链接形式注明文章来源(planckscale.info)、作者信息和本声明,否则将追究法律责任。

上一篇中谈到了编程模型中的Block等概念如何映射到硬件上执行,以及CUDA如何用并行来掩盖延迟。这一篇继续剖析SIMT,谈一谈控制流分叉,指令吞吐和线程间通讯机制。
虽然我们说warp中的线程类似于SIMD,但事实上它是真正的线程。warp中的每一个thread都有自己的指令地址寄存器,允许它们各自执行不同的任务(控制流分叉)。最简单的,比如一个

if(threadIdx &lt 10)
 {}
 else
 {}

语句,将threadIdx=0…31这一个warp划分成两个分支,各自做不同的事情。这个灵活性以性能为代价,当一个warp中控制流出现分叉时,不同分支的线程会被分组相继执行,直到各分支执行完毕后,控制流重新汇聚成一支(上例中即if语句的结束点)。这种情况下执行单元的利用率较低,因为每个分支执行时都需要关闭其他分支的线程,所以这时一些执行单元是用不到的。
为了尽可能高效的计算,需要约束控制流分叉的出现。除了减少流程控制语句外,还需要注意,并不是只要有流程控制语句就一定会带来控制流分叉。关键是,控制流分叉只是针对同一warp中的线程而言,不同warp的线程原本就是串行化执行的,分叉对其无影响。因此,只有流程控制语句的条件在
同一warp内不一致时,才会有控制流分叉。这样,诸如

if(threadIdx.x / WARPSIZE &lt n)
{...}
else
{...}

这样的语句是不会有分叉的。当然,更宽松的条件如

if(blockIdx.x &lt n)
{...}
else
{...}

也不会有分叉。依赖于输入数据的条件如

if(globalArray[threadIdx.x] &lt n)
{...}
else
{...}

则会带来分叉。

对于简单的指令如32位浮点数的加、乘,32位整数的加减等,通常CUDA Core在一个时钟周期内可以完成一次操作,因而一个周期内完成的操作数目就等于CUDA Cores数目;而对于一些较复杂的指令,如sin/cos等超越函数,执行单元并不能提供这么高的吞吐率。
我们可以用单位周期内进行的操作数目N除以32来计算指令的吞吐率。以GM204为例,它的SM中有32*4 = 128个CUDA Cores,32个SFU(特殊函数单元),在计算32位浮点加法时具有最高吞吐,一个周期内完成128次操作,单位周期内指令吞吐为128/32 = 4;而计算如sin/cos等超越函数时线程不再一一分配到CUDA Cores上,而是要在32个SFU上计算,单位周期内只能完成32次操作,指令吞吐为1条指令每周期.
指令的吞吐率数据可参考CUDA C Programming Guide中 5.4.1. Arithmetic Instructions,该小节以单位时钟周期每SM上能够进行的操作数的形式给出了各指令的吞吐率。
指令吞吐率是我们进行性能优化的有一个重要指标。通常,影响指令吞吐率的因素除了数值计算操作的复杂度、精确度之外,控制流分叉也是一个贡献因子。这里的原因不难理解,控制流分叉时执行单元的利用率下降,使得单位周期内执行的操作数目下降,从而降低了指令吞吐。

到这里,硬件图景下线程的执行就基本说完了,只剩下一个留到最后的话题:线程间交互。通常,不存在任何相互作用的线程,它们之间才能够以任意的顺序执行,像block。但对于warp这样的线程组,是可能与同一block中其他warp通讯或同步的,这时执行顺序就不能任意。所幸即便在block之内,线程间的交互仍然是较弱的,因而底层可以将block划分成warp来分组串行化执行,遇到交互时再另作处理。我们现在来看看这些交互机制。

线程间交互可以细分为通讯和同步两类。通讯主要由公共存储区域交换数据来实现,但也不排除像shuffle这样的特殊方式存在。
从通讯的粒度来看,可以分为warp内部线程间通讯,block内部线程间通讯,block间通讯,更粗的粒度这里不考虑。block之间的通讯则只能基于global memory,block内部的通讯主要基于shared memory/global memory,warp内部线程间除了可以利用上述所有方式,还有一种特殊的shuffle机制.下面我们以通讯的粒度分类陈述各种通讯的实现方式。

block间通讯通常基于两次kernel发射,一次将通讯数据写入global memory,另一次发射读global memory进行后续处理。这种通讯开销较大,主要来自于global memory访存和kernel发射,所以如果有可能,尽量把任务放在一次kernel发射中完成。
或许有人会问,同一个kernel发射中的两个block具有共同的global memory,是不是也可以利用这个特点来构造同一kernel下block间的通讯呢?通常的答案是no,因为block之间执行顺序不定,很难构造有意义的通讯;但如果要较真,答案是yes,我们真的可以构造一些特殊的block间通讯方式。一个例子如下所示,该实例来自于CUDA C Programming Guide B.5. Memory Fence Functions:

__device__ unsigned int count = 0;
__shared__ bool isLastBlockDone;
__global__ void sum(const float* array, unsigned int N,
	volatile float* result)
{
	// Each block sums a subset of the input array.
	float partialSum = calculatePartialSum(array, N);
	if (threadIdx.x == 0) {
		// Thread 0 of each block stores the partial sum
		// to global memory. The compiler will use
		// a store operation that bypasses the L1 cache
		// since the "result" variable is declared as
		// volatile. This ensures that the threads of
		// the last block will read the correct partial
		// sums computed by all other blocks.
		result[blockIdx.x] = partialSum;
		// Thread 0 makes sure that the incrementation
		// of the "count" variable is only performed after
		// the partial sum has been written to global memory.
		__threadfence();
		// Thread 0 signals that it is done.
		unsigned int value = atomicInc(& count, gridDim.x);
		// Thread 0 determines if its block is the last
		// block to be done.
		isLastBlockDone = (value == (gridDim.x - 1));
	}
	// Synchronize to make sure that each thread reads
	// the correct value of isLastBlockDone.
	__syncthreads();
	if (isLastBlockDone) {
		// The last block sums the partial sums
		// stored in result[0 .. gridDim.x-1]
		float totalSum = calculateTotalSum(result);
		if (threadIdx.x == 0) {
			// Thread 0 of last block stores the total sum
			// to global memory and resets the count
			// varialble, so that the next kernel call
			// works properly.
			result[0] = totalSum;
			count = 0;
		}
	}
}

代码 1. block间通讯实现数组求和
本代码摘录自 CUDA C Programming Guide B.5. Memory Fence Functions

该例实现一个数组的求和,首先各个block计算部分和,然后由最后一个完成部分和计算的block再把所有的部分和加和出最终结果。block间通过一个位于global memory的变量count通讯,它记录了目前已经完成计算的线程数。这样,最后一个完成部分和计算的block就会发现count的数值为最大线
程id,因此可以判定需要由它自己来完成最后从部分和向总和的计算。
不过,为了更好的软件结构,最好还是避免同一kernel的block间产生耦合。同一kernel中block的通讯还涉及到CUDA的weakly-ordered内存模型问题,一个线程中先后两次内存操作在另一个线程看来未必能够保持原有顺序,这产生了相当大的复杂性。我们在下文还会提到这一问题。

block内的线程通讯机制较为丰富,尤其是线程同属一个warp时的shuffle机制。shuffle在Kepler后出现,是一种相当快的线程间通讯方式,它允许同属一个warp的线程间可以互相引用彼此的寄存器,比如下例:

__global__ void bcast(int arg)
{
	int laneId = threadIdx.x & 0x1f;
	int value;
	if (laneId == 0) // Note unused variable for
		value = arg; // all threads except lane 0
	value = __shfl(value, 0); // Get "value" from lane 0
	if (value != arg)
		printf("Thread %d failed.\n", threadIdx.x);
}

代码 2. shuffle机制实现一个值向整个warp的广播
本代码摘录自 CUDA C Programming Guide B.14. Warp Shuffle Functions

laneId是warp中线程的一个index,有threadIdx对32取余得到。__shfl(value, 0)语句使得各线程能够访问laneId==0这一线程中value的值。

更常用的通讯机制自然是shared memory和global memory了。其中shared memory更快速,在大多数时候是构建高性能CUDA程序的必由之路。这些常识不再赘述。基于shared/global memory的线程间数据交换,一定要注意线程的同步。block中线程的同步由__syncthreads()实现。线程会等待同block中其他线程都执行到这一点,并且__syncthreads()语句之前的所有shared/global memory操作都尘埃落定,保证block内所有线程在__syncthreads()之后都能看到这些操作的结果。

最后谈一下CUDA采用的weakly-ordered内存模型。它导致一个线程中相继执行的两个存储器操作在另一个线程看来未必是一样的顺序。例如:

__device__ int X = 1, Y = 2;
//thread 0
__device__ void writeXY()
{
	X = 10;
	Y = 20;
}



//thread 1
__device__ void readXY()
{
	int B = Y;
	int A = X;
}

代码 3. weakly-ordered内存模型示例
本代码摘录自 CUDA C Programming Guide B.5. Memory Fence Functions

这段代码可能产生A=1,B=20这样的结果。原因是有多种可能的,要么thread 1看到的X、Y的写入顺序被颠倒,要么thread 1中读取顺序被颠倒。这种看似相当毁三观的事情确确实实发生在我们的代码背后。在一个线程里两个相继但无依赖的内存操作,其实际完成的顺序可能是不确定的。在这个线程
看来这并没有导致什么不同,因为两个操作无依赖,并不会破坏因果链;但在另一个线程的眼里,它就暴露出来了。
忍不住插句嘴,这简直就是狭义相对论的世界观在计算机世界的翻版:一个参考系的观察者所看到两个类空间隔事件(可以是相继发生但因距离遥远而无因果关联)在另一个参考系中看来是颠倒的,但有因果关联的两事件在所有观察者看来时序都不会改变。好玩吧?

所以,表面的秩序井然背后有着巨大的复杂性怪兽,为了关牢它的笼子,我们需要约束我们的代码,用合适的机制来实现线程间通讯。要保证另一个线程看起来,两组存储器操作具有我们所希望的顺序,需要用 Memory Fence Function. 这里不再涉及,对更多细节感兴趣的同学,请参考CUDA C Programming Guide B.5. Memory Fence Functions等章节。
(未完待续)

CUDA, 软件抽象的幻影背后 之二

先更新到这儿,稍后再回来抛光查错。CUDA比较杂,我一写起来容易满嘴跑火车弄出错误,欢迎拍砖。

**********************************************************************

版权声明:原创作品,欢迎转载,但转载请以超链接形式注明文章来源(planckscale.info)、作者信息和本声明,否则将追究法律责任。

上一篇里说到,有两点对CUDA的计算能力影响甚大:数据并行,以及用多线程掩盖延迟。接下来我们要深入到其硬件实现,看一看这些机制是如何运作的。

通常人们经常说某GPU有几百甚至数千的CUDA核心,这很容易让人联想到多核CPU。不过事实上两种“核心”是不一样的概念,GPU的CUDA核心只相当于处理器中的执行单元,负责执行指令进行运算,并不包含控制单元。可以类比到CPU核心的是流多处理器(Streaming Multiprocessor,简写为SM. Kepler中叫做SMX,Maxwell中叫做SMM),通常一个GPU中有数个SM,而每个SM中包含几十或者上百个CUDA核心,以及数个warp scheduler(相当于控制单元)。如下图GM204中有16个SM,每个SM中有128个CUDA核心,4个warp scheduler。

GeForce_GTX_980_SM_Diagram-545x1024

图 1.  GM204的SM结构图

每个SM中有大量的寄存器资源,在GM204的例子中,有总共64k 32-bit寄存器,可以养活成千上万的线程。SM中另外一个重要资源是Shared Memory,没错,它正是软件抽象中Shared Memory的对应物。在GM204中,每个SM有96KB的Shared Memory.

到这里,SM在软件抽象里的对应也呼之欲出了,没错,正是Block。我们不妨先摆出这个对应:
Block <-> SM
Thread执行 <-> CUDA Cores
Thread数据 <-> Register/Local Memory

同一Grid下的不同Block会被分发到不同的SM上执行。SM上可能同时存在多个Block被执行,它们不一定来自同一个kernel函数。每个Thread中的局域变量被映射到SM的寄存器上,而Thread的执行则由CUDA核心来完成。

SM上可以同时存在多少个Block?这由硬件资源的消耗决定:每个SM会占用一定数量的寄存器和Shared Memory,因此SM上同时存活的Block数目不应当超过这些硬件资源的限制。由于SM上可以同时有来自不同kernel的Block存在,因此有时候即便SM上剩余资源不足以再容纳一个kernel A的Block,但却仍可能容纳下一个kernel B的Block.

接下来一个很重要的问题是Block如何被执行。我们可以看到,SM上的CUDA核心是有限的,它们代表了能够在物理上真正并行的线程数——软件抽象里,Block中所有的线程是并行执行的,这只是个逻辑上无懈可击的抽象,事实上我们不可能对一个任意大小的Block都给出一个同等大小的CUDA核心阵列,来真正并行的执行它们。
因而有了Warp这个概念:物理上,Block被划分成一块块分别映射到CUDA核心阵列上执行,每一块就叫做一个Warp.目前,CUDA中的Warp都是从threadIdx = 0开始,以threadIdx连续的32个线程为一组划分得到,即便最后剩下的线程不足32个,也将其作为一个Warp.CUDA kernel的配置中,我们经常把Block的size设置为32的整数倍,正是为了让它能够精确划分为整数个Warp(更深刻的原因和存储器访问性能有关,但这种情况下仍然和Warp的size脱不了干系)。
在GM204的SM结构图里我们可以看到,SM被划分成四个相同的块,每一块中有单独的Warp Scheduler,以及32个CUDA核心。Warp正是在这里被执行。
Warp的执行非常类似于SIMD. Warp中的活动线程由Warp Scheduler驱动,同步执行。我们可以看到,GM204中32个CUDA核心共享一个Warp Scheduler. 关于Warp执行中可能出现的复杂些的问题,留到下文另外说。

现在可以整理一下这个世界的图景了。SM上存活着几个Block,每个Block中的变量占据着自己的寄存器和Shared Memory,Block被划分成32个线程组成的Warp. 这样,大量的Warp生存在SM上,等待被调度到CUDA核心阵列去执行。

Warp Scheduler正如其名,是这个Warp世界里的调度者。当一个Warp执行中出现等待(存储器读写延迟等)后,Warp Scheduler就迅速切换到下一个可执行的Warp,对其发送指令直到这个Warp又一次出现等待,周而复始。这就是上一篇所说“用多线程掩盖延迟”在硬件图景下的模样。

CPU_GPU_COMPARE

图 2.  GPU用多个Warp掩盖延迟 / 与CPU计算模式的对比

本图引用自PPT “CUDA Overview” from Cliff Woolley, NVIDIA.

如图,GPU用多个Warp快速切换来掩盖延迟,而CPU用快速的寄存器来减小延迟。两者的重要区别是寄存器数目,CPU的寄存器快但少,因此Context Switch代价高;GPU寄存器多而慢,但寄存器数量保证了线程Context Switch非常快。

多少线程才能够掩盖掉常见的延迟呢?对于GPU,最常见的延迟大概要数寄存器写后读依赖,即一个局域变量被赋值后接着不久又被读取,这时候会产生大约24个时钟周期的延迟。为了掩盖掉这个延迟,我们需要至少24个Warp轮流执行,一个Warp遇到延迟后的空闲时间里执行其余23个Warp,从而保持硬件的忙碌。在Compute Capability 2.0,SM中有32个CUDA核心,平均每周期发射一条指令的情况下,我们需要24*32 = 768个线程来掩盖延迟。
保持硬件忙碌,用CUDA的术语来说,就是保持充分的Occupancy,这是CUDA程序优化的一个重要指标。

(未完待续)

***********************************************

一些后续补充。

网友邵:

SM上可以同时存在多少个Block,除了受到资源的限制之外,还受到设备上限的限制,每个SM有一个Device Limit,warps和blocks不能超过对应的上限。

CUDA, 软件抽象的幻影背后

版权声明:原创作品,欢迎转载,但转载请以超链接形式注明文章来源(planckscale.info)、作者信息和本声明,否则将追究法律责任。

今天最酷炫的事情应该就是来自老黄的这条消息:1TFLOPS,P < 15W, ARM Cortex A57 * 4 + ARM Cortex A53 * 4 +  Maxwell 256 CUDA Cores,  Tegra X1.

tegrax1
图1.  Tegra X1

本想挖掘一下写篇博,但目前报道满天飞没太大必要了。于是又想起了这个命途多舛的话题:CUDA. 关于CUDA我写了两次,第一次不满意未发,第二次成文后保存失败灰飞烟灭在热力学第二定律决定的命运里。今天借X1的东风,我们再来聊聊CUDA.

**********************************************************

CUDA是个以性能为第一目标的语言,这也决定了CUDA开发者所要面对的复杂性远远要多于CUDA语言所抽象出来的编程模型本身。这大概会是软件抽象所要面对的永恒话题,我们可以去抽象出一组逻辑上漂亮完备的功能基元,却不能保证从性能的观点看它们同样也是小开销的基本操作。具体在CUDA里,最典型的例子是内存<->显存数据交换,一个简单的拷贝操作在性能上却是让人难以接受的,这背后是PCIE总线;对性能影响稍小些的例子比如Global Memory的读写需要考虑对齐,这是由于硬件层面warp和cache机制的体现;再者如过度臃肿的kernel或block过大导致寄存器耗尽,局域变量被吐到Local Memory导致的性能损失。

所有这些,都要求我们透过CUDA简洁干净的编程模型,看到软件抽象的美丽幻影背后那个不同的世界,它存在于抽象之下我们不熟悉的另一个层次,却透过性能这一个几乎是唯一的方式来影响着我们的软件。这颇类似万有引力与我们世界的关系:引力是唯一能透入额外维度的基本相互作用,如果世界有我们所不知道的维度存在,如何才能感受到那个世界对我们的影响?答案就是用引力。看过《星际穿越》的同学们想必对此有些印象。

在深入GPU的硬件架构之前,我们不妨先探讨一下这个问题:为什么GPU具有这么高的计算能力?我们试着归纳两条最主要的原因。

目前典型的计算模式有两种,CPU式的高速低延迟串行计算,和GPU式的高延迟高吞吐大规模并行计算。CPU是人们熟知的,它具有高速的内部寄存器和Cache,现代CPU又加入了多级流水线,猜测、乱序执行,超线程等技术加速其指令吞吐能力,具有快速的响应能力,但是对于大量数据的处理却相对不够用。这一点3D游戏应用就是典型的例子,当然,这就是GPU崛起的契机。
GPU天生为数据的批量处理而生,它擅长的是在大量数据上同时做同样或几乎一致(这点很重要)的计算。为什么要求一样的计算?这一点可以从很多个角度来回答。
最重要的一个回答是,多个线程同步执行一致的运算,使得我们可以用单路指令流对多个执行单元进行控制,大幅度减少了控制器的个数和系统的复杂度(设想成千上万的线程各自做不同的事情,如果再有线程间通讯/同步,将会是怎样的梦魇)。
另一方面,现实世界中应用在大规模数据上的计算,通常都涵盖在这一计算模式之中,因而考虑更复杂的模式本质上是不必要的。比如计算大气的流动,每一点的风速仅仅取决于该点邻域上的密度和压强分布;再如计算图像的卷积,每一个输出像素都仅是对应源点邻域和一个卷积核的内积。从这些例子中我们可以看到,除了各个数据单元上进行的计算是一样的,计算中数据之间的相互影响也具有某种“局域性”,一个数据单元上的计算最多需要它某个邻域上的数据。这一点意味着线程之间是弱耦合的,邻近线程之间会有一些共享数据(或者是计算结果),远距离的线程间则独立无关。
这个性质反映在CUDA里,就是Block划分的两重天地:Block内部具有Shared Memory,线程间可以共享数据、通讯和同步,Block外部则完全独立,Block间没有通讯机制,相互执行顺序不影响计算结果。这一划分使得我们既可以利用线程间通讯做一些复杂的应用和算法加速,又可以在Block的粒度上自由调度计算任务,在不同计算能力的硬件平台上自适应的调整任务安排。
现在我们把注意力放在“几乎一致”这里。最简单的并行计算方案是多路数据上同时进行完全一致的计算,即SIMD(单指令流多数据流)。这种方案是非常受限的。事实上我们可以看出,“完全一致”是不必要的。只要这些计算在大多数时候完全一致,就可以对它们做SIMD加速,而在计算分叉,各个线程不一致的特殊情况下,只需要分支内并行,分支间串行执行即可,毕竟这些只是很少出现的情况。这样,把“完全一致”这个限制稍微放松,就可以得到更广阔的应用范围和不输于SIMD的计算性能,即SIMT(单指令流多线程)的一个重要环节,这是GPU强大处理能力的第一个原因。

一个或许让每个初学者都惊讶的事实是这样一组数据:Global Memory访存延迟可以达到数百个时钟周期,即便是最快的Shared Memory和寄存器在有写后读依赖时也需要数十个时钟周期。这似乎和CUDA强大的处理能力完全相悖——如果连寄存器都这么慢,怎么会有高性能呢?难道这不会成为最大的瓶颈吗?
答案恰恰就出乎意料:不,这不是瓶颈,这个高延迟的开销被掩盖了,掩盖在大量线程之下。更清楚的说,当一组线程(同步执行,类似于SIMD的一个线程组,在CUDA里叫做warp)因为访存或其他原因出现等待时,就将其挂起,转而执行另一组线程,GPU的硬件体系允许同时有大量线程存活于GPU的SM(流多处理器)之中,控制单元在多组线程之间快速切换,从而保证资源的最大利用率——控制单元始终有指令可以发放,执行单元始终有任务可以执行,仍然可以保持最高的指令吞吐,每个单元基本都能保持充分的忙碌。
这就是GPU硬件设计中非常有特色的基本思想:用多线程掩盖延迟。这一设计区别于CPU的特点是,大量高延迟寄存器取代了少量低延迟寄存器,寄存器的数量保证了可以有大量线程同时存活,且可以在各组线程间快速切换。尽管每个线程是慢的,但庞大的线程数成就了GPU的数据吞吐能力。此为高性能的第二个原因。

这文又要写成未完待续了。接下来的日子,不填完旧坑不再开新话题。

算法描述与性能优化的解耦——Halide语言 (1)

版权声明:原创作品,欢迎转载,但转载请以超链接形式注明文章来源(planckscale.info)、作者信息和本声明,否则将追究法律责任。

程序的结构和运行效率常常被人们看作是难以调和的。这个事实源于我们把一个数学上结构清晰良好的(比如,用递归形式刻画的)算法映射到一个现实的不完美的计算模型上,这个模型计算是有代价的,要极小化这个代价就需要尽可能的重用中间计算结果,减少依赖增加指令级并行,充分利用空间和时间的Locality和存储体系的分级结构,等等…

于是长久以来,程序中描述算法要”做什么“的逻辑,掩盖在了用来优化性能,描述”怎么做“的芜杂逻辑之下。当然,我们也有能够只通过清晰刻画”做什么“来编程的语言(如函数式语言家族,尤其是以Haskell为代表的纯函数式语言),这些语言收起了让用户自己规划计算过程的权限,把计算过程的优化交给编译器自动完成。这导致了人们对其效率的不信任,事实上这个问题也确实普遍存在,聪明的编译器只是少数,而且他们也未必能保证给出一个高度优化的程序。所以至今,对高性能有要求的各种库仍然采用C等提供了底层操作的语言实现。

性能优化会破坏软件结构这一点带来了无尽的苦恼,和风险。这使得我们每做一次优化都需要格外慎重。性能优化后的代码基本被凝固,不再具有良好的可维护、可修改性质,这样我们跨平台移植或者改变优化方案时就会面临代码需要大片重写的风险。所以人们确立了”不成熟的优化是万恶之源“这样的信条,强迫自己把性能优化留到最后一步,万不得已时采用最成熟最保守的思路去优化。

真的只能接受这个现状吗?我们知道抽象和解耦是软件的灵魂,当两件不同目的的事情纠缠在同一段代码中时,意味着需要把两件事情各自抽象出来,解耦成简单独立的逻辑。这里我们正面对一个典型的解耦问题:算法描述(做什么)和性能优化(怎么做)需要被解耦。函数式语言虽然能做到这个解耦,但是它把优化工作交给不那么靠谱的编译器了,那能否把这个优化工作交还给设计者自己呢?设想我们如果能够独立构造两个逻辑,就可以利用如纯函数式语言这样具有强大描述力的工具,寥寥数笔刻画出算法;然后再针对具体的硬件平台,实现相应的优化计算方案。不同的平台间的移植只需要替换这个优化方案部分。更好的一点是,我们可以尝试更激进的优化方案,测试各种各样的方案,这不过是个替换而已,而且对算法的功能本身没有影响。

解耦工作的难度一定程度上取决于要解耦的两个概念是否能够清晰的区分开来。算法描述和性能优化的解耦是不容易的,因为一般说来这两个概念不易区分。但在图像处理这样的领域里,计算具有典型的模式(数据在pipeline上流动,被各个节点依次处理),我们仍然可以把二者很好地解耦。

Halide就是这样一门语言。

Halide是由MIT、Adobe和Stanford等机构合作实现的图像处理语言,它的核心思想即解耦算法和优化,事实也证明它是成功的,在各种实例中它均以几分之一的代码量实现出同等或者数倍于手工C++代码的效能,更不用提代码的可维护性和开发效率。

先上例子(来自Halide的文献”Decoupling Algorithms from Schedules for Easy Optimization of Image Processing Pipelines”),三个版本的图像模糊算法,以及他们各自的性能。




void box_filter_3x3(const Image & in, Image & blury) {
	Image blurx(in.width(), in.height()); // allocate blurx array
	for (int y = 0; y &lt in.height(); y++)
	for (int x = 0; x &lt in.width(); x++)
		blurx(x, y) = (in(x - 1, y) + in(x, y) + in(x + 1, y)) / 3;
	for (int y = 0; y &lt in.height(); y++)
	for (int x = 0; x &lt in.width(); x++)
		blury(x, y) = (blurx(x, y - 1) + blurx(x, y) + blurx(x, y + 1)) / 3;
}



9.96 ms/megapixel
(quad core x86)

代码1.  C++实现图像模糊,结构良好但效率差。




void box_filter_3x3(const Image & in, Image & blury) {
	__m128ione_third = _mm_set1_epi16(21846);
#pragmaomp parallel for
	for (int yTile = 0; yTile &lt in.height(); yTile += 32) {
		__m128ia, b, c, sum, avg;
		__m128i blurx[(256 / 8)*(32 + 2)]; // allocate tile blurx array
		for (int xTile = 0; xTile &lt in.width(); xTile += 256) {
			__m128i*blurxPtr = blurx;
			for (int y = -1; y &lt 32 + 1; y++) {
				const uint16_t *inPtr = & (in[yTile + y][xTile]);
				for (int x = 0; x &lt 256; x += 8) {
					a = _mm_loadu_si128((__m128i*)(inPtr - 1));
					b = _mm_loadu_si128((__m128i*)(inPtr + 1));
					c = _mm_load_si128((__m128i*)(inPtr));
					sum = _mm_add_epi16(_mm_add_epi16(a, b), c);
					avg = _mm_mulhi_epi16(sum, one_third);
					_mm_store_si128(blurxPtr++, avg);
					inPtr += 8;
				}
			}
			blurxPtr = blurx;
			for (int y = 0; y &lt 32; y++) {
				__m128i*outPtr = (__m128i*)(& (blury[yTile + y][xTile]));
				for (int x = 0; x &lt 256; x += 8) {
					a = _mm_load_si128(blurxPtr + (2 * 256) / 8);
					b = _mm_load_si128(blurxPtr + 256 / 8);
					c = _mm_load_si128(blurxPtr++);
					sum = _mm_add_epi16(_mm_add_epi16(a, b), c);
					avg = _mm_mulhi_epi16(sum, one_third);
					_mm_store_si128(outPtr++, avg);
				}
			}
		}
	}
}



11x fasterthan a
naïve implementation
0.9 ms/megapixel
(quad core x86)

代码2.  上一段代码的优化版本,效率好但结构性破坏。




Func halide_blur(Func in) {
	Func tmp, blurred;
	Var x, y, xi, yi;
	// The algorithm
	tmp(x, y) = (in(x - 1, y) + in(x, y) + in(x + 1, y)) / 3;
	blurred(x, y) = (tmp(x, y - 1) + tmp(x, y) + tmp(x, y + 1)) / 3;
	// The schedule
	blurred.tile(x, y, xi, yi, 256, 32)
		.vectorize(xi, 8).parallel(y);
	tmp.chunk(x).vectorize(x, 8);
	return blurred;
}



0.9 ms/megapixel

代码3.  Halide代码,清晰简短,且一样高效。

我们可以看到,在Halide所实现的版本中,代码分成两部分,一部分是描述算法的algorithm部分,采用典型的函数式风格定义出要计算什么;另一部分则是指定”如何计算“的schedule部分。Halide目前没有自己的语法解析器,它的前端直接嵌入在C++里,作为一个库来使用。我们构造出一个图像处理算法后,可以把它编译到诸如x86/SSE, ARM v7/NEON, CUDA, Native Client, OpenCL各种平台上。

用schedule这样一个概念来抽象各种各样的底层优化技巧是整个方案里最关键的一环。不同平台有不同的优化技巧,怎样才能用一个统一的观点去处理它们,使得我们能在一个足够简单、与底层细节无关的世界观里处理优化问题?Halide用这样的观点来统一各种性能优化方法:它们都是控制存储或计算顺序的手段。这样,Halide通过提供一系列控制计算过程中存储和计算顺序的工具而帮助我们描绘性能优化方案。

Halide目前并没有太多的考虑编译器自动优化的问题,但这是一个漂亮的开端。如果将来在手动优化的同时仍有强大的编译器优化做后盾,将会是一番什么景象?

本站将持续跟踪这方面的进展。关于Halide语言进一步的剖析,请等待下篇。

(未完待续)