Remove parabolic or curved trend
6 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Konvictus177
am 31 Aug. 2022
Kommentiert: Image Analyst
am 31 Aug. 2022
Hi,
I collected data from a laser sensor that measures a flat plate with two elevations. Due to the curved lense of the sensor the data has a parabolic trend.
How do I get rid of that parabolic nature of the graph so I can analyze the elevations based on a flat plate/signal? Is there a way to filter it out and make the signal "flat" with the two elevations?
I tried detrending(y_data,1) but it does not filter out a polynomial trend and it does not leave out the real elevations.
The data is attached.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1112960/image.png)
0 Kommentare
Akzeptierte Antwort
Star Strider
am 31 Aug. 2022
LD = load(websave('y_data','https://www.mathworks.com/matlabcentral/answers/uploaded_files/1112945/y_data.mat'));
y_data = LD.y_data;
L = numel(y_data)
x = linspace(0, L-1, L);
y_data_dt = detrend(y_data, 5)
figure
yyaxis left
plot(x, y_data)
ylabel('Original')
yyaxis right
plot(x, y_data_dt)
ylabel('Detrended')
grid
Experiment with different polynomial orders to get the result you want. Another option might be to use the bandpass function (ideally with 'ImpulseResponse','iir') depending on the result you want.
.
2 Kommentare
Star Strider
am 31 Aug. 2022
A highpass filter seems to work to eliminate the low-frequency trends without changing the rest of the data. Choose the cutoff frequency to get the result you want.
This uses the first ‘valley’ freqeuency —
LD = load(websave('y_data','https://www.mathworks.com/matlabcentral/answers/uploaded_files/1112945/y_data.mat'));
y_data = LD.y_data.';
L = numel(y_data)
x = linspace(0, L-1, L);
figure
plot(x,y_data)
grid
Fs = 1/(x(2)-x(1))
Fn = Fs/2;
NFFT = 2^nextpow2(L)
FTds = fft(y_data-mean(y_data))/L;
Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);
figure
plot(Fv, abs(FTds(Iv))*2, 'DisplayName','Fourier Transform')
grid
xlim([0 0.1])
[pks,plcs] = findpeaks(abs(FTds(Iv))*2, 'MinPeakProminence',0.00025);
[vys,vlcs] = findpeaks(-abs(FTds(Iv))*2, 'MinPeakProminence',0.00025);
hold on
plot(Fv(plcs), pks, '^r', 'DisplayName','Peaks')
plot(Fv(vlcs), -vys, 'vr', 'DisplayName','Valleys')
hold off
xlabel('Frequency')
ylabel('Magnitude')
legend('Location','best')
Fvpk = Fv(plcs); % Peak Frequencies
Fvvy = Fv(vlcs); % Valley Frequencies
% BPv = [Fvvy([1 end])] % Peak Frequencies ± Offset To Use In 'bandpass' Call
N = 50; % Padding Vector Length
% y_data_filt = bandpass([zeros(N,1)+mean(y_data); y_data], BPv, Fs, 'ImpulseResponse','iir'); % Filter Signal
y_data_filt = highpass([zeros(N,1)+y_data(1); y_data], Fvvy(1), Fs, 'ImpulseResponse','iir'); % Filter Signal
y_data_filt = y_data_filt(N+1:end); % Adding Vector Of Mean Values Eliminates Initial Filter Transient, Remove Those Values Later
figure
plot(x,y_data,'-b', 'DisplayName','Original')
hold on
plot(x,y_data_filt, '-r', 'DisplayName','Filtered')
hold off
grid
legend('Location','best')
This is the best I can do with this signal.
.
Weitere Antworten (1)
Image Analyst
am 31 Aug. 2022
Bearbeitet: Image Analyst
am 31 Aug. 2022
If you have the Computer Vision Toolbox you can fit a quadratic to the lower part of your curve, without explicitly telling it what elements those are, by using fitPolynomialRANSAC
s = load('y_data.mat')
y = s.y_data;
subplot(2, 1, 1);
plot(y, 'b.-', 'LineWidth', 2);
hold on;
xlabel('Index', 'FontSize',fontSize);
ylabel('y', 'FontSize',fontSize);
grid on;
x = 1 : length(y);
xy = [x(:), y(:)]
N = 2; % second-degree polynomial
maxDistance = 1; % maximum allowed distance for a point to be inlier
coefficients = fitPolynomialRANSAC(xy, N, maxDistance)
yFitted = polyval(coefficients, x);
plot(yFitted, 'r.-', 'LineWidth', 2);
% Subtract
baseLineCorrectedY = y - yFitted
subplot(2, 1, 2);
plot(baseLineCorrectedY, 'b.-', 'LineWidth', 2);
hold on;
grid on;
xlabel('Index', 'FontSize',fontSize);
ylabel('y', 'FontSize',fontSize);
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1113065/image.png)
2 Kommentare
Image Analyst
am 31 Aug. 2022
Yeah. Hmmm...I guess that's just because of the way RANSAC works by sampling random points. You could put it in a loop and call it, say 50 times and take the one set of y values where the mean absolute corrected y value is closest to 0.
% Demo by Image Analyst
% Initialization Steps.
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 = 18;
s = load('y_data.mat')
y = s.y_data;
subplot(2, 1, 1);
plot(y, 'b.-', 'LineWidth', 2);
hold on;
xlabel('Index', 'FontSize',fontSize);
ylabel('y', 'FontSize',fontSize);
grid on;
x = 1 : length(y);
xy = [x(:), y(:)];
N = 2; % second-degree polynomial
maxDistance = 1; % maximum allowed distance for a point to be inlier
numTrials = 50;
meanAbsValue = zeros(1, numTrials);
% Try RANSAC 50 times and we'll choose the best.
subplot(2, 1, 2);
for k = 1 : numTrials
coefficients{k} = fitPolynomialRANSAC(xy, N, maxDistance);
yFitted = polyval(coefficients{k}, x);
% Subtract
baseLineCorrectedY = y - yFitted;
plot(baseLineCorrectedY, '-', 'LineWidth', 1);
hold on;
grid on;
xlabel('Index', 'FontSize',fontSize);
ylabel('y', 'FontSize',fontSize);
% Compute mean value
meanAbsValue(k) = mean(abs(baseLineCorrectedY));
end
% Find the closest
[~, index] = min(meanAbsValue);
yFitted = polyval(coefficients{index}, x);
% Subtract
baseLineCorrectedY = y - yFitted;
subplot(2, 1, 1);
plot(yFitted, 'r.-', 'LineWidth', 2);
% Plot final baseline corrected signal on a separate figure.
figure;
plot(baseLineCorrectedY, 'r.-', 'LineWidth', 2);
grid on;
xlabel('Index', 'FontSize',fontSize);
ylabel('y', 'FontSize',fontSize);
% Put a line along 0
yline(0, 'LineWidth', 2)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1113320/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1113325/image.png)
Siehe auch
Kategorien
Mehr zu Bartlett finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!