Невязки в игрушечном видеоизмерителе
Ещё одна вещь, до которой не доходили руки до поры до времени. Есть видеоизмеритель параметров сближения, он получает изображение уголковых отражателей и находит такие 6 параметров взаимного расположения (3 поступательных и 3 угловых), при которых его МОДЕЛЬ ИЗМЕРЕНИЯ (то, как по его мнению, формируется изображение при конкретном положении уголковых отражателей в пространстве) выдаст изображение, максимально близкое к полученному в реальности.
Максимально близкое - но всё же без точного совпадения. Мы отдаём себе отчёт, что все измерения неизбежно "зашумлены". Узнать уровень шума АПРИОРИ очень сложно, поскольку он состоит из нескольких компонентов. Это и точность изготовления мишени, и аберрации оптической системы, и, самое главное, "геометрические шумы" матрицы, в первую очередь неравномерность чувствительности пикселей (PhotoResponse Non-Uniformity, или PRNU). Это КМОП-матрица, в каждом пикселе содержится свой полевой транзистор, на затворе которого накапливается сигнал, а на истоке он "повторяется", чтобы пойти в АЦП и не утеряться в процессе. Так что параметры пикселей могут прилично "гулять", причём на разных температурах, и просто в течение времени - по-разному. А потом ещё добавить сюда "износ" матрицы радиацией, всякие фокусы с пульсацией напряжения питания, потом ещё добавить шум квантования и прочие "пакости" от АЦП, и это всё будет меняться в зависимости от светотехнической обстановки, в которой мы будем работать...
Поэтому давным-давно я решил вообще не связываться с априорной оценкой шума, а определять его как раз-таки по неустранимой невязке, т.е по отличию реального изображения от наиболее близкого "модельного", соответствующего измеренным параметрам сближения. И сделал там одно допущение, которое толком не проверял.
Сейчас на " игрушечном видеоизмерителе" хочу понять, насколько это допущение было верно...
Мы снова поработаем в одной плоскости. Мишень на этот раз состоит из 3 точек, A1 (1;1), A2 (1;0) и A3 (1;-1).
Линза видеоизмерителя (точнее говоря, передняя главная точка) сидит в начале координат, (0;0), и проецирует эти точки на плоскость x=-1. Как обычно, игнорируем, что картинка переворачивается (поскольку считываем мы её как положено), и мы должны получить на "фотоприёмной линейке" пятна в координатах G1 = 1, G2 = 0, G3 = -1.
Но из-за разных шумов (то ли мишень плохо изготовили, то ли пятно "провалилось" между пикселями) мы получили результаты F1 = 1 + ξ, F2 = η, F3 = -1 + χ, где ξ,η,χ - независимые случайные величины с нулевым средним и одинаковой дисперсией, равной σ2. Модель немного упрощённая, но пока нас устроит.
И сейчас начнём варьировать число измеренных параметров.
Вырожденный случай: мы и не собирались измерять параметры, мы считаем, что поместили мишень во вполне конкретную точку, и просто "калибруем" нашу оптику, и заодно проверяем уровень шумов. Модельные координаты: G1 = 1, G2 = 0, G3 = -1. Разность между реальными и модельными: ξ,η и χ, соответственно.
Неустранимой невязкой назовём среднеквадратичное значение, т.е
^2+(F_2-G_2)^2+(F_3-G_3)^2}{3}}=\sqrt{\frac{\xi^2+\eta^2+\chi^2}{3}})
С вычислением математического ожидания этой невязки есть некоторые проблемы, а вот математическое ожидание выражения под корнем посчитается легко:
=\frac{E(\xi^2)+E(\eta^2)+E(\chi^2)}{3})
Поскольку у каждой из наших случайных величин нулевое математическое ожидание, наши слагаемые - не что иное, как дисперсии. Получаем:
=\frac{\sigma^2+\sigma^2+\sigma^2}{3}=\sigma^2)
Напишем так:
\approx\sigma)
Точного равенства не получается, но и не очень-то хотелось. Нас бы и дисперсия устроила, корень извлекаем просто для удобной размерности. Невязка в 1/10 пикселя - это понятно. А в 1/100 кв. пикселей - уже не очень...
[Пример, почему нельзя поставить знак равенства]
Обойдёмся одной случайной величиной ξ. И пусть она принимает лишь 3 значения: -1, 0 и +1, каждое с вероятностью 1/3. Найдём сначала корень из математического ожидания квадрата такой величины:
}=\sqrt{(-1)^2/3+0^2/3+1^2/3}=\sqrt{2/3})
А теперь наоборот, математическое ожидание от корня:
=E(|\xi|)=|-1|/3+|0|/3+|1|/3=2/3)
Впрочем, когда у нас складывается много случайных величин, результаты становятся уже вполне сравнимыми. Только что для очистки совести в экселе помучал исходный пример с 3 случайными величинами, взял их равномерно распределёнными между -1 и +1. Получил отличие в 3..5%, а для настоящего видеоизмерителя, где 8 точек и для каждой по 2 координаты, разница и того меньше.
Ожидание вполне ожидаемое: находя эту неустранимую невязку, мы как раз оцениваем среднеквадратичное значение шума.
Случай посложнее: измеряем только поперечное смещение, то есть смещение по оси Y. Считаем, что приближаться или отдаляться мишень не может, и поворачиваться тоже.
Вся хитрючая математика здесь практически "вырубается". Частные производные по одному-единственному параметру сближения оказываются единичными, поэтому ненормированная оценка "вектора параметров" 1×1 равна
+(F_2-G_2)+(F_3-G_3)=\xi+\eta+\chi)
Нормирующая матрица будет иметь размеры 1×1 и станет равна:
J=3.
Обратная к ней: J-1 = 1/3.
Итого, нормированная оценка сдвига:
Кто бы мог подумать: среднее арифметическое...
Теперь, получив такой "сдвиг", мы полагаем, что координаты точек на фотоприёмной линейке станут равны
Поэтому невязки составят:
Найдём их сумму квадратов:
=\frac{2}{3}(\xi^2+\eta^2+\chi^2-\xi\eta-\xi\chi-\eta\chi))
Когда будем искать математическое ожидание, последние 3 слагаемых дадут ноль, поскольку случайные величины имеют нулевое мат. ожидание и независимы друг от друга. Останется следующее:
=\frac{2}{3}(\sigma^2+\sigma^2+\sigma^2)=2\sigma^2)
И наконец, оцениваем математическое ожидание нашей "неустранимой невязки":
\approx\sqrt{\frac{E(h_1^2+h_2^2+h_3^2)}{3}}=\sqrt{\frac{2}{3}}\cdot\sigma)
Мы обнаруживаем, что неустранимая невязка уменьшилась!
Это произошло из-за того, что нам "удалось" отыскать такой сдвиг, при котором картинки лучше соответствуют друг другу. Кстати, именно поэтому в среднеквадратичном отклонении ДЛЯ ВЫБОРКИ в знаменателе стоит N-1. По сути, ровно этот механизм мы и увидели: ведь когда мы считаем это среднеквадратичное отклонение, мы не знаем истинного среднего, его мы тоже ОЦЕНИВАЕМ, и поэтому если у нас все слагаемые случайно сдвинулись в одну сторону, мы это уже не почувствуем, это у нас войдёт в оценку среднего значения, а не в среднеквадратичное отклонение.
Так что наша "неустранимая невязка", если её напрямую интерпретировать как шум каждого измерения, будет давать заниженные ответы. Правильнее будет считать её так (если мы по 3 точкам измеряем ОДИН параметр сближения):
Поехали дальше: ДВА параметра сближения, сдвиги по осям X (дальность) и Y (поперечное смещение).
Частные производные по поперечному смещению, как и раньше, остались единичными. А по продольному получаем -1, 0, +1 для точек (1;1), (1;0) и (1;-1) соответственно. То есть, если дальность увеличивается, мы обнаружим, как верхняя точка ползёт вниз, центральная остаётся на одном месте, и нижняя ползёт вверх, все они "стягиваются" в ноль, т.е "в горизонт".
Найдём ненормированную оценку сдвига:
\xi+\left(\begin{array}{c}0\\1\end{array}\right)\eta+\left(\begin{array}{c}1\\1\end{array}\right)\chi=\left(\begin{array}{c}\chi-\xi\\\xi+\eta+\chi\end{array}\right))
Матрица J получается довольно-таки простой, без перекрёстных связей, за счёт симметричности мишени и её расположения:
)
Наконец, находим нормированную оценку сдвига:
=\left(\begin{array}{cc}\frac{1}{2}&0\\0&\frac{1}{3}\end{array}\right)\left(\begin{array}{c}\chi-\xi\\\xi+\eta+\chi\end{array}\right)=\left(\begin{array}{c}\frac{\chi-\xi}{2}\\\frac{\xi+\eta+\chi}{3}\end{array}\right))
Чтобы посчитать неустранимую невязку, нужно "применить" оба этих сдвига. Будем считать, что шумы у нас всё-таки маленькие, так что можно здесь работать с "линеаризованной" задачей. Тогда получим такие "ожидаемые" положения пятен на фотоприёмной линейке:
В результате получим невязки:
Мы уже сообразили, что при возведении каждой невязки в квадрат можно сразу оставить лишь квадраты случайных величин, а где одна умножается на другую - вычёркиваем, т.к мат. ожидание равно нулю. Для начала посчитаем мат. ожидания квадратов отдельных невязок:
=E(h_3^2)=\frac{1}{36}E(\xi^2)+\frac{1}{9}E(\eta^2)+\frac{1}{36}E(\chi^2)=\frac{1}{6}\sigma^2,)
=\frac{4}{9}E(\eta^2)+\frac{1}{9}E(\xi^2)+\frac{1}{9}E(\chi^2)=\frac{2}{3}\sigma^2)
Разница в дисперсии в 4 раза между крайними и центральной точками! По среднеквадратичному значению - в 2 раза. Исходя из "теории оптимального оценивания", мы догадываемся, что правильнее было бы складывать квадраты невязок с разными весами, гораздо более "пристальное внимание" обратив на центральную точку.
Но допустим, что мы снова просто складываем все квадраты невязок. Тогда получим математическое ожидание от их суммы квадратов:
=\frac{1}{6}\sigma^2+\frac{2}{3}\sigma^2+\frac{1}{6}\sigma^2=\sigma^2)
А это мы ещё даже среднее не брали! Как мы первый раз делили на 3 (когда никаких параметров сближения не находили), во второй раз делили на 2 (когда находили один параметр), ну а теперь уже и не надо...
То есть, чтобы оценить шум измерений, эту самую сигму (которую мы на самом деле не знаем изначально, априори), для нынешней задачи (3 точки, 2 параметра сближения) предлагается такое выражение:
Наконец, будем находить 3 параметра сближения, т.е добавим ещё и поворот в плоскости. И тут оказывается, что никакой "неустранимой невязки" больше нет, оптимизационная задача превратилась в задачу прямую: 3 уравнения с 3 неизвестными.
Пока что моя догадка подтверждается. Я в реальном видеоизмерителе, где каждая точка даёт 2 значения (по горизонтали и вертикали), и всего точек N штук, стал оценивать шум измерений вот так:
Интуитивно вроде всё правильно: если каждая точка на фотоприёмной матрице может смещаться в двух направлениях, то весь шум можно изобразить в виде вектора из 2N значений. Нахождение параметров сближения - это взятие проекций этого вектора на 6 направлений с последующим вычитанием этих проекций. В итоге размерность остающегося вектора уменьшается до 2N-6.
Но всё же пока я не успокоился. Возможно, в нынешнем примере просто "повезло", надо попробовать это "абстрактно" доказать. И также надо попробовать решить ещё одну оптимизационную задачу - с какими весами брать неустранимые невязки по точкам, чтобы оценка шума имела наименьшую возможную дисперсию...
Эх, а ещё надо с терминологией разобраться, каждой участвующей величине назначив своё имя. А то у меня куда ни плюнь - невязки...
Максимально близкое - но всё же без точного совпадения. Мы отдаём себе отчёт, что все измерения неизбежно "зашумлены". Узнать уровень шума АПРИОРИ очень сложно, поскольку он состоит из нескольких компонентов. Это и точность изготовления мишени, и аберрации оптической системы, и, самое главное, "геометрические шумы" матрицы, в первую очередь неравномерность чувствительности пикселей (PhotoResponse Non-Uniformity, или PRNU). Это КМОП-матрица, в каждом пикселе содержится свой полевой транзистор, на затворе которого накапливается сигнал, а на истоке он "повторяется", чтобы пойти в АЦП и не утеряться в процессе. Так что параметры пикселей могут прилично "гулять", причём на разных температурах, и просто в течение времени - по-разному. А потом ещё добавить сюда "износ" матрицы радиацией, всякие фокусы с пульсацией напряжения питания, потом ещё добавить шум квантования и прочие "пакости" от АЦП, и это всё будет меняться в зависимости от светотехнической обстановки, в которой мы будем работать...
Поэтому давным-давно я решил вообще не связываться с априорной оценкой шума, а определять его как раз-таки по неустранимой невязке, т.е по отличию реального изображения от наиболее близкого "модельного", соответствующего измеренным параметрам сближения. И сделал там одно допущение, которое толком не проверял.
Сейчас на " игрушечном видеоизмерителе" хочу понять, насколько это допущение было верно...
Мы снова поработаем в одной плоскости. Мишень на этот раз состоит из 3 точек, A1 (1;1), A2 (1;0) и A3 (1;-1).
Линза видеоизмерителя (точнее говоря, передняя главная точка) сидит в начале координат, (0;0), и проецирует эти точки на плоскость x=-1. Как обычно, игнорируем, что картинка переворачивается (поскольку считываем мы её как положено), и мы должны получить на "фотоприёмной линейке" пятна в координатах G1 = 1, G2 = 0, G3 = -1.
Но из-за разных шумов (то ли мишень плохо изготовили, то ли пятно "провалилось" между пикселями) мы получили результаты F1 = 1 + ξ, F2 = η, F3 = -1 + χ, где ξ,η,χ - независимые случайные величины с нулевым средним и одинаковой дисперсией, равной σ2. Модель немного упрощённая, но пока нас устроит.
И сейчас начнём варьировать число измеренных параметров.
Вырожденный случай: мы и не собирались измерять параметры, мы считаем, что поместили мишень во вполне конкретную точку, и просто "калибруем" нашу оптику, и заодно проверяем уровень шумов. Модельные координаты: G1 = 1, G2 = 0, G3 = -1. Разность между реальными и модельными: ξ,η и χ, соответственно.
Неустранимой невязкой назовём среднеквадратичное значение, т.е
С вычислением математического ожидания этой невязки есть некоторые проблемы, а вот математическое ожидание выражения под корнем посчитается легко:
Поскольку у каждой из наших случайных величин нулевое математическое ожидание, наши слагаемые - не что иное, как дисперсии. Получаем:
Напишем так:
Точного равенства не получается, но и не очень-то хотелось. Нас бы и дисперсия устроила, корень извлекаем просто для удобной размерности. Невязка в 1/10 пикселя - это понятно. А в 1/100 кв. пикселей - уже не очень...
[Пример, почему нельзя поставить знак равенства]
Обойдёмся одной случайной величиной ξ. И пусть она принимает лишь 3 значения: -1, 0 и +1, каждое с вероятностью 1/3. Найдём сначала корень из математического ожидания квадрата такой величины:
А теперь наоборот, математическое ожидание от корня:
Впрочем, когда у нас складывается много случайных величин, результаты становятся уже вполне сравнимыми. Только что для очистки совести в экселе помучал исходный пример с 3 случайными величинами, взял их равномерно распределёнными между -1 и +1. Получил отличие в 3..5%, а для настоящего видеоизмерителя, где 8 точек и для каждой по 2 координаты, разница и того меньше.
Ожидание вполне ожидаемое: находя эту неустранимую невязку, мы как раз оцениваем среднеквадратичное значение шума.
Случай посложнее: измеряем только поперечное смещение, то есть смещение по оси Y. Считаем, что приближаться или отдаляться мишень не может, и поворачиваться тоже.
Вся хитрючая математика здесь практически "вырубается". Частные производные по одному-единственному параметру сближения оказываются единичными, поэтому ненормированная оценка "вектора параметров" 1×1 равна
Нормирующая матрица будет иметь размеры 1×1 и станет равна:
J=3.
Обратная к ней: J-1 = 1/3.
Итого, нормированная оценка сдвига:
Кто бы мог подумать: среднее арифметическое...
Теперь, получив такой "сдвиг", мы полагаем, что координаты точек на фотоприёмной линейке станут равны
Поэтому невязки составят:
Найдём их сумму квадратов:
Когда будем искать математическое ожидание, последние 3 слагаемых дадут ноль, поскольку случайные величины имеют нулевое мат. ожидание и независимы друг от друга. Останется следующее:
И наконец, оцениваем математическое ожидание нашей "неустранимой невязки":
Мы обнаруживаем, что неустранимая невязка уменьшилась!
Это произошло из-за того, что нам "удалось" отыскать такой сдвиг, при котором картинки лучше соответствуют друг другу. Кстати, именно поэтому в среднеквадратичном отклонении ДЛЯ ВЫБОРКИ в знаменателе стоит N-1. По сути, ровно этот механизм мы и увидели: ведь когда мы считаем это среднеквадратичное отклонение, мы не знаем истинного среднего, его мы тоже ОЦЕНИВАЕМ, и поэтому если у нас все слагаемые случайно сдвинулись в одну сторону, мы это уже не почувствуем, это у нас войдёт в оценку среднего значения, а не в среднеквадратичное отклонение.
Так что наша "неустранимая невязка", если её напрямую интерпретировать как шум каждого измерения, будет давать заниженные ответы. Правильнее будет считать её так (если мы по 3 точкам измеряем ОДИН параметр сближения):
Поехали дальше: ДВА параметра сближения, сдвиги по осям X (дальность) и Y (поперечное смещение).
Частные производные по поперечному смещению, как и раньше, остались единичными. А по продольному получаем -1, 0, +1 для точек (1;1), (1;0) и (1;-1) соответственно. То есть, если дальность увеличивается, мы обнаружим, как верхняя точка ползёт вниз, центральная остаётся на одном месте, и нижняя ползёт вверх, все они "стягиваются" в ноль, т.е "в горизонт".
Найдём ненормированную оценку сдвига:
Матрица J получается довольно-таки простой, без перекрёстных связей, за счёт симметричности мишени и её расположения:
Наконец, находим нормированную оценку сдвига:
Чтобы посчитать неустранимую невязку, нужно "применить" оба этих сдвига. Будем считать, что шумы у нас всё-таки маленькие, так что можно здесь работать с "линеаризованной" задачей. Тогда получим такие "ожидаемые" положения пятен на фотоприёмной линейке:
В результате получим невязки:
Мы уже сообразили, что при возведении каждой невязки в квадрат можно сразу оставить лишь квадраты случайных величин, а где одна умножается на другую - вычёркиваем, т.к мат. ожидание равно нулю. Для начала посчитаем мат. ожидания квадратов отдельных невязок:
Разница в дисперсии в 4 раза между крайними и центральной точками! По среднеквадратичному значению - в 2 раза. Исходя из "теории оптимального оценивания", мы догадываемся, что правильнее было бы складывать квадраты невязок с разными весами, гораздо более "пристальное внимание" обратив на центральную точку.
Но допустим, что мы снова просто складываем все квадраты невязок. Тогда получим математическое ожидание от их суммы квадратов:
А это мы ещё даже среднее не брали! Как мы первый раз делили на 3 (когда никаких параметров сближения не находили), во второй раз делили на 2 (когда находили один параметр), ну а теперь уже и не надо...
То есть, чтобы оценить шум измерений, эту самую сигму (которую мы на самом деле не знаем изначально, априори), для нынешней задачи (3 точки, 2 параметра сближения) предлагается такое выражение:
Наконец, будем находить 3 параметра сближения, т.е добавим ещё и поворот в плоскости. И тут оказывается, что никакой "неустранимой невязки" больше нет, оптимизационная задача превратилась в задачу прямую: 3 уравнения с 3 неизвестными.
Пока что моя догадка подтверждается. Я в реальном видеоизмерителе, где каждая точка даёт 2 значения (по горизонтали и вертикали), и всего точек N штук, стал оценивать шум измерений вот так:
Интуитивно вроде всё правильно: если каждая точка на фотоприёмной матрице может смещаться в двух направлениях, то весь шум можно изобразить в виде вектора из 2N значений. Нахождение параметров сближения - это взятие проекций этого вектора на 6 направлений с последующим вычитанием этих проекций. В итоге размерность остающегося вектора уменьшается до 2N-6.
Но всё же пока я не успокоился. Возможно, в нынешнем примере просто "повезло", надо попробовать это "абстрактно" доказать. И также надо попробовать решить ещё одну оптимизационную задачу - с какими весами брать неустранимые невязки по точкам, чтобы оценка шума имела наименьшую возможную дисперсию...
Эх, а ещё надо с терминологией разобраться, каждой участвующей величине назначив своё имя. А то у меня куда ни плюнь - невязки...