diff --git a/projects/braidingTightBinding/bin/plot.py b/projects/braidingTightBinding/bin/plot.py
new file mode 100644
index 0000000000000000000000000000000000000000..fff45edbbd01dbcb59266045baf756d375dff456
--- /dev/null
+++ b/projects/braidingTightBinding/bin/plot.py
@@ -0,0 +1,67 @@
+import matplotlib.pyplot as plt
+import numpy as np
+import sys
+
+
+if len(sys.argv) < 4:
+    print(f"Usage: python3.9 {argv[0]} ../path/to/data/ ../path/to/plots/ name [modes_range_start modes_range_end]")
+    sys.exit()
+
+readpath = sys.argv[1]
+writepath = sys.argv[2]
+name = sys.argv[3]
+
+#eigenvalues
+ev = np.loadtxt(readpath+"eigenvalues.txt")
+
+starting_mode = 1110
+end_mode = 1180
+
+if len(sys.argv) == 6:
+    starting_mode = int(sys.argv[4])
+    end_mode = int(sys.argv[5])
+
+if end_mode - starting_mode > 70:
+    print(f"Maximum number of modes: 70 (entered: {end_mode - starting_mode})")
+    end_mode = starting_mode + 70
+    print(f"Using Modes {starting_mode}-{end_mode}")
+
+plt.plot(np.array(list(range(starting_mode,end_mode))),ev[starting_mode:end_mode],'o', ms=1)
+plt.title(f"Frequency levels under {name}")
+plt.xlabel("Level index")
+plt.ylabel("Frequency [Hz]")
+plt.savefig(writepath + "Frequencies_zoomed.pdf", dpi = 500)
+
+plt.clf()
+plt.plot(ev[:],'o', ms=1)
+plt.title(f"Frequency levels under {name}")
+plt.xlabel("Level index")
+plt.ylabel("Frequency [Hz]")
+plt.savefig(writepath + "Frequencies.pdf", dpi = 500)
+
+
+
+
+#eigenvectors
+data = np.loadtxt(readpath+"eigenvectors.txt")
+pos_data = np.loadtxt(readpath+"PARAMS.txt", dtype=str) #positions are rows 4 and 5
+desired_length = len(data[0,:])
+pos_x = pos_data[0:desired_length,4].astype('float')
+pos_y = pos_data[:desired_length,5].astype('float')
+
+plt.figure(figsize=(8,5*(end_mode - starting_mode + 1)))
+for i in range(starting_mode,end_mode):
+    plt.subplot(int(end_mode-starting_mode),1,i-starting_mode+1)
+    cmax = np.amax(data[i,:])
+    cmin = np.amin(data[i,:])
+    if(-cmin > cmax):
+        cmax = -cmin
+    elif(-cmin <= cmax):
+        cmin = -cmax
+    plt.scatter(pos_x, pos_y, c=data[i,:], edgecolors='k', linewidth=0.1, cmap='seismic', s=16, vmin=cmin, vmax=cmax)
+    plt.colorbar()
+    plt.title(f"Mode {i+1}, ({ev[i]} Hz)")
+    plt.xlabel("x")
+    plt.ylabel("y")
+
+plt.savefig(writepath + f"Modes.pdf", dpi = 40, bbox_inches='tight', pad_inches=0)