![]() |
|
![]() ![]() ![]() |
|
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
Работаем с модулем FFT !!! Можно 2 раза - т.к. одну точку ты уже нашёл через FFT. Если алгоритм расчёта параболы рассчитан на равномерную сетку, тогда три раза интеграл считать. Но это лучше перебора с шагом 0.5. И лучше с 44100 работать. На увеличение быстродействия особо не надейся ![]() Прикинь сколько вообще ты можешь итераций себе позволить. Учти, что на коротком сигнале ты не получишь хорошую точность! И ничего не поможет! Я думаю точность 0.1 Гц для 6 струны получишь при длительности сигнала 1 секунда. |
|||
|
||||
Pavia |
|
|||
Опытный ![]() ![]() Профиль Группа: Участник Сообщений: 418 Регистрация: 6.12.2008 Репутация: 11 Всего: 12 |
Самый простой способ увеличить в 10-20 раз это взять библиотеку FFTW, она быстрее обычной самописной реализации.
Берем 2 буфера 1) N=2048 отсчётов 2) NewN=64*1024 отсчётов В первый пишем сигнал пишем сигнал потом умножаем на оконную функцию лучше Ланцоша (sinc) Копируем результат во второй остававшийся часть заполняем нулями. После FFT и ищем максимум. i*Fd/(NewN); частоту определяем как. Точность будет 0.1 Гц. Больше ничего не требуется. Без FFTW по скорости не отличается от предыдущего алгоритма.
|
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
Сейчас проверил ещё один интересный способ уточнения частоты.
Алгоритм такой: 1. FFT и находим основную частоту 2. Определяем сигнал с линейно меняющейся частотой Fmax +- 0.5 Гц (44100) 3. FFT этого сигнала 4. Перемножаем вектора 5. Реверсное FFT--> X(i) 5. Отношение Х(N/2)/X(1) даст поправку к определённой в п.1 частоте Проверил на тестовом сигнале. Точность 0.001 Гц !!! (Fd=44100) |
|||
|
||||
Kpeved |
|
|||
Новичок Профиль Группа: Участник Сообщений: 42 Регистрация: 26.2.2012 Репутация: нет Всего: нет |
Santik, а как по производительности? Сколько итераций это все приблизительно займет ? И что значит "Сигнал с линейно меняющейся частотой"?
|
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
А итераций нет совсем. Грубо говоря - всё забито в прямое и обратное преобразование Фурье.
Алгоритм я уже приводил на Фортране. Там нужно только ввести сигнал с переменной частотой: f=f0+(0.5)/T*t Все примеры для 44100, Т=1сек. Пусть мы из FFT определили fmax= 328 Гц . Если значение при 329 больше чем в 327 то нам надо двигаться в +. f(t)=328+(0.5)/T*t , т.е. частота меняется линейно во времени от 328 до 328.5 Гц Интеграл от 2*Pi*f(t) даст фазу сигнала Ф(t). cos(Ф(t)) +i sin(Ф(t)) - это и будет наш сигнал с линейно меняющейся частотой. Далее делаем свёртку этого сигнала с исходным сигналом (см. программу) (умножаем вектора результатов FFT и делаем обратное FFT Максимум модуля свёртки (или точная частота) гарантированно будет в этом диапазоне. Завтра приведу данные по производительности. Программу полноценную мне конечно лень писать. Но для примера текст отладочной программы приведу, чтобы можно разобраться было. Ну а дальше - сам. ![]() |
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
Данные по быстродействию:
Fd=44100 N=44100 (1 сек) Тестовый сигнал 329.123456 Гц с шумом. 600 значений частот рассчитаны за 53 секунды. (P-4 3ГГц) Погрешность 0.001 Гц. |
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
|
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
||||
|
||||
Kpeved |
|
||||
Новичок Профиль Группа: Участник Сообщений: 42 Регистрация: 26.2.2012 Репутация: нет Всего: нет |
Наконецто набрася храбрости прочитать код на ещё одном неизвестном мне языке .) По коду мне все ещё непонятно следующее :
Ещё не знаю - есть ли в фортране перегрузка операторов - но как именно он умножает операнды ? Да, и что делает DCONJG тоже не ясно )
Вроде бы по коду все . Ещё по быстродействию - Если мы будем использовать ЧД 8К при 1024 значениях - будет ли у нас нужная точность ? Получается что на 1 отсчет придется по 8Гц и +- надо брать намного больше . |
||||
|
|||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
Даже это уже положительный момент!!! Да и что может быть проще Фортрана? Приивлекает богатейшая библиотека IMSL. БПФ, к примеру, - с любым количеством точек и очень быстро! DCONJG -функция комплексного сопряжения (комплексного вектора FX - в данном случае в FX - БПФ исходного сигнала ) CDABS - -функция возвращае модуль комплексного числа (D- двойная точность) По быстродействию сказать ничего не могу. При изменении N "магические числа" могут быть другими! По точности см. http://dsplib.ru/forum/viewtopic.php?f=7&a...42&start=10 Выбираешь нужную погрешность и находишь М и потом N - число точек (для 6 струны - обязательно, т.к. НЧ хуже по точности) Я пока этот алгоритм не "обкатыва". Надеюсь совместными усилиями .... Это сообщение отредактировал(а) Santik - 26.3.2012, 07:39 |
|||
|
||||
Kpeved |
|
|||
Новичок Профиль Группа: Участник Сообщений: 42 Регистрация: 26.2.2012 Репутация: нет Всего: нет |
Santik, ну а что вообще за числа в 73 ? Как их рассчитывать ?
Как я понял М и N рассчитываются по последней формуле ? Тогда если Т=0.128 ,f= 70Гц, то M =8.96. Погрешность около 0.1 Гц . И ??) |
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
Определение коэффициентов: Снимаем характеристику er=CDABS(X1(N/2))/CDABS(X1(1)) при расстройке частоты от 0 до 0.5 с шагом 0.01 Гц. Находим ошибку. Находим тренд функции Ошибка(er). (я для начала в Excel делал) Получаем df=-0.861*er+0.85504 По полученному значению er находим ошибку df. У меня сейчас: df=0.11356244E+01*er**7-0.56608063E+01*er**6+0.11240384E+02*er**5-0.11300708E+02*er**4+ 0.64610102E+01*er**3-0.25383018E+01*er**2 -0.60435858E-01*er+ 0.72327669E+00 ![]() ... И ??) Это говорит о том, что при данной частоте и времени сигнала ты определишь частоту с погрешностью не менее 0.1 Гц . Если тебе нужна меньшая погрешность - придется увеличивать Т и другого выхода нет! |
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
Для Fd=8000 Гц и N=1024
df=-0.69019186E+01*er+ 0.68520498E+01 При частоте 78 Гц метод даёт погрешность 0.2 Гц |
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
В файле некоторые результаты вычисления при Fd=8000 Гц.
Присоединённый файл ( Кол-во скачиваний: 10 ) ![]() |
|||
|
||||
Santik |
|
|||
Новичок Профиль Группа: Участник Сообщений: 45 Регистрация: 13.3.2012 Где: Мирный (Якутия) Репутация: нет Всего: нет |
Исследуем полифонический сигнал (файл 44k-all.txt)!
Это сообщение отредактировал(а) Santik - 29.3.2012, 19:09 Присоединённый файл ( Кол-во скачиваний: 10 ) ![]() |
|||
|
||||
![]() ![]() ![]() |
Правила форума "Алгоритмы" | |
|
Форум "Алгоритмы" предназначен для обсуждения вопросов, связанных только с алгоритмами и структурами данных, без привязки к конкретному языку программирования и/или программному продукту.
Если Вам понравилась атмосфера форума, заходите к нам чаще! С уважением, maxim1000. |
1 Пользователей читают эту тему (1 Гостей и 0 Скрытых Пользователей) | |
0 Пользователей: | |
« Предыдущая тема | Алгоритмы | Следующая тема » |
|
По вопросам размещения рекламы пишите на vladimir(sobaka)vingrad.ru
Отказ от ответственности Powered by Invision Power Board(R) 1.3 © 2003 IPS, Inc. |