xpbias.m

GPS TOOLBOX包含以下内容： 1、GPS相关常量和转换因子； 2、角度变换； 3、坐标系转换： &#61656 点变换； &#61656 矩阵变换； &#61656 向量变换

%                             xpbias.m
%  Scope:   This MATLAB program computes the value of the parity bias (pbiasb); 
%           for the 5 satellites case the Gaussian distribution is used, while 
%           for 6 up to 14 satellites case the Chi-square distribution is used 
%           [1].
%  Usage:   xpbias
%  Inputs:  - probability of missed detection (pmd)
%           - td_table, table containing fault detection normalized thresholds 
%             td for degree-of-freedom (dof) 1 to 10 and the corresponding 
%             probability of false alarm (pfa); a default table is given for 
%             pfa = 3.3333e-7
%           - for degree-of-freedom  = 1 (Gausssian ditribution)  
%                                         --> initial guess for pbiasb
%             for degree-of-freedom >= 2  (Chi-square distribution) 
%                                         --> initial guess for lambda
%  Outputs: - input data (normalized thresholds, pfa and pmd) for dof 1 to 10
%           - parity bias (pbiasb) for a specified set of input data for dof
%             1 to 10 
%  Reference: 
%           [1] Brown, R. G., Chin, G. V., GPS RAIM: Calculation of thresholds 
%               and protection radius using Chi-square methods - A geometric 
%               approach. RTCA Paper No. 491-94/SC159-584, Nov. 1994
%  Remarks:  1) Table td_table can be generated by using the program XFDNT.m.
%            2) The parity bias values can be used in the computation of horizon-
%               tal protection limit, and fault detection availability in the
%               case of RAIM with constant false alarm 
%  External Matlab macros used: fbias_c.m, fbias_g.m, ncchi2_2.m, ncchi2_3.m,
%               ncchi2_4.m, ncchi2_5.m , ncchi2_6.m, ncchi2_7.m, , ncchi2_8.m 
%               ncchi2_9.m, ncchi2_10.m 
%  Last update: 11/23/00
%  Copyright (C) 1997-00 by LL Consulting. All Rights Reserved.

clear

global   a dof pmd td 

yes = 'y';

%  Enter input data

disp('  ');
disp('Enter the data for the precomputed threshold table  td ');
answer = input('Do you want to use the default data for pfa = 3.3333e-7? (y/n) --> ','s');
if (strcmp(answer,yes) == 1)
   %  default input data for  pfa = 3.3333e-7 (probability of false alarm)
   td_table = [5.10383 5.46152 5.73854 5.97522 6.18628 6.37900 6.55774 6.72514 6.88354 7.03383];
   pfa = 3.3333e-7;
else
   td_table(1) = input('Enter value for td_table(1) --> ');
   td_table(2) = input('Enter value for td_table(2) --> ');
   td_table(3) = input('Enter value for td_table(3) --> ');
   td_table(4) = input('Enter value for td_table(4) --> ');
   td_table(5) = input('Enter value for td_table(5) --> ');
   td_table(6) = input('Enter value for td_table(6) --> ');
   td_table(7) = input('Enter value for td_table(7) --> ');
   td_table(8) = input('Enter value for td_table(8) --> ');
   td_table(9) = input('Enter value for td_table(9) --> ');
   td_table(10) = input('Enter value for td_table(10) --> ');
   pfa = input('Enter corresponding probability of false alarm (pfa) --> ');
end

disp('  ');
disp('Enter value for the probability of missed detection (pmd) ');
answer = input('Use the default data, pmd = 0.001 ? (y/n) --> ','s');
if (strcmp(answer,yes) == 1)
   pmd = 0.001;        %  pmd - probability of missed detection
else
   pmd = input('Enter probability of missed detection --> ');
end
if (pmd < 0) | (pmd > 1)
   error('Error - XPBIAS ; check the input data value for pmd ');
end

%  Save computed results into an external file or display on screen

disp('  ');
disp('Enter the output desired (output file or on screen)');
answer = input('Do you want to display the results on screen? (y/n) --> ','s');
disp('  ');
if (strcmp(answer,yes) == 1)
   f2 = 1;          %  on screen
else
   f2 = input('Specify the output filename -->  ','s');
end

%  Computation for dof = 1 with Gaussian distribution

dof = 1;
disp('  '); 
disp('*****   For  dof = 1   *****'); 
answer = input('Use default for initial guess for pbiasb? (y/n) --> ','s');
if (strcmp(answer,yes) == 1)
   xinit = 7;
else
   xinit = input('Enter the initial guess for the pbiasb -->  ');
end
td = td_table(dof);               %  td = sqrt(a)
 
pbiasb(1) = fzero('fpbias_g',xinit,1.e-6);

for dof = 2:10

%  Computation for dof >= 2 with Chi-square distribution

   disp('  '); 
   disp(['*****   For  dof = ' int2str(dof) '   *****']);
   answer = input('Use default for initial guess for lambda? (y/n) --> ','s');
   if (strcmp(answer,yes) == 1)
      lam_init = 70;
   else
      lam_init = input('Enter the initial guess for the lambda -->  ');
   end
 
   a = td_table(dof)* td_table(dof); 

   lambda = fzero('fpbias_c',lam_init,1.e-6);
   pbiasb(dof) = sqrt(lambda);                    %  normalized pbiasb

end

%  Display/save the results

fprintf(f2,'\n');
fprintf(f2,'****************************************************************\n');
fprintf(f2,'\n***  COMPUTATION OF FAULT DETECTION PARITY BIAS THRESHOLDS  ****\n');
fprintf(f2,'                  for RAIM with constant false alarm          \n\n');
fprintf(f2,'*****                      Input data                      *****\n\n');
for dof = 1:10	
   fprintf(f2,'When  dof = %2.0f , td = %14.7f \n',dof,td_table(dof));
end
fprintf(f2,'\nProbability of false alarm = %e \n\n',pfa);
fprintf(f2,'Probability of missed detection = %e \n\n',pmd);

fprintf(f2,'*****                Results by using XPBIAS               *****\n\n')
for dof = 1:10	
   fprintf(f2,'When  dof = %2.0f , Normalized pbiasb = %14.7f \n', ...
           dof, pbiasb(dof));
end
fprintf(f2,'\n');
fprintf(f2,'****************************************************************\n');

disp('  ');
disp('End of the program  XPBIAS');
disp('  ');