function [lon,lat,delta]=pole(lon1,lat1,lon2,lat2)
%Определение положени? полюса вращени? и угол поворота
%дл? перемещени? одной точки в другую с известными
%сферическими координатами.
%Используютс? векторное и скал?рное произведени?
%
%Вызов:   [lon,lat,delta]=pole(lon1,lat1,lon2,lat2);
%
%   lon1,lat1 - долгота и широта точки 1;
%   lon2,lat2 - долгота и широта точки 2;
%
%     lon,lat - долгота и широта полюса вращени?;
%       delta - угол поворота.
%
%Все углы - в радианах !!!
%
%Вектор [1,2] показывает на точку, противоположную полюсу,
%определ?емую правилом правого винта при вращении от 1 к 2
%на угол DELTA.


%...........Пересчет в декартовые координаты на единичной сфере
[x1,y1,z1] = sph2car(lon1,lat1);
[x2,y2,z2] = sph2car(lon2,lat2);
%...........Угол DELTA между точками 1 и 2 (в рад.)
delta = acos( x1.*x2 + y1.*y2 + z1.*z2 );
%...........Вычисление компонент векторного произведени?
c = abs(sin(delta));   % величина вект. произведени? (нормир. множитель)
x3 = (y1.*z2 - z1.*y2)./c;
y3 = (z1.*x2 - x1.*z2)./c;
z3 = (x1.*y2 - y1.*x2)./c;
disp([x3 y3 z3])
%...........Перевод в сферическую СК
[lon,lat] = car2sph(x3,y3,z3);
lon = lon+pi;
lat = -lat;
