% [CL_de]=CL_de(S_w,S_h,AR_h,ce_ch,eta_oe,eta_ie,beta,kappa_h,lambda_h,Lambda_c2_h,tc_h,delta_e,Cl_alpha_h)
%
%  CL_de - Variation of airplane lift coefficient with  
%           elevator deflection angle
%
% Version 3.0 2/15/02. 
% This version corrects Cm_alpha, rudder control power, units, and comments.
% AAE565 Spring 2002
%
%
%	[Ref.2]	Jan Roskam
%		Methods for Estimating
%		Stability and Control Derivatives
%		of Conventional Subsonic Airplanes  -  1977
%
%		S_w			- Reference area - main wing [ft^2]
%		S_h 			- Surface area - horizontal tail [ft^2]
%		AR_h			- Aspect Ratio of horizontal tail    
%		ce_ch			- Elevator flap chord to horizontal tail chord ratio 
%		eta_oe		- Percent span position of outboard edge of elevator
%		eta_ie		- Percent span position of inboard edge of elevator
%		beta			- Compressibility correction factor (sqrt(1-Mach^2))
%		kappa_h		- Ratio of 2D lift curve slope to 2pi
%		Lambda_c2_h - Horizontal tail sweep angle at 1/2chord [rad] 
%		lambda_h		- Horizontal tail taper ratio
%		tc_h			- Thickness to chord ratio of horizontal tail
%		delta_e 		- Elevator deflection [rad]
%

function [CL_de]=CL_de(S_w,S_h,AR_h,ce_ch,eta_oe,eta_ie,beta,kappa_h, ...
                       lambda_h,Lambda_c2_h,tc_h,delta_e,Cl_alpha_h)
                    
CL_alpha_h = 2*pi*AR_h/(2+sqrt((AR_h*beta/kappa_h)^2*(1+(tan(Lambda_c2_h))^2/beta^2)+4 ));%eqn 3.8 Ref 2

% Figure 10.2  -  Reference 2
cf_c_fig = [0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1];
alpha_delta_cl_fig = [0	-0.394 -0.553 -0.667 -0.758 -0.829 -0.884 -0.93 -0.968 -0.994 -1];
alpha_delta_cl=interlim1(cf_c_fig,alpha_delta_cl_fig,ce_ch);
AR_fig = [1 2 4 6 8 10];
alpha_delta_cl_fig = [-0.1 -0.2 -0.3 -0.4 -0.5 -0.6 -0.7 -0.8 -0.9 -1.0];
alpha_delta_ratio_fig = [1.2	1.25	1.33	1.46	1.78	2.2	;...
		1.14	1.16	1.2	1.3	1.49	1.71	;...	
		1.1	1.11	1.14	1.2	1.36	1.52	;...
		1.06	1.08	1.1	1.14	1.25	1.39	;...
		1.05	1.06	1.07	1.1	1.18	1.28	;...
		1.03	1.04	1.05	1.07	1.13	1.2	;...
		1.02	1.03	1.03	1.05	1.09	1.14	;...
		1.01	1.01	1.02	1.03	1.06	1.09	;...
		1.01	1.01	1.01	1.02	1.03	1.05	;...
		1	1	1	1	1	1	];
alpha_delta_ratio=interlim2(AR_fig,alpha_delta_cl_fig,alpha_delta_ratio_fig,AR_h,alpha_delta_cl);

% Figure 10.3  -  Reference 2
eta_fig = [0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1];
lambda_fig = [0 0.5 1.0];
K_b_fig = [0	0.149	0.3	0.433	0.555	0.669	0.768	0.857	0.933	0.986	1	;...	
		0	0.132	0.263	0.393	0.514	0.624	0.729	0.829	0.916	0.975	1	;...
      0	0.12	0.246	0.368	0.483	0.598	0.701	0.799	0.888	0.957	1] ;	
K_bo=interlim2(eta_fig,lambda_fig,K_b_fig,eta_oe,lambda_h);
K_bi=interlim2(eta_fig,lambda_fig,K_b_fig,eta_ie,lambda_h);
K_b=K_bo-K_bi;

% Figure 10.5  -  Reference 2
cf_c_fig = [0.05 0.1 0.15 0.2 0.25 0.3 0.4 0.5];
t_c_fig = [0.15 0.12 0.10 0.08 0.06 0.04 0.02 0];
Cl_delta_theory_fig=[1.81	2.57	3.09	3.44	3.82	4.61	5.33	5.96 	;...
		1.81	2.55	3.11	3.48	3.88	4.54	5.21	5.8	;...
		1.8	2.53	3.13	3.57	3.97	4.46	5.12	5.7	;...
		1.8	2.51	3.16	3.52	4.03	4.41	5.03	5.57	;...
		1.79	2.5	3.19	3.6	4.08	4.35	4.95	5.47	;...
		1.79	2.49	3		3.63	4.14	4.29	4.86	5.36	;...
		1.78	2.48	3.04	3.73	4.19	4.23	4.78	5.24	;...
      1.78	2.47	3.07	3.69	3.93	4.16	4.7	5.14	];
Cl_delta_theory=interlim2(cf_c_fig,t_c_fig,Cl_delta_theory_fig,ce_ch,tc_h);

% Figure 10.6  -  Reference 2
cf_c_fig = [0.05 0.1 0.15 0.2 0.25 0.3 0.4 0.5];
Cl_alpha_ratio_fig = [1 .98 .96 .94 .92 .9 .88 .86 .84 .82 .8 .78 .76 .74 .72 .7];
Cl_delta_ratio_fig=[1	1	1	1	1	1	1	1	;...	
		0.952	0.961	0.965	0.965	0.966	0.968	0.971	0.973	;...	
		0.91	0.92	0.928	0.932	0.935	0.939	0.944	0.95	;...
		0.863	0.882	0.894	0.9	0.905	0.909	0.917	0.923	;...
		0.825	0.846	0.858	0.865	0.871	0.878	0.891	0.901	;...
		0.782	0.809	0.824	0.835	0.846	0.85	0.863	0.875	;...
		0.739	0.771	0.786	0.8	0.81	0.819	0.835	0.848	;...
		0.699	0.733	0.751	0.765	0.778	0.789	0.805	0.818	;...
		0.656	0.695	0.714	0.731	0.744	0.756	0.777	0.793	;...
		0.616	0.654	0.673	0.691	0.706	0.72	0.742	0.762	;...
		0.574	0.613	0.636	0.655	0.671	0.684	0.709	0.729	;...
		0.527	0.566	0.59	0.611	0.627	0.642	0.67	0.693	;...
		0.483	0.52	0.547	0.567	0.587	0.602	0.631	0.656	;...
		0.441	0.471	0.497	0.521	0.543	0.562	0.593	0.616	;...
		0.398	0.427	0.452	0.476	0.497	0.515	0.549	0.582	;...
      0.354	0.384	0.408	0.431	0.453	0.471	0.509	0.545	];
Cl_delta_ratio=interlim2(cf_c_fig,Cl_alpha_ratio_fig,Cl_delta_ratio_fig,ce_ch,kappa_h);
 
% Figure 10.7  -  Reference 2
delta_f_fig = [0 10 15 20 25 30 35 40 45 50 55 60]*pi/180;
%changed to radians as input will be in radians - Murthy
cf_c_fig = [0.1 0.15 0.2 0.25 0.3 0.4 0.5];
K_prime_fig = [1	1	0.982	0.901	0.787	0.716	0.674	0.642	0.618	0.596	0.578	0.563	;...	
		1	1	0.982	0.889	0.769	0.697	0.653	0.62	0.591	0.569	0.549	0.531	;...	
		1	1	0.982	0.878	0.743	0.674	0.63	0.597	0.568	0.544	0.522	0.502	;...
		1	1	0.982	0.859	0.704	0.644	0.602	0.57	0.541	0.518	0.498	0.481	;...
		1	1	0.981	0.817	0.669	0.611	0.573	0.543	0.517	0.496	0.477	0.46	;...
		1	1	0.951	0.75	0.633	0.578	0.543	0.515	0.492	0.473	0.455	0.439	;...
		1	1	0.918	0.696	0.597	0.546	0.514	0.49	0.469	0.453	0.437	0.423	];	
K_prime=interlim2(delta_f_fig,cf_c_fig,K_prime_fig,delta_e,ce_ch);

Cl_alpha_M = Cl_alpha_h/beta;				                  		% Eqn 10.2 Ref2
Cl_df = Cl_delta_ratio*Cl_delta_theory*K_prime;  					% Eqn 10.3 Ref2
CL_df = Cl_df*CL_alpha_h/Cl_alpha_M*alpha_delta_ratio*K_b; 		% Eqn 10.1 Ref2
CL_de = CL_df*S_h/S_w;                     							% Eqn 10.6 Ref2

return
