% pplane2f.m 

function dz=pplane2f(F,ng,ppxlabel,ppylabel) 
% Plots the phase plane portion of SYS3PLOT.M 

global y 

cla 
ax=[min(y(:,1)), max(y(:,1)), min(y(:,2)), max(y(:,2))]; 
hold on 
x1g=linspace(ax(1),ax(2),ng); 
x2g=linspace(ax(3),ax(4),ng); 
dx1=x1g(2)-x1g(1); 
dx2=x2g(2)-x2g(1); 
x=[x1g(1); x2g(1)]; 
t=0; 

for j=1:ng 
 x(1)=x1g(j); 
 for k=1:ng 
  x(2)=x2g(k); 
  V=feval(F,t,x); 
  v=V(2)/V(1); 
  d=.3/max(1/dx1,abs(v)/dx2); 
  plot([x(1)-d,x(1)+d],[x(2)-v*d,x(2)+v*d],'g') 
 end 
end 
axis(ax) 
xlabel(ppxlabel) 
ylabel(ppylabel)