diff --git a/scripts/python3/data_manip/formats/mascaretgeo_file.py b/scripts/python3/data_manip/formats/mascaretgeo_file.py
index 8d9582f10634551ef2c5ec6eaf5a166b859a680b..249b1228e10c1f7449fe750c78c9d1f3cfaae143 100644
--- a/scripts/python3/data_manip/formats/mascaretgeo_file.py
+++ b/scripts/python3/data_manip/formats/mascaretgeo_file.py
@@ -5,7 +5,8 @@ from data_manip.formats.mascaret_file import Reach, Section
 from utils.exceptions import MascaretException
 
 
-class MascaretGeoFile:
+
+class MascaretGeoFile():
     """
     Parse Mascaret geometry file (geo/geoC/georef/georefC)
     Handles multiple reaches
@@ -23,9 +24,10 @@ class MascaretGeoFile:
     """
     OUTPUT_FLOAT_FMT = '%.6f'
 
-    def __init__(self, file_name, fformat=None, mode='read'):
+    def __init__(self, file_name, load_file=True, fformat=None, mode='read'):
         """
         @param file_name (str) file name
+        @param load_file (bool) option to activate loading of the file
         @param fformat (str) file format ('opt' or 'rub')
         @param mode (str) define the mode for the class,
             'read' by default to read a file,
@@ -44,18 +46,19 @@ class MascaretGeoFile:
                 self.fformat = os.path.splitext(file_name)[1][1:]
             else:
                 self.fformat = fformat.lower().strip()
-            if self.fformat not in ('geo', 'georef'):
+            if self.fformat not in ('geo', 'geoC', 'georef', 'georefC'):
                 raise NotImplementedError(
                     'Format `%s` not supported,\
-                     only geo and georef formats are supported as input' %
-                    self.fformat)
+                     only geo, geoC, georef and georefC formats are supported\
+                     as input' % self.fformat)
             self.has_ref = 'ref' in self.fformat
 
             # Layers for sediments (Courlis)
             self.has_layers = self.fformat.endswith('C')
 
-            # Load file content
-            self.load()
+            if load_file:
+                # Load file content
+                self.load()
 
     def load(self):
         """
@@ -67,8 +70,9 @@ class MascaretGeoFile:
             section_id = 0
             section_name = ''
             section_pk = -1.0
-            dist, x_list, y_list, z_list = [], [], [], []
-            xa, ya = None, None
+            dist, x_list, y_list, z_list, topo_bath_list, layers_elev_list = \
+                [], [], [], [], [], []
+            x_axis, y_axis = None, None
 
             for line in filein:
                 if line.startswith('#'):
@@ -78,19 +82,25 @@ class MascaretGeoFile:
                     if dist:
                         # Add previous Section
                         section = Section(section_id, section_pk, section_name)
-                        if xa is not None and ya is not None:
-                            section.set_axis(xa, ya)
+                        if x_axis is not None and y_axis is not None:
+                            section.set_axis(x_axis, y_axis)
                         if self.has_ref:
-                            section.set_points_from_xyz(x_list, y_list, z_list)
+                            section.set_points_from_xyz(x_list, y_list, z_list,
+                                                        topo_bath_list)
                         else:
-                            section.set_points_from_trans(dist, z_list)
+                            section.set_points_from_trans(dist, z_list,
+                                                          topo_bath_list)
+                        section.distances = dist
+                        if self.has_layers:
+                            section.add_layers_from_elevations(
+                                    layers_elev_list)
                         reach.add_section(section)
 
                     if self.has_ref:
                         _, reach_name, section_name, pk_str, x1, y1, x2, y2,\
-                         _, xa, ya = line.split()
-                        xa = float(xa)
-                        ya = float(ya)
+                         _, x_axis, y_axis = line.split()
+                        x_axis = float(x_axis)
+                        y_axis = float(y_axis)
                     else:
                         _, reach_name, section_name, pk_str = line.split()
 
@@ -108,28 +118,89 @@ class MascaretGeoFile:
 
                     # Reset variables to store section
                     section_pk = float(pk_str)
-                    dist, x_list, y_list, z_list = [], [], [], []
+                    dist, x_list, y_list, z_list, topo_bath_list, \
+                        layers_elev_list = [], [], [], [], [], []
                     section_id += 1
                 else:
-                    if self.has_ref:
-                        dist_str, z_str, _, x1, y1 = line.split()
-                        x_list.append(float(x1))
-                        y_list.append(float(y1))
+                    list_line = line.split()
+
+                    if self.has_ref and self.has_layers is False:
+                        dist_str = list_line[0]
+                        z_str = list_line[1]
+                        if "T" in list_line or "B" in list_line:
+                            shift = 0
+                            topo_bath = list_line[2]
+                        else:
+                            shift = -1
+                            topo_bath = "B"
+
+                        x_list.append(float(list_line[shift + 3]))
+                        y_list.append(float(list_line[shift + 4]))
+
+                    elif self.has_ref and self.has_layers:
+                        dist_str = list_line[0]
+                        z_str = list_line[1]
+                        if "T" in list_line or "B" in list_line:
+                            shift = 0
+                            topo_bath = list_line[-3]
+                        else:
+                            shift = 1
+                            topo_bath = "B"
+
+                        layers_elev = list_line[2:shift - 3]
+
+                        x_list.append(float(list_line[-2]))
+                        y_list.append(float(list_line[-1]))
+
+                        layers_elev = [float(elev) for elev in layers_elev]
+
+                        if self.nlayers == 0:
+                            self.nlayers = len(layers_elev)
+
+                    elif self.has_ref is False and self.has_layers:
+                        dist_str = list_line[0]
+                        z_str = list_line[1]
+                        if "T" in list_line or "B" in list_line:
+                            layers_elev = list_line[2:-1]
+                            topo_bath = list_line[-1]
+                        else:
+                            layers_elev = list_line[2:]
+                            topo_bath = "B"
+
+                        layers_elev = [float(elev) for elev in layers_elev]
+
+                        if self.nlayers == 0:
+                            self.nlayers = len(layers_elev)
+
                     else:
-                        dist_str, z_str, _ = line.split()
+                        dist_str = list_line[0]
+                        z_str = list_line[1]
+
+                        if "T" in list_line or "B" in list_line:
+                            topo_bath = list_line[2]
+                        else:
+                            topo_bath = "B"
 
                     # Add new point to current section
                     dist.append(float(dist_str))
                     z_list.append(float(z_str))
+                    topo_bath_list.append(topo_bath)
+                    if self.has_layers:
+                        layers_elev_list.append(layers_elev)
 
             # Add last section
             section = Section(section_id, section_pk, section_name)
-            if xa is not None and ya is not None:
-                section.set_axis(xa, ya)
+            if x_axis is not None and y_axis is not None:
+                section.set_axis(x_axis, y_axis)
             if self.has_ref:
-                section.set_points_from_xyz(x_list, y_list, z_list)
+                section.set_points_from_xyz(x_list, y_list, z_list,
+                                            topo_bath_list)
             else:
-                section.set_points_from_trans(dist, z_list)
+                section.set_points_from_trans(dist, z_list, topo_bath_list)
+            section.distances = dist
+            if self.has_layers:
+                section.add_layers_from_elevations(layers_elev_list)
+                self.layer_names = section.layer_names
             reach.add_section(section)
 
     def save(self, output_file_name):
@@ -149,60 +220,79 @@ class MascaretGeoFile:
 
         if ref and not self.has_ref:
             raise MascaretException('Could not write `%s` format without\
-             any geo-referenced data' % fformat)
+                    any geo-referenced data' % fformat)
 
-        with open(output_file_name, 'w') as fileout:
+        with open(output_file_name, 'w', encoding="utf-8") as fileout:
             for _, reach in self.reaches.items():
 
                 for sec in reach:
                     positions_str = ''
                     if ref:
                         # Get river_banks and `AXE` coordinates if necessary
-                        xa, ya = sec.axis
-                        positions_str += ' %f %f %f %f' %\
+                        x_axis, y_axis = sec.axis
+                        positions_str += '%f %f %f %f' %\
                             (sec.x[0], sec.y[0], sec.x[-1], sec.y[-1])
-                        positions_str += ' AXE %f %f' % (xa, ya)
+                        positions_str += ' AXE %f %f' % (x_axis, y_axis)
 
                     # Write profile header
                     fileout.write(
-                        'Profil %s %s %f%s\n' %
+                        'Profil %s %s %f %s\n' %
                         (reach.name, sec.name, sec.pk, positions_str))
 
                     # Write points and layers if necessary
                     if not ref and not layers:
-                        for dist, z in zip(sec.distances, sec.z):
-                            fileout.write('%f %f B\n' % (dist, z))
+                        for dist, z, topo_bath in zip(sec.distances, sec.z,
+                                                      sec.topo_bath):
+                            fileout.write('%f %f %s\n' % (dist, z, topo_bath))
 
                     elif ref and not layers:
-                        for dist, x, y, z in zip(sec.distances,
-                                                 sec.x, sec.y, sec.z):
-                            fileout.write('%f %f B %f %f\n' % (dist, z, x, y))
+                        for dist, x, y, z, topo_bath in zip(sec.distances,
+                                                            sec.x,
+                                                            sec.y,
+                                                            sec.z,
+                                                            sec.topo_bath):
+                            fileout.write('%f %f %s %f %f\n'
+                                          % (dist, z, topo_bath, x, y))
 
                     elif not ref and layers:
-                        for i, (dist, z) in enumerate(zip(sec.distances,
-                                                          sec.z)):
+                        for i, (dist, z, topo_bath) in \
+                                enumerate(zip(sec.distances, sec.z,
+                                              sec.topo_bath)):
                             if self.nlayers == 0:
                                 layers_str = ''
                             else:
                                 layers_str = ' ' +\
-                                    ' '.join(
-                                        [MascaretGeoFile.OUTPUT_FLOAT_FMT % zl
-                                            for zl in sec.layers_elev[:, i]])
-                            fileout.write('%f %f%s B\n' %
-                                          (dist, z, layers_str))
+                                        ' '.join(
+                                                [MascaretGeoFile.OUTPUT_FLOAT_FMT % zl
+                                                    for zl in sec.layers_elev[:, i]])
+
+                            try:
+                                fileout.write('%f %f %s %s\n' %
+                                            (dist, z, layers_str, topo_bath))
+                            except Exception as e:
+                                raise e
 
                     elif ref and layers:
-                        for i, (dist, x, y, z) in enumerate(zip(sec.distances,
-                                                            sec.x, sec.y,
-                                                            sec.z)):
+                        for i, (dist, x, y, z, topo_bat) in enumerate(zip(sec.distances,
+                            sec.x,
+                            sec.y,
+                            sec.z,
+                            sec.topo_bath)):
+
                             if self.nlayers == 0:
                                 layers_str = ''
                             else:
                                 layers_str = ' ' + ' '\
-                                    .join([MascaretGeoFile.OUTPUT_FLOAT_FMT %
-                                          zl for zl in sec.layers_elev[:, i]])
-                            fileout.write('%f %f%s B %f %f\n' %
-                                          (dist, z, layers_str, x, y))
+                                        .join([MascaretGeoFile.OUTPUT_FLOAT_FMT %
+                                            zl for zl in sec.layers_elev[:, i]])
+
+                            try:
+                                fileout.write('%f %f %s %s %f %f\n'
+                                          % (dist, z, layers_str,
+                                             topo_bath, x, y))
+                            except Exception as e:
+                                raise e
+
 
     def __repr__(self):
         return 'MascaretGeoFile: %s' % self.file_name
@@ -241,8 +331,8 @@ if __name__ == '__main__':
                                    'examples', 'mascaret'), '*.geo')
         geo_files += recursive_glob(os.path.join(os.environ['HOMETEL'],
                                     'examples', 'mascaret'), 'geometrie')
-        for file_name in sorted(geo_files):
-            geo_file = MascaretGeoFile(file_name, 'geo')
+        for filename in sorted(geo_files):
+            geo_file = MascaretGeoFile(filename, 'geo')
             print(geo_file.summary())
     except MascaretException as e:
         print(str(e))