Заметки по поводу фильтрации FITS-ов
Mar. 30th, 2012 09:45 am![[personal profile]](https://www.dreamwidth.org/img/silk/identity/user.png){kind=small}
eddy_emВ связи с тем, что я уже пятый день мучаю тех. задание на систему регистрации ИК-спектрофотометра, у меня возникли кое-какие соображения касательно фильтрации шумов на изображении.
Итак, принцип работы КМОП-мультиплексоров таков, что позволяет производить в процессе накопления недеструктивные считывания сигнала. Пошукав по литературе и пообсуждав этот факт, было решено добавить в систему регистрации фильтрацию изображения посредством линейной аппроксимации сигнала в каждом пикселе изображения. А вот как эту аппроксимацию сделать - уже особый разговор.
Допустим, что мы делаем экспозицию в 1 час с квантом времени в 5с на считывание промежуточных изображений. В итоге получаем 3600/5=720 16-битных картинок размером 1024х1024, на хранение которых понадобится 1024*1024*2*720=1.5ГБ (плюс еще какое-то место на «шапки»)! Теоретически, все их можно «запихать» в оперативку, для этого рассчитаем предельный размер оперативной памяти: если экспозиция будет длиться всю ночь (самая длинная ночь у нас — где-то 15 часов, вычтя время на сумерки и всякие подготовки, получим в пределе часов 13), квант времени будет равен 1с, то максимальное количество файлов составит 13*3600=46800 и займут эти файлы как минимум 93ГБ! Есть материнки, поддерживающие 128ГБ оперативки, т.е. теоретически такой компьютер собрать можно, практически же стоить он будет слишком дорого. Поэтому идем другим путем.
Другой путь — складывать все временные файлы на диске (лучше всего - на SSD, для увеличения скорости доступа). Т.е. каждое промежуточное изображение с соответствующей шапкой мы сохраняем на диск (всю шапку заполнять не надо - только меняющиеся данные, вроде температур и т.п. параметров, а также - времени от предыдущего считывания, т.к. выдержать требуемый квант с точностью в микросекунды вряд ли получится; постоянные параметры достаточно занести в шапку только первого файла). В итоге по окончании экспозиции у нас на диске получается n-е количество FITS-файлов, которые можно в фоновом режиме, не мешая наблюдениям, обработать и отдать под утро наблюдателю все его файлы в «чистеньком» виде.
Для обработки нам надо для каждого пикселя проанализировать кривую зависимости уровня сигнала от времени экспозиции (в основном - ее линейную часть, аппроксимирующуюся выражением I=I₀+ΔIt), выдернув оттуда такие параметры, как bias (I₀), произведение gain на освещенность в данном пикселе (ΔI), предел линейности (тот участок графика, где наша кривая начинает «загибаться» вправо) и уровень насыщения (предельный уровень сигнала в пикселе, не изменяющийся при повышении времени экспозиции). Нелинейные параметры, в принципе, наблюдателю неинтересны, поэтому расчет их можно сделать опциональным. Помимо этих характеристик обработка позволит выявить «горячие» пиксели (находящиеся в насыщении вне зависимости от уровня сигнала), «битые» пиксели (вроде «горячих», но с низким уровнем сигнала) и «неправильные» пиксели (с нелинейной характеристикой или прочими проблемами).
Теория, конечно, хороша, но опять вернемся к нашим баранам: современные компьютеры пока еще не позволяют за разумные деньги сделать такие вычисления без ухищрений. Во-первых, cfitsio не позволяет открыть такое количество файлов (а если бы библиотека и позволяла, ее особенности в том, что она, по сути, mmap'ит файлы в память - а то и просто копирует их туда, так что, оперативки не хватит), следовательно, если захотеть открывать эти файлы «разом» и поочередно считывать, придется писать свой FITS-парсер (чего делать не хочется — это во-вторых). В-третьих, можно придумать еще уйму «отмазок», чтобы закончить эти размышления и послать все нафиг. Но работать-то надо!
Итак, единственным более-менее приличным (и легко реализуемым) вариантом предварительной обработки данных мне предоставляется следующий:
До обработки создаем итоговый FITS и заполняем его шапку: заносим всю статику, обрабатываем динамику и заносим туда же экстремальные, средние и медианные значения, а также стандартное отклонение. Далее начинаем заполнять данные изображения:
Да, забыл добавить: благодаря тому, что этот способ очистки компенсирует «перекоп» от слишком ярких объектов (т.к. анализирует лишь линейный участок, а затем умножает крутизну на время экспозиции), данные надо будет сохранять уже как минимум в формате 32-битных целых, а то и 32-битных числах с плавающей точкой. FITS-файлы с дополнительной информацией позволят упростить предварительную обработку изображений, да и всякие калибровки будут уже значительно точнее.
Итак, принцип работы КМОП-мультиплексоров таков, что позволяет производить в процессе накопления недеструктивные считывания сигнала. Пошукав по литературе и пообсуждав этот факт, было решено добавить в систему регистрации фильтрацию изображения посредством линейной аппроксимации сигнала в каждом пикселе изображения. А вот как эту аппроксимацию сделать - уже особый разговор.
Допустим, что мы делаем экспозицию в 1 час с квантом времени в 5с на считывание промежуточных изображений. В итоге получаем 3600/5=720 16-битных картинок размером 1024х1024, на хранение которых понадобится 1024*1024*2*720=1.5ГБ (плюс еще какое-то место на «шапки»)! Теоретически, все их можно «запихать» в оперативку, для этого рассчитаем предельный размер оперативной памяти: если экспозиция будет длиться всю ночь (самая длинная ночь у нас — где-то 15 часов, вычтя время на сумерки и всякие подготовки, получим в пределе часов 13), квант времени будет равен 1с, то максимальное количество файлов составит 13*3600=46800 и займут эти файлы как минимум 93ГБ! Есть материнки, поддерживающие 128ГБ оперативки, т.е. теоретически такой компьютер собрать можно, практически же стоить он будет слишком дорого. Поэтому идем другим путем.
Другой путь — складывать все временные файлы на диске (лучше всего - на SSD, для увеличения скорости доступа). Т.е. каждое промежуточное изображение с соответствующей шапкой мы сохраняем на диск (всю шапку заполнять не надо - только меняющиеся данные, вроде температур и т.п. параметров, а также - времени от предыдущего считывания, т.к. выдержать требуемый квант с точностью в микросекунды вряд ли получится; постоянные параметры достаточно занести в шапку только первого файла). В итоге по окончании экспозиции у нас на диске получается n-е количество FITS-файлов, которые можно в фоновом режиме, не мешая наблюдениям, обработать и отдать под утро наблюдателю все его файлы в «чистеньком» виде.
Для обработки нам надо для каждого пикселя проанализировать кривую зависимости уровня сигнала от времени экспозиции (в основном - ее линейную часть, аппроксимирующуюся выражением I=I₀+ΔIt), выдернув оттуда такие параметры, как bias (I₀), произведение gain на освещенность в данном пикселе (ΔI), предел линейности (тот участок графика, где наша кривая начинает «загибаться» вправо) и уровень насыщения (предельный уровень сигнала в пикселе, не изменяющийся при повышении времени экспозиции). Нелинейные параметры, в принципе, наблюдателю неинтересны, поэтому расчет их можно сделать опциональным. Помимо этих характеристик обработка позволит выявить «горячие» пиксели (находящиеся в насыщении вне зависимости от уровня сигнала), «битые» пиксели (вроде «горячих», но с низким уровнем сигнала) и «неправильные» пиксели (с нелинейной характеристикой или прочими проблемами).
Теория, конечно, хороша, но опять вернемся к нашим баранам: современные компьютеры пока еще не позволяют за разумные деньги сделать такие вычисления без ухищрений. Во-первых, cfitsio не позволяет открыть такое количество файлов (а если бы библиотека и позволяла, ее особенности в том, что она, по сути, mmap'ит файлы в память - а то и просто копирует их туда, так что, оперативки не хватит), следовательно, если захотеть открывать эти файлы «разом» и поочередно считывать, придется писать свой FITS-парсер (чего делать не хочется — это во-вторых). В-третьих, можно придумать еще уйму «отмазок», чтобы закончить эти размышления и послать все нафиг. Но работать-то надо!
Итак, единственным более-менее приличным (и легко реализуемым) вариантом предварительной обработки данных мне предоставляется следующий:
- создать во временной директории (назовем ее пока /tmp) 1024 (в общем, в зависимости от количества строк) поддиректории, в каждой из которых завести по 1024 (в зависимости от количества столбцов) файла;
- открыть очередной FITS-файл;
- прочитать «шапку» и сохранить в оперативке динамические данные (статику возьмем из нулевого изображения);
- перейти к сегменту данных изображения и начать их попиксельное считывание;
- значения интенсивности в пикселе с координатами (X, Y) дописать в файл /tmp/Y/X;
- считать так все пиксели изображения;
- перейти к следующему изображению;
До обработки создаем итоговый FITS и заполняем его шапку: заносим всю статику, обрабатываем динамику и заносим туда же экстремальные, средние и медианные значения, а также стандартное отклонение. Далее начинаем заполнять данные изображения:
- открываем очередной файл /tmp/Y/X;
- считываем все данные в оперативку (или просто mmap'им его);
- получаем графики первой и второй производных;
- из анализа графиков производных находим точки перегиба (определяя таким образом участок «разгона», линейный участок, переход в насыщение и собственно насыщение);
- если мы не имеем в пикселе линейного участка — определяем, «горячий» он или «битый», увеличивая соответствующий счетчик и занося в соответствующий пиксель 65535 или 0 (плюс можно создавать дополнительные изображения-маски, в которых помечать все «косяки») и переходим к следующему пикселю;
- аппроксимируем линейный участок, вычисляем требуемые параметры (I₀, ΔI и т.п.) и заносим в очередной пиксель итогового изображения значение, равное ΔI·T (T — время экспозиции), производя таким образом сразу же очистку от bias'а; сопутствующую информацию можно записывать в дополнительные файлы;
- переходим к следующему пикселю.
Да, забыл добавить: благодаря тому, что этот способ очистки компенсирует «перекоп» от слишком ярких объектов (т.к. анализирует лишь линейный участок, а затем умножает крутизну на время экспозиции), данные надо будет сохранять уже как минимум в формате 32-битных целых, а то и 32-битных числах с плавающей точкой. FITS-файлы с дополнительной информацией позволят упростить предварительную обработку изображений, да и всякие калибровки будут уже значительно точнее.
