#function Test()
E = wgs84Ellipsoid();
%% reference inputs
az = 33; el=70; srange = 1e3;
lat = 42; lon= -82; alt = 200;
%% reference outputs
er = 186.277521; nr = 286.84222; ur = 939.69262; % aer2enu
xl = 660.93019e3; yl = -4701.42422e3; zl = 4246.5796e3; % aer2ecef
lat1 = 42.0026; lon1 = -81.9978; alt1 = 1.1397e3; % aer2geodetic
x0 = 660.67525e3; y0 = -4700.94868e3; z0 = 4245.73766e3; % geodetic2ecef, ecef2geodetic
angleUnit="d";
%% 
    %% aer2ecef contains:
    [x1,y1,z1] = geodetic2ecef(E,lat,lon,alt, angleUnit);
    assert([x1,y1,z1],[x0,y0,z0],10e-2)

    [e1,n1,u1] = aer2enu(az, el, srange, angleUnit);
    assert([e1,n1,u1], [er,nr,ur],10e-2)

    [x2,y2,z2] = aer2ecef(az,el,srange,lat,lon,alt,E, angleUnit);
    assert([x2,y2,z2], [xl,yl,zl],10e-2)

    %% ecef2geodetic is self-contained iterative algorithm.
    [lat2, lon2, alt2] = ecef2geodetic(E, x1, y1, z1, angleUnit); % round-trip
    assert([lat2, lon2, alt2], [lat, lon, alt],10e-2)

    [az2, el2, rng2] = enu2aer(e1,n1,u1, angleUnit); % round-trip
    assert([az2,el2,rng2],[az,el,srange],10e-2)

    [az3, el3, rng3] = ecef2aer(x2,y2,z2, lat,lon,alt, E, angleUnit); % round-trip
    assert([az3,el3,rng3], [az,el,srange],10e-2)


    [lat3,lon3,alt3] = aer2geodetic(az,el,srange,lat,lon,alt, E, angleUnit);
    assert([lat3,lon3,alt3], [lat1, lon1, alt1],10e-2)

    [e2, n2, u2] = geodetic2enu(lat3, lon3, alt3, lat, lon, alt, E, angleUnit);
    assert([e2,n2,u2],[e1,n1,u1],10e-2)

    [az4, el4, rng4] = geodetic2aer(lat3,lon3,alt3,lat,lon,alt, E, angleUnit); % round-trip
    assert([az4,el4,rng4], [az,el,srange],10e-2)
    %%
    [x3, y3, z3] = enu2ecef(e1,n1,u1,lat,lon,alt, E, angleUnit);
    assert([x3,y3,z3],[x2,y2,z2],10e-2)

    [lat4, lon4, alt4] = enu2geodetic(e2,n2,u2,lat,lon,alt, E, angleUnit); % round-trip
    assert([lat4, lon4, alt4],[lat3, lon3, alt3],10e-2)

    [e3,n3,u3] = ecef2enu(x3,y3,z3,lat,lon,alt, E, angleUnit); % round-trip
    assert([e3,n3,u3],[e1,n1,u1],10e-2)
