Compare commits

..

No commits in common. "main" and "minor_fixes" have entirely different histories.

5 changed files with 59 additions and 332 deletions

View file

@ -1,28 +1,13 @@
#! /usr/bin/make -f
MUSS_OPT=-h
MUSC_OPT=-h -c2
%.MUSS.bz2: %.dat.xz
xzcat $< \
| ./nm64file -n 32 \
| timeout -v 300 python3 mus_sync.py $(MUSS_OPT) \
| nm64file -n 32 \
| timeout -v 300 python3 mus_sync.py -h \
| bzip2 -9c > $@
%.MUSS.bz2: %.dat
cat $< \
| ./nm64file -n 32 \
| timeout -v 300 python3 mus_sync.py $(MUSS_OPT) \
| nm64file -n 32 \
| timeout -v 300 python3 mus_sync.py -h \
| bzip2 -9c > $@
%.MUSC: %.dat.xz
xzcat $< \
| ./nm64file -n 32 \
| timeout -v 300 python3 mus_sync.py $(MUSC_OPT) \
> $@
%.MUSC: %.dat
cat $< \
| ./nm64file -n 32 \
| timeout -v 300 python3 mus_sync.py $(MUSC_OPT) \
> $@

View file

@ -1,82 +0,0 @@
#! /usr/bin/python3
import sys
import numpy as np
from matplotlib import pyplot as plt
verbose=True
f = sys.stdin
l=[]
for line in f:
if 'EOF' in line:
continue
if line[0] == 'E':
l.append(list(map(int,line.split()[5:])))
data=np.array(l)
if verbose:
print("data.shape: ", data.shape, file=sys.stderr)
data=data.reshape(data.shape[0],18,3)
print(data[0])
lmus=[]
lnichtmus=[]
for line in data:
if line[5,1]<2 or line[6,1]<2 or line[9,1]<2 or line[10,1]<2:
continue
else:
if line[0,1]>=2:
lmus.append(line)
else:
lnichtmus.append(line)
lnichtmus = np.array(lnichtmus)
i=0
for line in lnichtmus[:,:,0]:
print (line)
i+=1
if i>200:
break
l=[]
for line in lnichtmus[:,5,0]:
l.append(line)
for line in lnichtmus[:,6,0]:
l.append(line)
for line in lnichtmus[:,9,0]:
l.append(line)
for line in lnichtmus[:,10,0]:
l.append(line)
l=np.array(l)
plt.hist(l/14000,bins=100, range=(0.000000001,73))
plt.xlabel('Pulseheight in mV')
plt.ylabel('Counts')
plt.show()
'''
tab = f.read()
f.close
#print(tab)
tab = tab.split('\n')
h = [x for x in tab if 'E' in x]
#print(h)
l = []
for i in h :
l.append([i])
print(l)
'''

View file

@ -1,21 +0,0 @@
import numpy as np
import getopt
import fileinput
import sys
opt,files = getopt.getopt(sys.argv[1:], 'v')
verbose=False
for o,v in opt:
if o=='-v':
verbose=True
data=[]
for line in fileinput.input(files):
if line[0:2]=='EI':
l=line.split()
if int(l[6])>0:
print(l[5],l[6],l[7])

View file

@ -8,9 +8,6 @@ All lines after the last C64 line in the .MUS file are discarded, if the script
If executed with otion -h this script starts writing at the first C64 line in the hour that startet in the .MUS file. The last
line printed is the last EMUS line, before the first C64 line in the hour that started approx one hour after the start of the file.
Thus there are no redundat lines, and no missing lines, when concatenating multiple outputfiles of this script.
Option -v writes verbosity in sys.stderr
Option -c1 adds Coincidences to the outputfile after the housekeeping blocks
Option -c2 has the ame output as -c1 but all EMUS lines are not in the outputfile
'''
import numpy as np
import getopt
@ -19,22 +16,14 @@ import sys
import datetime
import time
opt,files = getopt.getopt(sys.argv[1:], 'hvc:')
opt,files = getopt.getopt(sys.argv[1:], 'h')
hourly=False
verbose=False
coincopt=False
for o,v in opt:
if o=='-h':
hourly=True
if o=='-v':
verbose=True
if o=='-c':
coincopt=int(v)
print(f'FMUS {" ".join(sys.argv)}')
def fixclockflip(data):
clockflipch1=0
@ -46,17 +35,11 @@ def fixclockflip(data):
line[2]+=clockflipch1*4294967296
if line[2]<tb1-1000000:
clockflipch1+=1
line[2]+=4294967296
tb1=line[2]
continue
tb1=line[2]
else:
line[2]+=clockflipch2*4294967296
if line[2]<tb2-1000000:
clockflipch2+=1
line[2]+=4294967296
tb2=line[2]
continue
tb2=line[2]
return data
@ -74,40 +57,28 @@ def get_first_t(data,ch):
i+=1
return data[i][2],i
def get_last_t(data,ch):
i=1
if ch==1:
while True:
if data[-i][0]<16:
break
i+=1
if ch==2:
while True:
if data[-i][0]>15:
break
i+=1
return data[-i][2],i
def count_coinc(data):
lastline=np.zeros(3)
n=data.shape[0]-1
delta_times=np.zeros(n)
count=0
for line in data:
if -0.5<line[2]-lastline[2]<0.5 and ((lastline[0] <=15 and line[0]>=16) or (line[0] <=15 and lastline[0]>=16)):
for i in range (n):
delta_times[i]=data[i+1][2]-data[i][2]
for t in delta_times:
if t>=-1.5 and t<=1.5:
count+=1
lastline=line
return count
def approx_drift(data):
offset1=0
offset2=0
lastt=0
for t in data[:5000,2]:
for t in data[:2000,2]:
if t-lastt<offset1:
offset1=t-lastt
lastt=t
offset1=-offset1
lastt=0
for t in data[-5000:,2]:
for t in data[-2000:,2]:
if t-lastt<offset2:
offset2=t-lastt
lastt=t
@ -115,24 +86,18 @@ def approx_drift(data):
slope=(offset2-offset1)/(get_first_t(data,1)[0]-get_first_t(data[-2000:,:],1)[0])
if abs(slope)>4e-6 or abs(slope)<2e-7:
if found_weird_lines:
raise Exception('something went wrong, getting the approximate drift of the two clocks. There were Lines in the Data, that seem strange.'+str(slope))
raise Exception('something went wrong, getting the approximate drift of the two clocks. There were Lines in the Data, that seem strange.')
else:
raise Exception('something went wrong, getting the approximate drift of the two clocks.'+str(slope))
raise Exception('something went wrong, getting the approximate drift of the two clocks.')
return slope,offset1
def fixtimes(data,slope,offset,firsttime):
def fixtimes(data,slope,offset):
fixing_times=True
highest_count=0
additional=0
a=0
t1=get_first_t(data,unit)[0]
t2=get_first_t(data,runit)[0]
if abs(abs(t1-t2)-abs(offset))>1000000:
if found_weird_lines:
raise Exception('Could not find a Correction of the Times. There may be a problem in the data. No hits on one Unit for more than 1s less than 2000 lines after EMUS '+str(data[0][0])+' '+str(data[0][1])+' '+str(data[0][2])+'. There were Lines in the Data, that seem strange.')
else:
raise Exception('Could not find a Correction of the Times. There may be a problem in the data, No hits on one Unit for more than 1s less than 2000 lines after EMUS '+str(data[0][0])+' '+str(data[0][1])+' '+str(data[0][2])+'.')
te=get_last_t(data,runit)[0]
for line in data:
if unit==1:
if line[0]>15:
@ -148,29 +113,28 @@ def fixtimes(data,slope,offset,firsttime):
modi=a%2
highest_count=count
additional=0
offset=te-get_last_t(data,runit)[0]
offset=t2-get_first_t(data,runit)[0]
if unit==1:
offset=-offset
if firsttime and a>2000 and a%2==modi:
if highest_count<15:
if found_weird_lines:
raise Exception('A suspiciously low amount of coincidences between the two units were detected. The script did not find a better fit. There were Lines in the Data, that seem strange')
else:
raise Exception('A suspiciously low amount of coincidences between the two units were detected. The script did not find a better fit.')
firsttime=False
if count*1.1<highest_count and a%2==modi:
break
if (not firsttime) and a>=3 and a%2==modi:
break
a+=1
for line in data:
if unit==1:
if line[0]>15:
line[2]-=a*(-1)**(a%2)
line[2]+=a*(-1)**(a%2)
if unit==2:
if line[0]<16:
line[2]+=a*(-1)**(a%2)
line[2]-=a*(-1)**(a%2)
additional+=1
a+=1
if a>20000:
if found_weird_lines:
raise Exception('Could not find a Correction of the Times. There may be a problem in the data, the script could not work around. There were Lines in the Data, that seem strange')
else:
raise Exception('Could not find a Correction of the Times. There may be a problem in the data, the script could not work around.')
toomuch=additional/2
for line in data:
if unit==1:
if line[0]>15:
@ -179,7 +143,7 @@ def fixtimes(data,slope,offset,firsttime):
if line[0]<16:
line[2]-=(-1)**(a%2)*toomuch
return offset,data,firsttime
return offset,data
unit=None
@ -190,9 +154,7 @@ ctime=0
lasttu1=0
lasttu2=0
found_weird_lines=False
if verbose:
st=time.time()
startt=time.time()
#startt=time.time()
for line in fileinput.input(files):
if 'EOF' in line:
continue
@ -200,29 +162,17 @@ for line in fileinput.input(files):
numbers=line.split()[1:4]
if int(numbers[0])<16:
dt=int(numbers[2])-lasttu1
if 8<int(numbers[1])<11 and(lasttu1==0 or abs(dt)<1000000 or 2**32-1000000<int(abs(dt))<2**32+1000000):
if 8<int(numbers[1])<11 and(lasttu1==0 or abs(dt)<1000000 or 2**32-1000000<int(abs(dt))<2**32+1000000):#wahrscheinlichkeit 1s ohne event in der selben unit ca. 1e-300
datalist.append([int(numbers[0]), int(numbers[1]), float(numbers[2])])
lasttu1=int(numbers[2])
false1=0
elif lasttu1!=0 and not (abs(dt)<1000000 or 2**32-1000000<int(abs(dt))<2**32+1000000):
if verbose :
sys.stderr.write(line)
false1+=1
if false1>10000:
raise Exception('Unit 1 did not write data for more than one second.')
else:
found_weird_lines=True
else:
dt=int(numbers[2])-lasttu2
if 8<int(numbers[1])<11 and(lasttu2==0 or abs(dt)<1000000 or 2**32-1000000<abs(dt)<2**32+1000000):
datalist.append([int(numbers[0]), int(numbers[1]), float(numbers[2])])
lasttu2=int(numbers[2])
false2=0
elif lasttu2!=0 and not (abs(dt)<1000000 or 2**32-1000000<abs(dt)<2**32+1000000):
if verbose:
sys.stderr.write(line)
false2+=1
if false2>10000:
raise Exception('Unit 2 did not write data for more than one second.')
else:
found_weird_lines=True
elif any(char.isdigit() for char in line):
if line[0:3]=='C64':
@ -245,92 +195,47 @@ for line in fileinput.input(files):
otherlines[ctime].append(line)
data=np.array(datalist)
keys=list(otherlines.keys())
keys.sort()
if len(keys)<60:
sys.exit()
if verbose:
sys.stderr.write('reading data took '+str(time.time()-startt)+'s\n')
startt=time.time()
#sys.stderr.write('reading data took '+str(time.time()-startt)+'s\n')
#startt=time.time()
data=fixclockflip(data)
if verbose:
sys.stderr.write('fixing clockflip took '+str(time.time()-startt)+'s\n')
startt=time.time()
#sys.stderr.write('fixing clockflip took '+str(time.time()-startt)+'s\n')
#startt=time.time()
slope,offset=approx_drift(data)
if verbose:
sys.stderr.write(str((slope,offset)))
sys.stderr.write('getting approximate drift took '+str(time.time()-startt)+'s\n')
startt=time.time()
#sys.stderr.write('getting approximate drift took '+str(time.time()-startt)+'s\n')
#startt=time.time()
def fat(data,slope,offset,biggertu1):
fixing_times=True
i=0
firsttime=True
if biggertu1:
while fixing_times:
if i<data.shape[0]-4001:
if firsttime:
offset,data[i:i+10000],firsttime=fixtimes(data[i:i+10000],slope,offset,firsttime)
i+=10000
else:
offset,data[i:i+2000],firsttime=fixtimes(data[i:i+2000],slope,offset,firsttime)
i+=2000
if i<data.shape[0]-2001:
offset,data[i:i+2000]=fixtimes(data[i:i+2000],slope,offset-2)
i+=2000
else:
offset,data[i:],firsttime=fixtimes(data[i:],slope,offset,firsttime)
offset,data[i:]=fixtimes(data[i:],slope,offset)
fixing_times=False
else:
offset=-offset
while fixing_times:
if i<data.shape[0]-4001:
if firsttime:
offset,data[i:i+10000],firsttime=fixtimes(data[i:i+10000],-slope,offset,firsttime)
i+=10000
else:
offset,data[i:i+2000],firsttime=fixtimes(data[i:i+2000],-slope,offset,firsttime)
i+=2000
if i<data.shape[0]-2001:
offset,data[i:i+2000]=fixtimes(data[i:i+2000],-slope,offset-2)
i+=2000
else:
offset,data[i:],firsttime=fixtimes(data[i:],-slope,offset,firsttime)
offset,data[i:]=fixtimes(data[i:],-slope,offset-2)
fixing_times=False
biggertu1=get_first_t(data,1)[0]>get_first_t(data,2)[0]
fat(data,slope,offset,biggertu1)
if verbose:
sys.stderr.write('fixing times took '+str(time.time()-startt)+'s\n')
startt=time.time()
data=data[np.argsort(data[:,2])]
if verbose:
sys.stderr.write('sorting data took '+str(time.time()-startt)+'s\n')
startt=time.time()
if coincopt:
k=0
coinces=np.zeros((32,32))
lastline=data[0,:]
lastline[2]-=10
llastline=data[0,:]
llastline[2]-=10
lllastline=data[0,:]
lllastline[2]-=10
for line in data:
if k>=len(keys):
break
if line[2]>keys[k]:
for i in [8,9,10,11,12,13,14,15,24,25,26,27,28,29,30,31]:
for j in [0,1,2,3,4,5,6,7,16,17,18,19,20,21,22,23]:
c=coinces[i,j]+coinces[j,i]
otherlines[keys[k-1]].append('CMUS '+str(i)+' '+str(j)+' '+str(int(c))+'\n')
coinces=np.zeros((32,32))
k+=1
if lastline[2]-llastline[2]<2 and llastline[2]-lllastline[2]>2 and line[2]-lastline[2]>2:
coinces[int(llastline[0]),int(lastline[0])]+=1
lllastline=llastline
llastline=lastline
lastline=line
if verbose:
sys.stderr.write('finding coinces took '+str(time.time()-startt)+'s\n')
startt=time.time()
#sys.stderr.write('fixing times took '+str(time.time()-startt)+'s\n')
#startt=time.time()
data=data[np.argsort(data[:,2])]
#sys.stderr.write('sorting data took '+str(time.time()-startt)+'s\n')
#startt=time.time()
keys=list(otherlines.keys())
keys.sort()
def printhour(data,keys):
i=0
@ -352,7 +257,7 @@ def printhour(data,keys):
if writing:
print(''.join(otherlines[keys[i]]))
i+=1
if writing and coincopt!=2:
if writing:
print('EMUS',int(line[0]),int(line[1]),f'{line[2]:.2f}')
return
@ -366,9 +271,6 @@ else:
if line[2]>keys[i]:
print(''.join(otherlines[keys[i]]))
i+=1
if coincopt!=2:
print('EMUS',int(line[0]),int(line[1]),f'{line[2]:.2f}')
if verbose:
sys.stderr.write('printing fixed data took '+str(time.time()-startt)+'s\n')
dt=time.time()-st
sys.stderr.write(str(files)+' done in '+f'{dt:.2f}'+'s\n')
print('EMUS',int(line[0]),int(line[1]),f'{line[2]:.2f}')
#sys.stderr.write('printing fixed data took '+str(time.time()-startt)+'s\n')
sys.stderr.write(str(files)+'done\n')

View file

@ -1,57 +0,0 @@
import numpy as np
import getopt
import fileinput
import sys
import matplotlib.pyplot as plt
opt,files = getopt.getopt(sys.argv[1:], 'b:c:t:a')
bins=None
col=0
allmus=False
for o,v in opt:
if o=='-b':
bins=int(v)
if o=='-c':
col=int(v)
if o=='-a':
allmus=True
data=[]
with open(files[0], "r", encoding="us-ascii", errors="ignore") as file:
for line in file:
if line[0:2]=='EI':
try:
data.append(list(map(int,line.split()[5:])))
except:
continue
data=np.array(data)
muonmus=[]
for line in data:
if (line[16]>2 and line[19]>2) and (line[28]>2 and line[31]>2) and line[15]>200000 and line[18]>245000 and line[27]>185000 and line[30]>208000:
muonmus.append(line[col])
muonmus=np.array(muonmus)
if allmus:
hist,bins,_=plt.hist(data[:,0],bins=bins,histtype='step',density=True)
else:
if col==15 or col==18 or col==27 or col==30:
hist,bins,_=plt.hist(muonmus,bins=bins,range=[0,0.75e6],histtype='step',density=True)
else:
hist,bins,_=plt.hist(muonmus,bins=bins,histtype='step',density=True)
dateiname=files[0]+'col'+str(col)+'.hist'
np.savez(dateiname,hist=hist,bins=bins)
#plt.xlabel('A')
#plt.ylabel('Density')
#plt.show()