def CCT_to_xy_illuminant_D(CCT):
    """
    Converts from the correlated colour temperature :math:`T_{cp}` of a
    *CIE Illuminant D Series* to the chromaticity of that
    *CIE Illuminant D Series* using *D. B. Judd, D. L. Macadam, G. Wyszecki,
    H. W. Budde, H. R. Condit, S. T. Henderson and J. L. Simonds* method.

    Parameters
    ----------
    CCT : numeric
        Correlated colour temperature :math:`T_{cp}`.

    Returns
    -------
    tuple
        *xy* chromaticity coordinates.

    Raises
    ------
    ValueError
        If the correlated colour temperature is not in appropriate domain.

    References
    ----------
    .. [12] **Wyszecki & Stiles**,
            *Color Science - Concepts and Methods Data and Formulae -
            Second Edition*,
            Wiley Classics Library Edition, published 2000,
            ISBN-10: 0-471-39918-3,
            page  145.

    Examples
    --------
    >>> CCT_to_xy_illuminant_D(6504.38938305)  # doctest: +ELLIPSIS
    (0.3127077..., 0.3291128...)
    """

    if 4000 <= CCT <= 7000:
        x = (-4.607 * 10 ** 9 / CCT ** 3 +
             2.9678 * 10 ** 6 / CCT ** 2 +
             0.09911 * 10 ** 3 / CCT +
             0.244063)
    elif 7000 < CCT <= 25000:
        x = (-2.0064 * 10 ** 9 / CCT ** 3 +
             1.9018 * 10 ** 6 / CCT ** 2 +
             0.24748 * 10 ** 3 / CCT +
             0.23704)
    else:
        raise ValueError(
            'Correlated colour temperature must be in domain [4000, 25000]!')

    y = -3 * x ** 2 + 2.87 * x - 0.275

    return x, y