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