Compare output to p_pod and opt3d

lasthdr.txt

@@ -0,0 +1 @@
+D:\data\anc\1anc.hdr

mne/io/boxy/boxy.py

@@ -539,8 +539,14 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult):
 
                                 # Save our data based on data type.
                                 data_[index_loc, :] = boxy_array[:, channel]
+
+                    ###need to fix the first few bad points in the recording###
+                    n_bad_points = 12
+                    for i_point in range(n_bad_points):
+                        data_[:,i_point] = data_[:,n_bad_points] 
 
                     # Phase unwrapping.
+
                     if i_data == 'Ph':
                         print('Fixing phase wrap')
                         # Accounts for sharp, sudden changes in phase
@@ -558,7 +564,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult):
 
                         print('Detrending phase data')
                         # Remove trends and drifts that occur over time.
-                        y = np.linspace(0, np.size(data_, axis=1) - 1,
+                        y = np.linspace(1, np.size(data_, axis=1),
                                         np.size(data_, axis=1))
                         x = np.transpose(y)
                         for i_chan in range(np.size(data_, axis=0)):
@@ -580,7 +586,7 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult):
                         ph_out_thr = 3
 
                         # Set ddof to 1 to mimic matlab.
-                        sdph = np.std(data_, 1, ddof=1)
+                        sdph = np.std(data_[:, n_bad_points - 1:], 1, ddof=1)
                         n_ph_out = np.zeros(np.size(data_, axis=0),
                                             dtype=np.int8)
 
@@ -591,18 +597,27 @@ def _read_segment_file(self, data, idx, fi, start, stop, cals, mult):
                             if len(outliers) > 0:
                                 if outliers[0] == 0:
                                     outliers = outliers[1:]
-                                if len(outliers) > 0:
-                                    if (outliers[-1] == np.size(data_,
-                                                                axis=1) - 1):
-                                        outliers = outliers[:-1]
-                                    n_ph_out[i_chan] = int(len(outliers))
-                                    for i_pt in range(n_ph_out[i_chan]):
-                                        j_pt = outliers[i_pt]
-                                        data_[i_chan, j_pt] = (
-                                            (data_[i_chan, j_pt - 1] +
-                                             data_[i_chan, j_pt + 1]) / 2)
-
-                        # Convert phase to pico seconds.
+                                if (outliers[-1] == np.size(data_,
+                                                            axis=1) - 1):
+                                    outliers = outliers[:-1]
+                                n_ph_out[i_chan] = int(len(outliers))
+                                for i_pt in range(n_ph_out[i_chan]):
+                                    j_pt = outliers[i_pt]
+                                    data_[i_chan, j_pt] = (
+                                        (data_[i_chan, j_pt - 1] +
+                                         data_[i_chan, j_pt + 1]) / 2)
+
+                    # now let's normalise our data #
+                    for i_chan in range(len(data_)):
+                        if i_data == 'Ph':
+                            data_[i_chan,:] = (data_[i_chan,:] - 
+                                               np.mean(data_[i_chan,:]))
+                        else:
+                            data_[i_chan,:] = (data_[i_chan,:] /
+                                               np.mean(data_[i_chan,:]) - 1)
+
+                    # Convert phase to pico seconds.
+                    if i_data == 'Ph':
                         for i_chan in range(np.size(data_, axis=0)):
                             data_[i_chan, :] = ((1e12 * data_[i_chan, :]) /
                                                 (360 * mtg_mdf[i_mtg][i_chan]))