Matlab и след това c реализация: Гладки данни

Имам някои данни, които са малко шумни. Искам да го почистя, както и производното. Много прост "филтър" като Y(n)=0.2*X(n)+0.8*Y(n-1) ми подхожда добре: с други думи, изходът е достатъчно гладък. Ето пример както с данни, така и с техните производни:

Отгоре: данни, отдолу: производни на данни

Въпреки това, тъй като е в реално време, не мога да приема голямо забавяне на данните. За производното забавянето е около 10 итерации: 5 итерации (за получаване на филтрираните данни) и след това още 5 итерации (за получаване на филтрираните производни данни).

Как мога да внедря алгоритъм в Matlab - добре, имам нужда от c реализация - която произвежда подобен резултат, но с много малко забавяне, като максимум 3 итерации ДОРИ и ОСОБЕНО за производното.

Ето кода на Matlab:

%Here is data input;
data = [5554 5767 5737 5828 6011 6011 6103 6225 6408 6286 6347 6439 6591 6683 6530 6652 6805 7019 7080 7019 6988 6988 6958 6958 6927 6896 6866 6805 6835 6713 6591 6713 6713 6652 6713 6713 6835 6896 6744 6622 6683 6561 6591 6591 6622 6561 6500 6652 6774 6805 6958 7049 7171 7232 7293 7293 7415 7476 7690 8087 8056 7965 8178 8270 8148 8178 8392 8239 8056 8026 7873 7873 7476 7446 7476 7293 7171 7141 7263 7019 7141 6774 6652 6530 6530 6469 6561 6622 6652 6713 6866 6958 6958 7110 6958 7049 6805 6866 6866 6744 7049 6805 6652 6530 6622 6744 6408 6652 6683 6713 6713 6622 6744 6805 6713 6683 6805 6835 6927 6866 6805 6744 7049 7232 7171 7141 7202 7354 7324 7415 7415 7415 7965 8178 8178 8178 8178 7995 7995 7965 7843 7904 7904 7812 7812 7720 7690 7843 7873 7904 8056 7904 7843 7934 7873 7873 7781 7781 7873 7873 7873 7690 7659 7690 7720 7751 7659 7598 7629 7659 7568 7537 7629 7537 7659 7568 7568 7476 7568 7598 7354 7324 7476 7354 7385 7202 7293 7293 6958 6866 6866 6683 6927 6927 6927 7019 7019 7049 7141 7202 7293 7385 7568 7446 7446 7415 7354 7263 7293 7202 7202 7202 7232 7202 7263 7141 7171 7019 7141 7232 7080 7232 7141 7263 6958 7019 7019 6866 6896 6744 6774 6805 6835 6866 6774 6744 6774 6744 6744 6866 6958 6958 7019 6988 6988 6896 6958 7080 7080 7171 7232 7202 7354 7385 7385 7476 7720 7690 7598 7659 7659 7690 7720 7720 7690 7690 7781 7751 7659 7812 7873 7659 7629 7446 7446 7629 7354 7415 7476 7568 7415 7385 7385 7324 7202 7324 7110 6958 7019 7019 7049 6774 7202 6561 6408 6195 5920 5920 5920 5859 5615 5645 5584 5676 5432 5401 5126 5096 5004 5004 4791 4852 4791 4791 4730 4760 4730 4730 4699 4760 4669 4730 4760 4760 4791 4882 4882 4882 4974 4943 4974 5065 5096 5157 5218 5279 5279]';

%Define smooth coefficient
c=0.2;

%Filter data
data_filtered=zeros(size(data),1);
for n=drange(2:size(data))
    data_filtered(n,1)=c*data(n)+(1-c)*data_filtered(n-1,1);
end

%Graph both data and filtered data
subplot(2,1,1);
plot([1:size(data)], data, [1:size(data)], data_filtered);

%Calculate derivative data. Delta time is constant.
for n=drange(2:size(data_filtered))
    data_derivative(n,1)=data(n)-data(n-1,1);
end

%Calculate derivative based on filtered data. Delta time is constant.
for n=drange(2:size(data_filtered))
    data_filtered_derivative(n,1)=data_filtered(n)-data_filtered(n-1,1);
end

%Filter derivative (which is based on filtered data)
data_filtered_derivative_filtered=data_filtered_derivative;
for n=drange(2:size(data_filtered_derivative))
    data_filtered_derivative_filtered(n,1)=c*data_filtered_derivative(n)+(1-c)*data_filtered_derivative_filtered(n-1,1);
end

%Graph both derivative data and filtered derivative data
subplot(2,1,2);
plot([1:size(data_derivative)], data_derivative, [1:size(data_derivative)], data_filtered_derivative_filtered);

person gregoiregentil    schedule 09.09.2013    source източник
comment
Какво имаш предвид с итерации? Брой времеви стъпки на вашия алгоритъм в реално време? Искам да кажа, ограничени ли сте от времето за изчисление или от „груповото забавяне“, ако вашите цифрови филтри?   -  person Bas Swinckels    schedule 10.09.2013
comment
Тук не говорим за ограничение на изчислението. На първата графика можете ясно да видите, че зелената крива е изместена надясно в сравнение със синята крива. Това е, което искам да избегна, доколкото е възможно. Предполагам, че някои хора ще предложат филтър на Калман, но все пак моят експеримент показва, че има голямо забавяне.   -  person gregoiregentil    schedule 10.09.2013
comment
Не съм експерт по филтрите, но мисля, че това забавяне е нещо фундаментално. Вие правите нещо в реално време, така че не можете да правите трикове като filtfilt, като прилагате един и същ филтър напред и назад, за да избегнете забавяния. Закъснението вероятно се мащабира обратно пропорционално на граничната честота на вашия нискочестотен, така че опитайте да увеличите граничната честота колкото е възможно повече. Накрая имайте предвид, че нискочестотен филтър от първи ред, последван от производна, е ефективно високочестотен филтър от първи ред. Може би бихте могли да намалите забавянето, като приложите това в единичен IIR или така, вместо да го правите на две стъпки.   -  person Bas Swinckels    schedule 10.09.2013
comment
@Bas Swinckels: Благодаря. Точно това е смисълът: не можете да използвате данни в бъдеще с филтър напред/назад, тъй като е в реално време. Просто се опитвам да минимизирам колкото е възможно повече забавянето и все пак да получа гладки данни.   -  person gregoiregentil    schedule 10.09.2013


Отговори (1)


Ако искате отличен филтър с ниско забавяне, изглежда нещо като Jurik JMA може да е за Вие. Вероятно е доста скъпо (и определено не е безплатно), но статията може да ви даде идеи за други (вероятно безплатни) алгоритми.

Не съм тествал този алгоритъм, просто го намерих и харесах снимките.

person DasKrümelmonster    schedule 10.09.2013