linse/frappy
3
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

take over changes from ultrasound PC

Browse Source 
Change-Id: I1eae717a5963e618d87ddf52db991d428a046d24
This commit is contained in:
zolliker 2025-01-09 10:47:23 +01:00
parent 4397d8db1a
commit a25a368491
4 changed files with 184 additions and 72 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

67
cfg/PEUS.py Normal file
View File

@ -0,0 +1,67 @@
Node(equipment_id = 'pe_ultrasound.psi.ch',
description = 'pulse echo ultra sound setup',
interface = 'tcp://5000',
)
Mod('f',
cls = 'frappy_psi.ultrasound.Frequency',
description = 'ultrasound frequency and acquisition loop',
uri = 'serial:///dev/ttyS1',
pars = 'pars',
pollinterval = 0.1,
time = 900, # start time
size = 5000,
freq = 1.17568e+06,
basefreq = 4.14902e+07,
control = False,
rusmode = False,
amp = 5.0,
nr = 1000, #500 #300 #100 #50 #30 #10 #5 #3 #1 #1000 #500 #300 #100 #50 #30 #10 #5 #3 #1 #500
sr = 32768, #16384
plot = True,
maxstep = 100000,
bw = 10E6, #butter worth filter bandwidth
maxy = 0.7, # y scale for plot
curves = 'curves', # module to transmit curves:
)
Mod('curves',
cls = 'frappy_psi.ultrasound.Curves',
description = 't, I, Q and pulse arrays for plot',
)
Mod('delay',
cls = 'frappy__psi.dg645.Delay',
description = 'delay line with 2 channels',
uri = 'serial:///dev/ttyS2',
on1 = 1e-9,
on2 = 1E-9,
off1 = 400e-9,
off2 = 600e-9,
)
Mod('pars',
cls = 'frappy_psi.ultrasound.Pars',
description = 'SEA parameters',
)
def roi(nr, time=None, size=300):
Mod(f'roi{nr}',
cls = 'frappy_psi.ultrasound.Roi',
description = f'I/Q of region {nr}',
main = 'f',
time=time or 4000,
size=size,
enable=time is not None,
)
roi(0, 2450) # you may add size as argument if not default
roi(1, 5950)
roi(2, 9475)
roi(3, 12900)
roi(4, 16100)
roi(5) # disabled
roi(6)
roi(7)
roi(8)
roi(9)

62
cfg/RUS.py Normal file
View File

@ -0,0 +1,62 @@
Node(equipment_id = 'r_ultrasound.psi.ch',
description = 'resonant ultra sound setup',
interface = 'tcp://5000',
)
Mod('f',
cls = 'frappy_psi.ultrasound.Frequency',
description = 'ultrasound frequency and acquisition loop',
uri = 'serial:///dev/ttyS1',
pars = 'pars',
pollinterval = 0.1,
time = 900, # start time
size = 5000,
freq = 1.e+03,
basefreq = 1.E+3,
control = False,
rusmode = False,
amp = 2.5,
nr = 1, #500 #300 #100 #50 #30 #10 #5 #3 #1 #1000 #500 #300 #100 #50 #30 #10 #5 #3 #1 #500
sr = 1E8, #16384
plot = True,
maxstep = 100000,
bw = 10E6, #butter worth filter bandwidth
maxy = 0.7, # y scale for plot
curves = 'curves', # module to transmit curves:
)
Mod('curves',
cls = 'frappy_psi.ultrasound.Curves',
description = 't, I, Q and pulse arrays for plot',
)
Mod('roi0',
cls = 'frappy_psi.ultrasound.Roi',
description = 'I/Q of region in the control loop',
time = 300, # this is the center of roi:
size = 5000,
main = f,
)
Mod('roi1',
cls = 'frappy_psi.ultrasound.Roi',
description = 'I/Q of region 1',
time = 100, # this is the center of roi:
size = 300,
main = f,
)
Mod('delay',
cls = 'frappy__psi.dg645.Delay',
description = 'delay line with 2 channels',
uri = 'serial:///dev/ttyS2',
on1 = 1e-9,
on2 = 1E-9,
off1 = 400e-9,
off2 = 600e-9,
)
Mod('pars',
cls = 'frappy_psi.ultrasound.Pars',
description = 'SEA parameters',
)

99
frappy_psi/adq_mr.py
View File

@ -5,15 +5,17 @@ Created on Tue Nov 26 15:42:43 2019
"""
import sys
import atexit
import signal
import time
import numpy as np
import ctypes as ct
import time
from numpy import sqrt, arctan2, sin, cos
import scipy.signal
#from pylab import *
from scipy import signal
#ADQAPI = ct.cdll.LoadLibrary("ADQAPI.dll")
ADQAPI = ct.cdll.LoadLibrary("libadq.so.0")
@ -36,7 +38,7 @@ ADQ_CHANNELS_MASK = 0x3
def butter_lowpass(cutoff, sr, order=5):
nyq = 0.5 * sr
normal_cutoff = cutoff / nyq
b, a = signal.butter(order, normal_cutoff, btype = 'low', analog = False)
b, a = scipy.signal.butter(order, normal_cutoff, btype = 'low', analog = False)
return b, a
@ -82,12 +84,12 @@ class Adq(object):
else :
print('SVN Updated')
print('')
ADQ_CLOCK_INT_INTREF = 0 #internal clock source
ADQ_CLOCK_EXT_REF = 1 #internal clock source, external reference
ADQ_CLOCK_EXT_CLOCK = 2 #External clock source
ADQAPI.ADQ_SetClockSource(self.adq_cu, self.adq_num, ADQ_CLOCK_EXT_REF);
##########################
# Test pattern
#ADQAPI.ADQ_SetTestPatternMode(self.adq_cu, self.adq_num, 4)
@ -95,7 +97,7 @@ class Adq(object):
# Sample skip
#ADQAPI.ADQ_SetSampleSkip(self.adq_cu, self.adq_num, 1)
##########################
# Set trig mode
self.trigger = EXT_TRIG_1
#trigger = LVL_TRIG
@ -105,13 +107,13 @@ class Adq(object):
if (self.trigger == LVL_TRIG):
success = ADQAPI.ADQ_SetLvlTrigLevel(self.adq_cu, self.adq_num, -100)
if (success == 0):
print('ADQ_SetLvlTrigLevel failed.')
print('ADQ_SetLvlTrigLevel failed.')
success = ADQAPI.ADQ_SetTrigLevelResetValue(self.adq_cu, self.adq_num, 1000)
if (success == 0):
print('ADQ_SetTrigLevelResetValue failed.')
print('ADQ_SetTrigLevelResetValue failed.')
success = ADQAPI.ADQ_SetLvlTrigChannel(self.adq_cu, self.adq_num, 1)
if (success == 0):
print('ADQ_SetLvlTrigChannel failed.')
print('ADQ_SetLvlTrigChannel failed.')
success = ADQAPI.ADQ_SetLvlTrigEdge(self.adq_cu, self.adq_num, LVL_RISING)
if (success == 0):
print('ADQ_SetLvlTrigEdge failed.')
@ -124,35 +126,38 @@ class Adq(object):
# print('SetTriggerThresholdVoltage failed.')
print("CHANNEL:"+str(ct.c_int(ADQAPI.ADQ_GetLvlTrigChannel(self.adq_cu, self.adq_num))))
self.setup_target_buffers()
atexit.register(self.deletecu)
signal.signal(signal.SIGTERM, lambda *_: sys.exit(0))
def setup_target_buffers(self):
# Setup target buffers for data
self.target_buffers=(ct.POINTER(ct.c_int16 * self.samples_per_record * self.number_of_records)
* self.max_number_of_channels)()
for bufp in self.target_buffers:
bufp.contents = (ct.c_int16 * self.samples_per_record * self.number_of_records)()
def deletecu(self):
# Only disarm trigger after data is collected
ADQAPI.ADQ_DisarmTrigger(self.adq_cu, self.adq_num)
ADQAPI.ADQ_MultiRecordClose(self.adq_cu, self.adq_num);
# Delete ADQControlunit
ADQAPI.DeleteADQControlUnit(self.adq_cu)
print('ADQ closed')
def start(self):
"""start datat acquisition"""
"""start data acquisition"""
# samples_per_records = samples_per_record/number_of_records
# Change number of pulses to be acquired acording to how many records are taken
# Start acquisition
ADQAPI.ADQ_MultiRecordSetup(self.adq_cu, self.adq_num,
self.number_of_records,
self.samples_per_record)
ADQAPI.ADQ_DisarmTrigger(self.adq_cu, self.adq_num)
ADQAPI.ADQ_ArmTrigger(self.adq_cu, self.adq_num)
def getdata(self):
"""wait for aquisition to be finished and get data"""
"""wait for aquisition to be finished and get data"""
#start = time.time()
while(ADQAPI.ADQ_GetAcquiredAll(self.adq_cu,self.adq_num) == 0):
time.sleep(0.001)
@ -163,7 +168,7 @@ class Adq(object):
self.samples_per_record * self.number_of_records, 2,
0, self.number_of_records, ADQ_CHANNELS_MASK,
0, self.samples_per_record, ADQ_TRANSFER_MODE_NORMAL);
#print(time.time()-mid,mid-start)
#print(time.time()-mid,mid-start)
if not status:
raise ValueError('no succesS from ADQ_GetDATA')
# Now this is an array with all records, but the time is artificial
@ -176,52 +181,7 @@ class Adq(object):
def acquire(self):
self.start()
return self.getdata()
'''
def average(self, data):
#Average over records
return [data[ch].sum(axis=0) / self.number_of_records for ch in range(2)]
def iq(self, channel, f_LO):
newx = np.linspace(0, self.samples_per_record /2, self.samples_per_record)
s0 = channel /((2**16)/2)*0.5*np.exp(1j*2*np.pi*f_LO/(1e3)*newx)
I0 = s0.real
Q0 = s0.imag
return I0, Q0
def fitting(self, data, f_LO, ti, tf):
# As long as data[0] is the pulse
si = 2*ti #Those are for fitting the pulse
sf = 2*tf
phase = np.zeros(self.number_of_records)
amplitude = np.zeros(self.number_of_records)
offset = np.zeros(self.number_of_records)
for i in range(self.number_of_records):
phase[i], amplitude[i] = sineW(data[0][i][si:sf],f_LO*1e-9,ti,tf)
offset[i] = np.average(data[0][i][si:sf])
return phase, amplitude, offset
def waveIQ(self, channel,ti,f_LO):
#channel is not the sample data
t = np.linspace(0, self.samples_per_record /2, self.samples_per_record + 1)[:-1]
si = 2*ti # Again that is where the wave pulse starts
cwi = np.zeros((self.number_of_records,self.samples_per_record))
cwq = np.zeros((self.number_of_records,self.samples_per_record))
iq = np.zeros((self.number_of_records,self.samples_per_record))
q = np.zeros((self.number_of_records,self.samples_per_record))
for i in range(self.number_of_records):
cwi[i] = np.zeros(self.samples_per_record)
cwq[i] = np.zeros(self.samples_per_record)
cwi[i] = amplitude[i]*sin(t[si:]*f_LO*1e-9*2*np.pi+phase[i]*np.pi/180)+bias[i]
cwq[i] = amplitude[i]*sin(t[si:]*f_LO*1e-9*(2*np.pi+(phase[i]+90)*np.pi/180))+bias[i]
iq[i] = channel[i]*cwi[i]
q[i] = channel[i]*cwq[i]
return iq,q
'''
def sinW(self,sig,freq,ti,tf):
# sig: signal array
# freq
@ -257,16 +217,16 @@ class Adq(object):
iorq = np.array([self.mix(data[0][i], data[1][i], freq, ti, tf) for i in range(self.number_of_records)])
# iorq = np.array([self.mix(data[0][:], data[1][:], freq, ti, tf)])
return iorq.sum(axis=0) / self.number_of_records
def filtro(self, iorq, cutoff):
b, a = butter_lowpass(cutoff, self.samp_freq*1e9)
#ifi = np.array(signal.filtfilt(b,a,iorq[0]))
#qf = np.array(signal.filtfilt(b,a,iorq[1]))
iqf = [signal.filtfilt(b,a,iorq[i]) for i in np.arange(len(iorq))]
#ifi = np.array(scipy.signal.filtfilt(b,a,iorq[0]))
#qf = np.array(scipy.signal.filtfilt(b,a,iorq[1]))
iqf = [scipy.signal.filtfilt(b,a,iorq[i]) for i in np.arange(len(iorq))]
return iqf
return iqf
def box(self, iorq, ti, tf):
si = int(self.samp_freq * ti)
sf = int(self.samp_freq * tf)
@ -275,12 +235,15 @@ class Adq(object):
def gates_and_curves(self, data, freq, pulse, roi):
"""return iq values of rois and prepare plottable curves for iq"""
self.ndecimate = int(round(2E9/freq))
times = []
times.append(('aviq', time.time()))
iq = self.averageiq(data,freq*1e-9,*pulse)
times.append(('filtro', time.time()))
iqf = self.filtro(iq,self.bw_cutoff)
m = len(iqf[0]) // self.ndecimate
ll = m * self.ndecimate
iqf = [iqfx[0:ll] for iqfx in iqf]
times.append(('iqdec', time.time()))
iqd = np.average(np.resize(iqf, (2, m, self.ndecimate)), axis=2)
t_axis = np.arange(m) * self.ndecimate / self.samp_freq

28
frappy_psi/ultrasound.py
View File

@ -28,7 +28,8 @@ import numpy as np
import frappy_psi.iqplot as iqplot
from frappy_psi.adq_mr import Adq
from frappy.core import Attached, BoolType, Done, FloatRange, HasIO, \
IntRange, Module, Parameter, Readable, StringIO, StringType
IntRange, Module, Parameter, Readable, StringIO, StringType, \
IDLE, DISABLED, TupleOf, ArrayOf
from frappy.properties import Property
@ -52,7 +53,6 @@ class Roi(Readable):
time = Parameter('start time', FloatRange(unit='nsec'), readonly=False)
size = Parameter('interval (symmetric around time)', FloatRange(unit='nsec'), readonly=False)
enable = Parameter('calculate this roi', BoolType(), readonly=False, default=True)
#status = Parameter(export=False)
pollinterval = Parameter(export=False)
interval = (0,0)
@ -65,6 +65,9 @@ class Roi(Readable):
def calc_interval(self):
self.interval = (self.time - 0.5 * self.size, self.time + 0.5 * self.size)
def read_status(self):
return (IDLE, '') if self.enable else (DISABLED, 'disabled')
def write_time(self, value):
self.time = value
self.calc_interval()
@ -82,7 +85,7 @@ class Pars(Module):
timestamp = Parameter('unix timestamp', StringType(), default='0', readonly=False)
temperature = Parameter('T', FloatRange(unit='K'), default=0, readonly=False)
mf = Parameter('field', FloatRange(unit='T'), default=0, readonly=False)
sr = Parameter('rotaion angle', FloatRange(unit='deg'), default=0, readonly=False)
sr = Parameter('rotation angle', FloatRange(unit='deg'), default=0, readonly=False)
class FreqStringIO(StringIO):
@ -91,16 +94,18 @@ class FreqStringIO(StringIO):
class Frequency(HasIO, Readable):
pars = Attached()
sr = Property('samples per record', datatype=IntRange(), default=16384)
curves = Attached(mandatory=False)
maxy = Property('plot y scale', datatype=FloatRange(), default=0.5)
value = Parameter('frequency@I,q', datatype=FloatRange(unit='Hz'), default=0)
basefreq = Parameter('base frequency', FloatRange(unit='Hz'), readonly=False)
nr = Parameter('number of records', datatype=IntRange(1,10000), default=500)
sr = Parameter('samples per record', datatype=IntRange(1,1E9), default=16384)
freq = Parameter('target frequency', FloatRange(unit='Hz'), readonly=False)
bw = Parameter('bandwidth lowpassfilter', datatype=FloatRange(unit='Hz'),default=10E6)
amp = Parameter('amplitude', FloatRange(unit='dBm'), readonly=False)
control = Parameter('control loop on?', BoolType(), readonly=False, default=True)
rusmode = Parameter('RUS mode on?', BoolType(), readonly=False, default=False)
time = Parameter('pulse start time', FloatRange(unit='nsec'),
readonly=False)
size = Parameter('pulse length (starting from time)', FloatRange(unit='nsec'),
@ -128,6 +133,7 @@ class Frequency(HasIO, Readable):
self.adq = Adq(self.nr, self.sr, self.bw)
self.roilist = []
self.write_nr(self.nr)
self.write_sr(self.sr)
self.skipctrl = 0
self.plotter = iqplot.Plot(self.maxy)
self.calc_interval()
@ -149,6 +155,9 @@ class Frequency(HasIO, Readable):
# self.pollinterval = value * 0.0001
return value
def write_sr(self, value):
return value
def register_roi(self, roi):
self.roilist.append(roi)
@ -174,7 +183,11 @@ class Frequency(HasIO, Readable):
"""main poll loop body"""
if self.lastfreq is None:
self.lastfreq = self.set_freq()
if self.rusmode:
self.sr = int(12e9/self.lastfreq) #picking up 12 period at the ith frequency in the time scale
# self.adq.samples_per_record = self.sr
self.adq.start()
if self.starttime is None:
self.starttime = time.time()
times = []
@ -198,6 +211,8 @@ class Frequency(HasIO, Readable):
gates = self.adq.gates_and_curves(data, freq, self.interval,
[r.interval for r in roilist])
if self.curves: # if attached Curves module is defined, update it
self.curves.value = self.adq.curves
if self.save:
times.append(('save',time.time()))
tdata, idata, qdata, pdata = self.adq.curves
@ -257,3 +272,8 @@ class Frequency(HasIO, Readable):
self.freq = sorted((self.freq - self.maxstep, newfreq, self.freq + self.maxstep))[1]
#print(times)
return Done
class Curves(Readable):
value = Parameter("t, i, q, pulse curves",
TupleOf(*[ArrayOf(FloatRange(), 0, 16283) for _ in range(4)]), default=[[]] * 4)