地理坐标转换行列号问题

Question

上次u给的这个函数也用不啊，我是2022 a版本，问题就是转换后出现负值，而且还超出矩阵

Answer 1

引自免费微信小程序：皆我百晓生

在MATLAB R2022a版本中，worldToPixel函数可能未包含在标准库中。你需要创建一个自定义函数来实现从经纬度到图像像素坐标的转换。这里提供一个基于UTM（Universal Transverse Mercator）投影的简单示例。首先确保你有proj4库支持，如果没有，请安装它。

安装 proj4 库:

!pip install pyproj

然后在MATLAB中使用Python交互式环境调用这个库。下面是一个完整的例子，包括如何将经纬度转换为像素坐标，并处理可能出现的负值和超出矩阵范围的问题。

function [px, py] = worldToPixel(info, lat, lon)
    % 使用Python交互式环境导入必要的库
    python.eval('import numpy as np')
    python.eval('from osgeo import ogr, osr')

    % 创建一个地理坐标系对象
    python.eval('src_srs = osr.SpatialReference()')
    python.eval('src_srs.ImportFromEPSG(4326)') % WGS84地理坐标系

    % 获取图像信息中的投影坐标系
    python.eval('dst_srs = osr.SpatialReference()')
    python.eval(['dst_srs.ImportFromWkt("' info.Projection '+")'])

    % 创建一个变换对象
    python.eval('transform = osr.CoordinateTransformation(src_srs, dst_srs)')

    % 将经纬度转换为UTM投影坐标
    utm_x, utm_y, _ = python.call('transform.TransformPoint', lon, lat)

    % 根据info.RefMatrix获取图像的大小
    width = size(info.RefMatrix, 2);
    height = size(info.RefMatrix, 1);

    % 计算像素坐标
    px = utm_x * width / info.RefMatrix.RasterXSize + 0.5;
    py = height - (utm_y * height / info.RefMatrix.RasterYSize + 0.5);

    % 处理可能出现的负值和超出矩阵范围的问题
    px(px < 1) = 1;
    px(px > width) = width;
    py(py < 1) = 1;
    py(py > height) = height;

end

接下来，你需要在你的代码中使用这个自定义的worldToPixel函数。确保你的info结构体包含了正确的投影信息，例如Projection、RefMatrix等。以下是使用这个函数的例子：

[x, y] = worldToPixel(info, lat, lon);

这样应该能解决你的问题。注意，这个方法依赖于Python环境中的pyproj库，如果你没有这个库，你需要先安装。另外，确保你的数据是在合适的投影系统下，这里是基于UTM投影进行转换的。如果你的数据使用了其他投影系统，你需要相应地修改代码中的投影设置。