%--------------------------------------------------------------------------
% KT or J Airfoil COORDINATES AND VELOCITY DISTRIBUTION
%--------------------------------------------------------------------------
function [x, y, V,Cu] = KATJO(N,Voo,Chord,F,G,m,AlfaD)

b   = Chord /(2*m)* (1 - (F/(1+F))^m);
mb  = m*b;
xLE = mb - Chord;

BetaR = atan(G/(1+F));
AlfaR = AlfaD*pi/180;
ab  = sqrt((1+F)^2+G^2);
db  = sqrt(F^2+G^2);
da  = db/ab;
dlt = atan(G/F);

dTeta= 2*pi/N;
for i=2:N

    teta =-(i-1)*dTeta;
    arg  = G*sin(teta) - F*cos(teta);
    r  = b * (arg + sqrt(1 + 2*F + arg^2));
    xc = r*cos(teta);
    yc = r*sin(teta);
    z = complex(xc, yc);

    zeta = mb * ((z+b)^m+(z-b)^m) / ((z+b)^m-(z-b)^m);
    x(i) = real(zeta) - xLE;
    y(i) = imag(zeta);

    teta1 = teta - AlfaR - asin(da * sin(teta+dlt));
    Uz = 2*Voo * abs(sin(teta1)+sin(AlfaR+BetaR));
    dZetadz = (2*mb)^2 * (z^2-b^2)^(m-1)/((z+b)^m-(z-b)^m)^2;
    AbsdZetadz = abs(dZetadz);
    V(i) = Uz / AbsdZetadz;
    Cu(i)= V(i)^2;

end

x(1) = Chord;
y(1) = 0;
V(1) = 0;
if m==2
    V(1) = Voo / ab * cos(Alfa+beta);
end
Cu(1) = V(1)^2;

x(N+1) = x(1);
y(N+1) = y(1);
V(N+1) = V(1);
Cu(N+1) = Cu(1);

%   Aerodynamic coefficients

Gamma = 4 * pi * ab * b * sin(AlfaR + BetaR)
CL    = 2 * Gamma / Chord
CY    = CL * cos(AlfaR)
Cbb   = Chord / b 
C1    = 4 * (m * m - 1) / 3 / Cbb ^ 2 * pi * sin(2 * AlfaR)
C2    = F * cos(AlfaR) - G * sin(AlfaR)
C2    = 8 * pi * ab / (Cbb * Cbb) * C2 * sin(AlfaR + BetaR)
CMO   = C1 + C2
CMhk  = CMO - CY * (abs(xLE) / Chord)
CMc4  = CMO - CY * (abs(xLE) / Chord - 0.25)
Alfa0 = -BetaR * 180 / pi
