Dilatation and Erosion as custom function for signal processing
10 次查看(过去 30 天)
显示 更早的评论
I want to define the two function for accomplish the morphological operation of erosion and dilatation on a 1-D signal (in this case an ecg signal). I know that I can use functions like imerode and imdilate, but I want to define them by myself, obtaining a result comparable to the functions mentioned before.
To design the two function I have to followed the example propose below:
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1060405/image.png)
Based on this knowledge I have defined two functions: dilatation and erosion (below). However when I apply the algorithm (show below in the code) the result is something that I do not understand. Based on the methods presented above I can't understand what I is missing or what is wrong in the two functions.
My Code:
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 15;
% Signal info.
G = 200; % Gain
Fs = 360; % [Hz]
L = 3600; % lenght of ECG signals
T = linspace(0,L/Fs,L); % time axis
F = linspace(-Fs/2, Fs/2, L); % Frequency axis
% Load Signal
load("100m.mat");
ecg = val/G;
% Plot the ecg signal
figure(1); plot(T, ecg); grid on;
title("ECG Signal","FontSize",fontSize);
xlabel("Time (sec)", "FontSize", fontSize);
ylabel("voltage [mV]", "FontSize", fontSize);
% Apply MMF algorithm
% Declaration of the two Structuring Elements
B1 = [0 1 5 1 0]; B2 = [1 1 1 1 1];
% modified morphological filtering algorithm
d1 = dilatation(ecg, B1); % d1 = imdilate(ecg, B1);
e1 = erosion(d1,B2); % e1 = imerode(ecg, B2);
e2 = erosion(ecg,B1); % e2 = imerode(ecg, B1);
d2 = dilatation(e2,B2); % d2 = imdilate(e2, B2);
mmf = (e1+d2)/2; % average
% Plot
figure(2); plot(T,mmf,"b-"); grid on;
title("MMF Denoised Ecg signal", "FontSize", fontSize);
xlabel("Time (sec)", "FontSize", fontSize);
ylabel("Voltage (Hz)", "FontSize", fontSize);
Dilatation Function:
function [output] = dilatation(inputSignal,structuringElement)
% Dilatation operation on a 1-D signal
% insputSignal = signal in input
% StructuringElement = structuring element in input
output = zeros(length(inputSignal));
m = length(structuringElement);
n = length(inputSignal);
for i = abs((m-1)/2):n-abs((m+1)/2)-1
for j=1:m
output(i) = max(output(i-1), inputSignal(i-abs((m-1)/2)+j)+structuringElement(j));
end
end
end
Erosion Function:
function [output] = erosion(inputSignal,structuringElement)
% Erosion operation on a 1-D signal
% insputSignal = signal in input
% StructuringElement = structuring element in input
output = zeros(length(inputSignal));
m = length(structuringElement);
n = length(inputSignal);
for i = abs((m-1)/2):n-abs((m+1)/2)-1
for j=1:m
output(i) = min(output(i-1), inputSignal(i-abs((m-1)/2)+j)-structuringElement(j));
end
end
end
0 个评论
采纳的回答
DGM
2022-7-10
编辑:DGM
2022-7-10
The loop part can be cleaned up, but a big part of why the output looks messed up starts here:
output = zeros(length(inputSignal));
This will create a square matrix, not a vector. So then this assignment inside the loop:
output(i) = min(...);
... will populate the first column. Then if you feed the result to plot(), it plots each row as a series.
So just make sure to preallocate a vector.
I'm going to assume that this is closer to what you're looking for.
insignal = [0 0 0 0 0 1 1 1 0 0 0 0 0 1 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0];
se = [1 1 1];
outsignal = erosion(insignal,se);
subplot(2,1,1)
plot(insignal)
subplot(2,1,2)
plot(outsignal)
function [output] = erosion(inputSignal,structuringElement)
% Erosion operation on a 1-D signal
% insputSignal = signal in input
% StructuringElement = structuring element in input
% get sizes
M = length(structuringElement);
N = length(inputSignal);
% this relies on the strel being odd width
if mod(M,2)==0
error('strel length must be odd')
end
% perform dilation
hw = floor(M/2); % half-width of strel
output = zeros(1,N);
for n = ceil(M/2):N-ceil(M/2)+1
output(n) = min(inputSignal(n-hw:n+hw) - structuringElement);
end
end
Bear in mind that the behavior at the ends is due to the way the operation is constrained there. The alternative to avoiding the ends of the signal vector is to add enough padding to them that you can filter the whole vector.
9 个评论
DGM
2022-7-13
编辑:DGM
2022-7-13
Augh. It's picking up the padding anyway. Try doing this instead of padding with zeros:
input = padarray(input,[0 hw],'replicate','both');
I am kind of confused by the behavior of this method of dilation/erosion. I'm used to doing this with images, where the strel is multiplied instead of added/subtracted. It's just not a convention I'm used to thinking about. At least with the above change, it seems to work okay. See what you think.
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Multirate Signal Processing 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!