function [CRLB,Aux] = XcvrTOALocCRLB(uo,t0o,vo,s,K,Q,OfstUnkn,OfstBias)
% function [CRLB,Aux] = XcvrTOALocCRLB(uo,t0o,vo,s,K,Q,OfstUnkn,OfstBias)
%
% CRLB for distant object localization using moving transceivers
%
% Input parameter list:
% uo: (dim x 1), object location
% t0o:((M dim) x M), initial positions of M transceivers
% vo: ((M dim) x M), velocities of M transceivers
% s:  (Dim x N), N sensor positions
% K:  K = num. of time instants
% Q:  (KMN x KMN), covariance matrix of distance measurements
% OfstUnkn:  flag, 1 if measurement offsets known and 0 if offsets unknown
% OfstBias:  residual bias offsets when OfstUnkn=0
%
% Output parameter list:
% CRLB: CRLB matrix,
%       ((L(2M+1)+M) x (L(2M+1)+M)) when OfstUnkn=1,
%       ((L(2M+1)) x (L(2M+1))) when OfstUnkn=0
% Aux:  auxiliary values related to bias when OfstUnkn=1
%
% The code can work with 2-D(Dim=2) or 3-D(Dim=3) localization.
%
% Reference:
% K. C. Ho, "Localization through transceivers in unknown constant velocity
% trajectories," IEEE Trans. Signal Process., vol. 70, pp. 3011-3028, 2022.
%
% K. C. Ho     08-28-2022
%
%       Copyright (C) 2022
%       Computational Intelligence Signal Processing Laboratory
%       University of Missouri
%       Columbia, MO 65211, USA
%       hod@missouri.edu
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[L,M]=size(t0o); N=size(s,2);

k=0:K-1;

for i=1:M,
    tio=repmat(t0o(:,i),1,K)+vo(:,i)*k;
    liuo(:,i)=sqrt(sum((tio-repmat(uo,1,K)).^2));
    rhoiuo=-(tio-repmat(uo,1,K))./repmat(liuo(:,i)',L,1);
    
    Gi=[];
    Hi=[];
    for j=1:N,
        
        lijo(:,i,j)=sqrt(sum((tio-repmat(s(:,j),1,K)).^2));
        rhoijo=(tio-repmat(s(:,j),1,K))./repmat(lijo(:,i,j)',L,1);
        if (OfstUnkn)
            Gij=[rhoiuo; (-rhoiuo+rhoijo); repmat(k,L,1).*(-rhoiuo+rhoijo); ones(1,length(rhoiuo))]';
        else
            Gij=[rhoiuo; (-rhoiuo+rhoijo); repmat(k,L,1).*(-rhoiuo+rhoijo)]';
        end;
        Gi=[Gi; Gij];
    end;
    
    if (i>1)
        G=[G, zeros((i-1)*N*K,size(Gi,2)-L); Gi(:,1:L), zeros(N*K,size(G,2)-L), Gi(:,L+1:end)];
    else
        G=Gi;
    end;
    
end;

% -- rearrange G and H with all time offsets at the end --
if (OfstUnkn)
    tpG=G;
    for i=1:M,
        G(:,L+(2*L)*(i-1)+1:L+(2*L)*(i))=tpG(:,L+(2*L+1)*(i-1)+1:L+(2*L+1)*(i-1)+2*L);
        G(:,end-M+i)=tpG(:,L+(2*L+1)*(i-1)+2*L+1);
    end;
end;

FIM=G'*inv(Q)*G;
CRLB=inv(FIM);

if (~OfstUnkn)
    BIAS=inv(FIM)*G'*inv(Q)*ones(M*N*K,1)*OfstBias;
    Aux.BIAS=BIAS;
else
    Aux.BIAS=0;
end;


