diff --git a/covjsonkit/encoder/TimeSeries.py b/covjsonkit/encoder/TimeSeries.py
index 310f630..4ed7395 100644
--- a/covjsonkit/encoder/TimeSeries.py
+++ b/covjsonkit/encoder/TimeSeries.py
@@ -143,6 +143,7 @@ def from_polytope(self, result, date_key: str = "date") -> dict:
         fields["step"] = 0
         fields["dates"] = []
         fields["levels"] = [0]
+        fields["identifiers"] = []
 
         start = time.time()
         logging.debug("Tree walking starts at: %s", start)  # noqa: E501
@@ -243,13 +244,22 @@ def from_polytope(self, result, date_key: str = "date") -> dict:
                                         f"Key {key} not found in range_dict. "
                                         f"Please ensure all axes are compressed in config"
                                     )
-                        mm = mars_metadata.copy()
-                        mm["number"] = num
-                        mm["Forecast date"] = date
-                        mm["levelist"] = level
-                        coordinates[date][i]["levelist"] = [level]
-                        del mm["step"]
-                        self.add_coverage(mm, coordinates[date][i], val_dict)
+                        # Determine identifiers for this point
+                        point_identifiers = [None]
+                        if fields["identifiers"] and i < len(fields["identifiers"]):
+                            point_identifiers = fields["identifiers"][i]
+
+                        # Emit one coverage per identifier (duplicates coverage for merged points)
+                        for identifier in point_identifiers:
+                            mm = mars_metadata.copy()
+                            mm["number"] = num
+                            mm["Forecast date"] = date
+                            mm["levelist"] = level
+                            if identifier is not None:
+                                mm["label"] = identifier
+                            coordinates[date][i]["levelist"] = [level]
+                            del mm["step"]
+                            self.add_coverage(mm, coordinates[date][i], val_dict)
 
         end = time.time()
         delta = end - start
diff --git a/covjsonkit/encoder/encoder.py b/covjsonkit/encoder/encoder.py
index 0332de2..2bcafcd 100644
--- a/covjsonkit/encoder/encoder.py
+++ b/covjsonkit/encoder/encoder.py
@@ -383,6 +383,17 @@ def append_composite_coords(dates, tree_values, lat, coords):
             for value in tree_values:
                 coords[dates]["composite"].append([lat, value])
 
+        def collect_tags(tree, fields, count):
+            """Collect tags from the latitude node into fields['identifiers']."""
+            if "identifiers" in fields:
+                tags = fields.get("_lat_tags", None) or getattr(tree, "tags", None)
+                if tags:
+                    tag_list = sorted(tags, key=str)
+                else:
+                    tag_list = [None]
+                for _ in range(count):
+                    fields["identifiers"].append(tag_list)
+
         if len(tree.children) != 0:
             for child in tree.children:
                 handle_non_leaf_node(child)
@@ -390,6 +401,7 @@ def append_composite_coords(dates, tree_values, lat, coords):
                 if result is not None:
                     if child.axis.name == "latitude":
                         fields["lat"] = result
+                        fields["_lat_tags"] = getattr(child, "tags", None)
                     elif child.axis.name == "levelist":
                         fields["levels"] = result
                         if "l" in fields:
@@ -426,6 +438,7 @@ def append_composite_coords(dates, tree_values, lat, coords):
                 step_len = para_len / len(fields["step"])
 
                 append_composite_coords(fields["dates"][-1], tree.values, fields["lat"], coords)
+                collect_tags(tree, fields, len(tree.values))
 
                 for l, level in enumerate(fields["levels"]):  # noqa: E741
                     for i, num in enumerate(fields["number"]):