def flat_ind_zp_so3(l, m, n, b):
    """
    The SO3 spectrum consists of matrices f_hat^l of size (2l+1, 2l+1) for l=0, ..., L_max = b - 1.
    These can be stored in a zero-padded array A of shape (b, 2b, 2b) with axes l, m, n with zero padding around
     the center of the last two axes. If we flatten this array A we get a vector v of size 4b^3.
    This function gives the flat index in this array v corresponding to element (l, m, n) 

    The zero-based 3D index of (l, m, n) in A is (l, b + m, b + n).
    The corresponding flat index is i = l * 4b^2 + (b + m) * 2b + b + n

    :param l, m, n: spectral indices
    :return: index of (l, m, n) in flat zero-padded vector 
    """
    return l * 4 * (b ** 2) + (b + m) * 2 * b + b + n