farm_modeling.py

import geopandas as gpd
import numpy as np
import matplotlib.pyplot as plt
from shapely.geometry import Point, LineString, Polygon
from shapely.ops import transform
import pyproj
import networkx as nx
import os
import time
import json
from functools import partial
import math
import utm
from scipy.spatial import cKDTree

class CommercialFarmPathPlanner:
    """商业化农业路径规划系统 - 完整坐标转换版"""
    
    def __init__(self, geojson_file, auto_detect_utm=True, utm_zone=None):
        """
        初始化路径规划系统
        
        :param geojson_file: 统一GeoJSON文件路径
        :param auto_detect_utm: 是否自动检测UTM分区
        :param utm_zone: 手动指定UTM分区 (如45)
        """
        self.geojson_file = geojson_file
        self.utm_zone = utm_zone
        self.auto_detect_utm = auto_detect_utm
        self.paths = {}
        self.transport_routes = {}
        
        # 加载数据
        self.load_data()
        
        # 设置坐标系转换
        self.setup_coordinate_system()
        
        # 创建空间索引
        self.create_spatial_index()
        
        # 提取要素
        self.extract_features()
    
    def load_data(self):
        """加载GeoJSON数据"""
        print(f"加载农场数据: {self.geojson_file}")
        start_time = time.time()
        
        # 加载GeoJSON
        self.gdf = gpd.read_file(self.geojson_file)
        
        # 验证坐标系
        if self.gdf.crs is None:
            self.gdf.crs = "EPSG:4326"  # 默认WGS84
            print("警告: 未检测到坐标系，已设置为WGS84 (EPSG:4326)")
        
        # 计算加载时间
        load_time = time.time() - start_time
        print(f"数据加载完成! 耗时: {load_time:.2f}s")
        print(f"总要素数量: {len(self.gdf)}")
        print(f"坐标系: {self.gdf.crs}")
    
    def setup_coordinate_system(self):
        """设置坐标系转换系统"""
        # 确定UTM分区
        if self.auto_detect_utm and not self.utm_zone:
            self.determine_utm_zone()
        
        # 创建转换器
        self.wgs84_crs = pyproj.CRS("EPSG:4326")
        self.utm_crs = pyproj.CRS(f"EPSG:326{self.utm_zone}")  # 北纬UTM
        
        # 创建转换函数
        self.project_to_utm = partial(
            pyproj.transform,
            pyproj.Proj(self.wgs84_crs),
            pyproj.Proj(self.utm_crs)
        )
        
        self.project_to_wgs84 = partial(
            pyproj.transform,
            pyproj.Proj(self.utm_crs),
            pyproj.Proj(self.wgs84_crs)
        )
        
        print(f"坐标系转换器已创建: UTM 分区 {self.utm_zone}N (EPSG:326{self.utm_zone})")
    
    def determine_utm_zone(self):
        """自动确定UTM分区"""
        # 计算农场中心点
        centroid = self.gdf.unary_union.centroid
        lon, lat = centroid.x, centroid.y
        
        # 计算UTM分区
        _, _, zone_number, zone_letter = utm.from_latlon(lat, lon)
        
        # 验证是否在北半球
        if zone_letter < 'N':
            raise ValueError("农场位于南半球，需要特殊处理")
        
        self.utm_zone = zone_number
        print(f"自动检测UTM分区: {self.utm_zone}N")
    
    def create_spatial_index(self):
        """创建空间索引"""
        print("创建空间索引...")
        
        # 创建UTM投影数据集
        self.gdf_utm = self.gdf.copy()
        self.gdf_utm['geometry'] = self.gdf_utm['geometry'].apply(
            lambda geom: transform(self.project_to_utm, geom) if geom else geom
        )
        
        # 创建KDTree索引
        points = self.gdf_utm.geometry.centroid
        self.coords = np.array([[p.x, p.y] for p in points if p])
        self.indices = np.arange(len(self.coords))
        
        self.spatial_index = cKDTree(self.coords)
        print("空间索引创建完成")
    
    def extract_features(self):
        """提取农场要素"""
        print("提取农场要素...")
        
        # 从UTM坐标系提取
        self.boundary = self.gdf_utm[self.gdf_utm['feature_type'] == 'boundary'].iloc[0]
        self.fields = self.gdf_utm[self.gdf_utm['feature_type'] == 'field']
        self.roads = self.gdf_utm[self.gdf_utm['feature_type'] == 'road']
        self.entries = self.gdf_utm[self.gdf_utm['feature_type'] == 'entry']
        self.warehouses = self.gdf_utm[self.gdf_utm['feature_type'] == 'warehouse']
        
        # 创建田块索引
        self.field_dict = {row['field_id']: row for _, row in self.fields.iterrows()}
        
        print(f"提取完成: {len(self.fields)}田块, {len(self.roads)}道路, {len(self.entries)}入口点")
    
    def generate_field_path(self, field_id, machine_width=3.0):
        """
        为田块生成作业路径 (UTM坐标系)
        
        :param field_id: 田块ID
        :param machine_width: 农机宽度 (米)
        :return: 作业路径GeoDataFrame (WGS84坐标系)
        """
        print(f"为田块 {field_id} 生成作业路径...")
        
        # 获取田块和入口点 (UTM坐标)
        field = self.field_dict.get(field_id)
        if field is None:
            raise ValueError(f"未找到ID为{field_id}的田块")
            
        field_entries = self.entries[self.entries['field_id'] == field_id]
        if field_entries.empty:
            raise ValueError(f"田块{field_id}缺少入口点")
        
        # 确定主入口点
        primary_entries = field_entries[field_entries['entry_type'] == 'primary']
        if not primary_entries.empty:
            entry = primary_entries.iloc[0]
        else:
            # 如果没有主入口点，使用第一个入口点
            entry = field_entries.iloc[0]
        
        # 计算最优方向
        minx, miny, maxx, maxy = field.geometry.bounds
        if (maxx - minx) > (maxy - miny):
            direction = 'east_west'
        else:
            direction = 'north_south'
        
        # 生成平行作业路径 (使用米制单位)
        work_lines = self.generate_parallel_lines(
            field.geometry, 
            direction, 
            machine_width
        )
        
        # 连接路径
        path_lines = self.connect_work_lines(work_lines, entry.geometry)
        
        # 创建路径数据框 (UTM坐标系)
        path_data = {
            'field_id': [field_id] * len(path_lines),
            'path_type': ['operation'] * len(path_lines),
            'sequence': range(len(path_lines)),
            'machine_width': [machine_width] * len(path_lines),
            'geometry': path_lines
        }
        
        path_gdf_utm = gpd.GeoDataFrame(path_data, crs=self.utm_crs)
        
        # 转换为WGS84坐标系
        path_gdf = path_gdf_utm.copy()
        path_gdf['geometry'] = path_gdf['geometry'].apply(
            lambda geom: transform(self.project_to_wgs84, geom)
        )
        path_gdf.crs = self.wgs84_crs
        
        # 保存路径
        self.paths[field_id] = path_gdf
        
        # 计算路径统计
        total_length = sum(line.length for line in path_lines)  # 米
        num_segments = len(path_lines)
        
        print(f"路径生成完成: {num_segments}段, 总长{total_length:.1f}米")
        
        return path_gdf
    
    def generate_parallel_lines(self, polygon, direction, spacing_m):
        """
        生成平行作业线 (UTM坐标系)
        
        :param polygon: 田块多边形 (UTM)
        :param direction: 方向 ('east_west' 或 'north_south')
        :param spacing_m: 间距 (米)
        :return: 作业线列表 (UTM)
        """
        minx, miny, maxx, maxy = polygon.bounds
        buffer = 50  # 50米缓冲区确保覆盖
        
        lines = []
        
        if direction == 'east_west':
            # 东西方向
            y = miny + spacing_m / 2
            while y < maxy:
                line = LineString([(minx - buffer, y), (maxx + buffer, y)])
                intersection = polygon.intersection(line)
                
                if not intersection.is_empty:
                    if intersection.geom_type == 'MultiLineString':
                        lines.extend(list(intersection.geoms))
                    else:
                        lines.append(intersection)
                y += spacing_m
        else:
            # 南北方向
            x = minx + spacing_m / 2
            while x < maxx:
                line = LineString([(x, miny - buffer), (x, maxy + buffer)])
                intersection = polygon.intersection(line)
                
                if not intersection.is_empty:
                    if intersection.geom_type == 'MultiLineString':
                        lines.extend(list(intersection.geoms))
                    else:
                        lines.append(intersection)
                x += spacing_m
        
        return lines
    
    def connect_work_lines(self, lines, entry_point):
        """
        连接作业线形成连续路径 (UTM坐标系)
        
        :param lines: 作业线列表
        :param entry_point: 入口点
        :return: 连接后的路径列表
        """
        if not lines:
            return []
        
        # 找到最近的作业线
        start_line = min(lines, key=lambda line: entry_point.distance(line))
        
        # 从入口点到第一条作业线
        path = [LineString([entry_point, start_line.interpolate(0.1)])]
        
        # 连接作业线
        remaining = [line for line in lines if line != start_line]
        current_line = start_line
        
        while remaining:
            # 找到最近的下一段
            next_line = min(remaining, key=lambda line: current_line.distance(line))
            
            # 添加连接线
            end_point = current_line.interpolate(1.0, normalized=True)
            start_point = next_line.interpolate(0.1, normalized=True)
            path.append(LineString([end_point, start_point]))
            
            # 添加作业线
            path.append(next_line)
            
            # 更新状态
            current_line = next_line
            remaining.remove(next_line)
        
        return path
    
    def generate_transport_network(self):
        """
        生成农场运输网络 (UTM坐标系)
        
        :return: 运输网络GeoDataFrame (WGS84坐标系)
        """
        print("生成运输网络...")
        
        # 创建道路网络图
        G = nx.Graph()
        
        # 添加道路节点
        node_counter = 0
        node_positions = {}
        road_segments = []
        
        for _, road in self.roads.iterrows():
            line = road.geometry
            coords = list(line.coords)
            
            # 添加节点
            for i, coord in enumerate(coords):
                node_id = f"R_{road.name}_{i}"
                G.add_node(node_id, pos=coord, type="road")
                node_positions[coord] = node_id
                
                # 连接连续节点
                if i > 0:
                    prev_coord = coords[i-1]
                    prev_node = node_positions[prev_coord]
                    distance = Point(coord).distance(Point(prev_coord))
                    G.add_edge(prev_node, node_id, length=distance, road_id=road['road_id'])
        
        # 添加入口点
        for _, entry in self.entries.iterrows():
            entry_point = entry.geometry.coords[0]
            closest_node = self.find_closest_road_node(entry_point, G)
            
            if closest_node:
                distance = Point(entry_point).distance(Point(G.nodes[closest_node]['pos']))
                G.add_node(f"E_{entry['entry_id']}", pos=entry_point, type="entry", field_id=entry['field_id'])
                G.add_edge(f"E_{entry['entry_id']}", closest_node, length=distance, type="access")
        
        # 添加仓库
        for _, warehouse in self.warehouses.iterrows():
            wh_point = warehouse.geometry.coords[0]
            closest_node = self.find_closest_road_node(wh_point, G)
            
            if closest_node:
                distance = Point(wh_point).distance(Point(G.nodes[closest_node]['pos']))
                G.add_node(f"W_{warehouse['warehouse_id']}", pos=wh_point, type="warehouse", name=warehouse['name'])
                G.add_edge(f"W_{warehouse['warehouse_id']}", closest_node, length=distance, type="access")
        
        # 创建网络GeoDataFrame
        edges = []
        for u, v, data in G.edges(data=True):
            line = LineString([G.nodes[u]['pos'], G.nodes[v]['pos']])
            edges.append({
                'from': u,
                'to': v,
                'length': data['length'],
                'type': data.get('type', 'road'),
                'geometry': line
            })
        
        network_gdf_utm = gpd.GeoDataFrame(edges, crs=self.utm_crs)
        
        # 转换为WGS84坐标系
        network_gdf = network_gdf_utm.copy()
        network_gdf['geometry'] = network_gdf['geometry'].apply(
            lambda geom: transform(self.project_to_wgs84, geom)
        )
        network_gdf.crs = self.wgs84_crs
        
        print(f"运输网络生成完成: {len(edges)}条路段")
        return network_gdf
    
    def find_closest_road_node(self, point, graph, max_distance=300):
        """
        在道路网络中查找最近节点
        
        :param point: 查询点 (UTM坐标)
        :param graph: 道路网络图
        :param max_distance: 最大搜索距离 (米)
        :return: 最近节点ID
        """
        min_distance = float('inf')
        closest_node = None
        
        for node, data in graph.nodes(data=True):
            if data['type'] == 'road':
                node_point = Point(data['pos'])
                distance = node_point.distance(Point(point))
                
                if distance < min_distance and distance < max_distance:
                    min_distance = distance
                    closest_node = node
        
        return closest_node
    
    def plan_transport_route(self, from_field, to_warehouse=None):
        """
        规划从田块到仓库的运输路径
        
        :param from_field: 起始田块ID
        :param to_warehouse: 目标仓库ID (可选)
        :return: 运输路径GeoDataFrame (WGS84坐标系)
        """
        print(f"规划运输路径: 田块 {from_field} → 仓库 {to_warehouse or '最近仓库'}")
        
        # 获取田块入口点
        field_entries = self.entries[self.entries['field_id'] == from_field]
        if field_entries.empty:
            raise ValueError(f"田块{from_field}缺少入口点")
        
        entry_point = field_entries.iloc[0].geometry.coords[0]
        
        # 确定目标仓库
        if to_warehouse:
            warehouse = self.warehouses[self.warehouses['warehouse_id'] == to_warehouse].iloc[0]
        else:
            warehouse = self.find_closest_warehouse(entry_point)
        
        wh_point = warehouse.geometry.coords[0]
        
        # 使用A*算法查找最短路径
        route = self.find_shortest_path(entry_point, wh_point)
        
        if not route:
            raise RuntimeError("无法找到有效路径")
        
        # 创建路径数据框 (UTM坐标系)
        path_lines = []
        for i in range(len(route) - 1):
            path_lines.append(LineString([route[i], route[i+1]]))
        
        route_gdf_utm = gpd.GeoDataFrame({
            'from': from_field,
            'to': warehouse['warehouse_id'],
            'distance': [line.length for line in path_lines],
            'geometry': path_lines
        }, crs=self.utm_crs)
        
        # 转换为WGS84坐标系
        route_gdf = route_gdf_utm.copy()
        route_gdf['geometry'] = route_gdf['geometry'].apply(
            lambda geom: transform(self.project_to_wgs84, geom)
        )
        route_gdf.crs = self.wgs84_crs
        
        # 保存路径
        self.transport_routes[(from_field, warehouse['warehouse_id'])] = route_gdf
        
        total_distance = sum(line.length for line in path_lines)
        print(f"运输路径规划完成: 总长 {total_distance:.1f}米")
        
        return route_gdf
    
    def find_closest_warehouse(self, point):
        """查找最近的仓库"""
        min_distance = float('inf')
        closest_warehouse = None
        
        for _, warehouse in self.warehouses.iterrows():
            distance = Point(point).distance(warehouse.geometry)
            if distance < min_distance:
                min_distance = distance
                closest_warehouse = warehouse
        
        return closest_warehouse
    
    def find_shortest_path(self, start_point, end_point):
        """
        使用A*算法查找最短路径 (UTM坐标系)
        
        :param start_point: 起点 (UTM坐标)
        :param end_point: 终点 (UTM坐标)
        :return: 路径点列表
        """
        # 创建简化道路网络
        G = nx.Graph()
        
        # 添加关键节点
        G.add_node("start", pos=start_point)
        G.add_node("end", pos=end_point)
        
        # 添加道路节点
        for _, road in self.roads.iterrows():
            coords = list(road.geometry.coords)
            for i, coord in enumerate(coords):
                node_id = f"{road.road_id}_{i}"
                G.add_node(node_id, pos=coord)
                
                # 连接到起点/终点
                if i == 0:
                    distance = Point(start_point).distance(Point(coord))
                    G.add_edge("start", node_id, weight=distance)
                
                if i == len(coords) - 1:
                    distance = Point(end_point).distance(Point(coord))
                    G.add_edge("end", node_id, weight=distance)
                
                # 连接连续节点
                if i > 0:
                    prev_coord = coords[i-1]
                    prev_node = f"{road.road_id}_{i-1}"
                    distance = Point(coord).distance(Point(prev_coord))
                    G.add_edge(prev_node, node_id, weight=distance)
        
        # 添加启发式函数 (直线距离)
        def heuristic(u, v):
            pos_u = G.nodes[u]['pos']
            pos_v = G.nodes[v]['pos']
            return Point(pos_u).distance(Point(pos_v))
        
        # 执行A*算法
        try:
            path = nx.astar_path(G, "start", "end", heuristic=heuristic, weight="weight")
            return [G.nodes[node]['pos'] for node in path]
        except nx.NetworkXNoPath:
            return None
    
    def visualize(self, show_paths=True, save_path=None):
        """可视化农场地图"""
        print("生成农场可视化...")
        
        # 创建绘图
        fig, ax = plt.subplots(figsize=(15, 12))
        
        # 绘制农场边界 (WGS84)
        boundary_wgs = transform(self.project_to_wgs84, self.boundary.geometry)
        gpd.GeoSeries([boundary_wgs]).plot(ax=ax, color='none', edgecolor='gray', linewidth=2)
        
        # 绘制田块 (WGS84)
        fields_wgs = self.fields.copy()
        fields_wgs['geometry'] = fields_wgs['geometry'].apply(
            lambda geom: transform(self.project_to_wgs84, geom)
        )
        fields_wgs.plot(ax=ax, column='crop_type', legend=True, alpha=0.5)
        
        # 绘制道路 (WGS84)
        roads_wgs = self.roads.copy()
        roads_wgs['geometry'] = roads_wgs['geometry'].apply(
            lambda geom: transform(self.project_to_wgs84, geom)
        )
        roads_wgs.plot(ax=ax, color='gray', linewidth=2)
        
        # 绘制入口点 (WGS84)
        entries_wgs = self.entries.copy()
        entries_wgs['geometry'] = entries_wgs['geometry'].apply(
            lambda geom: transform(self.project_to_wgs84, geom)
        )
        entries_wgs.plot(ax=ax, color='blue', markersize=50)
        
        # 绘制仓库 (WGS84)
        warehouses_wgs = self.warehouses.copy()
        warehouses_wgs['geometry'] = warehouses_wgs['geometry'].apply(
            lambda geom: transform(self.project_to_wgs84, geom)
        )
        warehouses_wgs.plot(ax=ax, color='red', markersize=100)
        
        # 绘制路径
        if show_paths:
            for field_id, path_gdf in self.paths.items():
                path_gdf.plot(ax=ax, color='green', linewidth=1.5)
            
            for _, route_gdf in self.transport_routes.items():
                route_gdf.plot(ax=ax, color='purple', linewidth=2.5)
        
        # 添加标题和标签
        plt.title(f"新疆华兴农场 - 路径规划系统 (UTM Zone {self.utm_zone}N)", fontsize=16)
        plt.xlabel("经度", fontsize=12)
        plt.ylabel("纬度", fontsize=12)
        
        # 添加比例尺
        self.add_scale_bar(ax)
        
        # 添加图例
        legend_elements = [
            plt.Line2D([0], [0], color='green', lw=1.5, label='作业路径'),
            plt.Line2D([0], [0], color='purple', lw=2.5, label='运输路径'),
            plt.Line2D([0], [0], marker='o', color='w', markerfacecolor='blue', markersize=8, label='田块入口'),
            plt.Line2D([0], [0], marker='o', color='w', markerfacecolor='red', markersize=10, label='仓库')
        ]
        ax.legend(handles=legend_elements, loc='upper right')
        
        # 保存或显示
        if save_path:
            plt.savefig(save_path, dpi=300, bbox_inches='tight')
            print(f"地图已保存至: {save_path}")
        else:
            plt.show()
    
    def add_scale_bar(self, ax):
        """添加比例尺 (UTM坐标系)"""
        # 获取当前坐标范围
        x_min, x_max = ax.get_xlim()
        y_min, y_max = ax.get_ylim()
        
        # 计算比例尺位置 (左下角)
        bar_x = x_min + 0.05 * (x_max - x_min)
        bar_y = y_min + 0.05 * (y_max - y_min)
        
        # 计算比例尺长度 (对应1公里)
        # 使用UTM坐标系计算真实距离
        start_point = Point(bar_x, bar_y)
        end_point = Point(bar_x + 0.01, bar_y)  # 初始值
        
        # 转换为UTM计算距离
        start_utm = transform(self.project_to_utm, start_point)
        end_utm = transform(self.project_to_utm, end_point)
        current_distance = start_utm.distance(end_utm)
        
        # 调整到精确1公里
        scale_factor = 1000 / current_distance
        end_point = Point(bar_x + 0.01 * scale_factor, bar_y)
        
        # 绘制比例尺
        ax.plot([bar_x, end_point.x], [bar_y, bar_y], color='black', linewidth=2)
        ax.text((bar_x + end_point.x)/2, bar_y - 0.001, '1 km', 
                ha='center', fontsize=10)
    
    def export_all_paths(self, output_dir="path_output"):
        """导出所有路径"""
        os.makedirs(output_dir, exist_ok=True)
        print(f"导出路径到目录: {output_dir}")
        
        # 导出田块作业路径
        for field_id, path_gdf in self.paths.items():
            output_file = os.path.join(output_dir, f"field_{field_id}_path.geojson")
            path_gdf.to_file(output_file, driver="GeoJSON")
            print(f"  田块 {field_id} 路径: {output_file}")
        
        # 导出运输路径
        for (field_id, wh_id), route_gdf in self.transport_routes.items():
            output_file = os.path.join(output_dir, f"transport_{field_id}_to_{wh_id}.geojson")
            route_gdf.to_file(output_file, driver="GeoJSON")
            print(f"  运输路径 {field_id}→{wh_id}: {output_file}")
        
        # 导出元数据
        meta = {
            "export_time": time.strftime("%Y-%m-%d %H:%M:%S"),
            "utm_zone": self.utm_zone,
            "total_fields": len(self.paths),
            "total_routes": len(self.transport_routes)
        }
        with open(os.path.join(output_dir, "metadata.json"), "w") as f:
            json.dump(meta, f, indent=2)
    
    def generate_report(self):
        """生成农场运营报告"""
        report = {
            "generated_time": time.strftime("%Y-%m-%d %H:%M:%S"),
            "farm_name": "新疆华兴农场",
            "utm_zone": self.utm_zone,
            "total_fields": len(self.fields),
            "fields": [],
            "warehouses": []
        }
        
        # 田块信息
        for field_id, field in self.field_dict.items():
            field_info = {
                "field_id": field_id,
                "name": field.get("name", ""),
                "crop_type": field.get("crop_type", ""),
                "area_mu": field.get("area_mu", 0),
                "path_generated": field_id in self.paths
            }
            
            if field_id in self.paths:
                path_gdf = self.paths[field_id]
                total_length = sum(path_gdf.geometry.length) * 111000  # 近似转换为米
                field_info["path_length"] = total_length
                field_info["path_segments"] = len(path_gdf)
            
            report["fields"].append(field_info)
        
        # 仓库信息
        for _, warehouse in self.warehouses.iterrows():
            wh_info = {
                "warehouse_id": warehouse["warehouse_id"],
                "name": warehouse["name"],
                "capacity": warehouse.get("capacity", 0),
                "connected_fields": []
            }
            
            # 查找连接的田块
            for (field_id, wh_id), _ in self.transport_routes.items():
                if wh_id == warehouse["warehouse_id"]:
                    wh_info["connected_fields"].append(field_id)
            
            report["warehouses"].append(wh_info)
        
        return report


# ======================================
# 系统使用示例 - 新疆华兴农场
# ======================================

if __name__ == "__main__":
    # 初始化路径规划系统 (自动检测UTM分区)
    farm_planner = CommercialFarmPathPlanner("huaxing_farm.geojson")
    
    # 为所有田块生成作业路径
    for field_id in farm_planner.field_dict.keys():
        try:
            print(f"\n{'='*40}")
            print(f"处理田块: {field_id}")
            print(f"{'='*40}")
            
            # 生成作业路径
            field_path = farm_planner.generate_field_path(field_id, machine_width=4.0)
            
            # 规划到最近仓库的运输路径
            transport_route = farm_planner.plan_transport_route(field_id)
            
            # 生成路径可视化
            farm_planner.visualize(
                show_paths=True,
                save_path=f"field_{field_id}_paths.png"
            )
            
        except Exception as e:
            print(f"处理田块 {field_id} 时出错: {str(e)}")
    
    # 导出所有路径
    farm_planner.export_all_paths()
    
    # 生成最终报告
    report = farm_planner.generate_report()
    with open("farm_operation_report.json", "w") as f:
        json.dump(report, f, indent=2)
    
    print("\n" + "="*50)
    print("  农场路径规划完成!")
    print("="*50)
    print(f"  总田块数: {len(report['fields'])}")
    print(f"  总路径数: {len(farm_planner.paths)}")
    print(f"  总运输路径: {len(farm_planner.transport_routes)}")