function Pd = rocpdcalc(Pfa,SNR,N,type)
% This function is for internal use and may be modified or removed in a future release.

%ROCPDCALC Calculate Pd value for ROC curves
%   Pd = ROCPDCALC(Pfa,SNR,N,TYPE) returns the probability of detection
%   Pd corresponding the given probability of false alarm Pfa, signal to
%   noise ratio SNR and number of pulses N used in pulse integration.
%   TYPE specifies the signal type and can be one of the following values:
%   [ 'Real' | 'NonfluctuatingNoncoherent' | {'NonfluctuatingCoherent'} | 
%   'Swerling1' | 'Swerling2' | 'Swerling3' | 'Swerling4' ].
%
%   Both Pfa and SNR can be vectors and Pd's dimension is given by
%   [length(Pfa) length(SNR)];

%   Copyright 2008-2018 The MathWorks, Inc.

%   Reference
%   [1] Mark Richards, Fundamentals of Radar Signal Processing

%#codegen

% Checking letters was enough and ok here because type was already validated and strcmp
% is auto extrinsic.

if type(1) == 'R' %Real
    Pd =localRealROC(Pfa,SNR,N);
elseif type(1) == 'N' 
    if type(15) == 'C' %NonfluctuatingCoherent
        Pd =localNonfluctuatingCoherentROC(Pfa,SNR,N);
    else               %NonfluctuatingNoncoherent
        Pd = phasedfusion.internal.pdNonfluctuatingNoncoherent(Pfa,SNR,N);
    end
else 
    if type(9) == '1'     %Swerling1
        Pd = phasedfusion.internal.pdSwerling1(Pfa,SNR,N);
    elseif type(9) == '2' %Swerling2
        Pd =localSwerling2ROC(Pfa,SNR,N);
    elseif type(9) == '3' %Swerling3
        Pd = phasedfusion.internal.pdSwerling3(Pfa,SNR,N);
    else %Swerling4
        Pd = phasedfusion.internal.pdSwerling4(Pfa,SNR,N);
    end
end

%-------------------------------------------------------------------------
function Pd = localNonfluctuatingCoherentROC(Pfa,SNR,N) 

SNRlen = numel(SNR);
Pd = zeros(numel(Pfa),SNRlen);  % preallocate

erfcinv2PFA = erfcinv(2*Pfa);
for k = 1:SNRlen
    Pd(:,k) = 0.5*erfc( erfcinv2PFA-(sqrt(SNR(k)*N)) );
end

%-------------------------------------------------------------------------
function Pd = localRealROC(Pfa,SNR,N)

SNRlen = numel(SNR);
Pd = zeros(numel(Pfa),SNRlen);  % preallocate

erfcinv2PFA = erfcinv(2*Pfa);
SQRT2 = sqrt(2);
for k = 1:SNRlen
    Pd(:,k) = 0.5*erfc( erfcinv2PFA-(sqrt(SNR(k)*N)/SQRT2) );
end

%-------------------------------------------------------------------------
function Pd = localSwerling2ROC(Pfa,SNR,N)

SNRlen = numel(SNR);
Pd = zeros(numel(Pfa),SNRlen);  % preallocate
T = real(gammaincinv(1-Pfa,N));
T = T(:).';

for k = 1:SNRlen
    Pd(:,k) = 1 - real(gammainc(T/(1+SNR(k)),N));
end
% [EOF]
