function [dga,dgp]=ffsdcal(a,b,f,na,nb,nf,nk,GC,OM,P,T)
% FFSDCAL	Auxiliary function to TH2FF
%
%	[dga,dgp]=ffsdcal(a,b,f,na,nb,nf,nk,GC,OM,P,T)


%	L. Ljung 4-8-90
%	Copyright (c), 1990 The MathWorks, Inc.
%	All Rights Reserved

% 
% *** Now compute the standard deviations ***
%
%	D3 = " dGC/dTHETA "
%
if T>0,
    D3=[((-GC./(a*OM(1:na+1,:)))'*ones(1,na)).*OM(2:na+1,:)',((GC./(b*OM(1:length(b),:)))'*ones(1,nb)).*OM(nk+1:nk+nb,:)',((-GC./(f*OM(1:nf+1,:)))'*ones(1,nf)).*OM(2:nf+1,:)'];
else
    D3=[((-GC./(a*OM(na+1:-1:1,:)))'*ones(1,na)).*OM(na:-1:1,:)',((GC./(b*OM(length(b):-1:1,:)))'*ones(1,nb)).*OM(nb:-1:1,:)',((-GC./(f*OM(nf+1:-1:1,:)))'*ones(1,nf)).*OM(nf:-1:1,:)'];
end


    D4=D3*P;
%
%   The matrix [C1 C3;conj(C3) C2] is the covariance matrix of [Re GC; Im GC]
%   according to Gauss' approximation formula
%    
    C1=sum((real(D4).*real(D3))')';
    C2=sum((imag(D4).*imag(D3))')';
    C3=sum((imag(D4).*real(D3))')';
%
%   Now translate these covariances to those of abs(GC) and arg(GC)
%
    dga=sqrt((real(GC').^2).*C1+2*((real(GC')).*(imag(GC'))).*C3+(imag(GC').^2).*C2)./abs(GC');
    dgp=(180/pi)*sqrt((imag(GC').^2).*C1-2*((real(GC')).*imag(GC')).*C3+(real(GC').^2).*C2)./(abs(GC').^2);
