Main Content

Receiver Operating Characteristic to Tracker Operating Characteristic

This example shows how to assess the probability of target track and the probability of false track. The example discusses different methods to perform these calculations with varying levels of fidelity and computation time.

In the assessment of tracker performance, there are four types of probabilities that are often desired:

  1. Probability of single target track in the absence of false alarms (Pfa=0)

  2. Probability of single false track in the absence of targets(Pd=0)

  3. Probability of single target track in the presence of false alarms

  4. Probability of single false track in the presence of targets

This example will first delve into the foundations of the Bernoulli sum, which will help us calculate the probability of single target track in the absence of false alarms. We will then discuss the Common Gate History (CGH) algorithm that can be used to calculate all 4 types of probabilities and introduce the concept of the Tracker Operating Characteristic (TOC). The CGH algorithm provides an estimate of system capacity and thus offers a means to assess end-to-end system performance. We will then present an automotive radar design scenario and demonstrate the power of the CGH algorithm for high-level system analysis.

Calculating Single Target Track Probability in the Absence of False Alarms

Bernoulli Sum

The Bernoulli sum allows for quick and easy performance analysis in the case of a single target in the absence of false alarms. The track detection probability Pdt can be defined in terms of the receiver detection probability Pd for a window period defined as

Tw=NT,

where T is the basic sampling period and N represents the number of opportunities for a detection.

For a confirmation threshold logic of M-of-N, the target track probability Pdt is defined as

Pdt=i=MNC(N,i)Pdi(1-Pd)(N-i),

where

C(N,i)=N!(N-i)!i!.

The confirmation threshold logic denoted as M-of-N or M/N is a one-stage logic where a track must associate to a detection, also known as a hit, at least M times out of N consecutive looks. For example, consider a 2-of-3 logic. In the figure below, the solid yellow represents a hit, which can be from either a target or a false alarm. The patterned blue blocks represent misses. The cases represented below are not intended to be exhaustive, but as can be seen, cases 2 and 3 satisfy the threshold, but case 1 does not.

Investigate the probability of target track versus probability of detection for a confirmation threshold of 2/3 using the Bernoulli sum method. Note how quickly the computation is performed. Perform the Bernoulli sum calculation assuming that one hit is required to initialize a track.

% Define probabilities for analysis 
Pd = linspace(0,1,100).';

% Define confirmation threshold M/N
M = 2; % Number of hits
N = 3; % Number of observations or opportunities

% Perform Bernoulli sum assuming 1 hit is required to initialize a track
tic 
PdtBernoulli = helperBernoulliSum(Pd,M,N);
elapsedTime = toc;
helperUpdate('Bernoulli',elapsedTime);
Bernoulli calculation completed. Total computation time is 0.0182 seconds.
% Plot the probability of detection versus the probability of target track 
hAxes = helperPlot(Pd,PdtBernoulli,'M/N = 2/3','P_D','P_{DT}', ...
    sprintf('Bernoulli Sum\nProbability of Target Track in the Absence of False Alarms'));

% Set desired probability of target track
yline(hAxes,0.9,'--','LineWidth',2,'DisplayName','Desired P_{DT}'); 

Figure Bernoulli Sum Probability of Target Track in the Absence of False Alarms contains an axes. The axes with title Bernoulli Sum Probability of Target Track in the Absence of False Alarms contains 2 objects of type line, constantline. These objects represent M/N = 2/3, Desired P_{DT}.

Assuming that a probability of target track of 0.9 was desired, the above plot informs an analyst that a detection probability of about 0.7 is required.

Calculating Probabilities for a Target in Clutter

The values obtained from a Bernoulli sum calculation are useful in quick analyses, but are not generally representative of real tracking environments where target tracks are affected by the presence of false alarms. Consider the scenario where a target is operating in the presence of clutter.

Assuming that false alarms occur on a per-look per-cell basis, the probability of false alarms in a tracking gate depend upon the number of cells in the gate. Assume that there are three types of events,

  • Miss: No detection

  • Hit: False alarm

  • Hit: Target detection

The number of cells in a gate depends upon the history of events and the order in which events occur. These factors dictate the gate growth sequence of the tracker.

The Bernoulli sum method assumes that there are no false alarms and that the order of detections does not matter. Thus, when using the Bernoulli sum in the target in clutter scenario, it produces overly optimistic results.

One approach that is used to analyze such scenarios is to evaluate every possible track sequence and determine which sequences satisfy the confirmation threshold logic. This brute-force approach of building the Markov chain is generally too computationally intensive.

Another approach is to utilize a Monte Carlo-type analysis. Rather than generating the full Markov chain, a Monte Carlo simulation manually generates random sequences of N events. The confirmation threshold is applied to each sequence, and statistics are aggregated. The Monte Carlo method is based upon the law of large numbers, so performance improves as the number of iterations increases. Monte Carlo analysis lends itself well to parallelization, but in the case of small probabilities of false alarm, the number of iterations can become quickly untenable. Thus, we need to consider alternative methods to quickly calculate track probability measures.

The Common Gate History (CGH) Algorithm

The Common Gate History algorithm greatly reduces computation times and memory requirements. The algorithm prevents the need for manual generation of sequences as in the case of the Monte Carlo analysis, which can be costly for low-probability events.

The algorithm begins by making the assumption that there are three types of tracks, which can contain:

  1. Detections from targets

  2. Detections from targets and false alarms

  3. Detections from false alarms only

A target track is defined as any track that contains at least one target detection and satisfies the M/N confirmation threshold. Thus, track types 1 and 2 are considered to be target tracks, whereas 3 is considered to be a false track.

Given the previously defined track types and with a confirmation threshold logic M/N, a statistic is proposed to identify a track state. A track state is defined as

ω=[ωl,ωlt,λ]

where ωl is the number of time steps since the last detection (target or false alarm), ωltis the number of time steps since the last target detection, and λ is the total count of detections (targets or false alarms). As the algorithm proceeds, the track state vector evolves according to a Markov chain.

The algorithm assumes that a track can be started given only two types of events:

  • Target detection

  • False alarm

Once a track is initiated, the following four types of events continue a track:

  • No detection

  • Target detection

  • False alarm

  • Target detection and false alarm

The track probability at look m is multiplied by the probability of the event that continues the track at look m+1. The tracks are then lumped by adding the track probabilities of the track files with a common gate history vector. This keeps the number of track states in the Markov chain within reasonable limits.

The assumptions of the Common Gate History algorithm are as follows:

  • The probability of more than one false alarm in a gate is low, which is true when the probability of false alarm is low (10-3 or less).

  • Assumes a uniform spatial distribution of the location of a target in a gate.

  • Assumes a track splitting algorithm is used.

The Common Gate History algorithm can be used to calculate all four probability types:

  1. Probability of single target track in the absence of false alarms (Pfa=0)

  2. Probability of single false track in the absence of targets(Pd=0)

  3. Probability of single target track in the presence of false alarms

  4. Probability of single false track in the presence of targets

The probability of target track and the probability of false track form the basis of the Tracker Operating Characteristic (TOC). The (TOC) complements the Receiver Operating Characteristic (ROC), which is commonly used in the analysis and performance prediction of receivers. Combining the ROC and the TOC provides an end-to-end system analysis tool.

Calculate and plot the TOC using ROC curves from rocsnr as inputs. Assume a signal-to-noise ratio (SNR) of 8 dB. Continue to use the 2/3 confirmation threshold logic that was used in the Bernoulli sum example. Use the toccgh built-in tracker. Note how quickly the computation is performed despite the low probabilities of false alarms.

% Receiver operating characteristic (ROC)
snrdB = 8; % SNR (dB) 
[Pd,Pfa] = rocsnr(snrdB,'MaxPfa',1e-3,'MinPfa',1e-12,'NumPoints',20);

% Plot ROC
helperPlotLog(Pfa,Pd,snrdB, ...
    'Probability of False Alarm (P_{FA})', ...
    'Probability of Detection (P_D)', ...
    'Receiver Operating Characteristic');

Figure Receiver Operating Characteristic contains an axes. The axes with title Receiver Operating Characteristic contains an object of type line. This object represents SNR (dB) = 8.00.

% Common gate history algorithm
tic
[PdtCGH,PftCGH] = toccgh(Pd,Pfa,'ConfirmationThreshold',[M N]);
elapsedTime = toc;
helperUpdate('Common Gate History',elapsedTime);
Common Gate History calculation completed. Total computation time is 0.5302 seconds.
% Plot Common Gate History results 
hAxes = helperPlotLog(PftCGH,PdtCGH,'CGH','P_{FT}','P_{DT}', ...
    'Tracker Operating Characteristic (TOC) Curve');

Figure Tracker Operating Characteristic (TOC) Curve contains an axes. The axes with title Tracker Operating Characteristic (TOC) Curve contains an object of type line. This object represents CGH.

The Common Gate History algorithm permits the assessment of tracker performance similar to a Monte Carlo analysis but with acceptable computation times despite low-probability events. The Common Gate History algorithm thus permits high-level investigation and selection of options prior to more intensive, detailed simulations.

Using CGH with Custom Trackers

Consider developing a tracker for an automotive application. Define a custom one-dimensional, nearly constant velocity (NCV) tracker using trackingKF. Assume that the update rate Δt is 1 second. Assume that the state transition matrix is of the form

A=[1Δt01]

and the process noise is of the form

Q=[Δt44Δt32Δt32Δt2]q,

where q is a tuning factor defined as

q=amax2Δt.

The input amax2is the maximum target acceleration expected. Assume a maximum accelaration of 4 m/s2 is expected for vehicles.

% Define the state transition matrix
dt = 1; % Update rate (sec)
A = [1 dt; 0 1];

% Define the process noise 
Q = [dt^4/4 dt^3/2; dt^3/2 dt^2];

% Tune the process noise
amax = 4; % Maximum target acceleration (m/s^2)
q = amax^2*dt; 

% Update the process noise 
Q = Q.*q; 

% Initialize the Kalman filter 
trkfilt = trackingKF('MotionModel','Custom', ...
    'StateCovariance', [0 0; 0 0], ...
    'StateTransitionModel',A, ...
    'ProcessNoise',Q, ...
    'MeasurementNoise',0, ...
    'MeasurementModel',[1 0]);

An error ellipse is used to model track uncertainty. From this uncertainty ellipse, the gate growth sequence can be calculated.

The 1-σ values of the error ellipse are calculated as the square root of the eigenvalues λ of the predicted state covariance Pk+1|k:

[λ1,λ2]=eig(Pk+1|k).

The area of the error ellipse is then calculated as

ErrorEllipseArea=πλ1λ2.

The area of the bins is calculated as

BinArea=(ΔRange)(ΔRangeRate).

Finally, the gate size in bins is

GateSize=ErrorEllipseAreaBinArea.

The gate size is thus dependent upon the tracker, the event sequence, and the resolution of the bins. Calculate the gate growth sequence assuming a confirmation threshold N equal to 3. Assume the range and range rate resolutions for the automotive radar are equal to 1 m and 1 m/s, respectively.

% Calculate gate growth sequence
res = [1, 1]; % Bin resolutions [range (m), range-rate (m/s)]
gs = helperCalculateGateSize(N,trkfilt,res) 
gs = 1×3

     1    51   124

% Common gate history algorithm
tic
[PdtCGHcustom,PftCGHcustom] = toccgh(Pd,Pfa,'ConfirmationThreshold',[M N],'GateGrowthSequence',gs);
elapsedTime = toc;
helperUpdate('Common Gate History',elapsedTime);
Common Gate History calculation completed. Total computation time is 0.2272 seconds.
% Add plot to previous plot
helperAddPlotLog(hAxes,PftCGHcustom,PdtCGHcustom,'CGH with Custom Gate Growth Sequence');

Figure Tracker Operating Characteristic (TOC) Curve contains an axes. The axes with title Tracker Operating Characteristic (TOC) Curve contains 2 objects of type line. These objects represent CGH, CGH with Custom Gate Growth Sequence.

Tracker Performance Assessment for an Automotive Radar System

Probability of False Alarm and Probability of Target Track Requirements

By using the ROC and TOC in conjunction, a system analyst can select a detector's operation point that satisfies the overall system requirements. Consider an automotive radar case. Due to the nature of the application, it is desired that false alarms remain very low probability events. Additionally, the probability of target track should be high for safety purposes. Consider the following two requirements.

  • Requirement 1: Probability of false alarm must be less than 10-6

  • Requirement 2: Probability of target track must be equal to 0.9 or above

Calculate the ROC curves for SNRs equal to 6, 8, 10, and 12 dB.

% Calculate ROC curves
snrdB = [6 8 10 12]; % SNR (dB)
[Pd,Pfa] = rocsnr(snrdB,'MaxPfa',1e-3,'MinPfa',1e-10,'NumPoints',10);

Use the calculated ROC curves as inputs to the toccgh function to generate the associated TOC curves. Use the same confirmation threshold and gate growth sequence previously generated.

% Generate TOC curves
tic
toccgh(Pd,Pfa,'ConfirmationThreshold',[M N],'GateGrowthSequence',gs);
elapsedTime = toc;
helperUpdate('Common Gate History',elapsedTime);
Common Gate History calculation completed. Total computation time is 0.6146 seconds.
% Requirement 1: Probability of false alarms must be less than 1e-6
hAxesROC = subplot(2,1,1);
xlim([1e-10 1e-2])
ylim([0 1])
reqPfa = 1e-6; 
helperColorZonesReqPfa(hAxesROC,reqPfa)
legend(hAxesROC,'Location','eastoutside')

% Requirement 2: Probability of target track must be equal to 0.9 or above
hAxesTOC = subplot(2,1,2);
xlim([1e-14 1e-4])
ylim([0 1])
reqPdt = 0.9; 
helperColorZonesReqPdt(hAxesTOC,reqPdt)

Figure Operating Characteristics contains 2 axes. Axes 1 with title Receiver Operating Characteristic (ROC) Curve contains 47 objects of type patch, line, text, constantline. These objects represent Requirement Not Met, Requirement Met, ROC 1, ROC 2, ROC 3, ROC 4, Requirement. Axes 2 with title Tracker Operating Characteristic (TOC) Curve contains 47 objects of type patch, line, text, constantline. These objects represent Requirement Not Met, Requirement Met, Requirement.

Requirement 1 dictates that only points 1 through 6 on the ROC curves may be included for further analysis since they satisfy a probability lower than 10-6. Points 7 through 10 do not meet this requirement.

Considering requirement 2, it can be seen that the ROC curve corresponding to 6 dB SNR does not meet the second requirement at any point. The only curves to continue considering are the 8, 10, and 12 dB curves. Requirement 2 is met only for point 10 on the 8 dB curve, points 9 and 10 on the 10 dB curve, and points 5 through 10 on the 12 dB curve.

Combining requirements 1 and 2, there are only two analysis points that satisfy both requirements: Points 5 and 6 on the 12 dB curve. Point 5 corresponds to a probability of target track of 0.90 and a probability of false track of 1.27×10-13. Translating these back to the ROC curve, this corresponds to a probability of detection of 0.68 and a probability of false alarm of 1.29×10-7. Similarly, point 6 corresponds to a probability of target track of 0.96, probability of false track of 1.64×10-12, a probability of detection of 0.80, and a probability of false alarm of 7.74×10-7. In this case, operating point 6 is selected. It is a tradeoff of improved probability of target track at the expense of a slightly higher but reasonable probability of false track.

The Common Gate History algorithm permits an estimation of the expected number of false tracks based on the number of targets anticipated in an environment and the number of cells in the radar data. The expected number of false tracks Eft is calculated as

Eft=Pft,ntNc+PftNt,

where Pft,nt is the probability of false track in the absence of targets, Nc is the number of cells, Pft is the probability of false track in the presence of targets, and Nt is the number of targets.

Consider an environment where the number of targets are expected to be equal to 10 and the number of cells is equal to

NumberofCells=(NumberofRangeCells)×(NumberofRangeRateCells)=1000×100=105.

% Calculate expected number of false tracks using toccgh
numCells = 1e5;         % Number of cells in radar data
numTargets = 10;        % Number of targets in scenario
selectedPd = Pd(6,4) ;  % Selected probability of detection
selectedPfa = Pfa(6);   % Selected probability of false alarm
[Pdt,Pft,Eft] = toccgh(selectedPd,selectedPfa, ...
    'ConfirmationThreshold',[M N],'GateGrowthSequence',gs, ...
    'NumCells',numCells,'NumTargets',10);

% Output results
helperPrintTrackProbabilities(Pdt,Pft,Eft);
Probability of Target Track in Presence of False Alarms = 0.9581
Probability of False Track in the Presence of Targets = 1.6410e-12
Expected Number of False Tracks = 5

Thus, based on the system parameters, an analyst could expect about five false tracks.

Analysis of Confirmation Thresholds

Consider the same automotive radar design case but investigate the effect of confirmation thresholds 2/4, 3/4, and 4/4. Assume the following:

  1. The probability of false alarm is equal to 10-6

  2. The objective probability of target track is 0.9 or better

First, calculate the ROC curve using the rocpfa function.

% Calculate ROC curve assuming a probability of false alarm of 1e-6
Pfa = 1e-6;
numPoints = 20;
[Pd,snrdB] = rocpfa(Pfa,'NumPoints',numPoints,'MaxSNR',15);

Update the gate growth sequence due to the longer number of observations.

% Update the gate growth sequence
N = 4; 
trkfilt = trackingKF('MotionModel','Custom', ...
    'StateCovariance', [0 0; 0 0], ...
    'StateTransitionModel',A, ...
    'ProcessNoise',Q, ...
    'MeasurementNoise',0, ...
    'MeasurementModel',[1 0]);
gs = helperCalculateGateSize(N,trkfilt,res)
gs = 1×4

     1    51   124   225

Calculate the TOC given the ROC curves as inputs. Analyze for confirmation thresholds equal to 2/4, 3/4, and 4/4.

% Calculate TOC
cp = [2 4; 3 4; 4 4];
numCp = size(cp,1);
PdtMat = zeros(numPoints,numCp);
PftMat = zeros(numPoints,numCp);
EftMat = zeros(numPoints,numCp);
for ii = 1:numCp
    [PdtMat(:,ii),PftMat(:,ii),EftMat(:,ii)] = toccgh(Pd.',Pfa, ...
        'ConfirmationThreshold',cp(ii,:),'GateGrowthSequence',gs, ...
        'NumCells',numCells,'NumTargets',10);
end

% Plot ROC and TOC
reqPdt = 0.9; 
helperPlotROCTOC(reqPdt,Pfa,Pd,snrdB,PdtMat,cp);

Figure contains 2 axes. Axes 1 with title Receiver Operating Characteristic (ROC) contains 21 objects of type line, text. This object represents 1.0e-06. Axes 2 with title SNR versus P_{DT} contains 66 objects of type patch, line, text, constantline. These objects represent Requirement Not Met, Requirement Met, 2/4, 3/4, 4/4, Requirement.

helperPrintReqValues(reqPdt,Pd,snrdB,PdtMat,EftMat,cp);
Confirmation Threshold            = 2/4
Required Probability of Detection = 0.55
Required SNR (dB)                 = 10.76
Expected Number of False Tracks   = 18

Confirmation Threshold            = 3/4
Required Probability of Detection = 0.81
Required SNR (dB)                 = 12.03
Expected Number of False Tracks   = 1

Confirmation Threshold            = 4/4
Required Probability of Detection = 0.97
Required SNR (dB)                 = 13.37
Expected Number of False Tracks   = 1

Reviewing the results, it can be seen that the more stringent the confirmation threshold, the higher the required SNR. However, the more stringent confirmation thresholds result in improved numbers of false tracks.

Summary

In the assessment of tracker performance, there are four types of probabilities that are often desired:

  1. Probability of single target track in the absence of false alarms (Pfa=0)

  2. Probability of single false track in the absence of targets(Pd=0)

  3. Probability of single target track in the presence of false alarms

  4. Probability of single false track in the presence of targets

For the calculation of 1, Bernoulli sums can be used. However, to obtain the other probabilities, a different method must be employed. While Monte Carlo analysis can be used for the computation of the latter three types of probabilities, the computational resources and time required can become untenable, which is particularly true for events with low probability. The Common Gate History (CGH) algorithm can be used to calculate all four quantities and greatly reduces computational resources needed.

The Common Gate History algorithm can be used to generate the Tracker Operating Characteristic (TOC). The TOC compliments the Receiver Operating Characteristic (ROC) and provides a means to assess overall system performance. The TOC and ROC curves can be used in a multitude of ways such as determining:

  • Which track confirmation threshold should be selected

  • Which detection threshold on the ROC curve should be used

  • Whether the system SNR can meet the track probablity requirements.

Finally, the CGH algorithm permits the calculation of an expected number of false tracks, which offers insights into system capacity. The expected number of false tracks can be used to ascertain computational load and to assist with decisions related to hardware and processing.

References

  1. Bar-Shalom, Y., L. J. Campo, and P. B. Luh. "From Receiver Operating Characteristic to System Operating Characteristic: Evaluation of a Track Formation System." IEEE Transactions on Automatic Control. Vol. 35, No. 2, February 1990.

  2. Bar-Shalom, Y., P. K. Willett, and X. Tian. Tracking and Data Fusion: A Handbook of Algorithms. Storrs, CT: YBS Publishing, 2011.

Helper Functions

function Pcnf = helperBernoulliSum(Pd,Mc,Nc)
% Calculate simple Bernoulli sum. Use the start TOC logic, which assumes
% that there is already one hit that initializes the track. 

% Update M and N for probability of deletion
Nd = Nc - 1;      % Need one hit to start counting. Assume first hit initializes track. 
Md = Nc - Mc + 1; % Need this many misses to delete

% Bernoulli sum. Probability of deletion calculation.
ii = Md:Nd;
C = arrayfun(@(k) nchoosek(Nd,k),ii);
P = (1 - Pd);
Pdel = sum(C.*P(:).^ii.*(1 - P(:)).^(Nd - ii),2);

% Probability of confirmation
Pcnf = 1 - Pdel;
end

function helperUpdate(calculationType,elapsedTime)
% Output elapsed time 
fprintf('%s calculation completed. Total computation time is %.4f seconds.\n', ...
        calculationType,elapsedTime);
end

function varargout = helperPlot(x,y,displayName,xAxisName,yAxisName,titleName,varargin)
% Create a plot with logarithmic scaling on the x-axis

% Create a figure
figure('Name',titleName)
hAxes = gca;

% Plot data
plot(hAxes,x,y,'LineWidth',2,'DisplayName',displayName,varargin{:})
hold(hAxes,'on')
grid(hAxes,'on')

% Update axes
hAxes.Title.String = titleName;
hAxes.XLabel.String = xAxisName;
hAxes.YLabel.String = yAxisName;

% Make sure legend is on and in best location
legend(hAxes,'Location','Best')

% Set axes as optional output 
if nargout == 1
    varargout{1} = hAxes;
end
end

function varargout = helperPlotLog(x,y,displayName,xAxisName,yAxisName,titleName,varargin)
% Create a plot with logarithmic scaling on the x-axis

% Create a figure
figure('Name',titleName)
hAxes = gca;

% Plot data
numCol = size(y,2);
for ii = 1:numCol
    idxX = min(ii,size(x,2));
    hLine = semilogx(hAxes,x(:,idxX),y(:,ii),'LineWidth',2,varargin{:});
    if ischar(displayName)
        hLine.DisplayName = displayName;
    else
        hLine.DisplayName = sprintf('SNR (dB) = %.2f',displayName(ii));
    end
    hold on
end
grid on

% Update axes
hAxes.Title.String = titleName;
hAxes.XLabel.String = xAxisName;
hAxes.YLabel.String = yAxisName;

% Make sure legend is on and in best location 
legend(hAxes,'Location','Best')

% Set axes as optional output 
if nargout == 1
    varargout{1} = hAxes;
end
end

function helperAddPlotLog(hAxes,x,y,displayName,varargin)
% Add an additional plot to the axes hAxes with logarithmic scaling on the
% x-axis

% Plot data
hold on
hLine = semilogx(hAxes,x,y,'LineWidth',2,varargin{:});
hLine.DisplayName = displayName;
end

function gs = helperCalculateGateSize(N,trkfilt,res)
% Calculate a gate growth sequence in bins

% Initialize tracker gate growth sequence
gs = zeros(1,N); % Gate growth sequence

% Calculate gate growth sequence by projecting state uncertainty using
% linear approximations.
areaBin = prod(res(:),1);
for n = 1:N
    [~,Ppred] = predict(trkfilt); % Predict
    
    % Calculate the products of the 1-sigma values
    E = eig(Ppred);
    E(E<0) = 0; % Remove negative values
    sigma1Prod = sqrt(prod(E(:),1));
    
    % Calculate error ellipse area
    areaErrorEllipse = pi*sigma1Prod; % Area of ellipse = pi*a*b
    
    % Translate to bins
    gs(n) = max(ceil(areaErrorEllipse/areaBin),1);
end
end

function helperColorZonesReqPfa(hAxes,req)
% Plot color zones for requirement type 1 

% Vertical requirement line
xline(req,'--','DisplayName','Requirement',...
    'HitTest','off'); 

% Get axes limits
xlims = get(hAxes,'XLim');
ylims = get(hAxes,'YLim');

% Green box
pos = [xlims(1) ylims(1) req ylims(2)];
x = [pos(1) pos(1) pos(3) pos(3) pos(1)];
y = [pos(1) pos(4) pos(4) pos(1) pos(1)];
hP = patch(hAxes,x,y,[0.4660 0.6740 0.1880], ...
    'FaceAlpha',0.3,'EdgeColor','none','DisplayName','Requirement Met');
uistack(hP,'bottom');

% Red box
pos = [req ylims(1) xlims(2) ylims(2)];
x = [pos(1) pos(1) pos(3) pos(3) pos(1)];
y = [pos(1) pos(4) pos(4) pos(1) pos(1)];
hP = patch(hAxes,x,y,[0.6350 0.0780 0.1840], ...
    'FaceAlpha',0.3,'EdgeColor','none','DisplayName','Requirement Not Met',...
    'HitTest','off');
uistack(hP,'bottom');
end

function helperColorZonesReqPdt(hAxes,req)
% Plot color zones for requirement type 2 

% Horizontal requirement line
yline(req,'--','DisplayName','Requirement',...
    'HitTest','off')

% Get axes limits
xlims = get(hAxes,'XLim');
ylims = get(hAxes,'YLim');

% Green box
pos = [xlims(1) req xlims(2) ylims(2)];
x = [pos(1) pos(1) pos(3) pos(3) pos(1)];
y = [pos(1) pos(4) pos(4) pos(1) pos(1)];
hP = patch(hAxes,x,y,[0.4660 0.6740 0.1880], ...
    'FaceAlpha',0.3,'EdgeColor','none','DisplayName','Requirement Met',...
    'HitTest','off');
uistack(hP,'bottom');

% Red box
pos = [xlims(1) req xlims(2) req];
x = [pos(1) pos(1) pos(3) pos(3) pos(1)];
y = [pos(1) pos(4) pos(4) pos(1) pos(1)];
hP = patch(hAxes,x,y,[0.6350 0.0780 0.1840], ...
    'FaceAlpha',0.3,'EdgeColor','none','DisplayName','Requirement Not Met',...
    'HitTest','off');
uistack(hP,'bottom');
end

function helperPrintTrackProbabilities(Pdt,Pft,Eft)
% Print out results 

fprintf('Probability of Target Track in Presence of False Alarms = %.4f\n',Pdt)
fprintf('Probability of False Track in the Presence of Targets = %.4e\n',Pft)
fprintf('Expected Number of False Tracks = %d\n',Eft)
end

function helperPlotROCTOC(reqPdt,Pfa,Pd,snrdB,PdtMat,cp)
% Plot ROC/TOC 

% Plot ROC curves
figure
hAxesROC = subplot(2,1,1);
plot(hAxesROC,snrdB,Pd,'-o')
title(hAxesROC,'Receiver Operating Characteristic (ROC)')
xlabel(hAxesROC,'SNR (dB)')
ylabel(hAxesROC,'P_D')
grid(hAxesROC,'on')
legend(hAxesROC,sprintf('%.1e',Pfa),'Location','Best')

% Plot SNR versus probability of target track
hAxesTOC = subplot(2,1,2);
numCp = size(cp,1);
for ii = 1:numCp
    plot(hAxesTOC,snrdB,PdtMat(:,ii),'-o', ...
        'DisplayName',sprintf('%d/%d',cp(ii,1),cp(ii,2)))
    hold(hAxesTOC,'on')
end
title(hAxesTOC,'SNR versus P_{DT}')
xlabel(hAxesTOC,'SNR (dB)')
ylabel(hAxesTOC,'P_{DT}')
grid(hAxesTOC,'on')
legend(hAxesTOC,'Location','Best')

% Label points
colorVec = get(hAxesROC,'ColorOrder'); 
numSnr = numel(snrdB);
textArray = arrayfun(@(x) sprintf('  %d',x),1:numSnr,'UniformOutput',false).';
xPosROC = snrdB;
colorFont = brighten(colorVec,-0.75);
numColors = size(colorVec,1);
idxC = mod(1:numCp,numColors); % Use only available default colors
idxC(idxC == 0) = numColors; % Do not let color index equal 0

% Label points ROC
yPosROC = Pd;
text(hAxesROC,xPosROC,yPosROC,textArray,'FontSize',6,'Color',colorFont(1,:),'Clipping','on')

% Label points TOC
xPosTOC = snrdB;
for ii = 1:numCp
    yPosTOC = PdtMat(:,ii);
    text(hAxesTOC,xPosTOC,yPosTOC,textArray,'FontSize',6,'Color',colorFont(idxC(ii),:),'Clipping','on')
end

% Add requirement zone color blocks
helperColorZonesReqPdt(hAxesTOC,reqPdt);
end

function helperPrintReqValues(reqPdt,Pd,snrdB,PdtMat,EftMat,cp)
% Output information about required values given a required probability of
% target track requirement

% Get values 
numCp = size(PdtMat,2);
reqPd = zeros(1,numCp);
reqSNRdB = zeros(1,numCp); 
expEft = zeros(1,numCp);
for ii = 1:numCp
    reqPd(ii) = interp1(PdtMat(:,ii),Pd,reqPdt);
    reqSNRdB(ii) = interp1(PdtMat(:,ii),snrdB,reqPdt);
    expEft(ii) = interp1(PdtMat(:,ii),EftMat(:,ii),reqPdt);
end

 % Display required probability of detection, SNR, and expected false
 % tracks
 for ii = 1:numCp
     fprintf('Confirmation Threshold            = %d/%d\n',cp(ii,1),cp(ii,2));
     fprintf('Required Probability of Detection = %.2f\n',reqPd(ii));
     fprintf('Required SNR (dB)                 = %.2f\n',reqSNRdB(ii));
     fprintf('Expected Number of False Tracks   = %d\n\n',expEft(ii));
 end
end