Учебное пособие для студентов высших учебных заведений
Скачать 5.41 Mb.
|
maxflat производит расчет обобщенного цифрового фильтра Баттерворта. Формы обращения к ней таковы: [ b, a] = maxflat(nb, na, Wc) [ b, a] = maxflat(nb, 'sym', Wc) [ b, a, b1, b2] = maxflat(nb,na, Wc) [ b, a] = maxflat(nb,na, Wc, 'design_flag'). Первое обращение позволяет вычислить коэффициенты b - числителя и а - знаменателя дискретной передаточной функции H(z) цифрового ФНЧ Баттерворта с частотой среза Wc, порядок числителя которой равен nb, а знаменателя - na. При обращении второго вида вычисляются коэффициенты цифрового сим- метричного КИХ-фильтра Баттерворта. В этом случае na принимается равным 0. Параметр nb должен быть четным. 5.4. Проектирование фильтров 235 Если обратиться к процедуре так, как указано в третьем обращении, т. е. указать в качестве выходных четыре величины, то дополнительные параметры b1 и b2 дадут коэффициенты двух полиномов, произведение которых является поли- номом числителя b искомой дискретной передаточной функции, причем все нули полинома b1 равны -1, а полином b2 содержит все остальные нули полинома b. Добавление в список входных параметров процедуры параметра 'design_flag' позволяет изменять характер выводимой на экран информации. Если значение этого параметра равно 'trace', на экране отображаются параметры, ис- пользуемые в процессе проектирования: >> nb=10; na = 2; w = 0.5; >> [b,a,b1,b2] = maxflat(nb,na,w,'trace') Table: L M N wo_min/pi wo_max/pi 10.0000 0 2.0000 0 0.2394 9.0000 1.0000 2.0000 0.2394 0.3259 8.0000 2.0000 2.0000 0.3259 0.3991 7.0000 3.0000 2.0000 0.3991 0.4669 6.0000 4.0000 2.0000 0.4669 0.5331 5.0000 5.0000 2.0000 0.5331 0.6009 4.0000 6.0000 2.0000 0.6009 0.6741 3.0000 7.0000 2.0000 0.6741 0.7606 2.0000 8.0000 2.0000 0.7606 1.0000 b = Columns 1 through 7 0.0414 0.2263 0.4932 0.5252 0.2494 0.0079 -0.0266 Columns 8 through 11 0.0008 0.0023 -0.0004 -0.0000 a = 1.0000 0 0.5195 b1 = 1 6 15 20 15 6 1 b2 = 0.0414 -0.0221 0.0049 -0.0004 -0.0000 Если же этому параметру задать значение 'plots', то на экран будут выведе- ны графики амплитудной характеристики, группового времени замедления, а также графическое изображение нулей и полюсов: >>[b,a,b1,b2] = maxflat(nb,na,w,'plots') b = Columns 1 through 7 0.0414 0.2263 0.4932 0.5252 0.2494 0.0079 -0.0266 Columns 8 through 11 0.0008 0.0023 -0.0004 -0.0000 a = 1.0000 0 0.5195 b1 = 1 6 15 20 15 6 1 b2 =-0.0221 0.0049 -0.0004 -0.0000 5.4. Проектирование фильтров 236 Рис. 5.46 Расчет БИХ-фильтра по заданной амплитудно-частотной характеристи- ке производится процедурой yulewalk. Если набрать в командном окне MatLAB строку [ b, a] = yulewalk(n, f, m), то будут вычислены коэффициенты b - числителя и a - знаменателя дискретной передаточной функции БИХ-фильтра порядка n, АЧХ которого задана векторами f (частоты в нормированных значениях) и m - соответствующих значений отноше- ний амплитуд выхода и входа. Первый элемент вектора f должен быть равен 0, а последний 1. Все остальные элементы должны быть расположены в неубываю- щем порядке. Частоты, при которых происходит скачок АЧХ, указываются два раза с разными значениями соответствующих им "амплитуд". Рис. 5.47 Приведем пример расчета ФНЧ 8-го порядка и построим желаемую АЧХ и АЧХ полученного фильтра. Результат приведен на рис. 5.47. f = [ 0 0.5 0.5 1]; m=[1 1 0 0]; [b,a] = yulewalk(8,f,m); [h,w] = freqz(b,a, 128); 5.4. Проектирование фильтров 237 plot(f,m, w/pi,abs(h)) grid, title('Пример применения процедуры YULEWALK') xlabel('Нормализованная частота'), ylabel('А Ч Х') 5.4.4. Проектирование КИХ-фильтров В отличие от БИХ-фильтров, которые характеризуются двумя векторами - коэффициентов b числителя и а - знаменателя своей дискретной передаточной функции, КИХ-фильтры описываются только одним вектором b. Знаменатель их дискретной передаточной функции тождественно равен 1. Группа функций fir1 предназначена для расчета коэффициентов b цифрово- го КИХ-фильтра с линейной фазой методом взвешивания с использованием окна. Общий вид обращения к этой процедуре имеет вид b = fir1(n, Wn, 'ftype', window). Процедура вычисляет вектор n+1 коэффициентов b КИХ-фильтра с норма- лизованной частотой среза Wn. Параметр 'ftype' задает желаемый тип фильтра (ФНЧ, ФВЧ, полосовой или режекторный). Он может отсутствовать - и тогда по умолчанию рассчитываются параметры ФНЧ c частотой среза Wn, если последняя задана как скаляр, или по- лосовой фильтр с полосой пропускания от W1 до W2, если Wn задан в виде век- тора из двух элементов [W1 W2], - или принимать одно из четырех значений : 'high', 'stop', 'DC-1' и 'DC-0'. В первом случае синтезируется ФВЧ с частотой среза Wn. Во втором, - режекторный фильтр (при этом Wn должен быть вектором из двух элементов, значения которых определяют границы полосы задерживания по отношению к частоте Найквиста). В третьем случае рассчитываются параметры многополосового фильтра, первая полоса которого является полосой пропуска- ния, а в четвертом, - тоже многополосовый фильтр, первая полоса которого явля- ется полосой задерживания. При расчете режекторных фильтров и ФВЧ порядок фильтра следует назна- чать четным числом. Параметр window позволяет задавать отсчеты окна в векторе-столбце window длины n+1. Если этот параметр не указан, то, по умолчанию, будет ис- пользовано окно Хемминга. Для вычисления окон различного типа в MatLAB предусмотрены следую- щие функции: bartlett(n) - создает вектор-столбец из n элементов окна Бартлетта; blackman(n) - создает вектор-столбец из n элементов окна Блэкмана; boxcar(n) - создает вектор-столбец из n элементов прямоугольного окна; chebwin(n,r) - создает вектор-столбец из n элементов окна Чебышева, где r - желаемый уровень допустимых пульсаций в полосе задерживания в децибелах; hamming(n) - создает вектор-столбец из n элементов окна Хэмминга; hanning(n) - создает вектор-столбец из n элементов окна Хэннинга; 5.4. Проектирование фильтров 238 kaizer(n, beta) - создает вектор-столбец из n элементов окна Кайзера, где параметр beta определяет затухание боковых лепестков преобразования Фурье окна; triang(n) - создает вектор-столбец из n элементов треугольного окна. Приведем пример. Произведем расчет полосового КИХ-фильтра 24 порядка с полосой пропускания 65 0 / 35 0 ≤ ≤ N ω ω : b = fir1(48,[0.35 0.65]); freqz(b,1,512) Результат показан на рис.5.48. Рис. 5.48 Группа процедур fir2 служит для расчета коэффициентов цифрового КИХ- фильтра с произвольной амплитудно-частотной характеристикой, задаваемой век- торами f частот и m - соответствующих желаемых значений АЧХ. Общий вид об- ращения к процедуре таков: b = fir2(n, f, m, npt,lap, window). Вектор f должен содержать значения нормализованной частоты в неубы- вающем порядке от 0 до 1. Вектор m должен быть той же длины, что и вектор f, и содержать желаемые значения АЧХ на соответствующих частотах. Параметр npt позволяет задать число точек, по которым выполняется ин- терполяция АЧХ. Параметр lap определяет размер (число точек) области около точек скачкообразного изменения АЧХ, в которой выполняется сглаживание. Ес- ли эти параметры не указаны, то, по умолчанию, принимается npt = 512 и lap = 25. Рассчитаем двухполосный фильтр 30-го порядка: f = [0 0.2 0.2 0.6 0.6 0.8 0.8 1]; m = [ 1 1 0 0 0.5 0.5 0 0]; b = fir2(30,f,m); 5.4. Проектирование фильтров 239 [h,w] = freqz(b,1,512); plot(f,m,w/pi,abs(h)), grid title('АЧХ КИХ-фильтра (процедура FIR2)') xlabel('Нормализованная частота'), ylabel('А Ч Х') На рис. 5.49 приведены желаемая АЧХ и АЧХ, полученная в результате расчета. Рис. 5.49 Следующая процедура - fircls - также рассчитывает многополосный фильтр, но в несколько другой форме - путем задания кусочно-постоянной желаемой АЧХ. Формат обращения к ней таков b = fircls(n, f, amp, up, lo, 'design_flag'). Здесь f, как и ранее, вектор значений нормализованных частот (от 0 до 1), определяющих границы полос фильтра. Вектор amp определяет кусочно- посто- янную желаемую АЧХ фильтра, количество его элементов равно числу полос фильтра и, следовательно, на 1 меньше числа элементов вектора f. Векторы up и lo определяют соответственно верхние и нижние допустимые отклонения АЧХ спроектированного фильтра от желаемой для каждой из полос. Размер их совпа- дает с размером вектора amp. Параметр 'design_flag' может принимать три значения: trace - для обеспечения вывода результатов в виде текстовой таблицы; plots - для графического отображения АЧХ, групповой задержки, нулей и полюсов; both - для отображения результатов как в текстовой, так и в графической форме. Приведем пример разработки прежнего двухполосного фильтра: n= 30; f = [0 0.2 0.6 0.8 1]; amp = [1 0 0.5 0]; up = [1.02 0.02 0.51 0.02]; lo = [0.98 -0.02 0.49 -0.02 ]; b = fircls(n,f,amp,up,lo,'both') Результат приведен ниже и на рис. 5.50. 5.4. Проектирование фильтров 240 Bound Violation = 0.0755112846369 Bound Violation = 0.0116144793011 Bound Violation = 0.0004154355279 Bound Violation = 0.0000905996658 Bound Violation = 0.0000214272508 Bound Violation = 0.0000009624286 Bound Violation = 0.0000002393147 Bound Violation = 0.0000000596813 Bound Violation = 0.0000000146532 Bound Violation = 0.0000000036610 b = Columns 1 through 7 -0.0001 -0.0031 0.0226 0.0101 0.0060 0.0011 -0.0105 Columns 8 through 14 -0.0231 -0.0626 0.0090 -0.0001 -0.0145 0.1775 0.1194 Columns 15 through 21 0.1272 0.3023 0.1272 0.1194 0.1775 -0.0145 -0.0001 Columns 22 through 28 0.0090 -0.0626 -0.0231 -0.0105 0.0011 0.0060 0.0101 Columns 29 through 31 0.226 -0.0031 -0.0001 Рис. 5.50 Для сравнения с результатами работы процедуры fir2 построим график по- лученной АЧХ, аналогичный приведенному на рис. 5.49: [h,w] = freqz(b,1,512); plot(w/pi,abs(h)), grid title('АЧХ КИХ-фильтра (процедура FIRСLS)') xlabel('Нормализованная частота'), ylabel('А Ч Х') Результат представлен на рис. 5.51. 5.4. Проектирование фильтров 241 Рис. 5.51 Процедура fircls1 предназначается для расчета параметров ФНЧ и ФВЧ с КИХ методом наименьших квадратов с учетом допусков на отклонения АЧХ. Предусмотрены следующие виды обращения к этой процедуре: b =fircls1(n, Wo, dp, ds) b =fircls1(n, Wo, dp, ds, 'high') b =fircls1(n, Wo, dp, ds, Wt) b =fircls1(n, Wo, dp, ds, Wt, 'high') b =fircls1(n, Wo, dp, ds, Wp, Ws, k) b =fircls1(n, Wo, dp, ds, Wp, Ws, k, 'high') b =fircls1(n, Wo, dp, ds, . . ., 'design_flag'). Параметр Wo представляет собой нормализованную частоту среза; dp оп- ределяет максимально допустимое отклонение АЧХ рассчитанного фильтра от 1 в полосе пропускания, а ds - максимальное отклонение АЧХ рассчитанного фильтра от 0 в полосе задерживания. Наличие флажка 'high' определяет, что рассчитываются параметры ФВЧ. Если этот флажок отсутствует, рассчитывается ФНЧ. Указание параметра Wt позволяет задать частоту Wt, выше которой при Wt > Wo или ниже которой при Wt < Wo гарантируется выполнение требований к АЧХ синтезируемого фильтра. Параметры Wp, Ws и k позволяют соответственно задать граничную часто- ту пропускания, граничную частоту задерживания и отношение ошибки в полосе пропускания к ошибке в полосе задерживания. Флаг 'design_flag' имеет тот же смысл и принимает те же значения, что и у предыдущей процедуры. Группа процедур remez осуществляет расчет коэффициентов цифрового КИХ-фильтра с линейной ФЧХ по алгоритму Паркса-МакКлелла, в котором ис- пользован обменный алгоритм Ремеза и метод аппроксимации Чебышева. При 5.4. Проектирование фильтров 242 этом минимизируется максимальное отклонение АЧХ спроектированного фильтра от желаемой АЧХ. Приведем наиболее полный вид обращения к процедуре: b = remez(n, f, а, W, 'ftype'). Вектор f должен состоять из последовательных, в возрастающем порядке записанных пар нормализованных (от 0 до 1) частот, определяющих соответст- венно нижнюю и верхнюю границы диапазона полосы пропускания или задержи- вания. Вектор 'a' должен содержать желаемые значения АЧХ на частотах, опре- деляемых соответствующими элементами вектора 'f'. Желаемая АЧХ в полосе частот от f(k) до f(k+1) при k нечетном представляет собой отрезок прямой от точки f(k), a(k) до точки f(k+1), a(k+1). В диапазонах от f(k) до f(k+1) при k - чет- ном значение желаемой АЧХ не определено (а, значит, при проектировании фильтра АЧХ в этих диапазонах может принимать любое значение). Следует за- метить, что f(1) должно всегда быть равным 0. Векторы f и a должны быть одина- ковой длины, причем общее количество элементов каждого вектора должно быть четным числом. Вектор W задает значения коэффициентов веса каждого из полос АЧХ, за- данных парами частот вектора f. Эти коэффициенты используются при аппрокси- мации АЧХ и определяют достигаемое при аппроксимации соотношение между реальным и желаемым значением АЧХ в каждом из диапазонов. Число элементов вектора W равно половине числа элементов вектора f. Флаг 'ftype' может принимать одно из двух значений: 'hilbert' - в этом случае процедура проектирует фильтры с нечетной сим- метрией и линейной фазой; 'differentiator' - синтезируется фильтр с использованием специальных ме- тодов взвешивания; при этом для ошибок задаются веса, пропорциональ- ные 1/f; поэтому ошибки аппроксимации на низких частотах меньше, чем на высоких; для дифференциаторов, АЧХ которых пропорциональна час- тоте, минимизируется максимальная относительная ошибка. Ниже приводится пример проектирования полосового фильтра 17-го поряд- ка: f = [0.0.3 0.4 0.6 0.7 1]; a = [0 0 1 1 0 0]; b = remez(17,f,a); [h, w] = freqz(b, 1, 512); plot(f,a,w/pi,abs(h)), grid title('АЧХ КИХ-фильтра (процедура REMEZ)') xlabel('Нормализованная частота'), ylabel('А Ч Х') Результат приведен на рис.5.52. 5.4. Проектирование фильтров 243 Рис. 5.52 Особенностью следующей процедуры сremez является то, что исходные данные по желаемой форме АЧХ фильтра задаются в виде функции, условно обо- значенной fresp. Формы обращения к этой процедуре приведены ниже: b = cremez(n, f, 'fresp') b = cremez(n, f, 'fresp', w) b = cremez(n, f, {'fresp', p1,p2,...},w) b = cremez(n, f, a, w) b = cremez(..., 'sym') b = cremez(..., 'debug') b = cremez(..., 'skip_stage2') [b, delta, opt] = cremez(...). Параметры n, f имеют тот же смысл и требования к их представлению такие же, как и при применении процедуры remez. В отличие от последней, вектор зна- чений желаемой АЧХ, соответствующих заданным значениям вектора f, опреде- ляется путем обращения к функции fresp. Функция fresp может принимать одно из следующих: значений lowpass, highpass, bandpass, bandstop (ФНЧ, ФВЧ, полосовой и режектор- ный фильтры); при этом рассчитываются параметры указанного типа фильтра; если к функции 'fresp' не указаны дополнительные параметры (обращения к процедуре первого и второго видов), то групповое время замедления (ГВЗ) принимается равным n/2; в случае же обращения к про- цедуре в третьей форме, где в качестве дополнительного параметра функции fresp указан один - d, ГВЗ = n/2+d; multiband (многополосовой фильтр); синтезируется фильтр, заданный вектором 'a' желаемой АЧХ при значениях частот, определенных векто- ром 'f'; при этом вектор 'a' указывается в качестве первого дополнитель- ного параметра к функции multiband (третья форма обращения); если, кроме этого вектора, не указаны другие дополнительные параметры, то 5.4. Проектирование фильтров 244 ГВЗ принимается равным n/2, если же указан еще один дополнительный параметр d, то ГВЗ = n/2+d; differentiator (дифференциатор); эта функция позволяет рассчитывать ко- эффициенты дифференцирующего фильтра с линейной фазой; при обра- щении к этой функции в качестве дополнительного параметра необходи- мо указать частоту дискретизации Fs; по умолчанию Fs=1; hilbfilt (фильтр Гильберта); в этом случае находятся коэффициенты фильтра Гильберта с линейной фазой. Обращение к процедуре четвертого вида эквивалентно обращению b = cremez(n, f, {'multiband', a},w). Параметр 'sym' позволяет задать тип симметрии импульсной характеристи- ки (ИХ) фильтра. Он может принимать следующие значения: none - в этом случае ИХ может быть произвольной; это значение пара- метра используется по умолчанию, если при определении желаемой АЧХ задаются отрицательные значения частот; even - АЧХ должна быть вещественной с четным типом симметрии; такое значение параметра используется по умолчанию при проектировании ФНЧ, ФВЧ, полосовых и режекторных фильтров; odd - АЧХ должна быть вещественной с нечетным типом симметрии; та- кое значение по умолчанию используется при проектировании фильтров Гильберта и дифференциаторов; real - АЧХ должна иметь сопряженный тип симметрии. Использование флага 'skip_stage2' (см. седьмой вид обращения к процеду- ре) позволяет не выполнять второй этап алгоритма оптимизации, который рассчи- тывает коэффициенты фильтра в тех случаях, когда этого нельзя сделать с помо- щью алгоритма Ремеза. Исключение второго этапа сокращает время расчетов, но может повлечь снижение точности. По умолчанию выполняются оба этапа опти- мизации. Параметр 'debug' (см. шестой вид вызова процедуры) определяет вид выводимых на екран результатов расчета фильтра и может принимать следующие значения: 'trace', 'plots', 'both' и 'off'. По умолчанию используется 'off' (т. е. на эк- ран не выводится информация). Использование дополнительного выходного параметра delta (см. восьмой вид обращения к процедуре) дает возможность использовать в дальнейших опе- рациях значение максимальной амплитуды пульсаций АЧХ. Выходной параметр opt содержит набор дополнительных характеристик: opt.grid - вектор отсчетов частоты, использованных при оптимизации; opt.H - вектор значений АЧХ, соответствующих значениям элементов в векторе opt.grid; |