дискретное преобразование фурье на AVR\Pic24F

[YS]

Так, я чую, что все-таки надо прояснить, что, откуда и зачем получается в ДПФ. А то тут некоторые собираются делать фильтры из преобразования Фурье, не переходя к временной области...

1. На вход ПФ подается некоторая функция времени на некотором интервале.
2. На выходе ПФ получается набор комплексных амплитуд гармонических функций, которые, если сложить их на заданном интервале, дадут исходную функцию времени.
3. Исходная функция обязательно должна быть периодичной. Если она непериодична, она полагается таковой. Отсюда растут ноги элайзинга (aliasing).

Если интересует только один частотный компонент - удобнее использовать алгоритм Гёрцеля. Это специальная реализация ДПФ, дающая комплексную амплитуду одной заданной гармоники.

Если стоит задача именно фильтровать сигнал, то лучше использовать реализации цифровых фильтров. В общем случае построение фильтра на ДПФ нецелесообразно, для этого есть гораздо более эффективные методы.

[стьюдент]

1. На вход ПФ подается некоторая функция времени на некотором интервале.
2. На выходе ПФ получается набор комплексных амплитуд гармонических функций, которые, если сложить их на заданном интервале, дадут исходную функцию времени.
3. Исходная функция обязательно должна быть периодичной. Если она непериодична, она полагается таковой. Отсюда растут ноги элайзинга (aliasing).

Как определить период неизвестного wav сигнала?

Если интересует только один частотный компонент - удобнее использовать алгоритм Гёрцеля. Это специальная реализация ДПФ, дающая комплексную амплитуду одной заданной гармоники.

Если стоит задача именно фильтровать сигнал, то лучше использовать реализации цифровых фильтров. В общем случае построение фильтра на ДПФ нецелесообразно, для этого есть гораздо более эффективные методы.

Да тут бы хотя бы с дпф разобраться.. От ДПФ мне нужен только спектр
Непонятен следующий момент:
«Если мы рассчитаем ДПФ для всех k, мы получим массив из N элементов, которые будут иметь нужные нам значения (амплитуды) в первой половине [0..N/2-1] массива. Вторая половина [N/2..N-1], является зеркальным отражением первой. Это получается за счет математической симметрии ДПФ. Из-за чего рассчитывать вторую половину нет необходимости».
У меня симметрии почему-то не получается

[Gudd-Head]

Как определить период неизвестного wav сигнала?

У него нет и не может быть периода. Мы полагаем его периодическим — его период равен Nfft (число точек БПФ), умноженное на период дискретизации.

[Alexeyslav]

Как определить период неизвестного wav сигнала?

Для ДПФ не существует такого понятия как период, функция оперирует с отсчетами. А собственно сам период отсчетов которые ты подаешь на вход задаются именно тобой посредством задания частоты выборки с которой массив отсчетов сохраняется в памяти и передается на обработку.

для самой функции абсолютно не важно с какой скоростью делалась выборка - что 10 секунд на отсчет что 1милиард отсчетов в секунду - ей все равно. От выбора частоты дискретизации будут зависеть лишь частоты соответствующие результатам разложения. 0-й результат - это постоянный ток, 1-й это f/N, k-й это k*f/N вроде бы так. f - это частота дискретизации.

[стьюдент]

Собрал все свои мысли (файл "Исходные данные.pdf"), если не трудно, поправьте что не так. А что то явно не так, по такой методике спектр получается несимметричным(файл "исправления.jpg").

[Alexeyslav]

Вполне симметричный. рассчитайте так же значение для точки k = 8.
а в середине графика значения порядка 10^-14 - это погрешности FPU.

Я как-то делал решение уравнений с 20-ю неизвестными, там было много вычислений и если брать числа одинарной точности (тип Real) то исходное уравнение по рассчитанному решению сходилось с точностью до одной десятой... так что все что меньше 10^-10 в таких вычислениях можно считать нулем.

[YS]

Как определить период неизвестного wav сигнала?

Мы делаем из сигнала выборку Nfft точек, и кагбэ считаем, что анализируемый сигнал повторяется перед и после выдранного куска. Понятное дело, это приводит к искажениям результата. Иногда это неважно, а когда важно - применяют специальные взвешивающие окна (оконные функции).

Да тут бы хотя бы с дпф разобраться..

С практической точки зрения что алгоритм Гёрцеля, что фильтры гораздо проще ДПФ. А фильтры проще еще и в плане понимания теории.

[стьюдент]

Кварц - 16 МГц, коэф. преддилителя АЦП=128. Преобразование - 13 тактов (второе и последующие преобразования).
Следовательно: 125кГц/13 = 9,615кГц - частота дискретизации.
Кварц - 10 МГц, коэф. преддилителя АЦП=128. 78,125кГц/13=6.009 кГц - частота дискретизации.
В итоге шаг расчета частот в первом случае: 9,615кГц/64=150,23 Гц.
Во втором: 6,009 кГц/64=93,89 Гц.

Есть ли способ уменьшить шаг, кроме как повышать число точек? Есть идея в том, чтобы заставить АЦП выполнять преобразование за 25 тактов, как при первом запуске, вопрос только как это сделать?

[Alexeyslav]

Шаг ты уменьшишь, только уменьшишь и диапазон частот. Если тебе надо выделить исключительно низкие частоты с мелким шагом - уменьшая частоту дискретизации. 25 тактов на преобразование можно получить запуская АЦП вручную! В прерывании по окончанию преобразования например.

[YS]

25 тактов на преобразование можно получить запуская АЦП вручную!

Эм... По-моему, 25 тактов работает только первое преобразование после включения АЦП (ADEN), вне зависимости, запускается он вручную (ADSC) или нет. Т.е., надо каждый раз дергать ADEN. Но это извращение.

Чтобы регулировать шаг, проще действительно запускать руками через ADSC по таймеру.

[Alexeyslav]

А не после переключения канала перед любым одиночным преобразованием?

[YS]

Вроде нет, но это надо глядеть даташит.

Вроде как это зависит от опоры.

UPD: посмотрел.

A normal conversion takes 13 ADC clock cycles. The first conversion after the ADC is switched
on (ADEN in ADCSRA is set) takes 25 ADC clock cycles in order to initialize the analog circuitry.

Так что все же ADEN.

[Alexeyslav]

Ладно, пропускай каждый второй результат, получится 26 тактов.

[YS]

Я не понимаю смысла извращений.

Нам нужен период выборки Tsmp? Ну так настроить АЦП так, чтобы он отрабатывал заведомо быстрее Tsmp - АЦП AVR может выдавать полное разрешение при тактовой до 200 кГц. И запускать его по таймеру, который будет обеспечивать Tsmp, а временем преобразования пренебречь.

[стьюдент]

И запускать его по таймеру, который будет обеспечивать Tsmp, а временем преобразования пренебречь.

Спасибо, мысль понял, но с тем, что временем преобразования можно пренебречь не согласен. Интересуют частоты до 1 кГц. Попробую пропускать каждый второй результат или каждый раз дергать ADEN.

[YS]

Интересуют частоты до 1 кГц.

Дискретизации? Тогда точно можно принебречь.

Пусть кварц 10 МГц. Ставим делитель на 64, получаем тактирование АЦП чуть больше 156 кГц. 13 тактов на измерение - 83 мкс. При этом собственно захват уровня будет занимать 1.5 такта - <10 мкс, что в 100 раз меньше желаемого периода измерений (1 мс).

[стьюдент]

Дискретизации? Тогда точно можно принебречь.

Анализ сигналов, которые имеют частоту <1 кГц. Вариант с таймером, кажется более замороченным, чем вариант с пропуском каждого второго измерения

[YS]

Анализ сигналов, которые имеют частоту <1 кГц.

А, это. Тогда по теореме Котельникова частота дискретизации должна быть не менее 2 кГц, что не меняет дела. Ну будет выборка в 50 раз короче, а не в 100...

Вариант с таймером, кажется более замороченным, чем вариант с пропуском каждого второго измерения

Вот как раз с пропуском замороченно.

С таймером система получается гораздо прозрачнее и гибче.

1. Настраиваем АЦП, включаем прерывание по окончанию преобразования.
2. Настраиваем таймер, конфигурируем прерывание по совпадению.

В прерывании таймера запускаем АЦП, в прервании АЦП вынимаем данные и складываем в буфер. Буфер заполнился - ставим флаг готовности, заодно запрещающий запись в буфер. Ну а в главном алгоритме по флагу вызываем преобразование Фурье, или чего там надо... Готовы, надо еще данных - в основном алгоритме сбрасываем флаг и ждем, пока он опять не будет взведен.

Можно завести два буфера, и пока один заполняется, манипулировать с другим...

В такой конфигурации период дискретизации меняется одним значением - содержимым регистра совпадения таймера. Причем меняется гораздо плавнее. Про остальные настройки можно забыть.

[Alexeyslav]

Это пример того как НАДО делать, а не как проще. Проще пропустить одно измерение чем городить огород с прерываниями и синхронизацией обмена между прерыванием и основной программой. Особенно, когда не нужна поточная обработка данных(по всей видимости тут она не нужна).

Оно-то правильно так делать, но мест где можно будет накосячить в 20 раз больше. После отработки линейного алгоритма при достатке свободного времени можно будет заморочится и с правильным методом... если это будет востребовано.

Я бы выбирал максимально простой алгоритм подходящий под условия задачи.

[YS]

Все-таки я не понимаю, почему всем тут кажется, что с таймером сложнее. Как по мне, оба варианта - практически "Hello world". Только второй еще и правильный/расширяемый.

Оно-то правильно так делать, но мест где можно будет накосячить в 20 раз больше.

Э-э-э... Например, где? Я теряюсь в догадках.

Вообще, это, похоже, какая-то инерция мышления программистов-прикладников, пришедших с PC. Мне доводилось натурально переучивать под эмбед пару товарищей, которые неплохо писали под ПК, но под МК рожали такие шедевры, что смотреть было страшно.

городить огород с прерываниями и синхронизацией обмена между прерыванием и основной программой

Где там огород? Одна переменная, в одном месте установили бит, в другом проверили/очистили. При пропуске тоже нужна будет переменная и ее тоже надо будет проверять.

UPD:

Пример с таймером.

Код: Выделить всё

#include <avr/io.h>
#include <avr/interrupt.h>
#include <stdint.h>

#define TIMER_PERIOD		0x04E2
#define BUFFER_RDY			0x01
#define BUF_SIZE			100

volatile uint8_t status=0;
volatile uint16_t buffer[BUF_SIZE],buf_ptr=0;

ISR (ADC_vect)
{
  if (!(status & BUFFER_RDY))
  {
    buffer[buf_ptr]=ADCL | (ADCH << 8);
	buf_ptr++;

	if (buf_ptr>BUF_SIZE)
	{
	  buf_ptr=0;
	  status|=BUFFER_RDY;
	}
  }
}

ISR (TIMER1_COMPB_vect)
{
  TCNT1=0;
  ADCSRA|=(1<<ADSC);
}

void main(void)
{
  //***ADC configuration***

  //REF=AVcc, input=ADC0
  ADMUX=(1<<REFS0);
  //ADC on, ADC int. on, prescaler=64
  ADCSRA=(1<<ADEN) | (1<<ADIE) | (1<<ADPS2) | (1<<ADPS1);

  //***T/C1 ch. B configuration***

  //Ch. B int enabled
  TIMSK=(1<<OCIE1B);
  //1kHz interrupt @ 10MHz CPU clk
  OCR1B=TIMER_PERIOD;
  //prescaler=8, starting timer
  TCCR1B=(1<<CS11);

  sei();

  while (1);
  {
    if (status & BUFFER_RDY)
	{

	  /*TODO: process buffer here*/

	  status&=~BUFFER_RDY;
	}
  }
}