Средства пакета расширения Calculus системы
Mathematica 4 содержат функции для вычисления
тригонометрического ряда Фурье для
периодических функций, представленных на
отрезке времени от -T/2 до +T/2. Это крайне неудобно
для практических вычислений, когда заданная
функция y(t) определена на интервале [0,T], причем
она может быть не периодической. Кроме того,
нередко такая функция задана таблицей своих
отсчетов. Если таких отсчетов Ni, то при
классическом спектральном анализе, на основе
теоремы отсчетов, для представления функции
можно использовать не более Ni/2 гармоник. Это
резко снижает точность представления y(t)
тригонометрическим рядом Фурье, что, в частности,
связано со значительным проявлением эффекта
Гиббса. Он проявляется в виде колебаний и
выбросов аппроксимирующей функции, амплитуда
которых достигает 18% от амплитуды сигнала.
Предлагаемый ниже метод обеспечивает
высокоточное представление сложных сигналов,
заданных в табличной форме и допускающих
линейную (или иную) аппроксимацию. Суть метода
очень проста - сигнал y(t) заменяется его
непрерывной на отрезке [0,T] аппроксимацией - в
частности кусочно линейной. Это возмоляет найти
произвольно большое число гармоник, что повышает
точность спектрального представления y(t) и
снижает проявление эффекта Гиббса,
обусловленного малым конечным числом гармоник.
Основы метода описаны в книгах:
Дьяконов В. П. Справочник по расчетам на
микрокалькуляторах. Изд. 3-е. М.: Наука. Физматлит.
Дьяконов В. П. Справочник по адгоритмам и
программам на языке Бейсик для персональных ЭВМ.
М.: Наука. Физматлит.- 1987.
Дьяконов В. П., Абраменкова И. В. Mathcad 7 в
математике, в физике и в Internet. М.: Нолидж.- 1999.
Пусть сигнал представлен вектором его
значений:
Найдем число элементов вектора Yi:
Выполним линейную интерполяцию без
трансформации масштаба времени:
График зависимости yy(i+1) при изменении
i от 0 до (Ni-1) представлен ниже:
Зададим интервал времени на котором
определен сигнал T и нужное нам число гармоник M:
Еще раз отметим, что при данном методе
число гармоник может значительно превышать N/2.
Выполним кусочную линейную интерполяцию сигнала
с преведением интервала интерполяции к
интервалу времени [0,T]:
Если нужна разрывно-ступенчатая
аппроксимация y[t] переформатируйте следущую
ячейку под входную.
y[t_]:=If[t<-dti+T,yy[1-dti+Floor[t/dti]],0]
Специфика интерполяции связана стем,
что функция Interpolation дает интерполяцию на отрезке
[1,Ni]. Введенная функция пользователя y[t]
определяет сигнал на отрезке [0,T]. Ниже
представлен график сигнала на этом отрезке.
Когда проивольная функция задана лишь
своими отсчетами, расчет косинусных и синусных
коэффиуиентовпрямо по их интегральным
представлениям нерационален, поскольку не дает
никаких преимуществ перед вычислением
коэффициентов по формуле интегрирования методом
прямоугольников. Ниже представлен такой расчет,
причем для ускорения вычислений в общем случае
нулевые отсчеты функции (если таковые есть) не
обрабатываются:
По найденным косинусным и синусным
коэффициентам ряда Фурье можно найти амплитуды
Mg[k] и фазы Phg[k] гармоник:
Показатель степени m позволяет
реализовать три метода расчета амплитуды
гармоник:
m=0 - приближенный (в этом случае
целесообразно отказаться от возводимого в
степень множителя);
m=1 - точный метод приямоугольников;
m=2 - точный метод трапеций (для кусочно-линейного
представления y(t) кусочно-линейная).
Построим амплитудно-частотную и
фазо-частотную характеристики спектра:
Зададим функцию пользователя,
вычисляющую мгновенное значение сигнала по
заданным M гармоникам:
Нетрудно получить аналитическое
выражение для ряда Фурье зависимости F(t):
Для ускорения
вычислений при вычислении F(t) исключаются
гармоники с нулевой амплитудой, если таковые
есть. Теперь можно построить график
синтезированной по гармоникам функции.
Нетрудно заметить, что при числе
гармоник порядка 40-50 график почти повторяет
исходную функцию при очень малой волнистости,
обусловленной эффектом Гиббса.
О том, насколько точно приближение
исходной функции тригонометрическим рядом для M
гармоник позволяет судить совмещенный график
этих функций:
Зададим амплитудно-частотную АЧХ и
фазо-частотную ФЧХ характеристику устройства
(усилителя, фильтра и т.д.) искажающего сигнал.
Построим график АЧХ (верхняя линия) и
ФЧХ (нижняя линия):
Запишем выражение для гармонического
синтеза сигнала, учитывая его искажения
(амплитуда каждой гармоники умножается на
значение Af на ее частоте k*f1, а к фазовому сдвигу
гармоники прибавляется фазовый сдвиг Pf на
частоте гармоники).
Теперь можно построить
исходную временную зависимость сигнала и
временнуб зависимость сигнала после прохождения
им искажающей цепи.
Эксперименты
1. Чтобы судить о том, насколько
описанный выше метод эффективен, побробуйте
вначале сделать число гармоник M равным
предельно допустимому по теореме отсчетов - в
нашем случае M=10 и запустите данной документ с
начала - командой Kernel-Evalution-Evalution Notebook. Затем
сделайте M=20 и 100.
2. Опробуйте данный метод используя
глобальную полиномиальную аппроксимацию
функции y(t). Для этого измените опцию в функции
Interpolation InterpolationOrder®1 на InterpolationOrder®5. Это означает
аппроксимацию полиномом пятой степени для всей
зависимости y(t). Убедитель в том, что результат
приближения в этом случае заметно хуже, из-за
погрешности полиномиальной аппроксимации
(однако для гладких функций он может оказаться
весьма неплохим - проверьте это задав вектор Yi
для более гладкой функции).
3. Опробуйте данный метод для
ступенчатого сигнала (с разрывами). Для этого переформатируйте ячейку с y(t) для
ступенчато-разрывной аппроксимации и заново
пустите документ.
4. Исследуйте, как влияет показатель m,
задающий метод расчета амплитуд гармоник, на
точность представления сигналов той или иной
формы, в частности насколько эффективно
уточнение расчета для подавления эффекта Гиббса.
5. Измените АЧХ и ФЧХ фильтра и оцените
степень искажений сигнала при этом.
6. Опробуйте данный "ноутбук" для сигналов
другой формы и оцените эффективность применения
спектрального подхода при их анализе и синтезе.
|