Введение

Статистическая радарная модель

Статистические радарные модели, такие как radarDataGenerator предложите ценные данные для молодой разработки радиолокационных систем. Система моделируется с точки зрения ее показателей производительности, и минимальный набор параметров используется наряду с установленной теорией определить, как эта система будет взаимодействовать со своей средой. Это включает спецификацию основных системных параметров, таких как рабочая частота, пропускная способность и импульсная частота повторения (PRF), наряду с важными метриками как угол и разрешение области значений, ложный сигнальный уровень и вероятность обнаружения. Радарная модель выделила режимы для моностатической, бистатической, и пассивной операции. Первичный выход этой модели является необработанными данными об обнаружении или последовательностью обновлений дорожки. В то время как сигналы временной области не сгенерированы, эмиттер может задать числовой ID для формы волны, которую он использует так, чтобы усиление сигнала approriate и отклонение могли быть рассмотрены. Обнаружения могут быть выведены в радарной системе координат, или в системе координат сценария, если радар получения имеет знание местоположения и ориентацию передатчика. Дорожки сгенерированы с помощью одного из разнообразия отслеживания и отслеживают алгоритмы управления. Скорость передачи данных обработана системной спецификацией частоты обновления, так, чтобы обнаружения были сгенерированы на соответствующем уровне.

Управление сценарием

Модель работает в тандеме с radarScenario генерировать обнаружение или данные о дорожке в зависимости от времени в динамической среде. Этот объект сценария не только управляет временем симуляции, но и может обработать передачу структур данных между объектами в сцене как требуется, когда система координат обнаружения или отслеживает обновление, требуется. Радарный объект смонтирован на платформе, которая обрабатывает динамику, такую как положение, скорость и ориентация в системе координат сценария. Эти платформы сценария также могут содержать информацию о подписи, которая описывает, как платформа появляется к различным типам датчиков. В наших целях, rcsSignature класс будет использоваться, чтобы обеспечить целевые платформы радарным сечением.

Сканирование

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

Сценарий: воздушное наблюдение

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

Радиолокационная система

Используйте центральную частоту на 1 ГГц и пропускную способность 1,5 МГц, чтобы дать нам разрешение области значений приблизительно 100 метров. Чтобы смоделировать один импульс на положение скана, установите частоту обновления, равную желаемому PRF. Установите размер неоднозначности области значений отражать наш выбор PRF и устанавливать наш верхний предел области значений к дважды размеру неоднозначности области значений, чтобы включить обнаружения объектов в первой неоднозначности области значений.

% Set the seed for repeatability
rng(2021);

% System parameters
c = physconst('lightspeed');
fc = 1e9; % Hz
bandwidth = 1.5e6; % Hz
prf = 4e3; % Hz
updateRate = prf;
rangeRes = c/(2*bandwidth); % m
rangeAmb = c/(2*prf); % m
rangeLims = [0 2*rangeAmb];

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

Мы отсканируем +/-20 степени в области азимута и +/-10 степени в области вертикального изменения. Установите FoV на 4 азимута степеней и 8 вертикальных изменений степеней. Наконец, задайте общее количество полных сканов, чтобы симулировать, который продиктует общее время симуляции.

% Scanning parameters
azLimits = [-20 20];
elLimits = [-10 10];
fov = [4;8]; % deg
numFullScans = 20; % number of full scans to simulate

Общее количество точек скана может быть найдено можно следующим образом. Пределы скана обозначают минимальные и максимальные углы сканирования, на которые ссылаются к вектору опорного направления антенны, таким образом, общая площадь, отсканированная в одном направлении, больше степени пределов скана до половины FoV в том направлении. Шаблон всегда является прямоугольной сеткой по az/el.

numScanPointsAz = floor(diff(azLimits)/fov(1)) + 1;
numScanPointsEl = floor(diff(elLimits)/fov(2)) + 1;
numScanPoints = numScanPointsAz*numScanPointsEl;

Спецификация углового разрешения, отдельного от FoV, будет использоваться, который позволяет нам моделировать что-то как моноимпульсная угловая оценка. Это угловое разрешение является для determing нашей способностью различить между целями. и используется, чтобы вывести фактическую угловую точность из нижней границы Крамера-Рао (CRLB) для моноимпульсной угловой оценки с тем разрешением. Точно так же разрешение области значений определяет минимальный интервал в области значений, требуемой отличать две цели, но фактическая точность измерения области значений прибывает из CRLB для оценки области значений. Точный measureables может использоваться путем выключения шума измерения с HasNoise свойство. Как правило, погрешность измерения не переходит в 0, когда ОСШ увеличивается без связанного. Это может быть получено свойствами "смещения", из которых существует один для каждого типа измеримых. В наших целях мы оставим их в их значениях по умолчанию.

Позвольте нашему угловому разрешению быть 1/4-м FoV в каждом направлении, которое позволило бы нам различать максимум между 16 целями точки в FoV в той же области значений.

angRes = fov/4;

Цель ссылки и радарное усиление цикла

Вместо того, чтобы непосредственно задать вещи как степень передачи, шумовая степень и специфические особенности используемого алгоритма обнаружения, эта радарная модель использует концепцию радарного усиления цикла, чтобы перевести целевой RCS и расположиться в ОСШ, который переводится непосредственно в вероятность обнаружения. Радарное усиление цикла является важным свойством любой радиолокационной системы и вычисляется с целью ссылки, которая состоит из диапазона ссылки и RCS и рабочей характеристики приемника (ROC) для данного обнаружения и ложных сигнальных вероятностей.

Мы будем использовать вероятность обнаружения по умолчанию 90% для цели ссылки на уровне 20 км и 0 dBsm. Мы будем использовать ложный сигнальный уровень по умолчанию 1e-6 ложных предупреждений на ячейку разрешения на обновление.

refTgtRange = 20e3;
refTgtRCS = 0;
detProb = 0.9;
faRate = 1e-6;

Построение радара

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

pitch = elLimits(1); % mount the antenna rotated upwards

Включите измерения угла возвышения, ложную сигнальную генерацию и неоднозначности области значений с соответствующим Has... свойство. Чтобы смоделировать радар, который только измеряет область значений и угол, устанавливает HasRangeRate ко лжи. Чтобы получить обнаружения в координатах сценария, наш радар должен знать об ориентации передатчика. Поскольку мы используем моностатический радар, это может быть достигнуто путем добавления функциональности инерционной системы навигации (INS) с HasINS свойство.

% Create radar object
sensorIndex = 1; % a unique identifier is required
radar = radarDataGenerator(sensorIndex,'UpdateRate',updateRate,...
    'DetectionMode','Monostatic','ScanMode','Electronic','TargetReportFormat','Detections','DetectionCoordinates','scenario',...
    'HasElevation',true,'HasINS',true,'HasRangeRate',false,'HasRangeAmbiguities',true,'HasFalseAlarms',true,...
    'CenterFrequency',fc,'Bandwidth',bandwidth,...
    'RangeResolution',rangeRes,'AzimuthResolution',angRes(1),'ElevationResolution',angRes(2),...
    'ReferenceRange',refTgtRange,'ReferenceRCS',refTgtRCS,'DetectionProbability',detProb,'FalseAlarmRate',faRate,...
    'RangeLimits',rangeLims,'MaxUnambiguousRange',rangeAmb,'ElectronicAzimuthLimits',azLimits,'ElectronicElevationLimits',elLimits,...
    'FieldOfView',fov,'MountingAngles',[0 pitch 0]);

Обеспеченный класс помощника может использоваться, чтобы визуализировать наш шаблон скана. Мы будем также использовать этот плоттер, чтобы анимировать шаблон скана в цикле симуляции.

scanplt = helperScanPatternDisplay(radar);
scanplt.makeOverviewPlot;

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

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

numResolutionCells = diff(rangeLims)*prod(fov)/(rangeRes*prod(angRes));
numFrames = numFullScans*numScanPoints; % total number of frames for detection
expNumFAs = faRate*numFrames*numResolutionCells % expected number of FAs

expNumFAs = 7.9200

Сценарий и цель

Мы будем использовать radarScenario управлять потоком симуляции верхнего уровня. Этот объект обеспечивает функциональность, чтобы добавить новые платформы и быстро сгенерировать обнаружения от всех объектов датчика в сценарии. Это также управляет течением времени в симуляции. Используйте общее количество систем координат, вычисленных выше, чтобы найти необходимый StopTime, и адаптивно обновить время симуляции на основе частоты обновления объектов в сцене, установите частоту обновления сценария на 0.

% Create scenario
stopTime = numFullScans*numScanPoints/prf;
scene = radarScenario('StopTime',stopTime,'UpdateRate',0);

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

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

% Create platforms
rdrPlat = platform(scene,'Position',[0 0 12]); % place 12 m above the origin
tgtPlat(1) = platform(scene,'Trajectory',kinematicTrajectory('Position',[38e3 6e3 10e3],'Velocity',[-380 -50 0]));
tgtPlat(2) = platform(scene,'Trajectory',kinematicTrajectory('Position',[25e3 -3e3 1e3],'Velocity',[-280 10 -10]));

Чтобы включить методы обнаружения удобства и автоматическое продвижение времени симуляции, сценарий должен знать о радарном объекте, созданном ранее. Для этого смонтируйте радар на платформу путем заполнения свойства Sensors платформы. Чтобы смонтировать больше чем один датчик на платформе, свойство Sensors может быть массивом ячеек объектов датчика.

% Mount radar to platform
rdrPlat.Sensors = radar;

rcsSignature класс может использоваться, чтобы задать платформу RCS в зависимости от угла обзора. Для простой постоянной-RCS подписи, просто устанавливает Pattern свойство к скалярному значению. Позвольте одной платформе быть 20 dBsm и другими 4 dBsm, чтобы продемонстрировать различие в обнаружительной способности.

% Set target platform RCS profiles
tgtRCS = [20, 4]; % dBsm
tgtPlat(1).Signatures = rcsSignature('Pattern',tgtRCS(1));
tgtPlat(2).Signatures = rcsSignature('Pattern',tgtRCS(2));

Обнаружительная способность

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

time = (0:numFrames-1).'/updateRate;
tgtPosTruth(:,:,1) = tgtPlat(1).Position + time*tgtPlat(1).Trajectory.Velocity;
tgtPosTruth(:,:,2) = tgtPlat(2).Position + time*tgtPlat(2).Trajectory.Velocity;
tgtRangeTruth(:,1) = sqrt(sum((tgtPosTruth(:,:,1)-rdrPlat.Position).^2,2));
tgtRangeTruth(:,2) = sqrt(sum((tgtPosTruth(:,:,2)-rdrPlat.Position).^2,2));

Мы можем затем использовать известный целевой RCS и радарное усиление цикла, вычисленное радарной моделью, чтобы получить ОСШ истины для каждой цели:

tgtSnr = tgtRCS - 40*log10(tgtRangeTruth) + radar.RadarLoopGain;

Минимальный ОСШ, требуемый для обнаружения, учитывая схему обнаружения только для величины, такую как CFAR, учитывая нашу вероятность обнаружения и ложного сигнального уровня,

minSnr = 20*log10(erfcinv(2*radar.FalseAlarmRate) - erfcinv(2*radar.DetectionProbability)); % dB

Идя только немного далее, мы можем вычислить нашу минимальную область значений для обнаружения:

Rd = 10.^((tgtRCS - minSnr + radar.RadarLoopGain)/40); % meters

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

snrMargin = mean(tgtSnr,1) - minSnr

snrMargin = 1×2

    8.0827    0.0019

rangeMargin = (Rd - mean(tgtRangeTruth,1))/1e3 % km

rangeMargin = 1×2

   23.5299    0.0028

Первая цель на приблизительно 8 дБ более ярка, чем необходимый для обнаружения данными параметрами CFAR, в то время как вторая цель рассматривается едва обнаруживаемой.

Симулируйте обнаружения

Первичный цикл симуляции запускается путем вызова advance на объекте сценария. Этот метод продвигается, сценарий передает следующему разу, когда в который было определено, что объект в сцене требует обновления и возвращает false, когда заданное время остановки достигнуто, выйдя из цикла. Обратите внимание на то, что первый вызов advance не продвигается время симуляции, так, чтобы первая система координат сбора данных могла произойти в 0 с. Как альтернатива использованию advance, scene.SimulationStatus может проверяться, чтобы определить, выполнение ли sim до завершения.

В цикле обнаружение и методы генерации дорожки могут быть названы для любого объекта датчика индивидуально (см. датчик step метод), но удобные методы существуют, чтобы сгенерировать обнаружения и дорожки для всех датчиков и от всех целей в сцене с одним вызовом функции. Поскольку у нас только есть один датчик, detect(scene) хорошее решение должно получить список обнаружений для всех целей в сцене. Путем инициирования генерации обнаружений в верхнем уровне как это сама сцена обработает требуемую синхронизацию передачи, INS и другие данные конфигурации к датчику. Класс визуализации использовал, ранее будет использоваться здесь также, чтобы анимировать шаблон скана.

allDets = []; % initialize list of all detections generated
lookAngs = zeros(2,numFrames);
simTime = zeros(1,numFrames);
frame = 0; % current frame
while advance(scene) % advance until StopTime is reached
    
    % Increment frame index and record simulation time
    frame = frame + 1;
    simTime(frame) = scene.SimulationTime;
    
    % Generate new detections
    dets = detect(scene);
    
    % Record look angle for inspection
    lookAngs(:,frame) = radar.LookAngle;
    
    % Update coverage plot
    scanplt.updatePlot(radar.LookAngle);
    
    % Compile any new detections
    allDets = [allDets;dets];
    
end

Отсканируйте поведение

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

figure;
plot(simTime(1:2*numScanPoints)*1e3,lookAngs(:,1:2*numScanPoints));
xlabel('Time (ms)');
ylabel('Look Angle (deg)');
legend('Azimuth','Elevation');
title('Look Angles');

Обнаружения

Давайте смотреть содержимое objectDetection как выведено нашим радаром.

allDets{1} % look at the first detection

ans = 
  objectDetection with properties:

                     Time: 7.5000e-04
              Measurement: [3x1 double]
         MeasurementNoise: [3x3 double]
              SensorIndex: 1
            ObjectClassID: 0
    MeasurementParameters: [1x1 struct]
         ObjectAttributes: {[1x1 struct]}

Time время симуляции, в котором было сгенерировано обнаружение. SensorIndex говорит нам, какой датчик это обнаружение было сгенерировано (который важен, когда мы имеем несколько датчиков и используем уровень сцены, обнаруживают метод). Measurement измеримое значение, сопоставленное с обнаружением. Формат зависит от выбора нашего режима обнаружения и выходных координат. MeasurementNoise дает отклонение или ковариацию измерения и используется моделями средства отслеживания. Поскольку мы выводим обнаружения в координатах сценария, поле измерения является просто предполагаемым радиус-вектором цели, и шум измерения дает ковариацию той оценки положения.

ObjectAttributes поле содержит некоторые полезные метаданные:

allDets{1}.ObjectAttributes{1}

ans = struct with fields:
    TargetIndex: 3
            SNR: 12.5894

Это говорит нам, что обнаружение прибыло из платформы с индексом 2, который является уникальным идентификатором, присвоенным сценарием нашей первой целевой платформе (вторая платформа, которая будет добавлена к сцене). Мы также видим, что ОСШ этого обнаружения о том, что мы ожидали для нашей первой цели. Если обнаружение будет ложным предупреждением, TargetIndex будет-1, который является недопустимым ID платформы. Скомпилируйте целевой индекс, ОСШ, и время каждого обнаружения.

detTgtIdx = cellfun(@(t) t.ObjectAttributes{1}.TargetIndex, allDets);
detTgtSnr = cellfun(@(t) t.ObjectAttributes{1}.SNR, allDets);
detTime   = cellfun(@(t) t.Time, allDets);
firstTgtDet  = find(detTgtIdx == tgtPlat(1).PlatformID);
secondTgtDet = find(detTgtIdx == tgtPlat(2).PlatformID);

Мы можем найти общее количество обнаружений FA с:

numFAs = numel(find(detTgtIdx==-1))

numFAs = 8

который является близко к ожидаемому количеству ФАСА, вычисленного ранее.

Соберите данные об измерении через обнаружения и график наряду с положениями истины для этих двух целевых платформ.

tgtPosObs = cell2mat(cellfun(@(t) t.Measurement,allDets,'UniformOutput',0).');
subplot(1,2,1);
helperPlotPositionError( tgtPosTruth(:,:,1),time,tgtPosObs(:,firstTgtDet),detTime(firstTgtDet),scene.SimulationTime );
title('Target 1 Position Error');
subplot(1,2,2);
helperPlotPositionError( tgtPosTruth(:,:,2),time,tgtPosObs(:,secondTgtDet),detTime(secondTgtDet),scene.SimulationTime );
title('Target 2 Position Error');
set(gcf,'Position',get(gcf,'Position')+[0 0 600 0]);

Большинство наших обнаружений происходит из первой цели, которая не удивляет рассмотрение анализа обнаружительной способности ранее. Однако это не было обнаружено на каждой системе координат, и вторая цель сгенерировала много обнаружений несмотря на наличие низкого ОСШ. Хотя первая цель легко обнаруживается, погрешность измерения большая, потому что цель находится в нашей второй неоднозначности области значений, и никакое разрешение неоднозначности не было выполнено, тогда как наша вторая цель, когда обнаружено, показывает ошибку положения, которая сопоставима с нашим 100-метровым разрешением области значений и относительно плохим угловым разрешением.

Давайте также посмотрим на общее отклонение обнаружений против ОСШ. Общее отклонение является суммой крайних отклонений в каждом декартовом направлении (X, Y, и Z). Это отклонение включает эффекты оценки пробел угла области значений наряду с преобразованием тех статистических данных к координатам сценария.

detTotVar = cellfun(@(t) trace(t.MeasurementNoise),allDets);
figure;
subplot(1,2,1);
helperPlotPositionTotVar( detTgtSnr(firstTgtDet),detTotVar(firstTgtDet) );
title('First Target');
subplot(1,2,2);
helperPlotPositionTotVar( detTgtSnr(secondTgtDet),detTotVar(secondTgtDet) );
title('Second Target');
set(gcf,'Position',get(gcf,'Position')+[0 0 600 0]);

Обратите внимание на то, что несмотря на то, что абсолютная ошибка положения для первой цели была больше из-за неоднозначности области значений, более низкое отклонение измерения отражает и больший ОСШ первой цели и ее по-видимому малую дальность.

Симулируйте дорожки

radarDataGenerator также способно к выполнению отслеживания разомкнутого контура. Вместо того, чтобы вывести объекты обнаружения, это может вывести обновления дорожки.

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

release(radar); % unlock for re-simulation
radar.TargetReportFormat = 'Tracks';
radar.TrackCoordinates = 'Scenario';

Отслеживание может быть выполнено многими различными алгоритмами. Существуют простые фильтры альфы - беты и Фильтры Калмана линейного, расширенного, и недушистого разнообразия, наряду с постоянной скоростью, постоянным ускорением и постоянными-turnrate моделями движения. Алгоритм, который будет использоваться, задан с FilterInitializationFcn свойство. Это - указатель на функцию или имя функции в векторе символов, т.е. @initcvekf или 'initcvekf', это берет начальное обнаружение, как введено и возвращает инициализированный объект средства отслеживания. Любая пользовательская функция с той же подписью может использоваться. Мы будем использовать расширенный Фильтр Калмана предоставленной постоянной скорости.

radar.FilterInitializationFcn = @initcvekf;

Полная структура sim цикла является тем же самым, однако мы должны вызвать step радара функционируйте вручную. Необходимые входные параметры являются целевым struct положений, полученным путем вызова targetPoses на нашей радарной платформе, struct INS, полученном путем вызова pose на радарной платформе и текущем времени симуляции.

restart(scene); % restart scenario

allTracks = []; % initialize list of all track data
while advance(scene)
    
    % Generate new track data
    tgtPoses = targetPoses(rdrPlat);
    ins = pose(rdrPlat);
    tracks = radar(tgtPoses,ins,scene.SimulationTime);
    
    % Compile any new track data
    allTracks = [allTracks;tracks];
    
end

Осмотр отслеживаемого объекта:

allTracks(1) % look at first track update

ans = 
  objectTrack with properties:

             TrackID: 1
            BranchID: 0
         SourceIndex: 1
          UpdateTime: 0.0093
                 Age: 34
               State: [6x1 double]
     StateCovariance: [6x6 double]
     StateParameters: [1x1 struct]
       ObjectClassID: 0
          TrackLogic: 'History'
     TrackLogicState: [1 1 0 0 0]
         IsConfirmed: 1
           IsCoasted: 0
      IsSelfReported: 1
    ObjectAttributes: [1x1 struct]

TrackID уникальный идентификатор файла дорожки и SourceIndex ID объекта средства отслеживания создание отчетов о дорожке, которая полезна, когда у вас есть несколько средств отслеживания в сцене. Когда средство отслеживания мультигипотезы используется, BranchID дает индекс используемой гипотезы. State поле содержит информацию состояния об этом обновлении. Поскольку мы используем модель постоянной скорости и координаты дорожки системы координат сценария, State состоит из радиус-векторов и векторов скорости. UpdateTime дает время, в которое было сгенерировано это обновление дорожки. Другое важное свойство, IsCoasted, говорит вам, если это обновление дорожки использовало новое целевое обнаружение, чтобы обновить фильтр, или если это просто распространило вперед вовремя от последнего целевого обнаружения (дорожка "курсируется").

Структура ObjectAttributes дает ID платформы, подпись которой обнаруживается, и ОСШ той подписи во время обновления.

allTracks(1).ObjectAttributes

ans = struct with fields:
    TargetIndex: 3
            SNR: 18.5910

Соберите целевой индекс и время обновления. Также проверяйте, сколько обнаружений использовалось, чтобы обновить наши дорожки.

trackTgtIdx = arrayfun(@(t) t.TargetIndex,[allTracks.ObjectAttributes]);
updateTime = [allTracks.UpdateTime];
firstTgtTrack = find(trackTgtIdx == tgtPlat(1).PlatformID);
secondTgtTrack = find(trackTgtIdx == tgtPlat(2).PlatformID);

numUpdates = numel(find(~[allTracks.IsCoasted])) % total number of track updates with new detection data

numUpdates = 31

Это - примерно количество обнаружений неFA, сгенерированных в первой части.

Давайте посмотрим на положения, о которых сообщают, в плоскости XY и их ковариациях, наряду с данными о положении истины. Поскольку первая цель находится в первой неоднозначности области значений, и никакое разрешение неоднозначности не выполняется, давайте сначала найдем "истину" неоднозначным положением для первой цели как таковой:

los = tgtPosTruth(:,:,1) - rdrPlat.Position;
R = sqrt(sum(los.^2,2));
los = los./R;
Ramb = mod(R,rangeAmb);
tgtPosTruthAmb = rdrPlat.Position + los.*Ramb;

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

figure;
subplot(1,2,1);
helperPlotPositionTrack( tgtPosTruthAmb,allTracks(firstTgtTrack) );
title('First Target Track');
subplot(1,2,2);
helperPlotPositionTrack( tgtPosTruth(:,:,2),allTracks(secondTgtTrack) );
title('Second Target Track');
set(gcf,'Position',get(gcf,'Position')+[0 0 600 0]);

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

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

figure;
subplot(1,2,1);
helperPlotPositionTrackError( tgtPosTruthAmb,time,allTracks(firstTgtTrack),updateTime(firstTgtTrack) );
title('First Target Ambiguous Position Error');
subplot(1,2,2);
helperPlotPositionTrackError( tgtPosTruth(:,:,2),time,allTracks(secondTgtTrack),updateTime(secondTgtTrack) );
title('Second Target Position Error');
set(gcf,'Position',get(gcf,'Position')+[0 0 600 0]);

Для первой цели ошибка положения уменьшается на выборках, которые имели новое обнаружение, чтобы обновить файл дорожки (обозначенный красным кругом). Поскольку вторая цель имела намного более низкий ОСШ, ошибка в некоторых отдельных обнаружениях была достаточно большой увеличить ошибку положения дорожки. Несмотря на это, вторая целевая ошибка положения дорожки имеет небольшую тенденцию к понижению, указывая, что предполагаемое состояние будет, вероятно, сходиться. Однако мы видим проблему с нашей способностью отследить первую цель. Поскольку это находится в первой неоднозначности области значений и имеет компонент перпендикуляра движения к углу обзора, очевидная скорость на самом деле изменяется, и Фильтр Калмана постоянной скорости не может сходиться.

Заключение

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

Функции помощника

function helperPlotPositionError( tgtPosTruth,time,tgtPosObs,detTime,T )

tgtPosTruthX = interp1(time,tgtPosTruth(:,1),detTime).';
tgtPosTruthY = interp1(time,tgtPosTruth(:,2),detTime).';
tgtPosTruthZ = interp1(time,tgtPosTruth(:,3),detTime).';

err = sqrt((tgtPosTruthX - tgtPosObs(1,:)).^2+(tgtPosTruthY - tgtPosObs(2,:)).^2+(tgtPosTruthZ - tgtPosObs(3,:)).^2);

plot(detTime*1e3,err/1e3,'--o');
grid on;
xlim([0 T*1e3]);
ylabel('RMS Error (km)');
xlabel('Time (ms)');

end

function helperPlotPositionTotVar( detTgtSnr,detTotVar )

plot(detTgtSnr,detTotVar,'--o');
grid on;
xlabel('SNR (dB)');
ylabel('Total Var (m^2)');

end

function helperPlotPositionTrack( tgtPosTruth,tracks )

plot(tgtPosTruth(:,1)/1e3,tgtPosTruth(:,2)/1e3);
hold on;

updateIdx = find(~[tracks.IsCoasted]);
theta = 0:0.01:2*pi;

for ind = 1:2:numel(updateIdx)
    
    % plot update
    T = tracks(updateIdx(ind));
    plot(T.State(1)/1e3,T.State(3)/1e3,'*black');
    sigma = T.StateCovariance([1 3],[1 3]);
    [evec,eval] = eigs(sigma);
    mag = sqrt(diag(eval)).';
    u = evec(:,1)*mag(1);
    v = evec(:,2)*mag(2);
    x = cos(theta)*u(1) + sin(theta)*v(1);
    y = cos(theta)*u(2) + sin(theta)*v(2);
    plot(x/1e3 + T.State(1)/1e3,y/1e3 + T.State(3)/1e3,'magenta');
    
end

hold off;
grid on;
xlabel('X (km)');
ylabel('Y (km)');
legend('Truth Position','Position Estimates','Covariance','Location','northwest');

end

function helperPlotPositionTrackError( tgtPosTruth,time,tracks,updateTime )

state = [tracks.State];

sx = state(1,:);
sy = state(3,:);
sz = state(5,:);

x = interp1(time,tgtPosTruth(:,1),updateTime);
y = interp1(time,tgtPosTruth(:,2),updateTime);
z = interp1(time,tgtPosTruth(:,3),updateTime);

err = sqrt((x-sx).^2+(y-sy).^2+(z-sz).^2);

isUpdate = ~[tracks.IsCoasted];

plot(updateTime*1e3,err);
hold on;
plot(updateTime(isUpdate)*1e3,err(isUpdate),'or');
hold off;
grid on;
xlabel('Update Time (ms)');
ylabel('Error (m)')
legend('Coasted','Update','Location','northwest');

end