diff --git a/microring_dynamics.py b/microring_dynamics.py
index 3eda9ca..fbb4c48 100644
--- a/microring_dynamics.py
+++ b/microring_dynamics.py
@@ -7,8 +7,7 @@ K = 0.1
 t = np.sqrt(T) + 0.0j
 k = 0.0 - 1j * np.sqrt(K)
 
-Coupler = np.array([[t, k],
-                    [k.conjugate(), -t.conjugate()]])
+Coupler = np.array([[t, k], [k.conjugate(), -t.conjugate()]])
 
 num_cycle = 1000
 a = np.zeros((num_cycle, 2), dtype=complex)
diff --git a/ref/simple_anim.py b/ref/simple_anim.py
index 59552cd..39c98a0 100644
--- a/ref/simple_anim.py
+++ b/ref/simple_anim.py
@@ -11,22 +11,23 @@ import matplotlib.animation as animation
 
 fig, ax = plt.subplots()
 
-x = np.arange(0, 2*np.pi, 0.01)
+x = np.arange(0, 2 * np.pi, 0.01)
 line, = ax.plot(x, np.sin(x))
 
 
 def init():  # only required for blitting to give a clean slate.
     line.set_ydata([np.nan] * len(x))
-    return line,
+    return (line,)
 
 
 def animate(i):
     line.set_ydata(np.sin(x + i / 100))  # update the data.
-    return line,
+    return (line,)
 
 
 ani = animation.FuncAnimation(
-    fig, animate, init_func=init, interval=2, blit=True, save_count=50)
+    fig, animate, init_func=init, interval=2, blit=True, save_count=50
+)
 
 # To save the animation, use e.g.
 #
diff --git a/ref/strip_chart.py b/ref/strip_chart.py
index 8e9bf2d..11c50fc 100644
--- a/ref/strip_chart.py
+++ b/ref/strip_chart.py
@@ -21,7 +21,7 @@ class Scope(object):
         self.ydata = [0]
         self.line = Line2D(self.tdata, self.ydata)
         self.ax.add_line(self.line)
-        self.ax.set_ylim(-.1, 1.1)
+        self.ax.set_ylim(-0.1, 1.1)
         self.ax.set_xlim(0, self.maxt)
 
     def update(self, y):
@@ -36,18 +36,19 @@ class Scope(object):
         self.tdata.append(t)
         self.ydata.append(y)
         self.line.set_data(self.tdata, self.ydata)
-        return self.line,
+        return (self.line,)
 
 
 def emitter(p=0.03):
-    'return a random value with probability p, else 0'
+    "return a random value with probability p, else 0"
     while True:
         v = np.random.rand(1)
         if v > p:
-            yield 0.
+            yield 0.0
         else:
             yield np.random.rand(1)
 
+
 # Fixing random state for reproducibility
 np.random.seed(19680801)
 
@@ -56,7 +57,6 @@ fig, ax = plt.subplots()
 scope = Scope(ax)
 
 # pass a generator in "emitter" to produce data for the update func
-ani = animation.FuncAnimation(fig, scope.update, emitter, interval=10,
-                              blit=True)
+ani = animation.FuncAnimation(fig, scope.update, emitter, interval=10, blit=True)
 
 plt.show()
diff --git a/split_step_fourier_method.py b/split_step_fourier_method.py
index 4435d86..eedff4b 100644
--- a/split_step_fourier_method.py
+++ b/split_step_fourier_method.py
@@ -12,58 +12,57 @@ import matplotlib.pyplot as plt
 import time
 
 # ---Specify input parameters
-distance = 150.  # Enter fiber length (in units of L_c)=
-# Normalized 2nd-dispersion: kappa=beta2*f^2*L/2): 
+distance = 150.0  # Enter fiber length (in units of L_c)=
+# Normalized 2nd-dispersion: kappa=beta2*f^2*L/2):
 # +ve for normal,-ve for anomalous*)
 kappa = -0.001
-sigma = 0.  # Normalized 3rd-dispersion: sigma=beta3*f^3*L/6
-G = 1.  # small signal gain of Ramam Amp: G=g*L
-Is = 1.  # gain saturation parameter
+sigma = 0.0  # Normalized 3rd-dispersion: sigma=beta3*f^3*L/6
+G = 1.0  # small signal gain of Ramam Amp: G=g*L
+Is = 1.0  # gain saturation parameter
 alpha = 0.4  # Normalized fiber amplitude absorption coeff: alpha=l*L
 
 # Nonlinear parameter n=')
 # sqrt(L_D/L_NL)=sqrt(gamma*P0*T0^2/|beta2|) or QT: n=kappa^0.5
-n = 2. ** 0.5
+n = 2.0 ** 0.5
 
 # ---Specify filter parameters
-bdwidth = 2. * np.pi * 6.
-delta = 2. * np.pi * 0.5
+bdwidth = 2.0 * np.pi * 6.0
+delta = 2.0 * np.pi * 0.5
 a = np.log(np.sqrt(0.95))
 perta = 0.3
 pertfsr = 0.17
 T = 0.2
 t = np.sqrt(T)
-r = 1.j * np.sqrt(1. - T)
+r = 1.0j * np.sqrt(1.0 - T)
 
 # ---Specify input parameters
-mshape = -1.  # m=0 for sech,m>0 for super-Gaussian=
-chirp0 = 0.  # % input pulse chirp (default value)
+mshape = -1.0  # m=0 for sech,m>0 for super-Gaussian=
+chirp0 = 0.0  # % input pulse chirp (default value)
 
 # P = 1/(gamma*L); (P is the ref peak power);
-# uu = A/sqrt(P); 
+# uu = A/sqrt(P);
 # z = z0/L_c; (z0 is the real length, L_c is the cavity length);
 # tau = f*t; (t is the time of reference traveling frame);
 
 
 # ---set simulation parameters
 nt = 2 ** 13  # % FFT points  (powers of 2)
-Tmax = 100.  # (half) window size
+Tmax = 100.0  # (half) window size
 stepno = 1 * round(20 * distance * n ** 2)  # No.of z steps to
 dz = distance / stepno  # step size in z
-dtau = (2. * Tmax) / nt  # step size in tau
+dtau = (2.0 * Tmax) / nt  # step size in tau
 
-Twin = 5.
-fmax = (1. / (2. * Tmax)) * nt / 2.
-fwin = 5.
+Twin = 5.0
+fmax = (1.0 / (2.0 * Tmax)) * nt / 2.0
+fwin = 5.0
 
 filterz = 0.5
 plotz = 5
 
 # ---tau and omega arrays
-tau = np.arange(-nt / 2., nt / 2.) * dtau  # temporal grid
+tau = np.arange(-nt / 2.0, nt / 2.0) * dtau  # temporal grid
 # [(0:nt/2-1) (-nt/2:-1)]
-omega = (np.pi / Tmax) * \
-        np.append(np.arange(0.0, nt / 2.), np.arange(-nt / 2., 0.0))
+omega = (np.pi / Tmax) * np.append(np.arange(0.0, nt / 2.0), np.arange(-nt / 2.0, 0.0))
 
 # frequency grid
 delaytau = dtau * np.arange(-round(Twin / dtau), round(Twin / dtau) + 1)
@@ -71,21 +70,22 @@ delaytau = dtau * np.arange(-round(Twin / dtau), round(Twin / dtau) + 1)
 # Input Field profile
 if mshape == 0:
     # ;% soliton
-    uu = np.exp(-0.5j * chirp0 * tau ** 2.) / np.cosh(tau)
+    uu = np.exp(-0.5j * chirp0 * tau ** 2.0) / np.cosh(tau)
 elif mshape > 0:
     # super-Gaussian
-    uu = np.exp(-0.5 * (1. + 1.j * chirp0) * tau ** (2. * mshape))
+    uu = np.exp(-0.5 * (1.0 + 1.0j * chirp0) * tau ** (2.0 * mshape))
 else:
     # White noise
-    uu = (np.random.randn(nt) + 1.j * np.random.randn(nt)) * np.sqrt(0.5)
+    uu = (np.random.randn(nt) + 1.0j * np.random.randn(nt)) * np.sqrt(0.5)
 
 # temp = np.fft.fftshift(np.fft.ifft(uu) * (nt * dtau) / np.sqrt(2. * np.pi))
 tempomega = np.fft.fftshift(omega)
 
 # ---store dispersive phase shifts to speedup code
 # % nonlinear phase factor
-dispersion = np.exp((-alpha + 1.j * kappa * omega ** 2. +
-                     1.j * sigma * omega ** 3.) * dz)
+dispersion = np.exp(
+    (-alpha + 1.0j * kappa * omega ** 2.0 + 1.0j * sigma * omega ** 3.0) * dz
+)
 
 # comb filter type
 # original comb filter + BPF
@@ -93,23 +93,27 @@ dispersion = np.exp((-alpha + 1.j * kappa * omega ** 2. +
 #           (1 - r ** 2 * np.exp(-1.j * (omega + delta) + a))
 
 # perturbated comb filter + BPF
-filtert = np.exp(-omega ** 2. / bdwidth ** 2. -
-                 perta * np.sin(0.5 * omega / pertfsr) ** 2) * \
-          (t ** 2) / (1.0 - r ** 2 * np.exp(-1.j * (omega + delta) + a))
+filtert = (
+    np.exp(-omega ** 2.0 / bdwidth ** 2.0 - perta * np.sin(0.5 * omega / pertfsr) ** 2)
+    * (t ** 2)
+    / (1.0 - r ** 2 * np.exp(-1.0j * (omega + delta) + a))
+)
 
 fig1 = plt.figure()
-plt.plot(tempomega / (2. * np.pi),
-         np.fft.fftshift(10. * np.log10(np.abs(filtert) ** 2.)))
+plt.plot(
+    tempomega / (2.0 * np.pi), np.fft.fftshift(10.0 * np.log10(np.abs(filtert) ** 2.0))
+)
 plt.title("Perturbated comb filter + BPF")
 plt.xlim(-fwin, fwin)
 plt.xlim(-fwin, fwin)
-plt.ylim(-30., 10.)
+plt.ylim(-30.0, 10.0)
 plt.show()
 
 # %*********[Beginning of MAIN Loop]***********
 # % scheme:1/2N\[Rule]D\[Rule]1/2N;first half step nonlinear
-temp = uu * np.exp((1.j * np.abs(uu) ** 2. +
-                    G / (1.0 + np.abs(uu) ** 2. / Is)) * dz / 2.)
+temp = uu * np.exp(
+    (1.0j * np.abs(uu) ** 2.0 + G / (1.0 + np.abs(uu) ** 2.0 / Is)) * dz / 2.0
+)
 # % note hhz/2
 
 start_time = time.time()
@@ -136,7 +140,9 @@ for i in range(stepno):
         ftemp = np.fft.ifft(temp) * dispersion
 
     uu = np.fft.fft(ftemp)
-    temp = uu * np.exp((1.j * np.abs(uu) ** 2. + G / (1.0 + np.abs(uu) ** 2 / Is)) * dz)
+    temp = uu * np.exp(
+        (1.0j * np.abs(uu) ** 2.0 + G / (1.0 + np.abs(uu) ** 2 / Is)) * dz
+    )
     z = z + dz
 
     if round((z % plotz) / dz) == 0 or round(((z % plotz) - plotz) / dz) == 0:
@@ -144,7 +150,7 @@ for i in range(stepno):
         print("Z: " + str(z))
         # fig2.suptitle("i = " + str(i) + ", z = " + str(z) + " Time: " + str(time_used))
 
-        line1.set_data(tau, np.abs(temp) ** 2.)
+        line1.set_data(tau, np.abs(temp) ** 2.0)
         ax1.relim()
         ax1.autoscale_view(True, True, True)
         ax1.set_xlim([-Twin, Twin])
@@ -152,7 +158,9 @@ for i in range(stepno):
 
         ftemp0 = np.fft.fftshift(ftemp * (nt * dtau) / np.sqrt(2 * np.pi))
 
-        line2.set_data(tempomega / (2. * np.pi), 10. * np.log10(np.abs(ftemp0) ** 2.))
+        line2.set_data(
+            tempomega / (2.0 * np.pi), 10.0 * np.log10(np.abs(ftemp0) ** 2.0)
+        )
         #
         ax2.relim()
         ax2.autoscale_view(True, True, True)
@@ -165,14 +173,20 @@ for i in range(stepno):
         ax3.set_xlim([-Tmax, Tmax])
         ax3.set_title("Time domain (Full Scale)")
 
-        line4.set_data(tempomega / (2. * np.pi), 10. * np.log10(np.abs(ftemp0) ** 2.))
+        line4.set_data(
+            tempomega / (2.0 * np.pi), 10.0 * np.log10(np.abs(ftemp0) ** 2.0)
+        )
         ax4.relim()
         ax4.autoscale_view(True, True, True)
         ax4.set_xlim([-fmax, fmax])
         ax4.set_title("Spectrum (Full Scale)")
 
         autocorr0 = np.fft.fftshift(
-            np.fft.ifft(np.fft.fft(np.abs(temp) ** 2) * np.conjugate(np.fft.fft(np.abs(temp) ** 2))))
+            np.fft.ifft(
+                np.fft.fft(np.abs(temp) ** 2)
+                * np.conjugate(np.fft.fft(np.abs(temp) ** 2))
+            )
+        )
 
         line5.set_data(tau, np.abs(autocorr0) / max(np.abs(autocorr0)))
         ax5.relim()
@@ -187,7 +201,23 @@ plt.show()
 
 # Exporting results
 fname = "result"  # file name
-np.savetxt(fname + "_time.csv", (tau, np.real(uu), np.imag(uu), np.abs(uu) ** 2), delimiter=",")
-np.savetxt(fname + "_freq.csv", (tempomega / (2. * np.pi), np.real(temp), np.imag(temp),
-                                 np.abs(temp) ** 2, np.fft.fftshift(np.abs(filtert) ** 2.)), delimiter=",")
-np.savetxt(fname + "_autocorr.csv", (tau, np.abs(autocorr0), np.abs(autocorr0) / max(np.abs(autocorr0))), delimiter=",")
+np.savetxt(
+    fname + "_time.csv", (tau, np.real(uu), np.imag(uu), np.abs(uu) ** 2), delimiter=","
+)
+np.savetxt(
+    fname + "_freq.csv",
+    (
+        tempomega / (2.0 * np.pi),
+        np.real(temp),
+        np.imag(temp),
+        np.abs(temp) ** 2,
+        np.fft.fftshift(np.abs(filtert) ** 2.0),
+    ),
+    delimiter=",",
+)
+np.savetxt(
+    fname + "_autocorr.csv",
+    (tau, np.abs(autocorr0), np.abs(autocorr0) / max(np.abs(autocorr0))),
+    delimiter=",",
+)
+
diff --git a/split_step_fourier_method_v0.py b/split_step_fourier_method_v0.py
index 7270a7c..f033651 100644
--- a/split_step_fourier_method_v0.py
+++ b/split_step_fourier_method_v0.py
@@ -12,63 +12,62 @@ import matplotlib.pyplot as plt
 import time
 
 # ---Specify input parameters
-distance = 150.  # Enter fiber length (in units of L_c)=
-# Normalized 2nd-dispersion: kappa=beta2*f^2*L/2): 
+distance = 150.0  # Enter fiber length (in units of L_c)=
+# Normalized 2nd-dispersion: kappa=beta2*f^2*L/2):
 # +ve for normal,-ve for anomalous*)
-kappa = -0.001  
-sigma = 0.  # Normalized 3rd-dispersion: sigma=beta3*f^3*L/6
-G = 1.  # small signal gain of Ramam Amp: G=g*L
-Is = 1.  # gain saturation parameter
+kappa = -0.001
+sigma = 0.0  # Normalized 3rd-dispersion: sigma=beta3*f^3*L/6
+G = 1.0  # small signal gain of Ramam Amp: G=g*L
+Is = 1.0  # gain saturation parameter
 alpha = 0.4  # Normalized fiber amplitude absorption coeff: alpha=l*L
 
 
 # Nonlinear parameter n=')
 # sqrt(L_D/L_NL)=sqrt(gamma*P0*T0^2/|beta2|) or QT: n=kappa^0.5
-n = 2. ** 0.5
+n = 2.0 ** 0.5
 
 
 # ---Specify filter parameters
-bdwidth = 2. * np.pi * 6.
-delta = 2. * np.pi * 0.5
+bdwidth = 2.0 * np.pi * 6.0
+delta = 2.0 * np.pi * 0.5
 a = np.log(np.sqrt(0.95))
 perta = 0.3
 pertfsr = 0.17
 T = 0.2
 t = np.sqrt(T)
-r = 1.j * np.sqrt(1. - T)
+r = 1.0j * np.sqrt(1.0 - T)
 
 
 # ---Specify input parameters
-mshape = -1.  # m=0 for sech,m>0 for super-Gaussian=
-chirp0 = 0.  # % input pulse chirp (default value)
+mshape = -1.0  # m=0 for sech,m>0 for super-Gaussian=
+chirp0 = 0.0  # % input pulse chirp (default value)
 
 
-# P = 1/(gamma*L); (P is the ref peak power);  
-# uu = A/sqrt(P); 
+# P = 1/(gamma*L); (P is the ref peak power);
+# uu = A/sqrt(P);
 # z = z0/L_c; (z0 is the real length, L_c is the cavity length);
 # tau = f*t; (t is the time of reference traveling frame);
 
 
 # ---set simulation parameters
 nt = 2 ** 13  # % FFT points  (powers of 2)
-Tmax = 100.  # (half) window size
+Tmax = 100.0  # (half) window size
 stepno = 1 * round(20 * distance * n ** 2)  # No.of z steps to
 dz = distance / stepno  # step size in z
-dtau = (2. * Tmax) / nt  # step size in tau
+dtau = (2.0 * Tmax) / nt  # step size in tau
 
 
-Twin = 5.
-fmax = (1. / (2. * Tmax)) * nt / 2.
-fwin = 5.
+Twin = 5.0
+fmax = (1.0 / (2.0 * Tmax)) * nt / 2.0
+fwin = 5.0
 
 filterz = 0.5
 plotz = 5
 
 # ---tau and omega arrays
-tau = np.arange(-nt / 2., nt / 2.) * dtau  # temporal grid
+tau = np.arange(-nt / 2.0, nt / 2.0) * dtau  # temporal grid
 # [(0:nt/2-1) (-nt/2:-1)]
-omega = (np.pi / Tmax) * \
-        np.append(np.arange(0.0 , nt / 2.), np.arange(-nt / 2., 0.0))
+omega = (np.pi / Tmax) * np.append(np.arange(0.0, nt / 2.0), np.arange(-nt / 2.0, 0.0))
 
 
 # frequency grid
@@ -78,23 +77,24 @@ delaytau = dtau * np.arange(-round(Twin / dtau), round(Twin / dtau) + 1)
 # Input Field profile
 if mshape == 0:
     # ;% soliton
-    uu = np.exp(-0.5j * chirp0 * tau ** 2.) / np.cosh(tau)
+    uu = np.exp(-0.5j * chirp0 * tau ** 2.0) / np.cosh(tau)
 elif mshape > 0:
     # super-Gaussian
-    uu = np.exp(-0.5 * (1. + 1.j * chirp0) * tau ** (2. * mshape))
+    uu = np.exp(-0.5 * (1.0 + 1.0j * chirp0) * tau ** (2.0 * mshape))
 else:
     # White noise
-    uu = (np.random.randn(nt) + 1.j * np.random.randn(nt)) * np.sqrt(0.5)
+    uu = (np.random.randn(nt) + 1.0j * np.random.randn(nt)) * np.sqrt(0.5)
 
 
-temp = np.fft.fftshift(np.fft.ifft(uu) * (nt * dtau) / np.sqrt(2.* np.pi))
+temp = np.fft.fftshift(np.fft.ifft(uu) * (nt * dtau) / np.sqrt(2.0 * np.pi))
 tempomega = np.fft.fftshift(omega)
 
 
 # ---store dispersive phase shifts to speedup code
 # % nonlinear phase factor
-dispersion = np.exp((-alpha + 1.j * kappa * omega **2. + \
-                     1.j * sigma * omega ** 3.) * dz)
+dispersion = np.exp(
+    (-alpha + 1.0j * kappa * omega ** 2.0 + 1.0j * sigma * omega ** 3.0) * dz
+)
 
 # comb filter type
 # original comb filter + BPF
@@ -102,22 +102,26 @@ dispersion = np.exp((-alpha + 1.j * kappa * omega **2. + \
 #           (1 - r ** 2 * np.exp(-1.j * (omega + delta) + a))
 
 # perturbated comb filter + BPF
-filtert = np.exp(-omega ** 2. / bdwidth ** 2. - \
-                 perta * np.sin(0.5 * omega / pertfsr) ** 2) * \
-          (t ** 2) / (1.0 - r ** 2 * np.exp(-1.j*(omega + delta) + a))
+filtert = (
+    np.exp(-omega ** 2.0 / bdwidth ** 2.0 - perta * np.sin(0.5 * omega / pertfsr) ** 2)
+    * (t ** 2)
+    / (1.0 - r ** 2 * np.exp(-1.0j * (omega + delta) + a))
+)
 
 plt.figure()
-plt.plot(tempomega / (2. * np.pi),
-         np.fft.fftshift(10. * np.log10(np.abs(filtert) ** 2.)))
+plt.plot(
+    tempomega / (2.0 * np.pi), np.fft.fftshift(10.0 * np.log10(np.abs(filtert) ** 2.0))
+)
 plt.xlim(-fwin, fwin)
 plt.xlim(-fwin, fwin)
-plt.ylim(-30., 10.)
+plt.ylim(-30.0, 10.0)
 plt.show()
 
 # %*********[Beginning of MAIN Loop]***********
 # % scheme:1/2N\[Rule]D\[Rule]1/2N;first half step nonlinear
-temp = uu * np.exp((1.j * np.abs(uu) ** 2. + 
-                   G / (1.0 + np.abs(uu) ** 2. / Is )) * dz / 2.)
+temp = uu * np.exp(
+    (1.0j * np.abs(uu) ** 2.0 + G / (1.0 + np.abs(uu) ** 2.0 / Is)) * dz / 2.0
+)
 # % note hhz/2
 
 start_time = time.time()
@@ -131,13 +135,15 @@ z = 0
 plt.ion()
 fig, ((ax1, ax2), (ax3, ax4), (ax5, ax6)) = plt.subplots(3, 2)
 for i in range(stepno):
-    if round((z % 1 - filterz)/dz) == 0:
+    if round((z % 1 - filterz) / dz) == 0:
         ftemp = np.fft.ifft(temp) * filtert * dispersion
     else:
         ftemp = np.fft.ifft(temp) * dispersion
 
     uu = np.fft.fft(ftemp)
-    temp = uu*np.exp((1.j * np.abs(uu)**2. + G / (1.0 + np.abs(uu)**2 / Is)) * dz)
+    temp = uu * np.exp(
+        (1.0j * np.abs(uu) ** 2.0 + G / (1.0 + np.abs(uu) ** 2 / Is)) * dz
+    )
     z = z + dz
 
     if round((z % plotz) / dz) == 0 or round(((z % plotz) - plotz) / dz) == 0:
@@ -145,29 +151,33 @@ for i in range(stepno):
         print(i)
 
         ax1.clear()
-        ax1.plot(tau, np.abs(temp)**2.)
+        ax1.plot(tau, np.abs(temp) ** 2.0)
         ax1.set_xlim([-Twin, Twin])
         ax1.set_title("i = " + str(i) + ", z = " + str(z))
 
         ftemp0 = np.fft.fftshift(ftemp * (nt * dtau) / np.sqrt(2 * np.pi))
 
         ax2.clear()
-        ax2.plot(tempomega / (2. * np.pi), 10. * np.log10(np.abs(ftemp0)**2.))
+        ax2.plot(tempomega / (2.0 * np.pi), 10.0 * np.log10(np.abs(ftemp0) ** 2.0))
         ax2.set_xlim([-fwin, fwin])
         ax2.set_title(" Time = " + str(time_used))
 
         ax3.clear()
-        ax3.plot(tau, np.abs(temp)**2)
+        ax3.plot(tau, np.abs(temp) ** 2)
         ax3.set_xlim([-Tmax, Tmax])
         ax3.set_title("i = " + str(i) + ", z = " + str(z))
 
         ax4.clear()
-        ax4.plot(tempomega / (2. * np.pi), 10. * np.log10(np.abs(ftemp0)**2.))
+        ax4.plot(tempomega / (2.0 * np.pi), 10.0 * np.log10(np.abs(ftemp0) ** 2.0))
         ax4.set_xlim([-fmax, fmax])
         ax4.set_title(" Time = " + str(time_used))
 
         autocorr0 = np.fft.fftshift(
-            np.fft.ifft(np.fft.fft(np.abs(temp)**2) * np.conjugate(np.fft.fft(np.abs(temp)**2))))
+            np.fft.ifft(
+                np.fft.fft(np.abs(temp) ** 2)
+                * np.conjugate(np.fft.fft(np.abs(temp) ** 2))
+            )
+        )
 
         ax5.clear()
         ax5.plot(tau, autocorr0 / max(autocorr0))
diff --git a/tmp.py b/tmp.py
index d480db8..27132d3 100644
--- a/tmp.py
+++ b/tmp.py
@@ -1 +1,2 @@
-print(1.3%1.0)
\ No newline at end of file
+print(1.3 % 1.0)
+