def average(da, dim=None, weights=None):
    """
        Returns the weighted average
        Author: Dougie Squire
        Date: 06/08/2018

        Parameters
        ----------
        da : xarray DataArray
            Array to be averaged
        dim : str or sequence of str, optional
            Dimension(s) over which to compute weighted average. If None, average is computed over all
            dimensions
        weights : xarray DataArray, optional
            Weights to apply during averaging. Shape of weights must be broadcastable to shape of da.
            If None, unity weighting is applied
            
        Returns
        -------
        weighted : xarray DataArray
            Weighted average of input array along specified dimensions
            
        Examples
        --------
        >>> A = xr.DataArray(np.random.normal(size=(4,4)), coords=[('lat', np.arange(-90,90,45)), 
        ...                                                        ('lon', np.arange(0,360,90))])
        >>> degtorad = doppyo.utils.constants().pi / 180
        >>> cos_lat = xr.ufuncs.cos(A['lat'] * degtorad) 
        >>> doppyo.utils.average(A, dim='lat', weights=cos_lat)
        <xarray.DataArray (lon: 4)>
        array([-0.473632, -0.241208, -0.954826,  0.498559])
        Coordinates:
          * lon      (lon) int64 0 90 180 270
    """
    if weights is None:
        return da.mean(dim)
    else:
        weights = (0 * da + 1) * weights
        return (da * weights).sum(dim) / weights.sum(dim)
