diff --git a/src/View/Results/PlotSedProfile.py b/src/View/Results/PlotSedProfile.py
index 8f77ee3b420fcffa18ac4fd009a300570c95e593..62f4e29059972bd35b0b70bf2dea9cccb433a12d 100644
--- a/src/View/Results/PlotSedProfile.py
+++ b/src/View/Results/PlotSedProfile.py
@@ -31,28 +31,7 @@ class PlotSedProfile(APlot):
         self._current_reach_id = reach_id
         self._current_profile_id = profile_id
 
-    @timer
-    def draw(self):
-        self.canvas.axes.cla()
-        self.canvas.axes.grid(color='grey', linestyle='--', linewidth=0.5)
-
-        if self.data is None:
-            return
-
-        reach = self._results.river.reach(self._current_reach_id)
-        profile = reach.profile(self._current_profile_id)
-        if profile.geometry.number_points == 0:
-            return
-
-        self.canvas.axes.set_xlabel(
-            _translate("MainWindow_reach", "X (m)"),
-            color='green', fontsize=12
-        )
-        self.canvas.axes.set_ylabel(
-            _translate("MainWindow_reach", "Height (m)"),
-            color='green', fontsize=12
-        )
-
+    def get_zsl(self, profile):
         x = profile.geometry.get_station()
         z = profile.geometry.z()
 
@@ -62,56 +41,101 @@ class PlotSedProfile(APlot):
                 profile.get_ts_key(self._current_timestamp, "sl")
             )
         )
+        profiles_sl_0 = list(
+            map(
+                lambda sl: sl[0],
+                profile.get_ts_key(0.0, "sl")
+            )
+        )
+
+        f = list(map(lambda p: 0, range(profile.geometry.number_points)))
 
         sl = []
-        for i in range(len(profiles_sl)):
+        sl_0 = []
+        for profile_sl, profile_sl_0 in zip(profiles_sl, profiles_sl_0):
             cur = []
-            for p in range(profile.geometry.number_points):
-                    cur.append(profiles_sl[i])
-            sl.append(cur)
+            cur_0 = []
 
-        # logger.info(sl)
-
-        self.canvas.axes.set_xlim(
-            left = min(x), right = max(x)
-        )
+            for p in range(profile.geometry.number_points):
+                 cur.append(profile_sl)
+                 cur_0.append(profile_sl_0)
 
-        # Dummy layer with height = 0
-        f = list(map(lambda p: 0, range(profile.geometry.number_points)))
+            sl.append(cur)
+            sl_0.append(cur_0)
 
-        # We compute Z sediment layer in reverse order, from last layer to
-        # fake river bottom
-        r_sl = list(reversed(sl))
+        # Compute sediment layer from initial data in function to
+        # profile z_min
         z_sl = reduce(
             lambda acc, v: acc + [
                 list(
-                    map(lambda x, y: y + x, v, acc[-1])
+                    map(lambda x, y: y - x, v, acc[-1])
                 )
             ],
-            r_sl,
-            [f]
+            sl_0,
+            [z]
         )
 
-        # We normalize Z coordinate to 0 (the maximum must be 0)
-        f_z_max = max(z_sl[-1])
-        z_sl = list(
+        # Diff between initial data and data att current timestamp
+        d_sl = list(
             map(
-                lambda p: list(map(lambda z: z - f_z_max, p)),
-                z_sl
+                lambda ln0, lni: list(
+                    map(
+                        lambda z0, zi: z0 - zi,
+                        ln0, lni
+                    )
+                ),
+                sl_0, sl
             )
         )
 
-        # We apply the river geometry bottom height at each layers to
-        # fond the new river geometry
+        # Apply diff for t0 for each layer Z
         z_sl = list(
             map(
-                lambda sl: list(
-                    map(lambda pz, m: pz + m, sl, z)
+                lambda z, d: list(
+                    map(
+                        lambda zn, dn: zn - dn,
+                        z, d
+                    )
                 ),
-                z_sl
+                z_sl,
+                d_sl + [f]      # HACK: Add dummy data for last layer
             )
         )
 
+        return list(reversed(z_sl))
+
+
+    @timer
+    def draw(self):
+        self.canvas.axes.cla()
+        self.canvas.axes.grid(color='grey', linestyle='--', linewidth=0.5)
+
+        if self.data is None:
+            return
+
+        reach = self._results.river.reach(self._current_reach_id)
+        profile = reach.profile(self._current_profile_id)
+        if profile.geometry.number_points == 0:
+            return
+
+        self.canvas.axes.set_xlabel(
+            _translate("MainWindow_reach", "X (m)"),
+            color='green', fontsize=12
+        )
+        self.canvas.axes.set_ylabel(
+            _translate("MainWindow_reach", "Height (m)"),
+            color='green', fontsize=12
+        )
+
+        x = profile.geometry.get_station()
+        z = profile.geometry.z()
+
+        self.canvas.axes.set_xlim(
+            left = min(x), right = max(x)
+        )
+
+        z_sl = self.get_zsl(profile)
+
         self.line_kp_sl = []
         for i, zsl in enumerate(z_sl):
             self.line_kp_sl.append(None)
diff --git a/src/View/Results/PlotSedReach.py b/src/View/Results/PlotSedReach.py
index d88e26df6b18fcc8c6b255ef04e894a763af71b7..657b70e855a8c10c393e028428fa45c33132a1b8 100644
--- a/src/View/Results/PlotSedReach.py
+++ b/src/View/Results/PlotSedReach.py
@@ -31,34 +31,93 @@ class PlotSedReach(APlot):
         self._current_reach_id = reach_id
         self._current_profile_id = profile_id
 
-    @timer
-    def draw(self):
-        self.canvas.axes.cla()
-        self.canvas.axes.grid(color='grey', linestyle='--', linewidth=0.5)
-
-        if self.data is None:
-            return
-
-        reach = self._results.river.reach(self._current_reach_id)
-        if reach.geometry.number_profiles == 0:
-            return
-
-        self.canvas.axes.set_xlabel(
-            _translate("MainWindow_reach", "Kp (m)"),
-            color='green', fontsize=12
-        )
-        self.canvas.axes.set_ylabel(
-            _translate("MainWindow_reach", "Height (m)"),
-            color='green', fontsize=12
-        )
-
+    # DEPRECATED version of sediment layser display
+    # def _get_zsl(self, reach):
+    #     kp = reach.geometry.get_kp()
+    #     z_min = reach.geometry.get_z_min()
+    #     z_max = reach.geometry.get_z_max()
+
+    #     profiles_sl = list(
+    #         map(
+    #             # Get SL list for profile p
+    #             lambda p: p.get_ts_key(self._current_timestamp, "sl"),
+    #             reach.profiles
+    #         )
+    #     )
+
+    #     max_sl_num = reduce(
+    #         lambda acc, sl: max(acc, len(sl)),
+    #         profiles_sl,
+    #         0
+    #     )
+
+    #     sl = []
+    #     for i in range(max_sl_num):
+    #         cur = []
+    #         for profile_sl in profiles_sl:
+    #             if i < len(profile_sl):
+    #                 cur.append(profile_sl[i][0])
+    #             else:
+    #                 cur.append(0)
+    #         sl.append(cur)
+
+    #     self.canvas.axes.set_xlim(
+    #         left = min(kp) - 10, right = max(kp) + 10
+    #     )
+
+    #     # Dummy layer with height = 0
+    #     f = list(map(lambda p: 0, reach.profiles))
+
+    #     # We compute Z sediment layer in reverse order, from last layer to
+    #     # fake river bottom
+    #     r_sl = list(reversed(sl))
+    #     z_sl = reduce(
+    #         lambda acc, v: acc + [
+    #             list(
+    #                 map(lambda x, y: y + x, v, acc[-1])
+    #             )
+    #         ],
+    #         r_sl,
+    #         [f]
+    #     )
+
+    #     # We normalize Z coordinate to 0 (the maximum must be 0)
+    #     f_z_max = max(z_sl[-1])
+    #     z_sl = list(
+    #         map(
+    #             lambda p: list(map(lambda z: z - f_z_max, p)),
+    #             z_sl
+    #         )
+    #     )
+
+    #     # We apply the river geometry bottom height at each layers to
+    #     # fond the new river geometry
+    #     z_sl = list(
+    #         map(
+    #             lambda sl: list(
+    #                 map(lambda z, m: z + m, sl, z_min)
+    #             ),
+    #             z_sl
+    #         )
+    #     )
+
+    #     return z_sl
+
+    def get_zsl(self, reach):
         kp = reach.geometry.get_kp()
         z_min = reach.geometry.get_z_min()
         z_max = reach.geometry.get_z_max()
 
+        profiles_sl_0 = list(
+            map(
+                # Get SL list for profile p at time 0 (initial data)
+                lambda p: p.get_ts_key(0.0, "sl"),
+                reach.profiles
+            )
+        )
         profiles_sl = list(
             map(
-                # Get SL list for profile p
+                # Get SL list for profile p at current time
                 lambda p: p.get_ts_key(self._current_timestamp, "sl"),
                 reach.profiles
             )
@@ -69,57 +128,91 @@ class PlotSedReach(APlot):
             profiles_sl,
             0
         )
+        f = list(map(lambda p: 0, reach.profiles))
 
         sl = []
+        sl_0 = []
         for i in range(max_sl_num):
             cur = []
-            for profile_sl in profiles_sl:
-                if i < len(profile_sl):
+            cur_0 = []
+            for profile_sl, profile_sl_0 in zip(profiles_sl, profiles_sl_0):
+                if i < len(profile_sl_0):
                     cur.append(profile_sl[i][0])
+                    cur_0.append(profile_sl_0[i][0])
                 else:
                     cur.append(0)
+                    cur_0.append(0)
             sl.append(cur)
+            sl_0.append(cur_0)
 
-        self.canvas.axes.set_xlim(
-            left = min(kp) - 10, right = max(kp) + 10
-        )
-
-        # Dummy layer with height = 0
-        f = list(map(lambda p: 0, reach.profiles))
-
-        # We compute Z sediment layer in reverse order, from last layer to
-        # fake river bottom
-        r_sl = list(reversed(sl))
+        # Compute sediment layer from initial data in function to
+        # profile z_min
         z_sl = reduce(
             lambda acc, v: acc + [
                 list(
-                    map(lambda x, y: y + x, v, acc[-1])
+                    map(lambda x, y: y - x, v, acc[-1])
                 )
             ],
-            r_sl,
-            [f]
+            sl_0,
+            [z_min]
         )
 
-        # We normalize Z coordinate to 0 (the maximum must be 0)
-        f_z_max = max(z_sl[-1])
-        z_sl = list(
+        # Diff between initial data and data att current timestamp
+        d_sl = list(
             map(
-                lambda p: list(map(lambda z: z - f_z_max, p)),
-                z_sl
+                lambda ln0, lni: list(
+                    map(
+                        lambda z0, zi: z0 - zi,
+                        ln0, lni
+                    )
+                ),
+                sl_0, sl
             )
         )
 
-        # We apply the river geometry bottom height at each layers to
-        # fond the new river geometry
+        # Apply diff for t0 for each layer Z
         z_sl = list(
             map(
-                lambda sl: list(
-                    map(lambda z, m: z + m, sl, z_min)
+                lambda z, d: list(
+                    map(
+                        lambda zn, dn: zn - dn,
+                        z, d
+                    )
                 ),
-                z_sl
+                z_sl,
+                d_sl + [f]      # HACK: Add dummy data for last layer
             )
         )
 
+        return list(reversed(z_sl))
+
+    @timer
+    def draw(self):
+        self.canvas.axes.cla()
+        self.canvas.axes.grid(color='grey', linestyle='--', linewidth=0.5)
+
+        if self.data is None:
+            return
+
+        reach = self._results.river.reach(self._current_reach_id)
+        if reach.geometry.number_profiles == 0:
+            return
+
+        self.canvas.axes.set_xlabel(
+            _translate("MainWindow_reach", "Kp (m)"),
+            color='green', fontsize=12
+        )
+        self.canvas.axes.set_ylabel(
+            _translate("MainWindow_reach", "Height (m)"),
+            color='green', fontsize=12
+        )
+
+        kp = reach.geometry.get_kp()
+        z_min = reach.geometry.get_z_min()
+        z_max = reach.geometry.get_z_max()
+
+        z_sl = self.get_zsl(reach)
+
         # Draw
         self.line_kp_sl = []
         for i, z in enumerate(z_sl):