
%NBBBBB CAMBIARE I SEGNI DI TUTTO IN MODO CHE SINAO CONCORDI CON XYZ
%STABILITA IN VERSO
format long
alldata_1=datcomimport('usato2.out')


alpha=alldata_1{1}.alpha;
alpha_rad=alpha.*pi./180;

%FLAP%%%%%
deltaflap= alldata_1{1}.delta;
deltacl= alldata_1{1}.dcl_sym;
deltacd=alldata_1{1}.dcdmin_sym;
cdindotta=alldata_1{1}.dcdi_sym;

%derivata Cmdf
deltacm=alldata_1{1}.dcm_sym;% incremnto cm dovuto flap
deltacm=[deltacm(1) deltacm(2) deltacm(3)];

%%%%%%%% VELIVOLO %%%%%%

clfinale=alldata_1{3}.cl;
cdfinale=alldata_1{3}.cd;
cm=alldata_1{3}.cm;

%componenti coeff. di forza in assi datum
cxpuro=alldata_1{3}.ca;
czpuro=alldata_1{3}.cn;

% coeff forza in assi bodi centrati in datum 
CX_velivolo=-cxpuro;
CZ_velivolo=-czpuro;

% derivata Clr media peche non unico valore
Clr=alldata_1{3}.clr;
Clr_fin=sum(Clr)/length(alpha)

% derivata Cnr media peche non unico valore
Cnr=alldata_1{3}.cnr;
Cnr_fin=sum(Cnr)/length(alpha)

%derivata Clp media peche non unico valore
Clp=alldata_1{3}.clp;
Clp_fin=sum(Clp)/length(alpha)


% derivata Cnp media peche non unico valore
Cnp=alldata_1{3}.cnp;
Cnp_fin=sum(Cnp)/length(alpha)


% derivata Cyp media peche non unico valore
Cyp=alldata_1{3}.cyp;
Cyp_fin=sum(Cyp)/length(alpha)

%derivata Clb media perche non unico valore
Clb=alldata_1{3}.clb;
Clb_fin=sum(Clb)/length(alpha)

%EQILIBRATORE%%%%%
deltae=(alldata_1{3}.delta).*pi./180;% deflessione equilibratore in rad
%derivata Cmde
delta_e_cm=alldata_1{3}.dcm_sym; %incremento cm dovuto equilibratore



%%%% aileron clroll delta a 
Clroll_aileron=alldata_1{2}.clroll; %incrementi sul cl roll dovuti a escursioni aileron 
delta_aileron=(alldata_1{2}.deltar).*pi./180;%escursioni aileron 

%%%% cn delta a 
cn_aileron=alldata_1{2}.cn_asy; % variazione di cn causato al variare di aileron e incidenza 

for i=1:1:length(alpha)
    if(alpha(i)==0)
        t=i;
    cn_aileron_aplha_null=cn_aileron(t,:);
  end 
end 

%%%% COMPONENTI DI FORZA CALCOLATI%%%%
%aumento cl dovuto flap e variazione alpha
cl_flap0=deltacl(1);
cl_flap1=deltacl(2);
cl_flap2=deltacl(3);

%aumento cd dovuto flap e variazione alpha e cd indotta
cd_flap0=cdindotta(:,1);
cd_flap1=cdindotta(:,2);
cd_flap2=cdindotta(:,3);




%cx e xz dovuti flapp 0
cl= cl_flap0;
cd= cd_flap0;
i=1;
for k=1:1:length(alpha)
    alpha1=alpha(k);
cx0(i)=-cl*sind(alpha1)+cd(i)*cosd(alpha1);
cz0(i)=cl*cosd(alpha1)+cd(i)*sind(alpha1);
i=i+1;
end



%cx e xz dovuti flapp 1
cl= cl_flap1;
cd= cd_flap1;
i=1;
for k=1:1:length(alpha)
    alpha1=alpha(k);
cx1(i)=-cl*sind(alpha1)+cd(i)*cosd(alpha1);
cz1(i)=cl*cosd(alpha1)+cd(i)*sind(alpha1);
i=i+1;
end

%cx e xz dovuti flapp 2
cl= cl_flap2;
cd= cd_flap2;
i=1;
for k=1:1:length(alpha)
    alpha1=alpha(k);
cx2(i)=-cl*sind(alpha1)+cd(i)*cosd(alpha1);
cz2(i)=cl*cosd(alpha1)+cd(i)*sind(alpha1);
i=i+1;
end


% cambaindo di segno ai coeff di forza si ha il segno che rispecchia xyz
% body

cx0=(-1)*cx0; 
cx1=(-1)*cx1; 
cx2=(-1)*cx2;
cz0=(-1)*cz0;
cz1=(-1)*cz1;
cz2=(-1)*cz2;



Xf=[-cxpuro cx0' cx1'  cx2' alpha']
Zf=[-czpuro cz0' cz1'  cz2' alpha']




%calcolo Cxdf e Xzdf
%NNNNBBB CX dleta flap Cz delta flap contempla il caso non contemplato Cx e
%Cz alpha*delta flap, ossia cosidero Cx e Cz valori a alpha=0

escusione_flap=[0 deltaflap(2)  deltaflap(3)];
escusione_flap_rad=escusione_flap.*pi./180;

escusione_flap_rad20= deltaflap.*pi./180;

%trovo valore per incidenza nulla
for i=1:1:length(alpha)
    if(alpha(i)==0)
     forza_x_alpha_null=[cx0(i) cx1(i)  cx2(i)];
     forza_z_alpha_null=  [ cz0(i) cz1(i)   cz2(i)];
  end 
end 



%calcolo Cxadf e Xzadf

adf0=0.*alpha_rad;% moltiplicazione delta flap0 per incidenze 
adf10=escusione_flap_rad(2).*alpha_rad;
adf30=escusione_flap_rad(3).*alpha_rad;

af=[ adf10 adf30 ];
CXaf=[ cx1 cx2];
CZAf=[ cz1 cz2 ];




%%%FORZE EQUILIBRATORE %%%
%Incremento dovuto deflessioen eq
delta_e_cl= alldata_1{3}.dcl_sym;
delta_e_cd=alldata_1{3}.dcdmin_sym;
cd_idotta_eq=alldata_1{3}.dcdi_sym

% considero solo 2 valori di escursione eq escluso zero e copreso max
% escursioni positivi e negative 



for j=1:1:length(deltae)
cl= delta_e_cl(j)
cd=cd_idotta_eq(:,j);
i=1;
for k=1:1:length(alpha)
    alpha1=alpha(k);
cze(j,i)=cl*cosd(alpha1)+cd(i)*sind(alpha1);
i=i+1;
end
end

%concorde con assi bodi 
cze=-cze



%NNNNBBB  Cz delta eq contempla il caso per cui volgio variazione di z
%dovuto a flap e non interessa l'incidenza, già tutto concorde con xyz stab


for i=1:1:length(alpha)
    if(alpha(i)==0)
    forza_zee_alpha_null= [cze(:,i)];
  end 
end 




