# Useful additional packages 
import matplotlib.pyplot as plt
%matplotlib inline
import numpy as np
import math
import random

from qiskit import QuantumCircuit, ClassicalRegister, QuantumRegister, execute
from qiskit.tools.visualization import circuit_drawer
from qiskit import BasicAer, Aer

qc = QuantumCircuit(65, 32)

# Z-type Stabilizers
Zstab = [
    [47, [0,3,5]],
    [48, [1,3,4,6]],
    [49, [2,4,7]],
    [50, [5,8,10]],
    [51, [6,8,9,11]],
    [52, [7,9,12]],
    [53, [10,13,15]],
    [54, [11,13,14,16]],
    [55, [12,14,17]],
    [56, [15,18,20]],
    [57, [16,18,19,21]],
    [58, [17,19,22]],
    [59, [20,23,25]],
    [60, [21,23,24,26]],
    [61, [22,24,27]],
    [62, [25,28,30]],
    [63, [26,28,29,31]],
    [64, [27,29,32]]
]

# X-type Stabilizers
Xstab = [
    [33, [0,1,3]],
    [34, [1,2,4]],
    [35, [3,5,6,8]],
    [36, [4,6,7,9]],
    [37, [8,10,11,13]],
    [38, [9,11,12,14]],
    [39, [13,15,16,18]],
    [40, [14,16,17,19]],
    [41, [18,20,21,23]],
    [42, [19,21,22,24]],
    [43, [23,25,26,28]],
    [44, [24,26,27,29]],
    [45, [28,30,31]],
    [46, [29,31,32]]
]

# Make stabilizer state
Xinitialize = [
    [0, 1, 3],
    [2, 1, 4],
    [5, 3, 6, 8],
    [7, 4, 6, 9],
    [10, 8, 11, 13],
    [12, 9, 11, 14],
    [15, 13, 16, 18],
    [17, 14, 16, 19],
    [20, 18, 21, 23],
    [22, 19, 21, 24],
    [25, 23, 26, 28],
    [27, 24, 26, 29],
    [30, 28, 31],
    [32, 29, 31],
]

for x in Xinitialize:
    q0 = x[0]
    qc.h(q0)
    for q1 in x[1:]:
        qc.cx(q0, q1)

# Place stabilizers
for z in Zstab:
    qc.h(z[0])
    for q in z[1]:
        qc.cz(z[0], q)
    qc.h(z[0])

for x in Xstab:
    qc.h(x[0])
    for q in x[1]:
        qc.cx(x[0], q)
    qc.h(x[0])

# Measurement stabilizer errors
for m, q in enumerate(Zstab + Xstab):
    qc.measure(q[0], m)


backend = Aer.get_backend('qasm_simulator')
job = execute(qc, backend, shots=1000)
job.result().get_counts(qc)

qc = QuantumCircuit(65, 32+1)  # classical bit-32 for logical bit test

# Z-type Stabilizers
Zstab = [
    [47, [0,3,5]],
    [48, [1,3,4,6]],
    [49, [2,4,7]],
    [50, [5,8,10]],
    # [51, [6,8,9,11]],  #  defect Z-cut hole 1
    [52, [7,9,12]],
    [53, [10,13,15]],
    [54, [11,13,14,16]],
    [55, [12,14,17]],
    [56, [15,18,20]],
    [57, [16,18,19,21]],
    [58, [17,19,22]],
    [59, [20,23,25]],
    #  [60, [21,23,24,26]],  # defect Z-cut hole 2
    [61, [22,24,27]],
    [62, [25,28,30]],
    [63, [26,28,29,31]],
    [64, [27,29,32]]
]

# X-type Stabilizers
Xstab = [
    [33, [0,1,3]],
    [34, [1,2,4]],
    [35, [3,5,6,8]],
    [36, [4,6,7,9]],
    [37, [8,10,11,13]],
    [38, [9,11,12,14]],
    [39, [13,15,16,18]],
    [40, [14,16,17,19]],
    [41, [18,20,21,23]],
    [42, [19,21,22,24]],
    [43, [23,25,26,28]],
    [44, [24,26,27,29]],
    [45, [28,30,31]],
    [46, [29,31,32]]
]

# Make stabilizer state
Xinitialize = [
    [0, 1, 3],
    [2, 1, 4],
    [5, 3, 6, 8],
    [7, 4, 6, 9],
    [10, 8, 11, 13],
    [12, 9, 11, 14],
    [15, 13, 16, 18],
    [17, 14, 16, 19],
    [20, 18, 21, 23],
    [22, 19, 21, 24],
    [25, 23, 26, 28],
    [27, 24, 26, 29],
    [30, 28, 31],
    [32, 29, 31],
]

for x in Xinitialize:
    q0 = x[0]
    qc.h(q0)
    for q1 in x[1:]:
        qc.cx(q0, q1)

# Place stabilizers
for z in Zstab:
    qc.h(z[0])
    for q in z[1]:
        qc.cz(z[0], q)
    qc.h(z[0])

for x in Xstab:
    qc.h(x[0])
    for q in x[1]:
        qc.cx(x[0], q)
    qc.h(x[0])

# logical X-operator
logical_x = True
if logical_x:
    qc.x(11)
    qc.x(16)
    qc.x(21)

# Measurement for logical bit of # [51, [6,8,9,11]],  #  defect Z-cut hole 1
qc.h(51)
qc.cz(51, 6)
qc.cz(51, 8)
qc.cz(51, 9)
qc.cz(51, 11)
qc.h(51)

# Measurement stabilizer errors
for m, q in enumerate(Zstab + Xstab):
    qc.measure(q[0], m)

# Z-Measurement of logical bit (Observe bit-32 change, with no stabilizer errors)
qc.measure(51, 32)

backend = Aer.get_backend('qasm_simulator')
job = execute(qc, backend, shots=1000)
job.result().get_counts(qc)

