Композитное объединение видимого и инфракрасного изображений получаемых с тепловизора
Как известно, большинство современных тепловизоров имеют два канала: TV – формирующий цветное изображение в видимом диапазоне и IR – визуализирующий собственно распределение температур. Во многих задачах тепловизионного обследования возникает необходимость анализа композитного изображения (когда на видимое изображение накладывается инфракрасное). Правда не все модели тепловизоров предоставляют такую возможность. Как правило, наличие функции синтеза композитного кадра – существенно увеличивает стоимость тепловизионного прибора. Поэтому возникает закономерный вопрос: а возможно ли на видимое изображение наложить инфракрасное, иными средствами, не используя функционал тепловизора? Оказывается можно. О том как это делается средствами Wolfram Mathematica и будет рассказано в настоящей заметке. Стоит сразу же заметить, что это не единственный способ, но один из самых простых.
Исходная задача. В нашем научно-исследовательском проекте возникла потребность объединить координатные системы TV и IR каналов тепловизора, а также иметь нестандартное композитное изображение (с отличными от типовых правилами наложения инфракрасного изображения на видимое). По условиям съёмки – сцена статичная, т.е. камера и объект не перемещаются относительно друг друга. Статичность существенно облегчает задачу, так как преобразование координат вычисленное единожды, затем применяется ко всему тепловизионному фильму. К сожалению, преобразование координат существенно меняется от сцены к сцене (при изменении взаимного положения тепловизора и объекта съёмки), поэтому его приходится перевычислять каждый раз заново.
Поставленную задачу мы решили при помощи реперного объекта, помещённого в поле зрения тепловизора. Существенных требований к реперному объекту несколько: измерительная (геометрическая) стабильность; контраст относительно объектов сцены, как в видимом, так и в инфракрасном диапазонах; стерилизуемость и биологическая безопасность – весьма жёсткое требование, учитывая что основной объект наблюдаемой сцены, это открытая кора живого головного мозга живого человека. Вычисление координатного преобразования естественно должно происходить автоматически – с использованием методов компьютерного зрения (computer vision).
Итак, вначале зададим пути к файлам с данными и типовые имена файлов:
$HistoryLength = 0;
FrameTVFileName = "measuring_frame_TV.jpg";
FrameIRFileName = "measuring_frame_IR.bmp";
FrameIRdFileName = "measuring_frame_IR.txt";
PathFiles = "...\\BrainIR";
Затем загрузим данные и проведём их первичные преобразования:
imgIR = ImageTake[Import[FileNameJoin[{PathFiles, FrameIRFileName}], "Image"],
{1, 240}, {1, 320}];
imgIRd = Image[ToExpression@StringReplace[#, "," -> "."] &
/@Rest@Import[FileNameJoin[{PathFiles, FrameIRdFileName}],
"Table"] // Transpose // Rescale, ColorSpace -> "Gray"];
imgTV = ImageTake[Import[FileNameJoin[{PathFiles, FrameTVFileName}], "Image"],
{701, 1400}, {601, 1500}];
Отобразим результат загрузки. На рисунке ниже изображена сцена и реперный объект в различных каналах (слева направо): imgIR – инфракрасное изображение; imgIRd – температурная карта; imgTV – цветное изображение в видимом диапазоне. И хотя на рисунке изображения почти одного размера, в действительности их размер весьма отличается: 320х240 – ИК канал и температурная карта; 900х700 – TB канал (в оригинале 2048х1536).
GraphicsGrid[{{imgIR, imgIRd, imgTV}}, ImageSize -> 693]
Далее преобразуем цветной TV канал в HSB модель и посмотрим на результат. Как видно из рисунка, наиболее контрастный и информативный канал – средний, это saturation (насыщенность).
imgHSBChannelTV = ColorSeparate@ColorConvert[imgTV, "HSB"];
GraphicsGrid[{imgHSBChannelTV}, ImageSize -> 693]
Здесь необходимо отметить, что информативность канала насыщенности мы конечно же оцениваем не «на глаз», а посредством вычисления и анализа гистограммы локальной вариации яркости (но это немного далеко от нашей темы, поэтому соответствующий код не приводим). Кроме того, возможно также использовать пространство XYZ и применить байесовский подход, аналогичный изложенному в этом посте.
Далее начинается магия Computer Vision.
Выделим контур интересующего нас реперного объекта и посмотрим на результат: IR – слева; TV – справа.
ContrastIRd = 1;
ThresholdBinIRd = 0.81;
KerrIRd = CrossMatrix[7];
ContrastTV = 0.45;
ThresholdBinTV = 0.19;
KerrTV = CrossMatrix[5];
imgDetectIRd = SkeletonTransform@TopHatTransform[
ColorNegate@Binarize[ImageAdjust[imgIRd, ContrastIRd], ThresholdBinIRd]
, KerrIRd];
imgDetectTV = SkeletonTransform@TopHatTransform[
ColorNegate@Binarize[ImageAdjust[imgHSBChannelTV[[2]], ContrastTV], ThresholdBinTV]
, KerrTV];
GraphicsGrid[{{imgDetectIRd, imgDetectTV}}, ImageSize -> 693]
Здесь необходимо отметить, что вышеуказанные параметры в заметке приведены в виде констант для простоты понимания сути идеи, в оригинальном алгоритме (который эксплуатируется нами) они автоматически адаптивно подстраиваются, для наилучшего выделения реперных контуров.
Идентифицируем опорные линии реперного контура и посмотрим на результат:
ThresholdLineIRd = 0.22;
DistinctnessLineIRd = 0.1;
ThresholdLineTV = 0.12;
DistinctnessLineTV = 0.1;
linesIRd = ImageLines[imgDetectIRd, ThresholdLineIRd, DistinctnessLineIRd];
linesTV = ImageLines[imgDetectTV, ThresholdLineTV, DistinctnessLineTV];
GraphicsGrid[{{
Show[imgIRd, Graphics[{Orange, Thick, Line /@ linesIRd}], ImageSize -> 320],
Show[imgTV, Graphics[{Yellow, Thick, Line /@ linesTV}], ImageSize -> 320]}}
, ImageSize -> 693]
Здесь вышеуказанные параметры также приведены в виде констант, а в оригинальном алгоритме они автоматически адаптивно изменяются, для наилучшего выделения опорных линий реперных контуров.
Выделим точки пересечений опорных линий и отфильтруем корректные точки:
LineParamTV = LineParam[linesTV];
LineParamIRd = LineParam[linesIRd];
IntersPointTVi = Cases[InterSectionPoint[LineParamTV], {x_, y_, \[Gamma]_} /;
1.36 < Abs[\[Gamma]] < 1.6 -> {x, y, \[Gamma]}];
IntersPointIRdi = Cases[InterSectionPoint[LineParamIRd], {x_, y_, \[Gamma]_} /;
1.45 < Abs[\[Gamma]] < 1.6 -> {x, y, \[Gamma]}];
CoordIntersPointTVi = Most/@IntersPointTVi;
CoordIntersPointIRdi = Most/@IntersPointIRdi;
IteratorIntersPoint = Range@4;
GraphicsGrid[{{
Show[imgIRd, Graphics[{Orange, Thick, Line@linesIRd,
Blue, PointSize@Large, Point/@CoordIntersPointIRdi,
Text[Style[#, FontSize -> 12], CoordIntersPointIRdi[[#]] - {12, 0} ] & /@IteratorIntersPoint}],
ImageSize -> 320],
Show[imgTV, Graphics[{Yellow, Thick, Line@linesTV,
White, PointSize@Large, Point/@CoordIntersPointTVi,
Text[Style[#, FontSize -> 12], CoordIntersPointTVi[[#]] - {21, 0} ] & /@ IteratorIntersPoint}],
ImageSize -> 320]}}
, ImageSize -> 693]
Код функций LineParam[...] и InterSectionPoint[...] приведён ниже, в не самой оптимальной форме, но зато понятно и читабельно:
LineParam[Lines_] := Module[{kl, bl},
kl = Function[coords,
(Subtract @@ Last@coords)/(Subtract @@ First@coords)][Reverse /@ Transpose@#] & /@ Lines;
bl = MapThread[Function[{coords, kll},
Last@coords - kll First@coords][First@#1, #2] &, {Lines, kl}];
Thread@{kl, bl}];
InterSectionPoint[LineParam_] := Module[{kl, bl, Numberlines},
Numberlines = Length@LineParam;
kl = LineParam[[All, 1]];
bl = LineParam[[All, 2]];
Flatten[Table[
If[kl[[i]] != kl[[j]],
{
(bl[[i]] - bl[[j]])/(kl[[j]] - kl[[i]]),
(kl[[j]] bl[[i]] - kl[[i]] bl[[j]])/(kl[[j]] - kl[[i]]),
ArcTan[(kl[[j]] - kl[[i]])/(1 + kl[[i]] kl[[j]])]
}
,
{\[Infinity], \[Infinity], 0}]
, {i, 1, Numberlines - 1}, {j, i + 1, Numberlines}], 1]
];
Как видно из рисунка выше, порядок нумерации точек в каналах IR и TV не совпадает. И он зачастую разный, так как детектирование и фильтрация паразитных линий и точек пересечения вносят своё «чёрное» дело. Поэтому нормализуем порядок следования точек пересечения:
IntersPointTV = ReorderingPoint@IntersPointTVi;
IntersPointIRd = ReorderingPoint@IntersPointIRdi;
CoordIntersPointTV = Most /@ IntersPointTV;
CoordIntersPointIRd = Most /@ IntersPointIRd;
Код функции ReorderingPoint[...] – не сложный, но интересный. Поэтому его реализацию оставляем читателям в качестве домашнего задания.
На этом магия компьютерного зрения заканчивается, и начинается магия проективной геометрии. По полученным ключевым точкам реперного объекта найдём матрицу преобразования координат из координатной системы IR канала в координатную систему TV канала. На Wolfram Mathematica это делается в одну функцию:
{error, ftr} = FindGeometricTransform[CoordIntersPointTV, CoordIntersPointIRd, Method ->"FindFit"]
Чтобы визуально оценить нетривиальность преобразований построим координатную сетку:
Show[{
Graphics[{Rectangle[{1, 1}, ImageDimensions@imgTV]}],
ParametricPlot[Evaluate@Array[ftr[{#, j}] &, 20, {1, 320}], {j, 1, 240}, PlotStyle -> Green],
ParametricPlot[Evaluate@Array[ftr[{i, #}] &, 20, {1, 240}], {i, 1, 320}, PlotStyle -> Green]
}, Axes -> None]
Точность совмещения каналов тепловизора можно оценить на рисунке ниже, где приведено простейшее композитное изображение с направляющими линиями реперного объекта. Из рисунка видно фактически попиксельное совпадение видимого изображения с инфракрасным.
Вместо заключения. Данная заметка демонстрирует некоторые приёмы Computer Vision при решении реальной задачи: объединения координатных систем TV и IR каналов тепловизора. И показана реализация этих приёмов в среде Wolfram Mathematica.
Обсудить данную заметку Вы можете на площадке Русскоязычной поддержки Wolfram Mathematica по ссылке.
Данный документ распространяется на условиях лицензии: Creative Commons «Attribution-NonCommercial-ShareAlike» («Атрибуция-Некоммерческое использование-На тех же условиях») 3.0 Непортированная.
16 ноября 2013 года.
Андрей Макаренко,
группа «Конструктивная Кибернетика».
Обсуждение: contact@rdcn.ru
Ключевые слова: computer vision, композитный кадр, тепловизор, проективные преобразования, Wolfram Mathematica.
