% 14-09-03, Some special plots

% area.m  -----------------------------------------------------------------------------------
  
  X=[0:100]';
  Y=exp(  -(X-50).^2/(sqrt(2*pi)*20) );
  area(X,Y)

  disp('  -------- area.m ---------') 
  disp('  X=[0:100]''; Y=exp(  -(X-50).^2/(sqrt(2*pi)*20) );')
  disp('  area(X,Y)')

  disp('  Pause - Press any key'); pause; disp('  ')
  
% bar.m ------------------------------------------------------------------------------------------   
 
  Year = [2000 2001 2002 2003];
  Cost = [0.1 0.2 0.3 ; 0.2 0.22 0.15 ; 0.11 0.12 0.23 ; 0.05 0.02 0.01 ];
  bar(Year,Cost,'stacked')
  legend('1','2','3')
  xlabel('Year')
  ylabel('Total cost')
  
  disp('  -------- bar.m ---------') 
  disp('  Year = [2000 2001 2002 2003];')
  disp('  Cost = [0.1 0.2 0.3 ; 0.2 0.22 0.15 ; 0.11 0.12 0.23 ; 0.05 0.02 0.01 ];')
  disp('  bar(Year,Cost,''stacked'')')
  disp('  legend(''1'',''2'',''3'')')
  disp('  xlabel(''Year'')')
  disp('  ylabel(''Total cost'')')

  disp('  Pause - Press any key'); pause; disp('  ')

  
% feather.m -----------------------------------------------------------------------------------------
 
  theta = (-90:10:90)*pi/180; r = 2*ones(size(theta));
  [u,v] = pol2cart(theta,r);
  feather(u,v)
  axis([0 20 -4 4])

  disp('  -------- feather.m ---------') 
  disp('  theta = (-90:10:90)*pi/180;')
  disp('  r = 2*ones(size(theta));')
  disp('  [u,v] = pol2cart(theta,r);')
  disp('  feather(u,v)')
  disp('  axis([0 20 -4 4])')

  disp('  Pause - Press any key'); pause; disp('  ')


% quiver.m -----------------------------------------------------------------------------------------
 
  [x,y] = meshgrid(-2:.2:2,-1:.15:1);
  z = x .* exp(-x.^2 - y.^2); [px,py] = gradient(z,.2,.15);
  contour(x,y,z), hold on
  quiver(x,y,px,py), hold off, axis image
         
  disp('  -------- quiver.m ---------') 
    
  disp('  [x,y] = meshgrid(-2:.2:2,-1:.15:1);')
  disp('  z = x .* exp(-x.^2 - y.^2); [px,py] = gradient(z,.2,.15);')
  disp('  contour(x,y,z), hold on')
  disp('  quiver(x,y,px,py), hold off, axis image')

  disp('  Pause - Press any key'); pause; disp('  ')

  
% quiver3.m -----------------------------------------------------------------------------------------
 
  [x,y] = meshgrid(-2:.2:2,-1:.15:1);
  z = x .* exp(-x.^2 - y.^2);
  [u,v,w] = surfnorm(x,y,z);
  quiver3(x,y,z,u,v,w); hold on, surf(x,y,z), hold off
  xlabel('X');
  ylabel('Y');
  axis([-2 2 -2 2 -0.5 0.7]);

  disp('  -------- quiver3.m ---------') 
    
  disp('  [x,y] = meshgrid(-2:.2:2,-1:.15:1);')
  disp('  z = x .* exp(-x.^2 - y.^2);')
  disp('  [u,v,w] = surfnorm(x,y,z);')
  disp('  quiver3(x,y,z,u,v,w); hold on, surf(x,y,z), hold off')
  disp('  xlabel(''X'');');
  disp('  ylabel(''Y'');');
  disp('  axis([-2 2 -2 2 -0.5 0.7]);')

  disp('  Pause - Press any key'); pause; disp('  ')

% scatter.m -----------------------------------------------------------------------------------------
 
  load seamount
  scatter(x,y,5,z)
 
  disp('  -------- scatter.m ---------') 
  disp('  load seamount')
  disp('  scatter(x,y,5,z)')

  disp('  Pause - Press any key'); pause; disp('  ')

% scatter3.m -----------------------------------------------------------------------------------------
 
  [x,y,z] = sphere(16);
  X = [x(:)*.5 x(:)*.75 x(:)];
  Y = [y(:)*.5 y(:)*.75 y(:)];
  Z = [z(:)*.5 z(:)*.75 z(:)];
  S = repmat([1 .75 .5]*10,prod(size(x)),1);
  C = repmat([1 2 3],prod(size(x)),1);
  scatter3(X(:),Y(:),Z(:),S(:),C(:),'filled'), view(-60,60)
       
  disp('  -------- scatter3.m ---------') 
   
  disp('  [x,y,z] = sphere(16);')
  disp('  X = [x(:)*.5 x(:)*.75 x(:)];')
  disp('  Y = [y(:)*.5 y(:)*.75 y(:)];')
  disp('  Z = [z(:)*.5 z(:)*.75 z(:)];')
  disp('  S = repmat([1 .75 .5]*10,prod(size(x)),1);')
  disp('  C = repmat([1 2 3],prod(size(x)),1);')
  disp('  scatter3(X(:),Y(:),Z(:),S(:),C(:),''filled''), view(-60,60)')
  
  disp('  Pause - Press any key'); pause; disp('  ')

% plotmatrix.m -----------------------------------------------------------------------------------------
  
  x = randn(50,3); y = x*[-1 2 1;2 0 1;1 -2 3;]';
  plotmatrix(y)
       
  disp('  -------- plotmatrix.m ---------') 
  
  disp('  x = randn(50,3); y = x*[-1 2 1;2 0 1;1 -2 3;]'';')
  disp('  plotmatrix(y)')

  disp('  Pause - Press any key'); pause; disp('  ')

% plot3.m -----------------------------------------------------------------------------------------
  
  t = 0:pi/50:10*pi;
  plot3(sin(t),cos(t),t);
  
  disp('  -------- plot3.m ---------') 

  disp('  t = 0:pi/50:10*pi;')
  disp('  plot3(sin(t),cos(t),t);')

  disp('  Pause - Press any key'); pause; disp('  ')

  
  
  
% subplots ---------------------------------------------------------------------------------------
  
% subplot dimension
n1 = 2; % number of rows
n2 = 3; % number of columns
nw = 0.9; % normalized width
nh = 0.9; % normalized height
for k1 = 1:n1
    for k2 = 1:n2
        subplot(n1,n2,(k1-1)*n2 + k2,...
            'position', [(1-nw)/n2/2 + (k2-1)/n2, (1-nh)/n1/2 + 1-k1/n1, nw/n2 nh/n1]);
        plot(rand(5));
        set(gca, 'XTick', []);    set(gca, 'YTick', []);
    end
end

disp('  % subplot dimension')
disp('  n1 = 2; % number of rows')
disp('  n2 = 3; % number of columns')
disp('  % These values would define the space between the graphs')
disp('  % if equal to 1 there will be no space between graphs')
disp('  nw = 0.9; % normalized width')
disp('  nh = 0.9; % normalized height')
disp('  for k1 = 1:n1')
disp('      for k2 = 1:n2')
disp('          subplot(n1,n2,(k1-1)*n2 + k2,...')
disp('              ''position'', [(1-nw)/n2/2 + (k2-1)/n2, (1-nh)/n1/2 + 1-k1/n1,...')
disp('              nw/n2 nh/n1]);')
disp('          plot(rand(5));')
disp('          % turn off the labels if you want')
disp('          set(gca, ''XTick'', []);')
disp('          set(gca, ''YTick'', []);')
disp('      end')
disp('  end')

disp('  Pause - Press any key'); pause; disp('  ')
clf


% contourslice.m ------------------------------------------------------------------------------------

  [x y z v] = flow;
  h=contourslice(x,y,z,v,[1:9],[],[0], linspace(-8,2,10));
  axis([0 10 -3 3 -3 3]); daspect([1 1 1])
  camva(24); camproj perspective;
  campos([-3 -15 5])
  set(gcf, 'Color', [.3 .3 .3], 'renderer', 'zbuffer')
  set(gca, 'Color', 'black' , 'XColor', 'white', ...
     'YColor', 'white' , 'ZColor', 'white')
  box on

  disp('  -------- contourslice.m  ---------') 
  disp('  [x y z v] = flow;') 
  disp('  h=contourslice(x,y,z,v,[1:9],[],[0], linspace(-8,2,10));')
  disp('  axis([0 10 -3 3 -3 3]); daspect([1 1 1])')
  disp('  camva(24); camproj perspective;')
  disp('  campos([-3 -15 5])')
  disp('  set(gcf, ''Color'', [.3 .3 .3], ''renderer'', ''zbuffer'')')
  disp('  set(gca, ''Color'', ''black'' , ''XColor'', ''white'', ...')
  disp('     ''YColor'', ''white'' , ''ZColor'', ''white'')')
  disp('   box on')

  disp('  Pause - Press any key'); pause; disp('  ')

% slice.m ------------------------------------------------------------------------------------

  
  disp('  -------- slice.m  ---------') 

  disp('  [x,y,z] = meshgrid(-2:.2:2, -2:.25:2, -2:.16:2);')
  disp('  v = x .* exp(-x.^2 - y.^2 - z.^2);')
  disp('  slice(x,y,z,v,[-1.2 .8 2],2,[-2 -.2])')

  [x,y,z] = meshgrid(-2:.2:2, -2:.25:2, -2:.16:2);
  v = x .* exp(-x.^2 - y.^2 - z.^2);
  slice(x,y,z,v,[-1.2 .8 2],2,[-2 -.2])

  disp('  Pause - Press any key'); pause; disp('  ')

  clf
% streamparticles.m ------------------------------------------------------------------------------------

  disp('  -------- streamparticles.m  Example 1 ---------') 

  load wind
  [sx sy sz] = meshgrid(80, 20:1:55, 5);
  verts = stream3(x,y,z,u,v,w,sx,sy,sz);
  sl = streamline(verts);
  iverts = interpstreamspeed(x,y,z,u,v,w,verts,.025);
  axis tight; view(30,30); daspect([1 1 .125])
  %set(gca, 'drawmode', 'fast')
  camproj perspective; box on
  camva(44); camlookat; camdolly(0,0,.4, 'f');
  h = line; 
  streamparticles(h, iverts, 35, 'animate', 10, ...
     'ParticleAlignment', 'on');
  

  
  disp('  Pause - Press any key'); pause; disp('  ')
  clf
% streamparticles.m ------------------------------------------------------------------------------------

  disp('  -------- streamparticles.m  Example 2 ---------') 

  load wind
  daspect([1 1 1]); view(2)
  [verts averts] = streamslice(x,y,z,u,v,w,[],[],[5]); 
  sl = streamline([verts averts]);
  axis tight off;
  set(sl, 'linewidth', 2, 'color', 'r', 'vis', 'off')
  iverts = interpstreamspeed(x,y,z,u,v,w,verts,.05);
  %set(gca,'drawmode','fast','position',...
  %   [0 0 1 1],'zlim',[4.9 5.1]);
  set(gca,'position',[0 0 1 1],'zlim',[4.9 5.1]);
  set(gcf, 'color', 'k')
  h = line; 
  streamparticles(h, iverts, 200, ...
     'animate', 10, 'framerate',40, ...
     'markers', 10, 'markerf', 'y');  
  

  disp('  Pause - Press any key'); pause; disp('  ')
  clf
  
% avifile ------------------------------------------------------------------------------------
%{
  fig=figure;
  set(fig,'DoubleBuffer','on');
  set(gca,'xlim',[-80 80],'ylim',[-80 80],...
     'NextPlot','replace','Visible','off')
  mov = avifile('example.avi')
  x = -pi:.1:pi;
  radius = [0:length(x)];
  for i=1:length(x)
     h = patch(sin(x)*radius(i),cos(x)*radius(i),[abs(cos(x(i))) 0 0]);
     set(h,'EraseMode','xor');
     F = getframe(gca);
     mov = addframe(mov,F);
  end
  mov = close(mov);
  
  disp('  -------- avifile.m ---------') 

  disp('  Pause - Press any key'); pause; disp('  ')
  clf
%}  
% avifile ------------------------------------------------------------------------------------

  
