function [pathloss, PLhata, testout] = hata_davidson2(f,d,h1,h2,terrain)
% parameters
% f frequenz (30-1500MHz)
% h1 senderhöhe (HAAT) (20-2500m)
% h2 empfängerhöhe (HAAT) (1-10m)
% d distanz sender-empfänger (1-300km)
% terrain (1 oder 2)

testout = A(h1,d);

PLhata = 69.55+26.16*log10(f)-13.82*log10(h1)-a(h2,f,terrain)...
    +(44.9-6.55*log10(h1))*log10(d);


pathloss = PLhata + A(h1,d) - S1(d) - S2(h1,d) - S3(f) - S4(f,d);


function out = a(h2,f,terrain)
if terrain==1
    out = (1.1*log10(f)-0.7)*h2-(1.56*log10(f)-0.8);
else
    out = 8.29*(log10(1.54*h2))^2-1.1;
end


function out = A(h1,d)
% if d<20
%     out=0;
% else
%     out = 0.62137*(d-20)*(0.5+0.15*log10(h1/121.92));
% end
out = zeros(size(d));
out(d>=20) = 0.62137*(d(d>=20)-20)*(0.5+0.15*log10(h1/121.92));


function out = S1(d)
% if d<64.38
%     out = 0;
% else
%     out = 0.174*(d-64.38);
% end
out = zeros(size(d));
out(d>=64.38) = 0.174*(d(d>=64.38)-64.38);


function out = S2(h1,d)
if h1>300
    out=0.00784*abs(log10(9.98./d)*(h1-300));
else
    out=zeros(size(d));
end


function out = S3(f)
out = f/250*log10(1500/f);


function out = S4(f,d)
% if d>64.38
%     out = (0.112*log10(1500/f))*(d-64.38);
% else
%     out=0;
% end
out = zeros(size(d));
out(d>64.38) = (0.112*log10(1500/f))*(d(d>64.38)-64.38);