Реализация операции scatter на GPU на примере построения гистограмм.

Ряд операций по обработке данных можно разделить как операции типа gather и операции типа scatter.

Типичным примером операции типа gather является применение матричного фильтра к изображению - каждый результирующий пиксел (для сглаживания с ядром k_ij) будет задаваться следующей формулой:

Из приведенной формулы легко видно, что для получение одного тексела результирующего изображения необходимо прочесть несколько текселов исходного и объединить их при помощи некоторой операции (в данном случае свертки).

Типичным примером операции scatter является построение гистограммы заданного изображения - весь диапазон возможных значений некоторого параметра (например, яркости) разбивается на ряд отрезков ("корзин", bins) и для каждого такого отрезка считается сколько исходных текселов попали в него.

Классическая реализация построения гистограмм на CPU инициализирует массив "корзин" и для каждого пиксела определяет номер соответствующей ему корзины и увеличивает счетчик по этому номеру (см. псевдокод ниже).

Однако если операции первого типа (gather) очень легко и эффективно реализуются на GPU, то операции второго типа (scatter) обычно крайне плохо ложатся на архитектуру современных GPU, так как требуют доступа к памяти типа random write.

Однако сравнительно недавно было замечено, что ряд операций данного типа (включая и построение гистограмм) можно довольно легко и эффективно реализовать на GPU.

Для этого сопоставим каждому текселу исходного изображения одну точку. Массиву "корзин" (bins) сопоставим текстуру в которую мы будем осуществлять рендеринг.

После этого выведем каждую точку исходного изображения именно как точку. Воспользовавшись операций чтения из текстуры в вершинном шейдере (поддерживаемой начиная с GPU GeForce 6xxx) можно легко прочесть цвет соответствующей точки, определить индекс корзины и изменить координаты данной точки таким образом, чтобы при рендеринге она попадала точно в нужную корзину.

Теперь если задать в качестве режима вывода alpha-blinding с параметрами (GL_ONE, GL_ONE), то вывод каждой точки будет увеличивать значение в соответствующей корзине и в результате вывода всех точек мы получим в массиве корзин интересующие нас количества точек (может быть только отмасштабированные некоторым множителем).

Ниже приводится вершинный шейдер, как раз осуществляющий перенаправление точек по корзинам в зависимости от яркости соответствующего тексела текстуры.

// // Fast GPU scatter via vertex texture fetch - vertex shader // uniform sampler2D srcMap; // map to compute histogram uniform float numBins; // number of bins (in [0,1] segment] void main(void) { vec4 color = texture2D ( srcMap, gl_MultiTexCoord0.xy ); // get color float lum = dot ( color.rgb, vec3 ( 0.3, 0.59, 0.11 ) ); float bin = floor ( lum * numBins ); gl_Position = gl_ModelViewProjectionMatrix * vec4 ( bin, 0.0, 0.0, 1.0 ); }

Соответствующий код на С++ приводится на листинге ниже.

Существует и еще один способ, позволяющий эффективно реализовать построение гистограмм на GPU, к тому же он часто оказывается быстрее. Этот способ заключается в использовании рендеринга в вершинный буфер для создания массива вершин, при этом каждой создаваемой вершине сразу же назначается правильные координаты.

После этого остается только вывести данный вершинный массив как набор точек, используя alpha-blinding.

Ниже приводятся соответствующие вершинный и фрагментный шейдеры, осуществляющие построение такого массива вершин.

// // Fragment shader GPU scatter via vertex texture fetch - fragment shader // uniform sampler2D srcMap; // map to compute histogram uniform float numBins; // number of bins (in [0,1] segment] void main () { vec3 clr = texture2D ( srcMap, gl_TexCoord [0].xy ).xyz; float lum = dot ( clr, vec3 ( 0.3, 0.59, 0.11 ) ); float bin = floor ( lum * numBins ); gl_FragColor = vec4 ( bin, 0.0, 0.0, 1.0 ); }

Ниже приводится фрагмент кода на С++, реализующий данный подход.

#include "libExt.h" #ifdef MACOSX #include <GLUT/glut.h> #else #include <glut.h> #endif #include <stdio.h> #include <stdlib.h> #include "libTexture.h" #include "Vector3D.h" #include "Vector2D.h" #include "Vector4D.h" #include "GlslProgram.h" #include "FrameBuffer.h" #include "utils.h" #include "fpTexture.h" #include "VertexBuffer.h" #define TEX_SIZE 256 // width & height of texture #define NUM_BINS 128 GLenum mapFormat; // texture format used for vertex map unsigned srcMap; // map with source image FrameBuffer buffer ( NUM_BINS, 1 ); FrameBuffer buffer2 ( TEX_SIZE, TEX_SIZE ); GlslProgram program; GlslProgram program2; VertexBuffer * vertexBuffer = NULL; // vertex coordinates float res [NUM_BINS*4]; float res2[TEX_SIZE*TEX_SIZE*4]; void reshape ( int w, int h ); void createVertexBuffer () { if ( vertexBuffer != NULL ) return; vertexBuffer = new VertexBuffer (); vertexBuffer -> bind ( GL_PIXEL_PACK_BUFFER_ARB ); vertexBuffer -> setData ( buffer2.getWidth () * buffer2.getHeight () * 4 * sizeof ( float ), NULL, GL_DYNAMIC_DRAW ); vertexBuffer -> unbind (); } // // create vertex coordinates using r2vb // void buildVertices () { createVertexBuffer (); // load source image map glActiveTextureARB ( GL_TEXTURE0_ARB ); glBindTexture ( GL_TEXTURE_2D, srcMap ); program.bind (); // perform rendering to FBO texture buffer2.bind (); // render to vertex buffer startOrtho ( TEX_SIZE, TEX_SIZE ); drawQuad ( TEX_SIZE, TEX_SIZE ); endOrtho (); program.unbind (); // copy vertex data from texture to vertex buffer vertexBuffer -> bind ( GL_PIXEL_PACK_BUFFER_ARB ); glReadPixels ( 0, 0, TEX_SIZE, TEX_SIZE, GL_RGBA, GL_FLOAT, NULL ); buffer2.unbind (); vertexBuffer -> getSubData ( 0, TEX_SIZE*TEX_SIZE*4*sizeof(float), res2 ); vertexBuffer -> unbind (); } void display () { float h = 1.0 / (float) TEX_SIZE; float h2 = h * 0.5; buildVertices (); startOrtho ( buffer.getWidth (), buffer.getHeight () ); glDisable ( GL_DEPTH_TEST ); glEnable ( GL_BLEND ); glBlendFunc ( GL_ONE, GL_ONE ); glPointSize ( 1.0 ); glDisable ( GL_POINT_SMOOTH ); glEnableClientState ( GL_VERTEX_ARRAY ); vertexBuffer -> bind ( GL_ARRAY_BUFFER ); glVertexPointer ( 4, GL_FLOAT, 0, NULL ); program2.bind (); buffer.bind (); glClear ( GL_COLOR_BUFFER_BIT ); glDrawArrays ( GL_POINTS, 0, TEX_SIZE*TEX_SIZE ); glBindTexture ( GL_TEXTURE_2D, 0 ); glBindTexture ( GL_TEXTURE_RECTANGLE_ARB, buffer.getColorBuffer () ); glReadPixels ( 0, 0, NUM_BINS, 1, GL_RGBA, GL_FLOAT, res ); glBindTexture ( GL_TEXTURE_RECTANGLE_ARB, 0 ); buffer.unbind (); program2.unbind (); endOrtho (); glutSwapBuffers (); float sum = 0; for ( int i = 0; i < NUM_BINS; i++ ) { sum += res [4*i]; printf ( "%4d ", (int)res [4*i] ); } printf ( "\nTotal: %d pixels.\n", (int) sum ); exit ( 0 ); } void reshape ( int w, int h ) { glViewport ( 0, 0, (GLsizei)w, (GLsizei)h ); glMatrixMode ( GL_PROJECTION ); glLoadIdentity (); glMatrixMode ( GL_MODELVIEW ); glLoadIdentity (); } void key ( unsigned char key, int x, int y ) { if ( key == 27 || key == 'q' || key == 'Q' ) // quit requested exit ( 0 ); } int main ( int argc, char * argv [] ) { // initialize glut glutInit ( &argc, argv ); glutInitDisplayMode ( GLUT_DOUBLE | GLUT_RGB | GLUT_DEPTH ); glutInitWindowSize ( 512, 512 ); // create window glutCreateWindow ( "Fast GPU scatter via r2vb" ); // register handlers glutDisplayFunc ( display ); glutReshapeFunc ( reshape ); glutKeyboardFunc ( key ); init (); initExtensions (); assertExtensionsSupported ( "ARB_pixel_buffer_object" ); if ( (mapFormat = fpFormatWithPrecision ( 32 ) ) == 0 ) { printf ( "Floating-point textures not supported !\n" ); return 1; } if ( !FrameBuffer :: isSupported () ) { printf ( "Fbo not supported\n" ); return 1; } if ( !GlslProgram :: isSupported () ) { printf ( "GLSL not supported\n" ); return 1; } srcMap = createTexture2D ( false, "../../Textures/floor.bmp" ); //256*256 pixels if ( srcMap == 0 ) { printf ( "Texture floor.bmp not found.\n" ); return 1; } unsigned tempMap = buffer.createColorRectTexture ( GL_RGBA, mapFormat ); buffer.create (); buffer.bind (); if ( !buffer.attachColorTexture ( GL_TEXTURE_RECTANGLE_ARB, tempMap ) ) printf ( "buffer error with color attachment\n"); if ( !buffer.isOk () ) printf ( "Error with framebuffer\n" ); buffer.unbind (); tempMap = buffer2.createColorRectTexture ( GL_RGBA, mapFormat ); buffer2.create (); buffer2.bind (); if ( !buffer2.attachColorTexture ( GL_TEXTURE_RECTANGLE_ARB, tempMap ) ) printf ( "buffer error with color attachment\n"); if ( !buffer2.isOk () ) printf ( "Error with framebuffer\n" ); buffer2.unbind (); if ( !program.loadShaders ( "r2vb.vsh", "r2vb.fsh" ) ) { printf ( "Error loading shaders:\n%s\n", program.getLog ().c_str () ); return 3; } if ( !program2.loadShaders ( "scatter2.vsh", "scatter-2.fsh" ) ) { printf ( "Error loading shaders:\n%s\n", program.getLog ().c_str () ); return 3; } program.bind (); program.setTexture ( "srcMap", 0 ); program.setUniformFloat ( "numBins", NUM_BINS ); program.unbind (); glutMainLoop (); return 0; }

При использовании данных подходов следует обратить внимание на точность представление числе в текстуре, в которой происходит накопление результатов - если использовать обычную текстуру типа GL_RBA8, то 8 бит дадут нам только 256 различных значений счетчика, после чего происходит переполнение.

Существует два варианта борьбы с этим два подхода. Первый вариант заключается в разбиении изображения на ряд блоков и одновременно строим гистограммы для каждого такого блока и в конце объединяем счетчики по блокам.

Второй подход заключается в использовании форматов текстур с большей разрешающей способностью - так использование 16-битовых floating point чисел позволяет получать до 2048 различных значения счетчика, а если использовать GPU серии GeForce 8xxx, поддерживающие alpha blending в 32-битовые целочисленные и floating point-форматы, то проблема фактически полностью решается (правда ценой более высоких требований к используемому GPU).

Данная статья написана по материалам статьи Efficient Histogram Generation Using Scattering on GPUs.

По этой ссылке можно скачать весь исходный код к этой статье. Также доступны для скачивания откомпилированные версии для M$ Windows, Linux и Mac OS X.