diff --git a/tracking/central_electrons.sh b/tracking/central_electrons.sh
index c3ad2cc99a77343b8bf4035d45b6e305474d9c1f..3b8366c84b018e7e1b2d6d98494a7873eaa576cf 100644
--- a/tracking/central_electrons.sh
+++ b/tracking/central_electrons.sh
@@ -71,7 +71,7 @@ mkdir -p results/tracking
 
 root -b -q "tracking/scripts/rec_central_electrons.cxx(\"${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}\")"
 if [[ "$?" -ne "0" ]] ; then
-  echo "ERROR running juggler"
+  echo "ERROR running root script"
   exit 1
 fi
 
diff --git a/tracking/scripts/rec_central_electrons.cxx b/tracking/scripts/rec_central_electrons.cxx
index fdcb0d6b18913bab627446a3875e00ba3c724228..d58cd1ba60337b3b12fc88477deea3c76d06f721 100644
--- a/tracking/scripts/rec_central_electrons.cxx
+++ b/tracking/scripts/rec_central_electrons.cxx
@@ -1,6 +1,12 @@
 #include "ROOT/RDataFrame.hxx"
+#include "TCanvas.h"
+#include "TLegend.h"
+#include "TH1D.h"
+
 #include <iostream>
 
+R__LOAD_LIBRARY(libeicd.so)
+R__LOAD_LIBRARY(libDD4pod.so)
 #include "dd4pod/Geant4ParticleCollection.h"
 #include "eicd/TrackParametersCollection.h"
 #include "eicd/ClusterCollection.h"
@@ -60,6 +66,16 @@ auto delta_p = [](const std::vector<double>& tracks, const std::vector<double>&
   return res;
 };
 
+auto delta_p_over_p = [](const std::vector<double>& tracks, const std::vector<double>& thrown) {
+  std::vector<double> res;
+  for (const auto& p1 : thrown) {
+    for (const auto& p2 : tracks) {
+      res.push_back((p1 - p2)/p1);
+    }
+  }
+  return res;
+};
+
 int rec_central_electrons(const char* fname = "topside/rec_central_electrons.root")
 {
 
@@ -74,16 +90,23 @@ int rec_central_electrons(const char* fname = "topside/rec_central_electrons.roo
                  .Define("p_track", p_track, {"outputTrackParameters"})
                  .Define("p_track1", p_track, {"outputTrackParameters1"})
                  .Define("p_track2", p_track, {"outputTrackParameters2"})
-                 .Define("delta_p",delta_p, {"p_track", "p_thrown"})
+                 .Define("delta_p0",delta_p, {"p_track", "p_thrown"})
                  .Define("delta_p1",delta_p, {"p_track1", "p_thrown"})
-                 .Define("delta_p2",delta_p, {"p_track2", "p_thrown"});
+                 .Define("delta_p2",delta_p, {"p_track2", "p_thrown"})
+                 .Define("delta_p_over_p0",delta_p_over_p, {"p_track", "p_thrown"})
+                 .Define("delta_p_over_p1",delta_p_over_p, {"p_track1", "p_thrown"})
+                 .Define("delta_p_over_p2",delta_p_over_p, {"p_track2", "p_thrown"});
 
   auto h_nTracks = df0.Histo1D({"h_nTracks", "; N tracks ", 10, 0, 10}, "nTracks");
   auto h_pTracks = df0.Histo1D({"h_pTracks", "; GeV/c ", 100, 0, 10}, "p_track");
 
-  auto h_delta_p  = df0.Histo1D({"h_delta_p", "; GeV/c ",  100, -10, 10}, "delta_p");
-  auto h_delta_p1 = df0.Histo1D({"h_delta_p1", "; GeV/c ", 100, -10, 10}, "delta_p1");
-  auto h_delta_p2 = df0.Histo1D({"h_delta_p2", "; GeV/c ", 100, -10, 10}, "delta_p2");
+  auto h_delta_p0  = df0.Histo1D({"h_delta_p0", "Truth Track Init; GeV/c ",  100, -10, 10}, "delta_p0");
+  auto h_delta_p1 = df0.Histo1D({"h_delta_p1", "Ecal Cluster Init; GeV/c ", 100, -10, 10}, "delta_p1");
+  auto h_delta_p2 = df0.Histo1D({"h_delta_p2", "Ecal Cluster, innner vtx hit Init; GeV/c ", 100, -10, 10}, "delta_p2");
+
+  auto h_delta_p0_over_p = df0.Histo1D({"h_delta_p0_over_p",  "Truth Track Init; delta p/p ",  100, -0.1, 0.1}, "delta_p_over_p0");
+  auto h_delta_p1_over_p = df0.Histo1D({"h_delta_p1_over_p", "Ecal Cluster Init; delta p/p ", 100, -0.1, 0.1}, "delta_p_over_p1");
+  auto h_delta_p2_over_p = df0.Histo1D({"h_delta_p2_over_p", "Ecal Cluster, innner vtx hit Init; delta p/p ", 100, -0.1, 0.1}, "delta_p_over_p2");
 
   auto c = new TCanvas();
 
@@ -95,8 +118,9 @@ int rec_central_electrons(const char* fname = "topside/rec_central_electrons.roo
   c->SaveAs("results/tracking/rec_central_electrons_pTracks.png");
   c->SaveAs("results/tracking/rec_central_electrons_pTracks.pdf");
 
-  THStack * hs = new THStack("hs_delta_p","; GeV/c "); 
-  TH1D* h1 = (TH1D*) h_delta_p->Clone();
+  c = new TCanvas();
+  THStack * hs = new THStack("hs_delta_p","; GeV/c ");
+  TH1D* h1 = (TH1D*) h_delta_p0->Clone();
   hs->Add(h1);
   h1 = (TH1D*) h_delta_p1->Clone();
   h1->SetLineColor(2);
@@ -107,8 +131,26 @@ int rec_central_electrons(const char* fname = "topside/rec_central_electrons.roo
   h1->SetFillColor(4);
   hs->Add(h1);
   hs->Draw("nostack");
+  c->BuildLegend();
   c->SaveAs("results/tracking/rec_central_electrons_delta_p.png");
   c->SaveAs("results/tracking/rec_central_electrons_delta_p.pdf");
 
+  c  = new TCanvas();
+  hs = new THStack("hs_delta_p_over_p","; delta p/p ");
+  h1 = (TH1D*) h_delta_p0_over_p->Clone();
+  hs->Add(h1);
+  h1 = (TH1D*) h_delta_p1_over_p->Clone();
+  h1->SetLineColor(2);
+  hs->Add(h1);
+  h1 = (TH1D*) h_delta_p2_over_p->Clone();
+  h1->SetLineColor(4);
+  h1->SetFillStyle(3001);
+  h1->SetFillColor(4);
+  hs->Add(h1);
+  hs->Draw("nostack");
+  c->BuildLegend();
+  c->SaveAs("results/tracking/rec_central_electrons_delta_p_over_p.png");
+  c->SaveAs("results/tracking/rec_central_electrons_delta_p_over_p.pdf");
+
   return 0;
 }