diff --git a/.DS_Store b/.DS_Store
index a74eb23937f1a297326f4be205ecc327dd7de2c3..0546b8ed14131c54188822cc34de5acd427c4dc8 100644
Binary files a/.DS_Store and b/.DS_Store differ
diff --git a/Analisys_With_Features.m b/Analisys_With_Features.m
index 138fa7f9f1c6662095e132d2624e53ffb22880a2..169fdf336bd5960ad5c63ca791bb7ae81a90d935 100644
--- a/Analisys_With_Features.m
+++ b/Analisys_With_Features.m
@@ -41,8 +41,8 @@ clc
data = '2026-03-16 00:00:00 (UTC)';
data = cspice_str2et(data);
day = 86400;
- pars.et0 = data;
- [Ph,pars] = Phobos_States_NewModel(data,pars);
+ pars.et0 = data;
+ [Ph,pars] = Phobos_States_NewModel(data,pars);
% Covariance analysis parsameters
[pars, units] = CovarianceAnalysisParameters(pars, units);
@@ -61,10 +61,10 @@ clc
St0 = [MMX0; Ph0(1:4); Ph0(7:8); pars.I2; pars.bias];
pars.X0_reference = [MMX0; Ph0; pars.I2; pars.bias];
- % Est0.dx = zeros(size(St0,1),1);
- % pars.flag_FILTERING = 0;
- Est0.dx = St0.*randn(size(St0))*1e-6;
- pars.flag_FILTERING = 1;
+ Est0.dx = zeros(size(St0,1),1);
+ pars.flag_FILTERING = 0;
+ % Est0.dx = St0.*randn(size(St0))*1e-6;
+ % pars.flag_FILTERING = 1;
Est0.X = St0+Est0.dx;
Est0.P0 = pars.P0;
Est0.t0 = data*units.tsf;
@@ -84,9 +84,9 @@ clc
% [Est] = UKF_Tool(Est0,@Dynamics_MPHandMMX_Inertia,...
% @Cov_Dynamics_Good, @Observables_model_with_Features,...
% pars.R,YObs,pars,units,file_features);
- % [Est] = UKF_Tool(Est0,@Dynamics_MPHandMMX_Inertia,...
- % @Sigma_dynamics_NoPhi, @Observables_model_with_Features,...
- % pars.R,YObs,pars,units,file_features);
+ [Est] = UKF_Tool(Est0,@Dynamics_MPHandMMX_Inertia,...
+ @Sigma_dynamics_NoPhi, @Observables_model_with_Features,...
+ pars.R,YObs,pars,units,file_features);
% Questo è l'UKF puro, funziona giocando con l'alpha! Ma non tutti gli
% stati sono visibilissimi. Vedo il libration ma non la sua velocita'
@@ -97,11 +97,11 @@ clc
% @Observables_model_with_Features,...
% pars.R,YObs,pars,units,file_features);
- pars.alpha = 1e-1;
- pars.beta = 2;
- [Est] = UKF_features(Est0, @SigmaPoints_Dynamics_Good_NoPhi,...
- @Observables_model_with_Features,...
- pars.R,YObs,pars,units,file_features);
+ % pars.alpha = 1e-1;
+ % pars.beta = 2;
+ % [Est] = UKF_features(Est0, @SigmaPoints_Dynamics_Good_NoPhi,...
+ % @Observables_model_with_Features,...
+ % pars.R,YObs,pars,units,file_features);
Results(Est, pars, units)
diff --git a/CovarianceAnalysisParameters.m b/CovarianceAnalysisParameters.m
index 31e7d6f34faf8fc58671fc49c4a98bee32a707d8..4c250fd4953419b854359dcba384a86fa0eaa748 100644
--- a/CovarianceAnalysisParameters.m
+++ b/CovarianceAnalysisParameters.m
@@ -55,7 +55,7 @@ pars.FOV = deg2rad(45);
pixels = 2472;
pars.ObsNoise.camera = pars.FOV/pixels; % [rad/pixel]
-pars.ObsNoise.pixel = .25; % [n. pixel]
+pars.ObsNoise.pixel = .2; % [n. pixel]
% Useful to adimensionalize the observables
diff --git a/Deimos_Pic.m b/Deimos_Pic.m
new file mode 100644
index 0000000000000000000000000000000000000000..7f0da830f9a6ecec1df184b48fead58a43ce3f4f
--- /dev/null
+++ b/Deimos_Pic.m
@@ -0,0 +1,100 @@
+function [Deimos_LOS, Deimos_Limb] = Deimos_Pic(MMX, Phobos, Deimos, Sun, pars, units)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% Y = Deimos_Pic(MMX, Phobos, Deimos, Sun, Phi, pars, units)
+%
+% Identify wether Deimos is in the FOV and how far is it
+%
+% Input:
+%
+% MMX MMX state vector in the MARSIAU reference frame
+% Phobos Phobos state vector in the MARSIAU reference frame
+% Deimos Deimos state vector in the MARSIAU reference frame
+% Sun Sun state vector in the MARSIAU reference frame
+% pars parameters structure
+% units units structure
+%
+% Output:
+%
+% Limb_Range Range derived from Mars limb size
+%
+%
+%
+% Author: E.Ciccarelli
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+% Rotation from MARSIAU to Phobos radial reference frame
+
+ i_feat = Phobos(1:3)/norm(Phobos(1:3));
+ v = Phobos(4:6)/norm(Phobos(4:6));
+ k = cross(i_feat, v)/norm(cross(i_feat, v));
+ j = cross(k, i_feat)/norm(cross(k, i_feat));
+
+ R2rot = [i_feat, j, k]';
+
+% MMX-Phobos direction in the rotating reference frame
+ r_sb = MMX(1:3) - Phobos(1:3);
+ r_sb_Ph = R2rot*r_sb;
+ R_sb_Ph = norm(r_sb_Ph);
+ I_sb = r_sb_Ph/R_sb_Ph;
+
+% MMX-Mars direction in the rotating reference frame
+ r_sM = R2rot*MMX(1:3);
+ I_sM = r_sM/norm(r_sM);
+
+% Mars pole position in the Phobos rotating frame
+ Deimos_center = Deimos(1:3) - Phobos(1:3);
+ R_ini = [0, -1, 0; 1, 0, 0; 0, 0, 1];
+ Deimos_center = R_ini*R2rot'*Deimos_center;
+ Deimos_pole = [0; 15; 0]+Deimos(1:3) - Phobos(1:3);
+ R_ini = [0, -1, 0; 1, 0, 0; 0, 0, 1];
+ Deimos_pole = R_ini*R2rot'*Deimos_pole;
+
+% Phobos pole position in the Phobos rotating frame
+ Phobos_pole = [0; 0; -pars.Phobos.gamma];
+ Phobos_pole = R_ini*Phobos_pole;
+ Phobos_centre = zeros(3,1);
+ Phobos_centre = R_ini*Phobos_centre;
+
+% Sunlight direction in the rotating reference frame
+ Sun = R2rot*Sun(1:3);
+ I_Sun = Sun/norm(Sun);
+
+ angle_MMXSun = acos(dot(I_sM,I_Sun));
+
+ [~, ~, PM] = ProjMatrix(-r_sb_Ph, pars, units);
+ DC = PM*[-Deimos_center; 1];
+ DC = DC./repmat(DC(3),3,1);
+ DP = PM*[-Deimos_pole; 1];
+ DP = DP./repmat(DP(3),3,1);
+
+ PhP = PM*[-Phobos_pole; 1];
+ PhP = PhP./repmat(PhP(3),3,1);
+ PhC = PM*[-Phobos_centre; 1];
+ PhC = PhC./repmat(PhC(3),3,1);
+
+ Phobos_range = norm(PhP - PhC);
+ Phobos_ellipsoid = @(x,y) (x-PhC(1))^2 + (y-PhC(2))^2 - Phobos_range^2;
+
+ if (any(DC<0))||(DC(1)>pars.cam.nPixX)||(DC(2)>pars.cam.nPixY)||...
+ (dot(r_sb_Ph, Deimos_center)>0)||(Phobos_ellipsoid(DC(1),DC(2))<0)
+ % if (Phobos_ellipsoid(DC(1),DC(2))<0)
+ % Plotfeatures_Pic([DC, PhC], [DP, PhP], pars);
+ % end
+ Deimos_LOS = [];
+ Deimos_Limb = [];
+ else
+ Deimos_LOS = DP;
+ Limb_radius = DP - DC;
+ Deimos_Limb = sqrt(Limb_radius(1)^2 + Limb_radius(2)^2) + random('Normal',0, pars.ObsNoise.pixel);
+ Deimos_Limb = Deimos_Limb - rem(Deimos_Limb,pars.ObsNoise.pixel);
+ if Deimos_Limb<1.5
+ Deimos_Limb = [];
+ end
+ % Plotfeatures_Pic([DC, PhC], [DP, PhP], pars);
+ end
+
+
+end
\ No newline at end of file
diff --git a/Mars_LimbRange.m b/Mars_LimbRange.m
index 3962de9f211bd9dd15d8adfa36f4f4ff7567efb6..9ed3711cde37d63b83a7e10ea97ff08eafd8eb82 100644
--- a/Mars_LimbRange.m
+++ b/Mars_LimbRange.m
@@ -73,15 +73,24 @@ function Limb_Range = Mars_LimbRange(MMX, Phobos, Sun, pars, units)
PhC = PM*[-Phobos_centre; 1];
PhC = PhC./repmat(PhC(3),3,1);
+ Phobos_range = norm(PhP - PhC);
+ Phobos_ellipsoid = @(x,y) (x-PhC(1))^2 + (y-PhC(2))^2 - Phobos_range^2;
+
if ((any(MC<0))||(any(MP<0))||((MC(1)>pars.cam.nPixX))||...
(MP(1)>pars.cam.nPixX)||(MC(2)>pars.cam.nPixY)||...
(MP(2)>pars.cam.nPixY))||(dot(r_sb_Ph, Mars_centre)>0)||...
((dot(I_sM,I_Sun)<0)&&((angle_MMXSun-pi/2)<pars.FOV))||...
- (PhP(2)<MP(2)*1.1)
+ (Phobos_ellipsoid(MP(1),MP(2))<0)
+
+ % if (Phobos_ellipsoid(MP(1),MP(2))<0)
+ % Plotfeatures_Pic([MC, PhC], [MP, PhP], pars);
+ % end
Limb_Range = [];
+
else
Limb_radius = MP - MC;
- Limb_Range = round(sqrt(Limb_radius(1)^2 + Limb_radius(2)^2) + random('Normal',0, pars.ObsNoise.pixel));
+ Limb_Range = sqrt(Limb_radius(1)^2 + Limb_radius(2)^2) + random('Normal',0, pars.ObsNoise.pixel);
+ Limb_Range = Limb_Range - rem(Limb_Range,pars.ObsNoise.pixel);
% Plotfeatures_Pic([MC, PhC], [MP, PhP], pars);
end
diff --git a/Observables_model_with_Features.m b/Observables_model_with_Features.m
index 085cfc21b213024a69072700f0eb69444f4d30f8..85ed41ae7a7c6d1bdd829347ca066edd48c2799b 100644
--- a/Observables_model_with_Features.m
+++ b/Observables_model_with_Features.m
@@ -61,10 +61,6 @@ function [G, Rm] = Observables_model_with_Features(et, X, St_ID_Range, St_ID_Rat
r_Phobos = pars.perifocal2MARSIAU*X(d+1)*[cos(X(d+3)); sin(X(d+3)); 0];
switch size(X,1)
- case 17
- bias = pars.bias;
- case 18
- bias = X(16:end);
case 20
bias = X(16:end);
case 22
@@ -212,7 +208,7 @@ function [G, Rm] = Observables_model_with_Features(et, X, St_ID_Range, St_ID_Rat
Phobos = [r_Phobos.*units.lsf; v_Phobos.*units.vsf];
Limb = Mars_LimbRange_model(X_MMX.*units.sfVec, Phobos, Sun, pars, units);
- Limb = Limb + sqrt(bias(4)^2 + bias(5)^2);
+ Limb = Limb + (bias(4) + bias(5));
Rm = [Rm; O(8,8)];
end
diff --git a/Observations_Generation.m b/Observations_Generation.m
index 82af073e5710265584866c70300b4b009f9a9a15..b8845221768bcd3127c48e65edb816304709ede0 100644
--- a/Observations_Generation.m
+++ b/Observations_Generation.m
@@ -30,7 +30,7 @@ clc
% Time of the analysis
data = '2026-03-16 00:00:00 (UTC)';
data = cspice_str2et(data);
- Ndays = 5;
+ Ndays = 7;
t = Ndays*86400;
date_end = data+t;
diff --git a/Observations_with_features.m b/Observations_with_features.m
index 5b6a9b3ce56dd507097ecd94ea6c308a4f378ffe..40143be958c2cd97ed5f80b74482deb01d54d44a 100644
--- a/Observations_with_features.m
+++ b/Observations_with_features.m
@@ -225,7 +225,6 @@ function YObs = Observations_with_features(date_0, date_end, file_features, pars
else
Y_LOS = [];
Y_pix = [];
- Limb = [];
end
if (rem(round(i-et0),interval_camera)==0)&&pars.flag_Limb
@@ -239,6 +238,17 @@ function YObs = Observations_with_features(date_0, date_end, file_features, pars
end
+ if (rem(round(i-et0),interval_camera)==0)&&pars.flag_Limb
+
+ Deimos = cspice_spkezr('402', i, 'MARSIAU', 'none', 'MARS');
+% Is Deimos in the picture?
+ [Deimos_LOS, Deimos_Limb] = Deimos_Pic(MMX, Phobos, Deimos, Sun, pars, units);
+
+ if ~isempty(Deimos_Limb)
+ got_it = 1;
+ end
+
+ end
%--------------------------------------------------------------------------
@@ -261,9 +271,10 @@ function YObs = Observations_with_features(date_0, date_end, file_features, pars
Obs.Rate = [];
end
- Obs.Lidar = Lidar;
- Obs.Features= Y_pix;
- Obs.Limb = Limb;
+ Obs.Lidar = Lidar;
+ Obs.Features = Y_pix;
+ Obs.Limb = Limb;
+ Obs.Deimos_Limb = Deimos_Limb;
YObs = [YObs; Obs];
end
end
diff --git a/Radial_Along_Cross_Covariance.m b/Radial_Along_Cross_Covariance.m
index d0e1ebd789679e39fe643b8d949bebe4d66ffee3..e571b6eb82564083aacdc0b260be20f3d0d466b8 100644
--- a/Radial_Along_Cross_Covariance.m
+++ b/Radial_Along_Cross_Covariance.m
@@ -31,8 +31,7 @@ function P_RAC = Radial_Along_Cross_Covariance(P_XYZ, X, pars, units)
K = cross(I,v)/norm(cross(I,v));
J = cross(K,I);
- R_NOPh = [I, J, K]';
- R_tot = R_NOPh;
+ R_tot = [I, J, K]';
% Covariance rotation
P_RAC = R_tot*P_XYZ(1:3,1:3)*R_tot';
diff --git a/Results.m b/Results.m
index df0739c0c4c3896befdddc1daaa2abfe008e2089..24ceda83567983f46f9efefc13c5cf4969790dee 100644
--- a/Results.m
+++ b/Results.m
@@ -531,65 +531,65 @@ function Results(Est, pars, units)
%--------------------------------------------------------------------------
- figure()
- subplot(1,2,1)
- semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
- hold on;
- grid on;
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,SQRT_X_Ph,'Color','k','LineWidth',2)
- xlabel('time $[hour]$','FontSize',30)
- ylabel('$[km]$','FontSize',30)
- title('Phobos position vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
- legend('$3\sigma_{x}$','$3\sigma_{y}$','$3\sigma_{z}$','$3 RMS$','Interpreter','latex','FontSize',26)
- subplot(1,2,2)
- semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhVXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
- hold on;
- grid on;
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,SQRT_V_Ph,'Color','k','LineWidth',2)
- xlabel('time $[hour]$','FontSize',30)
- ylabel('$[km]$','FontSize',30)
- title('Phobos velocity vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
- legend('$3\sigma_{\dot{x}}$','$3\sigma_{\dot{y}}$','$3\sigma_{\dot{z}}$','$3 RMS$','Interpreter','latex','FontSize',26)
-
-
- figure()
- subplot(1,3,1)
- semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
- hold on;
- grid on;
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,SQRT_X_Ph,'Color','k','LineWidth',2)
- xlabel('time $[hour]$','FontSize',30)
- ylabel('$[km]$','FontSize',30)
- title('Phobos position vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
- legend('$3\sigma_{x}$','$3\sigma_{y}$','$3\sigma_{z}$','$3 RMS$','Interpreter','latex','FontSize',26)
- subplot(1,3,2)
- semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhVXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
- hold on;
- grid on;
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,SQRT_V_Ph,'Color','k','LineWidth',2)
- xlabel('time $[hour]$','FontSize',30)
- ylabel('$[km]$','FontSize',30)
- title('Phobos velocity vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
- legend('$3\sigma_{\dot{x}}$','$3\sigma_{\dot{y}}$','$3\sigma_{\dot{z}}$','$3 RMS$','Interpreter','latex','FontSize',26)
- subplot(1,3,3)
- semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_rel_MMXPh(1,1,1:end-1)))),'Color','b','LineWidth',1)
- hold on;
- grid on;
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_rel_MMXPh(2,2,1:end-1)))),'Color','r','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_rel_MMXPh(3,3,1:end-1)))),'Color','g','LineWidth',1)
- semilogy(t_obs(1:end-1)/3600,SQRT_X_rel,'Color','k','LineWidth',2)
- xlabel('time $[hour]$','FontSize',30)
- ylabel('$[km]$','FontSize',30)
- title('MMX-Phobos relative distance $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
- legend('$3\sigma_{x}$','$3\sigma_{y}$','$3\sigma_{z}$','$3 RMS$','Interpreter','latex','FontSize',26)
+ % figure()
+ % subplot(1,2,1)
+ % semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
+ % hold on;
+ % grid on;
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,SQRT_X_Ph,'Color','k','LineWidth',2)
+ % xlabel('time $[hour]$','FontSize',30)
+ % ylabel('$[km]$','FontSize',30)
+ % title('Phobos position vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
+ % legend('$3\sigma_{x}$','$3\sigma_{y}$','$3\sigma_{z}$','$3 RMS$','Interpreter','latex','FontSize',26)
+ % subplot(1,2,2)
+ % semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhVXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
+ % hold on;
+ % grid on;
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,SQRT_V_Ph,'Color','k','LineWidth',2)
+ % xlabel('time $[hour]$','FontSize',30)
+ % ylabel('$[km]$','FontSize',30)
+ % title('Phobos velocity vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
+ % legend('$3\sigma_{\dot{x}}$','$3\sigma_{\dot{y}}$','$3\sigma_{\dot{z}}$','$3 RMS$','Interpreter','latex','FontSize',26)
+
+
+ % figure()
+ % subplot(1,3,1)
+ % semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
+ % hold on;
+ % grid on;
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,SQRT_X_Ph,'Color','k','LineWidth',2)
+ % xlabel('time $[hour]$','FontSize',30)
+ % ylabel('$[km]$','FontSize',30)
+ % title('Phobos position vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
+ % legend('$3\sigma_{x}$','$3\sigma_{y}$','$3\sigma_{z}$','$3 RMS$','Interpreter','latex','FontSize',26)
+ % subplot(1,3,2)
+ % semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_PhVXYZ(1,1,1:end-1)))),'Color','b','LineWidth',1)
+ % hold on;
+ % grid on;
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(2,2,1:end-1)))),'Color','r','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_PhVXYZ(3,3,1:end-1)))),'Color','g','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,SQRT_V_Ph,'Color','k','LineWidth',2)
+ % xlabel('time $[hour]$','FontSize',30)
+ % ylabel('$[km]$','FontSize',30)
+ % title('Phobos velocity vector $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
+ % legend('$3\sigma_{\dot{x}}$','$3\sigma_{\dot{y}}$','$3\sigma_{\dot{z}}$','$3 RMS$','Interpreter','latex','FontSize',26)
+ % subplot(1,3,3)
+ % semilogy(t_obs(1:end-1)/3600,3.*squeeze(real(sqrt(P_rel_MMXPh(1,1,1:end-1)))),'Color','b','LineWidth',1)
+ % hold on;
+ % grid on;
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_rel_MMXPh(2,2,1:end-1)))),'Color','r','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,3*squeeze(real(sqrt(P_rel_MMXPh(3,3,1:end-1)))),'Color','g','LineWidth',1)
+ % semilogy(t_obs(1:end-1)/3600,SQRT_X_rel,'Color','k','LineWidth',2)
+ % xlabel('time $[hour]$','FontSize',30)
+ % ylabel('$[km]$','FontSize',30)
+ % title('MMX-Phobos relative distance $3\sigma$ envelopes','Interpreter','latex','FontSize',30)
+ % legend('$3\sigma_{x}$','$3\sigma_{y}$','$3\sigma_{z}$','$3 RMS$','Interpreter','latex','FontSize',26)
figure()
diff --git a/YObs.mat b/YObs.mat
index c67140bb505a949fa5225d3554934badb37ebc20..8098dfb4c16bbb7ab3ecd00428b7a4c8f51bed77 100644
Binary files a/YObs.mat and b/YObs.mat differ
diff --git a/visible_features.m b/visible_features.m
index 2119bc4b17607595c608f8c0b1f6b96011a4187c..3df18910ba0a62559c80d3af76a25d197b67746f 100644
--- a/visible_features.m
+++ b/visible_features.m
@@ -115,8 +115,9 @@ function [Y_LOS, Y_pix] = visible_features(MMX, Phobos, Sun, Phi, file_features,
[~, ~, PM] = ProjMatrix(-r_sb_Ph, pars, units);
pix_feat = PM*[-candidate; 1];
pix_feat = pix_feat./repmat(pix_feat(3),3,1);
- Y_pix = [[(pix_feat(1:2)+ [random('Normal',0, pars.ObsNoise.pixel);...
- random('Normal',0, pars.ObsNoise.pixel)]); point_in_light(end,n)], Y_pix];
+ pix_feat = [(pix_feat(1:2)+ [random('Normal',0, pars.ObsNoise.pixel);...
+ random('Normal',0, pars.ObsNoise.pixel)]); point_in_light(end,n)];
+ Y_pix = [pix_feat - rem(pix_feat,pars.ObsNoise.pixel), Y_pix];
end
end