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 22 : TEST_CASE_TEMPLATE("EmissionLogLikelihood: Testing with only data no residual", TestType, float,
28 : double)
29 : {
30 : using Vector = Eigen::Matrix<TestType, Eigen::Dynamic, 1>;
31 :
32 20 : GIVEN("just data (no residual)")
33 : {
34 20 : IndexVector_t numCoeff(3);
35 10 : numCoeff << 9, 15, 19;
36 20 : VolumeDescriptor dd(numCoeff);
37 :
38 20 : Vector y(dd.getNumberOfCoefficients());
39 10 : y.setRandom();
40 20 : DataContainer<TestType> dcY(dd, y);
41 :
42 16 : WHEN("instantiating without r")
43 : {
44 12 : EmissionLogLikelihood func(dd, dcY);
45 :
46 8 : THEN("the functional is as expected")
47 : {
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 : }
55 :
56 8 : THEN("a clone behaves as expected")
57 : {
58 4 : auto emllClone = func.clone();
59 :
60 2 : REQUIRE_NE(emllClone.get(), &func);
61 2 : REQUIRE_EQ(*emllClone, func);
62 : }
63 :
64 8 : THEN("the evaluate, gradient and Hessian work as expected")
65 : {
66 4 : 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 2644 : dataVec[i] *= -1;
71 4 : DataContainer<TestType> x(dd, dataVec);
72 :
73 : // compute the "true" values
74 2 : TestType trueValue = 0;
75 4 : Vector trueGrad(dd.getNumberOfCoefficients());
76 4 : 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 : }
83 :
84 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
85 4 : 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 : }
91 : }
92 :
93 14 : WHEN("instantiating with r")
94 : {
95 8 : 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 5184 : r[i] *= -1;
100 8 : DataContainer<TestType> dcR(dd, r);
101 :
102 8 : EmissionLogLikelihood func(dd, dcY, dcR);
103 :
104 6 : THEN("a clone behaves as expected")
105 : {
106 4 : auto emllClone = func.clone();
107 :
108 2 : REQUIRE_NE(emllClone.get(), &func);
109 2 : REQUIRE_EQ(*emllClone, func);
110 : }
111 :
112 6 : THEN("the evaluate, gradient and Hessian work as expected")
113 : {
114 4 : 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 2510 : dataVec[i] *= -1;
119 4 : DataContainer<TestType> x(dd, dataVec);
120 :
121 : // compute the "true" values
122 2 : TestType trueValue = 0;
123 4 : Vector trueGrad(dd.getNumberOfCoefficients());
124 4 : 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 : }
131 :
132 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
133 4 : 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 : }
139 : }
140 : }
141 10 : }
142 :
143 22 : TEST_CASE_TEMPLATE("EmissionLogLikelihood: Testing with residual", TestType, float, double)
144 : {
145 : using Vector = Eigen::Matrix<TestType, Eigen::Dynamic, 1>;
146 :
147 20 : GIVEN("a residual with data")
148 : {
149 20 : IndexVector_t numCoeff(3);
150 10 : numCoeff << 3, 15, 21;
151 20 : VolumeDescriptor dd(numCoeff);
152 :
153 20 : Vector resData(dd.getNumberOfCoefficients());
154 10 : resData.setRandom();
155 20 : DataContainer<TestType> dcResData(dd, resData);
156 20 : Identity<TestType> idOp(dd);
157 20 : LinearResidual<TestType> linRes(idOp, dcResData);
158 :
159 20 : Vector y(dd.getNumberOfCoefficients());
160 10 : y.setRandom();
161 20 : DataContainer<TestType> dcY(dd, y);
162 :
163 16 : WHEN("instantiating without r")
164 : {
165 12 : EmissionLogLikelihood func(linRes, dcY);
166 :
167 8 : THEN("the functional is as expected")
168 : {
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 : }
175 :
176 8 : THEN("a clone behaves as expected")
177 : {
178 4 : auto emllClone = func.clone();
179 :
180 2 : REQUIRE_NE(emllClone.get(), &func);
181 2 : REQUIRE_EQ(*emllClone, func);
182 : }
183 :
184 8 : THEN("the evaluate, gradient and Hessian work as expected")
185 : {
186 4 : 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 998 : dataVec[i] -= 2 * (dataVec[i] - resData[i]);
192 : }
193 4 : DataContainer<TestType> x(dd, dataVec);
194 :
195 : // compute the "true" values
196 2 : TestType trueValue = 0;
197 4 : Vector trueGrad(dd.getNumberOfCoefficients());
198 4 : 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 : }
205 :
206 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
207 4 : DataContainer<TestType> dcTrueGrad(dd, trueGrad);
208 2 : REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
209 :
210 4 : auto hessian = func.getHessian(x);
211 4 : auto hx = hessian.apply(x);
212 1892 : for (index_t i = 0; i < hx.getSize(); ++i)
213 1890 : REQUIRE_UNARY(checkApproxEq(hx[i], dataVec[i] * trueScale[i]));
214 : }
215 : }
216 :
217 14 : WHEN("instantiating with r")
218 : {
219 8 : 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 1966 : r[i] *= -1;
224 8 : DataContainer<TestType> dcR(dd, r);
225 :
226 8 : EmissionLogLikelihood func(linRes, dcY, dcR);
227 :
228 6 : THEN("a clone behaves as expected")
229 : {
230 4 : auto emllClone = func.clone();
231 :
232 2 : REQUIRE_NE(emllClone.get(), &func);
233 2 : REQUIRE_EQ(*emllClone, func);
234 : }
235 :
236 6 : THEN("the evaluate, gradient and Hessian work as expected")
237 : {
238 4 : 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 950 : dataVec[i] -= 2 * (dataVec[i] - resData[i]);
244 : }
245 4 : DataContainer<TestType> x(dd, dataVec);
246 :
247 : // compute the "true" values
248 2 : TestType trueValue = 0;
249 4 : Vector trueGrad(dd.getNumberOfCoefficients());
250 4 : 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 : }
257 :
258 2 : REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
259 4 : DataContainer<TestType> dcTrueGrad(dd, trueGrad);
260 2 : REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
261 :
262 4 : auto hessian = func.getHessian(x);
263 4 : auto hx = hessian.apply(x);
264 1892 : for (index_t i = 0; i < hx.getSize(); ++i)
265 1890 : REQUIRE_EQ(hx[i], dataVec[i] * trueScale[i]);
266 : }
267 : }
268 : }
269 10 : }
270 :
271 : TEST_SUITE_END();
|
