xhpe.m

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

%                              xhpe.m
%  Scope:   This MATLAB program determines the horizontal position error (range)
%           when latitude and longitude for two points are given; WGS-84 constants
%           are used.
%  Usage:   xhpe
%  Inputs:  - selection of the output data file or dispay the results on screen
%           - input from keyboard, latitude of the first position point, 
%             in radians or degrees/minutes/seconds
%           - input from keyboard, longitude of the first position point, 
%             in radians or degrees/minutes/seconds
%           - input from keyboard, latitude of the second position point, 
%             in radians or degrees/minutes/seconds
%           - input from keyboard, longitude of the second position point,
%             in radians or degrees/minutes/seconds
%  Outputs: - save on a specified file or display on screen, the input data and
%             the horizontal position error (range) between the two defined 
%             position points
%  Remark:  The computation accuracy is limited due to algorithm approximations;
%           it is recommended to be used when the range is less than 100 km.
%  External Matlab macros used:  hpe, tadmsrad, wgs84con
%  Last update:  09/05/00
%  Copyright (C) 1996-00 by LL Consulting. All Rights Reserved.

clear

wgs84con;
% global constants used:  a_smaxis, eccentr2

yes = 'y';
answer = yes;
disp(' ');

%  Save the generated data into a specified file or displayed on screen

answer1 = input('Do you want to save the generated data? (y/n)[n] ','s');
if  isempty(answer1)
   answer1 = 'n';
end
disp('  ');
if strcmp(answer1,yes) == 1
   f2 = input('Specify the output filename (with extension) -->  ','s');
   disp('  ');
else
   f2 = 1;     %   output to the screen
end

while (strcmp(answer,yes) == 1)     
   
   disp('Enter:  1 - latitude/longitude in radians ');
   disp('        2 - latitude/longitude in degrees/minutes/seconds ');
   option = input('Make selection --> ');
   disp('  ');

%  Read input data

   if  (option == 1)
      lat1 = input('Enter latitude for the first point, in radians   -->  ');
      lon1 = input('Enter longitude for the first point, in radians  -->  ');
      lat2 = input('Enter latitude for the second point, in radians  -->  ');
      lon2 = input('Enter longitude for the second point, in radians -->  ');
   elseif (option == 2)
      disp('Enter latitude for the first point:');
      lat1deg = input('   degrees (including sign) -->  ');
      lat1min = input('   minutes -->  ');
      lat1sec = input('   seconds -->  ');
      lat1 = tadmsrad(lat1deg,lat1min,lat1sec);
      disp('Enter longitude for the first point:');
      lon1deg = input('   degrees (including sign) -->  ');
      lon1min = input('   minutes -->  ');
      lon1sec = input('   seconds -->  ');
      lon1 = tadmsrad(lon1deg,lon1min,lon1sec);
      disp('Enter latitude for the second point:');
      lat2deg = input('   degrees (including sign) -->  ');
      lat2min = input('   minutes -->  ');
      lat2sec = input('   seconds -->  ');
      lat2 = tadmsrad(lat2deg,lat2min,lat2sec);
      disp('Enter longitude for the second point:');
      lon2deg = input('   degrees (including sign) -->  ');
      lon2min = input('   minutes -->  ');
      lon2sec = input('   seconds -->  ');
      lon2 = tadmsrad(lon2deg,lon2min,lon2sec);
   else
      disp('Error XHPE; option is not in range');
      disp('  ');
      return
   end   
   
%  Determine the horizontal position error (range)

   hperror = hpe(lat1,lon1,lat2,lon2);
   
%  Save the data
   
   fprintf(f2,'\n*******************************************************************\n\n');
   fprintf(f2,'*****Initial conditions:\n');
   if  (option == 1)
      fprintf(f2,'First point  - Latitude   (radians) =  %12.8f\n',lat1);
      fprintf(f2,'             - Longitude  (radians) =  %12.8f\n',lon1);
      fprintf(f2,'Second point - Latitude   (radians) =  %12.8f\n',lat2);
      fprintf(f2,'             - Longitude  (radians) =  %12.8f\n',lon2);
   elseif (option == 2)
      fprintf(f2,'First point  - Latitude   (deg/min/sec) =  %4.0f %3.0f %7.3f\n',...
                  lat1deg,lat1min,lat1sec);
      fprintf(f2,'             - Longitude  (deg/min/sec) =  %4.0f %3.0f %7.3f\n',...
                  lon1deg,lon1min,lon1sec);
      fprintf(f2,'Second point - Latitude   (deg/min/sec) =  %4.0f %3.0f %7.3f\n',...
                  lat2deg,lat2min,lat2sec);
      fprintf(f2,'             - Longitude  (deg/min/sec) =  %4.0f %3.0f %7.3f\n',...
                  lon2deg,lon2min,lon2sec);
   end
   fprintf(f2,'\nHorizontal position error (range) =  %12.4f  meters\n',hperror); 
   fprintf(f2,'\n*******************************************************************\n');

%  Select another computation, if desired

   disp('  ');
   answer = input('Do you want another computation? (y/n)[n] --> ','s');
   disp('  ');
   
end

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