clear all;
cpt=0; % un ptit compteur ...

% Definition des angles pour les spheres
n=20;
theta=(-n:2:n)/n*pi;
phi=(-n:2:n)'/n*pi/2;
cosphi=cos(phi);cosphi(1)=0;cosphi(n+1)=0;
sintheta=sin(theta);sintheta(1)=0;sintheta(n+1)=0;

x=cosphi*cos(theta);
y=cosphi*sintheta;
z=sin(phi)*ones(1,n+1);

% Definition des angles pour l onde incidente
%  du satellite
thetao=((-n:2:n)/n+3)*pi/3;%(-n:2:n)/n*pi;
phio=((-n:2:n)'/(4*n)-0.75)*pi/2;%(-pi/3:1/n:0)';%(-n:2:n)'/n*pi/2;
cosphio=cos(phio);cosphio(1)=0;cosphio(n+1)=0;
sinthetao=sin(thetao);sinthetao(1)=0;sinthetao(n+1)=0;

xo=cosphio*cos(thetao);
yo=cosphio*sinthetao;
zo=sin(phio)*ones(1,n+1);


figure(1);clf;whitebg('w');hold on
view(3);

Rinit=[-10:0.4:-8,-6:0.4:-4,-2:0.4:0];
cpt=0;
for t=(0:1:120)/10;
    clf;hold on

C2=[0.5 -1 3.5];R2=0.1;
Sat2=surf(0.1*x+C2(1),0.1*y+C2(2),0.1*z+C2(3),1*ones(size(z)));
S2=surf(R2*x+C2(1),R2*y+C2(2),R2*z+C2(3),10*ones(size(z)));


shading interp
alpha(0.4)
alpha(Sat2,1)
%camlight;lighting gouraud;
axis square equal;axis off

    R=Rinit+t;
    Rs=R.*(R>0);
    tht=2*pi/3:0.01:4*pi/3;

    for R=Rs
        po=surf(R*xo+C2(1),R*yo+C2(2),R*zo+C2(3),ones(size(zo)));
        set(po,'linewidth',5)
        shading interp
        alpha(po,0.4);
        %camlight;lighting gouraud;
    end


%% Earthmap : 
load('topo.mat','topo','topomap1');
contour(0:359,-89:90,topo,[0 0],'b')

[x,y,z] = sphere(50);

props.AmbientStrength = 0.1;
props.DiffuseStrength = 1;
props.SpecularColorReflectance = .5; 
props.SpecularExponent = 20;
props.SpecularStrength = 1;
props.FaceColor= 'texture';
props.EdgeColor = 'none';
props.FaceLighting = 'phong';
props.Cdata = topo;
surface(x,y,z,props);
%light('position',[-1 0 1]);
light('position',[0 0 1]);

light('position',[-1.5 0.5 -0.5], 'color', [.6 .2 .2]);
view(3)

%axis([-1.5 1.5 -2 2 -1 4])
axis([-1. 1. -1 1 -1 4])

cpt=cpt+1;;
M(cpt)=getframe;
%pause(0.5)
end