PhotonStatsAnalysis/SlidingFrame.py at main · JCodeZero/PhotonStatsAnalysis · GitHub

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
import sys
from astropy.io import fits
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np

filename = sys.argv[1]
eventfile=fits.open(filename)
events=eventfile[1].data
shortdat=events[(events['TIME']>660920900) & (events['TIME']<660921200)]
cleanevts = shortdat[(shortdat['EVENT_FLAGS']==0) & (shortdat['ENERGY']>15)]
filterdat = shortdat[(shortdat['EVENT_FLAGS']==0) & (shortdat['ENERGY']>15)]

time = []

df = pd.DataFrame(filterdat, columns = ['TIME','DET_ID','EVENT_FLAGS','PHA','MASK_WEIGHT','DETX','DETY','PI','ENERGY'])
df['TIME'] = df['TIME'].apply(lambda x: round(x))
time.extend(df['TIME'].tolist())

photons = np.empty((0,1), int)
temp = 0

timearray = np.array(list(range(time[0], time[-1])))


for i in timearray:
    for j in range(len(time)):
        if time[j] == i:
            temp = temp + 1
    photons = np.append(photons, np.array([[temp]]))
    temp = 0

length = 5
edgesx = np.empty((0,1))
means = np.empty((0,1))
stds = np.empty((0,1))
topy = np.empty((0,1))
bottomy = np.empty((0,1))
for i in range(int((time[-1] - time[0])/length)):
    edgesx = np.append(edgesx, np.array([[time[0] + length*i]]))
    means = np.append(means, np.array([[np.mean(photons[length*i:length*(i+1)])]]))
    std = np.std(photons[length*i:length*(i+1)])
    stds = np.append(stds, np.array([[std]]))
    topy = np.append(topy, np.array([[np.mean(photons[length*i:length*(i+1)]) + 3 * std]]))
    bottomy = np.append(bottomy, np.array([[np.mean(photons[length*i:length*(i+1)]) - 3 * std]]))

edgesx = np.append(edgesx, np.array([[edgesx[-1] + length]]))
means = np.append(means, np.array([[means[-1]]]))
stds = np.append(stds, np.array([[stds[-1]]]))
topy = np.append(topy, np.array([[topy[-1]]]))
bottomy = np.append(bottomy, np.array([[bottomy[-1]]]))

intervals = pd.DataFrame(columns = ['START', 'END', 'COUNTMEAN'])
maxstart = []
maxend = []
maxmean = []
for i in range(len(timearray)):
    tempstartlist = []
    tempendlist = []
    tempmeanlist = []

    for j in range(0, time[-1] - time[0]):
        if j+i+1 <= 300:
            tempstartlist.append(j)
            tempendlist.append(j+i+1)
            tempmeanlist.append(np.mean(photons[j:j+i+1]))

    for i in range(len(tempstartlist)):
        if tempmeanlist[i] == np.max(tempmeanlist):
            maxstart.append(tempstartlist[i])
            maxend.append(tempendlist[i])
            maxmean.append(tempmeanlist[i])

intervals['START'] = maxstart
intervals['END'] = maxend
intervals['COUNTMEAN'] = maxmean

freqstart = np.argmax(np.bincount(intervals['START']))
freqend = 0

for i in range(len(maxstart)):
    if intervalranges['START'][i] == freqstart:
        if intervals['ErangeND'][i] > freqend:
            freqend = intervals['END'][i]

slidemethod = filterdat[(filterdat['TIME']>timearray[freqstart]) & (filterdat['TIME']<timearray[freqend])]

_=plt.hist(slidemethod['TIME'], bins=300)
plt.xlabel('Time (s)')
plt.ylabel('# of photons')
plt.show()