В MATLAB нет встроенной функции, которая будет делать то, что вы хотите. Скорее всего, вы могли бы построить свой собственный алгоритм, используя функцию INPOLYGON как , предложенную Джонасом , но вы можете самостоятельно создать более быструю реализацию, используя некоторые стандартные алгоритмы для определения, находится ли точка в пределах полигон.
Некоторое время назад я написал свой собственный код для нахождения точек пересечения между отрезком и набором треугольных поверхностей в 3-D и обнаружил, что эта ссылка для софтсерфера наиболее полезна для реализации. алгоритм. Мой случай был более сложным, чем ваш. Поскольку вы работаете в 2D, вы можете игнорировать первый раздел ссылки о поиске точки, где отрезок пересекает плоскость треугольника.
Я включил упрощенную версию моего кода MATLAB ниже для вас, чтобы использовать. Функция interpolate примет ваши матрицы faces, vertices и values в качестве входных данных и вернет дескриптор функции f, который можно оценить в заданной точке (x, y), чтобы получить кусочное значение в пределах ограничивающий треугольник. Вот несколько особенностей этого кода:
- Обработка, которая будет выполняться при оценке
f, содержится во вложенной функции evaluate_function. Эта функция имеет доступ к другим переменным в interpolate, поэтому ряд переменных, необходимых для теста в треугольнике, предварительно вычисляется, так что evaluate_function выполняется максимально быстро.
- В тех случаях, когда у вас много треугольников, проверка того, находится ли ваша точка внутри, может оказаться дорогостоящей. Следовательно, код сначала находит треугольники, которые находятся в пределах заданного радиуса (т.е. длины самого длинного участка треугольника) от вашей точки. Только эти соседние треугольники проверяются, чтобы определить, находится ли точка внутри них.
- Если точка не попадает в какую-либо треугольную область, значение NaN возвращается
f.
Есть некоторые вещи, которые не включены в код, которые вы можете добавить в зависимости от того, для чего вы его используете:
- Проверка ввода: В настоящее время в коде предполагается, что
faces является матрицей N-3, vertices является матрицей M-2, а values является вектором длины N , Возможно, вы захотите добавить проверку ошибок, чтобы убедиться, что входные данные соответствуют этим требованиям, и выдать ошибку, указывающую, когда один или несколько из них неверны.
- Проверка вырожденного треугольника: Возможно, что один или несколько треугольников, определенных вашими входами
faces и vertices, могут быть вырожденными (т.е. они могут иметь площадь 0). Это происходит, когда две или более вершин треугольников являются одной и той же точной точкой, или когда все три вершины треугольника лежат на одной прямой. Вы, вероятно, захотите добавить проверку, которая будет игнорировать такие треугольники, когда дело доходит до оценки f.
- Обработка краевых случаев: Возможно, что точка может оказаться на краю двух или более треугольных областей. Поэтому вам нужно будет решить, какое значение примет точка (то есть наибольшее из номиналов, среднее значение номиналов и т. Д.). Для таких случаев, как этот, приведенный ниже код автоматически выберет значение лица, которое ближе к началу списка лиц в вашей переменной
faces.
Наконец, вот код:
function f = interpolate(faces,vertices,values)
%# Precompute some data (helps increase speed):
triVertex = vertices(faces(:,2),:); %# Triangles main vertices
triLegLeft = vertices(faces(:,1),:)-triVertex; %# Triangles left legs
triLegRight = vertices(faces(:,3),:)-triVertex; %# Triangles right legs
C1 = sum(triLegLeft.*triLegRight,2); %# Dot product of legs
C2 = sum(triLegLeft.^2,2); %# Squared length of left leg
C3 = sum(triLegRight.^2,2); %# Squared length of right leg
triBoundary = max(C2,C3); %# Squared radius of triangle boundary
scale = C1.^2-C2.*C3;
C1 = C1./scale;
C2 = C2./scale;
C3 = C3./scale;
%# Return a function handle to the nested function:
f = @evaluate_function;
%# The nested evaluation function:
function val = evaluate_function(x,y)
w = [x-triVertex(:,1) y-triVertex(:,2)];
nearIndex = find(sum(w.^2,2) <= triBoundary); %# Find nearby triangles
if isempty(nearIndex)
val = NaN; %# Return NaN if no triangles are nearby
return
end
w = w(nearIndex,:);
wdotu = sum(w.*triLegLeft(nearIndex,:),2);
wdotv = sum(w.*triLegRight(nearIndex,:),2);
C = C1(nearIndex);
s = C.*wdotv-C3(nearIndex).*wdotu;
t = C.*wdotu-C2(nearIndex).*wdotv;
inIndex = find((s >= 0) & (t >= 0) & (s+t <= 1),1);
if isempty(inIndex)
val = NaN; %# Return NaN if point is outside all triangles
else
val = values(nearIndex(inIndex));
end
end
end