-
Notifications
You must be signed in to change notification settings - Fork 0
/
BayesRate.py
134 lines (129 loc) · 4.38 KB
/
BayesRate.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
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
from astropy.stats import bayesian_blocks as bb
from numpy import histogram as hist
from numpy import diff
from datetime import datetime
import numpy as np
import matplotlib.pyplot as plt
from os import system as sys
import pickle
ttedata = []
#store_open = datetime(2024, 5, 18, 10, 30, 0, 0)
try:
with open("ttedata.pkl", 'rb') as file:
ttedata = pickle.load(file)
except:
print("Could not load TTE file...\n")
try:
with open("resume.pkl", 'rb') as file:
begin = pickle.load(file)
except:
print("Could not load last TTE for resuming...\n")
print("Using current system time as beginning...\n")
begin = datetime.now()
event = ""
p0 = 0.05 # default
while event != "q":
sys('cls')
print("Change-Point and Event Rate Tracker\n")
print("Plot Rate History (p)\n")
print("Display Statistics (d)\n")
print("Save times to file (s)\n")
print("Load times from file (l)\n")
print("View Event Times (t)\n")
print("Reset Everything (R)\n")
print("Remove last event (r)\n")
print("Adjust p0 (n)\n")
print("Quit Program (q)\n")
event = input("Press Enter to log a TTE>>")
if event == "R":
ttedata = []
begin = datetime.now()
with open("ttedata.pkl", 'wb') as file:
pickle.dump(ttedata, file)
with open("resume.pkl", 'wb') as file:
pickle.dump(datetime.now(), file)
if event == "l":
try:
with open("ttedata.pkl", 'rb') as file:
ttedata = pickle.load(file)
except:
print("Could not load TTE file...\n")
input("Press Enter to continue...\n")
try:
with open("resume.pkl", 'wb') as file:
pickle.dump(datetime.now(), file)
except:
print("Could not load save time...\n")
input("Press Enter to continue...\n")
if event == "s":
with open("ttedata.pkl", 'wb') as file:
pickle.dump(ttedata, file)
with open("resume.pkl", 'wb') as file:
pickle.dump(datetime.now(), file)
if event == "n":
sys('cls')
print("p0 = ", str(p0), "\n")
p0 = float(input("Enter new value: "))
if event == "t":
sys('cls')
try:
for n, e in enumerate(ttedata):
print("Event ", str(n+1), ": ", str(round(e/60.0, 3)))
input("Press Enter to continue...\n")
except:
print("Not enough data...\n")
if event == "d":
sys('cls')
ttedata.append((datetime.now()-begin).total_seconds())
#ttedata.append((datetime.now()-store_open).total_seconds())
try:
change_points = bb(ttedata, fitness='events', p0 = p0)
block_sizes = diff(change_points)
counts, _ = hist(ttedata, bins = change_points)
rates = counts / block_sizes
current_rate = round(60*rates[-1], 2)
print("Change-Points: \n")
for n, cp in enumerate(change_points, 1):
print("CP", str(n+1), ": ", str(round(cp/60.0, 3)), " minutes")
print("\nCurrent Rate: \n")
print(str(current_rate), "events/min\n")
ttedata.pop()
input("Press Enter to continue...\n")
except:
print("\nBayesian Blocks failed... (not enough data?)\n")
input("Press Enter to continue...\n")
if event == "":
ttedata.append((datetime.now()-begin).total_seconds())
#ttedata.append((datetime.now()-store_open).total_seconds())
with open('ttedata.pkl', 'wb') as file:
pickle.dump(ttedata, file)
if event == "p":
try:
ttedata.append((datetime.now()-begin).total_seconds())
#ttedata.append((datetime.now()-store_open).total_seconds())
change_points = bb(ttedata, fitness='events', p0 = p0)
block_sizes = diff(change_points)
counts, _ = hist(ttedata, bins = change_points)
rates = counts / block_sizes
current_rate = round(60*rates[-1], 2)
for n, cp in enumerate(change_points, 1):
print("CP", str(n+1), ": ", str(round(cp/60.0, 3)), " minutes")
print("\nCurrent Rate: ", str(current_rate), " events / min.\n")
plt.bar(x = change_points[0:-1]/60.0, height = 60*rates, width = block_sizes/60.0, align = 'edge', alpha = 0.5)
plt.vlines(change_points/60.0, 0.0, 60*max(rates), color = 'r', alpha = 0.5)
plt.title("Event Rate History")
plt.xlabel("Time (minutes)")
plt.ylabel("Events/Min")
plt.show()
ttedata.pop()
#this was not a real event, but a query on current rate if an event happened right now
except:
print("\nBayesian Blocks failed... (not enough data?)\n")
input("Press Enter to continue...\n")
if event == "r":
if len(ttedata) < 1:
break
else:
ttedata.pop()
with open('ttedata.pkl', 'wb') as file:
pickle.dump(ttedata, file)