code: from __future__ import print_function
  # Authors: Denis Engemann <denis.engemann@gmail.com>
  #
  # License: BSD (3-clause)
  
  import os
  import os.path as op
  from functools import reduce, partial
  import warnings
  
  import numpy as np
  from numpy.testing import (assert_array_almost_equal, assert_array_equal,
                             assert_allclose)
  from nose.tools import assert_true, assert_raises, assert_equal
  
  from mne.datasets import testing
  from mne.io import read_raw_fif, read_raw_bti
  from mne.io.bti.bti import (_read_config, _process_bti_headshape,
                              _read_bti_header, _get_bti_dev_t,
                              _correct_trans, _get_bti_info)
  from mne.io.tests.test_raw import _test_raw_reader
  from mne.tests.common import assert_dig_allclose
  from mne.io.pick import pick_info
  from mne.io.constants import FIFF
  from mne import pick_types
  from mne.utils import run_tests_if_main, slow_test
  from mne.transforms import Transform, combine_transforms, invert_transform
  from mne.externals import six
  
  warnings.simplefilter('always')
  
  base_dir = op.join(op.abspath(op.dirname(__file__)), 'data')
  
  archs = 'linux', 'solaris'
  pdf_fnames = [op.join(base_dir, 'test_pdf_%s' % a) for a in archs]
  config_fnames = [op.join(base_dir, 'test_config_%s' % a) for a in archs]
  hs_fnames = [op.join(base_dir, 'test_hs_%s' % a) for a in archs]
  exported_fnames = [op.join(base_dir, 'exported4D_%s_raw.fif' % a)
                     for a in archs]
  tmp_raw_fname = op.join(base_dir, 'tmp_raw.fif')
  
  fname_2500 = op.join(testing.data_path(download=False), 'BTi', 'erm_HFH',
                       'c,rfDC')
  
  # the 4D exporter doesn't export all channels, so we confine our comparison
  NCH = 248
  
  
  @testing.requires_testing_data
  def test_read_2500():
      """Test reading data from 2500 system."""
      _test_raw_reader(read_raw_bti, pdf_fname=fname_2500, head_shape_fname=None)
  
  
  def test_read_config():
      """Test read bti config file."""
      # for config in config_fname, config_solaris_fname:
      for config in config_fnames:
          cfg = _read_config(config)
          assert_true(all('unknown' not in block.lower() and block != ''
                          for block in cfg['user_blocks']))
  
  
  def test_crop_append():
      """Test crop and append raw."""
      raw = _test_raw_reader(
          read_raw_bti, pdf_fname=pdf_fnames[0],
          config_fname=config_fnames[0], head_shape_fname=hs_fnames[0])
      y, t = raw[:]
      t0, t1 = 0.25 * t[-1], 0.75 * t[-1]
      mask = (t0 <= t) * (t <= t1)
      raw_ = raw.copy().crop(t0, t1)
      y_, _ = raw_[:]
      assert_true(y_.shape[1] == mask.sum())
      assert_true(y_.shape[0] == y.shape[0])
  
  
  def test_transforms():
      """Test transformations."""
      bti_trans = (0.0, 0.02, 0.11)
      bti_dev_t = Transform('ctf_meg', 'meg', _get_bti_dev_t(0.0, bti_trans))
      for pdf, config, hs, in zip(pdf_fnames, config_fnames, hs_fnames):
          raw = read_raw_bti(pdf, config, hs, preload=False)
          dev_ctf_t = raw.info['dev_ctf_t']
          dev_head_t_old = raw.info['dev_head_t']
          ctf_head_t = raw.info['ctf_head_t']
  
          # 1) get BTI->Neuromag
          bti_dev_t = Transform('ctf_meg', 'meg', _get_bti_dev_t(0.0, bti_trans))
  
          # 2) get Neuromag->BTI head
          t = combine_transforms(invert_transform(bti_dev_t), dev_ctf_t,
                                 'meg', 'ctf_head')
          # 3) get Neuromag->head
          dev_head_t_new = combine_transforms(t, ctf_head_t, 'meg', 'head')
  
          assert_array_equal(dev_head_t_new['trans'], dev_head_t_old['trans'])
  
  
  @slow_test
  def test_raw():
      """Test bti conversion to Raw object."""
      for pdf, config, hs, exported in zip(pdf_fnames, config_fnames, hs_fnames,
                                           exported_fnames):
          # rx = 2 if 'linux' in pdf else 0
          assert_raises(ValueError, read_raw_bti, pdf, 'eggs', preload=False)
          assert_raises(ValueError, read_raw_bti, pdf, config, 'spam',
                        preload=False)
          if op.exists(tmp_raw_fname):
              os.remove(tmp_raw_fname)
          ex = read_raw_fif(exported, preload=True)
          ra = read_raw_bti(pdf, config, hs, preload=False)
          assert_true('RawBTi' in repr(ra))
          assert_equal(ex.ch_names[:NCH], ra.ch_names[:NCH])
          assert_array_almost_equal(ex.info['dev_head_t']['trans'],
                                    ra.info['dev_head_t']['trans'], 7)
          assert_dig_allclose(ex.info, ra.info)
          coil1, coil2 = [np.concatenate([d['loc'].flatten()
                          for d in r_.info['chs'][:NCH]])
                          for r_ in (ra, ex)]
          assert_array_almost_equal(coil1, coil2, 7)
  
          loc1, loc2 = [np.concatenate([d['loc'].flatten()
                        for d in r_.info['chs'][:NCH]])
                        for r_ in (ra, ex)]
          assert_allclose(loc1, loc2)
  
          assert_allclose(ra[:NCH][0], ex[:NCH][0])
          assert_array_equal([c['range'] for c in ra.info['chs'][:NCH]],
                             [c['range'] for c in ex.info['chs'][:NCH]])
          assert_array_equal([c['cal'] for c in ra.info['chs'][:NCH]],
                             [c['cal'] for c in ex.info['chs'][:NCH]])
          assert_array_equal(ra._cals[:NCH], ex._cals[:NCH])
  
          # check our transforms
          for key in ('dev_head_t', 'dev_ctf_t', 'ctf_head_t'):
              if ex.info[key] is None:
                  pass
              else:
                  assert_true(ra.info[key] is not None)
                  for ent in ('to', 'from', 'trans'):
                      assert_allclose(ex.info[key][ent],
                                      ra.info[key][ent])
  
          ra.save(tmp_raw_fname)
          re = read_raw_fif(tmp_raw_fname)
          print(re)
          for key in ('dev_head_t', 'dev_ctf_t', 'ctf_head_t'):
              assert_true(isinstance(re.info[key], dict))
              this_t = re.info[key]['trans']
              assert_equal(this_t.shape, (4, 4))
              # cehck that matrix by is not identity
              assert_true(not np.allclose(this_t, np.eye(4)))
          os.remove(tmp_raw_fname)
  
  
  def test_info_no_rename_no_reorder_no_pdf():
      """Test private renaming, reordering and partial construction option."""
      for pdf, config, hs in zip(pdf_fnames, config_fnames, hs_fnames):
          info, bti_info = _get_bti_info(
              pdf_fname=pdf, config_fname=config, head_shape_fname=hs,
              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)
          info2, bti_info = _get_bti_info(
              pdf_fname=None, config_fname=config, head_shape_fname=hs,
              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)
  
          assert_equal(info['ch_names'],
                       [ch['ch_name'] for ch in info['chs']])
          assert_equal([n for n in info['ch_names'] if n.startswith('A')][:5],
                       ['A22', 'A2', 'A104', 'A241', 'A138'])
          assert_equal([n for n in info['ch_names'] if n.startswith('A')][-5:],
                       ['A133', 'A158', 'A44', 'A134', 'A216'])
  
          info = pick_info(info, pick_types(info, meg=True, stim=True,
                                            resp=True))
          info2 = pick_info(info2, pick_types(info2, meg=True, stim=True,
                                              resp=True))
  
          assert_true(info['sfreq'] is not None)
          assert_true(info['lowpass'] is not None)
          assert_true(info['highpass'] is not None)
          assert_true(info['meas_date'] is not None)
  
          assert_equal(info2['sfreq'], None)
          assert_equal(info2['lowpass'], None)
          assert_equal(info2['highpass'], None)
          assert_equal(info2['meas_date'], None)
  
          assert_equal(info['ch_names'], info2['ch_names'])
          assert_equal(info['ch_names'], info2['ch_names'])
          for key in ['dev_ctf_t', 'dev_head_t', 'ctf_head_t']:
              assert_array_equal(info[key]['trans'], info2[key]['trans'])
  
          assert_array_equal(
              np.array([ch['loc'] for ch in info['chs']]),
              np.array([ch['loc'] for ch in info2['chs']]))
  
      # just check reading data | corner case
      raw1 = read_raw_bti(
          pdf_fname=pdf, config_fname=config, head_shape_fname=None,
          sort_by_ch_name=False, preload=True)
      # just check reading data | corner case
      raw2 = read_raw_bti(
          pdf_fname=pdf, config_fname=config, head_shape_fname=None,
          rename_channels=False,
          sort_by_ch_name=True, preload=True)
  
      sort_idx = [raw1.bti_ch_labels.index(ch) for ch in raw2.bti_ch_labels]
      raw1._data = raw1._data[sort_idx]
      assert_array_equal(raw1._data, raw2._data)
      assert_array_equal(raw2.bti_ch_labels, raw2.ch_names)
  
  
  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):
          raw_info, _ = get_info(pdf, config, hs, convert=False)
          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_bytes_io():
      """Test bti bytes-io API."""
      for pdf, config, hs in zip(pdf_fnames, config_fnames, hs_fnames):
          raw = read_raw_bti(pdf, config, hs, convert=True, preload=False)
  
          with open(pdf, 'rb') as fid:
              pdf = six.BytesIO(fid.read())
          with open(config, 'rb') as fid:
              config = six.BytesIO(fid.read())
          with open(hs, 'rb') as fid:
              hs = six.BytesIO(fid.read())
  
          raw2 = read_raw_bti(pdf, config, hs, convert=True, preload=False)
          repr(raw2)
          assert_array_equal(raw[:][0], raw2[:][0])
  
  
  def test_setup_headshape():
      """Test reading bti headshape."""
      for hs in hs_fnames:
          dig, t = _process_bti_headshape(hs)
          expected = set(['kind', 'ident', 'r'])
          found = set(reduce(lambda x, y: list(x) + list(y),
                             [d.keys() for d in dig]))
          assert_true(not expected - found)
  
  run_tests_if_main()
   # -*- coding: utf-8 -*-
  """
  Created on Tue May 30 13:06:52 2017
  
  @author: Mostafa Hammoud
  This script contains various functions that will be used
  """
  
  import numpy as np
  import operator
  from os import listdir
  from os.path import isfile, join
  from matplotlib import pyplot as plt
  """
  This function takes the R waves that were detected and makes sure they are good 
  R waves
  """
  def fix_r(peaks,long_filt_ecg):
      # Loop through the peaks 
      for i in range (0,len(peaks)):
          # Select a range close to the peaks to check the maximum
          temp_range = np.arange(peaks[i] - 10,peaks[i] + 10)
          
          # If we have a negative R wave we don't want to confuze it with a 
          # significant S wave
          if long_filt_ecg[peaks[i]] < 0:
              temp_range = np.arange(peaks[i] - 15,peaks[i] + 10)
          # Make sure the range does not go out of bound
          if temp_range[-1] >= len(long_filt_ecg):
              temp_range = np.arange(peaks[i] - 6,len(long_filt_ecg) - 1)
          # Find the maximum
          max_index, max_value = max(enumerate(long_filt_ecg[temp_range]), 
                                     key=operator.itemgetter(1))
          # Set the peak to the maximum that was found
          peaks[i] = max_index + temp_range[0]
          
      return peaks
      
  """
  This function takes an uper limit and a lower limit and removes all values that
  are not in that range
  upper: Upper limit
  lower: Lower limit
  equal: Boolean to select if the limits are included. Default is True
  
  """    
  def remove_unwanted(array, upper, lower, equal = True):
      
      out = []
      
      for i in array:
          # i can be []
          if np.size(i) > 0:
              if equal is True:
                  if (i >= lower and i <= upper):
                      out.append(i)
              else:
                  if (i > lower and i < upper):
                      out.append(i)            
      
      # Select the peaks that are inside the range start to end
      array = out
      
      return array
  """
  This function returns the an array with the path to all three ECGs
  """
  def get_ecg_path(file_name):
      
      ecg_path = [f for f in listdir(file_name) if isfile(join(file_name,f))]
      
      j = 0
      for i in ecg_path:
          
          if 'ECG_I.csv' in i or 'ECG_II.csv' in i or 'ECG_III.csv' in i:
              ecg_path[j] = file_name + '\\' + i
              j += 1
          
      return ecg_path
  """
  This function calculates the derivative
  """
  def derivative(y, fs):
      
      if len(y) < 2:
          return 0
          
      dy = np.zeros(len(y),np.float)
      dx = 1/fs
      
      dy[0:-1] =  np.diff(y)/dx
        
      dy[-1] = (y[-1] - y[-2])/dx
      
      return dy
  
  """
  This function plots the points on the figure
  """
  def plot_points(long_ecg, index, long_time, start, end):
      
      r = remove_unwanted(index.r, end, start)  
      
      if (len(r)):
          
          p = remove_unwanted(index.p, end, start)  
          
          q = remove_unwanted(index.q, end, start)  
  
          s = remove_unwanted(index.s, end, start)  
  
          t = remove_unwanted(index.t, end, start)  
          
          pon  = remove_unwanted(index.pon , end, start) 
          poff = remove_unwanted(index.poff, end, start)  
          qon  = remove_unwanted(index.qon , end, start)  
          soff = remove_unwanted(index.soff, end, start)  
          ton  = remove_unwanted(index.ton , end, start)  
          toff = remove_unwanted(index.toff, end, start)  
          
          plt.scatter(long_time[r],long_ecg[r], color = 'red')
          plt.scatter(long_time[p],long_ecg[p], color = 'blue')
          plt.scatter(long_time[q],long_ecg[q], color = 'green')
          plt.scatter(long_time[s],long_ecg[s], color = 'magenta')
          plt.scatter(long_time[t],long_ecg[t], color = 'black')
          
          #plt.scatter(long_time[pon],long_ecg[pon])
          #plt.scatter(long_time[poff],long_ecg[poff])
          #plt.scatter(long_time[qon],long_ecg[qon])
          #plt.scatter(long_time[soff],long_ecg[soff])
          #plt.scatter(long_time[ton],long_ecg[ton])
          #plt.scatter(long_time[toff],long_ecg[toff])
  
  """
  
  This function loops thorugh the index from start to end and fixes the location 
  of the peaks or troughs. Since sometimes we might have a small shift in their 
  location
  
  If the start and end were not indicated then we loop through the whole index
  
  If 
  
  """
  def fix_peaks(long_ecg, wave_index, start_index, end_index, peak = False, 
                trough = False):
      # If we are dealing with a peak, then we keep the ecg as it is
      if peak:
          ecg = long_ecg.copy()
      # If we are dealng with a trough, then we inverse the ecg so the min is max
      elif trough:
          ecg = -long_ecg.copy()  
      else:
          print('A type must be selected')
          return
  
      for i in range(start_index, end_index):
          if np.size(wave_index[i]) > 0:
              temp_range = np.arange(wave_index[i] - 2,wave_index[i] + 3)
          
              if temp_range[-1] >= len(ecg):
                  temp_range = np.arange(wave_index[i] - 2,len(ecg) - 1)
          
              # Find the maximum
              max_index, max_value = max(enumerate(ecg[temp_range]), 
                                         key=operator.itemgetter(1))
          
              # Set the peak to the maximum that was found
              wave_index[i] = max_index + temp_range[0]
  
      return wave_index
apis: numpy.eye numpy.allclose numpy.testing.assert_array_equal numpy.testing.assert_array_almost_equal numpy.testing.assert_allclose numpy.array
apis: numpy.arange numpy.size matplotlib.pyplot.scatter numpy.diff