Нахождение симметричных областей / образцов в изображении

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

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

Конечная цель - найти непрерывную, скорее всего, изогнутую продольную линию, которая адаптивно разделяет спину на симметричные «половинки».

Образец изображения размещен ниже. Обратите внимание, что не все области являются симметричными (в частности, чуть выше центра изображения, красная вертикальная «полоса» отклоняется вправо). Этот регион должен получить плохую оценку или что-то еще, но тогда локальная симметрия будет определяться из симметричных точек, расположенных дальше. В любом случае, мне придется адаптировать любой алгоритм к своей области приложения, но мне нужна стратегия сом корреляции / свертки / сопоставления с образцом, я думаю, что уже должно быть что-то вокруг.

(РЕДАКТИРОВАТЬ: есть еще изображения ниже, и еще несколько объяснений)

РЕДАКТИРОВАТЬ: в соответствии с просьбой, я буду включать более типичные изображения, либо с хорошим поведением и проблемные. Но вместо цветных картинок они являются полутоновыми, поэтому цвет напрямую связан с величиной данных, чего не было с цветным изображением (предназначенным только для связи). Хотя серые изображения, по-видимому, лишены контраста по сравнению с цветными, градиенты данных присутствуют, и при желании их можно повысить с помощью некоторого адаптивного контраста.

1) Изображение очень симметричного предмета:

2) Изображение одного и того же объекта в другой момент. Хотя есть больше «особенностей» (больше градиентов), он не «чувствует» себя таким симметричным, как раньше:

3) Тонкий молодой субъект с выпуклостями (костные выступы, обозначенные более светлыми областями) по средней линии вместо более обычной вогнутой средней линии:

4) Молодой человек с отклонением позвоночника, подтвержденным рентгеном (обратите внимание на асимметрию):

5) Типичный «наклоненный» предмет (хотя в основном симметричный относительно изогнутой средней линии и как таковой не должным образом «деформированный»):

Любая помощь приветствуется!

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

Я взломал быстрый алгоритм в Mathematica, чтобы продемонстрировать это. Это не алгоритм, который вы должны использовать в медицинском приложении, а только демонстрация основных идей.

Сначала я загружаю ваше изображение, отражаю его и делю эти изображения на маленькие блоки:

src = ColorConvert[Import["http://i.stack.imgur.com/jf709.jpg"], 
   "Grayscale"];
mirror = ImageReflect[src, Left -> Right];
blockSize = 30;
partsS = ImagePartition[src, {blockSize, blockSize}];
partsM = ImagePartition[mirror, {blockSize, blockSize}];
GraphicsGrid[partsS]

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

Если мы посмотрим на один блок и его зеркальное отражение:

{partsS[[6, 10]], partsM[[6, 10]]}

Мы видим, что они похожи, но смещены. Количество и направление сдвига - это то, что мы пытаемся выяснить.

Чтобы оценить сходство совпадений, я могу использовать квадрат евклидова расстояния:

ListPlot3D[
  ImageData[
   ImageCorrelate[partsM[[6, 10]], partsS[[6, 10]], 
    SquaredEuclideanDistance]]]

к сожалению, использование этих данных напрямую оказалось сложнее, чем я думал, поэтому я использовал вместо этого приближение 2-го порядка:

fitTerms = {1, x, x^2, y, y^2, x*y};

fit = Fit[
   Flatten[MapIndexed[{#2[[1]] - blockSize/2, #2[[2]] - 
        blockSize/2, #1} &, 
     ImageData[
      ImageCorrelate[partsM[[6, 10]], partsS[[6, 10]], 
       SquaredEuclideanDistance]], {2}], 1], fitTerms, {x, y}];

Plot3D[fit, {x, -25, 25}, {y, -25, 25}]

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

distancesFit = MapThread[
   Function[{part, template},
    Fit[Flatten[
      MapIndexed[{#2[[2]] - blockSize/2, #2[[1]] - blockSize/2, #1} &,
        ImageData[
        ImageCorrelate[part, template, 
         SquaredEuclideanDistance]], {2}], 1], 
     fitTerms, {x, y}]], {partsM, partsS}, 2];

Это дает нам наш первый энергетический термин для оптимизации:

variablesX = Array[dx, Dimensions[partsS]];
variablesY = Array[dy, Dimensions[partsS]];

matchEnergyFit = 
  Total[MapThread[#1 /. {x -> #2, y -> #3} &, {distancesFit, 
     variablesX, variablesY}, 2], 3];

variablesX/Yсодержит смещения для каждого блока и matchEnergyFitприближает квадрат евклидовой разницы между исходным изображением и зеркальным отражением с примененными смещениями.

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

Итак, мы установили второй энергетический член для гладкости:

smoothnessEnergy = Total[Flatten[
    {
     Table[
      variablesX[[i, j - 1]] - 2 variablesX[[i, j]] + 
       variablesX[[i, j + 1]], {i, 1, Length[partsS]}, {j, 2, 
       Length[partsS[[1]]] - 1}],
     Table[
      variablesX[[i - 1, j]] - 2 variablesX[[i, j]] + 
       variablesX[[i + 1, j]], {i, 2, Length[partsS] - 1}, {j, 1, 
       Length[partsS[[1]]]}],
     Table[
      variablesY[[i, j - 1]] - 2 variablesY[[i, j]] + 
       variablesY[[i, j + 1]], {i, 1, Length[partsS]}, {j, 2, 
       Length[partsS[[1]]] - 1}],
     Table[
      variablesY[[i - 1, j]] - 2 variablesY[[i, j]] + 
       variablesY[[i + 1, j]], {i, 2, Length[partsS] - 1}, {j, 1, 
       Length[partsS[[1]]]}]
     }^2]];

К счастью, в Mathematica встроена ограниченная оптимизация:

allVariables = Flatten[{variablesX, variablesY}];
constraints = -blockSize/3. < # < blockSize/3. & /@ allVariables;
initialValues = {#, 0} & /@ allVariables;
solution = 
  FindMinimum[{matchEnergyFit + 0.1 smoothnessEnergy, constraints}, 
   initialValues];

Давайте посмотрим на результат:

grid = Table[{(j - 0.5)*blockSize - dx[i, j], (i - 0.5)*blockSize - 
      dy[i, j]}, {i, Length[partsS]}, {j, Length[partsS[[1]]]}] /. 
   solution[[2]];
Show[src, Graphics[
  {Red,
   Line /@ grid,
   Line /@ Transpose[grid]
   }]]

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

Возможные улучшения:

• Как я уже сказал, сначала выполните жесткую регистрацию. На белом фоне простая регистрация на основе моментов изображения должна работать нормально.
• Это только один шаг. Вы можете использовать смещения, которые вы нашли за один шаг, и улучшить их на втором шаге, возможно, с меньшим окном поиска или меньшим размером блока
• Я читал статьи, где они делают это вообще без блоков, но оптимизируют смещение на пиксель.
• Попробуйте разные функции сглаживания

Интересный вопрос. Во-первых, возможно, вам нужны подходы, основанные на детекторе ключевых точек интереса и сопоставлении. Это может включать SIFT (масштабно-инвариантное преобразование объектов), SURF, ORB и т. Д. ... или даже более простой подход, основанный исключительно на операторе Харриса (csce.uark.edu/~jgauch/library/Features/Harris.1988.pdf ). Из вашего поста непонятно, что вы пробовали, поэтому извините, если я здесь наивен.

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

Я взял ваше примерное изображение и преобразовал его в цветовое пространство L a b *, используя ImageMagick, и использовал только полосу L *:

convert x.jpg -colorspace Lab -separate %d.png

0.png соответствует полосе L *. Теперь я уверен, что у вас есть фактические данные изображения, но я имею дело с артефактами сжатия jpg, а с какими нет. Чтобы частично решить эту проблему, я выполнил морфологическое открытие с последующим морфологическим закрытием с плоским диском радиуса 5. Это основной способ уменьшения шума с ММ, и, учитывая радиус диска, незначительная часть изображения изменяется. Затем моя идея была основана на этом единственном изображении, которое имеет большой шанс потерпеть неудачу в других случаях. Ваш интересующий регион визуально выделяется тем, что он темнее («горячее» в цветном изображении), поэтому я предположил, что статистический бинаризатор может работать хорошо. Я использовал подход Оцу, который является автоматическим.

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

Если вы все еще следуете за мной, мы еще не нашли фактическую предполагаемую центральную область интереса. Вот мой взгляд на это. Независимо от того, насколько изогнут человек (на самом деле я вижу некоторые проблемные случаи), регион напоминает вертикальную линию. С этой целью я упрощаю текущее изображение, выполняя морфологическое раскрытие с вертикальной линией длиной 100. Эта длина является чисто произвольной, если у вас нет проблем с масштабированием, то это не сложно определить. Теперь мы снова отбрасываем компоненты, но я был немного более осторожен на этом этапе. Я использовал открытие по области с дополнением изображения, чтобы отбросить то, что я считал маленькими областями, это можно сделать более контролируемым образом, выполнив что-то в форме гранулометрического анализа (из ММ тоже).

Теперь у нас примерно три части: левая часть изображения, центральная часть и правая часть изображения. Ожидается, что центральная часть будет меньшим компонентом из трех, поэтому она получается тривиально.

Вот окончательный результат, нижнее правое изображение - это просто наложенное изображение слева от исходного. Отдельные цифры не все выровнены, извините за спешку.