forked from MelancholeyLemon/EasyNanoporeCode
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy patheventDetect.py
117 lines (102 loc) · 5.57 KB
/
eventDetect.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
# -*- coding : utf-8 -*-
# coding: utf-8
import numpy as np
import scipy.signal
def eventDownFast(rawSignal, startCoeff, endCoeff, filterCoeff, minDuration, maxDuration, fileName, sampleRate, fileNumber, iterNumber):
padLen = np.int64(sampleRate)
prepadded = np.ones(padLen) * np.mean(rawSignal[0:1000])
signalToFilter = np.concatenate((prepadded, rawSignal))
rawSignal = np.array(rawSignal)
mlTemp = scipy.signal.lfilter([1 - filterCoeff, 0], [1, -filterCoeff], signalToFilter)
vlTemp = scipy.signal.lfilter([1 - filterCoeff, 0], [1, -filterCoeff], np.square(signalToFilter - mlTemp))
ml = np.delete(mlTemp, np.arange(padLen))
vl = np.delete(vlTemp, np.arange(padLen))
sl = ml - startCoeff * np.sqrt(vl)
Ni = len(rawSignal)
points = np.array(np.where(rawSignal <= sl)[0])
to_pop = np.array([])
for i in range(1, len(points)):
if points[i] - points[i - 1] == 1:
to_pop = np.append(to_pop, i)
to_pop = np.int64(to_pop)
points = np.unique(np.delete(points, to_pop))
NumberOfEvents = 0;
RoughEventLocations = np.zeros((10000, 3))
event_current = np.zeros((10000, 3))
for i in points:
event_start_mean = ml[i-100]
if i >= Ni - 10:
break;
start = i
El = ml[i] - endCoeff * np.sqrt(vl[i])
while rawSignal[i + 1] < El and i <= Ni - 10:
i = i + 1
if ((minDuration * sampleRate / 1000) < (i + 1 - start)) and ((i + 1 - start) < (maxDuration * sampleRate / 1000)) and (event_start_mean - min(rawSignal[start:i + 1])) > 0:
NumberOfEvents = NumberOfEvents + 1
RoughEventLocations[NumberOfEvents - 1, 2] = i + 1 - start
RoughEventLocations[NumberOfEvents - 1, 0] = start
RoughEventLocations[NumberOfEvents - 1, 1] = i + 1
event_current[NumberOfEvents - 1, 0] = (event_start_mean - min(rawSignal[start:i + 1]))
event_current[NumberOfEvents - 1, 1] = event_start_mean
event_current[NumberOfEvents - 1, 2] = abs(
event_current[NumberOfEvents - 1, 0] / event_current[NumberOfEvents - 1, 1]) * 10000
event_statistic = np.zeros((NumberOfEvents, 6))
event_statistic[:, 0] = RoughEventLocations[0: NumberOfEvents, 2] / sampleRate * 1000
event_statistic[:, 1:4] = event_current[0: NumberOfEvents, 0:3]
if iterNumber == 0:
event_statistic[:, 4] = (RoughEventLocations[0: NumberOfEvents, 0] + iterNumber * sampleRate * 10) * 1000 / sampleRate
else:
event_statistic[:, 4] = (RoughEventLocations[0: NumberOfEvents,
0] + iterNumber * sampleRate * 10) * 1000 / sampleRate -10
event_statistic[:, 5] = fileNumber
with open(fileName, "a+") as fp:
np.savetxt(fp, event_statistic, fmt='%.3f', delimiter="\t")
def eventUpFast(rawSignal, startCoeff, endCoeff, filterCoeff, minDuration, maxDuration, fileName, sampleRate, fileNumber, iterNumber):
padLen = np.int64(sampleRate)
prepadded = np.ones(padLen) * np.mean(rawSignal[0:1000])
signalToFilter = np.concatenate((prepadded, rawSignal))
rawSignal = np.array(rawSignal)
mlTemp = scipy.signal.lfilter([1 - filterCoeff, 0], [1, -filterCoeff], signalToFilter)
vlTemp = scipy.signal.lfilter([1 - filterCoeff, 0], [1, -filterCoeff], np.square(signalToFilter - mlTemp))
ml = np.delete(mlTemp, np.arange(padLen))
vl = np.delete(vlTemp, np.arange(padLen))
sl = ml + startCoeff * np.sqrt(vl)
Ni = len(rawSignal)
points = np.array(np.where(rawSignal >= sl)[0])
to_pop = np.array([])
for i in range(1, len(points)):
if points[i] - points[i - 1] == 1:
to_pop = np.append(to_pop, i)
to_pop = np.int64(to_pop)
points = np.unique(np.delete(points, to_pop))
NumberOfEvents = 0
RoughEventLocations = np.zeros((10000, 3))
event_current = np.zeros((10000, 3))
for i in points:
event_start_mean = ml[i-100]
if i >= Ni - 10:
break;
start = i
El = ml[i] + endCoeff * np.sqrt(vl[i])
while rawSignal[i + 1] > El and i <= Ni - 10:
i = i + 1
if ((minDuration * sampleRate / 1000) < (i + 1 - start)) and ((i + 1 - start) < (maxDuration * sampleRate / 1000)) and (max(rawSignal[start:i + 1]) - event_start_mean) > 0:
NumberOfEvents = NumberOfEvents + 1
RoughEventLocations[NumberOfEvents - 1, 2] = i + 1 - start
RoughEventLocations[NumberOfEvents - 1, 0] = start
RoughEventLocations[NumberOfEvents - 1, 1] = i + 1
event_current[NumberOfEvents - 1, 0] = (max(rawSignal[start:i + 1]) - event_start_mean)
event_current[NumberOfEvents - 1, 1] = event_start_mean
event_current[NumberOfEvents - 1, 2] = abs(
event_current[NumberOfEvents - 1, 0] / event_current[NumberOfEvents - 1, 1]) * 10000
event_statistic = np.zeros((NumberOfEvents, 6))
event_statistic[:, 0] = RoughEventLocations[0: NumberOfEvents, 2] / sampleRate * 1000
event_statistic[:, 1:4] = event_current[0: NumberOfEvents, 0:3]
if iterNumber == 0:
event_statistic[:, 4] = (RoughEventLocations[0: NumberOfEvents, 0] + iterNumber * sampleRate * 10) * 1000 / sampleRate
else:
event_statistic[:, 4] = (RoughEventLocations[0: NumberOfEvents,
0] + iterNumber * sampleRate * 10) * 1000 / sampleRate -10
event_statistic[:, 5] = fileNumber
with open(fileName, "a+") as fp:
np.savetxt(fp, event_statistic, fmt='%.3f', delimiter="\t")