Contents
NAIF/MatLab script
Author: A. Lucas lucas@ipgp.fr version: Jul 2008 Description: This script reads data from NAIF kernels files so as to get spacecraft and camera model informations (position, velocity...) Requierments: Matlab, MICE toolkit, NAIF kernels files Tested under Linux/Solaris/Mac OS X
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ close all hidden clear all %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ %
Initialisation
%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ addpath('/net/crunch/data1/Lucas/PGMw/SOFTS/NAIF_MICE/mice/src/mice') addpath('/net/crunch/data1/Lucas/PGMw/SOFTS/NAIF_MICE/mice/src/cspice') addpath('/net/crunch/data1/Lucas/PGMw/SOFTS/NAIF_MICE/mice/lib/') %~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
User options
Define the number of divisions of the time interval and the time interval.
STEP = 1*12; nn=round(STEP/2); utc = strvcat( '00:02:59,Dec 8, 2006', '00:03:02, Dec 8, 2006' ); ramars = 3396.19 ; rbmars= ramars; rcmars = 3376.2 ; X= imread('mars_mdim.png'); X=X(end:-1:1,:,:); % Define some stuffs targ='MRO'; ref='IAU_MARS'; abcorr='NONE'; obs='MARS'; cam_ref='MRO_HIRISE_LOOK_DIRECTION'; method='near point'; txt = sprintf( 'NAIF/MAtLab script. A. Lucas'); disp( txt ) disp( ' ' ) txt = sprintf( 'This script will output and display :'); disp( txt ) txt = sprintf( 'The position of : %s', targ); disp( txt ) txt = sprintf( 'As observed from : %s', obs ); disp( txt ) txt = sprintf( 'In reference frame : %s', ref ); disp( txt ) txt = sprintf( 'and the transpose matrix of : %s', cam_ref ); disp( txt ) disp( ' ' )
NAIF/MAtLab script. A. Lucas This script will output and display : The position of : MRO As observed from : MARS In reference frame : IAU_MARS and the transpose matrix of : MRO_HIRISE_LOOK_DIRECTION
Load the needed kernels. Use a meta kernel "standard.ker" to load the kernels
%kernels path: k_path='/net/crunch/data1/Lucas/PGMw/SOFTS/ISIS/data'; %base suite: lsk_base=[k_path,'/base/kernels/lsk/naif0008.tls']; spk_base=[k_path,'/base/kernels/spk/de405.bsp']; pck_base=[k_path,'/base/kernels/pck/pck00008.tpc']; %Spacecraft suite: ck_sc=[k_path,'/mro/kernels/ck/mro_sc_psp_061205_061211.bc']; spk_sc=[k_path,'/mro/kernels/spk/mro_psp1.bsp']; sclk_sc=[k_path,'/mro/kernels/sclk/MRO_SCLKSCET.00026.65536.tsc']; fk_sc=[k_path,'/mro/kernels/fk/mro_v11.tf']; %loading using cspice_furnsh command cspice_furnsh(lsk_base) cspice_furnsh(spk_base) cspice_furnsh(pck_base) cspice_furnsh(spk_sc) cspice_furnsh(ck_sc) cspice_furnsh(sclk_sc) cspice_furnsh(fk_sc);
Processing
et = cspice_str2et( utc );
times = (0:STEP-1) * ( et(2) - et(1) )/STEP + et(1);
time = times(nn);
[pos,ltime ]=cspice_spkpos( targ, times, ref, abcorr, obs);
[state, lt] = cspice_spkezr(targ, times, ref, abcorr, obs);
[spoint, alt] = cspice_subpt( method, obs, time, abcorr, targ);
[radius, longitude, latitude] = cspice_reclat(spoint);
longitude = longitude * cspice_dpr;
latitude = latitude * cspice_dpr;
r2d = cspice_dpr;
mat = cspice_pxform( cam_ref,ref, time );
pole = mat(:,3,:);
[radius, ra, dec] = cspice_recrad( pole );
ra = ra * r2d;
dec = dec * r2d;
output = [ time; ra; dec ];
%fprintf( '%24.8f %16.8f %16.8f\n', output );
x = pos(1,:);
y = pos(2,:);
z = pos(3,:);
vx = state(4,nn);
vy = state(5,nn);
vz = state(6,nn);
Results
txt = sprintf( '1/* The position of: %s', targ); disp( txt ) txt = sprintf( ' X position (km): %s', num2str(x(nn))); disp( txt ) txt = sprintf( ' Y position (km): %s', num2str(y(nn))); disp( txt ) txt = sprintf( ' Z position (km): %s', num2str(z(nn))); disp( txt ) disp( ' ' ) txt = sprintf( ' Latitude: %s', num2str(latitude)); disp( txt ) txt = sprintf( ' Longitude: %s', num2str(longitude)); disp( txt ) txt = sprintf( ' Altitude (km): %s', num2str(alt)); disp( txt ) disp( ' ' ) txt = sprintf( '2/* The velocity of: %s', targ); disp( txt ) txt = sprintf( ' X velocity (km/s): %s', num2str(vx)); disp( txt ) txt = sprintf( ' Y velocity (km/s): %s', num2str(vy)); disp( txt ) txt = sprintf( ' Z velocity (km/s): %s', num2str(vz)); disp( txt ) disp( ' ' ) txt = sprintf( '3/* Transpose Matrix of the camera to the frame:'); disp( txt ) mat txt = sprintf( ' Azimuth and declinaison:'); disp( txt ) output = [ra; dec ]; fprintf( '%16.8f %16.8f\n', output );
1/* The position of: MRO
X position (km): -346.9951
Y position (km): 2832.5302
Z position (km): -2261.2288
Latitude: -38.3699
Longitude: 96.9841
Altitude (km): 252.5456
2/* The velocity of: MRO
X velocity (km/s): 0.14314
Y velocity (km/s): 2.1856
Z velocity (km/s): 2.6726
3/* Transpose Matrix of the camera to the frame:
mat =
0.0416 -0.9911 0.1261
0.6297 -0.0720 -0.7735
0.7757 0.1116 0.6212
Azimuth and declinaison:
279.26278211 38.40153191
Displaying Results
Plot the resulting trajectory if fig=1;
fig=1; if (fig==1) figure1=figure(1); clf; set(gcf, 'color', 'white') set(0,'defaultaxesfontname','Helvetica','defaultaxesfontsize',2); subplot1 = subplot(1,2,1,'ZColor',[0.2471 0.2471 0.2471],... 'YColor',[0.2471 0.2471 0.2471],... 'XColor',[0.2471 0.2471 0.2471],... 'PlotBoxAspectRatio',[1 1 1],... 'DataAspectRatio',[1 1 1],... 'Color',[0.9412 0.9412 0.9412],... 'CameraViewAngle',9.108); box('on'); grid('on'); hold('all'); axis([-4000 4000 -4000 4000 -4000 4000]); x = pos(1,:); y = pos(2,:); z = pos(3,:); scatter3(x,y,z,'+k') %plot3(x,y,z) xlabel('Distance in X-axis (km)') ylabel('Distance in Y-Axis (km)') zlabel('Distance in Z-Axis (km)') title(num2str(utc)) hold on [mx,my,mz]=ellipsoid(0,0,0,ramars,rbmars,rcmars,50); %mesh(mx,my,mz) %shading interp surf(mx,my,mz,'facecolor','texturemap',... 'cdata',X,'edgecolor','none'); %axes('Color',[0.9412 0.9412 0.9412]); axis equal axis vis3d % Uncomment the following line to preserve the X-limits of the axes % xlim([-3396 3396]); % Uncomment the following line to preserve the Y-limits of the axes % ylim([-3389 3389]); % Uncomment the following line to preserve the Z-limits of the axes % zlim([-3541 3575]); view([-128.5 28]); box('on'); grid('on'); % % Create figure % figure1 = figure('XVisual',... % '0x23 (TrueColor, depth 24, RGB mask 0xff0000 0xff00 0x00ff)',... % 'PaperType','a4letter',... % 'PaperSize',[20.98 29.68],... % 'Color',[1 1 1]); % Create subplot subplot2 = subplot(1,2,2); %subplot(1,2,2) vx = state(4,nn); vy = state(5,nn); vz = state(6,nn); %scatter3(vx,vy,vz) quiver3(vx,vy,vz,15) hold on %plot3(x(nn),y(nn),z(nn),'+k') xlabel('Velocity in X-axis (km/s)') ylabel('Velocity in Y-Axis (km/s)') zlabel('Velocity in Z-Axis (km/s)') [mx,my,mz]=ellipsoid(0,0,0,ramars,rbmars,rcmars,50); %mesh(mx,my,mz) %shading interp hold on %surf(mx,my,mz,'facecolor','texturemap',... % 'cdata',X,'edgecolor','none'); axis vis3d aminx=min(x); aminy=min(y); aminz=min(z); amaxx=max(x); amaxy=max(y); amaxz=max(z); %axis([aminx 2000 aminy 3000 aminz 2000]) title(num2str(utc)) xt = 0; yt = 0; zt=0; % Create text annotation(figure1,'textbox',... 'LineStyle','none',... 'String',{['Velocity (km/s): V_x =', num2str(vx),', V_y = ', num2str(vy),... ', V_z = ', num2str(vz)]}',... 'Position',[0 0 0.1 0.1]); % Create text annotation(figure1,'textbox',... 'String',{['Position of spacecraft (km): X =', num2str(x(nn)),... ', Y = ', num2str(y(nn)),', Z = ', num2str(z(nn))]}',... 'LineStyle','none',... 'Position',[0 0.08 .1 0.1]); % Create text annotation(figure1,'textbox',... 'LineStyle','none',... 'String',{['Longitude = ', num2str(longitude),'°E, Latitude = ',... num2str(latitude),'°N, Altitude = ', num2str(alt),'km']},... 'Position',[0 0.05 0.1 0.1]); % Create textbox annotation(figure1,'textbox',... 'String',{'MRO Spacecraft position and velocity for PSP\_001714\_1415'},... 'FontWeight','bold',... 'FontSize',14,... 'FitHeightToText','off',... 'EdgeColor',[1 1 1],... 'Position',[0.2012 0.9405 0.546 0.04415]); % Create textbox annotation(figure1,'textbox','String',{'A. Lucas <lucas@ipgp.fr>'},... 'FontSize',6,... 'FitHeightToText','off',... 'LineStyle','none',... 'Position',[0.8065 0.03113 0.3254 0.02689]); % Output into eps & png files %print 'spk.mro.psp_1714_1415.eps' -depsc %print 'spk.mro.psp_mro-1714-1415.png' -dpng -adobecset end cspice_kclear
