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 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86
| al=input('entre la valeur de alpha')
be=input('entre la valeur de betta')
b2=be*be;
bt=b2/(1+b2)/(1+b2);
v1=[ 0 10 20 30 40 50 60 70 75 80 85 90 95 100 110 120 130 140 150 160 170 180 ];
v2=[-180 -175 -170 -167 -150 -130 -110 -90 -70 -50 -30 -13 -10 -5 0 5 10 13 ...
30 50 70 90 110 130 150 167 170 175 ];
tet=0:pi/100:pi/2;
tetd=tet/pi*180;
phy=-pi:pi/100:pi;
tet=0:pi/500:pi/2;
tetd=tet/pi*180;
phy=-pi:pi/500:pi;
[p,t]=meshgrid(phy,tet);
[pp,tt]=meshgrid(phy,tet/pi*180);
[X,Y]=pol2cart(pp,tt);
kz=sin(t).*cos(p)/al;
kz2=kz.*kz;
ky=sin(t).*sin(p)/al;
ky2=ky.*ky;
kx2=(((1-kz2).*(1-b2*kz2))./((1+b2)*(1-(b2/(1+b2))*kz2)))-ky2;
kx=sqrt(kx2);
kxx=kx;
kx(kx2<0)=nan;
k2=kx2+ky2+kz2;
flag=(1-b2*kz2);
flag(flag>0)=1;
flag(flag<0)=-1;
flag(~imag(kxx)==0)=0;
flagf=flag;
flags=flag;
flagf(flag<0)=0;
flags(flag>0)=0;
% invkx2=1./kx2;
%value=abs(mean(mean(kx2)))+std(std(kx2));
%l=abs(kx2)>value;
%l2=1-l;
%kxp=kx2+invkx2;
theta = linspace(0,2*pi,100);
r_max = 90;
h_fake = polar(theta,r_max*ones(size(theta)));
hold on;
set(h_fake, 'Visible');
h=polar([-pi pi],[0 90])
delete (h)
hold on
contour(X,Y,kx2,[0 0],'b')
k2(~imag(kxx)==0)=nan;
cgx=((1+b2)/2)*kx.*(1-(k2-2*bt*kz2)./sqrt(k2.*k2-4*bt*k2.*kz2));
cgx(~imag(kxx)==0)=nan;
cgy=((1+b2)/2)*ky.*(1-(k2-2*bt*kz2)./sqrt(k2.*k2-4*bt*k2.*kz2));
cgy(~imag(kxx)==0)=nan;
cgz=((1+b2)/2)*kz.*(1-(k2-2*bt*(kz2+k2))./sqrt(k2.*k2-4*bt*k2.*kz2));
cgz(~imag(kxx)==0)=nan;
cgx1=((1+b2)/2)*kx.*(1+(k2-2*bt*kz2)./sqrt(k2.*k2-4*bt*k2.*kz2));
cgx1(~imag(kxx)==0)=nan;
cgy1=((1+b2)/2)*ky.*(1+(k2-2*bt*kz2)./sqrt(k2.*k2-4*bt*k2.*kz2));
cgy1(~imag(kxx)==0)=nan;
cgz1=((1+b2)/2)*kz.*(1+(k2-2*bt*(kz2+k2))./sqrt(k2.*k2-4*bt*k2.*kz2));
cgz1(~imag(kxx)==0)=nan;
phyy=atan2(cgy,cgz)*180/pi;
phyys=phyy;
phyyf=phyy;
phyyf(flagf==0)=nan;
phyys(flags==0)=nan;
phyy1=atan2(cgy1,cgz1)*180/pi;
phyys1=phyy1;
phyyf1=phyy1;
phyyf1(flagf==0)=nan;
phyys1(flags==0)=nan;
[C,h]=contour(X,Y,phyys,v2,'b');
[C,h]=contour(X,Y,phyyf1,v1,'b');
[C,h]=contour(X,Y,phyyf1,-v1,'b');
clear |
Partager