【问题标题】:Data evaluation [duplicate]数据评估 [重复]
【发布时间】:2019-11-05 09:20:35
【问题描述】:

我在MATLAB中有一些数据,想区分这些数据超过指定阈值时的起点和终点(例如-50),并保存它们然后计算-50下该部分的近似面积和如果它低于某个定义的值,请忽略这些点并检查接下来的两个点。见下图:

图中左侧的两点用红色标记x,所需区域用绿色显示。我想为整个图这样做。

有什么想法吗?谢谢。

【问题讨论】:

  • 所以你想计算曲线下低于给定阈值的面积?
  • 先区分出那些通过-50的点并保存,然后分别计算曲线中-50以下各部分的面积。
  • 你的问题和this one很像,已经有了答案,你应该在那里找到灵感。
  • 这能回答你的问题吗? Shade and calculate specific area

标签: matlab plot figure


【解决方案1】:

关于绘图,cmets 中提到了一些可能的方法,而我通常会使用patch 来绘制填充的多边形区域。对于面积近似,您可以使用trapz 函数进行梯形数值积分。

这将是我的解决方案,包括检测间隔,以及忽略面积不足的间隔(它有点长,并且充满了用于绘制所有间隔的循环;当然可以优化):

% Set up function, and parameter(s)
x = linspace(-0.125*pi, 4.125*pi, 10001);
y = linspace(60, 100, 10001) .* sin(x);
thr = -50;
thr_area = 30;

% Find y values lower than threshold
y_idx = find(y <= thr);

% Get start and end of intervals
idx_int = find(diff(y_idx) > 1);
n_int = numel(idx_int)+1;
s = zeros(n_int, 1);
e = zeros(n_int, 1);
s(1) = y_idx(1);
e(end) = y_idx(end);
for k = 1:n_int-1
  e(k) = y_idx(idx_int(k));
  s(k+1) = y_idx(idx_int(k)+1);
end

% Calculate areas
Q = zeros(n_int, 1);
for k = 1:n_int
  Q(k) = abs(trapz(x(s(k):e(k)), y(s(k):e(k))-thr));
end

% Visualization
figure(1);
hold on;
plot(x, y);
xlim([x(1), x(end)]);
ylim([min(y)-10, max(y)+10]);
plot([x(1), x(end)], [thr thr], 'k');
for k = 1:n_int
  patch(x(s(k):e(k)), y(s(k):e(k)), 'k');
  plot([x(s(k)), x(e(k))], [y(s(k)), y(e(k))], 'r.', 'MarkerSize', 15);
  text(x(s(k)), thr+20, num2str(Q(k)));
  if (Q(k) < thr_area)
    text(x(s(k)), thr+10, 'Area too low');
  else
    text(x(s(k)), thr+10, 'Area OK');
  end
end
hold off;

结果如下:

您现在应该拥有所有信息,可以进行任何进一步的计算、分析等。

希望有帮助!

免责声明:我使用 Octave 5.1.0 测试了代码,但我很确定它应该与 MATLAB 完全兼容。如果没有,请发表评论,我会尽力解决可能的问题。

【讨论】:

  • 感谢@HansHirse,此代码完美运行,但存在一个小问题,它也选择了 y = -53 或 y = -52。也许为此添加另一个阈值会有所帮助!
  • @FR 我不明白。所有低于阈值的y 值都会被考虑在内,是的。但是,为了找到正确的间隔开始和结束(y 跨越阈值),我采用了该“部分”的第一个和最后一个元素。要计算面积,我需要区间中的所有 y 值。
  • 是的,我注意到了。它的编码非常完美,但在我的数据中,它也选择了 -52 和 -53 的 y 值。我必须详细检查它们:) 谢谢。 @HansHirse
  • @FR 这正是我添加答案的原因。由于您的数据不包含确切的阈值,matlab 将在阈值之后选取最接近的值。
【解决方案2】:

除了@HansHirse 的回答之外,插入数据以找到阈值交叉点可能很有用。

例如,如果您的数据如下所示:

x = [ 1  2  3  4];
y = [47 49 51 53];

y 不包含确切的阈值 (50),因此我们可以对这些数据进行插值,以猜测根据x,我们将到达y = 50

x_interp = [ 1  2 2.5  3  4];
y_interp = [47 49 50  51 53];

没有插值的交叉点:

% Dummy data
x = 0:0.2:5*pi;
y = sin(x)*10;

% Threshold
T = 5;

% Crossing points
ind = find(abs(diff(sign(y-T)))==2)+1
xind = x(ind)
yind = y(ind)

% Plot
plot(x,y);
hold on
plot(xind,yind,'o','markersize',2,'color','r')

插值交叉点:

% Dummy data
x = 0:0.2:5*pi;
y = sin(x)*10;

% Threshold
T = 5;

%% Crossing points interpolation
% Index where intersection occurs
ind = [find(abs(diff(sign(y-T)))==2)+1].'+[-1,0]

% For example we could obtain:
%        [5;              [4, 5;  %We cross the threshold between y(4) and y(5)
% ind =  10;   + [-1,0] =  9,10;  %We cross the threshold between y(9) and y(10)
%        18]              17,18]  %...

xind = x(ind)
yind = y(ind)-T
% Linear interpolation
xint = xind(:,1)-yind(:,1)./(diff(yind,1,2)./diff(xind,1,2))
yint = T

% Plot
plot(x,y);
hold on
plot(xint,yint,'o','markersize',2,'color','r')

然后只需将这些新的插值添加到您的原始向量中:

[x,pos] = sort([x xint]);
y       = [y yint];
y       = y(pos);

% Now apply the @HansHirse's solution.
% ...

【讨论】:

  • 我没有安静地掌握你的代码中的 .'+[-1,0]。
  • @FR 我在我的代码中添加了一个示例,.' 是转置运算符:[column vector].' = [row vector]
  • 好的,我明白了,所以这段代码通常用于选择最接近阈值的值,对吧? @obchardon
猜你喜欢
  • 1970-01-01
  • 2012-07-25
  • 2021-03-02
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2020-03-27
  • 2017-02-19
  • 2018-01-05
相关资源
最近更新 更多