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             : TEST_CASE_TEMPLATE("TransmissionLogLikelihood: Testing with only data no residual", TestType, float,
      31             :                    double)
      32          10 : {
      33          10 :     using Vector = Eigen::Matrix<TestType, Eigen::Dynamic, 1>;
      34             : 
      35          10 :     GIVEN("just data (no residual)")
      36          10 :     {
      37          10 :         IndexVector_t numCoeff(3);
      38          10 :         numCoeff << 3, 7, 13;
      39          10 :         VolumeDescriptor dd(numCoeff);
      40             : 
      41          10 :         Vector y(dd.getNumberOfCoefficients());
      42          10 :         y.setRandom();
      43          10 :         DataContainer<TestType> dcY(dd, y);
      44             : 
      45          10 :         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        1366 :                 b[i] *= -1;
      51        2730 :             if (b[i] == 0)
      52           0 :                 b[i] += 1;
      53        2730 :         }
      54          10 :         DataContainer<TestType> dcB(dd, b);
      55             : 
      56          10 :         WHEN("instantiating without r")
      57          10 :         {
      58           6 :             TransmissionLogLikelihood func(dd, dcY, dcB);
      59             : 
      60           6 :             THEN("the functional is as expected")
      61           6 :             {
      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           2 :             }
      69             : 
      70           6 :             THEN("a clone behaves as expected")
      71           6 :             {
      72           2 :                 auto tmllClone = func.clone();
      73             : 
      74           2 :                 REQUIRE_NE(tmllClone.get(), &func);
      75           2 :                 REQUIRE_EQ(*tmllClone, func);
      76           2 :             }
      77             : 
      78           6 :             THEN("the evaluate, gradient and Hessian work as expected")
      79           6 :             {
      80           2 :                 Vector dataVec(dd.getNumberOfCoefficients());
      81           2 :                 dataVec.setRandom();
      82           2 :                 DataContainer<TestType> x(dd, dataVec);
      83             : 
      84             :                 // compute the "true" values
      85           2 :                 TestType trueValue = 0;
      86           2 :                 Vector trueGrad(dd.getNumberOfCoefficients());
      87           2 :                 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         546 :                 }
      94             : 
      95           2 :                 REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
      96           2 :                 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           2 :             }
     102           6 :         }
     103             : 
     104          10 :         WHEN("instantiating with r")
     105          10 :         {
     106           4 :             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         561 :                     r[i] *= -1;
     112           4 :             DataContainer<TestType> dcR(dd, r);
     113             : 
     114           4 :             TransmissionLogLikelihood func(dd, dcY, dcB, dcR);
     115             : 
     116           4 :             THEN("a clone behaves as expected")
     117           4 :             {
     118           2 :                 auto tmllClone = func.clone();
     119             : 
     120           2 :                 REQUIRE_NE(tmllClone.get(), &func);
     121           2 :                 REQUIRE_EQ(*tmllClone, func);
     122           2 :             }
     123             : 
     124           4 :             THEN("the evaluate, gradient and Hessian work as expected")
     125           4 :             {
     126           2 :                 Vector dataVec(dd.getNumberOfCoefficients());
     127           2 :                 dataVec.setRandom();
     128           2 :                 DataContainer<TestType> x(dd, dataVec);
     129             : 
     130             :                 // compute the true values
     131           2 :                 TestType trueValue = 0;
     132           2 :                 Vector trueGrad(dd.getNumberOfCoefficients());
     133           2 :                 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         546 :                 }
     141             : 
     142           2 :                 REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
     143           2 :                 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           2 :             }
     149           4 :         }
     150          10 :     }
     151          10 : }
     152             : 
     153             : TEST_CASE_TEMPLATE("TransmissionLogLikelihood: Testing with residual", TestType, float, double)
     154          10 : {
     155          10 :     using Vector = Eigen::Matrix<TestType, Eigen::Dynamic, 1>;
     156             : 
     157          10 :     GIVEN("a residual with data")
     158          10 :     {
     159          10 :         IndexVector_t numCoeff(2);
     160          10 :         numCoeff << 13, 17;
     161          10 :         VolumeDescriptor dd(numCoeff);
     162             : 
     163          10 :         Vector resData(dd.getNumberOfCoefficients());
     164          10 :         resData.setRandom();
     165          10 :         DataContainer<TestType> dcResData(dd, resData);
     166          10 :         Identity<TestType> idOp(dd);
     167          10 :         LinearResidual<TestType> linRes(idOp, dcResData);
     168             : 
     169          10 :         Vector y(dd.getNumberOfCoefficients());
     170          10 :         y.setRandom();
     171          10 :         DataContainer<TestType> dcY(dd, y);
     172             : 
     173          10 :         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        1135 :                 b[i] *= -1;
     179        2210 :             if (b[i] == 0)
     180           0 :                 b[i] += 1;
     181        2210 :         }
     182          10 :         DataContainer<TestType> dcB(dd, b);
     183             : 
     184          10 :         WHEN("instantiating without r")
     185          10 :         {
     186           6 :             TransmissionLogLikelihood func(linRes, dcY, dcB);
     187             : 
     188           6 :             THEN("the functional is as expected")
     189           6 :             {
     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           2 :             }
     196             : 
     197           6 :             THEN("a clone behaves as expected")
     198           6 :             {
     199           2 :                 auto tmllClone = func.clone();
     200             : 
     201           2 :                 REQUIRE_NE(tmllClone.get(), &func);
     202           2 :                 REQUIRE_EQ(*tmllClone, func);
     203           2 :             }
     204             : 
     205           6 :             THEN("the evaluate, gradient and Hessian work as expected")
     206           6 :             {
     207           2 :                 Vector dataVec(dd.getNumberOfCoefficients());
     208           2 :                 dataVec.setRandom();
     209           2 :                 DataContainer<TestType> x(dd, dataVec);
     210             : 
     211             :                 // compute the "true" values
     212           2 :                 TestType trueValue = 0;
     213           2 :                 Vector trueGrad(dd.getNumberOfCoefficients());
     214           2 :                 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         442 :                 }
     221             : 
     222           2 :                 REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
     223           2 :                 DataContainer<TestType> dcTrueGrad(dd, trueGrad);
     224           2 :                 REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
     225             : 
     226           2 :                 auto hessian = func.getHessian(x);
     227           2 :                 auto hx = hessian.apply(x);
     228         444 :                 for (index_t i = 0; i < hx.getSize(); ++i)
     229           2 :                     REQUIRE_EQ(hx[i], dataVec[i] * trueScale[i]);
     230           2 :             }
     231           6 :         }
     232             : 
     233          10 :         WHEN("instantiating with r")
     234          10 :         {
     235           4 :             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         429 :                     r[i] *= -1;
     241           4 :             DataContainer<TestType> dcR(dd, r);
     242             : 
     243           4 :             TransmissionLogLikelihood func(linRes, dcY, dcB, dcR);
     244             : 
     245           4 :             THEN("a clone behaves as expected")
     246           4 :             {
     247           2 :                 auto tmllClone = func.clone();
     248             : 
     249           2 :                 REQUIRE_NE(tmllClone.get(), &func);
     250           2 :                 REQUIRE_EQ(*tmllClone, func);
     251           2 :             }
     252             : 
     253           4 :             THEN("the evaluate, gradient and Hessian work as expected")
     254           4 :             {
     255           2 :                 Vector dataVec(dd.getNumberOfCoefficients());
     256           2 :                 dataVec.setRandom();
     257           2 :                 DataContainer<TestType> x(dd, dataVec);
     258             : 
     259             :                 // compute the true values
     260           2 :                 TestType trueValue = 0;
     261           2 :                 Vector trueGrad(dd.getNumberOfCoefficients());
     262           2 :                 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         442 :                 }
     270             : 
     271           2 :                 REQUIRE_UNARY(checkApproxEq(func.evaluate(x), trueValue));
     272           2 :                 DataContainer<TestType> dcTrueGrad(dd, trueGrad);
     273           2 :                 REQUIRE_UNARY(isApprox(func.getGradient(x), dcTrueGrad));
     274             : 
     275           2 :                 auto hessian = func.getHessian(x);
     276           2 :                 auto hx = hessian.apply(x);
     277         444 :                 for (index_t i = 0; i < hx.getSize(); ++i)
     278           2 :                     REQUIRE_UNARY(checkApproxEq(hx[i], dataVec[i] * trueScale[i]));
     279           2 :             }
     280           4 :         }
     281          10 :     }
     282          10 : }
     283             : 
     284             : TEST_SUITE_END();