Python maxwell.maxwell_filter函数代码示例

def test_head_translation():
    """Test Maxwell filter head translation"""
    with warnings.catch_warnings(record=True):  # maxshield
        raw = Raw(raw_fname, allow_maxshield=True).crop(0.0, 1.0, False)
    # First try with an unchanged destination
    raw_sss = maxwell_filter(raw, destination=raw_fname, origin=mf_head_origin, regularize=None, bad_condition="ignore")
    assert_meg_snr(raw_sss, Raw(sss_std_fname).crop(0.0, 1.0, False), 200.0)
    # Now with default
    with catch_logging() as log:
        raw_sss = maxwell_filter(
            raw,
            destination=mf_head_origin,
            origin=mf_head_origin,
            regularize=None,
            bad_condition="ignore",
            verbose="warning",
        )
    assert_true("over 25 mm" in log.getvalue())
    assert_meg_snr(raw_sss, Raw(sss_trans_default_fname), 125.0)
    # Now to sample's head pos
    with catch_logging() as log:
        raw_sss = maxwell_filter(
            raw,
            destination=sample_fname,
            origin=mf_head_origin,
            regularize=None,
            bad_condition="ignore",
            verbose="warning",
        )
    assert_true("= 25.6 mm" in log.getvalue())
    assert_meg_snr(raw_sss, Raw(sss_trans_sample_fname), 350.0)
    # Degenerate cases
    assert_raises(RuntimeError, maxwell_filter, raw, destination=mf_head_origin, coord_frame="meg")
    assert_raises(ValueError, maxwell_filter, raw, destination=[0.0] * 4)

def test_fine_calibration():
    """Test Maxwell filter fine calibration"""

    # Load testing data (raw, SSS std origin, SSS non-standard origin)
    with warnings.catch_warnings(record=True):  # maxshield
        raw = Raw(raw_fname, allow_maxshield=True).crop(0., 1., False)
    sss_fine_cal = Raw(sss_fine_cal_fname)

    # Test 1D SSS fine calibration
    raw_sss = maxwell_filter(raw, calibration=fine_cal_fname,
                             origin=mf_head_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, sss_fine_cal, 82, 611)
    py_cal = raw_sss.info['proc_history'][0]['max_info']['sss_cal']
    assert_true(py_cal is not None)
    assert_true(len(py_cal) > 0)
    mf_cal = sss_fine_cal.info['proc_history'][0]['max_info']['sss_cal']
    # we identify these differently
    mf_cal['cal_chans'][mf_cal['cal_chans'][:, 1] == 3022, 1] = 3024
    assert_allclose(py_cal['cal_chans'], mf_cal['cal_chans'])
    assert_allclose(py_cal['cal_corrs'], mf_cal['cal_corrs'],
                    rtol=1e-3, atol=1e-3)

    # Test 3D SSS fine calibration (no equivalent func in MaxFilter yet!)
    # very low SNR as proc differs, eventually we should add a better test
    raw_sss_3D = maxwell_filter(raw, calibration=fine_cal_fname_3d,
                                origin=mf_head_origin, regularize=None,
                                bad_condition='ignore')
    assert_meg_snr(raw_sss_3D, sss_fine_cal, 1.0, 6.)

def test_cross_talk():
    """Test Maxwell filter cross-talk cancellation"""
    raw = Raw(raw_fname, allow_maxshield='yes').crop(0., 1., False)
    raw.info['bads'] = bads
    sss_ctc = Raw(sss_ctc_fname)
    raw_sss = maxwell_filter(raw, cross_talk=ctc_fname,
                             origin=mf_head_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, sss_ctc, 275.)
    py_ctc = raw_sss.info['proc_history'][0]['max_info']['sss_ctc']
    assert_true(len(py_ctc) > 0)
    assert_raises(ValueError, maxwell_filter, raw, cross_talk=raw)
    assert_raises(ValueError, maxwell_filter, raw, cross_talk=raw_fname)
    mf_ctc = sss_ctc.info['proc_history'][0]['max_info']['sss_ctc']
    del mf_ctc['block_id']  # we don't write this
    assert_equal(object_diff(py_ctc, mf_ctc), '')
    raw_ctf = Raw(fname_ctf_raw)
    assert_raises(ValueError, maxwell_filter, raw_ctf)  # cannot fit headshape
    raw_sss = maxwell_filter(raw_ctf, origin=(0., 0., 0.04))
    _assert_n_free(raw_sss, 68)
    raw_sss = maxwell_filter(raw_ctf, origin=(0., 0., 0.04), ignore_ref=True)
    _assert_n_free(raw_sss, 70)
    raw_missing = raw.crop(0, 0.1, copy=True).load_data().pick_channels(
        [raw.ch_names[pi] for pi in pick_types(raw.info, meg=True,
                                               exclude=())[3:]])
    with warnings.catch_warnings(record=True) as w:
        maxwell_filter(raw_missing, cross_talk=ctc_fname)
    assert_equal(len(w), 1)
    assert_true('Not all cross-talk channels in raw' in str(w[0].message))
    # MEG channels not in cross-talk
    assert_raises(RuntimeError, maxwell_filter, raw_ctf, origin=(0., 0., 0.04),
                  cross_talk=ctc_fname)

def test_fine_calibration():
    """Test Maxwell filter fine calibration."""

    # Load testing data (raw, SSS std origin, SSS non-standard origin)
    raw = read_crop(raw_fname, (0., 1.))
    sss_fine_cal = read_crop(sss_fine_cal_fname)

    # Test 1D SSS fine calibration
    raw_sss = maxwell_filter(raw, calibration=fine_cal_fname,
                             origin=mf_head_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, sss_fine_cal, 82, 611)
    py_cal = raw_sss.info['proc_history'][0]['max_info']['sss_cal']
    assert_true(py_cal is not None)
    assert_true(len(py_cal) > 0)
    mf_cal = sss_fine_cal.info['proc_history'][0]['max_info']['sss_cal']
    # we identify these differently
    mf_cal['cal_chans'][mf_cal['cal_chans'][:, 1] == 3022, 1] = 3024
    assert_allclose(py_cal['cal_chans'], mf_cal['cal_chans'])
    assert_allclose(py_cal['cal_corrs'], mf_cal['cal_corrs'],
                    rtol=1e-3, atol=1e-3)

    # Test 3D SSS fine calibration (no equivalent func in MaxFilter yet!)
    # very low SNR as proc differs, eventually we should add a better test
    raw_sss_3D = maxwell_filter(raw, calibration=fine_cal_fname_3d,
                                origin=mf_head_origin, regularize=None,
                                bad_condition='ignore')
    assert_meg_snr(raw_sss_3D, sss_fine_cal, 1.0, 6.)
    raw_ctf = read_crop(fname_ctf_raw).apply_gradient_compensation(0)
    assert_raises(RuntimeError, maxwell_filter, raw_ctf, origin=(0., 0., 0.04),
                  calibration=fine_cal_fname)

def test_maxwell_noise_rejection():
    """Test Maxwell filter shielding factor using empty room"""
    with warnings.catch_warnings(record=True):  # maxshield
        raw_erm = Raw(erm_fname, allow_maxshield=True, preload=True)
    picks = pick_types(raw_erm.info, meg=True)
    erm_power = raw_erm[picks][0].ravel()
    erm_power = np.sqrt(np.sum(erm_power * erm_power))
    # Vanilla SSS
    _assert_shielding(Raw(sss_erm_std_fname), erm_power, 1.5)
    raw_sss = maxwell_filter(raw_erm, coord_frame='meg')
    _assert_shielding(raw_sss, erm_power, 1.5)
    # tSSS
    _assert_shielding(Raw(sss_erm_st_fname), erm_power, 5)
    raw_sss = maxwell_filter(raw_erm, st_duration=1., coord_frame='meg')
    _assert_shielding(raw_sss, erm_power, 5.)
    # Fine cal
    _assert_shielding(Raw(sss_erm_fine_cal_fname), erm_power, 2)
    raw_sss = maxwell_filter(raw_erm, calibration=fine_cal_fname,
                             coord_frame='meg')
    _assert_shielding(raw_sss, erm_power, 2.)
    # Crosstalk
    _assert_shielding(Raw(sss_erm_ctc_fname), erm_power, 2.1)
    raw_sss = maxwell_filter(raw_erm, cross_talk=ctc_fname, coord_frame='meg')
    _assert_shielding(raw_sss, erm_power, 2.1)
    # tSSS + fine cal + ctc
    _assert_shielding(Raw(sss_erm_st1FineCalCrossTalk_fname), erm_power, 6.)
    raw_sss = maxwell_filter(raw_erm, calibration=fine_cal_fname,
                             cross_talk=ctc_fname, st_duration=1.,
                             coord_frame='meg')
    _assert_shielding(raw_sss, erm_power, 100)  # somehow this is really high?

def test_head_translation():
    """Test Maxwell filter head translation."""
    raw = read_crop(raw_fname, (0., 1.))
    # First try with an unchanged destination
    with use_coil_def(elekta_def_fname):
        raw_sss = maxwell_filter(raw, destination=raw_fname,
                                 origin=mf_head_origin, regularize=None,
                                 bad_condition='ignore')
    assert_meg_snr(raw_sss, read_crop(sss_std_fname, (0., 1.)), 200.)
    # Now with default
    with use_coil_def(elekta_def_fname):
        with pytest.warns(RuntimeWarning, match='over 25 mm'):
            raw_sss = maxwell_filter(raw, destination=mf_head_origin,
                                     origin=mf_head_origin, regularize=None,
                                     bad_condition='ignore', verbose=True)
    assert_meg_snr(raw_sss, read_crop(sss_trans_default_fname), 125.)
    destination = np.eye(4)
    destination[2, 3] = 0.04
    assert_allclose(raw_sss.info['dev_head_t']['trans'], destination)
    # Now to sample's head pos
    with pytest.warns(RuntimeWarning, match='= 25.6 mm'):
        raw_sss = maxwell_filter(raw, destination=sample_fname,
                                 origin=mf_head_origin, regularize=None,
                                 bad_condition='ignore', verbose=True)
    assert_meg_snr(raw_sss, read_crop(sss_trans_sample_fname), 13., 100.)
    assert_allclose(raw_sss.info['dev_head_t']['trans'],
                    read_info(sample_fname)['dev_head_t']['trans'])
    # Degenerate cases
    pytest.raises(RuntimeError, maxwell_filter, raw,
                  destination=mf_head_origin, coord_frame='meg')
    pytest.raises(ValueError, maxwell_filter, raw, destination=[0.] * 4)

def test_head_translation():
    """Test Maxwell filter head translation"""
    with warnings.catch_warnings(record=True):  # maxshield
        raw = Raw(raw_fname, allow_maxshield=True).crop(0., 1., False)
    # First try with an unchanged destination
    raw_sss = maxwell_filter(raw, destination=raw_fname,
                             origin=mf_head_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, Raw(sss_std_fname).crop(0., 1., False), 200.)
    # Now with default
    with warnings.catch_warnings(record=True):
        with catch_logging() as log:
            raw_sss = maxwell_filter(raw, destination=mf_head_origin,
                                     origin=mf_head_origin, regularize=None,
                                     bad_condition='ignore', verbose='warning')
    assert_true('over 25 mm' in log.getvalue())
    assert_meg_snr(raw_sss, Raw(sss_trans_default_fname), 125.)
    destination = np.eye(4)
    destination[2, 3] = 0.04
    assert_allclose(raw_sss.info['dev_head_t']['trans'], destination)
    # Now to sample's head pos
    with warnings.catch_warnings(record=True):
        with catch_logging() as log:
            raw_sss = maxwell_filter(raw, destination=sample_fname,
                                     origin=mf_head_origin, regularize=None,
                                     bad_condition='ignore', verbose='warning')
    assert_true('= 25.6 mm' in log.getvalue())
    assert_meg_snr(raw_sss, Raw(sss_trans_sample_fname), 350.)
    assert_allclose(raw_sss.info['dev_head_t']['trans'],
                    read_info(sample_fname)['dev_head_t']['trans'])
    # Degenerate cases
    assert_raises(RuntimeError, maxwell_filter, raw,
                  destination=mf_head_origin, coord_frame='meg')
    assert_raises(ValueError, maxwell_filter, raw, destination=[0.] * 4)

def test_basic():
    """Test Maxwell filter basic version"""
    # Load testing data (raw, SSS std origin, SSS non-standard origin)
    with warnings.catch_warnings(record=True):  # maxshield
        raw = Raw(raw_fname, allow_maxshield=True).crop(0., 1., False)
        raw_err = Raw(raw_fname, proj=True, allow_maxshield=True)
        raw_erm = Raw(erm_fname, allow_maxshield=True)
    assert_raises(RuntimeError, maxwell_filter, raw_err)
    assert_raises(TypeError, maxwell_filter, 1.)  # not a raw
    assert_raises(ValueError, maxwell_filter, raw, int_order=20)  # too many

    n_int_bases = int_order ** 2 + 2 * int_order
    n_ext_bases = ext_order ** 2 + 2 * ext_order
    nbases = n_int_bases + n_ext_bases

    # Check number of bases computed correctly
    assert_equal(_get_n_moments([int_order, ext_order]).sum(), nbases)

    # Test SSS computation at the standard head origin
    raw_sss = maxwell_filter(raw, origin=mf_head_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, Raw(sss_std_fname), 200., 1000.)
    py_cal = raw_sss.info['proc_history'][0]['max_info']['sss_cal']
    assert_equal(len(py_cal), 0)
    py_ctc = raw_sss.info['proc_history'][0]['max_info']['sss_ctc']
    assert_equal(len(py_ctc), 0)
    py_st = raw_sss.info['proc_history'][0]['max_info']['max_st']
    assert_equal(len(py_st), 0)
    assert_raises(RuntimeError, maxwell_filter, raw_sss)

    # Test SSS computation at non-standard head origin
    raw_sss = maxwell_filter(raw, origin=[0., 0.02, 0.02], regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, Raw(sss_nonstd_fname), 250., 700.)

    # Test SSS computation at device origin
    sss_erm_std = Raw(sss_erm_std_fname)
    raw_sss = maxwell_filter(raw_erm, coord_frame='meg',
                             origin=mf_meg_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, sss_erm_std, 100., 900.)
    for key in ('job', 'frame'):
        vals = [x.info['proc_history'][0]['max_info']['sss_info'][key]
                for x in [raw_sss, sss_erm_std]]
        assert_equal(vals[0], vals[1])

    # Check against SSS functions from proc_history
    sss_info = raw_sss.info['proc_history'][0]['max_info']
    assert_equal(_get_n_moments(int_order),
                 proc_history._get_sss_rank(sss_info))

    # Degenerate cases
    raw_bad = raw.copy()
    raw_bad.comp = True
    assert_raises(RuntimeError, maxwell_filter, raw_bad)
    del raw_bad
    assert_raises(ValueError, maxwell_filter, raw, coord_frame='foo')
    assert_raises(ValueError, maxwell_filter, raw, origin='foo')
    assert_raises(ValueError, maxwell_filter, raw, origin=[0] * 4)

def test_basic():
    """Test Maxwell filter basic version."""
    # Load testing data (raw, SSS std origin, SSS non-standard origin)
    raw = read_crop(raw_fname, (0., 1.))
    raw_err = read_crop(raw_fname).apply_proj()
    raw_erm = read_crop(erm_fname)
    assert_raises(RuntimeError, maxwell_filter, raw_err)
    assert_raises(TypeError, maxwell_filter, 1.)  # not a raw
    assert_raises(ValueError, maxwell_filter, raw, int_order=20)  # too many

    n_int_bases = int_order ** 2 + 2 * int_order
    n_ext_bases = ext_order ** 2 + 2 * ext_order
    nbases = n_int_bases + n_ext_bases

    # Check number of bases computed correctly
    assert_equal(_get_n_moments([int_order, ext_order]).sum(), nbases)

    # Test SSS computation at the standard head origin
    assert_equal(len(raw.info['projs']), 12)  # 11 MEG projs + 1 AVG EEG
    raw_sss = maxwell_filter(raw, origin=mf_head_origin, regularize=None,
                             bad_condition='ignore')
    assert_equal(len(raw_sss.info['projs']), 1)  # avg EEG
    assert_equal(raw_sss.info['projs'][0]['desc'], 'Average EEG reference')
    assert_meg_snr(raw_sss, read_crop(sss_std_fname), 200., 1000.)
    py_cal = raw_sss.info['proc_history'][0]['max_info']['sss_cal']
    assert_equal(len(py_cal), 0)
    py_ctc = raw_sss.info['proc_history'][0]['max_info']['sss_ctc']
    assert_equal(len(py_ctc), 0)
    py_st = raw_sss.info['proc_history'][0]['max_info']['max_st']
    assert_equal(len(py_st), 0)
    assert_raises(RuntimeError, maxwell_filter, raw_sss)

    # Test SSS computation at non-standard head origin
    raw_sss = maxwell_filter(raw, origin=[0., 0.02, 0.02], regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, read_crop(sss_nonstd_fname), 250., 700.)

    # Test SSS computation at device origin
    sss_erm_std = read_crop(sss_erm_std_fname)
    raw_sss = maxwell_filter(raw_erm, coord_frame='meg',
                             origin=mf_meg_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, sss_erm_std, 100., 900.)
    for key in ('job', 'frame'):
        vals = [x.info['proc_history'][0]['max_info']['sss_info'][key]
                for x in [raw_sss, sss_erm_std]]
        assert_equal(vals[0], vals[1])

    # Check against SSS functions from proc_history
    sss_info = raw_sss.info['proc_history'][0]['max_info']
    assert_equal(_get_n_moments(int_order),
                 proc_history._get_sss_rank(sss_info))

    # Degenerate cases
    assert_raises(ValueError, maxwell_filter, raw, coord_frame='foo')
    assert_raises(ValueError, maxwell_filter, raw, origin='foo')
    assert_raises(ValueError, maxwell_filter, raw, origin=[0] * 4)
    assert_raises(ValueError, maxwell_filter, raw, mag_scale='foo')

def test_spatiotemporal_maxwell():
    """Test spatiotemporal (tSSS) processing"""
    # Load raw testing data
    with warnings.catch_warnings(record=True):  # maxshield
        raw = Raw(raw_fname, allow_maxshield=True)

    # Create coils
    picks = pick_types(raw.info)

    # Test that window is less than length of data
    assert_raises(ValueError, maxwell_filter, raw, st_dur=1000.)

    # Check both 4 and 10 seconds because Elekta handles them differently
    # This is to ensure that std/non-std tSSS windows are correctly handled
    st_durs = [4., 10.]
    for st_dur in st_durs:
        # Load tSSS data depending on st_dur and get data
        tSSS_fname = op.join(data_path, 'SSS', 'test_move_anon_raw_' +
                             'spatiotemporal_%0ds_sss.fif' % st_dur)

        with warnings.catch_warnings(record=True):  # maxshield, naming
            tsss_bench = Raw(tSSS_fname, allow_maxshield=True)
            # Because Elekta's tSSS sometimes(!) lumps the tail window of data
            # onto the previous buffer if it's shorter than st_dur, we have to
            # crop the data here to compensate for Elekta's tSSS behavior.
            if st_dur == 10.:
                tsss_bench.crop(0, st_dur, copy=False)
        tsss_bench_data = tsss_bench[picks, :][0]
        del tsss_bench

        # Test sss computation at the standard head origin. Same cropping issue
        # as mentioned above.
        if st_dur == 10.:
            raw_tsss = maxwell_filter(raw.crop(0, st_dur), st_dur=st_dur)
        else:
            raw_tsss = maxwell_filter(raw, st_dur=st_dur)
        assert_allclose(raw_tsss[picks][0], tsss_bench_data,
                        rtol=1e-12, atol=1e-4, err_msg='Spatiotemporal (tSSS) '
                        'maxwell filtered data at standard origin incorrect.')

        # Confirm SNR is above 500. Single precision is part of discrepancy
        bench_rms = np.sqrt(np.mean(tsss_bench_data * tsss_bench_data, axis=1))
        error = raw_tsss[picks][0] - tsss_bench_data
        error_rms = np.sqrt(np.mean(error * error, axis=1))
        assert_true(np.mean(bench_rms / error_rms) >= 500,
                    'SNR (%0.1f) < 500' % np.mean(bench_rms / error_rms))

    # Confirm we didn't modify other channels (like EEG chs)
    non_picks = np.setdiff1d(np.arange(len(raw.ch_names)), picks)
    assert_allclose(raw[non_picks, 0:raw_tsss.n_times][0],
                    raw_tsss[non_picks, 0:raw_tsss.n_times][0])

    # Degenerate cases
    assert_raises(ValueError, maxwell_filter, raw, st_dur=10., st_corr=0.)

def test_spatiotemporal_maxwell():
    """Test Maxwell filter (tSSS) spatiotemporal processing"""
    # Load raw testing data
    raw = Raw(raw_fname, allow_maxshield='yes')

    # Test that window is less than length of data
    assert_raises(ValueError, maxwell_filter, raw, st_duration=1000.)

    # Check both 4 and 10 seconds because Elekta handles them differently
    # This is to ensure that std/non-std tSSS windows are correctly handled
    st_durations = [4., 10.]
    tols = [325., 200.]
    for st_duration, tol in zip(st_durations, tols):
        # Load tSSS data depending on st_duration and get data
        tSSS_fname = op.join(sss_path,
                             'test_move_anon_st%0ds_raw_sss.fif' % st_duration)
        tsss_bench = Raw(tSSS_fname)
        # Because Elekta's tSSS sometimes(!) lumps the tail window of data
        # onto the previous buffer if it's shorter than st_duration, we have to
        # crop the data here to compensate for Elekta's tSSS behavior.
        if st_duration == 10.:
            tsss_bench.crop(0, st_duration, copy=False)

        # Test sss computation at the standard head origin. Same cropping issue
        # as mentioned above.
        if st_duration == 10.:
            raw_tsss = maxwell_filter(raw.crop(0, st_duration),
                                      origin=mf_head_origin,
                                      st_duration=st_duration, regularize=None,
                                      bad_condition='ignore')
        else:
            raw_tsss = maxwell_filter(raw, st_duration=st_duration,
                                      origin=mf_head_origin, regularize=None,
                                      bad_condition='ignore', verbose=True)
            raw_tsss_2 = maxwell_filter(raw, st_duration=st_duration,
                                        origin=mf_head_origin, regularize=None,
                                        bad_condition='ignore', st_fixed=False,
                                        verbose=True)
            assert_meg_snr(raw_tsss, raw_tsss_2, 100., 1000.)
            assert_equal(raw_tsss.estimate_rank(), 140)
            assert_equal(raw_tsss_2.estimate_rank(), 140)
        assert_meg_snr(raw_tsss, tsss_bench, tol)
        py_st = raw_tsss.info['proc_history'][0]['max_info']['max_st']
        assert_true(len(py_st) > 0)
        assert_equal(py_st['buflen'], st_duration)
        assert_equal(py_st['subspcorr'], 0.98)

    # Degenerate cases
    assert_raises(ValueError, maxwell_filter, raw, st_duration=10.,
                  st_correlation=0.)

def test_spatiotemporal_only():
    """Test tSSS-only processing."""
    # Load raw testing data
    tmax = 0.5
    raw = read_crop(raw_fname, (0, tmax)).load_data()
    picks = pick_types(raw.info, meg=True, exclude='bads')[::2]
    raw.pick_channels([raw.ch_names[pick] for pick in picks])
    mag_picks = pick_types(raw.info, meg='mag', exclude=())
    power = np.sqrt(np.sum(raw[mag_picks][0] ** 2))
    # basics
    raw_tsss = maxwell_filter(raw, st_duration=tmax / 2., st_only=True)
    assert_equal(len(raw.info['projs']), len(raw_tsss.info['projs']))
    assert_equal(raw_tsss.estimate_rank(), len(picks))
    _assert_shielding(raw_tsss, power, 9)
    # with movement
    head_pos = read_head_pos(pos_fname)
    raw_tsss = maxwell_filter(raw, st_duration=tmax / 2., st_only=True,
                              head_pos=head_pos)
    assert_equal(raw_tsss.estimate_rank(), len(picks))
    _assert_shielding(raw_tsss, power, 9)
    with pytest.warns(RuntimeWarning, match='st_fixed'):
        raw_tsss = maxwell_filter(raw, st_duration=tmax / 2., st_only=True,
                                  head_pos=head_pos, st_fixed=False)
    assert_equal(raw_tsss.estimate_rank(), len(picks))
    _assert_shielding(raw_tsss, power, 9)
    # should do nothing
    raw_tsss = maxwell_filter(raw, st_duration=tmax, st_correlation=1.,
                              st_only=True)
    assert_allclose(raw[:][0], raw_tsss[:][0])
    # degenerate
    pytest.raises(ValueError, maxwell_filter, raw, st_only=True)  # no ST
    # two-step process equivalent to single-step process
    raw_tsss = maxwell_filter(raw, st_duration=tmax, st_only=True)
    raw_tsss = maxwell_filter(raw_tsss)
    raw_tsss_2 = maxwell_filter(raw, st_duration=tmax)
    assert_meg_snr(raw_tsss, raw_tsss_2, 1e5)
    # now also with head movement, and a bad MEG channel
    assert_equal(len(raw.info['bads']), 0)
    bads = [raw.ch_names[0]]
    raw.info['bads'] = list(bads)
    raw_tsss = maxwell_filter(raw, st_duration=tmax, st_only=True,
                              head_pos=head_pos)
    assert_equal(raw.info['bads'], bads)
    assert_equal(raw_tsss.info['bads'], bads)  # don't reset
    raw_tsss = maxwell_filter(raw_tsss, head_pos=head_pos)
    assert_equal(raw_tsss.info['bads'], [])  # do reset MEG bads
    raw_tsss_2 = maxwell_filter(raw, st_duration=tmax, head_pos=head_pos)
    assert_equal(raw_tsss_2.info['bads'], [])
    assert_meg_snr(raw_tsss, raw_tsss_2, 1e5)

def test_bads_reconstruction():
    """Test Maxwell filter reconstruction of bad channels"""
    raw = Raw(raw_fname, allow_maxshield='yes').crop(0., 1., False)
    raw.info['bads'] = bads
    raw_sss = maxwell_filter(raw, origin=mf_head_origin, regularize=None,
                             bad_condition='ignore')
    assert_meg_snr(raw_sss, Raw(sss_bad_recon_fname), 300.)

def test_bads_reconstruction():
    """Test Maxwell filter reconstruction of bad channels"""
    with warnings.catch_warnings(record=True):  # maxshield
        raw = Raw(raw_fname, allow_maxshield=True).crop(0.0, 1.0, False)
    raw.info["bads"] = bads
    raw_sss = maxwell_filter(raw, origin=mf_head_origin, regularize=None, bad_condition="ignore")
    assert_meg_snr(raw_sss, Raw(sss_bad_recon_fname), 300.0)

def test_bads_reconstruction():
    """Test Maxwell filter reconstruction of bad channels"""
    with warnings.catch_warnings(record=True):  # maxshield
        raw = Raw(raw_fname, allow_maxshield=True).crop(0., 1., False)
    raw.info['bads'] = bads
    raw_sss = maxwell_filter(raw)
    _assert_snr(raw_sss, Raw(sss_bad_recon_fname), 300.)