diff --git a/full-requirements.txt b/full-requirements.txt
index bfef1edbce296deca2b9e5695b2db39da0ad9b6f..124d97f83832b091d4dab9ddea7ac45f8a939cc9 100644
--- a/full-requirements.txt
+++ b/full-requirements.txt
@@ -11,3 +11,4 @@ pyinstaller>=5.11.0
 shapely>=2.0.1
 lxml>=4.9.3
 platformdirs>=4.2.0
+pyshp>=2.3.1
diff --git a/requirements.txt b/requirements.txt
index 1f0a31942f0f27f41b0c3f1602751ecc13196341..a439384228ffb7c508423e1e552d94a02182d4eb 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -11,3 +11,4 @@ pyinstaller>=5.11.0
 shapely>=2.0.1
 lxml>=4.9.3
 platformdirs>=4.2.0
+pyshp>=2.3.1
diff --git a/src/Model/Geometry/PointXYZ.py b/src/Model/Geometry/PointXYZ.py
index 2e806687cbc5bb9f5bfeac1031e06cb02c93ac04..cf678c24f26e6f25e9ff728406971c94e405c2e9 100644
--- a/src/Model/Geometry/PointXYZ.py
+++ b/src/Model/Geometry/PointXYZ.py
@@ -192,6 +192,22 @@ class PointXYZ(Point, SQLSubModel):
     def dist(self, p2):
         return PointXYZ.distance(self, p2)
 
+    def dist_2d(self, p2):
+        return PointXYZ.distance_2d(self, p2)
+
+    @staticmethod
+    def distance_2d(p1, p2):
+        """Euclidean distance between p1 and p2.
+
+        Args:
+            p1: A XYZ Point
+            p2: A XYZ Point
+
+        Returns:
+            Euclidean 2D distance between the two points
+        """
+        return dist((p1.x, p1.y), (p2.x, p2.y))
+
     @staticmethod
     def distance(p1, p2):
         """Euclidean distance between p1 and p2.
@@ -201,6 +217,6 @@ class PointXYZ(Point, SQLSubModel):
             p2: A XYZ Point
 
         Returns:
-            Euclidean distance between the two points
+            Euclidean 3D distance between the two points
         """
         return dist((p1.x, p1.y, p1.z), (p2.x, p2.y, p2.z))
diff --git a/src/Model/Geometry/Reach.py b/src/Model/Geometry/Reach.py
index 448bc6220f8c761b67bffb339a21b098eda3234a..e294a087d630429cf0b73fcf6a3af7928cc9c4c7 100644
--- a/src/Model/Geometry/Reach.py
+++ b/src/Model/Geometry/Reach.py
@@ -16,7 +16,9 @@
 
 # -*- coding: utf-8 -*-
 
+import os
 import logging
+import shapefile
 import numpy as np
 
 from time import time
@@ -533,8 +535,58 @@ class Reach(SQLSubModel):
 
     # Import/Export
 
-    @timer
     def import_geometry(self, file_path_name: str):
+        if file_path_name.endswith(".st"):
+            return self.import_geometry_st(file_path_name)
+        elif file_path_name.endswith(".shp"):
+            return self.import_geometry_shp(file_path_name)
+        else:
+            return []
+
+    @timer
+    def import_geometry_shp(self, file_path_name: str):
+        sf = shapefile.Reader(
+            os.path.abspath(file_path_name)
+        )
+
+        if sf.shapeType != shapefile.POLYLINEZ:
+            raise FileFormatError(
+                file_path_name,
+                "Shapefile is not a POLYLINEZ"
+            )
+
+        profiles = []
+
+        for i in range(len(sf)):
+            s_profile = sf.shape(i)
+            z = s_profile.z
+
+            ind = 0
+            points = []
+            for point in s_profile.points:
+                points.append([
+                    point[0],
+                    point[1],
+                    z[ind]
+                ])
+
+                ind += 1
+
+            prof = ProfileXYZ(
+                reach=self, status=self._status
+            )
+            prof.import_points(points)
+
+            profiles.append(prof)
+
+        self.profiles = profiles + self.profiles
+        self._update_profile_numbers()
+        self._recompute_kp()
+
+        return profiles
+
+    @timer
+    def import_geometry_st(self, file_path_name: str):
         """Import a geometry from file (.ST or .st)
 
         Args:
@@ -762,3 +814,24 @@ class Reach(SQLSubModel):
         except Exception as e:
             logger_exception(e)
             return 0
+
+    def _recompute_kp(self, offset=0.0):
+        self._recompute_kp_no_gl(offset=offset)
+
+    def _recompute_kp_no_gl(self, offset=0.0):
+        profiles = iter(self.profiles)
+
+        previous = next(profiles)
+        previous.kp = offset
+
+        for profile in profiles:
+            prev_points = previous.points
+            curr_points = profile.points
+
+            dist = (
+                prev_points[0].dist_2d(curr_points[0]) +
+                prev_points[-1].dist_2d(curr_points[-1])
+            ) / 2.0
+
+            profile.kp = previous.kp + dist
+            previous = profile
diff --git a/src/View/Geometry/Translate.py b/src/View/Geometry/Translate.py
index 7b607d83c537adb783a9242fef6ae383131b0b9b..23034fd51475dd711e868528cc381cc8963ecb90 100644
--- a/src/View/Geometry/Translate.py
+++ b/src/View/Geometry/Translate.py
@@ -36,6 +36,8 @@ class GeometryTranslate(MainTranslate):
             "Geometry", "File mage geometry (*.ST *.st)")
         self._dict["file_m"] = _translate(
             "Geometry", "File mage meshed geometry (*.M *.m)")
+        self._dict["file_shp"] = _translate(
+            "Geometry", "Shapefile (*.SHP *.shp)")
         self._dict["file_all"] = _translate("Geometry", "All file (*)")
 
         self._dict["cross_section"] = _translate("Geometry", "cross-section")
diff --git a/src/View/Geometry/Window.py b/src/View/Geometry/Window.py
index d8bcc06b47c8e913659c5b26ccb13637197f40c0..3631ff1003cd7dbc2541993ac2b92ce48e43376f 100644
--- a/src/View/Geometry/Window.py
+++ b/src/View/Geometry/Window.py
@@ -241,6 +241,7 @@ class GeometryWindow(PamhyrWindow):
         file_types = [
             self._trad["file_st"],
             self._trad["file_m"],
+            self._trad["file_shp"],
             self._trad["file_all"],
         ]