Е.С. Семенистая
КОРРЕКЦИЯ ФАЗОВЫХ ИСКАЖЕНИЙ ПРИ ОБРАБОТКЕ ФОТОПЛЕТИЗМОГРАФИЧЕСКОГО СИГНАЛА
Создание неинвазивных диагностических медицинских приборов и систем является перспективным и активно развивающимся направлением. Для получения стабильных, воспроизводимых результатов измерения физических параметров биообъекта необходимо стандартизировать условия измерения и устранить систематические ошибки, связанные с амплитудно-фазочастотными характеристиками (АФЧХ) измерительных преобразователей биомедицинских датчиков. В рассчитываемых параметрах состояния пациента эти искажения могут давать, как показали исследования, погрешность до 50 % (например, при контурном анализе импе-дансной реовазограммы). Амплитудно-фазочастотные искажения можно устранять апостериорной обработкой полученного сигнала. Такая обработка предполагает большой объем вычислений, поэтому ее использование возможно в рамках диагностического комплекса, включающего ЭВМ.
Для оценки параметров состояния пациента по пульсовым кривым необходимо быть уверенным, что полученные результаты достоверны, сопоставимы и могут быть воспроизведены многократно. Таким образом, естественно возникает требование стандартизации условий получения данных для анализа.
При решении проблемы фазочастотных искажений существуют 2 пути:
1) зафиксировать какую-либо электрическую схему, что фактически останавливает техническое совершенствование;
2) использовать алгоритм приведения выходной кривой (или рассчитываемых параметров) к какой-либо стандартной ФЧХ прибора, при этом естественным является использование в качестве стандартной ФЧХ, тождественно равной нулю для всех частот (т.е. прибор не искажает фазу сигнала).
Поскольку в оптике хорошо разработаны алгоритмы математического восстановления искаженных изображений при известных характеристиках передающих систем, то при использовании данного подхода мы можем опираться на уже готовый математический аппарат. Для решения задач такого рода предлагается использовать методы спектрального анализа. Исходными данными при этом являются получаемые на выходе прибора сигналы, АЧХ и ФЧХ прибора и допущения о линейности преобразования сигнала прибором и известности полосы частот, в которой лежат информативные гармоники спектра сигнала. Сравнительно низкочастотные диагностические сигналы, тем не менее, приводят к необходимости обрабатывать большие объемы данных, что заставляет принимать во внимание время выполнения программ и применять эффективные алгоритмы дискретного преобразования Фурье.
Одним из способов устранения искажений является апостериорная обработка сигнала с целью уменьшения влияния реальных характеристик прибора (сведение их к "идеальному" виду). Вопрос о редукции к идеальному прибору, не искажающему полезный входной сигнал, был поставлен еще Рэлеем в 1871 г. в связи с некоторыми задачами спектроскопии и состоял в следующем: можно ли после опыта так обработать функцию, наблюдаемую на выходе некоторого реального прибора, чтобы полностью восстановить сигнал, поступивший на его вход. В такой постановке задача относится к классу некорректно поставленных, и в общем случае не может быть решена, однако на практике при соблюдении некоторых требований (ограниченность спектра рассматриваемого сигнала по частоте, отсутствие в рассматриваемом диапазоне частот нулей в передаточной функции прибо-
ра и т.д.) она может быть решена с достаточной точностью и иметь практическое значение.
Уровень современного развития элементной базы позволяет реализовать достаточно сложные схемотехнические решения и устранить большинство искажений входного сигнала и помех. Применение высокоточных фотодатчиков даёт возможность получать сигнал в широком динамическом диапазоне. Построение схем на основе современных микропроцессоров цифровой обработки сигналов позволяет реализовать алгоритмы коррекции амплитудных и фазовых искажений на этапе съема сигнала и его предварительной обработки, что существенно повышает качество и точность получаемых сигналов и повышает качество диагностического комплекса в целом.
Решить проблему фазовых искажений позволяют вычислительные алгоритмы.
Коррекция исходного сигнала есть не что иное, как фильтрация с передаточной функцией 1/И('^). Синтез подобного фильтра аналоговыми методами в общем случае невозможен, то можно реализовать его двумя способами: либо программно в микроконтроллере, либо апостериорно с использованием алгоритмов цифровой обработки сигналов в ПК.
В теории цифровой обработки сигналов часто используется метод квантованного БПФ. Суть его в том, что значения синуса и косинуса квантуются, т.е. принимают дискретные значения, например -1, 0, 1 при трех уровнях и -1, -1/2, 0, 1/2, 1 при пяти уровнях (вообще говоря, это разложение сигнала по другим базисным функциям). В этом случае умножения либо вообще отсутствуют, либо заменяются сдвигами. Для фильтрации изображений этот метод дает вполне приемлемые результаты.
Существует еще ряд базисов (функции Уолша, Хаара, и т.д.), разложение по которым привлекает внимание тем, что операции можно вести с целыми числами.
Наиболее удобна следующая схема вычислений:
1. Предварительно измеряем или вычисляем передаточную функцию прибора ЩУ). Вычисляем передаточную функцию восстанавливающего фильтра 1/И('м?) (если предполагается прямая свертка во временной области, то затем вычисляем импульсную характеристику (обратное ДПФ от передаточной характеристики), с учетом приведенных выше замечаний).
2. Исходный сигнал умножаем на функцию окна.
3. Проводим фильтрацию-восстановление (во временной или в частотной области).
4. Отфильтрованный сигнал делим на функцию окна.
Вопрос о том, в какой области проводить фильтрацию-восстановление - во
временной или в частотной, решается исходя из требуемой точности и скорости
обработки.
Как правило, реальные приборы с достаточной степенью точности являются линейными, и сигнал ОД на их выходе описывается интегральным уравнением:
р(^о = • Б(^о , (1)
где И(1) - импульсная переходная характеристика прибора, 8(1) - сигнал, поступающий на вход прибора. Это уравнение является уравнением свертки, и может быть записано как
Б(^0 = Б^) / И(^0 , (2)
где Р(^0, И(^0, Б(^0 - Фурье-образы функций ОД, И© и 8(1), т.е. соответственно спектр отклика, передаточная функция системы и спектр входного сигнала. Из уравнения (2) следует, что? по крайней мере формально, можно записать
¥
8(1) = 1/(2р) |F(w)/H(w)exp(jwt)dw , (3)
— ¥
откуда обратным преобразованием Фурье находим
¥
ст2 = 1/(2р) | Яп (w)/1 H(w) |2 dw . (4)
—¥
На практике применение формулы (4) сопряжено с рядом сложностей:
1. В действительности Н^) отлична от нуля лишь в конечном интервале частот. Радиотехнические приборы разрабатываются таким образом, что Н(^ находится вне частотного диапазона сигнала, поэтому естественной становится коррекция в диапазоне частот, ограниченном частотными характеристиками сигнала, и зануление всех частот, где спектральная плотность мощности шума превышает сигнал.
2. Присутствие в сигнале аддитивной помехи. Выходной сигнал прибора ОД известен всегда приближенно, и содержит в себе случайную помеху у(1):
ВД = йОД + уф,
тогда согласно выражению (3) спектр восстановленного сигнала равен
Б^) = Бт^) + У^)/Н^),
где У^) - спектр помехи. Спектральная плотность случайной ошибки восстановления (дисперсия решения задачи) равна
0(1) = ехр(-2(а1/(2Т))2) , (5)
где Яу^)= |У(^)|2 - спектральная плотность мощности помехи. Реальная помеха содержит компоненту белого шума, и ее спектральная мощность стремится к конечному пределу при w . Поэтому интеграл (5) может не сходиться либо быть недопустимо большим. В этом случае можно рекомендовать ограничение частотного диапазона коррекции с подавлением ненужных частот спектра.
3. Погрешности, связанные с квантованием уровня сигнала, округлением при машинном представлении действительных чисел, и т.д. Оценка ошибок восстановления, обусловленных данными погрешностями, рассматривается в литературе и должна в каждом случае оцениваться отдельно. В наших задачах сколько-нибудь заметного появления этих факторов замечено не было.
4. Граничные эффекты. Наиболее неприятным фактом при коррекции искажений сигналов являются краевые эффекты. На практике сигнал известен не на всей временной оси, а в диапазоне (-Т, Т), как показано на рис. 1,а. Поэтому для формулы фильтрации в пространственной области (1) возникает вопрос о доопределении сигнала за границами рассматриваемого интервала. В простейшем случае считают, что за границами (-Т, Т) сигнал равен своему математическому ожиданию как показано на рис. 1,б. Также используется различная экстраполяция (например, считают, что сигнал вне (-Т, Т) равен своему значению на границах интервала).
в) г)
Рис. 1. Варианты доопределения сигналов за границами выборки
При цифровой фильтрации (особенно при использовании аппарата преобразований Фурье) часто исходный сигнал считают периодически продолженным с периодом, равным длине исходной последовательности. Так как начало и окончание записи сигнала расположены в произвольные (независимые от сигнала) моменты времени, то длина исследуемого интервала не равна целому числу периодов сигнала, что приводит к несовпадению амплитудного значения (и значения производных амплитуды) сигнала в начальной и конечной точках и, как следствие, к появлению в спектре высокочастотных гармоник, а после фильтрации - к возникновению осцилляций на краях последовательности (эффект Гиббса). Некоторые исследователи рекомендуют в этом случае использование линейной аппроксимации между начальной и конечной точками, как на рис. 1,в, или четного продолжения исходной последовательности, как на рис. 1,г. Проведенные численные эксперименты показали, что хотя эти способы уменьшают высокочастотные выбросы спектра, однако сильно искажают его низкочастотную область. Более удачные результаты при использовании дискретного преобразования Фурье (ДПФ) получаются при использовании так называемых "весовых окон". Обычно они применяются при оценке спектра по малой выборке: предварительно сигнал умножают на быстро спадающую функцию, например на кривую Г аусса:
ва) = ехр(-2(а1/(2Т))2) , (6)
где 2Т - длительность выборки, |1|<Т, а - параметр (наилучшие результаты получаются при 2,5<а<3), затем проводят Фурье-преобразование для оценки спектра. Другие весовые окна рассматриваются, например, в работе [4]. Их применимость в нашем случае определяют два фактора: скорость спада боковых лепестков и малые вычислительные погрешности. Окно Ханна в(1) = ехр(-2(а^(2Т))2) дает лучшие результаты в смысле отсутствия боковых лепестков, но на границах его значения близки к нулю, поэтому после фильтрации часть информации на концах интервала будет потеряна.
Исходный сигнал можно рассматривать как бесконечный, умноженный на прямоугольное весовое окно, что в частотной области эквивалентно свертке с Фу-рье-образом этого окна, функцией бшс(х) = бш(х)/х. В результате истинный спектр как бы "размывается" по частоте (боковые лепестки имеют уровень -13.3 дБ), что приводит:
• к неограниченности полученного спектра (даже если в качестве исходного сигнала взят бш^));
• к тому, что гармоника вносит ненулевой вклад во весь спектр сигнала (т.е. их невозможно разделить при фильтрации).
Предварительное умножение на функцию Гаусса приводит к тому, что свертка производится с ее Фурье-образом, а не с бшс(х), таким образом, за счет
некоторого уширения самой спектральной полосы сильно подавляются ее боковые лепестки (до уровня - 45.. - 55 дБ, в зависимости от параметра а), и спектр более точно аппроксимирует спектр реального бесконечного сигнала.
Д.А. Беспалов, В.Е. Золотовский, Т.А. Головченко, Е.В. Ляпунцова
АППРОКСИМАЦИОННОЕ СЖАТИЕ СЕЙСМОДАННЫХ ПОЛИНОМАМИ НАИЛУЧШЕГО ПРИБЛИЖЕНИЯ
Введение. Передача и хранение сейсмоданных в современных системах обработки подразумевает использование процедур сжатия информации. В настоящее время наибольший интерес вызывают алгоритмы, основанные на разложении сигналов по некоторому ортогональному базису [1,2]. Однако ввиду огромной информационной емкости обрабатываемых данных, остро стоит вопрос о создании процедур, минимизирующих временные и вычислительные затраты на этапах анализа, кодирования и декодирования.
Такие алгоритмы должны быть «быстрыми» и обеспечивать максимальную скорость работы с минимальными искажениями и потерями информации. Для существующих многочисленных базисов можно выделить ограниченный спектр преобразований, допускающих быструю реализацию. Это различные варианты дискретного преобразования Фурье, разложения по базису разрывных кусочнопостоянных функций (Уолша и Хаара), а также вейвлет-преобразование. В последнее время выполнен ряд работ, расширяющий этот список разложениями по базису классических ортогональных полиномов Эрмита, Якоби, Лежандра и Чебышева [3, 4, 5]. Такие разложения (по результатам экспериментов) позволяют достичь максимального коэффициента сжатия с минимальными потерями и допускают применение «быстрых» алгоритмов.
В данной работе предлагается метод компрессии сейсмоданных с использованием «быстрых» алгоритмов кодирования анализируемого сигнала полиномами наилучшего приближения, спектрального сжатия, а также алгоритмов вторичного статистического сжатия без потерь.
Эксперименты показали, что специфика обрабатываемых данных не позволяет использовать классические методы компрессии «в чистом» виде. Так, например, представление данных в числах с плавающей запятой ограничивает применение алгоритмов словарного или статистического сжатия без потерь, а наличие шумов и высокая волотильность данных (быстрые изменения амплитуды сигнала) значительно снижают применение таких алгоритмов сжатия с потерями, как, например, JPEG.
Исходя из этого, был разработан адаптивный алгоритм компрессии сейсмоданных, включающий в себя следующие этапы:
1) пространственное разделение сигнала на подобласти высокой или низкой волотильности;
2) полиномиальное сжатие данных;
3) вторичное кодирование коэффициентов методами компрессии без потерь.
Рассмотрим последовательно каждый этап разработанного алгоритма.
Пространственное разделение сигнала необходимо для определения областей с плавными или резкими изменениями амплитуды. Участки малой воло-тильности легко аппроксимируются малым количеством коэффициентов, в то вре-