如何从非网格 XZY 点生成 DEM？ （八度/MATLAB）

Question

我有相当大的激光雷达数据块（约 1600 万个点云），我已将其缩减为 X、Y、Z 数据（例如）：

436003.73  8758000.16  330.97  
436007.17  8758000.40  330.49  
436009.98  8758000.25  330.26  
436011.71  8758000.16  330.10 

我想为此数据创建一个 DEM，但由于数据集庞大，我的计算机发现这很困难。我不是特别精通 MATLAB/Octave，但我正在学习。 目前我正在使用griddata 和meshgrid 命令，例如：

cd ~/Octave
dem=load('436003_ed.txt');
x=dem(:,1);
y=dem(:,2);
z=dem(:,3);
minx=min(x);
miny=min(y);
minz=min(z);
maxz=max(z);
maxx=max(x);
maxy=max(y);
rijx=minx:2:maxx;
rijy=miny:2:maxy;
[XI,YI]=meshgrid(rijx,rijy);
ZI=griddata(x,y,z,XI,YI);

...这就是我在网格化数据时内存不足的情况。

我认为问题的一部分是因为我将整个数据集作为变量 (dem=load('436003_ed.txt');) 打开为 ~500MB 文件，这当然是非常消耗内存的。有人告诉我可以使用fopen 命令的另一种方法：

fid=fopen('436003_ed.txt');
XYZ=fscanf(fid,'%f %f  %f');  

但这似乎提供了一种奇怪的格式，因为它本身没有加载到 RAM 中？

有没有一种方法可以逐行加载它，将每个 Z 值放入适当的 DEM 单元格，同时计算该单元格的总高程和命中数（这样我就可以平均高程之后）生成 DEM？这将否定meshgrid 和griddata 命令的使用。
在这方面完全是新手，我不知道要查找和使用哪些命令！

提前非常感谢。我有 MATLAB R2013b，如果有帮助的话，我也一直在通过 Octave 做一些事情。

更新 我使用一些名为 GDAL 的命令行工具库（特别是http://www.gdal.org/ GDAL 网格）解决了这个问题。这有助于我对数据进行网格化。虽然这是对我查询的答案，但它不是我正在寻找的确切答案。 虽然我的原始方法可能会起作用，但如果我的计算机有足够的内存，还有其他方法可以在 MATLAB 中编写脚本 - 使用 fscanf 和 fopen 命令 - 我还没有发现如何正确执行此操作。

Answer 1

Octave 目前使用 GraphicsMagick 或 ImageMagick 进行图像 IO。但是，这两个项目都只处理整数，并且您的数据是浮点数，因此您将无法使用典型的 imwrite 函数将其保存为图像。如果您真的想将其保存为图像，例如，您可以编写自己的 libtiff 接口。我相信 Matlab 可以处理浮点图像。

关于您的代码，我不知道您可能会如何处理内存不足，但您有一些不必要的行。您可以替换：

x=dem(:,1);
y=dem(:,2);
z=dem(:,3);
minx=min(x);
miny=min(y);
minz=min(z);
maxz=max(z);
maxx=max(x);
maxy=max(y);
rijx=minx:2:maxx;
rijy=miny:2:maxy;

与：

mins = min (dem);
maxs = max (dem);
rijx = mins(1):2:maxs(1);
rijy = mins(2):2:maxs(2);