이번 포스팅에선 Matlab으로 구(Sphere)를 그리기 위해 몇가지 방법으로 접근해 보려한다. 첫번째 방법으로는 Matlab 내장 함수, 두번 째 매개변수 방정식을 이용하여 세번 째 DCM을 이용하여 그리는 방법으로 진행해볼 것이다.
흔히 알고있는 구의 방정식은 다음과 같다.
$$(x-x_c)^2+(y-y_c)^2+(z-z_c)^2=r^2$$
여기서 $x_c$,$y_c$, 그리고 $z_c$는 구의 중심을 의미하고, $r$은 구의 반지름을 의미한다.
이 식을 이용해 직관적으로 구를 한점한점 그릴 수도있지만 조금은 어려워 보인다. 다른 형태로서 매개변수 방정식인 다음 형태로 주어진다.
$$x = rsin\theta cos\phi + x_c$$
$$y = rsin\theta cos\phi + y_c$$
$$z = rcos\theta + z_c$$
여기서 $\theta$는 $\phi$는 매개변수이고 다음과 같은 범위를 가진다.
$$0\leq \theta\leq\pi$$
$$0\leq \phi \leq 2\pi$$
이 매개변수 방정식은 비교적 그리기 쉬워보이는 식이다.
Wave 1) Matlab 내장 함수로 구 그려보기
Matlab에서 지원되는 함수 중 sphere(n)이라는 함수가 있다. 이 함수는 "반지름이 1이고 n×n 면이 있는 구의 x, y, z 좌표를 반환합니다"설명한다. 즉 n이 커질수록 해상도가 높아져서 덜 각지고 완벽한 구체가 된다(당연히 그 만큼 데이타 양도 많아짐).
Matlab 내장함수를 이용한 구를 그리는 코드는 다음과 같다.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
|
close all;
clear all;
clc;
% 반지름이 7인 구를 그려보자
r = 7;
[X,Y,Z] = sphere(50);% 반지름이 1인 sphere에 대한 포인트(점)을
% 50에 해당하는 조밀함으로 return 해줌
Xr = X * r;
Yr = Y * r;
Zr = Z * r;
figure;
surf(Xr,Yr,Zr);
title('Draw Sphere using Matlab Function');
axis equal;
axis([-r r -r r -r r]);
xlabel('X-pos[m]');
ylabel('Y-pos[m]');
zlabel('Z-pos[m]');
|
cs |
이에대한 결과는 다음 그림과 같다.
후에 그려질 구와 달리 Surface가 이쁘게 덮혀져 있는 걸 볼 수 있는데 이는 surf함수를 통해
그렸기 때문이다.
Wave 2) 매개변수 방정식을 이용
매개변수 방정식을 이용해 그리는 코드는 다음과 같다.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
|
clear all;
close all;
clc;
deg2rad = pi/180;
r = 10;
dd = 2;
phi = (0:dd:180)'*deg2rad;
Theta = linspace(0,360,length(phi))'*deg2rad;
% Theta = (-180:10:180)'*deg2rad;
x = r*sin(phi).*cos(Theta);
y = r*sin(phi).*sin(Theta);
z = r*cos(phi);
cnt = 1;
for i=1:length(phi)
for j=1:length(Theta)
x = r*sin(phi(i,1)).*cos(Theta(j,1));
y = r*sin(phi(i,1)).*sin(Theta(j,1));
z = r*cos(phi(i,1));
sphere_pt(cnt,:) = [x,y,z];
cnt = cnt + 1;
end
end
figure;
plot3(0,0,0,'k.');
hold on; grid on;
title('Draw Sphere');
plot3(sphere_pt(:,1),sphere_pt(:,2),sphere_pt(:,3),'k.');
xlabel('X-pos[m]');
ylabel('Y-pos[m]');
zlabel('Z-pos[m]');
axis equal;
axis([-r r -r r -r r]);
|
cs |
Line 20~28 까지 2중 for문으로 구성되있는 것을 알 수 있는데 이는 구의 좌표 하나하나는 n by 3인 column vector로서 저장하기 위함이다. 이렇게 하면 각 데이터를 확인하는데 좋고 점 하나하나를 찍어볼 수 있다는 장점이 있다. 물론 변수 dd가 작아질수록 배열의 크기가 커진다.
Wave 3) DCM을 이용하여 원을 회전시켜 그리는 방법
앞서 그렸던 부분들은 어떻게든 원을 그리는데 초점이 맞춰져있었다.
하지만 인공위성의 궤도처럼 그리고싶다면 어떻게 접근해야 비슷하게 그려볼 수 있을 까?라는 생각에서 진행해본다. 먼저 y-z 평면상에서 중심이 0,0이고 어떠한 반지름 r을 가지는 원을 z축 방향으로 조금씩 회전시켜서 그려보면 어느정도 비슷할것 같고 구도 그려질거 같다.
방법은 다음과 같다.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
|
변수 초기화및 설정
구의 반지름(r), 구의 중심(sphere_center)
원을 그리기위한 각도 배열 생성
angle...
z축 회전을 위한 각도 배열 생성
z_angle...
for...
DCM(z_angle)
for...
DCM행렬을 이용한 원의 점 변환
DCM*circle_pt...
변환된 포인트 저장
res...
end
end
3차원 위치 이동
res = res + sphere_center...
결과 plotting...
|
cs |
이에대한 코드는 다음과 같다.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
|
clear all;%변수 초기화
close all;%모든 plot 화면 삭제
clc;%cmd-line clear
R = 10; % 구의 반지름 정의
shere_center = [0 0 0]'; % 구의 중심점 정의
dd_1 = 0.1;%회전될 원의 조밀함을 결정하는 변수
angle = (0:dd_1:2*pi)';%회전될 원을 그리기위한 각도 저장 변수
%회전될 원은 y-z plane에 존재
y1 = (c_x + R*cos(angle));%회전될 원의 y축 point 저장
z1= (c_y + R*sin(angle));%회전될 원의 z축 point 저장
circle_pt = [y1, z1];%회전될 원의 y,z point 통합및 저장
dd_2 = 0.1;%구를만들기 위한 원이 회전할 각도의 조밀함을 결정하는 변수
z_angle = (0:dd_2:2*pi)';%원이 회전할 각도 저장
cnt = 1;%회전변환된 원의 점을 저장하기위한 counter변수
for i=1:1:length(z_angle)
dcm_ = angle2dcm(z_angle(i),0,0,'ZYX');%z축 회전 각도에 대한 DCM 얻음
for j=1:1:length(circle_pt)%기준이되는 y-z plane에 그려진 원을 회전변환 하는 for문
sphere_pt = dcm_*[0;circle_pt(j,:)'];%원의 한 점마다 회전변환을 한다.
sphere_pt_save(cnt,:) = sphere_pt';%변환된 원의 점을 저장한다.
cnt = cnt + 1;%변환된 점이늘때마다 1씩 증가
end
end
%저장된 중심이 모두 0인 원을 원하는 원의 중심으로 3차원 위치 단순 이동
sphere_pt_save = sphere_pt_save+shere_center';
%plotting....
figure;
plot3(shere_center(1,1),sphere_pt_save(2,1),sphere_pt_save(3,1),'r.');
hold on; grid on;
plot3(sphere_pt_save(:,1),sphere_pt_save(:,2),sphere_pt_save(:,3),'k');
plot3(sphere_pt_save(:,1),sphere_pt_save(:,2),sphere_pt_save(:,3),'k.');
axis equal;
axis([-R R -R R -R R]);
title('Drawing Sphere');
xlabel('X-pos[m]');ylabel('Y-pos[m]');zlabel('Z-pos[m]');
|
cs |
결과 그림은 다음과 같다.
Extra Wave)
이번에는 이미 그려진 원의 궤적을 어떠한 점이 따라가는 걸 해볼것이다.
기본적으로 Wave3를 기반으로 하고 handler이용하여 살짝만 바꿔주면되고 코드는 다음과 같다.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
|
clear all;
close all;
clc;
R = 10;
shere_center = [0 0 0]';
dd_1 = 0.01;
angle = (0:dd_1:2*pi)';
y1 = R*cos(angle);
z1= R*sin(angle);
circle_pt = [y1, z1];
dd_2 = 0.2;
z_angle = (0:dd_2:pi)';
cnt = 1;
for i=1:1:length(z_angle)
dcm_ = angle2dcm(z_angle(i),0,0,'ZYX');
for j=1:1:length(circle_pt)
sphere_pt = dcm_*[0;circle_pt(j,:)'];
sphere_pt_save(cnt,:) = sphere_pt';
cnt = cnt + 1;
end
end
sphere_pt_save = sphere_pt_save+shere_center';
figure;
pt_handler = plot3(sphere_pt_save(1,1),sphere_pt_save(1,2),sphere_pt_save(1,3),'r.','MarkerSize',20);;
hold on; grid on;
plot3(sphere_pt_save(:,1),sphere_pt_save(:,2),sphere_pt_save(:,3),'k');
axis equal;
axis([-R R -R R -R R]);
xlabel('X-pos[m]');ylabel('Y-pos[m]');zlabel('Z-pos[m]');
title('Sphere Point Moving');
for i=2:length(sphere_pt_save)
pt_handler.XData = sphere_pt_save(i,1);
pt_handler.YData = sphere_pt_save(i,2);
pt_handler.ZData = sphere_pt_save(i,3);
drawnow();
pause(0.001);
end
|
cs |
크게 바뀐 점은 pt_handler뿐이다. plot3도 마찬가지로 단순히 점만 찍는 명령어로 쓸 수 있지만 이처럼 움직이는 점 즉 잔상없는 움직임을 보여주고 싶을 때 위처럼 handler로서 정의해서 사용하면 path와 같이 움직임을 보여줄 수 있다.
그려진 원을 회전변환하기 때문에 한 궤도를 그리고 다음 궤도에 jump가 생겨 다소 어색한 느낌이 있 목표가 궤도처럼 똑같이 그리는것이 목표는 아니라 어디까지나 Matlab의 plot, data저장 등등을 재밌게(‘•̀ ▽ •́ )φ 익히기 위한것이다. 다음 포스팅에서는 좀더 Real하게 그려보도록 하자 (*ᴗ͈ˬᴗ͈)ꕤ*.ﾟ벌써 봄이왔다..
Comments
Post a Comment