In [1]:
%pylab inline
In [2]:
import DCRC_controls
import Trigger_Simulation.trigger_simulation as trig_sim
import trig_controls
import picoscope
import ipywidgets as widgets
from test_utilities import *
import scipy.signal
import scipy.optimize
import pymc3
import pickle
In [3]:
# Set up some plotting defaults
matplotlib.rc('figure', figsize=(1.5*3.375, 1.5*3.375*(4/6)), dpi=140)
#matplotlib.rc('text', usetex=True)
matplotlib.rc('font', size=10, family='serif')
matplotlib.rc('legend', scatterpoints=1)
In [4]:
# Connect to the picoscope
scope = picoscope.picoscope()
scope.open()
In [5]:
DCRC = DCRC_controls.DCRC('165.91.180.60', 5001)
DCRC.open() # Connect to the DCRC
DCRC.write_reg(0x20, 0x10) # Make this DCRC the bus master
DCRC.enable_ADC() # Turn on digitization
# Check serial number etc.
IDinfo = DCRC.ID()
for key in IDinfo.keys():
print('% 13s: % 15s'%(key, IDinfo[key]))
print()
print("FPGA Firmware version string:", DCRC.firmware_version())
In [6]:
DCRC.DCRC.send(b"SOCK\r")
time.sleep(1)
print(DCRC.DCRC.recv(32768).decode())
In [7]:
DCRC.DCRC.send(b"drec 0\r")
time.sleep(1)
DCRC.DCRC.recv(32768)
pass
Check for phonon channel inversion¶
Read out all 12 phonon channels, with the test signal on and connected to all twelve, simultaneously. The twelve sine waves should be in phase. If any of the sine waves are upside down (180 degrees out of phase), then that is a problem and you should consult an expert.
In [8]:
%%time
DCRC.write_reg(0x403, 2) # sine wave
DCRC.write_reg(0x405, 0x080)
DCRC.write_reg(0x410, 0) # upper frequency bits
DCRC.write_reg(0x411, 0x4325) # 100 Hz (0x4325)
DCRC.write_reg(0x402, 1) # Turn test signal on
DCRC.write_reg(0x2D4, 1) # Turn on signal-generator triggers
N=4
for i in range(12):
DCRC.write_reg(0x420+i, 1) # Connect the test signal (1)
DCRC.write_reg(0x700+i, 0) # Set driver gain to x1 (0)
DCRC.write_reg(0x710+i, 0x7f) # Set input gain to x5 (0x7f)
DCRC.write_reg(0x730+i, 1) # Set signal source to preamp (1)
DCRC.reset_ADC()
time.sleep(0.25)
trigs = DCRC.read_trig_FIFO()
DCRC.disable_ADC()
sigtrigs = trigs[(trigs['amplitude'] == 4) &
(trigs['trigword'] == 0) &
(trigs['trigbits'] == 255)]
t0 = sigtrigs[0]['timestamp']
dat = array([DCRC.readout(i, 6250*N, t0) for i in range(12)])
plot(dat.T)
tight_layout()
fdat = fft.fft(dat, axis=1)
obs_amp = (abs(fdat[:,N]) * sign(angle(fdat[:,N])) * 2 / (6250*N))
assert all(obs_amp > 0)
In [9]:
def par(*Rs):
"Compute the resistance of several parallel resistances"
return 1/sum(1/array(Rs))
class Phonon_channel:
def __init__(self, R1, R2, R3, R4, R5, R7, R8, R9, R10, R12, R13,
R14, R15, R16, R17, R18, R19, R20, R21, R22, R23):
self.R1 = R1
self.R2 = R2
self.R3 = R3
self.R4 = R4
self.R5 = R5
self.R7 = R7
self.R8 = R8
self.R9 = R9
self.R10 = R10
self.R12 = R12
self.R13 = R13
self.R14 = R14
self.R15 = R15
self.R16 = R16
self.R17 = R17
self.R18 = R18
self.R19 = R19
self.R20 = R20
self.R21 = R21
self.R22 = R22
self.R23 = R23
self.G1 = 1 + self.R3 / par(self.R1, self.R2)
self.G2 = -self.R5 / self.R4
self.G3 = 1 + self.R8 / self.R7
self.G4 = -self.R10 / self.R9
self.G5 = -self.R14 / self.R12
self.GL = -self.R3 / self.R1
self.GT = -self.R14 / self.R13
self.d1S = 1 - self.R17 / (self.R17 + par(self.R16 + self.R15,
self.R18 + par(self.R20, self.R19 + par(self.R21, self.R22))))
self.d15 = 1 - ((self.R15 + self.R16) /
(self.R15 + self.R16 +
par(self.R17, self.R18 + par(self.R20, self.R19 + par(self.R21, self.R22)))))
self.d3i = 1 - (self.R22 /
(self.R22 +
par(self.R21, self.R19 + par(self.R20, self.R18 + par(self.R18, self.R16 + self.R15)))))
self.d23 = 1 - (self.R19 / (self.R19 + par(self.R20, self.R18 + par(self.R17, self.R16 + self.R15))))
self.d12 = 1 - (self.R18 / (self.R18 + par(self.R17, self.R16 + self.R15)))
self.d1i = self.d3i * self.d23 * self.d12
def ADC(self, input_gain, driver_gain, VS_dac, VL_dac, Vi, VT_dac, Voff_dac, mode, scopemode):
"""
mode 1 is preamp
mode 0 is feedback
driver_gain should be 0, 1, 2, or 3, to select x1, x2, x4, or x8, respectively
input_gain should be an integer in the range [0, 255], same as the DCRC setting
VS_dac, VL_dac, VT_dac, and Voff_dac should be the DCRC DAC settings
Vi should be the picoscope signal-generator voltage
scopemode should be "none", "drive", or "read"
"""
w = uint8(127 - input_gain) / 255
G2w = (1-w) + w*self.G2
G123 = self.G1*G2w*self.G3
GL23 = self.GL*G2w*self.G3
G12345 = self.G1*G2w*self.G3*self.G4*self.G5
GL2345 = self.GL*G2w*self.G3*self.G4*self.G5
if scopemode == 'none':
R18_19_20_21 = self.R18 + par(self.R20, self.R19 + self.R21)
elif scopemode == 'drive':
R18_19_20_21 = self.R18 + par(self.R20, self.R19 + par(self.R21, self.R22))
elif scopemode == 'read':
R18_19_20_21 = self.R18 + par(self.R20, self.R19 + par(self.R21, self.R23))
d1S = 1 - self.R17 / (self.R17 + par(self.R16 + self.R15, R18_19_20_21))
d15 = 1 - (self.R15 + self.R16) / (self.R15 + self.R16 + par(self.R17, R18_19_20_21))
R18_17_16_15 = self.R18 + par(self.R17, self.R16 + self.R15)
if scopemode == "drive":
d3i = 1 - self.R22 / (self.R22 + par(self.R21, self.R19 + par(self.R20, R18_17_16_15)))
else:
d3i = 0
d23 = 1 - self.R19 / (self.R19 + par(self.R20, R18_17_16_15))
d12 = 1 - self.R18 / R18_17_16_15
d1i = self.d3i * self.d23 * self.d12
dgain = [1,2,4,8][driver_gain]
V_to_ADC = dgain * 4 * (32768 / 4.096)
# The factor of 4 here is the amplifier after the driver
# programmable-gain-amplifier (PGA) and before the actual
# ADC ("THS4522" on the block diagram). According to the
# block diagram, this is supposed to have a gain of 2. But,
# according to the Rev E schematic, this should have a
# gain of 4 ("U103" on the schematic). However, on DCRC
# serial number 8, I think the actual gain was about 1.8 for
# no apparent reason. I asked Sten to look in to this, but
# never heard back about it. On other boards, the gain will
# probably be approximately 4, as it is supposed to be.
VS = (VS_dac-0x8000) * 125e-6
VL = (VL_dac-0x8000) * 125e-6
Voff = (Voff_dac-0x8000) * 125e-6
VT = (VT_dac / 0x400) * (0xfff/0x1000) * (32 * 200 / 6800)
if mode == 0: # feedback mode
f1 = (self.R15/G12345 + self.R16) / (self.R15 + self.R16)
f2 = (d15*self.R15 + self.R16) / (self.R15 + self.R16)
den = (1 - d15*G12345)
dFS = f1 * (d1S*G12345) / den
dFi = 2 * f1 * (d1i*G12345) / den
dFL = f2 * GL2345 / den
dFT = f2 * self.GT / den
return V_to_ADC * array([VS*dFS + VL*dFL - Voff, VT*dFT, Vi*dFi])
elif mode == 1: # preamp mode
return V_to_ADC * array([VS*G123*d1S + VL*GL23 - Voff, VT*d15*self.GT*G123, 2*Vi*G123*d1i])
In [10]:
pc = Phonon_channel(R1=1.5e3,
R2=10,
R3=390,
R4=1e3,
R5=1e3,
R7=1e3,
R8=4.99e3,
R9=1.5e3,
R10=30e3,
R12=2.49e3,
R13=4.99e3,
R14=4.99e3,
R15=24.9,
R16=10e3,
R17=20e3,
R18=92,
R19=909,
R20=10,
R21=52.3,
R22=600,
R23=1e6)
In [11]:
def test(DCRC, input_gain, driver_gain,
VS_dac, VL_dac, Vi, VT_dac, Voff_dac,
mode, scope=None, scopemode='none',
sleeptime=10):
"""
mode 1 is preamp
mode 0 is feedback
driver_gain should be 0, 1, 2, or 3, to select x1, x2, x4, or x8, respectively
input_gain should be an integer in the range [0, 255], same as the DCRC setting
VS_dac, VL_dac, VT_dac, and Voff_dac should be the DCRC DAC settings
Vi should be the picoscope signal-generator voltage
scopemode should be "none", "drive", or "read"
"""
DCRC.write_reg(0x403, 2) # sine wave
DCRC.write_reg(0x405, VT_dac)
DCRC.write_reg(0x410, 0) # upper frequency bits
DCRC.write_reg(0x411, 0x4325) # 100 Hz (0x4325)
DCRC.write_reg(0x402, 1) # Turn test signal on
DCRC.write_reg(0x2D4, 1) # Turn on signal-generator triggers
N=16
need_sleep = ((DCRC.read_reg(0x760) != Voff_dac) or
(DCRC.read_reg(0x770) != VL_dac) or
(DCRC.read_reg(0x780) != VS_dac))
for i in range(12):
DCRC.write_reg(0x420+i, 1) # Connect the test signal (1)
DCRC.write_reg(0x700+i, driver_gain) # Set driver gain to x1 (0)
DCRC.write_reg(0x710+i, input_gain) # Set input gain to x5 (0x7f)
DCRC.write_reg(0x730+i, mode) # Set signal source to preamp (1)
DCRC.write_reg(0x760+i, Voff_dac)
DCRC.write_reg(0x770+i, VL_dac)
DCRC.write_reg(0x780+i, VS_dac)
if need_sleep:
print("Sleeping: 0x%04x 0x%04x 0x%04x"%(VS_dac, VL_dac, Voff_dac), flush=True)
time.sleep(sleeptime)
if (scope is not None) and (scopemode == "drive"):
pass # Do some stuff to set up the scope to send the input signal
DCRC.reset_ADC()
time.sleep(0.25)
trigs = DCRC.read_trig_FIFO()
DCRC.disable_ADC()
sigtrigs = trigs[(trigs['amplitude'] == 4) &
(trigs['trigword'] == 0) &
(trigs['trigbits'] == 255)]
t0 = sigtrigs[0]['timestamp']
dat = array([DCRC.readout(i, 6250*N, t0) for i in range(12)])
fdat = fft.fft(dat, axis=1) * 2 / float(6250*N)
AmpDC = real(fdat[:,0]) / 2 # DC component gets a factor of 2 relative to other components
AmpT = abs(fdat[:,N]) * sign(angle(fdat[:,N]))
AmpScope = abs(fdat[:,5*N]) * 2*(1 - (angle(fdat[:,10*N])/pi - 2*angle(fdat[:,5*N])/pi) % 2)
return (array([AmpDC, AmpT, AmpScope]),
{'input_gain':input_gain, 'driver_gain':driver_gain,
'VS_dac':VS_dac, 'VL_dac':VL_dac, 'Vi':Vi, 'VT_dac':VT_dac, 'Voff_dac':Voff_dac,
'mode':mode, 'scopemode':scopemode})
In [12]:
feedback_gains = [0x80, 0xc0, 0xe0, 0xf0, 0xf8, 0xfc, 0xfe, 0xff, 0x00]
preamp_gains = [0x80, 0xc0, 0xf0, 0xff, 0x00, 0x0f, 0x3f, 0x7f]
VS_dacs = [0x7000, 0x8000, 0x9000]
VL_dacs = [0x7000, 0x8000, 0x9000]
Voff_dacs = [0x7000, 0x8000, 0x9000]
VT_dacs = [0x080, 0x100, 0x200]
driver_gains = [0,1,2,3]
#driver_gains = [0]
In [13]:
%%time
tests = []
print("SQUID bias tests", flush=True)
for VS_dac in VS_dacs:
for driver_gain in driver_gains:
for input_gain in feedback_gains:
tests.append(test(DCRC, input_gain, driver_gain, VS_dac, 0x8000, 0, 0x000, 0x8000, 0, None, 'none'))
for input_gain in preamp_gains:
tests.append(test(DCRC, input_gain, driver_gain, VS_dac, 0x8000, 0, 0x000, 0x8000, 1, None, 'none'))
print("SQUID lockpoint tests", flush=True)
for VL_dac in VL_dacs:
for driver_gain in driver_gains:
for input_gain in feedback_gains:
tests.append(test(DCRC, input_gain, driver_gain, 0x8000, VL_dac, 0, 0x000, 0x8000, 0, None, 'none'))
for input_gain in preamp_gains:
tests.append(test(DCRC, input_gain, driver_gain, 0x8000, VL_dac, 0, 0x000, 0x8000, 1, None, 'none'))
print("Test signal tests", flush=True)
for VT_dac in VT_dacs:
for driver_gain in driver_gains:
for input_gain in feedback_gains:
tests.append(test(DCRC, input_gain, driver_gain, 0x8000, 0x8000, 0, VT_dac, 0x8000, 0, None, 'none'))
for input_gain in preamp_gains:
tests.append(test(DCRC, input_gain, driver_gain, 0x8000, 0x8000, 0, VT_dac, 0x8000, 1, None, 'none'))
print("ADC offset tests", flush=True)
for Voff_dac in Voff_dacs:
for driver_gain in driver_gains:
tests.append(test(DCRC, 0x80, driver_gain, 0x8000, 0x8000, 0, 0x000, Voff_dac, 0, None, 'none'))
tests.append(test(DCRC, 0x00, driver_gain, 0x8000, 0x8000, 0, 0x000, Voff_dac, 1, None, 'none'))
all_tests = tests
In [15]:
'''
with open("all_tests.pkl", "rb") as f:
all_tests = pickle.load(f)
'''
In [14]:
%%time
tests = [test for test in all_tests if test[1]['driver_gain'] in (0,1)]
with pymc3.Model() as model:
Rs = array([1.5e3, 10, 390, 1e3, 1e3, 1e3, 4.99e3, 1.5e3, 30e3, 2.49e3, 4.99e3, 4.99e3, 24.9,
10e3, 20e3, 92, 909, 10, 52.3])
R = pymc3.Normal("R", mu=Rs, sd=0.01*Rs, shape=19, testval=Rs)
off = pymc3.Normal("off", mu=0, sd=3000, shape=len(feedback_gains)+len(preamp_gains))
ifFeedback = pymc3.math.constant(array([test_parms['mode']==0 for result, test_parms in tests]))
dct = dict()
j = 0
for gain in feedback_gains:
dct[(0, gain)] = j
j+= 1
for gain in preamp_gains:
dct[(1, gain)] = j
j+= 1
offinds = pymc3.math.constant([dct[(test_parms['mode'], test_parms['input_gain'])]
for result, test_parms in tests])
offs = off[offinds]
G1 = pymc3.Deterministic("G1", 1 + R[2] * (1/R[0] + 1/R[1]))
G2 = pymc3.Deterministic("G2", - R[4] / R[3])
G3 = pymc3.Deterministic("G3", 1 + R[6] / R[5])
G4 = pymc3.Deterministic("G4", - R[8] / R[7])
G5 = pymc3.Deterministic("G5", - R[11] / R[9])
GL = pymc3.Deterministic("GL", - R[2] / R[0])
GT = pymc3.Deterministic("GT", R[11] / R[10])
w = pymc3.math.constant(array([uint8(127 - test_parms['input_gain'])/255 for result, test_parms in tests]))
G2w = (1-w) + w*G2
R18_19_20_21 = R[15] + 1/(1/R[17] + 1/(R[16] + R[18]))
d1S = pymc3.Deterministic("d1S", 1 - R[14] / (R[14] + 1/(1/(R[13] + R[12]) + 1/R18_19_20_21)))
d15 = pymc3.Deterministic("d15", 1 - (R[12] + R[13]) / (R[12] + R[13] + 1/(1/R[14] + 1/R18_19_20_21)))
f1 = (R[12]/(G1*G2w*G3*G4*G5) + R[13]) / (R[12] + R[13])
f2 = (d15*R[12] + R[13]) / (R[12] + R[13])
#V_to_ADC = pymc3.math.constant(array([[1,2,4,8][test_parms['driver_gain']] * 4 * (32768 / 4.096)
# for result, test_parms in tests]))
GADC = pymc3.Normal("GADC", mu=1.823, sd=0.05, testval=1.823)
# Related to the earlier comment about the gain of the amplifier
# between the driver PGA and the ADCs. On DCRC serial number 8, this
# gain appeared to be about 1.8. On other Rev E DCRCs it should be very
# close to 4. So you'll want to change the "1.823" in the line above this to 4.
DGain = pymc3.math.constant(array([[1,2,4,8][test_parms['driver_gain']] for result, test_parms in tests]))
V_to_ADC = pymc3.Deterministic("V_to_ADC", DGain * GADC * (32768 / 4.096))
VS0 = pymc3.Normal("VS0", mu=0x8000, sd=0x100, testval=0x8000)
VL0 = pymc3.Normal("VL0", mu=0x8000, sd=0x100, testval=0x8000)
Voff0 = pymc3.Normal("Voff0", mu=0x8000, sd=0x100, testval=0x8000)
VSm = pymc3.Normal("VSm", mu=125e-6, sd=10e-6, testval=125e-6)
VLm = pymc3.Normal("VLm", mu=125e-6, sd=10e-6, testval=125e-6)
Voffm = pymc3.Normal("Voffm", mu=125e-6, sd=10e-6, testval=125e-6)
# This is R422 on the Rev E "Analog section" schematic
RL = pymc3.Normal("RL", mu=200, sd=2, testval=200)
# This is R469 and R470 on the Rev E "Analog section" schematic
Ra = pymc3.Normal("Ra", mu=6800, sd=68, testval=6800)
VTm = pymc3.Deterministic("VTm", (1/0x400) * (0xfff/0x1000) * 32 * RL / Ra)
VS_dac = pymc3.math.constant(array([test_parms['VS_dac'] for result, test_parms in tests]))
VS = (VS_dac - VS0) * VSm
VL_dac = pymc3.math.constant(array([test_parms['VL_dac'] for result, test_parms in tests]))
VL = (VL_dac - VL0) * VLm
Voff_dac = pymc3.math.constant(array([test_parms['Voff_dac'] for result, test_parms in tests]))
Voff = (Voff_dac - Voff0) * Voffm
VT_dac = pymc3.math.constant(array([test_parms['VT_dac'] for result, test_parms in tests]))
VT = VT_dac * VTm
# DC (Direct Current, or 0 Hz) component (feedback mode)
VF0 = V_to_ADC * (((f1*VS*d1S*G1 + f2*VL*GL)*G2w*G3*G4*G5) / (1 - d15*G1*G2w*G3*G4*G5) - Voff) - offs
# 100 Hz component (feedback mode)
VF1 = V_to_ADC * (f2*VT*GT) / (1 - d15*G1*G2w*G3*G4*G5)
VF = pymc3.math.stack(VF0, VF1)
# DC (Direct Current, or 0 Hz) component (preamp mode)
Vp0 = V_to_ADC * ((VS*d1S*G1 + VL*GL)*G2w*G3 - Voff) - offs
# 100 Hz component (preamp mode)
Vp1 = V_to_ADC * (VT*d15*GT*G1*G2w*G3)
Vp = pymc3.math.stack(Vp0, Vp1)
Vout = pymc3.Deterministic("Vout", pymc3.math.where(ifFeedback, VF, Vp))
fit_channel = 3 # Fit one channel at a time
result = pymc3.math.constant(array([result[:2, fit_channel] for result, test_parms in tests]).T)
obs = pymc3.Normal("obs", mu=Vout, sd=50, observed=result)
trace = pymc3.sample(2000, tune=2000, init='adapt_diag')
In [10]:
with open("trace.pkl", "rb") as f:
trace = pickle.load(f)
In [13]:
with open("alltrace.pkl", "rb") as f:
alltrace = pickle.load(f)
In [15]:
i=0
hist(trace['R'][:,i], bins='auto', histtype='step')
axvline(Rs[i], c='k', ls='--')
tight_layout()
In [16]:
hist(-trace['GL'], bins='auto', histtype='step')
axvline(-pc.GL)
Out[16]:
In [17]:
hist(trace['G1'], bins='auto', histtype='step')
axvline(pc.G1)
Out[17]:
In [18]:
hist(-trace['G2'], bins='auto', histtype='step')
axvline(-pc.G2)
Out[18]:
In [19]:
std(trace['G2']) / abs(mean(trace['G2']))
Out[19]:
In [20]:
hist(trace['G3'], bins='auto', histtype='step')
axvline(pc.G3)
Out[20]:
In [21]:
100 * std(trace['G3']) / abs(mean(trace['G3']))
Out[21]:
In [22]:
hist(-trace['G4'], bins='auto', histtype='step')
axvline(-pc.G4)
100 * std(trace['G4']) / abs(mean(trace['G4']))
Out[22]:
In [23]:
hist(-trace['G5'], bins='auto', histtype='step')
axvline(-pc.G5)
100 * std(trace['G5']) / abs(mean(trace['G5']))
Out[23]:
In [24]:
hist(trace['GT'], bins='auto', histtype='step')
axvline(-pc.GT)
100 * std(trace['GT']) / abs(mean(trace['GT']))
Out[24]:
In [25]:
hist(trace['G1']*trace['G3']*trace['G4']*trace['G5'], bins='auto', histtype='step')
axvline(pc.G1 * pc.G3 * pc.G4 * pc.G5, c='k', ls='--')
tight_layout()
In [26]:
hist(-4096 * trace['GL'] / trace['G1'], bins='auto', histtype='step')
axvline(4096 * 0.26 / 40.26, c='k', ls='--')
tight_layout()
In [27]:
hist(trace['d1S']*100, bins='auto', histtype='step')
axvline(pc.d1S*100)
Out[27]:
In [28]:
hist(trace['d15']*100, bins='auto', histtype='step')
axvline(pc.d15*100)
Out[28]:
In [29]:
hist(trace['VSm']*1e6, bins='auto', histtype='step')
axvline(125)
Out[29]:
In [30]:
hist(trace['VLm']*1e6, bins='auto', histtype='step')
axvline(125)
Out[30]:
In [31]:
hist(trace['Voffm']*1e6, bins='auto', histtype='step')
axvline(125)
Out[31]:
In [32]:
hist(trace['VTm']*1e6, bins='auto', histtype='step')
axvline(941, c='k', ls='--')
Out[32]:
In [33]:
hist(trace['GADC'], bins='auto', histtype='step')
tight_layout()
In [34]:
scatter(trace['VS0'], trace['VL0'], s=0.2, lw=0, c='k')
xs = linspace(32000, 33750)
plot(xs, 0x8000-260 - (xs-0x8000)*pc.d1S*pc.G1/pc.GL)
plot(xs, 0x8000-260 - (xs-0x8000)*mean(trace['d1S']*trace['G1']/trace['GL']))
axhline(0x8000)
axvline(0x8000)
Out[34]:
In [35]:
trace.varnames
Out[35]:
In [36]:
alltrace = concatenate([array([trace[name] for name in ['VSm', 'VS0', 'VLm', 'VL0', 'Voffm', 'Voff0',
'G1', 'G2', 'G3', 'G4', 'G5', 'GL', 'GT',
'd1S', 'd15']]).T,
trace['R'],
trace['off']], axis=1)
names = ['VSm', 'VS0', 'VLm', 'VL0', 'Voffm', 'Voff0',
'G1', 'G2', 'G3', 'G4', 'G5', 'GL', 'GT',
'd1S', 'd15',
'R1 [0]', 'R2 [1]', 'R3 [2]', 'R4 [3]', 'R5 [4]',
'R7 [5]', 'R8 [6]', 'R9 [7]', 'R10 [8]', 'R12 [9]',
'R13 [10]', 'R14 [11]', 'R15 [12]', 'R16 [13]', 'R17 [14]',
'R18 [15]', 'R19 [16]', 'R20 [17]', 'R21 [18]',
'offF0', 'offF1', 'offF2', 'offF3', 'offF4', 'offF5', 'offF6', 'offF7', 'offF8',
'offP0', 'offP1', 'offP2', 'offP3', 'offP4', 'offP5', 'offP6', 'offP7']
In [37]:
len(names)
Out[37]:
In [38]:
alltrace.shape
Out[38]:
In [39]:
corrs = (abs(corrcoef(alltrace.T)) - diag(ones(51)))
pcolor(corrs, vmax=1, vmin=0)
colorbar()
Out[39]:
In [40]:
hist(trace['G1']*trace['G3'], bins='auto', histtype='step')
axvline(pc.G1*pc.G3)
tight_layout()
In [41]:
hist(trace['G4']*trace['G5'], bins='auto', histtype='step')
axvline(pc.G4*pc.G5)
tight_layout()
In [42]:
mostcorr = argsort(corrs.flatten())[::-1]
In [43]:
k=2
i1, i2 = unravel_index(mostcorr[2*k], (51,51))
print(names[i1])
print(names[i2])
In [44]:
scatter(abs(trace['G5']), abs(trace['G4']), s=1, lw=0, c='k')
xlabel('G5')
ylabel('G4')
G5s = linspace(1.82, 1.98)
plot(G5s, mean(trace['G5']*trace['G4'])/G5s)
tight_layout()
In [45]:
scatter(trace['G3'], Rs[6]/1000 * trace['R'][:,6], s=1, lw=0, c='k')
xlabel('G3')
ylabel('R8')
# G3 with R8 makes good sense
tight_layout()
In [46]:
scatter(trace['off'][:,16], trace['off'][:,8], s=1, lw=0, c='k')
xlabel('preamp offset 7')
ylabel('feedback offset 8')
# Not sure why these would be correlated
tight_layout()
In [47]:
scatter(Rs[13]/1000 * trace['R'][:,13], trace['d15']*1000, s=1, lw=0, c='k')
xlabel('R16')
ylabel('d15')
# d15 with R16 makes sense
tight_layout()
In [48]:
scatter(trace['GL'], Rs[2]/1000 * trace['R'][:,2], s=1, lw=0, c='k')
xlabel('GL')
ylabel('R3')
# GL with R3 makes sense
tight_layout()
In [49]:
scatter(Rs[14]/1000 * trace['R'][:,14], trace['VSm']*1e6, s=1, lw=0, c='k')
xlabel('R17')
ylabel('VSm')
# VS slope with R17 makes sense
tight_layout()
In [50]:
scatter(Rs[0]/1000 * trace['R'][:,0], trace['GL'], s=1, lw=0, c='k')
xlabel('R1')
ylabel('GL')
# GL with R1 makes really good sense
tight_layout()
In [51]:
scatter(trace['VLm']*1e6, trace['GL'], s=1, lw=0, c='k')
xlabel('VL slope')
ylabel('GL')
# GL with VL slope makes really good sense
tight_layout()
In [52]:
scatter(trace['G5'], Rs[9]*trace['R'][:,9]/1000, s=1, lw=0, c='k')
xlabel('G5')
ylabel('R12')
# G5 with R12 makes sense
tight_layout()
In [53]:
scatter(trace['d1S']*1000, Rs[14]*trace['R'][:,14]/1000, s=1, lw=0, c='k')
xlabel('d1S')
ylabel('R17')
# d1S with R17 makes good sense
tight_layout()
In [54]:
scatter(trace['VS0'], trace['VL0'], s=0.2, lw=0, c='k', zorder=10)
xlabel('Zero of $V_S$ [DAC setting]')
ylabel('Zero of $V_L$ [DAC setting]')
axhline(0x8000, c='k', ls='-', alpha=0.4)
axvline(0x8000, c='k', ls='-', alpha=0.4)
grid(True)
def to_hex(x, pos):
return '0x%x' % int(x)
fmt = matplotlib.ticker.FuncFormatter(to_hex)
gca().get_xaxis().set_major_locator(matplotlib.ticker.MultipleLocator(0x200))
gca().get_xaxis().set_major_formatter(fmt)
gca().get_yaxis().set_major_locator(matplotlib.ticker.MultipleLocator(0x200))
gca().get_yaxis().set_major_formatter(fmt)
gca().set_aspect('equal', 'datalim')
tight_layout()
In [55]:
scatter(Rs[3]/1000 * trace['R'][:,3], Rs[4]/1000 * trace['R'][:,4], s=1, lw=0, c='k', zorder=2)
xlabel('R4')
ylabel('R5')
axhline(1000, c='k', alpha=0.5)
axvline(1000, c='k', alpha=0.5)
grid(True)
# R4 with R5 makes sense iff G2 is tightly constrained, which it is
tight_layout()
In [56]:
pymc3.plots.pairplot(trace, varnames=['G1', 'G2', 'G3', 'G4', 'G5', 'GL', 'GT', 'd1S', 'd15'])
Out[56]:
In [57]:
hist(trace['VSm']*1e6, bins='auto', histtype='step')
axvline(125)
tight_layout()
In [58]:
hist(trace['VLm']*1e6, bins='auto', histtype='step')
axvline(125)
tight_layout()
In [59]:
hist(trace['Voffm']*1e6, bins='auto', histtype='step')
axvline(125)
tight_layout()
In [60]:
scatter(trace['VS0'], trace['VL0'], s=0.2, lw=0, c='k', zorder=10)
xlabel('Zero of $V_S$ [DAC setting]')
ylabel('Zero of $V_L$ [DAC setting]')
axhline(0x8000, c='k', ls='-', alpha=0.4)
axvline(0x8000, c='k', ls='-', alpha=0.4)
grid(True)
def to_hex(x, pos):
return '0x%x' % int(x)
fmt = matplotlib.ticker.FuncFormatter(to_hex)
gca().get_xaxis().set_major_locator(matplotlib.ticker.MultipleLocator(0x200))
gca().get_xaxis().set_major_formatter(fmt)
gca().get_yaxis().set_major_locator(matplotlib.ticker.MultipleLocator(0x200))
gca().get_yaxis().set_major_formatter(fmt)
gca().set_aspect('equal', 'datalim')
tight_layout()
In [61]:
scatter(Rs[3]*trace['R'][:,3]/1000, Rs[4]*trace['R'][:,4]/1000, s=0.2, lw=0, c='k', zorder=10)
xlabel('$R_4$ [$\Omega$]')
ylabel('$R_5$ [$\Omega$]')
axhline(pc.R5, c='k', alpha=0.4)
axvline(pc.R4, c='k', alpha=0.4)
grid(True)
gca().set_aspect('equal', 'datalim')
tight_layout()
In [62]:
errorbar(x=arange(result.data.shape[1]),
y=mean(trace['Vout'][:,0,:], axis=0),
yerr=sqrt(var(trace['Vout'][:,0,:], axis=0) + 50**2),
ls='none', marker='.', markersize=1, capsize=1.3)
#plot(Vout.eval({R:Rs, Voff0:0x8000, VS0:0x8000, VL0:0x8000, off:zeros(17)})[0,:], c='k')
plot(result.data[0,:], c='r')
xlim(17*0, 17*6)
ylim(-20000, 20000)
Out[62]:
In [63]:
errorbar(x=arange(result.data.shape[1]),
y=mean(trace['Vout'][:,1,:], axis=0),
yerr=std(trace['Vout'][:,1,:], axis=0),
ls='none', marker='.', markersize=1, capsize=1.3)
#plot(Vout.eval({R:Rs, Voff0:0x8000, VS0:0x8000, VL0:0x8000, off:zeros(17)})[1,:], c='k')
plot(result.data[1,:], c='r')
xlim(200, 320)
Out[63]:
In [64]:
k=10 + 17*4
hist(trace['Vout'][:,0,k], bins='auto', histtype='step')
axvline(result.data[0,k])
tight_layout()
In [65]:
errorbar(x=arange(17), y=trace['off'].mean(axis=0), yerr=trace['off'].std(axis=0), ls='none', marker='.', capsize=2)
Out[65]:
In [ ]: