Main Content

Track Heterogeneous Aircraft using Heterogeneous Sensors

Since R2024b

This example shows how to track heterogeneous aircraft using data collected by heterogeneous sensors. In the example, you track commercial passenger aircraft, general aviation aircraft, and helicopters flying in the airspace near the Logan International Airport in Boston, MA. You use active and passive radars and their outputs to track the aircraft using a JIPDATracker.

Define Scenario

This example shows tracking of heterogeneous aircraft using heterogeneous sensors. The helper function helperCreateScenario defines a scenario around Logan International Airport in Boston, MA metropolitan area. The scenario contains two ground-based sensors mounted on towers and seven aircraft.

scenario = helperCreateScenario;

The first sensor, which is mounted on the scenario platform number 1, is an airport primary surveillance radar, which is a monostatic, active, rotating radar. As it scans, the radar reports measurements that contain azimuth, elevation, range, and range rate. The second sensor is a passive secondary receiver, which scans the sky for signals emitted by the traveling aircraft. It measures the direction from which the signals arrive, in terms of azimuth and elevation angles.

The seven aircraft are listed in the table below.

ID

Type

Position

Speed

1

Primary surveillance radar tower

[0 0 -50]

Stationary

2

Secondary radar tower

[-5000 -3000 -50]

Stationary

3

Passenger aircraft

Traveling east towards Logan at 3000m altitude

900 km/hr

4

Passenger aircraft

Traveling north away from Logan at 4000m altitude

700 km/hr

5

Passenger aircraft

Traveling east tangential to Logan at 4000m altitude

600 km/hr

6

General aviation aircraft

Traveling south tangential to Logan at 2000m altitude

300 km/hr

7

General aviation aircraft

Traveling east tangential to Logan at 1000m altitude

200 km/hr

8

Helicopter

Traveling north climbing from 700m to 1000m.

50 km/hr

9

Helicopter

Traveling east descending from 600m to 320m

10 km/hr

The following code creates a trackingGlobeViewer and visualizes the scenario.

mapOrigin = [42.39423231362 -70.95934958874 0];
mapViewer = trackingGlobeViewer('ReferenceLocation',mapOrigin,'Basemap','streets-dark');
campos(mapViewer, mapOrigin + [0 0 1e5]);
drawnow;
plotScenario(mapViewer,scenario);
snapshot(mapViewer);

Figure contains an axes object. The hidden axes object contains an object of type image.

Define Tracker

Tracking heterogeneous aircraft using heterogeneous sensors is not trivial. You need to consider the varying data formats and sensor capabilities. Additionally, the system must be robust enough to handle real-time data integration and ensure accurate tracking across different aircraft types.

The tracker must be a multi-sensor data fusion algorithm capable of integrating inputs from all available sensors to provide a unified and accurate representation of each aircraft's state.

As new sensor data become available, the tracker takes into consideration existing tracks it maintains, associates the new sensor data to the existing tracks, updates the tracks using common filtering techniques, like Kalman filtering, and manages the track existence probability, taking into account the sensor instantaneous field of view. All these steps are enabled by defining the target and sensor specifications, which describe how the targets move and how the sensors observe the scenario. Therefore, the steps to define the tracker are shown in the following process:

Specify Aircraft Types

The first step in defining the tracker requires specifying the types of objects you want to track. In this example, these objects are passenger aircraft, general aviation aircraft, and helicopters. Use the trackerTargetSpec function to create a target specification for a passenger aircraft and observe the properties of the specification.

passengerSpec = trackerTargetSpec("aerospace","aircraft","passenger");
disp(passengerSpec)
  PassengerAircraft with properties:

           MaxHorizontalSpeed: 250    m/s 
             MaxVerticalSpeed: 20     m/s 
    MaxHorizontalAcceleration: 10     m/s²
      MaxVerticalAcceleration: 1      m/s²

The passenger aircraft target specification defines typical speeds and accelerations in the horizontal and vertical planes for a passenger aircraft. The speeds are used to define a prior distribution for a constant velocity model while the accelerations are used to define the process noise of the constant velocity model. This is all taken care of by the spec, when it integrates with the tracker shown below. In the passenger aircraft case, the speeds and accelerations represent a typical passenger jet flying at common cruise speeds. You can use this specification for tracking any passenger jet-propelled aircraft, including business jets. The speeds and accelerations are tunable and can be modified to fit specific aircraft type.

To add the general aviation and helicopter specifications, use the same process. Observe the differences between the relatively fast jet passenger aircraft, the slower general aviation fixed-wing aircraft, and the even slower helicopter aircraft.

generalAviationSpec = trackerTargetSpec("aerospace","aircraft","general-aviation")
generalAviationSpec = 
  GeneralAviation with properties:

           MaxHorizontalSpeed: 150    m/s 
             MaxVerticalSpeed: 15     m/s 
    MaxHorizontalAcceleration: 1      m/s²
      MaxVerticalAcceleration: 5      m/s²
helicopterSpec = trackerTargetSpec("aerospace","aircraft","helicopter")
helicopterSpec = 
  Helicopter with properties:

           MaxHorizontalSpeed: 80    m/s 
             MaxVerticalSpeed: 10    m/s 
    MaxHorizontalAcceleration: 5     m/s²
      MaxVerticalAcceleration: 1     m/s²

Specify Sensor Types

After specifying what you want to track, you must specify the sensors you use for tracking. In this example, you use the airport primary surveillance radar and a secondary radar. The airport surveillance radar is a monostatic rotating radar that scans the sky as it rotates. The secondary radar is a passive radar that looks for RF returns from the tracked objects. First, use the trackerSensorSpec to create the monostatic radar specification for the airport surveillance radar and observe its properties.

activeRadarSpec = trackerSensorSpec("aerospace","radar","monostatic")
activeRadarSpec = 
  AerospaceMonostaticRadar with properties:

           MaxNumLooksPerUpdate: 30    
    MaxNumMeasurementsPerUpdate: 10    
    IsPlatformStationary: 1                  
        PlatformPosition: [0 0 0]         m  
     PlatformOrientation: [3⨯3 double]       
        MountingLocation: [0 0 0]         m  
          MountingAngles: [0 0 0]         deg
       HasElevation: 1                
       HasRangeRate: 1                
        FieldOfView: [60 20]       deg
        RangeLimits: [0 1e+05]     m  
    RangeRateLimits: [-500 500]    m/s
      AzimuthResolution: 1      deg
        RangeResolution: 100    m  
    ElevationResolution: 5      deg
    RangeRateResolution: 10     m/s
    DetectionProbability: 0.9      
          FalseAlarmRate: 1e-06    

You can divide the radar properties into the following broad categories:

  1. Properties related to how the radar operates and how it reports to the tracker, for example MaxNumLooksPerUpdate, IsPlatformStationary, PlatformPosition, etc. These properties must be specified using information from the radar operator, for example Logan airport.

  2. Properties related to the radar measurement information, for example HasElevation, FieldOfView, RangeLimits, etc. These properties are specified by the radar manufacturer for a specific radar mode in use and can be found on a radar spec sheet.

Use the following code to define the Logan airport surveillance radar. The radar reports its measurements to the tracker once every two seconds.

towerPosition = [0 0 -50];
activeRadarSpec.MaxNumLooksPerUpdate = ceil(75*2/1.4); % Defines the number of looks in two seconds
activeRadarSpec.MaxNumMeasurementsPerUpdate = 20;
activeRadarSpec.PlatformPosition = towerPosition;
activeRadarSpec.PlatformOrientation = eye(3);

You follow a similar process to define the scanning passive radar.

passiveRadarSpec = trackerSensorSpec("aerospace","radar","direction-finder");
passiveRadarSpec.MaxNumMeasurementsPerUpdate = 10;
passiveRadarSpec.MaxNumLooksPerUpdate = ceil(150/1.4);
passiveRadarSpec.PlatformPosition = [-5000 -3000 -50];
passiveRadarSpec.PlatformOrientation = eye(3);
passiveRadarSpec.RangeLimits = [0 1e5];

The following two lines of code ensure that the properties of the sensor specifications match the equivalent properties of the scenario sensors.

activeRadarSpec = helperSyncSensor2spec(scenario.Platforms{1}.Sensors{1},activeRadarSpec);
passiveRadarSpec = helperSyncSensor2spec(scenario.Platforms{2}.Sensors{1},passiveRadarSpec);

Configure Tracker to Use Target and Sensor Specifications

The last step in specifying the tracker is to configure it to use the heterogeneous target and sensor specifications you have defined.

tracker = multiSensorTargetTracker({passengerSpec,generalAviationSpec,helicopterSpec},{activeRadarSpec,passiveRadarSpec},"jipda");

Using the target specifications, the tracker configures the right motion models. Using the sensor specifications, the tracker configures the right observability, measurement, clutter, birth, and state initialization models. This is all done within the tracker and does not require additional tuning. Note that all the target and sensor specification properties are tunable.

disp(tracker)
  fusion.tracker.JIPDATracker with properties:

                TargetSpecifications: {[1×1 PassengerAircraft]  [1×1 GeneralAviation]  [1×1 Helicopter]}
                SensorSpecifications: {[1×1 AerospaceMonostaticRadar]  [1×1 AerospaceESMRadar]}
              MaxMahalanobisDistance: 5
    ConfirmationExistenceProbability: 0.9000
        DeletionExistenceProbability: 0.1000

In addition to the target and sensor specifications, the tracker exposes three tuning parameters to help control data association, track confirmation, and track deletion.

The MaxMahalanobisDistance property controls the size of the association gate. Generally, this distance should be small to limit the association gate, but if you observe that the tracker fails to associate a measurement to the track, you can increase this value slightly.

Any multi-object tracker needs to handle the possibility of false alarms reported by the sensors. As a result, the tracker creates tracks that are initially tentative, meaning they are not confirmed to be real objects. The ConfirmationExistenceProbability property controls how quickly the tracker confirms a track. A higher value of this property means that the track must have a higher probability of existence to be confirmed as a real object. A lower value means that tracks with lower probability can be confirmed and it may lead to a higher rate of false tracks.

Similarly, any multi-object tracker needs to handle the possibility that an object may not exist anymore. The DeletionExistenceProbability property defines how quickly such tracks are removed from the tracker. A higher value means more aggressive deletion and can help improve the tracker performance at the cost of a higher likelihood of erroneously deleting a real track.

Understand Sensor Data Format

The sensor specifications allow you to pass sensor data into the tracker and act as an interpreter between the sensor data and the tracker. To understand this point better, display the data format for the active radar spec.

dataFormat(activeRadarSpec)
ans = struct with fields:
             LookTime: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
          LookAzimuth: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
        LookElevation: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
        DetectionTime: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
              Azimuth: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
            Elevation: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
                Range: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
            RangeRate: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
      AzimuthAccuracy: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
    ElevationAccuracy: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
        RangeAccuracy: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
    RangeRateAccuracy: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]

The active radar data format includes all the information needed to understand the sensor scanning and measurements. The LookTime, LookAzimuth, and LookElevation fields provide information about the sensor beam azimuth and elevation angles as a function of time. Each of these fields contains as many zeros as the MaxNumLooksPerUpdate property of the active radar specification. When populating the struct with your sensor data, these fields are variable sized and you can pass only the number of elements that you need, for example, [0 0.1 0.2] for time.

The rest of the fields provide information about the sensor detections: their time, measurement, and measurement accuracy. The measurement and measurement accuracy are further split into azimuth, elevation, range, and range rate. Furthermore, the number of zeros in each of these fields is equal to the MaxNumMeasurementsPerUpdate property and you use them as variable sized.

A similar pattern is used for the passive radar data format.

The next line loads sensor data that was recorded in the active and passive radar data formats. The first cell in the sensorData cell array contains the monostatic sensor data from the beginning to the end of the scenario. Similarly, the second cell contains the passive radar data. To recreate the data, uncomment the line below, but note that it may take a considerable amount of time to run.

load("sensorData.mat","sensorData");
% sensorData = helperRecordSensorData(scenario);

Run Tracker

The following controls allow tuning the tracker.

release(tracker);
tracker.DeletionExistenceProbability = 0.3;
tracker.ConfirmationExistenceProbability = 0.9;
tracker.MaxMahalanobisDistance = 5;

You can define how much time passes between tracker updates using the control below. The helperSplitSensorData and helperSavePlatformPoses functions split the sensor and truth poses data based on this interval.

trackerUpdateInterval = 1; % seconds
[activeRadarData,passiveRadarData] = helperSplitSensorData(sensorData,scenario.StopTime,trackerUpdateInterval);
platformPoses = helperSavePlatformPoses(scenario,trackerUpdateInterval);

The following code prepares the map display.

clear(mapViewer);
plotPlatform(mapViewer,[scenario.Platforms{3:end}],TrajectoryMode="Full");

Update the tracker and display in a loop.

for ind = 1:numel(activeRadarData)
    % Update tracker
    tracks = tracker(activeRadarData(ind), passiveRadarData(ind));
    
    % Visualize data and tracks (helper)
    plotActiveRadarData(mapViewer, activeRadarSpec,activeRadarData(ind));
    plotESMRadarData(mapViewer, passiveRadarSpec,passiveRadarData(ind));

    % Visualize tracks and platforms
    plotTrack(mapViewer, tracks, LabelStyle="ATC",LineWidth=6);
    plotPlatform(mapViewer,platformPoses(:,ind));
end

Take a final snapshot.

snapshot(mapViewer);

Figure contains an axes object. The hidden axes object contains an object of type image.

Summary

In this example, you learned how to track heterogeneous aircraft using a diverse set of radars. You saw how to specify the target types and the sensor types using the trackerTargetSpec and trackerSensorSpec functions, respectively. These functions allow you to define target and sensor specifications that are common in the aerospace and automotive domains. You also saw how to configure the tracker to use these specifications with the multiSensorTargetTracker function.

This example was set up in a way that allows you to explore the tracker capabilities. You can use the controls in the "Run the Tracker" section to modify the tracker parameters and update rate. The tracker parameters, including those parameters of each specification, are tunable during tracker's operation.

To explore further, you can modify the scenario. Ensure that you use the correct target and sensor specifications to match the new scenario, and update the sensor specifications accordingly.

Utility Functions

helperCreateScenario - A helper to create the scenario

function scenario = helperCreateScenario()
scenarioDuration = 60; % s
scenario = trackingScenario(UpdateRate=0,StopTime=scenarioDuration);
radarTower = platform(scenario, Position=[0 0 -50]);
rpm = 12.5;
fov = [1.4;10];
scanrate = rpm*360/60;  % deg/s
updaterate = scanrate/fov(1); % Hz

activeRadar = fusionRadarSensor(1,"Rotator", ...
    UpdateRate=updaterate, ...           % Hz
    FieldOfView=fov, ...                 % [az;el] deg
    MaxAzimuthScanRate=scanrate, ...     % deg/sec
    AzimuthResolution=fov(1), ...        % deg
    ReferenceRange=111e3, ...            % m
    ReferenceRCS=0, ...                  % dBsm
    RangeResolution=135, ...             % m
    HasElevation=true, ...
    HasRangeRate=true, ...
    RangeRateLimits=[-300 300], ...
    MountingLocation=[0 0 -15], ...
    MountingAngles=[0 0 0], ...
    DetectionCoordinates="Sensor spherical"); % Report spherical measurements

% Mount radar at the top of the tower
radarTower.Sensors = activeRadar;

% Set mechanical elevation scan to begin at 2 degrees above the horizon
elFov = fov(2);
tilt = 2; % deg
activeRadar.MechanicalElevationLimits = [-fov(2) 0]-tilt; % deg
activeRadar.FieldOfView(2) = elFov+1e-3;

esmTower = platform(scenario, Position=[-5000 -3000 -50]);
esmfov = [1.4 15];
passiveRadar = fusionRadarSensor(2,"Rotator", ...
    DetectionMode="ESM", ...
    MechanicalElevationLimits=[-esmfov(2) 0]-tilt, ... % deg
    UpdateRate=updaterate, ...           % Hz
    FieldOfView=esmfov, ...                 % [az;el] deg
    MaxAzimuthScanRate=scanrate, ...     % deg/sec
    AzimuthResolution=esmfov(1), ...        % deg
    HasElevation=true, ...
    MountingLocation=[0 0 -15], ...      % m
    MountingAngles=[-90 0 0], ...    % deg. Make sure it rotates after the active radar
    Sensitivity=-200 ...                 % dB
    ); % Reports spherical measurements

% Set mechanical elevation scan to begin at 2 degrees above the horizon
elFov = esmfov(2);
passiveRadar.FieldOfView(2) = elFov+1e-3;

% Mount radar at the top of the tower
esmTower.Sensors = passiveRadar;

emitter = radarEmitter(1,"No scanning",UpdateRate=updaterate,FieldOfView=[360 180]);

% Inbound passenger aircraft
ht = 3e3;
spd = 900*1e3/3600; % m/s
wp = [-5e3 -40e3 -ht;-5e3 -40e3+spd*scenarioDuration -ht];
traj = waypointTrajectory(wp,[0 scenarioDuration]);
emitter1 = clone(emitter);
emitter1.EmitterIndex = 3;
platform(scenario,Trajectory=traj,Emitters=emitter1);

% Outbound passenger aircraft
ht = 4e3;
spd = 700*1e3/3600; % m/s
wp = [20e3 10e3 -ht;20e3+spd*scenarioDuration 10e3 -ht];
traj = waypointTrajectory(wp,[0 scenarioDuration]);
emitter2 = clone(emitter);
emitter2.EmitterIndex = 4;
platform(scenario,Trajectory=traj,Emitters=emitter2);

% Tangential passenger aircraft
ht = 4e3;
spd = 600*1e3/3600; % m/s
wp = [-20e3 -spd*scenarioDuration/2 -ht;-20e3 spd*scenarioDuration/2 -ht];
traj = waypointTrajectory(wp,[0 scenarioDuration]);
emitter3 = clone(emitter);
emitter3.EmitterIndex = 5;
platform(scenario,Trajectory=traj,Emitters=emitter3);

% General aviation aircraft
ht = 2e3;
spd = 300*1e3/3600; % m/s
wp = [spd*scenarioDuration/2 10e3 -ht;-spd*scenarioDuration/2 10e3 -ht];
traj = waypointTrajectory(wp,[0 scenarioDuration]);
emitter4 = clone(emitter);
emitter4.EmitterIndex = 6;
platform(scenario,Trajectory=traj,Emitters=emitter4);

% General aviation aircraft
ht = 1e3;
spd = 200*1e3/3600; % m/s
wp = [-12e3 -4e3-spd*scenarioDuration/2 -ht;-12e3 -4e3+spd*scenarioDuration/2 -ht];
traj = waypointTrajectory(wp,[0 scenarioDuration]);
emitter5 = clone(emitter);
emitter5.EmitterIndex = 7;
platform(scenario,Trajectory=traj,Emitters=emitter5);

% % Helicopter taking off
ht = 700;
spd = 50*1e3/3600; % m/s
wp = [-spd*scenarioDuration/2 -13e3 -ht;spd*scenarioDuration/2 -13e3 -ht-300];
traj = waypointTrajectory(wp,[0 scenarioDuration]);
emitter6 = clone(emitter);
emitter6.EmitterIndex = 8;
platform(scenario,Trajectory=traj,Emitters=emitter6);
% 
% Helicopter landing
ht = 600;
spd = 10*1e3/3600; % m/s
wp = [-1e3 -5e3+spd*scenarioDuration/2 -ht;-1e3 -5e3+-spd*scenarioDuration/2 -ht+280];
traj = waypointTrajectory(wp,[0 scenarioDuration]);
emitter7 = clone(emitter);
emitter7.EmitterIndex = 9;
platform(scenario,Trajectory=traj,Emitters=emitter7);
end

helperSyncSensor2spec - A helper to synchronize sensor properties to the sensor specification

function spec = helperSyncSensor2spec(sensor,spec)
if isa(sensor,'fusionRadarSensor') && strcmpi(sensor.DetectionMode,'monostatic')
    spec = syncRadar2spec(sensor,spec);
elseif isa(sensor,'fusionRadarSensor') && strcmpi(sensor.DetectionMode,'esm')
    spec = syncESM2spec(sensor,spec);
end
end

function spec = syncRadar2spec(sensor,spec)
spec.MountingLocation = sensor.MountingLocation;
spec.MountingAngles = sensor.MountingAngles;
spec.HasElevation = sensor.HasElevation;
spec.HasRangeRate = sensor.HasRangeRate;
spec.FieldOfView = sensor.FieldOfView;
spec.RangeLimits = sensor.RangeLimits;
spec.RangeRateLimits = sensor.RangeRateLimits;
spec.AzimuthResolution = sensor.AzimuthResolution;
spec.RangeResolution = sensor.RangeResolution;
spec.ElevationResolution = sensor.ElevationResolution;
spec.RangeRateResolution = sensor.RangeRateResolution;
spec.DetectionProbability = sensor.DetectionProbability;
spec.FalseAlarmRate = sensor.FalseAlarmRate;
end

function spec = syncESM2spec(sensor,spec)
spec.MountingLocation = sensor.MountingLocation;
spec.MountingAngles = sensor.MountingAngles;
spec.HasElevation = sensor.HasElevation;
spec.FieldOfView = sensor.FieldOfView;
spec.AzimuthResolution = sensor.AzimuthResolution;
spec.ElevationResolution = sensor.ElevationResolution;
spec.DetectionProbability = sensor.DetectionProbability;
spec.FalseAlarmRate = sensor.FalseAlarmRate;
end

helperSavePlatformPoses - Record platform poses at tracker update timestamps

function platformPoses = helperSavePlatformPoses(scenario,trackerInterval)
endind = ceil(scenario.StopTime/trackerInterval);
samplePose = struct(PlatformID=1,ClassID=0,Position=zeros(1,3),Velocity=zeros(1,3),Acceleration=zeros(1,3),Orientation=eye(3),AngularVelocity=zeros(1,3));
platformPoses = repmat(samplePose,numel(scenario.Platforms),endind);
sampleTimes = trackerInterval*(1:endind);

for p = 1:numel(scenario.Platforms)
    thisPlatform = scenario.Platforms{p};
    trajectory = thisPlatform.Trajectory;
    if isa(trajectory,"waypointTrajectory")
        [position,orientation,velocity,acceleration,angularVelocity] = lookupPose(trajectory,sampleTimes);
        orientation = rotmat(orientation,"frame");
        for i = 1:endind
            platformPoses(p,i) = struct(PlatformID=thisPlatform.PlatformID,ClassID=thisPlatform.ClassID, ...
                Position=position(i,:),Velocity=velocity(i,:),Acceleration=acceleration(i,:), ...
                Orientation=orientation(:,:,i),AngularVelocity=angularVelocity(i,:));
        end
    else % Assumes that non waypointTrajectory means stationary object
        platformPoses(p,:) = repmat(struct(PlatformID=thisPlatform.PlatformID,ClassID=thisPlatform.ClassID, ...
            Position=thisPlatform.Position,Velocity=zeros(1,3),Acceleration=zeros(1,3),Orientation=thisPlatform.Orientation, ...
            AngularVelocity=zeros(1,3)),1,endind);
    end
end
end

helperRecordSensorData - A helper to record sensor data from the scenario

function sensorData = helperRecordSensorData(scenario) %#ok<DEFNU>
if isa(scenario,"trackingScenarioRecording")
    recording = scenario.RecordedData; % We only need the struct
elseif isa(scenario,"struct")
    recording = scenario;
elseif isa(scenario,"trackingScenario") || isa(scenario,"radarScenario")
    disp("Recording the scenario to get sensor data. This may take several minutes.");
    recording = record(scenario,"Rotmat",IncludeSensors=true,RecordingFormat="struct",InitialSeed=seed);
else
    error("This example utility function can only be used with a trackingScenario a trackingScenarioRecording");
end

if ~isfield(recording,"SensorConfigurations") || ~isfield(recording,"Detections")
    error("Scenario recording must contain sensor configurations and detections to record sensor data");
end
% Assume all steps have the same sensors.
sensors = recording(1).SensorConfigurations;
detBuffer = vertcat(recording.Detections);
for i = 1:numel(recording)
    for j = 1:numel(recording(i).SensorConfigurations)
        recording(i).SensorConfigurations(j).Time = recording(i).SimulationTime;
    end
end
configBuffer = reshape([recording.SensorConfigurations],[],1);
sensorData = cell(1,numel(sensors));
for k = 1:numel(sensors)
    sensor = sensors(k);
    sensorIdx = sensor.SensorIndex;
    sensorMountingAngles =  rotmat(recording(1).CoverageConfig(k).Orientation,"frame");
    
    detSensorIndex = cellfun(@(x)x.SensorIndex,detBuffer);
    dets = detBuffer(detSensorIndex == sensorIdx);

    cfgSensorIndex = arrayfun(@(x)x.SensorIndex, configBuffer);
    cfgs = configBuffer(cfgSensorIndex == sensorIdx);
    if sensor.MeasurementParameters.HasRange
        sensorData{k} = recordMonostaticData(sensorMountingAngles, dets, cfgs);
    else
        sensorData{k} = recordESMData(sensorMountingAngles, dets, cfgs);
    end
end
end

function activeRadarData = recordMonostaticData(sensorMountingAngles, detections, cfgs)
% Reported config information
lookMountRot = cell(numel(cfgs),1);
for i = 1:numel(cfgs)
    if cfgs(i).MeasurementParameters(1).IsParentToChild
        lookMountRot{i} = cfgs(i).MeasurementParameters(1).Orientation;
    else
        lookMountRot{i} = cfgs(i).MeasurementParameters(1).Orientation';
    end
end
beamRot = cellfun(@(x)x*sensorMountingAngles',lookMountRot,UniformOutput=false);
ypr = cellfun(@(x)fusion.internal.frames.rotmat2ypr(x,"degrees"),beamRot,UniformOutput=false);
lookTime = arrayfun(@(x)x.Time,cfgs);
lookAzimuth = cellfun(@(x)x(1),ypr);
lookElevation = cellfun(@(x)-x(2),ypr);

% Reported measurement information
detectionTime = cellfun(@(x)x.Time,detections);
azimuth = cellfun(@(x)x.Measurement(1),detections);
elevation = cellfun(@(x)x.Measurement(2),detections);
range = cellfun(@(x)x.Measurement(3),detections);
rangeRate = cellfun(@(x)x.Measurement(4),detections);
azimuthAccuracy =  cellfun(@(x)x.MeasurementNoise(1,1),detections);
elevationAccuracy =  cellfun(@(x)x.MeasurementNoise(2,2),detections);
rangeAccuracy =  cellfun(@(x)x.MeasurementNoise(3,3),detections);
rangeRateAccuracy =  cellfun(@(x)x.MeasurementNoise(4,4),detections);

activeRadarData = struct(LookTime=lookTime(:)',LookAzimuth=lookAzimuth(:)', ...
    LookElevation=lookElevation(:)', DetectionTime=detectionTime(:)', ...
    Azimuth=azimuth(:)',Elevation=elevation(:)', Range=range(:)', ...
    RangeRate=rangeRate(:)',AzimuthAccuracy=sqrt(azimuthAccuracy(:)'), ...
    ElevationAccuracy=sqrt(elevationAccuracy(:)'),RangeAccuracy=sqrt(rangeAccuracy(:)'), ...
    RangeRateAccuracy=sqrt(rangeRateAccuracy(:)'));
end

function esmData = recordESMData(sensorMountingAngles, detections, cfgs)
% Reported config information
lookMountRot = cell(numel(cfgs),1);
for i = 1:numel(cfgs)
    if cfgs(i).MeasurementParameters(1).IsParentToChild
        lookMountRot{i} = cfgs(i).MeasurementParameters(1).Orientation;
    else
        lookMountRot{i} = cfgs(i).MeasurementParameters(1).Orientation';
    end
end
beamRot = cellfun(@(x)x*sensorMountingAngles',lookMountRot,UniformOutput=false);
ypr = cellfun(@(x)fusion.internal.frames.rotmat2ypr(x,"degrees"),beamRot,UniformOutput=false);
lookAzimuth = cellfun(@(x)x(1),ypr);
lookElevation = cellfun(@(x)-x(2),ypr);
lookTime = arrayfun(@(x)x.Time,cfgs);

detectionTime = cellfun(@(x)x.Time,detections);
azimuth = cellfun(@(x)x.Measurement(1),detections);
elevation = cellfun(@(x)x.Measurement(2),detections);
azimuthAccuracy =  cellfun(@(x)x.MeasurementNoise(1,1),detections);
elevationAccuracy =  cellfun(@(x)x.MeasurementNoise(2,2),detections);

esmData = struct(LookTime=lookTime(:)',LookAzimuth=lookAzimuth(:)', ...
    LookElevation=lookElevation(:)',DetectionTime=detectionTime(:)', ...
    Azimuth=azimuth(:)',Elevation=elevation(:)',AzimuthAccuracy=sqrt(azimuthAccuracy(:)'), ...
    ElevationAccuracy=sqrt(elevationAccuracy(:)'));
end

helperSplitSensorData - A helper to split recorded sensor data to the tracker update interval

function varargout=helperSplitSensorData(sensorData,stopTime,trackerInterval)
varargout = cell(1,nargout);
if iscell(sensorData)
    varargout = cellfun(@(s) splitOneSensorData(s,stopTime,trackerInterval), sensorData, UniformOutput=false);
else
    varargout{1} = splitOneSensorData(sensorData,stopTime,trackerInterval);
end
end

function data=splitOneSensorData(sensorData,stopTime,trackerInterval)
frameTime = 0;
ind = 1;
numFrames = ceil(stopTime/trackerInterval);
data = repmat(sensorData,numFrames,1);
while frameTime < stopTime
    withinInterval = sensorData.LookTime >= frameTime & sensorData.LookTime < frameTime+trackerInterval;
    data(ind).LookTime = sensorData.LookTime(withinInterval);
    data(ind).LookAzimuth = sensorData.LookAzimuth(withinInterval);
    data(ind).LookElevation = sensorData.LookElevation(withinInterval);
    withinInterval = sensorData.DetectionTime >=frameTime & sensorData.DetectionTime < frameTime+trackerInterval;
    data(ind).DetectionTime = sensorData.DetectionTime(withinInterval);
    if isfield(sensorData,"Azimuth")
        data(ind).Azimuth = sensorData.Azimuth(1,withinInterval);
        data(ind).AzimuthAccuracy = sensorData.AzimuthAccuracy(1,withinInterval);
    end
    if isfield(sensorData,"Elevation")
        data(ind).Elevation = sensorData.Elevation(1,withinInterval);
        data(ind).ElevationAccuracy = sensorData.ElevationAccuracy(1,withinInterval);
    end
    if isfield(sensorData,"Range")
        data(ind).Range = sensorData.Range(1,withinInterval);
        data(ind).RangeAccuracy = sensorData.RangeAccuracy(1,withinInterval);
    end
    if isfield(sensorData,"RangeRate")
        data(ind).RangeRate = sensorData.RangeRate(1,withinInterval);
        data(ind).RangeRateAccuracy = sensorData.RangeRateAccuracy(1,withinInterval);
    end
    if isfield(sensorData,"TargetIndex")
        data(ind).TargetIndex = sensorData.TargetIndex(1,withinInterval);
    end
    frameTime = frameTime + trackerInterval;
    ind = ind + 1;
end
end

Related Topics