python 东北区域：雷达回波图PPI（平面位置显示）

用的雷达基数据：

Z_RADR_I_Z9026_20260427010755_O_DOR_CAD_CAP_FMT_DPCTEST.bin.bz2
仰角：0.5°

雷达是绕着天线一圈一圈扫的（像切圆锥）。PPI 就是把某一个固定仰角（比如 0.5°、1.5°）扫描的那一层圆面，直接平摊在地图上展现出来。

• 特点：它是“扇面”或者“圆锥面”。因为你离雷达越远，波束飞得越高（受地球曲率影响）。

• 缺点：所以在 PPI 图上，离雷达近的地方看到的是低空，离雷达远的地方看到的是高空。它并不是一个真正高度一致的水平面。

• 

   

#!usr/bin/env python
# -*- coding:utf-8 -*-
"""
@author: Suyue
@file: leida.py
@time: 2024/07/30
@desc: 白底，全地图显示，盟市边界
"""
import os
import time
import warnings
import numpy as np
import cinrad
import matplotlib
matplotlib.use('TkAgg')
from cinrad.visualize import PPI

warnings.filterwarnings("ignore")

# ================= 定义路径区域 =================
DATA_DIR = "E:/东北区域雷达数据/20260427/Z9026"
SAVE_DIR = "E:/东北区域雷达数据/images2"
CITY_BOUNDARY_SHP = "D:/标准地图/内蒙古shp干旱业务机/盟市界.shp"  # 替换为你真实的shp路径
# ================================================

file_name = "Z_RADR_I_Z9026_20260427010755_O_DOR_CAD_CAP_FMT_DPCTEST.bin.bz2"
nFiles = os.path.join(DATA_DIR, file_name)

if not os.path.exists(SAVE_DIR):
    os.makedirs(SAVE_DIR)
    print(f"已自动创建保存目录: {SAVE_DIR}")

# 读取数据
f = cinrad.io.read_auto(nFiles)

# 获取REF数据
data = f.get_data(0, 230, "REF")

print(data)
print(type(f).__name__)

# 数据预处理
data["REF"].values = np.ma.masked_less(data["REF"].values, 0)

# ================= 关键修改：地图缩放 =================
# 雷达站点在 119.13E, 43.99N 附近
# extent 格式: [西经(最小), 东经(最大), 南纬(最小), 北纬(最大)]
# 比如：东经 115 到 124，北纬 40 到 48，这样整个内蒙古东部及周边都能看到
map_extent = [115, 124, 40, 48]

# 创建PPI对象
# section=False：不裁剪成扇形，让图变成完整的矩形地图
# extent=map_extent：强制扩大地图的显示范围
fig = cinrad.visualize.PPI(data,
                           style="white",
                           section=False,
                           extent=map_extent,
                           add_city_names=True)
# ===================================================

# 叠加自定义的盟市边界shp文件
try:
    fig.add_custom_shp(CITY_BOUNDARY_SHP, encoding='gbk', color='red', linewidth=1)
    print("已成功加载自定义盟市边界。")
except Exception as e:
    print(f"加载shp文件失败，请检查路径或文件内容。报错信息: {e}")

print("绘图完成...")

# 自动生成带时间戳的文件名
timestamp = time.strftime("%Y%m%d_%H%M%S")
output_path = os.path.join(SAVE_DIR, f"radar_{timestamp}.png")

# 保存图片
fig(output_path)
print(f"图片已成功保存至: {output_path}")

 叠加飞机作业路线

#!usr/bin/env python
# -*- coding:utf-8 -*-
"""
@author: Suyue
@file: leida.py
@time: 2026/09/09
@desc: 白底，全地图显示，盟市边界，叠加飞机作业路线
"""
import os
import time
import warnings
import numpy as np
import pandas as pd
import cinrad
import matplotlib

matplotlib.use('TkAgg')
from cinrad.visualize import PPI
import matplotlib.pyplot as plt

warnings.filterwarnings("ignore")

# ================= 定义路径区域 =================
DATA_DIR = "E:/东北区域雷达数据/20260427/Z9026"
SAVE_DIR = "E:/东北区域雷达数据/images2"
CITY_BOUNDARY_SHP = "D:/标准地图/内蒙古shp干旱业务机/盟市界.shp"
FLIGHT_CSV = "E:/东北区域雷达数据/NMG_B3857_202604271350.csv"  # 你的CSV路径
# ================================================

file_name = "Z_RADR_I_Z9026_20260427010755_O_DOR_CAD_CAP_FMT_DPCTEST.bin.bz2"
nFiles = os.path.join(DATA_DIR, file_name)

if not os.path.exists(SAVE_DIR):
    os.makedirs(SAVE_DIR)
    print(f"已自动创建保存目录: {SAVE_DIR}")

# 读取数据
f = cinrad.io.read_auto(nFiles)
data = f.get_data(0, 230, "REF")

print(data)

# 数据预处理
data["REF"].values = np.ma.masked_less(data["REF"].values, 0)

# 地图范围
map_extent = [115, 124, 40, 47]

# 创建PPI对象
fig = cinrad.visualize.PPI(data,
                           style="white",
                           section=False,
                           extent=map_extent,
                           add_city_names=True)

# 叠加盟市边界
try:
    fig.add_custom_shp(CITY_BOUNDARY_SHP, encoding='gbk', color='red', linewidth=1)
    print("已成功加载自定义盟市边界。")
except Exception as e:
    print(f"加载shp文件失败，请检查路径或文件内容。报错信息: {e}")

# ================= 叠加高亮飞机作业路线 =================
try:
    # 1. 读取CSV文件
    flight_data = pd.read_csv(FLIGHT_CSV, encoding='gb18030')
    lon = flight_data['经度'].values
    lat = flight_data['纬度'].values

    # 2. 【核心修改】获取真正的 Cartopy GeoAxes
    import cartopy.crs as ccrs

    ax = None
    # 遍历所有的轴，寻找哪个是 GeoAxes（地图轴）
    for axis in fig.fig.axes:
        if hasattr(axis, 'projection'):  # Cartopy 的 GeoAxes 都有 projection 属性
            ax = axis
            break

    if ax is None:
        raise ValueError("未找到地图坐标轴")

    # 3. 画“描边”线（黑色粗线，让飞机线在彩色回波上极其清晰）
    ax.plot(lon, lat, transform=ccrs.PlateCarree(), color='black', linewidth=1.5, zorder=100, solid_capstyle='round')

    # 4. 画主体航线（蓝色粗线，带上圆点，绝对在最顶层）
    ax.plot(lon, lat, transform=ccrs.PlateCarree(), color='blue', linewidth=1.5, marker='o', markersize=1.5, zorder=101,
            label='road')

    # 5. 标记起点和终点（同样必须加上 transform）
    ax.scatter(lon[0], lat[0], transform=ccrs.PlateCarree(), color='green', marker='*', s=150, zorder=102, label='start')
    # ax.scatter(lon[-1], lat[-1], transform=ccrs.PlateCarree(), color='red', marker='X', s=100, zorder=102, label='end')

    # 6. 添加图例
    ax.legend(loc='upper left', fontsize=10, framealpha=0.9)

    print("已成功叠加飞机作业路线。")
except Exception as e:
    print(f"叠加飞机路线失败，请检查CSV文件路径或内容。报错信息: {e}")
# ===================================================

print("绘图完成...")

# 自动生成带时间戳的文件名
timestamp = time.strftime("%Y%m%d_%H%M%S")
output_path = os.path.join(SAVE_DIR, f"radar_{timestamp}.png")

# 保存图片
fig(output_path)
print(f"图片已成功保存至: {output_path}")