Проверьте пожалуйста алгоритмы IIR фильтров

[thermit 0]

В данных из гита постоянь есть. У вас нет. У самурая спектр более адекватный.

[p_v 0]

17 minutes ago, thermit said:

В данных из гита постоянь есть. У вас нет. У самурая спектр более адекватный.

Не понял. Я в спектр раскладываю предварительно фильтрованный сигнал. На последних графиках отрезано ниже 800 и выше 15000 герц. Иначе сетевые пульсации все забьют нафик. Соответственно, из за high-pass 800hz постоянная составляющая уходит.

 

У бормашинки минимально возможные обороты 4500 rpm. При восьми полюсах это 1500 герц. Вроде параметры пред-фильтра адекватные, режут мусор но лишнего не задевают.

Edited October 6, 2021 by p_v

[Самурай 0]

1 час назад, p_v сказал:

Я в спектр раскладываю предварительно фильтрованный сигнал. На последних графиках отрезано ниже 800 и выше 15000 герц. Иначе сетевые пульсации все забьют нафик.

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

1 час назад, p_v сказал:

У бормашинки минимально возможные обороты 4500 rpm. При восьми полюсах это 1500 герц.

4500*8/60 = 600Гц. Разве нет? Ниже часть спектра для hilda_250kHz_rpm_low, 100к отсчетов, по оси х частота в Гц. Хорошо виден пик как раз на 600Гц. Это хорошая новость. Плохая новость в том, что этих пиков там как.... Много в общем. А самое скверное, что по амплитуде они не очень то сильно и отличаются. Совершенно не очевидно, что пик на 600Гц, это та самая искомая частота...

 

И еще пара нюансов:

1. Если сигнал не выпрямлять перед оцифровкой, то количество пиков будет существенно меньше.

2. Явно не хватает разрядности АЦП.

[thermit 0]

3 часа назад, thermit сказал:

В данных из гита постоянь есть. У вас нет. У самурая спектр более адекватный.

Вам нужно внятно сформулировать задачу для начала. Заниматься хиромантией лично у меня никакого желания нет.

[p_v 0]

1 hour ago, Самурай said:

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

Ну нема в контроллере флоатов, и тем паче даблов. А скриптами я эмулирую то, что потом в контроллере реально повторить.

 

1 hour ago, Самурай said:

4500*8/60 = 600Гц. Разве нет?

Хм... промахнулсо. Поправлю. Правда я дрюкал дамп _high (максимальные обороты), там ничего не задето. Но для low не прокатит, вы правы. Откачу на 400 герц, чтобы хотя бы первые 3 гармоники прибрать. Максимальная частота тогда получается 4-6 килогерц (30000-45000rpm).

 

И еще, на дампах забыли тахометром обороты померить, так что я не знаю, какие на _low. Знакомый питал хильду с автотрансформатора через выпрямитель. Мог меньше выставить.

 

1 hour ago, Самурай said:

Если сигнал не выпрямлять перед оцифровкой, то количество пиков будет существенно меньше.

Я плохо представляю, как сделать высокочастотный ШИМ без выпрямления. А если делать фазовую регулировку симистором, там спектр должен еще сильнее убиваться. Диаграммы снимали, мне с них сигнал совсем не понравился, поэтому даже не стал в репозиторий заливать.

 

Переделал спектрограммы, без фильтрации. 250000/16384 => 15,2588 - шаг в герцах на 1 точку.

 

hilda_250kHz_rpm_low.txt

 

Виден пичок на ~ 600 герцах => 4500 rpm. охоже на правду

 

hilda_250kHz_rpm_high.txt

 

Группа пиков около 3600 Гц (возможно обороты слегка дрожат) => 27000 rpm. Тоже похоже на правду.

 

Еще какой-то паразитный пик на 2800 Гц, фик знает почему. Может обмотки не симметричные. Я могу понять как это даст пики на 1/8 частоты от основной, но не понимаю откуда 3/4 берутся.

 

@Самурай, проверьте, может ваш спектроанализатор хочет в качестве разделителя точку а не запятую. У вас очень странный результат выходит.

1 hour ago, thermit said:

Вам нужно внятно сформулировать задачу для начала. Заниматься хиромантией лично у меня никакого желания нет.

Хочется определять обороты универсального двигаетеля (в бормашинке) по частоте колебаний тока. Через противо-ЭДС уже сделано, работает приемлимо, но хочется еще лучше.

 

Судя по прикрепленным сонограммам, задача имеет решение. Но надо разобраться, возможно ли оно на чахлом контроллере. Текущая гипотеза - попробовать почистить сигнал и померить период через zero cross.

[p_v 0]

Вот графики из предыдущего поста, но с вырезанными 30 отсчетами в начале (отбрасываем 0-450 герц, чтобы не херили масштаб)

 

hilda_250kHz_rpm_low.txt

 

hilda_250kHz_rpm_high.txt

 

Вроде нужные пики хорошо видны. На втором графике, правда, двоятся, но это не критично.

 

А можно FFT как-то "ступенчато" считать? Сначала 32 точки, чтобы отбросить 0-400 герц и найти диапазон с пиком, а потом уже найденный диапазон просчитать точнее?

[stealthisname 0]

по условию задачи и по обсуждению выше, в цифровой обработке допускается использовать БИХ фильтры второго порядка (только если нет высоких требований к точности коэффициентов и к точности вычислений)

я отфильтровал сигналы из hilda_50kHz_rpm_high и hilda_50kHz_rpm_low в ВЧ фильтре с коэффициентами 

b = [1,-1]; a = [1,-0.9375]; 

с такой АЧХ:

и затем отфильтровал в группе полосовых фильтров на разных частотах с такими коэффициентами:

Скрытый текст

f =   400; b = [  0.03281759059,0, -0.03281759059]; a = [1,     -1.9340458,         0.9375];
f =   900; b = [  0.03156366737,0, -0.03156366737]; a = [1,   -1.924120029,         0.9375];
f =  1400; b = [   0.0313784307,0,  -0.0313784307]; a = [1,   -1.906600635,         0.9375];
f =  1900; b = [  0.03131856954,0, -0.03131856954]; a = [1,    -1.88155676,         0.9375];
f =  2400; b = [  0.03129199717,0, -0.03129199717]; a = [1,    -1.84908724,         0.9375];
f =  2900; b = [  0.03127793598,0, -0.03127793598]; a = [1,   -1.809320217,         0.9375];
f =  3400; b = [  0.03126961201,0, -0.03126961201]; a = [1,   -1.762412634,         0.9375];
f =  3900; b = [  0.03126428538,0, -0.03126428538]; a = [1,   -1.708549613,         0.9375];
f =  4400; b = [  0.03126067664,0, -0.03126067664]; a = [1,   -1.647943728,         0.9375];
f =  4900; b = [  0.03125812333,0, -0.03125812333]; a = [1,   -1.580834161,         0.9375];

 

АЧХ фильтров:

от результата фильтрации в группе полосовых фильтров взял по модулю и накопил

процесс такого накопления в hilda_50kHz_rpm_high:

результат такого накопления:

результат накопления для hilda_50kHz_rpm_low:

 

в сигнале hilda_50kHz_rpm_high больше всего накопил фильтр на частоте 3900. в сигнале hilda_50kHz_rpm_low больше всего накопил фильтр на частоте 400.

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

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

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

 

[p_v 0]

Надо 1% точности, да по корявому сигналу. На низко-порядковых фильтрах откровенно стремно.

Есть вопросы

• По БПФ. Где можно посмотреть, как выбирать "разрядность" (Q15?) и сколько циклов займет на Cortex-M0+ 256 точек?
• По ресамплингу. Что использовать на МК для понижения частоты в 2^n раз? Децимация там тривиальная, но еще ведь фильтр нужен.

По поводу точности. Если правильно понимаю, то для погрешности меньше 1% во всем диапазоне мне надо Constant-Q transform. Чтобы не мудрить с имплементацией, можно сделать в 2 фазы:

1. Делаем БПФ на 32 точки, чтобы грубо определить искомый диапазон.
2. Включаем нужный ресамплинг (отрезаем выше искомого диапазона), набираем новые данные и делаем БПФ по 128-256 точкам. Т.к. после ресамплинга искомый диапазон всегда окажется в верхней половине, точность определения частоты будет достаточной.

Я в правильном направлении думаю?

[Ruslan1 0]

Я тоже сначала без плавающей точки бодался (На майкрочипе pic24). Потом плюнул на все и сделал честный float32. Памяти не хватало- добавил внешнюю последовательную RAM. И всем пофиг что он 10 секунд считает, зато с нужной точностью (4096 точек, по 16-битный данным АЦП).

Вот повозитесь Вы с этими всем (диапазон, ресемплинг, набираем новые данные...) и придете к классике.