diff --git a/nonos/_readers/binary.py b/nonos/_readers/binary.py
index 9ea3258c..f270749a 100644
--- a/nonos/_readers/binary.py
+++ b/nonos/_readers/binary.py
@@ -29,6 +29,8 @@
 
 @final
 class VTKReader(ReaderMixin):
+    NATIVE_COORDINATE_REGEXP = re.compile(r"X(1|2|3)(L|C)_NATIVE_COORDINATES")
+
     @staticmethod
     def parse_output_number_and_filename(
         file_or_number: Union[PathT, int],
@@ -81,6 +83,8 @@ def read(file, /, **meta) -> BinData:
         if meta["geometry"] is not None:
             V["geometry"] = Geometry(meta["geometry"])
 
+        native_coordinates = {}
+
         # datatype we read
         dt = np.dtype(">f")  # Big endian single precision floats
         dint = np.dtype(">i4")  # Big endian integer
@@ -136,6 +140,11 @@ def read(file, /, **meta) -> BinData:
                 elif entry == "PERIODICITY":
                     # skip
                     fid.seek(dint.itemsize * 3, os.SEEK_CUR)
+                elif VTKReader.NATIVE_COORDINATE_REGEXP.match(entry):
+                    native_name, _ncomp, native_dim, _dtype = entry.split()
+                    native_coordinates[native_name] = np.fromfile(
+                        fid, dtype=dt, count=int(native_dim)
+                    )
                 else:  # pragma: no cover
                     fid.close()
                     raise ValueError(f"Received unknown field: {entry!r}")
@@ -470,6 +479,35 @@ def read(file, /, **meta) -> BinData:
                 fid.readline()  # extra line feed
         fid.close()
 
+        geom: Geometry = V["geometry"]
+        if geom is Geometry.SPHERICAL:
+            native2attr = {
+                "X1L_NATIVE_COORDINATES": "rl",
+                "X1C_NATIVE_COORDINATES": "r",
+                "X2L_NATIVE_COORDINATES": "thetal",
+                "X2C_NATIVE_COORDINATES": "theta",
+                "X3L_NATIVE_COORDINATES": "phil",
+                "X3C_NATIVE_COORDINATES": "phi",
+            }
+        elif (
+            geom is Geometry.CARTESIAN
+            or geom is Geometry.CYLINDRICAL
+            or geom is Geometry.POLAR
+        ):
+            native2attr = {
+                "X1L_NATIVE_COORDINATES": "xl",
+                "X1C_NATIVE_COORDINATES": "x",
+                "X2L_NATIVE_COORDINATES": "yl",
+                "X2C_NATIVE_COORDINATES": "y",
+                "X3L_NATIVE_COORDINATES": "zl",
+                "X3C_NATIVE_COORDINATES": "z",
+            }
+        else:
+            assert_never(geom)
+
+        for native_field, attr in native2attr.items():
+            V[attr][native_field] = native_coordinates[native_field]
+
         return BinData(**V)