필터 지우기
필터 지우기

Remove parabolic or curved trend

조회 수: 6 (최근 30일)
Konvictus177
Konvictus177 2022년 8월 31일
댓글: Image Analyst 2022년 8월 31일
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.

채택된 답변

Star Strider
Star Strider 2022년 8월 31일
Use the detrend function —
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)
L = 799
x = linspace(0, L-1, L);
y_data_dt = detrend(y_data, 5)
y_data_dt = 1×799
0.0028 0.0027 0.0025 0.0023 0.0022 0.0020 0.0018 0.0016 0.0014 0.0013 0.0011 0.0010 0.0008 0.0006 0.0004 0.0003 0.0001 -0.0001 -0.0002 -0.0002 -0.0001 -0.0001 -0.0000 -0.0000 0.0000 0.0001 0.0000 -0.0001 -0.0001 -0.0001
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
Konvictus177
Konvictus177 2022년 8월 31일
편집: Konvictus177 2022년 8월 31일
Thanks for the quick reply.
I want to measure the distance to both elevations from my base profile which should be around 10 microns.
Simple detrending seems not to be the right thing to do in my opinion since detrending also removes the trend caused by the elevations and higher detrending orders will detrend more then they should.
Star Strider
Star Strider 2022년 8월 31일
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)
L = 799
x = linspace(0, L-1, L);
figure
plot(x,y_data)
grid
Fs = 1/(x(2)-x(1))
Fs = 1
Fn = Fs/2;
NFFT = 2^nextpow2(L)
NFFT = 1024
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.
.

댓글을 달려면 로그인하십시오.

추가 답변 (1개)

Image Analyst
Image Analyst 2022년 8월 31일
편집: Image Analyst 2022년 8월 31일
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);
  댓글 수: 2
Konvictus177
Konvictus177 2022년 8월 31일
편집: Konvictus177 2022년 8월 31일
Looks promising. Just tried it out quickly but my results looks a bit different than yours. Used exactly the code in your answer.
Actually each time I run the code I will get a different result.
Image Analyst
Image Analyst 2022년 8월 31일
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)

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Bartlett에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by