diff --git a/.DS_Store b/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..f9f4fed87a611331bed99dcf016a61ae889c7c07
Binary files /dev/null and b/.DS_Store differ
diff --git a/MMX_Fcn_CovarianceAnalyses/Cov_Dynamics_Good.m b/MMX_Fcn_CovarianceAnalyses/Cov_Dynamics_Good.m
index c85e66e358ca03eb73e021aeda18e314efbfdca4..8397310095a54ee298e06b55ded41de4561eeb3d 100644
--- a/MMX_Fcn_CovarianceAnalyses/Cov_Dynamics_Good.m
+++ b/MMX_Fcn_CovarianceAnalyses/Cov_Dynamics_Good.m
@@ -69,7 +69,7 @@ function dx = Cov_Dynamics_Good(~,x,par,~)
IPhz_bar = x(17+i*size_St);
% Position of MMX wrt Phobos
- rPh = par.perifocal2MARSIAU*[RPh*cos(theta); RPh*sin(theta); 0];
+ rPh = par.perifocal2MARSIAU*RPh*[cos(theta); sin(theta); 0];
rsb = r - rPh;
% Rotation matrix from orbit's to Phobos fixed coordinates frame
diff --git a/MMX_Fcn_CovarianceAnalyses/Dynamics_MPHandMMX_Inertia.m b/MMX_Fcn_CovarianceAnalyses/Dynamics_MPHandMMX_Inertia.m
index f7ea0e9ffca908473b01d6dcbcf653430e0b9482..3b14c94daa63aea4af137d99e8b200694d160aa2 100644
--- a/MMX_Fcn_CovarianceAnalyses/Dynamics_MPHandMMX_Inertia.m
+++ b/MMX_Fcn_CovarianceAnalyses/Dynamics_MPHandMMX_Inertia.m
@@ -105,9 +105,10 @@ function dx = Dynamics_MPHandMMX_Inertia(~,x,par,~)
C_20_Ph = -.5*(2*IPhz_bar - IPhx_bar - IPhy_bar);
par.SH.Clm = [1, 0, 0; 0, 0, 0; C_20_Ph, 0, C_22_Ph]./par.SH.Norm(1:3,1:3);
- [gb,~,~,~] = SHGM_Full(BO*ON*rsb, par.SH);
+% [gb,~,~,~] = SHGM_Full(BO*ON*rsb, par.SH);
+ [gb,~,G,~] = CGM(BO*ON*rsb, par.SH);
gPh = NO*OB*gb;
-
+
% Trascinamento
G = par.G;
@@ -143,7 +144,7 @@ function dx = Dynamics_MPHandMMX_Inertia(~,x,par,~)
dbetady = -rPh(3)*rPh(2)/((1 - rPh(3)^2/RPh^2)^(1/2)*RPh^3);
dbetadz = (1/RPh - rPh(3)^2/RPh^3)/(1 - rPh(3)^2/RPh^2)^(1/2);
- dVdtheta = 3*G/(2*norm(rPh(1:3))^5)*rPh(1:3)'*(M2*(dA1dtheta'*I1*A1 + A1'*I1*dA1dtheta) + M1*(dA2dtheta'*I2*A2 + A2'*I2*dA2dtheta))*rPh(1:3);
+ dVdtheta = 3*G/(2*norm(rPh(1:3))^5)*rPh'*(M2*(dA1dtheta'*I1*A1 + A1'*I1*dA1dtheta) + M1*(dA2dtheta'*I2*A2 + A2'*I2*dA2dtheta))*rPh;
% Gravitational acceleration in x,y,z
g_M_on_PH = -[dRPhdx, dthetadx, dbetadx; dRPhdy, dthetady, dbetady; dRPhdz, dthetadz, dbetadz]*[dVdRPh; dVdtheta; 0];
diff --git a/MMX_Fcn_CovarianceAnalyses/NewCov_ResultsPlot.m b/MMX_Fcn_CovarianceAnalyses/NewCov_ResultsPlot.m
index 7ccaa4050b8b0993192fd7604a02bb3f6edff680..1efbc7a75d8368a2704f833239f3c8710cb41670 100644
--- a/MMX_Fcn_CovarianceAnalyses/NewCov_ResultsPlot.m
+++ b/MMX_Fcn_CovarianceAnalyses/NewCov_ResultsPlot.m
@@ -255,28 +255,41 @@ function NewCov_ResultsPlot(Est, YObs_Full, pars, units)
% corr_label(end+1) = {'$\dot{\rho}$'};
% corr_label(end+1) = {'$LIDAR_b$'};
-% figure()
-% subplot(1,pars.nBias,1)
-% semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+1,place0+1,1:end-1)))),'LineWidth',1,...
-% 'DisplayName',('$\rho$'),'Color','b');
-% grid on;
-% hold on;
-% xlabel('time $[hour]$')
-% ylabel('$[km]$')
-% subplot(1,pars.nBias,2)
-% semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+2,place0+2,1:end-1)))),'LineWidth',1,...
-% 'DisplayName',('$\dot{\rho}$'),'Color','r');
-% grid on;
-% hold on;
-% xlabel('time $[hour]$')
-% ylabel('$[km/s]$')
-% subplot(1,pars.nBias,3)
-% semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+3,place0+3,1:end-1)))),'LineWidth',1,...
-% 'DisplayName',('$LIDAR_b$'),'Color','g');
-% grid on;
-% hold on;
-% xlabel('time $[hour]$')
-% ylabel('$[km]$')
+ figure()
+ subplot(1,pars.nBias,1)
+ semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+1,place0+1,1:end-1)))),'LineWidth',1,...
+ 'DisplayName',('$\rho$'),'Color','b');
+ grid on;
+ hold on;
+ xlabel('time $[hour]$')
+ ylabel('$[km]$')
+ subplot(1,pars.nBias,2)
+ semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+2,place0+2,1:end-1)))),'LineWidth',1,...
+ 'DisplayName',('$\dot{\rho}$'),'Color','r');
+ grid on;
+ hold on;
+ xlabel('time $[hour]$')
+ ylabel('$[km/s]$')
+ subplot(1,pars.nBias,3)
+ semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+3,place0+3,1:end-1)))),'LineWidth',1,...
+ 'DisplayName',('$LIDAR_b$'),'Color','g');
+ grid on;
+ hold on;
+ xlabel('time $[hour]$')
+ ylabel('$[km]$')
+
+ subplot(1,pars.nBias,4)
+ semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+4,place0+4,1:end-1)))),'LineWidth',1,...
+ 'DisplayName',('$cam_x$'),'Color','m');
+ grid on;
+ hold on;
+ xlabel('time $[hour]$')
+ subplot(1,pars.nBias,5)
+ semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(Est.P_t(place0+5,place0+5,1:end-1)))),'LineWidth',1,...
+ 'DisplayName',('$cam_y$'),'Color','b');
+ grid on;
+ hold on;
+ xlabel('time $[hour]$')
%
diff --git a/MMX_Fcn_CovarianceAnalyses/New_MMX_CovarianceAnalysisParameters.m b/MMX_Fcn_CovarianceAnalyses/New_MMX_CovarianceAnalysisParameters.m
index 711573167df0486fe8f2e481d69454d1fb05b7ec..f28b1a462131f8856b35ea31ba34e041c77587a6 100644
--- a/MMX_Fcn_CovarianceAnalyses/New_MMX_CovarianceAnalysisParameters.m
+++ b/MMX_Fcn_CovarianceAnalyses/New_MMX_CovarianceAnalysisParameters.m
@@ -60,12 +60,8 @@ pars.R(9,9) = pars.ObsNoise.camera^2;
% Seconds between observation
pars.interval_Range = 3600; % s
pars.interval_Range_rate = 600; % s
-pars.interval_lidar = 100; % s
+pars.interval_lidar = 600; % s
pars.interval_camera = 3600; % s
-% pars.interval_Range = 100e100; % s
-% pars.interval_Range_rate = 100e100; % s
-% pars.interval_lidar = 100e100; % s
-% pars.interval_camera = 100e100; % s
pars.cutoffLIDAR = 200; % km
@@ -75,28 +71,36 @@ pars.elevation = 80; % deg
pars.I2 = [pars.IPhx_bar; pars.IPhy_bar; pars.IPhz_bar];
% Measurements biases
-bias = [0; 0; 0]...
- ./[units.lsf; units.vsf; units.lsf];
+% bias = [0; 0; 0]...
+% ./[units.lsf; units.vsf; units.lsf];
+% pars.bias = bias;
+
+bias = [0; 0; 0; 0; 0]...
+ ./[units.lsf; units.vsf; units.lsf; 1; 1];
pars.bias = bias;
% Unit vector and parameters number useful for later use
pars.nCoeff = size(pars.I2,1);
pars.nBias = size(bias,1);
% pars.nBias = 0;
+% units.Dim = [units.sfVec; units.lsf; units.vsf; 1; units.tsf;...
+% 1; units.tsf; 1; units.tsf; ones(pars.nCoeff,1); units.lsf; units.vsf;...
+% units.lsf];
units.Dim = [units.sfVec; units.lsf; units.vsf; 1; units.tsf;...
1; units.tsf; 1; units.tsf; ones(pars.nCoeff,1); units.lsf; units.vsf;...
- units.lsf];
-% units.Dim = [units.sfVec; units.lsf; units.vsf; 1; units.tsf;...
-% 1; units.tsf; 1; units.tsf; ones(pars.nCoeff,1)];
+ units.lsf; 1; 1];
units.CovDim = (repmat(units.Dim, 1, length(units.Dim)).*...
repmat(units.Dim', length(units.Dim), 1));
% A-priori covariance
k = 1e-1;
+% pars.P0 = diag(([.3*ones(3,1); .3e-3*ones(3,1);...
+% 1e-1; 1e-4; 1e-5; 1e-7; 3e-1; 1e-3; 3e-1; 1e-3; 1e-2*ones(pars.nCoeff,1);...
+% 1; 1e-3; 1]./units.Dim).^2);
pars.P0 = diag(([.3*ones(3,1); .3e-3*ones(3,1);...
1e-1; 1e-4; 1e-5; 1e-7; 3e-1; 1e-3; 3e-1; 1e-3; 1e-2*ones(pars.nCoeff,1);...
- 1; 1e-3; 1]./units.Dim).^2);
+ 1; 1e-3; 1; 1e-1; 1e-1]./units.Dim).^2);
% pars.P0 = diag(([.3*ones(3,1); .3e-3*ones(3,1);...
% 1e-1; 1e-4; 1e-6; 1e-8; 3e-1; 1e-3; 3e-1; 1e-3; 1e-2*ones(pars.nCoeff,1)]./units.Dim).^2);
diff --git a/MMX_Fcn_CovarianceAnalyses/New_Observables_model.m b/MMX_Fcn_CovarianceAnalyses/New_Observables_model.m
index 8b35543c7a900c901f727edc542fabb986f70262..3f3ced3acd81acc33071ad9064ca478099c1df04 100644
--- a/MMX_Fcn_CovarianceAnalyses/New_Observables_model.m
+++ b/MMX_Fcn_CovarianceAnalyses/New_Observables_model.m
@@ -64,6 +64,8 @@ function G = New_Observables_model(et,Obs,X,par,units)
bias = X(16:end);
case 20
bias = X(18:end);
+ case 22
+ bias = X(18:end);
end
% Definition of Mars position WRT Earth J2000
@@ -176,7 +178,13 @@ function G = New_Observables_model(et,Obs,X,par,units)
k = cross(I, j)/norm(cross(I, j));
T = [I, j, k]';
- camera = ([0 1 0; 0 0 1]*T*I)';
+
+ switch size(X,1)
+ case 20
+ camera = ([0 1 0; 0 0 1]*T*I)';
+ case 22
+ camera = ([0 1 0; 0 0 1]*T*I)' + [bias(4), bias(5)];
+ end
end
diff --git a/MMX_Fcn_CovarianceAnalyses/Phobos_States_NewModel.m b/MMX_Fcn_CovarianceAnalyses/Phobos_States_NewModel.m
index 4769fc3cfa8efdbdcc7640b03abda09f1a7785f1..dc21e010bcee4f66d49219986d2e3461224ddd0f 100644
--- a/MMX_Fcn_CovarianceAnalyses/Phobos_States_NewModel.m
+++ b/MMX_Fcn_CovarianceAnalyses/Phobos_States_NewModel.m
@@ -13,8 +13,8 @@ function [Ph, par] = Phobos_States_NewModel(data,par)
ey = cross(v,h)/par.mu - r/norm(r);
ey = ey/norm(ey);
Y = cross(eh,ey)/norm(cross(eh,ey));
- par.MARSIAU2perifocal = [ey,Y,eh]';
- par.perifocal2MARSIAU = par.MARSIAU2perifocal';
+ par.perifocal2MARSIAU = [ey,Y,eh];
+ par.MARSIAU2perifocal = par.perifocal2MARSIAU';
Ph_0 = par.MARSIAU2perifocal*Phobos_0(1:3);
theta = atan2(Ph_0(2),Ph_0(1));
diff --git a/MMX_Fcn_CovarianceAnalyses/UKF_CovarianceAnalysis.m b/MMX_Fcn_CovarianceAnalyses/UKF_CovarianceAnalysis.m
index 5865ebb38ef7d2b7da86e5292e1b4fbe6cbd528a..4683531b2f64b8e392c66d36ddefb79744fc146b 100644
--- a/MMX_Fcn_CovarianceAnalyses/UKF_CovarianceAnalysis.m
+++ b/MMX_Fcn_CovarianceAnalyses/UKF_CovarianceAnalysis.m
@@ -97,6 +97,8 @@ function [Est] = UKF_CovarianceAnalysis(Est0,f,fsigma,G,O,Measures,par,units)
St0 = par.X0_reference;
case 20
St0 = Xold(1:17);
+ case 22
+ St0 = Xold(1:17);
end
tspan = (0:min_interval:tobs(end))*units.tsf;
@@ -133,6 +135,8 @@ function [Est] = UKF_CovarianceAnalysis(Est0,f,fsigma,G,O,Measures,par,units)
% For-loop on the observations
for i = 1:m
+ fprintf('\nProcesso osservazione n.: %d', i)
+
% Look for the observation to process
if tobs(i) == told
idx = 1;
@@ -146,6 +150,8 @@ function [Est] = UKF_CovarianceAnalysis(Est0,f,fsigma,G,O,Measures,par,units)
par.PhiM = X(idx,11);
case 20
Xmean_bar = [X(idx,:)'; 0; 0; 0];
+ case 22
+ Xmean_bar = [X(idx,:)'; 0; 0; 0; 0; 0];
end
elseif tobs(i) > told
idx = find(round(t/units.tsf)==tobs(i));
@@ -158,6 +164,8 @@ function [Est] = UKF_CovarianceAnalysis(Est0,f,fsigma,G,O,Measures,par,units)
Xmean_bar = [X(idx,1:10)'; X(idx,13:end)'; 0; 0; 0];
case 20
Xmean_bar = [X(idx,:)'; 0; 0; 0];
+ case 22
+ Xmean_bar = [X(idx,:)'; 0; 0; 0; 0; 0];
end
end
@@ -208,7 +216,6 @@ function [Est] = UKF_CovarianceAnalysis(Est0,f,fsigma,G,O,Measures,par,units)
Xbar = reshape(X_sigma(end,1:n*(2*n))',[n,(2*n)]);
end
- i
% Process noise if needed
@@ -231,6 +238,10 @@ function [Est] = UKF_CovarianceAnalysis(Est0,f,fsigma,G,O,Measures,par,units)
Gamma = [eye(q)*deltaT^2/2; eye(q)*deltaT; zeros(n-6,q)];
Q = Gamma*Q*Gamma' + [zeros(6,n); zeros(n-6,6), par.sigmaPh^2*...
diag([deltaT^2/2, deltaT, deltaT^2/2, deltaT, deltaT^2/2, deltaT, deltaT^2/2, deltaT, zeros(1,n-14)])];
+ case 22
+ Gamma = [eye(q)*deltaT^2/2; eye(q)*deltaT; zeros(n-6,q)];
+ Q = Gamma*Q*Gamma' + [zeros(6,n); zeros(n-6,6), par.sigmaPh^2*...
+ diag([deltaT^2/2, deltaT, deltaT^2/2, deltaT, deltaT^2/2, deltaT, deltaT^2/2, deltaT, zeros(1,n-14)])];
end
% Covariance matrix's time update
diff --git a/MMX_Product/.DS_Store b/MMX_Product/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..889ff5794b329496a4322246213cc9c1ae3258b0
Binary files /dev/null and b/MMX_Product/.DS_Store differ
diff --git a/MMX_Product/MMX_BSP_Files_GravLib/.DS_Store b/MMX_Product/MMX_BSP_Files_GravLib/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..418b13c2cf4e440f631bf902850e834fef6f6693
Binary files /dev/null and b/MMX_Product/MMX_BSP_Files_GravLib/.DS_Store differ
diff --git a/MMX_Product/MMX_BSP_Files_GravLib/MARPHOFRZ.txt b/MMX_Product/MMX_BSP_Files_GravLib/MARPHOFRZ.txt
new file mode 100644
index 0000000000000000000000000000000000000000..db722d0e2c3839c223ddb77c50f2c9a0a51feef1
--- /dev/null
+++ b/MMX_Product/MMX_BSP_Files_GravLib/MARPHOFRZ.txt
@@ -0,0 +1,22 @@
+\begindata
+ FRAME_MARPHOFRZ = 14499401
+ FRAME_14499401_NAME = 'MARPHOFROZEN'
+ FRAME_14499401_CLASS = 5
+ FRAME_14499401_CLASS_ID = 14499401
+ FRAME_14499401_CENTER = 499
+ FRAME_14499401_RELATIVE = 'J2000'
+ FRAME_14499401_DEF_STYLE = 'PARAMETERIZED'
+ FRAME_14499401_FAMILY = 'TWO-VECTOR'
+ FRAME_14499401_PRI_AXIS = 'X'
+ FRAME_14499401_PRI_VECTOR_DEF = 'OBSERVER_TARGET_POSITION'
+ FRAME_14499401_PRI_OBSERVER = 499
+ FRAME_14499401_PRI_TARGET = 401
+ FRAME_14499401_PRI_ABCORR = 'NONE'
+ FRAME_14499401_SEC_AXIS = 'Y'
+ FRAME_14499401_SEC_VECTOR_DEF = 'OBSERVER_TARGET_VELOCITY'
+ FRAME_14499401_SEC_OBSERVER = 499
+ FRAME_14499401_SEC_TARGET = 401
+ FRAME_14499401_SEC_ABCORR = 'NONE'
+ FRAME_14499401_SEC_FRAME = 'J2000'
+ FRAME_14499401_FREEZE_EPOCH = @2026-MAR-16/00:00:00.00
+\begintext
diff --git a/Model_Check/TestIntegrCovAnalMMXPhonly.m b/Model_Check/TestIntegrCovAnalMMXPhonly.m
index b5bd9c98ee0050dfa66e52341cfa0177df845083..b47f34f3b3c48c52591b85de5c4c31e9c6ed3ec8 100644
--- a/Model_Check/TestIntegrCovAnalMMXPhonly.m
+++ b/Model_Check/TestIntegrCovAnalMMXPhonly.m
@@ -42,7 +42,7 @@ clc
Ph0 = Ph./units.sfVec2;
% Initial MMX's State vector
- MMX0 = cspice_spkezr('-33', data, 'MARSIAU', 'none', '499')./units.sfVec;
+ MMX0 = cspice_spkezr('-34', data, 'MARSIAU', 'none', '499')./units.sfVec;
% Initial state vector
St0 = [MMX0; Ph0];
diff --git a/Model_Check/TestIntegrWithMMX.m b/Model_Check/TestIntegrWithMMX.m
index 38ed5ac8f67ebbd2dcd8cf762344102a830e2d8d..5e212d0177f4b480940a2a869dd4978d498d980a 100644
--- a/Model_Check/TestIntegrWithMMX.m
+++ b/Model_Check/TestIntegrWithMMX.m
@@ -18,6 +18,7 @@ clc
addpath(genpath('../MMX_Product/MMX_BSP_Files_GravLib/'));
MMX_InitializeSPICE
cspice_furnsh(which('mar097.bsp'));
+ cspice_furnsh(which('MARPHOFRZ.txt'));
cspice_furnsh(which('MMX_QSO_031_2x2_826891269_828619269.bsp'));
cspice_furnsh(which('Phobos_826891269_828619269.bsp'));
@@ -46,10 +47,10 @@ clc
% Initial state vector
St0 = [MMX0; Ph0];
-
+
% Integration
day = 24*3600;
- tspan = (0:150:10*day)*units.tsf;
+ tspan = (0:150:day)*units.tsf;
RelTol = 1e-13;
AbsTol = 1e-16;
opt = odeset('RelTol',RelTol,'AbsTol',AbsTol);
diff --git a/Model_Check/TestIntegrWithMMXInertia.m b/Model_Check/TestIntegrWithMMXInertia.m
index 999834fccbd629b1f281310a30f991abbee9b66a..d9e69855ad1951d783b435506a3e51ccfa621a1e 100644
--- a/Model_Check/TestIntegrWithMMXInertia.m
+++ b/Model_Check/TestIntegrWithMMXInertia.m
@@ -19,7 +19,7 @@ clc
MMX_InitializeSPICE
cspice_furnsh('MARPHOFRZ.txt');
cspice_furnsh(which('mar097.bsp'));
- cspice_furnsh(which('MMX_QSO_049_2x2_826891269_831211269.bsp'));
+ cspice_furnsh(which('MMX_QSO_031_2x2_826891269_828619269.bsp'));
% cspice_furnsh(which('MMX_3DQSO_027_009_2x2_826891269_828619269.bsp'));
% cspice_furnsh(which('MMX_SwingQSO_027_007_2x2_826891269_831211269.bsp'));
cspice_furnsh(which('Phobos_826891269_828619269.bsp'));
@@ -52,10 +52,11 @@ clc
% Initial state vector
St0 = [MMX0; Ph0];
St02 = [MMX0; Ph0; par.IPhx_bar; par.IPhy_bar; par.IPhz_bar];
+ dx = Dynamics_MPHandMMX_Inertia(0,St02,par,units);
% Integration
day = 24*3600;
- tspan = (0:100:5*day)*units.tsf;
+ tspan = (0:100:2*day)*units.tsf;
RelTol = 1e-13;
AbsTol = 1e-16;
opt = odeset('RelTol',RelTol,'AbsTol',AbsTol,'event', @(t,X) landing_Phobos(t,X,par,units));
@@ -212,4 +213,27 @@ clc
end
fclose(file);
- movefile('MMX_QSOLa_Blender.txt','../../computer-vision/MMX/');
\ No newline at end of file
+ movefile('MMX_QSOLa_Blender.txt','../../computer-vision/MMX/');
+
+
+%% Difference wrt SPICE
+
+
+ MMX_t = cspice_spkezr('-34', data+t', 'MARSIAU', 'none', '499')./units.sfVec;
+
+
+
+ figure()
+ plot(t/3600, X1(:,1) - MMX_t(1,:)')
+ grid on
+ hold on
+ plot(t/3600, X1(:,2) - MMX_t(2,:)')
+ plot(t/3600, X1(:,3) - MMX_t(3,:)')
+
+ figure()
+ plot(t/3600, X1(:,4) - MMX_t(4,:)')
+ grid on
+ hold on
+ plot(t/3600, X1(:,5) - MMX_t(5,:)')
+ plot(t/3600, X1(:,6) - MMX_t(6,:)')
+
\ No newline at end of file
diff --git a/Model_Check/TestIntegratorNew.m b/Model_Check/TestIntegratorNew.m
index 87d4736aae2cd01cbd0925991e85d2fdbc3e69e2..43416f0f79336f0fe101e665a1006a2c5a3043f5 100644
--- a/Model_Check/TestIntegratorNew.m
+++ b/Model_Check/TestIntegratorNew.m
@@ -35,7 +35,7 @@ clc
startdate = '2026-03-16 00:00:00 (UTC)';
data = cspice_str2et(startdate);
- [Ph, par] = Phobos_States_NewModel(data,par);
+ [Ph, par] = Phobos_States_NewModelPiu(data,par);
Ph = Ph./units.sfVec2;
day = 24*3600;
@@ -43,7 +43,7 @@ clc
St0 = [Ph; reshape(eye(par.d2),[par.d2^2,1])];
% Integration
- tspan = [0, 10*day]*units.tsf;
+ tspan = [0, 2*day]*units.tsf;
RelTol = 1e-13;
AbsTol = 1e-16;
opt = odeset('RelTol',RelTol,'AbsTol',AbsTol);
@@ -283,13 +283,13 @@ clc
grid on;
hold on;
plot(t/3600,mod(Phi2_t,2*pi));
- ylabel('$\dot{\Phi}_{Ph}$ [rad]','Fontsize',26)
+ ylabel('$\Phi_{Ph}$ [rad]','Fontsize',26)
xlabel('[hours]','Fontsize',26)
legend('Integrated','Spice','Fontsize',22)
subplot(2,4,7)
plot(t/3600,angdiff(X(:,7)',Phi2_t));
grid on;
- ylabel('$\Delta \dot{\Phi}_{Ph}$ [rad]','Fontsize',26)
+ ylabel('$\Delta \Phi_{Ph}$ [rad]','Fontsize',26)
xlabel('[hours]','Fontsize',26)
% Phi2_dot Spice and integration
diff --git a/New_MMXCov_Analysis.m b/New_MMXCov_Analysis.m
index fdb21e6655192656fba72da5a3d38f4591e7f04f..45f4e21c2de48ea6ec4e6374bbf7e13d4157d969 100644
--- a/New_MMXCov_Analysis.m
+++ b/New_MMXCov_Analysis.m
@@ -20,17 +20,17 @@ clc
MMX_InitializeSPICE
cspice_furnsh(which('mar097.bsp'));
cspice_furnsh(which('MARPHOFRZ.txt'));
-% cspice_furnsh(which('MMX_QSO_027_2x2_826891269_828619269.bsp'));
- cspice_furnsh(which('MMX_3DQSO_031_009_2x2_826891269_828619269.bsp'));
-% cspice_furnsh(which('MMX_SwingQSO_031_011_2x2_826891269_831211269.bsp'));
+% cspice_furnsh(which('MMX_QSO_031_2x2_826891269_828619269.bsp'));
+% cspice_furnsh(which('MMX_3DQSO_031_009_2x2_826891269_828619269.bsp'));
+ cspice_furnsh(which('MMX_SwingQSO_031_011_2x2_826891269_831211269.bsp'));
cspice_furnsh(which('Phobos_826891269_828619269.bsp'));
%% Load observations and path for syntethic pictures
- load("./ObservationsHPC/Observations-3DQSOLb_16.mat");
-% load("./Observations.mat");
+% load("./ObservationsHPC/Observations-3DQSOLb_16.mat");
+ load("./Observations.mat");
%% Initial conditions for the analysis
@@ -67,15 +67,13 @@ clc
%% Analysis
-% YObs_Full(2,:) = NaN;
-% YObs_Full(4,:) = NaN;
-% YObs_Full(6,:) = NaN;
+% YObs_Full(8,1) = NaN;
par.alpha = 1;
par.beta = 2;
- [Est] = UKF_CovarianceAnalysis(Est0,@Dynamics_MPHandMMX_Inertia,...
+ [Est] = UKF_CovarianceAnalysis(Est0, @Dynamics_MPHandMMX_Inertia,...
@Cov_Dynamics_Good, @New_Observables_model,...
- par.R,YObs_Full,par,units);
+ par.R, YObs_Full, par,units);
% save('SwingQSOLb_PN10_alpha51_16.mat','Est');
% movefile('SwingQSOLb_PN10_alpha51_16.mat','./ResultsHPC/')
diff --git a/New_MMXCov_ObsGen.m b/New_MMXCov_ObsGen.m
index 04f78970fc6f5f71985b3782e313c3a28c99c195..18d9f7d2a9250c669bde324930eed6c6b1065f9b 100644
--- a/New_MMXCov_ObsGen.m
+++ b/New_MMXCov_ObsGen.m
@@ -18,16 +18,16 @@ clc
MMX_InitializeSPICE
cspice_furnsh(which('mar097.bsp'));
cspice_furnsh(which('MARPHOFRZ.txt'));
- cspice_furnsh(which('MMX_QSO_094_2x2_826891269_828619269.bsp'));
-% cspice_furnsh(which('MMX_3DQSO_094_009_2x2_826891269_828619269.bsp'));
-% cspice_furnsh(which('MMX_SwingQSO_031_011_2x2_826891269_828619269.bsp'));
+% cspice_furnsh(which('MMX_QSO_031_2x2_826891269_828619269.bsp'));
+% cspice_furnsh(which('MMX_3DQSO_031_009_2x2_826891269_828619269.bsp'));
+ cspice_furnsh(which('MMX_SwingQSO_031_011_2x2_826891269_828619269.bsp'));
cspice_furnsh(which('Phobos_826891269_828619269.bsp'));
% Time of the analysis
data = '2026-03-16 00:00:00 (UTC)';
data = cspice_str2et(data);
day = 86400;
- PropTime = 7*day;
+ PropTime = 10*day;
date_end = data+PropTime;
% Model parameters
@@ -47,5 +47,5 @@ clc
% YObs_Full(7,:) = NaN;
% YObs_Full(end-2:end,:) = NaN;
- save('./ObservationsHPC/Observations-QSOM_16_Lidar5Min.mat','YObs_Full');
-% save('Observations.mat','YObs_Full');
+% save('./ObservationsHPC/Observations-QSOM_16_Lidar5Min.mat','YObs_Full');
+ save('Observations.mat','YObs_Full');
diff --git a/Observations.mat b/Observations.mat
index c7e8cef8a60edcb382b07c9351f2571d04f615eb..937e61d65b1f4cef677935340708541fcc579288 100644
Binary files a/Observations.mat and b/Observations.mat differ
diff --git a/PlotFinali.m b/PlotFinali.m
index 28595be058f9e4157fbae3e3d68a9551313f7352..0be4b29333b0ac98dab2a381243f9f42afac63d8 100644
--- a/PlotFinali.m
+++ b/PlotFinali.m
@@ -549,8 +549,8 @@ fineCss = round(coeff*size(t_obs,2));
[par,~] = MMX_DefineNewModelParametersAndUnits;
-C_20 = par.SH.Clm(3,1);
-C_22 = par.SH.Clm(3,3);
+C_22 = .25*(par.IPhy_bar - par.IPhx_bar);
+C_20 = -.5*(2*par.IPhz_bar - par.IPhx_bar - par.IPhy_bar);
t_obs = EstQSO.t(1,:);
@@ -590,25 +590,25 @@ xlabel('Time $[hour]$','Fontsize',26)
ylabel('$3\sigma_{C_{22}}$','Fontsize',30)
legend('QSO','3D-QSO','Swing-QSO','Interpreter','latex','Fontsize',26)
-% figure()
-% subplot(1,2,1)
-% semilogy(t_obs(1:end-1)/3600,C_20_QSO_perc,'LineWidth',1.2);
-% grid on;
-% hold on;
-% semilogy(t_obs(1:end-1)/3600,C_20_3D_perc,'LineWidth',1.2);
-% semilogy(t_obs(1:end-1)/3600,C_20_Swing_perc,'LineWidth',1.2);
-% xlabel('Time $[hour]$','Fontsize',26)
-% ylabel('%','Fontsize',26)
-% legend('QSO','3D-QSO','Swing-QSO','Interpreter','latex','Fontsize',26)
-% subplot(1,2,2)
-% semilogy(t_obs(1:end-1)/3600,C_22_QSO_perc,'LineWidth',1.2);
-% grid on;
-% hold on;
-% semilogy(t_obs(1:end-1)/3600,C_22_3D_perc,'LineWidth',1.2);
-% semilogy(t_obs(1:end-1)/3600,C_22_Swing_perc,'LineWidth',1.2);
-% xlabel('Time $[hour]$','Fontsize',26)
-% ylabel('%','Fontsize',26)
-% legend('QSO','3D-QSO','Swing-QSO','Interpreter','latex','Fontsize',26)
+figure()
+subplot(1,2,1)
+semilogy(t_obs(1:end-1)/3600,C_20_QSO_perc,'LineWidth',1.2);
+grid on;
+hold on;
+semilogy(t_obs(1:end-1)/3600,C_20_3D_perc,'LineWidth',1.2);
+semilogy(t_obs(1:end-1)/3600,C_20_Swing_perc,'LineWidth',1.2);
+xlabel('Time $[hour]$','Fontsize',26)
+ylabel('\% $|C_{20}|$','Fontsize',26,'Interpreter','latex')
+legend('QSO','3D-QSO','Swing-QSO','Interpreter','latex','Fontsize',26)
+subplot(1,2,2)
+semilogy(t_obs(1:end-1)/3600,C_22_QSO_perc,'LineWidth',1.2);
+grid on;
+hold on;
+semilogy(t_obs(1:end-1)/3600,C_22_3D_perc,'LineWidth',1.2);
+semilogy(t_obs(1:end-1)/3600,C_22_Swing_perc,'LineWidth',1.2);
+xlabel('Time $[hour]$','Fontsize',26)
+ylabel('\% $|C_{22}|$','Fontsize',26,'Interpreter','latex')
+legend('QSO','3D-QSO','Swing-QSO','Interpreter','latex','Fontsize',26)
diff --git a/Reference_Definition/BSP_generator.m b/Reference_Definition/BSP_generator.m
index ce74f4a86794d501481872ac1efcad26fe9b7f2d..5ee32809506ed32332331a943620beb4508d6b67 100644
--- a/Reference_Definition/BSP_generator.m
+++ b/Reference_Definition/BSP_generator.m
@@ -51,7 +51,7 @@ clc
% X_MMX = X_good(1:6,1);
% If MMX comes from a .bsp file
- cspice_furnsh(which('MMX_QSO_198_2x2_826891269_828619269.bsp'));
+ cspice_furnsh(which('MMX_QSO_031_2x2_826891269_828619269.bsp'));
X_MMX = cspice_spkezr('-34', data, 'MARSIAU', 'none', '499')./units.sfVec;