평면 표면에서 점 투영 좌표를 찾는 방법
잘 지내길 바랍니다. 두 개의 numpy 배열이 있으며 둘 다 공간의 일부 지점입니다. 파이썬을 사용하여 먼저 첫 번째 데이터 세트 ( surface_maker)를 전달하는 표면 을 찾은 다음 contact_maker생성 된 표면에서 두 번째 배열 ( ) 의 opoint에 인접한 투영의 x, y 및 z를 찾고 싶습니다 . surface_maker항상 평면형 표면을 작성했습니다. 프로젝션의 경우 인접한 지점에서 표면을 향해 수직으로 만 가고 싶습니다. 실제로 두 세트 모두에 많은 포인트가 있지만 여기에 간단한 사례를 복사합니다.
surface_maker=np.array([[50., 15., 46.04750574],
[50., 5., 45.56400925],
[44.83018398, 5., 25.],
[44.76296902, 15., 25.],
[50., 25., 45.56400925],
[44.83018398, 25., 25.],
[59.8336792, 5., 75.],
[59.71483707, 15., 75.],
[59.8336792, 25., 75.]])
contact_maker=np.array([[10., 5., 70.00014782],
[10., 15., 70.00018358],
[10., 25., 70.0001955 ],
[30., 5., 69.99981105],
[30., 15., 69.99982297],
[30., 25., 69.99985874],
[70., 5., 50.00000298],
[70., 15., 50.00002682],
[70., 25., 50.00005066],
[90., 5., 49.99996871],
[90., 15., 49.99999255],
[90., 25., 50.00001788]])
1 , 2 등과 같은 여러 솔루션을 시도했습니다 . 그러나 나는 내 문제를 성공적으로 해결했습니다. 나에게는 투영 위치를 x, y 및 z로 지정하는 것이 중요합니다. 그림은 또한 내가 원하는 것을 보여줍니다 (보여지는 것처럼으로 contact_maker생성 된 표면에 투영 된 6 개의 인접한 지점 만 필요합니다 surface_maker).
미리 도움을 주셔서 진심으로 감사드립니다.
답변
두 가지 문제를 해결해야한다는 것을 이해합니다.
- 점 집합에 맞는 평면 찾기
- 특정 방향을 따라 해당 평면에 두 번째 점 집합을 투영합니다.
두 번째 문제는 다른 답변에서 완전히 해결되었으므로 첫 번째 문제에 대한보다 일반적인 접근 방식을 제공하고 있습니다.
모든 포인트가 평면에 있다는 것을 확실하게 알고 있으면 정렬되지 않은 세 개의 포인트를 선택하고 평면을 계산할 수 있습니다. 그러나 포인트는 약간의 노이즈가있는 실제 측정에서 나올 수 있으며 포인트에 가장 잘 맞는 평면을 찾는 것이 좋습니다.
다음 함수는 점 집합에 가장 잘 맞는 평면을 찾는 일반적인 문제를 해결합니다. 주석의 설명을 참조하십시오.
import numpy as np
PRECISION = 1e-8 # Arbitrary zero for real-world purposes
def plane_from_points(points):
# The adjusted plane crosses the centroid of the point collection
centroid = np.mean(points, axis=0)
# Use SVD to calculate the principal axes of the point collection
# (eigenvectors) and their relative size (eigenvalues)
_, eigenvalues, eigenvectors = np.linalg.svd(points - centroid)
# Each eigenvalue is paired with its eigenvector and they are sorted from
# largest to smallest eigenvalue.
# The adjusted plane plane must contain the eigenvectors corresponding to
# the two largest eigenvalues. If only one eigenvector is different
# from zero, then points are aligned and they don't define a plane.
if eigenvalues[1] < PRECISION:
raise ValueError("Points are aligned, can't define a plane")
# So the plane normal is the eigenvector with the smallest eigenvalue
normal = eigenvectors[2]
# Calculate the coefficients (a,b,c,d) of the plane's equation ax+by+cz+d=0.
# The first three coefficients are given by the normal, and the fourth
# one (d) is the plane's signed distance to the origin of coordinates
d = -np.dot(centroid, normal)
plane = np.append(normal, d)
# If the smallest eigenvector is close to zero, the collection of
# points is perfectly flat. The larger the eigenvector, the less flat.
# You may wish to know this.
thickness = eigenvalues[2]
return plane, thickness
다음을 확인할 수 있습니다.
>>> surface_maker=np.array([[50., 15., 46.04750574], [50., 5., 45.56400925], [44.83018398, 5., 25.], [44.76296902, 15., 25.], [50., 25., 45.56400925], [44.83018398, 25., 25.], [59.8336792, 5., 75.], [59.71483707, 15., 75.], [59.8336792, 25., 75.]])
>>> plane, thickness = plane_from_points(surface_maker)
>>> print(plane)
[-0.95725318 0. 0.28925136 35.2806339 ]
>>> print(thickness)
1.3825669490602308
따라서 실제로 포인트 분포는 평평하지 않으며 (두께는 0과 명확하게 다름) 문제를 해결하기 위해 세 개의 임의의 포인트를 선택할 수 없습니다.
첫 번째 세트의 동일 선상이 아닌 세 점을 사용하여 평면을 만들 수 있습니다.
포인트를 A, B, C로 둡니다. 처음에는 벡터를 계산합니다.
AB = B - A (ab.x = b.x - a.x and so on)
AC = C - A
이제 외적을 사용하여 정규 벡터를 계산합니다.
N = AB x AC
N이 벡터가 0이면 점은 동일 선상에 있고 다른 삼중 항을 선택해야합니다.
(numpy에는 이러한 모든 벡터 작업에 대해 즉시 사용할 수있는 함수가 포함되어 있다고 확신합니다.)
이제 우리는 평면 방정식의 세 가지 구성 요소 (일반 구성 요소)를 가지고 있습니다.
N.x * x + N.y * y + N.z * z + D = 0
네 번째 성분 D를 얻으려면이 방정식에 A 포인트를 대입하십시오.
D = - (N.x * A.x + Ny * A.y + Nz * A.z)
투영이 OX 축을 따라있는 것 같습니다. 이 경우 모든 점 Q에 대해 평면 해석에 대한 투영을 쉽게 찾을 수 있습니다.
N.x * x + N.y * Q.y + N.z * Q.z + D = 0
x = -(N.y * Q.y + N.z * Q.z + D) / N.x
알 수없는 x의 경우 투영의 y 및 z 좌표는 Qy 및 Qz와 같습니다.
import numpy as np
S = np.array([[50., 15., 46.04750574], [50., 5., 45.56400925], [44.83018398, 5., 25.]])
AB = S[1] - S[0]
AC = S[2] - S[0]
N = np.cross(AB, AC)
D = - (N[0] * S[0][0] + N[1] * S[0][1] + N[2] * S[0][2])
Q = np.array([10., 5., 70.00014782])
x = -(N[1] * Q[1] + N[2] * Q[2] + D) / N[0]
print(x,Q[1],Q[2])
>>> 56.143273867965505 5.0 70.00014782