在各种多边形上迭代时设置恒定的子地块大小

2024-05-15 12:44:35 发布

您现在位置:Python中文网/ 问答频道 /正文

所以我在那之前也问过类似的问题。但我没能按计划完成

我正在尝试(有些成功地)自动化每月更新工作流程,其中我们制作地图,突出显示各个分区的施工进度。我遇到的问题是,我无法保持子图(?)背景图(?)的恒定物理大小。基本上,我尝试在背景中使用航空图像,当程序在多边形(细分边界)上迭代时,每次生成的地图的缩放级别和范围都会发生变化,而背景中的航空图像总是占据一张纸上相同的物理尺寸

下面是一个用GIS制作的地图示例,我正试图用Python复制它

GIS Produced Map

背景图像被拉伸以适合页面,而不会扭曲地图的焦点,以便平移和缩放到细分边界的范围

这是一个我能够使用Python为相同的细分生成的地图

Python Produced Map

我尝试了边界框方法和setxlim/set_ylim方法,但这总是裁剪背景航空图像,我希望它每次都通过具有不同位置和大小的多个细分填充页面

下面是我用来创建此输出的代码:

import pandas as pd # tabular data in Python
import geopandas as gpd # extension to Pandas to work with geospatial data
import urllib.request, json # download from web
import os.path # work with the local file system
import matplotlib as mpl # plotting
from matplotlib import pyplot as plt # matplotlib convenience functions
from matplotlib.patches import Patch # build legend items
from matplotlib.lines import Line2D # build legend items
from matplotlib import figure # change size of figures
from matplotlib.font_manager import FontProperties
import matplotlib.patheffects as pe # add halos to text
import contextily as ctx # basemaps
import csv # work with .csv files
import operator # store lists alphabetically
import time # time tracking for processes
import warnings # ignore warnings

warnings.filterwarnings("ignore") # ignore warnings of deprecated functions
start_time = time.time() # Keep track of execution time
skipped_phases = [] # empty list of skipped phases to be printed at the end
maps_made = 0 # empty number of maps made
completed_phases = 0 # empty list for number of completed phases found
completed_phase_names = [] # empty list for the names any completed phases found
plt.rcParams['legend.title_fontsize'] = 'xx-small'

# Read differences table
fname = table_path + "diff_con.csv"
diff_con = pd.read_csv(fname)
print("Number of maps to be made for this update: " + str(len(diff_con)) + "\n")

# Build Legend Elements
legend_elements = [Patch(facecolor="green", edgecolor="black", alpha=0.5, 
                         label="COMPLETED LOTS"),
                   Patch(facecolor="red", edgecolor="black", alpha=0.5, 
                         label="LOTS UNDER CONSTRUCTION"),
                   Patch(facecolor="none", edgecolor="r",
                         label="SUBDIVISION BOUNDARY")]

for index, row in diff_con.iterrows():
    name = row["Subdivision"]
    print(name) # Print name of current selected phase
    
    lots_complete_select = lots_complete[(lots_complete["SUBDIVISION"] 
                                          == row["Subdivision"])] # Select completed lots w/ same phase name
    print("Completed Lots: " + str(len(lots_complete_select)))
    
    lots_uc_select = lots_uc[(lots_uc["SUBDIVISION"] 
                              == row["Subdivision"])] # Select incomplete lots w/ same phase name
    print("Lots Under Construction: " + str(len(lots_uc_select)))
    
    phase_select = phases[(phases["NAME"] == row["Subdivision"])] # Select phase w/ same name from loop
    
    lots_select = lots[(lots["SUBDIVISION"] == row["Subdivision"])] # Select all lots w/ same name as phase
    print("Lots in subdivision: " + str(len(lots_select)))
    
    minx, miny, maxx, maxy = (lots_select).total_bounds # Set zoom to lots for area of interest

    # Skip phases with zero lots found so error is not thrown and stops script
    lots_filter = row["Total Lots"]
    if lots_filter == 0:
        print("Subdivision phase skipped due to 0 lots found.\n")
        skipped_phases.append(row["Subdivision"]) # store phase name for later reference
        pass
    else:
        if row["Lots Completed"] == row["Proposed Lots"]: # Check to see if phase is completed
            print(name + " has been completely built out!")
            save_name = "COMPLETED_" + name # Change name of map file produced
            completed_phases += 1 # Update number of phases completed this update
            completed_phase_names.append(row["Subdivision"])
        else:
            save_name = name

        map_time = time.time()
        print("Building map...")
        
        # Map Building Time
        fig, ax = plt.subplots() # Create a figure with multiple plots within it
        #ax.set_aspect("equal")
        ax.axis("off") # Turn off axes

        margin = 0.0005 # Set margin to remove excess white space from final output
        ax.set_xlim(minx - margin, maxx + margin) # Apply x zoom level
        ax.set_ylim(miny - margin, maxy + margin) # Apply y zoom level
        # Unfortunately, setting zoom also crops the output, still working on a way around that

        lots_complete_select.plot(ax=ax, linewidth=0.5, color="green", 
                                  edgecolor="black", alpha=0.5) # Plot completed lots

        lots_uc_select.plot(ax=ax, linewidth=0.5, color="red", edgecolor="black", alpha=0.5) # Plot U.C. lots

        lots_select.plot(ax=ax, linewidth=0.5, color="none", edgecolor="black") # Plot lot lines
        
        # Phase Boundary
        phase_select.plot(ax=ax, linewidth=3, color="none", edgecolor="black") # Plot Active Subdivision Phase
        phase_select.plot(ax=ax, linewidth=1, color="none", edgecolor="red") # Plot Active Subdivision Phase
        phase_select.plot(ax=ax, linewidth=0.5, color="none", edgecolor="white") # Plot Active Subdivision Phase
        
        # Label address points for each figure
        phase_buffer = phase_select["geometry"].buffer(0) # Create a 'buffered' polygon of the selected phase
        df = gpd.GeoDataFrame(gpd.GeoSeries(phase_buffer), columns=["geometry"]) # Save as new dataframe
        addr_select = gpd.overlay(addresses, df, how="intersection") # Use overlay to select addr points in phase

        # Label each address point located in selected phase
        for x, y, label in zip(addr_select.geometry.x, addr_select.geometry.y, addr_select.ADDR_NUM):
            ax.annotate(label, xy=(x, y), xytext=(0.2, 0.2), color="black", 
                        path_effects=[pe.withStroke(linewidth=1, foreground="white")], ha="center", va="center", 
                        fontsize=2, fontweight ="bold", textcoords="offset points")

        ctx.add_basemap(ax, crs=lots.crs.to_string(), source=ctx.providers.Esri.WorldStreetMap, 
                        attribution="") # Add ESRI Streets Basemap
        
        ctx.add_basemap(ax, crs=lots.crs.to_string(), source=ctx.providers.Esri.WorldImagery, 
                        attribution="",alpha=0.5) # Add ESRI Imagery Basemap

        ax.set_title((name + " Residential Construction"), fontsize=8, fontweight ="bold") # Add title

        ax.legend(handles=legend_elements, prop={"size": 3}, title="Legend", framealpha=1, fancybox=False, 
                  edgecolor="black", loc="upper left") # Add legend

        filePath = map_path + save_name + ".jpg" # Build filepath for saving maps
        
        lims = plt.axis('equal')
        
        plt.savefig(filePath, edgecolor="black", dpi=300, bbox_inches="tight") # Save map
        
        print("Map saved.")
        print("Done building map.")
        print("--- %s Seconds ---" % (time.time() - map_time) + "\n")

        plt.close("all") # Close figures to save memory

        maps_made += 1 # Update number of maps made
        
        print("Number of maps remaining: " + str(len(diff_con) - maps_made) + "\n") # Print how many maps are left

# Print update statistics
print("Done configuring maps.")
print("--- %s Minutes ---\n" % ((time.time() - start_time)/60))
print("Number of maps made for this month's update: " + str(maps_made) + "\n")
print("List of subdivision phases skipped due to 0 lots found:")
print(skipped_phases)
print("\n" + str(completed_phases) + " phase(s) completed this update! See names below:")
print(completed_phase_names)

是否有任何方法可以尝试使背景航空图像在每一细分阶段的范围和缩放变化时生成的纸张jpg文件上具有恒定的物理大小

多谢各位


Tags: oftonameimportfortimeaxselect