PyOpenCLでアプリケーション(PIL & PyOpenGL)

つい最近も書いた気がするんですがGPGPU Advent Calendarの15日目です。もうゴールしてもいいよね？

PyOpenCLを使えばGPGPUのプログラムを作るのがちょっと簡単になるということで、Pythonの他のライブラリも使ってGPGPUを使ったアプリケーションを簡単に作る方法を紹介します。OpenCLのImageオブジェクトとPythonの画像処理ライブラリPILを使った例と、OpenCLのOpenGL InteroperabilityとPyOpenGLを使った例の2つを取り上げます。

PyOpenCL + PIL

PyOpenCLでImageオブジェクトを扱うにはpyopencl.Imageを使います。そのまんまですね。ImageFormatオブジェクトで画像の表現形式を与えます。今回は1ピクセルあたりRGBA各8ビットの32ビットにするので cl.ImageFormat(cl.channel_order.RGBA, cl.channel_type.UNSIGNED_INT8) という感じです。ホスト側メモリを受け持つnumpyのndarrayをnumpy.empty((height, width, 4), dtype=numpy.uint8)という風に作るとImageオブジェクトから直接ピクセルの情報を受け取れます。PILを使えばndarrayで持っている画像情報を簡単にファイルに保存したりができます。

ソースコード全体は下の通りです。ちなみに、PyOpenCLでPythonのdocstringの中にカーネル関数を書くときにはVim用のシンタックスハイライトを使っておくと、OpenCL C言語の部分も色づけしてくれるので便利ですよ。

polarchecker.py

# vim: filetype=pyopencl.python
import numpy as np
from PIL import Image
import pyopencl as cl
 
SRC = '''//CL//
__constant const int CHECKER_SIZE = 10;
 
float2 polar(const int2);
int disc(const int2);
 
 
float2 polar(const int2 pos)
{
    const float x = pos.x;
    const float y = pos.y;
    return (float2)(sqrt(hypot(x, y)), atan2(x, y));
}
 
 
int disc(const int2 pos)
{
    const float2 ppos = polar(pos);
    return ((int)(ppos.x) + (int)(ppos.y * CHECKER_SIZE / M_PI_F)) % 2;
}
 
 
__kernel void gen(write_only image2d_t dest)
{
    const int2 pos = (int2)(get_global_id(0), get_global_id(1));
    uint4 pix;
    if (disc(pos)) {
        pix = (uint4)(0, 0, 255, 255);
    } else {
        pix = (uint4)(255, 255, 255, 255);
    }
    write_imageui(dest, pos, pix);
}
'''
 
ctx = cl.create_some_context()
queue = cl.CommandQueue(ctx)
prg = cl.Program(ctx, SRC).build()
 
width = 640
height = 480
fmt = cl.ImageFormat(cl.channel_order.RGBA, cl.channel_type.UNSIGNED_INT8)
buf = cl.Image(ctx, cl.mem_flags.WRITE_ONLY, fmt, shape=(width, height))
 
prg.gen(queue, (width, height), None, buf)
 
result = np.empty((height, width, 4), dtype=np.uint8)
cl.enqueue_copy(queue, result, buf, origin=(0, 0), region=(width, height))
 
image = Image.fromarray(result)
image.save('checker.png')

カーネル関数部分でOpenCL C言語独特な書き方をしてると思います。この辺はOpenCL API 1.2 Reference Cardとか眺めてみるといろいろ新たな発見があるかもしれませんよ。

実行するとこんな感じの画像が生成されます。

カーネル関数がピクセル単位の処理になってプログラムの見通しもよくなるので良いですね。

PyOpenCL + PyOpenGL

PyOpenCLとPyOpenGLの連携を真面目にやると手順が多くてちょっと説明が大変なので、PyOpenCLに同梱のサンプルの解説をすることにします。このプログラムはOpenGLでsin波のグラフを出力するだけのプログラムで、初期化時にOpenCLを使って頂点位置の計算を行っています。

プログラムのメイン部分が以下のようになっています。

if __name__ == '__main__':
    import sys
    glutInit(sys.argv)
    if len(sys.argv) > 1:
        n_vertices = int(sys.argv[1])
    glutInitWindowSize(800, 160)
    glutInitWindowPosition(0, 0)
    glutCreateWindow('OpenCL/OpenGL Interop Tutorial: Sin Generator')
    glutDisplayFunc(display)
    glutReshapeFunc(reshape)
    initialize()
    glutMainLoop()

OpenCL/OpenGL連携で重要なのはinitializeの中身だけですので順番に見ていきましょう。

まず始めにOpenGL連携を有効にしたOpenCL Contextの作成です。

    platform = cl.get_platforms()[0]

    from pyopencl.tools import get_gl_sharing_context_properties
    import sys
    if sys.platform == "darwin":
        ctx = cl.Context(properties=get_gl_sharing_context_properties(),
                devices=[])
    else:
        # Some OSs prefer clCreateContextFromType, some prefer
        # clCreateContext. Try both.
        try:
            ctx = cl.Context(properties=[
                (cl.context_properties.PLATFORM, platform)]
                + get_gl_sharing_context_properties())
        except:
            ctx = cl.Context(properties=[
                (cl.context_properties.PLATFORM, platform)]
                + get_gl_sharing_context_properties(),
                devices = [platform.get_devices()[0]])

ややこしそうなことを書いてますが、多くの場合はtry節の中の部分だけで事足りると思います。私の環境はMacで最初のifの条件に入るのですが、逆にエラーで落ちるのでコメントアウトしましたが普通に動いてそうな気がします。

次にOpenGLのBufferを初期化します。注目するところはglGenBuffersでvboを取得するところです。ここで取得した値はPyOpenCLのGLBufferオブジェクトを生成するのに必要です。

    glClearColor(1, 1, 1, 1)
    glColor(0, 0, 1)
    vbo = glGenBuffers(1)
    glBindBuffer(GL_ARRAY_BUFFER, vbo)
    rawGlBufferData(GL_ARRAY_BUFFER, n_vertices * 2 * 4, None, GL_STATIC_DRAW)
    glEnableClientState(GL_VERTEX_ARRAY)
    glVertexPointer(2, GL_FLOAT, 0, None)

続くGLBufferオブジェクトの生成が以下の一行です。

    coords_dev = cl.GLBuffer(ctx, cl.mem_flags.READ_WRITE, int(vbo))

最後にOpenCLからOpenGLのBufferを操作します。OpenCL側でOpenGLのBufferを操作するときはenqueue_acquire_gl_objectsとenqueue_release_gl_objectsではさみます。

    cl.enqueue_acquire_gl_objects(queue, [coords_dev])
    prog.generate_sin(queue, (n_vertices,), None, coords_dev)
    cl.enqueue_release_gl_objects(queue, [coords_dev])

基本的なところはこれだけです。プログラム全体を実行すると以下のような画面が表示されます。

C言語などで書くよりはだいぶ楽になります。ただしOpenGLの知識はある程度必要かもしれませんね。

もうちょっと動いたりするデモはこの辺とかこの辺にあります。