Оптический поток с CUDA

1. Оптический поток с CUDA Михаил Смирнов (ВМК МГУ)

2. План • Что такое оптический поток • Классический метод Horn-Schunck • Современный подход к реализации модели Horn-Schunck

7. Оптический поток • Набор векторов смещения • Плотный (для каждого пикселя) • Разреженный (для некоторого подмножества пикселей) • Описывает движение в кадре

11. Как найти оптический поток? • Яркость объекта в кадре не меняется • Если никто не меняет освещения

12. Как найти оптический поток? • Яркость объекта в кадре не меняется • Если никто не меняет освещения • Будем сопоставлять пиксели одинаковой яркости

14. Не самый удачный вариант • Нужна дополнительная информация

15. Aperture Problem

16. Дополнительная информация • Посмотрим на движение в целом

17. Дополнительная информация • Посмотрим на движение в целом • Соседние точки скорее всего двигаются одинаково

18. Гладкость поля смещений

19. Обозначения • Смещение по вертикали • Смещение по горизонтали

20. Гладкость поля смещений

21. Гладкость поля смещений

22. Гладкость поля смещений • «Хорошее» поле • Небольшие абсолютные значения производных

23. «Хорошее» поле смещений • Сохраняет яркость • Гладкое

24. Обозначения • Прямоугольная область изображения • Яркость изображения с номером в точке

25. Классическая модель Horn-Schunck • Сохранениеяркости • data term • Гладкость поля • smoothness term

26. Точка минимума • Линеаризуем подинтегральное выражение в • Уравнения Эйлера-Лагранжа для линеаризованной модели Граничные условия: отражение

27. Метод Якоби для дискретных уравнений Эйлера-Лагранжа где

28. РеализацияОбщая схема • Копирование изображений в память GPU • Вычисление производных изображения • Решение СЛАУ • Копирование результат в основную память

29. Производные изображения • Пятиточечный шаблон: • Усреднение по времени

30. Производные изображенияРеализация • Чтение через текстуры • Нормализованные координаты • cudaAddressModeMirror

31. Производные изображения float t0, t1; // x derivative t0 = tex2D(texSource, x - 2.0f * dx, y); t0 -= tex2D(texSource, x - 1.0f * dx, y) * 8.0f; t0 += tex2D(texSource, x + 1.0f * dx, y) * 8.0f; t0 -= tex2D(texSource, x + 2.0f * dx, y); t0 /= 12.0f; t1 = tex2D(texTarget, x - 2.0f * dx, y); t1 -= tex2D(texTarget, x - 1.0f * dx, y) * 8.0f; t1 += tex2D(texTarget, x + 1.0f * dx, y) * 8.0f; t1 -= tex2D(texTarget, x + 2.0f * dx, y); t1 /= 12.0f; Ix[pos] = (t0 + t1) * 0.5f; // t derivative Iz[pos] = tex2D(texTarget, x, y) - tex2D(texSource, x, y);

32. Решение уравнений Эйлера-Лагранжа • Ядро выполняют одну итерацию • Обновляются и , и • Каждый поток обрабатывает вектор для одного пикселя • Исходные данные общие для всех потоков • Возможен выигрыш от использования кеша

33. Одна итерация // handle borders if (ix != 0) left = pos - 1; else left = pos; … floatsumU = (u0[left] + u0[right] + u0[up] + u0[down]) * 0.25f; floatsumV = (v0[left] + v0[right] + v0[up] + v0[down]) * 0.25f; floatfrac = (Ix[pos] * sumU + Iy[pos] * sumV + Iz[pos]) / (Ix[pos] * Ix[pos] + Iy[pos] * Iy[pos] + alpha); u1[pos] = sumU - Ix[pos] * frac; v1[pos] = sumV - Iy[pos] * frac;

34. Современный подход • Классический метод хорошо работает только для небольших смещений • Линеаризация

35. Современный подход • Уменьшим изображения • «Большие» смещения станут «маленькими» • Решим задачу для уменьшенных изображений • Классический метод Horn-Schunck • Используем найденное решение как стартовую точку в исходной задаче

36. Современный подход • Можно рассмотреть пирамиду изображений • Каждое следующее меньше предыдущего • Используем решение задачи для текущего уровня как начальное приближение для следующего • Стартуем с самых маленьких изображений

37. Современный подход • Линеаризация работает для небольших смещений • Разобьем решение на текущем уровне на два слагаемыхрешение с предущего уровня

38. Деформация • На текущем уровне вместо будем работать с • Задача сводится к нахождению небольших смещений • Классический Horn-Schunck

39. Общая схема • Копирование изображений в память GPU • Построение пирамиды (коэффициент 0.5) • Начинаяснаименьшегоразрешения • Замена на его деформированную версию • Вычисление производных (для и ) • Нахождение • Масштабирование (переход к большему разрешению) • Копирование результатов в основную память

40. Деформация constint ix = threadIdx.x + blockIdx.x * blockDim.x; constintiy = threadIdx.y + blockIdx.y * blockDim.y; constintpos = ix + iy * stride; if (ix >= width || iy >= height) return; float x = ((float)ix + u[pos] + 0.5f) / (float)width; float y = ((float)iy + v[pos] + 0.5f) / (float)height; out[pos] = tex2D(texToWarp, x, y);

41. Современный подход • На каждом уровне изображение можно деформировать несколько раз • warping iterations

42. Современный подход • Создание пирамиды • Для всех уровней пирамиды начиная с наименьшего • Деформация • Вычисление производных • Нахождение • Перейти к следующему уровню Выполнить несколько итераций(warping iterations)

43. Warping iterations решение Размер изображения Warping iterations Продолжитьu, v Warping iterations Продолжитьu, v

44. Решение СЛАУОптимизация доступа к памяти • Выравнивание адреса начала строк матрицы • Использование текстур • Использование общей памяти мультипроцессора (shared memory)

45. Решение СЛАУОптимизация доступа к памяти • 640x480, 500 итераций • G92 (GeForce 8800GTX)

46. Решение СЛАУОптимизация доступа к памяти • 640x480, 500 итераций • GF104 (GeForce GTX 460)

50. Ресурсы нашего курса • Steps3d.Narod.Ru • Google Site CUDA.CS.MSU.SU • Google Group CUDA.CS.MSU.SU • Google Mail CS.MSU.SU • Google SVN • Tesla.Parallel.Ru • Twirpx.Com • Nvidia.Ru