import numpy as np
from sympy import *
from sympy.solvers import solve
import random
from qiskit import *
import scipy as sp
from numpy import linalg as la
from scipy import linalg as sla
from IPython.display import clear_output
from numpy import linalg as la
from qiskit.tools.visualization import plot_histogram
from qiskit.circuit import QuantumCircuit, ParameterVector
from qiskit.compiler import transpile, assemble
from qiskit.providers.ibmq.managed import IBMQJobManager
from qiskit.visualization import *
from operator import itemgetter
from functools import reduce

# Pauli Matrices
Id = np.matrix([[1,0],[0,1]],dtype=complex)
X = np.matrix([[0,1],[1,0]],dtype=complex)
Y = np.matrix([[0,-1j],[1j,0]],dtype=complex)
Z = np.matrix([[1,0],[0,-1]],dtype=complex)



# Generate random values. 
a,b,c,d = [random.uniform(-1,1) for n in range(4)]
print('a=',a,',b=',b,',c=',c,',d=',d)

# The values used in our paper. 
a= -0.08496099317164885 
b= -0.8913449722923625 
c= 0.2653567873681286 
d= 0.5720506205204723

a = 2.0641811586466403
b = -0.8795521785527299
c = 0.2618460278326864
d = 0.5644821984323585

# Generate random values. 
c_i,c_x,c_y,c_z = [random.uniform(-1,1) for n in range(4)]
print('c_i=', c_i,',c_x=',c_x,',c_y=',c_y,',c_z=',c_z)

# The values used in our paper. 
c_i= -0.8453688208879107 
c_x= 0.006726766902514836 
c_y= -0.29354208134440296 
c_z= 0.18477436355081123

# The Hamiltonian. 
H = a*Id+b*X+c*Y+d*Z
print(H)

# Find the exact eigenstates and eigenvalues. 
Awh, Avh = la.eig(H)
Aeval0h = float(re(Awh[0]))
Aeval1h = float(re(Awh[1]))
Aevec0h = Avh[:,0]
Aevec1h = Avh[:,1]
print('The energy eigenvalues of H are:',Aeval0h, 'and', Aeval1h)
print('The corresponding eigenvectors are:', '\n', Aevec0h,'\n', 'and','\n', Aevec1h)

O = c_i*Id+c_x*X+c_y*Y+c_z*Z
print(O)

Aw, Av = la.eig(O)
Aeval0 = float(re(Aw[0]))
Aeval1 = float(re(Aw[1]))
Aevec0 = Av[:,0]
Aevec1 = Av[:,1]
print('The energy eigenvalues are:',Aeval0, 'and', Aeval1)
print('The corresponding eigenvectors are:', '\n', Aevec0,'\n', 'and','\n', Aevec1)



Zero = np.matrix([[1],[0]],dtype=complex)
One = np.matrix([[0],[1]],dtype=complex)
Ze = np.squeeze(np.asarray(Zero))
On = np.squeeze(np.asarray(One))

def H_obj(phi):  # All the information needed for H_obj. 
    # Hamiltonian.
    Hobj = H + phi*O
    wo, vo = la.eig(Hobj)
    # Eigenvalues.
    eval0o = float(re(wo[0]))
    eval1o = float(re(wo[1]))
    # Eigenvectors. 
    evec0o = vo[:,0]
    evec1o = vo[:,1]
    # Direct calculation of the expectation value. 
    Expe0 = evec0o.getH() @ O @ evec0o
    Expe1 = evec1o.getH() @ O @ evec1o
    return eval0o, eval1o, evec0o, evec1o, Expe0, Expe1

def GH(phimax,num,sigma):
    plt.figure(figsize=(16.18,10)) 
    philist = np.linspace(-phimax,phimax,num)
 #   epslist = np.concatenate((epslist[epslist<0], [0], epslist[epslist>0]))
    x = np.arange(-1.5,1.5,0.005)
    
    cm = LinearSegmentedColormap.from_list('defcol', ["#979A9A", "#17F514"])
    trace_color = cm(np.linspace(0,1,len(philist)))
    evalf = []
    evals = []
    expfi = []
    expse = []
    PS0 = []
    PS1 = []
    VEC0 = []
    VEC1 = []
    for phi, color in zip(philist, trace_color):
        val0,val1,vec0,vec1,expe0,expe1 = H_obj(phi)
        evalf.append(val0)
        evals.append(val1)
        expfi.append(expe0)
        expse.append(expe1)
        EV0 = np.squeeze(np.asarray(vec0))
        EV1 = np.squeeze(np.asarray(vec1))
        Pini0 = simplify(np.inner(Ze,EV0)*conjugate(np.inner(Ze,EV0)))
        Pini1 = simplify(np.inner(Ze,EV1)*conjugate(np.inner(Ze,EV1)))
        PS0.append(Pini0)
        PS1.append(Pini1)
        Probs = [Pini0,Pini1]
        Evals = [val0,val1]
        VEC0.append(vec0)
        VEC1.append(vec1)
        Yt = 0
        for k in range (2):
            y = Probs[k]
            for i in range (3):
                y *= 0.5*(1+np.exp((-0.5)*(x-Evals[k])**2*sigma**2)) 
            Yt += y  
        if phi == 0:
            expc = plt.plot(x,Yt,color = 'red',linewidth=3.6)

        else:
            expc = plt.plot(x,Yt,color = color)
    plt.xlabel(r'$E_{target}$',fontsize = 25)
    plt.ylabel(r'$P_{success}$',fontsize = 25)
    plt.grid()
    plt.show()
    
    dx = 4*phimax/(num-1)
    drf = [(evalf[i+1]-evalf[i-1])/dx for i in range(1,num-1)]
    drs = [(evals[i+1]-evals[i-1])/dx for i in range(1,num-1)]
    Df = dict(zip(philist[1:-1],drf))
    Ds = dict(zip(philist[1:-1],drs))
    Expf = dict(zip(philist[1:-1],expfi[1:-1]))
    Exps = dict(zip(philist[1:-1],expse[1:-1]))
    Energy1 = dict(zip(philist[1:-1],evalf[1:-1])) # Energy Eigenvalues
    Energy2 = dict(zip(philist[1:-1],evals[1:-1]))
    Pros0 = dict(zip(philist[1:-1],PS0[1:-1])) # Energy Eigenvalues
    Pros1 = dict(zip(philist[1:-1],PS1[1:-1]))

    return Df, Ds, Expf, Exps, Energy1, Energy2,Pros0,Pros1,VEC0,VEC1
    

Dc01, Dc02, Dc03, Dc04, Espec1, Espec2, PR0,PR1,VL0,VL1 = GH(0.1,21,12)



N = 3       # Number of cycles. 
obj = 3     # Numbering of the physical (object) qubit on quantum device.
arena = 1   # Numbering of the physical (arena) qubit on quantum device.



def qc1(E_targ,tstep,phi):                      
    
    ao = a + phi*c_i
    bo = b + phi*c_x
    co = c + phi*c_y
    do = d + phi*c_z
    
    sigo = Matrix([[bo],[co],[do]])
    n_hato = np.dot(1/sqrt(bo**2+co**2+do**2),sigo)
    nxo, nyo, nzo = [float(i.item()) for i in n_hato]    
    
    theta = 2*tstep*sqrt(bo**2+co**2+do**2) % (4*pi)   
    delta = float(re(atan(nzo*tan(theta/2))+atan(nxo/nyo)))
    beta = float(re(atan(nzo*tan(theta/2))-atan(nxo/nyo)))
    gamma = float(2*re(acos(cos(theta/2)/cos((delta+beta)/2)))) 
    if nyo > 0 and theta > 2*pi:
        gamma *= (-1)
    elif nyo < 0 and theta <= 2*pi:
        gamma *= (-1)   
    
    quanc1t = QuantumCircuit(2)
    quanc1t.h(1)
    quanc1t.cu(gamma,beta,delta,0,1,0)
    quanc1t.p(tstep*(E_targ-ao)-0.5*(delta+beta),1)               # The Energy Sensor.
    quanc1t.h(1) 
    quanc1t.barrier()
    return quanc1t

def qct(E_targ,sigma,phi):              
    qrz = QuantumRegister(2,'q')
    crz = ClassicalRegister(N,'c')
    qc = QuantumCircuit(qrz,crz)
    t = np.random.normal(0, sigma, N).tolist() 
    init_gate = Initialize([i.item() for i in Aevec0h])
    qc.append(init_gate, [1])
 #   qc.initialize([i.item() for i in Aevec0h],1)
   # qc.append(init_gate, [1])
    for i in range (N):
   #     sub_inst = qc1(E_targ,t[i],epsilon).to_instruction()
   #     qc.append(sub_inst,[obj,arena])
        qc.append(qc1(E_targ,t[i],phi),[0,1])   # 0:Obj, 1:Arena.
        qc.measure(1,i)    # Mid-circuit measurement, arena is measured. 
        qc.barrier()

    qc = transpile(qc,backend=backend,optimization_level=2,initial_layout = [obj,arena])
  #  qc = transpile(qc,backend=backend,optimization_level=2,layout_method = "noise_adaptive")
   # qc.barrier()
    return qc

def numsli_expc (L,rang1,rang2,nuni,nshots,sigma,k,phi):   
    
    val0, val1,_,_,_,_ = H_obj(phi)
    EW = np.linspace(rang1,rang2,L)   # The energy spectrum we're interested in.
    circuits = []                     # The 'circuit matrix'. 
    
    for i in range(L):
        circuitsi = []
        for j in range (nuni):
            circuit =  qct(EW[i],sigma,phi)      
            circuitsi.append(circuit)
        circuits.append(circuitsi)
        
    all_circuits = []
    for i in range(L):
        all_circuits += circuits[i][:]
    
    
    MExperiments = job_manager.run(all_circuits, backend=backend, shots = nshots)
    results = MExperiments.results() 
      
  #  cresults = results.combine_results()             # We are not doing read-out error mitigation here.   
  #  mitigated_results = meas_filter.apply(cresults)
        
    probsu = []
    succs = '0'*N
    for i in range (L):
        probsui = []
        for j in range (nuni):
       #     counts = mitigated_results.get_counts((nuni*i)+j)
            counts = results.get_counts((nuni*i)+j)
            
            if succs in counts:
                prob = counts[succs]/sum(counts.values())
            else:
                prob = 0
            probsui.append(prob)
        probsu.append(np.average(probsui))
            
    spec = dict(zip(EW,probsu))
    
    
    plt.figure(figsize=(16.18,10)) 
    plt.xlabel(r'$Energy\ Spectrum\ (E_{target})$',fontsize = 15)
    plt.ylabel('Probability to measure the success state',fontsize = 15)
    hist = plt.bar(spec.keys(), spec.values(), width=np.abs((rang2-rang1)/L), color='lightseagreen',alpha = 0.7)  
    plt.legend([hist],['Results from qasm'])
    plt.grid()
    plt.show()
    
    # Weighted avg based on the highest bins. 
    #Larg = dict(sorted(spec.items(), key = itemgetter(1), reverse = True)[:Hisk])
    #Wavg = np.average(list(Larg.keys()), weights=list(Larg.values()))
    
    # Weighted avg of the highest consecutive bins. 
    keys = np.fromiter(spec.keys(), float, count=len(spec))
    vals = np.fromiter(spec.values(), float, count=len(spec))
    idx = np.convolve(vals, np.ones(k), 'valid').argmax()
    Wavg = np.average(keys[idx:idx + k], weights=vals[idx:idx + k])
    
    print('Range for the next scan', keys[idx:idx + k])
    print('The estimated Eigenvalue is', Wavg)
    print('The precise values are',val0,'and',val1)
    return spec

# (L,rang1,rang2,nuni,nshots,sigma,k,phi) 9/26 qasm
firstscan = numsli_expc (18,0,np.pi+1,25,2500,2,6,0)
print(firstscan)

firstscan = numsli_expc (18,0,np.pi+1,15,2500,2,6,0)
print(firstscan)

# (L,rang1,rang2,nuni,nshots,sigma,k,phi) 7/22
firstscan = numsli_expc (18,-1-np.sqrt(3),1+np.sqrt(3),25,2500,2,6,0)

print(firstscan) #7/22

secondscan_a = numsli_expc (18, -1.768, -0.161 ,25,2500,7,6,0)     # Small Peak. 7/22
print(secondscan_a)

secondscan_b = numsli_expc (18, 0.161, 1.768 ,25,2500,7,6,0)     # Large Peak. 
print(secondscan_b)



thirdscan_a = numsli_expc (18, -1.390, -0.917 ,25,2500,12,6,0)     # Small Peak. 
print(thirdscan_a)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,0)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,0)     # Large Peak. 
print(thirdscan_b)



firstscan_722 = {-2.732050807568877: 0.11027200000000001, -2.4106330655019503: 0.118928, -2.089215323435024: 0.12784, -1.767797581368097: 0.184432, -1.4463798393011702: 0.264032, -1.1249620972342433: 0.300608, -0.8035443551673167: 0.254656, -0.48212661310039007: 0.200192, -0.1607088710334632: 0.14988800000000002, 0.16070887103346365: 0.19345600000000002, 0.4821266131003905: 0.37512, 0.8035443551673169: 0.5875360000000001, 1.1249620972342438: 0.6517439999999999, 1.4463798393011706: 0.499424, 1.767797581368097: 0.22110400000000002, 2.0892153234350244: 0.101648, 2.410633065501951: 0.11871999999999998, 2.732050807568877: 0.17247999999999997}
secondscan_722a = {-1.768: 0.131792, -1.673470588235294: 0.115936, -1.5789411764705883: 0.122, -1.4844117647058823: 0.09486399999999998, -1.3898823529411766: 0.14913600000000002, -1.2953529411764706: 0.22596799999999997, -1.2008235294117648: 0.239616, -1.1062941176470589: 0.247136, -1.011764705882353: 0.223568, -0.9172352941176471: 0.16216, -0.8227058823529413: 0.08616000000000001, -0.7281764705882354: 0.09582399999999999, -0.6336470588235295: 0.10512, -0.5391176470588237: 0.10432000000000001, -0.44458823529411773: 0.09491199999999997, -0.350058823529412: 0.08795199999999999, -0.255529411764706: 0.08761599999999997, -0.161: 0.07192}
secondscan_722b = {0.161: 0.182816, 0.2555294117647059: 0.15763200000000002, 0.35005882352941176: 0.10408, 0.4445882352941176: 0.13350399999999998, 0.5391176470588235: 0.08361600000000001, 0.6336470588235293: 0.12872, 0.7281764705882353: 0.15812800000000002, 0.8227058823529412: 0.24163199999999999, 0.917235294117647: 0.554208, 1.011764705882353: 0.695728, 1.1062941176470587: 0.509696, 1.2008235294117646: 0.296912, 1.2953529411764706: 0.184352, 1.3898823529411763: 0.116848, 1.4844117647058823: 0.10819200000000001, 1.578941176470588: 0.145056, 1.673470588235294: 0.09142399999999999, 1.768: 0.09601599999999999}
thirdscan_722a = {-1.39: 0.105148, -1.362176470588235: 0.132608, -1.3343529411764705: 0.163096, -1.3065294117647057: 0.18905999999999998, -1.2787058823529411: 0.188512, -1.2508823529411763: 0.227852, -1.2230588235294118: 0.238048, -1.195235294117647: 0.311544, -1.1674117647058824: 0.288476, -1.1395882352941176: 0.27693999999999996, -1.111764705882353: 0.22833199999999998, -1.0839411764705882: 0.21373999999999999, -1.0561176470588234: 0.20185599999999998, -1.0282941176470588: 0.18544799999999997, -1.0004705882352942: 0.1288, -0.9726470588235294: 0.14307599999999998, -0.9448235294117647: 0.12436799999999998, -0.917: 0.121528}
thirdscan_722b = {0.823: 0.221904, 0.8507647058823529: 0.16145600000000002, 0.8785294117647058: 0.23401600000000003, 0.9062941176470588: 0.307184, 0.9340588235294117: 0.42928000000000005, 0.9618235294117646: 0.533504, 0.9895882352941177: 0.6538719999999999, 1.0173529411764706: 0.6308799999999999, 1.0451176470588235: 0.577024, 1.0728823529411764: 0.438288, 1.1006470588235293: 0.30886399999999997, 1.1284117647058824: 0.261488, 1.1561764705882354: 0.21272, 1.1839411764705883: 0.207696, 1.2117058823529412: 0.144432, 1.239470588235294: 0.117392, 1.267235294117647: 0.12107200000000003, 1.295: 0.13726400000000002}
firstscan_926qasm = {-2.732050807568877: 0.095152, -2.4106330655019503: 0.109664, -2.089215323435024: 0.147952, -1.767797581368097: 0.19204800000000002, -1.4463798393011702: 0.315872, -1.1249620972342433: 0.33392000000000005, -0.8035443551673167: 0.301248, -0.48212661310039007: 0.22039999999999998, -0.1607088710334632: 0.14299199999999998, 0.16070887103346365: 0.27012800000000003, 0.4821266131003905: 0.41830399999999995, 0.8035443551673169: 0.6923039999999999, 1.1249620972342438: 0.7712320000000001, 1.4463798393011706: 0.47627200000000003, 1.767797581368097: 0.21656000000000003, 2.0892153234350244: 0.157872, 2.410633065501951: 0.134464, 2.732050807568877: 0.154464}



#plt.figure(figsize=(5*(1+np.sqrt(5)),10)) 
plt.figure(figsize=(14*1.2,12)) 
hist1 = plt.bar(firstscan_722.keys(), firstscan_722.values(), width = 2*(1+sqrt(3))/18 , color='limegreen',alpha = 0.32)  
plt.bar(secondscan_722a.keys(), secondscan_722a.values(), width= 1.608/18, color='deepskyblue',alpha = 0.49) 
hist2 = plt.bar(secondscan_722b.keys(), secondscan_722b.values(), width= 1.607/18, color='deepskyblue',alpha = 0.49) 
plt.bar(thirdscan_722a.keys(), thirdscan_722a.values(), width=0.5/18, color='darkorchid',alpha = 0.65) 
hist3 = plt.bar(thirdscan_722b.keys(), thirdscan_722b.values(), width=0.5/18, color='darkorchid',alpha = 0.65) 
hist_ref = plt.bar(firstscan_926qasm.keys(), firstscan_926qasm.values(),width = 2*(1+sqrt(3))/18,fill = False , color='grey',alpha = 0.49) 

Zero = np.matrix([[1],[0]],dtype=complex)
One = np.matrix([[0],[1]],dtype=complex)
Ze = np.squeeze(np.asarray(Zero))
On = np.squeeze(np.asarray(One))

x = np.arange(-1-np.sqrt(3),1+np.sqrt(3),0.005)
EV0 = np.squeeze(np.asarray(Aevec0h))
EV1 = np.squeeze(np.asarray(Aevec1h))
Pini0 = np.inner(Ze,EV0)*conjugate(np.inner(Ze,EV0))
Pini1 = np.inner(Ze,EV1)*conjugate(np.inner(Ze,EV1))
Probs = [Pini0,Pini1]
Evals = [Aeval0h,Aeval1h]
Yt = 0
for k in range (2):
    y = Probs[k]
    for i in range (N):
        y *= 0.5*(1+np.exp((-0.5)*(x-Evals[k])**2*2**2)) 
    Yt += y  
    
Ytl = 0
for k in range (2):
    yl = Probs[k]
    for i in range (N):
        yl *= 0.5*(1+np.exp((-0.5)*(x-Evals[k])**2*7**2)) 
    Ytl += yl 
    
Yt2 = 0
for k in range (2):
    y2 = Probs[k]
    for i in range (N):
        y2 *= 0.5*(1+np.exp((-0.5)*(x-Evals[k])**2*12**2)) 
    Yt2 += y2 

expc, = plt.plot(x,Yt,color = 'royalblue',linewidth=2.7,linestyle='--',alpha=1)
expc1, = plt.plot(x,Ytl,color = 'red',linewidth=2.7,linestyle='--',alpha=1)
expc2, = plt.plot(x,Yt2,color = 'black',linewidth=2.7,linestyle='--',alpha=1)
#plt.xlabel(r'$E$',fontsize = 35)
#plt.ylabel(r'$P_{0^N}$',fontsize = 35)


plt.legend([hist1,hist2,hist3,hist_ref,expc,expc1,expc2],
           ['$\sigma=2$ (1st scan)','$\sigma=7$ (2nd scan)','$\sigma=12$ (3rd scan)','Noiseless simulation (1st scan)','$\sigma=2$ (expected)','$\sigma=7$ (expected)','$\sigma=12$ (expected)'], 
           ncol = 2,loc="upper left",fontsize=20,columnspacing = -4.5)
#plt.grid(alpha = 0.6)
plt.xticks(fontsize=25)
plt.yticks(fontsize=25)
plt.text(0.087, 0.9, r'$P_{0^N}$', fontsize=35, transform=plt.gcf().transFigure)
plt.text(0.9, 0.087, r'$E$', fontsize=35, transform=plt.gcf().transFigure)
plt.savefig("SeqScan430_0610.png",dpi=430,bbox_inches='tight',pad_inches = 0)
plt.show()


thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,0.01)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,0.01)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.01)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.01)     # Large Peak. 
print(thirdscan_b)



thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,0.02)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,0.02)     # Large Peak. 
print(thirdscan_b)



thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.02)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.02)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,0.03)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,0.03)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.03)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.03)     # Large Peak. 
print(thirdscan_b)





thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.04)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.04)     # Large Peak. 
print(thirdscan_b)



thirdscan_a = numsli_expc (18, -1.490, -1.017 ,50,5000,12,6,0.04)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.723, 1.195 ,25,2500,12,6,0.04)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.490, -1.017 ,50,5000,12,6,0.05)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.723, 1.195 ,25,2500,12,6,0.05)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.05)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.05)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.490, -1.017 ,50,5000,12,6,0.06)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.723, 1.195 ,25,2500,12,6,0.06)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.06)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.06)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.490, -1.017 ,50,5000,12,6,0.07)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.723, 1.195 ,25,2500,12,6,0.07)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.07)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.07)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.490, -1.017 ,50,5000,12,6,0.08)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.723, 1.195 ,25,2500,12,6,0.08)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.390, -0.917 ,50,5000,12,6,-0.08)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.08)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.490, -1.017 ,50,5000,12,6,0.09)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.723, 1.195 ,25,2500,12,6,0.09)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.290, -0.817 ,50,5000,12,6,-0.09)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.09)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.490, -1.017 ,50,5000,12,6,0.1)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.723, 1.195 ,25,2500,12,6,0.1)     # Large Peak. 
print(thirdscan_b)

thirdscan_a = numsli_expc (18, -1.290, -0.817 ,50,5000,12,6,-0.1)     # Small Peak. 
print(thirdscan_a)

thirdscan_b = numsli_expc (18, 0.823, 1.295 ,25,2500,12,6,-0.1)     # Large Peak. 
print(thirdscan_b)



