Данные потенциометра фильтрации неисправностей (шум, высокие пики)
Я собрал некоторые данные с потенциометра, используя микроконтроллер Arduino. Вот данные, которые были выбраны с частотой 500 Гц (это много данных):
http://cl.ly/3D3s1U3m1R1T?_ga=1.178935463.2093327149.1426657579
Если вы увеличите масштаб, вы увидите, что по сути у меня есть банк, который просто вращается взад и вперед, т.е. я должен увидеть линейное увеличение, а затем линейное уменьшение. Хотя общая форма данных подтверждает это, почти каждый раз, когда появляются действительно ужасные раздражающие (иногда удивительно широкие) всплески, которые мешают действительно хорошей форме. Есть ли способ, которым я могу сделать какой-то тип алгоритма или фильтра, который исправляет это? Я пробовал медианный фильтр и использовал процентили, но ни один не работал. Я имею в виду, я чувствую, что это не должно быть самым сложным, потому что я могу ясно видеть, как это должно выглядеть - в основном минимум того, где происходят всплески - но по какой-то причине все, что я пытаюсь, с треском проваливается или, по крайней мере, теряет целостность исходные данные.
Буду признателен за любую помощь, которую я могу получить с этим.
3 ответа
Есть много способов решить вашу проблему. Однако ни один из них никогда не будет идеальным. Я дам вам 2 подхода здесь.
Скользящая средняя (фильтр нижних частот)
В Matlab одним из простых способов фильтрации низких частот без необходимости явного использования БПФ является использование filter function` (доступно в базовом пакете, вам не нужно никаких специальных инструментов).
Вы создаете ядро для фильтра и применяете его дважды (по одному в каждом направлении), чтобы отменить введенный фазовый сдвиг. По сути, это фильтр "Скользящее среднее" с нулевым фазовым сдвигом. Размер (длина) ядра будет контролировать, насколько тяжелым будет процесс усреднения.
Так, например, 2 разных длины фильтра:
n = 100 ; %// length of the filter
kernel = ones(1,n)./n ;
q1 = filter( kernel , 1 , fliplr(p) ) ; %// apply the filter in one direction
q1 = filter( kernel , 1 , fliplr(q1) ) ; %// re-apply in the opposite direction to cancel phase shift
n = 500 ; %// length of the filter
kernel = ones(1,n)./n ;
q2 = filter( kernel , 1 , fliplr(filter( kernel , 1 , fliplr(p) )) ) ; %// same than above in one line
По вашим данным выдадим:
Как видите, каждый размер фильтра имеет свои плюсы и минусы. Чем больше вы фильтруете, тем больше ваших пиков вы отменяете, но тем больше вы искажаете свой исходный сигнал. Это зависит от вас, чтобы найти ваши оптимальные настройки.
2) Поиск производных аномалий
Это другой подход. Вы можете наблюдать по своему сигналу, что выбросы в большинстве случаев внезапны, это означает, что изменение значения вашего сигнала происходит быстро и, к счастью, быстрее, чем "нормальная" скорость изменения желаемого сигнала. Это означает, что вы можете вычислить производную вашего сигнала и идентифицировать все пики (производная будет намного выше, чем для остальной части кривой).
Поскольку это идентифицирует только "начало" и "конец" шипов (а не случайное плато в середине), нам потребуется немного расширить зону, определенную как неисправную этим методом.
Когда идентификация ошибочных данных завершена, вы просто отбрасываете эти точки данных и повторно интерполируете свою кривую в течение исходного интервала (принимая поддержку по оставленным точкам).
%% // Method 2 - Reinterpolation of cancelled data
%// OPTIONAL slightly smooth the initial data to get a cleaner derivative
n = 10 ; kernel = ones(1,n)./n ;
ps = filter( kernel , 1 , fliplr(filter( kernel , 1 , fliplr(p) )) ) ;
%// Identify the derivative anomalies (too high or too low)
dp = [0 diff(ps)] ; %// calculate the simplest form of derivative (just the difference between consecutive points)
dpos = dp >= (std(dp)/2) ; %// identify positive derivative above a certain threshold (I choose the STD but you could choose something else)
dneg = dp <= -(std(dp)/2) ; %// identify negative derivative above the threshold
ixbad = dpos | dneg ; %// prepare a global vector of indices to cancel
%// This will cancel "nPtsOut" on the RIGHT of each POSITIVE derivative
%// point identified, and "nPtsOut" on the LEFT of each NEGATIVE derivative
nPtsOut = 100 ; %// decide how many points after/before spikes we are going to cancel
for ii=1:nPtsOut
ixbad = ixbad | circshift( dpos , [0 ii]) | circshift( dneg , [0 -ii]) ;
end
%// now we just reinterpolate the missing gaps
xp = 1:length(p) ; %// prepare a base for reinterpolation
pi = interp1( xp(~ixbad) , p(~ixbad) , xp ) ; %// do the reinterpolation
Это даст:
Красный сигнал - результат вышеуказанного скользящего среднего, зеленый сигнал - результат производного подхода.
Есть также настройки, которые вы можете изменить, чтобы скорректировать этот результат (порог для производной, nPtsOut и даже начальное сглаживание данных).
Как вы можете видеть, при той же величине подавления всплесков, что и при использовании метода скользящего среднего, он несколько более уважает целостность исходных данных. Однако это также не идеально, и некоторый интервал все еще будет деформирован. Но, как я сказал в начале, ни один метод не является идеальным.
Пики шума вызваны отскоком стеклоочистителя POT вдоль резистивной дорожки при повороте ручки. Это общая проблема с ними. В будущем вы должны смотреть на добавление конденсатора 0.1 мкФ к выходу POT, и это должно решить проблему.
С вашими текущими данными самый простой вариант - просто сделать простое скользящее среднее и визуально настроить количество выборок, усредненных до тех пор, пока выбросы не будут подавлены в достаточной степени, не влияя на базовые данные. Обратите внимание, что скользящее среднее - это просто фильтр нижних частот с частотной характеристикой sinc.
Обычный способ постобработки данных такого рода - это выполнить БПФ (используя соответствующую оконную функцию), обнулить значения шума над сигналом, представляющим интерес, а затем принять обратное БПФ. Это также просто низкочастотная фильтрация (с взвешенной скользящей средней функции sinc *), но вы используете информацию, предоставленную БПФ, для выбора частоты среза. Если вас не устраивает математика, вовлеченная в это, то просто используйте простой фильтр скользящих средних. Это должно быть хорошо для ваших нужд.
Кажется, у вас есть большие шипы вблизи максимальных и минимальных очков вашего банка. Например, вы можете ограничить диапазон ваших действительных данных от 200 до 300.
Другой вариант - фильтр нижних частот 1-го порядка, такой как этот
alpha = 0.01 %parameter to tune!
p_filtered(1) = p(1);
for i=2:length(p)
p_filtered(i) = alpha*p(i) + (1-alpha)* p_filtered(i-1);
end