In [4]:

    

%matplotlib inline
import pandas as pd
import numpy as np
import scipy as sp 
import os, sys, glob
import matplotlib as mpl
from matplotlib import pylab as plt
import wave
import pyaudio  
from datetime import timedelta, datetime


    

In [5]:

    

def timecodematch(x, base_date):
    delta = []
    b_d = datetime.strptime(base_date, "%H:%M:%S:%f")
    for y in range(len(x)):
        d = datetime.strptime(x[y], "%H:%M:%S:%f")
        delta.append(abs(d - b_d))
    return delta

def getSec(minutes,sec):
    print str(minutes*60 + sec) + ' sec'

def getMin(sec):
    print str(int(np.floor(sec/60))) + ' min ' + str(np.mod(sec,60)) + ' sec'

def playMusic(s1,s2):
    f = wave.open('music/Science_of_Disco_Set_16bit.wav')
    chunk = 1024
    p = pyaudio.PyAudio()  
    stream = p.open(format = p.get_format_from_width(f.getsampwidth()),  
                    channels = f.getnchannels(),  
                    rate = f.getframerate(),  
                    output = True)  
    count = 0
    d = f.readframes(chunk)
    while d != '':
        d = f.readframes(chunk)
        count += 1  
        if count >= s1*(f.getframerate()/chunk) and not count > s2*(f.getframerate()/chunk):
            stream.write(d)

    stream.stop_stream()  
    stream.close()  
    p.terminate()


    

In [6]:

    

# Read in data for each subject between start and end times

sampling_freq = 1000 # in milliseconds
beginTS = '20:29:33:113' # from 'data/76.csv'
endTS   = '21:00:20:802' # from 'data/76.csv'
da = datetime.strptime(endTS, "%H:%M:%S:%f") - datetime.strptime(beginTS, "%H:%M:%S:%f")
totMS = ((da.seconds * 1000000) + da.microseconds) / 1000

gforce_all = []
sub = []

count = 1
countall = 1

plotting = False
zeroGravity = False 

for file in glob.glob("data/*.csv"):
    
    s = os.path.splitext(os.path.split(file)[1])[0]
                  
    try:        

        df = pd.read_csv(file, delimiter=';')        

        print '                      subject ' + s + ', start time: ' + df['time'].get_values()[0] + ', end time: ' + df['time'].get_values()[-1]

        startTime = np.argmin(timecodematch(df['time'].get_values(), beginTS))
        endTime = np.argmin(timecodematch(df['time'].get_values(), endTS))

        # instead, try interpolating to specific time intervals:
        time_s = df['time'].get_values()[startTime:endTime]
        gforce_s = df['gforce'].get_values()[startTime:endTime]

        if plotting:
            fig, (ax1, ax2) = plt.subplots(2)
            ax1.plot(gforce_s)

        ts = timecodematch(time_s, beginTS)
        for t in range(len(time_s)):
            time_s[t] = ((ts[t].seconds * 1000000) + ts[t].microseconds) / 1000
        
        # !!! check interpolation method. Probably better ways to do this:
        gforce_interp = np.interp(range(0,totMS,sampling_freq), 
                                    np.array(time_s, dtype=np.float64), 
                                    np.array(gforce_s, dtype=np.float64))

        if zeroGravity: # take absolute value to account for decrease in gforce
            gforce_interp = abs(gforce_interp - 1)                

        if plotting:
            ax2.plot(gforce_interp)

        if s not in ['15','35','37','40','42','47','48','66','67','76']: # removed from visual inspection of plots & '76' 

            gforce_all.append(gforce_interp)
            sub.append(s)
            print "good: " + str(count) + " of total: " + str(countall)
            count += 1
        
        else:
            print "removed: subject " + str(s)
        countall += 1

    except:
        print "bad: subject " + str(s)
        countall += 1
        pass

gforce_mean = np.mean(gforce_all, axis=0)

print "percentage good: " + str(count/float(countall))


    

    




                    subject 10, start time: 19:57:58:3, end time: 20:13:00:6
bad: subject 10
                    subject 11, start time: 20:14:32:301, end time: 21:25:59:245
good: 1 of total: 2
                    subject 13, start time: 20:29:47:906, end time: 21:35:43:960
good: 2 of total: 3
                    subject 14, start time: 20:28:38:990, end time: 21:08:15:422
good: 3 of total: 4
                    subject 15, start time: 20:46:01:9, end time: 21:03:26:6
removed: subject 15
                    subject 16, start time: 19:55:57:462, end time: 21:10:16:567
good: 4 of total: 6
                    subject 17, start time: 19:36:37:265, end time: 21:23:29:78
good: 5 of total: 7
                    subject 18, start time: 19:39:36:7, end time: 19:54:56:4
bad: subject 18
                    subject 19, start time: 19:45:43:8, end time: 20:00:44:7
bad: subject 19
                    subject 2, start time: 19:58:57:323, end time: 21:09:13:424
good: 6 of total: 10
                    subject 20, start time: 19:32:28:421, end time: 21:22:20:477
good: 7 of total: 11
                    subject 21, start time: 20:03:16:4, end time: 20:18:16:6
bad: subject 21
                    subject 22, start time: 20:06:08:9, end time: 20:21:08:9
bad: subject 22
                    subject 23, start time: 19:46:32:7, end time: 20:01:38:2
bad: subject 23
                    subject 24, start time: 19:44:50:331, end time: 21:16:29:631
good: 8 of total: 15
                    subject 26, start time: 19:40:20:2, end time: 19:55:31:2
bad: subject 26
                    subject 27, start time: 19:34:11:6, end time: 19:49:15:6
bad: subject 27
                    subject 28, start time: 19:48:08:601, end time: 02:49:52:90
good: 9 of total: 18
                    subject 29, start time: 19:54:54:5, end time: 20:09:59:7
bad: subject 29
                    subject 3, start time: 19:43:37:613, end time: 21:24:11:467
good: 10 of total: 20
                    subject 30, start time: 19:53:22:979, end time: 21:23:53:360
good: 11 of total: 21
                    subject 31, start time: 20:10:33:75, end time: 02:05:23:951
good: 12 of total: 22
                    subject 32, start time: 20:00:09:476, end time: 21:09:46:592
good: 13 of total: 23
                    subject 33, start time: 19:41:30:682, end time: 21:27:49:991
good: 14 of total: 24
                    subject 34, start time: 19:29:20:6, end time: 19:44:25:0
bad: subject 34
                    subject 35, start time: 20:18:29:37, end time: 20:53:33:629
removed: subject 35
                    subject 36, start time: 19:26:54:662, end time: 19:44:46:626
bad: subject 36
                    subject 37, start time: 20:32:13:4, end time: 20:47:19:2
removed: subject 37
                    subject 38, start time: 20:11:31:62, end time: 21:33:19:769
good: 15 of total: 29
                    subject 4, start time: 20:21:48:746, end time: 21:33:44:779
good: 16 of total: 30
                    subject 40, start time: 20:46:06:89, end time: 21:14:00:432
removed: subject 40
                    subject 41, start time: 19:28:06:933, end time: 01:37:06:561
good: 17 of total: 32
                    subject 42, start time: 20:35:21:5, end time: 20:50:21:1
removed: subject 42
                    subject 43, start time: 20:09:17:1, end time: 20:24:28:9
bad: subject 43
                    subject 44, start time: 20:13:13:839, end time: 21:32:34:870
good: 18 of total: 35
                    subject 45, start time: 19:18:51:8, end time: 19:36:30:4
bad: subject 45
                    subject 46, start time: 20:07:14:6, end time: 20:24:58:2
bad: subject 46
                    subject 47, start time: 20:42:42:371, end time: 02:49:29:54
removed: subject 47
                    subject 48, start time: 20:25:15:3, end time: 20:40:12:9
removed: subject 48
                    subject 49, start time: 19:14:29:634, end time: 21:18:41:511
good: 19 of total: 40
                    subject 5, start time: 20:15:33:620, end time: 23:53:53:958
good: 20 of total: 41
                    subject 50, start time: 20:08:16:689, end time: 21:39:40:349
good: 21 of total: 42
                    subject 6, start time: 20:02:08:7, end time: 20:17:11:2
bad: subject 6
                    subject 61, start time: 20:12:42:98, end time: 21:11:22:120
good: 22 of total: 44
                    subject 63, start time: 20:14:22:691, end time: 21:08:50:824
good: 23 of total: 45
                    subject 64, start time: 20:16:30:911, end time: 21:27:06:497
good: 24 of total: 46
                    subject 65, start time: 20:20:18:542, end time: 21:28:36:300
good: 25 of total: 47
                    subject 66, start time: 20:21:29:9, end time: 20:36:34:5
removed: subject 66
                    subject 67, start time: 20:49:08:723, end time: 21:05:32:96
removed: subject 67
                    subject 7, start time: 20:23:42:480, end time: 21:28:55:801
good: 26 of total: 50
                    subject 71, start time: 19:56:49:4, end time: 20:11:56
bad: subject 71
                    subject 72, start time: 19:57:56:8, end time: 20:13:02:4
bad: subject 72
                    subject 73, start time: 19:59:07:769, end time: 21:18:15:688
good: 27 of total: 53
                    subject 74, start time: 20:00:12:476, end time: 21:42:01:407
good: 28 of total: 54
                    subject 75, start time: 20:02:11:834, end time: 21:31:49:125
good: 29 of total: 55
                    subject 76, start time: 20:25:59:175, end time: 21:04:58:26
removed: subject 76
                    subject 77, start time: 20:19:01:12, end time: 21:28:11:175
good: 30 of total: 57
                    subject 78, start time: 20:04:41:514, end time: 21:32:11:781
good: 31 of total: 58
                    subject 79, start time: 20:23:16:86, end time: 21:34:22:323
good: 32 of total: 59
                    subject 8, start time: 20:32:52:681, end time: 21:40:28:865
good: 33 of total: 60
                    subject 80, start time: 20:10:20:6, end time: 20:25:24:4
bad: subject 80
                    subject 81, start time: 19:36:42:933, end time: 21:23:05:514
good: 34 of total: 62
                    subject 82, start time: 19:38:20:842, end time: 21:24:43:314
good: 35 of total: 63
                    subject 83, start time: 19:41:52:6, end time: 19:57:00:4
bad: subject 83
                    subject 84, start time: 19:43:27:335, end time: 21:36:40:507
good: 36 of total: 65
                    subject 85, start time: 19:44:58:4, end time: 20:00:12:2
bad: subject 85
                    subject 86, start time: 19:46:43:610, end time: 21:20:20:202
good: 37 of total: 67
                    subject 87, start time: 19:48:45:825, end time: 21:13:10:635
good: 38 of total: 68
                    subject 88, start time: 19:52:20:86, end time: 21:32:56:72
good: 39 of total: 69
                    subject 89, start time: 19:54:11:782, end time: 21:15:39:871
good: 40 of total: 70
                    subject 9, start time: 20:17:39:998, end time: 21:17:17:455
good: 41 of total: 71
                    subject 90, start time: 19:55:33:807, end time: 21:37:02:652
good: 42 of total: 72
                    subject 91, start time: 19:14:50:5, end time: 19:30:00:6
bad: subject 91
                    subject 92, start time: 19:23:51:31, end time: 21:22:03:765
good: 43 of total: 74
                    subject 93, start time: 19:27:19:6, end time: 19:42:28:1
bad: subject 93
                    subject 94, start time: 19:29:28:768, end time: 21:14:52:340
good: 44 of total: 76
                    subject 95, start time: 19:30:57:580, end time: 21:19:22:762
good: 45 of total: 77
                    subject 96, start time: 19:33:05:166, end time: 21:30:02:476
good: 46 of total: 78
                    subject 97, start time: 19:34:57:371, end time: 21:38:12:546
good: 47 of total: 79
percentage good: 0.6

In [17]:

    

# Smoothing:
# !!! might be better way to smooth data. Currently using Gaussian kernel with:

smoothingSize = 20 # in seconds

from scipy.ndimage.filters import gaussian_filter

# normalize data
gforce_all_normed = []
for i in range(len(gforce_all)):
    gforce_all_normed.append((gforce_all[i] - np.mean(gforce_all[i])) / np.std(gforce_all[i]))

# smooth data
gforce_all_smoothed = []
for i in range(len(gforce_all)):
    gforce_all_smoothed.append(gaussian_filter(gforce_all_normed[i], sigma=smoothingSize))

# get mean gforce across individuals:
gforce_mean = np.mean(gforce_all_smoothed, axis=0)


    

In [9]:

    

# Visually check timeseries of each 'good' participant:

for i in range(len(gforce_all_smoothed)):
    fig, ax1 = plt.subplots(1,figsize=(15,2))
    ax1.plot(gforce_all_smoothed[i])
    plt.title(str(sub[i]))


    

    




/Applications/miniconda3/envs/topography/lib/python2.7/site-packages/matplotlib/pyplot.py:516: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    













    

In [18]:

    

# Plotting group:

fig, ax1 = plt.subplots(1,figsize=(20,7))
for i in range(len(gforce_all_smoothed)):
    ax1.plot(gforce_all_smoothed[i])

ax1.plot(gforce_mean, 'k',linewidth=5,  linestyle='--')
ax1.plot(np.linspace(0,0,len(gforce_mean)), 'r', linewidth=2, linestyle='--')
plt.xticks(range(0,len(gforce_mean), 50), np.array(range(0,totMS,sampling_freq*50))/1000, rotation='0')
ax1.set_xlabel('seconds')
ax1.set_ylabel('gforce')
ax1.set_xlim(0,1870)
ax1.grid(True)
plt.show()


    

In [12]:

    

# functions for converting between sec and [min sec]
getMin(930)
#getSec(19,49)


    

    




15 min 30 sec

In [31]:

    

# function for playing music between start and end seconds
playMusic(1200,1210)


    

In [37]:

    

#import IPython
#IPython.display.Audio('music/Science_of_Disco_Set_16bit.wav')


    

In [32]:

    

# For phase clustering, instead try: https://github.com/alexminnaar/time-series-classification-and-clustering


    

In [33]:

    

# !!! the following isn't correct

import nitime.timeseries as ts
from nitime.analysis.coherence import CoherenceAnalyzer

np.set_printoptions(precision=4)  # for doctesting
t1 = ts.TimeSeries(data = np.array(gforce_all, dtype=np.float32), sampling_rate=1, time_unit='s')
c1 = CoherenceAnalyzer(t1)

print np.shape(np.array(gforce_all, dtype=np.float32))
print np.shape(c1.coherence)
print np.shape(c1.phase)


    

    




(47, 1848)
(47, 47, 33)
(47, 47, 33)

In [ ]: