Repository: demiangomez/Parallel.GAMIT (language: python)

Issue:
Is checking for regional vs global still needed?
It appears that we are still checking for two different solution types

https://github.com/demiangomez/Parallel.GAMIT/blob/99a7a4dc85ce9dd34af5bbe716733ef0a00c6205/pgamit/pyNetwork.py#L167

but the calls to make_clusters is the same. If the regional branch is not needed anymore, please remove it or duplicate the if statement that checks for the number of stations (the the else branch of the above referenced if). We just run a test and it failed:

 >> 2024-10-17 16:37:38 g2regain 2024 268 -> Creating network clusters
 --  Processing type is regional with 7 active stations
Traceback (most recent call last):
  File "/fs/project/gomez.124/opt/conda/pgamit/bin/ParallelGamit.py", line 756, in <module>
    main()
  File "/fs/project/gomez.124/opt/conda/pgamit/bin/ParallelGamit.py", line 399, in main
    sessions = ExecuteGamit(cnn, JobServer, GamitConfig, stations, check_stations, args.ignore_missing, dates,
  File "/fs/project/gomez.124/opt/conda/pgamit/bin/ParallelGamit.py", line 700, in ExecuteGamit
    net_object = Network(cnn, archive, GamitConfig, stations, date, check_stations, ignore_missing)
  File "/fs/project/gomez.124/opt/conda/pgamit/lib/python3.10/site-packages/pgamit/pyNetwork.py", line 170, in __init__
    clusters, ties = self.make_clusters(stations.get_active_coordinates(date), stn_active)
  File "/fs/project/gomez.124/opt/conda/pgamit/lib/python3.10/site-packages/pgamit/pyNetwork.py", line 208, in make_clusters
    OC = over_cluster(qmean.labels_, points, metric='euclidean',
  File "/fs/project/gomez.124/opt/conda/pgamit/lib/python3.10/site-packages/pgamit/NetClusters.py", line 174, in over_cluster
    D, _ = nbrs.kneighbors(coordinates[nonmembers, :])
  File "/fs/project/gomez.124/opt/conda/pgamit/lib/python3.10/site-packages/sklearn/neighbors/_base.py", line 825, in kneighbors
    X = self._validate_data(X, accept_sparse="csr", reset=False, order="C")
  File "/fs/project/gomez.124/opt/conda/pgamit/lib/python3.10/site-packages/sklearn/base.py", line 633, in _validate_data
    out = check_array(X, input_name="X", **check_params)
  File "/fs/project/gomez.124/opt/conda/pgamit/lib/python3.10/site-packages/sklearn/utils/validation.py", line 1087, in check_array
    raise ValueError(
ValueError: Found array with 0 sample(s) (shape=(0, 3)) while a minimum of 1 is required by NearestNeighbors.

the reason is because this day only has 7 active stations (it is a recent day, not all RINEX available yet).



Discussion / hints:

Fix patch:
diff --git a/com/ParallelGamit.py b/com/ParallelGamit.py
index fcaa57ae..6c7bc931 100644
--- a/com/ParallelGamit.py
+++ b/com/ParallelGamit.py
@@ -34,7 +34,7 @@
 from pgamit import pyParseZTD
 from pgamit import pyArchiveStruct
 from pgamit.pyETM import pyETMException
-from pgamit.pyNetwork import Network
+from pgamit.network import Network
 from pgamit.pyStation import Station, StationCollection
 from pgamit.Utils import (process_date,
                           process_stnlist,
diff --git a/pgamit/__init__.py b/pgamit/__init__.py
index 74d5aa1e..65571a83 100644
--- a/pgamit/__init__.py
+++ b/pgamit/__init__.py
@@ -1,5 +1,7 @@
 __version__ = "1.2.3"
 __all__ = [
+    'cluster',
+    'network',
     'pyRinexName',
     'Utils',
     'pyJobServer',
@@ -27,12 +29,10 @@
     'pyProducts',
     'pyZTD',
     'pyStack',
-    'NetClusters',
     'pyGamitConfig',
     'pyGamitSession',
     'pyGamitTask',
     'pyGlobkTask',
-    'pyNetwork',
     'pyParseZTD',
     'pyStation'
 ]
diff --git a/pgamit/NetClusters.py b/pgamit/cluster.py
similarity index 96%
rename from pgamit/NetClusters.py
rename to pgamit/cluster.py
index 6233cda2..13e24d0c 100644
--- a/pgamit/NetClusters.py
+++ b/pgamit/cluster.py
@@ -202,67 +202,14 @@ def over_cluster(labels, coordinates, metric='haversine', neighborhood=5,
 # Modifications by: Shane Grigsby <refuge@rocktalus.com>
 
 
-class _BisectingTree:
-    """Tree structure representing the hierarchical clusters from bisecting"""
-
-    def __init__(self, center, indices, score):
-        """Create a new cluster node in the tree.
-
-        The node holds the center of this cluster and the indices of the data
-        points that belong to it.
-        """
-        self.center = center
-        self.indices = indices
-        self.score = score
-
-        self.left = None
-        self.right = None
-
-    def split(self, labels, centers, scores):
-        """Split the cluster node into two subclusters."""
-        self.left = _BisectingTree(indices=self.indices[labels == 0],
-                                   center=centers[0], score=scores[0])
-        self.right = _BisectingTree(indices=self.indices[labels == 1],
-                                    center=centers[1], score=scores[1])
-
-        # reset the indices attribute to save memory
-        self.indices = None
-
-    def get_cluster_to_bisect(self):
-        """Return the cluster node to bisect next.
-
-        It's based on the score of the cluster, which can be either the number
-        of data points assigned to that cluster or the inertia of that cluster.
-        """
-        max_score = None
-
-        for cluster_leaf in self.iter_leaves():
-            if max_score is None or cluster_leaf.score > max_score:
-                max_score = cluster_leaf.score
-                best_cluster_leaf = cluster_leaf
-
-        if np.isneginf(max_score):
-            self.bisect = False
-        else:
-            return best_cluster_leaf
-
-    def iter_leaves(self):
-        """Iterate over all the cluster leaves in the tree."""
-        if self.left is None:
-            yield self
-        else:
-            yield from self.left.iter_leaves()
-            yield from self.right.iter_leaves()
-
-
 class BisectingQMeans(_BaseKMeans):
     """Bisecting Q-Means clustering; modified from sklearn Bisecting K-means.
 
     In contrast to Bisecting K-Means, Bisecting Q-Means clustering will infer
     the number of clusters based on a termination condition. For this
     implementation the bisecting termination occurs according to the minimum
-    and optimum cluster sizes, which are set by the `opt_clust_size` and
-    `min_clust_size` parameters respectively. The child cluster of the bisected
+    and optimum cluster sizes, which are set by the `opt_size` and
+    `min_size` parameters respectively. The child cluster of the bisected
     root cluster with the biggest inertia as determined by SSE (Sum of Squared
     Errors) will be selected bisection-- provided that the child cluster
     exceeds the set `*clust_size` boundary conditions. Cluster bisection
@@ -271,19 +218,19 @@ class BisectingQMeans(_BaseKMeans):
 
     Parameters
     ----------
-    min_clust_size : int, default=4
+    min_size : int, default=4
         The minimum acceptable cluster size. Clusters of size <= to this
         parameter will **not** be produced by this algorithm.
 
-    opt_clust_size : int, default=20
+    opt_size : int, default=12
         Target optimum cluster size. If the sum membership of a proposed
         cluster bisection is less than this value, the cluster will not be
-        bisected. When combined with the `min_clust_size` parameter above,
+        bisected. When combined with the `min_size` parameter above,
         these conditions together mean that clusters of sizes smaller than
-        (`opt_clust_size` - `min_clust_size`) are *a priori* ineligible to be
+        (`opt_size` - `min_size`) are *a priori* ineligible to be
         bisected.
 
-    max_clust_size: int, default=35
+    max_size: int, default=25
         Hard cutoff to bypass the heuristic when bisecting clusters; no
         clusters greater than this size will be produced.
 
@@ -383,9 +330,9 @@ class BisectingQMeans(_BaseKMeans):
 
     def __init__(
         self,
-        min_clust_size=4,
-        opt_clust_size=20,
-        max_clust_size=35,
+        min_size=4,
+        opt_size=12,
+        max_size=25,
         *,
         init="random",
         n_init=1,
@@ -407,9 +354,9 @@ def __init__(
             n_init=n_init,
         )
 
-        self.min_clust_size = min_clust_size
-        self.opt_clust_size = opt_clust_size
-        self.max_clust_size = max_clust_size
+        self.min_size = min_size
+        self.opt_size = opt_size
+        self.max_size = max_size
         self.copy_x = copy_x
         self.algorithm = algorithm
         self.bisect = True
@@ -516,17 +463,17 @@ def _bisect(self, X, x_squared_norms, sample_weight, cluster_to_bisect):
                                            sample_weight)
         counts = np.bincount(best_labels, minlength=2)
         scores[np.where(counts <
-                        (self.opt_clust_size - self.min_clust_size))] = -np.inf
+                        (self.opt_size - self.min_size))] = -np.inf
         # case where bisecting is not optimum
-        if (counts[0] + counts[1]) < self.opt_clust_size:
+        if (counts[0] + counts[1]) < self.opt_size:
             cluster_to_bisect.score = -np.inf
         # bisect as long as the smallest child meets membership constraints
-        elif ((counts[0] >= self.min_clust_size) and
-              (counts[1] >= self.min_clust_size)):
+        elif ((counts[0] >= self.min_size) and
+              (counts[1] >= self.min_size)):
             cluster_to_bisect.split(best_labels, best_centers, scores)
         # one child will have membership of 3 or less; don't split
         else:
-            if (counts[0] + counts[1] >= self.max_clust_size):
+            if (counts[0] + counts[1] >= self.max_size):
                 cluster_to_bisect.split(best_labels, best_centers, scores)
             else:
                 cluster_to_bisect.score = -np.inf
@@ -637,3 +584,57 @@ def fit(self, X, y=None, sample_weight=None):
 
     def _more_tags(self):
         return {"preserves_dtype": [np.float64, np.float32]}
+
+
+class _BisectingTree:
+    """Tree structure representing the hierarchical clusters from bisecting"""
+
+    def __init__(self, center, indices, score):
+        """Create a new cluster node in the tree.
+
+        The node holds the center of this cluster and the indices of the data
+        points that belong to it.
+        """
+        self.center = center
+        self.indices = indices
+        self.score = score
+
+        self.left = None
+        self.right = None
+
+    def split(self, labels, centers, scores):
+        """Split the cluster node into two subclusters."""
+        self.left = _BisectingTree(indices=self.indices[labels == 0],
+                                   center=centers[0], score=scores[0])
+        self.right = _BisectingTree(indices=self.indices[labels == 1],
+                                    center=centers[1], score=scores[1])
+
+        # reset the indices attribute to save memory
+        self.indices = None
+
+    def get_cluster_to_bisect(self):
+        """Return the cluster node to bisect next.
+
+        It's based on the score of the cluster, which can be either the number
+        of data points assigned to that cluster or the inertia of that cluster.
+        """
+        max_score = None
+
+        for cluster_leaf in self.iter_leaves():
+            if max_score is None or cluster_leaf.score > max_score:
+                max_score = cluster_leaf.score
+                best_cluster_leaf = cluster_leaf
+
+        if np.isneginf(max_score):
+            self.bisect = False
+        else:
+            return best_cluster_leaf
+
+    def iter_leaves(self):
+        """Iterate over all the cluster leaves in the tree."""
+        if self.left is None:
+            yield self
+        else:
+            yield from self.left.iter_leaves()
+            yield from self.right.iter_leaves()
+
diff --git a/pgamit/pyNetwork.py b/pgamit/network.py
similarity index 90%
rename from pgamit/pyNetwork.py
rename to pgamit/network.py
index 695edfcc..47b901d1 100644
--- a/pgamit/pyNetwork.py
+++ b/pgamit/network.py
@@ -37,7 +37,7 @@
 # app
 from pgamit.pyGamitSession import GamitSession
 from pgamit.pyStation import StationCollection
-from pgamit.NetClusters import over_cluster, select_central_point, BisectingQMeans
+from pgamit.cluster import over_cluster, select_central_point, BisectingQMeans
 from pgamit.NetPlots import plot_global_network
 
 BACKBONE_NET = 45
@@ -164,42 +164,21 @@ def __init__(self, cnn, archive, GamitConfig, stations, date,
             tqdm.write(' --  Processing type is %s with %i active stations'
                        % (GamitConfig.NetworkConfig['type'], len(stn_active)))
 
-            if GamitConfig.NetworkConfig['type'] == 'regional':
-                # create station clusters
-                # cluster centroids will be used later to tie the networks
+            if len(stn_active) > BACKBONE_NET + 5:
+                backbone = self.backbone_delauney(stations.get_active_coordinates(date), stn_active)
                 clusters, ties = self.make_clusters(stations.get_active_coordinates(date), stn_active)
-
-                if len(clusters['stations']) > 1:
-                    # build the backbone network
-                    backbone = self.backbone_delauney(stations.get_active_coordinates(date), stn_active)
-                else:
-                    backbone = []
             else:
-                # DDG: if active stations is greater than BACKBONE_NET + 5,
-                # then we need to split the processing into smaller subnets.
-                # If not, then we just use all stations as the processing
-                # network. We add 5 to BACKBONE_NET to create a hysteresis
-                # behavior in subnets_delaunay. In other words, the backbone
-                # net will have BACKBONE_NET stations and 6 stations will be
-                # available to create some subnets. A single network solution
-                # will have a max size of BACKBONE_NET + 5
-                if len(stn_active) > BACKBONE_NET + 5:
-                    backbone = self.backbone_delauney(stations.get_active_coordinates(date), stn_active)
-                    clusters, ties = self.make_clusters(stations.get_active_coordinates(date), stn_active)
-                else:
-                    # no need to create a set of clusters, just use them all
-                    clusters = {'stations': [stn_active]}
-                    backbone = []
-                    ties = []
+                # no need to create a set of clusters, just use them all
+                clusters = {'stations': [stn_active]}
+                backbone = []
+                ties = []
 
         self.sessions = self.create_gamit_sessions(cnn, archive, clusters,
                                                    backbone, ties, date)
 
     def make_clusters(self, points, stations, net_limit=NET_LIMIT):
         # Run initial clustering using bisecting 'q-means'
-        qmean = BisectingQMeans(min_clust_size=4, opt_clust_size=20,
-                                init='random', n_init=50, algorithm='lloyd',
-                                max_iter=8000, random_state=42, n_clusters=2)
+        qmean = BisectingQMeans(min_size=2, random_state=42)
         qmean.fit(points)
         # snap centroids to closest station coordinate
         central_points = select_central_point(qmean.labels_, points,


Test patch:
diff --git a/pgamit/tests/test_make_clusters.py b/pgamit/tests/test_make_clusters.py
index 72cf7f7f..451eeebb 100644
--- a/pgamit/tests/test_make_clusters.py
+++ b/pgamit/tests/test_make_clusters.py
@@ -5,11 +5,11 @@
 import numpy as np
 
 from .common import gen_variable_density_clusters, generate_clustered_data
-from ..NetClusters import BisectingQMeans, over_cluster
+from ..cluster import BisectingQMeans, over_cluster
 
 
 @pytest.mark.parametrize(
-    ("max_clust_size", "clust_size"),
+    ("max_size", "clust_size"),
     [
         [5, 10],
         [7, 10],
@@ -20,20 +20,20 @@
         [30, 250],
     ],
 )
-def test_ceiling_variable_density(max_clust_size, clust_size):
+def test_ceiling_variable_density(max_size, clust_size):
     """Test algorithmic guarantee of BisectingQMeans on variable density data
     
-    Verify that when `min_clust_size=2`, that the max per cluster membership is
+    Verify that when `min_size=2`, that the max per cluster membership is
     under (<, less than) what parameter `opt_cluster_size` is set to"""
 
     data = gen_variable_density_clusters(clust_size)
-    clust = BisectingQMeans(min_clust_size=1, opt_clust_size=max_clust_size,
+    clust = BisectingQMeans(min_size=1, opt_size=max_size,
                             init='random', n_init=50, algorithm='lloyd',
                             max_iter=8000, random_state=42)
     clust.fit(data)
 
     _, counts = np.unique(clust.labels_, return_counts=True)
-    assert np.max(counts) < max_clust_size
+    assert np.max(counts) < max_size
 
 @pytest.mark.parametrize(
     ("min_clust", "max_clust", "neighbors", "overlap"),
@@ -53,7 +53,7 @@ def test_max_clust_expansion(min_clust, max_clust, neighbors, overlap):
     [initial cluster size + (neighbors * overlap)]"""
     
     data = gen_variable_density_clusters()
-    clust = BisectingQMeans(min_clust_size=min_clust, opt_clust_size=max_clust,
+    clust = BisectingQMeans(min_size=min_clust, opt_size=max_clust,
                             init='random', n_init=50, algorithm='lloyd',
                             max_iter=8000, random_state=42)
     clust.fit(data)