How to use the pyriemann.utils.base.expm function in pyriemann

def test_expm():
    """Test matrix exponential"""
    C = 2*np.eye(3)
    Ctrue = np.exp(2)*np.eye(3)
    assert_array_almost_equal(expm(C), Ctrue)

def geodesic_logeuclid(A, B, alpha=0.5):
    """Return the matrix at the position alpha on the log euclidean geodesic between A and B  :

    .. math::
            \mathbf{C} =  \exp \left( (1-\\alpha) \log(\mathbf{A}) + \\alpha \log(\mathbf{B}) \\right)

    C is equal to A if alpha = 0 and B if alpha = 1

    :param A: the first coavriance matrix
    :param B: the second coavriance matrix
    :param alpha: the position on the geodesic
    :returns: the covariance matrix

    """
    return expm((1 - alpha) * logm(A) + alpha * logm(B))

tmp = logm(numpy.dot(numpy.dot(B.T, Ci), B))
            J += sample_weight[index] * tmp

        update = numpy.diag(numpy.diag(expm(J)))
        B = numpy.dot(B, invsqrtm(update))

        crit = distance_riemann(numpy.eye(Ne), update)

    A = numpy.linalg.inv(B)

    J = numpy.zeros((Ne, Ne))
    for index, Ci in enumerate(covmats):
        tmp = logm(numpy.dot(numpy.dot(B.T, Ci), B))
        J += sample_weight[index] * tmp

    C = numpy.dot(numpy.dot(A.T, expm(J)), A)
    return C

tau = numpy.finfo(numpy.float64).max
    crit = numpy.finfo(numpy.float64).max
    # stop when J<10^-9 or max iteration = 50
    while (crit > tol) and (k < maxiter) and (nu > tol):
        k = k + 1
        C12 = sqrtm(C)
        Cm12 = invsqrtm(C)
        J = numpy.zeros((Ne, Ne))

        for index in range(Nt):
            tmp = numpy.dot(numpy.dot(Cm12, covmats[index, :, :]), Cm12)
            J += sample_weight[index] * logm(tmp)

        crit = numpy.linalg.norm(J, ord='fro')
        h = nu * crit
        C = numpy.dot(numpy.dot(C12, expm(nu * J)), C12)
        if h < tau:
            nu = 0.95 * nu
            tau = h
        else:
            nu = 0.5 * nu

    return C

sample_weight = _get_sample_weight(sample_weight, covmats)
    Nt, Ne, Ne = covmats.shape
    crit = numpy.inf
    k = 0

    # init with AJD
    B, _ = ajd_pham(covmats)
    while (crit > tol) and (k < maxiter):
        k += 1
        J = numpy.zeros((Ne, Ne))

        for index, Ci in enumerate(covmats):
            tmp = logm(numpy.dot(numpy.dot(B.T, Ci), B))
            J += sample_weight[index] * tmp

        update = numpy.diag(numpy.diag(expm(J)))
        B = numpy.dot(B, invsqrtm(update))

        crit = distance_riemann(numpy.eye(Ne), update)

    A = numpy.linalg.inv(B)

    J = numpy.zeros((Ne, Ne))
    for index, Ci in enumerate(covmats):
        tmp = logm(numpy.dot(numpy.dot(B.T, Ci), B))
        J += sample_weight[index] * tmp

    C = numpy.dot(numpy.dot(A.T, expm(J)), A)
    return C

.. math::
            \mathbf{C} = \exp{(\\frac{1}{N} \sum_i \log{\mathbf{C}_i})}

    :param covmats: Covariance matrices set, Ntrials X Nchannels X Nchannels
    :param sample_weight: the weight of each sample

    :returns: the mean covariance matrix

    """
    sample_weight = _get_sample_weight(sample_weight, covmats)
    Nt, Ne, Ne = covmats.shape
    T = numpy.zeros((Ne, Ne))
    for index in range(Nt):
        T += sample_weight[index] * logm(covmats[index, :, :])
    C = expm(T)

    return C

The reference covariance matrix
    :returns: np.ndarray
        A set of Covariance matrix, Ntrials X Nchannels X Nchannels

    """
    Nt, Nd = T.shape
    Ne = int((numpy.sqrt(1 + 8 * Nd) - 1) / 2)
    C12 = sqrtm(Cref)

    idx = numpy.triu_indices_from(Cref)
    covmats = numpy.empty((Nt, Ne, Ne))
    covmats[:, idx[0], idx[1]] = T
    for i in range(Nt):
        triuc = numpy.triu(covmats[i], 1) / numpy.sqrt(2)
        covmats[i] = (numpy.diag(numpy.diag(covmats[i])) + triuc + triuc.T)
        covmats[i] = expm(covmats[i])
        covmats[i] = numpy.dot(numpy.dot(C12, covmats[i]), C12)

    return covmats