Создайте цифровую модель вертикального изменения (DEM) данных об облаке точек
elevModel = pc2dem(ptCloudIn,gridResolution)
[___] = pc2dem(___, задает опции с помощью одних или нескольких аргументов name-value. Например, Name,Value)'CornerFillMethod','min' задает для функции, чтобы вычислить обобщенные значения вертикального изменения для углов сетки в DEM как минимальное вертикальное изменение в поисковом радиусе по умолчанию каждого угла сетки.
Создайте lasFileReader возразите, чтобы считать данные об облаке точек, хранимые в aerialLidarData.laz.
fileName = fullfile(toolboxdir("lidar"),"lidardata","las", ... "aerialLidarData.laz"); lasReader = lasFileReader(fileName);
Считайте данные об облаке точек из файла с помощью readPointCloud функция.
ptCloud = readPointCloud(lasReader);
Визуализируйте данные об облаке точек.
figure
pcshow(ptCloud.Location)
title("Point Cloud")
Сегментируйте наземные точки из данных об облаке точек.
groundPtsIdx = segmentGroundSMRF(ptCloud);
Извлеките наземные точки.
ptCloudWithGround = select(ptCloud,groundPtsIdx);
Визуализируйте наземные точки.
figure
pcshow(ptCloudWithGround.Location)
title("Ground Points")
Создайте цифровую модель ландшафта (DTM) из сегментированных наземных точек.
terrainModel = pc2dem(ptCloudWithGround);
Визуализируйте DTM.
figure
imagesc(terrainModel)
colormap(gray)
title("Digital Terrain Model")
Создайте lasFileReader возразите, чтобы считать воздушные данные об облаке точек из "aerialLidarData.laz".
fileName = fullfile(toolboxdir("lidar"),"lidardata","las", ... "aerialLidarData.laz"); lasReader = lasFileReader(fileName);
Считайте данные об облаке точек первого возврата датчика лидара из файла LAS с помощью readPointCloud функция.
ptCloud = readPointCloud(lasReader,"LaserReturns",1);Создайте цифровую поверхностную модель (DSM) облака точек с разрешением элемента сетки 1,1 метров.
gridRes = 1.1; surfaceModel = pc2dem(ptCloud,gridRes,"CornerFillMethod","max");
Задайте местоположение источника освещения.
azimuthAng = 135; zenithAng = 45;
Вычислите направленные градиенты DSM использование imgradientxy функция.
[gx,gy] = imgradientxy(surfaceModel,"sobel");Нормируйте градиенты с помощью разрешения элемента сетки.
gx = gx/(8*gridRes); gy = gy/(8*gridRes);
Вычислите наклон и аспект DSM.
slopeAngle = atand(sqrt(gx.^2 + gy.^2)); aspectAngle = atan2d(gy,-gx); aspectAngle(aspectAngle < 0) = aspectAngle(aspectAngle < 0) + 360;
Вычислите hillshade использование алгоритма от Esri®. hillshade является 3-D полутоновым представлением поверхности с относительным положением источника освещения, учтенного при штриховке изображения.
h = 255.0*((cosd(zenithAng).*cosd(slopeAngle)) ...
+ (sind(zenithAng).*sind(slopeAngle).*cosd(azimuthAng - aspectAngle)));
h(h < 0) = 0;Визуализируйте hillshade DSM.
figure imagesc(h) colormap(gray)
Функция использует локальный алгоритм раскладывания, чтобы создать цифровую модель вертикального изменения (DEM) данных об облаке точек. Алгоритм принимает, что вертикальное изменение точек приезжает z - ось.
Локальный алгоритм раскладывания:
Разделите облако точек на сетку по измерениям xy- (вид с высоты птичьего полета). Задайте размерности сетки с помощью gridResolution аргумент.
Используйте информацию о вертикальном изменении всех точек в круговой области вокруг каждого угла сетки, чтобы вычислить обобщенные значения сетки. Можно задать поисковый радиус и метод расчета с помощью 'SearchRadius' и 'CornerFillMethod' аргументы name-value, соответственно.
Если нет никаких точек в круговой области, алгоритм не вычисляет значение, и те углы сетки остаются незаполненными. Функция представляет их как NaN. Алгоритм использует интерполяцию обратного расстояния взвешивается (IDW), чтобы заполнить незаполненные углы сетки. Чтобы задать размер фильтра для метода интерполяции IDW, используйте 'FilterSize' аргумент значения имени.
pointCloud | velodyneFileReader | lasFileReader | segmentGroundSMRF | pcread