%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
# Importing standard Qiskit libraries and configuring account
from qiskit import *
from qiskit.visualization import *
from qiskit.circuit import Parameter

# Simple example

# Define two parameters, t1 and t2
theta1 = Parameter('t1')
theta2 = Parameter('t2')

# Build a 1-qubit circuit
qc = QuantumCircuit(1, 1)

# First parameter, t1, is used for a single qubit rotation of a controlled qubit
qc.ry(theta1,0)
qc.rz(theta2,0)

qc.barrier()

qc.measure(0, 0)

qc.draw(output='mpl')

theta1_range = np.linspace(0, 2 * np.pi, 20)
theta2_range = np.linspace(0, np.pi, 2)

circuits = [qc.bind_parameters({theta1: theta_val1, theta2: theta_val2})
            for theta_val2 in theta2_range for theta_val1 in theta1_range ]

# Visualize several circuits to check that correct circuits are generated correctly.
display(circuits[0].draw(output='mpl'))
display(circuits[1].draw(output='mpl'))
display(circuits[20].draw(output='mpl'))

# Execute multiple circuits
job = execute(qc, backend=Aer.get_backend('qasm_simulator'), shots = 8192,
              parameter_binds=[{theta1: theta_val1, theta2: theta_val2} 
                               for theta_val2 in theta2_range for theta_val1 in theta1_range])

# Store all counts
counts = [job.result().get_counts(i) for i in range(len(job.result().results))]

# Plot to visualize the result
plt.figure(figsize=(12,6))
plt.plot(range(len(theta1_range)*len(theta2_range)),
         list(map(lambda counts: (counts.get('0',0)-counts.get('1',1))/8192,counts)))
plt.show()

# IBMQ.disable_account()
provider = IBMQ.enable_account('TOKEN')

from qiskit.tools.monitor import backend_overview, backend_monitor, job_monitor
from qiskit.tools.visualization import plot_gate_map, plot_error_map

# Retrieve IBM Quantum device information
backend_overview()

# Execute multiple circuits
job_exp = execute(qc, backend=provider.get_backend('ibmq_essex'), shots = 8192,
              parameter_binds=[{theta1: theta_val1, theta2: theta_val2} 
                               for theta_val2 in theta2_range for theta_val1 in theta1_range])

# Monitor job status
job_monitor(job_exp)

# Store all counts
counts_exp = [job_exp.result().get_counts(i) for i in range(len(job_exp.result().results))]

# Plot to visualize the result
plt.figure(figsize=(12,6))
plt.rcParams.update({'font.size': 16})

plt.plot(range(len(theta1_range)*len(theta2_range)),
         list(map(lambda counts: (counts.get('0',0)-counts.get('1',1))/8192,counts)),'r',label='Simulation')
plt.plot(range(len(theta1_range)*len(theta2_range)),
         list(map(lambda counts_exp: (counts_exp.get('0',0)-counts_exp.get('1',1))/8192,counts_exp)),
         'b',label='Experiment')
plt.legend(loc='best')
plt.show()

def swaptest(qc,ancilla,qubit1,qubit2):
    qc.h(ancilla)
    qc.cswap(ancilla,qubit1,qubit2)
    qc.h(ancilla)

qr = QuantumRegister(5,'q')
cr = ClassicalRegister(2,'c')
qc = QuantumCircuit(qr,cr,name='qclassifier')

# Put equal weights on the training data
qc.h(1)

# Prepare the test data quantum state
theta1 = Parameter('t1')
qc.rx(theta1,4)

# Prepare the training data quantum state
qc.h(2)
qc.rz(-np.pi,2)
qc.s(2)
qc.cz(1,2)

# Put the label
qc.cx(1,3)
qc.barrier()

# Perform swap-test
swaptest(qc,0,2,4)
qc.barrier()

# Measurement
qc.measure(0,0)
qc.measure(3,1)

qc.draw(output='mpl')

# Select the range of parameters
theta1_range = np.linspace(0, 2 * np.pi, 20)

# Execute multiple circuits
job = execute(qc, backend=Aer.get_backend('qasm_simulator'), shots = 8192,
              parameter_binds=[{theta1: theta_val1} for theta_val1 in theta1_range])

# Store all counts
counts = [job.result().get_counts(i) for i in range(len(job.result().results))]

# Let's see how the 'counts' looks like
counts

# Plot to visualize the result
plt.figure(figsize=(12,6))
plt.rcParams.update({'font.size': 20})
# How to extract <ZZ>? Answer: <ZZ> = P(00) - P(01) - P(10) + P(11)
plt.plot(theta1_range,
         list(map(lambda counts: (counts.get('00')-counts.get('10')-counts.get('01')+counts.get('11'))/8192,
                  counts)))
plt.rc('text', usetex=True)
plt.xlabel(r'$\theta$ (Parameter for test data)')
plt.ylabel(r'$ \langle Z\otimes Z \rangle $')
plt.show()

# Let's also try the same experiment on the 5-qubit device.
# Select the range of parameters
theta1_range = np.linspace(0, 2 * np.pi, 20)

# Execute multiple circuits
job_classifier = execute(qc, backend=provider.get_backend('ibmq_essex'), shots = 8192,
              parameter_binds=[{theta1: theta_val1} for theta_val1 in theta1_range])

# Monitor job status
job_monitor(job_classifier)

# Store all counts
counts_exp = [job_classifier.result().get_counts(i) for i in range(len(job_classifier.result().results))]

# Plot to visualize the result
plt.figure(figsize=(12,6))
plt.rcParams.update({'font.size': 20})

# Plot theory
plt.plot(theta1_range,
         list(map(lambda counts: (counts.get('00')-counts.get('10')-counts.get('01')+counts.get('11'))/8192,
                  counts)),'r',label='Theory')

# Plot experimental result
plt.plot(theta1_range,
         list(map(lambda counts_exp: 4*(counts_exp.get('00')-counts_exp.get('10')
                                      -counts_exp.get('01')+counts_exp.get('11'))/8192,counts_exp)),'b',
         label='Exp. (scaled by a factor of 4)')

plt.rc('text', usetex=True)
plt.xlabel(r'$\theta$ (Parameter for test data)')
plt.ylabel(r'$ \langle Z\otimes Z \rangle $')
plt.legend(loc='best')
plt.show()

