% угол между проекциями планет на эклиптику (их эклиптикальных долгот)

      clear all % очищаем среду МатЛаб от предыдущих, возможных, остатков в памяти после прерываний
      cspice_furnsh( 'd:\Documents and Settings\user\Мои документы\MATLAB\astro\mice\data\standard.tm' );
      MAXWIN  = 1000; nintvls= MAXWIN ; abcorr = 'LT+S'; 
      TFMT='DD.MM.YYYY ERA HR:MN:SC ::RND::MCAL::TDT'; TFMTH='HR:MN:SC.## ::RND::UTC' ;

      et1 = cspice_str2et({ 'JD 1721415.13', 'JD 1721418.63'}); 
      et2 = cspice_str2et({ 'JD 2141813.48', 'JD 2141818.96'}); 
%       cnfine1 = cspice_wninsd( et1(1), et1(2) ); cnfine2 = cspice_wninsd( et1(1), et1(2) );
      
      targ1   = 'Sun';
      targ2   = 'Mercury';
      observer = 'DYN'; frame = strcat(observer,'_TOPO'); 
      lat=44.6 ; lon=33.5333333; alt=0.0; % start topo place %obsrvr  = 'EARTH';
            
      s_step=(et1(2)-et1(1))/(cspice_spd/2);
      ii=1;
      for ii=1:s_step
      [pos1, lt1] = cspice_spkpos ( targ1, et1(1)+(ii-1)*(cspice_spd/2), 'ECLIPJ2000', abcorr, observer);  % позиция планеты
      [r, ra, lat1] = cspice_recrad(pos1); % картезианские координаты в радиальные небесной сферы
      lat1=cspice_convrt(lat1,'radians','degrees'); 
      %r = cspice_convrt( r, 'km' , 'au'); ra =cspice_convrt(ra,'radians','degrees'); 
      [pos2, lt2] = cspice_spkpos ( targ2, et1(1)+(ii-1)*(cspice_spd/2), 'ECLIPJ2000', abcorr, observer);  % позиция планеты
      [r, ra, lat2] = cspice_recrad(pos2); % картезианские координаты в радиальные небесной сферы
      lat2=cspice_convrt(lat2,'radians','degrees'); 
      %r = cspice_convrt( r, 'km' , 'au'); ra =cspice_convrt(ra,'radians','degrees');  abcorr, observer);  % позиция планеты
      sep1(ii) = lat1-lat2;
      end
      
      s_step=(et2(2)-et2(1))/(cspice_spd/2);
      ii=1;
      for ii=1:s_step
      [pos1, lt1] = cspice_spkpos ( targ1, et2(1)+(ii-1)*(cspice_spd/2), 'ECLIPJ2000', abcorr, observer);  % позиция планеты
      [r, ra, lat1] = cspice_recrad(pos1); % картезианские координаты в радиальные небесной сферы
      lat1=cspice_convrt(lat1,'radians','degrees'); 
      %r = cspice_convrt( r, 'km' , 'au'); ra =cspice_convrt(ra,'radians','degrees'); 
      [pos2, lt2] = cspice_spkpos ( targ2, et2(1)+(ii-1)*(cspice_spd/2), 'ECLIPJ2000', abcorr, observer);  % позиция планеты
      [r, ra, lat2] = cspice_recrad(pos2); % картезианские координаты в радиальные небесной сферы
      lat2=cspice_convrt(lat2,'radians','degrees'); 
      %r = cspice_convrt( r, 'km' , 'au'); ra =cspice_convrt(ra,'radians','degrees');  abcorr, observer);  % позиция планеты
      sep2(ii) = lat1-lat2;
      end
      
      p_name1=horzcat('Разница эклиптикальных долгот ', targ1, ' и ', targ2, ' с ',  ... 
          cspice_timout( et1(1), 'DD.MM.YYYY ERA HR:MN:SC ::RND::JCAL::UTC'), ' по ',  cspice_timout( et1(2), 'DD.MM.YYYY ERA HR:MN:SC ::RND::JCAL::UTC'), ' JCAL');
      p_name2=horzcat('Разница эклиптикальных долгот ', targ1, ' и ', targ2, ' с ',  ... 
          cspice_timout( et2(1), 'DD.MM.YYYY ERA HR:MN:SC ::RND::JCAL::UTC'), ' по ',  cspice_timout( et2(2), 'DD.MM.YYYY ERA HR:MN:SC ::RND::JCAL::UTC'), ' JCAL');
      
      p=figure('Name','Видимый угол между планетами (центрами) на определенный период','NumberTitle','off');      
      subplot(2,1,1); p=plot(sep1(:));
      set(p,'Color','red','LineWidth',2)
      title(p_name1)
      subplot(2,1,2); p=plot(sep2(:));
      title(p_name2)
      set(p,'Color','green','LineWidth',2)

        cspice_kclear