commit 2830c112f30715abb5984088500c37e447e46d8e Author: adamr Date: Thu Apr 13 14:57:04 2023 +0200 Main initial version diff --git a/.Rhistory b/.Rhistory new file mode 100644 index 0000000..e69de29 diff --git a/EStiMo_GUI_0123.py b/EStiMo_GUI_0123.py new file mode 100644 index 0000000..43b493b --- /dev/null +++ b/EStiMo_GUI_0123.py @@ -0,0 +1,1460 @@ +#!/usr/bin/env python3 +# -*- coding: utf-8 -*- +""" +Created on Tue Nov 23 11:11:41 2021 + +@author: adamr +""" +import NeurOne +import RDA + +import time +import matplotlib +import mne +#import inspect +import sys, os +import numpy as np +import scipy.signal as ss +import matplotlib.pyplot as plt +import matplotlib as mpl +import matplotlib.cm as cm +import PyQt5.uic as uic +import pandas as pd +import json +#import scipy.stats as st +import csv +import datetime +import ctypes +import scipy +import pywt +import queue + +from cycler import cycler +from matplotlib.backend_bases import MouseButton +from PyQt5.QtCore import QTimer, Qt +from PyQt5.QtGui import QImage, QPixmap, QIcon, QFont +from PyQt5.QtWidgets import (QMainWindow, QFileDialog, QMessageBox, QCheckBox, QLineEdit, QWidget, QPushButton, + QLabel, QHBoxLayout, QGridLayout, QAction, QApplication, QDialog, QDialogButtonBox, + QVBoxLayout, QFrame) +from matplotlib.backends.backend_qt5agg import FigureCanvasQTAgg as FigureCanvas +from matplotlib.backends.backend_qt5agg import NavigationToolbar2QT as NavigationToolbar +from matplotlib.figure import Figure +from mne.channels.layout import _find_topomap_coords as get_pos +if sys.platform=='darwin': + import multiprocessing as mp + from multiprocessing import Process, Value + #from multiprocessing import Queue as QueueOld #Queue doesn't work on MacOS +else: + from multiprocessing import Process, Queue, Value +from Waiting import Waiting_window +from FirstWindow import First_window +from Functions import (apply_montage, eye_reg, most_frequent, connect_sig, set_to_gray, update_stem, + adjust_hist, min_zero, pentropy, entropy, check_integrity, doubleMADsfromMedian) + +if sys.platform=='darwin': + from multiprocessing.queues import Queue as QueueOld + class SharedCounter(object): + #source: https://gist.github.com/FanchenBao/d8577599c46eab1238a81857bb7277c9 + """A synchronized shared counter. + The locking done by multiprocessing.Value ensures that only a single + process or thread may read or write the in-memory ctypes object. However, + in order to do n += 1, Python performs a read followed by a write, so a + second process may read the old value before the new one is written by the + first process. The solution is to use a multiprocessing.Lock to guarantee + the atomicity of the modifications to Value. + This class comes almost entirely from Eli Bendersky's blog: + http://eli.thegreenplace.net/2012/01/04/shared-counter-with-pythons-multiprocessing/ + """ + + def __init__(self, n=0): + self.count = mp.Value('i', n) + + def increment(self, n=1): + """ Increment the counter by n (default = 1) """ + with self.count.get_lock(): + self.count.value += n + + @property + def value(self): + """ Return the value of the counter """ + return self.count.value + + + class Queue(QueueOld): + #source: https://gist.github.com/FanchenBao/d8577599c46eab1238a81857bb7277c9 + """ A portable implementation of multiprocessing.Queue. + Because of multithreading / multiprocessing semantics, Queue.qsize() may + raise the NotImplementedError exception on Unix platforms like Mac OS X + where sem_getvalue() is not implemented. This subclass addresses this + problem by using a synchronized shared counter (initialized to zero) and + increasing / decreasing its value every time the put() and get() methods + are called, respectively. This not only prevents NotImplementedError from + being raised, but also allows us to implement a reliable version of both + qsize() and empty(). + Note the implementation of __getstate__ and __setstate__ which help to + serialize MyQueue when it is passed between processes. If these functions + are not defined, MyQueue cannot be serialized, which will lead to the error + of "AttributeError: 'MyQueue' object has no attribute 'size'". + See the answer provided here: https://stackoverflow.com/a/65513291/9723036 + + For documentation of using __getstate__ and __setstate__ to serialize objects, + refer to here: https://docs.python.org/3/library/pickle.html#pickling-class-instances + """ + + def __init__(self): + super().__init__(ctx=mp.get_context()) + self.size = SharedCounter(0) + + def __getstate__(self): + """Help to make MyQueue instance serializable. + Note that we record the parent class state, which is the state of the + actual queue, and the size of the queue, which is the state of MyQueue. + self.size is a SharedCounter instance. It is itself serializable. + """ + return { + 'parent_state': super().__getstate__(), + 'size': self.size, + } + + def __setstate__(self, state): + super().__setstate__(state['parent_state']) + self.size = state['size'] + + def put(self, *args, **kwargs): + super().put(*args, **kwargs) + self.size.increment(1) + + def get(self, *args, **kwargs): + item = super().get(*args, **kwargs) + self.size.increment(-1) + return item + + def qsize(self): + """ Reliable implementation of multiprocessing.Queue.qsize() """ + return self.size.value + + def empty(self): + """ Reliable implementation of multiprocessing.Queue.empty() """ + return not self.qsize() + + +class NeurOneOffline(): + def __init__(self): + self.data = None + + def NO(self,ip,port=50000,buffersize=2**15,ringbuffersize=2000, + sendqueue=False,ringbuf_factor=2,dump=False,avgPackets=1): + self.ringbuffersize = ringbuffersize + tmp_path = '/mnt/projects/P_BCT_EEG/DLPFCM1_iTBS/DLPFC/nobeep/subj_8/X47851_iTBS.vhdr'# + #'/mnt/projects/P_BCT_EEG/DLPFCM1_iTBS/DLPFC/beep/subj_14/X13193_adam.vhdr' #'/mnt/projects/P_BCT_EEG/DLPFCM1_iTBS/DLPFC/beep/subj_6/X77384_iTBS.vhdr' + num_electr = 18 + eeg_chn = np.arange(0,num_electr,1) + hdr = mne.io.read_raw_brainvision(tmp_path) + + # hdr.set_channel_types({'EMGleft': 'emg', 'EOGright': 'eog'}) + # hdr.set_montage(mne.channels.read_custom_montage('easycap-M10_63_NO.txt')) + # mrk_fullpath = tmp_path[:-4]+'vmrk' + # eeg_fullpath = tmp_path[:-4]+'eeg' #this two are made by hand instead of function. + #Maybe there is some func for this + # Annotations returns all events - stimA, stimB, stopA, stopB, start of experiment etc... We chose only stim + stim = hdr.annotations.onset[np.logical_or(hdr.annotations.description=="Stimulus/A", + hdr.annotations.description=="Stimulus/B")] + # Separate stimA and stimB + stimA = hdr.annotations.onset[hdr.annotations.description=="Stimulus/A"] + stimB = hdr.annotations.onset[hdr.annotations.description=="Stimulus/B"] + #divide for stim A and B + #stimR = hdr.annotations.onset[hdr.annotations.description=='Response/R 16'] + npts = hdr.n_times + nfft = int(hdr.info['sfreq']) # Sampling rate [Hz] + fs = int(hdr.info['sfreq']) # Sampling rate [Hz] + endsample = npts + begsample = 0 + print('Sampling rate [Hz]: ' + str(fs)) + #[stim, resp, segment, timezero] = mne.read_annotations(mrk_fullpath) thinkkkkk + eeg_raw = hdr.get_data(start=begsample, stop=endsample); + eeg = eeg_raw[eeg_chn, :] # Select all electrodes + + data_raw = (hdr.get_data()*1e7) + # data_raw = data_raw[:,:]#int(data_raw.shape[1]%5)] + print(data_raw.shape) + stimB_arr = (stimB*fs).astype(int) + dividing_factor = int(fs/1000) + data_raw_new = np.zeros([data_raw.shape[0], int(data_raw.shape[1]//dividing_factor)-1]) + for i in range(data_raw.shape[0]): + a = data_raw[i, :-int(data_raw.shape[1]%dividing_factor)-dividing_factor] + R = dividing_factor + data_raw_new[i,:] = a.reshape(-1, R).mean(1) + data_raw = data_raw_new + fs=1000 + stimB_arr = (stimB_arr/dividing_factor).astype(int) + stim_sig = np.zeros(data_raw.shape[1]).reshape([1,data_raw.shape[1]]) + for i in stimB_arr: + stim_sig[0,i] = 1 + self.fs = fs + self.data = np.concatenate((data_raw, stim_sig)) + self.buffer = self.data[:,:ringbuffersize] + + def start(self): + self.time = time.time() + + def getBuffer(self): + time_now = time.time() + time_diff = int((time_now - self.time)*self.fs) + return self.data[:,time_diff-self.ringbuffersize:time_diff] + + + +def acquire_data(q, size, run, speed, downsample, sleep_time, ip = '192.168.200.201', port = 5000, offline=False): + """Acquires data from NeurOne_v3 and pass it to the queue. Function is supposed to work + in separated process. + + Parameters + ------------ + q: Queue class from multiprocessing library + size: number of samples to acquire + run: Value class from multiprocessing library. That value can be changed in main process + downsample: boolean value. Says if data will be downsampled to 1000 Hz + sleep_time: int, set how often function should refresh. Usually it takes a bit more that that""" + offline = 'offline' + #import NeurOne_v3 + if offline=="offline": + NO = NeurOneOffline() + NO.NO(ip=ip,port=50000,buffersize=2**15,ringbuffersize=size,\ + sendqueue=False,ringbuf_factor=2,dump=False,avgPackets=downsample) #offline + elif offline=="NeurOne": + # ip='192.168.200.201' + # port = 50000 + NO = NeurOne.NO(ip=ip,port=port,buffersize=2**15,ringbuffersize=size,\ + sendqueue=False,ringbuf_factor=2, dump=None,avgPackets=downsample) #online + elif offline=="BrainProducts": + # ip = "169.254.252.66" + # port = 51244 + NO = RDA.RDA(ip=ip, port=port, buffersize=2**15, ringbuffersize=size,\ + sendqueue=False, avgPackets=downsample) + NO.start() + while True: + startt = time.time() + if run.value: + data = NO.getBuffer()#[channel_list] + try: + q.put(data) + except: print('ERROR') + endd=time.time() + time.sleep(sleep_time*speed.value-(endd-startt)) + +class AppForm(QMainWindow): + """TMS GUI script.""" + def __init__(self, passed_params = None, parent=None): + super().__init__() + #self.setStyleSheet("background: whitesmoke") + # self.offline = 'NeurOne' #"BrainProducts"#False + self.time_start = time.time() + QMainWindow.__init__(self, parent) + self.montage_file_path = 'montage_18ch.csv' + self.features() + self.speed_general = 250 # fs has to be divisible by this value + self.all_params(passed_params) + self.run = Value('i', 0) # Value that can be changed inside separated process + # To control data flow speed we can initialize another variable used inside the process + self.speed = Value(ctypes.c_float, self.speed_general/1000) + self.q = Queue() # Queue for data from separated process (data acquire) + # Run acquire_data function in separated process so it doesn't freeze when other operations are going + self.p = Process(target=acquire_data, args=(self.q, self.size_of_up, + self.run, self.speed, True, + self.update_time, self.ip, self.port, + self.offline)) + self.p.start() #start it + self.timer = QTimer(timerType = Qt.PreciseTimer) #preparing timer + self.timer.setInterval(self.speed_general) #every 1000 ms = 1 sec + self.timer.timeout.connect(self.update) #use update fuction every 1 sec + #self.timer.start() + + #set name and icon (for some reason the icon doesn't appear on the taskbar) + self.setWindowTitle('EStiMo') + self.setWindowIcon(QIcon('icon.jpg')) + + def features(self): + """Collection of implemented features. Should be moved to separated .py + file in the future, to make adding new features easier. For now it's part + of the class, what also optimizes it slightly, because some of parameters are + shared within a class + + TO ADD NEW FEATURES: check the bottom of this function""" + + def calculate_theta(): #1 + """Theta band calculation based on previously calculated FFT/Welch""" + theta = self.S[:, self.theta_band[0]:self.theta_band[1]].mean(1) + return theta.mean() + + def calculate_alpha(): #2 + """Alpha band calculation based on previously calculated FFT/Welch""" + alpha = self.S[:, self.alpha_band[0]:self.alpha_band[1]].mean(1) + return alpha.mean() + + def calculate_beta(): #3 + """Beta band calculation based on previously calculated FFT/Welch""" + beta = self.S[:,self.beta_band[0]:self.beta_band[1]].mean(1) + return beta.mean() + + def calculate_h_gamma(): #4 + """High Gamma band calculation based on previously calculated FFT/Welch""" + high = self.S[:,80:120].mean(1) + return high.mean() + + def calculate_temp_entropy(): #5 + """Temporal Entropy calculation""" + return np.mean(np.mean([entropy(self.second_to_analyze[e_ch]) for e_ch + in range(self.second_to_analyze.shape[0])])) + + def calculate_spect_entropy(): #6 + """Spectral Entropy calculation based on previously calculated FFT/Welch""" + return np.mean(np.mean([entropy(self.S[e_ch]) for e_ch in range(self.S.shape[0])])) + + def calculate_line_length(): #7 + """Line length calculation. The function scipy.spatial.distance.cdist makes it + very fast compared to other methods""" + ch_sum = 0 + for i in range(self.second_to_analyze.shape[0]): + a_seg = np.zeros([2, self.second_to_analyze.shape[1]]) + a_seg[0,:] = self.second_to_analyze[0]#*1e6 + a_seg[1,:] = np.linspace(0,1,self.Fs) + matrix_dist = scipy.spatial.distance.cdist(a_seg, a_seg, 'euclidean') + ch_sum += sum(np.diagonal(matrix_dist,1)) + return ch_sum + + def wavelets(band): #8,9,10,11 + """Power in given band based on DWT calculated before. Works correctly + only for fs=1000 Hz!!! With different Fs it will give power for different band. + Ranges for 1000 Hz are: 1000, 500, 250, 125, 62.5, 31.25, 15.625, 7.8125, 3.9062, 0. + We are interested in lasts of them: + 8: 0 - 3.906 Hz + 9: 3.906 - 7.812 Hz + 10: 7.812 - 15.625 Hz + 11: 15.625 - 31.25 Hz + """ + band = band - 9 #weird way, but then call from function works well for particular band + #but works correctly only for this set of features! + dwt_pw1, dwt_pw2, dwt_pw3, dwt_pw4 = self.dwt[:4] + #sum of squares + dwt_pw1 = np.sum(dwt_pw1**2) + dwt_pw2 = np.sum(dwt_pw2**2) + dwt_pw3 = np.sum(dwt_pw3**2) + dwt_pw4 = np.sum(dwt_pw4**2) + dwts = dwt_pw1, dwt_pw2, dwt_pw3, dwt_pw4 + return dwts[band] + + def calculate_variance(): #12 + """Sum of Variance calculation""" + return sum(np.var(self.second_to_analyze,1)) + + def calculate_corr(): #13 + """Temporal Correlation between channels calculation. Returnes mean of absoulte + values of correlations between channels""" + R = np.corrcoef(self.second_to_analyze) + if type(R)==np.ndarray: + R = abs(np.tril(R, -1)) #tril zeroes one triangle of array + R_sum = R[R!=0] #we take only what is left + if type(R)==np.ndarray: #in case it is just only one value (that means we analyze only one channels) + R_sum = R_sum[~np.isnan(R_sum)] + return np.mean(np.abs(R_sum)) + ############################################### + # You can add new functions now: + + + + ############################################### + # Add all new functions to this list + self.functions = [None, calculate_theta, calculate_alpha, calculate_beta, calculate_h_gamma, + calculate_spect_entropy, calculate_temp_entropy, calculate_line_length, + wavelets, wavelets, wavelets, wavelets, calculate_variance, calculate_corr] + + # Set the unit for the output of each function (they correspond to functions based on index) + self.unit_label = ['None', "Power", "Power", "Power", "Power", "Entropy", "Entropy", "Length", + "Power", "Power", "Power", "Power", "Variance", "Correlation"] + + # Add name of the added feature to the end of this list + self.features_names = ['None', 'Theta FFT Power', 'Alpha FFT Power', 'Beta FFT Power', + 'High Gamma FFT Power', 'Spectral entropy', 'Temporal entropy', + 'Line length', 'DWT Power 0-4 Hz', 'DWT 4-8 Hz', + 'DWT 8-16 Hz', 'DWT 16-31 Hz','Variance','Correlation'] + # One last thing to use your new function is to add string with its name + # to the other First_window.py: variable features_names in class First_window + + def all_params(self, passed_params = None, restarted=''): + #reads values from TMS_protocol.txt file + + #montage matrix is a matrix that is multiplied with the signal + #identity matrix multiply all channels by 1, so we have the same signal as an output + settings_file = pd.read_csv('TMS_protocol.txt',sep=':', header=None) #read file with settings + print(restarted) + #open the file to save values during the experiment + self.log_file = open('logs//TMS_log_'+str(f"{datetime.datetime.now():%Y-%m-%d-%H-%M}"+restarted+'.csv'), 'w') + self.log_file_writer = csv.writer(self.log_file) + + self.log_file.flush() #flushing is applying changes we made to the file + self.Fs = 1000 #5000 + self.size_of_up = 2*self.Fs #5000 how much data we get from NeurOne function + + #file with settings and assigning variables to them + settings_file = pd.read_csv('TMS_protocol.txt',sep=':', header=None) + + #params can be set in GUI, so then no need for using ones from the file + if passed_params is not None: + self.montage_file_path = passed_params['montage_path'] + self.time_between_bursts = int(passed_params['time_between']) + self.breaktime = int(passed_params['cut_time']) + self.included_ch = passed_params['included_ch'] + self.eog_ch_l = list(passed_params['eog_ch']) + self.emg_ch_l = list(passed_params['emg_ch']) + self.num_of_ch = int(passed_params['num_of_ch']) + self.num_of_lines = int(passed_params['num_of_lines']) + self.ch_names = passed_params['ch_names'] + self.slow_mode = passed_params['slow_mode'] + self.used_features = passed_params['features'] + self.use_notch = passed_params['notch'] + self.use_regression = passed_params['eye_reg'] + self.remove_outliers = passed_params['remove_outliers'] + self.ip = passed_params['ip'] + self.port = passed_params['port'] + self.offline = passed_params['offline'] + self.exp_trig = passed_params['exp_trig'] + self.exp_time = passed_params['exp_time'] + else: + self.montage_file_path = 'montage_18ch.csv' + self.time_between_bursts = int(settings_file[settings_file[0]=='time_between_trains'].values[0][1]) + self.breaktime = int(settings_file[settings_file[0]=='cut_time'].values[0][1]) + self.included_ch = settings_file[settings_file[0]=='included_channels'].values[0][1].strip() + self.eog_ch_l = [int(settings_file[settings_file[0]=='eog_channel'].values[0][1])] + self.emg_ch_l = [int(settings_file[settings_file[0]=='emg_channel'].values[0][1])] + self.num_of_ch = int(settings_file[settings_file[0]=='number_of_channels'].values[0][1]) + self.num_of_lines = int(settings_file[settings_file[0]=='number_of_lines'].values[0][1]) + self.exp_trig_loaded = int(settings_file[settings_file[0]=='expected_triggers'].values[0][1]) + self.exp_time_loaded = int(settings_file[settings_file[0]=='expected_time'].values[0][1]) + self.slow_mode = False + self.used_features = [0,1,2,3,4,5] + self.use_notch = True + self.use_regression = True + self.remove_outliers = True + self.ip = '192.168.200.201' + self.port = 50000 + self.offline = False + + self.unit_label = np.array(self.unit_label)[self.used_features] + self.log_file_writer.writerow(['time', 'state', self.used_features]) + + if len(self.eog_ch_l)==0: + self.eog_ch = '' + if len(self.emg_ch_l)==0: + self.emg_ch = '' + + self.eog_ch = self.eog_ch_l[0] + self.emg_ch = self.emg_ch_l[0] + + if self.slow_mode: + self.speed_general = 1000 + else: + self.speed_general = 250 + + if self.included_ch in ['all', 'All']: + self.included_ch = np.arange(len(self.ch_names)) + elif type(self.included_ch)==list: + pass + else: + self.included_ch = json.loads(self.included_ch) + print(self.included_ch) + # self.montage_matrix = np.identity(self.num_of_ch) + if self.montage_file_path in [None, '']: + self.montage_file_path = 'montage_18ch.csv' + self.montage_matrix = np.array(pd.read_csv(self.montage_file_path, header=None)) + + + #initialization of variables. Some of them are loaded from file and cannot be changed in the GUI + self.no_eeg = len(self.eog_ch_l) + len(self.emg_ch_l) + # print(self.no_eeg) + # self.ch_names = json.loads(settings_file[settings_file[0]=='names'].values[0][1]) + self.update_time = 1 #how many seconds between updates of everything + #prepare empty arrays for data + self.feature1 = np.full([1000, self.time_between_bursts-self.breaktime], None) + self.feature2 = np.full([1000, self.time_between_bursts-self.breaktime], None) + self.feature3 = np.full([1000, self.time_between_bursts-self.breaktime], None) + self.feature4 = np.full([1000, self.time_between_bursts-self.breaktime], None) + self.feature5 = np.full([1000, self.time_between_bursts-self.breaktime], None) + self.feature6 = np.full([1000, self.time_between_bursts-self.breaktime], None) + #json.loads takes string and returns list, so we can read list from txt file + self.alpha_band = json.loads(settings_file[settings_file[0]=='alpha_range'].values[0][1]) + self.beta_band = json.loads(settings_file[settings_file[0]=='beta_range'].values[0][1]) + self.theta_band= json.loads(settings_file[settings_file[0]=='theta_range'].values[0][1]) + self.colors = ['b', 'm', 'r', 'k', 'c', ] #colors used for lines if less that 6 of them + self.data_len = 30*self.Fs #length of the data array in seconds + self.plot_len = 4 #length of data to plot (last seconds of data array) + self.plot_dividing_factor = 100 + self.previous_state = np.zeros(6) + if self.num_of_lines>5: #if more than 5 lines then colors of them from colormap + self.colors = cm.plasma(np.linspace(0,1,self.num_of_lines)) + #Example of thresholds for each measurment. Changes after first second of calibration. + self.thr_parameter = int(settings_file[settings_file[0]=='threshold_parameter'].values[0][1]) + self.thr_2 = [20,100] + self.thr_1 = [20,100] + self.thr_3 = [20,100] + self.thr_4 = [0, 0.5] + self.thr_4 = [0, 0.5] + self.thr_5 = [0.3, 0.8] + self.thr_6 = [1,3] + self.linetemps = np.zeros(6, dtype=object) + self.labels_size = 7 + self.lab_x_dist = 1. + self.lab_y_dist = 0.5 + self.measures_x_lims = np.array([[-1,1],[0,350],[0,350],[0,350],[0,2],[3,8]], dtype=float) + self.measures_y_lims = np.array([[0,70],[0,70],[0,70],[0,70],[0,70],[0,70]], dtype=float) + self.create_main_frame() #create plots, buttons, figures etc... + + #create data array + self.loaded = np.zeros([self.num_of_ch,self.data_len]) + self.loaded_full = np.zeros([self.num_of_ch+1,self.data_len]) + self.data = np.random.rand(self.num_of_ch,self.data_len) + self.trigg_data = np.zeros(self.data_len) #array to keep trigger data in + self.num=0 + self.doit=0 #to count number of seconds after last stimuli in train + self.plot_story = np.zeros([self.num_of_lines, 6], dtype=object) #thigs that are being plot + self.trial_num = 0 #how many train there were + self.auto_readout_lim = True + self.disp_max = True + self.do_calibration = False #if True then calibration process is going on + self.calibration_counter = 0 + self.calibration_counter_max = int(1000/self.speed_general) + #Lists to keep calibration values for each measurment + self.f1_cal = [] + self.f3_cal = [] + self.f2_cal = [] + self.f4_cal = [] + self.f6_cal = [] + self.f5_cal = [] + self.qmbx = None + self.red_dots = [] + + def data_processing(self, calibration=False): + """Pre-processing of the data used for extracting features""" + + #in case of readout period times have to be adjusted to the last pulse, during calibration it's just a last second. + if calibration: + self.last_sec = self.loaded_noeye[:,-self.Fs:].copy() + if self.use_regression: + eog = self.loaded_noeye[self.eog_ch, -self.Fs:] + else: + self.last_sec = self.loaded_noeye[:,self.id1:self.id2].copy() + if self.use_regression: + eog = self.loaded_noeye[self.eog_ch ,self.id1:self.id2] + + [A,B] = ss.butter(2, 0.1/(self.Fs/2), 'highpass') + self.last_sec = ss.filtfilt(A, B, self.last_sec) + if self.use_notch: + [A,B] = ss.butter(2, [49/(self.Fs/2), 51/(self.Fs/2)], 'bandstop') + self.last_sec = ss.filtfilt(A, B, self.last_sec) + self.last_sec = ss.detrend(self.last_sec, axis=1) + + # plt.figure() + # plt.plot(self.last_sec[self.included_ch].T) + # plt.plot(eog) + #that's stupid, move channel selection before!!!!!!!!!!!!! + if self.use_regression: + self.last_sec = eye_reg(self.last_sec[self.included_ch], eog) + return self.last_sec + + def update_methods(self): + """Calculates all features and plot their values""" + times = time.time() + self.second_to_analyze = self.data_processing() #Take data + + #calculate fft + self.S = abs(np.fft.rfft(self.second_to_analyze)) + + #do dwt only if any features requires it + if any(feature_num in self.used_features for feature_num in [9,10,11,12]): + self.dwt = pywt.wavedec(self.second_to_analyze, 'db1', level=8) #db1, 8 levels. Works well for 1000 Hz only + + self.results = np.zeros(len(self.used_features)) + + #Put all results to one array + for idx, feature in enumerate(self.used_features): + if feature==0: + self.results[idx] = None + elif feature in [8,9,10,11]: #to specify which band + self.results[idx] = self.functions[feature](feature) + else: + self.results[idx] = self.functions[feature]() + + #Checks how many fields were filled already, so we know what stage are we on + #and where to save the data + x = np.where(self.feature1==None)[0][0] + y = np.where(self.feature1==None)[1][0] + + if y==0: + self.previous_state = np.zeros(6) + + self.feature1[x,y] = self.results[0] + self.feature2[x,y] = self.results[1] + self.feature3[x,y] = self.results[2] + self.feature4[x,y] = self.results[3] + self.feature5[x,y] = self.results[4] + self.feature6[x,y] = self.results[5] + + #save in the log + self.log_file_writer.writerow([time.time() - self.time_start, 'between bursts', + self.results]) + self.log_file.flush() + times1 = time.time() + print("Calculation of features: {}".format(times1-times)) + + # #set xlim of plots + # self.axesMap[0,0].set_xlim(0.5,self.time_between_bursts-self.breaktime+0.5) + # self.axesMap[0,1].set_xlim(0.5,self.time_between_bursts-self.breaktime+0.5) + # self.axesMap[1,0].set_xlim(0.5,self.time_between_bursts-self.breaktime+0.5) + # self.axesMap[1,1].set_xlim(0.5,self.time_between_bursts-self.breaktime+0.5) + # self.axesMap[2,0].set_xlim(0.5,self.time_between_bursts-self.breaktime+0.5) + # self.axesMap[2,1].set_xlim(0.5,self.time_between_bursts-self.breaktime+0.5) + + #if value is not within threshold values then background color is red (salmon), otherwise green + + #checks if there is a need to change a color of the background + old_prv_state = self.previous_state.copy() + if not all(np.isnan(self.thr_1)): + if self.results[0]>=self.thr_1[0] and self.results[0]<=self.thr_1[1]: + if self.previous_state[0]==1: + self.axesMap[0,0].set_facecolor('whitesmoke') + self.previous_state[0]=0 + elif self.previous_state[0]==0: + self.red_dots.append(self.axesMap[0,0].plot(y+1, self.results[0], 'r*')) + self.axesMap[0,0].set_facecolor('xkcd:salmon') + self.previous_state[0]=1 + + if not all(np.isnan(self.thr_2)): + if self.results[1]>=self.thr_2[0] and self.results[1]<=self.thr_2[1]: + if self.previous_state[1]==1: + self.axesMap[0,1].set_facecolor('whitesmoke') + self.previous_state[1]=0 + elif self.previous_state[1]==0: + self.red_dots.append(self.axesMap[0,1].plot(y+1, self.results[1], 'r*')) + self.axesMap[0,1].set_facecolor('xkcd:salmon') + self.previous_state[1]=1 + + if not all(np.isnan(self.thr_3)): + if self.results[2]>=self.thr_3[0] and self.results[2]<=self.thr_3[1]: + if self.previous_state[2]==1: + self.axesMap[1,0].set_facecolor('whitesmoke') + self.previous_state[2]=0 + elif self.previous_state[2]==0: + self.red_dots.append(self.axesMap[1,0].plot(y+1, self.results[2], 'r*')) + self.axesMap[1,0].set_facecolor('xkcd:salmon') + self.previous_state[2]=1 + + if not all(np.isnan(self.thr_4)): + if self.results[3]>=self.thr_4[0] and self.results[3]<=self.thr_4[1]: + if self.previous_state[3]==1: + self.axesMap[1,1].set_facecolor('whitesmoke') + self.previous_state[3]=0 + elif self.previous_state[3]==0: + self.red_dots.append(self.axesMap[0,0].plot(y+1, self.results[3], 'r*')) + self.axesMap[1,1].set_facecolor('xkcd:salmon') + self.previous_state[3]=1 + + if not all(np.isnan(self.thr_5)): + if self.results[4]>=self.thr_5[0] and self.results[4]<=self.thr_5[1]: + if self.previous_state[4]==1: + self.axesMap[2,0].set_facecolor('whitesmoke') + self.previous_state[4]=0 + elif self.previous_state[4]==0: + self.red_dots.append(self.axesMap[2,0].plot(y+1, self.results[4], 'r*')) + self.axesMap[2,0].set_facecolor('xkcd:salmon') + self.previous_state[4]=1 + + if not all(np.isnan(self.thr_6)): + if self.results[5]>=self.thr_6[0] and self.results[5]<=self.thr_6[1]: + if self.previous_state[5]==1: + self.axesMap[2,1].set_facecolor('whitesmoke') + self.previous_state[5]=0 + elif self.previous_state[5]==0: + self.red_dots.append(self.axesMap[2,1].plot(y+1, self.results[5], 'r*')) + self.axesMap[2,1].set_facecolor('xkcd:salmon') + self.previous_state[5]=1 + #If there was no change in a color, not need to redraw figure + if all(self.previous_state==old_prv_state): + need_redraw = False + else: #but otherwise we need to redraw the figure + need_redraw = True + alfa = 1 + color = 'b' + #if first plot + if self.trial_num0: + self.timer.setInterval(int(1*self.speed_general*0.97)) + if self.q.qsize()<1 and self.timer.interval()!= int(self.speed_general*1.1): + self.timer.setInterval(int(self.speed_general*1.03)) + + times = time.time() + self.offline='offline' #remove this! + if self.offline=="offline": + incl = [0,2,6,7,8,10,13,16,18,22,25,28,31,34,41,43,-3,-2,-1] # For offline only + loaded_temp = self.q.get()[incl]/10 # Load data + try: + self.loaded, most_fr = connect_sig(self.loaded_full, loaded_temp, self.speed_general) + except ValueError: + print("ValueError, wait...") + return 0 + else: + data_divide = 1 + if self.offline=="NeurOne": + data_divide = 10 + loaded_temp = self.q.get()/data_divide + try: + self.loaded, most_fr = connect_sig(self.loaded_full, loaded_temp.T, self.speed_general) + except ValueError: + print("ValueError, wait...") + return 0 + + self.loaded_noeye = self.loaded.copy() + + step1 = time.time()-time_start + + #if connection point is too far or too close it modifies the acquisition rate. + try: + if most_fr>900: + self.speed.value=1.01*self.speed_general/1000 + elif most_fr<500: + self.speed.value=0.8*self.speed_general/1000 + elif most_fr<800 and most_fr>500: + self.speed.value = 0.98*self.speed_general/1000 + print('speed:',self.speed.value) + #TODO: check if that is still needed + except UnboundLocalError: + pass + + self.trigg = self.loaded[-1].copy() #keeps trigger in a different array + self.loaded_full = self.loaded.copy() #keeps it for the next second + step2 = time.time()-time_start + # ZeroDivisionError means that data is not yet loaded + try: + # [A,B] = ss.butter(2, 0.1/(self.Fs/2), 'highpass') + # self.loaded[:self.num_of_ch,-4*self.Fs:] = ss.filtfilt(A, B, self.loaded[:self.num_of_ch, -4*self.Fs:]) + # self.loaded[:self.num_of_ch,-4*self.Fs:] = self.loaded[:self.num_of_ch,-4*self.Fs:] - np.mean(self.loaded[:self.num_of_ch,-4*self.Fs:],1, keepdims=True) + self.loaded[:self.num_of_ch,-4*self.Fs:] = ss.detrend(self.loaded[:self.num_of_ch,-4*self.Fs:]) + except ZeroDivisionError: # This error means that buffer is still not full + if self.qmbx == None: + self.qmbx = Waiting_window() # Small window with a message to wait + print('Waiting for data (should take up to few seconds)') + return(0) # Stop function here, because in this case rest of update is not needed + + # Getting to this point of the code means there were no exception before + # So the waiting window can be destroyed + if self.qmbx!=None: + self.qmbx.setText('Data received!') + self.qmbx=None #then it dies + + step3 = time.time()-time_start + + # Eye regression + self.trigg_data = self.trigg.copy() + od = self.data_len-self.plot_len*self.Fs # First point + do = self.data_len # Last point + + # For interpolation of TMS artifact + int_from = 0.3 #first point for interpolation + int_to = 0.3 #last point for interpolation + + step4 = time.time()-time_start + + # If two stimuli (A and B) are sent together (shouldn't be) then it removes one to not create chaos + stim = check_integrity(np.where(self.trigg_data[od:do]>0)[0]) + + step5 = time.time()-time_start + + # If there are some pulses detected then it interpolates the data used for the visualisation + if len(stim)>0: + for ind in stim[::-1]: + size = self.data_len - (od+ind-int(int_from*self.Fs)) + # print('size:', size) + # print(len(self.loaded)) + # Interpolation - pretty long line, but basically it chooses ranges and + # assign boundary value as a baseline and does that in (I guess) more optimal way than using loops + self.loaded[:, od+ind-int(int_from*self.Fs):od+ind+int(int_to*self.Fs)] = np.outer( + self.loaded[:,min(od+ind+int(int_to*self.Fs), 30000-1)], np.ones(min(size, int((int_from+int_to)*self.Fs)))) + + # for i in range(self.loaded.shape[0]): + # self.loaded[i, od+ind-int(int_from*self.Fs):od+ind+int(int_to*self.Fs)] = np.linspace( + # self.loaded[i,min(od+ind-int(int_from*self.Fs),30000-1)], + # self.loaded[i,max(od+ind+int(int_to*self.Fs),30000-1)], + # od+ind+int(int_to*self.Fs)-(od+ind-int(int_from*self.Fs))) + + step6 = time.time()-time_start + # Eye regression for visualisation + #plt.figure() + #plt.plot(self.loaded[3, -4*self.Fs:]) + #plt.plot(self.loaded[self.eog_ch,-4*self.Fs:], '--b') + if self.use_regression: + self.loaded[:(self.num_of_ch-self.no_eeg),-4*self.Fs:] = eye_reg( + self.loaded[:self.num_of_ch-self.no_eeg, -4*self.Fs:], self.loaded[self.eog_ch,-4*self.Fs:]) + #plt.plot(self.loaded[3,-4*self.Fs:], '--r') + #plt.show() + step7 = time.time()-time_start + # Applies montage for both types of data (visualisation and processing) + self.loaded = apply_montage(self.loaded, self.montage_matrix) + self.loaded_noeye = apply_montage(self.loaded_noeye, self.montage_matrix) + step8 = time.time()-time_start + self.update_plot() #updates main plot with data + stim_where = check_integrity(np.where(self.trigg_data>0)[0]) #indexes of triggers + + step9 = time.time()-time_start + + # If do_calibration is True and there were trigger recently then stop calibration + # TODO: THAT'S CONDITION REQUIRED FOR STARTING MEASURMENTS. PROBABLY IT WON'T WORK FOR SOME MORE EXTREME SETTINGS + # if self.do_calibration and len(stim_where[stim_where>(self.data_len-max( + # 2000, round(self.Fs*self.time_between_bursts*0.7, -3)))])>0: + if self.do_calibration and len(stim_where[stim_where>(self.data_len-max(2000, self.exp_time+1500))])>0: + self.calibration() + + if self.do_calibration: + # If do_calibration is True then continue calibration process and exit this function + if self.calibration_counter%self.calibration_counter_max==0: + self.calibration_process() + self.calibration_counter+=1 + return 0 #ends run of this function so nothing else happens + step10 = time.time()-time_start + # If set number of stimuli is detected + # doit --> set length of the measurment between bursts + # if self.doit==0 and sum(stim_where>self.data_len-max( + # 2000, round(self.Fs*self.time_between_bursts*0.7, -3)))==10: + print(sum(stim_where>self.data_len-max(2000, self.exp_time+1500))) + if self.doit==0 and sum(stim_where>self.data_len- + max(2000, self.exp_time+1500))==self.exp_trig: + self.axesMap[0,0].set_facecolor('whitesmoke') + self.axesMap[0,1].set_facecolor('whitesmoke') + self.axesMap[1,0].set_facecolor('whitesmoke') + self.axesMap[1,1].set_facecolor('whitesmoke') + self.axesMap[2,0].set_facecolor('whitesmoke') + self.axesMap[2,1].set_facecolor('whitesmoke') + self.doit=self.time_between_bursts-self.breaktime + self.last_stim = stim_where[-1] + set_to_gray(self.plot_story) + for dot in self.red_dots: + for line in dot: line.remove() + self.red_dots = [] + self.canvasMap.draw() + + # If doit is more than 0 --> readout phase + elif self.doit>0: + print(self.doit) + self.last_stim = stim_where[-1] + # In some situations last stimuli might be already from another train, while + # First 100ms of loaded signal belong to previous one. That is why this exception exists + # EDIT: not sure if it's still needed after other changes I made + if any(np.diff(stim_where[stim_where>self.data_len-max( + 2000, round(self.Fs*self.time_between_bursts*1.2, -3))])>1000): + print('UWAGA NA TO') + self.last_stim = stim_where[np.where(np.diff(stim_where)>1000)[0]][-1] + #index of last stim + half of breaktime + seconds that were already read before + self.id1 = int(self.last_stim+0.5*(self.breaktime*self.Fs)+( + self.time_between_bursts-self.breaktime-self.doit)*self.Fs) + self.id2 = int(self.id1 + self.Fs) #one second later + print('id1',self.id1) + if self.id20)[0]) + + # dif = np.max(loaded[:,-self.plot_len*self.Fs:])-np.min( + # loaded[:,-self.plot_len*self.Fs:]) #difference between max and min + #include only plot + for i in range(self.num_of_ch): + self.data[i] = loaded[i].copy() #dif*i+dif #to fit into one plot + #set data, last plot_len seconds + plot_data = self.data[i,self.data_len-self.plot_len* + self.Fs:self.data_len] + #EMG has higher amplitude usually. A special case to make it smaller + #self.emg_ch+1 because self.num_of_ch doesn't include trigger + if i==np.arange(self.num_of_ch)[self.emg_ch+1] and self.emg_ch!='': + plot_data = plot_data/(self.plot_dividing_factor*0.1) + self.num_of_ch - i + + if i==np.arange(self.num_of_ch)[self.eog_ch+1]: + plot_data = plot_data/(self.plot_dividing_factor*2) + self.num_of_ch - i + + else: + plot_data = plot_data/self.plot_dividing_factor + self.num_of_ch - i + # plot_data = plot_data/(3.5*np.max(np.abs(plot_data))) + self.num_of_ch - i + self.line[i].set_data(np.arange(0,self.plot_len*self.Fs), plot_data) + self.axes.set_ylim(0, self.num_of_ch+1) + #self.axes.set_ylim(0,np.max(self.data[:,-self.plot_len*self.Fs:])+dif) #set ylim to fit everything on the plot + if len(stim)>0: + for ind in stim: + self.axes.axvline(ind) #plot vertical line for each trigger + self.num = len(stim) + # plt.figure() + # plt.plot(ss.detrend(self.data[10,self.data_len-self.plot_len* + # self.Fs:self.data_len])) + self.canvas.draw() #update canvas + for i in range(self.num): + self.axes.lines[-1].remove() #if trigger outside of plot range then remove vertical line + self.num=0 + print('main plot time:', time.time()-startt) + + def start_stop(self): + """Starts and stops timer and changes value of run value for data acquiring process""" + if self.timer.isActive(): + self.timer.stop() + else: + self.timer.start() + self.run.value = not self.run.value + + def calibration(self): + """Prepares the process of calibration, changes a button's text and color""" + self.do_calibration = not self.do_calibration + if self.do_calibration: + #This prepares the software for calibration phase + self.axesMap[0,0].set_facecolor('whitesmoke') + self.axesMap[0,1].set_facecolor('whitesmoke') + self.axesMap[1,0].set_facecolor('whitesmoke') + self.axesMap[1,1].set_facecolor('whitesmoke') + self.axesMap[2,0].set_facecolor('whitesmoke') + self.axesMap[2,1].set_facecolor('whitesmoke') + + self.button4.setText("Stop calibration") + self.button4.setStyleSheet('background-color : lawngreen') + #Need to clean the figure to prepare is for a different type of plot + for i in range(6): + self.axesMap[i//2,i%2].cla() + + self.f1_cal = [] + self.f2_cal = [] + self.f3_cal = [] + self.f4_cal = [] + self.f5_cal = [] + self.f6_cal = [] + + #Change labels and make sure labels are ok + for i in range(6): + #if not all(np.isnan(self.thrs[i])): + self.axesMap[i//2,i%2].set_ylabel("Counts", fontsize=self.labels_size) + self.axesMap[i//2,i%2].set_xlabel(self.unit_label[i], fontsize=self.labels_size) + self.axesMap[i//2,i%2].set_title(self.Titles[i], fontsize=9) + + else: + self.trial_num = 0 + self.plot_story = np.zeros([self.num_of_lines, 6], dtype=object) + #This is post-calibration phase + self.button4.setText("Start calibration") + self.button4.setStyleSheet('') + #par = self.thr_parameter + #Need to clean the figure to prepare is for a different type of plot + for i in range(6): + self.axesMap[i//2,i%2].cla() + + #Draw horizontal lines - out calculated thresholds + self.axhlines = [None for _ in range(12)] + num_temp = 0 + cals = [self.f1_cal, self.f2_cal, self.f3_cal, self.f4_cal, self.f5_cal, self.f6_cal] + print(len(self.f1_cal)) + #remove outliers! + if self.remove_outliers: + for cal_idx, cal in enumerate(cals): + cals[cal_idx] = np.array(cal)[~doubleMADsfromMedian(cal)] + self.f1_cal, self.f2_cal, self.f3_cal, self.f4_cal, self.f5_cal, self.f6_cal = cals + print(len(self.f1_cal)) + + self.thr_1 = [np.min(self.f1_cal)-0.1*(np.max(self.f1_cal) - np.min(self.f1_cal)), + np.max(self.f1_cal)+0.1*(np.max(self.f1_cal) - np.min(self.f1_cal))] + self.thr_2 = [np.min(self.f2_cal)-0.1*0.1*(np.max(self.f2_cal) - np.min(self.f2_cal)), + np.max(self.f2_cal)+0.1*0.1*(np.max(self.f2_cal) - np.min(self.f2_cal))] + self.thr_3 = [np.min(self.f3_cal)-0.1*0.1*(np.max(self.f3_cal) - np.min(self.f3_cal)), + np.max(self.f3_cal)+0.1*0.1*(np.max(self.f3_cal) - np.min(self.f3_cal))] + self.thr_4 = [np.min(self.f4_cal)-0.1*0.1*(np.max(self.f4_cal) - np.min(self.f4_cal)), + np.max(self.f4_cal)+0.1*0.1*(np.max(self.f4_cal) - np.min(self.f4_cal))] + self.thr_5 = [np.min(self.f5_cal)-0.1*0.1*(np.max(self.f5_cal) - np.min(self.f5_cal)), + np.max(self.f5_cal)+0.1*0.1*(np.max(self.f5_cal) - np.min(self.f5_cal))] + self.thr_6 = [np.min(self.f6_cal)-0.1*0.1*(np.max(self.f6_cal) - np.min(self.f6_cal)), + np.max(self.f6_cal)+0.1*0.1*(np.max(self.f6_cal) - np.min(self.f6_cal))] + self.thrs = [self.thr_1, self.thr_2, self.thr_3, self.thr_4, self.thr_5, self.thr_6] + + #put axhlines + for i in range(6): + self.axhlines[num_temp] = self.axesMap[i//2,i%2].axhline(self.thrs[i][0], color = 'r') + self.axhlines[num_temp+1] = self.axesMap[i//2,i%2].axhline(self.thrs[i][1], color = 'r') + num_temp+=2 + + + #Auto_readout_lim - automatic adjustment of the figure limits: + if self.auto_readout_lim: + for i in range(6): + if not all(np.isnan(self.thrs[i])): + self.axesMap[i//2,i%2].set_ylim(min_zero(np.array(self.thrs[i]) + + np.diff(np.array(self.thrs[i]))*np.array([-2,2]))) + + else: + for i in range(6): + if not all(np.isnan(self.thrs[i])): + self.axesMap[i//2,i%2].set_ylim(self.measures_x_lims[i,0],self.measures_x_lims[i,1]) + + + #set xlim, titles and ylabels + for i in range(6): + if not all(np.isnan(self.thrs[i])): + self.axesMap[i//2,i%2].set_xlim(0.5,self.time_between_bursts-self.breaktime+0.5) + self.axesMap[i//2,i%2].set_title(self.Titles[i], fontsize=9) + self.axesMap[i//2,i%2].set_ylabel(self.unit_label[i], fontsize=self.labels_size, rotation=270) + self.axesMap[i//2,i%2].set_xlabel("Time [s]", fontsize=self.labels_size) + + texts = [] + #print on the figure whgat was the max of each feature + for ind,val in enumerate([np.max(self.f1_cal), np.max(self.f2_cal), np.max(self.f3_cal), + np.max(self.f4_cal), np.max(self.f5_cal), np.max(self.f6_cal)]): + if not all(np.isnan(self.thrs[ind])): + texts.append(self.axesMap[ind//2,ind%2].text(0.8,0.9, round(val,2), + transform=self.axesMap[ind//2,ind%2].transAxes, color='r', fontsize=7)) + + self.canvasMap.draw() + + def settings(self): + self.options = Ui(self) + + def update_plot_lim(self): + if not self.do_calibration: + if self.auto_readout_lim: + for i in range(6): + self.axesMap[i//2,i%2].set_ylim(min_zero(np.array(self.thrs[i]) + + np.diff(np.array(self.thrs[i]))*np.array([-2,2]))) + # self.axesMap[0,0].set_ylim(-1,1) + # self.axesMap[0,1].set_ylim(min_zero(np.array(self.thr_2) + + # np.diff(np.array(self.thr_2))*np.array([-2,2]))) + # self.axesMap[1,0].set_ylim(min_zero(np.array(self.thr_1) + + # np.diff(np.array(self.thr_1))*np.array([-2,2]))) + # self.axesMap[1,1].set_ylim(min_zero(np.array(self.thr_3) + + # np.diff(np.array(self.thr_3))*np.array([-2,2]))) + # self.axesMap[2,0].set_ylim(min_zero(np.array(self.thr_5) + + # np.diff(np.array(self.thr_5))*np.array([-2,2]))) + # self.axesMap[2,1].set_ylim(min_zero(np.array(self.thr_6) + + # np.diff(np.array(self.thr_6))*np.array([-2,2]))) + else: + for i in range(6): + self.axesMap[i//2,i%2].set_ylim(self.measures_x_lims[i,0],self.measures_x_lims[i,1]) + # self.axesMap[0,0].set_ylim(self.measures_x_lims[0,0],self.measures_x_lims[0,1]) + # self.axesMap[0,1].set_ylim(self.measures_x_lims[1,0],self.measures_x_lims[1,1]) + # self.axesMap[1,0].set_ylim(self.measures_x_lims[2,0],self.measures_x_lims[2,1]) + # self.axesMap[1,1].set_ylim(self.measures_x_lims[3,0],self.measures_x_lims[3,1]) + # self.axesMap[2,0].set_ylim(self.measures_x_lims[4,0],self.measures_x_lims[4,1]) + # self.axesMap[2,1].set_ylim(self.measures_x_lims[5,0],self.measures_x_lims[5,1]) + else: + for i in range(6): + self.axesMap[i//2,i%2].set_xlim(self.measures_x_lims[i,0],self.measures_x_lims[i,1]) + # self.axesMap[0,0].set_xlim(self.measures_x_lims[0,0],self.measures_x_lims[0,1]) + # self.axesMap[0,1].set_xlim(self.measures_x_lims[1,0],self.measures_x_lims[1,1]) + # self.axesMap[1,0].set_xlim(self.measures_x_lims[2,0],self.measures_x_lims[2,1]) + # self.axesMap[1,1].set_xlim(self.measures_x_lims[3,0],self.measures_x_lims[3,1]) + # self.axesMap[2,0].set_xlim(self.measures_x_lims[4,0],self.measures_x_lims[4,1]) + # self.axesMap[2,1].set_xlim(self.measures_x_lims[5,0],self.measures_x_lims[5,1]) + self.canvasMap.draw() + + def calibration_process(self): + """Calibrate thresholds values""" + #same things as in update_methods + time_start = time.time() + self.second_to_analyze = self.data_processing(True) + + #calculate fft + self.S = abs(np.fft.rfft(self.second_to_analyze)) + + #this is a bit shady, but should work. check it out if doesn't! + #do dwt only if any features requires it + if any(feature_num in self.used_features for feature_num in [8,9,10,11]): + self.dwt = pywt.wavedec(self.second_to_analyze, 'db1', level=8) + + self.results = np.zeros(len(self.used_features)) + + #print(self.used_features) + for idx, feature in enumerate(self.used_features): + if feature==0: + self.results[idx] = None + continue + elif feature in [8,9,10,11]: #to specify which band + self.results[idx] = self.functions[feature](feature) + else: + self.results[idx] = self.functions[feature]() + + self.log_file_writer.writerow([time.time() - self.time_start, 'Calibration', + self.results]) + self.log_file.flush() + + #add mean to the list + self.f1_cal.append(self.results[0]) + self.f2_cal.append(self.results[1]) + self.f3_cal.append(self.results[2]) + self.f4_cal.append(self.results[3]) + self.f5_cal.append(self.results[4]) + self.f6_cal.append(self.results[5]) + + self.cals = [self.f1_cal, self.f2_cal, self.f3_cal, self.f4_cal, self.f5_cal, + self.f6_cal] + + #!!! Temporary, it's wrong but it's overwritten later. It's needed to check if all features are used, + #but there could be more optimal solution. Remove it at some point! + self.thr_1 = [np.min(self.f1_cal)-0.1*np.min(self.f1_cal), + np.max(self.f1_cal)+0.1*np.max(self.f1_cal)] + self.thr_2 = [np.min(self.f2_cal)-0.1*np.min(self.f2_cal), + np.max(self.f2_cal)+0.1*np.max(self.f2_cal)] + self.thr_3 = [np.min(self.f3_cal)-0.1*np.min(self.f3_cal), + np.max(self.f3_cal)+0.1*np.max(self.f3_cal)] + self.thr_4 = [np.min(self.f4_cal)-0.1*np.min(self.f4_cal), + np.max(self.f4_cal)+0.1*np.max(self.f4_cal)] + self.thr_5 = [np.min(self.f5_cal)-0.1*np.min(self.f5_cal), + np.max(self.f5_cal)+0.1*np.max(self.f5_cal)] + self.thr_6 = [np.min(self.f6_cal)-0.1*np.min(self.f6_cal), + np.max(self.f6_cal)+0.1*np.max(self.f6_cal)] + self.thrs = [self.thr_1, self.thr_2, self.thr_3, self.thr_4, self.thr_5, self.thr_6] + + times1 = time.time() + print("Thr features calculation: {}".format(times1-time_start)) + num_bins=35 + if self.auto_readout_lim: + if len(self.f4_cal)==1: + for line_idx in range(len(self.linetemps)): + if not all(np.isnan(self.thrs[line_idx])): + y,x = np.histogram(self.cals[line_idx], num_bins, self.thrs[line_idx]) + self.linetemps[line_idx] = self.axesMap[line_idx//2,line_idx%2].plot(adjust_hist(x), y, '.-') + + else: + for line_idx in range(len(self.linetemps)): + if not all(np.isnan(self.thrs[line_idx])): + data = np.histogram(self.cals[line_idx], num_bins, self.thrs[line_idx]) + update_stem(self.linetemps[line_idx], data, self.axesMap[line_idx//2,line_idx%2]) + + for i in range(6): + if not all(np.isnan(self.thrs[i])): + self.axesMap[i//2,i%2].set_xlim(self.thrs[i][0], self.thrs[i][1]) + # self.axesMap[0,0].set_xlim(self.thr_1[0], self.thr_1[1]) + # self.axesMap[0,1].set_xlim(self.thr_2[0], self.thr_2[1]) + # self.axesMap[1,0].set_xlim(self.thr_3[0], self.thr_3[1]) + # self.axesMap[1,1].set_xlim(self.thr_4[0], self.thr_4[1]) + # self.axesMap[2,0].set_xlim(self.thr_5[0], self.thr_5[1]) + # self.axesMap[2,1].set_xlim(self.thr_6[0], self.thr_6[1]) + + elif len(self.f4_cal)==1: + for line_idx in range(len(self.linetemps)): + if not all(np.isnan(self.thrs[line_idx])): + y,x = np.histogram(self.cals[line_idx], num_bins, self.measures_x_lims[line_idx]) + self.linetemps[line_idx] = self.axesMap[line_idx//2,line_idx%2].plot(adjust_hist(x), y, '.-') + + else: + for line_idx in range(len(self.linetemps)): + if not all(np.isnan(self.thrs[line_idx])): + data = np.histogram(self.cals[line_idx], num_bins, self.thrs[line_idx]) + update_stem(self.linetemps[line_idx], data, self.axesMap[line_idx//2,line_idx%2]) + for i in range(6): + if not all(np.isnan(self.thrs[i])): + self.axesMap[i//2,i%2].set_xlim(self.measures_x_lims[i,0],self.measures_x_lims[i,1]) + # self.axesMap[0,0].set_xlim(self.measures_x_lims[0,0],self.measures_x_lims[0,1]) + # self.axesMap[0,1].set_xlim(self.measures_x_lims[1,0],self.measures_x_lims[1,1]) + # self.axesMap[1,0].set_xlim(self.measures_x_lims[2,0],self.measures_x_lims[2,1]) + # self.axesMap[1,1].set_xlim(self.measures_x_lims[3,0],self.measures_x_lims[3,1]) + # self.axesMap[2,0].set_xlim(self.measures_x_lims[4,0],self.measures_x_lims[4,1]) + # self.axesMap[2,1].set_xlim(self.measures_x_lims[5,0],self.measures_x_lims[5,1]) + + for i in range(6): + if not all(np.isnan(self.thrs[i])): + self.axesMap[i//2,i%2].set_ylim(self.measures_y_lims[i,0],self.measures_y_lims[i,1]) + # self.axesMap[0,0].set_ylim(self.measures_y_lims[0,0],self.measures_y_lims[0,1]) + # self.axesMap[0,1].set_ylim(self.measures_y_lims[1,0],self.measures_y_lims[1,1]) + # self.axesMap[1,0].set_ylim(self.measures_y_lims[2,0],self.measures_y_lims[2,1]) + # self.axesMap[1,1].set_ylim(self.measures_y_lims[3,0],self.measures_y_lims[3,1]) + # self.axesMap[2,0].set_ylim(self.measures_y_lims[4,0],self.measures_y_lims[4,1]) + # self.axesMap[2,1].set_ylim(self.measures_y_lims[5,0],self.measures_y_lims[5,1]) + + # self.linetemps = [self.linetemp1, self.linetemp2, self.linetemp3, self.linetemp4, + # self.linetemp5, self.linetemp6] + + # for i in range(3): + # for j in range(2): + # self.axesMap[i,j].draw_artist(self.linetemps[i*2+j][0]) + # if len(self.f1_cal)>0: + # self.canvasMap.update() + # self.canvasMap.flush_events() + # else: + self.canvasMap.draw() #update canvas + print('Thr histogram plots', time.time()-time_start) + + def zoom_out(self): + self.plot_dividing_factor = self.plot_dividing_factor*2 + + def zoom_in(self): + self.plot_dividing_factor = self.plot_dividing_factor*0.5 + + def create_main_frame(self): + """Prepare figure, buttons, plots, etc...""" + #set colors using for plotting, I THINK I DONT NEED THAT ANYMORE + mpl.rcParams['axes.prop_cycle'] = cycler(color=cm.nipy_spectral( + np.linspace(0,0.2,6))) + self.main_frame = QWidget() + self.dpi = 100 + self.fig = Figure((11, 10), dpi=self.dpi)#, facecolor='whitesmoke') #figure for signal plot + + self.figMap = Figure((8, 10), dpi=self.dpi) #figure for measurments + self.axes = self.fig.subplots(1,1) + #self.axes.set_position([0, 0, 1, 1]) + self.line = [None for _ in range(self.num_of_ch)] #list for lines + #self.axes.set_axis_off() #turn all axes off + for tick in self.axes.xaxis.get_major_ticks(): + tick.tick1line.set_visible(False) + tick.tick2line.set_visible(False) + tick.label1.set_visible(False) + tick.label2.set_visible(False) + + # for tick in self.axes.yaxis.get_major_ticks(): + # tick.tick1line.set_visible(False) + # tick.tick2line.set_visible(False) + # tick.label1.set_visible(False) + # tick.label2.set_visible(False) + + self.axes.set_yticks(np.arange(1, (self.num_of_ch)*1.01, 1)) + self.axes.set_yticklabels(self.ch_names[::-1]) + + self.axes.set_xticks(np.arange(self.Fs,self.plot_len*self.Fs, self.Fs)) + self.axes.grid(True) + + self.canvas = FigureCanvas(self.fig) + self.canvas.setParent(self.main_frame) + self.canvasMap = FigureCanvas(self.figMap) # + self.canvasMap.setParent(self.main_frame) + self.axesMap = self.figMap.subplots(3,2) + + self.Titles = [self.features_names[idx] for idx in self.used_features] + + for i in range(6): + self.axesMap[i//2,i%2].set_title(self.Titles[i], fontsize=9) + self.axesMap[i//2,i%2].set_ylabel(self.unit_label[i], fontsize=self.labels_size, rotation=270) + + # self.axesMap[0,0].set_title(self.Titles[0], fontsize=9) + # self.axesMap[0,1].set_title(self.Titles[1], fontsize=9) + # self.axesMap[1,0].set_title(self.Titles[2], fontsize=9) + # self.axesMap[1,1].set_title(self.Titles[3], fontsize=9) + # self.axesMap[2,0].set_title(self.Titles[4], fontsize=9) + # self.axesMap[2,1].set_title(self.Titles[5], fontsize=9) + + # self.axesMap[0,0].set_ylabel(self.unit_label[0], fontsize=self.labels_size, rotation=270) + # self.axesMap[0,1].set_ylabel(self.unit_label[1], fontsize=self.labels_size, rotation=270) + # self.axesMap[1,0].set_ylabel(self.unit_label[2], fontsize=self.labels_size, rotation=270) + # self.axesMap[1,1].set_ylabel(self.unit_label[3], fontsize=self.labels_size, rotation=270) + # self.axesMap[2,0].set_ylabel(self.unit_label[4], fontsize=self.labels_size, rotation=270) + # self.axesMap[2,1].set_ylabel(self.unit_label[5], fontsize=self.labels_size, rotation=270) + + label_size = 7 + for i in range(3): + for j in range(2): + self.axesMap[i,j].set_xlabel("Time [s]", fontsize=self.labels_size) + self.axesMap[i,j].yaxis.set_label_coords(self.lab_x_dist, self.lab_y_dist) + self.axesMap[i,j].xaxis.set_label_coords(self.lab_y_dist, self.lab_x_dist) + self.axesMap[i,j].xaxis.set_tick_params(labelsize = label_size) + self.axesMap[i,j].yaxis.set_tick_params(labelsize = label_size) + + #using list with lines we plot all channels + for i in range(self.num_of_ch): + if i in self.included_ch or self.ch_names[i] in ["EMG", "EOG"]: + self.line[i], = self.axes.plot([] , color = 'black', linewidth=0.4) + else: + self.line[i], = self.axes.plot([] , color = 'silver', linewidth=0.3) + self.axes.set_xlim(0, self.plot_len*self.Fs) + self.axes.set_ylim(0, (self.num_of_ch+1)*1) + #self.axes.axvspan((self.plot_len-1)*self.Fs, + # self.plot_len*self.size_of_up, alpha=0.3, color='lightcoral') + self.fig.tight_layout(pad=1, rect=(0.02,-0.02,1.02,1.02)) #tighten the layout by removing free spaces + self.figMap.tight_layout() + # self.fig.subplots_adjust(left=0.02,right=1, + # bottom=0,top=1, + # hspace=0,wspace=0) + + self.axes.spines['left'].set_visible(True) + self.axes.spines['bottom'].set_visible(False) + self.axes.spines['top'].set_visible(False) + self.axes.spines['right'].set_visible(False) + + self.canvas.draw() + self.canvasMap.draw() #update canvas + #NOT NEEDED + for i in range(self.num_of_ch): + self.axbackground = self.canvas.copy_from_bbox(self.axes.bbox) + + # texts = [] + # for ind,name in enumerate(self.ch_names): + # texts.append(self.axes.text(-0.03,1-(ind+1)/(len(self.ch_names)+1), + # name, transform=self.axes.transAxes, color='r', fontsize=13)) + + self.button1 = QPushButton("&Settings") + self.button1.setCheckable(False) + self.button1.clicked.connect(self.settings) + + self.button2 = QPushButton("&Start/Stop") + self.button2.setCheckable(True) + self.button2.clicked.connect(self.start_stop) + + self.button4 = QPushButton("&Start calibration") + self.button4.setCheckable(False) + self.button4.clicked.connect(self.calibration) + + self.button_zoom_in = QPushButton("+") + self.button_zoom_in.setCheckable(False) + self.button_zoom_in.clicked.connect(self.zoom_in) + self.button_zoom_in.setMaximumWidth(20) + + self.button_zoom_out = QPushButton("-") + self.button_zoom_out.setCheckable(False) + self.button_zoom_out.clicked.connect(self.zoom_out) + self.button_zoom_out.setMaximumWidth(20) + + hbox_zoom = QHBoxLayout() + hbox_zoom.addWidget(self.button_zoom_in, Qt.AlignLeft) + hbox_zoom.addWidget(self.button_zoom_out, Qt.AlignLeft) + hbox_zoom.addStretch(1) + # hbox_zoom.setAlignment(self.button_zoom_in, Qt.AlignLeft) + # hbox_zoom.setAlignment(self.button_zoom_out, Qt.AlignLeft) + + #set layout of buttons + hbox = QHBoxLayout() + for w in [self.button1, self.button2, #self.button3, + self.button4]: + hbox.addWidget(w) + hbox.setAlignment(w, Qt.AlignVCenter) + #set genera layout + layout = QGridLayout() + layout.addWidget(self.canvas, 0, 0, 1, 1) + layout.addWidget(self.canvasMap, 0, 1, 1, 1) + layout.addLayout(hbox, 1, 1, 1, 1) + layout.addLayout(hbox_zoom, 1,0,1,1) + + + self.main_frame.setLayout(layout) + self.setCentralWidget(self.main_frame) + +class Ui(QMainWindow): + """""" + def __init__(self, main_file): + self.main_file = main_file + super(Ui, self).__init__() + uic.loadUi('soft2.ui', self) + self.show() + self.load_button = self.findChild(QPushButton, 'pushButton') + self.load_button.clicked.connect(self.get_file) + self.auto_ylim = self.findChild(QCheckBox, 'checkBox') + self.auto_ylim.setChecked(self.main_file.auto_readout_lim) + self.disp_max_box = self.findChild(QCheckBox, 'checkBox_2') + self.disp_max_box.setChecked(self.main_file.disp_max) + self.apply_lims_button = self.findChild(QPushButton, 'pushButton_2') + self.apply_lims_button.clicked.connect(self.apply_lims) + self.apply_restart = self.findChild(QPushButton, 'pushButton_3') + self.apply_restart.clicked.connect(lambda: self.main_file.all_params('restarted')) + self.file_path = self.findChild(QLabel, 'label_12') + names_xlim = ['lineEdit_2', 'lineEdit_3', 'lineEdit_4', 'lineEdit_5', 'lineEdit_6', 'lineEdit_7'] + names_ylim = ['lineEdit_9', 'lineEdit_13', 'lineEdit_10', 'lineEdit_11', 'lineEdit_8', 'lineEdit_12'] + self.xlim_values = np.zeros(6, dtype=object) + self.ylim_values = np.zeros(6, dtype=object) + for i in range(6): + self.xlim_values[i] = self.findChild(QLineEdit, names_xlim[i]) + self.xlim_values[i].setText(str(list(self.main_file.measures_x_lims[i]))) + self.ylim_values[i] = self.findChild(QLineEdit, names_ylim[i]) + self.ylim_values[i].setText(str(list(self.main_file.measures_y_lims[i]))) + + def get_file(self): + self.path = QFileDialog.getOpenFileName(self, 'Open File', os.path.dirname(os.getcwd()), 'Text files (*.txt *.csv)') + print(self.path) + self.file_path.setText(self.path[0]) + self.main_file.montage_matrix = np.array(pd.read_csv(self.path[0], header=None)) + self.main_file.montage_file_path = self.path[0] + + def apply_lims(self): + for i in range(6): + try: + self.main_file.measures_x_lims[i] = json.loads(self.xlim_values[i].text()) + self.main_file.measures_y_lims[i] = json.loads(self.ylim_values[i].text()) + self.main_file.auto_readout_lim = self.auto_ylim.isChecked() + self.main_file.disp_max = self.disp_max_box.isChecked() + except: + print('Fix line:', self.xlim_values[i].text()) + print(self.main_file.measures_x_lims) + print(self.main_file.measures_y_lims) + + self.main_file.update_plot_lim() + + + +if __name__ == '__main__': + app = QApplication(sys.argv) + form = First_window(AppForm) #AppForm() + form.show() + app.exec_() + +# cut time different from both sides +# some deafult settings. Maybe remember last configuration? +# EMG and EOG - none, more than one? +# change names to final names \ No newline at end of file diff --git a/Electrode_selection.txt b/Electrode_selection.txt new file mode 100644 index 0000000..e49df4d --- /dev/null +++ b/Electrode_selection.txt @@ -0,0 +1,2 @@ +full_cap_file_path: /home/adamr/Documents/PYTHON/TMS TV/Up to date version Feb 23/easycap-M10_63_NO.txt +cap_file_path: /home/adamr/Documents/PYTHON/TMS TV/Up to date version Feb 23/easycap-M10_16_NO.txt \ No newline at end of file diff --git a/FirstWindow.py b/FirstWindow.py new file mode 100644 index 0000000..73f556b --- /dev/null +++ b/FirstWindow.py @@ -0,0 +1,992 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Jan 24 11:26:23 2022 + +@author: Basics +""" +import mne +import os, sys, traceback, time +from PyQt5.QtWidgets import (QMainWindow, QFileDialog, QMessageBox, QCheckBox, QLineEdit, QWidget, QPushButton, + QLabel, QHBoxLayout, QGridLayout, QAction, QApplication, QDialog, QDialogButtonBox, + QVBoxLayout, QFrame, QTabWidget, QComboBox, QScrollArea) +from PyQt5.QtCore import QTimer, Qt +from PyQt5 import QtCore +from matplotlib.backends.backend_qt5agg import FigureCanvasQTAgg as FigureCanvas +from matplotlib.backends.backend_qt5agg import NavigationToolbar2QT as NavigationToolbar +from matplotlib.figure import Figure +from matplotlib.widgets import Button as MatplotlibButton +import multiprocessing as mp +if sys.platform=='darwin': + from multiprocessing import Process, Value + #from multiprocessing import Queue as StupidNotWorkingQueue +else: + from multiprocessing import Process, Queue, Value +import NeurOne +import RDA +import ctypes + +import pandas as pd +from mne.channels.layout import _find_topomap_coords as get_pos +import json +import numpy as np + +if sys.platform=='darwin': + from multiprocessing.queues import Queue as QueueOld + class SharedCounter(object): + #source: https://gist.github.com/FanchenBao/d8577599c46eab1238a81857bb7277c9 + """A synchronized shared counter. + The locking done by multiprocessing.Value ensures that only a single + process or thread may read or write the in-memory ctypes object. However, + in order to do n += 1, Python performs a read followed by a write, so a + second process may read the old value before the new one is written by the + first process. The solution is to use a multiprocessing.Lock to guarantee + the atomicity of the modifications to Value. + This class comes almost entirely from Eli Bendersky's blog: + http://eli.thegreenplace.net/2012/01/04/shared-counter-with-pythons-multiprocessing/ + """ + + def __init__(self, n=0): + self.count = mp.Value('i', n) + + def increment(self, n=1): + """ Increment the counter by n (default = 1) """ + with self.count.get_lock(): + self.count.value += n + + @property + def value(self): + """ Return the value of the counter """ + return self.count.value + + + class Queue(QueueOld): + #source: https://gist.github.com/FanchenBao/d8577599c46eab1238a81857bb7277c9 + """ A portable implementation of multiprocessing.Queue. + Because of multithreading / multiprocessing semantics, Queue.qsize() may + raise the NotImplementedError exception on Unix platforms like Mac OS X + where sem_getvalue() is not implemented. This subclass addresses this + problem by using a synchronized shared counter (initialized to zero) and + increasing / decreasing its value every time the put() and get() methods + are called, respectively. This not only prevents NotImplementedError from + being raised, but also allows us to implement a reliable version of both + qsize() and empty(). + Note the implementation of __getstate__ and __setstate__ which help to + serialize MyQueue when it is passed between processes. If these functions + are not defined, MyQueue cannot be serialized, which will lead to the error + of "AttributeError: 'MyQueue' object has no attribute 'size'". + See the answer provided here: https://stackoverflow.com/a/65513291/9723036 + + For documentation of using __getstate__ and __setstate__ to serialize objects, + refer to here: https://docs.python.org/3/library/pickle.html#pickling-class-instances + """ + + def __init__(self): + super().__init__(ctx=mp.get_context()) + self.size = SharedCounter(0) + + def __getstate__(self): + """Help to make MyQueue instance serializable. + Note that we record the parent class state, which is the state of the + actual queue, and the size of the queue, which is the state of MyQueue. + self.size is a SharedCounter instance. It is itself serializable. + """ + return { + 'parent_state': super().__getstate__(), + 'size': self.size, + } + + def __setstate__(self, state): + super().__setstate__(state['parent_state']) + self.size = state['size'] + + def put(self, *args, **kwargs): + super().put(*args, **kwargs) + self.size.increment(1) + + def get(self, *args, **kwargs): + item = super().get(*args, **kwargs) + self.size.increment(-1) + return item + + def qsize(self): + """ Reliable implementation of multiprocessing.Queue.qsize() """ + return self.size.value + + def empty(self): + """ Reliable implementation of multiprocessing.Queue.empty() """ + return not self.qsize() + + +class NeurOneOffline(): + def __init__(self): + self.data = None + + def NO(self,ip,port=50000,buffersize=2**15,ringbuffersize=2000, + sendqueue=False,ringbuf_factor=2,dump=False,avgPackets=1): + self.ringbuffersize = ringbuffersize + tmp_path = '/mnt/projects/P_BCT_EEG/DLPFCM1_iTBS/DLPFC/beep/subj_14/X13193_adam.vhdr' + num_electr = 18 + eeg_chn = np.arange(0,num_electr,1) + hdr = mne.io.read_raw_brainvision(tmp_path) + # hdr.set_channel_types({'EMGleft': 'emg', 'EOGright': 'eog'}) + # hdr.set_montage(mne.channels.read_custom_montage('easycap-M10_63_NO.txt')) + mrk_fullpath = tmp_path[:-4]+'vmrk' + eeg_fullpath = tmp_path[:-4]+'eeg' #this two are made by hand instead of function. + #Maybe there is some func for this + #annotations returns all events - stimA, stimB, stopA, stopB, start of experiment etc... We chose only stim + stim = hdr.annotations.onset[np.logical_or(hdr.annotations.description=="Stimulus/A", + hdr.annotations.description=="Stimulus/B")] + #and here we separate stimA and stimB + stimA = hdr.annotations.onset[hdr.annotations.description=="Stimulus/A"] + stimB = hdr.annotations.onset[hdr.annotations.description=="Stimulus/B"] + #divide for stim A and B + #stimR = hdr.annotations.onset[hdr.annotations.description=='Response/R 16'] + npts = hdr.n_times + nfft = int(hdr.info['sfreq']) # Sampling rate [Hz] + fs = int(hdr.info['sfreq']) # Sampling rate [Hz] + endsample = npts + begsample = 0 + print('Sampling rate [Hz]: ' + str(fs)) + #[stim, resp, segment, timezero] = mne.read_annotations(mrk_fullpath) thinkkkkk + eeg_raw = hdr.get_data(start=begsample, stop=endsample); + eeg = eeg_raw[eeg_chn, :] # Select all electrodes + + data_raw = (hdr.get_data()*1e7) + data_raw = data_raw[:,:-5]#int(data_raw.shape[1]%5)] + stimB_arr = (stimB*fs).astype(int) + data_raw_new = np.zeros([data_raw.shape[0], int(data_raw.shape[1]//5)]) + for i in range(data_raw.shape[0]): + a = data_raw[i] + R = 5 + data_raw_new[i,:] = a.reshape(-1, R).mean(1) + data_raw = data_raw_new + fs=1000 + stimB_arr = (stimB_arr/5).astype(int) + stim_sig = np.zeros(data_raw.shape[1]).reshape([1,data_raw.shape[1]]) + for i in stimB_arr: + stim_sig[0,i] = 1 + self.fs = fs + self.data = np.concatenate((data_raw, stim_sig)) + self.buffer = self.data[:,:ringbuffersize] + + def start(self): + self.time = time.time() + + def getBuffer(self): + time_now = time.time() + time_diff = int((time_now - self.time)*self.fs) + return self.data[:,time_diff-self.ringbuffersize:time_diff] + +def acquire_data(q, size, run, speed, downsample, sleep_time, ip = '192.168.200.201', port = 5000, offline=False): + """Acquires data from NeurOne_v3 and pass it to the queue. Function is supposed to work + in separated process. + + Parameters + ------------ + q: Queue class from multiprocessing library + size: number of samples to acquire + run: Value class from multiprocessing library. That value can be changed in main process + downsample: boolean value. Says if data will be downsampled to 1000 Hz + sleep_time: int, set how often function should refresh. Usually it takes a bit more that that""" + #import NeurOne_v3 + #channel_list = [1,2,3,4,5,6,7,8,9,10,11,12,13,14,15, 16,-1] #just a temp solution to have less channels + if offline=="offline": + NO = NeurOneOffline() + NO.NO(ip=ip,port=50000,buffersize=2**15,ringbuffersize=size,\ + sendqueue=False,ringbuf_factor=2,dump=False,avgPackets=downsample) #offline + elif offline=="NeurOne": + # ip='192.168.200.201' + # port = 50000 + NO = NeurOne.NO(ip=ip,port=port,buffersize=2**15,ringbuffersize=size,\ + sendqueue=False,ringbuf_factor=2, dump=None,avgPackets=downsample) #online + elif offline=="BrainProducts": + # ip = "169.254.252.66" + # port = 51244 + NO = RDA.RDA(ip=ip, port=port, buffersize=2**15, ringbuffersize=size,\ + sendqueue=False, avgPackets=downsample) + NO.start() + while True: + startt = time.time() + if run.value: + data = NO.getBuffer()#[channel_list] + try: + q.put(data) + except: print('ERROR') + endd=time.time() + time.sleep(sleep_time*speed.value-(endd-startt)) + + +class First_window(QMainWindow): + """This window opens at the beggining before the main program. It provides + a GUI to choose a lot of parameters without modyfying any files.""" + def __init__(self, AppForm): + self.AppForm = AppForm + super(First_window, self).__init__() + self.setWindowTitle('EStiMo Configuration') + + try: + cap_loc_file = pd.read_csv('Electrode_selection.txt', sep=':', header=None) + self.cap_file_path = cap_loc_file[cap_loc_file[0]=='cap_file_path'].values[0][1].strip() + self.full_cap_file_path = cap_loc_file[cap_loc_file[0]=='full_cap_file_path'].values[0][1].strip() #'easycap-M10_63_NO.txt' + except: + print("CAP FILE EXCEPTION") + self.cap_file_path = 'easycap-M10_16_NO.txt' + self.full_cap_file_path = 'easycap-M10_63_NO.txt' + + self.conf_path = 'TMS_protocol.txt' + montage = mne.channels.read_custom_montage(self.cap_file_path) + montage_file = pd.read_csv(self.cap_file_path, sep='\t') + + montage_full = mne.channels.read_custom_montage(self.full_cap_file_path) + montage_full_file = pd.read_csv(self.full_cap_file_path, sep='\t') + + self.ch_names_loaded = montage_file.to_numpy()[:,0] + self.ch_names_loaded_full = montage_full_file.to_numpy()[:,0] + self.montage_file_path = None + self.all_names = montage_file.iloc[:,0] + self.full_all_names = montage_full_file.iloc[:,0] + ch_info = mne.create_info(montage.ch_names, sfreq=1000, ch_types='eeg') + ch_info.set_montage(montage) #2d projection of montage (positions) + self.pos2d = get_pos(ch_info, None) + ch_info_full = mne.create_info(montage_full.ch_names, sfreq=1000, ch_types='eeg') + ch_info_full.set_montage(montage_full) #2d projection of montage (positions) + self.pos2d_full = get_pos(ch_info_full, None) + self.set_values() + self.create_figure() + + + + def set_values(self, update_plot=False): + """ + Parameters + ---------- + update_plot : bool, optional + Update plot if needed. The default is False. + + Returns + ------- + bool + returns False if error occured. + + """ + #Try to load all data from the file. + self.ip = '192.168.200.201' + self.port = 50000 + try: + self.included_ch = [] + self.included_ch_full = [] + settings_file = pd.read_csv(self.conf_path,sep=':', header=None) + temp_names = settings_file[settings_file[0]=='names'].values[0][1].strip() + self.num_of_ch_loaded = int(settings_file[settings_file[0]=='number_of_channels'].values[0][1]) + self.num_of_lines_loaded = int(settings_file[settings_file[0]=='number_of_lines'].values[0][1]) + self.time_between_bursts_loaded = int(settings_file[settings_file[0]=='time_between_trains'].values[0][1]) + self.breaktime_loaded = int(settings_file[settings_file[0]=='cut_time'].values[0][1]) + self.included_ch_loaded = json.loads(settings_file[settings_file[0]=='included_channels'].values[0][1]) + self.eog_ch_loaded = int(settings_file[settings_file[0]=='eog_channel'].values[0][1]) + self.emg_ch_loaded = int(settings_file[settings_file[0]=='emg_channel'].values[0][1]) + self.exp_trig_loaded = int(settings_file[settings_file[0]=='expected_triggers'].values[0][1]) + self.exp_time_loaded = int(settings_file[settings_file[0]=='expected_time'].values[0][1]) + self.BoxChecked = False + except Exception as e: + ex_type, ex_value, ex_traceback = sys.exc_info() + + # Extract unformatter stack traces as tuples + trace_back = traceback.extract_tb(ex_traceback) + + # Format stacktrace + stack_trace = list() + + for trace in trace_back: + stack_trace.append("File : %s , Line : %d, Func.Name : %s, Message : %s" % (trace[0], trace[1], trace[2], trace[3])) + self.text_conf_file.setText("No configuration file selected") + print('exception') + + ### Shows error message if error occured. Does not break the program + self.msgBox = QMessageBox(self) + self.msgBox.setIcon(QMessageBox.Critical) + self.msgBox.setDetailedText("Detailed error: \nException type: {} \nMessage: {} \nTrace: {}".format(ex_type.__name__,ex_value, stack_trace)) + self.msgBox.setText("Error occured while loading the file. Check the structure of the file and try again.") + self.msgBox.setWindowTitle("Protocol file error") + self.msgBox.setStandardButtons(QMessageBox.Ok) + self.msgBox.show() + return False + + #some code-words for default settings + if temp_names not in ['Def', 'def', 'None', 'none', 'default', 'Default']: + self.ch_names_loaded = json.loads(temp_names) + #else use settings provided + else: + if self.eog_ch_loaded < self.emg_ch_loaded: + self.ch_names_loaded = np.append(self.ch_names_loaded, ["EOG", "EMG"]) + self.ch_names_loaded_full = np.append(self.ch_names_loaded_full, ["EOG", "EMG"]) + else: + self.ch_names_loaded = np.append(self.ch_names_loaded, ["EMG", "EOG"]) + self.ch_names_loaded_full = np.append(self.ch_names_loaded_full, ["EMG", "EOG"]) + + if update_plot: + self.color = np.zeros([len(self.pos2d[:,0]),3]) + for idx in self.ch_names_loaded: + num = np.where(self.all_names.astype(str)==idx) + self.color[num] = [0,0,1] #int(self.color[event.ind]!=1) + self.ddd.set_color(self.color) + self.canvas.draw() + lista = [self.time_between_bursts_loaded, self.breaktime_loaded, self.num_of_ch_loaded, + self.num_of_lines_loaded, self.eog_ch_loaded, self.emg_ch_loaded] + for idx, line in enumerate([self.line_time_between, self.line_cut_time, self.line_num_ch, + self.line_num_lines, self.line_eog_ch, self.line_emg_ch]): + line.setText(str(lista[idx])) + + def clickBox(self, state): + #change state of the BoxChecked variable + if state == QtCore.Qt.Checked: + self.BoxChecked = True + else: + self.BoxChecked = False + + def systemChanged(self): + if self.combobox_system.currentIndex()==0: + self.line_ip.setText('192.168.200.201') + self.line_port.setText('50000') + elif self.combobox_system.currentIndex()==1: + self.line_ip.setText("169.254.252.66") + self.line_port.setText('51244') + + def check_connection(self): + self.ip = str(self.line_ip.text()) + self.port = int(self.line_port.text()) + + self.run = Value('i', 1) #value that can be changed inside separated process + self.speed = Value(ctypes.c_float, 250/1000) #to control data flow speed we can initialize another variable used inside the function + self.q = Queue() #queue for data from separated process (data acquire) + #Run acquire_data function in separated process so it doesn't freeze when other operations are going + + if self.combobox_system.currentIndex()==0: + self.offline = "NeurOne" #"BrainProducts" #False + else: self.offline = "BrainProducts" + self.checking_connection = QMessageBox(self) + self.checking_connection.setText('Checking connection, please wait!') + self.checking_connection.setWindowTitle("Checking connection") + # self.checking_connection.setStandardButtons(QMessageBox.Ok) + self.checking_connection.show() + + self.p = Process(target=acquire_data, args=(self.q, 100, + self.run, self.speed, True, + 1, self.ip, self.port, self.offline)) + done=False + self.p.start() #start it + tries = 0 + while tries<50: + print(self.q.qsize()) + if self.q.qsize()>0: + while not done: + try: + shape = self.q.get().shape + print(shape) + if self.offline==False or self.offline=="NeurOne": + shape=tuple(reversed(shape)) + done=True + except AttributeError: + print('exception') + continue + + if shape[0]>0 and shape[1]>0: + self.p.terminate() + self.p.join(0.1) + self.q.close() + self.checking_connection.done(1) + self.confirmBoxCheck = QMessageBox(self) + self.confirmBoxCheck.setText("{} channels detected. {}/100 samples received.".format( + shape[0], shape[1])) + self.confirmBoxCheck.setWindowTitle("Signal detected!") + self.confirmBoxCheck.setStandardButtons(QMessageBox.Ok) + self.confirmBoxCheck.show() + return None + tries+=1 + time.sleep(0.2) + if tries>=50: + self.p.terminate() + self.p.join(0.1) + self.q.close() + self.confirmBoxCheck = QMessageBox(self) + self.confirmBoxCheck.setText("Signal was not detected. Check the connection!") + self.confirmBoxCheck.setWindowTitle("No signal detected!") + self.confirmBoxCheck.setStandardButtons(QMessageBox.Ok) + self.confirmBoxCheck.show() + return None + + def create_figure(self): + """ + Creates figure, defines layout, plots figures and defines connections + between them and events. + """ + def line_picker(line, mouseevent): + """ + NOT USED ATM + + Find the points within a certain distance from the mouseclick in + data coords and attach some extra attributes, pickx and picky + which are the data points that were picked. + """ + if mouseevent.xdata is None: + return False, dict() + xdata = line.get_xdata() + ydata = line.get_ydata() + maxd = 0.005 + d = np.sqrt( + (xdata - mouseevent.xdata)**2 + (ydata - mouseevent.ydata)**2) + + ind, = np.nonzero(d <= maxd) + if len(ind): + pickx = xdata[ind] + picky = ydata[ind] + props = dict(ind=ind, pickx=pickx, picky=picky) + return True, props + else: + return False, dict() + + def reset_included_ch_choice(event): + if np.sum(self.color)==0: + for i in range(len(self.color)): + self.color[i] = [0,0,1] + else: + for i in range(len(self.color)): + self.color[i] = [0,0,0] + self.included_ch = [] + for idx, col in enumerate(self.color): + if list(col) == [0,0,1]: + self.included_ch.append(idx) + self.ddd.set_color(self.color) + print(self.included_ch) + self.canvas.draw() + + def onpick3(event): + """ + Event handler for electrode selection plot + + Parameters + ---------- + event : event type from PyQt5 event handler. + + Returns + ------- + int + returns 0 if the wrong plot clicked. + + """ + if event.artist!=self.ddd: + return 0 + ind = event.ind + if list(self.color[event.ind][0]) == [0,0,0]: + self.color[event.ind] = [0,0,1] + else: + self.color[event.ind] = [0,0,0] #int(self.color[event.ind]!=1) + self.included_ch = [] + for idx, col in enumerate(self.color): + if list(col) == [0,0,1]: + self.included_ch.append(idx) + print(self.included_ch) + self.ddd.set_color(self.color) + self.canvas.draw() + print('onpick3 scatter:', ind, self.pos2d[ind,0], self.pos2d[ind,1]) + + def onpick_full(event): + """ + Event handler for full cap plot + + Parameters + ---------- + event : event type from PyQt5 event handler. + + Returns + ------- + None. + """ + + ind = event.ind + if list(self.color_full[event.ind][0]) == [0,0,0]: + self.color_full[event.ind] = [0,0,1] + else: + self.color_full[event.ind] = [0,0,0] #int(self.color[event.ind]!=1) + self.included_ch_full = [] + for idx, col in enumerate(self.color_full): + if list(col) == [0,0,1]: + self.included_ch_full.append(idx) + print(self.included_ch_full) + self.ddd_full.set_color(self.color_full) + self.pos2d = self.pos2d_full[self.included_ch_full] + #self.ddd.set_offsets(np.c_[self.pos2d[:,0], self.pos2d[:,1]]) + #Clear the plot before replottig. Computationally very inefficient + #but at that stage doesn't matter. + self.axes.clear() + self.canvas2.draw() + self.color = np.zeros([len(self.pos2d[:,0]),3]) + self.all_names = self.full_all_names[self.included_ch_full].reset_index(drop=True) + draw_main() + #self.canvas.draw() + print('onpick_full scatter:', ind, self.pos2d_full[ind,0], self.pos2d_full[ind,1]) + + self.main_frame = QWidget() + self.dpi = 100 + self.fig = Figure((11, 10), dpi=self.dpi)#, facecolor='whitesmoke') #figure for signal plot + self.fig2 = Figure((11, 10), dpi=self.dpi) + self.canvas = FigureCanvas(self.fig) + self.canvas2 = FigureCanvas(self.fig2) + self.axes = self.fig.subplots() + self.axes2 = self.fig2.subplots() + self.color = np.zeros([len(self.pos2d[:,0]),3]) + self.color_full = np.zeros([len(self.pos2d_full[:,0]),3]) + + for idx in np.array(self.ch_names_loaded)[self.included_ch_loaded]: + num = np.where(self.all_names.astype(str)==idx) + self.color[num] = [0,0,1] + + for idx in np.array(self.ch_names_loaded)[self.included_ch_loaded]: + num = np.where(self.full_all_names.astype(str)==idx) + self.color_full[num] = [0,0,1] + + for idx, col in enumerate(self.color): + if list(col) == [0,0,1]: + self.included_ch.append(idx) + + for idx, col in enumerate(self.color_full): + if list(col) == [0,0,1]: + self.included_ch_full.append(idx) + + ### input electrodes plot ### + def draw_main(): + """ + Draws the main plot. Made as a function because it's being redrawn + in case of electrode selection change (second plot). + """ + self.ddd_bck = self.axes.scatter(self.pos2d_full[:,0], self.pos2d_full[:,1], + s=600, picker=True, c='silver')#, picker=line_picker) + self.ddd = self.axes.scatter(self.pos2d[:,0], self.pos2d[:,1], s=600, + picker=True, c=self.color)#, picker=line_picker) + self.axes.set_axis_off() + self.canvas.mpl_connect('pick_event', onpick3) + self.axes.text(0.5, 1.05, 'Anterior', transform=self.axes.transAxes, ha='center', va='center') + self.axes.text(0.5, -0.05, 'Posterior', transform=self.axes.transAxes, ha='center', va='center') + texts = [] + for i in range(len(self.pos2d[:,0])): + texts.append(self.axes.text(self.pos2d[i,0], self.pos2d[i,1], + self.all_names[i], c='w', ha='center', + va='center', fontsize=15)) + self.info_text1 = self.axes.text(0,1, "Click on the electrode to select or deselect it", + transform=self.axes.transAxes, ha='left', va='center', fontsize=7) + self.axnext = self.fig.add_axes([0.81, 0.05, 0.1, 0.075]) + self.bnext = MatplotlibButton(self.axnext, 'Reset selection') + self.bnext.on_clicked(reset_included_ch_choice) + self.canvas.draw() + + draw_main() + ### full cap plot ### + self.ddd_full = self.axes2.scatter(self.pos2d_full[:,0], self.pos2d_full[:,1], + s=600, picker=True, c=self.color_full)#, picker=line_picker) + self.axes2.set_axis_off() + self.canvas2.mpl_connect('pick_event', onpick_full) + self.axes2.text(0.5, 1.05, 'Anterior', transform=self.axes2.transAxes, ha='center', va='center') + self.axes2.text(0.5, -0.05, 'Posterior', transform=self.axes2.transAxes, ha='center', va='center') + texts2 = [] + for i in range(len(self.pos2d_full[:,0])): + texts2.append(self.axes2.text(self.pos2d_full[i,0], self.pos2d_full[i,1], + self.full_all_names[i], c='w', ha='center', + va='center', fontsize=15)) + self.info_text2 = self.axes2.text(0,1, "Select channels you defined in the EEG system's protocol as a real-time output, not\n"\ + "channels you want to use for the feature calculation! Number of channels have to match \n"\ + "number of the EEG channels from the EEG system", + transform=self.axes2.transAxes, ha='left', va='center', fontsize=7) + self.canvas2.draw() + ### Buttons setup ### + + self.button_change_electrodes = QPushButton("&Load electrode locations") + self.button_change_electrodes.setCheckable(False) + self.button_change_electrodes.setMaximumWidth(250) + self.button_change_electrodes.clicked.connect(self.get_electrodes_file) + + self.button1 = QPushButton("&Select montage") + self.button1.setCheckable(False) + self.button1.setMaximumWidth(370) + self.button1.clicked.connect(self.get_file) + + self.button2 = QPushButton("&Load configuration") + self.button2.setCheckable(True) + self.button2.setMaximumWidth(370) + self.button2.clicked.connect(self.get_configuration) + + self.button3 = QPushButton("&Check connection") + self.button3.setCheckable(False) + self.button3.setMaximumWidth(370) + self.button3.clicked.connect(self.check_connection) + + self.button4 = QPushButton("&Run program") + self.button4.setCheckable(False) + self.button4.clicked.connect(self.run_main_program) + + def add_thing(self, text, settext, LineMaxWidth = 70): + """add single row (label+lineedit) as a horizontal box layout""" + layout = QHBoxLayout() + lab = QLabel(self) + lab.setText(text) + lab.setMaximumWidth(300) + line = QLineEdit(self) + line.setMaximumWidth(LineMaxWidth) + line.setAlignment(Qt.AlignRight) + layout.addWidget(lab) + layout.addWidget(line) + if settext!=None: + line.setText(str(settext)) + return lab, line, layout + + #Creating label+lineedit as a layouts, to be able to add them to the main + #layout later and keep the layout. + self.time_between_label, self.line_time_between, time_between_layout = add_thing(self, "Time between trains [s]:", self.time_between_bursts_loaded) + self.cut_time_lab, self.line_cut_time, cut_time_layout = add_thing(self, "Cut time [s]:", self.breaktime_loaded) + self.num_ch_lab, self.line_num_ch, num_ch_layout = add_thing(self, "Number of channels:", self.num_of_ch_loaded) + self.num_lines_lab, self.line_num_lines, num_lines_layout = add_thing(self, "Number of lines:", self.num_of_lines_loaded) + self.eog_ch_lab, self.line_eog_ch, eog_ch_layout = add_thing(self, "EOG channel number:", self.eog_ch_loaded) + self.emg_ch_lab, self.line_emg_ch, emg_ch_layout = add_thing(self, "EMG channel number:", self.emg_ch_loaded) + self.exp_trig_lab, self.line_exp_trig, exp_trig_layout = add_thing(self, "Number of bursts within the train:", self.exp_trig_loaded) + self.exp_time_lab, self.line_exp_time, exp_time_layout = add_thing(self, "Expected time of a single train:", self.exp_time_loaded) + + # You can add feature name if function was added to the function "features" in the main file + features_names = ['None', 'Theta FFT Power', 'Alpha FFT Power', 'Beta FFT Power', + 'High Gamma FFT Power', 'Spectral entropy', 'Temporal entropy', + 'Line length', 'DWT Power 0-4 Hz', 'DWT 4-8 Hz', + 'DWT 8-16 Hz', 'DWT 16-31 Hz','Variance','Correlation'] + + self.combobox1 = QComboBox() + self.combobox1.addItems(features_names) + + self.combobox2 = QComboBox() + self.combobox2.addItems(features_names) + + self.combobox3 = QComboBox() + self.combobox3.addItems(features_names) + + self.combobox4 = QComboBox() + self.combobox4.addItems(features_names) + + self.combobox5 = QComboBox() + self.combobox5.addItems(features_names) + + self.combobox6 = QComboBox() + self.combobox6.addItems(features_names) + + #Montage file path text + self.file_path = QLabel(self) + self.file_path.setText("No montage selected") + self.file_path.setMaximumWidth(370) + self.file_path.setWordWrap(True) + + #Path for a electrode location file text + self.text_right = QLabel(self) + self.text_right.setText("Electrodes names and locations file: " + self.cap_file_path) + text_left = QLabel(self) + text_left.setText("General settings:") + + self.box = QCheckBox("Slow mode",self) + self.box.stateChanged.connect(self.clickBox) + self.eye_reg_box = QCheckBox("Eye regression",self) + self.notch_box = QCheckBox("Notch filter",self) + self.outliers_box = QCheckBox("Remove outliers",self) + + box_layout = QHBoxLayout() + box_layout.addWidget(self.box) + box_layout.addWidget(self.eye_reg_box) + box_layout.addWidget(self.notch_box) + # box_layout.addWidget(self.outliers_box) + + #Some description of settings + text_last_ch = QLabel(self) + text_last_ch.setWordWrap(True) + text_last_ch.setFixedWidth(370) + text_last_ch.setText("Channels can be set using Python notation (negative numbers as indices from the end). "\ + "Number of channels and timings have to correspond to these in EEG system protocol. \n\n"\ + "Channels chosen on the right screen will be used for feature calculation. "\ + "Visualization will include all channels that are defined in EEG system protocol. "\ + "To change used channels: \n"\ + "\u2022 change protocol in EEG system \n"\ + "\u2022 change TMS_protocol.txt file. Rows \"included channels\" and \"names\". " \ + "If \"names\" will be set as \"def\" it will get data automatically from the location file. \n"\ + "\u2022 change electrode location file to contain proper number of elements. \n"\ + "\n\n"\ + "Last channel (-1) is always reserved for TMS pulse marking.\n\n"\ + "Slow mode will reduce the refresh rate of the raw signal plot. It will not affect features readout or calibration.") + + #configuration file path text + self.text_conf_file = QLabel(self) + self.text_conf_file.setMaximumWidth(370) + self.text_conf_file.setWordWrap(True) + self.text_conf_file.setText("No configuration file selected") + + #Tab widget for both plots as tabs + tabwidget = QTabWidget() + + #add every widget and layout to the main layout + vbox = QVBoxLayout() + vbox.addWidget(self.button3) + vbox.addWidget(self.text_conf_file) + vbox.addWidget(self.button2) + vbox.addStretch(1) + vbox.addLayout(time_between_layout) + vbox.addLayout(cut_time_layout) + vbox.addLayout(num_ch_layout) + vbox.addLayout(num_lines_layout) + vbox.addLayout(eog_ch_layout) + vbox.addLayout(emg_ch_layout) + vbox.addLayout(exp_trig_layout) + vbox.addLayout(exp_time_layout) + + scroll = QScrollArea() + scroll.setWidget(text_last_ch) + scroll.setWidgetResizable(True) + scroll.setFixedWidth(390) + + vbox.addWidget(scroll) + vbox.addLayout(box_layout) + vbox.addWidget(self.outliers_box) + vbox.addStretch(5) + vbox.addWidget(self.button1) + vbox.addWidget(self.file_path) + + feature_choice_layout = QVBoxLayout() + text_for_combo = QLabel(self) + #text_for_combo.setMaximumWidth(370) + text_for_combo.setWordWrap(True) + text_for_combo.setText("You can choose up to 6 different measurements that will be "\ + "calculated during the intervention.") + + connection_settings_text = QLabel(self) + connection_settings_text.setWordWrap(True) + connection_settings_text.setText('Set type of EEG system, ip address, and port') + + self.combobox_system = QComboBox() + self.combobox_system.addItems(['NeurOne', 'Brain Products']) + self.combobox_system.currentIndexChanged.connect(self.systemChanged) + + self.ip_box, self.line_ip, ip_layout = add_thing(self, 'IP: ', self.ip, 100) + self.port_box, self.line_port, port_layout = add_thing(self, 'Port: ', self.port, 100) + + comboboxes = [text_for_combo, self.combobox1, self.combobox2, self.combobox3, + self.combobox4, self.combobox5, self.combobox6] + + for combobox in comboboxes: + feature_choice_layout.addWidget(combobox) + feature_choice_layout.addStretch(1) + feature_choice_layout.addWidget(connection_settings_text) + feature_choice_layout.addWidget(self.combobox_system) + feature_choice_layout.addLayout(ip_layout) + feature_choice_layout.addLayout(port_layout) + self.feature_choice_widget = QWidget() + self.feature_choice_widget.setLayout(feature_choice_layout) + feature_choice_layout.addStretch(2) + + hbox = QHBoxLayout() + hbox.addWidget(self.button_change_electrodes) + hbox.addWidget(self.text_right) + #set genera layout + layout = QGridLayout() + # layout.addWidget(self.button2, 0,1,1,1) + layout.addLayout(hbox, 0, 1, 1, 1) + # layout.addWidget(self.button_change_electrodes, 0, 1, 1, 1) + # layout.addWidget(self.text_right, 0, 1, 1, 1) + layout.addWidget(text_left, 0, 0, 1, 1) + + tabwidget.addTab(self.canvas, 'Select electrodes for features') + tabwidget.addTab(self.canvas2, 'See the whole cap') + tabwidget.addTab(self.feature_choice_widget, 'Features and connection') + + #layout.addWidget(self.canvas, 1, 1, 1, 1) + layout.addLayout(vbox, 1, 0, 1, 1) + layout.addWidget(self.button4, 2,1,1,1) + layout.addWidget(tabwidget, 1,1,1,1) + + self.main_frame.setLayout(layout) #set main layout + + #Tight layout to remove padding on the sides + self.fig.tight_layout(pad=0.95)#, rect=(0.02,-0.02,1.02,1.02)) + self.fig2.tight_layout(pad=0.95) + self.setCentralWidget(self.main_frame) + self.show() + + def get_file(self): + """ + Load file path and set it in the main program. + + """ + err = False + self.path = QFileDialog.getOpenFileName(self, 'Open File', os.path.dirname(os.getcwd()), 'Text files (*.txt *.csv)') + print(self.path) + self.montage_file_path = self.path[0] + montage_matrix = np.array([]) + try: + montage_matrix = np.array(pd.read_csv(self.montage_file_path, header=None)) + except: + err = True + if len(montage_matrix.shape)!=2 or montage_matrix.shape[0]!=montage_matrix.shape[1] or err==True: + self.montage_file_path = None + self.bad_montage_msg = QMessageBox(self) + self.bad_montage_msg.setIcon(QMessageBox.Warning) + self.bad_montage_msg.setText("The montage file is incorrect. It has to be csv file resulting in (x, x) size. Choose the correct file!") + self.bad_montage_msg.setWindowTitle("Montage file error") + self.bad_montage_msg.setStandardButtons(QMessageBox.Ok) + self.bad_montage_msg.show() + else: + self.file_path.setText(self.path[0]) + + def restart_kernel(i): + os._exit(00) + + def get_electrodes_file(self): + """ + Load file path for electrode placement and set it in the main program. + + """ + + temp_path = QFileDialog.getOpenFileName(self, 'Full cap file', os.path.dirname(os.getcwd()), 'Text files (*.txt *.csv)') + print(temp_path) + if temp_path[0]=="": + return 0 + + temp_path2 = QFileDialog.getOpenFileName(self, 'Selected electrodes', os.path.dirname(os.getcwd()), 'Text files (*.txt *.csv)') + print(temp_path2) + if temp_path2[0]=="": + return 0 + + with open('Electrode_selection.txt', 'w') as f: + f.write('full_cap_file_path: {}'.format(temp_path[0])) + f.write('\n') + f.write('cap_file_path: {}'.format(temp_path2[0])) + + self.msg = QMessageBox() + self.msg.setIcon(QMessageBox.Information) + self.msg.setText("Python Kernel will restart now. Run the program again.") + # msg.setInformativeText("This is additional information") + self.msg.setWindowTitle("Restart") + # msg.setDetailedText("The details are as follows:") + self.msg.setStandardButtons(QMessageBox.Ok) + self.msg.buttonClicked.connect(self.restart_kernel) + self.msg.show() + + def get_configuration(self): + """ + Load file path for configuration file and re-set values in the window. + + Returns + ------- + bool + Used only as a function breaking exception. Returns False when no + path was chosen. + + """ + self.conf_path = QFileDialog.getOpenFileName(self, 'Open File', os.path.dirname(os.getcwd()), 'Text files (*.txt *.csv)')[0] + print(self.conf_path) + if self.conf_path=="": + return False + self.text_conf_file.setText("Configuration file: " + self.conf_path) + self.set_values(True) + + def run_main_program(self): + """ + Runs main program and passes all parameters. + + """ + def confirmation(msg): + """ + If OK option was chosen passes all parameters and runs the main program. + + Parameters + ---------- + msg : msg type from QMessageBox + """ + print(str(msg.text())) + if msg.text()=="&OK" or msg.text()=="OK": + self.main_app = self.AppForm(self.params_to_pass) + self.main_app.show() + self.confirmBox.close() + self.showMinimized() + #self.close() + + + #add EOG and EMG name to make channel names list complete + + + arr_temp = np.arange(int(self.line_num_ch.text())) + if self.line_eog_ch.text().strip()=="" and self.line_emg_ch.text().strip()=="": + all_names_temp = self.all_names + elif self.line_eog_ch.text().strip()=="" and self.line_emg_ch.text().strip()!="": + all_names_temp = np.append(self.all_names, ["EMG"]) + elif self.line_eog_ch.text().strip()!="" and self.line_emg_ch.text().strip()=="": + all_names_temp = np.append(self.all_names, ["EOG"]) + elif arr_temp[int(self.line_eog_ch.text())] < arr_temp[int(self.line_emg_ch.text())]: + all_names_temp = np.append(self.all_names, ["EOG", "EMG"]) + else: + all_names_temp = np.append(self.all_names, ["EMG", "EOG"]) + print(all_names_temp) + chosen_features = [self.combobox1.currentIndex(), self.combobox2.currentIndex(), + self.combobox3.currentIndex(), self.combobox4.currentIndex(), + self.combobox5.currentIndex(), self.combobox6.currentIndex()] + + if self.combobox_system.currentIndex()==0: + self.offline = "NeurOne" #"BrainProducts" #False + else: self.offline = "BrainProducts" + + #Create a dictionary with structures to pass + self.params_to_pass = {'time_between': int(self.line_time_between.text()), + 'cut_time': int(self.line_cut_time.text()), + 'included_ch': self.included_ch, + 'eog_ch': json.loads('[{}]'.format(self.line_eog_ch.text())), #that will turn it into list (also '2,3' --> [2, 3]) + 'emg_ch': json.loads('[{}]'.format(self.line_emg_ch.text())), #but don't use more than one channel for now! + 'num_of_ch': int(self.line_num_ch.text()), + 'num_of_lines': int(self.line_num_lines.text()), + 'montage_path': self.montage_file_path, + 'ch_names': all_names_temp, + 'slow_mode': self.box.isChecked(), + 'features': chosen_features, + 'notch': self.notch_box.isChecked(), + 'eye_reg': self.eye_reg_box.isChecked(), + 'remove_outliers': self.outliers_box.isChecked(), + 'ip': str(self.line_ip.text()), + 'port': int(self.line_port.text()), + 'offline': self.offline, + 'exp_trig': int(self.line_exp_trig.text()), + 'exp_time': int(self.line_exp_time.text()), + } + print(self.params_to_pass) + + #ConfirmationBox to check and inform the user if everything is fine. + #If not, pass the comment about what is wrong + self.confirmBox = QMessageBox(self) + self.confirmBox.setDetailedText("Initial choice\nNames: {}\nIncluded:{}\n\n"\ + "Final choice\nNames: {}\nIncluded: {}".format( + self.ch_names_loaded, self.included_ch_loaded, + self.params_to_pass["ch_names"], + self.params_to_pass["included_ch"])) + if self.params_to_pass['num_of_ch'] == len(self.params_to_pass["ch_names"]): + msg_ending = "Set number of channels is equal with number of chosen electrodes :)\n\n"\ + "channels selected (parameter): {}\nchosen channels: {}".format( + self.params_to_pass['num_of_ch'], len(self.params_to_pass["ch_names"])) + else: + msg_ending = "Set number of channels is NOT equal with number of chosen electrodes :(\n\n"\ + "channels selected (parameter): {}\nchosen channels: {}".format( + self.params_to_pass['num_of_ch'], len(self.params_to_pass["ch_names"])) + if np.array_equiv(all_names_temp, self.ch_names_loaded): + self.confirmBox.setIcon(QMessageBox.Information) + self.confirmBox.setText("Confirm the choice of electrodes. Click details to check it.") + else: + self.confirmBox.setIcon(QMessageBox.Warning) + self.confirmBox.setText("You've changed the electrode selection in the software. Please check and confirm your choice. "\ + "Remember that the selection of electrodes has to correspond to the protocol from the EEG system! "\ + "If the number of input channels will be the same as the number of labels, but labels will be incorrect, "\ + "software will run showing incorrect labels.\n\n{}".format(msg_ending)) + self.confirmBox.setWindowTitle("Confirm changes") + self.confirmBox.setStandardButtons(QMessageBox.Ok | QMessageBox.Cancel) + self.confirmBox.buttonClicked.connect(confirmation) + self.confirmBox.show() + + +if __name__ == '__main__': + app = QApplication(sys.argv) + form = First_window('test') #AppForm() + form.show() + app.exec_() \ No newline at end of file diff --git a/Functions.py b/Functions.py new file mode 100644 index 0000000..efa0ee9 --- /dev/null +++ b/Functions.py @@ -0,0 +1,205 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Jan 24 11:48:34 2022 + +@author: Basics +""" +import numpy as np +import time +import scipy.signal as ss +import scipy + +def apply_montage(data, matrix): + #print(matrix) + if data.shape[0]!=matrix.shape[0]: + if data.shape[0]==matrix.shape[0]+1: + data = data[:-1] + times = time.time() + new_data = np.zeros(data.shape) + # for i in range(matrix.shape[0]): + # for j in range(matrix.shape[1]): + # if float(abs(matrix[i,j]))>0: + # new_data[i,:]+=data[j]*float(matrix[i,j]) + # print(time.time()-times) + print(data.shape, matrix.shape) + new_data = (data.T@matrix.T).T + return new_data + +def eye_reg(eeg, eog, regg = True, Fs=1000): + if regg: + print(eeg.shape) + eeg = eeg.T + #eeg = ss.detrend(eeg.T) + print(eog.shape) + [A,B] = ss.butter(2, 40/(Fs/2), 'lowpass') + eog = ss.filtfilt(A, B, eog) + eogt = ss.detrend(eog.reshape([len(eog),1]),axis=0) + #eogt = eog.reshape([len(eog),1]) + data_reg_eog = eeg - np.dot(eogt, np.linalg.lstsq(eogt,eeg, rcond=None)[0]) + print('EYE REG ZROBIONEEEE hmmm') + print(data_reg_eog.shape, eogt.shape,eeg.shape) + # plt.figure() + # plt.plot(data_reg_eog[:,-5],c='b') + # plt.plot(eeg[:,-5], c='orange') + # plt.plot(eogt, c='green') + # plt.show() + return(data_reg_eog).T + else: + return(eeg) + +def most_frequent(arr): + """returns most frequent item in the array + + Parameters + ------------ + arr: array type + + Returns + ------------ + most frequent item in the list or 0 + """ + try: + counts = np.bincount(arr) + return np.argmax(counts) + #return max(set(List), key = List.count) + except: + return 0 + +def connect_sig(data1, data2, fs): + """Knowing that signal doesn't updates perfectly every x seconds, I was looking for a way + to connect it in proper time, but finally other approach have been used + + Parameters + ------------- + data1, data2: numpy array types with MxN size + fs: int + sampling rate + """ + print(data2.shape) + data2 = data2 + startt = time.time() + num = data1.shape[0] + size = data1.shape[1] + pts = list() + data_ret = np.zeros(data1.shape) + data_ret[:, :-fs] = data1[:,fs:] + if data2.shape[0]==0 or data2.shape[1]==0: + return data2 + if fs<2000: + for i in range(num): + try: + pts.extend(np.where(data1[i,-1]==data2[i])[0].tolist()) + except: + print("ARBEJDE IKKEEE") + # data_ret = np.concatenate((data1, data2[:,-int(size):]),1) + data_ret[:,-fs:] = data2[:, -fs:] + return data_ret, None + #print('hehe', time.time()-startt) + most_fr = most_frequent(np.array(pts)) + #print('hehe', time.time()-startt) + print(list(pts).count(most_fr)>num//2, most_fr) + if fs<2000 and list(pts).count(most_fr)>num//2: + most_fr = most_fr+1 + #print(data2[:,int(pts[i]):int(pts[i]+size)].shape) + print('hehe', time.time()-startt) + #data_ret = np.concatenate((data1, data2[:,int(most_fr)+1:int(most_fr+size)+1]),1) + data_ret[:,-fs:] = data2[:, most_fr:most_fr+fs] + print('connect:', time.time()-startt) + return data_ret, most_fr + else: + data_ret[:,-fs:] = data2[:, -fs:] + print("NEEEJJJJJJJ") + print('connect:', time.time()-startt) + return data_ret, 800 + +def set_to_gray(lines): + for i in range(lines.shape[0]): + for j in range(lines.shape[1]): + if type(lines[i,j]) != int: + lines[i,j][0].set_color("gray") + +def update_stem(line, data, ax, relim=False): + x = adjust_hist(data[1]) + y = data[0] + + line[0].set_ydata(y) + line[0].set_xdata(x) # not necessary for constant x + + # stemlines + # line[1].set_paths([np.array([[xx, 0], [xx, yy]]) for (xx, yy) in zip(x, y)]) + # line[2].set_xdata([np.min(x), np.max(x)]) + # line[2].set_ydata([0, 0]) # not necessary for constant bottom + if relim: + ax.relim() + # update ax.viewLim using the new dataLim + ax.autoscale_view(scalex=True) + + +def adjust_hist(data): + new_data = np.zeros(len(data)-1) + for i in range(len(new_data)): + new_data[i] = np.mean([data[i], data[i+1]]) + return new_data + +def min_zero(thr): + if thr[0]<0: + return [0, thr[1]] + else: + return thr + +def pentropy(signal, Fs, nperseg=None, fmin=None, fmax=None): + #I think it's good to put own nperseg value, because default one is not always good in this case + f, time, S = ss.spectrogram(signal, Fs, nperseg=nperseg) #spectrogram of signal + if fmin and fmax: + idxs = np.where((ffmin)) #choosing only bands that we are interested in + S = S[idxs] + P = np.zeros(S.shape) + H = np.zeros(S.shape[1]) + for t in range(S.shape[1]): + for m in range(S.shape[0]): + P[m,t] = S[m,t]/np.sum(S[:,t],0) #according to matlab instruction + H[t] += -P[m,t]*np.log2(P[m,t]+0.000001)/np.log2(S.shape[0]) + return f, S,time, H + +def entropy(S, bins=100): + histo = np.histogram(S, bins)[0] + p = scipy.special.entr(histo/histo.sum()) + ent = sum(p) + return ent + +def check_integrity(stim): + """If the list contains two consecutive values, the function leaves only the first. + It was needed because NeurOne sometimes returns two stimuli pulses78. + + Parameters + ----------------- + stim: list or array type with numbers, in our case indexes of stimuli + + Returns + ----------------- + ent: array with removed values + + """ + #copy the list, I need two, because from one values are removed, and other gives + #proper index values + stim_c = list(stim.copy()) + stim = list(stim) + for ind in range(len(stim_c)-1): + if stim_c[ind]==(stim_c[ind+1]-1): + stim.remove(stim_c[ind+1]) + return np.array(stim) + +def doubleMADsfromMedian(y,thresh=3.5): + # warning: this function does not check for NAs + # nor does it address issues when + # more than 50% of your data have identical values + m = np.median(y) + abs_dev = np.abs(y - m) + left_mad = np.median(abs_dev[y <= m]) + right_mad = np.median(abs_dev[y >= m]) + y_mad = left_mad * np.ones(len(y)) + y_mad[y > m] = right_mad + modified_z_score = 0.6745 * abs_dev / y_mad + modified_z_score[y == m] = 0 + return modified_z_score > thresh + diff --git a/NeurOne.py b/NeurOne.py new file mode 100644 index 0000000..b56bfd3 --- /dev/null +++ b/NeurOne.py @@ -0,0 +1,199 @@ +# -*- coding: utf-8 -*- +""" +Created on Tue Aug 16 09:15:49 2016 + +@author: KHM +""" + +import socket,time,sys,datetime +import numpy as np +#import scipy.io +if sys.version_info[0]>=3: + import queue +else: + import Queue as queue + +import threading,platform +if platform.system()=='Windows': + tfunc=time.time +else: + tfunc=time.time + +def bytes_to_int32(databuf): + in_data=np.frombuffer(databuf,dtype=np.dtype([('1','i1'),('2','>u2')])) + data=in_data['1'].astype('i4') + data <<= 16 + data |= in_data['2'] + return data + +def sampleLoop(no): + firstpackage=True + # databuf = np.empty((no.bufsiz,), dtype=np.uint8) + databuf = bytearray(b' ' * no.bufsiz) + readdump=False + if not no.readdump is None: + fnin=open(no.readdump,'rb') + readdump=True + tsamp=tfunc()+no.si + sampidx=0 + firstsamp=0 + oldsampidx=0 + droppeds=0 + timeout=False + if not no.dump is None: + if no.dump==True: + fn=open('dump_'+datetime.datetime.now().strftime('%Y%m%d_%H-%M-%S')+'.raw','wb') + else: + fn=open(no.dump,'wb') + while not no.stop: + try: + if readdump: + while tsamp>tfunc(): + #time.sleep(tsamp-tfunc()) + time.sleep(0.) + tsamp+=no.si + tsin=np.frombuffer(fnin.read(10),dtype=np.uint16) + nbytes=tsin[-1] + tsin=tsin[:-1].view(np.float64)[0] + databuf[:nbytes]=np.frombuffer(fnin.read(nbytes),dtype=np.uint8) + else: + nbytes=no.sock.recv_into(databuf) + t0=tfunc() + if timeout: + print('Connection re-established.') + timeout=False + except: + nbytes=0 + if not timeout: + print('Timeout - package too small, will keep retrying every second.') + timeout=True + time.sleep(1.) + if nbytes>28: + #spacket = databuf[0] + #mainunit = databuf[1] + packetno = np.frombuffer(databuf[4:8],'>u4').copy() + nch = np.frombuffer(databuf[8:10],'>u2').copy() + nsamp = np.frombuffer(databuf[10:12],'>u2').copy() + sampidx = np.frombuffer(databuf[12:20],'>u8').copy() + tstamp = np.frombuffer(databuf[20:28],'>u8').copy() + if oldsampidx!=0: + ds=sampidx-oldsampidx + if ds>nsamp: + droppeds += ds-nsamp + print('Dropped %i samples'%(ds-nsamp,)) + elif ds!=nsamp: + print('delta samp %i, samples %i'%(ds,nsamp)) + else: + firstsamp=sampidx + oldsampidx=sampidx + data=bytes_to_int32(databuf[28:28 + 3 * nch[0] * nsamp[0]]).reshape((nsamp[0],nch[0])) + t1=tfunc() + if not no.ringbuffer is None: + if no.avgPackets: + data1=data.mean(axis=0, keepdims=True) + data1[:,-1]=np.max(data[:,-1]) + no.updateRingBuffer(data1,sampidx,(tstamp,t1)) + else: + no.updateRingBuffer(data,sampidx,(tstamp,t1)) + if no.sendqueue: + no.queue.put((data,(packetno,sampidx,tstamp,t0,t1))) + if firstpackage: + no.tstamp0=(tstamp,t1) + if no.dump: + fn.write(np.array(t1).tostring()) + fn.write(np.array(nbytes,dtype=np.uint16).tostring()) + fn.write(databuf[:nbytes]) + try: + fn.close() + except: + pass + try: + fnin.close() + except: + pass + totals=sampidx-firstsamp + if totals>0: + if droppeds>0: + print('Dropped %i out of %i samples (%.1f%%)'%(droppeds,totals,droppeds/totals*100.)) + else: + print('Acquired %i samples none were dropped.'%(totals,)) + else: + print('No samples acquired.') + +class NO(): + def __init__(self,ip='127.0.0.1',port=50000,buffersize=2**10,ringbuffersize=None,sendqueue=False,\ + ringbuf_factor=2,dump=None,readdump=None,si=1./1000.,avgPackets=False): + if readdump is None: + self.sock=socket.socket(socket.AF_INET, # Internet + socket.SOCK_DGRAM) #UDP + self.sock.bind((ip, port)) + self.sock.settimeout(2.) + self.avgPackets = avgPackets + self.bufsiz=buffersize + self.ip=ip + self.port=port + self.sampidx=0 + self.tstamp=None + self.tstamp0=None + self.queue=queue.Queue() + self.A=None + self.stop=False + self.dump=dump + self.readdump=readdump + if ringbuffersize is None: + self.ringbuffer=None + else: + self.ringbuffer=True + self.idx=0 + self.ringbufferinit=True + self.ringbuffersize=ringbuffersize + self.ringbuf_factor=ringbuf_factor + self.sendqueue=sendqueue + self.lock=threading.RLock() + self.si=si + + def updateRingBuffer(self,data,i=None,tstamp=None): + if self.ringbufferinit: + self.ringbuffer=np.zeros((self.ringbuffersize*self.ringbuf_factor,data.shape[1]),dtype=np.float32) + self.ringbufferinit=False + ringbuf=self.ringbuffer + wlen=self.ringbuffersize + self.lock.acquire() + if (self.idx+data.shape[0])<=ringbuf.shape[0]: + ringbuf[self.idx:self.idx+data.shape[0],:]=data + self.idx+=data.shape[0] + else: + ringbuf[0:wlen-data.shape[0],:]=ringbuf[self.idx-wlen+data.shape[0]:self.idx,:] + self.idx=wlen + ringbuf[wlen-data.shape[0]:wlen,:]=data + self.datawindow=ringbuf[self.idx-wlen:self.idx] + if not i is None: + self.sampidx=i + if not tstamp is None: + self.tstamp=tstamp + self.lock.release() + + def getBuffer(self,returnIdx=False): + self.lock.acquire() + try: + out=self.datawindow.copy() + except: + out=None + #print self.sampleno + if returnIdx: + out=(out,self.sampidx) + self.lock.release() + return out + + def start(self): + self.thread=threading.Thread(target=sampleLoop,args=(self,)) + self.thread.start() + + def stopit(self): + self.stop=True + self.thread.join() + + try: + self.sock.close() + except: + pass \ No newline at end of file diff --git a/RDA.py b/RDA.py new file mode 100644 index 0000000..6edb91d --- /dev/null +++ b/RDA.py @@ -0,0 +1,282 @@ +# -*- coding: utf-8 -*- +""" +Created on Wed Nov 16 15:00:23 2022 + +@author: s202442 +""" + + +# needs socket and struct library +from socket import socket, AF_INET, SOCK_STREAM +from struct import unpack +import sys +import numpy as np +import threading +import queue +import time + + +# Marker class for storing marker information +class Marker: + def __init__(self): + self.position = 0 + self.points = 0 + self.channel = -1 + self.type = "" + self.description = "" + +# Helper function for receiving whole message +def RecvData(socket, requestedSize): + returnStream = bytes() + while len(returnStream) < requestedSize: + databytes = socket.recv(requestedSize - len(returnStream)) + if databytes == '': + raise RuntimeError + # print(databytes) + returnStream += databytes + + return returnStream + + +# Helper function for splitting a raw array of +# zero terminated strings (C) into an array of python strings +def SplitString(raw): + stringlist = [] + s = bytes() + for i in range(len(raw)): + if raw[i] != 0: #'\x00': + s = s + raw[i].to_bytes(1, sys.byteorder) + else: + stringlist.append(s.decode()) + s = bytes() + + return stringlist + + +# Helper function for extracting eeg properties from a raw data array +# read from tcpip socket +def GetProperties(rawdata): + + # Extract numerical data + (channelCount, samplingInterval) = unpack(' lastBlock + 1: + print("*** Overflow with " + str(block - lastBlock) + " datablocks ***" ) + lastBlock = block + + data1s.extend(data) + data1s = np.array(data1s) + + # Print markers, if there are some in actual block + marker_sig = np.zeros([1, int(len(data1s)/channelCount)]) + if markerCount > 0: + for m in range(markerCount): + print("Marker " + markers[m].description + " of type " + markers[m].type) + marker_sig[0][markers[m].position] = 1 + t1 = time.time() + # Put data at the end of actual buffer + data_array = data1s.reshape([int(len(data1s)/channelCount), channelCount]) * np.array(resolutions) + data_array = np.vstack([data_array.T, marker_sig]).T #isn't that too slow? + obj.updateRingBuffer(data_array,block) + data1s = [] + + + elif msgtype == 3: + # Stop message, terminate program + print("Stop") + finish = True + obj.sock.close() + + +############################################################################################## +# +# Main RDA routine +# +############################################################################################## + + +class RDA(): + def __init__(self,ip='127.0.0.1', port=51244, buffersize=2**10, sendqueue=False, + si=1/1000, ringbuffersize = 2**12, avgPackets=False): + # Create a tcpip socket + #con = socket(AF_INET, SOCK_STREAM) + # Connect to recorder host via 32Bit RDA-port + # adapt to your host, if recorder is not running on local machine + # change port to 51234 to connect to 16Bit RDA-port + + #ip_client = "169.254.200.198 "#.96.224" + # ip_server = "169.254.252.66" + # port = 51244 + # con.connect((ip_server, port)) + self.sock=socket(AF_INET, # Internet + SOCK_STREAM) #UDP + #self.sock.bind((ip_server, port)) + self.sock.connect((ip, port)) + self.sock.settimeout(2.) + # s = socket(AF_INET, SOCK_DGRAM) + # s.bind((ip_client, port)) + # s.settimeout(5) + # print(s.recvfrom(1024)) + # con.settimeout(5) + # Flag for main loop + #finish = False + + self.avgPackets = avgPackets + self.bufsiz=buffersize + self.ip=ip + self.port=port + self.sampidx=0 + self.tstamp=None + self.tstamp0=None + self.queue=queue.Queue() + self.A=None + self.stop=False + self.idx=0 + self.ringbufferinit=True + self.ringbuffersize=ringbuffersize + self.sendqueue=sendqueue + self.lock=threading.RLock() + self.si=si + + def updateRingBuffer(self,data,i=None,tstamp=None): + if self.ringbufferinit: + self.ringbuffer=np.zeros((self.ringbuffersize ,data.shape[1]),dtype=np.float32) + self.ringbufferinit=False + ringbuf=self.ringbuffer + wlen=self.ringbuffersize + self.lock.acquire() + if (self.idx+data.shape[0])<=ringbuf.shape[0]: + ringbuf[self.idx:self.idx+data.shape[0],:]=data + self.idx+=data.shape[0] + else: + ringbuf[0:wlen-data.shape[0],:]=ringbuf[self.idx-wlen+data.shape[0]:self.idx,:] + self.idx=wlen + ringbuf[wlen-data.shape[0]:wlen,:]=data + self.datawindow=ringbuf[self.idx-wlen:self.idx] + if not i is None: + self.sampidx=i + if not tstamp is None: + self.tstamp=tstamp + self.lock.release() + + def getBuffer(self,returnIdx=False): + self.lock.acquire() + try: + out=self.datawindow.copy() + except: + out=None + if returnIdx: + out=(out,self.sampidx) + self.lock.release() + return out + + def start(self): + self.thread=threading.Thread(target=sampleLoop,args=(self,)) + self.thread.start() + + def stopit(self): + self.stop=True + self.thread.join() + + try: + self.sock.close() + except: + pass + diff --git a/README.md b/README.md new file mode 100644 index 0000000..0b8a6b9 --- /dev/null +++ b/README.md @@ -0,0 +1,92 @@ +# EStiMo + + + +## Getting started + +To make it easy for you to get started with GitLab, here's a list of recommended next steps. + +Already a pro? Just edit this README.md and make it your own. Want to make it easy? [Use the template at the bottom](#editing-this-readme)! + +## Add your files + +- [ ] [Create](https://docs.gitlab.com/ee/user/project/repository/web_editor.html#create-a-file) or [upload](https://docs.gitlab.com/ee/user/project/repository/web_editor.html#upload-a-file) files +- [ ] [Add files using the command line](https://docs.gitlab.com/ee/gitlab-basics/add-file.html#add-a-file-using-the-command-line) or push an existing Git repository with the following command: + +``` +cd existing_repo +git remote add origin https://git.drcmr.dk/adamr/estimo.git +git branch -M main +git push -uf origin main +``` + +## Integrate with your tools + +- [ ] [Set up project integrations](https://git.drcmr.dk/adamr/estimo/-/settings/integrations) + +## Collaborate with your team + +- [ ] [Invite team members and collaborators](https://docs.gitlab.com/ee/user/project/members/) +- [ ] [Create a new merge request](https://docs.gitlab.com/ee/user/project/merge_requests/creating_merge_requests.html) +- [ ] [Automatically close issues from merge requests](https://docs.gitlab.com/ee/user/project/issues/managing_issues.html#closing-issues-automatically) +- [ ] [Enable merge request approvals](https://docs.gitlab.com/ee/user/project/merge_requests/approvals/) +- [ ] [Automatically merge when pipeline succeeds](https://docs.gitlab.com/ee/user/project/merge_requests/merge_when_pipeline_succeeds.html) + +## Test and Deploy + +Use the built-in continuous integration in GitLab. + +- [ ] [Get started with GitLab CI/CD](https://docs.gitlab.com/ee/ci/quick_start/index.html) +- [ ] [Analyze your code for known vulnerabilities with Static Application Security Testing(SAST)](https://docs.gitlab.com/ee/user/application_security/sast/) +- [ ] [Deploy to Kubernetes, Amazon EC2, or Amazon ECS using Auto Deploy](https://docs.gitlab.com/ee/topics/autodevops/requirements.html) +- [ ] [Use pull-based deployments for improved Kubernetes management](https://docs.gitlab.com/ee/user/clusters/agent/) +- [ ] [Set up protected environments](https://docs.gitlab.com/ee/ci/environments/protected_environments.html) + +*** + +# Editing this README + +When you're ready to make this README your own, just edit this file and use the handy template below (or feel free to structure it however you want - this is just a starting point!). Thank you to [makeareadme.com](https://www.makeareadme.com/) for this template. + +## Suggestions for a good README +Every project is different, so consider which of these sections apply to yours. The sections used in the template are suggestions for most open source projects. Also keep in mind that while a README can be too long and detailed, too long is better than too short. If you think your README is too long, consider utilizing another form of documentation rather than cutting out information. + +## Name +Choose a self-explaining name for your project. + +## Description +Let people know what your project can do specifically. Provide context and add a link to any reference visitors might be unfamiliar with. A list of Features or a Background subsection can also be added here. If there are alternatives to your project, this is a good place to list differentiating factors. + +## Badges +On some READMEs, you may see small images that convey metadata, such as whether or not all the tests are passing for the project. You can use Shields to add some to your README. Many services also have instructions for adding a badge. + +## Visuals +Depending on what you are making, it can be a good idea to include screenshots or even a video (you'll frequently see GIFs rather than actual videos). Tools like ttygif can help, but check out Asciinema for a more sophisticated method. + +## Installation +Within a particular ecosystem, there may be a common way of installing things, such as using Yarn, NuGet, or Homebrew. However, consider the possibility that whoever is reading your README is a novice and would like more guidance. Listing specific steps helps remove ambiguity and gets people to using your project as quickly as possible. If it only runs in a specific context like a particular programming language version or operating system or has dependencies that have to be installed manually, also add a Requirements subsection. + +## Usage +Use examples liberally, and show the expected output if you can. It's helpful to have inline the smallest example of usage that you can demonstrate, while providing links to more sophisticated examples if they are too long to reasonably include in the README. + +## Support +Tell people where they can go to for help. It can be any combination of an issue tracker, a chat room, an email address, etc. + +## Roadmap +If you have ideas for releases in the future, it is a good idea to list them in the README. + +## Contributing +State if you are open to contributions and what your requirements are for accepting them. + +For people who want to make changes to your project, it's helpful to have some documentation on how to get started. Perhaps there is a script that they should run or some environment variables that they need to set. Make these steps explicit. These instructions could also be useful to your future self. + +You can also document commands to lint the code or run tests. These steps help to ensure high code quality and reduce the likelihood that the changes inadvertently break something. Having instructions for running tests is especially helpful if it requires external setup, such as starting a Selenium server for testing in a browser. + +## Authors and acknowledgment +Show your appreciation to those who have contributed to the project. + +## License +For open source projects, say how it is licensed. + +## Project status +If you have run out of energy or time for your project, put a note at the top of the README saying that development has slowed down or stopped completely. Someone may choose to fork your project or volunteer to step in as a maintainer or owner, allowing your project to keep going. You can also make an explicit request for maintainers. diff --git a/TMS_protocol.txt b/TMS_protocol.txt new file mode 100644 index 0000000..dca2794 --- /dev/null +++ b/TMS_protocol.txt @@ -0,0 +1,14 @@ +time_between_trains: 8 +cut_time: 1 +number_of_channels: 18 +number_of_lines: 4 +eog_channel: -3 +emg_channel: -2 +included_channels: [0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15] +names: [1,3,7,8,9,11,14,17,19,23,26,29,32,43,50,52, "EOG", "EMG"] +alpha_range: [8,15] +beta_range: [16,30] +theta_range: [4,8] +threshold_parameter: 2 +expected_triggers: 10 +expected_time: 2000 diff --git a/Waiting.py b/Waiting.py new file mode 100644 index 0000000..f82ceb8 --- /dev/null +++ b/Waiting.py @@ -0,0 +1,28 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Jan 24 11:26:23 2022 + +@author: Basics +""" +from PyQt5.QtWidgets import (QLabel, QDialog, QDialogButtonBox, + QVBoxLayout) + +class Waiting_window(QDialog): + def __init__(self, parent=None): + super().__init__(parent) + + self.setWindowTitle("Wait...") + + QBtn = QDialogButtonBox.Ok + + self.buttonBox = QDialogButtonBox(QBtn) + self.buttonBox.accepted.connect(self.accept) + + self.layout = QVBoxLayout() + self.message = QLabel('Waiting for the data (should take up to 30 sec)') + self.layout.addWidget(self.message) + self.layout.addWidget(self.buttonBox) + self.setLayout(self.layout) + self.show() + def setText(self,text): + self.message.setText(text) \ No newline at end of file diff --git a/easycap-M10_16_NO.txt b/easycap-M10_16_NO.txt new file mode 100644 index 0000000..f8e6bdd --- /dev/null +++ b/easycap-M10_16_NO.txt @@ -0,0 +1,18 @@ +Site Theta Phi +1 0 0 +3 23 30 +7 -23 -30 +8 46 90 +9 46 66 +11 46 0 +14 -46 90 +17 -46 0 +19 -46 -66 +23 69 18 +26 69 -54 +29 -69 54 +32 -69 -18 +43 92 90 +50 92 -68 +52 -92 68 + diff --git a/easycap-M10_63_NO.txt b/easycap-M10_63_NO.txt new file mode 100644 index 0000000..66a76fa --- /dev/null +++ b/easycap-M10_63_NO.txt @@ -0,0 +1,64 @@ +Site Theta Phi +1 0 0 +2 23 90 +3 23 30 +4 23 -30 +5 -23 90 +6 -23 30 +7 -23 -30 +8 46 90 +9 46 66 +10 46 33 +11 46 0 +12 46 -33 +13 46 -66 +14 -46 90 +15 -46 66 +16 -46 33 +17 -46 0 +18 -46 -33 +19 -46 -66 +20 69 90 +21 69 66 +22 69 42 +23 69 18 +24 69 -6 +25 69 -30 +26 69 -54 +27 69 -78 +28 -69 78 +29 -69 54 +30 -69 30 +31 -69 6 +32 -69 -18 +41 -69 -42 +42 -69 -66 +43 92 90 +44 92 68 +45 92 45 +46 92 22 +47 92 0 +48 92 -22 +49 92 -45 +50 92 -68 +51 -92 90 +52 -92 68 +53 -92 45 +54 -92 22 +55 -92 0 +56 -92 -22 +57 -92 -45 +58 -92 -68 +59 115.000000000000 35 +60 115.000000000000 10 +61 115.000000000000 -15 +62 115.000000000000 -40 +63 115.000000000000 -65 +64 -115.000000000000 90 +65 -115.000000000000 65 +66 -115.000000000000 40 +67 -115.000000000000 15 +68 -115.000000000000 -10 +69 -115.000000000000 -35 +70 -135 -24 +71 135 24 diff --git a/features.py b/features.py new file mode 100755 index 0000000..0c84b51 --- /dev/null +++ b/features.py @@ -0,0 +1,7 @@ +#!/usr/bin/env python3 +# -*- coding: utf-8 -*- +""" +Created on Wed Jan 18 17:01:32 2023 + +@author: adamr +""" diff --git a/montage_18ch.csv b/montage_18ch.csv new file mode 100644 index 0000000..f7938b1 --- /dev/null +++ b/montage_18ch.csv @@ -0,0 +1,18 @@ +1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0 diff --git a/soft2.ui b/soft2.ui new file mode 100644 index 0000000..81c2ef4 --- /dev/null +++ b/soft2.ui @@ -0,0 +1,326 @@ + + + Settings + + + + 0 + 0 + 680 + 331 + + + + Settings + + + + + 11 + 11 + 330 + 309 + + + + + + + + + + measure 2 ax limits + + + + + + + + 7 + + + + y-axis caiibration + + + + + + + Apply + + + + + + + + + + measure 3 ax limits + + + + + + + + + + measure 4 ax limits + + + + + + + + + + measure 6 ax limits + + + + + + + + + + + + + + 7 + + + + x-axis calibration + y-axis readout + + + + + + + measure 1 ax limits + + + + + + + + + + + + + + + + + 11 + + + + plot limits + + + + + + + + + + + + + measure 5 ax limits + + + + + + + + + + + + + + 6 + + + + false + + + All values should be set like: [min, max] + + + true + + + + + + + + 7 + + + + Use auto y-axis limits for readout + + + + + + + + + 348 + 11 + 321 + 151 + + + + + 321 + 0 + + + + + + + + + + + 7 + + + + false + + + No file loaded + + + true + + + + + + + + 11 + + + + Montage + + + + + + + + 9 + + + + false + + + Select file with montage + + + true + + + + + + + + 9 + + + + false + + + File path: + + + true + + + + + + + Load + + + + + + + + + 348 + 169 + 321 + 130 + + + + + 321 + 0 + + + + + + + + + + Restart without applying + + + + + + + Apply montage and restart + + + + + + + Display max of each measurments after +calibration + + + + + + groupBox_3 + groupBox_2 + groupBox + + + +