Оценка данных

У меня есть некоторые данные в MATLAB, и я хочу различать начальную и конечную точки, когда эти данные пересекают указанный порог (например, -50), и сохранять их, а затем вычислять приблизительную площадь этого раздела под -50 и если она была ниже определенного значения пренебречь этими пунктами и проверить следующие два пункта. См. следующее изображение:

Рисунок

Две точки в левой части рисунка отмечены x красным, а необходимая область показана зеленым. Я хочу сделать это для всей фигуры.

Есть идеи? Спасибо.


person F R    schedule 05.11.2019    source источник
comment
Итак, вы хотите рассчитать площадь под кривой, которая ниже заданного порога?   -  person obchardon    schedule 05.11.2019
comment
Сначала выделить те точки, которые переходят -50, и сохранить их, а затем вычислить площадь каждого участка ниже -50 отдельно на кривой.   -  person F R    schedule 05.11.2019
comment
Ваш вопрос очень похож на этот вопрос, на который уже есть ответ, вы должны найти вдохновение там.   -  person Hoki    schedule 05.11.2019
comment
Отвечает ли это на ваш вопрос? Затенить и вычислить определенную область   -  person Hoki    schedule 05.11.2019


Ответы (2)


Что касается построения графиков, в комментариях упоминались возможные способы, тогда как я обычно использовал 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. Если нет, оставьте комментарий, и я постараюсь исправить возможные проблемы.

person HansHirse    schedule 05.11.2019
comment
Спасибо @HansHirse, этот код отлично работает, но есть небольшая проблема с тем, что он также выбирает y = -53 или y = -52. Возможно, добавление другого порога для этого поможет! - person F R; 05.11.2019
comment
@FR Я не понимаю. Все значения y ниже порога учитываются, да. Но чтобы найти правильное начало и конец интервала (y пересечение порога), я беру первый и последний элемент этого раздела. Для вычисления площади мне нужны все y значения из интервала. - person HansHirse; 05.11.2019
comment
Да, я это заметил. Он отлично закодирован, но в моих данных он также выбирает те значения y, которые равны -52 и -53. Я должен проверить их подробно :) спасибо. @ХансХирс - person F R; 05.11.2019
comment
@FR именно поэтому я добавил свой ответ. Поскольку ваши данные не содержат точного порогового значения, Matlab подберет ближайшее значение после порога. - person obchardon; 05.11.2019

В дополнение к ответу @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.
% ...
person obchardon    schedule 05.11.2019
comment
Я не понял .'+[-1,0] в вашем коде. - person F R; 05.11.2019
comment
@FR Я добавил пример в свой код, .' - это оператор транспонирования: [column vector].' = [row vector]. - person obchardon; 05.11.2019
comment
Хорошо, я вижу, так что этот код в целом предназначен для выбора ближайших значений к порогу, верно? @obchardon - person F R; 05.11.2019