Многоплатформенный радарный Fusion обнаружения

В этом примере показано, как плавить радарные обнаружения от многоплатформенной радарной сети. Сеть включает два бортовых и наземные радарные платформы дальние. Смотрите Multiplatform Radar Detection Generation пример для деталей. Центральное средство отслеживания обрабатывает обнаружения со всех платформ в фиксированном интервале обновления. Это позволяет вам оценить производительность сети против различных целевых типов и маневров, а также конфигураций платформы и местоположений.

Загрузите симулированные данные о датчике и основную истину

MultiplatformRadarDetectionGeneration.mat файл содержит три переменные.

  • truthLog struct, который содержит основную истину кинематическая информация всех семи целей и платформ на каждом шаге симуляции (сохраненный в Truth поле . Configurations поле содержит ориентацию, положение луча и систему координат координаты датчика для всех четырех радарных датчиков на каждом шаге симуляции.

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

  • scene содержит trackingScenario который используется в целях визуализации.

path = fullfile(matlabroot,'examples','fusion','main');
addpath(path)
load('MultiplatformRadarDetectionGeneration.mat');

Задайте центральное средство отслеживания

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

Средство отслеживания использует initFilter поддерживание функции, чтобы инициализировать постоянную скорость расширило Фильтр Калмана для каждого нового трека. initFilter изменяет фильтр, возвращенный initcvekf совпадать с высокими целевыми скоростями. Шум процесса фильтра собирается в 1-G позволить отследить маневрирования целей в сценарии.

AssignmentThreshold средства отслеживания увеличен с его значения по умолчанию 30 - 50, чтобы позволить обнаружениям со смещениями большого спектра из-за атмосферных эффектов преломления в длинных областях значений обнаружения быть сопоставленными средством отслеживания. DeletionThreshold уменьшен с его значения по умолчанию от 5 до 3, чтобы удалить избыточные дорожки более быстро.

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

tracker = trackerGNN('FilterInitializationFcn', @initFilter, ...
        'AssignmentThreshold', 50, 'DeletionThreshold', 3, ...
        'HasDetectableTrackIDsInput', true);

Отследите цели путем плавления обнаружений в центральном средстве отслеживания

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

trackUpdateRate = 1/2;   % Update the tracker every 2 seconds
stop = scene.StopTime;
trackUpdateTimes = (1/trackUpdateRate:1/trackUpdateRate:stop)';

% Values used for tracker update logic.
allTracks = [];

% Create a display to show the true, measured, and tracked positions of the
% detected targets and platforms.
theaterDisplay = helperMultiPlatFusionDisplay(scene, 'XLim', [-45 45], ...
    'YLim', [-45 45], 'PlotAssignmentMetrics', true);

detectableTrackIDs = uint32([]);

% construct an object to analyze assignment and error metrics
tam = trackAssignmentMetrics('DistanceFunctionFormat','custom',...
    'AssignmentDistanceFcn',@truthAssignmentDistance,...
    'DivergenceDistanceFcn',@truthAssignmentDistance);
assignmentTable = [];

for i=1:numel(trackUpdateTimes)
    time = trackUpdateTimes(i);

    % Buffer all detections since last track update
    isDetReady = (dataLog.DetectionTime>time-1/trackUpdateRate) & (dataLog.DetectionTime<=time);
    detBuffer = dataLog.Detections(isDetReady);

    % Get all beam positions since last track update
    configSelect = (truthLog.SimulationTime>time-1/trackUpdateRate) & (truthLog.SimulationTime<=time);
    configsBuffer = truthLog.Configurations(configSelect,:);

    lastDetectionTime = max(dataLog.DetectionTime(isDetReady));
    if isempty(lastDetectionTime)
        lastDetectionTime = time; % There was no detections since last update
    end

    % Collect list of tracks which fell within at least one radar's field
    % of view since last tracker update
    if isLocked(tracker)
        predictedtracks = predictTracksToTime(tracker, 'all', lastDetectionTime);
        detectableTrackIDs = detectableTracks(allTracks,predictedtracks,configsBuffer);
    end

    % Update tracker, only run track logic on tracks that fell within at
    % least one radar's field of view since last tracker update
    [confirmedTracks, ~, allTracks] = tracker(detBuffer,lastDetectionTime,detectableTrackIDs);

    % Analyze and retrieve the current track-to-truth assignment metrics
    truths = truthLog.Truth(find(configSelect,1,'last'),:);

    [~,~] = tam(confirmedTracks, truths);

    % Store assignment metrics in a table
    currentAssignmentTable = trackMetricsTable(tam);
    rowTimes = seconds(time*ones(size(currentAssignmentTable,1),1));
    assignmentTable = cat(1,assignmentTable,table2timetable(currentAssignmentTable,'RowTimes',rowTimes));

    % Update display with detections and tracks
    theaterDisplay(detBuffer,configsBuffer(end,:),confirmedTracks,assignmentTable,truths);

end

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

endTime = assignmentTable.Time(end);
assignmentTable(endTime,{'TrackID','AssignedTruthID','TotalLength','DivergenceCount','RedundancyCount','RedundancyLength'})
ans =

  9x6 timetable

     Time     TrackID    AssignedTruthID    TotalLength    DivergenceCount    RedundancyCount    RedundancyLength
    ______    _______    _______________    ___________    _______________    _______________    ________________

    60 sec       1               2              27                0                  0                  0        
    60 sec       7               4              27                0                  0                  0        
    60 sec       8               5              26                0                  0                  0        
    60 sec       9               1              22                0                  0                  0        
    60 sec      11             NaN              10                0                  0                  0        
    60 sec      12               6              20                0                  0                  0        
    60 sec      23               7              20                0                  1                  5        
    60 sec      32             NaN               7                0                  1                  7        
    60 sec      41               3              10                0                  0                  0        

Заметьте, что платформы, которые испытывают большинство затруднений, обеспечивающих дорожки (платформы 4 и 7) являются также платформами, самое большее лежит в диапазоне от радаров. Эта плохая производительность отслеживания приписана эффектам моделирования шума измерения как Распределение Гаусса. Это работает хорошо на цели с близких расстояний, но с больших расстояний, неопределенность измерения отклоняется от Распределения Гаусса. Следующая фигура сравнивает эллипс с 1 сигмой, соответствующий Гауссовой неопределенности измерения в цели на уровне 30 км, измеренных радаром с 5 степенями угловое разрешение фактической неопределенности измерения. Фактическая неопределенность измерения имеет вогнутую форму, которая следует из системы координат сферической координаты, в которой радар оценивает положение цели.

maxCondNum = 300;
figure;
helperPlotLongRangeCorrection(maxCondNum)

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

Симулируйте с коррекцией ковариации дальней

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

correctedDdets = longRangeCorrection(dataLog.Detections,maxCondNum);
[confirmedTracks,correctedAssignmentTable,ctheaterDisplay] = ...
    runMultiPlatFusionSim(scene,tracker,correctedDdets,dataLog,truthLog);

endTime = correctedAssignmentTable.Time(end);
correctedAssignmentTable(endTime,{'TrackID','AssignedTruthID','TotalLength','DivergenceCount','RedundancyCount','RedundancyLength'})
ans =

  7x6 timetable

     Time     TrackID    AssignedTruthID    TotalLength    DivergenceCount    RedundancyCount    RedundancyLength
    ______    _______    _______________    ___________    _______________    _______________    ________________

    60 sec       1              2               27                0                  0                  0        
    60 sec       7              4               27                0                  0                  0        
    60 sec       8              5               26                0                  0                  0        
    60 sec       9              1               22                0                  0                  0        
    60 sec      11              7               25                0                  0                  0        
    60 sec      12              6               20                0                  0                  0        
    60 sec      39              3               10                0                  0                  0        

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

ctheaterDisplay(dataLog.Detections,[],confirmedTracks,correctedAssignmentTable,truths)
axes(ctheaterDisplay.TheaterPlot.Parent)
legend('off')
xlim([-1 5]); ylim([-41 -36]); zlim([-5 0]);
view([-90 90])
axis square
title('Jet Executing Horizontal Turn')

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

view([-60 25])

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

xlim([-25 -9]); ylim([-31 -19]); zlim([-9 -2]);
view([-45 10])
title('Crossing Airliners')

Переключая точку зрения, чтобы фокусироваться на пересекающихся авиалайнерах, те же неточные высотные измерения изображены. Отметьте, как красные обнаружения сосредоточены на высоте 8 км, в то время как эти два авиалайнера летят на высотах 3 и 4 км. Использование очень большой неопределенности в высоте позволяет средству отслеживания игнорировать ошибочную высоту, читающую из красных обнаружений и отслеживать высоту с помощью других двух радаров. Смотря на ковариацию неопределенности дорожек T07 и T08, вы видите, что они обеспечивают сопоставимую оценку платформ P04 и P05, соответственно.

xlim([-10 10]); ylim([-0 10]); zlim([-12 -5]);
view([-15 10])
title('Airborne Radar Platforms')
rmpath(path)

Последний график фокусируется на двух бортовых радарных платформах. Каждая платформа обнаруживается другой платформой, а также наземным радаром. Траектории платформы пересекают друг друга, разделенного на 1 000 м в высоте, и их дорожки сопоставимы с основной истиной.

Сводные данные

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

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

initFilter Эта функция изменяет функциональный initcvekf обрабатывать более высокие скоростные цели, такие как авиалайнеры в сценарии.

function filter = initFilter(detection)
filter = initcvekf(detection);
classToUse = class(filter.StateCovariance);

% Airliners can move at speeds around 900 km/h. The velocity is initialized
% to 0, but will need to be able to quickly adapt to aircraft moving at
% these speeds. Use 900 km/h as 1 standard deviation for the velocity
% noise.
spd = 900*1e3/3600; % m/s
velCov = cast(spd^2,classToUse);
cov = filter.StateCovariance;
cov(2,2) = velCov;
cov(4,4) = velCov;
filter.StateCovariance = cov;

% Set filter's process noise to allow for some horizontal maneuver
scaleAccel = cast(10,classToUse);
Q = blkdiag(scaleAccel^2, scaleAccel^2, 1);
filter.ProcessNoise = Q;
end

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

function trackIDs = detectableTracks(tracks,predictedtracks,configs)
% Identify which tracks fell within a sensor's field of view

numTrks = size(tracks,1);
[numsteps, numSensors] = size(configs);
allposTrack = zeros(3,numsteps);
isDetectable = false(numTrks,1);
for iTrk = 1:numTrks
    % Interpolate track positions between current position and predicted
    % positions for each simulation step
    posTrack = tracks(iTrk).State(1:2:end);
    posPreditedTrack = predictedtracks(iTrk).State(1:2:end);
    for iPos = 1:3
        allposTrack(iPos,:) = linspace(posTrack(iPos),posPreditedTrack(iPos),numsteps);
    end
    for iSensor = 1:numSensors
        thisConfig = configs(:,iSensor);
        for k = 1:numsteps
            if thisConfig(k).IsValidTime
                pos = trackToSensor(allposTrack(:,k),thisConfig(k));
                % Check if predicted track position is in sensor field of
                % view
                [az,el] = cart2sph(pos(1),pos(2),pos(3));
                az = az*180/pi;
                el = el*180/pi;
                inFov = abs(az)<thisConfig(k).FieldOfView(1)/2 && abs(el) < thisConfig(k).FieldOfView(2)/2;
                if inFov
                    isDetectable(iTrk) = inFov;
                    k = numsteps; %#ok<FXSET>
                    iSensor = numSensors; %#ok<FXSET>
                end
            end
        end
    end

end

trackIDs = [tracks.TrackID]';
trackIDs = trackIDs(isDetectable);
end

trackToSensor Эта функция возвращает положение дорожки в координатной системе координат датчика. Структура дорожки возвращена trackerGNN и структура config, задающая ориентацию датчика относительно координатной системы координат дорожки, возвращена monostaticRadarSensor объект.

function pos = trackToSensor(pos,config)

frames = config.MeasurementParameters;
for m = numel(frames):-1:1
    rotmat = frames(m).Orientation;
    if ~frames(m).IsParentToChild
        rotmat = rotmat';
    end
    offset = frames(m).OriginPosition;
    pos = bsxfun(@minus,pos,offset);
    pos = rotmat*pos;
end
end

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

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

function dets = longRangeCorrection(dets,maxCond)
for m = 1:numel(dets)
    R = dets{m}.MeasurementNoise;
    [Q,D] = eig(R);
    Q = real(Q);
    d = real(diag(D));
    dMax = max(d);
    condNums = dMax./d;
    iFix = condNums>maxCond;
    d(iFix) = dMax/maxCond;
    R = Q*diag(d)*Q';
    dets{m}.MeasurementNoise = R;
end
end