Hough transform, оптимизация

Подозреваю что много времени уходит на sin + cos. Попробуй заменить на один вызов sincos (FPU может вычислить синус и косинус одновременно по инструкции fsincos). Или вообще на лукап по таблице.

Приведение типов из целых в вещественные внутри цикла тоже шустрости не способствует.

Спасибо, что-то о таблице я сразу и не подумал. Попробую затабулировать. Для шага в 1° надо будет хранить всего лишь 180 значений (если вычислять cos(x)=sin(90°-x)).

В конечной функции планируется вычислять углы с шагом хотя бы 0.25°. Думаю, здесь табулирование будет намного быстрей постоянного вычисления тригонометрических функций.

А для вычисления

N = (double)j*cos(theta) + (double)i*sin(theta);

без преобразований типов не обойтись.

Можно попробовать в циклах сразу оперировать числами с плавающей точкой.

Как вариант - да.

Мне так кажется, что здесь самое дорогое - это вычисление косинуса и синуса. Можно не вычислять на каждом шаге синус и косинус, а вычислять через синус и косинус на предыдущем шаге (через формулу косинуса и синуса суммы аргументов):

cos(x+dx) = cos(x)cos(dx) - sin(x)sin(dx)

sin(x+dx) = sin(x)cos(dx) + cos(x)sin(dx)

Тригонометрию вычислять по таблице. Первый цикл на SSE сделать. Давай компилирующийся фрагмент, пооптимизируем :)

А ещё, если есть возможность, то лучше преобразование делать в YUV на слое яркости.

Буду делать табулирование синуса, а вот работающий код. Изображение берется с веб-камеры («микроскоп»).

А оно у тебя не с JPEG с камеры умеет работать?

Разобрался, умеет.

Чуть кода подправить надо было, чтобы без jpeg работало :)

В общем, заменил я вычисление тригонометрических функций на таблицы. Получил значительный прирост производительности! Теперь вычисление происходит за 10-20мс!

Правда, если опорных точек больше 20тыс., все равно ~70мс вычисляет...

Если не будешь менять формат на планарный, то хотя бы выравнивай пикселы на границу 4-х байт. У тебя сейчас 3, такое переварит только Nehalem (Core i7). На обычной Корке доступ к невыровненным данным ажно в 4 раза медленней.

А что именно у тебя на мультиэкране в 1-й, 3-й и 4-й четвертях выводится?

Попробуй каким-нибудь профайлером посмотреть, на чем оно тормозит.
Но, как уже сказали, sin и cos, особенно при высокой точности, тормозят.

В смысле, выравнивать двумерный массив? А разве это не забота компилятора?

А что именно у тебя на мультиэкране в 1-й, 3-й и 4-й четвертях выводится?

Вейвлет-преобразование первого уровня. С 1-й по 4-ю четверти выводятся, соответственно, HL, LL, LH и HH компоненты. Суммируя HL и LH и проводя обратное вейвлет-преобразование, я получаю картину краев на изображении (это быстрее, чем при помощи дифференциальных фильтров вроде лапласиана). По ширине гистограммы HL/LH изображений получаю информацию об общей фокусировке (по полю). Если строить гистограммы выделенных областей HL/LH, получим оценку фокусировки в этих областях.

sin/cos уже табулирую:

	for(i = 0, ii = 0.; i < GRAB_HEIGHT; i++, ii += 1.){
		ptr[0] = imageR[i]; ptr[1] = imageG[i]; ptr[2] = imageB[i];
		for(j = 0, jj = 0.; j < GRAB_WIDTH; j++, jj += 1.){
			N = 0;
			for(k = 0; k < 3; k++){
				if((wdth = *ptr[k]) > Max){
					*GLptr++ = (GLubyte) (*ptr[k] * scale );
					N++;
				}
				else *GLptr++ = 0;
				ptr[k]++;
				}
			if(N == 3){
				cntr++;
				cosptr = CosBuf; sinptr = SinBuf;
				for(k=0; k<SinBufSize; k++){
					N = (int)(jj * (*cosptr++) + ii * (*sinptr++));
					if(N > 0 && N < Rmax)
						Hough[k][N]++;					
				}
			}
		}
	}

В смысле, выравнивать двумерный массив? А разве это не забота компилятора?

Пиксели. Если RGB24 (настоящий, 3 байта на пиксель), то добивать пустым байтом.

Вейвлет-преобразование первого уровня. С 1-й по 4-ю четверти выводятся, соответственно, HL, LL, LH и HH компоненты. Суммируя HL и LH и проводя обратное вейвлет-преобразование, я получаю картину краев на изображении (это быстрее, чем при помощи дифференциальных фильтров вроде лапласиана).

Что-то типа Robert's cross не дешевле будет для выявления граней?

Сама картинка у меня в линейном массиве. В двумерных - только RGB компоненты (по-отдельности).

Что-то типа Robert's cross не дешевле будет для выявления граней?

Мне вейвлеты все равно нужны (для фокусировки). А вычисление вейвлет-преобразования дешевле, чем двух БПФ (для свертки). А делать свертку напрямую - вообще долгая песня.

Сама картинка у меня в линейном массиве. В двумерных - только RGB компоненты (по-отдельности).

Как хранятся пиксели в этом линейном массиве?

В дальнейшем для ускорения кода думаю перейти к черно-белому изображению (чтобы не считать по шесть вейвлет-преобразований на каждую картинку). Но все равно, хочется хотя бы немного ускорить код (хотя, здесь, конечно, хватит и 1 кадра в секунду, но есть подобные задачи, где производительность важна).

Возможно ли как-то еще оптимизировать? А ведь хочется еще и рассчитывать положение линий более точно - скажем, с шагом по 0.25° по тета и 0.1 по R... Уменьшение шага приведет к тому, что придется искать не максимум в точке на диаграмме тета-R, а пятно с максимальной интенсивностью (т.е. нужно будет вовлекать еще и БПФ для фильтра лапласиана гауссианы)...

По порядку, конечно (RGB). Какой смысл их выравнивать, если это - текстура OpenGL?

По порядку, конечно (RGB). Какой смысл их выравнивать, если это - текстура OpenGL?

Выровняй, посмотри ;) GL поддерживает 32-битный RGBA с незначащим альфа-каналом.

Время, за которое выполняется отрисовка текстуры, незначительно, по сравнению с остальным. Так что вообще никакого смысла делать RGBA не вижу (наоборот, процедура заполнения массива изображения станет выполняться дольше, т.к. этот массив станет в 4/3 раза больше, а данные из исходных массивов R, G и B в него все равно надо будет переписывать).

RedrawWindow можно разогнать немного. Там, фактически, нужно сделать bind новой текстуре и glTexSubImage2D. Всё остальное можно сделать один раз при инициализации.

Время, за которое выполняется отрисовка текстуры, незначительно, по сравнению с остальным.

Ты эту текстуру в готовом виде с потолка берёшь или ручками заполняешь? Хотя, что это я... У тебя побайтно же всё обрабатывается, до выравнивания ещё есть, что оптимизировать :)

Там, фактически, нужно сделать bind новой текстуре и glTexSubImage2D.

Какие-то у меня были с этим проблемы. Да я еще и никак не могу нормальной книжки найти по openGL (есть в электронном виде, но это же надо печатать...), а гугление привело к тому, что есть.

Вот я пока и не обращаю внимания на функции openGL: начинаю с самых узких мест. Оптимизировал, как мог, вейвлет-преобразования (хотя, там, возможно, тоже что-то можно выгадать), сейчас хочу преобразования Хоу/Хафа оптимизировать... А результирующие текстуры - это либо изображение краев объекта с нарисованной поверх линией, проходящей через наибольшее кол-во пикселей, либо изображение вейвлет-образа, либо изображение гистограммы выбранного компонента вейвлет-образа.

В память пиксел лучше (быстрее) писать за раз. Атомарно 3 байта записать нельзя, можно 4. На самом деле, не стоит пренебрегать обращениями к памяти, когда работаешь с графикой. Графика жирная, копейки, потерянные на каждом пикселе, в итоге даже не в рубли, в сотни вырастают :) По этой же причине в фильтрах нужно минимизировать обращения к памяти, где только можно.

> библиотеки, включающие в себя уйму дополнительных ненужных мне функций.

одна из них Intel IPP. Вообще, быстрее чем интел довольно трудно сделать. Если нет проблем с лицензией я бы взял IPP.

Если нет проблем с лицензией я бы взял IPP.

Не хочется мне с интеловскими проприетарными библиотеками связываться. Одно дело - украсть функцию из numerical recipes, другое - целую библиотеку у интела :)