xpelazs.m

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

%                                xpelazs.m
%  Scope:  This MATLAB program plots graphs related to the number of visible
%          satellites, elevation and azimuth angles, and azimuth-elevation 
%          configuration for a specified satellite and selected user; the input
%          file can be generated by the program xelaza.m.
%  Usage:  xpelazs
%  Inputs: - name of the input file; the default is xpgelaz1.dat (i.e. 
%            xelaz1.out generated by the program xelaza.m).       
%            The file contains for each time step and for each visible
%            satellite the following:
%            1) time step number, 2) user's position index, 3) number of
%            visible satellites, 4) id of visible satellite, 5) elevation
%            angle, 6) azimuth angle 
%          - selected time steps in seconds (see also program xelaza.m) - from
%            keyboard; the default is 120 seconds 
%          - selected satellite index to be analyzed - from keyboard
%          - selected user's location index to be analyzed - from keyboard 
%          - selected graph type (from the menu) for a selected satellite and
%            specified user, namely: 
%            1) number of visible satellites versus Time,
%            2) elevation angle versus Time,
%            3) azimuth angle versus time,
%            4) elevation and azimuth angles versus Time,
%            5) elevation versus azimuth (polar plot)
%  Outputs:  - selected graphs 
%  Remarks:  For a specific interval, manual scaling and statistics are
%            available by using several other programs, for example xyp1s.
%  External Matlab macros used:  selectd      
%  Last update:  01/04/00
%  Copyright (C) 1996-00 by LL Consulting. All Rights Reserved.

clear
yes = 'y';
answer1 = 'y';
close all

rad_deg = 180./pi;

%  Read the input data 

disp('  ');
disp('Enter the name of the input file ');
answer = input('Use the default --> xelaz1.out (y/n)[n]?  ','s');
if (strcmp(answer,yes) == 1)
   fname = 'xelaz1.out';
else
   fname = input('Specify the input filename --> ','s');
end
hh = load(fname);
[nrow,ncol] = size(hh);

disp('  ');
disp('Enter time step - the default value is 120 seconds ');
answer = input('Do you want to use the default value (y/n)[n]?  ','s');
if (strcmp(answer,yes) == 1)
   t_incr = 120.;         %  120 seconds - time increment
else
   disp('  ');
   t_incr = input('Specify the time step in seconds -->  ');
end
disp('  ');

zu = hh(:,2);
zsv = hh(:,4);
[u_sel,u_nr] = selectd(zu);
[sv_sel,sv_nr] = selectd(zsv);
disp('User"s indices = ');
for ku = 1:u_nr
   fprintf('  %5.0f ',u_sel(ku));
end
disp('  ');
disp('Satellite"s indices = ');
for ks = 1:sv_nr
   fprintf('%3.0f',sv_sel(ks));
end
disp('  ');

%  Initialization

mtitle  = 'Choose a plot for the selected user';
select1 = 'Number of visible satellites versus Time';
select2 = 'Elevation angle versus Time for a selected satellite';
select3 = 'Azimuth angle versus Time for a selected satellite';
select4 = 'Elevation and azimuth angles versus Time for a selected satellite';
select5 = 'Elevation versus azimuth plot for a specified satellite';
select6 = 'Exit program';

%  Enter the main loop

while (strcmp(answer1,yes) == 1)
   clear table
   disp('  ');
   u_index = input('Specify the user"s location index to be analyzed --> ');
   sv_index = input('Specify the satellite number to be analyzed --> ');
   kk = 1;
   kkk = 1;
   for i = 1:nrow
      if hh(i,2) == u_index           %  only records for selected user
         if hh(i,1) == kkk
            y1(kkk) = hh(i,3);
            kkk = kkk + 1;
         end
         if hh(i,4) == sv_index       %  only records for selected satellite
            if kk == 1
               table = hh(i,:);
            else
               table = [table' hh(i,:)']';
            end
            kk = kk + 1;
         end
      end
   end
   if (kk == 1) | (kkk == 1)
      disp('  ');
      disp('Error XPELAZS; out of range selection for user"s index or SV #');
      disp('  ');
      return
   end
   nsteps = kkk - 1;            % number of time steps
   nrecords = kk - 1;           % number of records
   x = 1:1:nsteps;              % x-axis            
   k = 1;
   for i = 1:nsteps             % extract data for elevation and azimuth angles
      if k <= nrecords
         if i == table(k,1)
            y2(i) = table(k,5) * rad_deg;       % elevation angle in degrees
            y3(i) = table(k,6) * rad_deg;       % azimuth angle in degrees
            if y3(i) < 0.
               y3(i) = y3(i) + 360.;
            end
            k = k + 1;
         else                   % insert value 0. when not available
            y2(i) = 0.;            
            y3(i) = 0.;
         end
      else                      % insert value 0. when not available at the end
         y2(i) = 0.;               
         y3(i) = 0.;
      end
   end
   answer2 = 'y';
   while (strcmp(answer2,yes) == 1)
      comp = menu(mtitle,select1,select2,select3,select4,select5,select6);

%     Execute the selected plot  

      atemp = num2str(sv_index);
      aa = ['SV #  ' atemp];
      btemp = num2str(u_index);
      bb = ['User # ' btemp];
      ff = 'Time sequence  (Time step = ';
      fftemp = num2str(t_incr);
      ff = [ff fftemp ' seconds)'];

      if comp == 1
         cc = ['Number of visible satellites for ' bb];
         ee = 'Number of visible satellites';
         figure(1)
         plot(x,y1),...
         title(cc), ylabel(ee), xlabel(ff), grid
      elseif comp == 2
         cc = ['Elevation angle for ' aa '  and ' bb];
         ee = 'Elevation angle, in degrees';
         figure(2)
         plot(x,y2),...
         title(cc), ylabel(ee), xlabel(ff), grid
      elseif comp == 3
         cc = ['Azimuth angle for ' aa '  and ' bb];
         ee = 'Azimuth angle, in degrees';
         figure(3)
         plot(x,y3),...
         title(cc), ylabel(ee), xlabel(ff), grid
      elseif comp == 4
         cc = ['Elevation and Azimuth angles for ' aa ' and ' bb];
         ee = 'Elevation and azimuth angles, in degrees';
         figure(4)
         plot(x,y2,x,y3),...
         title(cc), ylabel(ee), xlabel(ff), grid
      elseif comp == 5
         cc = ['Azimuth-Elevation plot for ' aa '  and ' bb];
         h0 = [0  0]; 
         h9 = [-0.9  0.9];
         azrad = y3*pi/180;
         dist = (90.*ones(1,nsteps) - y2)/100.;
         for i = 1:nsteps
            svx(i) = dist(i)*cos(azrad(i));
            svy(i) = dist(i)*sin(azrad(i));
         end
         an1 = [0:120]*pi/60;
         an2 = [0:100]*pi/50;
         an3 = [0:80]*pi/40;
         an4 = [0:60]*pi/30;
         an5 = [0:32]*pi/16;
         an6 = [0:24]*pi/12;
         figure(5)
         plot(0.9*sin(an1),0.9*cos(an1),'.',0.8*sin(an1),0.8*cos(an1),'.',...
              0.7*sin(an2),0.7*cos(an2),'.',0.6*sin(an2),0.6*cos(an2),'.',...
              0.5*sin(an3),0.5*cos(an3),'.',0.4*sin(an3),0.4*cos(an3),'.',...
              0.3*sin(an4),0.3*cos(an4),'.',0.2*sin(an5),0.2*cos(an5),'.',...
              0.1*sin(an6),0.1*cos(an6),'.');
         hold on
         plot(h9,h0,'-k',h0,h9,'-k',svx,svy,'.k')
         title(cc)
         axis('square');
         axis('off');
         text(svx(1),svy(1),'X');
         text(1.,0,'East')
         text(-1.25,0,'West')
         text(-0.12,1,'North')
         text(-0.12,-1,'South')
         text(0.13,-0.22,'60')
         text(0.31,-0.46,'30')
         text(0.55,-0.68,'0')
     elseif comp == 6
         clear xpgelaz;
         return
     else
         clear xpgelaz;
         error('Error XPELAS;  check the index of the selection');
         return
     end

%     Check if another plot is needed

      disp('  ');
      nswer1 = input('Do you want to make another plot? (y/n)[y] ','s');
      if  isempty(answer1)
         answer1 = 'y';
      end
      disp('  ');
      if (strcmp(answer1,yes) == 1)
         disp('Do you want to use the same satellite and user? ');
         answer2 = input('Make the selection (y/n)[y] -->  ','s');
         if  isempty(answer2)
            answer2 = 'y';
         end
         if ~strcmp(answer2,yes)
            close all
         end   
      else
         answer2 = 'n';
         close all
      end
   
   end

end
disp('  ');
disp('End of the program  XPELAZS ');
disp('  ');