def test_legendre_val():
    """Test Legendre polynomial (derivative) equivalence
    """
    rng = np.random.RandomState(0)
    # check table equiv
    xs = np.linspace(-1., 1., 1000)
    n_terms = 100

    # True, numpy
    vals_np = legendre.legvander(xs, n_terms - 1)

    # Table approximation
    for fun, nc in zip([_get_legen_lut_fast, _get_legen_lut_accurate],
                       [100, 50]):
        lut, n_fact = _get_legen_table('eeg', n_coeff=nc, force_calc=True)
        vals_i = fun(xs, lut)
        # Need a "1:" here because we omit the first coefficient in our table!
        assert_allclose(vals_np[:, 1:vals_i.shape[1] + 1], vals_i,
                        rtol=1e-2, atol=5e-3)

        # Now let's look at our sums
        ctheta = rng.rand(20, 30) * 2.0 - 1.0
        beta = rng.rand(20, 30) * 0.8
        lut_fun = partial(fun, lut=lut)
        c1 = _comp_sum_eeg(beta.flatten(), ctheta.flatten(), lut_fun, n_fact)
        c1.shape = beta.shape

        # compare to numpy
        n = np.arange(1, n_terms, dtype=float)[:, np.newaxis, np.newaxis]
        coeffs = np.zeros((n_terms,) + beta.shape)
        coeffs[1:] = (np.cumprod([beta] * (n_terms - 1), axis=0) *
                      (2.0 * n + 1.0) * (2.0 * n + 1.0) / n)
        # can't use tensor=False here b/c it isn't in old numpy
        c2 = np.empty((20, 30))
        for ci1 in range(20):
            for ci2 in range(30):
                c2[ci1, ci2] = legendre.legval(ctheta[ci1, ci2],
                                               coeffs[:, ci1, ci2])
        assert_allclose(c1, c2, 1e-2, 1e-3)  # close enough...

    # compare fast and slow for MEG
    ctheta = rng.rand(20 * 30) * 2.0 - 1.0
    beta = rng.rand(20 * 30) * 0.8
    lut, n_fact = _get_legen_table('meg', n_coeff=10, force_calc=True)
    fun = partial(_get_legen_lut_fast, lut=lut)
    coeffs = _comp_sums_meg(beta, ctheta, fun, n_fact, False)
    lut, n_fact = _get_legen_table('meg', n_coeff=20, force_calc=True)
    fun = partial(_get_legen_lut_accurate, lut=lut)
    coeffs = _comp_sums_meg(beta, ctheta, fun, n_fact, False)

def test_cache_dir():
    """Test use of cache dir
    """
    tempdir = _TempDir()
    orig_dir = os.getenv("MNE_CACHE_DIR", None)
    orig_size = os.getenv("MNE_MEMMAP_MIN_SIZE", None)
    rng = np.random.RandomState(0)
    X = rng.randn(9, 2, 10)
    log_file = op.join(tempdir, "log.txt")
    try:
        os.environ["MNE_MEMMAP_MIN_SIZE"] = "1K"
        os.environ["MNE_CACHE_DIR"] = tempdir
        # Fix error for #1507: in-place when memmapping
        permutation_cluster_1samp_test(
            X, buffer_size=None, n_jobs=2, n_permutations=1, seed=0, stat_fun=ttest_1samp_no_p, verbose=False
        )
        # ensure that non-independence yields warning
        stat_fun = partial(ttest_1samp_no_p, sigma=1e-3)
        set_log_file(log_file)
        permutation_cluster_1samp_test(
            X, buffer_size=10, n_jobs=2, n_permutations=1, seed=0, stat_fun=stat_fun, verbose=False
        )
        with open(log_file, "r") as fid:
            assert_true("independently" in "".join(fid.readlines()))
    finally:
        if orig_dir is not None:
            os.environ["MNE_CACHE_DIR"] = orig_dir
        else:
            del os.environ["MNE_CACHE_DIR"]
        if orig_size is not None:
            os.environ["MNE_MEMMAP_MIN_SIZE"] = orig_size
        else:
            del os.environ["MNE_MEMMAP_MIN_SIZE"]
        set_log_file(None)

def test_read_ch_connectivity():
    "Test reading channel connectivity templates"
    a = partial(np.array, dtype='<U7')
    # no pep8
    nbh = np.array([[(['MEG0111'], [[a(['MEG0131'])]]),
                     (['MEG0121'], [[a(['MEG0111'])],
                                    [a(['MEG0131'])]]),
                     (['MEG0131'], [[a(['MEG0111'])],
                                    [a(['MEG0121'])]])]],
                   dtype=[('label', 'O'), ('neighblabel', 'O')])
    mat = dict(neighbours=nbh)
    mat_fname = op.join(tempdir, 'test_mat.mat')
    savemat(mat_fname, mat)

    ch_connectivity = read_ch_connectivity(mat_fname)
    x = ch_connectivity
    assert_equal(x.shape, (3, 3))
    assert_equal(x[0, 1], False)
    assert_equal(x[0, 2], True)
    assert_true(np.all(x.diagonal()))
    assert_raises(ValueError, read_ch_connectivity, mat_fname, [0, 3])
    ch_connectivity = read_ch_connectivity(mat_fname, picks=[0, 2])
    assert_equal(ch_connectivity.shape[0], 2)

    ch_names = ['EEG01', 'EEG02', 'EEG03']
    neighbors = [['EEG02'], ['EEG04'], ['EEG02']]
    assert_raises(ValueError, ch_neighbor_connectivity, ch_names, neighbors)
    neighbors = [['EEG02'], ['EEG01', 'EEG03'], ['EEG 02']]
    assert_raises(ValueError, ch_neighbor_connectivity, ch_names[:2],
                  neighbors)
    neighbors = [['EEG02'], 'EEG01', ['EEG 02']]
    assert_raises(ValueError, ch_neighbor_connectivity, ch_names, neighbors)

def test_cache_dir():
    """Test use of cache dir
    """
    tempdir = _TempDir()
    orig_dir = os.getenv('MNE_CACHE_DIR', None)
    orig_size = os.getenv('MNE_MEMMAP_MIN_SIZE', None)
    rng = np.random.RandomState(0)
    X = rng.randn(9, 2, 10)
    try:
        os.environ['MNE_MEMMAP_MIN_SIZE'] = '1K'
        os.environ['MNE_CACHE_DIR'] = tempdir
        # Fix error for #1507: in-place when memmapping
        with catch_logging() as log_file:
            permutation_cluster_1samp_test(
                X, buffer_size=None, n_jobs=2, n_permutations=1,
                seed=0, stat_fun=ttest_1samp_no_p, verbose=False)
            # ensure that non-independence yields warning
            stat_fun = partial(ttest_1samp_no_p, sigma=1e-3)
            assert_true('independently' not in log_file.getvalue())
            with warnings.catch_warnings(record=True):  # independently
                permutation_cluster_1samp_test(
                    X, buffer_size=10, n_jobs=2, n_permutations=1,
                    seed=0, stat_fun=stat_fun, verbose=False)
            assert_true('independently' in log_file.getvalue())
    finally:
        if orig_dir is not None:
            os.environ['MNE_CACHE_DIR'] = orig_dir
        else:
            del os.environ['MNE_CACHE_DIR']
        if orig_size is not None:
            os.environ['MNE_MEMMAP_MIN_SIZE'] = orig_size
        else:
            del os.environ['MNE_MEMMAP_MIN_SIZE']

def test_no_conversion():
    """ Test bti no-conversion option """

    get_info = partial(
        _get_bti_info,
        rotation_x=0.0, translation=(0.0, 0.02, 0.11), convert=False,
        ecg_ch='E31', eog_ch=('E63', 'E64'),
        rename_channels=False, sort_by_ch_name=False)

    for pdf, config, hs in zip(pdf_fnames, config_fnames, hs_fnames):
        with warnings.catch_warnings(record=True):  # weight tables
            raw_info, _ = get_info(pdf, config, hs, convert=False)
        with warnings.catch_warnings(record=True):  # weight tables
            raw_info_con = read_raw_bti(
                pdf_fname=pdf, config_fname=config, head_shape_fname=hs,
                convert=True, preload=False).info

        pick_info(raw_info_con,
                  pick_types(raw_info_con, meg=True, ref_meg=True),
                  copy=False)
        pick_info(raw_info,
                  pick_types(raw_info, meg=True, ref_meg=True), copy=False)
        bti_info = _read_bti_header(pdf, config)
        dev_ctf_t = _correct_trans(bti_info['bti_transform'][0])
        assert_array_equal(dev_ctf_t, raw_info['dev_ctf_t']['trans'])
        assert_array_equal(raw_info['dev_head_t']['trans'], np.eye(4))
        assert_array_equal(raw_info['ctf_head_t']['trans'], np.eye(4))
        dig, t = _process_bti_headshape(hs, convert=False, use_hpi=False)
        assert_array_equal(t['trans'], np.eye(4))

        for ii, (old, new, con) in enumerate(zip(
                dig, raw_info['dig'], raw_info_con['dig'])):
            assert_equal(old['ident'], new['ident'])
            assert_array_equal(old['r'], new['r'])
            assert_true(not np.allclose(old['r'], con['r']))

            if ii > 10:
                break

        ch_map = dict((ch['chan_label'],
                       ch['loc']) for ch in bti_info['chs'])

        for ii, ch_label in enumerate(raw_info['ch_names']):
            if not ch_label.startswith('A'):
                continue
            t1 = ch_map[ch_label]  # correction already performed in bti_info
            t2 = raw_info['chs'][ii]['loc']
            t3 = raw_info_con['chs'][ii]['loc']
            assert_allclose(t1, t2, atol=1e-15)
            assert_true(not np.allclose(t1, t3))
            idx_a = raw_info_con['ch_names'].index('MEG 001')
            idx_b = raw_info['ch_names'].index('A22')
            assert_equal(
                raw_info_con['chs'][idx_a]['coord_frame'],
                FIFF.FIFFV_COORD_DEVICE)
            assert_equal(
                raw_info['chs'][idx_b]['coord_frame'],
                FIFF.FIFFV_MNE_COORD_4D_HEAD)

def test_cluster_permutation_t_test():
    """Test cluster level permutations T-test
    """
    condition1_1d, condition2_1d, condition1_2d, condition2_2d = \
        _get_conditions()

    # use a very large sigma to make sure Ts are not independent
    stat_funs = [ttest_1samp_no_p,
                 partial(ttest_1samp_no_p, sigma=1e-1)]

    for stat_fun in stat_funs:
        for condition1 in (condition1_1d, condition1_2d):
            # these are so significant we can get away with fewer perms
            T_obs, clusters, cluster_p_values, hist =\
                permutation_cluster_1samp_test(condition1, n_permutations=100,
                                               tail=0, seed=1,
                                               buffer_size=None)
            assert_equal(np.sum(cluster_p_values < 0.05), 1)

            T_obs_pos, c_1, cluster_p_values_pos, _ =\
                permutation_cluster_1samp_test(condition1, n_permutations=100,
                                               tail=1, threshold=1.67, seed=1,
                                               stat_fun=stat_fun,
                                               buffer_size=None)

            T_obs_neg, _, cluster_p_values_neg, _ =\
                permutation_cluster_1samp_test(-condition1, n_permutations=100,
                                               tail=-1, threshold=-1.67,
                                               seed=1, stat_fun=stat_fun,
                                               buffer_size=None)
            assert_array_equal(T_obs_pos, -T_obs_neg)
            assert_array_equal(cluster_p_values_pos < 0.05,
                               cluster_p_values_neg < 0.05)

            # test with 2 jobs and buffer_size enabled
            buffer_size = condition1.shape[1] // 10
            T_obs_neg_buff, _, cluster_p_values_neg_buff, _ = \
                permutation_cluster_1samp_test(-condition1, n_permutations=100,
                                               tail=-1, threshold=-1.67,
                                               seed=1, n_jobs=2,
                                               stat_fun=stat_fun,
                                               buffer_size=buffer_size)

            assert_array_equal(T_obs_neg, T_obs_neg_buff)
            assert_array_equal(cluster_p_values_neg, cluster_p_values_neg_buff)

def test_permutation_connectivity_equiv():
    """Test cluster level permutations with and without connectivity
    """
    try:
        try:
            from sklearn.feature_extraction.image import grid_to_graph
        except ImportError:
            from scikits.learn.feature_extraction.image import grid_to_graph
    except ImportError:
        return
    rng = np.random.RandomState(0)
    # subjects, time points, spatial points
    X = rng.randn(7, 2, 10)
    # add some significant points
    X[:, 0:2, 0:2] += 10  # span two time points and two spatial points
    X[:, 1, 5:9] += 10  # span four time points
    max_steps = [1, 1, 1, 2]
    # This will run full algorithm in two ways, then the ST-algorithm in 2 ways
    # All of these should give the same results
    conns = [None, grid_to_graph(2, 10),
             grid_to_graph(1, 10), grid_to_graph(1, 10)]
    stat_map = None
    thresholds = [2, dict(start=0.5, step=0.5)]
    sig_counts = [2, 8]
    sdps = [0, 0.05, 0.05]
    ots = ['mask', 'mask', 'indices']
    for thresh, count in zip(thresholds, sig_counts):
        cs = None
        ps = None
        for max_step, conn in zip(max_steps, conns):
            for stat_fun in [ttest_1samp_no_p,
                             partial(ttest_1samp_no_p, sigma=1e-3)]:
                for sdp, ot in zip(sdps, ots):
                    t, clusters, p, H0 = \
                            permutation_cluster_1samp_test(X,
                                                           threshold=thresh,
                                                           connectivity=conn,
                                                           n_jobs=2,
                                                           max_step=max_step,
                                                           stat_fun=stat_fun,
                                                           step_down_p=sdp,
                                                           out_type=ot)
                    # make sure our output datatype is correct
                    if ot == 'mask':
                        assert_true(isinstance(clusters[0], np.ndarray))
                        assert_true(clusters[0].dtype == bool)
                        assert_array_equal(clusters[0].shape, X.shape[1:])
                    else:  # ot == 'indices'
                        assert_true(isinstance(clusters[0], tuple))

                    # make sure all comparisons were done; for TFCE, no perm
                    # should come up empty
                    if count == 8:
                        assert_true(not np.any(H0 == 0))
                    inds = np.where(p < 0.05)[0]
                    assert_true(len(inds) == count)
                    this_cs = [clusters[ii] for ii in inds]
                    this_ps = p[inds]
                    this_stat_map = np.zeros((2, 10), dtype=bool)
                    for ci, c in enumerate(this_cs):
                        if isinstance(c, tuple):
                            this_c = np.zeros((2, 10), bool)
                            for x, y in zip(c[0], c[1]):
                                this_stat_map[x, y] = True
                                this_c[x, y] = True
                            this_cs[ci] = this_c
                            c = this_c
                        this_stat_map[c] = True
                    if cs is None:
                        ps = this_ps
                        cs = this_cs
                    if stat_map is None:
                        stat_map = this_stat_map
                    assert_array_equal(ps, this_ps)
                    assert_true(len(cs) == len(this_cs))
                    for c1, c2 in zip(cs, this_cs):
                        assert_array_equal(c1, c2)
                    assert_array_equal(stat_map, this_stat_map)

                 make_ad_hoc_cov)
from mne.io import Raw
from mne.utils import _TempDir, slow_test, requires_module
from mne.io.proc_history import _get_sss_rank
from mne.io.pick import channel_type, _picks_by_type
from mne.fixes import partial

_recent_sklearn_call = """
required_version = '0.15'
import sklearn
version = LooseVersion(sklearn.__version__)
if version < required_version:
    raise ImportError
"""

requires_sklearn_0_15 = partial(requires_module, name='sklearn',
                                call=_recent_sklearn_call)
warnings.simplefilter('always')  # enable b/c these tests throw warnings

base_dir = op.join(op.dirname(__file__), '..', 'io', 'tests', 'data')
cov_fname = op.join(base_dir, 'test-cov.fif')
cov_gz_fname = op.join(base_dir, 'test-cov.fif.gz')
cov_km_fname = op.join(base_dir, 'test-km-cov.fif')
raw_fname = op.join(base_dir, 'test_raw.fif')
ave_fname = op.join(base_dir, 'test-ave.fif')
erm_cov_fname = op.join(base_dir, 'test_erm-cov.fif')
hp_fif_fname = op.join(base_dir, 'test_chpi_raw_sss.fif')


def test_ad_hoc_cov():
    """Test ad hoc cov creation and I/O"""
    tempdir = _TempDir()