필터 지우기
필터 지우기

How to plot 3D plot from polar plot

조회 수: 21 (최근 30일)
David Harra
David Harra 2022년 8월 9일
댓글: William Rose 2022년 8월 17일
I am looking at plotting the the slowness profiles of Titinium (Inverser of velocity) by solving the christoffell equation.
I have managed to get a polar plot of this but would like to visualise this in 3d. For my polar plot I have done.
c11 = 162000000000;
c12=92000000000;
c44=47000000000;
c13= 69000000000;
c33= 181000000000;
c66= 35000000000
rho = 4430;
C=[c11,c12,c13,0,0,0;
c12,c11,c13,0,0,0;
c13,c13,c33,0,0,0;
0,0,0,c44,0,0;
0,0,0,0,c44,0;
0,0,0,0,0,c66];
phase_vel3=zeros(1,181);
r=0;
for theta=-pi:pi/180:pi
r=r+1;
n=[cos(theta),0,sin(theta)];
L_11=(C(1,1)*n(1)^2+C(6,6)*n(2)^2+C(5,5)*n(3)^2+2*C(1,6)*n(1)*n(2)+2*C(1,5)*n(1)*n(3)+2*C(5,6)*n(2)*n(3));
L_12=(C(1,6)*n(1)^2+C(2,6)*n(2)^2+C(4,5)*n(3)^2+(C(1,2)+C(6,6))*n(1)*n(2)+(C(1,4)+C(5,6))*n(1)*n(3)+(C(4,6)+C(2,5))*n(2)*n(3));
L_13=(C(1,5)*n(1)^2+C(4,6)*n(2)^2+C(3,5)*n(3)^2+(C(1,4)+C(5,6))*n(1)*n(2)+(C(1,3)+C(5,5))*n(1)*n(3)+(C(3,6)+C(4,5))*n(2)*n(3));
L_22=(C(6,6)*n(1)^2+C(2,2)*n(2)^2+C(4,4)*n(3)^2+2*C(2,6)*n(1)*n(2)+2*C(4,6)*n(1)*n(3)+2*C(2,4)*n(2)*n(3));
L_23=(C(5,6)*n(1)^2+C(2,4)*n(2)^2+C(3,4)*n(3)^2+(C(4,6)+C(2,5))*n(1)*n(2)+(C(3,6)+C(4,5))*n(1)*n(3)+(C(2,3)+C(4,4))*n(2)*n(3));
L_33=(C(5,5)*n(1)^2+C(4,4)*n(2)^2+C(3,3)*n(3)^2+2*C(4,5)*n(1)*n(2)+2*C(3,5)*n(1)*n(3)+2*C(3,4)*n(2)*n(3));
Christoffel_mat=[L_11, L_12, L_13;
L_12,L_22,L_23;
L_13,L_23,L_33];
[ev,d]=eig(Christoffel_mat);
% eigvecs = polarisation vectors
pl=ev(:,3);
% eigvals ~ slowness (phase)
phase_vel3(r)=sqrt(rho./d(3,3)); % quasi-longitudinal
r=0;
for theta=-pi:pi/180:pi
r=r+1;
n=[cos(theta),0,sin(theta)];
L_11=(C(1,1)*n(1)^2+C(6,6)*n(2)^2+C(5,5)*n(3)^2+2*C(1,6)*n(1)*n(2)+2*C(1,5)*n(1)*n(3)+2*C(5,6)*n(2)*n(3));
L_12=(C(1,6)*n(1)^2+C(2,6)*n(2)^2+C(4,5)*n(3)^2+(C(1,2)+C(6,6))*n(1)*n(2)+(C(1,4)+C(5,6))*n(1)*n(3)+(C(4,6)+C(2,5))*n(2)*n(3));
L_13=(C(1,5)*n(1)^2+C(4,6)*n(2)^2+C(3,5)*n(3)^2+(C(1,4)+C(5,6))*n(1)*n(2)+(C(1,3)+C(5,5))*n(1)*n(3)+(C(3,6)+C(4,5))*n(2)*n(3));
L_22=(C(6,6)*n(1)^2+C(2,2)*n(2)^2+C(4,4)*n(3)^2+2*C(2,6)*n(1)*n(2)+2*C(4,6)*n(1)*n(3)+2*C(2,4)*n(2)*n(3));
L_23=(C(5,6)*n(1)^2+C(2,4)*n(2)^2+C(3,4)*n(3)^2+(C(4,6)+C(2,5))*n(1)*n(2)+(C(3,6)+C(4,5))*n(1)*n(3)+(C(2,3)+C(4,4))*n(2)*n(3));
L_33=(C(5,5)*n(1)^2+C(4,4)*n(2)^2+C(3,3)*n(3)^2+2*C(4,5)*n(1)*n(2)+2*C(3,5)*n(1)*n(3)+2*C(3,4)*n(2)*n(3));
Christoffel_mat=[L_11, L_12, L_13;
L_12,L_22,L_23;
L_13,L_23,L_33];
[ev,d]=eig(Christoffel_mat);
% eigvecs = polarisation vectors
pl=ev(:,3);
% eigvals ~ slowness (phase)
phase_vel3(r)=sqrt(rho./d(3,3)); % quasi-longitudinal
figure
polar(0:pi/180:2*pi,1./phase_vel3)
Any assistance how to visualise this in 3d would be appreciated.
Thanks
Dave :)

채택된 답변

William Rose
William Rose 2022년 8월 9일
[moved this from Comment to Answer as I had originally intended]
Please post the simplest possible code fragment that clearly illustrates what you want to do. This is always a good idea, in order to raise the odss of getting a useful answer.
I cannot tell what you want to plot on the third dimension.
Matlab cannot directly plot in cylindrical or spherical coordinates, as far as I know. But it does have a built-in function to convert from cylindrical to Cartesian, which I use below.
z=0:200; theta=z*pi/25; r=1-z/400;
[x,y,z]=pol2cart(theta,r,z); %convert cylindrical to Cartesian coords
figure;
subplot(121); polar(theta,r); %polar plot
subplot(122); plot3(x,y,z); %3D plot
I hope that helps.
  댓글 수: 4
David Harra
David Harra 2022년 8월 11일
Hi William
I think that is my issue/ lack of understanding.
What you siad makes perfect sense, but I actually thought 1/v would be plotted in the z-direction as this shows the variation in velocity. I might need to revisit some of the maths etc.
Thanks for the input- its massively appreciated.
Thanks
Dave :)
William Rose
William Rose 2022년 8월 17일
@David Harra, you're welcome.

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Polar Plots에 대해 자세히 알아보기

태그

제품


릴리스

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by