2014-06-08 5 views
0

안녕하세요. 내 첫 번째 게시물입니다. Delaunay 삼각 측량을위한 matlab 스크립트를 작성하고 싶습니다. 여기 내 스크립트가 있습니다 :Delaunay 삼각 측량을 matlab에 구현했습니다.

clear all;clc 
%% Delaunay 
x=[ 160.1671 366.9226 430.7894 540.1208 660.2771 508.7287 252.1787]; 
y=[ 223.9615 259.5000 120.5769 245.5000 283.1923 472.7308 469.5000]; 

% 
x=x'; 
y=y'; 

%orginal plot 
dd=delaunay(x,y); 
dt=TriRep(dd,x,y); 
triplot(dt); 
z=[x.^2+y.^2] 
i=1:length(x); 
ptk=[i' x y] 
%% main loop 
l=0; 
for i=1:length(x)-2 
for j=1+i:length(x) 
    for k=1+i:length(x) 
     if (j ~= k) 
      l=l+1; 
      xn = (y(j)-y(i))*(z(k)-z(i)) - (y(k)-y(i))*(z(j)-z(i)); 
      yn = (x(k)-x(i))*(z(j)-z(i)) - (x(j)-x(i))*(z(k)-z(i)); 
      zn = (x(j)-x(i))*(y(k)-y(i)) - (x(k)-x(i))*(y(j)-y(i)); 
       if (zn < 0) 
        border=zn; 
         for m=1:length(x) 
          border = (border) & ... 
           ((x(m)-x(i))*xn +... 
           (y(m)-y(i))*yn +... 
           (z(m)-z(i))*zn <= 0); 
          if (border) 
          ii(m)=[i]; 
          jj(m)=[j]; 
          kk(m)=[k]; 
          end 
         end 
       end 
     end 
    end 
end 
end 
wart=[ii' jj' kk'] 
dd 
figure(2) 
triplot(wart,x,y) 

이것은이 스크립트에서 얻은 것입니다. 이 행렬은 들로네() MATLAB 함수로 생성됩니다

dd = 
6  7  2 
7  1  2 
4  6  2 
1  3  2 
4  3  5 
6  4  5 
2  3  4 

이 내가 구현에서 무엇을 얻을 수 있습니다 :

wart = 
4  7  6 
4  7  5 
4  7  5 
4  7  5 
4  7  5 
4  6  5 
4  6  5 

당신의 사람이 잘못이 무엇인지 말해 줄 수 있을까요? 실수는 어디에서 왔습니까?

jils.

+0

설명은 의미 론적 수준의 코드를 이해하는 데 필수적입니다. 예 : border를 음의 double 값으로 초기화하지만 논리 값으로 사용됩니다. 'border = false'를 초기화하는 것과 같은 효과가 있습니다. 무엇을 의도합니까? – Daniel

+0

나는 그렇게 생각한다. 나는이 C 프로그램의 matlab 스크립트를 작성 중이다 : http://fatcat.ftj.agh.edu.pl/~kaprzyk/Geomkomp/geomkomp/node58.html 그러나이 스크립트는 지연 일 함수와 C 프로그램 사이에 차이가있다. 유사한 구현이 Java로 제공됩니다. http://stackoverflow.com/questions/5825089/how-does-this-code-for-delaunay-triangulation-work – jilsu

+0

@jiksu :'border = (zn <0) ; '는 C 코드의 올바른 번역입니다. – Daniel

답변

1

문제는 여기에, 당신의 가장 안쪽의 루프에 있습니다

if (zn < 0) 
       border=zn; 
        for m=1:length(x) 
         border = (border) & ... 
          ((x(m)-x(i))*xn +... 
          (y(m)-y(i))*yn +... 
          (z(m)-z(i))*zn <= 0); 
         if (border) 
         ii(m)=[i]; 
         jj(m)=[j]; 
         kk(m)=[k]; 
         end 
        end 
      end 

당신은 점에 의해 정의 된 삼각형 [난, J, K]가 유효한지 확인합니다. 외곽선 안쪽에 점이없는 m 인 경우에만 3 점을 출력에 저장하고 싶습니까? 현재 체크하는 첫 번째 점이 circumcircle 바깥에 있으면, 세 점은 무엇이든 관계없이 저장됩니다. 또한 각 가능한 삼각형에 대해 동일한 m을 반복하므로 올바른 값을 찾더라도 나중에 겹쳐 쓸 가능성이 있습니다. 이것은 출력에서 ​​같은 값을 반복하는 이유이기도합니다.

이러한 경우에는 무엇이 일어나는지 루프를 통해 (정신적으로 수동으로 명령 행에서 또는 디버그 메소드를 사용하여) 단계별로 수행 할 가치가 있습니다. 첫 번째 출력 4 7 6은 inbuilt 함수 결과에 나타나지 않습니다. 따라서 i,j,k을 해당 값으로 설정하고 해당 내부 루프에서 어떤 일이 발생하는지 확인하십시오.

덧붙여서 실제로 거기에는 루프가 필요 없습니다. 같은 작업을 수행하여 한 번에 모든 값을 확인 :

border = (x-x(i)).*xn + (y-y(i)).*yn + (z-z(i)).*zn; 
if all(border<0) 
    % then store coordinates 
end 
당신은 빈 출력 ([])로 시작하고 ( end+1 사용) 추가, 또는 삼각형의 최대 수를 계산하고 그 크기에 출력을 미리 할당 할 수

, 카운터 변수를 사용하여 몇 개를 찾고 출력 배열의 올바른 위치에 놓았는지 확인한 다음 끝에 출력을 오른쪽 크기로 자릅니다. 더 큰 입력 데이터 세트를 사용하려는 경우 미리 할당하는 것이 좋습니다.

+0

그것은 매우 도움이되었습니다. 나는 당신의 지시에 따라 삼각형을 얻었습니다. 두 번째 질문은 필자가 필요로하지 않는 내부 라인을 제거하는 방법입니다. 여기에 나는 아이디어가 없다. – jilsu