
% Generating patterns by plane waves superposition

% wave vectors
k1=[1 1 ; 1 0 ; 1 0 ;  1 0 ] ;
k2=[0 0 ; 0 1 ; -0.5 1/sqrt(2) ; -0.5 1/sqrt(2)];
k3=[0 0 ; 0 0 ; -0.5 -1/sqrt(2) ; -0.5 -1/sqrt(2) ];

% their phases
phi1=[0 0 0 1*pi/3];
phi2=[0 0 0 1*pi/3];
phi3=[0 0 0 1*pi/3];

% coeff for linear combo
coef=[1 0 0 ; 1 1 0 ; 1 1 1 ; 1 1 1];

% level line to make countour
cutplane=[0.99 0 -1 1];

% domain
R=12;
L=-R;
Npts=100;
x=linspace(L,R,Npts);
y=linspace(L,R,Npts);
[X,Y]=meshgrid(x,y);

fs=16;
figure(1)
clf

N=length(coef);

for jj=1:N
  U=real(coef(jj,1)*exp(i*(k1(jj,1)*X+k1(jj,2)*Y+phi1(jj)))+...
         coef(jj,2)*exp(i*(k2(jj,1)*X+k2(jj,2)*Y+phi2(jj)))+...
         coef(jj,3)*exp(i*(k3(jj,1)*X+k3(jj,2)*Y+phi3(jj))));

  subplot(4,N,jj)
  set(gca,'FontSize',[fs]);
  surf(X,Y,U);
  shading interp
  %colormap(gray)
  xlabel('x'); ylabel('y'); zlabel('z');axis tight
  view([12 46])

  subplot(4,N,N+jj)
  set(gca,'FontSize',[fs]);
  surf(X,Y,U);
  shading interp
  colormap(jet)
  view([0 90])
  xlabel('x'); ylabel('y'); zlabel('z');axis tight

  subplot(4,N,2*N+jj)
  set(gca,'FontSize',[fs]);
%  v=-1;
%  v=max(max(U))-0.0025
  v=cutplane(jj);
  contour(X,Y,U,[v v],'b');
  xlabel('x'); ylabel('y'); zlabel('z');axis tight

  subplot(4,N,3*N+jj)
  set(gca,'FontSize',[fs]);
  vects=[k1(jj,:);k2(jj,:);k3(jj,:)]';
  plotv(vects,'-')
  xlabel('x'); ylabel('y'); zlabel('z');axis tight;axis equal
  ax=axis; 
  axis([ax(1)-0.1*(ax(2)-ax(1)) ax(2)+0.1*(ax(2)-ax(1)) ax(3)-0.1*(ax(4)-ax(3)) ax(4)+0.1*(ax(4)-ax(3))]);

  if(jj<N)
  fprintf('Press any key to continue...\n');
  pause
  end
end

