Текст
                    оо оо
Уральский
федеральный
университет
имени первого Президента
России Б.Н. Ельцина
Институт радиоэлектроники
и информационных
технологий — РТФ
В. Г. КОБЕРНИЧЕНКО
ОСНОВЫ
ЦИФРОВОЙ ОБРАБОТКИ
СИГНАЛОВ
Учебное пособие


МИНИСТЕРСТВО НАУКИ И ВЫСШЕГО ОБРАЗОВАНИЯ РОССИЙСКОЙ ФЕДЕРАЦИИ УРАЛЬСКИЙ ФЕДЕРАЛЬНЫЙ УНИВЕРСИТЕТ ИМЕНИ ПЕРВОГО ПРЕЗИДЕНТА РОССИИ Б. Н. ЕЛЬЦИНА В. Г. Коберниченко ОСНОВЫ ЦИФРОВОЙ ОБРАБОТКИ СИГНАЛОВ Учебное пособие Рекомендовано методическим советом Уральского федерального университета в качестве учебного пособия для студентов вуза, обучающихся по направлению подготовки 11.03.01 «Радиотехника», по специальностям 11.05.01 «Радиоэлектронные системы и комплексы», 11.05.02 «Информационная безопасность телекоммуникационных систем» Екатеринбург Издательство Уральского университета 2018
УДК 621.372.083.92(075.8) ББК 32-013я73 К55 Рецензенты: кафедра теоретических основ радиотехники и связи Поволжского государственного университета телекоммуникаций и информатики (заведующий кафедрой доктор технических наук, доцент О. В. Го р я ч к и н); В. Д. Захарченко, доктор технических наук, профессор, профессор кафедры радиофизики Волгоградского государственного университета, федеральный эксперт научно-технической сферы Коберниченко, В. Г. К55 Основы цифровой обработки сигналов : учеб, пособие / В. Г. Коберниченко ; М-во науки и высш. образования Рос. Феде¬ рации, Урал, федер. ун-т. — Екатеринбург : Изд-во Урал, ун-та, 2018.— 150 с. ISBN 978-5-7996-2464-4 Приведено краткое введение в необходимый математический аппа¬ рат и основные понятия цифровой обработки сигналов. Рассмотрены модели и преобразования дискретных и цифровых сигналов. Особое вни¬ мание уделено алгоритмам быстрого преобразования Фурье и их приме¬ нению при цифровом спектральном анализе. Рассмотрены методы опи¬ сания, классификации и методика проектирования линейных цифровых фильтров, анализируются эффекты квантования и округления в цифровых фильтрах. Для студентов, обучающихся по направлению «Электроника, радио¬ техника и системы связи». Может быть полезным для студентов, изучаю¬ щих информационные системы и технологии, системы связи, прикладную информатику. УДК 621.372.083.92(075.8) ББК 32-013я73 ISBN 978-5-7996-2464-4 © Уральский федеральный университет, 2018
ОГЛАВЛЕНИЕ Основные сокращения 5 Предисловие 6 Введение. Аналоговые, дискретные и цифровые сигналы и системы 9 1. Модели и преобразования дискретных и цифровых сигналов 1.1. Математическое описание дискретных сигналов. Теорема Уиттекера — Котельникова — Шеннона 17 1.1.1. Математическая модель дискретного сигнала в непрерывном времени 17 1.1.2. Спектральная плотность модулированной импульсной последовательности 19 1.1.3. Теорема отсчетов 21 1.1.4. Дискретное по времени преобразование Фурье 26 1.2. Дискретное преобразование Фурье 27 1.3. Алгоритмы быстрого преобразования Фурье 35 1.3.1. Идея быстрого преобразования Фурье 35 1.3.2. БПФ с основанием 2 37 1.4. Алгоритмы БПФ с произвольным основанием 45 1.5. Основы теории Z-преобразования 50 Контрольные вопросы 57 2. Дискретные и цифровые фильтры 2.1. Линейные дискретные фильтры и их характеристики 58 2.2. Формы реализации линейных цифровых фильтров 65 2.3. Реализация линейных цифровых фильтров в частотной области с помощью алгоритмов БПФ 69 2.4. Цифровой спектральный анализ 72 2.5. Проектирование цифровых фильтров с конечной импульсной характеристикой 76 2.5.1. Этапы проектирования цифрового фильтра 76 2.5.2. Синтез нерекурсивных фильтров методом «окна» 77 3
2.6. Синтез рекурсивных фильтров по аналоговому прототипу 86 2.7. Метод билинейного Z-преобразования 94 Контрольные вопросы 101 3. Эффекты конечной разрядности при представлении чисел в цифровых фильтрах 3.1. Шум квантования 103 3.2. Эффекты округления промежуточных результатов 114 3.3. Анализ влияния квантования коэффициентов 119 Контрольные вопросы 125 4. Специальные алгоритмы цифровой обработки сигналов в радиотехнических и телекоммуникационных системах 4.1. Изменение частоты дискретизации в линейных цифровых фильтрах 126 4.2. Демодуляция узкополосных сигналов. Цифровые преобразователи Гильберта 134 4.2.1. Узкополосные сигналы, комплексная огибающая, аналитический сигнал 134 4.2.2. Дискретное преобразование Гильберта 140 4.2.3. Реализация дискретного преобразователя Гильберта на основе КИХ-фильтра 141 Контрольные вопросы 145 Библиографические ссылки 147
ОСНОВНЫЕ СОКРАЩЕНИЯ АЦП - АЧХ - БИХ - БПФ - ДВПФ - ДПФ - Дф - ИХ - ких - кчх - мил - ПФ - РФ ФВЧ - ФНЧ - ФЧХ - ЦАП - ЦОС - цпг - ЦФ - - аналого-цифровой преобразователь - амплитудно-частотная характеристика - бесконечная импульсная характеристика (тип фильтра) - быстрое преобразование Фурье - дискретное по времени преобразование Фурье - дискретное преобразование Фурье - дискретный фильтр - импульсная характеристика конечная импульсная характеристика (тип фильтра) комплексная частотная характеристика - модулированная импульсная последовательность - полосовой фильтр режекторный фильтр - фильтр верхних частот - фильтр нижних частот фазочастотная характеристика цифроаналоговый преобразователь цифровая обработка сигналов цифровой преобразователь Гильберта цифровой фильтр
ПРЕДИСЛОВИЕ Цифровая обработка сигналов (ЦОС) уже давно перестала быть только разделом радиотехники и теории связи. Методы циф¬ ровой обработки сигналов используются в самых разных областях: от медицинской диагностики (компьютерная томография) до кос¬ мического мониторинга (обработка данных дистанционного зон¬ дирования Земли), от фото- и видеосъемки (обработка изображе¬ ний) до обеспечения информационной и физической безопасности. Чтобы грамотно использовать методы цифровой обработки сигналов в этих областях и правильно интерпретировать их резуль¬ таты, специалист должен представлять «а что там внутри?». А для этого, прежде всего, нужно овладеть «азбукой» этой предметной области, которая включает модели цифровых сигналов, аппарат анализа и базовые алгоритмы. Изучению этой «азбуки» и посвя¬ щено настоящее учебное пособие, написанное на основе лекций, читаемых по одноименной дисциплине в Уральском федеральном университете имени первого Президента России Б. Н. Ельцина в течение многих лет. За эти годы на русском языке вышло несколько фундаменталь¬ ных монографий, посвященных цифровой обработке сигналов [1, 5, 10, 12, 15], и два учебника, ставших классическими [17, 21]. Вместе с тем потребность в компактном изложении основ ЦОС для студен¬ тов различных направлений подготовки по-прежнему существует. Это подтверждает большое количество учебных пособий, изданных за последние годы во многих университетах [4, 9, 23, 25], непре¬ рывно растущее количество подробных конспектов лекций, разме¬ щенных в сети Internet, и их популярность среди студентов. Курс «Основы цифровой обработки сигналов», наряду с такими курсами, как «Теория электрических цепей», «Радиотехниче¬ ские цепи и сигналы», «Теория электрической связи», «Теория 6
информации и кодирования», обеспечивает базовую подготовку специалистов в области электроники, радиотехники, телекоммуни¬ каций, вычислительной техники и информационной безопасности. С другой стороны, этот курс является вводным в обширную область предметной деятельности, объединенной под названием «Цифровая обработка сигналов», которая включает алгоритмы, спецпроцессоры, средства моделирования и проектирования, а также особенности применения методов ЦОС для решения самых различных прикладных задач [10, 16, 22]. В этой связи важно с самого начала определить, что будет рассмот¬ рено в данном учебном пособии, а что останется за его пределами. Предметом изучения в курсе являются фундаментальные поло¬ жения ЦОС, такие как: — методы описания дискретных и цифровых сигналов и систем (математические модели); — базовые алгоритмы преобразований этих сигналов в цифро¬ вых устройствах (цифровых фильтрах); — способы реализации и методы расчета цифровых фильтров; — методы учета влияния эффектов изменения частоты дискрети¬ зации, квантования и округления результатов в процессе обработки. Получение базовой теоретической подготовки необходимо для изучения принципов функционирования и методов проектирова¬ ния цифровых устройств, используемых в инфотелекоммуника- ционных системах и системах управления. Глубокое понимание основ цифровой обработки сигналов необходимо также для спе¬ циалистов, чья работа связана с обработкой сигналов в различных прикладных областях. Знание фундаментальных основ ЦОС необ¬ ходимо также для грамотного математического моделирования процессов формирования и обработки сигналов в широко распро¬ страненных профессиональных программных комплексах, таких как MATLAB с пакетами расширений, Lab VIEW и др. [7, 20, 24]. За пределами курса остаются разделы, касающиеся таких спе¬ циальных вопросов, как: — алгоритмы вейвлет-анализа и нелинейного спектрального анализа; 7
— адаптивная фильтрация и алгоритмы линейного предсказа¬ ния; — цифровая обработка изображений; — архитектура сигнальных процессоров; — аппаратная реализация алгоритмов ЦОС на базе програм¬ мируемых логических интегральных схем (ПЛИС) и сигнальных процессоров; — инструментальные средства проектирования и моделирова¬ ния систем ЦОС Курс «Основы цифровой обработки сигналов» формирует базу для изучения этих разделов в специальных дисциплинах. Полученные в курсе знания, умения и навыки позволят при¬ ступить к освоению таких дисциплин, как «Цифровые устройства и микропроцессоры», «Аппаратные средства телекоммуникаци¬ онных систем», «Измерения в телекоммуникационных системах», «Информационно-измерительные системы», «Инфокоммуника- ционные системы и сети», «Прием и обработка радиосигналов», «Цифровая обработка изображений», «Основы беспроводной радиосвязи», «Техническая защита информации». Создание пособия стало возможным благодаря поддержке кол¬ лег по кафедре теоретических основ радиотехники и департамента радиоэлектроники и связи ИРИТ-РТФ, с которыми проводилось обсуждение как содержания дисциплины, так и методики ее изуче¬ ния. Автор глубоко признателен рецензентам учебного пособия: профессору, доктору технических наук Владимиру Дмитриевичу Захарченко и заведующему кафедрой теоретических основ радио¬ техники и связи ПГУТИ доктору технических наук Олегу Вале¬ рьевичу Горячкину за внимательное чтение рукописи и полезные рекомендации. С выходом в свет учебного пособия работа над ним не закан¬ чивается. Автор будет благодарен за все замечания и предложения по улучшению его содержания и формы, которые можно направ¬ лять по адресу: v.g. kobernichenko@urfu.ru.
ВВЕДЕНИЕ. АНАЛОГОВЫЕ, ДИСКРЕТНЫЕ И ЦИФРОВЫЕ СИГНАЛЫ И СИСТЕМЫ Сигналом называют физический носитель сообщения, т. е. информации, предназначенной для передачи. В качестве сигналов могут выступать параметры самых различных физических про¬ цессов (давления, температуры, освещенности и т. п.). В радиотех¬ нике сигнал, поступающий от первичного источника информации в виде изменения во времени или пространстве указанных пара¬ метров, преобразуется в электрическое колебание, описываемое законом изменения напряжения или тока. Для того, чтобы сделать сигнал объектом теоретического изу¬ чения, вводят математическую модель сигнала — способ его математического описания, представляющий собой функциональную зависимость, аргументом которой, как правило, является время — s(t), x(t), u{f). Математическая модель позволяет абстрагироваться от физи¬ ческой природы носителя сообщения и описывает наиболее суще¬ ственные свойства сигнала. Введение математической модели позволяет провести классификацию сигналов. Аналоговым или континуальным называют сигнал, про¬ извольный по величине и непрерывный во времени. Аналоговый сигнал х(і) описывается непрерывной или кусочно-непре¬ рывной функцией времени. Аргумент и сама функция принимают любые значения на интервале: х <х<х ; min — — max? Например: t < t < t min — — max' x(t) = ехр(-Щ), 0 < t < oo. 9
Дискретный или импульсный сигнал может принимать произвольные по величине значения в дискретные моменты вре¬ мени. Дискретный сигнал xjj) описывается решетчатой функцией — последовательностью выборочных значений (отсчетов) в соответствующие моменты времени: x0 = x(t0), xx=x(tx), ...,xn = x(tn). При постоянном интервале дискретизации M = ti-ti_l = ti_l-ti_2= ... = Тл. Величину Тд называют периодом дискретизации, а величину, обратную ей — частотой дискретизации: F = — = —. Д Тд At В этом случае значения решетчатой функции записываются как х(пТд), х(п) или просто хп. Таким образом, дискретный сигнал задается как *д(0 = {(х(пТ)}. Переход от аналогового сигнала к дискретному — операция дискретизации — состоит в том, что заданному аналого¬ вому сигналу ставится в соответствие дискретный сигнал: 40 -> xa(tn), причем хл(пТ,) = х(пТ,)- Для приведенного выше примера хд{п) = хд(лТд) = ехр(-апТд) = а", а = еа7Ч п = 0, 1, 2,... . Обратный переход — операция восстановления — состоит в том, что заданному дискретному сигналу ставится в соответствие аналоговый сигнал: 10 ХЮ -*■ 40,
причем х(пТа) = хд(пТд). Эти операции являются взаимно обратными при выполне¬ нии условий теоремы отсчетов (теоремы Уиттекера — Котельни¬ кова — Шеннона). Цифровой сигнал — это квантованный по уровню дискрет¬ ный сигнал. Он описывается квантованными решетча¬ тыми функциями (квантованными последовательностями отсчетных значений), принимающими конечный ряд дискретных значенийd0, dv ... dk, называемых уровнями квантования. Связь между решетчатой функцией хд(пТд) и квантованной решетчатой функцией х (пТ^ определяется нелинейной функ¬ цией — амплитудной характеристикой квантования 0{х): Q(x) = d0+iK1(x-av)' (В.1) ѵ=1 где \ = dv - dv_j — шаг квантования; аѵ — порог квантования; 1(х) — единичная ступенчатая функция; К + 1 — число уровней квантования. Общий вид амплитудной характеристики квантования при квантовании с постоянным шагом Дѵ = А приведен на рис. В.1. Рис. В.1. Общий вид амплитудной характеристики квантования 11
Каждый уровень квантуется кодом, чаще всего двоичным. В этом случае число разрядов кода, описывающего цифровой сиг¬ нал, определяется как т = int [log2(.K + 1)]. (В.2) Здесь функция int(«) означает определение наименьшего целого числа, не менее заданного. Таким образом, переход от дискретного сигнала к цифровому хд(пТд) —> хц(пТ) осуществляется путем применения операций квантования и кодирования (рис. В.2). Соответственно переход от аналогового сигнала к цифровому x(t) —> хц(пТд) проводится путем осуществления операций дискре¬ тизации, квантования и кодирования, составляющих аналого- цифровое преобразование (АЦП) сигнала (рис. В.З). 0 Та 2 тл 3 Гд 1ч N Ь \ \ 4N- в > ^ ~ ^ t 0 12: т ід Рис. В.2. Аналоговый (а), дискретный (б) и цифровой (в) сигналы 12
Та Рис. В.З. Структура аналого-цифрового преобразования При этом возможны два типа искажений сигнала — за счет дискретизации и за счет конечного числа уровней квантования. Выбирая достаточно большое число разрядов, можно увеличить точность представления сигнала, но это приводит к усложнению и удорожанию обработки. Обратный переход — операция цифроаналогового преобразования (ЦАП) состоит в построении сигнала х(1.) по заданному цифровому сигналу: ха(пТл) -> x(t). Эти операции не являются взаимно обратными из-за необрати¬ мых погрешностей при квантовании. Примеры дискретных сигналов: 1. Дискретная дельта-функция (единичный импульс): ~ / \ ІО, пФ к' »к-*г.)={ІЯ=і; <вз> Вид этого дискретного сигнала приведен на рис. В.4. Цп) и 1 ■■ -* * # # * »- п -1012 3 Рис. В.4. Дискретная дельта-функция 13
Введение дискретной дельта-функции позволяет описать любую дискретную последовательность в виде следующего соот¬ ношения, называемого «динамическим представлением сигнала»: I ( ПТ, )= 2 X ( И; )5 ( „Тд - Яд ^ k=-да 2. Единичная последовательность (рис. в.5): 0, n < k; 1, n > k. (в.4) (в.5) рис. в.5. единичная последовательность связь между единичным импульсом и единичной последова¬ тельностью устанавливается следующими соотношениями: Ц(„Тд)=и(пТд)-и(пТд-Т,у, (в.6) и („Т, )=2 Ч„Т. - kT, ). (в.7) k = 0 3. периодическая последовательность (рис. В.6): х(пТд) = х(пТд + mNTд), m = 1, 2, 3, ..., (В.8) при этом величина NT называется периодом последовательности. -2-1012 N 2 N Рис. В.6. Дискретная периодическая последовательность 14
Сдвиг периодической последовательности на к] N отсчетов нельзя отличить от сдвига на к < N отсчетов, где к = (T^mocLV — остаток от деления к1 на N. Такой сдвиг называется цикличес¬ ким. 4. Гармоническая последовательность: ч ( f\ (2 пп Л F х(пТ\ = Acos (2nJhTa} = Acos 2пп— =^4cos , р=—. (В.9) I FJ \ Р ) / Если этот параметр представим в виде р = а/p, где аир — прос¬ тые положительные числа, то гармоническая последовательность является периодической с периодом а отсчетов. В противном слу¬ чае последовательность, образованная путем дискретизации гар¬ монического сигнала, не является периодической. Цифровые фильтры Под цифровым фильтром (ЦФ) в общем случае пони¬ мают систему, преобразующую один цифровой сигнал в другой. Цифровой фильтр может быть реализован как программа на цифро¬ вом сигнальном процессоре или аппаратным способом в виде циф¬ ровой схемы, содержащей регистры, сумматоры, умножители и т. п. Эта схема может быть реализована и на универсальном кристалле. Как указано ниже, в настоящем пособии в соответствии со сло¬ жившейся методикой изучения курса «Основы цифровой обра¬ ботки сигналов» сначала рассматриваются цифровые фильтры без учета эффектов квантования. То есть вместо понятий «дискретный фильтр» и «цифровой фильтр» используется ниже только термин «цифровой». Затем отдельно анализируются эффекты, связанные с цифровыми представлениями сигнала. Структура учебного пособия Дискретные сигналы, так же как и аналоговые, образуют линейное пространство сигналов. В связи с этим аппарат теории дискретных сигналов и дискретных систем разработан так же 15
подробно, как и аппарат теории аналоговых сигналов и систем, и во многом ему аналогичен. Цифровые сигналы, представляемые кодами с ограниченным числом разрядов, не образуют линейного пространства вследствие возможного переполнения при выполнении операций сложения и умножения. Поэтому при использовании теории дискретных линейных систем для описания обработки цифровых сигналов необходимо вводить модели, учитывающие влияние представле¬ ния чисел ограниченным числом разрядов. В этой связи изложение материала в учебном пособии ведется в следующей последовательности. В первой главе рассматрива¬ ются свойства дискретных сигналов, связь с соответствующими характеристиками (прежде всего спектральными) аналоговых сиг¬ налов, а также аппарат и алгоритмы для их описания и преобразо¬ вания. Особое внимание уделяется алгоритмам быстрого преобра¬ зования Фурье (БПФ) Во второй главе рассматриваются характеристики и методы синтеза линейных дискретных фильтров. Третья глава посвящена анализу влияния эффектов, связанных с цифровым представлением сигнала (квантование, ограничение разрядности коэффициентов и промежуточных данных), на резуль¬ тат обработки в линейных цифровых фильтрах. В заключительной, четвертой главе рассматриваются некото¬ рые специальные алгоритмы цифровой обработки сигналов, нахо¬ дящие применение в радиотехнических и инфокоммуникацион- ных системах. Здесь анализируются особенности преобразования сигналов в дискретных системах с понижением или повышением частоты дискретизации в процессе обработки, а также вопросы цифровой демодуляции узкополосных сигналов на основе приме¬ нения цифровых преобразователей Гильберта.
1. МОДЕЛИ И ПРЕОБРАЗОВАНИЯ ДИСКРЕТНЫХ И ЦИФРОВЫХ СИГНАЛОВ 1 Л. Математическое описание дискретных сигналов. Теорема Уиттекера — Котельникова — Шеннона 1ЛЛ. Математическая модель дискретного сигнала в непрерывном времени Важнейшей характеристикой аналогового сигнала является его спектральная плотность. Для анализа трансформации спек¬ тра, происходящей при дискретизации произвольного аналогового сигнала, необходимо задать модель дискретного сигнала в непре¬ рывном времени. Ибо значения дискретного сигнала определены не для всех моментов времени, а лишь в счетном множестве точек. Вследствие чего для математической модели дискретного сиг¬ нала в виде решетчатой функции х(пТа) нельзя записать обычное (непрерывное) преобразование Фурье. В качестве такой математической модели дискретного сигнала используют модулированную импульсную после¬ довательность (МИЛ) [2], получаемую на выходе идеаль¬ ного импульсного дискретизатора (рис. 1.1, о) как результат пере¬ множения аналогового сигнала и идеальной дисретизирующей последовательности г|(/): x„(t) = x(t)r](t); (1.1) Л(0= £ 5((-”Тл)- (12) к=-со откуда *и(0= £ хі”тМ(-”тЛ (!-з) 17
Обратим внимание, что согласно данной модели значения дискретного сигнала в паузах считаются равными нулю, а сама модель представляет собой последовательность 8-импульсов, площади которых равны значениям дискретного сигнала (рис. 1.1, г) [18]. х(і) Идеальный (0 дискретизатор п(0 а Генератор дискретизирующей по следовательно сти Рис. 1.1. Получение модулированной импульсной последовательности: а — блок-схема идеального дискретизатора; б — аналоговый сигнал; в — дискретизирующая последовательность; г — МИП; д — дискретный сигнал
Импульсный сигнал можно превратить в дискретный, заменяя каждый импульс *К)5(*-иГд) на число х(пТа), равное его пло¬ щади. 1.1.2. Спектральная плотность модулированной импульсной последовательности Сначала определим спектральную плотность дискретизирую¬ щей последовательности г|(?). Как любая периодическая функция, она может быть представлена комплексным рядом Фурье: ' .271, Л /—kt Т„ л(0= S С*ехр к=—оо Комплексные коэффициенты этого ряда вычисляются как 1 (14) С* = Zk 1 - дг I §(ОехР 2 ,2к . -j—kt dt = - (1.5) Связь между комплексными коэффициентами ряда Фурье и спектральной плотностью периодического сигнала выражается следующим образом [2]: со .Ѵп (/со) = 2тг X С,5 ( ~ Л CD-k^ V T,J (1.6) Откуда спектральная плотность идеальной дискретизирующей последовательности (рис. 1.2) приобретает вид: (1.7) Спектральная плотность МИП, задаваемой произведением двух функций времени (1.1), представляет собой свертку двух спектральных плотностей — аналогового сигнала Х(/со) и дискре¬ тизирующей последовательности S (/со): 2% 00 0 ( ,2п] — Е 5 со -к— Гд k=-m Т \ л J I со *иМ = — КО;)Х(о-^. Z,7Z —со (1.8) 19
Рис. 1.2. Спектральная плотность идеальной дискретизирующей по следовательно сти Подставим (1.7) в (1.8): 1 00 97Г *,(» = — Z-J5 Z 7Г к=-со / _оо о 9 т>; X(n-$)d$. (1.9) Откуда с использованием «фильтрующего свойства» 8-функции окончательно получаем: (1.10) 1 ■*> ( 2тг Л -XI (0 - —к Тд к=-т V Т д / Таким образом, спектр модели дискретного сигнала — МИП представляет собой сумму бесконечного числа спектров исходного аналогового сигнала, сдвинутых по частоте на величину, кратную частоте дискретизации (рис. 1.3, б). Отметим, что идеальная периодическая дискретизация во вре¬ мени приводит к периодизации спектра по частоте. Спектр дискретного сигнала является периодическим с перио¬ дом, равным частоте дискретизации. Именно поэтому спектраль¬ ный анализ дискретных сигналов имеет смысл проводить только в пределах частотного интервала, равного частоте дискретиза¬ ции. 20
ад а ті/Гд 2тг/Гд -сйд/2 0 сйд/2 Рис. 1.3. Преобразование спектра аналогового сигнала (а-г) при идеальной дискретизации 1.1.3. Теорема отсчетов Теорема отсчетов определяет условия, при выполнении кото¬ рых возможно сколь угодно точное восстановление непрерыв¬ ного сигнала по его дискретным значениям (отсчетам). Как один 21
из частных результатов теории интерполяции функций эта тео¬ рема была известна еще в начале XX в. из трудов английских математиков Э. Т. Уиттекера и Дж. М. Уиттекера. Однако, несмо¬ тря на ее очевидную важность, она не приводилась в литературе по теории связи. Применительно к этим задачам она впервые была сформулирована и доказана В. А. Котельниковым в его докладе «О пропускной способности эфира и проволоки», опубликован¬ ном в книге «Материалы к 1-му Всесоюзному съезду по вопро¬ сам технической реконструкции связи», изданной малым тира¬ жом Всесоюзным энергетическим комитетом в 1933 г. В 1949 г. эту теорему независимо установил К. Шеннон в известной работе «Связь при наличии шума». Именно благодаря Котельникову и Шеннону результаты и выводы этой теоремы стали основой для использования цифровых фильтров при обработке последователь¬ ностей выборок, полученных из аналоговых сигналов, и, в конеч¬ ном итоге, создания цифровых систем связи. Поэтому в мировой технической литературе к настоящему времени за этой теоремой закрепилось название: теорема Уиттекера — Котельникова — Шеннона (УКШ). Она может быть сформулирована следующим образом: непре- рывный сигнал х(І). имеющий спектр {преобразование Фурье), огра- ниченный частотой сов, может бытъ полностью описан последо¬ вательностью его дискретных значений, отстоящих друг от друга на интервал TR = 2л/сод, причем сод > 2 сов, и представлен рядом'. x(t)= Z x(nAt)—Ц (1.11) »!=-« tt>B(t-nAt) Ограничениесод>2совназывается критерием Найквиста, а частота дискретизации сод = 2совтакже называется частотой Найквиста Критерий Найквиста — это условие, при котором теоретически возможно полное восстановление аналогового сиг¬ нала из последовательности равномерно распределенных выбо¬ рочных значений (отсчетов). Одно из доказательств этой теоремы выглядит следующим образом. 22
Пусть спектральная плотность (результат прямого преобразо¬ вания Фурье) непрерывной функции x(t) удовлетворяет условию: Д/co) = 0 при I со I > сов. Использовав обратное преобразование Фурье, запишем: х(7) = — J X(ja)exp(jat)da>. (112) 2% -Юв Для моментов взятия отсчетов tn nAt жение принимает вид: пТ = пк/а\ это выра- I СОВ улгтг = — J ^(усо)ехр(усо—)сіох ZlZ -сов со (1.13) Рассмотрим периодическое продолжение XJjcn) функции Д/'со) с периодом 2сов. Подобная функция приведена на рис. 1.3, б. Эту периодическую функцию можно представить комплексным рядом Фурье по частоте аналогично формуле (1.4) и с учетом того, что роль периода теперь играет 2со,. : f \ 71 ] «СО ДО)= І Аиехр га„ (114) Комплексные коэффициенты этого ряда вычисляются по фор¬ муле, аналогичной (1.5): ( \ А. =■ 1 2со„ J А(усо)ехр к -j—«га га„ d со. (1.15) Здесь зачтено, что на интервале от -сов до сов периодическая XJj(o) и непериодическая Д/со) функции совпадают. Из сопоставления выражений для Ая и х(«л/сов) следует, что д ТС , . ч тс, тс, Ап = —x(-nAt) = —х(-п—). гав гав сов Тогда спектральная плотность периодического продолжения будет представлена как ( \ тг со тг Д (./«) =— 2 х(п—)ехр С0„ и=—со со„ 71 -J—«га га„ 23
Подставим это значение в выражение для x(t) через интеграл Фурье (1.12), учтя, что на интервале от - оу, до сов периодическая A"n(/co) и непериодическая X(jm) функции совпадают. В результате получаем: |[ сов со Л x(t) = — J —{ 2 х(п—)ехр(-jn—)}exp(j(di)d<a. 2п -ив сов я=-*> сов сов Изменим порядок суммирования и интегрирования: j со "Л ®в x(t) = 2 х(п—) J expf/ccK/1 -nAt)]da. 2сОв п=-т С0В -шв Результат вычисления интеграла имеет вид: j ехр[усо(/-нА?)]й?со -exp[yco(f-«AO|l!;,B = 2smaB(t-nAt) t-nAt j(t-nAt) С учетом этого окончательно получаем выражение непрерыв¬ ной функции через ее дискретные значения, взятые в моменты вре¬ мени tr = нА t = пТд = wi/coB: n=-oo <j)B(t-nAt) Результат доказательства не изменится, если частота дискре¬ тизации ®д > 2сов. Представление функции в виде ряда Котельникова представ¬ ляет собой частный случай разложения в обобщенный ряд Фурье по системе базисных функций (рА.(0: x(t)= £ с4ф*(0- (116) к=-со Коэффициенты разложения являются отсчетами x(kAt), а базис¬ ные функции _sincoB(r-MQ сoB(t-kAt) ортогональны между собой на бесконечном интервале времени. Отметим, что теорема отсчетов определяет как условия дискретизации непрерывного сигнала, при которых воз¬ можно его неискаженное восстановление (ограниченный спектр, 24
частота дискретизации не ниже ширины спектра), так и сам спо¬ соб восстановления. Соотношение (1.11) можно интерпре¬ тировать как результат прохождения импульсного сигнала — МИП (1.3) через фильтр с импульсной характеристикой: sin(co//2) b(t) = - со//2 (1.18) Такой импульсной характеристикой обладает идеальный фильтр нижних частот с комплексной частотной характеристикой (КЧХ) вида (рис. 1.3, в): _\ТЛ |ФЮд/2 = jo |со| > сод/2. (1-19) Идеальный ФНЧ физически нереализуем (так как его импульс¬ ная характеристика является опережающей), однако его характе¬ ристику можно аппроксимировать тем или иным образом. Как следует из доказательства теоремы отсчетов, выборки сигнала полностью определяют лишь спектр, полученный путем периодического продолжения исходного спектра аналогового сигнала (рис. 1.3, а) с периодом, равным сод, т. е. спектр МИП (рис. 1.3, б). Если (од < 2<ов, то спектр МИП не совпадает в основ¬ ной полосе (-(Од/2, <од/2) со спектром непрерывного сигнала за счет эффекта, называемого наложением. Характер этого явления показан на рис. 1.3, г. Степень допустимого наложения спектров при выборе частоты дискретизации определяется отношением мощности помехового сигнала Р№ возникающего в основной полосе за счет всех «сдвину¬ тых» составляющих спектра МИП, к мощности полезного сигнала в той же полосе Рс. Это отношение задается в децибелах и вычис¬ ляется через спектральную плотность амплитуд Х(оо) исходного аналогового сигнала следующим образом: ІГіх2(®¥® f р ^ 1 X р V е У дБ §|0и1Х2 (со)с/со (1.20) где со 1 - 7і/д = 7і/7'д — половина частоты дискретизации. 25
1.1.4. Дискретное по времени преобразование Фурье Установим теперь, как и при каких условиях можно опреде¬ лить спектральную плотность непрерывного сигнала -У(/со) через его отсчеты в дискретные моменты времени, т. е. через значения дискретного сигнала хд(и) = хд(пТд). Для этого запишем общее выражение для спектральной плот¬ ности МИП — модели дискретного сигнала в непрерывном вре¬ мени: со Хп0'(а)= I Хя(Г)е~^Л- —СО Подставляя выражение для МИП (1.3), получаем: ОО 00 СО 0О *„(./«)= 1 I х(яГд)(!О-»•/;)e "■"dt= I х(пТд) \Ь(і-ПТл)е**Л. -СО П =-00 П — со -со Откуда, с учетом «фильтрующего свойства» 8-функции, окон¬ чательно приходим к соотношению: *и (■/«>)= Z хЛпУ~ІЫПТл- (1-21) п =—СО Это выражение называется прямым дискретным по времени преобразованием Фурье (ДВПФ). Раскладывая экспоненты по формуле Эйлера, нетрудно пока¬ зать, что для вещественных последовательностей: п*.(-н=|*.0Ч а22) [arg[Хя (>)] = -arg[Хв (>)], т. е. спектральная плотность амплитуд вещественной последо¬ вательности является четной, а спектральная плотность фаз — нечетной функцией частоты. При выполнении условий теоремы отсчетов спектральная плотность XJJ(o) в основной полосе (-л/Гд < со < л/Гд) совпадает с точностью до множителя 1/ГД со спектром непрерывного сиг¬ нала, что позволяет рассчитать последний через отсчеты: со хО®) = ТА Е хд(")е п=~ СО (1.23) 26
Обратное ДВПФ получается из формулы обратного преобразо¬ вания Фурье подстановкой t = пТд: 1 00 х(пТл) = Х~1 (1.24) Z7T —со При выполнении условий теоремы Уиттекера — Котельни¬ кова — Шеннона бесконечные пределы интегрирования можно а спектральную плотность заменить на конечные 71 71 Т Т \ д д У непрерывного сигнала — на спектральную плотность МИП В результате получим выражение для обратного ДВПФ в виде хАп) Z7Z л/ (1.25) Обратим внимание на то, что в выражениях ДВПФ (1.21) и (1.25) спектральная плотность дискретного сигнала остается непрерывной функцией частоты и может быть вычислена через дискретные отсчеты. Отметим также, что выражения (1.23) и (1.25) справедливы только для сигналов с ограниченным спектром, в то время как (1.21) справедливо всегда, в том числе и при наличии эффекта наложения. 1.2. Дискретное преобразование Фурье Рассмотрим теперь особенности спектрального представления сигнала, заданного на конечном интервале наблюдения, длитель¬ ностью Т. После дискретизации такой сигнал хд(п) на отрезке [0, 7] представляется конечным числом отсчетов: х(0), х(1), ..., xN_l, взятых через интервал дискретизации Тд. Полное число отсчетов N = Т/Тд. Считается, что никакой другой информации о спектральных свойствах сигналаха(п), кроме этих отсчетов, нет (рис. 1.4, а). Методика изучения таких дискретных сигналов — искус¬ ственная периодизация с последующим разложением в ряд Фурье 27
дискретного периодического сигнала, точнее, его модели — моду¬ лированной импульсной последовательности (рис. 1.4, б). МИП на интервале наблюдения описывается выражением (0 = Z "Та )• (1 26) и=О Соответствующая ему искусственно периодизированная МИП хи||(/), как любой периодический сигнал, может быть представлена комплексным рядом Фурье: СО .271 ■ <127> к=-со Коэффициенты этого ряда вычисляются по общей для преоб¬ разования Фурье формуле Ck=~]xn{ty~jfk,dt. (128) Рис. 1.4. Дискретный сигнал на конечном интервале, его искусственная периодизация (а) и соответствующие им модулированные импульсные последовательности (б) 28
В этом выражении на интервале интегрирования периоди¬ ческая МИП совпадает с непериодической. Поэтому, подставляя в (1.28) выражение (1.26), получим: с» =7ІІІ х(пЩі-пТ,)е^“л. (1.29) Г о я=0 Откуда, с учетом «фильтрующего» свойства 8-функции, полу¬ чаем: .2*. _ Ск = тЬх\п)* ■ (1.30) 1 п = 0 В этом выражении связь с временным масштабом определяет только сомножитель ИТ. Пр ямым дискретным преобразованием Фурье (ДПФ) называют выражение Х(к) = £х(п)е N . (1.31) п = 0 Сравнивая два последних выражения, приходим к выводу, что коэффициенты ДПФ Х(£) — это коэффициенты разложения в ряд Фурье периодического импульсного сигнала (МИП), пло¬ щади импульсов которого равны х(п). Обратное дискретное преобразование Фурье (ОДПФ) определяется выражением Й") = 77ІХ(^, (1.32) М к=0 .2л, J—кп где для краткости ядро преобразования обозначено как Wv = e N Справедливость (1.32) доказывается путем прямой подста¬ новки в него выражения прямого ДПФ: •(")=£ к=0 jV w=o km i N-1 N-1 m=0 k(n-m) N (1.33) *=о Внутренняя сумма представляет собой сумму N членов геоме¬ трической прогрессии, поэтому Т,щ к=0 к(п-т) N 1-е j2nN(n-m) 2п, \ /—(п-т) I-е'N Г N, при п = т; [О, при п Ф т, что и превращает выражение (1.33) в тождество. 29
Отметим основные свойства ДПФ: 1. Линейность преобразования. Коэффициенты ДПФ дискретной последовательности у(п) = = ах,(п) + Ьух{п) определяются суммой ДПФ коэффициентов: У (к) = аХ,(к) + ЬХ2(к). 2. Периодичность коэффициентов ДПФ Число различных коэффициентов Х(к) равно числу отсче¬ тов дискретного сигнала за период N, ибо функция 'WNkn перио¬ дична по к с периодом N. Поэтому X(N) = Х(0), X(N + 1) = Х(1), X(k1+N) = X(k]), к] < N. При этом коэффициент Х(0) — постоянная составляющая, представляет собой среднее значение по всем отсчетам: х(о) N N N-1 п=О 3. Симметрия коэффициентов ДПФ Если х(п) — вещественная последовательность, то коэффици¬ енты ДПФ, номера которых симметричны относительно N12, обра¬ зуют комплексно-сопряженные пары. Для доказательства этого свойства запишем выражение для коэффициента ДПФ с номером N-k: N-1 Х(Л'-і) = £х(„)\ѵу-п п=О С учетом того, что = 1, получаем: //-i Х(М-к) = ^х(п)^”=Х\к). (1.34) п=О То есть Г|Х (JV — Аг)| = |Х (А-)|; {ф(іѴ-Аг) = -ф(Аг). (L35) 4. ДПФ сдвинутой последовательности. Если у(п) — последовательность, образованная путем сдвига периодической (с периодом N) последовательности х(п) на т отсчетов (т < N), то ее ДПФ-коэффициенты Y(k) = X(k)W%”. (1.36) 30
Действительно: \(к) = ^х(п + т)Щкп. п=0 Произведя замену переменных п1 = п + т и учитывая, что x(AJ1)Wj^H1 — периодическая последовательность с периодом N, ползаем: m+N-l N-1 Y(*)= Z x(«1)w;foI-wr = S^(«i)w^w^=x(^w^' (1.37) nl = m «j = 0 5. ДПФ симметричной последовательности Если х(п) = x(N — п), т. е. дискретная последовательность обра¬ зована путем дискретизации четной функции времени, то ее ДПФ- коэффициенты являются вещественными. (Доказывается путем представления + W^”no формуле Эйлера.) 6. Сдвиг коэффициентов ДПФ. Определим, какой дискретной последовательности соответ¬ ствуют коэффициенты ДПФ, сдвинутые по частоте на / отсчетов: ¥(*) = Х(* + /) = U х(п) \Ѵ/+,)я = |>(и)W/W^. п = 0 п=О Таким образом, коэффициенты ДПФ, сдвинутые по частоте на I отсчетов, соответствуют дискретной последовательности, умноженной на WjJ". Такая последовательность образуется в результате операции цифрового гетеродинирования: у(п) = х{п) \Уд/,г. (1.38) 7. ДПФ круговой, или циклической свертки двух последова¬ тельностей. Круговой, или циклической сверткой двух перио¬ дических (с периодом N) дискретных последовательностей х(п) и h{n) называется последовательность у{п). образованная следую¬ щим образом: N-1 у(п) = 'YJx{tn)h{n-m)\ т, п = 0,1,..., N-L (1.39) т= О Последовательность у(п) — периодическая, с тем же периодом. 31
Выражение для коэффициентов ДПФ циклической свертки найдем через соответствующие обратные ДПФ: i N-1 *M = T7lxmw;*; k=0 (1.40) 1 N-1 /= 0 (1.41) Подставив (1.40) и (1.41) в (1.39), получим: у{п)=И т= 0 N-1 N-1 i N-1 —Zx(t)wi к= О ЛГ /=о 1 ІѴ-1ІѴ-1 1 =^ZZXMH(')W *4 к=01=0 ZW, т= 0 т{к~1) N Поскольку то ЛГ-1 /и= О Г7ѴД = /; [0, кФІ, У{п) тЕх«н(*Х Откуда следует, что Y(£) = Х(£)Н(£), (1.42) (143) т. е. ДПФ циклической свертки равно произведению соответству¬ ющих коэффициентов ДПФ Последнее свойство является очень важным для цифровой обработки сигналов, поскольку позволяет использовать ДПФ для вычисления реакции на выходе цифрового фильтра, которая, как будет показано ниже, описывается линейной сверткой входной последовательности и импульсной характеристики фильтра. Связь круговой (циклической) свертки с линейной (апериоди¬ ческой) сверткой устанавливается в разд. 2.3. 32
В заключение установим связь между спектральной плотно¬ стью непериодического аналогового сигнала X(j'o)), спектральной плотностью МИП XJjcn) и коэффициентами Д11Ф Х(к). Из основ теории спектрального представления аналоговых сигналов известно, что коэффициенты разложения в ряд Фурье периодического сигнала определяются значениями спектральной плотности соответствующего непериодического сигнала на часто¬ тах, кратных частоте повторения [2]. В частности, для модели дис¬ кретного сигнала — периодической МИП: с tXx„U^k)(1.44) Как показано при определении ДПФ, коэффициенты разло¬ жения в ряд Фурье периодического импульсного сигнала (искусственно периодизированной МИП) вычисляются по фор¬ муле С,=^І>(и)е " • (!-45) ■* п = О Из сопоставления этого выражения с выражением для коэффи¬ циентов ДПФ (1.31) следует, что коэффициенты ДПФ представ¬ ляют собой отсчеты спектральной плотности непериодического МИП на частотах, кратных 2пИ: О -гг Х(к) = Хии—к). (1.46) Откуда с учетом выражения (1.10) получаем связь между коэф¬ фициентами ДПФ и спектральной плотностью аналогового сиг¬ нала: 1 Д п=-со Приведенные соотношения проиллюстрированы на рис. 1.5. Лук~у») (1.47) 33
0 x„ 0 Хи (/'“) д д Хи(/'®) д д X(k) /Т\ /Т\71\ ІШІ 1/г k -1 01 N N 2 Рис. 1.5. Соотношения между спектрами аналогового сигнала, импульсной последовательности и дИФ: а — аналоговый сигнал; б — МИИ; в — периодическая МИИ; г — дискретный сигнал на конечном интервале а t x и в t д x д 34
1.3. Алгоритмы быстрого преобразования Фурье 1.3.1. Идея быстрого преобразования Фурье ДПФ является эффективным инструментом спектрального анализа дискретных сигналов. Коэффициенты ДПФ представ¬ ляют собой отсчеты спектра дискретного сигнала на частотах, кратных частоте дискретизации, деленной на число отсчетов. При выполнении условий теоремы отсчетов (Уиттекера — Котельни¬ кова — Шеннона), т. е. при отсутствии наложения спектров, эти коэффициенты являются отсчетами спектра исходного аналого¬ вого сигнала. Алгоритмы цифрового спектрального анализа являются базо¬ выми при решении таких задач ЦОС, как реализация ЦФ в частот¬ ной области, согласованная фильтрация, обработка изображений и т. и. Поэтому на ранней стадии развития теории цифровой обра¬ ботки сигналов столь большое внимание уделялось разработке эффективных в вычислительном отношении алгоритмов ДПФ. Непосредственно вычисление ДПФ в соответствии с базовым соотношением требует выполнения N комплексных умножений и N-1 комплекс¬ ных сложений для каждого коэффициента Х(£). Общее число вычислений составляет N2 комплексных умножений и N(N 1) ~ N2 комплексных сложений. Реализация такого объема вычислений при обработке больших массивов сигналов в реальном времени сопряжена с определенными трудностями. Поиски более эффективных путей решения этой задачи при¬ вели к созданию быстрых алгоритмов, под которыми понимают описание вычислительной процедуры, которое не явля¬ ется очевидным способом вычисления в соответствии с прямой записью алгоритма. Р. Блейхут отмечал, что «как правило, быстрый алгоритм жертвует концептуальной ясностью вычислений в пользу их 35
эффективности» [3]. Он же иллюстрирует это положение на следу¬ ющем простом примере. Пусть требуется произвести вычисления по формуле А = ас + ad + be + bd. Для их реализации требуется 4 умножения и 3 сложения. Однако тот же результат можно получить, проведя вычисления, преобразовав формулу к виду: А = (а + b)(c + d), и затратить на это только одно умножение и два сложения. Таким образом, быстрые алгоритмы можно представить себе как «хитроумную расстановку скобок в вычислениях». Однако для сложных задач быстрые алгоритмы не удается получить простым просмотром вычислений, их построение бази¬ руется на теории чисел. Быстрое преобразование Фурье (БПФ) основывается на воз¬ можности представления размерности массива сигналов N в виде произведения сомножителей rtи выполнении ДПФ для более корот¬ ких последовательностей, число членов в которых определяется соответствующими сомножителями. Коэффициенты ДПФ исход¬ ной «длинной» последовательности получаются путем комбина¬ ции коэффициентов ДПФ коротких последовательностей. Сомно¬ жители ^называются при этом «основанием» БПФ. Оказывается, что, если N = г/2, ...,гр, (1.48) то Х(к) могут быть найдены интерактивно путем расчета суммы р слагаемых следующего типа: дискретных преобразований Фурье размерности г,, общим количеством Nh\ (по /у комплексных умножений в каждом); дискретных преобразований Фурье размерности г2 общим количеством N/r2 (по г~: комплексных умножений в каждом); - дискретных преобразований Фурье размерности г , общим количеством N/r (по гр комплексных умножений в каждом). 36
Таким образом, общее число операций комплексного умноже¬ ния составит: N , N N ляет: К + — К + ... —r;=NZr,. Г1 Г2 ГР '=l Коэффициент ускорения вычислений (КУВ) при этом состав- TV2 N КУВ = р р i=1 1=1 (1.49) В частности, если TV = 2Р, то Zr, = 21og2TV і=і TV2 КУВ = (1.50) 21og2TV Дополнительное увеличение скорости вычислений происходит за счет того факта, что при г = 2 .271, , /—КП W2 = е 2 = ±1 и соответствующие умножения заменяются на сложения. 1.3.2. БПФ с основанием 2 В этом случае длина последовательности TV = 2Р. Методика получения быстрого алгоритма для последовательности такой длины позволяет наглядно продемонстрировать, как и за счет чего получается сокращение вычислительных операций, однако не дает общих правил получения быстрых алгоритмов для последователь¬ ностей произвольной длины. Такая методика рассматривается в разделе 1.4. Алгоритмы БПФ с основанием 2 разделяются на две группы. Если при реализации алгоритма требуется перестановка отсче¬ тов входной последовательности х(п), то его называют алгорит¬ мом с прореживанием по времени. Если при реали¬ зации алгоритма осуществляется перестановка отсчетов выходной последовательности, т. е. коэффициентов ДПФ Х(к), то его назы¬ вают алгоритмом с прореживанием по частоте. 37
По требуемому количеству комплексных умножений и сложе¬ ний эти две разновидности алгоритмов БПФ эквивалентны. Алгоритм БПФ с прореживанием по времени получается следующим образом. Разобьем входную последовательность х(п) на две части — с четными и нечетными номерами: хч» = х(2п): xjn)=x(2n + 1), и = 0, 1, •••, — - 1- 2 (1.51) Процедура этого разбиения для вычисления восьмиточечного БПФ приведена на рис. 1.6. xjli) б п xJn) 0 12 3 п Рис. 1.6. Процедура разбиения входной последовательности (а-в) для восьмиточечного ДПФ 38
Тогда ДПФ исходной последовательности, определенное как X(A:)=Ix(«)WJvfel, „ = 0, 1,...,У-1, £=0, 1,...,7V-1, n=0 можно разбить на две части следующим образом: Х(к) = [£ x(/i)W-fa + Е х(и)\ѵ;Ч п=О п=О где п в первом слагаемом четные, а во втором нечетные. Или: f-1 fi X(*) = |X x(2«)WJvfel + W/ Хх(2и + l)Ww2fa|, -fo? , wj-к ' п=О £ = 0, 1, 1. Учтем, что Тогда \у 2кп = w кп (1.52) (1.53) 7V X(£) = X4T(£) + W/XH4(*), £ = 0,1,...,_-1. (1.54) Таким образом, первая половина коэффициентов ДПФ исход¬ ного последовательности вычисляется через коэффициенты ДПФ дет последовательностей половинного длины, полученных из исходной путем прореживания. Вторую половину коэффициентов можно получить, учтя, что N' Хг|і(Аг) и Хич(к) — периодические функции с периодом — : N N -(*+—) N Х(к + -) = ХЧТ(к + -) + Щ, 2 Хнч(£ + —). Окончательно: Х(* + у) = Хчт(£) - W^XH4(£). (1.55) Соотношения (1.54) и (1.55) являются основой алгоритма БПФ с прореживанием по времени и поэтому носят название базовой операции. 39
Ее удобно представлять направленным графом (рис. 1.7, а). По его виду базовую операцию БпФ с основанием 2 называют «бабочкой». для сокращенного обозначения умножение на 1 опус¬ кают и договариваются всегда в правом верхнем углу записывать сумму, а в нижнем — разность. стрелки означают умножение на число, записанное над ней (рис. 1.7, б). Хчт(к) Х(к) ХчТ(к) X(k) а б рис. 1.7. графическое представление базовой операции БпФ с прореживанием по времени: а — полное; б — сокращенное таким образом, коэффициенты восьмиточечного дпФ вычис¬ ляются через коэффициенты двух четырехточечных дпФ (рис. 1.8). Х(0) Х(1) Х(2) Х(3) Х(4) Х(5) Х(6) Х(7) 40 рис. 1.8. вычисление коэффициентов восьмиточечного дпФ через коэффициенты четырехточечных дпФ
Дальнейшие вычисления строят по итерационному принципу: последовательности отсчетов с четными и нечетными номерами вновь разбивают на две части и продолжают процесс до тех пор, пока не получится последовательность из двух элементов. Так, N — -точечные ДПФ могут быть представлены как комбинации двух —-точечных ДПФ: Хчт (к) = А (к) + W/B 0,1,..., ~- 1; о Хчт (k + N I А) = А(&)-W^2*B(A:). (1.56) N Здесь А(£) и В(£) — коэффициенты —точечного ДПФ после- 4 довательностей, составленных из четных и нечетных членов последовательности хчт(п). Двухточечное ДПФ последовательности /(0), Д1) может быть рассчитана без умножений: F(0)=X0)+/1); F(l)=/0)-Xl). В качестве примера на рис. 1.9 приведен полученный таким образом граф восьмиточечного БПФ. Х(0) Х(1) Х(2) Х(3) Х(4) Х(5) Х(6) Х(7) Рис. 1.9. Граф алгоритма восьмиточечного БПФ с прореживанием по времени 41
Отметим некоторые характерные свойства этого алгоритма: 1. При его реализации на каждом этапе входная (временная) последовательность разделяется на две последовательности поло¬ винной длины, т. е. происходит прореживание, откуда и сле¬ дует название алгоритма. 2. Число этапов равно logJV. 3. Базовая операция каждого этапа — «бабочка» ДПФ: X = А + W^B (первая половина коэффициентов ДПФ); Y = А - W^B (вторая половина коэффициентов ДПФ). 4. Для каждой базовой операции требуется только одно комп¬ лексное умножение для вычисления произведения W^B. Базо¬ вая операция позволяет толковать алгоритм БПФ как комбинацию ДПФ коротких последовательностей с умножением на поворачи- -і—к вающиеся множители (множители поворота) \\ /' = с ' ѵ . _ „ N 5. Всего на каждом этапе выполняется — комплексных умно¬ жений. N 6. Общее число комплексных умножений составляет — log-,A (в это число входят и тривиальные умножения на ±1, [/). Число нетривиальных комплексных умножений еще меньше. Так, в восьми¬ точечном БПФ только 2 нетривиальных умножения (на Wx' и \УХ3). Коэффициент ускорения вычислений составляет: 2 N2 2 N КУВ = Mog2iV log2iV Алгоритм БПФ с прореживанием по частоте получают, раз¬ бивая исходную последовательность х(п). п = 0,N -1 на две после¬ довательности по — отсчетов следующим образом: х,(«) = х(п); N N х2(п) = х(п + —),п = 0, 1,1. Тогда ДПФ исходной последовательности можно представить как 42
N X(k) = ±x(vW~"k + Ntx(n)W1 пк N > к = 0, N-\ или "-i 2 AT Х(Аг) = X *i(« W* + X х2(п)щ о 1 , AT „ > 4__ -f п+—)к -пк х ХГ /„\Х\Т 2 N > п=0 А: = 0,іѴ - 1 • Учитывая, что к , Г1 для четных к: W 2 = е--7 = , , 'v I -1 для нечетных к, получим, что четные и нечетные коэффициенты ДПФ вычисля¬ ются соответственно как "-i 2 Х(2А-) = 2 [Хі(п) + х2(п)]Щм = 2/(л)W я=0 к= оД-1; 2 "-і 2 (1.58) Х(2* + 1) = Z [х,(и)-х2(«№"№1) =Zg(«)W/, А- = 0,^-1. /7=0 п. 2 Таким образом, четные и нечетные коэффициенты ДПФ исход¬ ной последовательности являются коэффициентами ДПФ двух вспомогательных последовательностей половинной длины f(n) и g{n), образованных из первой и второй половины исходной последовательности следующим образом: f(ti) = xl(n) + x2{n); (L59> g(n) = [xl(n)-x2(n)]WNn, n = 0,--l. Нетрудно видеть, что соотношения (1.59) также описы¬ вают базовую операцию — «бабочку», граф которой приведен 43
Соотношения (1.58) и (1.59) также позволяют толковать алго¬ ритм БПФ как сочетание умножения на множители поворота ДПФ над последовательностями половинной длины. Только в алгоритме с прореживанием по частоте умножение на множители поворота Ww” предшествует выполнению короткого БПФ. гг JV N Переходя далее от —точечных ДПФ к —точечным и т. д., 2 4 приходим к двухточечному ДПФ, которое вычисляется без ком¬ плексных умножений. Проиллюстрируем построение алгоритма БПФ с прорежива¬ нием по частоте на примере последовательности из 8 отсчетов. На первом этапе представим восьмиточечные ДПФ через четырех¬ точечные ДПФ (рис. 1.11). Рис. 1.11. Схема вычислений восьмиточечного БПФ с прореживанием по частоте Каждое четырехточечное ДПФ, в свою очередь, может быть представлено через комбинацию: получение вспомогательных последовательностей половинной длины — умножение на множи¬ тели поворота — двухточечное ДПФ. 44
Отметим, что коэффициенты БПФ будут формироваться в пере¬ становленном порядке, т. е. «прорежены» по частоте. Количество комплексных умножений и сложений для этого алгоритма такое же, как и для алгоритма с прореживанием по вре¬ мени. 1.4. Алгоритмы БПФ с произвольным основанием Рассмотрение алгоритма БПФ с основанием 2 позволило про¬ демонстрировать, как и за счет чего получается сокращение вычис¬ лительных операций, однако эта методика не дает общих правил получения быстрых алгоритмов для последовательностей произ¬ вольной длины. Различные алгоритмы БПФ могут быть получены с помощью последовательного применения следующих операций: представле¬ ние одномерного массива чисел {х(п)} двумерным и вычисление соответствующего двумерного БПФ, сводящегося к одномерным БПФ меньшей размерности. Для этого необходимо, чтобы размерность массива была пред¬ ставима в виде произведения (1.48). Если размерность одномер¬ ного массива чисел — простое число, то для такой дискретной последовательности алгоритма БПФ не существует. Формы БПФ различаются в зависимости от количества сомно¬ жителей р в формуле (1.48) и порядка их расположения. Сомножи¬ тель, как отмечалось выше, называют «основанием» БПФ. Таким образом, под алгоритмом БПФ со смешанным основанием пони¬ мают такой алгоритм, когда не все сомножители г; одинаковы. Установим, как, оперируя с двумерным массивом, можно полу¬ чить ДПФ исходного одномерного. Пусть N = LxM. (1.60) Представим одномерный массив х(п) двумерным х(/, т), как показано на рис. 1.12. При этом связь номера одного и того же 45
отсчета в одномерном и двумерном массивах выражается следую¬ щим образом: n = l х M + m. (1.61) здесь l — номер строки; m — номер столбца: 0 < l < L - 1; 0 < m < M - 1. 0 n N - 1 0 < n < N - 1 0 m M - 1 0 s L - 1 r M - 1 0 к N - 1 0 < к < N - 1 Рис. 1.12. Преобразование одномерного массива в двумерный Одномерный массив коэффициентов дискретного преобразова¬ ния Фурье Х(к), 0 < к < N - 1, также представим двумерным масси¬ вом X(r s), где r — номер строки, s — номер столбца: 0 < r < M - 1; 0 < s < L - 1. l L - 1 При этом к = r х L + s. (1.62) 46
Выразим коэффициенты одномерного ДПФ \(к) через эле¬ менты двумерного массива х(/, т). С учетом соотношений (1.61) и (1.62) получаем: (ч ЛГ-1 L-1 М-1 Г, sU^x^W^ =Х£^(^^-(г™+т).(1.63) ' и = 0 1=0 т = 0 Поскольку WjV = ехр(-/2л/7Ѵ), то Ww 'ow = \Vv-'w = exp(-/2nrl) = 1; 'WNrmL = cxp(-j2nrmL / TV) = ехр(-/2тгт?/M) = WA/™; = cxp(-j2nlsM / TV) = ехр(-/2тг/.у / L) = W“b, что позволяет представить выражение (1.63) в следующем виде: М-1 L-1 X(^) = XWM-"W^^(^)W;b (1.64) т = О /=0 В этом выражении внутренняя сумма q m(5) = yx(/,m)wL“fc (1.65) 1 = 0 представляет собой /.-точечное ДПФ /и-го столбца двумерного массива х(/, т). Затем полученные в результате преобразования каждого из М столбцов коэффициенты ДПФ умножаются на мно¬ жители поворота, образуя вспомогательный массив: h(s,iH) = qm(*)WJTM. (1.66) При этом внешняя сумма в выражении (1.64) является M-точечным ДПФ 5-й строки двумерного массива h(v. т): М-1 X(r,5) = £h(5,W)WM-™. (1.67) т=О Вычисления коэффициентов ДПФ в соответствии с выражениями (1.65), (1.66), (1.67) составляют суть обобщенного алгоритма БПФ (алгоритма Кули — Тьюки с произвольным основанием и множите¬ лями поворота). Блок-схема алгоритма приведена на рис. 1.13. 47
Рис. 1.13. Блок-схема алгоритма БПФ с произвольным основанием Вычисления в соответствии с блок-схемой на рис. 1.13 состоят из: — Z-точечных ДПФ (М преобразований, по L2 комплексных умножений в каждом); 48
- L комплексных умножений на множители поворота (М раз); - М-точечных ДПФ (L преобразований, по М2 комплексных умножений) Общее число комплексных умножений составляет: ML2 + ML +M2L = N(L +M+ 1), а коэффициент ускорения вычислений определяется как N2 N КУВ= , - = . (1.68) N(L+M +1) L+M + 1 Эффективность вычислений возрастает, если описанную про¬ цедуру можно применять рекурсивно, т. е. когда N разлагается на большое число сомножителей. Изменение порядка суммирования в исходной формуле (1.63) приводит к алгоритму, в котором операция умножения на повора¬ чивающие множители предшествует вычислению первого ДПФ (аналогично алгоритму с «прореживанием по частоте» для БПФ с основанием 2): X(r,s) = !>,(/, s)WL-b, (1.69) где /=0 М-1 q1(/,5)=Xh1(/,5,«)WM'm т = О (1.70) hl(l,s,m) = x(l,m)WNsm. (і.7і) Этот алгоритм обеспечивает такое же ускорение вычислений, но требует при своей реализации дополнительной памяти для про¬ межуточных результатов. Более эффективные алгоритмы БПФ, не содержащие мно¬ жители поворота, строятся на основе отображения одномерной последовательности в /\-мерную, в соответствии с так называемой китайской теоремой об остатках для целых чисел. Такое отображе¬ ние требует попарной взаимной простоты сомножителей N в про¬ изведении N = TVjAT ..., Nk. Этот способ отображения приводит к более сложной перестановке входной и выходной последователь¬ ности по сравнению с правилами (1.61) и (1.62), использованными в алгоритме Кули — Тьюки. 49
1.5. Основы теории Z-преобразования В анализе и синтезе дискретных и цифровых систем Z-преобразование играет такую же роль, как преобразование Лапласа для непрерывных систем. Это объясняется следующими причинами: — свойства дискретных последовательностей можно изучать, исследуя их Z-преобразования (обычными методами математиче¬ ского анализа); — при Z-преобразовании разностные уравнения, описываю¬ щие линейные дискретные фильтры, преобразуются в алгебраиче¬ ские, таким образом, Z-преобразование является удобным аппара¬ том для решения разностных уравнений; — свойства линейных дискретных фильтров полностью опи¬ сываются расположением нулей и полюсов системной функции на комплексной z-плоскости. Пр ямым односторонним Z-n реобразованием дискретной последовательности х(п) (конечной или бесконечной) называют ряд по степеням комплексной переменной z = а +ур: со A'(z) = Z{x(n)} = ^x(w)z“". (1.72) П-О Множество значений z, где ряд сходится, называется обла¬ стью сходимости. В этой области X(z) — аналитическая функция, не имеющая особых точек. Для равномерной сходимости доста¬ точно, чтобы \x(n)\<Mal, (1.73) где М и аи — положительные вещественные числа. Область сходимости определяется кругом радиуса R в z-плоскости, вне которого ряд сходится. Примеры Z-преобразований тестовых дискретных последова¬ тельностей: 1. Единичный импульс (дискретная дельта-функция): |х(«)} = 5(«) = {1,0,0,...}, X(z) = 1. (1.74) 50
2. Единичный скачок (функция включения): {x(w)} = t/(w) = {l,l,l,...}, Х(г) = ^г-Я =—(1.75) п=0 Z -I Область сходимости: |z| > 1. 3. Экспоненциальная последовательность: оэ I х(п)} = а”, п> 0, X(z) = y(az4) . П = О При \az~l\ < 1 этот ряд представляет собой бесконечно убыва¬ ющую геометрическую прогрессию, сумма членов которой равна: X(z) = —. (1.76) z-a Таким образом, область сходимости \z\ > |a| лежит вне окружно¬ сти радиуса а на комплексной z-плоскости (рис. 1.14). Im z Рис. 1.14. Область сходимости Обратное Z-npeобразование ставит в соответствие функции комплексной переменной X(z) решетчатую функцию (дискретную последовательность) х(п), определяемую следующим образом: х(п) = Z 1 {x(z)} = — j>X(z)zn-ldz, 2 ту" (1.77) 51
где с — замкнутый контур в z-плоскости, охватывающий все осо¬ бые точки (полюсы) функции X(z); интегрирование по контуру с производится в направлении против часовой стрелки. Обратное Z-преобразование удобно вычислять при помощи теоремы о вычетах: функция х(п) равна сумме вычетов подынте¬ гральной функции (1.77) в полюсах, расположенных внутри кон¬ тура интегрирования. Если подынтегральная функция в (1.77) может быть представ¬ лена в виде X(z)z-1=Xn(z)= t F(Z) , (1.78) YKz-zf' 1=1 где z,. — полюс подынтегральной функции; тІ кратность полюса; к — количество полюсов, то к x(”)=Z5®^o^- (L79) i=i " 1 Причем Re^o(z)] = -J-T77lun^r[^o (z-zj*]. (1.80) \ші — 1)! dz^ В частном случае, если X(z) — рациональная функция, имею¬ щая простые полюса, то ее можно разложить на простые дроби: Ш = ±Т^Ч, (1.81) 7~il- a,.z Тогда обратное Z-преобразование: к х(«) = ^Д<. (1.82) і=і Другой часто применяемый способ вычисления обратного Z-преобразования — использование таблиц соответствия несколь¬ ких базовых пар преобразований (табл. 1.1) и их комбинации на основе использования основных свойств Z-преобразования. 52
Таблица 1.1 Базовые пары преобразований Последовательность х(п) Z-преобразование X(z) 8(и) 1 U(n) 1 1 — z 1 ап 1 i i i г, \a\ < 1 1 -az 1 1 1 -1 z п О-')’ -1 z rut (i-oz1)2 QІ<йпТ 1 (l-eJ"V) sinlcowT) sin(co T^z-1 z“2 -2cos(cor)z_1 +1 cos(co«7) (1 -cos(co T))z_1 z-2 -2cos(ror)z_1 +1 Основные свойства Z-преобразования: 1. Линейность. Если у{п) = а1х1 (л) + аух2(л), где а] и а2 — посто¬ янные коэффициенты, не зависящие от и, то: Y(z) = a^X^z) + aJC^z). (1.83) 2. Сдвиг последовательности. Если_у(и) = х(п - m)U{n - т), то Y(z) = X](z)zm. (1.84) Z-преобразование задержанной на т отсчетов дискретной последовательности равно произведению Z-преобразования 53
незадержанной последовательности на множитель z m, называемый оператором запаздывания. В частности, задержке на один период дискретизации соответствует умножение на z-1. Поэтому опера¬ тор z-1 применяют для обозначения цифрового элемента задержки на такт в структурных схемах устройств цифровой обработки сиг¬ налов. 3. Умножение на экспоненту: г[а”х(п)] = Х 'z' \а) 4. Умножение на п: Z[nx(n)] = -zd^. (1.85) (1.86) 5. Свертка последовательностей. п п Если 3;(") = ZX](m)x2(n-m) = YJxl(n-m)x2(m)^ причем т=0 т=0 Xj (т) = 0, х2(т) = 0 при т< 0, то f(z) = X](z)X2(z). (1.87) 6. Перемножение последовательностей. Еслиу{п) = хх{п)хп(гі), то 2щ с f z\dv \v J (1.88) Здесь контур интегрирования с лежит внутри пересекающихся областей сходимости У, (ѵ) и У, (z/v). Соотношение (1.88) называется комплексной сверткой. Из этого соотношения можно получить выражение для спек¬ тральной плотности произведения двух дискретных после¬ довательностей. Поскольку при z = е70 и ѵ = е7® соответствую¬ щие Z-преобразования Х] (e,n j и Х2 (Ѵѳ j представляют собой ДВПФ, то из (1.88) следует: f(e,n) = -1-J Х,Ѵ®)Х2[ЛП-®М®, (1.89) 54
т. е. ДВПФ произведения последовательностей х,(и) и х2(п) есть свертка ДВПФ сомножителей. Эта свертка является периоди¬ ческой (циклической) в силу того, что X^JQ.) и %,(/Q) являются периодическими функциями частоты, поскольку представляют собой спектры дискретных последовательностей. В заключение установим взаимосвязь между ДПФ и Z-преобразованием. Рассмотрим периодическую после¬ довательность хп(п) = хп(п + mN). Эта последовательность не может быть представлена через Z-преобразование, так как ряд (1.72) рас¬ ходится. Представление ее через дискретный ряд Фурье описывается коэффициентами ДПФ: xw=2;1x(H)w'it. п-О Z-преобразование непериодической конечной последователь¬ ности х{п), образованной из одного периода хп(п), равно: X(z) = 2 x(n)z~" = 2 x(n)z~". п=0 п=О Поскольку в пределах периода повторения из N отсчетов х(п)=хп(п). то мы приходим к равенству: ,2л& Х(&) = X(z) при z = (19°) .2 лк J— Уравнению z = е соответствуют точки, равномерно рас¬ положенные на окружности единичного радиуса в комплексной z-плоскости (рис. 1.15, б). Поэтому можно говорить, что спек¬ тральная плотность сигнала — это сечение его z-образа по единич¬ ной окружности (рис. 1.15, а). Соотношения между непрерывными и дискретными сигна¬ лами и их преобразованиями обобщены в схеме, представленной на рис. 1.16. 55
Рис. 1.15. Спектральная плотность сигнала и его Z-преобразование (пояснения в тексте) Рис. 1.16. соотношение между непрерывными и дискретными сигналами и их преобразованиями 56
КОНТРОЛЬНЫЕ ВОПРОСЫ 1. Дискретизация аналогового сигнала x(t) производится с перио¬ дом Тя. Чему равно значение решетчатой функции х(п), описывающей дискретный сигнал, на интервале пТ < t< (п + 1)Т ? 2. При каком условии последовательность, полученная путем дис¬ кретизации гармонического сигнала, не является периодической? 3. Какие операции осуществляются при переходе от дискретного сигнала к цифровому? 4. Запишите математическую модель идеальной дискретизации. 5. Максимальная частота в спектре звукового сигнала равна 20 кГц. Каков должен быть минимальный период дискретизации в АЦП, чтобы эффект наложения отсутствовал? 6. В чем заключается способ восстановления непрерывного сиг¬ нала по дискретным отсчетам, непосредственно вытекающий из теоремы УКШ? 7. Дискретная последовательность образована путем дискретизации одного периода гармонического колебания. Частота дискретизации равна со/8. Чему равен коэффициент ДПФДО)? 8. Что понимают под термином «алгоритм БПФ»? 9. Каково общее количество комплексных умножений при реализа¬ ции базовой операции «бабочка»? 10. Какая операция лежит в основе построения алгоритма БПФ с про¬ извольным основанием (алгоритма Кули — Тьюки)? 11. Последовательность у (и) образуется как результат свертки двух последовательностей х(п) = (-0,9)" и h(ri). Определите Y(z), если H(z) = 1/(1 -bz ').
2. ДИСКРЕТНЫЕ И ЦИФРОВЫЕ ФИЛЬТРЫ 2Л. Линейные дискретные фильтры и их характеристики Под дискретным фильтром (ДФ) в общем случае понимают систему, преобразующую одну дискретную последовательность в другую (рис. 2.1). Соответственно — цифровой фильтр (ЦФ) есть система, преобразующая один цифровой сигнал в другой. Цифро¬ вой фильтр реализуется как программа на ЦВМ или аппаратным способом в виде цифровой схемы, содержащей регистры, сумма¬ торы, умножители и другие вспомогательные элементы. х(гі) Дискретный К”) Вход (цифровой) фильтр Выход Рис. 2.1. К определению дискретного и цифрового фильтра Как отмечалось во введении, в настоящем учебном пособии сначала рассматриваются цифровые фильтры без учета эффек¬ тов квантования. То есть вместо понятий «дискретный фильтр» и «цифровой фильтр» используется ниже только понятие «цифро¬ вой». В гл. 3 отдельно анализируются эффекты, связанные с кван¬ тованием сигнала, ограничением разрядности коэффициентов фильтра и округлением промежуточных результатов. В настоящем разделе мы ограничим класс преобразований входных последовательностей и будем рассматривать только линейные стационарные, физически реализуемые цифровые фильтры 58
ЦФ называется линейным, если выходная последователь¬ ность у{п) при нулевых начальных условиях при воздействии вида х(п) = аххх(п) + аус2{п) описывается как у(п) = axyx(n) + ау\ф), где ух{п) иу2(п) соответственно отклики ЦФ на хх(п) и х2(п). То есть при преобразовании входных последовательностей в линейном фильтре выполняется принцип суперпозиции. Линейный ЦФ называется стационарным или инва¬ риантным во времени (ЛИВ), если откликом на воздей¬ ствие х{п) =хх{п- ла) является у(п) =ух(п - Ио). Иными словами, сдвиг входной последовательности приводит к такому же сдвигу выходной последовательности без изменения ее формы. Связь между входной х(п) и выходной у{п) последовательно¬ стями в стационарном ЦФ описывается линейным разностным уравнением с постоянными коэффициентами вида М N ZamT(«-W) = ZM"-0- (2.1) т= 0 /=0 Откуда, полагая а0 = 1 (нормируя относительно а0), получаем алгоритм работы цифрового фильтра во временной области: 1V м y(n) = YJbix(n~i)-'Lamy(,1-m)’ (2.2) і=0 т=\ Текущий отсчет выходной последовательности определяется текущим и N предыдущими отсчетами входного сигнала, а также М предыдущими отсчетами выходной последовательности. В частном случае, когда текущий отсчет выходной последо¬ вательности определяется только отсчетами входного сигнала 59
(когда все ат = 0), алгоритм работы цифрового фильтра прини¬ мает вид: N у{п) = ^Ь,х{п-і). (2.3) і = О Уравнение (2.2) описывает так называемый рекурсивный фильтр, а уравнение (2.3) — нерекурсивный фильтр. Разност¬ ные уравнения непосредственно определяют способ построения ЦФ. В дальнейшем, чтобы отличать коэффициенты рекурсивного фильтра от коэффициентов нерекурсивного фильтра, для обозна¬ чения последних будем использовать символы с,. Импульсной характеристикой ЦФ h(n) называют реакцию нулевого состояния на воздействие в виде дискретной 5-функции. Из уравнений (2.2) и (2.3) следует, что рекурсивный фильтр имеет бесконечную импульсную характеристику (БИХ), а нерекурсивный фильтр — конечную (КИХ). С помощью импульсной характеристики можно получить опи¬ сание выходной последовательности при любом входном воздей¬ ствии. Пусть в стационарном ЦФ откликом нулевого состояния на 8(п) является последовательность h(n). Тогда, на основе инвари¬ антности, откликом на Ъ{п - т) будет h(n - т). Входную последова¬ тельность представим в виде динамического соотношения °о ^ х(п)= Y х(т)8(п-т). (2.4) т=- оо Тогда, на основании свойства линейности, можно утверждать, что отклик нулевого состояния стационарного ЦФ: оо оо у(п)= Y x{rn)h{n-m)= Y x(n-m)h(m), (2.5) m=-cfj rn--<x> T. e. выходной сигнал представляет собой дискретную свертку входного сигнала и импульсной характеристики ЦФ h(n). Линейный ЦФ называется устойчивым, если воздействие любой ограниченной входной последовательности дает ограни¬ ченную выходную последовательность, такую, что X |у(ц)|' <00. (2'6) 60
ЦФ устойчив, если его импульсная характеристика удовлетво¬ ряет условию и Х|й(и)|<0°. (2.7) п=—со ЦФ называется физически реализуемым, если отклик не появляется раньше воздействия, т. е., если его импульсная характеристика h{n) = 0 при п < 0. (2.8) Для физически реализуемого ЦФ соотношение (2.5) переписы¬ вается в виде я у(п)= X x(m)h(n-m). (2.9) т=—со Если входная последовательность имеет начало, т. е. х(т) = 0 при т < 0, тогда выходная последовательность п у(п)= Y,x(m)h(n-m). (2.10) т=0 То есть выходная последовательность физически реализуе¬ мого стационарного линейного фильтра представляет собой апе¬ риодическую дискретную свертку входной последовательности и импульсной характеристики. По аналогии с системной (передаточной) функцией и частот¬ ными характеристиками аналогового фильтра определяются пере¬ даточная функция и частотные характеристики ЦФ. При этом используется аппарат Z-преобразования. Системная (передаточная) функция ЦФ — отношение Z-преобразования выходной последовательности к Z-преобразованию входной последовательности при нулевых начальных условиях: Ш = Z{y(n)} _ Y(z) Z{x(n)} X(z) Передаточная функция рекурсивного фильтра имеет вид: (211) X b,z 1 = • (2.12) l+z т=1 61
Передаточная функция нерекурсивного фильтра описывается выражением: H(z) = ^c;z/. (2.13) 1 = 0 Из свойств Z-преобразования следует, что свертке после¬ довательностей х(п) и h(n) соответствует произведение их Z-преобразований: Y(z) = X(z)H(z), (2.14) откуда следует, что передаточная функция есть Z-преобразование импульсной характеристики ЦФ: оо H(z) = '£th{n)z~". (2.15) п=0 Условие устойчивости ЦФ Поскольку передаточная функция ЦФ есть Z-преобразование импульсной характеристики, то ее модуль со п-О Если \z~l\ < 1, то со \Н(?)\<^\Кп)\- (2.16) п= О Таким образом, в устойчивом ЦФ передаточная функция H(z) конечна во всех точках z-плоскости, где \z\ > 1, т. e. H(z) может иметь особые точки только внутри единичной окружности. Для описания стационарных цифровых фильтров в частотной области используется специальный класс входных воздействий — дискретные комплексные гармонические последовательности: {-Ф)} = МехрІДсойГ, + (р)]}; Яе{ф)} = {ZcosIXcohT^ + ф)]}; (2.17) Іт{х(и)і = {А$т\]{шТл + ф)]}. Значения этих последовательностей не изменяются при замене со на со + 2л(он/Гд, поскольку 62
(to+ 271П )пТд + ф = СОНГд +Ф+271772. (2.18) Если на вход ЦФ поступает такая последовательность, то на выходе в соответствии с (2.5) получим: СО со у(п) = '^Ь(т)е’ю('п~т)Тд =е'’ю"Тд^Г1к(т)е~'1атТд = х(п)Й(е’аТд). m=Q m=Q Таким образом, для выбранного класса входных последова¬ тельностей отклик совпадает со входной последовательностью с точностью до комплексного множителя: СО Й(е]юТд) = ^h(m)e~jamT-. (2.19) m= О Последнее соотношение представляет собой дискретное по времени преобразование Фурье (ДВПФ) импульсной характе¬ ристики. Эта функция называется комплексной частотной характеристикой (КЧХ) дискретного (цифрового) фильтра, аее модуль и фаза — соответственно амплитудно-частот¬ ной (АЧХ) и фазочастотной (ФЧХ)характеристиками. Комплексная частотная характеристика ЦФ формально полу¬ чается из передаточной функции (2.11) путем подстановки z = е7®7». КЧХ ЦФ является периодической функцией частоты со с перио¬ дом 2л/Та. Для действительных h(n) модуль КЧХ //(со) — четная функция, а аргумент ф(со) — нечетная. Эти свойства проиллюстрированы на рис. 2.2. Требования к фильтрам могут задаваться как во временной, так и в частотной областях, что определяется назначением фильтра и удобством его описания. Так, согласованные фильтры чаще зада¬ ются импульсной характеристикой, а избирательные фильтры — частотными характеристиками. В частотной области требования предъявляются к КЧХ либо к ее составляющим — АЧХ и ФЧХ. Требования к амплитудно-частотной характеристике филь¬ тра, в первую очередь, включают параметры полосы про¬ пускания, полосы подавления и переходной полосы. Под 63
амплитудно-частотной характеристикой затухания А (go) пони¬ мают функцию, обратную //(со). Рис. 2.2. Амлитудно-частотная (а) и фазочастотная (б) характеристики цифрового фильтра Диапазон частот, в котором затухание фильтра минимально (для идеального фильтра равно нулю), называется полосой пропускания. Обычно это диапазон частот, занимаемый пре¬ имущественно полезным сигналом. Диапазон частот, в котором затухание фильтра максимально (для идеального фильтра равно бесконечности), называется полосой подавления (задерживания). Обычно это диапа¬ зон частот, занимаемый преимущественно помехой. Диапазон частот, лежащий между полосой пропускания и поло¬ сой подавления, называют переходной полосой. В зависимости от взаимного расположения полос подавления и пропускания различают следующие типы фильтров: 1. Фильтр нижних частот (ФНЧ) — фильтр с полосой пропус¬ кания от 0 до частоты шс и полосой подавления от со3 до бесконеч¬ ности (оос < со3). 2. Фильтр верхних частот (ФВЧ) — фильтр с полосой пропус¬ кания от частоты оос до бесконечности и полосой подавления от О до ш3(сос>со3). 64
3. Полосовой фильтр (ПФ) — обе границы полосы пропуска¬ ния представляют собой ненулевые частоты сосн, шсв, а с каждой из сторон от полосы пропускания имеется по одной полосе подав¬ ления (от 0 до охн и от созв до да). 4. Режекторный (заграждающий) фильтр (РФ) — фильтр с двумя полосами пропускания (от 0 до соон и от сосв до оо) и одной полосой подавления. 5. Гребенчатый фильтр (ГФ) — фильтр с несколькими поло¬ сами подавления и несколькими полосам пропускания. 6. Всепропускающий фильтр постоянного затухания (ФПЗ) — фильтр с единичной (постоянной) передачей для всех частот (т. е. с полосой пропускания от 0 до оо). Используется для обеспечения требуемой фазовой коррекции и фазового сдвига. Рассмотрим структурные схемы цифровых фильтров. Струк¬ турная схема по сути отображает последовательность вычислений при формировании реакции (отклика) фильтра на входное воздей¬ ствие. Структурные схемы могут быть построены как на основа¬ нии описания алгоритма функционирования фильтра во времен¬ ной области (разностное уравнение), так и на основании описания на комплексной z-плоскости (передаточная функция). Структурная схема нерекурсивного фильтра следует из раз¬ ностного уравнения, которое принято записывать в следующих обозначениях: где L 1 — порядок фильтра. Взяв Z-преобразование от левой и правой частей, получим: 2.2. Формы реализации линейных цифровых фильтров L-1 ѵ(ц) = ^сгх(н-/), (2.20) 1=0 L-1 (2.21) / = 0 65
Откуда следует, что передаточная функция нерекурсивного или трансверсального фильтра описывается полиномом относительно переменной z1: L_x 6(Z) = P'Z"-(2.22) Уравнениям (2.21) и (2.22) соответствует функциональная схема, приведенная на рис. 2.3. Элементами структурной схемы линейного дискретного филь¬ тра являются элемент задержки на период дискретизации, умно¬ житель и сумматор. Рис. 2.3. Функциональная схема нерекурсивного фильтра Для получения системной функции рекурсивного фильтра возьмем Z-преобразование от левой и правой частей разностного уравнения, описывающего этот тип фильтра во временной обла¬ сти: N М у(п) = YsbAn -0 - X а,Мп 0, (2.23) /' = 0 т = 1 в результате получим: N М Y(z) = - X (2.24) / = 0 т = О Откуда следует выражение для передаточной функции рекур¬ сивного фильтра: 66
(2.25) H(z\ = ^(z) _ bo + b\z 1 + b2z 2 + ♦♦♦ + bNz N X(z) 1 + axz~x + a2z~2 +... + aMz~M Прямая форма реализации рекурсивного ЦФ следует непосредственно из разностного уравнения (2.23) или уравнения (2.24). Структурная схема приведена на рис. 2.4. Недостатком реа¬ лизации рекурсивного ЦФ по прямой форме принято считать боль¬ шое число элементов задержки (большой объем памяти регистров). Рис. 2.4. Функциональная схема рекурсивного фильтра (прямая форма реализации) Минимальное число элементов задержки, равное порядку пере¬ даточной функции М- 1, требует для своей реализации канони¬ ческая форма. Для ее получения запишем уравнение Y(z) = X(z)H(z) в виде Y(z) = F(z)B(z), 67
где B(z) = 1 = 0 (2.26) f (z)= тг—ад- (2.27) m= 1 Этому уравнению, заданному на z-плоскости, соответствует разностное уравнение м f{n) = х{п) - £ ат/(« - т). т=1 (2.28) Из выражения (2.26) следует, что у{п) = ИЬіЯп- 0 і-1 (2.29) Этому разностному уравнению соответствует структура, при- веденная на рис. 2.5. Рис. 2.5. Каноническая форма реализации рекурсивного фильтра Каскадная форма реализации соответствует представ¬ лению передаточной функции фильтра в виде произведения H(z)=UHk(z). (2.30) к=1 68
Структурная схема ЦФ при каскадной реализации имеет вид, приведенный на рис. 2.6. Рис. 2.6. Каскадная форма реализации цифрового фильтра На практике обычно используют однотипные звенья второго и первого порядка, передаточная функция которых имеет вид: нЛг)-ьлі±МІ±МІ l + anz- +a2kz -2 (2.31) Эти звенья называются биквадратными блоками. Биквадратный блок является универсальным звеном, пригодным для построения любых ЦФ. Звено первого порядка может быть получено из биквадратного блока при а2к = 0, Ь2к = 0. Каскадная форма реализации позволяет снизить влияние эффектов, связанных с конечной разрядностью представления коэффициентов и округлением промежуточных результатов. 2.3. Реализация линейных цифровых фильтров в частотной области с помощью алгоритмов БПФ Дискретный фильтр в переменных вход-выход описывается уравнением у(п) = ^х{т)к{п-т),п = 0,\,2,...,N^+N2-2, (2.32) т = О представляющим собой апериодическую (линейную) свертку входного сигнала и импульсной характеристики. Линейная свертка у(п) двух непериодических последовательностей х(п) и И(п), содер¬ жащих соответственно Nx и N2 отсчетов, представляет собой конеч¬ ную последовательность длиной /V, i /V, - 1 отсчетов. Другим типом свертки двух последовательностей явля¬ ется круговая (циклическая) свертка. Круговой сверткой двух 69
периодических последовательностей х(п) и h{n) п 0. 1.2. .... N - 1 называют последовательность у(п), образованную в соответствии с выражением (1.39): ы-\ у{п) = ^ х(/и)А(и - т), т,п = 0,1,TV-1. т= 0 Последовательность у(п) также периодическая с тем же перио¬ дом TV. Как было показано при рассмотрении свойств ДПФ, опера¬ ции вычисления круговой (циклической) свертки соответствует в частотной области произведение коэффициентов ДПФ (1.43). Тогда реакция фильтра у{п) может быть найдена как обратное ДПФ: 1 ы-\ Я") = Т7ІХ(£)Н(£)\Ѵ;\ n = 0,...,N-l. (2.33) ™ к = О Для применения этих алгоритмов к описанию процесса преоб¬ разования сигнала в линейном ЦФ, а следовательно, и для приме¬ нения БПФ для реализации ЦФ, необходимо свести апериодиче¬ скую свертку к эквивалентной циклической. С этой целью сформируем вспомогательные периодические последовательности х,(я) и А,(«) длиной по TV, + TV, - 1 отсчетов путем дополнения нулями последовательностей х(п) и h(n) следу¬ ющим образом: (и) /2і(«) Г х{п), и = 0,1,..., ТѴ, -1; {о,я = ТѴ„...,ТѴ,+ТѴ2-2; Гй(и) прин = 0, 1, ..., N2 -1; [О прип = N2, N2+1, ..., TV, + TV, - 2. (2.34) (2.35) Тогда линейная свертка последовательностей х(п) и h(n) на интервале [0, TV, + N2 - 2] будет равна TV, + TV, - 1 — точечной круговой свертке последовательностей х,(н) и А,(и): Nx+N2-2 У(п) = X N(m)hi(n~m),n = 0,\,...,Nl +N2 -2, (2.36) т = 0 и поэтому может быть вычислена с использованием алгоритмов ДПФ. Вывод проиллюстрирован на рис. 2.7. 70
x(m) x1(m) h(m) 0 12 N1 - 1 m T T f ■ i ■ 0 1 2 h1(m) N1 - 1 N1 + N2 - 2 -m ■ illi - m .riii ~m h(n - m) 2 L h1(n - m) г Tit.... ( 1 2 n y(n)> y(n) 0 12 N2 - 1 _Ll_l N + N2 - 2 ' T r ' 1 I T 0 12 n -m 0 12 0 12 N2 - 1 N1 - 1 N1 + N2 - 2 N2 - 1 N1 - 1 N1 + N2 - 2 а б Рис. 2.7. Иллюстрация процесса обработки сигнала в дискретном фильтре: а — апериодическая свертка; б — соответствующая ей циклическая свертка переход к вспомогательным последовательностям, описываю¬ щим как входной сигнал, так и импульсную характеристику цФ, позволяет реализовать последний в частотной области. схема такой реализации приведена на рис. 2.8. рис. 2.8. схема реализации линейного цифрового фильтра в частотной области 71
При вычислении прямого и обратного преобразования целе¬ сообразно использовать алгоритмы БПФ. При таком способе реа¬ лизации ЦФ для вычисления всех 7V, + /Ѵ2 - 1 отсчетов выходной последовательности у(п) требуется значительно меньше операций комплексного умножения и сложения. Повышение вычислитель¬ ной эффективности достигается за счет применения алгоритмов БПФ, а также за счет того, что коэффициенты Н(£), описывающие комплексную частотную характеристику фильтра, могут быть рас¬ считаны заранее и храниться в памяти. 2.4. Цифровой спектральный анализ Как мы установили в гл. 1 при выполнении условий теоремы Котельникова, коэффициенты ДПФ Х(&) с точностью до множи¬ теля 1/Т представляют собой отсчеты спектральной плотности аналогового сигнала на частотах, кратных 2л/7'д. Таким образом, устройство (или алгоритм), реализующее БПФ, можно рассматри¬ вать как спектроанализатор, осуществляющий анализ спектра наТѴ частотах в диапазоне от 0 до частоты дискретизации Fa. Обычно спектроанализатор представляет собой набор фильтров или один перестраиваемый фильтр. Важнейшей характеристикой спектро¬ анализатора является разрешающая способность, определяемая формой частотной характеристики фильтра. Для определения этой характеристики для БПФ установим связь между спектральным измерением и фильтрацией. Рассмотрим нерекурсивный фильтр, схема которого приведена на рис. 2.9. Его импульсная характеристика (отклик на дискретную 8-функцию) представляет собой последовательность коэффициен¬ тов сп. h(n) = с0, сІ5..., cN_i, 0, 0, 0... (2.37) Алгоритм работы такого фильтра во временной области описы¬ вается следующим уравнением: п у(п)= £ X(m)cn_m. (238) m=n-N+\ 72
x{n-N+ 1) x(ri) + У(п) Рис. 2.9. Структурная схема нерекурсивного фильтра Потребуем, чтобы в момент п = N — 1, т. е. после обработки в фильтре N отсчетов, сигнал на его выходе совпадал со значением А>го коэффициента ДПФ: y(N -1) = X cN_x_ = Х(к). (2.39) т = О Для выполнения этого условия, как это следует из сопоставле¬ ния выражений (1.31) и (2.39), необходимо, чтобы коэффициенты фильтра удовлетворяли соотношению cN-,-m= W/\ (2.40) Таким образом, процессор БПФ, рассматриваемый как фильтр с к-м отводом в качестве выхода, описывается структурной схе¬ мой, изображенной на рис. 2.10. Импульсная характеристика этого фильтра, как это следует из выражений (2.38) и (2.40), задается соотношением hk(n)= w;t(JV-1-") = e N ,0 {IA'> 73
Рис. 2.10. Структурная схема нерекурсивного фильтра, реализующего вычисление к-го коэффициента ДПФ Для определения комплексной частотной характеристики (КЧХ) этого фильтра вычислим Z-преобразование его импульсной характеристики: Hk(z) = Yihk(n)z-n =^ехр[]2пк(п + 1)/Ы]-г-п. (2.42) п=0 п=0 Последнее выражение представляет собой сумму N членов гео¬ метрической прогрессии, у которой первый член равен ехр[/27гА;ЛѴ], а знаменатель — ехр|J2nk/N]z~l. В связи с чем ^(z) = exp(-/2Tt%)^ l-z 1ехр(у27г^) ' ехр(у 2л/N) Комплексная частотная характеристика получается из этого выражения путем подстановки z = ехр(/соГ ): Йкт = е*Р(]Ш/м) 1 - ехр[ - jN(соГд + 2 )] 1-ехр[-у(соГд +2 (2.43) Модуль этого выражения — амплитудно-частотная характери¬ стика (АЧХ) БПФ-фильтра по к-му отводу описывается выраже¬ нием 74
(2.44) ~ , ч An[N{o)T+2nk Нt(C!)) = ; . * sin[(oo7;+27t^/iV)/2] Таким образом, процессор БПФ с к-м отводом, рассматрива¬ емым в качестве выхода, представляет собой дискретный фильтр с амплитудно-частотной характеристикой (в функции от нормиро¬ ванной частоты Q = соГд) вида Н( О): sin NC1/2 sinfi/2 сдвинутой по частоте на величину (2.45) Ц = со кТл = 2тг k/N, (2.46) с шириной главного лепестка (по нулевому уровню), равной 4n/N (рис. 2.11). Рис. 2.11. Амплитудно-частотная характеристика БПФ фильтра по к-му (сплошная линия) и к - 1-му (пунктир) отводам Полоса пропускания такого фильтра по уровню половинной мощности, определяющая разрешающую способность спектро¬ анализатора по частоте, составляет: AQ = Acor =0,89— «—. д N N Амплитудно-частотные характеристики фильтров, соответ¬ ствующих соседним отсчетам БПФ, перекрываются не только 75
боковыми лепестками, но и главными практически на уровне поло¬ винной мощности (см. рис. 2.11). Это объясняет наличие ненуле¬ вого отклика на всех отводах процессора БПФ при вычислении ДПФ гармонического сигнала с частотой, некратной 2n/NTa. 2.5. Проектирование цифровых фильтров с конечной импульсной характеристикой 2.5.1. Этапы проектирования цифрового фильтра Процесс проектирования любого ЦФ включает следующие этапы: • анализ требований к ЦФ; • синтез дискретного фильтра (выбор формы реализации, рас¬ чет порядка и определение коэффициентов фильтра, расчет АЧХ и ФЧХ); • определение разрядности коэффициентов фильтра; • квантование входных данных и промежуточных результатов (определение разрядности регистров памяти); • моделирование ЦФ с учетом квантования входных данных и ограничения разрядности коэффициентов и регистров памяти; • выбор элементной базы (типа цифрового сигнального про¬ цессора), аппаратная или программная реализация ЦФ. Как отмечалось ранее, требования к ЦФ могут задаваться как во временной, так и в частотной областях. Так, согласованные фильтры чаще задаются импульсной характеристикой, а избира¬ тельные фильтры — частотными характеристиками. В частотной области требования предъявляются к КЧХ либо к ее составляю¬ щим — АЧХ и ФЧХ. В частности, при задании требований к фильтру нижних частот (ФНЧ) фигурируют следующие параметры АЧХ (рис. 2.12): сос — частота среза, определяющая полосу пропускания фильтра; ш3— граница области затухания; 76
Нс — уровень АЧХ, определяющий неравномерность передачи в полосе пропускания, одной из границ которой является частота среза; Н3 — уровень АЧХ, определяемый гарантированным затуха¬ нием в полосе подавления. Рис. 2.12. Определение требований к фильтру Поведение АЧХ в полосе пропускания и области затухания в задании на расчет не регламентируется, кроме выполнения един¬ ственного, заранее обговариваемого условия — допускается или не допускается наличие пульсаций на этих участках АЧХ. В зави¬ симости от того, как формулируется заданное условие, возможны четыре основных типа аппроксимаций АЧХ: Баттерворта, Чебы¬ шева (первого и второго рода), Кауэра. 2.5.2. Синтез нерекурсивных фильтров методом «окна» Передаточная функция дискретного нерекурсивного фильтра описывается выражением H(z) = X c,z~l, 1=0 где L — 1 — порядок фильтра. 77
Синтез КИХ-фильтра сводится к определению числа отсчетов ИХ L и расчету коэффициентов сг. Из свойств КИХ-фильтра сле¬ дует, что эти коэффициенты — суть отсчеты импульсной характе¬ ристики h(n). Дискретные фильтры с конечной импульсной характеристикой (КИХ) обладают рядом положительных свойств, главное из кото¬ рых — они всегда устойчивы. Кроме того, они позволяют обес¬ печить совершенно линейную фазочастотную характеристику (постоянное групповое время запаздывания). Условия линейности фазочастотной характеристики нерекурсивного фильтра заключа¬ ются в наличии определенного типа симметрии импульсной харак¬ теристики h(n): т. е. импульсная характеристика должна быть симметрична или антисимметрична относительно среднего отсчета с номером При выполнении этого условия ФЧХ КИХ-фильтра будет иметь вид: а групповое время запаздывания будет постоянно и равно: Пример такой характеристики приведен на рис. 2.13. h(n) = ± h(L - 1- п\ п = 0, 1,..., L — 1, (2.48) (Z - 1)/2. (2.49) h(n) п О (L - 1)/2 L- 1 78 Рис. 2.13. Пример симметричной импульсной характеристики нерекурсивного фильтра с линейной ФЧХ
Одним из наиболее часто применяемых методов синтеза филь¬ тров с КИХ является метод «временных окон», или метод «взве¬ шивания». Его суть заключается в получении отсчетов импульсной характеристики конечной длины путем усечения последователь¬ ности, описывающей импульсную характеристику бесконечной длины. Исходной при синтезе является комплексная частотная характеристика идеального фильтра H(j£l), имеющая для дискрет¬ ного фильтра периодический характер по частоте (с периодом, рав¬ ным 2л). Поэтому она может быть представлена рядом Фурье (пря¬ мое ДВПФ): со Щ]П) = (2.50) п = — оо Коэффициенты этого разложения определяются по общим правилам через обратное дискретное по времени преобразование Фурье: Л(и) = J H(jQ)exp(jnQ)dQ. (2 51) Эти коэффициенты и используются в качестве импульсной характеристики фильтра. При этом, однако, возникает трудность, заключающаяся в том, что полученная таким образом импульсная характеристика имеет в общем случае бесконечное число отсче¬ тов, поскольку суммирование в (2.50) осуществляется в бесконеч¬ ных пределах. Это означает, что синтезированный таким образом фильтр является фильтром с БИХ. Кроме того, он физически нереа¬ лизуем, ибо при любом конечном сдвиге L не выполняется условие h{n - L) = 0 при п = L. Для выполнения условий физической реализуемости необ¬ ходимо сделать усечение импульсной характеристики, ограни¬ чив ее L отсчетами, и сдвинуть ее на половину длительности, т. е. на (L - 1)/2. При этом КЧХ фильтра будет аппроксимироваться усеченным рядом Фурье. Однако такое простое усечение вызывает явление Гиббса, проявляющееся в появлении выбросов и колеба¬ ний частотной характеристики вблизи точек разрыва. Величина 79
этих выбросов и пульсаций не уменьшается с увеличением длины импульсной характеристики при условии сохранения ее конеч¬ ности. Пульсации только локализуются в более узком диапазоне частот. Это означает, что простое усечение ряда Фурье для получе¬ ния аппроксимации КЧХ фильтра с КИХ не обеспечивает хороших результатов. Проиллюстрируем этот эффект на примере синтеза идеаль¬ ного дискретного фильтра нижних частот (ФНЧ), КЧХ которого (рис. 2.14) равна единице на интервале [-£2С, QJ и нулю на интер¬ валах [—7г/2, -QJ и [£2с, 7г/2]. Импульсную характеристику ФНЧ получим путем вычисления обратного ДВПФ от КЧХ в соответствии с выражением (2.51): г / ч 1 Qfc /.^4 7^ siiwQ. h{n)-— \ exp(jnQ)dQ = — -. (2.52) 2п -Q л nQr нт i 1—I—I ►Q -71 Q 0 Qc 71 271 Рис. 2.14. Комплексная частотная характеристика идеального дискретного фильтра нижних частот Эта характеристика содержит бесконечное число отсчетов. Замена ее на характеристику конечной длины, т. е. усечение ряда Фурье, приводит к частотной характеристике, характерными осо¬ бенностями которой является конечная переходная полоса и про¬ явление эффекта Гиббса (пульсации, достигающие максимальной величины порядка 9 % вблизи частоты среза). При увеличении длины импульсной характеристики максимальная величина этих пульсаций не изменяется. 80
Это явление можно объяснить следующим образом. Усечение ряда Фурье можно трактовать как умножение импульсной харак¬ теристики с бесконечным числом отсчетов h(n) на прямоугольное окно wn(«), содержащее L отсчетов: К(п) = КФ„(п); (2.53) wu(n) = 1 при L-1 2 <п < L-1 2 ; wu(ri) = 0 при п < L-1 2 ,п> L-1 2 Прямоугольная последовательность ип(и) описывается комп¬ лексной частотной характеристикой (спектральной плотностью) вида ^=(7«) sin(QL/2) sin(Q/2) Эта функция имеет главный лепесток шириной 4n/L, боковые лепестки шириной 2л//.. максимальный уровень бокового лепестка (первого) по отношению к главному -13 дБ. При увеличении L ширина лепестков уменьшается, но их уровень остается без изме¬ нения. Умножению двух последовательностей во временной обла¬ сти соответствует в частотной области циклическая свертка двух комплексных частотных характеристик (идеального фильтра с «неусеченной» импульсной характеристикой и частотной харак¬ теристики «окна») W(/'Q): HJjSl) = H(jQ) ® W(jQ). (2.54) Результат вычисления свертки в каждой частотной точке пред¬ ставляет собой величину площади перекрытия КЧХ идеального фильтра и сдвинутого зеркального отображения КЧХ «окна». Поэтому результирующая КЧХ КИХ-фильтра HJjCl) будет иметь пульсации, максимальные вблизи частоты среза идеального ФНЧ. Указанные преобразования во временной и частотной областях иллюстрируются графиками, приведенными на рис. 2.15-2.18. 81
h(n) 1 Qc/n 1 -2n/fic -n/Q.c 0 л/Пс 2л/Пс Рис. 2.15. Импульсная характеристика идеального фильтра нижних частот W» ЧL - 1)/2 0 (L - 1)/2 Рис. 2.16. Прямоугольная весовая функция WJJCI) Q Рис. 2.17. Комплексная частотная характеристика прямоугольного окна 1 n 82
Из рассмотрения алгебры свертки следует, что полосу филь¬ тра (частоту среза) определяет число отсчетов в главном лепестке импульсной характеристики (2.15): 71 _ 71 _ 1 _ FR ^~2^fJa~2fJa=Af' Избирательность фильтра, т. е. ширину его переходной полосы, определяет величина «окна» (для прямоугольного — 4л/L). Максимальный и интегральный уровни боковых лепестков частотной характеристики окна определяют неравномерность АЧХ фильтра в полосе пропускания и уровень гарантированного затухания в полосе подавления. Откуда следует, что спектральная плотность опти¬ мального «окна» Щв,п) должна обладать: — минимальной шириной главного лепестка, содержащего большую часть общей энергии (для обеспечения минимальной переходной полосы); — минимальным уровнем первого бокового лепестка и мини¬ мальной площадью под боковыми лепестками (для обеспечения минимума пульсаций). К сожалению, эти требования несовместимы. «Взвешивание» с подходящей весовой функцией обеспечивает сглаживание выбро¬ сов первоначальной частотной характеристики, т. е. подавление отклонений в полосе пропускания и уровня боковых лепестков в полосе задерживания, и, таким образом, достижения требуемого затухания в полосе задерживания ФНЧ. Платой за это является более пологий склон АЧХ, т. е. расширение переходной полосы. Поэтому отыскание подходящих весовых последовательностей («окон»), основывающееся на классических работах по ускорению сходимости рядов Фурье, является искусством нахождения разум¬ ного компромисса между этими требованиями. Таким образом, метод «окна» (или «взвешивания») заклю¬ чается в модификации коэффициентов фильтра (отсчетов беско¬ нечной импульсной характеристики h(n)), полученных в соответ¬ ствии с обратным дискретным по времени преобразованием Фурье 83
от заданной комплексной частотной характеристики, для получе¬ ния требуемой импульсной характеристики конечной длительно¬ сти Іік(ѵ) следующим образом: h¥(n) = h(n)w(n), (2.55) где w(n) — конечная весовая последовательность, называемая «окном», причем: , ч Л L-1 L-1 w(ti) = 0 при п < —, п > . Для выполнения условий физической реализуемости полу¬ ченная ИХ сдвигается вправо на половину ее длительности, т. е. на (L - 1)/2. Таким образом, процедура расчета дискретного фильтра с КИХ по методу взвешивания состоит из следующих этапов: 1. Задается требуемая «идеальная» комплексная частотная характеристика H(jQ). 2. Находится соответствующая импульсная характеристика «идеального» фильтра И(п) путем вычисления обратного дискрет¬ ного во времени преобразования Фурье функции H(jCl). 3. Подбирается подходящая функция окна w(n) для модифика¬ ции последовательности h(n) и получения импульсной характерис¬ тики конечной длительности hK(n), обеспечивающей требуемый вид КЧХ фильтра. В настоящее время известно несколько десятков функций окон¬ ного взвешивания. С наиболее полным перечнем и описанием их свойств можно познакомиться в монографиях [1, 12]. В табл. 2.1 приводятся некоторые наиболее часто используемые оконные функции. В качестве примера приведем описание и параметры обоб¬ щенной весовой функции Хэмминга, представляющей собой сумму прямоугольника и одного периода косинусоидальной функции: ѵгн (и) = а + (1 — cx)cos(27t———). (2.56) L 1 84
При а = 0,5 это выражение описывает весовую функцию Ханна (окно Хэннинга), а при а = 0,54 — весовую функцию Хэмминга. Комплексную частотную характеристику обобщенного окна Хэмминга можно представить в виде суммы трех КЧХ прямо- угольных окон с центральными частотами Ц, = 0 и Q0 = ± ——- : Wu(jQ) = aWu(jn) + ^(l-a)WM^±j^)\ (257) Боковые лепестки первого слагаемого находятся «в противо¬ фазе» с главным и боковыми лепестками двух последних сла¬ гаемых. За счет этого существенно снижается уровень боковых лепестков частотной характеристики оконной функции. Для окна Хэмминга уровень первого бокового лепестка на 43 дБ ниже глав¬ ного. Однако при этом расширяется главный лепесток вдвое. Это соответствует расширению переходной полосы между полосами пропускания и подавления проектируемого филь¬ тра, тогда как уменьшение уровня боковых лепестков частот¬ ной характеристики окна соответствует уменьшению пульсаций в полосе пропускания и лучшему подавлению в полосе задержи¬ вания. Как отмечалось выше, оптимальная оконная функция, пред¬ ставляющая собой последовательность конечной длины, должна иметь минимум энергии спектра за пределами некоторой задан¬ ной частоты. Одним из наилучших приближений к оптимальному окну является окно Кайзера, основанное на относительно простой аппроксимации так называемых вытянутых сфероидальных вол¬ новых функций: wK(n) = I0{fiyj\-[2n/(N-1)2]}//0(P), где /0(P) — модифицированная функция Бесселя первого рода нулевого порядка; Р — параметр, определяющий вид окна и уро¬ вень пульсаций. 85
Таблица 2.1 Основные весовые функции Временное окно Весовая функция w(n) Ширина главного лепестка W(Q) Естественное (прямоугольное) 1 2n/L Бартлетта (треугольное) 1—2|n|/(L - 1) 4n/L Вельша (параболическое) 1-4|n2|/(L - 1)2 4n/L Хана 0,5{1 + cos[2nn/(L - 1)]} 4n/L Хемминга 0,54 + 0,46cos[2nn/(L - 1)] 4n/L Блекмана 0,42 + 0,5cos[2nn/(X - 1) + + 0,08cos[4nn/(X - 1)] 6n/L Примечание: все весовые функции задаются на интервале -(L - 1)/2 < n < (L -1)/2; для прочих n значение весовой функции равно 0. Другим приближением к оптимальной форме окна является окно Дольф-Чебышева, обеспечивающее минимальную ширину главного лепестка кЧХ окна при фиксированном уровне боковых лепестков. При этом все боковые лепестки КЧХ этого окна имеют одинаковый уровень. 2.6. Синтез рекурсивных фильтров по аналоговому прототипу Задача синтеза цифрового фильтра с БИХ заключается в отыс¬ кании реализуемой передаточной функции вида H (z) = N 2 bz i = 0 M 1 + 2 amz- m = 1 (2.58) удовлетворяющей заданным требованиям. m 86
Первую группу методов расчета ЦФ с БИХ образуют прямые методы расчета в z-плоскости. С их помощью часто удается найти такое расположение полюсов и нулей фильтра, при котором обес¬ печивается некоторая аппроксимация непосредственно заданной характеристики ЦФ Однако вместо того, чтобы заново создавать теорию расчета ЦФ, можно использовать простые методы отображения, позволяющие преобразовать характеристики аналоговых фильтров из одной ком¬ плексной области (плоскость преобразований Лапласа) в другую (z-плоскость). Такие методы расчета ЦФ получили название «син¬ тез по аналоговому прототипу» и наиболее широко используются. Заметим, что при замене z нар выражение (2.58) представляет собой передаточную функцию аналогового фильтра. Проектирование ЦФ по аналоговому прототипу содержит два этапа: 1. Получение подходящей передаточной функции Н(р) анало¬ гового фильтра (проектирование аналогового прототипа). 2. Создание процедур перехода от Н(р) к H(z) (переход от ана¬ логового фильтра к цифровому). Проектирование аналогового прототипа ЦФ заключается в выборе вида аппроксимации АЧХ, расчете порядка фильтра, выборе структур¬ ной схемы и определении коэффициентов его передаточной функции. Проектирование завершается построением амплитудно-частотной и фазочастотной характеристик фильтра. Первая задача построения фильтра — аппроксимация идеаль¬ ной прямоугольной характеристики передаточной функцией, удов¬ летворяющей условиям физической реализуемости. Эта задача имеет многочисленные решения, доведенные для аналоговых фильтров до ряда стандартных таблиц и графиков. В зависимости от того, допускаются или нет пульсации АЧХ в полосе пропускания и области подавления, используются следующие типы фильтров, отличающиеся видом аппроксимирующей функции: 1. Фильтр Баттерворта, имеющий максимально плоскую АЧХ в полосе пропускания и монотонно возрастающее затухание в полосе задерживания (рис. 2.19, а). 87
2. Фильтр Чебышева с равноволновой АЧХ в полосе пропуска¬ ния и монотонно возрастающим затуханием в полосе подавления (рис. 2.19, б). 3. Инверсный фильтр Чебышева с монотонно возрастающим в полосе пропускания затуханием и равноволновой АЧХ в полосе подавления (рис. 2.19, в). 4. Эллиптический фильтр (фильтр Золотарева — Кауэра) с рав¬ новолновой АЧХ как в полосе пропускания, так и в полосе подав¬ ления (рис. 2.19, г). 5. Фильтр Бесселя (фильтр с максимально плоской характерис¬ тикой группового времени запаздывания) с аппроксимацией ФЧХ рядом Тейлора Рис. 2.19. Виды аппроксимации амплитудно-частотной характеристики (а-г) фильтра нижних частот Для использования на этапе расчета фильтра графиков и таблиц, помещенных в справочниках, либо стандартных программ рас¬ чета, т. е. для обращения к «каталогу фильтров», необходимо про¬ ектируемый фильтр привести к каноническому виду. Это приведе¬ ние осуществляется за счет двух процедур: преобразования частоты и нормирования частоты. Подробное описа¬ ние этих процедур приведено в монографии Г. Лэм [13]. Преобразование частоты представляет собой процедуру, с помо¬ щью которой требования к фильтру верхних частот (ФВЧ), полосо¬ вому фильтру (ПФ), заграждающему фильтру (ЗФ) преобразуются в требования к ФНЧ, называемому фильтром-прототипом. 88
Эта же процедура после расчета фильтра-прототипа дает простой способ перехода от ФНЧ к более сложным типам фильтров. При задании требований к ФНЧ фигурируют параметры АЧХ, приведенные ранее на рис. 2.12. АЧХ фильтра нижних частот (ФНЧ) Баттерворта описывается монотонно убывающей зависимостью вида Я(со) = 1 + ( \ со у «о у (2.59) и однозначно определяется двумя параметрами: характерной частотой со0 и порядком фильтра N [4]. Из анализа разложения квадрата этой функции в ряд Макло- рена следует, что первые (2N - 1) производные ее равны нулю вблизи со = 0. По этой причине фильтры Баттерворта также называ¬ ются фильтрами с максимально плоскими (гладкими) АЧХ. Частота со0 и порядок фильтра находятся из решения системы двух уравнений. Они составляются для заданных значений параметров, определяющих требования к уровню передачи Я. на частоте среза сос и гарантированному уровню подавления Я на границе полосы подавления сог Система этих уравнений согласно выражению (2.59) имеет вид: ЯГ2 = 1 + я;2 = 1 + / \2ІѴ у® о j ( Л1 со (2.60) у ®о j Решение этих уравнений относительно со0 и N дает: ®0=®с(Яс“2-1) lg (н;г-1 N =- кНс~ (2.61) 21g ( \ со 89
Из выражения (2.59) следует, что на характерной частоте со0 АЧХ фильтров Баттерворта любого порядка пересекаются на уровне 1Л/2 (т. е. на уровне -3 дБ). Поэтому чаще всего в качестве харак¬ терной частоты выбирают частоту среза. Значение А, определяемое выражением (2.61), округляется до ближайшего большего целого числа. В фильтрах Чебышева отличия аппроксимации от идеально прямоугольной АЧХ представляются равновеликими пульсациями. В зависимости от того, где допускаются эти пульсации — в полосе пропускания или в полосе затухания, — различают фильтры Чебы¬ шева I и II типов (рис. 2.20). Рис. 2.20. Амплитудно-частотные характеристики прямого и инверсного фильтров Чебышева АЧХ прямого фильтра Чебышева (I типа) описывается выра¬ жением #(©) = (2.62) где Гд^О) — полином Чебышева А-го порядка от аргумента Q = со/сос: Г cos(/V arccosQ), Q < 1; [ ch(yV arc hQ), £2 > 1. (2.63) 90
Параметр s в соотношении (2.62) характеризует неравномер¬ ность АЧХ в полосе пропускания. Порядок фильтра определяется из соотношения (2.62) при со = со,, тогда 1 1 + е2Т^) = Н2 С0„ откуда с использованием формулы (2.63) получим: = ch 7Varch(—) со„ Тогда соотношение для определения порядка N запишется так: N = - arch( -) arch(co3 /coc) (2.64) Используя свойства обратной функции от гиперболического косинуса, это выражение можно заменить на более удобное: N> lg(<7 + V<72 -1) і§(о3+Л/о^Т)' Здесь использованы обозначения: 4 = 10OJ4 -1 (2.65) (2.66) ! ю ' -1 и 03= со, /со0. При расчетах в выражении (2.66) абсолютные значения зату¬ хания в полосе пропускания Ас и полосе подавления А3 берутся в децибелах. АЧХ аналогового нормированного ФНЧ Чебышева II типа (инверсного) описывается следующим выражением: 1 Я(<в) = \А + s2T^(co3 /сос)/7^(со3 /со) (2.67) 91
Для определения порядка фильтра Чебышева II типа (инверс¬ ного) используются выражения, аналогичные выражениям (2.65) и (2.66), которые использовались для фильтра I типа. Фильтр Кауэра обладает АЧХ, отличительной особенностью которой является наличие пульсаций как в полосе пропускания, так и в области затухания. Выражение для АЧХ фильтра Кауэра имеет следующий вид: Н( со) 1 (2.68) где і?Л(Ос, L) — эллиптическая функция Якоби; L — параметр, характеризующий пульсации функции ЛЛ(Ос, L): , в Ѵяс-2-1 Ѵяз2-і Vя;2-1' Присутствие функции RNв формуле (2.68) определило и другое название фильтров этого типа— эллиптические фильтры. Порядок фильтра Кауэра определяется по формуле N = x(Q;1)x(>/rZ) к(Ь)к(^ 1-аг2)’ (2.69) где К— символ полного эллиптического интеграла первого рода. Более подробные сведения об эллиптических фильтрах можно найти в учебном пособии [6]. В пакете MATLAB (приложение Signal Processing) имеются функции выбора порядка фильтров buttord, cheblord, cheb2ord, ellipord, позволяющие рассчитывать на основе соотношений (2.61)-(2.69) минимально необходимый порядок как аналоговых, так и дискретных фильтров. После расчета порядка фильтра определение коэффициентов передаточной функции аналогового фильтра-прототипа осущест¬ вляется либо по справочникам, либо с использованием стандарт¬ ных функций пакета MATLAB 92
При синтезе цифровых фильтров по аналоговому прототипу главным вопросом является создание процедуры отображения ком¬ плексной /^-плоскости на комплексную г-плоскость. Или, иными словами, процедуры перехода от передаточной функции аналого¬ вого фильтра-прототипа На(р) к передаточной функции цифрового фильтра H(z). Передаточная функция аналогового фильтра HJp>) представ¬ ляет собой дробно-рациональную функцию комплексной перемен¬ ной р. Чтобы получить передаточную функцию дискретного фильтра H(z), необходимо перейти из комплексной /7-области в комплексную г-область, причем дробно-рациональный характер функции должен сохраниться. Эта процедура должна удовлетворять двум условиям: 1. Мнимая ось /7-плоскости должна отображаться в единичную окружность Г-ПЛОСКОСТИ 1/7 =/С0, -00 < СО < ос] —> [- = еі<лТ>\ -п < 0)ТД < 7і]. Это требование обеспечивает сохранение вида частотной характеристики аналогового фильтра (рис. 2.21, а). 2. Левая половина /7-плоскости Re|/7] < 0 отображается в часть z-плоскости внутри единичного круга [р\ Re[/7] < 0 } —» (z| |z| < 1}. Это условие необходимо для сохранения свойств устойчиво¬ сти, т. е. чтобы процедура перехода переводила устойчивый анало¬ говый фильтр в устойчивый ЦФ (рис. 2.21, б). Рис. 2.21. Требования к процедуре перехода (пояснения в тексте) 93
Процедуры перехода строятся с использованием метода инва¬ риантности импульсной характеристики, метода билинейного пре¬ образования, метода согласованного (прямого) Z-преобразования и ряда других. 2.7. Метод билинейного Z-преобразования При использовании метода инвариантности импульсной характеристики импульсная характе¬ ристика h(n) рекурсивного цифрового фильтра (ЦФ) получается путем дискретизации импульсной характеристики аналогового фильтра-прототипа ha(t) с периодом Тп: И{П) = ИМ.-пт =НпТаУ (2.70) В этом смысле сохраняется инвариантность, т. е. неизменность импульсной характеристики. Затем находят передаточную функцию ЦФ как Z-преобра- зование импульсной характеристики: tf(z) = Z{A>7;)}; hy \ ѵ1 и ( \ -я (2.71) H{z) = 2^K(n)z . п = 0 Поскольку в данном случае импульсная характеристика — дис¬ кретный аналог ha(t), то связь между комплексными частотными характеристиками (КЧХ) аналогового и цифрового фильтров имеет такой же вид, как связь между спектром аналогового и дис¬ кретного сигнала — периодическое повторение с частотой дискре¬ тизации: Й(е”г-)=^Х НМѵ+ЩкЯ -*Д *=-°° -*Д (2.72) То есть при использовании этого метода неизбежно возни¬ кает эффект наложения (рис. 2.22). Для того, чтобы этот эффект был мал, значение АЧХ аналогового фильтра-прототипа должно 94
быть пренебрежимо малым на частотах, превышающих половину частоты дискретизации (частота Найквиста). Рис. 2.22. Трансформация КЧХ при переходе от аналогового (а) к цифровому (б) фильтру методом инвариантности импульсной характеристики Проиллюстрируем переход на следующем примере. Переда¬ точная функция аналогового фильтра-прототипа представляет собой дробно-рациональную функцию вида N _ *.(р) = -4г—. \ + ІаіРг і=і Разложим ее на простые дроби. В частности, когда все полюсы Д — различные, это разложение принимает следующий вид: М л вд=Е^г- tip-p, 95
Тогда импульсная характеристика аналогового фильтра пред¬ ставляет собой сумму: и 4(0 = 14^1(0- і = 1 Импульсная характеристика ЦФ получается путем дискрети¬ зации: м h(n) = \(nTR)=YA,e»nTm. 1 = 1 Передаточная функция ЦФ — Z-преобразование от импульс¬ ной характеристики: H(z) = Z{h(n)}- оо H(z)=YKn)z-n и = 0 со М М со „ Z Z =Е4 Z (ей п=01 = 1 і = 1 и = О М Z l-z-V*' Таким образом, соответствие между передаточными функциями аналогового фильтра-прототипа и ЦФ следует из сопоставления выражений _4 , 4 4 />-/>, 1-гѴЛ l-z^-1’ где zi = еА7д — полюс ЦФ, соответствующий полюсу Д аналогового фильтра. Переход от /7-плоскости к z-плоскости при этом описывается базовым соотношением z = е?т». (2.73) Источник эффекта наложения в том, что этот переход не обес¬ печивает однозначного отображения мнимой оси /7-плоскости на окружность единичного радиуса в z-плоскости. Так, точки/7 = 0, и=/ —, p = -j — отображаются в точку z = 1, а полоса шириной 2л; — отображается на всю единичную окружность в z-плоскости 4 при формировании передаточной функции цифрового фильтра (рис. 2.23): 96
п . тс . / < р < —i, eT,It < z < е±,%. Рис. 2.23. Преобразование КЧХ при дискретизации импульсной характеристики Из-за эффекта наложения метод инвариантности импульсной характеристики применим только для фильтров с существенно ограниченной частотной характеристикой фильтра-прототипа: |Щ/св)| ~ 0, |со| > Асо, т. е. для фильтров нижних частот или полосовых фильтров. Одним из преобразований, обеспечивающих однозначное ото¬ бражение мнимой оси /^-плоскости (ось частот) на единичную окружность z-плоскости, является билинейное преобра¬ зование, которое определяется следующим образом: 1_ 1-z-1 _ 2 z-l Уж'\^~Уж'7Ѵ\ (2.74) Этот метод обеспечивает такое построение ЦФ, при кото¬ ром приближенно получают соответствие реакции этого фильтра и аналогового фильтра для любых воздействий и сохраняется общая форма АЧХ. 97
Выражение (2.74) может быть получено из преобразованного базового соотношения (2.73): P=(\rQ\nz. (2.75) Представим логарифмическую функцию в последнем выраже¬ нии рядом: In z = 2 z — 1 lf z —О 1 z + 1 3 z + 1 и + - 5 z-1 z + 1 Взяв из этого ряда только первый член разложения, после его подстановки в выражение (2.75) получим (2.74). Для прояснения физической сути основного соотношения билинейного преобразования (2.74), описывающего процедуру перехода от Нл{р) к //(z), рассмотрим связь между передаточ¬ ными функциями идеального аналогового интегратора и циф¬ рового. Идеальный аналоговый интегратор описывается передаточной функцией #» = -. (2.76) Р Его импульсная характеристика: h(t) = 1 О t > 0; t < 0. Реакция такого интегратора на произвольное воздействие x(f): t y(t) = Jx(x)h(t - x)dx. 0 Если 0 < < t2, то h к y(t2) -y(0 = ]x(x)h(t2 -x)dx-jx(x)h(t, -x)dx. о 0 При 0 < t < tx, t2 h{t2 - x) = h(tx - t) = 1. h y(t2) ~y(0 = \x(x)dx. к Приближенное вычисление этого интеграла методом трапеций проиллюстрировано на рис. 2.24. 98
Рис. 2.24. Приближенное вычисление интеграла При tx —► t2 y(h) -y(h)+ x(t2)\ Если tx-t2 = Гд, то получаем уравнение цифрового интегратора: у(пТл) -у[(п-\)Тд] = ^-[х(п) + х(п- (2.77) Применив Z-преобразование к обеим частям, получим: y(z) -f(Z)Z-i=E[Je(Z)+x(Z)z-1]. Откуда следует, что передаточная функция цифрового интегра¬ тора имеет вид: я(г) = ЕЕ=Е.£±і X(z) 2 z — \ (2.78) Сравнивая H(z) с Яа(р), получим правило замены: 2 z-1 Поскольку аналоговый фильтр представляет собой совокуп¬ ность аналоговых сумматоров, умножителей и интеграторов, то, заменяя каждый аналоговый элемент соответствующим цифро¬ вым, получим ЦФ. Таким образом, 2 z — 1 Щг)=Н,(р) при Р ■ 99
Обратное соотношение 2 + РТа z = - 2~рт* (2.79) показывает, что мнимая ось р = /со отображается на единую окруж¬ ность в z-плоскости однозначно. Откуда следует, что эффект нало¬ жения отсутствует. Однако, поскольку билинейное преобразование — нелинейная функция, то частотные характеристики аналогового и дискретного фильтров связаны друг с другом трансформацией (нелинейной деформацией) частотной оси. Поясним эффект деформации шкалы частот. Пусть со и соц — частотные переменные аналогового и цифрового фильтров соот¬ ветственно. Тогда из (2.74) для точек на оси частот получаем: 2 е>л-1 '/Ю~ 7^Ѵ'Л +1 и “ = ftg(coqT;/2)- 1 л (2.80) На низких частотах, когда соТд«1 (соц < 0,3/Тд), тангенс примерно равен своему аргументу: )«со. Поэтому в области низких частот частотные характеристики ана¬ логового и дискретного фильтров почти совпадают. Далее, по мере ускорения роста функции тангенса, частотная характеристика дис¬ кретного фильтра все сильнее сжимается по сравнению с аналого¬ вым прототипом и на частоте, равной половине частоты дискре¬ тизации, достигает значения, которое частотная характеристика аналогового фильтра имела бы на бесконечной частоте (рис. 2.25). Для получения дискретного фильтра с заданными частотами среза необходимо скорректировать частоты среза аналогового прототипа, чтобы компенсировать искажения частотной оси. 100
Q Рис. 2.25. Деформация шкалы частот (о-в) при билинейном преобразовании Билинейное преобразование обеспечивает простую процедуру перехода от аналогового фильтра к ЦФ и сохраняет вид ЧХ при пре¬ образовании. Широкополосные аналоговые фильтры преобразуются в широкополосные ЦФ без эффекта наложения. Однако нелиней¬ ность соотношения между соа и соц приводит к искажению частот¬ ных характеристик ЦФ по отношению к аналоговому фильтру- прототипу. В частности, этот факт не позволяет синтезировать этим методом ЦФ с линейной фазочастотной характеристикой. При этом преобразовании также не сохраняется точный вид импульсной характеристики. КОНТРОЛЬНЫЕ ВОПРОСЫ 1. Что называется передаточной функцией цифрового фильтра? 2. Что понимают под определением «физически реализуемый циф¬ ровой фильтр»? 3. Как связаны между собой комплексная частотная и импульсная характеристики дискретного фильтра? 4. Каков период комплексной частотной характеристики дискретного фильтра? 101
5. Фильтр какого типа и порядка описывает передаточная функция вида Й(г) = с0 + сх2гх + с2г2 + съгг. 6. Какая форма реализации цифрового фильтра называется канони¬ ческой? 7. Каков порядок и форма реализации фильтра, описываемого пере¬ даточной функцией вида 1 +z~' -0,5z-2 H(z) = : : . (l-0,2z^)(l+0,2z_1-0,2 z-) 8. Нерекурсивный фильтр, представляющий собой равновесный сум¬ матор 4-х отсчетов, предназначен для обработки дискретного сигнала, задаваемого 16 отсчетами. Какова должна быть размерность (количество точек) дискретного преобразования Фурье при реализации обработки в таком фильтре в частотной области? 9. С какой дискретизацией по частоте анализирует спектр цифровой спектроанализатор, использующий алгоритм 256-точечного БПФ, если частота дискретизации АЦП равна 512 кГц? 10. Какова разрешающая способность по частоте (по уровню -3 дБ) цифрового спектроанализатора, использующего алгоритм 256-точечного БПФ, если период дискретизации равен 0,1 мкс9 11. В чем заключается синтез нерекурсивных ЦФ методом «окна»? 12. Какой тип фильтра (при одинаковом порядке) обладает макси¬ мальной избирательностью9 13. Какому типу фильтра соответствует комплексная частотная харак¬ теристика: H(j(a) = ——. а + 7Ш 14. Чем объясняется деформация шкалы частот при синтезе рекур¬ сивных цифровых фильтров методом билинейного преобразования?
3. ЭФФЕКТЫ КОНЕЧНОЙ РАЗРЯДНОСТИ ПРИ ПРЕДСТАВЛЕНИИ ЧИСЕЛ В ЦИФРОВЫХ ФИЛЬТРАХ ЗЛ. Шум квантования Ранее при рассмотрении цифровых фильтров предполага¬ лось, что переменные и коэффициенты фильтров представляются с неограниченной точностью. То есть, по сути дела, рассматрива¬ лись дискретные, а не цифровые фильтры. Теперь пришло время оценить влияние эффектов, вызванных конечной разрядностью при описании цифровых сигналов и реализации процедур обра¬ ботки в цифровых фильтрах. К этим эффектам, в первую очередь, относятся: 1. Ошибки за счет квантования входного сигнала (шум кванто¬ вания сигнала). 2. Погрешности за счет округления промежуточных результа¬ тов выполнения арифметических операций (шум округления). 3. Погрешности характеристик фильтров, обусловленные кван¬ тованием коэффициентов. Результаты проявления указанных эффектов зависят от следу¬ ющих факторов: — формы представления и кодирования чисел (с фиксирован¬ ной или плавающей запятой; прямой, обратный или дополнитель¬ ный код); — способа квантования; — структуры фильтра. Поскольку число возможных комбинаций этих условий весьма велико, сквозной анализ проводится только для каждого из эффек¬ тов в отдельности, а полный анализ осуществляют методом моде¬ лирования. Прежде всего рассмотрим влияние квантования сигнала, которое осуществляется в аналого-цифровом преобразователе 103
и представляет собой безынерционное нелинейное преобразова¬ ние, описываемое ступенчатой функцией F(x) — характеристи¬ кой квантователя. При квантовании обычно используются два спо¬ соба: округление и усечение. Квантование приводит к ошибке: екв = *кв -x = F(x)-x.(3.1) Характеристика равношагового квантователя при использова¬ нии округления приведена на рис. 3.1. Шаг квантования Q в этом случае равен весовому коэффициенту младшего числового раз¬ ряда. При кодировании двоичным кодом он равен 2~ъ, где Ъ — число разрядов. Квантование с округлением соответствует выбору ближайшего уровня квантования. Уровень Рис. 3.1. Амплитудная характеристика квантователя при квантовании с округлением За исключением случаев превышения предельных значений Jtminи *тах ошибка округления е(п) лежит в пределах: -QI2 <e<Q/2 104
или 2~b 2~b <e< . 2 2 (3.2) При достаточно общих предположениях можно считать, что распределение ошибки является равномерным (рис. 3.2). w(e) Q Q 2 2 Рис. 3.2. Распределение ошибки квантования при округлении При усечении в качестве квантованного отсчета используется ближайший меньший уровень квантования (рис. 3.3). Для цифро¬ вого сигнала это эквивалентно отбрасыванию младших разрядов. Поскольку результат усечения равен результату округления, уменьшенному на половину шага квантования, то график плот¬ ности вероятности сигнала имеет вид, изображенный на рис. 3.4. Среднее значение ошибки (математическое ожидание) равно Q/2. Дисперсия равна Q2/12. Рис. 3.3. Характеристика квантователя при использовании усечения 105
w(e) 2-ъ < e < О Q Q e Рис. 3.4. Распределение ошибки квантования при усечении Процесс квантования можно рассматривать как наложение на сигнал, заданный точно в каждый дискретный момент времени, шума квантования е(п): Модель шумов квантования основана на следующих предпо¬ ложениях: — е(п) — выборка из стационарного случайного процесса; — е{п) — не коррелирован с х — е(п) — отсчеты не коррелированы между собой (дискрет¬ ный белый шум); — распределение ошибки равномерное во всем диапазоне зна¬ чений ошибок квантования. Есть случаи, когда эти предположения неверны (постоянный сигнал, синусоида с частотой, кратной FJ. Модель можно использовать, когда поведение сигнала таково, что при переходе от одной выборки к другой его значение меняется на несколько уровней квантования и в то же время число уровней (разрядов) не очень мало. Вклад этой ошибки можно условно измерять в виде отношения мощности сигнала к мощности шума квантования (С/Ш): xjn)=x(n) + e(n). (3.3) (3.4) 12 106
Если 0 = 2 b С/Ш = 12 • 22Ч2. Используя логарифмическую меру, получим: ^(дБ) = 101gl2 + 6201g2 + 101ga3r, так как 201g2 ~ 6; ^г(дБ) = 10,8 + lOlgc^ + 6Ь. (3.5) Таким образом, добавление одного разряда увеличивает отно¬ шение сигнал/шум квантования на 6 дБ. В общем случае при выборе числа уровней квантования глав¬ ным критерием является точность воспроизведения выходного сигнала. Использование понятия «шум квантования» позволяет ввести линейную модель учета влияния процесса квантования на выход¬ ной сигнал в цифровом фильтре (рис. 3.5). е(гі) Гг < h(ri) х(п) V J xjn) У(п) Рис. 3.5. Линейная модель учета квантования сигнала в цифровом фильтре Модель линейна, когда шаг квантования достаточно мал. При этом различные источники ошибок в линейном ЦФ можно изучать отдельно. В предположении, что обработка выполняется с абсо¬ лютной точностью, т. е. единственной ошибкой является ошибка квантования входного сигнала, выходной сигнал определится как свертка квантованного сигнала и импульсной характеристики фильтра следующим образом: у(п)= Цк(к)х..(п-к) = Y. h{k)x(n - к) + £h(k)e(n-k). (3.6) к=0 к=0 к=О 107
Первое слагаемое в этом выражении представляет собой истин¬ ное значение, т. е. выходной сигнал в дискретном фильтре, а вто¬ рое — ошибку выходного сигнала, обусловленную квантованием входного сигнала в ЦФ (шумом квантования): y(n) =y(n)ucT + e(n); (3.7) в (л) = ^ h(k)e(n - к). (3.8) к= О Оценку величины этой ошибки обычно производят, используя два подхода: детерминированный и статистический. При детерминированной оценке ошибки е определяют ее мак¬ симально возможное значение. Если разрядность отсчетов вход¬ ного сигнала (после запятой) равна Ьвх, то максимальная ошибка квантования входного сигнала (при квантовании с округлением) равна: ^ Евх = тах|<?(и)| = 2^ = (3.9) При этом максимальная ошибка выходного сигнала ЦФ, обу¬ словленная квантованием входного сигнала, оценивается следую¬ щим образом: Евш = тах|б(и)| < тах|ф)| Z \h(k)\ = Цг Z \Щ)\ (ЗЛО) к= О L к= О При вероятностной оценке ошибки сигнала на выходе ЦФ вычисляют математическое ожидание (среднее значение) и дис¬ персию в. Среднее значение ошибки на выходе ЦФ: те(п) = м| Y, е(к)h(n - к) 1 = те X h(n - к). (3.11) Ц=о J к=о Если среднее значение ошибки квантования входного сигнала те = 0, то и тг = 0. Дисперсия ошибки на выходе за счет шума квантования: ст;(«) = м{е2(«)}-/яе2; М{е2(я)}=МІ ^h(m)e(n-m)- £ h(k)e(n-k)\. 1 ’ I w=0 k = 0 I 108
Поскольку в модели шумов квантования считается, что отсчеты е(п) не коррелированъ! (дискретный белый шум), то м{е2(и)}= £ X h(m)h(k)6(k- £ h2(m) + m2e. (3.12) Таким образом, в установившемся режиме дисперсия ошибки на выходе ЦФ с конечной импульсной характеристикой: В установившемся режиме для устойчивого БИХ-фильтра h{m) —» оо, если т —>• ос. Поэтому дисперсия ошибки на выходе ЦФ равна: Г)2 00 °1=—ЦЬ2(т). (3.14) 1Z т = О Определим выражение для дисперсии ошибки выходного сиг¬ нала, обусловленной квантованием сигнала на входе для рекурсив¬ ного фильтра первого порядка (рис. 3.6), алгоритм работы кото¬ рого описывается выражением т=0к=0 (3.13) у {п) = хігі) + ay in - 1). УІп) Рис. 3.6. Структурная схема рекурсивного цифрового фильтра первого порядка 109
Импульсная характеристика этого фильтра имеет вид: h(n) = ап при п > 0. Тогда дисперсия ошибки выходного сигнала: = ^|;л2(/і)=ст^-Т_. (3.15) п=0 I £/ Поскольку в устойчивом рекурсивном ЦФ первого порядка \а\ < 1, то <з\ ><з2е. То есть в ЦФ происходит усиление шумов квантования. Этот эффект усиливается по мере приближения параметра а к единице, или иными словами, при приближении координаты полюса передаточной функции к окружности единич¬ ного радиуса на комплексной z-плоскости. Проанализируем влияние шумов квантования в ЦФ второго порядка, схема которого изображена на рис. 3.7. Рис. 3.7. Структурная схема рекурсивного цифрового фильтра второго порядка Алгоритм работы такого фильтра описывается разностным уравнением: у(п) = х(п) - ахух(п - 1) - ау2(п - 2). ПО
Передаточная функция фильтра имеет вид: 1 H(z) = Y(z) 1 + axz 1 + a9z 2 Импульсную характеристику этого ЦФ найдем как обратное Z-преобразование передаточной функции, предварительно преоб¬ разовав ее, используя теорему разложения: tf(z) = 1 —a,z 1 1 — a2z 2 ’ H(z) = — z' +a,z + a2 H(z) = z A - z" + axz + a2 = z 4, , A У z — a, z — cl2j Полюсы передаточной функции a, и a, найдем из решения уравнения: z2 + ajZ + а2 = 0; ai,2 =-y±yJa2- V ^ / Коэффициенты разложения при этом равны: B\z) a, a, a, 2at + aj V / a, + a9 a. Поскольку aj — a2 = 2. 4 = l|V| Ы A(z) B'(z) -a9,T0 a-, =«2 a9 +aj Таким образом, выражение для передаточной функции приоб¬ ретает вид: H(z) = - a, a9 - a2 z - 0Ц - a2 z - a2 111
Из табл. 1.1 для Z-преобразований следует, что 1 ►а", п> 0. z-ax i _ z Таким образом, общее выражение для импульсной характери¬ стики рассматриваемого ЦФ имеет вид: h(n) = - OL 0Со а, -ои Рассмотрим случай комплексно-сопряженных корней: «u=-f±. г<02 V ^ У -а0. Их расположение на комплексной z-плоскости показано на рис. 3.8. Связь координат полюсов с коэффициентами ЦФ задается соотношениями: г - yfa^, cosB = а\2= rQ±je- Для этого случая импульсную характеристику ЦФ можно пред¬ ставить в компактной форме: <+1 -а"+1 _ rn+l(QKn+l)Q -Q-j(n+l)e) _ _n sin2 [(^2 + 1)0] ^ ^ - / _ /Ѳ -7ЙЧ r -2/Ч- h(n)■ оц - a r(e7U-e-j0) sin 0 112 Рис. 3.8. Расположение полюсов передаточной функции цифрового фильтра на комплексной плоскости
Зависимость отсчетов импульсной характеристики от угла расположения комплексно-сопряженных полюсов передаточной функции иллюстрирует рис. 3.9. h Рис. 3.9. Зависимость отсчетов импульсной характеристики от углового расположения полюсов передаточной функции цифрового фильтра на комплексной плоскости Найдем дисперсию шума на выходе, обусловленную квантова¬ нием входного сигнала: 12 и=о (3.17) и=о iz- и=о sinG Интенсивность шумов растет по мере приближения Ѳ —> О и Ѳ —> я. В этих случаях выражение (3.17) приобретает простой вид: ‘ , ais» 12 „Го smѲ 12 1-г Из анализа этого выражения следует, что интенсивность шумов на выходе растет по мере увеличения г, т. е. с приближением полю¬ сов к единичной окружности. Дисперсия ошибки выходного сигнала за счет шума кванто¬ вания также может быть вычислена через амплитудно-частотную характеристику ЦФ Щсо). Согласно равенству Парсеваля: УА!=^7ие'"г')ГА». <зл9> ш=0 71 о 113
Тогда Г)- Т я/гд Ов=тг-— J \H(j<o)\-da>. (3.20) Iz 71 о Таким образом, по допустимой величине ошибки на выходе о: и известной АЧХ или импульсной характеристике ЦФ можно опре¬ делить допустимую величину дисперсии ошибки входного сиг¬ нала, а тем самым и шаг квантования (разрядность) АЦП. 3.2. Эффекты округления промежуточных результатов Процедура линейной цифровой фильтрации состоит из опера¬ ций сложения, умножения на постоянные числа (коэффициенты фильтра) и сдвига (запоминания). Погрешности за счет округления промежуточных результатов выполнения арифметических опера¬ ций вызываются конечной разрядностью используемых в цифро¬ вом фильтре регистров. При реализации вычислений с фиксированной запятой сло¬ жение двух чисел с разрядностью b и разрядности сумматора не менее b не приводит к ошибкам округления. Возможно лишь переполнение регистра сумматора, для исключения которого вво¬ дится масштабирование. Произведение двух чисел, представленных в формате с фиксиро¬ ванной запятой Ь1 и Ь2 разрядами, может содержать Ьх + Ь2 разрядов. Как правило, его надо разместить в регистр, содержащий b < ъх + ъ2 разрядов (иначе, в рекурсивном фильтре длина регист¬ ров будет увеличиваться до бесконечности). Такое преобразование после каждого умножения равносильно квантованию промежуточ¬ ного результата с округлением. Возникает ошибка округления. Эту ошибку можно оценить, если использовать рассмотренную выше статистическую модель шумов квантования. То есть процедуру округления можно предста¬ вить как наложение на точный результат умножения шума округле¬ ния. Модель шумов округления базируется на тех же допущениях, 114
что модель шумов квантования сигнала. То есть шум округления еок(п) представляет собой дискретный стационарный случайный процесс, значения которого не коррелированы с сигналом и между собой. Распределение значений шума округления можно считать равномерным. Модель можно применять, если число разрядов Ъх и Ъ2 не очень мало. Максимальное значение шума округления: таХКк(>)|=2б=22"А> (3.21) а дисперсия: (3.22) Поскольку b > Ьѵ то ошибка округления может быть меньше ошибки квантования входного сигнала. В этом случае умножитель с ограниченным числом разрядов может быть представлен линейной моделью в виде идеального умножителя и сумматора, на второй вход которого поступает шум округления (рис. 3.10). Когда такое представление можно считать верным, диспер¬ сию шума в выходном сигнале, обусловленную і-м умножителем, можно вычислить с помощью метода, аналогичного анализу влия¬ ния шумов квантования. Ошибка выходного сигнала, обусловлен¬ ная шумом округления после /-го умножителя, определяется сверт¬ кой шума округления с импульсной характеристикой части ЦФ от выхода /-го умножителя до выхода фильтра //,(£): S0K, 00 = I h, ( (п - (3.23) к = 0 Рис. ЗЛО. Линейная модель умножителя с ограниченным числом разрядов выходного регистра 115
Шум округления, обусловленный всеми L источниками шума: £окОО = ІХк <(«)• (3.24) і = 0 Детерминированная оценка выходного шума — максимальная ошибка, обусловленная і-м умножителем: 00 О) max|E0KI-(w)| ^ max|e0Kj(«)| £ К(*)| ^ ^ X К(^)|- (3-25) я п к=о z А-=0 Оценка максимального значения выходного шума от всех L умножителей (при одинаковой разрядности): £окш„ £|л,■(*)!• (3.26) Z і=1А- = 0 Дисперсия шума на выходе ЦФ, обусловленная і-м умножите¬ лем в установившемся режиме: <4, =-у^1>,2(£). 11 к = О Дисперсия полной ошибки выходного сигнала: L П2 L л <4=Е<4 ,=7tZI№)- і=1 iz J=1 Л- = 0 (3.27) (3.28) Влияние шума округления промежуточных результатов на выходной сигнал, в отличие от шума квантования, зависит от структуры фильтра. В КИХ-фильтре шум квантования произведения просто добав¬ ляется к выходному сигналу, поэтому ст^к = La2OKj и зависит только от порядка фильтра. В БИХ-фильтре не все источники шума влияют одинаково. В качестве примера рассмотрим фильтр второго порядка, реализо¬ ванный по канонической схеме (рис. 3.11). Источники шумов ех{п), е2(п), е3(п) непосредственно добавляют ошибку в выходной сигнал фильтра, тогда как источники е4(п), е5(п) введены в цепь обратной связи, поэтому существует возмож¬ ность их усиления полюсами фильтра. Полюс дает резонанс в АЧХ фильтра, что может привести к существенному усилению шума квантования 116
Найдем дисперсию выходного шума, обусловленного округле¬ нием произведений в устойчивом звене первого порядка (рис. 3.12), описываемом уравнением y(n) = x(n) - ay(n - 1), \a\ < 1. Рис. 3.12. Модель цифрового фильтра первого порядка с учетом шума округления произведения 117
Ошибка округления произведения проходит через ту же цепь, что и входной сигнал. Поэтому импульсная характеристика, соот¬ ветствующая точке подключения источника шума округления, совпадает с импульсной характеристикой фильтра: х(п) = ап. Дисперсия ошибки на выходе: вых = °ок XЬ2("> = тт- Ёа2п = ^ (3.29) Здесь g — шаг квантования произведения, весовой коэффици¬ ент младшего разряда регистра. Проанализируем теперь влияние шума квантования произведе¬ ний на выходе звена второго порядка, схема которого приведена на рис. 3.7. Влияние округления результатов умножения учитыва¬ ется введением двух шумовых последовательностей ех(п) и е2(п). Как видно из рис. 3.13, эти шумовые последовательности прохо¬ дят по той же цепи, что и входной сигнал. Импульсные характери¬ стики этих цепей совпадают с импульсной характеристикой всего фильтра, описываемой выражением (3.16). Поэтому выражение для дисперсии ошибки выходного сигнала, обусловленной округ¬ лением результатов умножения, имеет вид: Рис. 3.13. Модель цифрового фильтра второго порядка с учетом шума округления произведения 118
Общая ошибка квантования, обусловленная квантованием входного сигнала и результатов арифметических операций, опре¬ деляется суммой дисперсий соответствующих ошибок. При задан¬ ной допустимой ошибке выходного сигнала уже нетрудно опреде¬ лить требуемую разрядность кодов для АЦП и умножителей. 3.3. Анализ влияния квантования коэффициентов Требования к фильтру в частотной области задаются величи¬ ной неравномерности АЧХ в полосе пропускания и гарантиро¬ ванным уровнем затухания в полосе подавления. На первых эта¬ пах расчета ЦФ решается задача аппроксимации передаточной функции (и соответственно АЧХ) одной из функций (Баттерворта, Чебышева, эллиптической, Бесселя и др.). При равномерном критерии аппроксимации АЧХ фильтра Я(П) должна удовлетворять следующим требованиям (рис. 3.14): Здесь Я'(0) — идеальная АЧХ, соответствующая кусочно¬ линейной аппроксимации; e(Q) — допустимое отклонение от иде¬ альной АЧХ в полосе пропускания и полосе подавления; “ - ю - 2 л д — нормированная частота. В результате решения аппроксимационной задачи рассчитыва¬ ются коэффициенты а , bt или с„ описывающие АЧХ H(Q). Для дис¬ кретного фильтра они задаются с бесконечной точностью. Коэффициенты цифрового фильтра квантованы и представля¬ ются в системе счисления с конечной разрядностью, определяе¬ мой разрядностью памяти микропроцессора, на котором реализу¬ ется цифровой фильтр. |Я'(П)-Я(П)|<е(П); О < Q < 0,5, (3.31) где H' 8 119
Рис. 3.14. Аппроксимация амплитудно-частотной характеристики дискретного фильтра Процесс округления коэффициентов вносит погрешность, в результате которой нули и полюсы передаточной функции изме¬ няют свое расположение на комплексной плоскости, что приводит к искажению АЧХ. Вместо Н(С1) получаем 7^ (fi) — АЧХ фильтра, рассчитан¬ ную при округлении коэффициентов. Требования к АЧХ остаются прежними: |//'(Q)-//1(Q)|<8(Q); О < fi < 0,5. Деформация частотных характеристик при переходе к ЦФ не должна приводить к их выходу за пределы, установленные в задании на проектирование ЦФ. Это может быть в том случае, когда при расчете дискретного фильтра величины неравномерно¬ сти затухания в полосе пропускания и гарантированного затухания в полосе подавления взяты с некоторым запасом. Поэтому проце¬ дура проектирования ЦФ представляет собой итеративный про¬ цесс ввиду сложности и трудоемкости оценки влияния разрядно¬ сти на искажения комплексной частотной характеристики. 120
АЧХ фильтра обычно нормируется так, чтобы \Н(е£,т)\* 2 * * * * ~ 1 в полосе пропускания и 1-Ще7®1)!2 < 1 в полосе задерживания (рис. 3.15). і|Я(е^ЧІ2 Рис. 3.15. Нормирование амплитудно-частотной характеристики цифрового фильтра Поскольку по теореме Парсеваля: то при такой нормировке площадь под АЧХ в пределах половины частоты дискретизации меньше тт/У) и правая часть этого равенства меньше единицы. Для КИХ-фильтров коэффициенты фильтра — суть отсчеты импульсной характеристики: со 7Г 2 (3.33) h(n) = сп. Это означает, что 2;Ѵ(и)< 1 и іѴ|<1. Поэтому 121
Двоичный код коэффициентов КИХ-фильтра содержит лишь знаковый разряд и дробную часть и не содержит целой части. При этом АЧХ фильтра с КИХ в каждой частотной точке рас¬ считывается по формуле Затем осуществляется проверка на выполнение условий нера¬ венства (3.32). Если в каждой частотной точке Q неравенство выполняется, то это число разрядов считается допустимым. При проектировании ЦФ влияние квантования коэффи¬ циентов рассматривается для нескольких значений разрядно¬ сти. При этом квантованные значения коэффициентов bzq, azq получаются путем применения следующей процедуры в пакете MATLAB: bzq = round(bz • М)/М; azq = round(az • М)/М. Здесь М — число уровней квантования при заданной разряд¬ ности. Для рекурсивных фильтров степень искажения частот¬ ных характеристик зависит не только от его типа и порядка, но и от формы реализации (прямая, каскадная, параллельная). Исследования показали, что для фильтров с крутыми скатами АЧХ наименее чувствительной к округлению коэффициентов является каскадная форма реализации, использующая биквадрат¬ ные блоки. В качестве примера на рис. 3.16 приведены АЧХ ФВЧ Баттер- ворта с параметрами: /с = 100 кГц; /э = 50 кГц; А = 20 дБ. Пунк¬ тирной линией показана АЧХ дискретного фильтра, сплошной — цифрового при различных значениях разрядности коэффициентов. При прямой форме реализации удовлетворительные результаты получаются только при 14-разряд ном квантовании, в то время как каскадная форма обеспечивает хорошие результаты уже при 8 раз¬ рядах. 122
Рис. 3.16. Амплитудно-частотная характеристика ФВЧ Баттерворта: а — прямая форма реализации, 10 разрядов; 6 — прямая форма реализации, 14 разрядов; в — каскадная форма реализации, 8 разрядов Еще больше это различие проявляется при реализации высоко¬ добротных полосовых фильтров. На рис. 3.17 приведены АЧХ поло¬ сового эллиптического фильтра четвертого порядка с центральной частотой 2 МГц и полосой 10 кГц. Пунктиром показана АЧХ циф¬ рового фильтра с округлением коэффициентов. Прямая форма даже при 32 разрядах дает неприемлемый результат, в то время как при каскадной форме достаточным оказывается 8 разрядов. 123
X Ю6 Рис. 3.17. Амплитудно-частотная характеристика полосового эллиптического фильтра: а — прямая форма реализации, 32 разряда; 6 — каскадная форма реализации, 8 разрядов В заключение отметим, что в пакете расширения MATLAB Filter Design реализованы различные типы квантователей, а также 124
алгоритмы фильтрации, использующие представление чисел (в том числе и промежуточных результатов вычислений) в задан¬ ных форматах, что позволяет моделировать влияние квантования коэффициентов при проектировании самых различных ЦФ. КОНТРОЛЬНЫЕ ВОПРОСЫ 1. Опишите математическую модель шума квантования. В каких случаях процесс квантования входного сигнала не удовлетворяет этой модели? 2. На сколько разрядов надо увеличить АЦП сигнала, чтобы отно¬ шение дисперсии сигнала к дисперсии шумов квантования увеличилось на 18 дБ? 3. Биполярный аналоговый сигнал на входе пятиразрядного АЦП имеет максимальную амплитуду 200 мВ. Вычислите дисперсию шумов квантования сигнала на выходе цифрового фильтра первого порядка, описываемого разностным уравнением с комплексным коэффициентом: у(п) = х(п) + ау(п - 1), ах = 0,8 ехр(/ШЦ. 4. Биполярный аналоговый сигнал на входе восьмиразрядного АЦП имеет максимальную амплитуду 1 В. Вычислите дисперсию шумов кван¬ тования сигнала на выходе цифрового фильтра, описываемого уравне¬ нием у(п) = 0,8 х(п) - 0,5х(п - 1) + 0,4х(я - 2). 5. Биполярный аналоговый сигнал, имеющий максимальную ампли¬ туду 2 В после преобразования в четырехразрядном АЦП, поступает на вход рекурсивного цифрового фильтра первого порядка с коэффи¬ циентом а = 0,5. После умножения на коэффициент результат записы¬ вается в четырехразрядный регистр. Какова дисперсия ошибки сигнала на выходе, обусловленной округлением промежуточных результатов? 6. Определите дисперсию ошибки выходного сигнала на выходе нерекурсивного ЦФ четвертого порядка, обусловленную округлением результатов умножения до 4-х разрядов. Входной сигнал преобразован в четырехразрядном АЦП с шагом квантования 0,5 В. Коэффициенты фильтра задаются четырехразрядным двоичным кодом. 7. Чем объясняются искажения АЧХ и ФЧХ ЦФ при квантовании коэффициентов фильтра? 8. Какая форма реализации рекурсивных фильтров наименее чув¬ ствительна к квантованию коэффициентов?
4. СПЕЦИАЛЬНЫЕ АЛГОРИТМЫ ЦИФРОВОЙ ОБРАБОТКИ СИГНАЛОВ В РАДИОТЕХНИЧЕСКИХ И ТЕЛЕКОММУНИКАЦИОННЫХ СИСТЕМАХ 4Л. Изменение частоты дискретизации в линейных цифровых фильтрах В процессе обработки дискретных и цифровых сигналов часто возникает необходимость уменьшения или увеличения частоты дискретизации. Уменьшение частоты используется для снижения объема вычислений, когда в процессе обработки полоса частот, занима¬ емая сигналом, уменьшается вследствие низкочастотной филь¬ трации. Преобразование частоты дискретизации приходится при¬ менять при обмене данными между двумя процессорами ЦОС, работающими с разными тактовыми частотами, при объединении изображений, полученных в разных спектральных каналах с раз¬ личной разрешающей способностью. Преобразование частоты дискретизации позволяет также упростить реализации некоторых типов узкополосных цифровых фильтров. Процедуры понижения и повышения частоты дискретизации носят название децимации и интерполяции сигналов. Следует отметить, что изначально термин «децимация» означал умень¬ шение в десять раз. В настоящее время этот термин использу¬ ется в цифровой обработке сигналов для обозначения понижения частоты дискретизации в любое целое число раз. Процедура интерполяции реализует получение оценок значе¬ ний промежуточных отсчетов. Системы, реализующие процедуры децимации и интерполя¬ ции сигналов, являются частным случаем так называемых нисхо¬ дящих и восходящих дискретных систем. Н исходящей дискретной системой (НДС) назы¬ вают систему, частота дискретизации сигнала на выходе которой Ffll ниже частоты дискретизации сигнала на входе Fa. 126
Восходящей дискретной системой (ВДС) назы¬ вают систему, частота дискретизации сигнала на выходе которой выше частоты дискретизации сигнала на входе Fл. Простейшая НДС (рис. 4.1) состоит из дискретного фильтра, осуществляющего предварительную обработку входного сигнала с частотой дискретизации Fд, и элемента, уменьшающего частоту дискретизации в т раз — компрессора частоты дис¬ кретизации (КЧД), стоящего на выходе системы. Предвари¬ тельный фильтр (ПФ) необходим для предотвращения или умень¬ шения наложения спектров, возникающего в процессе децимации сигнала. Простейшая ВДС (рис. 4.2) состоит из элемента, увеличиваю¬ щего частоту дискретизации в т раз, — экспандера частоты дискретизации (ЭЧД), и дискретного фильтра, осуществля¬ ющего обработку сигнала с частотой дискретизации Рд1 (интерпо¬ ляционного фильтра). ДиГд) у(пТд) Уі(уТДІ) Рис. 4.1. Структура простейшей нисходящей дискретной системы фТд) х(пТд1) у(пТд і) Рис. 4.2. Структура простейшей восходящей дискретной системы При аппаратурной или программной реализации НДС и ВДС операции, выполняемые компрессором и экспандером частоты дискретизации, обычно совмещаются с операциями, выполня¬ емыми дискретным фильтром. Однако при анализе работы НДС и ВДС их целесообразно выделять в отдельный блок. 127
Проанализируем простейшую нисходящую дискретную систему. Связь между входной и выходной последовательностями простейшей НДС (дециматора) во временной области может быть представлена следующим образом. Процедура предварительной фильтрации описывается выражением п у(п) = ^Щ)х(п-1), (4.1) 1=0 где h(n) — импульсная характеристика ПФ. КЧД представляет собой ключ, который замыкается в моменты времени t = птТл. п = 0, 1, 2, ... . Таким образом, из входного сигнала, описываемого решетчатой функцией у(пТд) с «малым» интервалом дискретизации Тд, выбирается только каждый т-й отсчет. В результате чего формируется выходной сигнал, представ¬ ляющий собой отсчеты исходного аналогового сигнала, следую¬ щие с «большим» интервалом дискретизации ТДІ = тТд (рис. 4.3). Тем самым КЧД прореживает последовательность у{п), пропуская на выход системы каждый т-й отсчет и формируя, таким образом, выходную последовательность уДѵ) в соответствии с правилом: Ті(ѵ^д1) = }’І(ѵтТд) =y(\’mTJ =у(пТ). (4.2) Здесь ѵ = 0, 1, 2,. .. — номер выходного отсчета; п = ѵт — номер входного отсчета. Объединяя алгоритм преобразования дискретного сигнала в КЧД (4.2) и уравнение работы предварительного фильтра (4.1), получаем уравнение работы простейшей НДС: ѵт У\(ѵ) = У(ѵт)= 'YJh{l)x{vm-l). (4.3) /=о Отметим одну особенность алгоритма фильтрации с пониже¬ нием частоты дискретизации. Если обработку в обычном линей¬ ном ЦФ (4.1) можно трактовать как суммирование входной после¬ довательности х(п) в скользящем через один отсчет окне h(n), то обработка с понижением частоты дискретизации (4.3) представ¬ ляет собой суммирование в «прыгающем» через т отсчетов окне h(n). 128
ѵ(иТ’д) t = ѵтТл Рис. 4.3. Сигнал на входе (а) и выходе (б) компрессора частоты дискретизации Из выражения (4.3) также следует, что системы с пониже¬ нием частоты дискретизации (дециматоры) не инвариантны к вре¬ менному сдвигу и имеют m различных импульсных характери¬ стик (реакций на входную последовательность в виде дискретной 5-функции). Для определения преобразования спектра сигнала в КЧД уста¬ новим связь между Z-преобразованиями входной у(п) и прорежен¬ ной у, (ѵ) последовательностей: ~ 00 Y(z)= 2 y(n)z~n; (4.4) п- О ~ 00 Yi{zi) = 2>’і(ѵ)г1-ѵ. (4.5) ѵ = О Спектры находятся подстановкой в (4.4) и (4.5): z = е JwTn j«>T ДІ _ jcomT (4.6) 129
Из (4.6) следует, что zx = z” и Y,(zm)= Z Л (v)z“ (4.7) Для установления связи между F(z) и Fj(Zj) рассмотрим сумму вида m i 2 f(zej2%qlm) q= О Ѵ 7 (4.8) и преобразуем ее с учетом (4.4): "Xy( zQjlliqlm\ = £ f; y(n)z-"e-J2ngn/m 9=0 <7=0 n=0 oo m I S(Se-;2ngn/m)y(n)z-. (4.9) n=0 <7=1 Внутренняя сумма (в скобках) представляет собой сумму ш членов геометрической прогрессии с первым членом, равным единице, и знаменателем, равным ет,2™/т. Поэтому m 1 Ее -jlnqnjm q=о ГW при И = V7W, V = 0, 1, 2, ..., I 0 при других п. (4.10) Подставив (4.10) в (4.9) и заменив п = ѵт, получим: X Y(zej2%q/m) = m E y(ym)z =т^Уі(v)z,v = /wfj (zj). (4.11) # = 0 v=0 v=0 Откуда следует, что искомая связь между Z-преобразованиями входной и прореженной последовательностей имеет вид: 1 т— 1 ^ (z) = — E F(zeJ2ltq/m). (4.12) m q=o Связь между спектрами входной и прореженной последова¬ тельностей находится подстановкой (4.6) в (4.12): Y, (е 'ю7д1) = - E F{exp[/T> + ^%)]}. (4.13) ' ' m q=o mla Рассматривая последнее выражение как функцию частоты, мы можем сделать вывод, что спектр выходного сигнала КЧД пред¬ ставляет собой сумму спектров входного сигнала, сдвинутых один относительно другого по оси частот на величину Ъі/тТу. 1 т-1 Fi О) = — YJY\j(a> + 2%q/mT )\. (4.14) m q=0 130
Наложение спектров при уменьшении частоты дискретизации отсутствует, если спектр сигнала на входе компрессора занимает полосу частот: -л/тТа < со < п/пгТд. Теперь становится ясным назначение ПФ: он подавляет спек¬ тральные составляющие сигнала за пределами этой полосы и сни¬ жает тем самым эффект наложения спектров за счет понижения частоты дискретизации Спектр выходной последовательности ПФ представляет собой произведение спектра входной последовательности НДС и КЧХ ПФ: Т(/оо) X(jo>)Hm,(jo>). (4.15) Подставляя (4.15) в (4.14), получаем выражение, связывающее спектры входного и выходного сигналов в простейшей НДС (деци- маторе): 1 m-l ГіО'ю) ^ — Z *[/(“ + 2щ / тТд)]Нпф [/(со + 2жq /тТд)]. (4.16) Л1 q = О Первое слагаемое в правой части этого выражения (соответ¬ ствующее q = 0) представляет собой спектр полезного сигнала. Прочие слагаемые можно рассматривать как спектры помех, иска¬ жающих полезный сигнал в основной полосе, если АЧХ ПФ не обеспечивает требуемого подавления за пределами основной полосы сигнала. Указанный эффект проиллюстрирован на рис. 4.4. В простейшей ВДС (интерполяторе) экспандер частоты дис¬ кретизации преобразует входную дискретную последовательность х(ѵГд) с периодом дискретизации Тд в выходную последователь¬ ность х^пТд) с периодом Т* = TJm по следующему алгоритму: ( т \ х(уТя) 4)11" = Vm’ ѵ = °’ 2’ •••; (4.17) [ 0 при n ^ vm. 131
Х(со) Рис. 4.4. Спектральная плотность сигнала на входе (а) и выходе (б) простейшей нисходящей дискретной системы (т = 2) Таким образом, выходная последовательность в экспандере частоты дискретизации получается путем ввода между двумя соседними входными отсчетами (т - 1)-го нулевого отсчета. При таком преобразовании сигнала (дополнении нулями) выходной сигнал имеет тот же спектр, что и входной. При этом спектр выходного сигнала остается периодическим со «старой» (низкой) частотой дискретизации (рис. 4.5, а). Если этот сигнал проходит затем через идеальный ФНЧ (его АЧХ показана на рис. 4.5, б тонкой линией), то «лишние» спек¬ тральные полосы будут удалены и выходной сигнал на выходе такого фильтра будет точно соответствовать сигналу х(пТл). т. е. частота дискретизации его увеличится в m раз, а промежуточные отсчеты будут восстановлены (рис. 4.5, в). 132
Рис. 4.5. Преобразование спектра дискретного сигнала в интерполяторе: а — вход экспандера; б — выход экспандера, в — выход интерполятора а б в поскольку практически идеальный ФЫЧ нереализуем, полу¬ ченная последовательность не будет являться точной интерполя¬ цией х(пТд1). Точность интерполяции зависит от уровня подавле¬ ния в полосе задерживания интерполирующего ФЫЧ. Чем больше подавление, тем точнее интерполяция. 133
4.2. Демодуляция узкополосных сигналов. Цифровые преобразователи Гильберта 4.2.1. Узкополосные сигналы, комплексная огибающая, аналитический сигнал В радиотехнических и телекоммуникационных системах часто приходится иметь дело с сигналами, спектр которых сосредо¬ точен вблизи некоторой частоты. К ним, в частности, относятся радиосигналы, полученные в результате одновременной модуля¬ ции по амплитуде и частоте. Такие колебания называются «узкопо¬ лосными». По определению сигнал считается узкополосным, если его спектральная плотность Л'(/ю) отлична от нуля лишь в преде¬ лах частотных интервалов Лю, лежащих в окрестностях частот ±ю0, причем ширина спектра Лю « ю0. (рис. 4.6). S(ja) А I -®0 0 Рис. 4.6. Спектр узкополосного сигнала Математическая модель узкополосного процесса представляет собой описание квазигармонического колебания: s(t) = A(t) cos ¥(0 = A(t) cos [go,/ + \|/(0] • (4.18) Здесь A(t) огибающая узкополосного процесса (закон ампли¬ тудной модуляции); *Р(0 = ю,/ + \|/(0 — полная фаза; \|/(0 — закон фазовой модуляции; ю(t)= ю0+ d\\i(t)/c/l — мгновенная частота. Частоту ю0 называют опорной частотой, для симметричной спектральной плотности это центральная частота, однако в общем случае ее выбор произволен. 134
Узкополосный сигнал (4.18) может быть представлен в виде суммы синфазной и квадратурной составляющих: s(t) = Ac(i) cos со/ - A^t) sin со0t, (4.19) где синфазная и квадратурная амплитуды: Ac(t) = A(t) cos \\f(t); As(t) = A(t) sin i|/(f) являются медленно меняющимися (по сравнению с cos со0^ и sin со0t) функциями времени и содержат всю информацию об амплитудной и фазовой модуляции. Получить квадратурные амплитуды можно с помощью квадратурного детектора, структурная схема которого приведена на рис. 4.7. Рис. 4.7. Структурная схема квадратурного детектора Закон амплитудной и фазовой модуляции узкополосного сиг¬ нала описывает комплексная огибающая — комплексная функция времени, действительная часть которой представляет собой син¬ фазную амплитуду, а мнимая — квадратурную: A(t) = A(t)e^=Ac(t) + j As(t). (4.21) Модель узкополосного процесса в виде (4.18) предложил вели¬ кий венгерский физик, основоположник голографии, лауреат 135
Нобелевской премии Денеш Габор (Gabor) в 1946 г. (Gabor D. Theory of communication // J. IEE. 1946. Pt. III. Vol. 93. P. 429-457). Введение понятия «комплексная огибающая» позволяет упро¬ стить решение многих задач, связанных со сложными видами модуляции, демодуляции и другими преобразованиями узкополос¬ ных сигналов. Однако представление узкополосного процесса в виде (4.18- 4.20) несет в себе некоторую неоднозначность, поскольку при любой заданной функции \)>(t) можно подобрать соответствую¬ щую функцию A(t), чтобы эти равенства выполнялись. Эту неод¬ нозначность выбора A(t) и ѵ)/(t) можно сформулировать следу¬ ющим образом: «А что такое огибающая и фаза узкополосного процесса?» Огибающая A(t), в частности, должна удовлетворять вполне естественным требованиям, вытекающим из ее наименования: — в любой момент времени A(t) > s(t); — в точках соприкосновения графиков, т. е. при A(t) = s(t), dAldt = ds/dt. Установить однозначную связь между огибающей и фазой узкополосного процесса позволяет понятие аналитического (ком¬ плексного) сигнала, действительная часть которого представляет собой исходное узкополосное колебание: Мнимая часть аналитического сигнала называется сопря¬ женным сигналом или квадратурным дополнением. Она связана с исходным сигналом преобразованием Гильберта: S(t) = s(t) +jS(f); s(t) = Re[*S(0] =A(t) cos [a></ + \|/ (0]. (4.22) 71 —со t X При этом аналитический сигнал: S(t) = s(t) + js(t) =A(t)e>rvh. (4.23) (4.24) 136
Из этого представления следуют однозначные способы опреде¬ ления огибающей, фазы и частоты узкополосного колебания: A(t) = ИО + s2(t)],/2; (4.25) Т(0 = arctg^-; s(t) со (t) = cW(f) dt d г x s(t) п —[arctg——] = oo0 + dt s(t) d\\r(t) dt (4.26) (4.27) Таким образом, для определения закона амплитудной, фазо¬ вой и частотной модуляции узкополосного процесса надо уметь вычислять преобразование Гильберта. Из формулы прямого преобразования Гильберта (4.2) сле¬ дует, что сопряженный сигнал представляет собой свертку исход¬ ного узкополосного сигнала s(t) и функции 1/л/. Это означает, что преобразователь Гильберта (ПГ) представляет собой линейный фильтр с импульсной характеристике Гг , /а — П r(t) = \%t Oj = 0. (4.28) Комплексная частотная характеристика такого фильтра — ПГ при этом имеет вид: j, со < 0; со 1 HT(Ja)= J —exp(-jat)dt = < 0, со = 0; (4.29) —со ТС t -і, со > 0. ѵ. J Обратим внимание, что комплексная частотная характеристика чисто мнимая. Таким образом, ПГ представляет собой идеальный фазовращатель. Его АЧХ равна 1 на всех частотах, за исключением нулевой частоты (рис. 4.8). Устройство с такими характеристи¬ ками называют «квадратурным фильтром». Идеальный квадратур¬ ный фильтр физически нереализуем, поэтому в реальных системах речь может идти о той или иной степени приближения. 137
Я(ш) 1 а со ф(ш) 71 2 71 2 Рис. 4.8. Амплитудно-частотная (а) и фазочастотная (б) характеристики преобразователя Гильберта Использование понятий «аналитический сигнал» и «комплекс¬ ная огибающая» является очень удобным инструментом для ана¬ лиза сложных видов модуляции и преобразования узкополосных сигналов благодаря следующим свойствам: 1. Произведение аналитического сигнала на комплексно сопря¬ женный ему сигнал равно квадрату огибающей исходного (физи¬ ческого) сигнала: S{t)S\t)=A\t\ (4.30) т. е. модуль аналитического сигнала равен огибающей физиче¬ ского сигнала. 2. Спектр аналитического сигнала содержит только положи¬ тельные частоты: |25(/со),... со > 0; [ 0, ...со<0. (4.31) 138
3. Спектральная плотность комплексной огибающей совпадает со смещенной на величину ю0 влево спектральной плотностью аналитического сигнала: Л(/ш) = 5а[/(ш + со0)]. (4.32) Эти свойства проиллюстрированы на рис. 4.9. а -ю0 S(j&) >(0 £а(/ю) A(Jw) Рис. 4.9. Спектры действительного узкополосного сигнала (а), соответствующего ему аналитического сигнала ( и комплексной огибающей (в) 139
4.2.2. Дискретное преобразование Гильберта Спектр дискретного узкополосного процесса периодичен с периодом, равным частоте дискретизации. Поэтому в дискрет¬ ном преобразовании Гильберта спектральные соотношения (4.29) должны выполняться в основной полосе частот, т. е. в диапазоне от -л/Гд до п/Тд. Тогда комплексная частотная характеристика дис¬ кретного преобразователя Гильберта (ДПГ) должна иметь вид (рис. 4.10): 271 271 j,(k-l/2)—<(o<k—, 0,1,2,...; гд тл ЯдГ(/со) = - 0, оо = к - (4.33) 7 2тг к—< со < (к +1/2)—. Для определения импульсной характеристики фильтра, реали¬ зующего дискретное преобразование Гильберта, вычислим обрат¬ ное ДВПФ на интервале ее периодичности. Введя нормированную частоту Q. = оТд, получим: КЛ”) = 1 HaT(JCi)ejQndQ = 2-[ 1 уе^оЮ]. (4.34) 271 -л 271 -тс 0 ЯдгО'Я) Рис. 4.10. Комплексная частотная характеристика дискретного преобразования Гильберта 140
Откуда следует, что импульсная характеристика дискретного фильтра — преобразователя Гильберта описывается выражением J sm'( — ) hT(n) = \—(1-COSTO?) = -2—, пФ 0; R П 71 71?? / 2 0,..., п = 0. (4.35) Все отсчеты этой последовательности с четными номерами равны нулю, а с нечетными 21 пи. Вид импульсной характеристики ДПГ приведен на рис. 4.11 Рис. 4.11. Импульсная характеристика дискретного преобразователя Гильберта 4.2.3. Реализация дискретного преобразователя Гильберта на основе КИХ-фильтра Как показано выше, преобразователь Гильберта имеет линей¬ ную фазочастотную характеристику. Такой характеристикой 141
обладают нерекурсивные цифровые фильтры. Поэтому ДПГ может быть реализован на основе КИХ-фильтра с антисимметричной импульсной характеристикой. Однако такой фильтр с импульсной характеристикой вида (4.35) физически нереализуем, поскольку она является бесконечно протяженной в обоих направлениях (рис. 4.11). Простое усечение и сдвиг ИХ приводят к уже рассмот¬ ренному эффекту Гиббса. Поэтому дискретное преобразование Гильберта реализуется лишь приближенно. Структурная схема ДПГ на базе КИХ-фильтра 10-го порядка (L= 11) приведена на рис. 4.12. Рис. 4.12. Структурная схема преобразователя Гильберта на основе КИХ-фильтра 10-го порядка Импульсная характеристика такого фильтра получается из выражения (4.35) путем усечения последовательности кдТ(п) в пределах п ± (L - 1)/2 и сдвига вправо на (L - 1)/2 отсчет. Такой фильтр относится к так называемой третьей форме КИХ-фильтров с линейной ФЧХ (.L нечетное, антисимметричная импульсная 142
характеристика). Его комплексная частотная характеристика опи¬ сывается выражением ^ Я(е;П) = е;7І/2 ^ с(п) sin(Q^). (4.36) п = \ Здесь коэффициенты фильтра с{п) представляют собой удвоен¬ ные значения сдвинутой импульсной характеристики (4.35): с(п) = 2йдГ(и - -),... п = 1, 2, ...^—(4.37) Комплексная частотная характеристика Н(о>Сі) такого филь¬ тра является чисто мнимой функцией частоты, причем на нулевой частоте она равна нулю, а АЧХ описывается выражением tf(Q) = L-\ Л ^c(n)sin(Qn) п = 1 для -л < Q < л. (4.38) АЧХ и ФЧХ цифрового преобразователя Гильберта на базе КИХ-фильтра 10-го порядка приведены на рис. 4.13, а. Для сравне¬ ния на рис. 4.13,6 показана АЧХ для КИХ-фильтра 22-го порядка. Для уменьшения пульсаций АЧХ при проектировании преоб¬ разователя Гильберта задают переходные полосы вблизи нуле¬ вой частоты и частоты, равной половине частоты дискретизации. Считается, что АЧХ фильтра, реализующего ДПГ, в этих поло¬ сах может меняться произвольно, а основные соотношения (4.35) должны выполняться только за пределами переходных полос (в рабочей полосе). При таких ограничениях цифровой преобра¬ зователь Гильберта может быть реализован и на базе рекурсив¬ ного ЦФ Применение ДПГ позволяет сформировать дискретный аналог аналитического сигнала и затем организовать вычисление огибаю¬ щей А(п) и фазы \\і(п) дискретной последовательности s(n), полу¬ ченной путем дискретизации узкополосного сигнала (4.18). Струк¬ турная схема формирователя аналитического сигнала приведена на рис. 4.14. Согласующая цифровая линия задержки компенсирует 143
групповое время запаздывания в ДПГ и обеспечивает временное согласование последовательностей s'(//)h s(n)- -3,142 3,142 Q -3,142 Q 3,142 40 Рис. 4.13. Амплитудно-частотная и фазочастотная характеристики преобразователя Гильберта на основе КИХ-фильтра: а — 10-й порядок; б — 22-й порядок s(n) Согласующая цифровая линия задержки Ап) s(n) = s\n) +js(n) Цифровой преобразователь Гильберта s(n) Рис. 4.14. Формирование дискретного аналога аналитического сигнала 144
Поскольку групповая задержка в КИХ-фильтре третьей формы с линейной ФЧХ равна: Ѵ=(£-1)Гд/2, то для согласования во времени последовательностей, описыва¬ ющих действительную и мнимую части аналитического сигнала, последовательность х(п) должна быть задержана на (L - 1)/2 отсче¬ тов. Поэтому задержанную копию входной последовательности s’(n) можно взять со среднего отвода цифровой линии задержки ДПГ на базе КИХ-фильтра, как это показано на рис. 4.12. Применение преобразования Гильберта позволяет однозначно определить законы амплитудной и фазовой модуляции узкопо¬ лосного сигнала. Поэтому методы цифровой обработки сигналов, использующие преобразование Гильберта, применяются в радио¬ технических и телекоммуникационных системах при решении сле¬ дующих задач: — оценка мгновенной частоты (частотная демодуляция); — демодуляция однополосных AM сигналов; — квадратурная модуляция и демодуляция; — перенос спектра по частоте (цифровое гетеродинирование); — формирование и обработка сложных сигналов; — формирование радиолокационных изображений высокого разрешения в системах с синтезированной апертурой антенны. Преобразование Гильберта используется также при анализе сейсмических данных, в компьютерной томографии, при сжатии аудиосигналов и цветных изображений, формировании изображе¬ ний в приемниках телевидения высокой четкости и других обла¬ стях. Этим объясняется особое внимание к рассмотрению таких специфических ЦФ. КОНТРОЛЬНЫЕ ВОПРОСЫ 1. Каково назначение предварительного фильтра в простейшей системе с понижением частоты дискретизации (дециматоре)? 145
2. Какой должна быть КЧХ идеального интерполяционного фильтра в простейшей системе с повышением частоты дискретизации (интерпо¬ ляторе)? 3. Запишите алгоритм работы однородного трансверсального фильтра третьего порядка при понижении частоты дискретизации на его выходе в 2 раза. 4. Запишите алгоритм работы рекурсивного фильтра первого порядка с коэффициентом а при понижении частоты дискретизации на его выходе в 2 раза. 5. Как связан спектр выходного сигнала компрессора частоты дис¬ кретизации со спектром входного сигнала? 6. Какой сигнал называется узкополосным? 7. Определите сопряженный сигнал и огибающую для веществен¬ ного сигнала — гармонического колебания вида s(t) = cos (2тгft + Ѳ). 8. Какой тип дискретного фильтра приближенно реализует дискрет¬ ное преобразование Гильберта?
БИБЛИОГРАФИЧЕСКИЕ ССЫЛКИ 1. Айфичер Э„ Джервис Б. Цифровая обработка сигналов. Практиче¬ ский подход / [пер. с англ. И. Ю. Дорошенко, А. В. Назаренко ; под ред. А. В. Назаренко]. 2-е изд. М.; СПб.; Киев : Вильямс, 2004. 992 с. 2. Баскаков С. II. Радиотехнические цепи и сигналы : учебник. 4-е изд., испр. и доп. М. : ЛЕНАНД/URSS, 2016. 528 с. 3. Блейхут Р. Быстрые алгоритмы цифровой обработки сигналов : пер. с англ. М. : Мир, 1989. 448 с. 4. Васильев В. П., Мѵро Э. Л., Смольский С. М. Основы теории и рас¬ чета цифровых фильтров : учеб, пособие для высш. учеб, заведений / под ред. С. М. Смольского. М. : Изд. центр «Академия», 2007. 272 с. 5. Гадзиковский В. II. Цифровая обработка сигналов. М. : Солон- Пресс, 2013. 766 с. 6. Гадзиковский В. И., Калмыков А. А. Теория и проектирование устройств цифровой фильтрации : учеб, пособие. Екатеринбург : ГОУ ВПО «УГТУ-УПИ», 2006. 433 с. 7. Дьяконов В., Абраменкова II. MATLAB. Обработка сигналов и изо¬ бражений : спец, справочник. СПб. ; М. ; Харьков ; Минск : Питер, 2002. 608 с. : ил. 8. Заалъ Р. Справочник по расчету фильтров : пер. с нем. М. : Радио и связь, 1985. 9. Коберниченко В. Г. Расчет и проектирование цифровых фильтров : учеб.-метод. пособие. Екатеринбург : Изд-во Урал, ун-та, 2013. 64 с. 10. Куприянов М. С., Матюшкин Б. Д. Цифровая обработка сигна¬ лов: процессоры, алгоритмы, средства проектирования. 2-е изд., перераб. и доп. СПб. : Политехника, 2002. 592 с. 11. Курячий М. И. Цифровая обработка сигналов : учеб, пособие для вузов. Томск : Томск, гос. ун-т систем упр. и радиоэлектроники, 2009. 190 с. 12. Лайонс Р. Цифровая обработка сигналов. 2-е изд.: пер. с англ. М.: ООО «Бином-Пресс», 2006. 656 с. 13. Лэм Г. Аналоговые и цифровые фильтры. Расчет и реализация : пер. с англ. М.: Мир, 1982. 592 с. 147
14. Маккеллан Дж. X., Рейд Ч. М. Применение теории чисел в цифро¬ вой обработке сигналов / пер. с англ, под ред. Ю. И. Минина. М. : Радио и связь. 1983. 264 с. 15. Оппенгейм А., Шафер Р. Цифровая обработка сигналов / пер. с англ. С. А. Кулешова под ред. А. Б. Сергиенко. 2-е изд., испр. М. : Тех¬ носфера. 2007. 856 с. 16. Рабинер Р, Гоулд Б. Теория и применение цифровой обработки сигналов. М. : Мир. 1978. 848 с. 17. Сергиенко А. Б. Цифровая обработка сигналов : учеб, пособие для студ. вузов, обучающихся по направлению подгот. «Информатика и вычисл. техника». 3-е изд. СПб.: БХВ-Петербург, 2011. 768 с. 18. Сиберт У. М. Цепи, сигналы, системы : в 2 ч. : пер. с англ. М. : Мир. 1988. Ч. 1.336 с. 19. Сиберт У. М. Цепи, сигналы, системы : в 2 ч. : пер. с англ. М. : Мир, 1988. Ч. 2. 360 с. 20. Солонина А. II, Арбузов С. II. Цифровая обработка сигналов. Моделирование в MATLAB : учеб, пособие для студ. вузов, обучающихся по направлению подгот. дипломир. специалистов 210400 «Телекоммуни¬ кации». СПб. : БХВ-Петербург, 2008. 816 с. 21. Солонина А. И., Улахович Д. А., Арбузов С. М. и др. Основы циф¬ ровой обработки сигналов : курс лекций. 2-е изд., испр. и перераб. СПб. : БХВ-Петербург, 2013. 768 с. 22. Солонина А. II., Улахович Д. А., Яковлев Л. Алгоритмы и про¬ цессоры цифровой обработки сигналов : учеб, пособие для студ., обу¬ чающихся по направлению 654400 «Телекоммуникации». СПб. : БХВ- Петербург, 2002. 454 с. 23. Тропченко А. Ю., Тропченко А. А. Цифровая обработка сигналов. Методы предварительной обработки : учеб, пособие. СПб. : СПбГУ ИТМО, 2009. 100 с. 24. Федосов В. П., Нестеренко А. К. Цифровая обработка сигналов в Lab VIEW : учеб, пособие / под ред. В. П. Федосова. М. : ДМК Пресс, 2007. 456 с. 25. Якимов Е. В., Вавилова Г. В., Клубоеич II. А. Цифровая обработка сигналов : учеб, пособие. Томск : Изд-во Томск, политехи, ун-та, 2008. 307 с.
Учебное издание Коберниченко Виктор Григорьевич основы ЦИФРОВОЙ ОБРАБОТКИ СИГНАЛОВ Учебное пособие Заведующий редакцией М. А. Овечкина редактор С. Г. Галинова корректор С. Г. Галинова Оригинал-макет Л. А. Хухаревой
Подписано в печать 15.11.2018 г. Формат 60 х 841/16. Бумага офсетная. цифровая печать. усл. печ. л. 8,83. уч.-изд. л. 7,5. тираж 50 экз. заказ 310. издательство уральского университета редакционно-издательский отдел ипц урФу 620083, Екатеринбург, ул. тургенева, 4. тел.: +7 (343) 389-94-79, 350-43-28 E-mail: rio.marina.ovechkina@mail.ru отпечатано в издательско-полиграфическом центре урФу 620083, екатеринбург, ул. тургенева, 4. тел.: +7 (343) 358-93-06, 350-58-20, 350-90-13 Факс +7 (343) 358-93-06 http://print.urfu.ru
для заметок
для заметок
ѵн Ив КОБЕРНИЧЕНКО ВИКТОР ГРИГОРЬЕВИЧ Профессор департамента радиоэлектроники и связи Института радиоэлектроники и информа¬ ционных технологий Уральского федерального университета, научный руководитель Центра космического мониторинга. Закончил радиотехнический факультет Уральского политехни¬ ческого института по специальности «Радиоэлектронные устройства» (1965). Кандидат технических наук (1974), доцент (1981). Заведующий кафедрой теоретических основ радиотехники (1999-2009). Автор более 200 научных работ. Лауреат премии Правительства Российской Федерации в области образования (2006). Область научных интересов: радиолокационные системы с синтезированной апертурой, цифровая обработка информа¬ ции в радиоэлектронных системах космического дистанционного зондирования Земли, применение информационных технологий в задачах обеспечения природной и техногенной безопасности.