<jupyter_start><jupyter_text>Template (Place Useful Header Here) This template file serves as our "portable library." Whenever you want to run a simulation, just make a copy of this file and adjust it to your needs in order to run an experiment. The structure is as follows:1. Imports at top2. Your Experiment in Middle3. All Functions in function_list hosted together at bottomNote: Be sure to actually run the notebook before working on your own experiments to have access to the functions. Oh, and in the future remove this note and add description of experiment, etc., etc.<jupyter_code># these are all the imports needed across all experiments... some may be unnecessary for any given test
import numpy as np
import math
import datetime
import random
import matplotlib.pyplot as plt
import matplotlib
from qiskit import IBMQ, execute, Aer, QuantumCircuit, transpile, schedule as build_schedule, assemble
from qiskit import BasicAer
from qiskit.visualization import plot_histogram
from qiskit.providers.aer.noise import NoiseModel
from qiskit.pulse import Schedule, Play, DriveChannel, SamplePulse, ConstantPulse, Delay
from qiskit.pulse import InstructionScheduleMap
from qiskit.tools.monitor import job_monitor
%config InlineBackend.figure_format = 'svg' # makes the figures look nice
matplotlib.rcParams['text.usetex'] = False
%matplotlib inline<jupyter_output><empty_output><jupyter_text>Place your experiment runs here!!! The functions I wrote are broken up into three sections, but since this update strategy sucks, just look at function_list to see them organized better.<jupyter_code># some useful 'globals' or future commands to use
backend_list = ['ibmq_armonk', 'ibmq_ourense']
# inst_map = backend.defaults().instruction_schedule_map

def load_backend(name, info_dump=True):
    """
    Loads backend called (name).
    """
    provider = IBMQ.load_account()
    backend = provider.get_backend(name)
    now = datetime.datetime.now()
    if info_dump is True:
        print("Backend loaded at date/time : ")
        print(now.strftime("%Y-%m-%d %H:%M:%S"))
        print("---------------------------------")
        print("Backend Configuration Information")
        print("---------------------------------")
        print(backend.configuration())
        print("Backend Properties Information")
        print(backend.properties())
    
    return backend

def load_noisemodel(backend):
    """
    Given a backend, loads in the necessary information to run a noisy simulation emulating noise on actual
    device.
    """
    # set-up noise model for simulator
    noise_model = NoiseModel.from_backend(backend)
    # Get coupling map from backend
    coupling_map = backend.configuration().coupling_map
    # Get basis gates from noise model
    basis_gates = noise_model.basis_gates
    noise_info = {'noise_model': noise_model, 'coupling_map': coupling_map, 'basis_gates': basis_gates}
    
    return noise_info

def get_max_runs(backend):
    """
    Given the backend, returns the max # experiments and max # shots can be queued in a single job.
    """
    
    max_experiments = backend.configuration().max_experiments
    max_shots = backend.configuration().max_shots
    max_runs = {'max_experiments': max_experiments, 'max_shots': max_shots}
    
    return max_runs

def frev(circ, qubits, num_ids=1):
    """
    Appends free evolution to [qubits] over [n] identity gates.

    Inputs:
    * qubits -- list, qubits to add identities to
    * num_ids -- int, number of identities to append

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_ids):
        ch_depth += 1
        circ.barrier(qubits)
        circ.id(qubits)
        
    # create channel tag
    ch_tag = 'frev_{}'.format(num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}


def xy4(circ, qubits, num_rep=1, num_ids=0):
    """
    Appends standard xy4 sequence to [qubits] for [num_rep] times with [num_ids] in between
    each DD pulse gate.

    Inputs:
    * qubits -- list, qubits to append to
    * num_rep -- int, desired number of times to repeat xy4 sequence
    * num_ids -- int, number of identities to pad in between DD pulses

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_rep):
        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.y(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.x(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.y(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.x(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']
        
    # add channel tag
    ch_tag = 'xy4_{}r_{}i'.format(num_rep, num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def rga2x(circ, qubits, num_rep=1, num_ids=0):
    """
    Appends RGA2 sequence to [qubits] for [num_rep] times with [num_ids] in between
    each DD pulse gate.
    Sequence: Xb X
    NOTE: Xb means X-bar which is pi phase flip on all axes of Y gate.
    
    Inputs:
    * qubits -- list, qubits to append to
    * num_rep -- int, desired number of times to repeat xy4 sequence
    * num_ids -- int, number of identities to pad in between DD pulses

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_rep):
        # add Xb gate
        ch_depth += 1
        circ.barrier(qubits)
        add_xb(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.x(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']
        
    # add channel tag
    ch_tag = 'rga2x_{}r_{}i'.format(num_rep, num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def rga2y(circ, qubits, num_rep=1, num_ids=0):
    """
    Appends RGA2 sequence to [qubits] for [num_rep] times with [num_ids] in between
    each DD pulse gate.
    Sequence: Yb Y
    NOTE: Xb means X-bar which is pi phase flip on all axes of Y gate.
    
    Inputs:
    * qubits -- list, qubits to append to
    * num_rep -- int, desired number of times to repeat xy4 sequence
    * num_ids -- int, number of identities to pad in between DD pulses

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_rep):
        # add Yb gate
        ch_depth += 1
        circ.barrier(qubits)
        add_yb(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.y(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']
        
    # add channel tag
    ch_tag = 'rga2y_{}r_{}i'.format(num_rep, num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def rga2z(circ, qubits, num_rep=1, num_ids=0):
    """
    Appends RGA2z sequence to [qubits] for [num_rep] times with [num_ids] in between
    each DD pulse gate.
    Sequence: Zb Z
    NOTE: Xb means X-bar which is pi phase flip on all axes of Y gate.
    
    Inputs:
    * qubits -- list, qubits to append to
    * num_rep -- int, desired number of times to repeat xy4 sequence
    * num_ids -- int, number of identities to pad in between DD pulses

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_rep):
        # add zb gate
        ch_depth += 1
        circ.barrier(qubits)
        add_zb(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add z gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.z(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']
        
    # add channel tag
    ch_tag = 'rga2z_{}r_{}i'.format(num_rep, num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}


def rga4(circ, qubits, num_rep=1, num_ids=0):
    """
    Appends RGA4 sequence to [qubits] for [num_rep] times with [num_ids] in between
    each DD pulse gate.
    NOTE: We choose P2 = X and P1 = Y for no particular reason.
    NOTE: Xb means X-bar which is pi phase flip on all axes of Y gate.


    Inputs:
    * qubits -- list, qubits to append to
    * num_rep -- int, desired number of times to repeat xy4 sequence
    * num_ids -- int, number of identities to pad in between DD pulses

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_rep):
        # add Xb gate
        ch_depth += 1
        circ.barrier(qubits)
        add_xb(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.y(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Xb gate
        ch_depth += 1
        circ.barrier(qubits)
        add_xb(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.y(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']
        
    # add channel tag
    ch_tag = 'rga4_{}r_{}i'.format(num_rep, num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def rga8c(circ, qubits, num_rep=1, num_ids=0):
    """
    Appends RGA8c sequence to [qubits] for [num_rep] times with [num_ids] in between
    each DD pulse gate.
    NOTE: We choose P2 = X and P1 = Y for no particular reason.
    NOTE: Xb means X-bar which is pi phase flip on all axes of Y gate.
    That is -- RGA8c: Y X Y X X Y X Y


    Inputs:
    * qubits -- list, qubits to append to
    * num_rep -- int, desired number of times to repeat xy4 sequence
    * num_ids -- int, number of identities to pad in between DD pulses

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_rep):
        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        add_y(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.x(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        add_y(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.x(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.x(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.y(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.x(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        circ.y(qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']
        
    # add channel tag
    ch_tag = 'rga8c_{}r_{}i'.format(num_rep, num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def prep_theta_state(circ, qubits, theta=0):
    """
    Appends u3(theta, 0, 0) gate to [qubits] in order to prepare each qubit in
    1 / sqrt(2) * (cos(theta/2)|0> + sin(theta/2)|1>) state.

    Inputs:
    * qubits -- list, qubits to append to
    * theta -- float, superposition theta angle

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """
    ch_depth = 1
    circ.u3(theta, 0, 0, qubits)
    circ.barrier(qubits)
    # add channel tag
    ch_tag = 'ptheta_{}'.format(theta)
    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def decode_theta_state(circ, qubits, theta=0):
    """
    Appends u3(-theta, 0, 0) gate to [qubits] in order to decode qubits in
    1 / sqrt(2) * (cos(theta/2)|0> + sin(theta/2)|1>) state back to |0> state.

    Inputs:
    * qubits -- list, qubits to append to
    * theta -- float, superposition theta angle

    Output:
    * ch_depth -- int, depth added to circuit
    """
    ch_depth = 1
    circ.barrier(qubits)
    circ.u3(-1*theta, 0, 0, qubits)
    circ.barrier(qubits)
    # add channel tag
    ch_tag = 'dtheta_{}'.format(theta)
    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def xz4(circ, qubits, num_rep=1, num_ids=0):
        """
        Appends xz4 sequence to [qubits] for [num_rep] times with [num_ids] in between
        each DD pulse gate.

        Inputs:
        * qubits -- list, qubits to append to
        * num_rep -- int, desired number of times to repeat xy4 sequence
        * num_ids -- int, number of identities to pad in between DD pulses

        Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
        """

        ch_depth = 0
        for i in range(num_rep):
            # add Z gate
            ch_depth += 1
            circ.barrier(qubits)
            circ.z(qubits)

            # add [num_ids] I gates
            ch_depth += frev(circ, qubits, num_ids)['ch_depth']

            # add X gate
            ch_depth += 1
            circ.barrier(qubits)
            circ.x(qubits)

            # add [num_ids] I gates
            ch_depth += frev(circ, qubits, num_ids)['ch_depth']

            # add Z gate
            ch_depth += 1
            circ.barrier(qubits)
            circ.z(qubits)

            # add [num_ids] I gates
            ch_depth += frev(circ, qubits, num_ids)['ch_depth']

            # add X gate
            ch_depth += 1
            circ.barrier(qubits)
            circ.x(qubits)

            # add [num_ids] I gates
            ch_depth += frev(circ, qubits, num_ids)['ch_depth']
            
        # add channel tag
        ch_tag = 'xz4_{}r_{}i'.format(num_rep, num_ids)

        return {'ch_tag': ch_tag, 'ch_depth': ch_depth}

def native_xy4(circ, qubits, num_rep=1, num_ids=0):
    """
    Appends standard xy4 sequence to [qubits] for [num_rep] times with [num_ids] in between
    each DD pulse gate using native gates.

    Inputs:
    * qubits -- list, qubits to append to
    * num_rep -- int, desired number of times to repeat xy4 sequence
    * num_ids -- int, number of identities to pad in between DD pulses

    Output:
    * {'ch_tag': ch_tag, 'ch_depth': ch_depth)
    --> ch_tag -- str, channel tag with identifies this evolution
    --> ch_depth -- int, depth added to circuit with this channel
    """

    ch_depth = 0
    for i in range(num_rep):
        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        add_y(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        add_x(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add Y gate
        ch_depth += 1
        circ.barrier(qubits)
        add_y(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']

        # add X gate
        ch_depth += 1
        circ.barrier(qubits)
        add_x(circ, qubits)

        # add [num_ids] I gates
        ch_depth += frev(circ, qubits, num_ids)['ch_depth']
        
    # add channel tag
    ch_tag = 'n_xy4_{}r_{}i'.format(num_rep, num_ids)

    return {'ch_tag': ch_tag, 'ch_depth': ch_depth}


# Adding functions which implement all the gates we need in terms of the
# IBM native gates. These were done with armonk as of 1 May 2020 using
# transpile(circ, backend) for some basic gates.
def add_x(circ, qubits):
    """
    Appends x gate to circuit on [qubits] in terms of native u3 gate.
    """
    circ.u3(np.pi, 0, np.pi, qubits)
    return

def add_y(circ, qubits):
    """
    Appends y gate to circuit on [qubits] in terms of native u3 gate.
    """
    circ.u3(np.pi, np.pi/2, np.pi/2, qubits)
    return

def add_z(circ, qubits):
    """
    Appends z gate to circuit on [qubits] in terms of native u1 gate.
    """
    #NOTE:  u1(pi) = u3(0, pi, 0)
    circ.u1(np.pi, qubits)
    return

def add_xb(circ, qubits):
    """
    Appends x bar to circ on [qubits] in terms of native gates.
    xbar comes from RGAn sequences and is pi rotation of x around
    all axes.
    """
    #NOTE: u1(pi) = u3(0, pi, 0) = z
    circ.u1(np.pi, qubits)
    return

def add_yb(circ, qubits):
    """
    Appends y bar to circ on [qubits] in terms of native gates.
    ybar comes from RGAn sequences and is pi rotation of y around
    all axes.
    """
    #NOTE: u1(3pi) = u1(pi) = z
    circ.u1(3*np.pi, qubits)
    return 

def add_zb(circ, qubits):
    """
    Appends z bar to circ on [qubits] in terms of native gates.
    zbar comes from RGAx sequences and is pi rotation of z around
    all axes.
    """
    #NOTE: u3(pi, 0, pi) = x
    circ.u3(np.pi, 0, np.pi, qubits)
    return

def add_H(circ, qubits):
    """
    Appends hadamard to circ on [qubits] in terms of native u2 gate.
    """
    #NOTE: u2(0, pi) = u3(pi/2, 0, pi)
    circ.u2(0, np.pi, qubits)
    return

def submit_job(experiments, backend, qobj_id, num_shots='max'):
    """
    Submit [experiments] to [backend] and run [num_shots] for each of the [experiments]. Results get a
    qobj_id tag labelled with [qobj_id].
    """
    # set runs to max runs allowed by hardware if set to do so
    max_runs = get_max_runs(backend)
    if str(num_shots).lower() == 'max':
        num_shots = max_runs['max_shots']

    # submit the job in the background and output status information
    program = assemble(experiments, backend=backend, shots=num_shots, qobj_id=qobj_id)
    job = backend.run(program)
    print("job id: {}".format(job.job_id()))
    print(job_monitor(job))

    return job.result(timeout=120)

def submit_test(experiments, backend, qobj_id, num_shots = 1000):
    """
    Submits [experiments] to simulator with noisemodel of [backend] and runs [num_shots] for each
    """
    # set up noise model of backend
    noise_info = load_noisemodel(backend)
    coupling_map = noise_info['coupling_map']
    basis_gates = noise_info['basis_gates']
    noise_model = noise_info['noise_model']
    
    job = execute(experiments, Aer.get_backend('qasm_simulator'),
                     coupling_map=coupling_map,
                     basis_gates=basis_gates,
                     noise_model=noise_model,
                     shots=num_shots, qobj_id=qobj_id)
    
    return job.result()

def combine_experiments(results):
    """
    Combines the results over all experiments in a submitted job.
    """
    counts = results.get_counts() # results of all experiments
    totals = {} # data to combine all results

    # if only one experiment, qiskit resturns dict which is already what we want
    if type(counts) is dict:
        totals = counts

    # otherwise, need to iterate through experiments and create one single dictionary
    elif type(counts) is list:
        for ic in counts:
            print(ic)
            # iterate over key/value pairs in results of ith count
            for key, value in ic.items():
                if key in totals:
                    totals[key] += value
                else:
                    totals[key] = value
                    
    return totals

def save_raw_data(result):
    """
    Saves data from [result] into a raw data file named after the unique job id.
    """
    # save the raw data
    fname = (result.job_id + ".txt")
    dict_result = result.to_dict()
    with open(fname, 'w') as f:
        print(dict_result, file=f)
    print("Successfully saved raw data to file: {}".format(fname))

    return fname

def save_combined_experiments_data(result):
    """
    Saves data from [result] which combines all experiments into one data set named after job id with tag 'comb_sum.'
    """

    # create summary of data through counts and save summary file
    fname = (result.job_id + "_comb_sum.txt")
    totals = combine_experiments(result)
    with open(fname, 'w') as f:
        print(totals, file=f)
    print("Successfully saved combined experiments data to file: {}".format(fname))

    return fname

def save_dict_data(data, fname):
    """
    Saves [data] which is assumed to be a dictionary.
    """
    with open(fname, 'w') as f:
        print(data, file=f)
    print("Successfully saved combined experiments data to file: {}".format(fname))

    return

def read_dict(fname):
    """
    Reads data file formatted as dictionary back into python as a dictionary.
    NOTE: This uses eval()--which is NOT a secure option, but this is default python so requires no 
    imports and will be necessary since Bruno cannot run on local envrionment. Just don't use this
    in your own codes too often...
    """
    with open(fname, 'r') as f:
        data_dict = eval(f.read())

    return data_dict

def hex_to_bin(hexstr):
    """
    Converts a hex string into binary string.
    """
    binstr = "{0:08b}".format(int(hexstr, 16))
    return binstr

def hex_dict_to_bin_dict(hexdict):
    """
    Converts a dictionary with hexadecimal (str) keys into binary (str) keys.
    """
    return {hex_to_bin(key): hexdict[key] for key in hexdict}

def count_0s(istr):
    """
    Counts the number of 0s in [istr] and returns that integer.
    """
    num_0s = 0
    for s in istr:
        if s == '0':
            num_0s += 1
    return num_0s

def dict_data_to_n0_array(data_dict):
    """
    Converts a dictionary of the form {'101': 2, '000': 1} --> [101, 101, 000] --> [1, 1, 3] where last
    list specifies how many zeros each bitstring has. (We don't actually compute second list, but helps to see
    it this way, I think.)
    """
    # populate data_list with relative frequencies from dict
    data_list = []
    for key, value in data_dict.items():
        for i in range(value):
            # get number of zeros in this string and append to list
            data_list.append(count_0s(key))
            
    return np.array(data_list)

# my old bootstrapping code
def bootstrapci(data, n=1000, func=np.mean, ci=.95):
    """
    Generate `n` bootstrap samples, evaluating `func`
    at each resampling. Then computes error bar with
    'ci' confidence interval on data.
    """    
    simulations = np.zeros(n)
    sample_size = len(data)
    xbar = func(data)
    for c in range(n):
        itersample = np.random.choice(data, size=sample_size, replace=True)
        simulations[c] = func(itersample)
    diff = np.sort(simulations - xbar)
    lefterr_idx = int(np.floor(ci*n))
    righterr_idx = int(np.ceil((1 - ci)*n))
    
    conf_int = (xbar - diff[lefterr_idx], xbar - diff[righterr_idx])
    
    return conf_int

def per_err(mu, errbar):
    """
    Computes percent error with respect to mean. In particular, we calculate difference between
    (mu - errbar[0]) / mu and (mu - errbar[1]) / mu and takes the larger of the two in abs value then
    just multiply by 100.
    """
    diff1 = abs(mu - errbar[0])
    diff2 = abs(mu - errbar[1])
    
    perdif1 = (diff1 / mu) * 100
    perdif2 = (diff2 / mu) * 100
    
    if perdif1 > perdif2:
        return perdif1
    else:
        return perdif2
    
def calc_exper_fid(exper):
    """
    Given an single experiment's results from IBMQ, calculates the fidelity and computes bootstrapped error
    bars. By "fidelity" we mean the total number of zeros--that is, we assume experiment is designed to output
    0 at the end. (If you encode state--decode it!)
    For example, if 3 qubit experiment with 10 shots gives results {'101': 3, '000': 5, '111': 2}, then the
    total number of zeros is 3(1) + 5(3) + 2(0) = 18, and fid = 18 / (10 * 3).
    
    Input
    ----------------
    exper -- "results" of IBMQ experiemnt run, i.e. if you run result = submit_job(...).to_dict(), then
    this expects results['results'][j], i.e. the results of jth experiment.
    tol -- float, desired tolerance for data (i.e. 0.05 is 5% abs err mag tolerance)
    
    Output
    ----------------
    fid -- tuple, (fidelity, -.95ci, +.95ci, p_err) -(+).95ci is lower(upper) bound of 95% confidence
    interval via bootstrapping and err_mag is magnitude of error w.r.t. mean.
    
    Also prints if abs mag of error is less than input tol.
    """
    # results have hex keys--convert them to binary keys first
    counts = hex_dict_to_bin_dict(exper['data']['counts'])
    # get arbitary key of counts dict to counts number of qubits (i.e. length of bitstring)
    key0 = list(counts.keys())[0]
    n_qubits = len(key0)
    
    # obtain the number 0s from each shot of experiment and create a (degenerate) array
    # That is, if '001' shows up 50 times, populate array with 50 instances of 2
    num_0s = dict_data_to_n0_array(counts)
    # turn num0s to fideltiy by dividing by number of qubits
    fids = num_0s / n_qubits
    # calculate mean
    mean_fid = np.mean(fids)
    # calculate confidence interval with 1000 bootstrap samples and .95 confidence interval
    ci = bootstrapci(fids, 1000, np.mean, .95)
    # calculate percent error
    p_err = per_err(mean_fid, ci)
    
    return (mean_fid, ci[0], ci[1], p_err)

def theta_wrangle(result_dict):
    """
    wrangles thata data for real this time
    """
    # get global info about all experiments in this result (i.e. type of experiment run)
    batch_name = result_dict['qobj_id']
    job_id = result_dict['job_id']
    # iterate over the experimental results in results dict
    thetas = []
    fids = []
    l_ci = []
    u_ci = []
    p_errs = []
    for exper in result_dict['results']:
        # extract theta value from experiment label
        label = exper['header']['name']
        x = 'theta_'
        try:
            theta = label[label.find(x)+len(x):]
        except ValueError:
            raise ValueError("theta tag is not present in circuit names")
        thetas.append(float(theta))
        # get statistics on data
        stats = calc_exper_fid(exper)
        # append stats to relevant lists
        fids.append(stats[0])
        l_ci.append(stats[1])
        u_ci.append(stats[2])
        p_errs.append(stats[3])
        
    data = (thetas, fids, l_ci, u_ci, p_errs)
    # save summary of data to file with name as (job_id)_(sim_tag)_(theta_sweep).txt
    fname = job_id + "_" + batch_name + "_" + "theta_sweep"
    with open(fname + '.txt', 'w') as f:
        header = "theta, fidelity, lower confidence internval, upper ci, percent error\n"
        f.write(header)
        for i in range(len(thetas)):
            line = str(thetas[i]) + "," + str(fids[i]) + "," + str(l_ci[i])\
            + "," + str(u_ci[i]) + "," + str(p_errs[i]) + "\n"
            f.write(line)
    print("Successfully saved fid_data to file: {}".format(fname + '.txt'))
    
    # plot the data with error bars
    title = job_id + '_' + batch_name + '_theta_sweep'
    plt = plot_theta_sweep(data, title)
    
    return data

def plot_theta_sweep(data, title=''):
    """
    plots theta sweep given data as tuple of lists
    """
    # unpack data
    thetas, fids, l_ci, u_ci, p_errs = data
    # convert to arrays for easier manipulation
    thetas = np.array(thetas)
    fids = np.array(fids)
    l_ci= np.array(l_ci)
    u_ci = np.array(u_ci)
    p_errs = np.array(p_errs)
    # get difference between upper confidence internval and mean value for errbar formatting
    fid_up = u_ci - fids
    fid_low = fids - l_ci
    
    # make thetas into more readable format
    thetas = thetas / np.pi
    # finally plot it and add labels
    plt.errorbar(thetas, fids, yerr=(fid_low, fid_up))
    plt.xlabel('theta / pi')
    plt.ylabel('fidelity')
    plt.title(title)
    
    return plt<jupyter_output><empty_output>