# Authors: Denis A. Engemann  <denis.engemann@gmail.com>
#          Martin Luessi <mluessi@nmr.mgh.harvard.edu>
#          Alexandre Gramfort <alexandre.gramfort@telecom-paristech.fr>
#          Matti Hamalainen <msh@nmr.mgh.harvard.edu>
#          Yuval Harpaz <yuvharpaz@gmail.com>
#
#          simplified BSD-3 license

import os.path as op
from itertools import count
import numpy as np

from ...utils import logger, verbose, sum_squared
from ..constants import FIFF
from ..base import _BaseRaw
from .constants import BTI
from .read import (read_int32, read_int16, read_str, read_float, read_double,
                   read_transform, read_char, read_int64, read_uint16,
                   read_uint32, read_double_matrix, read_float_matrix,
                   read_int16_matrix)
from .transforms import (bti_identity_trans, bti_to_vv_trans,
                         bti_to_vv_coil_trans, inverse_trans, merge_trans)
from ..meas_info import Info
from ...externals import six

FIFF_INFO_CHS_FIELDS = ('loc', 'ch_name', 'unit_mul', 'coil_trans',
                        'coord_frame', 'coil_type', 'range', 'unit', 'cal',
                        'eeg_loc', 'scanno', 'kind', 'logno')

FIFF_INFO_CHS_DEFAULTS = (np.array([0, 0, 0, 1] * 3, dtype='f4'),
                          None, 0, None, 0, 0, 1.0,
                          107, 1.0, None, None, 402, None)

FIFF_INFO_DIG_FIELDS = ('kind', 'ident', 'r', 'coord_frame')
FIFF_INFO_DIG_DEFAULTS = (None, None, None, FIFF.FIFFV_COORD_HEAD)

BTI_WH2500_REF_MAG = ['MxA', 'MyA', 'MzA', 'MxaA', 'MyaA', 'MzaA']
BTI_WH2500_REF_GRAD = ['GxxA', 'GyyA', 'GyxA', 'GzaA', 'GzyA']

dtypes = zip(list(range(1, 5)), ('>i2', '>i4', '>f4', '>f8'))
DTYPES = dict((i, np.dtype(t)) for i, t in dtypes)

RAW_INFO_FIELDS = ['dev_head_t', 'nchan', 'bads', 'projs', 'dev_ctf_t',
                   'meas_date', 'meas_id', 'dig', 'sfreq', 'highpass',
                   'comps', 'chs', 'ch_names', 'file_id',
                   'lowpass', 'acq_pars', 'acq_stim', 'filename',
                   'ctf_head_t']


def _rename_channels(names, ecg_ch='E31', eog_ch=('E63', 'E64')):
    """Renames appropriately ordered list of channel names

    Parameters
    ----------
    names : list of str
        Lists of 4-D channel names in ascending order

    Returns
    -------
    new : list
        List of names, channel names in Neuromag style
    """
    new = list()
    ref_mag, ref_grad, eog, eeg, ext = [count(1) for _ in range(5)]
    for i, name in enumerate(names, 1):
        if name.startswith('A'):
            name = 'MEG %3.3d' % i
        elif name == 'RESPONSE':
            name = 'STI 013'
        elif name == 'TRIGGER':
            name = 'STI 014'
        elif any([name == k for k in eog_ch]):
            name = 'EOG %3.3d' % six.advance_iterator(eog)
        elif name == ecg_ch:
            name = 'ECG 001'
        elif name.startswith('E'):
            name = 'EEG %3.3d' % six.advance_iterator(eeg)
        elif name == 'UACurrent':
            name = 'UTL 001'
        elif name.startswith('M'):
            name = 'RFM %3.3d' % six.advance_iterator(ref_mag)
        elif name.startswith('G'):
            name = 'RFG %3.3d' % six.advance_iterator(ref_grad)
        elif name.startswith('X'):
            name = 'EXT %3.3d' % six.advance_iterator(ext)

        new += [name]

    return new


def _read_head_shape(fname):
    """ Helper Function """
    with open(fname, 'rb') as fid:
        fid.seek(BTI.FILE_HS_N_DIGPOINTS)
        _n_dig_points = read_int32(fid)
        idx_points = read_double_matrix(fid, BTI.DATA_N_IDX_POINTS, 3)
        dig_points = read_double_matrix(fid, _n_dig_points, 3)

    return idx_points, dig_points


def _convert_head_shape(idx_points, dig_points):
    """ Helper function """

    fp = idx_points.astype('>f8')
    dp = np.sum(fp[2] * (fp[0] - fp[1]))
    tmp1, tmp2 = sum_squared(fp[2]), sum_squared(fp[0] - fp[1])
    dcos = -dp / np.sqrt(tmp1 * tmp2)
    dsin = np.sqrt(1. - dcos * dcos)
    dt = dp / np.sqrt(tmp2)

    idx_points_nm = np.ones((len(fp), 3), dtype='>f8')
    for idx, f in enumerate(fp):
        idx_points_nm[idx, 0] = dcos * f[0] - dsin * f[1] + dt
        idx_points_nm[idx, 1] = dsin * f[0] + dcos * f[1]
        idx_points_nm[idx, 2] = f[2]

    # adjust order of fiducials to Neuromag
    idx_points_nm[[1, 2]] = idx_points_nm[[2, 1]]

    t = bti_identity_trans('>f8')
    t[0, 0] = dcos
    t[0, 1] = -dsin
    t[1, 0] = dsin
    t[1, 1] = dcos
    t[0, 3] = dt

    dig_points_nm = np.dot(t[BTI.T_ROT_IX], dig_points.T).T
    dig_points_nm += t[BTI.T_TRANS_IX].T

    return idx_points_nm, dig_points_nm, t


def _setup_head_shape(fname, use_hpi=True):
    """Read index points and dig points from BTi head shape file

    Parameters
    ----------
    fname : str
        The absolute path to the head shape file

    Returns
    -------
    dig : list of dicts
        The list of dig point info structures needed for the fiff info
        structure.
    use_hpi : bool
        Whether to treat additional hpi coils as digitization points or not.
        If False, hpi coils will be discarded.
    """
    idx_points, dig_points = _read_head_shape(fname)
    idx_points, dig_points, t = _convert_head_shape(idx_points, dig_points)
    all_points = np.r_[idx_points, dig_points].astype('>f4')

    idx_idents = list(range(1, 4)) + list(range(1, (len(idx_points) + 1) - 3))
    dig = []
    for idx in range(all_points.shape[0]):
        point_info = dict(zip(FIFF_INFO_DIG_FIELDS, FIFF_INFO_DIG_DEFAULTS))
        point_info['r'] = all_points[idx]
        if idx < 3:
            point_info['kind'] = FIFF.FIFFV_POINT_CARDINAL
            point_info['ident'] = idx_idents[idx]
        if 2 < idx < len(idx_points) and use_hpi:
            point_info['kind'] = FIFF.FIFFV_POINT_HPI
            point_info['ident'] = idx_idents[idx]
        elif idx > 4:
            point_info['kind'] = FIFF.FIFFV_POINT_EXTRA
            point_info['ident'] = (idx + 1) - len(idx_idents)

        if 2 < idx < len(idx_points) and not use_hpi:
            pass
        else:
            dig += [point_info]

    return dig, t


def _convert_coil_trans(coil_trans, bti_trans, bti_to_nm):
    """ Helper Function """
    t = bti_to_vv_coil_trans(coil_trans, bti_trans, bti_to_nm)
    loc = np.roll(t.copy().T, 1, 0)[:, :3].flatten()

    return t, loc


def _convert_dev_head_t(bti_trans, bti_to_nm, m_h_nm_h):
    """ Helper Function """
    nm_to_m_sensor = inverse_trans(bti_identity_trans(), bti_to_nm)
    nm_sensor_m_head = merge_trans(bti_trans, nm_to_m_sensor)

    nm_dev_head_t = merge_trans(m_h_nm_h, nm_sensor_m_head)
    nm_dev_head_t[3, :3] = 0.

    return nm_dev_head_t


def _correct_offset(fid):
    """ Align fid pointer """
    current = fid.tell()
    if ((current % BTI.FILE_CURPOS) != 0):
        offset = current % BTI.FILE_CURPOS
        fid.seek(BTI.FILE_CURPOS - (offset), 1)


def _read_config(fname):
    """Read BTi system config file

    Parameters
    ----------
    fname : str
        The absolute path to the config file

    Returns
    -------
    cfg : dict
        The config blocks found.

    """
    with open(fname, 'rb') as fid:
        cfg = dict()
        cfg['hdr'] = {'version': read_int16(fid),
                      'site_name': read_str(fid, 32),
                      'dap_hostname': read_str(fid, 16),
                      'sys_type': read_int16(fid),
                      'sys_options': read_int32(fid),
                      'supply_freq': read_int16(fid),
                      'total_chans': read_int16(fid),
                      'system_fixed_gain': read_float(fid),
                      'volts_per_bit': read_float(fid),
                      'total_sensors': read_int16(fid),
                      'total_user_blocks': read_int16(fid),
                      'next_der_chan_no': read_int16(fid)}

        fid.seek(2, 1)

        cfg['checksum'] = read_uint32(fid)
        cfg['reserved'] = read_char(fid, 32)
        cfg['transforms'] = [read_transform(fid) for t in
                             range(cfg['hdr']['total_sensors'])]

        cfg['user_blocks'] = dict()
        for block in range(cfg['hdr']['total_user_blocks']):
            ub = dict()

            ub['hdr'] = {'nbytes': read_int32(fid),
                         'kind': read_str(fid, 20),
                         'checksum': read_int32(fid),
                         'username': read_str(fid, 32),
                         'timestamp': read_int32(fid),
                         'user_space_size': read_int32(fid),
                         'reserved': read_char(fid, 32)}

            _correct_offset(fid)
            kind = ub['hdr'].pop('kind')
            if not kind:  # make sure reading goes right. Should never be empty
                raise RuntimeError('Could not read user block. Probably you '
                                   'acquired data using a BTi version '
                                   'currently not supported. Please contact '
                                   'the mne-python developers.')
            dta, cfg['user_blocks'][kind] = dict(), ub
            if kind in [v for k, v in BTI.items() if k[:5] == 'UB_B_']:
                if kind == BTI.UB_B_MAG_INFO:
                    dta['version'] = read_int32(fid)
                    fid.seek(20, 1)
                    dta['headers'] = list()
                    for hdr in range(6):
                        d = {'name': read_str(fid, 16),
                             'transform': read_transform(fid),
                             'units_per_bit': read_float(fid)}
                        dta['headers'] += [d]
                        fid.seek(20, 1)

                elif kind == BTI.UB_B_COH_POINTS:
                    dta['n_points'] = read_int32(fid)
                    dta['status'] = read_int32(fid)
                    dta['points'] = []
                    for pnt in range(16):
                        d = {'pos': read_double_matrix(fid, 1, 3),
                             'direction': read_double_matrix(fid, 1, 3),
                             'error': read_double(fid)}
                        dta['points'] += [d]

                elif kind == BTI.UB_B_CCP_XFM_BLOCK:
                    dta['method'] = read_int32(fid)
                    # handle difference btw/ linux (0) and solaris (4)
                    size = 0 if ub['hdr']['user_space_size'] == 132 else 4
                    fid.seek(size, 1)
                    dta['transform'] = read_transform(fid)

                elif kind == BTI.UB_B_EEG_LOCS:
                    dta['electrodes'] = []
                    while True:
                        d = {'label': read_str(fid, 16),
                             'location': read_double_matrix(fid, 1, 3)}
                        if not d['label']:
                            break
                        dta['electrodes'] += [d]

                elif kind in [BTI.UB_B_WHC_CHAN_MAP_VER,
                              BTI.UB_B_WHS_SUBSYS_VER]:
                    dta['version'] = read_int16(fid)
                    dta['struct_size'] = read_int16(fid)
                    dta['entries'] = read_int16(fid)

                    fid.seek(8, 1)

                elif kind == BTI.UB_B_WHC_CHAN_MAP:
                    num_channels = None
                    for name, data in cfg['user_blocks'].items():
                        if name == BTI.UB_B_WHC_CHAN_MAP_VER:
                            num_channels = data['entries']
                            break

                    if num_channels is None:
                        raise ValueError('Cannot find block %s to determine '
                                         'number of channels'
                                         % BTI.UB_B_WHC_CHAN_MAP_VER)

                    dta['channels'] = list()
                    for i in range(num_channels):
                        d = {'subsys_type': read_int16(fid),
                             'subsys_num': read_int16(fid),
                             'card_num': read_int16(fid),
                             'chan_num': read_int16(fid),
                             'recdspnum': read_int16(fid)}
                        dta['channels'] += [d]
                        fid.seek(8, 1)

                elif kind == BTI.UB_B_WHS_SUBSYS:
                    num_subsys = None
                    for name, data in cfg['user_blocks'].items():
                        if name == BTI.UB_B_WHS_SUBSYS_VER:
                            num_subsys = data['entries']
                            break

                    if num_subsys is None:
                        raise ValueError('Cannot find block %s to determine'
                                         ' number of subsystems'
                                         % BTI.UB_B_WHS_SUBSYS_VER)

                    dta['subsys'] = list()
                    for sub_key in range(num_subsys):
                        d = {'subsys_type': read_int16(fid),
                             'subsys_num': read_int16(fid),
                             'cards_per_sys': read_int16(fid),
                             'channels_per_card': read_int16(fid),
                             'card_version': read_int16(fid)}

                        fid.seek(2, 1)

                        d.update({'offsetdacgain': read_float(fid),
                                  'squid_type': read_int32(fid),
                                  'timesliceoffset': read_int16(fid),
                                  'padding': read_int16(fid),
                                  'volts_per_bit': read_float(fid)})

                        dta['subsys'] += [d]

                elif kind == BTI.UB_B_CH_LABELS:
                    dta['version'] = read_int32(fid)
                    dta['entries'] = read_int32(fid)
                    fid.seek(16, 1)

                    dta['labels'] = list()
                    for label in range(dta['entries']):
                        dta['labels'] += [read_str(fid, 16)]

                elif kind == BTI.UB_B_CALIBRATION:
                    dta['sensor_no'] = read_int16(fid)
                    fid.seek(2, 1)
                    dta['timestamp'] = read_int32(fid)
                    dta['logdir'] = read_str(fid, 256)

                elif kind == BTI.UB_B_SYS_CONFIG_TIME:
                    # handle difference btw/ linux (256) and solaris (512)
                    size = 256 if ub['hdr']['user_space_size'] == 260 else 512
                    dta['sysconfig_name'] = read_str(fid, size)
                    dta['timestamp'] = read_int32(fid)

                elif kind == BTI.UB_B_DELTA_ENABLED:
                    dta['delta_enabled'] = read_int16(fid)

                elif kind in [BTI.UB_B_E_TABLE_USED, BTI.UB_B_E_TABLE]:
                    dta['hdr'] = {'version': read_int32(fid),
                                  'entry_size': read_int32(fid),
                                  'n_entries': read_int32(fid),
                                  'filtername': read_str(fid, 16),
                                  'n_e_values': read_int32(fid),
                                  'reserved': read_str(fid, 28)}

                    if dta['hdr']['version'] == 2:
                        size = 16
                        dta['ch_names'] = [read_str(fid, size) for ch in
                                           range(dta['hdr']['n_entries'])]
                        dta['e_ch_names'] = [read_str(fid, size) for ch in
                                             range(dta['hdr']['n_e_values'])]

                        rows = dta['hdr']['n_entries']
                        cols = dta['hdr']['n_e_values']
                        dta['etable'] = read_float_matrix(fid, rows, cols)
                    else:  # handle MAGNES2500 naming scheme
                        dta['ch_names'] = ['WH2500'] * dta['hdr']['n_e_values']
                        dta['hdr']['n_e_values'] = 6
                        dta['e_ch_names'] = BTI_WH2500_REF_MAG
                        rows = dta['hdr']['n_entries']
                        cols = dta['hdr']['n_e_values']
                        dta['etable'] = read_float_matrix(fid, rows, cols)

                        _correct_offset(fid)

                elif any([kind == BTI.UB_B_WEIGHTS_USED,
                          kind[:4] == BTI.UB_B_WEIGHT_TABLE]):
                    dta['hdr'] = {'version': read_int32(fid),
                                  'entry_size': read_int32(fid),
                                  'n_entries': read_int32(fid),
                                  'name': read_str(fid, 32),
                                  'description': read_str(fid, 80),
                                  'n_anlg': read_int32(fid),
                                  'n_dsp': read_int32(fid),
                                  'reserved': read_str(fid, 72)}

                    if dta['hdr']['version'] == 2:
                        dta['ch_names'] = [read_str(fid, 16) for ch in
                                           range(dta['hdr']['n_entries'])]
                        dta['anlg_ch_names'] = [read_str(fid, 16) for ch in
                                                range(dta['hdr']['n_anlg'])]

                        dta['dsp_ch_names'] = [read_str(fid, 16) for ch in
                                               range(dta['hdr']['n_dsp'])]

                        rows = dta['hdr']['n_entries']
                        cols = dta['hdr']['n_dsp']
                        dta['dsp_wts'] = read_float_matrix(fid, rows, cols)
                        cols = dta['hdr']['n_anlg']
                        dta['anlg_wts'] = read_int16_matrix(fid, rows, cols)

                    else:  # handle MAGNES2500 naming scheme
                        dta['ch_names'] = ['WH2500'] * dta['hdr']['n_entries']
                        dta['anlg_ch_names'] = BTI_WH2500_REF_MAG[:3]
                        dta['hdr']['n_anlg'] = len(dta['anlg_ch_names'])
                        dta['dsp_ch_names'] = BTI_WH2500_REF_GRAD
                        dta['hdr.n_dsp'] = len(dta['dsp_ch_names'])
                        dta['anlg_wts'] = np.zeros((dta['hdr']['n_entries'],
                                                    dta['hdr']['n_anlg']),
                                                   dtype='i2')
                        dta['dsp_wts'] = np.zeros((dta['hdr']['n_entries'],
                                                   dta['hdr']['n_dsp']),
                                                  dtype='f4')
                        for n in range(dta['hdr']['n_entries']):
                            dta['anlg_wts'][d] = read_int16_matrix(fid, 1,
                                                       dta['hdr']['n_anlg'])
                            read_int16(fid)
                            dta['dsp_wts'][d] = read_float_matrix(fid, 1,
                                                        dta['hdr']['n_dsp'])

                        _correct_offset(fid)

                elif kind == BTI.UB_B_TRIG_MASK:
                    dta['version'] = read_int32(fid)
                    dta['entries'] = read_int32(fid)
                    fid.seek(16, 1)

                    dta['masks'] = []
                    for entry in range(dta['entries']):
                        d = {'name': read_str(fid, 20),
                             'nbits': read_uint16(fid),
                             'shift': read_uint16(fid),
                             'mask': read_uint32(fid)}
                        dta['masks'] += [d]
                        fid.seek(8, 1)

            else:
                dta['unknown'] = {'hdr': read_char(fid,
                                  ub['hdr']['user_space_size'])}

            ub.update(dta)  # finally update the userblock data
            _correct_offset(fid)  # after reading.

        cfg['chs'] = list()

        # prepare reading channels
        dev_header = lambda x: {'size': read_int32(x),
                                'checksum': read_int32(x),
                                'reserved': read_str(x, 32)}

        for channel in range(cfg['hdr']['total_chans']):
            ch = {'name': read_str(fid, 16),
                  'chan_no': read_int16(fid),
                  'ch_type': read_uint16(fid),
                  'sensor_no': read_int16(fid),
                  'data': dict()}

            fid.seek(2, 1)
            ch.update({'gain': read_float(fid),
                       'units_per_bit': read_float(fid),
                       'yaxis_label': read_str(fid, 16),
                       'aar_val': read_double(fid),
                       'checksum': read_int32(fid),
                       'reserved': read_str(fid, 32)})

            cfg['chs'] += [ch]
            _correct_offset(fid)  # before and after
            dta = dict()
            if ch['ch_type'] in [BTI.CHTYPE_MEG, BTI.CHTYPE_REFERENCE]:
                dev = {'device_info': dev_header(fid),
                       'inductance': read_float(fid),
                       'padding': read_str(fid, 4),
                       'transform': read_transform(fid),
                       'xform_flag': read_int16(fid),
                       'total_loops': read_int16(fid)}

                fid.seek(4, 1)
                dev['reserved'] = read_str(fid, 32)
                dta.update({'dev': dev, 'loops': []})
                for loop in range(dev['total_loops']):
                    d = {'position': read_double_matrix(fid, 1, 3),
                         'orientation': read_double_matrix(fid, 1, 3),
                         'radius': read_double(fid),
                         'wire_radius': read_double(fid),
                         'turns': read_int16(fid)}
                    fid.seek(2, 1)
                    d['checksum'] = read_int32(fid)
                    d['reserved'] = read_str(fid, 32)
                    dta['loops'] += [d]

            elif ch['ch_type'] == BTI.CHTYPE_EEG:
                dta = {'device_info': dev_header(fid),
                       'impedance': read_float(fid),
                       'padding': read_str(fid, 4),
                       'transform': read_transform(fid),
                       'reserved': read_char(fid, 32)}

            elif ch['ch_type'] == BTI.CHTYPE_EXTERNAL:
                dta = {'device_info': dev_header(fid),
                       'user_space_size': read_int32(fid),
                       'reserved': read_str(fid, 32)}

            elif ch['ch_type'] == BTI.CHTYPE_TRIGGER:
                dta = {'device_info': dev_header(fid),
                       'user_space_size': read_int32(fid)}
                fid.seek(2, 1)
                dta['reserved'] = read_str(fid, 32)

            elif ch['ch_type'] in [BTI.CHTYPE_UTILITY, BTI.CHTYPE_DERIVED]:
                dta = {'device_info': dev_header(fid),
                       'user_space_size': read_int32(fid),
                       'reserved': read_str(fid, 32)}

            elif ch['ch_type'] == BTI.CHTYPE_SHORTED:
                dta = {'device_info': dev_header(fid),
                       'reserved': read_str(fid, 32)}

            ch.update(dta)  # add data collected
            _correct_offset(fid)  # after each reading

    return cfg


def _read_epoch(fid):
    """Read BTi PDF epoch"""
    out = {'pts_in_epoch': read_int32(fid),
           'epoch_duration': read_float(fid),
           'expected_iti': read_float(fid),
           'actual_iti': read_float(fid),
           'total_var_events': read_int32(fid),
           'checksum': read_int32(fid),
           'epoch_timestamp': read_int32(fid)}

    fid.seek(28, 1)

    return out


def _read_channel(fid):
    """Read BTi PDF channel"""
    out = {'chan_label': read_str(fid, 16),
           'chan_no': read_int16(fid),
           'attributes': read_int16(fid),
           'scale': read_float(fid),
           'yaxis_label': read_str(fid, 16),
           'valid_min_max': read_int16(fid)}

    fid.seek(6, 1)
    out.update({'ymin': read_double(fid),
                'ymax': read_double(fid),
                'index': read_int32(fid),
                'checksum': read_int32(fid),
                'off_flag': read_str(fid, 16),
                'offset': read_float(fid)})

    fid.seek(12, 1)

    return out


def _read_event(fid):
    """Read BTi PDF event"""
    out = {'event_name': read_str(fid, 16),
           'start_lat': read_float(fid),
           'end_lat': read_float(fid),
           'step_size': read_float(fid),
           'fixed_event': read_int16(fid),
           'checksum': read_int32(fid)}

    fid.seek(32, 1)
    _correct_offset(fid)

    return out


def _read_process(fid):
    """Read BTi PDF process"""

    out = {'nbytes': read_int32(fid),
           'process_type': read_str(fid, 20),
           'checksum': read_int32(fid),
           'user': read_str(fid, 32),
           'timestamp': read_int32(fid),
           'filename': read_str(fid, 256),
           'total_steps': read_int32(fid)}

    fid.seek(32, 1)
    _correct_offset(fid)
    out['processing_steps'] = list()
    for step in range(out['total_steps']):
        this_step = {'nbytes': read_int32(fid),
                     'process_type': read_str(fid, 20),
                     'checksum': read_int32(fid)}
        ptype = this_step['process_type']
        if ptype == BTI.PROC_DEFAULTS:
            this_step['scale_option'] = read_int32(fid)

            fid.seek(4, 1)
            this_step['scale'] = read_double(fid)
            this_step['dtype'] = read_int32(fid)
            this_step['selected'] = read_int16(fid)
            this_step['color_display'] = read_int16(fid)

            fid.seek(32, 1)
        elif ptype in BTI.PROC_FILTER:
            this_step['freq'] = read_float(fid)
            fid.seek(32, 1)
        elif ptype in BTI.PROC_BPFILTER:
            this_step['high_freq'] = read_float(fid)
            this_step['low_frew'] = read_float(fid)
        else:
            jump = this_step['user_space_size'] = read_int32(fid)
            fid.seek(32, 1)
            fid.seek(jump, 1)

        out['processing_steps'] += [this_step]
        _correct_offset(fid)

    return out


def _read_assoc_file(fid):
    """Read BTi PDF assocfile"""

    out = {'file_id': read_int16(fid),
           'length': read_int16(fid)}

    fid.seek(32, 1)
    out['checksum'] = read_int32(fid)

    return out


def _read_pfid_ed(fid):
    """Read PDF ed file"""

    out = {'comment_size': read_int32(fid),
           'name': read_str(fid, 17)}

    fid.seek(9, 1)
    out.update({'pdf_number': read_int16(fid),
                'total_events': read_int32(fid),
                'timestamp': read_int32(fid),
                'flags': read_int32(fid),
                'de_process': read_int32(fid),
                'checksum': read_int32(fid),
                'ed_id': read_int32(fid),
                'win_width': read_float(fid),
                'win_offset': read_float(fid)})

    fid.seek(8, 1)

    return out


def _read_coil_def(fid):
    """ Read coil definition """
    coildef = {'position': read_double_matrix(fid, 1, 3),
               'orientation': read_double_matrix(fid, 1, 3),
               'radius': read_double(fid),
               'wire_radius': read_double(fid),
               'turns': read_int16(fid)}

    fid.seek(fid, 2, 1)
    coildef['checksum'] = read_int32(fid)
    coildef['reserved'] = read_str(fid, 32)


def _read_ch_config(fid):
    """Read BTi channel config"""

    cfg = {'name': read_str(fid, BTI.FILE_CONF_CH_NAME),
           'chan_no': read_int16(fid),
           'ch_type': read_uint16(fid),
           'sensor_no': read_int16(fid)}

    fid.seek(fid, BTI.FILE_CONF_CH_NEXT, 1)

    cfg.update({'gain': read_float(fid),
                'units_per_bit': read_float(fid),
                'yaxis_label': read_str(fid, BTI.FILE_CONF_CH_YLABEL),
                'aar_val': read_double(fid),
                'checksum': read_int32(fid),
                'reserved': read_str(fid, BTI.FILE_CONF_CH_RESERVED)})

    _correct_offset(fid)

    # Then the channel info
    ch_type, chan = cfg['ch_type'], dict()
    chan['dev'] = {'size': read_int32(fid),
                   'checksum': read_int32(fid),
                   'reserved': read_str(fid, 32)}
    if ch_type in [BTI.CHTYPE_MEG, BTI.CHTYPE_REF]:
        chan['loops'] = [_read_coil_def(fid) for d in
                         range(chan['dev']['total_loops'])]

    elif ch_type == BTI.CHTYPE_EEG:
        chan['impedance'] = read_float(fid)
        chan['padding'] = read_str(fid, BTI.FILE_CONF_CH_PADDING)
        chan['transform'] = read_transform(fid)
        chan['reserved'] = read_char(fid, BTI.FILE_CONF_CH_RESERVED)

    elif ch_type in [BTI.CHTYPE_TRIGGER,  BTI.CHTYPE_EXTERNAL,
                     BTI.CHTYPE_UTILITY, BTI.CHTYPE_DERIVED]:
        chan['user_space_size'] = read_int32(fid)
        if ch_type == BTI.CHTYPE_TRIGGER:
            fid.seek(2, 1)
        chan['reserved'] = read_str(fid, BTI.FILE_CONF_CH_RESERVED)

    elif ch_type == BTI.CHTYPE_SHORTED:
        chan['reserved'] = read_str(fid, BTI.FILE_CONF_CH_RESERVED)

    cfg['chan'] = chan

    _correct_offset(fid)

    return cfg


def _read_bti_header(pdf_fname, config_fname):
    """ Read bti PDF header
    """
    with open(pdf_fname, 'rb') as fid:
        fid.seek(-8, 2)
        start = fid.tell()
        header_position = read_int64(fid)
        check_value = header_position & BTI.FILE_MASK

        if ((start + BTI.FILE_CURPOS - check_value) <= BTI.FILE_MASK):
            header_position = check_value

        # Check header position for alignment issues
        if ((header_position % 8) != 0):
            header_position += (8 - (header_position % 8))

        fid.seek(header_position, 0)

        # actual header starts here
        info = {'version': read_int16(fid),
                'file_type': read_str(fid, 5),
                'hdr_size': start - header_position,  # add to info for convenience
                'start': start}

        fid.seek(1, 1)

        info.update({'data_format': read_int16(fid),
                    'acq_mode': read_int16(fid),
                    'total_epochs': read_int32(fid),
                    'input_epochs': read_int32(fid),
                    'total_events': read_int32(fid),
                    'total_fixed_events': read_int32(fid),
                    'sample_period': read_float(fid),
                    'xaxis_label': read_str(fid, 16),
                    'total_processes': read_int32(fid),
                    'total_chans': read_int16(fid)})

        fid.seek(2, 1)
        info.update({'checksum': read_int32(fid),
                    'total_ed_classes': read_int32(fid),
                    'total_associated_files': read_int16(fid),
                    'last_file_index': read_int16(fid),
                    'timestamp': read_int32(fid)})

        fid.seek(20, 1)
        _correct_offset(fid)

        # actual header ends here, so dar seems ok.

        info['epochs'] = [_read_epoch(fid) for epoch in
                          range(info['total_epochs'])]

        info['chs'] = [_read_channel(fid) for ch in
                       range(info['total_chans'])]

        info['events'] = [_read_event(fid) for event in
                          range(info['total_events'])]

        info['processes'] = [_read_process(fid) for process in
                             range(info['total_processes'])]

        info['assocfiles'] = [_read_assoc_file(fid) for af in
                              range(info['total_associated_files'])]

        info['edclasses'] = [_read_pfid_ed(fid) for ed_class in
                             range(info['total_ed_classes'])]

        info['extra_data'] = fid.read(start - fid.tell())
        info['pdf_fname'] = pdf_fname

    info['total_slices'] = sum(e['pts_in_epoch'] for e in
                               info['epochs'])

    info['dtype'] = DTYPES[info['data_format']]
    bps = info['dtype'].itemsize * info['total_chans']
    info['bytes_per_slice'] = bps

    cfg = _read_config(config_fname)
    info['bti_transform'] = cfg['transforms']

    # augment channel list by according info from config.
    # get channels from config present in PDF
    chans = info['chs']
    chans_cfg = [c for c in cfg['chs'] if c['chan_no']
                 in [c_['chan_no'] for c_ in chans]]

    # check all pdf chanels are present in config
    match = [c['chan_no'] for c in chans_cfg] == \
            [c['chan_no'] for c in chans]

    if not match:
        raise RuntimeError('Could not match raw data channels with'
                           ' config channels. Some of the channels'
                           ' found are not described in config.')

    # transfer channel info from config to channel info
    for ch, ch_cfg in zip(chans, chans_cfg):
        ch['upb'] = ch_cfg['units_per_bit']
        ch['gain'] = ch_cfg['gain']
        ch['name'] = ch_cfg['name']
        ch['coil_trans'] = (ch_cfg['dev'].get('transform', None)
                            if 'dev' in ch_cfg else None)
        if info['data_format'] <= 2:
            ch['cal'] = ch['scale'] * ch['upb'] * (ch['gain'] ** -1)
        else:
            ch['cal'] = ch['scale'] * ch['gain']

    by_index = [(i, d['index']) for i, d in enumerate(chans)]
    by_index.sort(key=lambda c: c[1])
    by_index = [idx[0] for idx in by_index]
    info['chs'] = [chans[pos] for pos in by_index]

    by_name = [(i, d['name']) for i, d in enumerate(info['chs'])]
    a_chs = filter(lambda c: c[1].startswith('A'), by_name)
    other_chs = filter(lambda c: not c[1].startswith('A'), by_name)
    by_name = sorted(a_chs, key=lambda c: int(c[1][1:])) + sorted(other_chs)

    by_name = [idx[0] for idx in by_name]
    info['chs'] = [chans[pos] for pos in by_name]
    info['order'] = by_name

    # finally add some important fields from the config
    info['e_table'] = cfg['user_blocks'][BTI.UB_B_E_TABLE_USED]
    info['weights'] = cfg['user_blocks'][BTI.UB_B_WEIGHTS_USED]

    return info


def _read_data(info, start=None, stop=None):
    """ Helper function: read Bti processed data file (PDF)

    Parameters
    ----------
    info : dict
        The measurement info.
    start : int | None
        The number of the first time slice to read. If None, all data will
        be read from the beginning.
    stop : int | None
        The number of the last time slice to read. If None, all data will
        be read to the end.
    dtype : str | dtype object
        The type the data are casted to.

    Returns
    -------
    data : ndarray
        The measurement data, a channels x time slices array.
    """

    total_slices = info['total_slices']
    if start is None:
        start = 0
    if stop is None:
        stop = total_slices

    if any([start < 0, stop > total_slices, start >= stop]):
        raise RuntimeError('Invalid data range supplied:'
                           ' %d, %d' % (start, stop))

    with open(info['pdf_fname'], 'rb') as fid:
        fid.seek(info['bytes_per_slice'] * start, 0)
        cnt = (stop - start) * info['total_chans']
        shape = [stop - start, info['total_chans']]
        data = np.fromfile(fid, dtype=info['dtype'],
                           count=cnt).astype('f4').reshape(shape)

    for ch in info['chs']:
        data[:, ch['index']] *= ch['cal']

    return data[:, info['order']].T


class RawBTi(_BaseRaw):
    """ Raw object from 4D Neuroimaging MagnesWH3600 data

    Parameters
    ----------
    pdf_fname : str | None
        absolute path to the processed data file (PDF)
    config_fname : str | None
        absolute path to system config file. If None, it is assumed to be in
        the same directory.
    head_shape_fname : str
        absolute path to the head shape file. If None, it is assumed to be in
        the same directory.
    rotation_x : float | int | None
        Degrees to tilt x-axis for sensor frame misalignment.
        If None, no adjustment will be applied.
    translation : array-like
        The translation to place the origin of coordinate system
        to the center of the head.
    ecg_ch: str | None
      The 4D name of the ECG channel. If None, the channel will be treated
      as regular EEG channel.
    eog_ch: tuple of str | None
      The 4D names of the EOG channels. If None, the channels will be treated
      as regular EEG channels.
    verbose : bool, str, int, or None
        If not None, override default verbose level (see mne.verbose).

    Attributes & Methods
    --------------------
    See documentation for mne.io.Raw

    """
    @verbose
    def __init__(self, pdf_fname, config_fname='config',
                 head_shape_fname='hs_file', rotation_x=None,
                 translation=(0.0, 0.02, 0.11), ecg_ch='E31',
                 eog_ch=('E63', 'E64'), verbose=None):

        if not op.isabs(pdf_fname):
            pdf_fname = op.abspath(pdf_fname)

        if not op.isabs(config_fname):
            config_fname = op.join(op.dirname(pdf_fname), config_fname)

        if not op.exists(config_fname):
            raise ValueError('Could not find the config file %s. Please check'
                             ' whether you are in the right directory '
                             'or pass the full name' % config_fname)

        if not op.isabs(head_shape_fname):
            head_shape_fname = op.join(op.dirname(pdf_fname), head_shape_fname)

        if not op.exists(head_shape_fname):
            raise ValueError('Could not find the head_shape file %s. You shoul'
                             'd check whether you are in the right directory o'
                             'r pass the full file name.' % head_shape_fname)

        logger.info('Reading 4D PDF file %s...' % pdf_fname)
        bti_info = _read_bti_header(pdf_fname, config_fname)

         # XXX indx is informed guess. Normally only one transform is stored.
        dev_ctf_t = bti_info['bti_transform'][0].astype('>f8')
        bti_to_nm = bti_to_vv_trans(adjust=rotation_x,
                                    translation=translation, dtype='>f8')

        use_hpi = False  # hard coded, but marked as later option.
        logger.info('Creating Neuromag info structure ...')
        info = Info()
        info['bads'] = []
        info['meas_id'] = None
        info['file_id'] = None
        info['projs'] = list()
        info['comps'] = list()
        date = bti_info['processes'][0]['timestamp']
        info['meas_date'] = [date, 0]
        info['sfreq'] = 1e3 / bti_info['sample_period'] * 1e-3
        info['nchan'] = len(bti_info['chs'])

        # browse processing info for filter specs.
        hp, lp = 0.0, info['sfreq'] * 0.4  # find better default
        for proc in bti_info['processes']:
            if 'filt' in proc['process_type']:
                for step in proc['processing_steps']:
                    if 'high_freq' in step:
                        hp, lp = step['high_freq'], step['low_freq']
                    elif 'hp' in step['process_type']:
                        hp = step['freq']
                    elif 'lp' in step['process_type']:
                        lp = step['freq']

        info['highpass'] = hp
        info['lowpass'] = lp
        info['acq_pars'], info['acq_stim'] = None, None
        info['filename'] = None
        chs = []

        ch_names = [ch['name'] for ch in bti_info['chs']]
        info['ch_names'] = _rename_channels(ch_names)
        ch_mapping = zip(ch_names, info['ch_names'])
        logger.info('... Setting channel info structure.')
        for idx, (chan_4d, chan_vv) in enumerate(ch_mapping):
            chan_info = dict(zip(FIFF_INFO_CHS_FIELDS, FIFF_INFO_CHS_DEFAULTS))
            chan_info['ch_name'] = chan_vv
            chan_info['logno'] = idx + BTI.FIFF_LOGNO
            chan_info['scanno'] = idx + 1
            chan_info['cal'] = bti_info['chs'][idx]['scale']

            if any([chan_vv.startswith(k) for k in ('MEG', 'RFG', 'RFM')]):
                t, loc = bti_info['chs'][idx]['coil_trans'], None
                if t is not None:
                    t, loc = _convert_coil_trans(t.astype('>f8'), dev_ctf_t,
                                                 bti_to_nm)
                    if idx == 1:
                        logger.info('... putting coil transforms in Neuromag '
                                    'coordinates')
                chan_info['coil_trans'] = t
                if loc is not None:
                    chan_info['loc'] = loc.astype('>f4')

            if chan_vv.startswith('MEG'):
                chan_info['kind'] = FIFF.FIFFV_MEG_CH
                chan_info['coil_type'] = FIFF.FIFFV_COIL_MAGNES_MAG
                chan_info['coord_frame'] = FIFF.FIFFV_COORD_DEVICE
                chan_info['unit'] = FIFF.FIFF_UNIT_T

            elif chan_vv.startswith('RFM'):
                chan_info['kind'] = FIFF.FIFFV_REF_MEG_CH
                chan_info['coil_type'] = FIFF.FIFFV_COIL_MAGNES_R_MAG
                chan_info['coord_frame'] = FIFF.FIFFV_COORD_DEVICE
                chan_info['unit'] = FIFF.FIFF_UNIT_T

            elif chan_vv.startswith('RFG'):
                chan_info['kind'] = FIFF.FIFFV_REF_MEG_CH
                chan_info['coord_frame'] = FIFF.FIFFV_COORD_DEVICE
                chan_info['unit'] = FIFF.FIFF_UNIT_T_M
                if chan_4d in ('GxxA', 'GyyA'):
                    chan_info['coil_type'] = FIFF.FIFFV_COIL_MAGNES_R_GRAD_DIA
                elif chan_4d in ('GyxA', 'GzxA', 'GzyA'):
                    chan_info['coil_type'] = FIFF.FIFFV_COIL_MAGNES_R_GRAD_OFF

            elif chan_vv.startswith('EEG'):
                chan_info['kind'] = FIFF.FIFFV_EEG_CH
                chan_info['coil_type'] = FIFF.FIFFV_COIL_EEG
                chan_info['coord_frame'] = FIFF.FIFFV_COORD_HEAD
                chan_info['unit'] = FIFF.FIFF_UNIT_V

            elif chan_vv == 'STI 013':
                chan_info['kind'] = FIFF.FIFFV_RESP_CH
            elif chan_vv == 'STI 014':
                chan_info['kind'] = FIFF.FIFFV_STIM_CH
            elif chan_vv.startswith('EOG'):
                chan_info['kind'] = FIFF.FIFFV_EOG_CH
            elif chan_vv == 'ECG 001':
                chan_info['kind'] = FIFF.FIFFV_ECG_CH
            elif chan_vv.startswith('EXT'):
                chan_info['kind'] = FIFF.FIFFV_MISC_CH
            elif chan_vv.startswith('UTL'):
                chan_info['kind'] = FIFF.FIFFV_MISC_CH

            chs.append(chan_info)

        info['chs'] = chs

        logger.info('... Reading digitization points from %s' %
                    head_shape_fname)
        logger.info('... putting digitization points in Neuromag c'
                    'oordinates')
        info['dig'], ctf_head_t = _setup_head_shape(head_shape_fname, use_hpi)
        logger.info('... Computing new device to head transform.')
        dev_head_t = _convert_dev_head_t(dev_ctf_t, bti_to_nm,
                                         ctf_head_t)

        info['dev_head_t'] = dict()
        info['dev_head_t']['from'] = FIFF.FIFFV_COORD_DEVICE
        info['dev_head_t']['to'] = FIFF.FIFFV_COORD_HEAD
        info['dev_head_t']['trans'] = dev_head_t
        info['dev_ctf_t'] = dict()
        info['dev_ctf_t']['from'] = FIFF.FIFFV_MNE_COORD_CTF_DEVICE
        info['dev_ctf_t']['to'] = FIFF.FIFFV_COORD_HEAD
        info['dev_ctf_t']['trans'] = dev_ctf_t
        info['ctf_head_t'] = dict()
        info['ctf_head_t']['from'] = FIFF.FIFFV_MNE_COORD_CTF_HEAD
        info['ctf_head_t']['to'] = FIFF.FIFFV_COORD_HEAD
        info['ctf_head_t']['trans'] = ctf_head_t
        logger.info('Done.')

        if False:  # XXX : reminds us to support this as we go
            # include digital weights from reference channel
            comps = info['comps'] = list()
            weights = bti_info['weights']
            by_name = lambda x: x[1]
            chn = dict(ch_mapping)
            columns = [chn[k] for k in weights['dsp_ch_names']]
            rows = [chn[k] for k in weights['ch_names']]
            col_order, col_names = zip(*sorted(enumerate(columns),
                                               key=by_name))
            row_order, row_names = zip(*sorted(enumerate(rows), key=by_name))
            # for some reason the C code would invert the signs, so we follow.
            mat = -weights['dsp_wts'][row_order, :][:, col_order]
            comp_data = dict(data=mat,
                             col_names=col_names,
                             row_names=row_names,
                             nrow=mat.shape[0], ncol=mat.shape[1])
            comps += [dict(data=comp_data, ctfkind=101,
                           #  no idea how to calibrate, just ones.
                           rowcals=np.ones(mat.shape[0], dtype='>f4'),
                           colcals=np.ones(mat.shape[1], dtype='>f4'),
                           save_calibrated=0)]
        else:
            logger.warning('Warning. Currently direct inclusion of 4D weight t'
                           'ables is not supported. For critical use cases '
                           '\nplease take into account the MNE command '
                           '\'mne_create_comp_data\' to include weights as '
                           'printed out \nby the 4D \'print_table\' routine.')

        # check that the info is complete
        assert not set(RAW_INFO_FIELDS) - set(info.keys())

        # check nchan is correct
        assert len(info['ch_names']) == info['nchan']

        cals = np.zeros(info['nchan'])
        for k in range(info['nchan']):
            cals[k] = info['chs'][k]['range'] * info['chs'][k]['cal']

        self.verbose = verbose
        self.cals = cals
        self.rawdir = None
        self.proj = None
        self.comp = None
        self._filenames = list()
        self.preload = True
        self._projector_hashes = [None]
        self.info = info

        logger.info('Reading raw data from %s...' % pdf_fname)
        self._data = _read_data(bti_info)
        self.first_samp, self.last_samp = 0, self._data.shape[1] - 1
        self._raw_lengths = np.array([self._data.shape[1]])
        self._first_samps = np.array([0])
        self._last_samps = self._raw_lengths - 1
        self.rawdirs = [[]]

        assert len(self._data) == len(self.info['ch_names'])
        self._times = np.arange(self.first_samp,
                                self.last_samp + 1) / info['sfreq']
        self._projectors = [None]
        logger.info('    Range : %d ... %d =  %9.3f ... %9.3f secs' % (
                    self.first_samp, self.last_samp,
                    float(self.first_samp) / info['sfreq'],
                    float(self.last_samp) / info['sfreq']))

        logger.info('Ready.')


@verbose
def read_raw_bti(pdf_fname, config_fname='config',
                 head_shape_fname='hs_file', rotation_x=None,
                 translation=(0.0, 0.02, 0.11), ecg_ch='E31',
                 eog_ch=('E63', 'E64'), verbose=True):
    """ Raw object from 4D Neuroimaging MagnesWH3600 data

    Note.
    1) Currently direct inclusion of reference channel weights
    is not supported. Please use 'mne_create_comp_data' to include
    the weights or use the low level functions from this module to
    include them by yourself.
    2) The informed guess for the 4D name is E31 for the ECG channel and
    E63, E63 for the EOG channels. Pleas check and adjust if those channels
    are present in your dataset but 'ECG 01' and 'EOG 01', 'EOG 02' don't
    appear in the channel names of the raw object.

    Parameters
    ----------
    pdf_fname : str | None
        absolute path to the processed data file (PDF)
    config_fname : str | None
        absolute path to system confnig file. If None, it is assumed to be in
        the same directory.
    head_shape_fname : str
        absolute path to the head shape file. If None, it is assumed to be in
        the same directory.
    rotation_x : float | int | None
        Degrees to tilt x-axis for sensor frame misalignment.
        If None, no adjustment will be applied.
    translation : array-like
        The translation to place the origin of coordinate system
        to the center of the head.
    ecg_ch: str | None
      The 4D name of the ECG channel. If None, the channel will be treated
      as regular EEG channel.
    eog_ch: tuple of str | None
      The 4D names of the EOG channels. If None, the channels will be treated
      as regular EEG channels.
    verbose : bool, str, int, or None
        If not None, override default verbose level (see mne.verbose).
    """
    return RawBTi(pdf_fname, config_fname=config_fname,
                  head_shape_fname=head_shape_fname,
                  rotation_x=rotation_x, translation=translation,
                  verbose=verbose)