From d5893cdcdc5282ad7aeca3fc446c2611fe0753a6 Mon Sep 17 00:00:00 2001 From: Alireza Adli Date: Mon, 26 Feb 2024 12:06:49 -0500 Subject: [PATCH] Add spatial joins with NRCan and GeoIndex into the main workflow --- mtl_buildings_workflow.py | 78 ++++++++++++++++++++++++++++++++++----- 1 file changed, 68 insertions(+), 10 deletions(-) diff --git a/mtl_buildings_workflow.py b/mtl_buildings_workflow.py index bcee528..6a85b4e 100644 --- a/mtl_buildings_workflow.py +++ b/mtl_buildings_workflow.py @@ -124,9 +124,9 @@ for each in range(num_layers): processing.run("native:extractbyexpression", splitting_params) -# Creating spatian index for each new layer - create_spatial_indext_params = {'INPUT': output_layer_path, 'OUTPUT': 'Output'} - processing.run("native:createspatialindex", create_spatial_indext_params) +# Creating spatial index for each new layer + create_index_params = {'INPUT': output_layer_path, 'OUTPUT': 'Output'} + processing.run("native:createspatialindex", create_index_params) # Putting the remaining features in one last layer (So, overall we're going to have [num_layers + 1] layers). remaining_features = num_layers @@ -139,8 +139,8 @@ last_splitting_params = {'INPUT': clipped_nrcan_layer, processing.run("native:extractbyexpression", last_splitting_params) # Create spatial index for the last partition (layer) -create_spatial_indext_params = {'INPUT': output_layer_path, 'OUTPUT': 'Output'} -processing.run("native:createspatialindex", create_spatial_indext_params) +create_index_params = {'INPUT': output_layer_path, 'OUTPUT': 'Output'} +processing.run("native:createspatialindex", create_index_params) # Now we clip the property assessment with each of the partitions (layers) # from the last process. The process outputs the pairwise clipped version of @@ -159,9 +159,9 @@ for each in range(num_clipped_layers): pairwise_clipping_params = {'INPUT': property_assessment_layer, 'OVERLAY': clipping_layer, 'FILTER_EXPRESSION': '', 'FILTER_EXTENT': None, 'OUTPUT': clipped_layer} processing.run("native:clip", pairwise_clipping_params) - # Creating spatian index for each new layer - create_spatial_indext_params = {'INPUT': clipped_layer, 'OUTPUT': 'Output'} - processing.run("native:createspatialindex", create_spatial_indext_params) + # Creating spatial index for each new layer + create_index_params = {'INPUT': clipped_layer, 'OUTPUT': 'Output'} + processing.run("native:createspatialindex", create_index_params) print("Pairwise Clipping is completed.") @@ -171,12 +171,70 @@ print("Pairwise Clipping is completed.") # First we extract all the .shp files (clipped layers) from each folder pairwise_clipped_layers_path = find_shp_files(clipped_layers_folder) -pairwise_clipped_property_assessment = 'C:/Users/a_adli/PycharmProjects/hydroquebec_archetype_gispy/data/output_data/pairwise_clipped_property_assessment/pairwise_clipped.shp' +pairwise_clipped_property_assessment_layer = 'C:/Users/a_adli/PycharmProjects/hydroquebec_archetype_gispy/data/output_data/pairwise_clipped_property_assessment/pairwise_clipped.shp' merging_params = {'LAYERS': pairwise_clipped_layers_path, 'CRS': None, - 'OUTPUT': pairwise_clipped_property_assessment} + 'OUTPUT': pairwise_clipped_property_assessment_layer} processing.run("native:mergevectorlayers", merging_params) +# Create spatial index for the pairwise clipped Property Assessment +create_index_params = {'INPUT': pairwise_clipped_property_assessment_layer, 'OUTPUT': 'Output'} +processing.run("native:createspatialindex", create_index_params) + +pairwise_clipped_property_assessment, pairwise_clipped_property_assessment_name = \ + load_layer(pairwise_clipped_property_assessment_layer, 'Pairwise Clipped Property Assessment Layer') +print(f'{pairwise_clipped_property_assessment_name} data count: {pairwise_clipped_property_assessment.featureCount()}') + +# ----- +# Delete path and layer fields (This is going to be added after testing the whole program). +# ----- + +# In QGIS spatial join is named Join Attributes by Location +# Spatial join Pairwise Clipped Property Assessment with Clipped NRCan +property_assessment_nrcan_layer = 'C:/Users/a_adli/PycharmProjects/hydroquebec_archetype_gispy/data/output_data/property_join_nrcan/property_join_nrcan.shp' + +spatial_join_nrcan_params = {'INPUT': pairwise_clipped_property_assessment_layer, + 'PREDICATE': [0], + 'JOIN': clipped_nrcan_layer, + 'JOIN_FIELDS': [], + 'METHOD': 0, + 'DISCARD_NONMATCHING': False, + 'PREFIX': '', + 'OUTPUT': property_assessment_nrcan_layer} + +feedback = QgsProcessingFeedback() +processing.run('native:joinattributesbylocation', spatial_join_nrcan_params, feedback=feedback) + +create_index_params = {'INPUT': property_assessment_nrcan_layer, 'OUTPUT': 'Output'} +processing.run("native:createspatialindex", create_index_params) + +property_assessment_nrcan, property_assessment_nrcan_name = \ + load_layer(property_assessment_nrcan_layer, 'Property Assessment Layer Joined with NRCan') +print(f'{property_assessment_nrcan_name} data count: {property_assessment_nrcan.featureCount()}') + +# Spatial join Pairwise Clipped (joined with NRCan) Property Assessment layer with Clipped GeoIndex layer +property_assessment_nrcan_geo_layer = 'C:/Users/a_adli/PycharmProjects/hydroquebec_archetype_gispy/data/output_data/propertynrcan_join_geo/propertynrcan_join_geo.shp' + + +spatial_join_nrcan_params = {'INPUT': property_assessment_nrcan_layer, + 'PREDICATE': [0], + 'JOIN': clipped_geoindex_layer, + 'JOIN_FIELDS': [], + 'METHOD': 0, + 'DISCARD_NONMATCHING': False, + 'PREFIX': '', + 'OUTPUT': property_assessment_nrcan_geo_layer} + +feedback = QgsProcessingFeedback() +processing.run('native:joinattributesbylocation', spatial_join_nrcan_params, feedback=feedback) + +create_index_params = {'INPUT': property_assessment_nrcan_geo_layer, 'OUTPUT': 'Output'} +processing.run("native:createspatialindex", create_index_params) + +property_assessment_nrcan_geo, property_assessment_nrcan_geo_name = \ + load_layer(property_assessment_nrcan_geo_layer, 'Property Assessment Layer Joined with NRCan') +print(f'{property_assessment_nrcan_geo_name} data count: {property_assessment_nrcan_geo.featureCount()}') + qgs.exitQgis() \ No newline at end of file