3

I have address points connected to a network of roads that I have built using GRASS v.net.connect. I would like to cluster addresses into fixed size groups along the network.

enter image description here

So far, I have written an algorithm that:

  • Uses the southern most address point as the "start point" to begin.
  • Finds the shortest path between the start point and all other address points using the native:shortestpathpointtopoint algorithm
  • Associates these two points (the start point and nearest end point) to a cluster
  • Set the previous nearest end point as the next start point
  • Find the nearest address point using the shortest path from this start point to get the next address point to add to the cluster
  • repeat this process until all address points have been clustered incrementing the cluster number every N points

I could not tell you if this method works because the algorithm I have written is estimated to take about 55 hours using 416 address points. I have included the code below and would like any suggestions you might have on how I can improve its speed. Or if there is a better solution to doing what I am trying to achieve.

def get_layer_by_name_quick(namepattern):
    layers = []
    for lyr in QgsProject.instance().mapLayers().values():
        if namepattern in str(lyr.name()):
            layers.append(lyr)
    if len(layers) == 0:
        print(f"No Layer named  : {namepattern}")
        return None
    elif len(layers) > 1:
        print(f"multiple layers returned: \nnamepattern: {namepattern}, \nlayers: {layers}\n")
        return layers[0]
    else:
        return layers[0]

class ClusterAddressesV2(object):
    
    def __init__(self, addresses, network, clustersize = 30):
        
        self.addresses = addresses
        self.network = network
        self.clustersize = clustersize
        self.cluster = 0
        self.clusters = {}
        self.collected_points = []
        self.startPoint, self.start_point = self.get_min_y_point()
        self.all_points = {feat.id() : feat.geometry().asPoint() for feat in self.addresses.getFeatures()}
        
    def get_min_y_point(self):
        pnts = [(addr.id(), addr.geometry().asPoint()) for addr in self.addresses.getFeatures() if addr.id() not in self.collected_points]
        xys = [(pnt[1].x(), pnt[1].y()) for pnt in pnts]
        xysorted = sorted(xys, key = lambda item: item[1]) # sort by y
        first_point = xysorted[0]
        pnt0 = [(p[0], p[1]) for p in pnts if p[1].y() == first_point[1]] [0]
        pid, point = pnt0
        return pid, point
    
    def get_all_points(self):            
        the_points = {feat.id() : feat.geometry().asPoint() for feat in self.addresses.getFeatures() if feat.id() not in self.collected_points}
        return the_points
    
    def get_shortest_paths(self):
        
                  
        the_points = self.get_all_points() # get all points        
        shortest_paths = []
        
        for k, v in the_points.items():              
            if v != self.start_point:                 
                end_point = v
                endPoint = k            
            else:                
                # k is startPoint
                self.startPoint = k
                endPoint = k
                end_point = v
                
            
            if self.startPoint not in self.collected_points:
                self.collected_points.append(self.startPoint)
            
            
            feedback = QgsProcessingFeedback()    
            parameters = {'INPUT': self.network,
                'STRATEGY': 0,
                'DIRECTION_FIELD': '',
                'VALUE_FORWARD': '',
                'VALUE_BACKWARD': '',
                'VALUE_BOTH': '',
                'DEFAULT_DIRECTION': 2,
                'SPEED_FIELD': '',
                'DEFAULT_SPEED': 1,
                'TOLERANCE': 15, #tolerance
                'START_POINT': self.start_point,
                'END_POINT': end_point,
                'OUTPUT': 'memory:shortestpath'}
            
            try:
                branch = processing.run('native:shortestpathpointtopoint', parameters, feedback = feedback)
                shortest_paths.append((self.startPoint, endPoint, branch['TRAVEL_COST']))
                
            except Exception as e:
                print(f'Exception: {e} startPoint : {self.startPoint}, endPoint: {endPoint}')
                
        # after running through all the points
        s = sorted( shortest_paths, key = lambda item: item[2])
        
        s1 = s[0][1] #endPoint
        s0 = s[0][0] #startPoint
        
        if s1 not in list(self.clusters.keys()):
            self.clusters[s1] = self.cluster
        if s0 not in list(self.clusters.keys()):
            self.clusters[s0] = self.cluster
        if s1 not in self.collected_points:
            self.collected_points.append(s1)
        if s0 not in self.collected_points:
            self.collected_points.append(s0)  
            
        self.startPoint = s1
        self.start_point = self.all_points[self.startPoint]
        
        if len(list(self.clusters.keys())) % self.clustersize == 0:
            self.cluster += 1
            
        return
    
    def write_attributes(self):
        all_features = [f for f in self.addresses.getFeatures()] 
        layer_provider = self.addresses.dataProvider()
        self.addresses.startEditing()
        fields = self.addresses.fields()
        clusterIdx = fields.indexFromName('Cluster')
        for f in all_features:
            id = f.id()
            if id in list(self.clusters.keys()):
                cluster = self.clusters[id]
                #print(cluster)
                attr_value={clusterIdx : cluster}
                layer_provider.changeAttributeValues({id:attr_value})
        self.addresses.commitChanges()
        return
        
                
    def run_algo(self):
        
        all_features = [f for f in self.addresses.getFeatures()]         
        while len(list(self.clusters.keys())) < len(all_features):             
            self.get_shortest_paths()
            if len(list(self.clusters.keys())) >= len(all_features):                
                break                
        self.write_attributes()        


addresses = get_layer_by_name_quick('repr_addresses')
network = get_layer_by_name_quick('Network')

cav = ClusterAddressesV2(addresses, network, 25)

cav.run_algo()

 

1 Answer 1

4

Try this algorithm, you need some preparation to imitate stream network. Flip edges where travel to discharge point for FROM node is less than one for TO node:

enter image description here

Keep 1st shortest edge coming out from each node, delete others:

enter image description here

Now your network can be presented as streams network. Picture shows (flow) accumulation of address points in such network:

enter image description here

Find upstream critical edges where flow equal or jumps above 10. If count of nodes with weight (address) upstream is greater than 10, unselect nearest n = N -10 upstream. Memorize each group. Picture shows flow accumulation in edges.

enter image description here

Remove assigned nodes from graph and run flow accumulation again. Repeat until remaining total total of nodes with weight (address point) is less than 10.

It might produce ugly results if multiple addresses 'discharge' into same node:

enter image description here

My script took 25 minutes to group 16447 address points.

1
  • extremely creative
    – unbutu
    Commented Apr 18, 2023 at 0:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.