public static ObservedPosition jauAtioq(double ri, double di, Astrom astrom)
    {
        /* Minimum cos(alt) and sin(alt) for refraction purposes */
        final double CELMIN = 1e-6;
        final double SELMIN = 0.05;

        double v[] = new double[3], x, y, z, xhd, yhd, zhd, f, xhdt, yhdt, zhdt,
                xaet, yaet, zaet, azobs, r, tz, w, del, cosdel,
                xaeo, yaeo, zaeo, zdobs, hmobs, dcobs, raobs;

        /*--------------------------------------------------------------------*/

        /* CIRS RA,Dec to Cartesian -HA,Dec. */
        v = jauS2c(ri-astrom.eral, di);
        x = v[0];
        y = v[1];
        z = v[2];

        /* Polar motion. */
        xhd = x + astrom.xpl*z;
        yhd = y - astrom.ypl*z;
        zhd = z - astrom.xpl*x + astrom.ypl*y;

        /* Diurnal aberration. */
        f = ( 1.0 - astrom.diurab*yhd );
        xhdt = f * xhd;
        yhdt = f * ( yhd + astrom.diurab );
        zhdt = f * zhd;

        /* Cartesian -HA,Dec to Cartesian Az,El (S=0,E=90). */
        xaet = astrom.sphi*xhdt - astrom.cphi*zhdt;
        yaet = yhdt;
        zaet = astrom.cphi*xhdt + astrom.sphi*zhdt;

        /* Azimuth (N=0,E=90). */
        azobs = ( xaet != 0.0 || yaet != 0.0 ) ? atan2(yaet,-xaet) : 0.0;

        /* ---------- */
        /* Refraction */
        /* ---------- */

        /* Cosine and sine of altitude, with precautions. */
        r = sqrt(xaet*xaet + yaet*yaet);
        r = r > CELMIN ? r : CELMIN;
        z = zaet > SELMIN ? zaet : SELMIN;

        /* A*tan(z)+B*tan^3(z) model, with Newton-Raphson correction. */
        tz = r/z;
        w = astrom.refb*tz*tz;
        del = ( astrom.refa + w ) * tz /
                ( 1.0 + ( astrom.refa + 3.0*w ) / ( z*z ) );

        /* Apply the change, giving observed vector. */
        cosdel = 1.0 - del*del/2.0;
        f = cosdel - del*z/r;
        xaeo = xaet*f;
        yaeo = yaet*f;
        zaeo = cosdel*zaet + del*r;

        /* Observed ZD. */
        zdobs = atan2(sqrt(xaeo*xaeo+yaeo*yaeo), zaeo);

        /* Az/El vector to HA,Dec vector (both right-handed). */
        v[0] = astrom.sphi*xaeo + astrom.cphi*zaeo;
        v[1] = yaeo;
        v[2] = - astrom.cphi*xaeo + astrom.sphi*zaeo;

        /* To spherical -HA,Dec. */
        SphericalCoordinate co = jauC2s ( v);
        hmobs = co.alpha;
        dcobs = co.delta;
        /* Right ascension (with respect to CIO). */
        raobs = astrom.eral + hmobs;

        /* Return the results. */
        return new ObservedPosition(
        jauAnp(azobs),
        zdobs,
        -hmobs,
        dcobs,
        jauAnp(raobs));

        /* Finished. */


    }