<reponame>Kotzly/EMG_AG<gh_stars>0
# -*- coding: utf-8 -*-
"""
Created on Fri Aug 24 09:40:28 2018

@author: <NAME>
"""

import numpy as np
#from numpy import fft
import matplotlib.pyplot as plt
#import scipy.signal as sig
import os
import random 
import emgReaderClass_v2 as erc
import threading
import multiprocessing
#import dataPlotter
import snkbrain
import scipy.stats as st
import cProfile
import re


                    # This script is compatible with 'emgReaderClass_v2', that
                    # reads the .csv files generated by 'movementSaver.m', from
                    # the folder './csv/'
profiling=False
    
bias=0              # If bias = 1, every cromossome will have a non frequency dependant DNA
maxGen=2000      # The max number of generations
startOver=True# If True, the code will not consider the last simulation
tamPop=30      # Population number
maxFreq=240        # This is the max Frequency to consider #240
freqStep=3         # For freqStep=3 -> The code will consider [1,2,3],[3,4,5], etc# 3

sourceType='ninapro'
ninaprofolders=['csv1','csv2','csv3','csv6','csv7','csv8']
fs=2000
##############################################################################
guid=0                          # Individual ID (logging variable)
log_val,log_train,log_test=[],[],[]
param=[]
real=[]                         # DATA
origin=[]                       # DATA
fv=[]                           # DATA
frv=[]                          # DATA

max_freq=240
min_freq=6


curtose = st.kurtosis
obliquidade = st.skew
variancia = np.var
media = np.mean
desvio_padrao = np.std

def rms(v):
    temp=sum([pow(i,2) for i in v])
    return np.sqrt(temp/len(v))
def centroide_espectral(v):
    temp=v
    n = sum([i*value for i,value in zip(range(1,1+len(temp)),temp)])
    mass_center = float(n)/sum(temp)
    return mass_center
def feature_scaling(v):
    mean = np.mean(v)
    temp = [a-mean for a in v]
    return temp/desvio_padrao(v)

#def get_parameters(timeValues,freqValues):
#    global max_freq,min_freq
#    max_i=int(max_freq*len(freqValues)/fs)
#    min_i=int(min_freq*len(freqValues)/fs)
#    mf=np.max([abs(a) for a in freqValues[min_i:max_i]])
#    mt=np.max([abs(a) for a in timeValues[min_i:max_i]])
#    imp_freq=[a*2/(len(freqValues)) for a in freqValues[min_i:max_i]]
#    imp_freq=freqValues
#    tyme=[a/(1000) for a in timeValues]
#
#    temp=[curtose(tyme),
#          obliquidade(tyme),
#          variancia(tyme),
#          media(tyme),
#          desvio_padrao(tyme),
#          rms(tyme)]
#    return temp
def get_parameters(timeValues,freqValues):
#    maxt=np.max(timeValues)
#    maxf=np.max(freqValues)
#    freq=[a*2/(len(freqValues)) for a in freqValues]
    freq=[a/(1000) for a in freqValues]
    tyme=[a/(1000) for a in timeValues]

#    freq=[a/(maxf) for a in freqValues]
#    tyme=[a/(maxt) for a in timeValues]
    temp=[curtose(tyme),
          obliquidade(tyme),
          variancia(tyme),
          media(tyme),
          desvio_padrao(tyme),
#          rms(tyme),
          centroide_espectral(freq)/10,
          curtose(freq),
          obliquidade(freq),
          variancia(freq),
          media(freq),
          desvio_padrao(freq),
          rms(freq)]
    return temp
    

# Individual class
class ind:
    def __init__(self,layers,biases):
        global guid
        self.uid=guid
        guid+=1
        self.fit=-1000
        self.brain=snkbrain.brain(layers,biases)

def getParameters():
    global param
    param=[[get_parameters(realV,freqV) for realV,freqV in zip(sr,sf)] for sr,sf in zip(real,frv)]
#    param=[[get_parameters(realV,freqV) for realV,freqV in zip(sr,sf)] for sr,sf in zip(real,frv)]

def feature_scaling_all():
    global flx_train,flx_test,ext_train,ext_test
    data_sets=[flx_train,flx_test,ext_train,ext_test]
    
    all_data=[]
    for data_set in data_sets:
        for arq in data_set:
            for data in arq:
                for data_piece in data:
                    all_data.append(data_piece)
    dp=desvio_padrao(all_data)
    mn=media(all_data)
    print dp,mn
    for i in range(0,len(data_sets)):
        for j in range(0,len(data_sets[i])):
            for k in range(0,len(data_sets[i][j])):
                data_sets[i][j][k]=(data_sets[i][j][k]-mn)/dp
                    
# This function takes the fft data od an signal, and returns a similar vector,
# but instead of getting one element per frequency it take a number of freqStep
# frequencies, sum it and divide by freqStep
def getFreqVector(fv):
    x=[]
    tam=float(len(fv))
    norm=int(np.ceil(tam*1/fs))
    step=freqStep*norm
    for j in range(0,norm*maxFreq,step):
        x.append(sum(fv[j:j+step])*2/tam)
    return x

def adjust(v):
    max_freq=240
    min_freq=6
    max_i=int(max_freq*len(v)/fs)
    min_i=int(min_freq*len(v)/fs)
    v=v[min_i:max_i]
    
    
# Read the data archives. The original signal is stored in origin. Each signal
# Is stored in real. real[arq][5] will contain the 5th signal of the arq'th file
# (as read by getArqs). The fft data will be stored at "fv" (indexes works the
# the same as for "real"). The frequency vector as got by getFrequencyVector
# is stored at frv
def readArqs(source,muscle,interval):
    it=interval
    reader=erc.emgReader()
    global real,fv,frv

            
    if source=='ninapro':
        global ninaprofolders
        realt,fvt=[],[]
        for folder in ninaprofolders:
            realt.append([])
            fvt.append([])
            realt[-1],fvt[-1]=reader.getCsvData(muscle,folder)
        
        for arq in range(0,len(realt[0])):
            real.append([])
            fv.append([])
           
            for r,f in zip(realt,fvt):                        
                real[arq].extend(r[arq][ it[0]:it[1] ])   
                fv[arq].extend(f[arq][ it[0]:it[1] ])
 
        training=[18-1,21-1,22-1,25-1,26-1,31-1]
        real=[real[i] for i in training]
        fv=[fv[i] for i in training]
    
    frv=[[getFreqVector(rep) for rep in arq]for arq in fv]
#    for arq in range(0,len(fv)):      
#        frv.append([])
#        for i in range(0,len(fv[arq])):
#            frv[arq].append(getFreqVector(fv[arq][i]))

# Fitness method. Each signal frequency vector is multiplied by indiv
# chromossome. The numbers got are reconized as the score of each archive.
# Let's say that the 0th element gets the largest number. That mean this 
# individual "thinks" that that signal belongs to archive 4 (getArqs()[0])
# The fitness is then calculated by the number of right guesses of each
# individual
            
def fitness(indiv,dp=False,indivs=None):
    
    global nArq
    score=0
    for arq in range(0,len(param)):
      
        for i in range(0,len(param[arq])):
            tam=len(param[arq][i])
            
#            inpt=frv[arq][i]
            if not dp:
                inpt=param[arq][i]
            else:
                l,o1=indivs[0].brain.run(ext_train[arq][i])
                l,o2=indivs[1].brain.run(flx_train[arq][i])
                inpt=[a for a in o1]
                inpt.extend([a for a in o2])
                
            l,pont= indiv.brain.run(inpt)#np.array(frv[arq][i])*indiv.cromo.freqFactor
            
            def error(pont,ref):
                score=0
                for i in range(0,len(pont)):
                    if i==ref:
                        t=1
                    else:
                        t=0
                    score+= t * np.log(pont[i]) + (1-t) * np.log(1-pont[i])
#                    score+= t*np.log((pont[i]+1)/2)+(1-t)*np.log(1-(pont[i]+1)/2)

                return score
            score+=error(pont,arq)
    return score

def trainSomeone(indiv,number,learning_rate):
    global param
    count=0
    while count<number:
        for arq in range(0,len(param)):
            for i in range(0,len(param[arq])):
                target=[0,0,0,0,0,0]
                target[arq]=1
                indiv.brain.train(param[arq][i],target,learning_rate)
        count+=1                

def trainSomeone_2(rn,indiv1,indiv2,number,learning_rate):
    global param
    count=0
    while count<number:
        for arq in range(0,len(param)):
            for i in range(0,len(param[arq])):
                target=[0,0,0,0,0,0]
                target[arq]=1

                l,o1=indiv1.brain.run(ext_train[arq][i])
                l,o2=indiv1.brain.run(flx_train[arq][i])
                f=[a for a in o1]
                f.extend([a for a in o2])
                
#                f=[0,0,0, 0,0,0, 0,0,0, 0,0,0]
#                f[arq],f[arq+6]=1,1
                    
                rn.brain.train(f,target,learning_rate)                
        count+=1  
# Ver acertos
def get_score(ind,data):
    score=0
    total=0
    for arq in range(0,len(data)):
        for i in range(0,len(data[arq])):
            total+=1
            l,o=ind.brain.run(data[arq][i])
            if np.argmax(o)==arq:
                score+=1
    return float(score)/total
        
#score=0
#total=0
#for arq in range(0,len(real)):
#    for i in range(0,len(fv[arq])):
#        total+=1
#        l,o1=a.brain.run(ext_test[arq][i])
#        l,o2=b.brain.run(flx_test[arq][i])
#        f=[j for j in o1]
#        f.extend([j for j in o2])
#        l,o=c.brain.run(f)
#        if np.argmax(o)==arq:
#            score+=1
#        else:
#            print arq,i
#print score,' de ',total

def treinar(ind,init_lr,goal,dp=False,indivs=None):
    global bestAll,tamPop
    last=-1000
    count=0
    lr=init_lr
    errorcount=0
    flag=0
    f=-1000
    lastbrain=[a.copy() for a in ind.brain.brainMatrixes]
    while f<goal:
        f=fitness(ind,dp=dp,indivs=indivs)

        if last>f:
            lr/=1.1
            errorcount+=1
            ind.brain.brainMatrixes=lastbrain
        else:
            lastbrain=[a.copy() for a in ind.brain.brainMatrixes]
            errorcount=0
            count+=1
        if count==11:
            lr*=1.03
            count=0
        if errorcount==3:
            if flag>=3:
                print 'to many errors'
                break
            else:
                flag+=1
                lr=init_lr
                errorcount=0
        if dp:
            trainSomeone_2(ind, indivs[0],indivs[1],10,lr)
        else:
            trainSomeone(ind,100,lr)
        last = f    
        print f,lr#,get_score(ind,ext_train),get_score(ind,ext_test)
    return lr


def treinar_v(a,init,reset_log=False):
    global ext_train,flx_train,ext_test,flx_test,param,param_test,param_vald,log_vald,log_train,log_test,ind_count
    if reset_log:
        log_vald,log_train,log_test=[],[],[]
    concatenate(ext_test,flx_test)
#    param_vald=[[[data for data in rep] for rep in arq[::2]] for arq in param]
    param_test=[[[data for data in rep] for rep in arq] for arq in param]
#    log_vald.append([])
    log_train.append([])
    log_test.append([])
    concatenate(ext_train,flx_train)
    count=0.0
    while count<3000:
        count+=1
        if count%100==0:
            print 'Log numero ',len(log_train),' ',count/30,'% completo'
        log_train[-1].append(get_score(a,param))
#        log_vald[-1].append(get_score(a,param_vald))
        log_test[-1].append(get_score(a,param_test))
#        print get_score(a,param),get_score(a,param_vald),get_score(a,param_test)
        trainSomeone(a,5,init)

def search_candidates():
    global inds
    inds=[]
    inds.append(ind([24,24,6],[1,1,0]))    
    inds.append(ind([24,20,6],[1,1,0]))    
    inds.append(ind([24,16,6],[1,1,0]))
    inds.append(ind([24,12,6],[1,1,0]))
    inds.append(ind([24,16,6],[1,0,0]))    
    inds.append(ind([24,24,20,6],[1,1,1,0]))
    inds.append(ind([24,20,16,6],[1,1,1,0]))
    inds.append(ind([24,16,12,6],[1,1,1,0]))
    inds.append(ind([24,20,16,6],[1,0,0,0]))
    
    for indiv in inds:
        treinar_v(indiv,0.005)
        ind_count+=1

#def get_all_parameters():
#    
#    global flx_train,flx_test,ext_train,ext_test
#   
#    def a1():
#        global ext_test
#        r1,f1=readArqs('ninapro','ext',[3,6])
#        ext_test=getParameters(r1,f1)
#        print 'done'
#
#    def a2():
#        global flx_test
#        r2,f2=readArqs('ninapro','flx',[3,6])
#        flx_test=getParameters(r2,f2)
#        print 'done'
#
#    def a3():
#        global ext_train
#        r3,f3=readArqs('ninapro','ext',[0,3])
#        ext_train=getParameters(r3,f3)
#        print 'done'
#
#    def a4():
#        global flx_train
#        r4,f4=readArqs('ninapro','flx',[0,3])
#        flx_train=getParameters(r4,f4)
#        print 'done'
#
#    threading.Thread(target=a1).start()
#    threading.Thread(target=a2).start()
#    threading.Thread(target=a3).start()
#    threading.Thread(target=a4).start()

def concatenate(v1,v2):
    global param
    param=[[[data for data in rep] for rep in arq] for arq in v1]
    for i in range(0,len(v2)):
        for j in range(0,len(v2[i])):
            param[i][j].extend(v2[i][j])

def get_all_parameters():
    
    global flx_train,flx_test,ext_train,ext_test
    
        
    readArqs('ninapro','flx',[0,6])
    getParameters()
    flx_test=   [[[data for data in rep] for rep in arq if     arq.index(rep)%3==2] for arq in param]
    flx_train=  [[[data for data in rep] for rep in arq if not arq.index(rep)%3==2] for arq in param]

#    readArqs('ninapro','flx',[3,6])
#    getParameters()
#    flx_test=[[[data for data in rep] for rep in arq] for arq in param]

    readArqs('ninapro','ext',[0,6])
    getParameters()
    ext_test=[[[data for data in rep] for rep in arq if arq.index(rep)%3==2 ] for arq in param]
    ext_train= [[[data for data in rep] for rep in arq if not arq.index(rep)%3==2] for arq in param]
    
    feature_scaling_all()
#    readArqs('ninapro','ext',[3,6])
#    getParameters()
#    ext_test=[[[data for data in rep] for rep in arq] for arq in param]
    
def print_all():
    global log_train,log_test
    for i in range(0,len(log_test)):
        plt.figure()
        plt.plot(range(0,len(log_train[i])*5,5),log_train[i],'r')
        plt.plot(range(0,len(log_test[i])*5,5),log_test[i],'b')
        plt.grid()
        plt.title('NN Topology '+chr(ord('a')-32+i))
        plt.legend(['Training','Testing'])
        plt.grid(True)
        plt.xlabel('Epoch')
        plt.ylabel('Percentage')
        movs=[18,21,22,25,26,31]
        m=''
        for n in movs:
            m+=str(n)+'_'
        plt.savefig(m+chr(ord('a')-32+i)+'.png')
        
if profiling:
    cProfile.run('main()')
else:
    1
#    main()
