写一个高斯投影正反算的python代码

高斯投影是一种常见的地图投影方式，用于将地球表面的经纬度坐标转换为平面直角坐标。在这个例子中，我们将编写Python代码来实现简单的高斯投影的正向（从地理坐标到投影坐标）和反向（从投影坐标到地理坐标）计算。这里假设我们使用的是Mercator投影作为高斯投影的一种，它适用于中纬度地区。

首先安装必要的地理库`pyproj`（如果你还没有安装，可以使用pip install pyproj）：

pip install pyproj

然后，我们可以编写以下代码：

from pyproj import Proj, transform

def geodetic_to_mercator(lat, lon):
    """
    正向高斯投影，将地理坐标转换为Mercator投影坐标
    """
    # 创建一个 Mercator 投影对象
    proj = Proj(proj='merc', lat_ts=0)
    
    # 将经度和纬度转换为米制
    x, y = transform(proj, Proj(init='epsg:4326'), lon * 1e7, lat * 1e7)  # EPSG 4326 是 WGS84 纬度和经度
    return x, y

def mercator_to_geodetic(x, y):
    """
    反向高斯投影，将Mercator投影坐标转换回地理坐标
    """
    # 创建一个地理坐标系的对象
    inv_proj = Proj(proj='merc', lat_ts=0, inverse=True)
    
    # 将投影坐标转换回经度和纬度
    lon, lat = transform(inv_proj, Proj(init='epsg:4326'), x / 1e7, y / 1e7)  # 回转至WGS84
    return lat, lon

# 示例
geographic_lat, geographic_lon = 51.5074, -0.1278  # 伦敦坐标
projected_x, projected_y = geodetic_to_mercator(geographic_lat, geographic_lon)
print(f"Geographic: ({geographic_lon:.2f}, {geographic_lat:.2f}) -> Mercator: ({projected_x:.2f}, {projected_y:.2f})")

reversed_lat, reversed_lon = mercator_to_geodetic(projected_x, projected_y)
print(f"Mercator: ({projected_x:.2f}, {projected_y:.2f}) -> Geographic: ({reversed_lon:.2f}, {reversed_lat:.2f})")