Всем привет! В этой статье собран материал по дифференциальным методам оптического потока, а также попытка, наконец, разобраться: что это такое, почему так, как так.

Если у вас тоже есть непонимание, как с этим работать и что стоит внутри непонятного потока, то давайте разбираться вместе.

База

В первую очереь посмотрим на определение:

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

Для каждой точки находится сдвиг, таким образом чтобы исходной точке I(x, y) на 1 изображении соответствовала точка на втором изображении I(x+dx, y+dy).

Оптический поток может считаться для:

  • всех точек (плотный поток)

  • заданных точек (выборочный, разреженный). Заданные точки можно, например, получить с помощью feature detector, например, SIFT, детектор Харриса. Для тех, кто не знает, что такое feature detector (мы говорим про классическое компьютерное зрение), это некий детектор (ого, так неожиданно), который позволяет выделять интересные точки на изображении. Это могут быть углы зданий, перекрестки дорог, точки гже происходит резкое изменение яркости, края.

Как определить соответствие точек? Методы для определения оптического потока (определение движения объекта) делятся на следующие категории:

  1. Корреляция фаз

  2. Блочные методы

  3. Дифференциальные методы

  4. Общие вариационные методы

  5. Методы дискретной оптимизации

Толка от перечисления методов нет, если мы не понимаем их смысл. Обо всем попорядку. Было бы слишком громоздно рассматривать все методы и постараться копнуть вглубь. Поэтому в этой статье хотелось бы затронутьтолько пункт 3. А для остальных подготовить отдельные статьи.

Дифференциальные методы

Метод Лукаса-Канаде

Начальное условие: один кадр переходит в другой без изменений (очень реалистично), то есть интенсивность изображения (освещение) не меняется.

Переходим к великой и могучей матетматике. То что мы сказали выше представим как:

Считаем, что перемещение небольшое, поэтому мы можем применить разложение в ряд Тейлора!

I(x, y, t) переносим в левую часть, получаем:

Отсюда можно получить следующее

Это уравнение можно представить следующим образом:

Таким образом:

Или как скалярное произведение градиента по T и вектора смещения V (Vx, Vy)

Полученное уравнение содержит две неизвестных (Vx, Vy) и не может быть однозначно разрешено. Что будем делать? Давайте сделаем предположение, что соседние пиксели смещаются на одинаковое расстояние.

Теперь мы можем взять некий фрагмент, например, 22, 33, 44, 55 и для каждого пикселя внутри этого фрагмента - (Vx, Vy) будет равно. Благодаря этому из одного уравнения мы можем получить гораздо больше. Например, для фрагмента 2*2, получим целых 4 уравнения!

Без этого предположения однозначно решить уравнение мы не можем, но зато можем минимизировать ошибку.

Здесь g — функция, определяющая весовые коэффициенты для пикселей. Обычно применяют двухмерная гауссиана, которая дает наибольший вес центральному пикселю и все меньший по мере удаления от центра.

Для нахождения минимума можно использовать метод наименьших квадратов. Можно было бы все расписать, но зачем, когда есть эта замечательная статья.

У этого метода много недостатоков (допущений). С допущением о достаточности первой производной, можно побороться с помощью итеративного вычисления. Для каждой пары изображений итеративно применяем алгоритм. То есть использовали алгоритм для изображений х1, х2, а на следующей итерации вместо х2, х3 используем не х2 оригинальный, а скомпенсированный смещением, вычисленный на 1 итерации.

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

Как будем решать проблемы работы алгоритма только на маленьких смещениях? 2 слова - Многомасштабный анализ (multi-scaling). Построим пирамиду изображений разных масштабов.

  1. Построение пирамиды изображений: Создается набор изображений разного размера, где каждое изображение вдвое меньше предыдущего.

  2. Расчет оптического потока по пирамиде:

    • Начинается с самого маленького изображения, где смещение пикселей минимально.

    • Результат на маленьком изображении масштабируется и используется как начальное значение для следующего, более крупного изображения.

    • Процесс повторяется до достижения исходного размера изображения.

Немного кода

Построим оптический поток на видео с котом.

# Выбор точек для отслеживания
points = cv2.goodFeaturesToTrack(prev_gray, 100, 0.01, 10)

def optical_flow_pyr_lk(prev_frame, curr_frame, points):
  """
  Вычисляет оптический поток с помощью метода Лукаса-Канаде с пирамидальным подходом.

  Args:
    prev_frame: Предыдущий кадр.
    curr_frame: Текущий кадр.
    points: Список точек, для которых нужно найти оптический поток.

  Returns:
    Список точек с обновленными координатами,  и список статусов (0 - не найдено, 1 - найдено)
  """

  # Параметры метода
  lk_params = dict( winSize  = (15,15),
                    maxLevel = 2,
                    criteria = (cv2.TERM_CRITERIA_EPS | cv2.TERM_CRITERIA_COUNT, 10, 0.03))

  # Преобразование точек в формат NumPy
  points = np.array(points, dtype=np.float32)

  # Вычисление оптического потока
  next_points, status, error = cv2.calcOpticalFlowPyrLK(
      prev_frame, curr_frame, points, None, **lk_params)

  # Возврат результата
  return next_points, status

Метод Хорна-Шунка:

У метода Лукаса-Канаде как мы убедились много недосатков, рассмотрим другой алгоритм. Метод опирается на предположение о том, что на всем изображении оптический поток будет достаточно гладким. В результате этого предположения мы добавляем требование на отсутствие резкого изменения сдвигов с весовым коэффициентом. Обратимся к математике. Мы хотим перейти от этого:

к этому

Выглядит страшно, но это только кажется. На самом деле все довольно просто, мы просто добавляем сглаживающий коэффициент. Обо все по порядку. Мы добавили условие гладкости, которое выглядит следующим образом:

Сглаживающий член в методе Хорна-Шунка выражается как сумма квадратов производных, потому что такая форма наказывает резкие изменения. Квадрат производной больше для больших изменений, а мал для небольших. Таким образом, функционал будет стремиться к более гладким решениям, где изменения оптического потока плавные.

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

Метод решения

Для минимизации функционала используется вариационное исчисление, приводящее к системе уравнений Эйлера-Лагранжа. Полученные уравнения можно решать итеративно:

Стоп! Ничего не понятно

Для начала, что такое уравнение Эйлера-Лангранда? Уравнение Эйлера-Лангранжа выглядит следующим образом:

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

Окей, поехали. Вспомним, что мы минимизируем:

А теперь:

  1. Возьмем частную производную ланграндиана по

  2. Возьмем частную производную по ux

  3. Возьмем частную производную по uy

  4. Возьмем производные по х, у

А теперь подставим это в уравнение Эйлера-Лангранжа и получим:

Получаем систему уравнений, которую записываем для каждого пикселя и решаем общую систему итеративно. Как в предыдущем методе. В данном алгоритме тоже предлагают использовать multi-scaling, причем рекомендуют масштабировать изображения не в 2 раза, а с коэффициентом 0.65

Farneback

Алгоритм был придуман в 2003 году. Основная идея метода состоит в аппроксимации локального окрестности изображения второй степенью полинома (до этого мы строили с одной степенью). Для каждого пикселя изображения строится локальная модель в виде квадратичной поверхности I = xAx + bx + c с симметричной матрицей A, описывающей изменения интенсивности.

Если изображение сдвинулось, то квадратичная аппроксимация изменится. Подставляя новые координаты в уравнение, можно раскрыть скобки и выразить новое положение пикселей через старое: I2(x)=I1(x-d), решаем систему уравнений и получаем:

Коэффициенты A, b и c вычисляются для обоих изображений (до и после сдвига). Затем можно решить уравнения для определения сдвига (оптического потока) d. В реальных условиях изображение часто содержит шум, поэтому предполагается, что все пиксели в небольшой окрестности имеют похожий сдвиг. Это позволяет сгладить вычисления и улучшить точность оценки.

Для компенсации шума вычисляется ошибка аппроксимации для каждого пикселя, которая затем усредняется по окну с гауссовыми весами (более близкие пиксели имеют больший вес)

import cv2
import numpy as np

def calculate_optical_flow_farneback(prev_frame, curr_frame):
    """
    Вычисляет оптический поток методом Farneback.

    Args:
        prev_frame: Предыдущий кадр.
        curr_frame: Текущий кадр.

    Returns:
        Векторное поле оптического потока (u, v) в виде NumPy массива.
    """

    # Преобразуем кадры в оттенки серого
    prev_gray = cv2.cvtColor(prev_frame, cv2.COLOR_BGR2GRAY)
    curr_gray = cv2.cvtColor(curr_frame, cv2.COLOR_BGR2GRAY)

    # Вычисляем оптический поток методом Farneback
    flow = cv2.calcOpticalFlowFarneback(
        prev_gray, curr_gray, None, 
        0.5, 3, 15, 3, 5, 1.2, 0 
    )

    return flow

# Загрузка видеофайла
cap = cv2.VideoCapture("video.mp4")
visual = []
# Получение первого кадра
ret, prev_frame = cap.read()

while(True):
    # Получение текущего кадра
    ret, curr_frame = cap.read()
    if not ret:
        break

    # Вычисление оптического потока
    flow = calculate_optical_flow_farneback(prev_frame, curr_frame)

    # Визуализация потока 
    for y in range(0, flow.shape[0], 10):
        for x in range(0, flow.shape[1], 10):
            fx, fy = flow[y, x]
            cv2.arrowedLine(curr_frame, (x, y), (x + int(fx), y + int(fy)), (0, 255, 0), 1)

    # Обновление предыдущего кадра
    visual.append(curr_frame)
    prev_frame = curr_frame

# Закрытие окна и освобождение ресурсов
cv2.destroyAllWindows()
cap.release()

Возьмем изображения по крупнее

Применение

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

Также оптический поток можно использовать для стабилизации видео, для сжатия. Оптический поток может использоваться для анализа движения глазных яблок и диагностики заболеваний глаз, для отслеживания спорстменов и объектов, если речь про приминение в спорте. Для анализа поведения животных и многое другое.

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

Спасибо за прочтения статьи!

Литература:

Комментарии (7)


  1. morheus9
    04.06.2024 09:03

    Я один прочитал как "оптический лоток"?


  1. Jury_78
    04.06.2024 09:03

    Не проще использовать Event camera?


    1. Kit_Cat Автор
      04.06.2024 09:03

      Да, проще, статья была написана с целью помочь разобраться или узнать новое для тех, кто хочет понять, как оптический поток работает


      1. Jury_78
        04.06.2024 09:03
        +1

        Да, статья нужная, сам сталкивался с этим. Но по мне название "оптический поток", что по-русски, что по-английски звучит неподходяще.


  1. tminnigaliev
    04.06.2024 09:03
    +1

    Спасибо за статью! Добавил в закладки.

    Заметил следующие опечатки:

    Диффиренциальные методы --> Дифференциальные методы

    Для нахождения минимум можно использовать метод наименьших квадрат. --> квадратов

    Как будем решать проблемы работы алгоритма только на маленьких смещениях? 2 слова - Многомасштабного анализа (multi-scaling). Пропущено "методом".


    1. tminnigaliev
      04.06.2024 09:03
      +1

      Во втором фрагменте ещё "минимум" --> "минимума"


      1. Kit_Cat Автор
        04.06.2024 09:03

        Спасибо за замечания! Опечатки моя вечная проблема