From 9b44e4859fc1883676866afa4ffc014731d632d0 Mon Sep 17 00:00:00 2001
From: Lionel GUEZ <guez@lmd.ens.fr>
Date: Thu, 5 Oct 2023 21:42:20 +0200
Subject: [PATCH] Give special treatment to phantom patterns
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit
This is case 2 correction in the Matlab program of Rémi Laxenaire.
---
Trajectories/trajectories.py | 51 ++++++++++++++++++++++++++++++++++++
1 file changed, 51 insertions(+)
diff --git a/Trajectories/trajectories.py b/Trajectories/trajectories.py
index 0ed32247..252a31c2 100755
--- a/Trajectories/trajectories.py
+++ b/Trajectories/trajectories.py
@@ -16,6 +16,8 @@ import graph_tool
from graph_tool import topology
import json
import sys
+import report_graph
+import numpy as np
def new_traj(ind_traj, traj_prop, n, traj_vert_ind):
"""Assign new trajectory to vertex index n."""
@@ -61,6 +63,55 @@ traj_vert_ind = []
# of vertex indices.
for n in topology.topological_sort(g):
+ # First, we check whether we have a phantom pattern. Phantom
+ # eddies may appear in altimetric data because of
+ # interpolation. In a phantom pattern, we have splitting
+ # immediately followed by merging. We do not want to let phantom
+ # patterns interrupt trajectories.
+
+ if g.vertex(n).in_degree() == 2:
+ pred = g.get_in_neighbors(n)
+
+ if np.all(g.get_in_degrees(pred) == 1) \
+ and np.all(g.get_out_degrees(pred) == 1):
+ n1 = g.get_in_neighbors(pred[0])[0]
+ n2 = g.get_in_neighbors(pred[1])[0]
+
+ if n1 == n2 and g.vertex(n1).out_degree() == 2:
+ # We have a splitting event at n1 immediately followed
+ # by a merging at n. Check the time interval between
+ # splitting and merging. The date of merging is the
+ # date of the first eddy in segment n and the date of
+ # splitting is the date of the last eddy in segment
+ # n1.
+ date_merging = report_graph.node_to_date_eddy(
+ g.vp.name[n], g.gp.e_overestim, only_date = True)
+ date_splitting = report_graph.node_to_date_eddy(
+ g.vp.inst_eddies[n1][- 1], g.gp.e_overestim,
+ only_date = True)
+
+ if date_merging - date_splitting <= 6:
+ # We have a phantom pattern. Get the vertex on the
+ # side of the shortest path:
+ v_short = topology.shortest_path(
+ g, n1, n, weights = g.ep.cost_function,
+ dag = True)[0][1]
+
+ if traj_prop[v_short] != traj_prop[n1]:
+ # Correct the trajectories that were assigned
+ # according to closest successor, that is:
+ # swap the trajectories of pred[0] and
+ # pred[1]:
+ traj_prop[pred[0]], traj_prop[pred[1]] \
+ = traj_prop[pred[1]], traj_prop[pred[0]]
+ traj_vert_ind[traj_prop[pred[0]]][- 1], \
+ traj_vert_ind[traj_prop[pred[1]]][- 1] \
+ = traj_vert_ind[traj_prop[pred[1]]][- 1], \
+ traj_vert_ind[traj_prop[pred[0]]][- 1]
+
+ traj_prop[n] = traj_prop[n1]
+ traj_vert_ind[traj_prop[n1]].append(n)
+
if traj_prop[n] == - 1:
if g.vertex(n).in_degree() >= 1:
# Find the index of the closest direct predecessor of a node:
--
GitLab