convert5.py

import math
import numpy
import mdp
import sys
import array
import matplotlib.pyplot as plt

def main():
    file_name =  sys.argv[1]
    fs = int(sys.argv[2])
    input_data = read_file(file_name)
    normalized_input_data = normalize_samples(input_data)
    #g_modes = ['pow3', 'tanh', 'gaus', 'skew']
    g_modes = ['gaus']
    for g in g_modes:
        channel2 = apply_ica(normalized_input_data, 1, g)
        x = calc_heart_rate(channel2, fs)
        #values = [i[1] for i in x]
        #plt.plot([i[0] for i in x], values)
        #plt.show()
        print x
        print "Heart beat: ", x*60

def read_file(file_name):
    
    '''
        Receive rgb data from a file in this format r,g,b, on each line        
    '''
    a = open(file_name).read().split('\n')
    def my_split(x):
        _float = float
        return map(_float, x.split(','))
    return map(my_split, a[:-1])

def normalize_sample(x):
    '''
        x is a 1-dimensional raw list of receive data
        Return a normalized list
    '''
    
    #calc mean value. This may be not correct
    N = len(x)
    mean = float(sum(x))/N
    x1 = [i**2 for i in x]
    mean_x1 = float(sum(x1))/N
    variance = mean_x1 - mean**2
    normalized_x = [float(i-mean)/variance for i in x]
    return normalized_x
    
def normalize_samples(x):
    '''
        Normalize all the 3 channels
        x is a 3-dimensional list of raw data. Just normalize each of them
    '''

    f_n = normalize_sample
    return map(f_n, x)


def apply_ica(x, i, g):
    '''
        Seperate the original rgb sources into 3 different channels
        i is between 0 to 2, which is the channel
        g is the mode to converge
    '''

    #x1 = numpy.transpose(x)
    x1 = numpy.array(x)
    #channels = mdp.fastica(x1, input_dim=3)
    #channels = mdp.fastica(x1, input_dim=3, verbose=True, g=g)
    channels = mdp.fastica(x1, input_dim=3, g=g, verbose=True)
    tmp = channels.transpose()
    return tmp[i]


def calc_heart_rate(x, fs):
    '''
        x is a discrete function, list of values. 
        Returns a DFT of these samples.
        fs is the sampling frequency        
        N is the total number of samples
    '''

    N = len(x)
    #use discrete frequency transform for the samplings
    dft_x = numpy.fft.fft(x)
    amp_max = 0
    tar_freq = 0
    def my_compare(x,y):
        return cmp(y[1], x[1])

    tmp_list = [(fs*i/float(N),math.sqrt(dft_x[i].real**2 + dft_x[i].imag**2)) for i in range(N) if fs*i/float(N)>0.75 and fs*i/float(N)<fs/2]
    tmp_result = []
    for i in tmp_list:
        
    #for i in range(1, N):
    #    tmp = dft_x[i]
    #    f = fs * i / float(N)
    #    if f<0.75 or f>fs/2: continue
    #    amp = math.sqrt(tmp.real**2 + tmp.imag**2)
    #    if amp_max < amp:
    #        amp_max = amp
    #        tar_freq = f
    #return tar_freq

if __name__ == '__main__':
    main()