function [lon_new,lat_new,R]=rot(lon_old,lat_old,lon_p,lat_p,omega)
%Преобразование сферических координат (на сфере с радиусом 1)
%при повороте относительно полюса с координатами (lon_p,lat_p)
%на угол OMEGA
%
%Вызов:  [lon_new,lat_new,R]=rot(lon_old,lat_old,lon_p,lat_p,omega);
%
%    lon_old,lat_old - старые координаты (долгота и широта);
%        lon_p,lat_p - координаты полюса вращени? (долгота и широта);
%              omega - угол поворота;
%
%    lon_new,lat_new - новые координаты (долгота и широта);
%                  R - матрица вращени? (3*3);
%
%Все углы - в радианах !!!
%

% Хомутов С.Ю. 30.09.96  Новосибирск

% Кокс А., Харт Р. Тектоника плит. - М: Мир, 1989. - с.265-267.

% Вычисление элементов матрицы поворота в декартовой СК
% ....... положение полюса в декартовых координатах
[x,y,z] = sph2car(lon_p,lat_p);
tmp = 1-cos(omega); s=sin(omega); c=cos(omega);
R = zeros(3,3);
R(1,1) = tmp*(x*x) +   c;
R(1,2) = tmp*(x*y) - z*s;
R(1,3) = tmp*(x*z) + y*s;
R(2,1) = tmp*(y*x) + z*s;
R(2,2) = tmp*(y*y) +   c;
R(2,3) = tmp*(y*z) - x*s;
R(3,1) = tmp*(z*x) - y*s;
R(3,2) = tmp*(z*y) + x*s;
R(3,3) = tmp*(z*z) +   c;

% ........ старые координаты в декартовой системе
m = length(lon_old); % размер векторов старых координат
tmp=zeros(m,3);
[tmp(:,1),tmp(:,2),tmp(:,3)] = sph2car(lon_old,lat_old);
% ........ новые координаты в декартовой системе
tmp2=zeros(3,m);
for k=1:m
   tmp2(:,k) = R*tmp(k,:)';
end
tmp2=tmp2';
% ........ новые координаты в сферической системе
[lon_new,lat_new,r] = car2sph(tmp2(:,1),tmp2(:,2),tmp2(:,3));
