function [Pop,Vop,EtaAircraft,Pre,Ve,ifig]=DesignAircraft2(W,rho,S,CD0,A,e,V,ifig,R)
% function [Pop,Vop,EtaAircraft,Pre,Ve,ifig]=DesignAircraft2(W,rho,S,CD0,A,e,V,ifig,R)
% This version designs for a steady turn with radius of turn equal to 
% R=50 feet.
% OUPUTS
% Pop is operating power of the aircraft (ft-lbf/sec)
% Vop is the operating velocity of the aircraft (ft/sec)
% Pre is the minimum poewr required (ft-lbf/sec)
% Ve is the  velocity of the aircraft for minimum power required (ft/sec)
%
% INPUTS
% CD0= drag coefficient when CL=0.
% A= aspect ratio span squared divided by reference area
% e= Oswalds efficiency factor
% V= velocity vector in ft/sec
% rho= air density in slugs/ft^3
% W= aircraft weight in pounds (lbf)
% S= wing area
% ifig is the figure number of the plot
% R aircraft turn radius (feet)

tanphi=(V.*V)/(32.17*R); %tangent of bank angle Roskam eqn 4.89
phirad=atan(tanphi); % bank angle, radians
n=1./cos(phirad); % load factor, g, Roskam eqn 4.92
k=1/(pi*A*e); % k from the drag polar CD=CD0+k*CL^2
CL=2*n*W./(rho*S*V.*V); % Lift cofficient as a function of velocity
CD=CD0+k*CL.*CL; % Drag coefficient as a function of velocity

ifig=ifig+1;
figure(ifig); 
subplot(311); plot(V,phirad*180/pi)
ylabel('Bank angle, phi (deg)')
string40=['Steady turn with turn radius of ',num2str(R),' ft'];
title(string40)

subplot(312); plot(V,n)
ylabel('Load factor n (g)')

subplot(313); plot(V,CL)
ylabel('Lift coefficient')
xlabel('Velocity (ft/sec)')

Pr=.5*rho*S*V.^3.*CD; % Power required ft/lbf/sec
ifig=ifig+1; 
figure(ifig); plot(V,Pr)
xlabel('Velocity (ft/sec)')
ylabel('Power required (ft-lbf/sec)')
string41=['Aircraft Power required vs trim speed for a turn radius of ',num2str(R),' ft'];
title(string41)
%CLe=sqrt(3*CD0/k); % Lift coefficient for maximum endurance (minimum power required)
%CDe=4*CD0; % Drag coefficient for maximum endurance (minimum power required)
%Ve= sqrt(2*W/(rho*S*CLe));   % Speed for maximum endurance (minimum power required)
%Pre=interp1(V,Pr,Ve);         % YI = INTERP1(X,Y,XI) % minimum power required
ipre=find(Pr==min(Pr));
Pre=Pr(ipre);
Ve=V(ipre);
tanphie=(Ve*Ve)/(32.17*R); %tangent of bank angle Roskam eqn 4.89
phirade=atan(tanphie); % bank angle, radians
ne=1/cos(phirade); % load factor, g, Roskam eqn 4.92
CLe=2*ne*W/(rho*S*Ve*Ve);
hold on
plot([Ve,Ve],[0,Pre],':'); plot([0,Ve],[Pre,Pre],':')
string1=['Minimum power of ',num2str(Pre),' ft-lbf/sec is achieved at a speed of ',num2str(Ve),' ft/sec.'];
string1a=['The lift coefficient at this minimum power condition is ',num2str(CLe),'.'];
disp(' '); disp(string1);
text(5,Pre-1,string1)
text(5,Pre-4,string1a)
hold off 
VopInput=input('At what speed would you like to operate the aircraft? Default=Ve (speed for min. power required) >>>>');
if isempty(VopInput)
	Vop=Ve;
	Pop=Pre;
    CLop=CLe;
else
	Vop=VopInput;
	
	tanphiop=(Vop*Vop)/(32.17*R); % bank angle Roskam eqn 4.89
	phiradop=atan(tanphiop); % bank angle
	nop=1./cos(phiradop); % load factor Roskam eqn 4.92
	CLop=2*nop*W./(rho*S*Vop*Vop); % Lift cofficient as a function of velocity
	CDop=CD0+k*CLop.*CLop; % Drag coefficient as a function of velocity
	Pop=.5*rho*S*Vop^3.*CDop; % Power required ft/lbf/sec
	hold on
	plot([Vop,Vop],[0,Pop],':'); plot([0,Vop],[Pop,Pop],':');
	string3=['Operating power of ',num2str(Pop),' ft-lbf/sec is required at an operating speed of ',num2str(Vop),' ft/sec'];
	text(5,Pop+1,string3)
	hold off
end

disp(' '); disp('AIRCRAFT SUMMARY');
string2=['  You are asking the Propeller to deliver ',num2str(Pop),' ft-lbf/sec of power at a speed of ',num2str(Vop),' ft/sec.'];
EtaAircraft=Pre/Pop;
string4=['The operating efficiency of your aircraft is ',num2str(EtaAircraft)];
disp(string2); disp(string4);
hold on; text(5,Pop+3,string4); 
string1c=['The operating lift coefficient at this condition is ',num2str(CLop),'.'];
text(5,Pop+5,string1c)
hold off
disp(string1c)
disp('Is your wing capable of generating this much lift coefficient, and do you have adequate stall margin?')
disp(' ');