function projektioner(arg)
% Funktionen projektioner behandlar projektioner och basbyten i tre dimensioner.
% Funktionen anropas utan argument eller med ett tex. :
%      projektioner
%      projektioner(1)
% Med argument 
%      Ritar funktionen i figurfönster(2) ut x/y-planet i ett
%      x/y/z-rum, som demo hur ett plan kan ritas ut.
% Med eller utan argumet 
%     Låter funktionen Dig mata in två 3D-vektorer och
%     ritar ut i figurfönster(1) planet, som spänns av vektorerna.
%     Låter Dig upprepade gånger mata in 3D-vektorer och ritar ur
%     vektorerna med deras projektioner på planet samt residualerna.
% I funktionens kommentarer kan Du se beräkningsstegen.

% Först har vi "de vanliga basvektorerna e1, e2, e3.
% Välj ut två punkter i x-led: -2 och +2 och samma punketr i y-led.
% Gör ett rutnät, och sätt z=0 i alla punkter i rutnätet.
e1_led=linspace(-2,2,2);
[E1,E2]=meshgrid(e1_led,e1_led);dim=size(E1);
E3=zeros(dim);

% Om Du vill se en 3D-bild av planet z=0:
if nargin~=0
    figure(2)
    clf
    % Ytan med gul ytfärg, "halvt genomskinlig" och utan kantmarkering:
    surf(E1,E2,E3,'facecolor','y','facealpha',0.5,'edgecolor','none');
    xlabel('E1');ylabel('E2');zlabel('E3=0');title('E1-E2-planet');
    shg
end

v1=input('Mata in ena vektorn som en kolonnvektor: ');
v2=input('Mata in andra vektorn som en kolonnvektor: ');
% Kolla att vektorerna är linjärt oberoende:
if rank([v1,v2])~=2
    disp('Felaktiga indata: vektorerna spänner inte ett plan i 3D')
    return
end
% Normera vektorerna:
v1=v1/norm(v1);v2=v2/norm(v2);
% Beräkna en tredje vektor ortogonal mot de inmatade:
A=[v1';v2'];v3=null(A);
% Dessa tre vektorer bildar en bas i rummet
% Bilda P: "Basbytarmatrisen från basen {v1,v2,v3} till basen {e1,e2,e3}:
P=[A',v3];

% Om nu ytan [E1,E2,E3]relateras till basen {v1,v2,v3} är den en bild 
% av det plan, som spänns v1 och v2.
% Hur skildras den i basen {e1,e2,e3}? Jo skriv varje punkt på ytan som en
% kolonnvektor och multiplicera alla punkterna
% med basbytarmatrisen, se Lay kap 4.

% Skriv alla punkterna som kolonnvektorer: 
x=reshape(E1,1,dim(1)*dim(2));
y=reshape(E2,1,dim(1)*dim(2));
z=reshape(E3,1,dim(1)*dim(2));
Z=[x;y;z];

% multiplicera med basbytarmatrisen:
V=P*Z;

% Gör om dessa vektorer till punkter på en yta: 
V1=reshape(V(1,:),dim(1),dim(2)); % Alla e1-koordinater
V2=reshape(V(2,:),dim(1),dim(2)); % Alla e2-koordinater
V3=reshape(V(3,:),dim(1),dim(2)); % Alla e3-koordinater
figure(1)
clf
% Rita ut planet:
surf(V1,V2,V3,'facecolor','y','facealpha',0.5,'edgecolor','none');
hold on
% Rita ut vektorerna, som spänner planet:
quiver3(0,0,0,v1(1),v1(2),v1(3),2,'k')
quiver3(0,0,0,v2(1),v2(2),v2(3),2,'k')
% Rita ut en normal till planet:
quiver3(0,0,0,v3(1),v3(2),v3(3),2,'--g')
axis equal
xlabel('e1=[1;0;0]');ylabel('e2=[0;1;0]');zlabel('e3=[0;0;1]')
title('Planet, som spänns av Dina inmatade vektorer')
% Vill Du nu rita ut projektioner?
while 1
    val=menu('Vill Du skriva in en vektor för projektion på planet?','Ja','Nej');
    if val==2
        break
    end
    u=input('Mata in Din vektor, som skall projiceras som en kolonnvektor: ');
    % Ekvationssystemet i "Studio 1 materialet sidan 2:
    x=[dot(v1,v1),dot(v1,v2);dot(v2,v1),dot(v2,v2)]\[dot(u,v1);dot(u,v2)];
    proj=x(1)*v1+x(2)*v2;
    % Rita ut vektor, projektion och residual:
    plot3([0,proj(1),u(1),0],[0,proj(2),u(2), 0],[0,proj(3),u(3),0],'color',abs(u)/norm(u))
    
    axis auto;axis equal;grid off; shg
end
hold off