def solar_irradiation(latitude, longitude, Z, moment, surface_tilt, 
                      surface_azimuth, T=None, P=None, solar_constant=1366.1,
                      atmos_refract=0.5667, albedo=0.25, linke_turbidity=None, 
                      extraradiation_method='spencer',
                      airmass_model='kastenyoung1989',
                      cache=None):
    r'''Calculates the amount of solar radiation and radiation reflected back
    the atmosphere which hits a surface at a specified tilt, and facing a
    specified azimuth. 
    
    This functions is a wrapper for the incredibly 
    comprehensive `pvlib library <https://github.com/pvlib/pvlib-python>`_, 
    and requires it to be installed.
    
    Parameters
    ----------
    latitude : float
        Latitude, between -90 and 90 [degrees]
    longitude : float
        Longitude, between -180 and 180, [degrees]
    Z : float, optional
        Elevation above sea level for the position, [m]
    moment : datetime
        Time and date for the calculation, in local UTC time (not daylight 
        savings time), [-]
    surface_tilt : float
        The angle above the horizontal of the object being hit by radiation,
        [degrees]
    surface_azimuth : float
        The angle the object is facing (positive North eastwards 0° to 360°),
        [degrees]
    T : float, optional
        Temperature of atmosphere at ground level, [K]
    P : float, optional
        Pressure of atmosphere at ground level, [Pa]
    solar_constant : float, optional
        The amount of solar radiation which reaches earth's disk (at a 
        standardized distance of 1 AU); this constant is independent of 
        activity or conditions on earth, but will vary throughout the sun's
        lifetime and may increase or decrease slightly due to solar activity,
        [W/m^2]
    atmos_refract : float, optional
        Atmospheric refractivity at sunrise/sunset (0.5667 deg is an often used
        value; this varies substantially and has an impact of a few minutes on 
        when sunrise and sunset is), [degrees]
    albedo : float, optional
        The average amount of reflection of the terrain surrounding the object
        at quite a distance; this impacts how much sunlight reflected off the
        ground, gest reflected back off clouds, [-]
    linke_turbidity : float, optional
        The amount of pollution/water in the sky versus a perfect clear sky;
        If not specified, this will be retrieved from a historical grid;
        typical values are 3 for cloudy, and 7 for severe pollution around a
        city, [-]
    extraradiation_method : str, optional
        The specified method to calculate the effect of earth's position on the
        amount of radiation which reaches earth according to the methods 
        available in the `pvlib` library, [-]
    airmass_model : str, optional
        The specified method to calculate the amount of air the sunlight
        needs to travel through to reach the earth according to the methods 
        available in the `pvlib` library, [-]
    cache : dict, optional
        Dictionary to to check for values to use to skip some calculations;
        `apparent_zenith`, `zenith`, `azimuth` supported, [-]

    Returns
    -------
    poa_global : float
        The total irradiance in the plane of the surface, [W/m^2]
    poa_direct : float
        The total beam irradiance in the plane of the surface, [W/m^2]
    poa_diffuse : float
        The total diffuse irradiance in the plane of the surface, [W/m^2]
    poa_sky_diffuse : float
        The sky component of the diffuse irradiance, excluding the impact 
        from the ground, [W/m^2]
    poa_ground_diffuse : float
        The ground-sky diffuse irradiance component, [W/m^2]

    Examples
    --------
    >>> solar_irradiation(Z=1100.0, latitude=51.0486, longitude=-114.07, 
    ... moment=datetime(2018, 4, 15, 13, 43, 5), surface_tilt=41.0, 
    ... surface_azimuth=180.0)
    (1065.7621896280812, 945.2656564506323, 120.49653317744884, 95.31535344213178, 25.181179735317063)
    
    >>> cache = {'apparent_zenith': 41.099082295767545, 'zenith': 41.11285376417578, 'azimuth': 182.5631874250523}
    >>> solar_irradiation(Z=1100.0, latitude=51.0486, longitude=-114.07, 
    ... moment=datetime(2018, 4, 15, 13, 43, 5), surface_tilt=41.0, 
    ... linke_turbidity=3, T=300, P=1E5,
    ... surface_azimuth=180.0, cache=cache)
    (1042.5677703677097, 918.2377548545295, 124.33001551318027, 99.6228657378363, 24.70714977534396)

    At night, there is no solar radiation and this function returns zeros:
        
    >>> solar_irradiation(Z=1100.0, latitude=51.0486, longitude=-114.07, 
    ... moment=datetime(2018, 4, 15, 2, 43, 5), surface_tilt=41.0, 
    ... surface_azimuth=180.0)
    (0.0, -0.0, 0.0, 0.0, 0.0)
    

    Notes
    -----    
    The retrieval of `linke_turbidity` requires the pytables library (and 
    Pandas); if it is not installed, specify a value of `linke_turbidity` to 
    avoid the dependency.
    
    There is some redundancy of the calculated results, according to the 
    following relations. The total irradiance is normally that desired for 
    engineering calculations.
    
    poa_diffuse = poa_ground_diffuse + poa_sky_diffuse
    
    poa_global = poa_direct + poa_diffuse
    
    FOr a surface such as a pipe or vessel, an approach would be to split it
    into a number of rectangles and sum up the radiation absorbed by each.
    
    This calculation is fairly slow. 

    References
    ----------
    .. [1] Will Holmgren, Calama-Consulting, Tony Lorenzo, Uwe Krien, bmu, 
       DaCoEx, mayudong, et al. Pvlib/Pvlib-Python: 0.5.1. Zenodo, 2017. 
       https://doi.org/10.5281/zenodo.1016425.
    '''
    # Atmospheric refraction at sunrise/sunset (0.5667 deg is an often used value)
    from fluids.optional import spa
    from fluids.optional.irradiance import (get_relative_airmass, get_absolute_airmass,
                                            ineichen, get_relative_airmass, 
                                            get_absolute_airmass, get_total_irradiance)
#    try:
#        import pvlib
#    except:
#        raise ImportError(PVLIB_MISSING_MSG)

    moment_timetuple = moment.timetuple()
    moment_arg_dni = (moment_timetuple.tm_yday if 
                      extraradiation_method == 'spencer' else moment)

    dni_extra = _get_extra_radiation_shim(moment_arg_dni, solar_constant=solar_constant, 
                               method=extraradiation_method, 
                               epoch_year=moment.year)
    
    if T is None or P is None:
        atmosphere = ATMOSPHERE_NRLMSISE00(Z=Z, latitude=latitude, 
                                           longitude=longitude, 
                                           day=moment_timetuple.tm_yday)
        if T is None:
            T = atmosphere.T
        if P is None:
            P = atmosphere.P
    
    if cache is not None and 'zenith' in cache:
        zenith = cache['zenith']
        apparent_zenith = cache['apparent_zenith']
        azimuth = cache['azimuth']
    else:
        apparent_zenith, zenith, _, _, azimuth, _ = solar_position(moment=moment,
                                                                   latitude=latitude, 
                                                                   longitude=longitude,
                                                                   Z=Z, T=T, P=P, 
                                                                   atmos_refract=atmos_refract)
    
    if linke_turbidity is None:
        from pvlib.clearsky import lookup_linke_turbidity
        import pandas as pd
        linke_turbidity = float(lookup_linke_turbidity(
            pd.DatetimeIndex([moment]), latitude, longitude).values)

        
    if airmass_model in apparent_zenith_airmass_models:
        used_zenith = apparent_zenith
    elif airmass_model in true_zenith_airmass_models:
        used_zenith = zenith
    else:
        raise Exception('Unrecognized airmass model')
    
    relative_airmass = get_relative_airmass(used_zenith, model=airmass_model)
    airmass_absolute = get_absolute_airmass(relative_airmass, pressure=P)


    ans = ineichen(apparent_zenith=apparent_zenith,
                   airmass_absolute=airmass_absolute, 
                   linke_turbidity=linke_turbidity,
                   altitude=Z, dni_extra=solar_constant, perez_enhancement=True)
    ghi = ans['ghi']
    dni = ans['dni']
    dhi = ans['dhi']
    
    
#    from pvlib.irradiance import get_total_irradiance
    ans = get_total_irradiance(surface_tilt=surface_tilt, 
                      surface_azimuth=surface_azimuth,
                      solar_zenith=apparent_zenith, solar_azimuth=azimuth,
                      dni=dni, ghi=ghi, dhi=dhi, dni_extra=dni_extra, 
                      airmass=airmass_absolute, albedo=albedo)
    poa_global = float(ans['poa_global'])
    poa_direct = float(ans['poa_direct'])
    poa_diffuse = float(ans['poa_diffuse'])
    poa_sky_diffuse = float(ans['poa_sky_diffuse'])
    poa_ground_diffuse = float(ans['poa_ground_diffuse'])
    return (poa_global, poa_direct, poa_diffuse, poa_sky_diffuse, 
            poa_ground_diffuse)
r'''Calculates the amount of solar radiation and radiation reflected back
    the atmosphere which hits a surface at a specified tilt, and facing a
    specified azimuth. 
    
    This functions is a wrapper for the incredibly 
    comprehensive `pvlib library <https://github.com/pvlib/pvlib-python>`_, 
    and requires it to be installed.
    
    Parameters
    ----------
    latitude : float
        Latitude, between -90 and 90 [degrees]
    longitude : float
        Longitude, between -180 and 180, [degrees]
    Z : float, optional
        Elevation above sea level for the position, [m]
    moment : datetime
        Time and date for the calculation, in local UTC time (not daylight 
        savings time), [-]
    surface_tilt : float
        The angle above the horizontal of the object being hit by radiation,
        [degrees]
    surface_azimuth : float
        The angle the object is facing (positive North eastwards 0° to 360°),
        [degrees]
    T : float, optional
        Temperature of atmosphere at ground level, [K]
    P : float, optional
        Pressure of atmosphere at ground level, [Pa]
    solar_constant : float, optional
        The amount of solar radiation which reaches earth's disk (at a 
        standardized distance of 1 AU); this constant is independent of 
        activity or conditions on earth, but will vary throughout the sun's
        lifetime and may increase or decrease slightly due to solar activity,
        [W/m^2]
    atmos_refract : float, optional
        Atmospheric refractivity at sunrise/sunset (0.5667 deg is an often used
        value; this varies substantially and has an impact of a few minutes on 
        when sunrise and sunset is), [degrees]
    albedo : float, optional
        The average amount of reflection of the terrain surrounding the object
        at quite a distance; this impacts how much sunlight reflected off the
        ground, gest reflected back off clouds, [-]
    linke_turbidity : float, optional
        The amount of pollution/water in the sky versus a perfect clear sky;
        If not specified, this will be retrieved from a historical grid;
        typical values are 3 for cloudy, and 7 for severe pollution around a
        city, [-]
    extraradiation_method : str, optional
        The specified method to calculate the effect of earth's position on the
        amount of radiation which reaches earth according to the methods 
        available in the `pvlib` library, [-]
    airmass_model : str, optional
        The specified method to calculate the amount of air the sunlight
        needs to travel through to reach the earth according to the methods 
        available in the `pvlib` library, [-]
    cache : dict, optional
        Dictionary to to check for values to use to skip some calculations;
        `apparent_zenith`, `zenith`, `azimuth` supported, [-]

    Returns
    -------
    poa_global : float
        The total irradiance in the plane of the surface, [W/m^2]
    poa_direct : float
        The total beam irradiance in the plane of the surface, [W/m^2]
    poa_diffuse : float
        The total diffuse irradiance in the plane of the surface, [W/m^2]
    poa_sky_diffuse : float
        The sky component of the diffuse irradiance, excluding the impact 
        from the ground, [W/m^2]
    poa_ground_diffuse : float
        The ground-sky diffuse irradiance component, [W/m^2]

    Examples
    --------
    >>> solar_irradiation(Z=1100.0, latitude=51.0486, longitude=-114.07, 
    ... moment=datetime(2018, 4, 15, 13, 43, 5), surface_tilt=41.0, 
    ... surface_azimuth=180.0)
    (1065.7621896280812, 945.2656564506323, 120.49653317744884, 95.31535344213178, 25.181179735317063)
    
    >>> cache = {'apparent_zenith': 41.099082295767545, 'zenith': 41.11285376417578, 'azimuth': 182.5631874250523}
    >>> solar_irradiation(Z=1100.0, latitude=51.0486, longitude=-114.07, 
    ... moment=datetime(2018, 4, 15, 13, 43, 5), surface_tilt=41.0, 
    ... linke_turbidity=3, T=300, P=1E5,
    ... surface_azimuth=180.0, cache=cache)
    (1042.5677703677097, 918.2377548545295, 124.33001551318027, 99.6228657378363, 24.70714977534396)

    At night, there is no solar radiation and this function returns zeros:
        
    >>> solar_irradiation(Z=1100.0, latitude=51.0486, longitude=-114.07, 
    ... moment=datetime(2018, 4, 15, 2, 43, 5), surface_tilt=41.0, 
    ... surface_azimuth=180.0)
    (0.0, -0.0, 0.0, 0.0, 0.0)
    

    Notes
    -----    
    The retrieval of `linke_turbidity` requires the pytables library (and 
    Pandas); if it is not installed, specify a value of `linke_turbidity` to 
    avoid the dependency.
    
    There is some redundancy of the calculated results, according to the 
    following relations. The total irradiance is normally that desired for 
    engineering calculations.
    
    poa_diffuse = poa_ground_diffuse + poa_sky_diffuse
    
    poa_global = poa_direct + poa_diffuse
    
    FOr a surface such as a pipe or vessel, an approach would be to split it
    into a number of rectangles and sum up the radiation absorbed by each.
    
    This calculation is fairly slow. 

    References
    ----------
    .. [1] Will Holmgren, Calama-Consulting, Tony Lorenzo, Uwe Krien, bmu, 
       DaCoEx, mayudong, et al. Pvlib/Pvlib-Python: 0.5.1. Zenodo, 2017. 
       https://doi.org/10.5281/zenodo.1016425.