Skip to content

Commit c3d8165

Browse files
committed
🚨 Improve tests of MRP
1 parent db26059 commit c3d8165

1 file changed

Lines changed: 83 additions & 94 deletions

File tree

test/mrp.jl

Lines changed: 83 additions & 94 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,8 @@
66

77
# == File: ./src/mrp.jl ====================================================================
88

9+
# -- Functions: MRP ------------------------------------------------------------------------
10+
911
@testset "MRP Constructors" begin
1012
m = MRP(1.0, 2.0, 3.0)
1113
@test m.q1 == 1.0
@@ -29,29 +31,48 @@
2931
@test m.q2 == 0
3032
@test m.q3 == 0
3133

32-
# Test indexing and iteration
3334
@test m[1] == 0
3435
@test m[2] == 0
3536
@test m[3] == 0
3637
@test collect(m) == [0, 0, 0]
3738
@test length(m) == 3
3839
end
3940

41+
# -- Functions: show -----------------------------------------------------------------------
42+
4043
@testset "MRP Show" begin
41-
m = MRP(1.0, 2.0, 3.0)
42-
io = IOBuffer()
44+
buf = IOBuffer()
45+
io = IOContext(buf)
46+
m = MRP(1.0, -2.0, 3.0)
47+
48+
# Extended printing.
49+
show(io, MIME"text/plain"(), m)
50+
expected = """
51+
MRP{Float64}:
52+
X : + 1.0
53+
Y : - 2.0
54+
Z : + 3.0"""
55+
@test String(take!(io.io)) == expected
56+
57+
# Compact printing.
4358
show(io, m)
44-
s = String(take!(io))
45-
@test occursin("MRP{Float64}", s)
46-
@test occursin("1.0", s)
47-
@test occursin("2.0", s)
48-
@test occursin("3.0", s)
59+
expected = "MRP{Float64}: [1.0, 2.0, 3.0]"
60+
@test String(take!(io.io)) == expected
61+
62+
# Colors.
63+
io = IOContext(buf, :color => true)
64+
show(io, MIME"text/plain"(), m)
65+
expected = """
66+
MRP{Float64}:
67+
\e[1mX : \e[0m+ 1.0
68+
\e[1mY : \e[0m- 2.0
69+
\e[1mZ : \e[0m+ 3.0"""
70+
@test String(take!(io.io)) == expected
4971
end
5072

73+
# -- Operators: * --------------------------------------------------------------------------
5174

5275
@testset "MRP Composition" begin
53-
# m3 = m2 * m1
54-
5576
eul1 = EulerAngles(0.3, 0.2, 0.1, :ZYX)
5677
eul2 = EulerAngles(-0.2, 0.4, -0.3, :ZYX)
5778

@@ -60,130 +81,114 @@ end
6081

6182
m3 = m2 * m1
6283

63-
# Verify with DCM
64-
dcm1 = angle_to_dcm(eul1)
65-
dcm2 = angle_to_dcm(eul2)
66-
dcm3 = dcm2 * dcm1
67-
84+
# Verify with DCM.
85+
dcm1 = angle_to_dcm(eul1)
86+
dcm2 = angle_to_dcm(eul2)
87+
dcm3 = dcm2 * dcm1
6888
dcm_m3 = mrp_to_dcm(m3)
6989
@test maximum(abs.(dcm3 - dcm_m3)) < 1e-12
7090

71-
# Verify compose_rotation
72-
m_comp = compose_rotation(m1, m2) # Apply m1 then m2
73-
@test isapprox(m_comp, m3; atol=1e-12)
91+
# Verify compose_rotation.
92+
m_comp = compose_rotation(m1, m2)
93+
@test isapprox(m_comp, m3; atol = 1e-12)
7494
end
7595

96+
# -- Operators: +, -, *, / ----------------------------------------------------------------
97+
7698
@testset "MRP Arithmetic" begin
7799
m1 = MRP(0.1, 0.2, 0.3)
78100
m2 = MRP(0.05, -0.05, 0.1)
79101

80-
# +
81102
m3 = m1 + m2
82-
@test isapprox(m3.q1, 0.15; atol=1e-12)
83-
@test isapprox(m3.q2, 0.15; atol=1e-12)
84-
@test isapprox(m3.q3, 0.4; atol=1e-12)
103+
@test isapprox(m3.q1, 0.15; atol = 1e-12)
104+
@test isapprox(m3.q2, 0.15; atol = 1e-12)
105+
@test isapprox(m3.q3, 0.4; atol = 1e-12)
85106

86-
# - (binary)
87107
m4 = m1 - m2
88-
@test isapprox(m4.q1, 0.05; atol=1e-12)
89-
@test isapprox(m4.q2, 0.25; atol=1e-12)
90-
@test isapprox(m4.q3, 0.2; atol=1e-12)
108+
@test isapprox(m4.q1, 0.05; atol = 1e-12)
109+
@test isapprox(m4.q2, 0.25; atol = 1e-12)
110+
@test isapprox(m4.q3, 0.2; atol = 1e-12)
91111

92-
# - (unary)
93112
m_neg = -m1
94113
@test m_neg.q1 == -0.1
95114
@test m_neg.q2 == -0.2
96115
@test m_neg.q3 == -0.3
97116

98-
# * (scalar)
99117
m_scaled = 2.0 * m1
100118
@test m_scaled.q1 == 0.2
101119
@test m_scaled.q2 == 0.4
102120
@test m_scaled.q3 == 0.6
121+
@test m1 * 2.0 == m_scaled
103122

104-
m_scaled2 = m1 * 2.0
105-
@test m_scaled2 == m_scaled
106-
107-
# / (scalar)
108123
m_div = m1 / 2.0
109124
@test m_div.q1 == 0.05
110125
@test m_div.q2 == 0.1
111126
@test m_div.q3 == 0.15
112127
end
113128

129+
# -- Operators: inv, /, \ -----------------------------------------------------------------
130+
114131
@testset "MRP Inverse and Division" begin
115132
m1 = MRP(0.1, 0.2, 0.3)
116133
m2 = MRP(0.05, -0.05, 0.1)
117134

118-
# inv
119135
m_inv = inv(m1)
120-
# The inverse of a MRP q is -q
121136
@test m_inv == -m1
122137

123-
# Identity
124-
# m * inv(m) should be identity (0, 0, 0)
125138
m_identity = m1 * m_inv
126-
@test isapprox(m_identity.q1, 0.0; atol=1e-12)
127-
@test isapprox(m_identity.q2, 0.0; atol=1e-12)
128-
@test isapprox(m_identity.q3, 0.0; atol=1e-12)
129-
130-
# / (m1 / m2 = m1 * inv(m2))
131-
m_div = m1 / m2
132-
m_mult = m1 * inv(m2)
133-
@test isapprox(m_div, m_mult)
134-
135-
# \ (m1 \ m2 = inv(m1) * m2)
136-
m_backdiv = m1 \ m2
137-
m_mult_back = inv(m1) * m2
138-
@test isapprox(m_backdiv, m_mult_back)
139+
@test isapprox(m_identity.q1, 0.0; atol = 1e-12)
140+
@test isapprox(m_identity.q2, 0.0; atol = 1e-12)
141+
@test isapprox(m_identity.q3, 0.0; atol = 1e-12)
142+
143+
@test isapprox(m1 / m2, m1 * inv(m2))
144+
@test isapprox(m1 \ m2, inv(m1) * m2)
139145
end
140146

147+
# -- Functions: norm, one, zero, copy, vect -----------------------------------------------
148+
141149
@testset "MRP Utils" begin
142-
# 3-4-5 triangle? no, just simple numbers
143150
m = MRP(3.0, 0.0, 4.0)
144151
@test norm(m) == 5.0
145152

146-
@test one(MRP) == MRP(0.0, 0.0, 0.0)
147-
@test one(m) == MRP(0.0, 0.0, 0.0)
148-
153+
@test one(MRP) == MRP(0.0, 0.0, 0.0)
154+
@test one(m) == MRP(0.0, 0.0, 0.0)
149155
@test zero(MRP) == MRP(0.0, 0.0, 0.0)
150-
@test zero(m) == MRP(0.0, 0.0, 0.0)
156+
@test zero(m) == MRP(0.0, 0.0, 0.0)
151157

152-
# copy
153158
m_copy = copy(m)
154159
@test m_copy == m
155160

156-
# vect
157161
v = vect(m)
158162
@test v isa SVector{3, Float64}
159163
@test v[1] == m.q1
160164
@test v[2] == m.q2
161165
@test v[3] == m.q3
162166

163-
# UniformScaling
164-
@test I * m == m
165-
@test m * I == m
166-
@test I / m == inv(m)
167-
@test m / I == m
168-
@test I \ m == m
169-
@test m \ I == inv(m)
167+
@test I * m == m
168+
@test m * I == m
169+
@test I / m == inv(m)
170+
@test m / I == m
171+
@test I \ m == m
172+
@test m \ I == inv(m)
170173
end
171174

175+
# -- Functions: shadow_rotation -----------------------------------------------------------
176+
172177
@testset "MRP Shadow Rotation" begin
173178
m = MRP(0.1, 0.2, 0.3)
174179
m_shadow = shadow_rotation(m)
175180

176-
# Check if the shadow rotation is indeed the shadow
177181
@test norm(m_shadow) > 1.0
178-
@test isapprox(norm(m_shadow), 1/norm(m); atol=1e-12)
182+
@test isapprox(norm(m_shadow), 1 / norm(m); atol = 1e-12)
179183

180-
# Check if they represent the same rotation
181-
dcm = mrp_to_dcm(m)
184+
# The shadow rotation must represent the same rotation.
185+
dcm = mrp_to_dcm(m)
182186
dcm_shadow = mrp_to_dcm(m_shadow)
183-
184187
@test maximum(abs.(dcm - dcm_shadow)) < 1e-12
185188
end
186189

190+
# -- Functions: rand -----------------------------------------------------------------------
191+
187192
@testset "MRP Random" begin
188193
Random.seed!(123)
189194
m = rand(MRP)
@@ -193,47 +198,31 @@ end
193198
@test m_float32 isa MRP{Float32}
194199
end
195200

201+
# -- Functions: dmrp -----------------------------------------------------------------------
202+
196203
@testset "MRP Kinematics" begin
197-
# Test analytical derivative against numerical derivative
198204
for T in (Float64,)
199-
# We need a seed to ensure reproducibility
200205
Random.seed!(123)
201206

202-
m = rand(MRP{T})
203-
w = @SVector randn(T, 3)
204-
207+
m = rand(MRP{T})
208+
w = @SVector randn(T, 3)
205209
dm = dmrp(m, w)
206210

207211
dt = 1e-8
208212

209-
# Use DCM mapping for propagation to avoid circular dependency if we were to use dmrp for integration,
210-
# but here we want to verify dmrp.
211-
D = mrp_to_dcm(m)
212-
213-
# In the existing tests (kinematics.jl), they use:
214-
# dDba = ddcm(Dba, Dba * wba_a)
215-
# Dba = Dba + dDba * Δ
216-
# Dba = orthonormalize(Dba)
217-
218-
dD = ddcm(D, w) # w is in body frame? Make sure. ddcm says wba_b. dmrp says wba_b. Matches.
219-
D_next = D + dD * dt
220-
D_next = orthonormalize(D_next)
213+
D = mrp_to_dcm(m)
214+
dD = ddcm(D, w)
215+
D_next = orthonormalize(D + dD * dt)
221216

222217
m_next = dcm_to_mrp(D_next)
223218

224-
# Check if we switched to the shadow set (norm > 1 or just large jump)
225-
# If m_next is far from m, try shadow
219+
# Switch to the shadow set if necessary to ensure continuity.
226220
if norm(m_next - m) > 0.1
227-
m_next = shadow_rotation(m_next)
221+
m_next = shadow_rotation(m_next)
228222
end
229223

230224
dm_num = (m_next - m) / dt
231225

232-
# Check alignment
233-
@test isapprox(dm, dm_num; rtol = 1e-4, atol=1e-6)
226+
@test isapprox(dm, dm_num; rtol = 1e-4, atol = 1e-6)
234227
end
235228
end
236-
237-
238-
239-

0 commit comments

Comments
 (0)