[p_v 0]

Задача: хочется определять обороты коллекторного двигателя по пульсациям тока, которые вызывает переключение обмоток. И чтобы математика на Cortex M0+ вписывалась в 0.1 сек за итерацию.

• Искомая частота - 600...6000 герц (450...45000 rpm для 8-полюсного коллекторника бормашинки).
• Желаемая точность измерений - 1% (относительно измеряемой частоты)

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

 

Что имеем:

• https://github.com/speedcontrols/dc_sc_grinder/tree/scripts/scripts - исходники, с которыми экспериментирую (в бранчах есть полные дампы).
• Написал честный FIR-дециматор, и перегнал дампы на 16 килогерц (примерно как в реальном девайсе)
• Дампы и прочие выклопы нормализовал в целочисленные.

Дамп на высоких оборотах

 

Дамп на низких оборотах

 

Спектр на высоких оборотах, 512 точек, обнулена постоянная составляющая (чтобы не портила масштаб). ~ 32 герца на точку. Нужен второй большой пик (который "троится"). ~ 4000 герц, похоже на правду (30000 rpm).

 

Спектр на низких оборотах, 512 точек, обнулена постоянная составляющая. Тоже второй пик, ~ 600 герц (4500 rpm), сходится.

 

С точностью на низах примерно понятно. Если видим что пик внизу - включаем ресамплер, и делаем FFT еще раз. Сетевые пульсации и гармоники в стороне от нужного диапазона - можно спокойно игнорировать. Правда, критично вырастает время выборки. Получается 0.8 сек при точности 1.3%. Терпимо. В крайнем случае похерю точность до 2%, приподняв частоту самплинга.

 

Остались вопросы:

• Что все-таки смотреть после FFT, реальную часть или корень из суммы квадратов?
• Просто искать максимум или как-то проверять соседние отсчеты?
• Как определять "отсутствие" оборотов (или когда ниже допустимых, что актуально для старта). Ведь если PID/ADRC подхватит ложный пик - может не дать разогнаться. Мерить при калибровке амплитуду частоты на низах и ставить порог в половину?
• Есть ли смысл накладывать перед FFT оконную функцию?
• Может еще чего-то не знаю и пропустил?

Изменено 10 октября, 2021 пользователем p_v

[blackfin 14]

33 minutes ago, p_v said:

Остались вопросы:

• Просто искать максимум или как-то проверять соседние отсчеты?
• Есть ли смысл накладывать перед FFT оконную функцию?
• Может еще чего-то не знаю и пропустил?

 

И тд..

[gridinp 2]

Может быть вам поможет автокорреляция?  Попробуйте измерить время от пика до первого перехода через 0, как четверть периода частоты

[stealthisname 7]

10.10.2021 в 08:37, p_v сказал:

Что все-таки смотреть после FFT, реальную часть или корень из суммы квадратов?

если FFT выдает комплексные значения, то учитывать нужно и действительную и мнимую часть.

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

[Eddy_Em 1]

On 10/10/2021 at 8:37 AM, p_v said:

И чтобы математика на Cortex M0+ вписывалась в 0.1 сек за итерацию.

А сейчас сколько времени эти все вычисления выполняются?

И зачем вам, собственно, БПФ? Попробуйте дискретные косинусные преобразования. Правда, флеш сильно раздуется таблицами косинусов…

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

P.S. А еще мне непонятно, зачем мучиться с определением частоты оборотов по пульсациям, если у вас УЖЕ есть устройство, которое эту самую частоту задает? Ну так и берите данные о скорости вращения со своего привода!..

Изменено 11 октября, 2021 пользователем Eddy_Em

[stealthisname 7]

10.10.2021 в 08:37, p_v сказал:

Как определять "отсутствие" оборотов (или когда ниже допустимых, что актуально для старта).

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

в качестве примера можно попробовать такой алгоритм на сигнале из Вашего файла hilda_15625Hz_zero_to_middle.

я вычислил fft для 512 отсчетов на протяжении всего файла. для каждого результата убрал первые несколько отсчетов - в них мы ожидаем сильную низкочастотную помеху.

получилась такая спектрограмма.

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

в начале файла шум не превышает порог

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

после выхода искомой частоты на постоянную все выборки дают максимум, превышающий порог

таким образом отсутствие оборотов не прошло алгоритм, PID/ADRC точно не подхватит ложный пик. Единственное, чем пришлось пожертвовать - скоростью захвата, так как мы определили включение немного позже действительного включения устройства.

 

вот как выглядит спектрограмма с порогом

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

результат работы алгоритма - оценка частоты

аналогичные результаты для сигнала из файла hilda_15625Hz_zero_to_low  с тем же размером fft и с тем же множителем при вычислении порога

 

Изменено 12 октября, 2021 пользователем stealthisname

[_pv 52]

для "оценки" частоты то же самое и без фурье можно получить

просто пропустив через ФВЧ чтобы отрезать 50Гц, а от остального посчитать отношение среднеквадратичного значения производной к среднеквадратичному значению сигнала, и среднеквадратичное значение сигнала в качестве порога и использовать.

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

[Eddy_Em 1]

ФВЧ можно и аппаратный поставить перед АЦП, сильно упростив себе работу (да и МК меньше придется телодвижений совершать). Про использование более шустрых ДКП вместо БПФ я уже выше говорил. Но сдается мне, какой-нибудь простой вариант периодограммы будет еще быстрей.

[p_v 0]

On 10/12/2021 at 12:28 PM, _pv said:

просто пропустив через ФВЧ чтобы отрезать 50Гц, а от остального посчитать отношение среднеквадратичного значения производной к среднеквадратичному значению сигнала, и среднеквадратичное значение сигнала в качестве порога и использовать.

Я не понял, что это за магия, и почему она вдруг стала пропорциональна оборотам, а не погоде на марсе например.

 

On 10/12/2021 at 7:21 AM, stealthisname said:

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

Так тоже можно наверное. Но по-моему это не улучшает более простой вариант - померить максимум на самых низких оборотах, и рубить все что меньше 70% (в текущих исходниках так и сделано). Либо я мог упустить какую-то важную мысль, тогда напишите пожалуйста точную формулу, что вы впредлагаете для порога.

 

[_pv 52]

1 hour ago, p_v said:

Я не понял, что это за магия, и почему она вдруг стала пропорциональна оборотам, а не погоде на марсе например.

потому что синус отличается от производной синуса по амплитуде ровно в w=2*Pi*f раз. то есть их отношение поделённое на 2Pi и есть частота.

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

в каких-то пределах и для не совсем синусоидальных сигналов в качестве "оценки" частоты тоже годится.

[p_v 0]

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

 

Скажем так, скорость самого БПФ на 1024 точки не самый критичный кусок. Там больше прижимает время набора отсчетов (1000/16000 => 0.0625 сек). Я пока попробую на БПФ, т.к. там предсказуемый результат, но если что-то пойдет не так - буду думать над вашим.

 

PS. Ну и что скрывать, просто интересно БПФ поюзать.

Изменено 13 октября, 2021 пользователем p_v

[Сергей Борщ 119]

12.10.2021 в 14:46, Eddy_Em сказал:

ФВЧ можно и аппаратный поставить перед АЦП

А что, разве можно его не ставить? Наложение спектров (aliasing) придумали слабаки и трУсы?

[_pv 52]

9 hours ago, p_v said:

Плюс подсчет девиации - не самая дешевая штука.

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

БПФ ~100 тактов на отсчёт.

[_pv 52]

5 hours ago, Сергей Борщ said:

А что, разве можно его не ставить? Наложение спектров (aliasing) придумали слабаки и трУсы?

Речь про фВч, там в спектре корявые 50Гц с гармониками на которые "полезный" сигнал наложен, и если пытаться без Фурье, во временной области, частоту оценивать, то они будут мешать, и их надо убрать.

On 10/12/2021 at 1:46 PM, Eddy_Em said:

ФВЧ можно и аппаратный поставить перед АЦП, сильно упростив себе работу

Аналоговый фильтр перед АЦП не сильно поможет, так как фильтр там на самом деле довольно злой должен быть (а не скользящее среднее как в примере на картинке с частотой среза 1.5кГц), у ТС частоты от 600Гц начинаются и он 1% хочет, при этом Гц до 200-300 надо бы всё задавить довольно хорошо. Ну то есть можно конечно попробовать такой фильтр сгородить снаружи, но это совсем не упростит, а скорее наоборот.

[p_v 0]

2 hours ago, _pv said:

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

Если быть точнее, там выпрямленное напряжение, поэтому 100 либо 120 герц, плюс конская вторая гармоника, и возможно приличные третяя-четвертая.

Поэтому надо сурово резать от нуля до 400 - 480 герц. И учитывать от 600 герц. Чтобы такое разделить - это от десятого порядка фильтры. Без плавучки - либо каскады IIR либо очень длинный FIR. Можно, но пипец как грустно.

 

7 hours ago, _pv said:

БПФ ~100 тактов на отсчёт.

Часто встречал такую оценку на этом форуме. А для каких она условий? Случайно не для аппаратной плавучки?

 

https://community.st.com/s/question/0D50X00009XkhK9/fft-in-stm32f0 - тут в четветом ответе пишут для M0+ 600К циклов на 1024 точки для Q15 (мне возможно понадобится Q31, т.к. постоянная составляющая очень большая).