使用DEM进行直接定位

本教程的目的是使用Duvenhage算法通过与DEM（数字高程模型）的视线相交来计算直接定位网格。这个算法是Rugged中最高效的算法。

下图显示了在计算纬度、经度和海拔时考虑DEM的效果：

RuggedExplained.png

使用DEM瓦片提供Rugged数据

Rugged不解析DEM文件，而是将高程数据的缓冲区作为输入。调用应用程序负责读取DEM并将数据加载到缓冲区中。Rugged提供了一个具有缓存的瓦片机制，允许用户一次加载一个瓦片。这与全球覆盖DEM（如SRTM）的格式一致。Rugged提供了一个接口，用于使用回调函数更新缓存中的DEM瓦片，每当一个坐标落在当前区域之外时触发。

调用应用程序必须实现加载瓦片的回调函数。我们建议使用GDAL（http://www.gdal.org/）来解析DEM文件，因为它可以处理大多数格式（包括Aster DEM的geoTIFF或SRTM的DTED）。Rugged没有包含DEM的解析，这是出于设计考虑，尽管我们可以使用GDAL。我们希望Rugged保持一个低级别的库，不要引入太多第三方库。

实现DEM加载的TileUpdater接口。

在本教程中，我们不会包含真实的DEM数据。相反，我们将创建一个代表火山的虚拟DEM，它呈现为一个完美的圆锥体，类似于菲律宾的Mayon火山，只是我们将其定位在卫星的下方某个位置。这个示例已经是Rugged测试用例的一部分，源代码在包org.orekit.rugged.raster中，文件名为VolcanicConeElevationUpdater.java。

类VolcanicConeElevationUpdater实现了接口TileUpdater，并实现了其方法updateTile。该方法负责加载一个瓦片。瓦片的范围必须至少覆盖传递给该方法的地面点的坐标（纬度，经度）。瓦片是一个UpdatableTile类型的对象，它有两个方法：

  • setGeometry(minLatitude, minLongitude, latitudeStep, longitudeStep, latitudeRows, longitudeRows)：在加载之前初始化新瓦片的范围

  • setElevation(latIdx, longIdx, elevation)：用值elevation填充索引（latIdx，lonIdx）处的瓦片缓冲区

下面是类VolcanicConeElevationUpdater的源代码：

import org.hipparchus.util.FastMath;
import org.orekit.rugged.raster.TileUpdater;
import org.orekit.rugged.raster.UpdatableTile;
import org.orekit.bodies.GeodeticPoint;
import org.orekit.rugged.errors.RuggedException;
import org.orekit.utils.Constants;

public class VolcanicConeElevationUpdater implements TileUpdater {

    private GeodeticPoint summit;
    private double        slope;
    private double        base;
    private double        size;
    private int           n;

    public VolcanicConeElevationUpdater(GeodeticPoint summit, double slope, double base,
                                        double size, int n) {
        this.summit = summit;
        this.slope  = slope;
        this.base   = base;
        this.size   = size;
        this.n      = n;
    }

    public void updateTile(double latitude, double longitude, UpdatableTile tile)
        throws RuggedException {
        double step         = size / (n - 1);
        double minLatitude  = size * FastMath.floor(latitude  / size);
        double minLongitude = size * FastMath.floor(longitude / size);
        double sinSlope     = FastMath.sin(slope);
        tile.setGeometry(minLatitude, minLongitude, step, step, n, n);
        for (int i = 0; i < n; ++i) {
            double cellLatitude = minLatitude + i * step;
            for (int j = 0; j < n; ++j) {
                double cellLongitude = minLongitude + j * step;
                double distance       = Constants.WGS84_EARTH_EQUATORIAL_RADIUS *
                                        FastMath.hypot(cellLatitude  - summit.getLatitude(),
                                                       cellLongitude - summit.getLongitude());
                double altitude = FastMath.max(summit.getAltitude() - distance * sinSlope,
                                               base);
                tile.setElevation(i, j, altitude);
            }
        }
    }

}

DEM瓦片的重要注意事项：

  • 地面点的高程是通过4个相邻单元之间的双线性插值获得的。对于边界管理没有特定的算法。因此，落在瓦片边界上的点被认为是在外部。DEM瓦片在所有方向上必须至少重叠一行/列。

  • 在Rugged术语中，最小纬度和经度对应于DEM的最西南单元的中心。如果使用GDAL传递正确的信息，请注意gdalinfo中与左下角坐标有半个像素的偏移。

下图说明了邻近瓦片之间有一行/列重叠的正确DEM切片：

DEM-tiles-overlap.png

此图尝试表示定义瓦片时不同参数的含义：

tile-description.png

使用DEM初始化Rugged

初始化步骤与第一个教程直接定位略有不同，因为我们需要传递有关TileUpdater的信息。

实例化一个派生自TileUpdater的对象：

TileUpdater updater = new VolcanicConeElevationUpdater(summit, slope, base, FastMath.toRadians(1.0), 100);

int nbTiles = 8 ; // Rugged缓存中的最大瓦片数
AlgorithmId algoId = AlgorithmId.DUVENHAGE;

使用这些参数初始化Rugged：

Rugged rugged = new RuggedBuilder().
                setDigitalElevationModel(updater, nbTiles).
                setAlgorithm(algoId). 
                setEllipsoid(EllipsoidId.WGS84, BodyRotatingFrameId.ITRF).
                setTimeSpan(startDate, stopDate, 0.1, 10.0). 
                setTrajectory(InertialFrameId.EME2000,
                              satellitePVList, 6, CartesianDerivativesFilter.USE_P, 
                              satelliteQList, 8, AngularDerivativesFilter.USE_R).
                addLineSensor(lineSensor).
                build();

计算直接定位网格

与第一个教程直接定位类似，我们调用Rugged的直接定位方法。这次我们在循环中调用它，以便在磁盘上生成完整的网格。

DataOutputStream dos = new DataOutputStream(new FileOutputStream("demDirectLoc.c1"));
int lineStep = (maxLine - minLine) / nbLineStep;
int pxStep = (maxPx - minPx) / nbPxStep;

 for (int i = 0; i < nbLineStep; i++) {

    List<GeodeticPoint> pointList = new ArrayList<GeodeticPoint>(nbPxStep);
    int currentLine = minLine + i * lineStep;
    for (int j = 0; j < nbPxStep; j++) {
        int currentPx = minPx + j * pxStep;
        // 在当前行上调用直接定位
        Vector3D position = lineSensor.getPosition();
        AbsoluteDate currentLineDate = lineSensor.getDate(currentLine);
        Vector3D los = lineSensor.getLOS(absDate, currentPx);
        pointList.add(rugged.directLocation(currentLineDate, position, los));
    }
    for (GeodeticPoint point : pointList) {
        if (point != null) {
            dos.writeFloat((float) FastMath.toDegrees(point.getLatitude()));
        } else {
            dos.writeFloat((float) base);
        }
    }
    for (GeodeticPoint point : pointList) {
        if (point != null) {
            dos.writeFloat((float) FastMath.toDegrees(point.getLongitude()));
        } else {
            dos.writeFloat((float) base);
        }
    }
    for (GeodeticPoint point : pointList) {
        if (point != null) {
            dos.writeFloat((float) point.getAltitude());
        } else {
            dos.writeFloat((float) base);
        }
    }
}

源代码

源代码可以在DirectLocationWithDEM.java中找到（位于src/tutorials下的fr.cs.examples包中）

返回页面顶部