MEG MEG and EEG files in single fif file

merge_eeg_meg.py

import pyimpress
import mne
import numpy as np
from mne import create_info
from mne.io import RawArray


mff_fname='/cluster/transcend/sheraz/data/MEG_EEG/EEG/taskforce_1_vis_20160310_063855.mff'

# mff is the class supporting mff io
eeg = pyimpress.mff(mff_fname)
eeg_data = eeg.get_data()
eeg_events = eeg.get_events()

#usual mne python stuff
fiff_fname='/cluster/transcend/sheraz/data/MEG_EEG/taskforce_1_vis.fif'
meg = mne.io.read_raw_fif(fiff_fname,preload=True)
meg_events = mne.find_events(meg, stim_channel='STI101',
                         mask=15, mask_type='and',min_duration=0.002)

# Time vector not necessary just checking
#ts = (np.arange(0,  eeg.hdr['nSamples']) - eeg.hdr['nSamplesPre']) * (1000 / eeg.hdr['Fs'])


# making sure eeg and meg length are equal, zero padding or trucating eeg if needed
eegevents=[item['start_sample'] for item in eeg_events]
megevents=meg_events[:,0]
megevents=megevents-meg.first_samp
last_meg_sample=meg.last_samp - meg.first_samp +1
last_eeg_sample=eeg_data.shape[1]
dist_last_meg_sample=last_meg_sample-megevents[-1]
dist_last_eeg_sample=last_eeg_sample-eegevents[-1]
if megevents[0] > eegevents[0]:
    zero_padding=np.zeros((eeg_data.shape[0],megevents[0]-eegevents[0]))
    eeg_data=np.hstack((zero_padding,eeg_data))
else:
    eeg_data=eeg_data[:,eegevents[0]-megevents[0]:]

if dist_last_meg_sample > dist_last_eeg_sample:
    zero_padding=np.zeros((eeg_data.shape[0],last_meg_sample - last_eeg_sample))
    eeg_data=np.hstack((eeg_data,zero_padding))
else:
    eeg_data=eeg_data[:,:-(dist_last_eeg_sample - dist_last_meg_sample+1)]


#creating dummy header
info = create_info(ch_names=eeg.hdr['label'][:-2],sfreq=eeg.hdr['Fs'],ch_types='eeg')
eeg_data = eeg_data[:-2,:]*1e-3# converting to volts

#writing to mne python raw array
raw_eeg = RawArray(eeg_data,info)

# adding eeg channels to meg
fiff_meg_eeg_fname='/cluster/transcend/sheraz/data/MEG_EEG/taskforce_1_meg_eeg_vis.fif'
raw_meg_eeg = meg.copy().add_channels([raw_eeg],force_update_info=True)
raw_meg_eeg.save(fiff_meg_eeg_fname,overwrite=True)

#checking evoked
meg_eeg = mne.io.read_raw_fif(fiff_meg_eeg_fname,preload=True)
event_id = dict(Q1=1, Q2=2, Q3=3, Q4=4)
tmin, tmax = -0.2, 0.6
reject = dict(grad=4000e-13, mag=4e-12,eeg=500e-3)
epochs = mne.Epochs(meg_eeg, meg_events, event_id, tmin, tmax,
                    baseline=(None, 0), reject=reject, proj=False)

evoked_Q1 = epochs['Q1'].average()
picks = mne.pick_types(meg_eeg.info, meg=False, eeg=True, eog=False, stim=False,
                       exclude='bads')
evoked_Q1.plot(picks=np.arange(309,320))