Классификация модуляций с глубоким обучением

Этот пример использование:

В этом примере вы генерируете синтетический продукт, поврежденные каналом формы волны. Используя сгенерированные формы волны как данные тренировки, вы обучаете CNN классификации модуляций. Вы затем тестируете CNN с оборудованием программно определяемого радио (SDR) и беспроводными сигналами.

Предскажите тип модуляции Используя CNN

Обученный CNN в этом примере распознает эти цифровые восемь и три аналоговых типа модуляции:

Бинарное манипулирование сдвига фазы (BPSK)
Манипулирование сдвига фазы Quadrature (QPSK)
8-ary сдвиг фазы, включающий (8-PSK)
16-ary квадратурная амплитудная (16-QAM) модуляция
64-ary квадратурная амплитудная (64-QAM) модуляция
4-ary импульсная амплитудная модуляция (PAM4)
Гауссово манипулирование сдвига частоты (GFSK)
Непрерывное манипулирование сдвига частоты фазы (CPFSK)
Широковещательно передайте FM (B-FM)
Двойная амплитудная модуляция боковой полосы (DSB-AM)
Одна амплитудная модуляция боковой полосы (SSB-AM)

modulationTypes = categorical(["BPSK", "QPSK", "8PSK", ...
  "16QAM", "64QAM", "PAM4", "GFSK", "CPFSK", ...
  "B-FM", "DSB-AM", "SSB-AM"]);

Во-первых, загрузите обучивший сеть. Для получения дополнительной информации на сетевом обучении, смотрите раздел Training a CNN.

load trainedModulationClassificationNetwork
trainedNet

trainedNet = 
  SeriesNetwork with properties:

    Layers: [28x1 nnet.cnn.layer.Layer]

Обученный CNN берет 1 024 поврежденных каналом выборки и предсказывает тип модуляции каждого кадра. Сгенерируйте несколько кадров BPSK, которым повреждают с Rician многопутевое исчезновение, центральная частота и дрейф времени выборки и AWGN. Используйте функцию randi, чтобы сгенерировать случайные биты, функция pskmod к BPSK - модулируют биты, функция rcosdesign, чтобы разработать квадратный корень повысила формирующий фильтр импульса косинуса, и функция filter к импульсу формирует символы. Затем используйте CNN, чтобы предсказать тип модуляции кадров.

% Set the random number generator to a known state to be able to regenerate
% the same frames every time the simulation is run
rng(123456)
% Random bits
d = randi([0 1],1024,1);
% BPSK modulation
syms = pskmod(d,2);
% Square-root raised cosine filter
filterCoeffs = rcosdesign(0.35,4,8);
tx = filter(filterCoeffs,1,upsample(syms,8));
% Channel
channel = helperModClassTestChannel(...
  'SampleRate',200e3, ...
  'SNR',30, ...
  'PathDelays',[0 1.8 3.4] / 200e3, ...
  'AveragePathGains',[0 -2 -10], ...
  'KFactor',4, ...
  'MaximumDopplerShift',4, ...
  'MaximumClockOffset',5, ...
  'CenterFrequency',900e6);
rx = channel(tx);
% Plot transmitted and received signals
scope = dsp.TimeScope(2,200e3,'YLimits',[-1 1],'ShowGrid',true,...
  'LayoutDimensions',[2 1],'TimeSpan',45e-3);
scope(tx,rx)

% Frame generation for classification
unknownFrames = getNNFrames(rx,'Unknown');
% Classification
[prediction1,score1] = classify(trainedNet,unknownFrames);

Возвратите прогнозы классификатора, которые походят на трудные решения. Сеть правильно идентифицирует кадры как кадры BPSK. Для получения дополнительной информации на генерации модулируемых сигналов, см. Приложение: Модуляторы.

prediction1

prediction1 = 7x1 categorical array
     BPSK 
     BPSK 
     BPSK 
     BPSK 
     BPSK 
     BPSK 
     BPSK

Классификатор также возвращает вектор музыки к каждому кадру. Счет соответствует вероятности, что каждый кадр имеет предсказанный тип модуляции. Постройте очки.

plotScores(score1,modulationTypes)

Затем, используйте CNN, чтобы классифицировать кадры PAM4.

% Random bits
d = randi([0 3], 1024, 1);
% PAM4 modulation
syms = pammod(d,4);
% Square-root raised cosine filter
filterCoeffs = rcosdesign(0.35, 4, 8);
tx = filter(filterCoeffs, 1, upsample(syms,8));
% Channel
rx = channel(tx);
% Plot transmitted and received signals
scope = dsp.TimeScope(2,200e3,'YLimits',[-2 2],'ShowGrid',true,...
  'LayoutDimensions',[2 1],'TimeSpan',45e-3);
scope(tx,rx)

% Frame generation for classification
unknownFrames = getNNFrames(rx,'Unknown');
% Classification
[estimate2,score2] = classify(trainedNet,unknownFrames);
estimate2

estimate2 = 7x1 categorical array
     PAM4 
     PAM4 
     PAM4 
     PAM4 
     PAM4 
     PAM4 
     PAM4

plotScores(score2,modulationTypes)

Прежде чем мы сможем использовать CNN для классификации модуляций или любую другую задачу, мы сначала должны обучить сеть с известным (или маркированный) данные. Первая часть этого примера показывает, как использовать функции Communications Toolbox, такие как модуляторы, фильтры, и нарушения канала, чтобы сгенерировать синтетические данные тренировки. Вторая часть фокусируется на определении, обучении и тестировании CNN для задачи классификации модуляций. Третья часть проверяет производительность сети с по воздушным сигналам с помощью платформы программно определяемого радио (SDR) ADALM-PLUTO.

Генерация формы волны для обучения

Сгенерируйте 10 000 кадров для каждого типа модуляции, где 80% используются для обучения, 10% используется для валидации, и 10% используется для тестирования. Мы используем кадры обучения и валидации во время сетевой учебной фазы. Итоговая точность классификации получена с помощью тестовых кадров. Каждый кадр является 1 024 выборками долго и имеет частоту дискретизации 200 кГц. Для цифровых типов модуляции восемь выборок представляют символ. Сеть принимает каждое решение на основе одного кадров, а не на нескольких последовательных кадрах (как в видео). Примите центральную частоту 900 МГц и 100 МГц для цифровых и аналоговых типов модуляции, соответственно.

numFramesPerModType = 10000;
percentTrainingSamples = 80;
percentValidationSamples = 10;
percentTestSamples = 10;

sps = 8;                % Samples per symbol
spf = 1024;             % Samples per frame
symbolsPerFrame = spf / sps;
fs = 200e3;             % Sample rate
fc = [900e6 100e6];     % Center frequencies

Создайте нарушения канала

Передайте каждый кадр через канал с

AWGN
Rician многопутевое исчезновение
Смещение часов, приводящее к центральному смещению частоты и дрейфу времени выборки

Поскольку сеть в этом примере принимает решения на основе одного кадров, каждый кадр должен пройти через независимый канал.

AWGN

Канал добавляет AWGN с ОСШ 30 дБ. Поскольку кадры нормированы, шумовое стандартное отклонение, может быть вычислен как

SNR = 30;
std = sqrt(10.^(-SNR/10))

std = 0.0316

Реализуйте канал с помощью comm.AWGNChannel,

awgnChannel = comm.AWGNChannel(...
  'NoiseMethod', 'Signal to noise ratio (SNR)', ...
  'SignalPower', 1, ...
  'SNR', SNR)

awgnChannel = 
  comm.AWGNChannel with properties:

     NoiseMethod: 'Signal to noise ratio (SNR)'
             SNR: 30
     SignalPower: 1
    RandomStream: 'Global stream'

Многопутевой Rician

Канал передает сигналы через Rician многопутевой исчезающий канал с помощью Системного объекта comm.RicianChannel. Примите профиль задержки [0 1.8 3.4] выборки с соответствующими средними усилениями пути [0 - 2 - 10] дБ. K-фактор равняется 4, и максимальный эффект Доплера составляет 4 Гц, который эквивалентен скорости обхода на уровне 900 МГц. Реализуйте канал со следующими настройками.

multipathChannel = comm.RicianChannel(...
  'SampleRate', fs, ...
  'PathDelays', [0 1.8 3.4]/fs, ...
  'AveragePathGains', [0 -2 -10], ...
  'KFactor', 4, ...
  'MaximumDopplerShift', 4)

multipathChannel = 
  comm.RicianChannel with properties:

                SampleRate: 200000
                PathDelays: [0 9.0000e-06 1.7000e-05]
          AveragePathGains: [0 -2 -10]
        NormalizePathGains: true
                   KFactor: 4
    DirectPathDopplerShift: 0
    DirectPathInitialPhase: 0
       MaximumDopplerShift: 4
           DopplerSpectrum: [1x1 struct]

  Show all properties

Смещение часов

Смещение часов происходит из-за погрешностей внутренних источников часов передатчиков и получателей. Синхронизируйте причины смещения центральная частота, которая используется к downconvert сигнал к основной полосе и уровень выборки цифро-аналогового преобразователя, чтобы отличаться от идеальных значений. Средство моделирования канала использует фактор смещения часов $C$ , выраженный как $C = 1 + \frac{Δ_{часы}}{10^{6}}$ , где $Δ_{часы}$ смещение часов. Для каждого кадра канал генерирует случайное $Δ_{часы}$ значение от равномерно распределенного множества значений в области значений [ ${- \max Δ}_{часы}$ ${\max Δ}_{часы}$ ], где ${\max Δ}_{часы}$ максимальное смещение часов. Смещение часов измеряется в частях на миллион (ppm). В данном примере примите максимальное смещение часов 5 страниц в минуту.

maxDeltaOff = 5;
deltaOff = (rand()*2*maxDeltaOff) - maxDeltaOff;
C = 1 + (deltaOff/1e6);

Смещение частоты

Подвергните каждый кадр смещению частоты на основе фактора смещения часов $C$ и центральная частота. Реализованный канал с помощью comm.PhaseFrequencyOffset.

offset = -(C-1)*fc(1);
frequencyShifter = comm.PhaseFrequencyOffset(...
  'SampleRate', fs, ...
  'FrequencyOffset', offset)

frequencyShifter = 
  comm.PhaseFrequencyOffset with properties:

              PhaseOffset: 0
    FrequencyOffsetSource: 'Property'
          FrequencyOffset: -2.4332e+03
               SampleRate: 200000

Выборка смещения уровня

Подвергните каждый кадр смещению уровня выборки на основе фактора смещения часов $C$ . Реализуйте канал с помощью функции interp1, чтобы передискретизировать кадр на новом уровне $C \times f_{s}$ .

Объединенный канал

Используйте объект helperModClassTestChannel применить все три нарушения канала к кадрам.

channel = helperModClassTestChannel(...
  'SampleRate', fs, ...
  'SNR', SNR, ...
  'PathDelays', [0 1.8 3.4] / fs, ...
  'AveragePathGains', [0 -2 -10], ...
  'KFactor', 4, ...
  'MaximumDopplerShift', 4, ...
  'MaximumClockOffset', 5, ...
  'CenterFrequency', 900e6)

channel = 
  helperModClassTestChannel with properties:

                    SNR: 30
        CenterFrequency: 900000000
             SampleRate: 200000
             PathDelays: [0 9.0000e-06 1.7000e-05]
       AveragePathGains: [0 -2 -10]
                KFactor: 4
    MaximumDopplerShift: 4
     MaximumClockOffset: 5

Можно просмотреть основную информацию о канале с помощью информационной функции объекта.

chInfo = info(channel)

chInfo = struct with fields:
               ChannelDelay: 6
     MaximumFrequencyOffset: 4500
    MaximumSampleRateOffset: 1

Генерация формы волны

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

% Set the random number generator to a known state to be able to regenerate
% the same frames every time the simulation is run
rng(1235)
tic

numModulationTypes = length(modulationTypes);

channelInfo = info(channel);
frameStore = helperModClassFrameStore(...
  numFramesPerModType*numModulationTypes,spf,modulationTypes);
transDelay = 50;
for modType = 1:numModulationTypes
  fprintf('%s - Generating %s frames\n', ...
    datestr(toc/86400,'HH:MM:SS'), modulationTypes(modType))
  numSymbols = (numFramesPerModType / sps);
  dataSrc = getSource(modulationTypes(modType), sps, 2*spf, fs);
  modulator = getModulator(modulationTypes(modType), sps, fs);
  if contains(char(modulationTypes(modType)), {'B-FM','DSB-AM','SSB-AM'})
    % Analog modulation types use a center frequency of 100 MHz
    channel.CenterFrequency = 100e6;
  else
    % Digital modulation types use a center frequency of 900 MHz
    channel.CenterFrequency = 900e6;
  end
  
  for p=1:numFramesPerModType
    % Generate random data
    x = dataSrc();
    
    % Modulate
    y = modulator(x);
    
    % Pass through independent channels
    rxSamples = channel(y);
    
    % Remove transients from the beginning, trim to size, and normalize
    frame = helperModClassFrameGenerator(rxSamples, spf, spf, transDelay, sps);
    
    % Add to frame store
    add(frameStore, frame, modulationTypes(modType));
  end
end

00:00:00 - Generating BPSK frames
00:01:22 - Generating QPSK frames
00:02:33 - Generating 8PSK frames
00:03:45 - Generating 16QAM frames
00:05:21 - Generating 64QAM frames
00:06:43 - Generating PAM4 frames
00:07:36 - Generating GFSK frames
00:08:44 - Generating CPFSK frames
00:09:50 - Generating B-FM frames
00:10:50 - Generating DSB-AM frames
00:11:35 - Generating SSB-AM frames

Затем разделите кадры на обучение, валидацию и тестовые данные. По умолчанию frameStore помещает основополосные выборки I/Q в строки в выходных кадрах. Выходные кадры имеют размер [2xspfx1xN], где первая строка является синфазными выборками, и вторая строка является квадратурными выборками.

[mcfsTraining,mcfsValidation,mcfsTest] = splitData(frameStore,...
  [percentTrainingSamples,percentValidationSamples,percentTestSamples]);
[rxTraining,rxTrainingLabel] = get(mcfsTraining);
[rxValidation,rxValidationLabel] = get(mcfsValidation);
[rxTest,rxTestLabel] = get(mcfsTest);

% Plot the amplitude of the real and imaginary parts of the example frames
% against the sample number
plotTimeDomain(rxTest,rxTestLabel,modulationTypes,fs)

% Plot a spectrogram of the example frames
plotSpectrogram(rxTest,rxTestLabel,modulationTypes,fs,sps)

Избегайте неустойчивости класса в своих данных тренировки путем обеспечения равномерного распределения меток (типы модуляции). Постройте дистрибутивы метки, чтобы проверять, равномерно распределены ли сгенерированные метки.

% Plot the label distributions
figure
subplot(3,1,1)
histogram(rxTrainingLabel)
title("Training Label Distribution")
subplot(3,1,2)
histogram(rxValidationLabel)
title("Validation Label Distribution")
subplot(3,1,3)
histogram(rxTestLabel)
title("Test Label Distribution")

Обучите CNN

Этот пример использует CNN, который состоит из шести слоев свертки и одного полносвязного слоя. Каждый слой свертки кроме последнего сопровождается пакетным слоем нормализации, исправил линейный модуль (ReLU) слой активации и макс. объединение слоя. В последнем слое свертки макс. слой объединения заменяется средним слоем объединения. Выходной слой имеет softmax активацию. Для руководства проектирования сети смотрите Советы Глубокого обучения и Приемы.

dropoutRate = 0.5;
numModTypes = numel(modulationTypes);
netWidth = 1;
filterSize = [1 sps];
poolSize = [1 2];
modClassNet = [
  imageInputLayer([2 spf 1], 'Normalization', 'none', 'Name', 'Input Layer')
  
  convolution2dLayer(filterSize, 16*netWidth, 'Padding', 'same', 'Name', 'CNN1')
  batchNormalizationLayer('Name', 'BN1')
  reluLayer('Name', 'ReLU1')
  maxPooling2dLayer(poolSize, 'Stride', [1 2], 'Name', 'MaxPool1')
  
  convolution2dLayer(filterSize, 24*netWidth, 'Padding', 'same', 'Name', 'CNN2')
  batchNormalizationLayer('Name', 'BN2')
  reluLayer('Name', 'ReLU2')
  maxPooling2dLayer(poolSize, 'Stride', [1 2], 'Name', 'MaxPool2')
  
  convolution2dLayer(filterSize, 32*netWidth, 'Padding', 'same', 'Name', 'CNN3')
  batchNormalizationLayer('Name', 'BN3')
  reluLayer('Name', 'ReLU3')
  maxPooling2dLayer(poolSize, 'Stride', [1 2], 'Name', 'MaxPool3')
  
  convolution2dLayer(filterSize, 48*netWidth, 'Padding', 'same', 'Name', 'CNN4')
  batchNormalizationLayer('Name', 'BN4')
  reluLayer('Name', 'ReLU4')
  maxPooling2dLayer(poolSize, 'Stride', [1 2], 'Name', 'MaxPool4')
  
  convolution2dLayer(filterSize, 64*netWidth, 'Padding', 'same', 'Name', 'CNN5')
  batchNormalizationLayer('Name', 'BN5')
  reluLayer('Name', 'ReLU5')
  maxPooling2dLayer(poolSize, 'Stride', [1 2], 'Name', 'MaxPool5')
  
  convolution2dLayer(filterSize, 96*netWidth, 'Padding', 'same', 'Name', 'CNN6')
  batchNormalizationLayer('Name', 'BN6')
  reluLayer('Name', 'ReLU6')
  
  averagePooling2dLayer([1 ceil(spf/32)], 'Name', 'AP1')
  
  fullyConnectedLayer(numModTypes, 'Name', 'FC1')
  softmaxLayer('Name', 'SoftMax')
  
  classificationLayer('Name', 'Output') ]

modClassNet = 
  28x1 Layer array with layers:

     1   'Input Layer'   Image Input             2x1024x1 images
     2   'CNN1'          Convolution             16 1x8 convolutions with stride [1  1] and padding 'same'
     3   'BN1'           Batch Normalization     Batch normalization
     4   'ReLU1'         ReLU                    ReLU
     5   'MaxPool1'      Max Pooling             1x2 max pooling with stride [1  2] and padding [0  0  0  0]
     6   'CNN2'          Convolution             24 1x8 convolutions with stride [1  1] and padding 'same'
     7   'BN2'           Batch Normalization     Batch normalization
     8   'ReLU2'         ReLU                    ReLU
     9   'MaxPool2'      Max Pooling             1x2 max pooling with stride [1  2] and padding [0  0  0  0]
    10   'CNN3'          Convolution             32 1x8 convolutions with stride [1  1] and padding 'same'
    11   'BN3'           Batch Normalization     Batch normalization
    12   'ReLU3'         ReLU                    ReLU
    13   'MaxPool3'      Max Pooling             1x2 max pooling with stride [1  2] and padding [0  0  0  0]
    14   'CNN4'          Convolution             48 1x8 convolutions with stride [1  1] and padding 'same'
    15   'BN4'           Batch Normalization     Batch normalization
    16   'ReLU4'         ReLU                    ReLU
    17   'MaxPool4'      Max Pooling             1x2 max pooling with stride [1  2] and padding [0  0  0  0]
    18   'CNN5'          Convolution             64 1x8 convolutions with stride [1  1] and padding 'same'
    19   'BN5'           Batch Normalization     Batch normalization
    20   'ReLU5'         ReLU                    ReLU
    21   'MaxPool5'      Max Pooling             1x2 max pooling with stride [1  2] and padding [0  0  0  0]
    22   'CNN6'          Convolution             96 1x8 convolutions with stride [1  1] and padding 'same'
    23   'BN6'           Batch Normalization     Batch normalization
    24   'ReLU6'         ReLU                    ReLU
    25   'AP1'           Average Pooling         1x32 average pooling with stride [1  1] and padding [0  0  0  0]
    26   'FC1'           Fully Connected         11 fully connected layer
    27   'SoftMax'       Softmax                 softmax
    28   'Output'        Classification Output   crossentropyex

Используйте функцию analyzeNetwork, чтобы отобразить интерактивную визуализацию сетевой архитектуры, обнаружить ошибки и проблемы с сетью, и получить подробную информацию о сетевых слоях. Эта сеть имеет 98,323 learnables.

analyzeNetwork(modClassNet)

Затем сконфигурируйте TrainingOptionsSGDM, чтобы использовать решатель SGDM с мини-пакетным размером 256. Определите максимальный номер эпох к 12, поскольку большее число эпох не обеспечивает дальнейшего учебного преимущества. Обучите сеть на графическом процессоре путем установки среды выполнения на 'gpu'. Установите начальный темп обучения на $2 x 10^{- 2}$ . Уменьшайте темп обучения фактором 10 каждых 9 эпох. Установите 'Plots' на 'training-progress' строить учебный прогресс. На Титане NVIDIA Xp графический процессор сеть занимает приблизительно 25 минут, чтобы обучаться..

maxEpochs = 12;
miniBatchSize = 256;
validationFrequency = floor(numel(rxTrainingLabel)/miniBatchSize);
options = trainingOptions('sgdm', ...
  'InitialLearnRate',2e-2, ...
  'MaxEpochs',maxEpochs, ...
  'MiniBatchSize',miniBatchSize, ...
  'Shuffle','every-epoch', ...
  'Plots','training-progress', ...
  'Verbose',false, ...
  'ValidationData',{rxValidation,rxValidationLabel}, ...
  'ValidationFrequency',validationFrequency, ...
  'LearnRateSchedule', 'piecewise', ...
  'LearnRateDropPeriod', 9, ...
  'LearnRateDropFactor', 0.1, ...
  'ExecutionEnvironment', 'gpu');

Или обучите сеть или используйте, уже обучил сеть. По умолчанию этот пример использует обучивший сеть.

trainNow = false;
if trainNow == true
  fprintf('%s - Training the network\n', datestr(toc/86400,'HH:MM:SS'))
  trainedNet = trainNetwork(rxTraining,rxTrainingLabel,modClassNet,options);
else
  load trainedModulationClassificationNetwork
end

Когда график учебного прогресса показывает, сеть сходится приблизительно в 12 эпох почти с 90%-й точностью.

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

fprintf('%s - Classifying test frames\n', datestr(toc/86400,'HH:MM:SS'))

00:12:50 - Classifying test frames

rxTestPred = classify(trainedNet,rxTest);
testAccuracy = mean(rxTestPred == rxTestLabel);
disp("Test accuracy: " + testAccuracy*100 + "%")

Test accuracy: 90.6455%

Постройте матрицу беспорядка для тестовых кадров. Когда матрица показывает, сеть путает 16-QAM и 64-QAM кадры. Эта проблема ожидается, поскольку каждый кадр несет только 128 символов, и 16-QAM подмножество 64-QAM. Сеть также путает QPSK и кадры 8-PSK, поскольку совокупности этих типов модуляции выглядят подобными когда-то вращаемый фазой из-за исчезающего канала и смещения частоты.

figure
cm = confusionchart(rxTestLabel, rxTestPred);
cm.Title = 'Confusion Matrix for Test Data';
cm.RowSummary = 'row-normalized';
cm.Parent.Position = [cm.Parent.Position(1:2) 740 424];

I/Q как Страницы

По умолчанию frameStore помещает основополосные выборки I/Q в строки в 2D массиве. Поскольку фильтры свертки имеют размер [1xsps], сверточный синфазный процесс слоев и квадратурные компоненты независимо. Только в полносвязном слое информация от синфазных компонентов и квадратурных объединенных компонентов.

Альтернатива должна представлять выборки I/Q как трехмерный массив размера [1xSPFx2], который помещает синфазные компоненты и квадратурные компоненты в 3-й размерности (страницы). Этот подход смешивает информацию во мне и Q даже в сверточных слоях и делает лучше нас информации о фазе. Установите свойство 'OutputFormat' хранилища кадра к "IQAsPages" и размер входного слоя к [1xSPFx2].

% Put the data in [1xspfx2] format
mcfsTraining.OutputFormat = "IQAsPages";
[rxTraining,rxTrainingLabel] = get(mcfsTraining);
mcfsValidation.OutputFormat = "IQAsPages";
[rxValidation,rxValidationLabel] = get(mcfsValidation);
mcfsTest.OutputFormat = "IQAsPages";
[rxTest,rxTestLabel] = get(mcfsTest);

% Set the options
options = trainingOptions('sgdm', ...
  'InitialLearnRate',2e-2, ...
  'MaxEpochs',maxEpochs, ...
  'MiniBatchSize',miniBatchSize, ...
  'Shuffle','every-epoch', ...
  'Plots','training-progress', ...
  'Verbose',false, ...
  'ValidationData',{rxValidation,rxValidationLabel}, ...
  'ValidationFrequency',validationFrequency, ...
  'LearnRateSchedule', 'piecewise', ...
  'LearnRateDropPeriod', 9, ...
  'LearnRateDropFactor', 0.1, ...
  'ExecutionEnvironment', 'gpu');

% Set the input layer input size to [1xspfx2]
modClassNet(1) = ...
  imageInputLayer([1 spf 2], 'Normalization', 'none', 'Name', 'Input Layer');

Следующий анализ сети показывает, что каждый сверточный фильтр имеет размерность 1x8x2, который позволяет сверточному слою использовать и меня и данные Q в вычислении одного фильтра вывод.

analyzeNetwork(modClassNet)

% Train or load the pretrained modified network
trainNow = false;
if trainNow == true
  fprintf('%s - Training the network\n', datestr(toc/86400,'HH:MM:SS'))
  trainedNet = trainNetwork(rxTraining,rxTrainingLabel,modClassNet,options);
else
  load trainedModulationClassificationNetwork2
end

Когда график учебного прогресса показывает, сеть сходится приблизительно в 12 эпох больше чем с 95%-й точностью. При представлении компонентов I/Q, когда страницы вместо строк могут улучшить точность сети приблизительно на 5%.

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

fprintf('%s - Classifying test frames\n', datestr(toc/86400,'HH:MM:SS'))

00:13:07 - Classifying test frames

rxTestPred = classify(trainedNet,rxTest);
testAccuracy = mean(rxTestPred == rxTestLabel);
disp("Test accuracy: " + testAccuracy*100 + "%")

Test accuracy: 95.2545%

Постройте матрицу беспорядка для тестовых кадров. Когда матрица показывает, представляя компоненты I/Q, когда страницы вместо строк существенно увеличивают способность сети точно дифференцировать 16-QAM и 64-QAM кадры и кадры 8-PSK и QPSK.

figure
cm = confusionchart(rxTestLabel, rxTestPred);
cm.Title = 'Confusion Matrix for Test Data';
cm.RowSummary = 'row-normalized';
cm.Parent.Position = [cm.Parent.Position(1:2) 740 424];

Протестируйте с SDR

Проверьте производительность обучившего сеть с беспроводными сигналами с помощью функции sdrTest. Чтобы выполнить этот тест, у вас должно быть два радио ADALM-PLUTO и Пакет Поддержки Communications Toolbox для Радио ADALM-PLUTO. sdrTest использует те же функции модуляции что касается генерации учебных сигналов и передает их использующий радио ADALM-PLUTO. Вместо того, чтобы моделировать канал, получите поврежденные каналом сигналы с другим радио ADALM-PLUTO. Используйте обучивший сеть с той же функцией classify, чтобы предсказать тип модуляции. Сеть достигает 99%-й общей точности, где два радио являются стационарными и разделяются приблизительно на 2 фута.

if isPlutoSDRInstalled() == true
  radios = findPlutoRadio();
  if length(radios) >= 2
    sdrTest();
  end
end

Дальнейшее исследование

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

Коммуникационный Тулбокс обеспечивает намного больше типов модуляции и нарушений канала. Для получения дополнительной информации смотрите Модуляцию (Communications Toolbox) и Модели Канала Распространения (Communications Toolbox) разделы. Можно также добавить стандартные определенные сигналы с LTE Toolbox, WLAN Toolbox и 5G Toolbox. Можно также добавить радарные сигналы с Phased Array System Toolbox.

Приложение: раздел Modulators обеспечивает, функции MATLAB раньше генерировали модулируемые сигналы. Можно также исследовать следующие функции и Системные объекты для получения дополнительной информации:

Ссылки

О'Ши, T. J. Дж. Коргэн и Т. К. Клэнси. "Сверточные Радио-Сети Распознавания Модуляции". Предварительно распечатайте, представленный 10 июня 2016. https://arxiv.org/abs/1602.04105
О'Ши, T. J. Т. Рой и Т. К. Клэнси. "Беспроводная основанная на глубоком обучении Радио-Классификация Сигнала". Журнал IEEE Выбранных Тем в Обработке сигналов. Издание 12, Номер 1, 2018, стр 168–179.
Лю, X., Д. Янг и А. Э. Джамаль. "Архитектура Глубокой нейронной сети для Классификации Модуляций". Предварительно распечатайте, представленный 5 января 2018. https://arxiv.org/abs/1712.00443v3

Приложение: Функции помощника

function testAccuracy = sdrTest()
%sdrTest Test CNN performance with over-the-air signals
%   A = sdrTest sends test frames from one to another ADALM-PLUTO radio,
%   performs classification with the trained network, and returns the overall 
%   classification accuracy. Transmitting radio uses transmit-repeat capability
%   to send the same waveform repeatedly without loading the main loop.

modulationTypes = categorical(["BPSK", "QPSK", "8PSK", ...
  "16QAM", "64QAM", "PAM4", "GFSK", "CPFSK", "B-FM"]);
load trainedModulationClassificationNetwork2 trainedNet
numFramesPerModType = 100;

sps = 8;                % Samples per symbol
spf = 1024;             % Samples per frame
fs = 200e3;             % Sample rate

txRadio = sdrtx('Pluto');
txRadio.RadioID = 'usb:0';
txRadio.CenterFrequency = 900e6;
txRadio.BasebandSampleRate = fs;

rxRadio = sdrrx('Pluto');
rxRadio.RadioID = 'usb:1';
rxRadio.CenterFrequency = 900e6;
rxRadio.BasebandSampleRate = fs;
rxRadio.SamplesPerFrame = spf;
rxRadio.ShowAdvancedProperties = true;
rxRadio.EnableQuadratureCorrection = false;
rxRadio.OutputDataType = 'single';

% Set random number generator to a known state to be able to regenerate
% the same frames every time the simulation is run
rng(1235)
tic

numModulationTypes = length(modulationTypes);
txModType = repmat(modulationTypes(1),numModulationTypes*numFramesPerModType,1);
estimatedModType = repmat(modulationTypes(1),numModulationTypes*numFramesPerModType,1);
frameCnt = 1;
for modType = 1:numModulationTypes
  fprintf('%s - Testing %s frames\n', ...
    datestr(toc/86400,'HH:MM:SS'), modulationTypes(modType))
  dataSrc = getSource(modulationTypes(modType), sps, 2*spf, fs);
  modulator = getModulator(modulationTypes(modType), sps, fs);
  if contains(char(modulationTypes(modType)), {'B-FM'})
    % Analog modulation types use a center frequency of 100 MHz
    txRadio.CenterFrequency = 100e6;
    rxRadio.CenterFrequency = 100e6;
    rxRadio.GainSource = 'Manual';
    rxRadio.Gain = 60;
  else
    % Digital modulation types use a center frequency of 900 MHz
    txRadio.CenterFrequency = 900e6;
    rxRadio.CenterFrequency = 900e6;
    rxRadio.GainSource = 'AGC Fast Attack';
  end
  
  % Start transmitter
  disp('Starting transmitter')
  x = dataSrc();
  y = modulator(x);
  y = y(4*sps+1:end,1);
  maxVal = max(max(abs(real(y))), max(abs(imag(y))));
  y = y *0.8/maxVal;
  % Download waveform signal to radio and repeatedly transmit it over the air
  transmitRepeat(txRadio, complex(y));
  
  disp('Starting receiver and test')
  for p=1:numFramesPerModType
    for frame=1:16
      rx = rxRadio();
    end
    
    frameEnergy = sum(abs(rx).^2);
    rx = rx / sqrt(frameEnergy);
    reshapedRx(1,:,1,1) = real(rx);
    reshapedRx(1,:,2,1) = imag(rx);
    
    % Classify
    txModType(frameCnt) = modulationTypes(modType);
    estimatedModType(frameCnt) = classify(trainedNet, reshapedRx);
    
    frameCnt = frameCnt + 1;
    
    % Pause for 0.1 seconds to get an independent channel (assuming a
    % channel coherence time of less than 0.1 seconds)
    pause(0.1)
  end
  disp('Releasing radios')
  release(txRadio);
  release(rxRadio);
end
testAccuracy = mean(txModType == estimatedModType);
disp("Test accuracy: " + testAccuracy*100 + "%")

figure
cm = confusionchart(txModType, estimatedModType);
cm.Title = 'Confusion Matrix for Test Data';
cm.RowSummary = 'row-normalized';
cm.Parent.Position = [cm.Parent.Position(1:2) 740 424];
end

function modulator = getModulator(modType, sps, fs)
%getModulator Modulation function selector
%   MOD = getModulator(TYPE,SPS,FS) returns the modulator function handle
%   MOD based on TYPE. SPS is the number of samples per symbol and FS is 
%   the sample rate.

switch modType
  case "BPSK"
    modulator = @(x)bpskModulator(x,sps);
  case "QPSK"
    modulator = @(x)qpskModulator(x,sps);
  case "8PSK"
    modulator = @(x)psk8Modulator(x,sps);
  case "16QAM"
    modulator = @(x)qam16Modulator(x,sps);
  case "64QAM"
    modulator = @(x)qam64Modulator(x,sps);
  case "GFSK"
    modulator = @(x)gfskModulator(x,sps);
  case "CPFSK"
    modulator = @(x)cpfskModulator(x,sps);
  case "PAM4"
    modulator = @(x)pam4Modulator(x,sps);
  case "B-FM"
    modulator = @(x)bfmModulator(x, fs);
  case "DSB-AM"
    modulator = @(x)dsbamModulator(x, fs);
  case "SSB-AM"
    modulator = @(x)ssbamModulator(x, fs);
end
end

function src = getSource(modType, sps, spf, fs)
%getSource Source selector for modulation types
%    SRC = getSource(TYPE,SPS,SPF,FS) returns the data source
%    for the modulation type TYPE, with the number of samples 
%    per symbol SPS, the number of samples per frame SPF, and 
%    the sampling frequency FS.

switch modType
  case {"BPSK","GFSK","CPFSK"}
    M = 2;
    src = @()randi([0 M-1],spf/sps,1);
  case {"QPSK","PAM4"}
    M = 4;
    src = @()randi([0 M-1],spf/sps,1);
  case "8PSK"
    M = 8;
    src = @()randi([0 M-1],spf/sps,1);
  case "16QAM"
    M = 16;
    src = @()randi([0 M-1],spf/sps,1);
  case "64QAM"
    M = 64;
    src = @()randi([0 M-1],spf/sps,1);
  case {"B-FM","DSB-AM","SSB-AM"}
    src = @()getAudio(spf,fs);
end
end

function x = getAudio(spf,fs)
%getAudio Audio source for analog modulation types
%    A = getAudio(SPF,FS) returns the audio source A, with the 
%    number of samples per frame SPF, and the sample rate FS.

persistent audioSrc audioRC

if isempty(audioSrc)
  audioSrc = dsp.AudioFileReader('audio_mix_441.wav',...
    'SamplesPerFrame',spf,'PlayCount',inf);
  audioRC = dsp.SampleRateConverter('Bandwidth',30e3,...
    'InputSampleRate',audioSrc.SampleRate,...
    'OutputSampleRate',fs);
  [~,decimFactor] = getRateChangeFactors(audioRC);
  audioSrc.SamplesPerFrame = ceil(spf / fs * audioSrc.SampleRate / decimFactor) * decimFactor;
end

x = audioRC(audioSrc());
x = x(1:spf,1);
end

function frames = getNNFrames(rx,modType)
%getNNFrames Generate formatted frames for neural networks
%   F = getNNFrames(X,MODTYPE) formats the input X, into frames 
%   that can be used with the neural network designed in this 
%   example, and returns the frames in the output F.

frames = helperModClassFrameGenerator(rx,1024,1024,32,8);
frameStore = helperModClassFrameStore(10,1024,categorical({modType}));
add(frameStore,frames,modType);
frames = get(frameStore);
end

function plotScores(score,labels)
%plotScores Plot classification scores of frames
%   plotScores(SCR,LABELS) plots the classification scores SCR as a stacked 
%   bar for each frame. SCR is a matrix in which each row is the score for a 
%   frame.

co = [0.08 0.9 0.49;
  0.52 0.95 0.70;
  0.36 0.53 0.96;
  0.09 0.54 0.67;
  0.48 0.99 0.26;
  0.95 0.31 0.17;
  0.52 0.85 0.95;
  0.08 0.72 0.88;
  0.12 0.45 0.69;
  0.22 0.11 0.49;
  0.65 0.54 0.71];
figure; ax = axes('ColorOrder',co,'NextPlot','replacechildren');
bar(ax,[score; nan(2,11)],'stacked'); legend(categories(labels),'Location','best');
xlabel('Frame Number'); ylabel('Score'); title('Classification Scores')
end

function plotTimeDomain(rxTest,rxTestLabel,modulationTypes,fs)
%plotTimeDomain Time domain plots of frames

numRows = ceil(length(modulationTypes) / 4);
spf = size(rxTest,2);
t = 1000*(0:spf-1)/fs;
if size(rxTest,1) == 2
  IQAsRows = true;
else
  IQAsRows = false;
end
for modType=1:length(modulationTypes)
  subplot(numRows, 4, modType);
  idxOut = find(rxTestLabel == modulationTypes(modType), 1);
  if IQAsRows
    rxI = rxTest(1,:,1,idxOut);
    rxQ = rxTest(2,:,1,idxOut);
  else
    rxI = rxTest(1,:,1,idxOut);
    rxQ = rxTest(1,:,2,idxOut);
  end
  plot(t,squeeze(rxI), '-'); grid on; axis equal; axis square
  hold on
  plot(t,squeeze(rxQ), '-'); grid on; axis equal; axis square
  hold off
  title(string(modulationTypes(modType)));
  xlabel('Time (ms)'); ylabel('Amplitude')
end
end

function plotSpectrogram(rxTest,rxTestLabel,modulationTypes,fs,sps)
%plotSpectrogram Spectrogram of frames

if size(rxTest,1) == 2
  IQAsRows = true;
else
  IQAsRows = false;
end
numRows = ceil(length(modulationTypes) / 4);
for modType=1:length(modulationTypes)
  subplot(numRows, 4, modType);
  idxOut = find(rxTestLabel == modulationTypes(modType), 1);
  if IQAsRows
    rxI = rxTest(1,:,1,idxOut);
    rxQ = rxTest(2,:,1,idxOut);
  else
    rxI = rxTest(1,:,1,idxOut);
    rxQ = rxTest(1,:,2,idxOut);
  end
  rx = squeeze(rxI) + 1i*squeeze(rxQ);
  spectrogram(rx,kaiser(sps),0,1024,fs,'centered');
  title(string(modulationTypes(modType)));
end
h = gcf; delete(findall(h.Children, 'Type', 'ColorBar'))
end

function flag = isPlutoSDRInstalled
%isPlutoSDRInstalled Check if ADALM-PLUTO is installed

spkg = matlabshared.supportpkg.getInstalled;
flag = ~isempty(spkg) && any(contains({spkg.Name},'ADALM-PLUTO','IgnoreCase',true));
end

Приложение: модуляторы

function y = bpskModulator(x,sps)
%bpskModulator BPSK modulator with pulse shaping
%   Y = bpskModulator(X,SPS) BPSK modulates the input X, and returns the 
%   root-raised cosine pulse shaped signal Y. X must be a column vector 
%   of values in the set [0 1]. The root-raised cosine filter has a 
%   roll-off factor of 0.35 and spans four symbols. The output signal 
%   Y has unit power.

persistent filterCoeffs
if isempty(filterCoeffs)
  filterCoeffs = rcosdesign(0.35, 4, sps);
end
% Modulate
syms = pskmod(x,2);
% Pulse shape
y = filter(filterCoeffs, 1, upsample(syms,sps));
end

function y = qpskModulator(x,sps)
%qpskModulator QPSK modulator with pulse shaping
%   Y = qpskModulator(X,SPS) QPSK modulates the input X, and returns the 
%   root-raised cosine pulse shaped signal Y. X must be a column vector 
%   of values in the set [0 3]. The root-raised cosine filter has a 
%   roll-off factor of 0.35 and spans four symbols. The output signal 
%   Y has unit power.

persistent filterCoeffs
if isempty(filterCoeffs)
  filterCoeffs = rcosdesign(0.35, 4, sps);
end
% Modulate
syms = pskmod(x,4,pi/4);
% Pulse shape
y = filter(filterCoeffs, 1, upsample(syms,sps));
end

function y = psk8Modulator(x,sps)
%psk8Modulator 8-PSK modulator with pulse shaping
%   Y = psk8Modulator(X,SPS) 8-PSK modulates the input X, and returns the 
%   root-raised cosine pulse shaped signal Y. X must be a column vector 
%   of values in the set [0 7]. The root-raised cosine filter has a 
%   roll-off factor of 0.35 and spans four symbols. The output signal 
%   Y has unit power.

persistent filterCoeffs
if isempty(filterCoeffs)
  filterCoeffs = rcosdesign(0.35, 4, sps);
end
% Modulate
syms = pskmod(x,8);
% Pulse shape
y = filter(filterCoeffs, 1, upsample(syms,sps));
end

function y = qam16Modulator(x,sps)
%qam16Modulator 16-QAM modulator with pulse shaping
%   Y = qam16Modulator(X,SPS) 16-QAM modulates the input X, and returns the 
%   root-raised cosine pulse shaped signal Y. X must be a column vector 
%   of values in the set [0 15]. The root-raised cosine filter has a 
%   roll-off factor of 0.35 and spans four symbols. The output signal 
%   Y has unit power.

persistent filterCoeffs
if isempty(filterCoeffs)
  filterCoeffs = rcosdesign(0.35, 4, sps);
end
% Modulate and pulse shape
syms = qammod(x,16,'UnitAveragePower',true);
% Pulse shape
y = filter(filterCoeffs, 1, upsample(syms,sps));
end

function y = qam64Modulator(x,sps)
%qam64Modulator 64-QAM modulator with pulse shaping
%   Y = qam64Modulator(X,SPS) 64-QAM modulates the input X, and returns the 
%   root-raised cosine pulse shaped signal Y. X must be a column vector 
%   of values in the set [0 63]. The root-raised cosine filter has a 
%   roll-off factor of 0.35 and spans four symbols. The output signal 
%   Y has unit power.

persistent filterCoeffs
if isempty(filterCoeffs)
  filterCoeffs = rcosdesign(0.35, 4, sps);
end
% Modulate
syms = qammod(x,64,'UnitAveragePower',true);
% Pulse shape
y = filter(filterCoeffs, 1, upsample(syms,sps));
end

function y = pam4Modulator(x,sps)
%pam4Modulator PAM4 modulator with pulse shaping
%   Y = pam4Modulator(X,SPS) PAM4 modulates the input X, and returns the 
%   root-raised cosine pulse shaped signal Y. X must be a column vector 
%   of values in the set [0 3]. The root-raised cosine filter has a 
%   roll-off factor of 0.35 and spans four symbols. The output signal 
%   Y has unit power.

persistent filterCoeffs amp
if isempty(filterCoeffs)
  filterCoeffs = rcosdesign(0.35, 4, sps);
  amp = 1 / sqrt(mean(abs(pammod(0:3, 4)).^2));
end
% Modulate
syms = amp * pammod(x,4);
% Pulse shape
y = filter(filterCoeffs, 1, upsample(syms,sps));
end

function y = gfskModulator(x,sps)
%gfskModulator GFSK modulator
%   Y = gfskModulator(X,SPS) GFSK modulates the input X and returns the 
%   signal Y. X must be a column vector of values in the set [0 1]. The 
%   BT product is 0.35 and the modulation index is 1. The output signal 
%   Y has unit power.

persistent mod meanM
if isempty(mod)
  M = 2;
  mod = comm.CPMModulator(...
    'ModulationOrder', M, ...
    'FrequencyPulse', 'Gaussian', ...
    'BandwidthTimeProduct', 0.35, ...
    'ModulationIndex', 1, ...
    'SamplesPerSymbol', sps);
  meanM = mean(0:M-1);
end
% Modulate
y = mod(2*(x-meanM));
end

function y = cpfskModulator(x,sps)
%cpfskModulator CPFSK modulator
%   Y = cpfskModulator(X,SPS) CPFSK modulates the input X and returns 
%   the signal Y. X must be a column vector of values in the set [0 1]. 
%   the modulation index is 0.5. The output signal Y has unit power.

persistent mod meanM
if isempty(mod)
  M = 2;
  mod = comm.CPFSKModulator(...
    'ModulationOrder', M, ...
    'ModulationIndex', 0.5, ...
    'SamplesPerSymbol', sps);
  meanM = mean(0:M-1);
end
% Modulate
y = mod(2*(x-meanM));
end

function y = bfmModulator(x,fs)
%bfmModulator Broadcast FM modulator
%   Y = bfmModulator(X,FS) broadcast FM modulates the input X and returns
%   the signal Y at the sample rate FS. X must be a column vector of
%   audio samples at the sample rate FS. The frequency deviation is 75 kHz
%   and the pre-emphasis filter time constant is 75 microseconds.

persistent mod
if isempty(mod)
  mod = comm.FMBroadcastModulator(...
    'AudioSampleRate', fs, ...
    'SampleRate', fs);
end
y = mod(x);
end

function y = dsbamModulator(x,fs)
%dsbamModulator Double sideband AM modulator
%   Y = dsbamModulator(X,FS) double sideband AM modulates the input X and
%   returns the signal Y at the sample rate FS. X must be a column vector of
%   audio samples at the sample rate FS. The IF frequency is 50 kHz.

y = ammod(x,50e3,fs);
end

function y = ssbamModulator(x,fs)
%ssbamModulator Single sideband AM modulator
%   Y = ssbamModulator(X,FS) single sideband AM modulates the input X and
%   returns the signal Y at the sample rate FS. X must be a column vector of
%   audio samples at the sample rate FS. The IF frequency is 50 kHz.

y = ssbmod(x,50e3,fs);
end

Смотрите также

trainNetwork | trainingOptions

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