Qubit-CPW Stark Shift

The system described here is a transmon qubit coupled to a voltage-biased CPW. We are going to model the qubit as a multi-level artificial atom coupled to a single harmonic mode of a superconductor resonator by a Jaynes-Cummings Hamiltonian [1]:

$ \newcommand{\ket}[1]{\left|{#1}\right\rangle} \newcommand{\bra}[1]{\left\langle{#1}\right|} $

$\displaystyle H_{JC} = \hbar\omega_r (a^\dagger a + 1/2) + \hbar\Sigma\omega_j \ket{j}\bra{j} + \hbar\Sigma g_{j,j+1}(a^\dagger\ket{j}\bra{j+1}+a\ket{j+1}\bra{j}) $

where:

$\omega_r$ is the bare resonator frequency
$\omega_j$ is the frequency of the $|j>\rightarrow |j+1>$ transition
$g_{j,j+1}$ is the coupling strength

This Hamiltonian can be diagonalized in the dispersive limit, $\Delta_{j,j+1} \equiv \omega_{j,j+1} - \omega_r \gg g_{j,j+1}$, resuting in:

$\displaystyle \tilde{H}^{2}_{JC} = \hbar\tilde{\omega}_r (a^\dagger a + 1/2) + \dfrac{\hbar\tilde{\omega}_{ge}}{2}\sigma_z + \hbar\chi (a^\dagger a + 1/2) \sigma_z$

where:

$\chi \simeq \chi_{ge}-\dfrac{\chi_{ef}}{2} $
$\chi_{j,j+1} \equiv \dfrac{g^2_{j,j+1}}{\Delta_{j,j+1}}$
$\tilde{\omega}_{r} \simeq \omega_r - \dfrac{\chi_{ef}}{2}$

1.Jaynes, E. T. & Cummings, F. W. Comparison of quantum and semiclassical radiation theories with application to the beam maser. Proceedings of the IEEE 51, 89–109 (1963).



In [1]:

    
#import functions
%pylab inline
from MyUnits import *
from MyFunctions import *
from qutip import *
from MyQubit import *
import mpld3









    



Populating the interactive namespace from numpy and matplotlib



In [2]:

    
# %load_ext autoreload



In [3]:

    
# %aimport MyQubit
# %autoreload 1

Functions

Hamiltonian in the rotating frame of both drives.

$ \displaystyle H_{tot} \approx \hbar \tilde{\Delta}_d\left(a^\dagger a + 1/2\right) + \dfrac{\hbar}{2}\chi(\left(a^\dagger a + 1/2\right)\sigma_z + H_{Kerr}^4+H_{Drive}^2$

$H_{Kerr}^{(4)} \simeq \hbar \zeta \left(a^\dagger a + 1/2\right)^2\sigma_z + \hbar\zeta'\left(a^\dagger a + 1/2\right)^2$

$H_{drive}^{(2)} = \dfrac{\hbar}{2}\Omega_d (a^\dagger + a) + \dfrac{\hbar}{2}\Omega_s (\sigma^+ + \sigma^-) $

Function add to file MyQubit.py

Import data



In [4]:

    
cd /Users/rouxinol/Dropbox/Lab LaHaye/Hugo/2014/12042014









    



/Users/rouxinol/Dropbox/lab lahaye/Hugo/2014/12042014



In [5]:

    
!Head CavityAttn\ MWFreq\ swp\ Flux2pt8\ b

P,F,A,P= ImportXYZZ("/Users/rouxinol/Dropbox/Lab LaHaye/Hugo/2014/12042014/CavityAttn MWFreq swp Flux2pt8 b")

fig, axes = subplots(1,1, figsize=(16,6))
for i in range(6):
    plot(F,A[i,])

A0 = dBmtoV(A[0,])
A1 = dBmtoV(A[1,])
A2 = dBmtoV(A[2,])
A3 = dBmtoV(A[3,])
A4 = dBmtoV(A[4,])
A5 = dBmtoV(A[5,])

Estimated qubit parameters from fabrication

From the measurement of the resistence of one junction the maximum Josephson Energy for the qubit was estimated.



In [6]:

    
Ej = 11.55 # Maximum Josephson Energy

Ec = 0.22 # Capacitive Energy

w_ge_max = sqrt(8*Ec*Ej)-Ec
print('Qubit max frequency: ',w_ge_max, 'GHz')
f = 2e9 # GHz
T = 50e-3 # K
n_th_a = 1/(exp(h*f/(kB*T)-1))

print('<n> thermal: ',n_th_a)









    



Qubit max frequency:  4.28865833702 GHz
<n> thermal:  0.398639675353

System Parameters

From the measurements take in December 2014 the mains system parametes were determined

The resonator bare frequency is $\omega_r$ = 5.049 GHz, with a internal quality factor $Q_I \simeq 200k$, and a loaded quality facot $Q_L \simeq 20k$.

From the resistence of the Josephson resistence a maximum Josephson energy $E_{J,max}/h \approx 11.92$ GHz and a Coulomb charging energy of $E_c/h \approx $ 0.22 GHz (See comments above)

Using these parametes $\chi_{ef}$, $\chi_{ge}$ and $\omega_{ge}$ can be estimated using:

$\chi_{ef} \simeq 2\left(\omega_r - \tilde{\omega}_r\right) $

and

$\chi_{ge} \simeq \chi + \dfrac{\chi_{ef}}{2}$

and finally:

$\omega_{ge} \simeq \tilde{\omega}_{ge} - \chi_{ge}$

$\omega{ef}$ was estimated using two photon spectroscopy of que qubit. We observe $\dfrac{\tilde{\omega}_{gf}}{2} \simeq $ 4.245 GHz one gets $\omega_{gf} \simeq$ 8.49 GHz; using:

$\tilde{\omega}_{gf} \simeq \tilde{\omega}_{ge}+\tilde{\omega}_{ef}$

$\tilde{\omega}_{ef} \simeq \tilde{\omega}_{gf}-\tilde{\omega}_{ge} \approx \omega_{ef}$

Recovering data from measurements realized December 2014



In [7]:

    
# Measured data
w_r = (5.049) *2 * pi # bare resonator frequency

Q_l = 20000 # resonator loaded quality factor

w_r_tilde = (5.065) * 2 * pi # max cavity frequency change

chi = -(0.0054) * 2 * pi # Change in the spectroscopy peaks -> 2*n*Chi = 4.35 - 4.3

w_ge_tilde = (4.365) *2 * pi # Dressed qubit frequency |g> -> |e>


w_gf_tilde = 2 * (4.245) *2 * pi # Dressed qubit frequency |g> -> |f> 


# Derived data

chi_ef = 2*(w_r - w_r_tilde)
print('Chi_ef/2pi = ',chi_ef/2/pi, 'GHz')

chi_ge = chi + chi_ef/2
print('Chi_ge/2pi = ', chi_ge/2/pi, 'GHz\n')

w_ge = w_ge_tilde - chi_ge
print('qubit bare frequency: w_ge/2pi =', w_ge/2/pi, 'GHz')

w_ef = w_gf_tilde - w_ge_tilde + chi_ge - chi_ef  # qubit e -> f transition 
print('qubit |e> -> |f> frequency: w_ef/2pi =', w_ef/2/pi, 'GHz\n')

Delta_ge = w_ge - w_r # g -> e detunning
print('Detunning |g> -> |e>: Delta_ge/2pi =', Delta_ge/2/pi,'GHz')

Delta_ef = w_ef - w_r # g -> e detunning
print('Detunning |e> -> |f>: Delta_ef/2pi =', Delta_ef/2/pi,'GHz\n')

g_ge = sqrt(chi_ge * Delta_ge)
print('Coupling strength |g> -> |e>: g_ge/2pi =', (g_ge/2/pi),'GHz')

g_ef = sqrt(chi_ef * Delta_ef)
print('Coupling strength |e> -> |f>: g_ef/2pi =', (g_ef/2/pi),'GHz\n')









    



Chi_ef/2pi =  -0.03200000000000027 GHz
Chi_ge/2pi =  -0.021400000000000134 GHz

qubit bare frequency: w_ge/2pi = 4.3864 GHz
qubit |e> -> |f> frequency: w_ef/2pi = 4.1356 GHz

Detunning |g> -> |e>: Delta_ge/2pi = -0.6626 GHz
Detunning |e> -> |f>: Delta_ef/2pi = -0.9134000000000004 GHz

Coupling strength |g> -> |e>: g_ge/2pi = 0.119078293572 GHz
Coupling strength |e> -> |f>: g_ef/2pi = 0.170964323764 GHz

Pertubative expansion parameter for dispersive regime



In [8]:

    
# lambda expansion parameter
lambda_ge = g_ge/Delta_ge # expansion parameter lambda << 1
print('Expansion parameter |g> -> |e> << 1: ', lambda_ge, ' << 1')

lambda_ef = g_ef/Delta_ef # expansion parameter lambda << 1
print('Expansion parameter |e> -> |f> << 1: ', lambda_ef, ' << 1 \n')









    



Expansion parameter |g> -> |e> << 1:  -0.179713693891  << 1
Expansion parameter |e> -> |f> << 1:  -0.187173553497  << 1

Kerr-type nonlinearities

The general drive Hamiltonian $\displaystyle H_{drive}$ in the dispersive regime result in a Kerr-type nonlinearities Hamiltonian $H_{Kerr}^{(4)}$ and a Hamiltonian to the cavity and qubit tones $H_{drive}^{(2)}$, where:

$H_{Kerr}^{(4)} \simeq \hbar \zeta \left(a^\dagger a + 1/2\right)^2\sigma_z + \hbar\zeta'\left(a^\dagger a + 1/2\right)^2$

with

$\zeta \approx \chi_{ef}\lambda_{ef}^2 - 2\chi_{ge}\lambda_{ge}^2 + \dfrac{7\chi_{ef}\lambda_{ge}^2}{4} -\dfrac{5\chi_{ge}\lambda_{ef}^2}{4}$

and

$ \zeta' \approx \left(\chi_{ge}-\chi_{ef}\right)\left(\lambda_{ge}^2+\lambda_{ef}^2\right)$

and

$H_{drive}^{(2)} = \dfrac{\hbar}{2}\Omega_d (a^\dagger + a) + \dfrac{\hbar}{2}\Omega_s (\sigma^+ + \sigma^-) $

where

$\Omega_s/\Omega_d$ are the amplitudes of the qubit/resonator.



In [9]:

    
zeta_1 = chi_ef * lambda_ef**2 - 2 * chi_ge * lambda_ge**2 + (7 * chi_ef*lambda_ge**2 - 5 * chi_ge * lambda_ef**2)/4
print('cross-Kerr coeficient: Zeta_1/2pi =', zeta_1/2/pi, 'GHz')

zeta_2 = (chi_ge - chi_ef ) * (lambda_ge**2 + lambda_ef**2)
print('self-Kerr coeficient: Zeta_2/2pi =', zeta_2/2/pi, 'GHz\n')









    



cross-Kerr coeficient: Zeta_1/2pi = -0.000610248735813 GHz
self-Kerr coeficient: Zeta_2/2pi = 0.000713708079544 GHz

Qubit-Cavity characteristic times

Relaxation and Dephasing of the qubit and the resonator dissipation are accounted in the master equation through:

Resonator photon dissipation: $\kappa_n = \dfrac{\omega_r}{Q_l}$
Qubit relaxation: $\Gamma_n$
Qubit dephasing: $\gamma_{\varphi}$

The effects of the finite temperature are incorporated to the master equation through:

Resonator thermal excitations: $\kappa_p$
Qubit thermal excitations: $\Gamma_p$

where $\dfrac{\kappa_p}{\kappa_n}\simeq\exp\left({\dfrac{\hbar\omega_r}{k_BT}}\right)$ and $\dfrac{\Gamma_p}{\Gamma_n}\simeq\exp\left({\dfrac{\hbar\omega_{ge}}{k_BT}}\right)$

Qubit relaxation time was measured to be: $T_1 = 4.3 \mu$s

From Rabi measurements the Rabi decay time was measure to be: $T_{Rabi} = 5.8 \mu$s, and Rabi frequency of: $f_{Rabi} = 2.42$MHz with $\bar{n}_s \simeq 1500$ photons

Using

$T_1 = \dfrac{1}{\Gamma_+ + \Gamma_-}$
$\dfrac{1}{T_2^*}=\dfrac{1}{2T_1}+\dfrac{1}{T_{\varphi }}$



In [10]:

    
Temp = 20e-3 # temperature in Kelvin


T_1 = 4300 # qubit relaxation time
print('T_1=',T_1*1e-9, 's\n')

Delta_w_ge = (0.0002) *2 * pi # Spectroscopy HWHM lowest <n> Stark shift peak
T_2_star = 1/Delta_w_ge 
print('T_2_star =', T_2_star*1e-9, 's\n')

T_rabi = 5800 # Rabi decay time
print('T_rabi=',T_rabi*1e-9, 's')

w_rabi = (0.00248) * 2 * pi # Rabi frequency time
print('w_rabi/2pi =',w_rabi/2/pi, 'GHz\n')

T_dep = 2*T_2_star*T_1/(2*T_1 - T_2_star)
print('T_dephasing =', T_dep*1e-9 , 's\n')


kappa_n = w_r/Q_l # resonator dissipation
print('Resonator dissipation k-/2pi =',kappa_n/2/pi, 'GHz\n')

gamma_rel = 1/T_1
print('Qubit Relaxarion: Gamma_rel =',gamma_rel, 'GHz\n')



gamma_dep = 1/T_dep
print('Qubit dephasing: Gamma_def =',gamma_dep, 'GHz\n')









    



T_1= 4.3e-06 s

T_2_star = 7.957747154594768e-07 s

T_rabi= 5.8e-06 s
w_rabi/2pi = 0.00248 GHz

T_dephasing = 8.769176059676528e-07 s

Resonator dissipation k-/2pi = 0.00025245000000000004 GHz

Qubit Relaxarion: Gamma_rel = 0.00023255813953488373 GHz

Qubit dephasing: Gamma_def = 0.0011403579916684755 GHz



In [11]:

    
# basic parameters

w_d = 5.066 * 2 * pi # resonator drive frequency

w_s = 4.36 * 2 * pi # qubit drive frequency

# derivated parameters

# Delta tilde parameters

Delta_d_tilde = w_r_tilde - w_d
print('Detuning of the dressed resonator/2pi: ', Delta_d_tilde/2/pi)

Delta_s_tilde = w_ge_tilde - w_s
print('Detuning of the dressed qubit/2pi: ', Delta_s_tilde/2/pi, '\n')









    



Detuning of the dressed resonator/2pi:  -0.0009999999999993723
Detuning of the dressed qubit/2pi:  0.005000000000000254

Power vs. cavity occupancy

The cavity mean occupancy is:

$$ \bar{n}=\dfrac{P_{out}Q}{2\pi f hf} $$

Using Blais et al. (PRA 062320 (2004))

$$ \left<n\right> = \dfrac{2P}{\hbar\omega_r\kappa} $$

Also from Blais et al.

$$ \left< V_{out} \right> = \sqrt{R\hbar\omega_r\kappa}\dfrac{\left(a^{\dagger}+a\right)}{2} $$



In [12]:

    
p_out = linspace(-111,-151,9)

f = w_r/2/pi*1e9
Q_c = 20000 # External coupling # from arxiv 1308.2245v1 Suri et alli (??? add by hand)
n_mean = dBmtoW(p_out)*Q_l/(pi*f*h*f) # standard equation 
n_mean_2 = 2*Q_l*dBmtoW(p_out)/(hbar*w_r_tilde*1e9*kappa_n*1e9)/Q_c # from arxiv 1308.2245v1 Suri et alli
n_mean_3 = Q_c*kappa_n*1e9*dBmtoW(p_out)/(2*Q_l*hbar*w_r_tilde*1e9)/((kappa_n*1e9/2)**2+(Delta_d_tilde*1e9-chi*1e9)**2) # from arxiv 1308.2245v1 Suri et alli
n_mean_4 = 2*dBmtoW(p_out)/(hbar*w_r*1e9*kappa_n*1e9)



In [13]:

    
print('Power   =>     <n>    => <n2>   => <n3> => <n4> ')
for i in range(len(p_out)):
    print(p_out[i], ' => ',n_mean[i],' => ', n_mean_2[i] ,' => ', n_mean_3[i],' => ', n_mean_4[i])









    



Power   =>     <n>    => <n2>   => <n3> => <n4> 
-111.0  =>  2993.73524412  =>  2984.27823249  =>  2.45395984213  =>  2993.73524412
-116.0  =>  946.702208294  =>  943.711638633  =>  0.776010238772  =>  946.702208294
-121.0  =>  299.373524412  =>  298.427823249  =>  0.245395984213  =>  299.373524412
-126.0  =>  94.6702208294  =>  94.3711638633  =>  0.0776010238772  =>  94.6702208294
-131.0  =>  29.9373524412  =>  29.8427823249  =>  0.0245395984213  =>  29.9373524412
-136.0  =>  9.46702208294  =>  9.43711638633  =>  0.00776010238772  =>  9.46702208294
-141.0  =>  2.99373524412  =>  2.98427823249  =>  0.00245395984213  =>  2.99373524412
-146.0  =>  0.946702208294  =>  0.943711638633  =>  0.000776010238772  =>  0.946702208294
-151.0  =>  0.299373524412  =>  0.298427823249  =>  0.000245395984213  =>  0.299373524412

Drive amplitude qubit

$\Omega_d \simeq \left(\dfrac{Q_C\kappa_- P_{RF}}{2Q_L\hbar\tilde{\omega}_r}\right)^{1/2}$



In [14]:

    
P_qubit = -140 # Input power qubit (0-20-20-76 - cavity(?))

P_cavity = -(106 + 32) # Input power cavity

Gamma_t = sqrt(kappa_n*dBmtoW(P_cavity)/(hbar*w_r_tilde))/1e9
print(Gamma_t)

Gamma_d = sqrt(Q_c*kappa_n*1e9*dBmtoW(P_cavity)/(2*Q_l*hbar*w_r_tilde*1e9))/1e9  # Drive amplitude qubit
print('Drive amplitude cavity =',Gamma_d/2/pi, 'GHz')

Gamma_s = sqrt(Q_c*kappa_n*1e9*dBmtoW(P_qubit)/(2*Q_l*hbar*w_ge_tilde*1e9))/1e9 # Drive amplitude cavity ????
print('Drive amplitude qubit =',Gamma_s/2/pi,'GHz')









    



0.00273690612269
Drive amplitude cavity = 0.000308010154756 GHz
Drive amplitude qubit = 0.0002635497665 GHz

Simulation



In [15]:

    
Delta = w_r - w_ge

Plot $Tr[\rho \hat{o}]$, $\hat{o} : \sigma_z, n,x,a$ as a function Pump frequency



In [16]:

    
N = 20 # number of states cavity


P_cavity = -(106 + 32) # Input power cavity

Gamma_d = sqrt(Q_c*kappa_n*1e9*dBmtoW(P_cavity)/(2*Q_l*hbar*w_r_tilde*1e9))/1e9  # Drive amplitude qubit
print('Drive amplitude cavity =',Gamma_d/2/pi, 'GHz')

# w_s = 2 * pi * linspace(f_s_i,f_s_f,steps_s)

wr = (5.0596)  * 2* pi
Delta_d_tilde = w_r_tilde - wr 

# qubit freq scan
f_s_i = 4.2 # qubit initial freq
f_s_f = 4.4  # qubit final freq
steps_s = 200
w_s = 2 * pi * linspace(f_s_i,f_s_f,steps_s)

# Power scan
P_i = -100 # qubit tone start
P_f = -140 # qubit tone final
steps_p = 50 # 
P_qubit = linspace(P_i,P_f, steps_p)




amp_sz = []
amp_n = []
amp_x = []
amp_a = []
for Pq in P_qubit:
    
    for wq in w_s:
        
        
        Delta_s_tilde = w_ge_tilde - wq
        Gamma_s = sqrt(Q_c*kappa_n*1e9*dBmtoW(Pq)/(2*Q_l*hbar*w_ge_tilde*1e9))/1e9

        rho,temp, temp2, temp3, temp4 = calc_spectrum_2(N,Delta,w_ge_tilde,Delta_d_tilde,Delta_s_tilde,chi,Gamma_d,Gamma_s, zeta_1,zeta_2,kappa_n,gamma_rel,gamma_dep, n_th_a)
        amp_sz.append((temp.tr()))
        amp_n.append((temp2.tr()))
        amp_x.append((temp3.tr()))
        amp_a.append((temp4.tr()))









    



Drive amplitude cavity = 0.000308010154756 GHz



In [17]:

    
Tr_sz_1 = reshape(amp_sz,(-1,len(w_s+1)))
Tr_n_1 = reshape(amp_n,(-1,len(w_s+1)))
Tr_x_1 = reshape(amp_x,(-1,len(w_s+1)))
Tr_a_1 = reshape(amp_a,(-1,len(w_s+1)))



In [18]:

    
# qubit freq scan
f_s_i = 4.2 # qubit initial freq
f_s_f = 4.4  # qubit final freq
steps_s = 200
w_s = 2 * pi * linspace(f_s_i,f_s_f,steps_s)

# Power scan
P_i = -110 # qubit tone start
P_f = -160 # qubit tone final
steps_p = 50 # 
P_qubit = linspace(P_i,P_f, steps_p)



In [19]:

    
fig, ax = subplots(2,2, figsize=(16,10))

im = ax[0,0].pcolor( w_s/2/pi,P_qubit,abs(Tr_sz_1))#,vmin=0, vmax=1)
fig.colorbar(im, ax=ax[0,0])
ax[0,0].set_xlim(4.27,4.39)
ax[0,0].set_ylim(P_i,P_f)
ax[0,0].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[0,0].set_xlabel(r'Qubit Tone Frequency (GHz)',fontsize=10)
ax[0,0].set_title(r'$Tr[\rho\sigma_z]$',fontsize=20)


im2 = ax[0,1].pcolor(w_s/2/pi,P_qubit,abs(Tr_n_1))#,vmin=0, vmax=1)
fig.colorbar(im2, ax=ax[0,1])
ax[0,1].set_ylim(P_i,P_f)
ax[0,1].set_xlim(4.27,4.39)
ax[0,1].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[0,1].set_xlabel(r'Qubit Tone Frequency (GHz)', fontsize=10)
ax[0,1].set_title(r'$Tr[\rho {a^{\dagger}a}]$',fontsize=20)

im3 = ax[1,0].pcolor(w_s/2/pi,P_qubit,abs(Tr_x_1))#,vmin=0, vmax=1)
fig.colorbar(im2, ax=ax[1,0])
ax[1,0].set_ylim(P_i,P_f)
ax[1,0].set_xlim(4.27,4.39)
ax[1,0].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[1,0].set_xlabel(r'Qubit Tone Frequency (GHz)', fontsize=10)
ax[1,0].set_title(r'$Tr[\rho (a^{\dagger} + a)]$',fontsize=20)

im3 = ax[1,1].pcolor(w_s/2/pi,P_qubit,abs(Tr_a_1))#,vmin=0, vmax=1)
fig.colorbar(im2, ax=ax[1,1])
ax[1,1].set_ylim(P_i,P_f)
ax[1,1].set_xlim(4.27,4.39)
ax[1,1].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[1,1].set_xlabel(r'Qubit Tone Frequency (GHz)', fontsize=10)
ax[1,1].set_title(r'$Tr[\rho a]$',fontsize=20)









    Out[19]:





<matplotlib.text.Text at 0x11028acc0>



In [21]:

    
fig, ax = subplots(1,1, figsize=(20,12))

im = ax.pcolor(w_s/2/pi,P_qubit,log10(abs(Tr_a_1)))#,vmin=0, vmax=1)
cbar = fig.colorbar(im, ax=ax)
cbar.set_label('Arb. Units',size=30)
# access to cbar tick labels:
cbar.ax.tick_params(labelsize=20) 
# im.tick_params(labelsize=10) 
ax.set_xlim(4.27,4.39)
ax.set_ylim(-110,-138)
ax.tick_params(labelsize=20)
ax.tick_params(axis='z', labelsize=20)
# ax.set_ylim(P_i,P_f)
ax.set_ylabel(r'Qubit Tone Power(dBm)',fontsize=30)
ax.set_xlabel(r'Qubit Tone Frequency (GHz)',fontsize=30)

ax.set_title(r'$|Tr[\rho a]|$',fontsize=40)









    Out[21]:





<matplotlib.text.Text at 0x11ab82a90>

$Tr [\rho \hat{o}]$ as a function cavity power (cavity tone)



In [22]:

    
N = 20 # number of states cavity

P_qubit = -140 # Input power qubit (0-20-20-76 - cavity(?))

Gamma_s = sqrt(Q_c*kappa_n*1e9*dBmtoW(P_qubit)/(2*Q_l*hbar*w_ge_tilde*1e9))/1e9 # Drive amplitude cavity ????
print('Drive amplitude qubit =',Gamma_s/2/pi,'GHz')

wr = (5.0596)  * 2* pi
Delta_d_tilde = w_r_tilde - wr 

# qubit freq scan
f_s_i = 4.2 # qubit initial freq
f_s_f = 4.4  # qubit final freq
steps_s = 600
w_s = 2 * pi * linspace(f_s_i,f_s_f,steps_s)

# Power scan
P_i = -160 # qubit tone start
P_f = -120 # qubit tone final
steps_p = 60 # 
P_cavity = linspace(P_i,P_f, steps_p)



rho_matrix = []
amp_sz = []
amp_n = []
amp_x = []
amp_a = []
for Pc in P_cavity:
    
    for wq in w_s:
        
        
        Delta_s_tilde = w_ge_tilde - wq
        Gamma_d = sqrt(Q_c*kappa_n*1e9*dBmtoW(Pc)/(2*Q_l*hbar*w_r_tilde*1e9))/1e9  # Drive amplitude qubit

        rho,temp, temp2, temp3, temp4 = calc_spectrum_2(N,Delta,w_ge_tilde,Delta_d_tilde,Delta_s_tilde,chi,Gamma_d,Gamma_s, zeta_1,zeta_2,kappa_n,gamma_rel,gamma_dep, n_th_a)
#         rho_matrix.append(rho)
        amp_sz.append((temp.tr()))
        amp_n.append((temp2.tr()))
        amp_x.append((temp3.tr()))
        amp_a.append((temp4.tr()))









    



Drive amplitude qubit = 0.0002635497665 GHz



In [23]:

    
Tr_sz = reshape(amp_sz,(-1,len(w_s+1)))
Tr_n = reshape(amp_n,(-1,len(w_s+1)))
Tr_x = reshape(amp_x,(-1,len(w_s+1)))
Tr_a = reshape(amp_a,(-1,len(w_s+1)))



In [24]:

    
fig, ax = subplots(2,2, figsize=(16,10))

im = ax[0,0].pcolor(w_s/2/pi,P_cavity,abs(Tr_sz))#,vmin=0, vmax=1)
fig.colorbar(im, ax=ax[0,0])
ax[0,0].set_xlim(4.27,4.39)
ax[0,0].set_ylim(P_i,P_f)
ax[0,0].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[0,0].set_xlabel(r'Qubit Tone Frequency (GHz)',fontsize=10)
ax[0,0].set_title(r'$Tr[\rho\sigma_z]$',fontsize=20)


im2 = ax[0,1].pcolor(w_s/2/pi,P_cavity,abs(Tr_n))#,vmin=0, vmax=1)
fig.colorbar(im2, ax=ax[0,1])
ax[0,1].set_ylim(P_i,P_f)
ax[0,1].set_xlim(4.27,4.39)
ax[0,1].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[0,1].set_xlabel(r'Qubit Tone Frequency (GHz)', fontsize=10)
ax[0,1].set_title(r'$Tr[\rho {a^{\dagger}a}]$',fontsize=20)

im3 = ax[1,0].pcolor(w_s/2/pi,P_cavity,abs(Tr_x))#,vmin=0, vmax=1)
fig.colorbar(im2, ax=ax[1,0])
ax[1,0].set_ylim(P_i,P_f)
ax[1,0].set_xlim(4.27,4.39)
ax[1,0].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[1,0].set_xlabel(r'Qubit Tone Frequency (GHz)', fontsize=10)
ax[1,0].set_title(r'$Tr[\rho (a^{\dagger} + a)]$',fontsize=20)

im3 = ax[1,1].pcolor(w_s/2/pi,P_cavity,abs(Tr_a))#,vmin=0, vmax=1)
fig.colorbar(im2, ax=ax[1,1])
ax[1,1].set_ylim(P_i,P_f)
ax[1,1].set_xlim(4.27,4.39)
ax[1,1].set_ylabel(r'Qubit Tone Power(dBm)',fontsize=10)
ax[1,1].set_xlabel(r'Qubit Tone Frequency (GHz)', fontsize=10)
ax[1,1].set_title(r'$Tr[\rho a]$',fontsize=20)









    Out[24]:





<matplotlib.text.Text at 0x11c705fd0>



In [ ]:

    
f_s_i = 4.2 # qubit initial freq
f_s_f = 4.4  # qubit final freq
steps_s = 600
w_s = 2 * pi * linspace(f_s_i,f_s_f,steps_s)

# Power scan
P_i = -160 # qubit tone start
P_f = -130 # qubit tone final
steps_p = 60 # 
P_cavity = linspace(P_i,P_f, steps_p)



In [29]:

    
t = 16
print(P_cavity[t])
File = A0

fig, axes = subplots(1,4, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes[0].plot(1e9*w_s/2/pi+4e6,abs(Tr_a[t]),'k-',linewidth=3, label="Abs")
axes[1].plot(1e9*w_s/2/pi+4e6,(Tr_a[t].real),'r-',linewidth=3,label="real")
axes[2].plot(1e9*w_s/2/pi+4e6,(Tr_a[t].imag),'b-',linewidth=3,label="imag")
axes[3].plot(1e9*w_s/2/pi+4e6,angle(Tr_a[t]),'g-',linewidth=3,label="angle")
axes[0].legend(loc=0)
axes[1].legend(loc=0)
axes[2].legend(loc=0)
axes[3].legend(loc=0)
# axes.set_xlabel(r'f (GHz)', fontsize=18)
# axes.set_ylabel(r'$Tr[\rho \sigma_z]$ or $Tr[\rho a^{\dagger}+a]$ ', fontsize=18)
# axes.set_title(r'$Tr[\rho a^{\dagger}+a]$',fontsize=20)









    



-151.86440678






    Out[29]:





<matplotlib.legend.Legend at 0x119141da0>



In [30]:

    
# Test
t = 14
print(P_cavity[t])
File = A1

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi,abs(Tr_a[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \hat{x}]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \hat{x}]$',fontsize=20)

axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label=r'\muV' )


# axes.legend(loc=2)
#left plot
# axes.set_ylim(0.75,1.3) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-152.881355932






    Out[30]:





[<matplotlib.lines.Line2D at 0x11cdd9940>]

Absolute value of $Tr[\rho \hat{x}]$



In [42]:

    
# Test
t = 30
print(P_cavity[t])
File = A0

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_x[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \hat{x}]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \hat{x}]$',fontsize=20)
#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(0.22,.77) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-144.745762712






    Out[42]:





[<matplotlib.lines.Line2D at 0x11fbfbe48>]



In [43]:

    
# Test
t = 33
print(P_cavity[t])
File = A1

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_x[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \hat{x}]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \hat{x}]$',fontsize=20)

axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(0.3,.9) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-143.220338983






    Out[43]:





[<matplotlib.lines.Line2D at 0x11fed7cc0>]



In [44]:

    
# Test
t = 40
print(P_cavity[t])
File = A2

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_x[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \hat{x}]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \hat{x}]$',fontsize=20)

axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(0.75,1.3) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-139.661016949






    Out[44]:





[<matplotlib.lines.Line2D at 0x12037ceb8>]



In [112]:

    
# Test
t = 43
print(P_cavity[t])
File = A3

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_x[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \hat{x}]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \hat{x}]$',fontsize=20)

axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(1.05,1.54) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-130.847457627






    Out[112]:





[<matplotlib.lines.Line2D at 0x119032e80>]



In [46]:

    
# Test
t = 44
print(P_cavity[t])
File = A4

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+3e6,abs(Tr_x[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \hat{x}]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \hat{x}]$',fontsize=20)

#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(1.15,1.62) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-137.627118644






    Out[46]:





[<matplotlib.lines.Line2D at 0x120cbed30>]

Absolute value of $Tr[\rho a]$



In [86]:

    
# Test
t = 36
print(P_cavity[t])
File = A0

fig, axes = subplots(1,1, figsize=(16,10),sharex=True)
fig.tight_layout
axes.set_title(r'$-151$dBm',fontsize=30)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_a[t]),'b-',linewidth=3, label=r'-151 dBm')
axes.set_xlabel(r'probe frequency (GHz)', fontsize=30)
axes.set_ylabel(r'$|Tr[\rho \hat{a}]|$ ', fontsize=30)
# axes.set_title(r'$|Tr[\rho \hat{a}]$',fontsize=20)

#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='-151 dBm')
axes0.set_ylabel(r'$\mu V$ ', fontsize=30)

axes.set_ylim(0.35,0.85) # simulation 
# axes0.set_ylim(0,5) # data
axes0.legend(loc=0)
axes.tick_params(labelsize=25)
axes0.tick_params(labelsize=25)

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-141.694915254






    Out[86]:





[<matplotlib.lines.Line2D at 0x128325470>]



In [98]:

    
# Test
t = 39
print(P_cavity[t])
File = A1

fig, axes = subplots(1,1, figsize=(16,10),sharex=True)
fig.tight_layout
axes.set_title(r'$-147$dBm',fontsize=30)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_a[t]),'b-',linewidth=3, label=r'-147 dBm')
axes.set_xlabel(r'probe frequency (GHz)', fontsize=30)
axes.set_ylabel(r'$|Tr[\rho \hat{a}]|$ ', fontsize=30)
# axes.set_title(r'$|Tr[\rho \hat{a}]$',fontsize=20)

#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='-147 dBm')
axes0.set_ylabel(r'$\mu V$ ', fontsize=30)

axes.set_ylim(0.43,1.0) # simulation 
# axes0.set_ylim(0,5) # data
axes0.legend(loc=0)
axes.tick_params(labelsize=25)
axes0.tick_params(labelsize=25)

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-140.169491525






    Out[98]:





[<matplotlib.lines.Line2D at 0x1269c7630>]



In [49]:

    
# Test
t = 41
print(P_cavity[t])
File = A2

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_a[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \hat{a}]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \hat{a}]$',fontsize=20)

#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(1.15,1.75) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-139.152542373






    Out[49]:





[<matplotlib.lines.Line2D at 0x1210eab70>]



In [133]:

    
# Test
t = 39
print(P_cavity[t])
File = A3

fig, axes = subplots(1,1, figsize=(16,10),sharex=True)
fig.tight_layout
axes.set_title(r'$-134$dBm',fontsize=30)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_a[t]),'b-',linewidth=3, label=r'-134 dBm')
axes.set_xlabel(r'probe frequency (GHz)', fontsize=30)
axes.set_ylabel(r'$|Tr[\rho \hat{a}]|$ ', fontsize=30)
# axes.set_title(r'$|Tr[\rho \hat{a}]$',fontsize=20)

#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='-134 dBm')
axes0.set_ylabel(r'$\mu V$ ', fontsize=30)

axes.set_ylim(1.4,1.9) # simulation 
# axes0.set_ylim(0,5) # data
axes0.legend(loc=0)
axes.tick_params(labelsize=25)
axes0.tick_params(labelsize=25)

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-133.559322034






    Out[133]:





[<matplotlib.lines.Line2D at 0x135b80710>]



In [127]:

    
# Test
t = 41
print(P_cavity[t])
File = A4

fig, axes = subplots(1,1, figsize=(16,10),sharex=True)
fig.tight_layout
axes.set_title(r'$-147$dBm',fontsize=30)
axes.plot(1e9*w_s/2/pi+4e6,abs(Tr_a[t]),'b-',linewidth=3, label=r'-147 dBm')
axes.set_xlabel(r'probe frequency (GHz)', fontsize=30)
axes.set_ylabel(r'$|Tr[\rho \hat{a}]|$ ', fontsize=30)
# axes.set_title(r'$|Tr[\rho \hat{a}]$',fontsize=20)

#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='-147 dBm')
axes0.set_ylabel(r'$\mu V$ ', fontsize=30)

axes.set_ylim(1.68,2.2) # simulation 
# axes0.set_ylim(0,5) # data
axes0.legend(loc=0)
axes.tick_params(labelsize=25)
axes0.tick_params(labelsize=25)

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-132.203389831






    Out[127]:





[<matplotlib.lines.Line2D at 0x1332e2780>]

Angle $Tr[\rho a]$



In [52]:

    
# Test
t = 28
print(P_cavity[t])
File = A0

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

axes.plot(1e9*w_s/2/pi+4e6,-angle(Tr_a[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \sigma_z]$ or $Tr[\rho a^{\dagger}+a]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho a]$',fontsize=20)
#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

axes.set_ylim(1.84,2.05) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-145.762711864






    Out[52]:





[<matplotlib.lines.Line2D at 0x12195dfd0>]



In [53]:

    
# Test
t = 35
print(P_cavity[t])
File = A1

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,-angle(Tr_a[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \sigma_z]$ or $Tr[\rho a^{\dagger}+a]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho a^{\dagger}+a]$',fontsize=20)
#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(1.91,2.02) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-142.203389831






    Out[53]:





[<matplotlib.lines.Line2D at 0x121e7ff60>]



In [54]:

    
# Test
t = 41
print(P_cavity[t])
File = A2

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,-angle(Tr_a[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \sigma_z]$ or $Tr[\rho a^{\dagger}+a]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho a^{\dagger}+a]$',fontsize=20)
#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(1.952,1.987) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-139.152542373






    Out[54]:





[<matplotlib.lines.Line2D at 0x1210ffeb8>]



In [55]:

    
# Test
t = 43
print(P_cavity[t])
File = A3

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,-angle(Tr_a[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \sigma_z]$ or $Tr[\rho a^{\dagger}+a]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho a^{\dagger}+a]$',fontsize=20)
#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(1.954,1.983) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-138.13559322






    Out[55]:





[<matplotlib.lines.Line2D at 0x1220c26d8>]



In [56]:

    
# Test
t = 44
print(P_cavity[t])
File = A4

fig, axes = subplots(1,1, figsize=(16,6),sharex=True)

# fig.subplots_adjust(left=0, wspace=0.25

# axes.plot(1e9*w_s/2/pi+4e6,(Tr_n[t].real),'k-',linewidth=3)
# axes.plot(1e9*w_s/2/pi+4e6,-(Tr_x[t].real),'g-',linewidth=3)
axes.plot(1e9*w_s/2/pi+4e6,-angle(Tr_a[t]),'b-',linewidth=3)
axes.set_xlabel(r'f (GHz)', fontsize=18)
axes.set_ylabel(r'$Tr[\rho \sigma_z]$ or $Tr[\rho a^{\dagger}+a]$ ', fontsize=18)
axes.set_title(r'$Tr[\rho \left( a^{\dagger}+a\right)]$',fontsize=20)
#####
axes0 = axes.twinx()
axes0.plot(F,(File)*1e6,'ro', label='data' )

# axes.legend(loc=2)
#left plot
axes.set_ylim(1.957,1.982) # simulation 
# axes0.set_ylim(0,5) # data

axes0.set_xlim(4.25*1e9,4.40*1e9);
plot(w_s/2/pi,abs(Tr_sz[t]))









    



-137.627118644






    Out[56]:





[<matplotlib.lines.Line2D at 0x1225ac898>]

Bishop Thesis

$$ A = \sqrt{\left(2V_0Tr[\rho a]+ b\zeta\right)^2 + \left(2n_{amps}\right)^2} $$

where:

$V_0$ is a voltage related to the gain of the experimental amplification chain,
$b$ describes the amplitude and phase of the leakage of the drive bypassing the cavity
$\zeta$ the strength of the driving
$n$ is the measurement noise in each of the I and Q channels.



In [57]:

    
V0 = dBmtoV(-106-30)
b = 0
zeta = dBmtoV(-130)
noise = dBmtoV(-150)



In [58]:

    
A = sqrt((2*V0*Tr_a)**2 + (2*noise)**2)



In [59]:

    
plot(1e9*w_s/2/pi+4e6,abs(A[40])*1e6,'b-',linewidth=3)









    Out[59]:





[<matplotlib.lines.Line2D at 0x1228eada0>]



In [59]:



In [ ]: