Fraxinus  16.5.0-fx-rc9
An IGT application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
cxtestSyntheticVolumeComparer.cpp
Go to the documentation of this file.
1 /*=========================================================================
2 This file is part of CustusX, an Image Guided Therapy Application.
3 
4 Copyright (c) 2008-2014, SINTEF Department of Medical Technology
5 All rights reserved.
6 
7 Redistribution and use in source and binary forms, with or without
8 modification, are permitted provided that the following conditions are met:
9 
10 1. Redistributions of source code must retain the above copyright notice,
11  this list of conditions and the following disclaimer.
12 
13 2. Redistributions in binary form must reproduce the above copyright notice,
14  this list of conditions and the following disclaimer in the documentation
15  and/or other materials provided with the distribution.
16 
17 3. Neither the name of the copyright holder nor the names of its contributors
18  may be used to endorse or promote products derived from this software
19  without specific prior written permission.
20 
21 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
22 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
23 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
24 DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
25 FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
26 DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
27 SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28 CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
29 OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
30 OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 =========================================================================*/
32 
34 #include "cxVolumeHelpers.h"
35 #include "catch.hpp"
36 #include "cxImage.h"
37 #include "cxtestUtilities.h"
38 #include <vtkImageData.h>
40 #include "cxDataReaderWriter.h"
41 #include <QDir>
42 #include <QString>
43 #include "cxTypeConversions.h"
44 
45 namespace cxtest
46 {
47 
49 {
50 }
51 
53 {
54  mPhantom = phantom;
55  mNominalImage.reset();
56 }
57 
59 {
60  mTestImage = image;
61  mNominalImage.reset();
62 }
63 
64 //vtkImageDataPtr Utilities::create3DVtkImageData(Eigen::Array3i dim, const unsigned int voxelValue)
65 //{
66 // return cx::generateVtkImageData(dim, cx::Vector3D(1,1,1), voxelValue);
67 //}
68 
69 //cx::ImagePtr Utilities::create3DImage(Eigen::Array3i dim, const unsigned int voxelValue)
70 //{
71 // vtkImageDataPtr vtkImageData = create3DVtkImageData(dim, voxelValue);
72 // QString unique_string = qstring_cast(reinterpret_cast<long>(vtkImageData.GetPointer()));
73 // QString imagesUid = QString("TESTUID_%2_%1").arg(unique_string);
74 // cx::ImagePtr image(new cx::Image(imagesUid, vtkImageData));
75 
76 // return image;
77 //}
78 
79 
80 cx::ImagePtr SyntheticVolumeComparer::getNominalOutputImage() const
81 {
82  REQUIRE(mTestImage);
83 
84  if (!mNominalImage)
85  {
86  Eigen::Array3i dim(mTestImage->getBaseVtkImageData()->GetDimensions());
87  cx::Vector3D spacing(mTestImage->getBaseVtkImageData()->GetSpacing());
88  mNominalImage = cxtest::Utilities::create3DImage(dim, spacing, 0);
89  mNominalImage->get_rMd_History()->setRegistration(mTestImage->get_rMd());
90 
91  mPhantom->fillVolume(mNominalImage);
92  }
93 
94  return mNominalImage;
95 }
96 
98 {
99  float sse = this->getRMS();
100  if (this->getVerbose())
101  std::cout << "RMS value: " << sse << std::endl;
102  CHECK(sse < threshold);
103 }
104 
105 double SyntheticVolumeComparer::getRMS() const
106 {
107  double sse = cx::calculateRMSError(mTestImage->getBaseVtkImageData(), this->getNominalOutputImage()->getBaseVtkImageData());
108 // float sse = mPhantom->computeRMSError(mOutputData);
109 // std::cout << "RMS value: " << sse << std::endl;
110  return sse;
111 }
112 
114 {
115  cx::Vector3D c_n = calculateCentroid(this->getNominalOutputImage());
116  cx::Vector3D c_r = calculateCentroid(mTestImage);
117 
118  double difference = (c_n-c_r).norm();
119 
120  if (this->getVerbose())
121  std::cout << "c_n=[" << c_n << "] c_r=[" << c_r << "] diff=[" << difference << "]" << std::endl;
122 
123  CHECK(difference < val);
124 }
125 
127 {
128  double v_n = calculateMass(this->getNominalOutputImage());
129  double v_r = calculateMass(mTestImage);
130  double normalized_difference = fabs(v_n-v_r)/(v_n+v_r);
131 
132  if (this->getVerbose())
133  std::cout << "v_n=[" << v_n << "] v_r=[" << v_r << "] diff=[" << normalized_difference << "]" << std::endl;
134 
135  CHECK(normalized_difference<val);
136 }
137 
138 void SyntheticVolumeComparer::checkValueWithin(cx::Vector3D p_r, int lowerLimit, int upperLimit)
139 {
140  int val = this->getValue(mTestImage, p_r);
141  if (this->getVerbose())
142  std::cout << QString("p_r[%1]=%2, limit=[%3, %4]").arg(qstring_cast(p_r)).arg(val).arg(lowerLimit).arg(upperLimit) << std::endl;
143  CHECK(val>=lowerLimit);
144  CHECK(val<=upperLimit);
145 }
146 
147 double SyntheticVolumeComparer::getValue(cx::ImagePtr image, cx::Vector3D p_r)
148 {
149  vtkImageDataPtr raw = image->getBaseVtkImageData();
150  cx::Vector3D p_d = image->get_rMd().inv().coord(p_r);
151  Eigen::Array3d spacing(raw->GetSpacing());
152  Eigen::Array3i voxel = cx::round(p_r.array() / spacing).cast<int>();
153  double val = raw->GetScalarComponentAsDouble(voxel.data()[0],voxel.data()[1],voxel.data()[2], 0);
154  return val;
155 }
156 
158 {
159  cx::MetaImageReader().saveImage(this->getNominalOutputImage(), this->addFullPath(filename));
160 }
161 
163 {
164  cx::MetaImageReader().saveImage(mTestImage, this->addFullPath(filename));
165 }
166 
167 QString SyntheticVolumeComparer::addFullPath(QString filename)
168 {
169  QDir dir(QDir::homePath() + "/test/");
170  dir.mkdir(".");
171  QString fullFilename = dir.absoluteFilePath(filename);
172  return fullFilename;
173 }
174 
175 } // namespace cxtest
176 
177 
QString qstring_cast(const T &val)
double calculateMass(cx::ImagePtr image)
boost::shared_ptr< class Image > ImagePtr
Definition: cxDicomWidget.h:48
static cx::ImagePtr create3DImage(Eigen::Array3i dim=Eigen::Array3i(3, 3, 3), const unsigned int voxelValue=100)
void checkValueWithin(cx::Vector3D p_r, int lowerLimit, int upperLimit)
boost::shared_ptr< cxSyntheticVolume > cxSyntheticVolumePtr
cx::Vector3D calculateCentroid(cx::ImagePtr image)
void setPhantom(cx::cxSyntheticVolumePtr phantom)
double calculateRMSError(vtkImageDataPtr a, vtkImageDataPtr b)
void saveImage(ImagePtr image, const QString &filename)
Eigen::Vector3d Vector3D
Vector3D is a representation of a point or vector in 3D.
Definition: cxVector3D.h:63
vtkSmartPointer< class vtkImageData > vtkImageDataPtr
Reader for metaheader .mhd files.
Vector3D round(const Vector3D &a)
Definition: cxVector3D.cpp:96