From 74d284c15258a20be2320081cf4cbde8877649c4 Mon Sep 17 00:00:00 2001 From: Manuel Penschuck Date: Sat, 15 Nov 2014 20:53:19 +0100 Subject: [PATCH] CBMNet Dev. Tools --- cbmnet/calib.gp | 30 ++++++ cbmnet/calib.py | 43 ++++++++ cbmnet/dlm_jitter.gp | 15 +++ cbmnet/dlm_jitter.py | 92 +++++++++++++++++ cbmnet/dlm_jitter_block.gp | 9 ++ cbmnet/event_jitter.gp | 12 +++ cbmnet/event_jitter.py | 51 +++++++++ cbmnet/hldLib.py | 206 +++++++++++++++++++++++++++++++++++++ cbmnet/mbs.py | 58 +++++++++++ cbmnet/stats.py | 8 ++ 10 files changed, 524 insertions(+) create mode 100644 cbmnet/calib.gp create mode 100644 cbmnet/calib.py create mode 100644 cbmnet/dlm_jitter.gp create mode 100644 cbmnet/dlm_jitter.py create mode 100644 cbmnet/dlm_jitter_block.gp create mode 100644 cbmnet/event_jitter.gp create mode 100644 cbmnet/event_jitter.py create mode 100644 cbmnet/hldLib.py create mode 100644 cbmnet/mbs.py create mode 100644 cbmnet/stats.py diff --git a/cbmnet/calib.gp b/cbmnet/calib.gp new file mode 100644 index 0000000..8e586d0 --- /dev/null +++ b/cbmnet/calib.gp @@ -0,0 +1,30 @@ +set terminal pdf +set output "data/calib.pdf" + +set xlabel "Fine-Time" +set ylabel "Calibrated Fine-Time [ns]" +set key left; +set xrange [-1:512] +set yrange [-0.1 : 5.1] + +plot \ + 'data/calib.txt' using 1 with line title "Ch 0 Timing Trg", \ + 'data/calib.txt' using 2 with line title "Ch 1 DLM sensed", \ + 'data/calib.txt' using 3 with line title "Ch 2 Pulser", \ + 'data/calib.txt' using 4 with line title "Ch 3 Timing Trg (CBM)", \ + 'data/calib.txt' using 5 with line title "Ch 4 DLM Ref"; + +reset +set terminal pdf +set xrange [-1:512] + +set output "data/histogram.pdf" + +set ylabel "Frequency" + +plot \ + 'data/histogram.txt' using 1 with line title "Ch 0 Timing Trg", \ + 'data/histogram.txt' using 2 with line title "Ch 1 DLM sensed", \ + 'data/histogram.txt' using 3 with line title "Ch 2 Pulser", \ + 'data/histogram.txt' using 4 with line title "Ch 3 Timing Trg (CBM)", \ + 'data/histogram.txt' using 5 with line title "Ch 4 DLM Ref"; \ No newline at end of file diff --git a/cbmnet/calib.py b/cbmnet/calib.py new file mode 100644 index 0000000..818101b --- /dev/null +++ b/cbmnet/calib.py @@ -0,0 +1,43 @@ +import datetime +import glob +from hldLib import * +import numpy as np +from stats import * +import os + +noEvents = 10000000 +counter = 0 +channels = 5 +histogram = np.zeros( (channels, 1024), dtype="uint" ) + +def analyseCTSSSE(ssEvt, evt, counter): + global epochTimeInNs, eventTimeInCBM + global calibDataTDC + global idle + global netDLMChan, pinDLMChan + + cts, remain = extractCTS(ssEvt[1]) + sync, tdc = extractSync(remain) + assert(len(sync) in {1, 5, 11}) + + for w in tdc: + if w & 0x80000000: + data = tdcTimeData(w) + assert(data['channelNo'] < channels) + if data['fineTime'] != 0x3ff: + histogram[ data['channelNo'], data['fineTime'] ] += 1 + +files = glob.glob("/local/mpenschuck/hldfiles/*.hld") +files.sort() +print "Read file ", files[-2] +callbacks = {0x7005: analyseCTSSSE} +counter = iteratorOverHld(files[-1], callbacks, noEvents) + +print "Processed %d events, found %d fineTimes (%d in smallest channel)" % (counter, np.sum(histogram), np.min(np.sum(histogram, axis=1))) + +calib = 5.* np.cumsum(histogram, axis=1).astype('float') / np.sum(histogram, axis=1)[:,np.newaxis] +calib *= 5.0 / calib[:,-1].reshape(-1,1) + +np.save("calib.npy", calib) +np.savetxt("data/calib.txt", calib.T) +np.savetxt("data/histogram.txt", histogram.T) \ No newline at end of file diff --git a/cbmnet/dlm_jitter.gp b/cbmnet/dlm_jitter.gp new file mode 100644 index 0000000..3668a04 --- /dev/null +++ b/cbmnet/dlm_jitter.gp @@ -0,0 +1,15 @@ +set terminal png +set output "data/dlm_jitter.png" + +load "data/dlm_jitter.label" +set label labelText at graph 0.05,0.95 + +set xlabel "Deviation from average" +set ylabel "Frequency" + +set title "Jitter of DLM" + +binwidth=10 +bin(x,width)=width*floor(x/width) + +plot 'data/dlm_jitter.txt' using (bin($2*1000,binwidth)):(1.0) smooth freq with boxes notitle \ No newline at end of file diff --git a/cbmnet/dlm_jitter.py b/cbmnet/dlm_jitter.py new file mode 100644 index 0000000..e02af01 --- /dev/null +++ b/cbmnet/dlm_jitter.py @@ -0,0 +1,92 @@ +import datetime +import glob +from hldLib import * +import numpy as np +from stats import * +import os + +noEvents = 10000000 +counter = 0 +netDLMChan = 1 +pinDLMChan = 4 + +epoch = 0 +epochTimeInNs = 0.0 + +calibDataTDC = np.load("calib.npy") + +latencies = [[]] +netTimeBuf = [] +pinTimeBuf = [] +idle = 0 + +def analyseCTSSSE(ssEvt, evt, counter): + global epochTimeInNs, eventTimeInCBM + global calibDataTDC + global idle + global netDLMChan, pinDLMChan + global netTimeBuf, pinTimeBuf + + cts, remain = extractCTS(ssEvt[1]) + sync, tdc = extractSync(remain) + assert(len(sync) in {1, 5, 11}) + + channels, epochTimeInNs = interpretTDCData(tdc, calibDataTDC, epochTimeInNs) + + netTimeBuf += channels.get(netDLMChan, []) + pinTimeBuf += channels.get(pinDLMChan, []) + + if 0==len(channels.get(netDLMChan, [])) or 0==len(channels.get(pinDLMChan, [])): + if len(netTimeBuf) and len(pinTimeBuf): + idle += 1 + + if idle >= 1000: + netTime = np.array( netTimeBuf ) + + for pinTime in pinTimeBuf: + interval = np.min( np.abs( netTime - pinTime ) ) + if interval < 5000: + latencies[-1].append( interval ) + + idle = 0 + print len(pinTimeBuf), len(netTimeBuf) + print "start new block, old contained %d matches with avg %f" % (len(latencies[-1]), np.average(latencies[-1])) + latencies.append([]) + netTimeBuf = [] + pinTimeBuf = [] + + else: + idle = 0 + + +files = glob.glob("/local/mpenschuck/hldfiles/*.hld") +files.sort() +#files.pop() +#files.pop() +print "Read file ", files[-1] +callbacks = {0x7005: analyseCTSSSE} +counter = iteratorOverHld(files[-1], callbacks, noEvents) + +print "Processed %d events, found %d matches in %d blocks" % (counter, sum((len(x) for x in latencies)), len(latencies)) + +stats = np.zeros ( (len(latencies),4) ) +i = 0 +for block in latencies: + if len(block) < 500: continue + times = np.array( block ) + avgInterval = np.average(times) + stdInterval = np.std(times) + + stats[i] = [i, times.shape[0], avgInterval, stdInterval] + + print "Avg: %f ns, Std: %f ns" % (avgInterval, stdInterval) + np.savetxt("data/dlm_jitter.txt", np.vstack([times, times - avgInterval]).T ) + text_file = open("data/dlm_jitter.label", "w") + text_file.write('labelText = "Events: %d\\nAvg: %.2f ns\\nRMS: %.2f ps"' % (times.shape[0], avgInterval, stdInterval * 1000) ) + text_file.close() + + os.system("gnuplot dlm_jitter.gp") + i+=1 + os.system("mv data/dlm_jitter.png data/dlm_jitter%04d.png" % i) + +np.savetxt("data/dlm_jitter_blocks.txt", stats[:i]) diff --git a/cbmnet/dlm_jitter_block.gp b/cbmnet/dlm_jitter_block.gp new file mode 100644 index 0000000..13b381f --- /dev/null +++ b/cbmnet/dlm_jitter_block.gp @@ -0,0 +1,9 @@ +set terminal pdf +set output "data/dlm_jitter_blocks.pdf" + +set xlabel "Run" +set ylabel "Average DLM Latency [ns]" + +#set xrange [0:10] + +plot 'data/dlm_jitter_blocks.txt' using ($1+1):3:4 with errorbar title "DLM Latency" \ No newline at end of file diff --git a/cbmnet/event_jitter.gp b/cbmnet/event_jitter.gp new file mode 100644 index 0000000..0e0e395 --- /dev/null +++ b/cbmnet/event_jitter.gp @@ -0,0 +1,12 @@ +set terminal png +set output "data/event_jitter.png" + +load "data/event_jitter.label" +set label labelText at graph 0.05,0.95 + +set xlabel "Deviation from average [ps]" +set ylabel "Frequency" + +binwidth=10 +bin(x,width)=width*floor(x/width) +plot 'data/event_jitter.txt' using (bin($2*1000,binwidth)):(1.0) smooth freq with boxes notitle \ No newline at end of file diff --git a/cbmnet/event_jitter.py b/cbmnet/event_jitter.py new file mode 100644 index 0000000..6c48f75 --- /dev/null +++ b/cbmnet/event_jitter.py @@ -0,0 +1,51 @@ +from hldLib import * +import glob +import numpy as np + +# set-up constants and fetch calibration data +noEvents = 1000000 +calibDataTDC = np.load("calib.npy") +trbRefTimeTDCChan = 0 +cbmRefTimeTDCChan = 3 + +# interpret the HLD file +epochTimeInNs = 0.0 +eventTimeInCBM = [] +def analyseCTSSSE(ssEvt, evt, counter): + global epochTimeInNs, calibDataTDC, eventTimeInCBM + global trbRefTimeTDCChan, cbmRefTimeTDCChan + + cts, remain = extractCTS(ssEvt[1]) + sync, tdc = extractSync(remain) + assert(len(sync) in {1, 5, 11}) + + channels, epochTimeInNs = interpretTDCData(tdc, calibDataTDC, epochTimeInNs) + + if (not channels[trbRefTimeTDCChan] or \ + not channels[cbmRefTimeTDCChan] or \ + len(channels[trbRefTimeTDCChan]) != 1 or \ + len(channels[cbmRefTimeTDCChan]) != 1): + return + + eventTimeInCBM.append( 8.0 * sync[2] + channels[trbRefTimeTDCChan][0] - channels[cbmRefTimeTDCChan][0] ) + +# get file and iterate over it +files = glob.glob("/local/mpenschuck/hldfiles/ct14290170201.hld") +files.sort() +print "Read file ", files[-1] +callbacks = {0xf3c0: analyseCTSSSE} +counter = iteratorOverHld(files[-1], callbacks, noEvents) + +# calculate slope and jitter +print "Found %d events" % counter +eventTimeInCBM = np.array(eventTimeInCBM) +timeBetweenEvents = eventTimeInCBM[1:] - eventTimeInCBM[:-1] +avgInterval = np.average(timeBetweenEvents) +stdInterval = np.std(timeBetweenEvents) + +# dump values computed +print "Avg: %f ns, Std: %f ns" % (avgInterval, stdInterval) +np.savetxt("data/event_jitter.txt", np.vstack([timeBetweenEvents, timeBetweenEvents - avgInterval]).T ) +text_file = open("data/event_jitter.label", "w") +text_file.write('labelText = "Events: %d\\nAvg: %.2f ns\\nFreq: %.2f KHz\\nRMS: %.2f ps"\n' % (counter, avgInterval, 1e6 / avgInterval, stdInterval * 1000) ) +text_file.close() diff --git a/cbmnet/hldLib.py b/cbmnet/hldLib.py new file mode 100644 index 0000000..e3d6a49 --- /dev/null +++ b/cbmnet/hldLib.py @@ -0,0 +1,206 @@ +import struct +import datetime + +# iterates over a HLD files yieldling events of the following structure: +# evt = array( +# 0: dict(size, triggerType, evtNumFile, timestamp, runNum) +# 1: array( #of sub-events +# array( +# 0: dict(size, subEvtId, trgNum, trgCode) +# 1: array ( #of sub-sub-events +# array( +# 0: dict(size, subsubEvtId) +# 1: array( data-words ) +# ))))) +def eventIterator(f): + f.seek(0) + pos = 0 + + while(1): + # we might need to skip some padding + assert(f.tell() <= pos) + f.seek(pos) + + evtHdrBuf = f.read(8*4) + if (len(evtHdrBuf) != 8*4): + break #eof + + hdr = struct.unpack("<" + ("I" * 8), evtHdrBuf) + if (hdr[1] > 0xf0000): hdr = struct.unpack(">" + ("I" * 8), evtHdrBuf) + + evtSize = hdr[0] + + pos += evtSize + if (evtSize % 8): + pos += 8 - (evtSize % 8) + + # decode hdr + hdrDec = { + 'size': evtSize, + 'triggerType': hdr[2] & 0xffff, + 'evtNumFile': hdr[3], + 'timestamp': datetime.datetime( + (hdr[4] >> 16) & 0xff, (hdr[4] >> 8) & 0xff, (hdr[4] >> 0) & 0xff, + (hdr[5] >> 16) & 0xff, (hdr[5] >> 8) & 0xff, (hdr[5] >> 0) & 0xff ), + 'runNum': hdr[6] + } + + # load payload + subEvents = [] + innerPos = 8*4 + while(innerPos < evtSize): + subEvtHdrBuf = f.read(4 * 4) + endian = "<" + subEvtHdr = struct.unpack("<" + ("I" * 4), subEvtHdrBuf) + if (subEvtHdr[1] > 0xf0000): + subEvtHdr = struct.unpack(">" + ("I" * 4), subEvtHdrBuf) + endian = ">" + + subEvtSize = subEvtHdr[0] + subEvtHdrDec = { + 'size': subEvtSize, + 'subEvtId': subEvtHdr[2] & 0xffff, + 'trgNum': (subEvtHdr[3] >> 8) & 0xffffff, + 'trgCode': subEvtHdr[3] & 0xff + } + + subsubEvents = [] + ssPos = 4*4 + while(ssPos < subEvtSize): + sseHdrBuf = f.read(4) + sseHdr = struct.unpack(endian + "I", sseHdrBuf) + sseSize = ((sseHdr[0] >> 16) & 0xffff) * 4 + sseHdrDec = { + 'size': sseSize, + 'subsubEvtId': sseHdr[0] & 0xffff + } + + sseBuf = f.read(sseSize) + sseCont = struct.unpack(endian + ("I" * (sseSize/4)), sseBuf) + subsubEvents.append( (sseHdrDec, sseCont) ) + + ssPos += sseSize + 4 + + subEvents.append( (subEvtHdrDec, subsubEvents) ) + + innerPos += subEvtSize + + yield (hdrDec, subEvents) + +def dumpEvt(evt): + res = str(evt[0]) + "\n" + for subEvt in evt[1]: + res += " " + dumpSubEvt(subEvt).replace("\n", "\n ") + "\n" + return res + +def dumpSubEvt(sEvt): + h = sEvt[0] + res = "subEvtId: 0x%04x, trgNum: % 9d, trgCode: % 4d, size: % 4d\n" % (h['subEvtId'], h['trgNum'], h['trgCode'], h['size']) + for ssEvt in sEvt[1]: + res += " " + dumpSubSubEvt(ssEvt).replace("\n", "\n ") + "\n" + return res + +def dumpSubSubEvt(ssEvt): + res = "ID: 0x%04x, Size: %d\n" % (ssEvt[0]['subsubEvtId'], ssEvt[0]['size']) + res += dumpArray(ssEvt[1]) + return res + +def dumpArray(arr): + res = "" + for i in range(0, len(arr)): + if i != 0 and i % 8 == 0: res += "\n" + res += " [% 3d] 0x%08x" % (i+1, arr[i]) + + return res + +def iteratorOverHld(filename, cbs, eventNum = -1): + i = 0 + with open(filename, "rb") as hld: + for evt in eventIterator(hld): + if not len(evt[1]): continue + for sEvt in evt[1]: + for ssEvt in sEvt[1]: + ssid = ssEvt[0]['subsubEvtId'] + if ssid in cbs: + if callable( cbs[ssid] ): + cbs[ssid](ssEvt, evt, i) + else: + for func in cbs[ssid]: + func(ssEvt, evt, i) + + i += 1 + if (0 < eventNum <= i): + break + + return i + +def extractCTS(ssEvtData): + hdr = ssEvtData[0] + length = 1 + \ + ((hdr >> 16) & 0xf) * 2 + \ + ((hdr >> 20) & 0x1f) * 2 + \ + ((hdr >> 25) & 0x1) * 2 + \ + ((hdr >> 26) & 0x1) * 3 + \ + ((hdr >> 27) & 0x1) * 1 + + if (((hdr >> 28) & 0x3) == 0x1): length += 1 + if (((hdr >> 28) & 0x3) == 0x2): length += 4 + + return (ssEvtData[:length], ssEvtData[length:]) + +def extractSync(ssEvtData): + hdr = ssEvtData[0] + assert(hdr >> 28 == 0x1) + packType = (hdr >> 26) & 0x3 + length = 1 + if (packType == 1): length += 4 + if (packType == 3): length += 10 + + return (ssEvtData[:length], ssEvtData[length:]) + + +def extractTDC(ssEvtData): + length = ssEvtData[0] + assert(length >= len(ssEvtData)) + return (ssEvtData[:length+1], ssEvtData[length:]) + +def tdcTimeData(w): + if not (w & 0x80000000): + return None + + return { + "coarseTime": (w >> 0) & 0x7ff, + "edge": (w >> 11) & 1, + "fineTime": (w >> 12) & 0x3ff, + "channelNo": (w >> 22) & 0x7f + } + + +def interpretTDCData(data, calibDataTDC, epochTimeInNs = 0): + channels = {} + for w in data: + if w & 0x80000000: + tdcData = tdcTimeData(w) + chan = tdcData["channelNo"] + if tdcData["edge"] == 1:# and tdcData["fineTime"] != 0x3ff: + fineTimeInNs = calibDataTDC[tdcData["channelNo"], tdcData["fineTime"]] + assert( 0 <= fineTimeInNs <= 5.0 ) + coarseTimeInNs = tdcData["coarseTime"] * 5.0 + tdcTime = coarseTimeInNs - fineTimeInNs + epochTimeInNs + + if not chan in channels: channels[chan] = [] + channels[chan].append(tdcTime) + + + elif (w >> 29) == 3: # epoch counter + epoch = w & 0x0fffffff + epochTimeInNs = epoch * 10240.0 + + elif (w >> 29) == 1: # tdc header + pass + elif (w >> 29) == 2: # debug + pass + else: + print "Unknown TDC word type: 0x%08x" % w + + return channels, epochTimeInNs diff --git a/cbmnet/mbs.py b/cbmnet/mbs.py new file mode 100644 index 0000000..bb84118 --- /dev/null +++ b/cbmnet/mbs.py @@ -0,0 +1,58 @@ +import datetime +import glob +from hldLib import * +import numpy as np +from stats import * +import os + +noEvents = 1000000 +counter = 0 + +timestamps = [] # in us +lastId = None +trgOffset = 0 +lastRegNum = None + +def analyseMBS(ssEvt, evt, counter): + global timestamps + global lastRegNum, lastId, trgOffset + + regNum = ssEvt[1][0] & 0x00ffffff + + if lastRegNum == regNum: return + lastRegNum = regNum + + inclTime = ssEvt[1][0] & 0x01000000 + error = ssEvt[1][0] & 0x80000000 + + if (error): print "Error in regNum", regNum + assert(inclTime) + + trgId = trgOffset + (ssEvt[1][2]&0xffff) + if lastId != None and lastId > trgId: + trgOffset += 0x10000 + trgId += 0x10000 + + lastId = trgId + #print ssEvt[1][1]*5e-3 + + ts=20*trgId + ssEvt[1][1]*5e-3 + #print ts + + timestamps.append( ts ) + + +files = glob.glob("/local/mpenschuck/hldfiles/*.hld") +files.sort() +print "Read file ", files[-1] +callbacks = {0xf30a: analyseMBS} +counter = iteratorOverHld(files[-1], callbacks, noEvents) + +print "Processed %d events, found %d timestampts" % (counter, len(timestamps)) + +timestamps = np.array(timestamps) +sTimestamps = slope(timestamps) + +dev = sTimestamps - np.average(sTimestamps) + +print 1e3/30.52, np.average(sTimestamps), np.std(sTimestamps), np.max(dev) diff --git a/cbmnet/stats.py b/cbmnet/stats.py new file mode 100644 index 0000000..fd20254 --- /dev/null +++ b/cbmnet/stats.py @@ -0,0 +1,8 @@ +import numpy as np + +def slopeWithoutAvg(X): + slope = X[1:]-X[:-1] + return slope - np.average(slope) + +def slope(X): + return X[1:] - X[:-1] \ No newline at end of file -- 2.43.0