求助：matlab里球面坐標系下的非圓齒輪漸開線命令

愛無止境

我想根據一篇文獻里的關系式重現齒輪圓錐曲線，前面公式及命令如下：
% I, I& Z& ]$ M" p  b0 L" N' mfi1=[0:0.1:2*pi];
0 L9 z: x$ }8 |! yc=0.16;
R=43.6;
$ }3 d2 W; Q( _8 ~+ y4 v, qfi2=0.75*fi1+(c/3)*0.75*cos(3*fi1+3*pi/2);      %0.75=z1/z2=3/4,我直接換成0.75代入關系式了。
7 A; V* h) |  M' ~4 h" hdelta1=atan(0.75*(1-c*sin(3*fi1+3*pi/2)));
delta2=acot(0.75*(1-c*sin(3*fi1+3*pi/2)));
- A1 r( ?. i" @
$ H' {4 r1 w9 u: ?$ y球面方程如下：
x1=R*(sin(delta1).*cos(fi1));
6 G/ @) H+ j5 _7 |" T8 G4 h1 oy1=R*sin(delta1).*sin(fi1);
: ~2 L( F% O& G! Z" X0 ez1=R*cos(delta1);
r1=x1*i+y1*j+z1*k     %這里的i,j,k前的系數是節曲線上的點在三個坐標上的坐標。寫到程序里面的不是這樣的，下面也是一樣。這樣寫出來大家理解直觀一些。
$ i9 c2 A) F/ Q. F" v* ^x2=R*(sin(delta2).*cos(fi2));
y2=R*sin(delta2).*sin(fi2);
) b- C2 t1 y  H2 E% A9 v' nz2=R*cos(delta2);
r2=x2*i+y2*j+z2*k
: v" o# |% |5 N  m! @: [
r1和r2的值驗算出來的都是43.6，沒有問題。但是我是不知道用什么樣的命令得到下面的球面圖形。由于本人對matlab不很了解，到這里用命令老是出錯，各位英雄女俠幫幫忙。感激不盡！！！


% V5 H8 m, e0 B* w