Line data Source code
1 : /**
2 : * @file test_EmissionLogLikelihood.cpp
3 : *
4 : * @brief Tests for the EmissionLogLikelihood class
5 : *
6 : * @author Matthias Wieczorek - initial code
7 : * @author David Frank - rewrite
8 : * @author Tobias Lasser - rewrite
9 : */
10 :
11 : #include <doctest/doctest.h>
12 :
13 : #include <cmath>
14 : #include "testHelpers.h"
15 : #include "EmissionLogLikelihood.h"
16 : #include "LinearResidual.h"
17 : #include "Scaling.h"
18 : #include "Identity.h"
19 : #include "VolumeDescriptor.h"
20 : #include "TypeCasts.hpp"
21 :
22 : using namespace elsa;
23 : using namespace doctest;
24 :
25 : TEST_SUITE_BEGIN("functionals");
26 :
27 : TEST_CASE_TEMPLATE("EmissionLogLikelihood: Testing with only data no residual", TestType, float,
28 : double)
29 10 : {
30 10 : using Vector = Eigen::Matrix<TestType, Eigen::Dynamic, 1>;
31 :
32 10 : GIVEN("just data (no residual)")
33 10 : {
34 10 : IndexVector_t numCoeff(3);
35 10 : numCoeff << 9, 15, 19;
36 10 : VolumeDescriptor dd(numCoeff);
37 :
38 10 : Vector y(dd.getNumberOfCoefficients());
39 10 : y.setRandom();
40 10 : DataContainer<TestType> dcY(dd, y);
41 :
42 10 : WHEN("instantiating without r")
43 10 : {
44 6 : EmissionLogLikelihood func(dd, dcY);
45 :
46 6 : THEN("the functional is as expected")
47 6 : {
48 2 : REQUIRE_EQ(func.getDomainDescriptor(), dd);
49 :
50 2 : auto* linRes = downcast_safe<LinearResidual<TestType>>(&func.getResidual());
51 2 : REQUIRE_UNARY(linRes);
52 2 : REQUIRE_UNARY_FALSE(linRes->hasDataVector());
53 2 : REQUIRE_UNARY_FALSE(linRes->hasOperator());
54 2 : }
55 :
56 6 : THEN("a clone behaves as expected")
57 6 : {
58 2 : auto emllClone = func.clone();
59 :
60 2 : REQUIRE_NE(emllClone.get(), &func);
61 2 : REQUIRE_EQ(*emllClone, func);
62 2 : }
63 :
64 6 : THEN("the evaluate, gradient and Hessian work as expected")
65 6 : {
66 2 : Vector dataVec(dd.getNumberOfCoefficients());
67 2 : dataVec.setRandom();
68 5132 : for (index_t i = 0; i < dataVec.size(); ++i) // ensure non-negative numbers
69 5130 : if (dataVec[i] < 0)
70 2640 : dataVec[i] *= -1;
71 2 : DataContainer<TestType> x(dd, dataVec);
72 :
73 : // compute the "true" values
74 2 : TestType trueValue = 0;
75 2 : Vector trueGrad(dd.getNumberOfCoefficients());
76 2 : Vector trueScale(dd.getNumberOfCoefficients());
77 5132 : for (index_t i = 0; i < dataVec.size(); ++i) {
78 5130 : TestType temp = dataVec[i];
79 5130 : trueValue += temp - y[i] * std::log(temp);
80 5130 : trueGrad[i] = 1 - y[i] / temp;
81 5130 : trueScale[i] = y[i] / (temp * temp);
82 5130 : }
83 :
84 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
85 2 : DataContainer<TestType> dcTrueGrad(dd, trueGrad);
86 2 : REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
87 :
88 2 : DataContainer<TestType> dcTrueScale(dd, trueScale);
89 2 : REQUIRE_EQ(func.getHessian(x), leaf(Scaling(dd, dcTrueScale)));
90 2 : }
91 6 : }
92 :
93 10 : WHEN("instantiating with r")
94 10 : {
95 4 : Vector r(dd.getNumberOfCoefficients());
96 4 : r.setRandom();
97 10264 : for (index_t i = 0; i < r.size(); ++i)
98 10260 : if (r[i] < 0)
99 5153 : r[i] *= -1;
100 4 : DataContainer<TestType> dcR(dd, r);
101 :
102 4 : EmissionLogLikelihood func(dd, dcY, dcR);
103 :
104 4 : THEN("a clone behaves as expected")
105 4 : {
106 2 : auto emllClone = func.clone();
107 :
108 2 : REQUIRE_NE(emllClone.get(), &func);
109 2 : REQUIRE_EQ(*emllClone, func);
110 2 : }
111 :
112 4 : THEN("the evaluate, gradient and Hessian work as expected")
113 4 : {
114 2 : Vector dataVec(dd.getNumberOfCoefficients());
115 2 : dataVec.setRandom();
116 5132 : for (index_t i = 0; i < dataVec.size(); ++i) // ensure non-negative numbers
117 5130 : if (dataVec[i] < 0)
118 2575 : dataVec[i] *= -1;
119 2 : DataContainer<TestType> x(dd, dataVec);
120 :
121 : // compute the "true" values
122 2 : TestType trueValue = 0;
123 2 : Vector trueGrad(dd.getNumberOfCoefficients());
124 2 : Vector trueScale(dd.getNumberOfCoefficients());
125 5132 : for (index_t i = 0; i < dataVec.size(); ++i) {
126 5130 : auto temp = dataVec[i] + r[i];
127 5130 : trueValue += temp - y[i] * std::log(temp);
128 5130 : trueGrad[i] = 1 - y[i] / temp;
129 5130 : trueScale[i] = y[i] / (temp * temp);
130 5130 : }
131 :
132 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
133 2 : DataContainer<TestType> dcTrueGrad(dd, trueGrad);
134 2 : REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
135 :
136 2 : DataContainer<TestType> dcTrueScale(dd, trueScale);
137 2 : REQUIRE_EQ(func.getHessian(x), leaf(Scaling(dd, dcTrueScale)));
138 2 : }
139 4 : }
140 10 : }
141 10 : }
142 :
143 : TEST_CASE_TEMPLATE("EmissionLogLikelihood: Testing with residual", TestType, float, double)
144 10 : {
145 10 : using Vector = Eigen::Matrix<TestType, Eigen::Dynamic, 1>;
146 :
147 10 : GIVEN("a residual with data")
148 10 : {
149 10 : IndexVector_t numCoeff(3);
150 10 : numCoeff << 3, 15, 21;
151 10 : VolumeDescriptor dd(numCoeff);
152 :
153 10 : Vector resData(dd.getNumberOfCoefficients());
154 10 : resData.setRandom();
155 10 : DataContainer<TestType> dcResData(dd, resData);
156 10 : Identity<TestType> idOp(dd);
157 10 : LinearResidual<TestType> linRes(idOp, dcResData);
158 :
159 10 : Vector y(dd.getNumberOfCoefficients());
160 10 : y.setRandom();
161 10 : DataContainer<TestType> dcY(dd, y);
162 :
163 10 : WHEN("instantiating without r")
164 10 : {
165 6 : EmissionLogLikelihood func(linRes, dcY);
166 :
167 6 : THEN("the functional is as expected")
168 6 : {
169 2 : REQUIRE_EQ(func.getDomainDescriptor(), dd);
170 :
171 2 : auto* lRes = downcast_safe<LinearResidual<TestType>>(&func.getResidual());
172 2 : REQUIRE_UNARY(lRes);
173 2 : REQUIRE_EQ(*lRes, linRes);
174 2 : }
175 :
176 6 : THEN("a clone behaves as expected")
177 6 : {
178 2 : auto emllClone = func.clone();
179 :
180 2 : REQUIRE_NE(emllClone.get(), &func);
181 2 : REQUIRE_EQ(*emllClone, func);
182 2 : }
183 :
184 6 : THEN("the evaluate, gradient and Hessian work as expected")
185 6 : {
186 2 : Vector dataVec(dd.getNumberOfCoefficients());
187 2 : dataVec.setRandom();
188 : // ensure non-negative numbers
189 1892 : for (index_t i = 0; i < dataVec.size(); ++i) {
190 1890 : if (dataVec[i] - resData[i] < 0)
191 995 : dataVec[i] -= 2 * (dataVec[i] - resData[i]);
192 1890 : }
193 2 : DataContainer<TestType> x(dd, dataVec);
194 :
195 : // compute the "true" values
196 2 : TestType trueValue = 0;
197 2 : Vector trueGrad(dd.getNumberOfCoefficients());
198 2 : Vector trueScale(dd.getNumberOfCoefficients());
199 1892 : for (index_t i = 0; i < dataVec.size(); ++i) {
200 1890 : auto temp = dataVec[i] - resData[i];
201 1890 : trueValue += temp - y[i] * std::log(temp);
202 1890 : trueGrad[i] = 1 - y[i] / temp;
203 1890 : trueScale[i] = y[i] / (temp * temp);
204 1890 : }
205 :
206 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
207 2 : DataContainer<TestType> dcTrueGrad(dd, trueGrad);
208 2 : REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
209 :
210 2 : auto hessian = func.getHessian(x);
211 2 : auto hx = hessian.apply(x);
212 1892 : for (index_t i = 0; i < hx.getSize(); ++i)
213 2 : REQUIRE_UNARY(checkApproxEq(hx[i], dataVec[i] * trueScale[i]));
214 2 : }
215 6 : }
216 :
217 10 : WHEN("instantiating with r")
218 10 : {
219 4 : Vector r(dd.getNumberOfCoefficients());
220 4 : r.setRandom();
221 3784 : for (index_t i = 0; i < r.size(); ++i)
222 3780 : if (r[i] < 0)
223 1937 : r[i] *= -1;
224 4 : DataContainer<TestType> dcR(dd, r);
225 :
226 4 : EmissionLogLikelihood func(linRes, dcY, dcR);
227 :
228 4 : THEN("a clone behaves as expected")
229 4 : {
230 2 : auto emllClone = func.clone();
231 :
232 2 : REQUIRE_NE(emllClone.get(), &func);
233 2 : REQUIRE_EQ(*emllClone, func);
234 2 : }
235 :
236 4 : THEN("the evaluate, gradient and Hessian work as expected")
237 4 : {
238 2 : Vector dataVec(dd.getNumberOfCoefficients());
239 2 : dataVec.setRandom();
240 : // ensure non-negative numbers
241 1892 : for (index_t i = 0; i < dataVec.size(); ++i) {
242 1890 : if (dataVec[i] - resData[i] < 0)
243 965 : dataVec[i] -= 2 * (dataVec[i] - resData[i]);
244 1890 : }
245 2 : DataContainer<TestType> x(dd, dataVec);
246 :
247 : // compute the "true" values
248 2 : TestType trueValue = 0;
249 2 : Vector trueGrad(dd.getNumberOfCoefficients());
250 2 : Vector trueScale(dd.getNumberOfCoefficients());
251 1892 : for (index_t i = 0; i < dataVec.size(); ++i) {
252 1890 : auto temp = dataVec[i] - resData[i] + r[i];
253 1890 : trueValue += temp - y[i] * std::log(temp);
254 1890 : trueGrad[i] = 1 - y[i] / temp;
255 1890 : trueScale[i] = y[i] / (temp * temp);
256 1890 : }
257 :
258 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
259 2 : DataContainer<TestType> dcTrueGrad(dd, trueGrad);
260 2 : REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
261 :
262 2 : auto hessian = func.getHessian(x);
263 2 : auto hx = hessian.apply(x);
264 1892 : for (index_t i = 0; i < hx.getSize(); ++i)
265 2 : REQUIRE_EQ(hx[i], dataVec[i] * trueScale[i]);
266 2 : }
267 4 : }
268 10 : }
269 10 : }
270 :
271 : TEST_SUITE_END();
|
