def test_spm_hrf_octave():
    # Test SPM hrf against output from SPM code running in Octave
    my_path = dirname(__file__)
    hrfs_path = pjoin(my_path, 'spm_hrfs.mat')
    # mat file resulting from make_hrfs.m
    hrfs_mat = sio.loadmat(hrfs_path, squeeze_me=True)
    params = hrfs_mat['params']
    hrfs = hrfs_mat['hrfs']
    for i, pvec in enumerate(params):
        dt, ppk, upk, pdsp, udsp, rat = pvec
        t_vec = np.arange(0, 32.1, dt)
        our_hrf = spm_hrf_compat(t_vec,
                                 peak_delay=ppk,
                                 peak_disp=pdsp,
                                 under_delay=upk,
                                 under_disp=udsp,
                                 p_u_ratio=rat)
        # Normalize integral to match SPM
        assert_almost_equal(our_hrf, hrfs[i])
    # Test basis functions
    # mat file resulting from get_td_dd.m
    bases_path = pjoin(my_path, 'spm_bases.mat')
    bases_mat = sio.loadmat(bases_path, squeeze_me=True)
    dt = bases_mat['dt']
    t_vec = np.arange(0, 32 + dt, dt)
    # SPM function divides by sum of values - revert with dt
    assert_almost_equal(spmt(t_vec), bases_mat['hrf'] / dt, 4)
    assert_almost_equal(dspmt(t_vec), bases_mat['dhrf'] / dt, 4)
    assert_almost_equal(ddspmt(t_vec), bases_mat['ddhrf'] / dt, 4)