diff --git a/ThrustCurve.png b/ThrustCurve.png
index f004877..6b336c4 100644
Binary files a/ThrustCurve.png and b/ThrustCurve.png differ
diff --git a/thrust.jl b/thrust.jl
index 7668d59..06ca7b0 100644
--- a/thrust.jl
+++ b/thrust.jl
@@ -1,9 +1,11 @@
 using Unitful
 using DataFrames
 using Plots
+theme(:ggplot2); 
 using UnitfulRecipes
 using CSV
 using Measurements
+using Measurements: value, uncertainty
 using Printf
 
 
@@ -30,7 +32,7 @@ t = 0.0u"s"
 P = P0 |> u"Pa"
 M = V * (P / (R * T)) |> u"kg"
 ts = 1u"ms"
-df = DataFrame(Thrust = (0 ± 0)u"N", Pressure = P, Time = t, Mass = M)
+df = DataFrame(Thrust=(0 ± 0)u"N", Pressure=P, Time=t, Mass=M)
 
 while M > 0.005u"kg"
     # while t < 30u"s"
@@ -47,7 +49,7 @@ while M > 0.005u"kg"
     P = (M * R * T) / V |> u"Pa"
     t = t + ts
 
-    df_step = DataFrame(Thrust = Thrust, Pressure = P, Time = t, Mass = M)
+    df_step = DataFrame(Thrust=Thrust, Pressure=P, Time=t, Mass=M)
     append!(df, df_step)
 end
 final_time = t |> u"s"
@@ -61,7 +63,17 @@ println("Mass Total: ", Wtank + Wsolenoid + maximum(df.Mass))
 
 print(describe(df))
 
-plot(df.Time[1:350:end], df.Thrust[1:350:end], title = "Thrust Over Time")
+thrust_values = df.Thrust .|> ustrip .|> value
+thrust_uncertainties = df.Thrust .|> ustrip .|> uncertainty 
+
+plot(df.Time .|> ustrip, thrust_values, 
+    title="Thrust Over Time", 
+    ribbon=(thrust_uncertainties, thrust_uncertainties), 
+    fillalpha=.2,label="Thrust",
+    xlabel="Time (s)", 
+    ylabel="Thrust (N)",
+    palette=:leonardo,
+    )
 
 ### Save data to readme.md