Hello guys,
Tanks a lot for all tour help. But : - Image Analyst, your way of doing delete to much information around the artifact - Star Rider, your way is too "binary", the artifact does not necessarily have an absolute value (threshold) - Cedric, your way seems efficient, but I don't understand it :)
So I was inspired by all of your solution and made my own that is based on a moving average, that allows me to apply a threshold as a percentage :
function signal = deleteArtifactMax(signal, threshold, interpolation)
global samplingFrequency maxInterpolationDuration;
signal(isnan(signal)) = 0;
time = repmat([1:size(signal,1)]',1,size(signal,2));
movingAverage = (...
signal(1:end-9,:) + ...
signal(2:end-8,:) + ...
signal(3:end-7,:) + ...
signal(4:end-6,:) + ...
signal(5:end-5,:) + ...
signal(6:end-4,:) + ...
signal(7:end-3,:) + ...
signal(8:end-2,:) + ...
signal(9:end-1,:) + ...
signal(10:end,:))/10;
isValid = abs(signal(5:end-5,:) - movingAverage) <= movingAverage*threshold;
notValid = logical([zeros(4,size(signal,2)); ~isValid; zeros(5,size(signal,2))]);
signal(notValid) = NaN;
signal(signal == 0) = NaN;
timeCut = time .* notValid;
for i = 1:size(signal,2)
timeCutCell{i} = find(timeCut(:,i) ~= 0);
end
if interpolation == true
for i = 1:size(signal,2)
if ~isempty(timeCutCell{1,i})
% Delete too long interpolation
iNotToInterpolate = [];
dTimeCutCell = [0; diff(timeCutCell{:,i})];
k=0;
for j=1:size(dTimeCutCell)
if dTimeCutCell(j) == 1
k=k+1;
else
if k > maxInterpolationDuration*60/samplingFrequency
iNotToInterpolate = [iNotToInterpolate j-k:j-1];
end
k=0;
end
end
% Selection of points to interpolate
iToInterpolate = setdiff(timeCutCell{1,i},timeCutCell{1,i}(iNotToInterpolate));
interpTime = time(:,i);
interpTime(iToInterpolate) = [];
interpSignal = signal(:,i);
interpSignal(iToInterpolate) = [];
% Interpolation
signal(iToInterpolate,i) = interp1(interpTime, interpSignal, iToInterpolate, 'linear');
end
end
end
end
Thanks a lot for your help and for your time.
Max