Семантическая сегментация Aerial Lidar с использованием Deep Learning

В этом примере используются:

В этом примере показано, как обучить сеть глубокого обучения, предназначенную для выполнения семантической сегментации воздушных данных лидара.

Данные Лидара, полученные с помощью бортовых лазерных сканирующих систем, все чаще используются в таких областях, как топографическое картирование, моделирование зданий и городов, измерение биомассы и борьба со стихийными бедствиями. Извлечение значимой информации из этих данных требует семантической сегментации, процесса, в котором каждому пикселю в облаке точек назначается метка класса.

В этом примере вы обучаете сеть, предназначенную для выполнения семантической сегментации, с помощью набора данных DALES (Dayton annotated lidar earth scan) [1], который содержит сцены плотных, помеченных lidar воздушных данных из городских, пригородных, сельских и коммерческих настроек. Из 40 сцен 29 сцен используются для обучения, а остальные 11 последовательностей используются для тестирования. Набор данных предоставляет семантические метки сегментации для 8, таких как здания, легковые автомобили, грузовики, столбы, линии электропередачи, ограждения, земля и растительность.

Загрузить данные DALES

Набор данных DALES содержит 40 сцен воздушных данных лидара. Набор данных уже разделен на 29 сцен для обучения и 11 сцен для тестирования. Каждый пиксель в данных имеет метку класса. Следуйте инструкциям на веб-сайте DALES, чтобы загрузить набор данных в папку, указанную dataFolder переменная. Создание папок для хранения данных обучения и тестирования.

dataFolder = fullfile(tempdir,'DALES');
trainDataFolder = fullfile(dataFolder,'dales_las','train');
testDataFolder = fullfile(dataFolder,'dales_las','test');

Предварительный просмотр облака точек из учебных данных.

lasReader = lasFileReader(fullfile(trainDataFolder,'5080_54435.las'));
[pc,attr] = readPointCloud(lasReader,"Attributes", "Classification");
labels = attr.Classification + 1;
classNames = [
    "background"
    "ground"
    "vegetation"
    "cars"
    "trucks"
    "powerlines"
    "poles"
    "fences"
    "buildings"
    ];
cmap = labelColorMap();
figure;
colormap(cmap);
ax = pcshow(pc.Location,labels);
labelColorbar(ax,cmap,classNames);
title("PointCloud with overlaid semantic labels");

Данные предварительной обработки

Каждое точечное облако в наборе данных DALES занимает площадь 500 на 500 метров, что намного больше, чем типичная площадь, покрытая наземными вращающимися лидарными точечными облаками. Для эффективной обработки памяти разделите облако точек на небольшие неперекрывающиеся сетки. Кроме того, сосредоточьтесь на трех доминирующих классах - земле, зданиях и растительности - чтобы избежать проблем из-за дисбаланса классов.

Используйте helperCropPointCloudsAndMergeLabels функция, присоединенная к этому примеру в качестве вспомогательного файла, для:

Обрезать облака точек на неперекрывающиеся сетки размером 50 на 50 метров.
Объедините все классы, отличные от доминирующих (земля, здания и растительность), в единый фоновый класс, чтобы данные содержали только четыре класса.
Сохраните обрезанные сетки и семантические метки в файлах PCD и PNG соответственно.

Определите размеры сетки и задайте фиксированное количество точек на сетку, чтобы ускорить обучение.

gridSize = [50,50];
numPoints = 115000;

Если этот рабочий процесс используется для собственных данных, установите writeFiles «false», если данные обучения разделены на сетки. Обратите внимание, что данные обучения должны быть в формате, поддерживаемом pcread функция.

writeFiles = true;
[pcCropTrainPath,labelsCropTrainPath] = helperCropPointCloudsAndMergeLabels(gridSize,trainDataFolder,numPoints,writeFiles);
[pcCropTestPath,labelsCropTestPath] = helperCropPointCloudsAndMergeLabels(gridSize,testDataFolder,numPoints,writeFiles);

Обработка может занять некоторое время. Код приостанавливает выполнение MATLAB ® до завершения обработки.

Создание объектов хранилища данных для обучения

Создать fileDatastore объект для загрузки PCD-файлов с помощью pcread функция.

ldsTrain = fileDatastore(pcCropTrainPath,'ReadFcn',@(x) pcread(x));

Использовать pixelLabelDatastore объект для хранения пиксельных меток из изображений пиксельных меток. Объект сопоставляет каждую пиксельную метку с именем класса. В этом примере грунт, растительность и здания являются единственными интересующими классами; все остальные пикселы являются фоном. Укажите эти классы и присвойте каждому классу уникальный идентификатор метки.

classNames = classNames([1,2,3,9],:)

classNames = 4×1 string
    "background"
    "ground"
    "vegetation"
    "buildings"

numClasses = numel(classNames);

% Specify label IDs from 1 to the number of classes.
labelIDs = 1 : numClasses;
pxdsTrain = pixelLabelDatastore(labelsCropTrainPath,classNames,labelIDs);

Загрузите и отобразите облака точек.

ptcld = preview(ldsTrain);
labels = preview(pxdsTrain);
cmap = cmap([1,2,3,9],:);
figure;
ax1 = pcshow(ptcld.Location,uint8(labels));
labelColorbar(ax1,cmap,classNames);
title("Cropped point cloud with overlaid semantic labels");

Используйте normalizePointCloud функция, определенная в конце примера, для нормализации облака точек между 0 и 1.

ldsTransformed = transform(ldsTrain,@(x) normalizePointCloud(x));

Используйте onehotEncodeLabels функция, определенная в конце примера, для преобразования меток в однокодированное представление для обучения сети сегментации с использованием пользовательского обучающего цикла.

pxdsTransformed = transform(pxdsTrain,@(x) onehotEncodeLabels(x,classNames));

Используйте combine функция для объединения облаков точек и меток пикселей в единое хранилище данных для обучения.

dsTrain = combine(ldsTransformed,pxdsTransformed);

Определение модели «Net++»

PointNet ++ [2] модель сегментации состоит из двух главных компонентов:

Установка модулей абстракции
Модули распространения функций

Серия модулей абстракции множеств постепенно подбирает интересующие точки путем иерархического группирования точек и использует пользовательскую архитектуру, чтобы кодировать точки в векторы признаков. Поскольку задачи семантической сегментации требуют точечных элементов для всех исходных точек, в этом примере используется ряд модулей распространения элементов для иерархической интерполяции элементов в исходные точки с использованием схемы интерполяции на основе обратного расстояния.

Модуль абстракции множества реализован с помощью sampleAndGroup и sharedMLP вспомогательные функции и модуль распространения функций реализован с использованием featurePropagation и sharedMLP вспомогательные функции. sharedMLP модуль является единственным обучаемым модулем, который состоит из уровней свертки и нормализации экземпляра. Использовать вспомогательную функцию pointNetplusModel для определения всей модели сегментации, выполняемой с помощью ПО StartNet++. Вспомогательные функции перечислены в конце этого примера.

Определите PointNet ++ образцовые параметры

sampleAndGroup функция модуля абстракции набора параметризуется количеством кластеров, количеством выборок на кластер и радиусом каждого кластера, тогда как sharedMLP функция реализует пользовательскую архитектуру, которая параметризуется с использованием количества входных каналов и скрытых размеров каналов.

В этом примере используются значения гиперпараметров, настроенные в наборе данных DALES. Если Вы хотите применить PointNet ++ к другому набору данных, Вы должны выполнить дополнительную настройку гиперпараметра.

Задайте параметры первого модуля абстракции набора. Установите количество кластеров 1024, число выборок в каждом кластере 128 и радиус каждого кластера 0,1. Также установите размер входного канала модели, равный трем, а размер скрытого канала - 32, 32 и 64.

inputChannelSize = 3;
hiddenChannelSize = [32,32,64];
nontrainable.nClusters1 = 1024;
nontrainable.nSamples1 = 128;
nontrainable.radius1 = 0.1;
trainable.SharedMLP1 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Задайте параметры модуля абстракции второго набора. Установите число кластеров 256, число выборок в каждом кластере 64, а радиус каждого кластера 0,2. Также установите размер входного канала модели, равный 67, а размер скрытого канала - 64, 64 и 128.

inputChannelSize = 67;
hiddenChannelSize = [64,64,128];
nontrainable.nClusters2 = 256;
nontrainable.nSamples2 = 64;
nontrainable.radius2 = 0.2;
trainable.SharedMLP2 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Задайте параметры модуля абстракции третьего набора. Установите число кластеров равным 64, число выборок в каждом кластере - 32, а радиус каждого кластера - 0,4. Также установите размер входного канала модели, равный 131, и размер скрытого канала, равный 128, 128 и 256.

inputChannelSize = 131;
hiddenChannelSize = [128,128,256];
nontrainable.nClusters3 = 64;
nontrainable.nSamples3 = 32;
nontrainable.radius3 = 0.4;
trainable.SharedMLP3 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Задайте параметры модуля абстракции четвертого набора. Установите число кластеров в 16, число выборок в каждом кластере в 32 и радиус каждого кластера в 0,8. Кроме того, установите размер входного канала для модели, равный 259, а размер скрытого канала - 256, 256 и 512.

inputChannelSize = 259;
hiddenChannelSize = [256,256,512];
nontrainable.nClusters4 = 16;
nontrainable.nSamples4 = 32;
nontrainable.radius4 = 0.8;
trainable.SharedMLP4 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Задайте параметры первого модуля распространения функций. Задайте для пятой общей модели MLP размер входного канала 768, а для скрытого канала 256 и 256.

inputChannelSize = 768;
hiddenChannelSize = [256,256];
trainable.SharedMLP5 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Задайте параметры второго модуля распространения функций. Установите размер входного канала совместно используемой модели MLP на 384, а размер скрытого канала на 256 и 256.

inputChannelSize = 384;
hiddenChannelSize = [256,256];
trainable.SharedMLP6 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Задайте параметры третьего модуля распространения функций. Установите размер входного канала совместно используемой модели MLP на 320, а размер скрытого канала на 256 и 128.

inputChannelSize = 320;
hiddenChannelSize = [256,128];
trainable.SharedMLP7 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Задайте параметры четвертого модуля распространения функций. Установите размер входного канала совместно используемой модели MLP на 128, а размер скрытого канала на 128 и 128.

inputChannelSize = 128;
hiddenChannelSize = [128,128];
trainable.SharedMLP8 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);

Установите окончательный размер входного канала совместно используемой модели MLP на 128, а размер скрытого канала на 128 и numClasses.

inputChannelSize = 128;
hiddenChannelSize = [128,numClasses];
trainable.SharedMLP9 = initializeSharedMLP(inputChannelSize,hiddenChannelSize);
params.trainable = trainable;
params.nontrainable = nontrainable;

Укажите параметры обучения

Тренируйся на одну эпоху. Установите начальную скорость обучения 2e-4 и коэффициент регуляризации L2 0,01.

numEpochs = 1;
miniBatchSize = 1;
learningRate = 2e-4;
l2Regularization = 0.01;

Инициализируйте параметры оптимизации Adam.

gradientDecayFactor = 0.9;
squaredGradientDecayFactor = 0.999;

Инициализируйте скользящее среднее градиентов параметров и квадратов по элементам градиентов, используемых оптимизатором Адама.

trailingAvg = [];
trailingAvgSq = [];

Модель поезда

Обучение модели с помощью ЦП или графического процессора. Для использования графического процессора требуются параллельные вычислительные Toolbox™ и графический процессор NVIDIA ® с поддержкой CUDA ®. Дополнительные сведения см. в разделе Поддержка графического процессора по выпуску (Панель инструментов параллельных вычислений). Чтобы автоматически определить, доступен ли графический процессор, установитеexecutionEnvironment кому 'auto'. Если у вас нет графического процессора или вы не хотите использовать его для обучения, установите executionEnvironment кому 'cpu'. Для обеспечения использования графического процессора для обучения, установить executionEnvironment кому 'gpu'.

Далее создайте minibatchqueue (Deep Learning Toolbox) объект для загрузки данных в пакеты miniBatchSize во время обучения.

executionEnvironment = "auto";
if canUseParallelPool
    dispatchInBackground = true;
else
    dispatchInBackground = false;
end

mbq = minibatchqueue(dsTrain,2,...
                     "MiniBatchSize",miniBatchSize,...
                     "OutputEnvironment",executionEnvironment,...
                     "MiniBatchFormat",["SSCB","SSCB"],...
                     "DispatchInBackground",dispatchInBackground);

Обучение модели с помощью пользовательского цикла обучения. Для каждой итерации:

Прочитайте пакет данных.
Оцените градиенты модели.
Примените регуляризацию веса L2.
Использовать adamupdate для обновления параметров модели.
Обновление плоттера хода обучения с помощью функции помощника configureTrainingProgressPlotter, определяется в конце этого примера.

doTraining = false;
if doTraining
 
    % Initialize plot.
    fig = figure;
    lossPlotter = configureTrainingProgressPlotter(fig);    
    iteration = 0;
    
    % Custom training loop.
    for epoch = 1:numEpochs
         
        % Reset datastore.
        reset(mbq);
        
        while(hasdata(mbq))
            iteration = iteration + 1;
            
            % Read batch of data.
            [ptCld, ptLabels] = next(mbq);
                        
            % Evaluate the model gradients and loss using dlfeval and the modelGradients function.
            [gradients, loss] = dlfeval(@modelGradients, ptCld, ptLabels, params);
                                        
            % L2 regularization.
            gradients = dlupdate(@(g,p) g + l2Regularization*p,gradients,params.trainable);
            
            % Update the network parameters using the Adam optimizer.
            [params.trainable, trailingAvg, trailingAvgSq] = adamupdate(params.trainable, gradients, ...
                                                       trailingAvg, trailingAvgSq, iteration,...
                                                       learningRate,gradientDecayFactor, squaredGradientDecayFactor);
            
            % Update training plot with new points.         
            addpoints(lossPlotter, iteration,double(gather(extractdata(loss))));
            title("Training Epoch " + epoch +" of " + numEpochs);
            drawnow;
        end      
    end
else
    % Download pretrained model parameters.
    pretrainedNetURL = 'https://www.mathworks.com/supportfiles/lidar/data/trainedPointNetplusplusNet.zip';
    
    preTrainedZipFile = fullfile(dataFolder,'trainedPointNetplusplusNet.zip');
    preTrainedMATFile = fullfile(dataFolder,'trainedPointNetplusplusNet.mat');
     if ~exist(preTrainedMATFile,'file')
            if ~exist(preTrainedZipFile,'file')
                disp('Downloading pretrained detector (3.2 MB)...');
                websave(preTrainedZipFile, pretrainedNetURL);
            end
            unzip(preTrainedZipFile,dataFolder);   
     end
    
    % Load pretrained model.
    load(preTrainedMATFile,'params');
    
    % Move model parameters to the GPU if possible and convert to a dlarray.
    params.trainable = prepareForPrediction(params.trainable,@(x)dlarray(toDevice(x,canUseGPU)));   
end

Облако воздушных точек сегмента

Прочитайте полное облако контрольных точек.

lasReader = lasFileReader(fullfile(testDataFolder,'5175_54395.las'));
[pc,attr] = readPointCloud(lasReader,"Attributes", "Classification");

Вычислить количество неперекрывающихся сеток на основе gridSize, XLimits, и YLimits облака точек.

numGridsX = round(diff(pc.XLimits)/gridSize(1));
numGridsY = round(diff(pc.YLimits)/gridSize(2));
xRange = linspace(pc.XLimits(1),pc.XLimits(2),numGridsX+1);
yRange = linspace(pc.YLimits(1),pc.YLimits(2),numGridsY+1);

Выполните итерацию по всем неперекрывающимся сеткам и спрогнозируйте метки.

pcFinal = [];
for r = 1:numGridsX
    for c = 1:numGridsY
        idx = (pc.Location(:,1) >= xRange(r)) & (pc.Location(:,1) < xRange(r+1))...
             & (pc.Location(:,2) >= yRange(c)) & (pc.Location(:,2) < yRange(c+1));
        ptCloudGrid = select(pc,idx);
        
        % Select fixed number of points and labels from the given point cloud.
        ptCloudGrid = selectPoints(ptCloudGrid,numPoints);
        
        % Apply the preprocessing.
        ptCloud = normalizePointCloud(ptCloudGrid);
        loc = ptCloud{1,1};
        X1 = dlarray(loc,'SSCB');
        if (executionEnvironment == "auto" && canUseGPU) || executionEnvironment == "gpu"
            X1 = gpuArray(X1);
        end
        
        % Get the output predictions.
        op = pointNetplusModel(X1,params.trainable,params.nontrainable);
        [~,opLabels] = max(gather(extractdata(op)),[],3);
        color = cmap(opLabels,:);
        ptCloudGrid.Color = uint8(color); 
        
        pcFinal = [pcFinal;ptCloudGrid];
    end
end

ptCloudOut = pccat(pcFinal);
color = 255.*ptCloudOut.Color;
figure;
pcshow(ptCloudOut.Location,color);
title('PointCloud overlaid with detected semantic labels');

Оценка сети

Оцените производительность сети на тестовых данных. Используйте evaluateSemanticSegmentation функция для вычисления метрик семантической сегментации по результатам тестового набора.

confusionMatrix = {};

% Define the number of samples to evaluate the network. Set numSamplesToTest to
% 1100 to evaluate the model on the entire test data set.
numSamplesToTest = 50;

for i = 1:numSamplesToTest
    ptCldPath = fullfile(pcCropTestPath,sprintf("%03d.pcd",i));
    ptCloud = pcread(ptCldPath);
    ptCloud = normalizePointCloud(ptCloud);
    loc = ptCloud{1,1};
    
    X1 = dlarray(loc,'SSCB');
    if (executionEnvironment == "auto" && canUseGPU) || executionEnvironment == "gpu"
        X1 = gpuArray(X1);
    end
    
    % Get the output predictions.
    op = pointNetplusModel(X1,params.trainable,params.nontrainable);
    [~,predictionLabels] = max(gather(extractdata(op)),[],3);
  
    targetLabelPath = fullfile(labelsCropTestPath,sprintf("%03d.png",i));
    targetLabels = imread(targetLabelPath);
    targetLabels = squeeze(targetLabels);
    
    % Calculate the confusion matrix.
    confMatrix = segmentationConfusionMatrix(double(predictionLabels),...
                 double(targetLabels),"Classes",1:numClasses);
    if size(confMatrix,1) ~= numel(classNames)
        continue
    end
    confusionMatrix{i} = confMatrix;
end

confusionMatrix = confusionMatrix';
metrics = evaluateSemanticSegmentation(confusionMatrix,classNames)

Evaluating semantic segmentation results
----------------------------------------
* Selected metrics: global accuracy, class accuracy, IoU, weighted IoU.
* Processed 50 images.
* Finalizing... Done.
* Data set metrics:

    GlobalAccuracy    MeanAccuracy    MeanIoU    WeightedIoU
    ______________    ____________    _______    ___________

       0.76168          0.53687       0.42486      0.60851

metrics = 
  semanticSegmentationMetrics with properties:

              ConfusionMatrix: [4×4 table]
    NormalizedConfusionMatrix: [4×4 table]
               DataSetMetrics: [1×4 table]
                 ClassMetrics: [4×2 table]
                 ImageMetrics: [50×4 table]

Градиенты модели

modelGradients функция принимает в качестве входных данных мини-пакет данных X, соответствующая цель yTargetи обучаемые параметры, и возвращает градиенты потерь относительно обучаемых параметров и соответствующих потерь. Чтобы вычислить градиенты, вычислите modelGradients с использованием функции dlfeval в учебном цикле.

function [gradients, loss] = modelGradients(X,yTarget,params)

    % Execute the model function.
    yPred = pointNetplusModel(X,params.trainable,params.nontrainable);
    
    % Compute the loss.
    loss = focalCrossEntropy(yPred,yTarget,'TargetCategories','independent');
    
    % Compute the parameter gradients with respect to the loss. 
    gradients = dlgradient(loss, params.trainable);  
end

Функция модели «»/Net + +

Функция модели принимает в качестве входных параметры модели parameters and входные данные dlX. Сеть выводит прогнозы для меток.

function dlY = pointNetplusModel(dlX,trainable,nontrainable)

    % Set abstraction module 1
   [points1, pointFeatures1] = sampleAndGroup(dlX,[],nontrainable.nClusters1,...
                                nontrainable.nSamples1,nontrainable.radius1);
    dlY = sharedMLP(pointFeatures1,trainable.SharedMLP1.Perceptron);
    pointNetFeatures1 = max(dlY,[],2);
    
    % Set abstraction module 2
    [points2, pointFeatures2] = sampleAndGroup(points1,pointNetFeatures1,nontrainable.nClusters2,...
                                nontrainable.nSamples2,nontrainable.radius2);
    dlY = sharedMLP(pointFeatures2,trainable.SharedMLP2.Perceptron);
    pointNetFeatures2 = max(dlY,[],2);
    
    % Set abstraction module 3
    [points3, pointFeatures3] = sampleAndGroup(points2,pointNetFeatures2,nontrainable.nClusters3,...
                                nontrainable.nSamples3,nontrainable.radius3);  
    dlY = sharedMLP(pointFeatures3,trainable.SharedMLP3.Perceptron);
    pointNetFeatures3 = max(dlY,[],2);
    
    % Set abstraction module 4
   [points4, pointFeatures4] = sampleAndGroup(points3,pointNetFeatures3,nontrainable.nClusters4,...
                                nontrainable.nSamples4,nontrainable.radius4); 
    dlY = sharedMLP(pointFeatures4,trainable.SharedMLP4.Perceptron);
    pointNetFeatures4 = max(dlY,[],2);
    
    % Feature propagation module 1
    pointsFP1 = featurePropagation(points3,points4,pointNetFeatures3,pointNetFeatures4);
    pointNetFP1 = sharedMLP(pointsFP1,trainable.SharedMLP5.Perceptron);
    
    % Feature propagation module 2
    pointsFP2 = featurePropagation(points2,points3,pointNetFeatures2,pointNetFP1);
    pointNetFP2 = sharedMLP(pointsFP2,trainable.SharedMLP6.Perceptron);
   
    % Feature propagation module 3
    pointsFP3 = featurePropagation(points1,points2,pointNetFeatures1,pointNetFP2);
    pointNetFP3 = sharedMLP(pointsFP3,trainable.SharedMLP7.Perceptron);
    
    % Feature propagation module 4
    pointsFP4 = featurePropagation(dlX,points1,[],pointNetFP3);
    dlY = sharedMLP(pointsFP4,trainable.SharedMLP8.Perceptron);
   
    % Shared MLP
    dlY = sharedMLP(dlY,trainable.SharedMLP9.Perceptron);
    dlY = softmax(dlY);
    dlY = dlarray(dlY,'SSCB');
end

Функции инициализации параметров модели

Инициализация общей многослойной функции перцептрона

initializeSharedMLP функция принимает в качестве входного значения размер канала и размер скрытого канала и возвращает инициализированные параметры в структуре. Параметры инициализируются с помощью инициализации He weight.

function params = initializeSharedMLP(inputChannelSize,hiddenChannelSize)
weights = initializeWeightsHe([1 1 inputChannelSize hiddenChannelSize(1)]);
bias = zeros(hiddenChannelSize(1),1,"single");
p.Conv.Weights = dlarray(weights);
p.Conv.Bias = dlarray(bias);

p.InstanceNorm.Offset = dlarray(zeros(hiddenChannelSize(1),1,"single"));
p.InstanceNorm.Scale = dlarray(ones(hiddenChannelSize(1),1,"single"));

params.Perceptron(1) = p;

for k = 2:numel(hiddenChannelSize)
    weights = initializeWeightsHe([1 1 hiddenChannelSize(k-1) hiddenChannelSize(k)]);
    bias = zeros(hiddenChannelSize(k),1,"single");
    p.Conv.Weights = dlarray(weights);
    p.Conv.Bias = dlarray(bias);
    
    p.InstanceNorm.Offset = dlarray(zeros(hiddenChannelSize(k),1,"single"));
    p.InstanceNorm.Scale = dlarray(ones(hiddenChannelSize(k),1,"single"));
   
    params.Perceptron(k) = p;
end
end

function x = initializeWeightsHe(sz)
fanIn = prod(sz(1:3));
stddev = sqrt(2/fanIn);
x = stddev .* randn(sz);
end

function x = initializeWeightsGaussian(sz)
x = randn(sz,"single") .* 0.01;
end

Отбор проб и группирование

function [newClusters,newpointFeatures] = sampleAndGroup(points,pointFeatures,nClusters,nSamples,radius)
% The sampleAndGroup layer first samples the point cloud to a given number of
% clusters and then constructs local region sets by finding neighboring
% points around the centroids using the queryBallPoint function.

    points = extractdata(squeeze(points));
    if ~isempty(pointFeatures)
        pointFeatures = extractdata(squeeze(pointFeatures));
    end
        
    % Find the centroids for nClusters - nClusters*3.
    centroids = farthestPointSampling(points,nClusters);
    newClusters = points(centroids,:);
    
    % Find the neareset nSamples for nClusters - nClusters*nSamples*3.
    idx = queryBallPoint(points,newClusters,nClusters,nSamples,radius);
    newPoints = reshape(points(idx,:),[nClusters,nSamples,3]);
    
    % Normalize the points in relation to the cluster center.
    newpointFeatures = newPoints - reshape(newClusters,nClusters,1,3);
    
    if ~isempty(pointFeatures)
        groupFeatures = reshape(pointFeatures(idx,:),...
                       [nClusters,nSamples,size(pointFeatures,2)]);
        newpointFeatures = cat(3,newPoints,groupFeatures);
    end
    
    newpointFeatures = dlarray(newpointFeatures,'SSC');
    newClusters = dlarray(newClusters,'SSC');
end

Наиболее дальняя точка выборки

function centroids = farthestPointSampling(pointLocations,numPoints)
% The farthestPointSampling function selects a set of points from input
% points, which defines the centroids of local regions.
% pointLocations - PointCloud locations N-by-3.
% numPoints - Number of clusters to find.
% centroids - centroids of each farthest cluster.
    
    % Initialize initial indices as zeros.
    centroids = zeros(numPoints,1);
    
    % Distance from centroid to each point.
    distance = ones(size(pointLocations,1),1) .* 1e10; 
    
    % Random Initialization of the first point.
    farthest = randi([1,size(pointLocations,1)],1);
    
    for i = 1:numPoints
        centroids(i) = farthest;
        centroid = pointLocations(farthest,:);
        dist = sum(((pointLocations - centroid).^2),2);
        mask = dist < distance;
        distance(mask) = dist(mask);
        [~,farthest] = max(distance,[],1);
    end
end

Запрос точки шарика

function groupIdx = queryBallPoint(XYZ,newXYZ,nClusters,nSamples,radius)
% Given the cluster center, the queryBallPoint finds all points that are
% within a radius to the query point.

    N = size(XYZ,1);
    groupIdx = reshape(1:N,[1,N]);
    groupIdx = repmat(groupIdx,[nClusters,1]);
    
    % Find the distance between the centroids and given points.
    sqDist = squareDistance(newXYZ,XYZ);    
    
    % Find the points that are inside the given radius.
    groupIdx(sqDist > (radius)^2) = N;
    groupIdx = sort(groupIdx,2,"ascend");
    
    % Find the closest nSamples points within the given radius.
    groupIdx = groupIdx(:,1:nSamples);
    groupFirst = repmat(groupIdx(:,1),1,nSamples);
    mask = (groupIdx == N);
    groupIdx(mask) = groupFirst(mask);
end

Распространение элементов

function newPoints = featurePropagation(points1,points2,pointNetFeatures1,pointNetFeatures2)
    % Use the inverse distance weighted average based on the k nearest neighbors to
    % interpolate features.
    
    points1 = extractdata(squeeze(points1));
    points2 = extractdata(squeeze(points2));
    if ~isempty(pointNetFeatures1)
    pointNetFeatures1 = extractdata(squeeze(pointNetFeatures1));
    end
    pointNetFeatures2 = extractdata(squeeze(pointNetFeatures2));
    
    % Find the K nearest neighbors for each point.
    dists = squareDistance(points1,points2);
    [dists,idx] = sort(dists,2,"ascend");
    dists = dists(:,1:3);
    idx = idx(:,1:3);
    
    % Calculate the weights for interpolation.
    dist_recip = 1./(dists+1e-8);
    normFactor = sum(dist_recip,2);
    weights = dist_recip./normFactor;
    
    % Perform weighted interpolation.
    interpolatedPoints = pointNetFeatures2(idx,:);
    interpolatedPoints = reshape(interpolatedPoints,[size(idx,1),size(idx,2),size(pointNetFeatures2,2)]);
    
    interpolatedPoints = interpolatedPoints .* weights;
    interpolatedPoints = squeeze(sum(interpolatedPoints,2));
    
    if ~isempty(pointNetFeatures1)
        % Calculate the new points.
        newPoints = cat(2,pointNetFeatures1,interpolatedPoints);
    else
        newPoints = interpolatedPoints;
    end
    
     newPoints = dlarray(newPoints,'SCS');
end

% Find the squared distance.
function dist = squareDistance(src,dst)
    dist = -2 * (src*permute(dst,[2,1,3]));
    tmp1 = sum(src.^2,2);
    tmp1 = reshape(tmp1,[size(src,1),1]);
    tmp2 = sum(dst.^2,2);
    tmp2 = reshape(tmp2,[1,size(dst,1)]);
    dist = dist + tmp1 + tmp2;
end

Общая многослойная функция перцептрона

function dlY= sharedMLP(dlX,parameters)
% The shared multilayer perceptron function processes the input dlX using a
% series of perceptron operations and returns the result of the last
% perceptron.
    dlY = dlX;
    for k = 1:numel(parameters) 
        dlY = perceptron(dlY,parameters(k));
    end
end

Функция перцептрона

function dlY = perceptron(dlX,parameters)
% The perceptron function processes the input dlX using a convolution, a
% instance normalization, and a ReLU operation and returns the output of the
% ReLU operation.

    % Convolution.
    W = parameters.Conv.Weights;
    B = parameters.Conv.Bias;
    dlY = dlconv(dlX,W,B);
    
    % Instance normalization.
    offset = parameters.InstanceNorm.Offset;
    scale = parameters.InstanceNorm.Scale;
    dlY = instancenorm(dlY,offset,scale);
    
    % ReLU.
    dlY = relu(dlY);
end

Вспомогательные функции

normalizePointCloud функция извлекает данные точек X, Y, Z из входных данных и нормализует данные между 0 и 1. Функция возвращает нормализованные данные X, Y, Z.

function data = normalizePointCloud(data)
    if ~iscell(data)
        data = {data};
    end
    
    numObservations = size(data,1);
    for i = 1:numObservations
        % Scale points between 0 and 1.
        xlim = data{i,1}.XLimits;
        ylim = data{i,1}.YLimits;
        zlim = data{i,1}.ZLimits;
        
        xyzMin = [xlim(1) ylim(1) zlim(1)];
        xyzDiff = [diff(xlim) diff(ylim) diff(zlim)];
        
        tmp = (data{i,1}.Location - xyzMin) ./ xyzDiff;
        % Convert the data to SCSB format.
        data{i,1} = permute(tmp,[1,3,2]);
    end
end

function data = onehotEncodeLabels(data,classNames)
    numObservations = size(data,1);
    for i = 1:numObservations
       labels = data{i,1}';
       encodedLabels = onehotencode(labels,2,'ClassNames',classNames);
       data{i,1} = permute(encodedLabels,[1,3,2]);
    end
end

`prepareForPrediction` Функция

prepareForPrediction используется для применения определяемой пользователем функции к вложенным данным структуры. Эта функция используется для перемещения параметров модели и данных о состоянии в графический процессор.

function p = prepareForPrediction(p,fcn)

for i = 1:numel(p)
    p(i) = structfun(@(x)invoke(fcn,x),p(i),'UniformOutput',0);
end

    function data = invoke(fcn,data)
        if isstruct(data)
            data = prepareForPrediction(data,fcn);
        else
            data = fcn(data);
        end
    end
end

% Move data to the GPU.
function x = toDevice(x,useGPU)
if useGPU
    x = gpuArray(x);
end
end

`selectPoints` Функция

selectPoints функция выполняет выборку требуемого количества точек. Если облако точек содержит больше требуемого количества точек, функция использует pcdownsample для случайного выбора точек. В противном случае функция реплицирует данные для получения требуемого количества точек.

function ptCloudOut = selectPoints(ptCloud,numPoints) 
    if ptCloud.Count > numPoints
        ind = 1:ptCloud.Count;
    else    
        replicationFactor = ceil(numPoints/ptCloud.Count);
        ind = repmat(1:ptCloud.Count,1,replicationFactor);
    end 
    
     ptCloudOut = select(ptCloud,ind(1:numPoints));
end

labelColorbar функция добавляет панель цветов к текущей оси. Панель цветов отформатирована для отображения имен классов с цветом.

function labelColorbar(ax, cmap, classNames)
    colormap(ax, cmap);
    
    % Add colorbar to current figure.
    c = colorbar(ax);
    c.Color = 'w';
        
    % Center tick labels and use class names for tick marks.
    numClasses = size(classNames, 1);
    c.Ticks = 1:1:numClasses;
    c.TickLabels = classNames; 

    % Remove tick mark.
    c.TickLength = 0;
end

function cmap = labelColorMap()
% Colormap for the original classes.
cmap = [[10,10,10];
        [0,0,255];
        [0,255,0];
        [255,192,203];
        [255,255,0];
        [255,0,255];
        [255,165,0];
        [139,0,150];
        [255,0,0]];
cmap = cmap./255;
end

function lossPlotter = configureTrainingProgressPlotter(f)
% The configureTrainingProgressPlotter function configures training
% progress plots for various losses.
    figure(f);
    clf
    ylabel('Total Loss');
    xlabel('Iteration');
    lossPlotter = animatedline;
end

Ссылки

[1] Варни, Нина, Виджаян К. Асари и Куинн Грейлинг. «DALES: Крупномасштабный набор воздушных данных LiDAR для семантической сегментации». ArXiv:2004.11985 [Cs, Stat], 14 апреля 2020 года. https://arxiv.org/abs/2004.11985.

[2] Ци, Чарльз Р., Ли И, Хао Су и Леонидас Дж. Гибас. «StartNet++: Глубокое изучение иерархических функций на наборах точек в метрическом пространстве». ArXiv:1706.02413 [Cs], 7 июня 2017 года. https://arxiv.org/abs/1706.02413.

Документация

Семантическая сегментация Aerial Lidar с использованием Deep Learning

Загрузить данные DALES

Данные предварительной обработки

Создание объектов хранилища данных для обучения

Определение модели «Net++»

Определите PointNet ++ образцовые параметры

Укажите параметры обучения

Модель поезда

Облако воздушных точек сегмента

Оценка сети

Градиенты модели

Функция модели «»/Net + +

Функции инициализации параметров модели

Инициализация общей многослойной функции перцептрона

Отбор проб и группирование

Наиболее дальняя точка выборки

Запрос точки шарика

Распространение элементов

Общая многослойная функция перцептрона

Функция перцептрона

Вспомогательные функции

`prepareForPrediction` Функция

`selectPoints` Функция

Ссылки

Документация по панели инструментов Lidar

Поддержка

Документация

Семантическая сегментация Aerial Lidar с использованием Deep Learning

Загрузить данные DALES

Данные предварительной обработки

Создание объектов хранилища данных для обучения

Определение модели «Net++»

Определите PointNet ++ образцовые параметры

Укажите параметры обучения

Модель поезда

Облако воздушных точек сегмента

Оценка сети

Градиенты модели

Функция модели «»/Net + +

Функции инициализации параметров модели

Инициализация общей многослойной функции перцептрона

Отбор проб и группирование

Наиболее дальняя точка выборки

Запрос точки шарика

Распространение элементов

Общая многослойная функция перцептрона

Функция перцептрона

Вспомогательные функции

prepareForPrediction Функция

selectPoints Функция

Ссылки

Документация по панели инструментов Lidar

Поддержка

`prepareForPrediction` Функция

`selectPoints` Функция