Идея в общем виде простая. Решением возвратного уравнения

является сумма

синусоид (со сдвигом по фазе и, вообще говоря, затухающих или нарастающих по экспоненту). Оценив его параметры, можем решить вспомогательное уравнение

, комплексные корни которого дадут нам частоты (как аргумент комплексных корней) и декременты затухания (логарифмы амплитуд).
Разные алгоритмы оценки (авторегрессия, через собственные значения и т.п.) дают разные методы оценки частот (подробнее, например, Марпл-мл., "Цифровой спектральный анализ" и другие книги). Достоинство и недостаток такого подхода состоят в том, что число пиков спектра у нас априори задано порядком модели. Это вносит элемент произвола, но позволяет получить параметры отдельных пиков точнее, чем через БПФ и тому подобное.
В частном случае Вашей задачи затухание отсутствует.
Пусть

Тогда можно найти


и

Что очевидным образом позволяет найти

После этого можно построить

и

, затем регрессию y на

, откуда придти к амплитуде и фазе (или предпочесть оставить коэффициенты при косинусе и синусе, как удобнее).