Преобразование фурье

Цитата Alexei, 11.10.02, 18:31:51

Уважаемые, а как опредилить основную гармонику?

Не знаю что такое "основная гармоника", но если это то же самое что доминантная частота, то есть способ,который считает ее как положение "центра тяжести спектра".
Я когда-то считал.Получалось очень неплохо.Довольно точно.
Если нужно,пришлю исходник с-шный по почте.

Я несколько дней назад выводил дискретное ПФ из обычного.
И могу помочь в вашем вопросе.

Дискретное преобразование Фурье - это
    k-1    
g(u)=SUM (f(n)*cos(2*pi*n*u/k)-j*f(n)*sin(2*pi*n*u/k));
    n=0

k - число квантов и в сигнале и в спектре, оно одинаковое.

И максимум она достигает, когда подынтегральная функция имеет максимальный модуль, а её знак совпадает со знаком косинуса или синуса(смотря под что максимизируем, тут же и определим максимизованную частоту).
Косинус на его знак на максимальный модуль равно максимум умножить на косинус квадрат.
А интеграл по периоду от косинуса квадрат(с любой частотой) равно 4
В результате сумма получается равна 4*max*k

Этот результат сразу делится на k, но это уже на практике, а не по формуле

получается, что спектр превышает сигнал не более чем в четыре раза, да и этого значения она достигает в редких случаях.

Спасибо, наконец-то нашелся умный человек. А то тут  не форум, а базар  кокай-то был.

Извините, но в прошлом сообщении я допустил ошибку.
Дискретное преобразование, как я говорил, - это
      k-1
g(u)=SUM cos(2*pi*u*n/k)*f(n)-j*sin(2*pi*u*n/k)*f(n)
      n=k

сразу хочу привести формулу обратного преобразования фурье

              k-1
f(n)=1/k * SUM cos(2*pi*u*n/k)*g(u)+j*sin(2*pi*u*n/k)*g(u)
              u=k
(разница в замене n<->u g<->f, в изменении знака перед синусом и в ведении коэффициента)


Первая ошибка - в том, что чатоты равнозначны. Если частота равна половине от числа квантов, то значение в косинуса будут принимать не плавные значения, а лишь +- 1. Теперь, если сопоставить функцию принимающую значение +- u (где u-максимальное значение(примем за еденицу)) так, что в сумме умножая на синус всегда получается u. Тогда вся сумма будет равна k*u
Если затем спектр поделить на k, то получится u. Т.е. тогда размерность спектра будет такая же РАДУЙТЕСЬ!

Только тогда не сработает обратное ДПФ, приведённое выше, так как деление уже произведено.

Вторая ошибка была - в цифре 4. Она возникла из-за того, что я неправильно перешёл от интеграла к сумме - надо было ещё поделить на 2pi. Получилось бы вместо 4-х 0,63661977. Это - для первой гармоники.

Но так как я исправил и первую ошибку, и считал для k/2-й гармоники, то получилось ровно 1.

Я проверил в МАТЛАБе:
» a=0:999;                      -нумеруем кванты
» s=sign(cos(2*pi*500*a/1000)); -содержит +1 -1 +1 -1 +1 -1  и так 1000 квантов
» plot(a,s);                    -на графике сплошная стена - слишком частый меандр
» f=fft(s,1000);                -делаем ДПФ(в реализация - быстрое ПФ)
» plot(a,real(f));              -на графике все кванты примерно 0, кроме пятисотого, он - 1000
» hold                          -здесь же строим мнимую часть
Current plot held                  
» plot(a,imag(f),'-r');         -Красным цветом, чтоб отличить. Получается примерно ноль везде.
»        

Вся "энергия" меандра( sign(cos(t)) ) уровня U превращается в один импульс размером в число квантов умножить на U. Так что, разделив на число квантов, получим U

И последнее: если делаешь преобразование сам, то сохраняй битов сколько хочешь - всё равно информация теряется, так что можешь сохранять хоть больше, хоть, если припрёт, меньше.

Ещё совет:
Предельный спектр( соответствие частоте значения ПФ при сигнале меандр с этой же частотой
c(i+1)=SUM (cos(2*pi*i*(0:999)/1000)*sign(cos(2*pi*i*(0:999)/1000)))
)
который я считал для кванта 1 (=2/pi)и кванта 500 (=1)
кроме 0-го и k/2-квантов не превышает 0,65
На этом можно или сэкономить место, или сохранить больше битов.

У меня вопрос :
уважаемые кто нить занимался ФП под дельфи ? только не дискретным а полным (N) ?
если да поделитесь пожалуста алгоритмами и насколько медленнее БПФ это работает ?
есть ли смысл вообще делать ПФП или нет ?

Обычно преобразование - N^2 операций
БПФ - N*log2(N)
А какая разница где в Дельфи или где-то еще?
Алгоритм везде одинаков, синтаксис написания разный.

а чо есть бпф? или я что-либо пропустил в вашем обсуждении?

в смысле в чём заключается алгоритм, делающий пф за столь хорошее время?

Простите я наверное немного невнимательно прочитал вопрос от GEO.
Все написаное мною выше про число операций, касается ДПФ(Дискретного Преобразования Фурье), а БПФ есть ничто иное как алгоритм вычисления ДПФ используя меньшее число операций.
Смысл примерно такой БПФ: если полностью расписать 8-ми точечное ДПФ, то среди слагаемых будет ab+ac (3 операции), что можно вычислить a(b+c) (2 операции).
У меня сейчас нет под рукой конспекта лекций, но я думаю что это есть в Инете.
Вот собственно здесь:
http://www.dsp.sut.ru/book/lections/l4/l4_3.htm

2Cubloid

там слишком мало воды, непонятно, что какая переменная означает

Что я вам могу сказать на это:
Там все ясно по теории, но нет конечной реализации.
Только надо взять бумагу и ручку или пойти на Я.ру
Последний пункт я выбрал и на второй странице поиска нашел вот это.
http://home.tula.net/avm/fft/index.html
И там же ссылка на www.fftw.org
)) На последнем сайте представлена библиотека, которую можно использовать и большой набор ссылок на ресурсы по теме БПФ.
Что еще вам надо по этому вопросу? Все есть в Инете.

2Cubloid

ура! я понял! там не N log N умножений, а те же самые N^2, просто их в 2 раза меньше. ради этого и не стоило огород городить, в 2 раза я и сам могу уменьшить.....

8) Иди ты на хуй, придурок. Что значит N^2, только в два раза меньше. Ты долбоеб, не видишь разницы между NlogN и N^2 ????? Тогда вали отсюда, сука и больше здесь  не появляйся!!!!!! Нехуй людям попусту  голову морочить, пиздюк!

2melan

я так тоже могу :

2cubloid

я тут на досуге подумал и пришёл к следующим выводам:

1. твоя лекция существует для того, чтобы пугать людей. я из неё ничего не понял, однако, придя домой, я придумал подобный алгоритм, исходя из "свойства симметрии и свойства периодичности". на самом деле он не такой уж и страшный.

2. там в натуре n log n умножений, однако там такое же количество сложений, даже больше. а вот у тебя в лекции приведён пример, который со всей присущей ему наглядностью показывает, как замечательно этот алгоритм делает 500 000 умножений вместо миллиона, что, конечно же, неправда.

3. всё это чушь.

как видно из рисунка, если мы умножим am или am+1 на тот же коэффициент k, что и an, то мы получим неверный результат.

Сообщение отредактировано: wormball - 25.03.04, 15:51

0. Наиболее авторитетный, ИМХО, источник - Рабинер и Голд. Выложил его на dsp-book.narod.ru
1. При выполнении Фурье по-разному решают вопрос нормирования. Математики чаще делят на нормирующий множитель и прямое, и обратное преобразования, а при программировании зачастую, экономя операции, нормируют (деля уже не на sqrt(n), а на n) только обратное преобразование. Но даже при математическом варианте - возможно драматическое возрастание элементов Фурье-образа по сравнению с исходными. В крайнем случае, согласно теореме Парсеваля, возможно возрастание в sqrt(n/2) раз, т.е. для того, чтобы пользоваться только 16-битной арифметикой при векторе из 512 элементов, исходный вектор должен быть не более чем 12-разрядный, а если брать "программистскую" нормировку - то 8-разрядный. Поэтому используется 32-битная или же плавающая (кстати, нынче это может быть и дешевле фиксированной, но с дополнительными операциями по нормировке на каждом шаге) арифметика.
2. Вычислять ДПФ "в лоб" стоит только при ОЧЕНЬ малых размерностях исходного вектора. В точной арифметике они дают одинаковый результат, а в реальной - БПФ точнее. Что до числа операций - "в лоб" надо O(n^2), а через БПФ не O(n^2/2) - понимающие поймут шутку - а O(n log n). При длине вектора порядка 100 и более - выигрыш в разы и десятки раз. (Лично я, искусно применив БПФ, смог ускорить вычисление корреляции в 300 раз...) Но, скажем, в JPEG считают прямо - поскольку там матрица 8х8 и вектора, соответственно, по 8 элементов.
3. Во многих случаях можно обойтись без Фурье. Скажем, если нужно установить наличие в сигнале строго определенных частот, поможет алгоритм Герцеля, если мощность в полосе - фильтрация, спектры могут вычисляться через авторегрессию или методом Писаренко и т.п.
4. В качестве доминирующей частоты сигнала можно брать максимальный пик в спектре, среднюю частоту, медианную частоту. Зависит от конкретной задачи.