def eigen(X, P, NSIG=None, method='<STR_LIT>', threshold=None, NFFT=default_NFFT,<EOL>criteria='<STR_LIT>', verbose=False):
if method not in ['<STR_LIT>', '<STR_LIT>']:<EOL><INDENT>raise ValueError("<STR_LIT>")<EOL><DEDENT>if NSIG != None and threshold != None:<EOL><INDENT>raise ValueError("<STR_LIT>")<EOL><DEDENT>if NSIG is not None:<EOL><INDENT>if NSIG < <NUM_LIT:0>:<EOL><INDENT>raise ValueError('<STR_LIT>')<EOL><DEDENT>if NSIG >= P:<EOL><INDENT>raise ValueError("<STR_LIT>")<EOL><DEDENT><DEDENT>N = len(X)<EOL>NP = N - P<EOL>assert <NUM_LIT:2> * NP > P-<NUM_LIT:1>, '<STR_LIT>'<EOL>if NP > <NUM_LIT:100>:<EOL><INDENT>NP = <NUM_LIT:100><EOL><DEDENT>FB = np.zeros((<NUM_LIT:2>*NP, P), dtype=complex)<EOL>Z = np.zeros(NFFT, dtype=complex)<EOL>PSD = np.zeros(NFFT)<EOL>for I in range(<NUM_LIT:0>, NP):<EOL><INDENT>for K in range(<NUM_LIT:0>, P):<EOL><INDENT>FB[I, K] = X[I-K+P-<NUM_LIT:1>]<EOL>FB[I+NP, K] = X[I+K+<NUM_LIT:1>].conjugate()<EOL><DEDENT><DEDENT>_U, S, V = svd (FB)<EOL>V = -V.transpose()<EOL>NSIG  = _get_signal_space(S, <NUM_LIT:2>*NP,<EOL>verbose=verbose, threshold=threshold,<EOL>NSIG=NSIG, criteria=criteria)<EOL>for I in range(NSIG, P):<EOL><INDENT>Z[<NUM_LIT:0>:P] = V[<NUM_LIT:0>:P, I]<EOL>Z[P:NFFT] = <NUM_LIT:0><EOL>Z  = fft(Z, NFFT)<EOL>if method == '<STR_LIT>':<EOL><INDENT>PSD = PSD + abs(Z)**<NUM_LIT><EOL><DEDENT>elif method == '<STR_LIT>' :<EOL><INDENT>PSD = PSD + abs(Z)**<NUM_LIT> / S[I]<EOL><DEDENT><DEDENT>PSD = <NUM_LIT:1.>/PSD<EOL>nby2 = int(NFFT/<NUM_LIT:2>)<EOL>newpsd = np.append(PSD[nby2:<NUM_LIT:0>:-<NUM_LIT:1>], PSD[nby2*<NUM_LIT:2>-<NUM_LIT:1>:nby2-<NUM_LIT:1>:-<NUM_LIT:1>])<EOL>return newpsd, S<EOL>
r"""Pseudo spectrum using eigenvector method (EV or Music)

    This function computes either the Music or EigenValue (EV) noise
    subspace frequency estimator.

    First, an autocorrelation matrix of order `P` is computed from
    the data. Second, this matrix is separated into vector subspaces,
    one a signal subspace and the other a noise
    subspace using a SVD method to obtain the eigen values and vectors.
    From the eigen values :math:`\lambda_i`, and eigen vectors :math:`v_k`,
    the **pseudo spectrum** (see note below) is computed as follows:

    .. math:: P_{ev}(f) = \frac{1}{e^H(f)\left(\sum\limits_{k=M+1}^{p} \frac{1}{\lambda_k}v_kv_k^H\right)e(f)}

    The separation of the noise and signal subspaces requires expertise
    of the signal. However, AIC and MDL criteria may be used to automatically
    perform this task.

    You still need to provide the parameter `P` to indicate the maximum number
    of eigen values to be computed. The criteria will just select a subset
    to estimate the pseudo spectrum (see :func:`~spectrum.criteria.aic_eigen`
    and :func:`~spectrum.criteria.mdl_eigen` for details.

    .. note:: **pseudo spectrum**. func:`eigen` does not compute a PSD estimate.
        Indeed, the method does not preserve the measured process power.

    :param X: Array data samples
    :param int P: maximum number of eigen values to compute. NSIG (if
        specified) must therefore be less than P.
    :param str method: 'music' or 'ev'.
    :param int NSIG: If specified, the signal sub space uses NSIG eigen values.
    :param float threshold: If specified, the signal sub space is made of the
        eigen values larger than :math:`\rm{threshold} \times \lambda_{min}`,
        where :math:`\lambda_{min}` is the minimum eigen values.
    :param int NFFT: total length of the final data sets (padded with zero 
        if needed; default is 4096)

    :return:
        * PSD: Array of real frequency estimator values (two sided for
                complex data and one sided for real data)
        * S, the eigen values

    .. plot::
        :width: 80%
        :include-source:

        from spectrum import eigen, marple_data
        from pylab import plot, log10, linspace, legend, axis

        psd, ev = eigen(marple_data, 15, NSIG=11)
        f = linspace(-0.5, 0.5, len(psd))
        plot(f, 10 * log10(psd/max(psd)), label='User defined')

        psd, ev = eigen(marple_data, 15, threshold=2)
        plot(f, 10 * log10(psd/max(psd)), label='threshold method (100)')

        psd, ev = eigen(marple_data, 15)
        plot(f, 10 * log10(psd/max(psd)), label='AIC method (8)')

        legend()
        axis([-0.5, 0.5, -120, 0])

    .. seealso::
        :func:`pev`,
        :func:`pmusic`,
        :func:`~spectrum.criteria.aic_eigen`

    :References: [Marple]_, Chap 13

    .. todo:: for developers:

        * what should be the second argument of the criteria N, N-P, P...?
        * what should be the max value of NP