diff --git a/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidgetControls.ui b/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidgetControls.ui
index 6daf485458..659f57055c 100644
--- a/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidgetControls.ui
+++ b/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidgetControls.ui
@@ -1,463 +1,463 @@
QmitkInteractiveTransformationWidgetControls
0
0
435
381
Form
-
Interactive Translation
-
-
-
50
false
x-Direction
(Frontal):
false
-
0
0
ArrowCursor
-500
500
1
Qt::Horizontal
false
false
QSlider::TicksAbove
100
-
- 1
+ 4
-500.000000000000000
500.000000000000000
-
-
50
false
y-Direction
(Sagittal):
false
-
0
0
-500
500
Qt::Horizontal
false
QSlider::TicksAbove
100
-
- 1
+ 4
-500.000000000000000
500.000000000000000
-
-
50
false
z-Direction
(Transversal):
false
-
0
0
-500
500
0
Qt::Horizontal
false
QSlider::TicksAbove
100
-
- 1
+ 4
-500.000000000000000
500.000000000000000
-
Interactive Rotation
-
-
-
50
false
x-Axis
(Frontal):
false
-
0
0
ArrowCursor
-180
180
Qt::Horizontal
QSlider::TicksAbove
45
-
- 1
+ 4
-180.000000000000000
180.000000000000000
-
-
50
false
y-Axis
(Sagittal):
false
-
0
0
-180
180
Qt::Horizontal
QSlider::TicksAbove
45
-
- 1
+ 4
-180.000000000000000
180.000000000000000
-
-
50
false
z-Axis
(Transversal):
false
-
0
0
-180
180
Qt::Horizontal
QSlider::TicksAbove
45
-
- 1
+ 4
-180.000000000000000
180.000000000000000
-
-
Use Manipulated ToolTip
-
Reset To Identity
-
Revert Changes
-
Cancel
diff --git a/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/NavigationStepWidgets/QmitkUSNavigationStepCtUsRegistration.cpp b/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/NavigationStepWidgets/QmitkUSNavigationStepCtUsRegistration.cpp
index cdfe37e46f..6f4da18cab 100644
--- a/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/NavigationStepWidgets/QmitkUSNavigationStepCtUsRegistration.cpp
+++ b/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/NavigationStepWidgets/QmitkUSNavigationStepCtUsRegistration.cpp
@@ -1,2331 +1,2369 @@
/*===================================================================
The Medical Imaging Interaction Toolkit (MITK)
Copyright (c) German Cancer Research Center,
Division of Medical and Biological Informatics.
All rights reserved.
This software is distributed WITHOUT ANY WARRANTY; without
even the implied warranty of MERCHANTABILITY or FITNESS FOR
A PARTICULAR PURPOSE.
See LICENSE.txt or http://www.mitk.org for details.
===================================================================*/
#include "QmitkUSNavigationStepCtUsRegistration.h"
#include "ui_QmitkUSNavigationStepCtUsRegistration.h"
#include
#include "mitkNodeDisplacementFilter.h"
#include "../QmitkUSNavigationMarkerPlacement.h"
#include
#include
#include
#include "mitkProperties.h"
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
static const int NUMBER_FIDUCIALS_NEEDED = 8;
QmitkUSNavigationStepCtUsRegistration::QmitkUSNavigationStepCtUsRegistration(QWidget *parent) :
QmitkUSAbstractNavigationStep(parent),
ui(new Ui::QmitkUSNavigationStepCtUsRegistration),
m_PerformingGroundTruthProtocolEvaluation(false)
{
this->UnsetFloatingImageGeometry();
this->DefineDataStorageImageFilter();
this->CreateQtPartControl(this);
+ this->InitializeTransformationSensorCSToMarkerCS();
}
QmitkUSNavigationStepCtUsRegistration::~QmitkUSNavigationStepCtUsRegistration()
{
delete ui;
}
bool QmitkUSNavigationStepCtUsRegistration::OnStartStep()
{
MITK_INFO << "OnStartStep()";
return true;
}
bool QmitkUSNavigationStepCtUsRegistration::OnStopStep()
{
MITK_INFO << "OnStopStep()";
return true;
}
bool QmitkUSNavigationStepCtUsRegistration::OnFinishStep()
{
MITK_INFO << "OnFinishStep()";
return true;
}
bool QmitkUSNavigationStepCtUsRegistration::OnActivateStep()
{
MITK_INFO << "OnActivateStep()";
ui->floatingImageComboBox->SetDataStorage(this->GetDataStorage());
ui->groundTruthImageComboBox->SetDataStorage(this->GetDataStorage());
ui->ctImagesToChooseComboBox->SetDataStorage(this->GetDataStorage());
return true;
}
bool QmitkUSNavigationStepCtUsRegistration::OnDeactivateStep()
{
MITK_INFO << "OnDeactivateStep()";
return true;
}
void QmitkUSNavigationStepCtUsRegistration::OnUpdate()
{
if (m_NavigationDataSource.IsNull()) { return; }
m_NavigationDataSource->Update();
}
void QmitkUSNavigationStepCtUsRegistration::OnSettingsChanged(const itk::SmartPointer settingsNode)
{
}
QString QmitkUSNavigationStepCtUsRegistration::GetTitle()
{
return "CT-to-US registration";
}
QmitkUSAbstractNavigationStep::FilterVector QmitkUSNavigationStepCtUsRegistration::GetFilter()
{
return FilterVector();
}
void QmitkUSNavigationStepCtUsRegistration::OnSetCombinedModality()
{
mitk::AbstractUltrasoundTrackerDevice::Pointer combinedModality = this->GetCombinedModality(false);
if (combinedModality.IsNotNull())
{
m_NavigationDataSource = combinedModality->GetNavigationDataSource();
}
}
void QmitkUSNavigationStepCtUsRegistration::UnsetFloatingImageGeometry()
{
m_ImageDimension[0] = 0;
m_ImageDimension[1] = 0;
m_ImageDimension[2] = 0;
m_ImageSpacing[0] = 1;
m_ImageSpacing[1] = 1;
m_ImageSpacing[2] = 1;
}
void QmitkUSNavigationStepCtUsRegistration::SetFloatingImageGeometryInformation(mitk::Image * image)
{
m_ImageDimension[0] = image->GetDimension(0);
m_ImageDimension[1] = image->GetDimension(1);
m_ImageDimension[2] = image->GetDimension(2);
m_ImageSpacing[0] = image->GetGeometry()->GetSpacing()[0];
m_ImageSpacing[1] = image->GetGeometry()->GetSpacing()[1];
m_ImageSpacing[2] = image->GetGeometry()->GetSpacing()[2];
}
double QmitkUSNavigationStepCtUsRegistration::GetVoxelVolume()
{
if (m_FloatingImage.IsNull())
{
return 0.0;
}
MITK_INFO << "ImageSpacing = " << m_ImageSpacing;
return m_ImageSpacing[0] * m_ImageSpacing[1] * m_ImageSpacing[2];
}
double QmitkUSNavigationStepCtUsRegistration::GetFiducialVolume(double radius)
{
return 1.333333333 * 3.141592 * (radius * radius * radius);
}
bool QmitkUSNavigationStepCtUsRegistration::FilterFloatingImage()
{
if (m_FloatingImage.IsNull())
{
return false;
}
ImageType::Pointer itkImage1 = ImageType::New();
mitk::CastToItkImage(m_FloatingImage, itkImage1);
this->InitializeImageFilters();
m_ThresholdFilter->SetInput(itkImage1);
m_LaplacianFilter1->SetInput(m_ThresholdFilter->GetOutput());
m_LaplacianFilter2->SetInput(m_LaplacianFilter1->GetOutput());
m_BinaryThresholdFilter->SetInput(m_LaplacianFilter2->GetOutput());
m_HoleFillingFilter->SetInput(m_BinaryThresholdFilter->GetOutput());
m_BinaryImageToShapeLabelMapFilter->SetInput(m_HoleFillingFilter->GetOutput());
m_BinaryImageToShapeLabelMapFilter->Update();
ImageType::Pointer binaryImage = ImageType::New();
binaryImage = m_HoleFillingFilter->GetOutput();
this->EliminateTooSmallLabeledObjects(binaryImage);
//mitk::CastToMitkImage(binaryImage, m_FloatingImage);
return true;
}
void QmitkUSNavigationStepCtUsRegistration::InitializeImageFilters()
{
//Initialize threshold filters
m_ThresholdFilter = itk::ThresholdImageFilter::New();
m_ThresholdFilter->SetOutsideValue(0);
m_ThresholdFilter->SetLower(500);
m_ThresholdFilter->SetUpper(3200);
//Initialize binary threshold filter 1
m_BinaryThresholdFilter = BinaryThresholdImageFilterType::New();
m_BinaryThresholdFilter->SetOutsideValue(0);
m_BinaryThresholdFilter->SetInsideValue(1);
m_BinaryThresholdFilter->SetLowerThreshold(350);
m_BinaryThresholdFilter->SetUpperThreshold(10000);
//Initialize laplacian recursive gaussian image filter
m_LaplacianFilter1 = LaplacianRecursiveGaussianImageFilterType::New();
m_LaplacianFilter2 = LaplacianRecursiveGaussianImageFilterType::New();
//Initialize binary hole filling filter
m_HoleFillingFilter = VotingBinaryIterativeHoleFillingImageFilterType::New();
VotingBinaryIterativeHoleFillingImageFilterType::InputSizeType radius;
radius.Fill(1);
m_HoleFillingFilter->SetRadius(radius);
m_HoleFillingFilter->SetBackgroundValue(0);
m_HoleFillingFilter->SetForegroundValue(1);
m_HoleFillingFilter->SetMaximumNumberOfIterations(5);
//Initialize binary image to shape label map filter
m_BinaryImageToShapeLabelMapFilter = BinaryImageToShapeLabelMapFilterType::New();
m_BinaryImageToShapeLabelMapFilter->SetInputForegroundValue(1);
}
double QmitkUSNavigationStepCtUsRegistration::GetCharacteristicDistanceAWithUpperMargin()
{
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
return 12.07;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
return 18.105;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
return 24.14;
}
return 0.0;
}
double QmitkUSNavigationStepCtUsRegistration::GetCharacteristicDistanceBWithLowerMargin()
{
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
return 12.07;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
return 18.105;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
return 24.14;
}
return 0.0;
}
double QmitkUSNavigationStepCtUsRegistration::GetCharacteristicDistanceBWithUpperMargin()
{
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
return 15.73;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
return 23.595;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
return 31.46;
}
return 0.0;
}
double QmitkUSNavigationStepCtUsRegistration::GetMinimalFiducialConfigurationDistance()
{
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
return 10.0;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
return 15.0;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
return 20.0;
}
return 0.0;
}
void QmitkUSNavigationStepCtUsRegistration::CreateMarkerModelCoordinateSystemPointSet()
{
if (m_MarkerModelCoordinateSystemPointSet.IsNull())
{
m_MarkerModelCoordinateSystemPointSet = mitk::PointSet::New();
}
else
{
m_MarkerModelCoordinateSystemPointSet->Clear();
}
mitk::Point3D fiducial1;
mitk::Point3D fiducial2;
mitk::Point3D fiducial3;
mitk::Point3D fiducial4;
mitk::Point3D fiducial5;
mitk::Point3D fiducial6;
mitk::Point3D fiducial7;
mitk::Point3D fiducial8;
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
fiducial1[0] = 0;
fiducial1[1] = 0;
fiducial1[2] = 0;
fiducial2[0] = 0;
fiducial2[1] = 10;
fiducial2[2] = 0;
fiducial3[0] = 10;
fiducial3[1] = 0;
fiducial3[2] = 0;
fiducial4[0] = 20;
fiducial4[1] = 20;
fiducial4[2] = 0;
fiducial5[0] = 0;
fiducial5[1] = 20;
fiducial5[2] = 10;
fiducial6[0] = 10;
fiducial6[1] = 20;
fiducial6[2] = 10;
fiducial7[0] = 20;
fiducial7[1] = 10;
fiducial7[2] = 10;
fiducial8[0] = 20;
fiducial8[1] = 0;
fiducial8[2] = 10;
break;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
fiducial1[0] = 0;
fiducial1[1] = 0;
fiducial1[2] = 0;
fiducial2[0] = 0;
fiducial2[1] = 15;
fiducial2[2] = 0;
fiducial3[0] = 15;
fiducial3[1] = 0;
fiducial3[2] = 0;
fiducial4[0] = 30;
fiducial4[1] = 30;
fiducial4[2] = 0;
fiducial5[0] = 0;
fiducial5[1] = 30;
fiducial5[2] = 15;
fiducial6[0] = 15;
fiducial6[1] = 30;
fiducial6[2] = 15;
fiducial7[0] = 30;
fiducial7[1] = 15;
fiducial7[2] = 15;
fiducial8[0] = 30;
fiducial8[1] = 0;
fiducial8[2] = 15;
break;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
fiducial1[0] = 0;
fiducial1[1] = 0;
fiducial1[2] = 0;
fiducial2[0] = 0;
fiducial2[1] = 20;
fiducial2[2] = 0;
fiducial3[0] = 20;
fiducial3[1] = 0;
fiducial3[2] = 0;
fiducial4[0] = 40;
fiducial4[1] = 40;
fiducial4[2] = 0;
fiducial5[0] = 0;
fiducial5[1] = 40;
fiducial5[2] = 20;
fiducial6[0] = 20;
fiducial6[1] = 40;
fiducial6[2] = 20;
fiducial7[0] = 40;
fiducial7[1] = 20;
fiducial7[2] = 20;
fiducial8[0] = 40;
fiducial8[1] = 0;
fiducial8[2] = 20;
break;
}
m_MarkerModelCoordinateSystemPointSet->InsertPoint(0, fiducial1);
m_MarkerModelCoordinateSystemPointSet->InsertPoint(1, fiducial2);
m_MarkerModelCoordinateSystemPointSet->InsertPoint(2, fiducial3);
m_MarkerModelCoordinateSystemPointSet->InsertPoint(3, fiducial4);
m_MarkerModelCoordinateSystemPointSet->InsertPoint(4, fiducial5);
m_MarkerModelCoordinateSystemPointSet->InsertPoint(5, fiducial6);
m_MarkerModelCoordinateSystemPointSet->InsertPoint(6, fiducial7);
m_MarkerModelCoordinateSystemPointSet->InsertPoint(7, fiducial8);
/*mitk::DataNode::Pointer node = this->GetDataStorage()->GetNamedNode("Marker Model Coordinate System Point Set");
if (node == nullptr)
{
node = mitk::DataNode::New();
node->SetName("Marker Model Coordinate System Point Set");
node->SetData(m_MarkerModelCoordinateSystemPointSet);
this->GetDataStorage()->Add(node);
}
else
{
node->SetData(m_MarkerModelCoordinateSystemPointSet);
this->GetDataStorage()->Modified();
}*/
}
void QmitkUSNavigationStepCtUsRegistration::InitializeTransformationSensorCSToMarkerCS()
{
- //At first calculate the inverse transformation (marker-CS to sensor-CS)
+ //The following calculations are related to the 3mm | 15mm fiducial configuration
+ m_TransformSensorCSToMarkerCS = mitk::AffineTransform3D::New();
+
mitk::Vector3D translation;
+ translation[0] = -18.175;
+ translation[1] = 15.000;
+ translation[2] = 8.001; // considering a distance from the base plate of 0.315 inch and not 0.313 inch
+
+ m_TransformSensorCSToMarkerCS->SetOffset(translation);
+
+ // Quaternion (x, y, z, r) --> n = (1,0,0) --> q(sin(90°),0,0,cos(90°))
+ mitk::Quaternion q1(1, 0, 0, 0); // corresponding to a rotation of 180° around the normal x-axis.
+ // .transpose() is needed for changing the rows and the columns of the returned rotation_matrix_transpose
+ vnl_matrix_fixed vnl_rotation = q1.rotation_matrix_transpose().transpose(); // :-)
+ mitk::Matrix3D rotationMatrix;
+
+ for (int i = 0; i < 3; ++i) {
+ for (int j = 0; j < 3; ++j) {
+ rotationMatrix[i][j] = vnl_rotation[i][j];
+ }
+ }
+
+ m_TransformSensorCSToMarkerCS->SetMatrix(rotationMatrix);
+ //The transformation from the sensor-CS to the marker-CS is calculated now.
+ MITK_INFO << "TransformSensorCSToMarkerCS = " << m_TransformSensorCSToMarkerCS;
+
+ /*mitk::NavigationData::Pointer navData = mitk::NavigationData::New();
+ navData->SetOrientation(q1);
+ navData->SetPosition(translation);
+ navData->SetHasOrientation(true);
+ navData->SetHasPosition(true);*/
}
void QmitkUSNavigationStepCtUsRegistration::InitializePointsToTransformForGroundTruthProtocol()
{
m_PointsToTransformGroundTruthProtocol.clear();
mitk::Point3D point0mm;
mitk::Point3D point20mm;
mitk::Point3D point40mm;
mitk::Point3D point60mm;
mitk::Point3D point80mm;
mitk::Point3D point100mm;
point0mm[0] = 0.0;
point0mm[1] = 0.0;
point0mm[2] = 0.0;
point20mm[0] = 0.0;
point20mm[1] = 0.0;
point20mm[2] = 0.0;
point40mm[0] = 0.0;
point40mm[1] = 0.0;
point40mm[2] = 0.0;
point60mm[0] = 0.0;
point60mm[1] = 0.0;
point60mm[2] = 0.0;
point80mm[0] = 0.0;
point80mm[1] = 0.0;
point80mm[2] = 0.0;
point100mm[0] = 0.0;
point100mm[1] = 0.0;
point100mm[2] = 0.0;
m_PointsToTransformGroundTruthProtocol.insert(std::pair(0, point0mm));
m_PointsToTransformGroundTruthProtocol.insert(std::pair(20, point20mm));
m_PointsToTransformGroundTruthProtocol.insert(std::pair(40, point40mm));
m_PointsToTransformGroundTruthProtocol.insert(std::pair(60, point60mm));
m_PointsToTransformGroundTruthProtocol.insert(std::pair(80, point80mm));
m_PointsToTransformGroundTruthProtocol.insert(std::pair(100, point100mm));
}
void QmitkUSNavigationStepCtUsRegistration::CreatePointsToTransformForGroundTruthProtocol()
{
this->InitializePointsToTransformForGroundTruthProtocol();
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
MITK_WARN << "For this marker configuration (10mm) there does not exist a point to transform.";
break;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
m_PointsToTransformGroundTruthProtocol.at(0)[0] = 130; // = 30mm to end of clipping plate + 100 mm to middle axis of measurement plate
m_PointsToTransformGroundTruthProtocol.at(0)[1] = 15;
m_PointsToTransformGroundTruthProtocol.at(0)[2] = -7; // = 5mm distance to clipping plate + 2mm to base
m_PointsToTransformGroundTruthProtocol.at(20)[0] = 130;
m_PointsToTransformGroundTruthProtocol.at(20)[1] = 15;
m_PointsToTransformGroundTruthProtocol.at(20)[2] = -27; // = 5mm distance to clipping plate + 2mm to base + 20mm depth
m_PointsToTransformGroundTruthProtocol.at(40)[0] = 130;
m_PointsToTransformGroundTruthProtocol.at(40)[1] = 15;
m_PointsToTransformGroundTruthProtocol.at(40)[2] = -47; // = 5mm distance to clipping plate + 2mm to base + 40mm depth
m_PointsToTransformGroundTruthProtocol.at(60)[0] = 130;
m_PointsToTransformGroundTruthProtocol.at(60)[1] = 15;
m_PointsToTransformGroundTruthProtocol.at(60)[2] = -67; // = 5mm distance to clipping plate + 2mm to base + 60mm depth
m_PointsToTransformGroundTruthProtocol.at(80)[0] = 130;
m_PointsToTransformGroundTruthProtocol.at(80)[1] = 15;
m_PointsToTransformGroundTruthProtocol.at(80)[2] = -87; // = 5mm distance to clipping plate + 2mm to base + 80mm depth
m_PointsToTransformGroundTruthProtocol.at(100)[0] = 130;
m_PointsToTransformGroundTruthProtocol.at(100)[1] = 15;
m_PointsToTransformGroundTruthProtocol.at(100)[2] = -107; // = 5mm distance to clipping plate + 2mm to base + 100mm depth
break;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
m_PointsToTransformGroundTruthProtocol.at(0)[0] = 135; // = 20 + 15mm to end of clipping plate + 100 mm to middle axis of measurement plate
m_PointsToTransformGroundTruthProtocol.at(0)[1] = 20;
m_PointsToTransformGroundTruthProtocol.at(0)[2] = -9; // = 7mm distance to clipping plate + 2mm to base
m_PointsToTransformGroundTruthProtocol.at(20)[0] = 135;
m_PointsToTransformGroundTruthProtocol.at(20)[1] = 20;
m_PointsToTransformGroundTruthProtocol.at(20)[2] = -29; // = 7mm distance to clipping plate + 2mm to base + 20mm depth
m_PointsToTransformGroundTruthProtocol.at(40)[0] = 135;
m_PointsToTransformGroundTruthProtocol.at(40)[1] = 20;
m_PointsToTransformGroundTruthProtocol.at(40)[2] = -49; // = 7mm distance to clipping plate + 2mm to base + 40mm depth
m_PointsToTransformGroundTruthProtocol.at(60)[0] = 135;
m_PointsToTransformGroundTruthProtocol.at(60)[1] = 20;
m_PointsToTransformGroundTruthProtocol.at(60)[2] = -69; // = 7mm distance to clipping plate + 2mm to base + 60mm depth
m_PointsToTransformGroundTruthProtocol.at(80)[0] = 135;
m_PointsToTransformGroundTruthProtocol.at(80)[1] = 20;
m_PointsToTransformGroundTruthProtocol.at(80)[2] = -89; // = 7mm distance to clipping plate + 2mm to base + 80mm depth
m_PointsToTransformGroundTruthProtocol.at(100)[0] = 135;
m_PointsToTransformGroundTruthProtocol.at(100)[1] = 20;
m_PointsToTransformGroundTruthProtocol.at(100)[2] = -109; // = 7mm distance to clipping plate + 2mm to base + 100mm depth
break;
}
}
void QmitkUSNavigationStepCtUsRegistration::TransformPointsGroundTruthProtocol()
{
if (m_GroundTruthProtocolTransformedPoints.find(0) == m_GroundTruthProtocolTransformedPoints.end())
{
mitk::PointSet::Pointer pointSet = mitk::PointSet::New();
pointSet->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(0)));
m_GroundTruthProtocolTransformedPoints.insert(std::pair(0, pointSet));
}
else
{
m_GroundTruthProtocolTransformedPoints.at(0)->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(0)));
}
if (m_GroundTruthProtocolTransformedPoints.find(20) == m_GroundTruthProtocolTransformedPoints.end())
{
mitk::PointSet::Pointer pointSet = mitk::PointSet::New();
pointSet->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(20)));
m_GroundTruthProtocolTransformedPoints.insert(std::pair(20, pointSet));
}
else
{
m_GroundTruthProtocolTransformedPoints.at(20)->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(20)));
}
if (m_GroundTruthProtocolTransformedPoints.find(40) == m_GroundTruthProtocolTransformedPoints.end())
{
mitk::PointSet::Pointer pointSet = mitk::PointSet::New();
pointSet->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(40)));
m_GroundTruthProtocolTransformedPoints.insert(std::pair(40, pointSet));
}
else
{
m_GroundTruthProtocolTransformedPoints.at(40)->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(40)));
}
if (m_GroundTruthProtocolTransformedPoints.find(60) == m_GroundTruthProtocolTransformedPoints.end())
{
mitk::PointSet::Pointer pointSet = mitk::PointSet::New();
pointSet->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(60)));
m_GroundTruthProtocolTransformedPoints.insert(std::pair(60, pointSet));
}
else
{
m_GroundTruthProtocolTransformedPoints.at(60)->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(60)));
}
if (m_GroundTruthProtocolTransformedPoints.find(80) == m_GroundTruthProtocolTransformedPoints.end())
{
mitk::PointSet::Pointer pointSet = mitk::PointSet::New();
pointSet->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(80)));
m_GroundTruthProtocolTransformedPoints.insert(std::pair(80, pointSet));
}
else
{
m_GroundTruthProtocolTransformedPoints.at(80)->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(80)));
}
if (m_GroundTruthProtocolTransformedPoints.find(100) == m_GroundTruthProtocolTransformedPoints.end())
{
mitk::PointSet::Pointer pointSet = mitk::PointSet::New();
pointSet->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(100)));
m_GroundTruthProtocolTransformedPoints.insert(std::pair(100, pointSet));
}
else
{
m_GroundTruthProtocolTransformedPoints.at(100)->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(m_PointsToTransformGroundTruthProtocol.at(100)));
}
}
void QmitkUSNavigationStepCtUsRegistration::AddTransformedPointsToDataStorage()
{
if (m_GroundTruthProtocolTransformedPoints.find(0) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(20) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(40) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(60) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(80) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(100) == m_GroundTruthProtocolTransformedPoints.end())
{
QMessageBox msgBox;
msgBox.setText("Cannot add transformed Points to DataStorage because they do not exist.\
Stopping evaluation the protocol.");
msgBox.exec();
return;
}
std::string nameNode0mm = "GroundTruthProt-Depth0mm";
std::string nameNode20mm = "GroundTruthProt-Depth20mm";
std::string nameNode40mm = "GroundTruthProt-Depth40mm";
std::string nameNode60mm = "GroundTruthProt-Depth60mm";
std::string nameNode80mm = "GroundTruthProt-Depth80mm";
std::string nameNode100mm = "GroundTruthProt-Depth100mm";
//Add transformed points of depth 0mm to the data storage
mitk::DataNode::Pointer node0mm = this->GetDataStorage()->GetNamedNode(nameNode0mm);
if (node0mm.IsNull())
{
node0mm = mitk::DataNode::New();
node0mm->SetName(nameNode0mm);
node0mm->SetData(m_GroundTruthProtocolTransformedPoints.at(0));
this->GetDataStorage()->Add(node0mm);
}
else
{
node0mm->SetData(m_GroundTruthProtocolTransformedPoints.at(0));
this->GetDataStorage()->Modified();
}
if(ui->protocolEvaluationTypeComboBox->currentText().compare("PLANE") == 0 )
{
//Add transformed points of depth 20mm to the data storage
mitk::DataNode::Pointer node20mm = this->GetDataStorage()->GetNamedNode(nameNode20mm);
if (node20mm.IsNull())
{
node20mm = mitk::DataNode::New();
node20mm->SetName(nameNode20mm);
node20mm->SetData(m_GroundTruthProtocolTransformedPoints.at(20));
this->GetDataStorage()->Add(node20mm);
}
else
{
node20mm->SetData(m_GroundTruthProtocolTransformedPoints.at(20));
this->GetDataStorage()->Modified();
}
//Add transformed points of depth 40mm to the data storage
mitk::DataNode::Pointer node40mm = this->GetDataStorage()->GetNamedNode(nameNode40mm);
if (node40mm.IsNull())
{
node40mm = mitk::DataNode::New();
node40mm->SetName(nameNode40mm);
node40mm->SetData(m_GroundTruthProtocolTransformedPoints.at(40));
this->GetDataStorage()->Add(node40mm);
}
else
{
node40mm->SetData(m_GroundTruthProtocolTransformedPoints.at(40));
this->GetDataStorage()->Modified();
}
//Add transformed points of depth 60mm to the data storage
mitk::DataNode::Pointer node60mm = this->GetDataStorage()->GetNamedNode(nameNode60mm);
if (node60mm.IsNull())
{
node60mm = mitk::DataNode::New();
node60mm->SetName(nameNode60mm);
node60mm->SetData(m_GroundTruthProtocolTransformedPoints.at(60));
this->GetDataStorage()->Add(node60mm);
}
else
{
node60mm->SetData(m_GroundTruthProtocolTransformedPoints.at(60));
this->GetDataStorage()->Modified();
}
//Add transformed points of depth 80mm to the data storage
mitk::DataNode::Pointer node80mm = this->GetDataStorage()->GetNamedNode(nameNode80mm);
if (node80mm.IsNull())
{
node80mm = mitk::DataNode::New();
node80mm->SetName(nameNode80mm);
node80mm->SetData(m_GroundTruthProtocolTransformedPoints.at(80));
this->GetDataStorage()->Add(node80mm);
}
else
{
node80mm->SetData(m_GroundTruthProtocolTransformedPoints.at(80));
this->GetDataStorage()->Modified();
}
//Add transformed points of depth 100mm to the data storage
mitk::DataNode::Pointer node100mm = this->GetDataStorage()->GetNamedNode(nameNode100mm);
if (node100mm.IsNull())
{
node100mm = mitk::DataNode::New();
node100mm->SetName(nameNode100mm);
node100mm->SetData(m_GroundTruthProtocolTransformedPoints.at(100));
this->GetDataStorage()->Add(node100mm);
}
else
{
node100mm->SetData(m_GroundTruthProtocolTransformedPoints.at(100));
this->GetDataStorage()->Modified();
}
}
//Do a global reinit
mitk::RenderingManager::GetInstance()->InitializeViewsByBoundingObjects(this->GetDataStorage());
}
double QmitkUSNavigationStepCtUsRegistration::CalculateMeanFRE()
{
double meanFRE = 0.0;
for (int counter = 0; counter < m_GroundTruthProtocolFRE.size(); ++counter)
{
meanFRE += m_GroundTruthProtocolFRE[counter];
}
return meanFRE / m_GroundTruthProtocolFRE.size();
}
double QmitkUSNavigationStepCtUsRegistration::CalculateStandardDeviationOfFRE(double meanFRE)
{
double variance = 0.0;
for (int counter = 0; counter < m_GroundTruthProtocolFRE.size(); ++counter)
{
variance += ((meanFRE - m_GroundTruthProtocolFRE[counter]) * (meanFRE - m_GroundTruthProtocolFRE[counter]));
}
variance /= m_GroundTruthProtocolFRE.size(); // calculate the empirical variance (n) and not the sampling variance (n-1)
return sqrt(variance);
}
void QmitkUSNavigationStepCtUsRegistration::CalculateGroundTruthProtocolTRE()
{
if (m_GroundTruthProtocolTransformedPoints.find(0) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(20) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(40) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(60) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(80) == m_GroundTruthProtocolTransformedPoints.end() ||
m_GroundTruthProtocolTransformedPoints.find(100) == m_GroundTruthProtocolTransformedPoints.end())
{
QMessageBox msgBox;
msgBox.setText("Cannot calculate TRE of Ground-Truth-Protocol because points were not transformed.");
msgBox.exec();
return;
}
// clear the std::map containing possibly data of earlier TRE calculations
m_GroundTruthProtocolTRE.clear();
// loop through all existing point sets containing the transformed points
for (int counter = 0;
m_GroundTruthProtocolTransformedPoints.find(counter) != m_GroundTruthProtocolTransformedPoints.end();
counter += 20)
{
//calculate the middle point of the point set
mitk::PointSet::Pointer pointSet = m_GroundTruthProtocolTransformedPoints.at(counter);
mitk::Point3D middlePoint;
middlePoint[0] = 0.0;
middlePoint[1] = 0.0;
middlePoint[2] = 0.0;
for (int position = 0; position < pointSet->GetSize(); ++position)
{
middlePoint[0] += pointSet->GetPoint(position)[0];
middlePoint[1] += pointSet->GetPoint(position)[1];
middlePoint[2] += pointSet->GetPoint(position)[2];
}
middlePoint[0] /= pointSet->GetSize();
middlePoint[1] /= pointSet->GetSize();
middlePoint[2] /= pointSet->GetSize();
MITK_INFO << "Calculated MiddlePoint: " << middlePoint;
//sum up the euclidean distances between the middle point and each transformed point
double meanDistance = 0.0;
for (int position = 0; position < pointSet->GetSize(); ++position)
{
meanDistance += middlePoint.SquaredEuclideanDistanceTo(pointSet->GetPoint(position));
MITK_INFO << "SquaredEuclideanDistance: " << middlePoint.SquaredEuclideanDistanceTo(pointSet->GetPoint(position));
}
meanDistance /= pointSet->GetSize(); // this can be interpreted as empirical variance
// the root of the empirical variance can be interpreted as the protocols registration TRE
m_GroundTruthProtocolTRE.insert(std::pair(counter, sqrt(meanDistance)));
MITK_INFO << "Ground-Truth-Protocol TRE: " << sqrt(meanDistance);
}
}
void QmitkUSNavigationStepCtUsRegistration::EliminateTooSmallLabeledObjects(
ImageType::Pointer binaryImage)
{
BinaryImageToShapeLabelMapFilterType::OutputImageType::Pointer labelMap =
m_BinaryImageToShapeLabelMapFilter->GetOutput();
double voxelVolume = this->GetVoxelVolume();
double fiducialVolume;
unsigned int numberOfPixels;
if (ui->fiducialDiameter3mmRadioButton->isChecked())
{
fiducialVolume = this->GetFiducialVolume(1.5);
numberOfPixels = ceil(fiducialVolume / voxelVolume);
}
else
{
fiducialVolume = this->GetFiducialVolume(2.5);
numberOfPixels = ceil(fiducialVolume / voxelVolume);
}
MITK_INFO << "Voxel Volume = " << voxelVolume << "; Fiducial Volume = " << fiducialVolume;
MITK_INFO << "Number of pixels = " << numberOfPixels;
labelMap = m_BinaryImageToShapeLabelMapFilter->GetOutput();
// The output of this filter is an itk::LabelMap, which contains itk::LabelObject's
MITK_INFO << "There are " << labelMap->GetNumberOfLabelObjects() << " objects.";
// Loop over each region
for (int i = labelMap->GetNumberOfLabelObjects() - 1; i >= 0; --i)
{
// Get the ith region
BinaryImageToShapeLabelMapFilterType::OutputImageType::LabelObjectType* labelObject = labelMap->GetNthLabelObject(i);
MITK_INFO << "Object " << i << " contains " << labelObject->Size() << " pixel";
//TODO: Threshold-Wert evtl. experimentell besser abstimmen,
// um zu verhindern, dass durch Threshold wahre Fiducial-Kandidaten elimiert werden.
if (labelObject->Size() < numberOfPixels * 0.8)
{
for (unsigned int pixelId = 0; pixelId < labelObject->Size(); pixelId++)
{
binaryImage->SetPixel(labelObject->GetIndex(pixelId), 0);
}
labelMap->RemoveLabelObject(labelObject);
}
}
}
bool QmitkUSNavigationStepCtUsRegistration::EliminateFiducialCandidatesByEuclideanDistances()
{
if (m_CentroidsOfFiducialCandidates.size() < NUMBER_FIDUCIALS_NEEDED)
{
return false;
}
for (int counter = 0; counter < m_CentroidsOfFiducialCandidates.size(); ++counter)
{
int amountOfAcceptedFiducials = 0;
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(counter));
//Loop through all fiducial candidates and calculate the distance between the chosen fiducial
// candidate and the other candidates. For each candidate with a right distance between
// Configuration A: 7.93mm and 31.0mm (10 mm distance between fiducial centers) or
// Configuration B: 11.895mm and 45.0mm (15 mm distance between fiducial centers) or
// Configuration C: 15.86mm and 59.0mm (20 mm distance between fiducial centers)
//
// increase the amountOfAcceptedFiducials.
for (int position = 0; position < m_CentroidsOfFiducialCandidates.size(); ++position)
{
if (position == counter)
{
continue;
}
mitk::Point3D otherCentroid(m_CentroidsOfFiducialCandidates.at(position));
double distance = fiducialCentroid.EuclideanDistanceTo(otherCentroid);
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
if (distance > 7.93 && distance < 31.0)
{
++amountOfAcceptedFiducials;
}
break;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
if (distance > 11.895 && distance < 45.0)
{
++amountOfAcceptedFiducials;
}
break;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
if (distance > 15.86 && distance < 59.0)
{
++amountOfAcceptedFiducials;
}
break;
}
}
//The amountOfAcceptedFiducials must be at least 7. Otherwise delete the fiducial candidate
// from the list of candidates.
if (amountOfAcceptedFiducials < NUMBER_FIDUCIALS_NEEDED - 1)
{
MITK_INFO << "Deleting fiducial candidate at position: " <<
m_CentroidsOfFiducialCandidates.at(counter);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + counter);
if (m_CentroidsOfFiducialCandidates.size() < NUMBER_FIDUCIALS_NEEDED )
{
return false;
}
counter = -1;
}
}
//Classify the rested fiducial candidates by its characteristic Euclidean distances
// between the canidates and remove all candidates with a false distance configuration:
this->ClassifyFiducialCandidates();
return true;
}
void QmitkUSNavigationStepCtUsRegistration::ClassifyFiducialCandidates()
{
MITK_INFO << "ClassifyFiducialCandidates()";
std::vector fiducialCandidatesToBeRemoved;
std::vector> distanceVectorsFiducials;
this->CalculateDistancesBetweenFiducials(distanceVectorsFiducials);
for (int counter = 0; counter < distanceVectorsFiducials.size(); ++counter)
{
int distanceA = 0; // => 10,00mm distance
int distanceB = 0; // => 14,14mm distance
int distanceC = 0; // => 17,32mm distance
int distanceD = 0; // => 22,36mm distance
int distanceE = 0; // => 24,49mm distance
int distanceF = 0; // => 28,28mm distance
std::vector &distances = distanceVectorsFiducials.at(counter);
for (int number = 0; number < distances.size(); ++number)
{
double &distance = distances.at(number);
switch (ui->fiducialMarkerConfigurationComboBox->currentIndex())
{
// case 0 is equal to fiducial marker configuration A (10mm distance)
case 0:
if (distance > 7.93 && distance <= 12.07)
{
++distanceA;
}
else if (distance > 12.07 && distance <= 15.73)
{
++distanceB;
}
else if (distance > 15.73 && distance <= 19.84)
{
++distanceC;
}
else if (distance > 19.84 && distance <= 23.425)
{
++distanceD;
}
else if (distance > 23.425 && distance <= 26.385)
{
++distanceE;
}
else if (distance > 26.385 && distance <= 31.00)
{
++distanceF;
}
break;
// case 1 is equal to fiducial marker configuration B (15mm distance)
case 1:
if (distance > 11.895 && distance <= 18.105)
{
++distanceA;
}
else if (distance > 18.105 && distance <= 23.595)
{
++distanceB;
}
else if (distance > 23.595 && distance <= 29.76)
{
++distanceC;
}
else if (distance > 29.76 && distance <= 35.1375)
{
++distanceD;
+ if (distance > 33.54)
+ {
+ ++distanceE;
+ }
}
else if (distance > 35.1375 && distance <= 39.5775)
{
++distanceE;
+ if (distance < 36.735)
+ {
+ ++distanceD;
+ }
}
else if (distance > 39.5775 && distance <= 45.00)
{
++distanceF;
}
break;
// case 2 is equal to fiducial marker configuration C (20mm distance)
case 2:
if (distance > 15.86 && distance <= 24.14)
{
++distanceA;
}
else if (distance > 24.14 && distance <= 31.46)
{
++distanceB;
}
else if (distance > 31.46 && distance <= 39.68)
{
++distanceC;
}
else if (distance > 39.68 && distance <= 46.85)
{
++distanceD;
}
else if (distance > 46.85 && distance <= 52.77)
{
++distanceE;
}
else if (distance > 52.77 && distance <= 59.00)
{
++distanceF;
}
break;
}
}// End for-loop distances-vector
//Now, having looped through all distances of one fiducial candidate, check
// if the combination of different distances is known. The >= is due to the
// possible occurrence of other fiducial candidates that have an distance equal to
// one of the distances A - E. However, false fiducial candidates outside
// the fiducial marker does not have the right distance configuration:
if ((distanceA >= 2 && distanceD >= 2 && distanceE >= 2 && distanceF >= 1 ||
distanceA >= 1 && distanceB >= 2 && distanceC >= 1 && distanceD >= 2 && distanceE >= 1 ||
distanceB >= 2 && distanceD >= 4 && distanceF >= 1 ||
distanceA >= 1 && distanceB >= 1 && distanceD >= 3 && distanceE >= 1 && distanceF >= 1) == false)
{
MITK_INFO << "Detected fiducial candidate with unknown distance configuration.";
fiducialCandidatesToBeRemoved.push_back(counter);
}
}
for (int count = fiducialCandidatesToBeRemoved.size() - 1; count >= 0; --count)
{
MITK_INFO << "Removing fiducial candidate " << fiducialCandidatesToBeRemoved.at(count);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin()
+ fiducialCandidatesToBeRemoved.at(count));
}
}
void QmitkUSNavigationStepCtUsRegistration::GetCentroidsOfLabeledObjects()
{
MITK_INFO << "GetCentroidsOfLabeledObjects()";
BinaryImageToShapeLabelMapFilterType::OutputImageType::Pointer labelMap =
m_BinaryImageToShapeLabelMapFilter->GetOutput();
for (int i = labelMap->GetNumberOfLabelObjects() - 1; i >= 0; --i)
{
// Get the ith region
BinaryImageToShapeLabelMapFilterType::OutputImageType::LabelObjectType* labelObject = labelMap->GetNthLabelObject(i);
MITK_INFO << "Object " << i << " contains " << labelObject->Size() << " pixel";
mitk::Vector3D centroid;
centroid[0] = labelObject->GetCentroid()[0];
centroid[1] = labelObject->GetCentroid()[1];
centroid[2] = labelObject->GetCentroid()[2];
m_CentroidsOfFiducialCandidates.push_back(centroid);
}
//evtl. for later: itk::LabelMapOverlayImageFilter
}
void QmitkUSNavigationStepCtUsRegistration::CalculatePCA()
{
//Step 1: Construct data matrix
int columnSize = m_CentroidsOfFiducialCandidates.size();
if (columnSize == 0)
{
MITK_INFO << "Cannot calculate PCA. There are no fiducial candidates.";
return;
}
vnl_matrix pointSetMatrix(3, columnSize, 0.0);
for (int counter = 0; counter < columnSize; ++counter)
{
pointSetMatrix[0][counter] = m_CentroidsOfFiducialCandidates.at(counter)[0];
pointSetMatrix[1][counter] = m_CentroidsOfFiducialCandidates.at(counter)[1];
pointSetMatrix[2][counter] = m_CentroidsOfFiducialCandidates.at(counter)[2];
}
//Step 2: Remove average for each row (Mittelwertbefreiung)
for (int counter = 0; counter < columnSize; ++counter)
{
m_MeanCentroidFiducialCandidates += mitk::Vector3D(pointSetMatrix.get_column(counter));
}
//TODO: für später überprüfen, ob Division durch integer nicht zu Rechenproblemen führt.
m_MeanCentroidFiducialCandidates /= columnSize;
for (int counter = 0; counter < columnSize; ++counter)
{
pointSetMatrix.get_column(counter) -= m_MeanCentroidFiducialCandidates;
}
//Step 3: Compute covariance matrix
vnl_matrix covarianceMatrix = (1.0 / (columnSize - 1.0)) * pointSetMatrix * pointSetMatrix.transpose();
//Step 4: Singular value composition
vnl_svd svd(covarianceMatrix);
//Storing results:
for (int counter = 0; counter < 3; ++counter)
{
mitk::Vector3D eigenVector = svd.U().get_column(counter);
double eigenValue = sqrt(svd.W(counter));
m_EigenVectorsFiducialCandidates[eigenValue] = eigenVector;
m_EigenValuesFiducialCandidates.push_back(eigenValue);
}
std::sort( m_EigenValuesFiducialCandidates.begin(), m_EigenValuesFiducialCandidates.end() );
mitk::DataNode::Pointer axis1Node = mitk::DataNode::New();
axis1Node->SetName("Eigenvector 1");
mitk::PointSet::Pointer axis1 = mitk::PointSet::New();
axis1->InsertPoint(0, m_MeanCentroidFiducialCandidates);
axis1->InsertPoint(1, (m_MeanCentroidFiducialCandidates + m_EigenVectorsFiducialCandidates.at(m_EigenValuesFiducialCandidates.at(2))*m_EigenValuesFiducialCandidates.at(2)));
axis1Node->SetData(axis1);
axis1Node->SetBoolProperty("show contour", true);
axis1Node->SetColor(1, 0, 0);
this->GetDataStorage()->Add(axis1Node);
mitk::DataNode::Pointer axis2Node = mitk::DataNode::New();
axis2Node->SetName("Eigenvector 2");
mitk::PointSet::Pointer axis2 = mitk::PointSet::New();
axis2->InsertPoint(0, m_MeanCentroidFiducialCandidates);
axis2->InsertPoint(1, (m_MeanCentroidFiducialCandidates + m_EigenVectorsFiducialCandidates.at(m_EigenValuesFiducialCandidates.at(1))*m_EigenValuesFiducialCandidates.at(1)));
axis2Node->SetData(axis2);
axis2Node->SetBoolProperty("show contour", true);
axis2Node->SetColor(2, 0, 0);
this->GetDataStorage()->Add(axis2Node);
mitk::DataNode::Pointer axis3Node = mitk::DataNode::New();
axis3Node->SetName("Eigenvector 3");
mitk::PointSet::Pointer axis3 = mitk::PointSet::New();
axis3->InsertPoint(0, m_MeanCentroidFiducialCandidates);
axis3->InsertPoint(1, (m_MeanCentroidFiducialCandidates + m_EigenVectorsFiducialCandidates.at(m_EigenValuesFiducialCandidates.at(0))*m_EigenValuesFiducialCandidates.at(0)));
axis3Node->SetData(axis3);
axis3Node->SetBoolProperty("show contour", true);
axis3Node->SetColor(3, 0, 0);
this->GetDataStorage()->Add(axis3Node);
MITK_INFO << "Mean: " << m_MeanCentroidFiducialCandidates;
MITK_INFO << "Eigenvektor 1: " << m_EigenVectorsFiducialCandidates.at(m_EigenValuesFiducialCandidates.at(2));
MITK_INFO << "Eigenvektor 2: " << m_EigenVectorsFiducialCandidates.at(m_EigenValuesFiducialCandidates.at(1));
MITK_INFO << "Eigenvektor 3: " << m_EigenVectorsFiducialCandidates.at(m_EigenValuesFiducialCandidates.at(0));
MITK_INFO << "Eigenwert 1: " << m_EigenValuesFiducialCandidates.at(2);
MITK_INFO << "Eigenwert 2: " << m_EigenValuesFiducialCandidates.at(1);
MITK_INFO << "Eigenwert 3: " << m_EigenValuesFiducialCandidates.at(0);
}
void QmitkUSNavigationStepCtUsRegistration::NumerateFiducialMarks()
{
MITK_INFO << "NumerateFiducialMarks()";
bool successFiducialNo1;
bool successFiducialNo4;
bool successFiducialNo2And3;
bool successFiducialNo5;
bool successFiducialNo8;
bool successFiducialNo6;
bool successFiducialNo7;
std::vector> distanceVectorsFiducials;
this->CalculateDistancesBetweenFiducials(distanceVectorsFiducials);
successFiducialNo1 = this->FindFiducialNo1(distanceVectorsFiducials);
successFiducialNo4 = this->FindFiducialNo4(distanceVectorsFiducials);
successFiducialNo2And3 = this->FindFiducialNo2And3();
successFiducialNo5 = this->FindFiducialNo5();
successFiducialNo8 = this->FindFiducialNo8();
successFiducialNo6 = this->FindFiducialNo6();
successFiducialNo7 = this->FindFiducialNo7();
if (!successFiducialNo1 || !successFiducialNo4 || !successFiducialNo2And3 ||
!successFiducialNo5 || !successFiducialNo8 || !successFiducialNo6 || !successFiducialNo7)
{
QMessageBox msgBox;
msgBox.setText("Cannot numerate/localize all fiducials successfully.");
msgBox.exec();
return;
}
if( ui->optimizeFiducialPositionsCheckBox->isChecked())
{
this->OptimizeFiducialPositions();
}
if (m_MarkerFloatingImageCoordinateSystemPointSet.IsNull())
{
m_MarkerFloatingImageCoordinateSystemPointSet = mitk::PointSet::New();
}
else if (m_MarkerFloatingImageCoordinateSystemPointSet->GetSize() != 0)
{
m_MarkerFloatingImageCoordinateSystemPointSet->Clear();
}
for (int counter = 1; counter <= m_FiducialMarkerCentroids.size(); ++counter)
{
m_MarkerFloatingImageCoordinateSystemPointSet->InsertPoint(counter - 1, m_FiducialMarkerCentroids.at(counter));
}
if( !m_PerformingGroundTruthProtocolEvaluation )
{
mitk::DataNode::Pointer node = mitk::DataNode::New();
node->SetData(m_MarkerFloatingImageCoordinateSystemPointSet);
node->SetName("MarkerFloatingImageCSPointSet");
//node->SetFloatProperty("pointsize", 5.0);
this->GetDataStorage()->Add(node);
}
}
void QmitkUSNavigationStepCtUsRegistration::ShowGroundTruthMarkerEdges()
{
mitk::Point3D edge1;
mitk::Point3D edge2;
mitk::Point3D edge3;
mitk::Point3D edge4;
mitk::Point3D edge5;
mitk::Point3D edge6;
mitk::Point3D edge7;
mitk::Point3D edge8;
edge1[0] = -5;
edge1[1] = -5;
edge1[2] = -7.5;
edge2[0] = -5;
edge2[1] = 25;
edge2[2] = -7.5;
edge3[0] = 25;
edge3[1] = 25;
edge3[2] = -7.5;
edge4[0] = 25;
edge4[1] = -5;
edge4[2] = -7.5;
edge5[0] = -5;
edge5[1] = -5;
edge5[2] = 15;
edge6[0] = -5;
edge6[1] = 25;
edge6[2] = 15;
edge7[0] = 25;
edge7[1] = 25;
edge7[2] = 15;
edge8[0] = 25;
edge8[1] = -5;
edge8[2] = 15;
mitk::PointSet::Pointer edgePointSet = mitk::PointSet::New();
edgePointSet->InsertPoint(0, edge1);
edgePointSet->InsertPoint(1, edge2);
edgePointSet->InsertPoint(2, edge3);
edgePointSet->InsertPoint(3, edge4);
edgePointSet->InsertPoint(4, edge1);
edgePointSet->InsertPoint(5, edge5);
edgePointSet->InsertPoint(6, edge6);
edgePointSet->InsertPoint(7, edge7);
edgePointSet->InsertPoint(8, edge8);
edgePointSet->InsertPoint(9, edge5);
edgePointSet->InsertPoint(10, edge6);
edgePointSet->InsertPoint(11, edge2);
edgePointSet->InsertPoint(12, edge6);
edgePointSet->InsertPoint(13, edge7);
edgePointSet->InsertPoint(14, edge3);
edgePointSet->InsertPoint(15, edge7);
edgePointSet->InsertPoint(16, edge8);
edgePointSet->InsertPoint(17, edge4);
mitk::DataNode::Pointer node = mitk::DataNode::New();
node->SetData(edgePointSet);
node->SetName("Marker Edges PointSet");
node->SetBoolProperty("show contour", true);
node->SetFloatProperty("pointsize", 0.25);
this->GetDataStorage()->Add(node);
mitk::PointSet* pointSet_orig = edgePointSet;
mitk::PointSet::Pointer pointSet_moved = mitk::PointSet::New();
for (int i = 0; i < pointSet_orig->GetSize(); i++)
{
pointSet_moved->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(pointSet_orig->GetPoint(i)));
}
pointSet_orig->Clear();
for (int i = 0; i < pointSet_moved->GetSize(); i++)
pointSet_orig->InsertPoint(pointSet_moved->GetPoint(i));
//Do a global reinit
mitk::RenderingManager::GetInstance()->InitializeViewsByBoundingObjects(this->GetDataStorage());
}
void QmitkUSNavigationStepCtUsRegistration::CalculateDistancesBetweenFiducials(std::vector>& distanceVectorsFiducials)
{
std::vector distancesBetweenFiducials;
for (int i = 0; i < m_CentroidsOfFiducialCandidates.size(); ++i)
{
distancesBetweenFiducials.clear();
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(i));
for (int n = 0; n < m_CentroidsOfFiducialCandidates.size(); ++n)
{
mitk::Point3D otherCentroid(m_CentroidsOfFiducialCandidates.at(n));
distancesBetweenFiducials.push_back(fiducialCentroid.EuclideanDistanceTo(otherCentroid));
}
//Sort the distances from low to big numbers
std::sort(distancesBetweenFiducials.begin(), distancesBetweenFiducials.end());
//First entry of the distance vector must be 0, so erase it
if (distancesBetweenFiducials.at(0) == 0.0)
{
distancesBetweenFiducials.erase(distancesBetweenFiducials.begin());
}
//Add the distance vector to the collecting distances vector
distanceVectorsFiducials.push_back(distancesBetweenFiducials);
}
for (int i = 0; i < distanceVectorsFiducials.size(); ++i)
{
MITK_INFO << "Vector " << i << ":";
for (int k = 0; k < distanceVectorsFiducials.at(i).size(); ++k)
{
MITK_INFO << distanceVectorsFiducials.at(i).at(k);
}
}
}
bool QmitkUSNavigationStepCtUsRegistration::FindFiducialNo1(std::vector>& distanceVectorsFiducials)
{
for (int i = 0; i < distanceVectorsFiducials.size(); ++i)
{
std::vector &distances = distanceVectorsFiducials.at(i);
if (distances.size() < NUMBER_FIDUCIALS_NEEDED - 1 )
{
MITK_WARN << "Cannot find fiducial 1, there aren't found enough fiducial candidates.";
return false;
}
double characteristicDistanceAWithUpperMargin = this->GetCharacteristicDistanceAWithUpperMargin();
if (distances.at(0) <= characteristicDistanceAWithUpperMargin &&
distances.at(1) <= characteristicDistanceAWithUpperMargin)
{
MITK_INFO << "Found Fiducial 1 (PointSet number " << i << ")";
m_FiducialMarkerCentroids.insert( std::pair(1, m_CentroidsOfFiducialCandidates.at(i)));
distanceVectorsFiducials.erase(distanceVectorsFiducials.begin() + i);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + i);
return true;
}
}
return false;
}
bool QmitkUSNavigationStepCtUsRegistration::FindFiducialNo2And3()
{
if (m_FiducialMarkerCentroids.find(1) == m_FiducialMarkerCentroids.end() )
{
MITK_WARN << "Cannot find fiducial No 2 and 3. Before must be found fiducial No 1.";
return false;
}
mitk::Point3D fiducialNo1(m_FiducialMarkerCentroids.at(1));
mitk::Vector3D fiducialVectorA;
mitk::Vector3D fiducialVectorB;
mitk::Point3D fiducialPointA;
mitk::Point3D fiducialPointB;
bool foundFiducialA = false;
bool foundFiducialB = false;
mitk::Vector3D vectorFiducial1ToFiducialA;
mitk::Vector3D vectorFiducial1ToFiducialB;
for (int i = 0; i < m_CentroidsOfFiducialCandidates.size(); ++i)
{
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(i));
double distance = fiducialNo1.EuclideanDistanceTo(fiducialCentroid);
if (distance <= this->GetCharacteristicDistanceAWithUpperMargin())
{
fiducialVectorA = m_CentroidsOfFiducialCandidates.at(i);
fiducialPointA = fiducialCentroid;
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + i);
foundFiducialA = true;
break;
}
}
for (int i = 0; i < m_CentroidsOfFiducialCandidates.size(); ++i)
{
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(i));
double distance = fiducialNo1.EuclideanDistanceTo(fiducialCentroid);
if (distance <= this->GetCharacteristicDistanceAWithUpperMargin())
{
fiducialVectorB = m_CentroidsOfFiducialCandidates.at(i);
fiducialPointB = fiducialCentroid;
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + i);
foundFiducialB = true;
break;
}
}
if (!foundFiducialA || !foundFiducialB)
{
MITK_WARN << "Cannot identify fiducial candidates 2 and 3";
return false;
}
else if (m_CentroidsOfFiducialCandidates.size() == 0)
{
MITK_WARN << "Too less fiducials detected. Cannot identify fiducial candidates 2 and 3";
return false;
}
vectorFiducial1ToFiducialA = fiducialVectorA - m_FiducialMarkerCentroids.at(1);
vectorFiducial1ToFiducialB = fiducialVectorB - m_FiducialMarkerCentroids.at(1);
vnl_vector crossProductVnl = vnl_cross_3d(vectorFiducial1ToFiducialA.GetVnlVector(), vectorFiducial1ToFiducialB.GetVnlVector());
mitk::Vector3D crossProduct;
crossProduct.SetVnlVector(crossProductVnl);
mitk::Vector3D vectorFiducial1ToRandomLeftFiducial = m_CentroidsOfFiducialCandidates.at(0) - m_FiducialMarkerCentroids.at(1);
double scalarProduct = (crossProduct * vectorFiducial1ToRandomLeftFiducial) /
(crossProduct.GetNorm() * vectorFiducial1ToRandomLeftFiducial.GetNorm());
double alpha = acos(scalarProduct) * 57.29578; //Transform into degree
MITK_INFO << "Scalar Product = " << alpha;
if (alpha <= 90)
{
m_FiducialMarkerCentroids[3] = fiducialVectorA;
m_FiducialMarkerCentroids[2] = fiducialVectorB;
}
else
{
m_FiducialMarkerCentroids[2] = fiducialVectorA;
m_FiducialMarkerCentroids[3] = fiducialVectorB;
}
MITK_INFO << "Found Fiducial 2, PointSet: " << m_FiducialMarkerCentroids.at(2);
MITK_INFO << "Found Fiducial 3, PointSet: " << m_FiducialMarkerCentroids.at(3);
return true;
}
bool QmitkUSNavigationStepCtUsRegistration::FindFiducialNo4(std::vector>& distanceVectorsFiducials)
{
double characteristicDistanceBWithLowerMargin = this->GetCharacteristicDistanceBWithLowerMargin();
double characteristicDistanceBWithUpperMargin = this->GetCharacteristicDistanceBWithUpperMargin();
for (int i = 0; i < distanceVectorsFiducials.size(); ++i)
{
std::vector &distances = distanceVectorsFiducials.at(i);
if (distances.size() < NUMBER_FIDUCIALS_NEEDED - 1)
{
MITK_WARN << "Cannot find fiducial 4, there aren't found enough fiducial candidates.";
return false;
}
if (distances.at(0) > characteristicDistanceBWithLowerMargin &&
distances.at(0) <= characteristicDistanceBWithUpperMargin &&
distances.at(1) > characteristicDistanceBWithLowerMargin &&
distances.at(1) <= characteristicDistanceBWithUpperMargin)
{
MITK_INFO << "Found Fiducial 4 (PointSet number " << i << ")";
m_FiducialMarkerCentroids.insert(std::pair(4, m_CentroidsOfFiducialCandidates.at(i)));
distanceVectorsFiducials.erase(distanceVectorsFiducials.begin() + i);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + i);
return true;
}
}
return false;
}
bool QmitkUSNavigationStepCtUsRegistration::FindFiducialNo5()
{
if (m_FiducialMarkerCentroids.find(2) == m_FiducialMarkerCentroids.end())
{
MITK_WARN << "To find fiducial No 5, fiducial No 2 has to be found before.";
return false;
}
double characteristicDistanceBWithUpperMargin = this->GetCharacteristicDistanceBWithUpperMargin();
mitk::Point3D fiducialNo2(m_FiducialMarkerCentroids.at(2));
for (int counter = 0; counter < m_CentroidsOfFiducialCandidates.size(); ++counter)
{
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(counter));
double distance = fiducialNo2.EuclideanDistanceTo(fiducialCentroid);
if (distance <= characteristicDistanceBWithUpperMargin)
{
m_FiducialMarkerCentroids[5] = m_CentroidsOfFiducialCandidates.at(counter);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + counter);
MITK_INFO << "Found Fiducial No 5, PointSet: " << m_FiducialMarkerCentroids[5];
return true;
}
}
MITK_WARN << "Cannot find fiducial No 5.";
return false;
}
bool QmitkUSNavigationStepCtUsRegistration::FindFiducialNo6()
{
if (m_FiducialMarkerCentroids.find(5) == m_FiducialMarkerCentroids.end())
{
MITK_WARN << "To find fiducial No 6, fiducial No 5 has to be found before.";
return false;
}
double characteristicDistanceAWithUpperMargin = this->GetCharacteristicDistanceAWithUpperMargin();
mitk::Point3D fiducialNo5(m_FiducialMarkerCentroids.at(5));
for (int counter = 0; counter < m_CentroidsOfFiducialCandidates.size(); ++counter)
{
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(counter));
double distance = fiducialNo5.EuclideanDistanceTo(fiducialCentroid);
if (distance <= characteristicDistanceAWithUpperMargin)
{
m_FiducialMarkerCentroids[6] = m_CentroidsOfFiducialCandidates.at(counter);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + counter);
MITK_INFO << "Found Fiducial No 6, PointSet: " << m_FiducialMarkerCentroids[6];
return true;
}
}
MITK_WARN << "Cannot find fiducial No 6.";
return false;
}
bool QmitkUSNavigationStepCtUsRegistration::FindFiducialNo7()
{
if (m_FiducialMarkerCentroids.find(8) == m_FiducialMarkerCentroids.end())
{
MITK_WARN << "To find fiducial No 7, fiducial No 8 has to be found before.";
return false;
}
double characteristicDistanceAWithUpperMargin = this->GetCharacteristicDistanceAWithUpperMargin();
mitk::Point3D fiducialNo8(m_FiducialMarkerCentroids.at(8));
for (int counter = 0; counter < m_CentroidsOfFiducialCandidates.size(); ++counter)
{
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(counter));
double distance = fiducialNo8.EuclideanDistanceTo(fiducialCentroid);
if (distance <= characteristicDistanceAWithUpperMargin)
{
m_FiducialMarkerCentroids[7] = m_CentroidsOfFiducialCandidates.at(counter);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + counter);
MITK_INFO << "Found Fiducial No 7, PointSet: " << m_FiducialMarkerCentroids[7];
return true;
}
}
MITK_WARN << "Cannot find fiducial No 7.";
return false;
}
bool QmitkUSNavigationStepCtUsRegistration::FindFiducialNo8()
{
if (m_FiducialMarkerCentroids.find(3) == m_FiducialMarkerCentroids.end())
{
MITK_WARN << "To find fiducial No 8, fiducial No 3 has to be found before.";
return false;
}
double characteristicDistanceBWithUpperMargin = this->GetCharacteristicDistanceBWithUpperMargin();
mitk::Point3D fiducialNo3(m_FiducialMarkerCentroids.at(3));
for (int counter = 0; counter < m_CentroidsOfFiducialCandidates.size(); ++counter)
{
mitk::Point3D fiducialCentroid(m_CentroidsOfFiducialCandidates.at(counter));
double distance = fiducialNo3.EuclideanDistanceTo(fiducialCentroid);
if (distance <= characteristicDistanceBWithUpperMargin)
{
m_FiducialMarkerCentroids[8] = m_CentroidsOfFiducialCandidates.at(counter);
m_CentroidsOfFiducialCandidates.erase(m_CentroidsOfFiducialCandidates.begin() + counter);
MITK_INFO << "Found Fiducial No 8, PointSet: " << m_FiducialMarkerCentroids[8];
return true;
}
}
MITK_WARN << "Cannot find fiducial No 8.";
return false;
}
void QmitkUSNavigationStepCtUsRegistration::OptimizeFiducialPositions()
{
MITK_INFO << "OptimizeFiducialPositions()";
//Initialization for planes 1 - 2 - 3 - 4 and 5 - 6 - 7 - 8
mitk::PlaneFit::Pointer planeFit1234 = mitk::PlaneFit::New();
mitk::PlaneFit::Pointer planeFit5678 = mitk::PlaneFit::New();
mitk::PointSet::Pointer pointSet1234 = mitk::PointSet::New();
mitk::PointSet::Pointer pointSet5678 = mitk::PointSet::New();
mitk::PlaneGeometry::Pointer planeGeometry1234 = mitk::PlaneGeometry::New();
mitk::PlaneGeometry::Pointer planeGeometry5678 = mitk::PlaneGeometry::New();
for (int counter = 1; counter <= 4; ++counter)
{
pointSet1234->InsertPoint(counter - 1, m_FiducialMarkerCentroids.at(counter));
}
for (int counter = 5; counter <= 8; ++counter)
{
pointSet5678->InsertPoint(counter - 5, m_FiducialMarkerCentroids.at(counter));
}
//Make planes parallel
this->CreateParallelPlanes(planeFit1234, planeFit5678,
pointSet1234, pointSet5678,
planeGeometry1234, planeGeometry5678,
true);
this->MovePlanes(planeGeometry1234, planeGeometry5678, this->GetMinimalFiducialConfigurationDistance());
//Move the points into the parallel planes
for (int counter = 1; counter <= 4; ++counter)
{
this->MovePoint(planeGeometry1234, counter);
}
for (int counter = 5; counter <= 8; ++counter)
{
this->MovePoint(planeGeometry5678, counter);
}
MITK_INFO << "NormalVector plane 1234 = " << planeGeometry1234->GetNormal();
MITK_INFO << "NormalVector plane 5678 = " << planeGeometry5678->GetNormal();
MITK_INFO << "Are parallel: " << planeGeometry1234->IsParallel(planeGeometry5678);
//Optimize now the positions of Fiducials 1 - 2 - 5 and 4 - 7 - 8
//Initialization for planes 1 - 2 - 5 and 4 - 7 - 8
mitk::PlaneFit::Pointer planeFit125 = mitk::PlaneFit::New();
mitk::PlaneFit::Pointer planeFit478 = mitk::PlaneFit::New();
mitk::PointSet::Pointer pointSet125 = mitk::PointSet::New();
mitk::PointSet::Pointer pointSet478 = mitk::PointSet::New();
mitk::PlaneGeometry::Pointer planeGeometry125 = mitk::PlaneGeometry::New();
mitk::PlaneGeometry::Pointer planeGeometry478 = mitk::PlaneGeometry::New();
pointSet125->InsertPoint(0, m_FiducialMarkerCentroids.at(1));
pointSet125->InsertPoint(1, m_FiducialMarkerCentroids.at(2));
pointSet125->InsertPoint(2, m_FiducialMarkerCentroids.at(5));
//Add the points projected onto the opposite (parallel) plane to the pointset
// By means of these projected points the calculated plane
// is orthogonal to the plane 1234 and 5678.
pointSet125->InsertPoint(3, planeGeometry5678->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(1)));
pointSet125->InsertPoint(4, planeGeometry5678->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(2)));
pointSet125->InsertPoint(5, planeGeometry1234->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(5)));
pointSet478->InsertPoint(0, m_FiducialMarkerCentroids.at(4));
pointSet478->InsertPoint(1, m_FiducialMarkerCentroids.at(7));
pointSet478->InsertPoint(2, m_FiducialMarkerCentroids.at(8));
//Add the points projected onto the opposite (parallel) plane to the pointset
// By means of these projected points the calculated plane
// is orthogonal to the plane 1234 and 5678.
pointSet478->InsertPoint(3, planeGeometry5678->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(4)));
pointSet478->InsertPoint(4, planeGeometry1234->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(7)));
pointSet478->InsertPoint(5, planeGeometry1234->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(8)));
//Make planes parallel
this->CreateParallelPlanes(planeFit125, planeFit478,
pointSet125, pointSet478,
planeGeometry125, planeGeometry478,
false);
this->MovePlanes(planeGeometry125, planeGeometry478, (2 * this->GetMinimalFiducialConfigurationDistance()));
//Move the points into the parallel planes
this->MovePoint(planeGeometry125, 1);
this->MovePoint(planeGeometry125, 2);
this->MovePoint(planeGeometry125, 5);
this->MovePoint(planeGeometry478, 4);
this->MovePoint(planeGeometry478, 7);
this->MovePoint(planeGeometry478, 8);
MITK_INFO << "NormalVector plane 125 = " << planeGeometry125->GetNormal();
MITK_INFO << "NormalVector plane 478 = " << planeGeometry478->GetNormal();
MITK_INFO << "Are parallel: " << planeGeometry125->IsParallel(planeGeometry478);
//Optimize now the positions of Fiducials 1 - 3 - 8 and 4 - 5 - 6
//Initialization for planes 1 - 3 - 8 and 4 - 5 - 6
mitk::PlaneFit::Pointer planeFit138 = mitk::PlaneFit::New();
mitk::PlaneFit::Pointer planeFit456 = mitk::PlaneFit::New();
mitk::PointSet::Pointer pointSet138 = mitk::PointSet::New();
mitk::PointSet::Pointer pointSet456 = mitk::PointSet::New();
mitk::PlaneGeometry::Pointer planeGeometry138 = mitk::PlaneGeometry::New();
mitk::PlaneGeometry::Pointer planeGeometry456 = mitk::PlaneGeometry::New();
pointSet138->InsertPoint(0, m_FiducialMarkerCentroids.at(1));
pointSet138->InsertPoint(1, m_FiducialMarkerCentroids.at(3));
pointSet138->InsertPoint(2, m_FiducialMarkerCentroids.at(8));
//Add the points projected onto the opposite (parallel) plane to the pointset
// By means of these projected points the calculated plane
// is orthogonal to the plane 1234 and 5678.
pointSet138->InsertPoint(3, planeGeometry5678->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(1)));
pointSet138->InsertPoint(4, planeGeometry5678->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(3)));
pointSet138->InsertPoint(5, planeGeometry1234->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(8)));
pointSet456->InsertPoint(0, m_FiducialMarkerCentroids.at(4));
pointSet456->InsertPoint(1, m_FiducialMarkerCentroids.at(5));
pointSet456->InsertPoint(2, m_FiducialMarkerCentroids.at(6));
//Add the points projected onto the opposite (parallel) plane to the pointset
// By means of these projected points the calculated plane
// is orthogonal to the plane 1234 and 5678.
pointSet456->InsertPoint(3, planeGeometry5678->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(4)));
pointSet456->InsertPoint(4, planeGeometry1234->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(5)));
pointSet456->InsertPoint(5, planeGeometry1234->ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(6)));
//Make planes parallel
this->CreateParallelPlanes(planeFit138, planeFit456,
pointSet138, pointSet456,
planeGeometry138, planeGeometry456,
false);
this->MovePlanes(planeGeometry138, planeGeometry456, (2 * this->GetMinimalFiducialConfigurationDistance()));
//Move the points into the parallel planes
this->MovePoint(planeGeometry138, 1);
this->MovePoint(planeGeometry138, 3);
this->MovePoint(planeGeometry138, 8);
this->MovePoint(planeGeometry456, 4);
this->MovePoint(planeGeometry456, 5);
this->MovePoint(planeGeometry456, 6);
MITK_INFO << "NormalVector plane 138 = " << planeGeometry138->GetNormal();
MITK_INFO << "NormalVector plane 456 = " << planeGeometry456->GetNormal();
MITK_INFO << "Are parallel: " << planeGeometry138->IsParallel(planeGeometry456);
}
void QmitkUSNavigationStepCtUsRegistration::CreateParallelPlanes(
mitk::PlaneFit::Pointer planeA, mitk::PlaneFit::Pointer planeB,
mitk::PointSet::Pointer pointSetA, mitk::PointSet::Pointer pointSetB,
mitk::PlaneGeometry::Pointer planeGeometryA, mitk::PlaneGeometry::Pointer planeGeometryB,
bool minimizeInfluenceOutliers )
{
planeA->SetInput(pointSetA);
planeA->Update();
mitk::PlaneGeometry::Pointer geometryA = dynamic_cast(planeA->GetOutput()->GetGeometry());
planeB->SetInput(pointSetB);
planeB->Update();
mitk::PlaneGeometry::Pointer geometryB = dynamic_cast(planeB->GetOutput()->GetGeometry());
mitk::DataNode::Pointer node1 = mitk::DataNode::New();
MITK_INFO << "Angle before = " << geometryA->Angle(geometryB) * 57.29578; //Transform into degree
//Minimize influence of outliers concerning the inclination angle of the plane:
if (minimizeInfluenceOutliers)
{
bool minimizeA = true;
bool minimizeB = true;
while(minimizeA)
{
MITK_INFO << "Minimize A";
std::vector> distancesA;
distancesA.push_back(std::pair(geometryA->Distance(m_FiducialMarkerCentroids.at(1)), 1));
distancesA.push_back(std::pair(geometryA->Distance(m_FiducialMarkerCentroids.at(2)), 2));
distancesA.push_back(std::pair(geometryA->Distance(m_FiducialMarkerCentroids.at(3)), 3));
distancesA.push_back(std::pair(geometryA->Distance(m_FiducialMarkerCentroids.at(4)), 4));
std::sort(distancesA.begin(), distancesA.end());
for (std::vector>::iterator it = distancesA.begin(); it != distancesA.end(); ++it)
{
MITK_INFO << "First = " << (*it).first << " Second = " << (*it).second;
}
if ((*distancesA.rbegin()).first < 0.005)
{
minimizeA = false;
break;
}
std::vector>::reverse_iterator it = distancesA.rbegin();
int fiducialNoToBeMovedToPlane1 = (*it).second;
++it;
int fiducialNoToBeMovedToPlane2 = (*it).second;
this->MovePoint(geometryA, fiducialNoToBeMovedToPlane1);
this->MovePoint(geometryA, fiducialNoToBeMovedToPlane2);
pointSetA->Clear();
for (int counter = 1; counter <= 4; ++counter)
{
pointSetA->InsertPoint(counter - 1, m_FiducialMarkerCentroids.at(counter));
}
planeA = mitk::PlaneFit::New();
planeA->SetInput(pointSetA);
planeA->Update();
geometryA = dynamic_cast(planeA->GetOutput()->GetGeometry());
}
while (minimizeB)
{
MITK_INFO << "Minimize B";
std::vector> distancesB;
distancesB.push_back(std::pair(geometryB->Distance(m_FiducialMarkerCentroids.at(5)), 5));
distancesB.push_back(std::pair(geometryB->Distance(m_FiducialMarkerCentroids.at(6)), 6));
distancesB.push_back(std::pair(geometryB->Distance(m_FiducialMarkerCentroids.at(7)), 7));
distancesB.push_back(std::pair(geometryB->Distance(m_FiducialMarkerCentroids.at(8)), 8));
std::sort(distancesB.begin(), distancesB.end());
for (std::vector>::iterator it = distancesB.begin(); it != distancesB.end(); ++it)
{
MITK_INFO << "First = " << (*it).first << " Second = " << (*it).second;
}
if ((*distancesB.rbegin()).first < 0.005)
{
minimizeB = false;
break;
}
std::vector>::reverse_iterator it = distancesB.rbegin();
int fiducialNoToBeMovedToPlane1 = (*it).second;
++it;
int fiducialNoToBeMovedToPlane2 = (*it).second;
this->MovePoint(geometryB, fiducialNoToBeMovedToPlane1);
this->MovePoint(geometryB, fiducialNoToBeMovedToPlane2);
pointSetB->Clear();
for (int counter = 5; counter <= 8; ++counter)
{
pointSetB->InsertPoint(counter - 5, m_FiducialMarkerCentroids.at(counter));
}
planeB = mitk::PlaneFit::New();
planeB->SetInput(pointSetB);
planeB->Update();
geometryB = dynamic_cast(planeB->GetOutput()->GetGeometry());
}
MITK_INFO << "Angle after = " << geometryA->Angle(geometryB) * 57.29578; //Transform into degree
}
// End outliers minimization
mitk::Point3D originPlaneA = geometryA->GetOrigin();
mitk::Point3D originPlaneB = geometryB->GetOrigin();
mitk::Vector3D normalPlaneA = geometryA->GetNormal();
mitk::Vector3D normalPlaneB = geometryB->GetNormal();
double scalarProduct = (normalPlaneA * normalPlaneB) /
(normalPlaneA.GetNorm() * normalPlaneB.GetNorm());
double alpha = acos(scalarProduct) * 57.29578; //Transform into degree
if (alpha > 90)
{
normalPlaneB *= -1;
}
mitk::Vector3D combinedNormalPlaneA_B = 0.5 * (normalPlaneA + normalPlaneB);
combinedNormalPlaneA_B.Normalize();
planeGeometryA->InitializePlane(originPlaneA, combinedNormalPlaneA_B);
planeGeometryB->InitializePlane(originPlaneB, combinedNormalPlaneA_B);
}
void QmitkUSNavigationStepCtUsRegistration::MovePlanes(mitk::PlaneGeometry::Pointer planeA, mitk::PlaneGeometry::Pointer planeB, double referenceDistance)
{
const mitk::PlaneGeometry *constPlaneB = planeB;
double distanceBetweenPlanes = planeA->DistanceFromPlane(constPlaneB);
MITK_INFO << "Distance between planes before moving = " << distanceBetweenPlanes;
//If distanceBetweenPlanes is > referenceDistance, the result will be negative,
// otherwise it will be positive:
double distanceToMove = 0.5 * (referenceDistance - distanceBetweenPlanes);
const mitk::Vector3D movingVectorPositive = distanceToMove * planeA->GetNormal();
const mitk::Vector3D movingVectorNegative = distanceToMove * planeA->GetNormal() * -1;
//Check, whether the planeB is above planeA
if (planeA->IsAbove(planeB->GetOrigin()))
{
//planeB is above planeA
planeA->Translate(movingVectorNegative);
planeB->Translate(movingVectorPositive);
}
else
{
//planeB is below planeA
planeA->Translate(movingVectorPositive);
planeB->Translate(movingVectorNegative);
}
const mitk::PlaneGeometry *constPlaneBMoved = planeB;
MITK_INFO << "Distance between planes after moving = " << planeA->DistanceFromPlane(constPlaneBMoved);
}
void QmitkUSNavigationStepCtUsRegistration::MovePoint(
mitk::PlaneGeometry::Pointer planeGeometry, int fiducialNo)
{
MITK_INFO << "Moved Point from " << m_FiducialMarkerCentroids.at(fiducialNo);
//Move all points outside the plane to the plane (finally all points lie into the plane):
mitk::Point3D movedPoint = planeGeometry->
ProjectPointOntoPlane(m_FiducialMarkerCentroids.at(fiducialNo));
m_FiducialMarkerCentroids.at(fiducialNo)[0] = movedPoint[0];
m_FiducialMarkerCentroids.at(fiducialNo)[1] = movedPoint[1];
m_FiducialMarkerCentroids.at(fiducialNo)[2] = movedPoint[2];
MITK_INFO << " to " << m_FiducialMarkerCentroids.at(fiducialNo);
}
void QmitkUSNavigationStepCtUsRegistration::DefineDataStorageImageFilter()
{
m_IsAPointSetPredicate = mitk::TNodePredicateDataType::New();
mitk::TNodePredicateDataType::Pointer isImage = mitk::TNodePredicateDataType::New();
auto isSegmentation = mitk::NodePredicateDataType::New("Segment");
mitk::NodePredicateOr::Pointer validImages = mitk::NodePredicateOr::New();
validImages->AddPredicate(mitk::NodePredicateAnd::New(isImage, mitk::NodePredicateNot::New(isSegmentation)));
mitk::NodePredicateNot::Pointer isNotAHelperObject = mitk::NodePredicateNot::New(mitk::NodePredicateProperty::New("helper object", mitk::BoolProperty::New(true)));
m_IsOfTypeImagePredicate = mitk::NodePredicateAnd::New(validImages, isNotAHelperObject);
mitk::NodePredicateProperty::Pointer isBinaryPredicate = mitk::NodePredicateProperty::New("binary", mitk::BoolProperty::New(true));
mitk::NodePredicateNot::Pointer isNotBinaryPredicate = mitk::NodePredicateNot::New(isBinaryPredicate);
mitk::NodePredicateAnd::Pointer isABinaryImagePredicate = mitk::NodePredicateAnd::New(m_IsOfTypeImagePredicate, isBinaryPredicate);
mitk::NodePredicateAnd::Pointer isNotABinaryImagePredicate = mitk::NodePredicateAnd::New(m_IsOfTypeImagePredicate, isNotBinaryPredicate);
m_IsASegmentationImagePredicate = mitk::NodePredicateOr::New(isABinaryImagePredicate, mitk::TNodePredicateDataType::New());
m_IsAPatientImagePredicate = mitk::NodePredicateAnd::New(isNotABinaryImagePredicate, mitk::NodePredicateNot::New(mitk::TNodePredicateDataType::New()));
}
void QmitkUSNavigationStepCtUsRegistration::CreateQtPartControl(QWidget *parent)
{
ui->setupUi(parent);
ui->floatingImageComboBox->SetPredicate(m_IsAPatientImagePredicate);
ui->groundTruthImageComboBox->SetPredicate(m_IsAPatientImagePredicate);
ui->ctImagesToChooseComboBox->SetPredicate(m_IsAPatientImagePredicate);
// create signal/slot connections
connect(ui->floatingImageComboBox, SIGNAL(OnSelectionChanged(const mitk::DataNode*)),
this, SLOT(OnFloatingImageComboBoxSelectionChanged(const mitk::DataNode*)));
connect(ui->groundTruthImageComboBox, SIGNAL(OnSelectionChanged(const mitk::DataNode*)),
this, SLOT(OnGroundTruthImageComboBoxSelectionChanged(const mitk::DataNode*)));
connect(ui->doRegistrationMarkerToImagePushButton, SIGNAL(clicked()),
this, SLOT(OnRegisterMarkerToFloatingImageCS()));
connect(ui->localizeFiducialMarkerPushButton, SIGNAL(clicked()),
this, SLOT(OnLocalizeFiducials()));
connect(ui->filterImageGroundTruthEvaluationPushButton, SIGNAL(clicked()),
this, SLOT(OnFilterGroundTruthImage()));
connect(ui->addCtImagePushButton, SIGNAL(clicked()),
this, SLOT(OnAddCtImageClicked()));
connect(ui->removeCtImagePushButton, SIGNAL(clicked()),
this, SLOT(OnRemoveCtImageClicked()));
connect(ui->evaluateProtocolPushButton, SIGNAL(clicked()),
this, SLOT(OnEvaluateGroundTruthFiducialLocalizationProtocol()));
}
void QmitkUSNavigationStepCtUsRegistration::OnFloatingImageComboBoxSelectionChanged(const mitk::DataNode* node)
{
MITK_INFO << "OnFloatingImageComboBoxSelectionChanged()";
if (m_FloatingImage.IsNotNull())
{
//TODO: Define, what will happen if the imageCT is not null...
}
if (node == nullptr)
{
this->UnsetFloatingImageGeometry();
m_FloatingImage = nullptr;
return;
}
mitk::DataNode* selectedFloatingImage = ui->floatingImageComboBox->GetSelectedNode();
if (selectedFloatingImage == nullptr)
{
this->UnsetFloatingImageGeometry();
m_FloatingImage = nullptr;
return;
}
mitk::Image::Pointer floatingImage = dynamic_cast(selectedFloatingImage->GetData());
if (floatingImage.IsNull())
{
MITK_WARN << "Failed to cast selected segmentation node to mitk::Image*";
this->UnsetFloatingImageGeometry();
m_FloatingImage = nullptr;
return;
}
m_FloatingImage = floatingImage;
this->SetFloatingImageGeometryInformation(floatingImage.GetPointer());
}
void QmitkUSNavigationStepCtUsRegistration::OnGroundTruthImageComboBoxSelectionChanged(const mitk::DataNode* node)
{
MITK_INFO << "OnGroundTruthImageComboBoxSelectionChanged()";
if (m_GroundTruthImage.IsNotNull())
{
//TODO: Define, what will happen if the imageCT is not null...
}
if (node == nullptr)
{
m_GroundTruthImage = nullptr;
return;
}
mitk::DataNode* selectedGroundTruthImage = ui->groundTruthImageComboBox->GetSelectedNode();
if (selectedGroundTruthImage == nullptr)
{
m_GroundTruthImage = nullptr;
return;
}
mitk::Image::Pointer groundTruthImage = dynamic_cast(selectedGroundTruthImage->GetData());
if (groundTruthImage.IsNull())
{
MITK_WARN << "Failed to cast selected groundTruth image node to mitk::Image*";
m_GroundTruthImage = nullptr;
return;
}
m_GroundTruthImage = groundTruthImage;
}
void QmitkUSNavigationStepCtUsRegistration::OnRegisterMarkerToFloatingImageCS()
{
this->CreateMarkerModelCoordinateSystemPointSet();
//Check for initialization
if( m_MarkerModelCoordinateSystemPointSet.IsNull() ||
m_MarkerFloatingImageCoordinateSystemPointSet.IsNull() )
{
MITK_WARN << "Fiducial Landmarks are not initialized yet, cannot register";
return;
}
//Retrieve fiducials
if (m_MarkerFloatingImageCoordinateSystemPointSet->GetSize() != m_MarkerModelCoordinateSystemPointSet->GetSize())
{
MITK_WARN << "Not the same number of fiducials, cannot register";
return;
}
else if (m_MarkerFloatingImageCoordinateSystemPointSet->GetSize() < 3)
{
MITK_WARN << "Need at least 3 fiducials, cannot register";
return;
}
//############### conversion to vtk data types (we will use the vtk landmark based transform) ##########################
//convert point sets to vtk poly data
vtkSmartPointer sourcePoints = vtkSmartPointer::New();
vtkSmartPointer targetPoints = vtkSmartPointer::New();
for (int i = 0; iGetSize(); i++)
{
double point[3] = { m_MarkerModelCoordinateSystemPointSet->GetPoint(i)[0],
m_MarkerModelCoordinateSystemPointSet->GetPoint(i)[1],
m_MarkerModelCoordinateSystemPointSet->GetPoint(i)[2] };
sourcePoints->InsertNextPoint(point);
double point_targets[3] = { m_MarkerFloatingImageCoordinateSystemPointSet->GetPoint(i)[0],
m_MarkerFloatingImageCoordinateSystemPointSet->GetPoint(i)[1],
m_MarkerFloatingImageCoordinateSystemPointSet->GetPoint(i)[2] };
targetPoints->InsertNextPoint(point_targets);
}
//########################### here, the actual transform is computed ##########################
//compute transform
vtkSmartPointer transform = vtkSmartPointer::New();
transform->SetSourceLandmarks(sourcePoints);
transform->SetTargetLandmarks(targetPoints);
transform->SetModeToRigidBody();
transform->Modified();
transform->Update();
//compute FRE of transform
double FRE = mitk::StaticIGTHelperFunctions::ComputeFRE(m_MarkerModelCoordinateSystemPointSet, m_MarkerFloatingImageCoordinateSystemPointSet, transform);
MITK_INFO << "FRE: " << FRE << " mm";
if (m_PerformingGroundTruthProtocolEvaluation)
{
m_GroundTruthProtocolFRE.push_back(FRE);
}
//#############################################################################################
//############### conversion back to itk/mitk data types ##########################
//convert from vtk to itk data types
itk::Matrix rotationFloat = itk::Matrix();
itk::Vector translationFloat = itk::Vector();
itk::Matrix rotationDouble = itk::Matrix();
itk::Vector translationDouble = itk::Vector();
vtkSmartPointer m = transform->GetMatrix();
for (int k = 0; k<3; k++) for (int l = 0; l<3; l++)
{
rotationFloat[k][l] = m->GetElement(k, l);
rotationDouble[k][l] = m->GetElement(k, l);
}
for (int k = 0; k<3; k++)
{
translationFloat[k] = m->GetElement(k, 3);
translationDouble[k] = m->GetElement(k, 3);
}
//create mitk affine transform 3D and save it to the class member
m_TransformMarkerCSToFloatingImageCS = mitk::AffineTransform3D::New();
m_TransformMarkerCSToFloatingImageCS->SetMatrix(rotationDouble);
m_TransformMarkerCSToFloatingImageCS->SetOffset(translationDouble);
MITK_INFO << m_TransformMarkerCSToFloatingImageCS;
//################################################################
//############### object is transformed ##########################
//transform surface/image
//only move image if we have one. Sometimes, this widget is used just to register point sets without images.
/*if (m_ImageNode.IsNotNull())
{
//first we have to store the original ct image transform to compose it with the new transform later
mitk::AffineTransform3D::Pointer imageTransform = m_ImageNode->GetData()->GetGeometry()->GetIndexToWorldTransform();
imageTransform->Compose(mitkTransform);
mitk::AffineTransform3D::Pointer newImageTransform = mitk::AffineTransform3D::New(); //create new image transform... setting the composed directly leads to an error
itk::Matrix rotationFloatNew = imageTransform->GetMatrix();
itk::Vector translationFloatNew = imageTransform->GetOffset();
newImageTransform->SetMatrix(rotationFloatNew);
newImageTransform->SetOffset(translationFloatNew);
m_ImageNode->GetData()->GetGeometry()->SetIndexToWorldTransform(newImageTransform);
}*/
//If this option is set, each point will be transformed and the acutal coordinates of the points change.
if( !m_PerformingGroundTruthProtocolEvaluation )
{
mitk::PointSet* pointSet_orig = m_MarkerModelCoordinateSystemPointSet;
mitk::PointSet::Pointer pointSet_moved = mitk::PointSet::New();
for (int i = 0; i < pointSet_orig->GetSize(); i++)
{
pointSet_moved->InsertPoint(m_TransformMarkerCSToFloatingImageCS->TransformPoint(pointSet_orig->GetPoint(i)));
}
pointSet_orig->Clear();
for (int i = 0; i < pointSet_moved->GetSize(); i++)
pointSet_orig->InsertPoint(pointSet_moved->GetPoint(i));
//Do a global reinit
mitk::RenderingManager::GetInstance()->InitializeViewsByBoundingObjects(this->GetDataStorage());
}
}
void QmitkUSNavigationStepCtUsRegistration::OnLocalizeFiducials()
{
m_FiducialMarkerCentroids.clear();
m_CentroidsOfFiducialCandidates.clear();
if (m_MarkerFloatingImageCoordinateSystemPointSet.IsNotNull())
{
m_MarkerFloatingImageCoordinateSystemPointSet->Clear();
}
if (!this->FilterFloatingImage())
{
QMessageBox msgBox;
msgBox.setText("Cannot perform filtering of the image. The floating image = nullptr.");
msgBox.exec();
return;
}
this->GetCentroidsOfLabeledObjects();
if (!this->EliminateFiducialCandidatesByEuclideanDistances() ||
m_CentroidsOfFiducialCandidates.size() != NUMBER_FIDUCIALS_NEEDED)
{
QMessageBox msgBox;
QString text = QString("Have found %1 instead of 8 fiducial candidates.\
Cannot perform fiducial localization procedure.").arg(m_CentroidsOfFiducialCandidates.size());
msgBox.setText(text);
msgBox.exec();
return;
}
//Before calling NumerateFiducialMarks it must be sure,
// that there rested only 8 fiducial candidates.
this->NumerateFiducialMarks();
}
void QmitkUSNavigationStepCtUsRegistration::OnFilterGroundTruthImage()
{
/*if (m_GroundTruthImage.IsNull())
{
return;
}*/
this->ShowGroundTruthMarkerEdges();
}
void QmitkUSNavigationStepCtUsRegistration::OnAddCtImageClicked()
{
mitk::DataNode* selectedCtImage = ui->ctImagesToChooseComboBox->GetSelectedNode();
if (selectedCtImage == nullptr)
{
return;
}
mitk::Image::Pointer ctImage = dynamic_cast(selectedCtImage->GetData());
if (ctImage.IsNull())
{
MITK_WARN << "Failed to cast selected segmentation node to mitk::Image*";
return;
}
QString name = QString::fromStdString(selectedCtImage->GetName());
for( int counter = 0; counter < ui->chosenCtImagesListWidget->count(); ++counter)
{
MITK_INFO << ui->chosenCtImagesListWidget->item(counter)->text() << " - " << counter;
MITK_INFO << m_ImagesGroundTruthProtocol.at(counter).GetPointer();
if (ui->chosenCtImagesListWidget->item(counter)->text().compare(name) == 0)
{
MITK_INFO << "CT image already exist in list of chosen CT images. Do not add the image.";
return;
}
}
ui->chosenCtImagesListWidget->addItem(name);
m_ImagesGroundTruthProtocol.push_back(ctImage);
}
void QmitkUSNavigationStepCtUsRegistration::OnRemoveCtImageClicked()
{
int position = ui->chosenCtImagesListWidget->currentRow();
if (ui->chosenCtImagesListWidget->count() == 0 || position < 0)
{
return;
}
m_ImagesGroundTruthProtocol.erase(m_ImagesGroundTruthProtocol.begin() + position);
QListWidgetItem *item = ui->chosenCtImagesListWidget->currentItem();
ui->chosenCtImagesListWidget->removeItemWidget(item);
delete item;
}
void QmitkUSNavigationStepCtUsRegistration::OnEvaluateGroundTruthFiducialLocalizationProtocol()
{
m_GroundTruthProtocolFRE.clear();
if (m_ImagesGroundTruthProtocol.size() != 6)
{
QMessageBox msgBox;
msgBox.setText("For evaluating the Ground-Truth-Fiducial-Localization-Protocol there must be loaded 6 different CT images.");
msgBox.exec();
return;
}
m_PerformingGroundTruthProtocolEvaluation = true;
this->CreatePointsToTransformForGroundTruthProtocol();
m_GroundTruthProtocolTransformedPoints.clear();
for (int cycleNo = 0; cycleNo < m_ImagesGroundTruthProtocol.size(); ++cycleNo)
{
m_FloatingImage = m_ImagesGroundTruthProtocol.at(cycleNo);
this->SetFloatingImageGeometryInformation(m_FloatingImage.GetPointer());
this->OnLocalizeFiducials();
this->OnRegisterMarkerToFloatingImageCS();
this->TransformPointsGroundTruthProtocol();
}
this->AddTransformedPointsToDataStorage();
double meanFRE = this->CalculateMeanFRE();
double sdOfFRE = this->CalculateStandardDeviationOfFRE(meanFRE);
this->CalculateGroundTruthProtocolTRE();
ui->meanFREValue->setText(QString("%1").arg(meanFRE));
ui->sdFREValue->setText(QString("%1").arg(sdOfFRE));
if (ui->protocolEvaluationTypeComboBox->currentText().compare("ANGLE") == 0)
{
if (m_GroundTruthProtocolTRE.find(0) != m_GroundTruthProtocolTRE.end())
{
ui->TREValue->setText(QString("%1").arg(m_GroundTruthProtocolTRE.at(0)));
}
}
else if (ui->protocolEvaluationTypeComboBox->currentText().compare("PLANE") == 0)
{
if (m_GroundTruthProtocolTRE.find(0) != m_GroundTruthProtocolTRE.end() &&
m_GroundTruthProtocolTRE.find(20) != m_GroundTruthProtocolTRE.end() &&
m_GroundTruthProtocolTRE.find(40) != m_GroundTruthProtocolTRE.end() &&
m_GroundTruthProtocolTRE.find(60) != m_GroundTruthProtocolTRE.end() &&
m_GroundTruthProtocolTRE.find(80) != m_GroundTruthProtocolTRE.end() &&
m_GroundTruthProtocolTRE.find(100) != m_GroundTruthProtocolTRE.end())
{
ui->TREValue->setText(QString("Depth 0mm: %1\nDepth 20mm: %2\nDepth 40mm: %3\
\nDepth 60mm: %4\nDepth 80mm: %5\nDepth 100mm: %6")
.arg(m_GroundTruthProtocolTRE.at(0))
.arg(m_GroundTruthProtocolTRE.at(20))
.arg(m_GroundTruthProtocolTRE.at(40))
.arg(m_GroundTruthProtocolTRE.at(60))
.arg(m_GroundTruthProtocolTRE.at(80))
.arg(m_GroundTruthProtocolTRE.at(100)));
}
}
m_PerformingGroundTruthProtocolEvaluation = false;
}
diff --git a/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/QmitkUltrasoundCalibration.cpp b/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/QmitkUltrasoundCalibration.cpp
index 019828609b..cf6efdaeef 100644
--- a/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/QmitkUltrasoundCalibration.cpp
+++ b/Plugins/org.mitk.gui.qt.igt.app.echotrack/src/internal/QmitkUltrasoundCalibration.cpp
@@ -1,1158 +1,1164 @@
/*===================================================================
The Medical Imaging Interaction Toolkit (MITK)
Copyright (c) German Cancer Research Center,
Division of Medical and Biological Informatics.
All rights reserved.
This software is distributed WITHOUT ANY WARRANTY; without
even the implied warranty of MERCHANTABILITY or FITNESS FOR
A PARTICULAR PURPOSE.
See LICENSE.txt or http://www.mitk.org for details.
===================================================================*/
// Blueberry
#include
#include
// Qmitk
#include "QmitkUltrasoundCalibration.h"
#include
// Qt
#include
#include
#include
#include
// MITK
#include
//#include
#include
#include
#include
#include
#include
#include "mitkIRenderingManager.h"
// us
#include
//VTK
#include
#include
#include
#include
#include
#include
#include "internal/org_mbi_gui_qt_usnavigation_Activator.h"
//sleep headers
#include
#include
const std::string QmitkUltrasoundCalibration::VIEW_ID = "org.mitk.views.ultrasoundcalibration";
QmitkUltrasoundCalibration::QmitkUltrasoundCalibration() :
m_USDeviceChanged(this, &QmitkUltrasoundCalibration::OnUSDepthChanged)
{
ctkPluginContext* pluginContext = mitk::PluginActivator::GetContext();
if (pluginContext)
{
// to be notified about service event of an USDevice
pluginContext->connectServiceListener(this, "OnDeviceServiceEvent",
QString::fromStdString("(" + us::ServiceConstants::OBJECTCLASS() + "=" + us_service_interface_iid() + ")"));
}
}
QmitkUltrasoundCalibration::~QmitkUltrasoundCalibration()
{
m_Controls.m_CombinedModalityManagerWidget->blockSignals(true);
mitk::AbstractUltrasoundTrackerDevice::Pointer combinedModality;
combinedModality = m_Controls.m_CombinedModalityManagerWidget->GetSelectedCombinedModality();
if (combinedModality.IsNotNull())
{
combinedModality->GetUltrasoundDevice()->RemovePropertyChangedListener(m_USDeviceChanged);
}
m_Timer->stop();
// Sleep(500); //This might be problematic... seems like sometimes some ressources are still in use at calling time.
this->OnStopCalibrationProcess();
this->OnStopPlusCalibration();
/*mitk::DataNode::Pointer node = this->GetDataStorage()->GetNamedNode("Tool Calibration Points");
if (node.IsNotNull())this->GetDataStorage()->Remove(node);
node = this->GetDataStorage()->GetNamedNode("Image Calibration Points");
if (node.IsNotNull())this->GetDataStorage()->Remove(node);
node = this->GetDataStorage()->GetNamedNode("US Image Stream");
if (node.IsNotNull())this->GetDataStorage()->Remove(node);*/
mitk::DataNode::Pointer node = this->GetDataStorage()->GetNamedNode("Needle Path");
if (node.IsNotNull())this->GetDataStorage()->Remove(node);
this->GetDataStorage()->Remove(m_VerificationReferencePointsDataNode);
delete m_Timer;
}
void QmitkUltrasoundCalibration::SetFocus()
{
m_Controls.m_ToolBox->setFocus();
}
void QmitkUltrasoundCalibration::CreateQtPartControl(QWidget *parent)
{
// create GUI widgets from the Qt Designer's .ui file
m_Controls.setupUi(parent);
m_Controls.m_CombinedModalityManagerWidget->SetCalibrationLoadedNecessary(false);
m_Timer = new QTimer(this);
m_StreamingTimer = new QTimer(this);
m_Controls.m_SpacingBtnFreeze->setEnabled(true);
m_Controls.m_SpacingAddPoint->setEnabled(false);
m_Controls.m_CalculateSpacing->setEnabled(false);
m_SpacingPointsCount = 0;
m_SpacingPoints = mitk::PointSet::New();
m_SpacingNode = mitk::DataNode::New();
m_SpacingNode->SetName("Spacing Points");
m_SpacingNode->SetData(this->m_SpacingPoints);
this->GetDataStorage()->Add(m_SpacingNode);
// Pointset for Calibration Points
m_CalibPointsTool = mitk::PointSet::New();
// Pointset for Worldpoints
m_CalibPointsImage = mitk::PointSet::New();
m_CalibPointsCount = 0;
// Evaluation Pointsets (Non-Visualized)
m_EvalPointsImage = mitk::PointSet::New();
m_EvalPointsTool = mitk::PointSet::New();
m_EvalPointsProjected = mitk::PointSet::New();
// Neelde Projection Filter
m_NeedleProjectionFilter = mitk::NeedleProjectionFilter::New();
// Tracking Status Widgets
m_Controls.m_CalibTrackingStatus->ShowStatusLabels();
m_Controls.m_EvalTrackingStatus->ShowStatusLabels();
// General & Device Selection
connect(m_Timer, SIGNAL(timeout()), this, SLOT(Update()));
//connect(m_Controls.m_ToolBox, SIGNAL(currentChanged(int)), this, SLOT(OnTabSwitch(int)));
// Calibration
connect(m_Controls.m_CalibBtnFreeze, SIGNAL(clicked()), this, SLOT(SwitchFreeze())); // Freeze
connect(m_Controls.m_CalibBtnAddPoint, SIGNAL(clicked()), this, SLOT(OnAddCalibPoint())); // Tracking & Image Points (Calibration)
connect(m_Controls.m_CalibBtnCalibrate, SIGNAL(clicked()), this, SLOT(OnCalibration())); // Perform Calibration
// Evaluation
connect(m_Controls.m_EvalBtnStep1, SIGNAL(clicked()), this, SLOT(OnAddEvalProjectedPoint())); // Needle Projection
connect(m_Controls.m_EvalBtnStep2, SIGNAL(clicked()), this, SLOT(SwitchFreeze())); // Freeze
connect(m_Controls.m_EvalBtnStep3, SIGNAL(clicked()), this, SLOT(OnAddEvalTargetPoint())); // Tracking & Image Points (Evaluation)
connect(m_Controls.m_EvalBtnSave, SIGNAL(clicked()), this, SLOT(OnSaveEvaluation())); // Save Evaluation Results
connect(m_Controls.m_CalibBtnSaveCalibration, SIGNAL(clicked()), this, SLOT(OnSaveCalibration())); // Save Evaluation Results
connect(m_Controls.m_BtnReset, SIGNAL(clicked()), this, SLOT(OnReset())); // Reset Pointsets
// PLUS Calibration
connect(m_Controls.m_GetCalibrationFromPLUS, SIGNAL(clicked()), this, SLOT(OnGetPlusCalibration()));
connect(m_Controls.m_StartStreaming, SIGNAL(clicked()), this, SLOT(OnStartStreaming()));
connect(m_StreamingTimer, SIGNAL(timeout()), this, SLOT(OnStreamingTimerTimeout()));
connect(m_Controls.m_StopPlusCalibration, SIGNAL(clicked()), this, SLOT(OnStopPlusCalibration()));
connect(m_Controls.m_SavePlusCalibration, SIGNAL(clicked()), this, SLOT(OnSaveCalibration()));
connect(this, SIGNAL(NewConnectionSignal()), this, SLOT(OnNewConnection()));
//Determine Spacing for Calibration of USVideoDevice
connect(m_Controls.m_SpacingBtnFreeze, SIGNAL(clicked()), this, SLOT(OnFreezeClicked()));
connect(m_Controls.m_SpacingAddPoint, SIGNAL(clicked()), this, SLOT(OnAddSpacingPoint()));
connect(m_Controls.m_CalculateSpacing, SIGNAL(clicked()), this, SLOT(OnCalculateSpacing()));
//connect( m_Controls.m_CombinedModalityManagerWidget, SIGNAL(SignalCombinedModalitySelected(mitk::USCombinedModality::Pointer)),
// this, SLOT(OnSelectDevice(mitk::USCombinedModality::Pointer)) );
connect(m_Controls.m_CombinedModalityManagerWidget, SIGNAL(SignalReadyForNextStep()),
this, SLOT(OnDeviceSelected()));
connect(m_Controls.m_CombinedModalityManagerWidget, SIGNAL(SignalNoLongerReadyForNextStep()),
this, SLOT(OnDeviceDeselected()));
connect(m_Controls.m_StartCalibrationButton, SIGNAL(clicked()), this, SLOT(OnStartCalibrationProcess()));
connect(m_Controls.m_StartPlusCalibrationButton, SIGNAL(clicked()), this, SLOT(OnStartPlusCalibration()));
connect(m_Controls.m_CalibBtnRestartCalibration, SIGNAL(clicked()), this, SLOT(OnReset()));
connect(m_Controls.m_CalibBtnStopCalibration, SIGNAL(clicked()), this, SLOT(OnStopCalibrationProcess()));
connect(m_Controls.m_AddReferencePoints, SIGNAL(clicked()), this, SLOT(OnAddCurrentTipPositionToReferencePoints()));
connect(m_Controls.m_AddCurrentPointerTipForVerification, SIGNAL(clicked()), this, SLOT(OnAddCurrentTipPositionForVerification()));
connect(m_Controls.m_StartVerification, SIGNAL(clicked()), this, SLOT(OnStartVerification()));
//initialize data storage combo box
m_Controls.m_ReferencePointsComboBox->SetDataStorage(this->GetDataStorage());
m_Controls.m_ReferencePointsComboBox->SetAutoSelectNewItems(true);
m_Controls.m_ReferencePointsComboBox->SetPredicate(mitk::NodePredicateDataType::New("PointSet"));
//initialize point list widget
if (m_VerificationReferencePoints.IsNull()) { m_VerificationReferencePoints = mitk::PointSet::New(); }
if (m_VerificationReferencePointsDataNode.IsNull())
{
m_VerificationReferencePointsDataNode = mitk::DataNode::New();
m_VerificationReferencePointsDataNode->SetName("US Verification Reference Points");
m_VerificationReferencePointsDataNode->SetData(m_VerificationReferencePoints);
this->GetDataStorage()->Add(m_VerificationReferencePointsDataNode);
}
m_Controls.m_ReferencePointsPointListWidget->SetPointSetNode(m_VerificationReferencePointsDataNode);
m_Controls.m_ToolBox->setCurrentIndex(0);
}
void QmitkUltrasoundCalibration::OnSelectionChanged(berry::IWorkbenchPart::Pointer /*source*/,
const QList& /*nodes*/)
{
}
void QmitkUltrasoundCalibration::OnTabSwitch(int index)
{
switch (index)
{
case 0:
if (m_Controls.m_ToolBox->isItemEnabled(1) || m_Controls.m_ToolBox->isItemEnabled(2))
{
this->OnStopCalibrationProcess();
}
break;
default:
;
}
}
//void QmitkUltrasoundCalibration::OnSelectDevice(mitk::USCombinedModality::Pointer combinedModality)
void QmitkUltrasoundCalibration::OnDeviceSelected()
{
mitk::AbstractUltrasoundTrackerDevice::Pointer combinedModality;
combinedModality = m_Controls.m_CombinedModalityManagerWidget->GetSelectedCombinedModality();
if (combinedModality.IsNotNull())
{
//m_Tracker = m_CombinedModality->GetNavigationDataSource();
// Construct Pipeline
//this->m_NeedleProjectionFilter->SetInput(0, m_Tracker->GetOutput(0));
combinedModality->GetUltrasoundDevice()->AddPropertyChangedListener(m_USDeviceChanged);
m_Controls.m_StartCalibrationButton->setEnabled(true);
m_Controls.m_StartPlusCalibrationButton->setEnabled(true);
m_Controls.m_ToolBox->setItemEnabled(1, true);
m_Controls.m_ToolBox->setItemEnabled(2, true);
}
}
void QmitkUltrasoundCalibration::OnDeviceDeselected()
{
mitk::AbstractUltrasoundTrackerDevice::Pointer combinedModality;
combinedModality = m_Controls.m_CombinedModalityManagerWidget->GetSelectedCombinedModality();
if (combinedModality.IsNotNull())
{
combinedModality->GetUltrasoundDevice()->RemovePropertyChangedListener(m_USDeviceChanged);
}
m_Controls.m_StartCalibrationButton->setEnabled(false);
m_Controls.m_StartPlusCalibrationButton->setEnabled(false);
m_Controls.m_ToolBox->setCurrentIndex(0);
m_Controls.m_ToolBox->setItemEnabled(1, false);
m_Controls.m_ToolBox->setItemEnabled(2, false);
}
void QmitkUltrasoundCalibration::OnAddCurrentTipPositionToReferencePoints()
{
if (m_Controls.m_VerificationPointerChoser->GetSelectedNavigationDataSource().IsNull() ||
(m_Controls.m_VerificationPointerChoser->GetSelectedToolID() == -1))
{
MITK_WARN << "No tool selected, aborting";
return;
}
mitk::NavigationData::Pointer currentPointerData = m_Controls.m_VerificationPointerChoser->GetSelectedNavigationDataSource()->GetOutput(m_Controls.m_VerificationPointerChoser->GetSelectedToolID());
mitk::Point3D currentTipPosition = currentPointerData->GetPosition();
m_VerificationReferencePoints->InsertPoint(m_VerificationReferencePoints->GetSize(), currentTipPosition);
}
void QmitkUltrasoundCalibration::OnStartVerification()
{
m_currentPoint = 0;
mitk::PointSet::Pointer selectedPointSet = dynamic_cast(m_Controls.m_ReferencePointsComboBox->GetSelectedNode()->GetData());
m_Controls.m_CurrentPointLabel->setText("Point " + QString::number(m_currentPoint) + " of " + QString::number(selectedPointSet->GetSize()));
m_allErrors = std::vector();
m_allReferencePoints = std::vector();
for (int i = 0; i < selectedPointSet->GetSize(); i++)
{
m_allReferencePoints.push_back(selectedPointSet->GetPoint(i));
}
}
void QmitkUltrasoundCalibration::OnAddCurrentTipPositionForVerification()
{
if (m_currentPoint == -1) { MITK_WARN << "Cannot add point"; return; }
if (m_Controls.m_VerificationPointerChoser->GetSelectedNavigationDataSource().IsNull() ||
(m_Controls.m_VerificationPointerChoser->GetSelectedToolID() == -1))
{
MITK_WARN << "No tool selected, aborting";
return;
}
mitk::NavigationData::Pointer currentPointerData = m_Controls.m_VerificationPointerChoser->GetSelectedNavigationDataSource()->GetOutput(m_Controls.m_VerificationPointerChoser->GetSelectedToolID());
mitk::Point3D currentTipPosition = currentPointerData->GetPosition();
double currentError = m_allReferencePoints.at(m_currentPoint).EuclideanDistanceTo(currentTipPosition);
MITK_INFO << "Current Error: " << currentError << " mm";
m_allErrors.push_back(currentError);
if (++m_currentPoint < static_cast(m_allReferencePoints.size()))
{
m_Controls.m_CurrentPointLabel->setText("Point " + QString::number(m_currentPoint) + " of " + QString::number(m_allReferencePoints.size()));
}
else
{
m_currentPoint = -1;
double meanError = 0;
for (std::size_t i = 0; i < m_allErrors.size(); ++i)
{
meanError += m_allErrors[i];
}
meanError /= m_allErrors.size();
QString result = "Finished verification! \n Verification of " + QString::number(m_allErrors.size()) + " points, mean error: " + QString::number(meanError) + " mm";
m_Controls.m_ResultsTextEdit->setText(result);
MITK_INFO << result.toStdString();
}
}
void QmitkUltrasoundCalibration::OnStartCalibrationProcess()
{
// US Image Stream
m_Node = mitk::DataNode::New();
m_Node->SetName("US Calibration Viewing Stream");
//create a dummy image (gray values 0..255) for correct initialization of level window, etc.
mitk::Image::Pointer dummyImage = mitk::ImageGenerator::GenerateRandomImage(100, 100, 1, 1, 1, 1, 1, 255, 0);
m_Node->SetData(dummyImage);
this->GetDataStorage()->Add(m_Node);
// data node for calibration point set
m_CalibNode = mitk::DataNode::New();
m_CalibNode->SetName("Tool Calibration Points");
m_CalibNode->SetData(this->m_CalibPointsImage);
this->GetDataStorage()->Add(m_CalibNode);
// data node for world point set
m_WorldNode = mitk::DataNode::New();
m_WorldNode->SetName("Image Calibration Points");
m_WorldNode->SetData(this->m_CalibPointsTool);
this->GetDataStorage()->Add(m_WorldNode);
m_CombinedModality = m_Controls.m_CombinedModalityManagerWidget->GetSelectedCombinedModality();
m_CombinedModality->SetCalibration(mitk::AffineTransform3D::New()); //dummy calibration because without a calibration the comined modality was laggy (maybe a bug?)
if (m_CombinedModality.IsNull()) { return; }
m_Tracker = m_CombinedModality->GetNavigationDataSource();
//QString curDepth = service.getProperty(QString::fromStdString(mitk::USDevice::US_PROPKEY_BMODE_DEPTH)).toString();
// Construct Pipeline
this->m_NeedleProjectionFilter->SetInput(0, m_Tracker->GetOutput(0));
QApplication::setOverrideCursor(Qt::WaitCursor);
// make sure that the combined modality is in connected state before using it
if (m_CombinedModality->GetUltrasoundDevice()->GetDeviceState() < mitk::USDevice::State_Connected) { m_CombinedModality->GetUltrasoundDevice()->Connect(); }
if (m_CombinedModality->GetUltrasoundDevice()->GetDeviceState() < mitk::USDevice::State_Activated) { m_CombinedModality->GetUltrasoundDevice()->Activate(); }
QApplication::restoreOverrideCursor();
this->SwitchFreeze();
//Trigger the ProbeChanged method for initializing/updating the spacing of the ultrasound image correctly
std::string probeName = m_CombinedModality->GetUltrasoundDevice()->GetCurrentProbe()->GetName();
m_CombinedModality->GetUltrasoundDevice()->ProbeChanged(probeName);
+ mitk::DataNode::Pointer usNode = this->GetDataStorage()->GetNamedNode("US Viewing Stream - Image 0");
+ if (usNode.IsNotNull())
+ {
+ this->GetDataStorage()->Remove(usNode);
+ }
+
// Todo: Maybe display this elsewhere
this->ShowNeedlePath();
// Switch active tab to Calibration page
m_Controls.m_ToolBox->setItemEnabled(1, true);
m_Controls.m_ToolBox->setCurrentIndex(1);
}
void QmitkUltrasoundCalibration::OnStartPlusCalibration()
{
if (m_CombinedModality.IsNull()){
m_CombinedModality = m_Controls.m_CombinedModalityManagerWidget->GetSelectedCombinedModality();
if (m_CombinedModality.IsNull()) { return; } //something went wrong, there is no combined modality
}
//setup server to send UltrasoundImages to PLUS
mitk::IGTLServer::Pointer m_USServer = mitk::IGTLServer::New(true);
m_USServer->SetName("EchoTrack Image Source");
m_USServer->SetHostname("127.0.0.1");
m_USServer->SetPortNumber(18944);
m_USMessageProvider = mitk::IGTLMessageProvider::New();
m_USMessageProvider->SetIGTLDevice(m_USServer);
m_USMessageProvider->SetFPS(5);
m_USImageToIGTLMessageFilter = mitk::ImageToIGTLMessageFilter::New();
m_USImageToIGTLMessageFilter->ConnectTo(m_CombinedModality->GetUltrasoundDevice());
m_USImageToIGTLMessageFilter->SetName("USImage Filter");
//setup server to send TrackingData to PLUS
m_TrackingServer = mitk::IGTLServer::New(true);
m_TrackingServer->SetName("EchoTrack Tracking Source");
m_TrackingServer->SetHostname("127.0.0.1");
m_TrackingServer->SetPortNumber(18945);
m_TrackingMessageProvider = mitk::IGTLMessageProvider::New();
m_TrackingMessageProvider->SetIGTLDevice(m_TrackingServer);
m_TrackingMessageProvider->SetFPS(5);
m_TrackingToIGTLMessageFilter = mitk::NavigationDataToIGTLMessageFilter::New();
m_TrackingToIGTLMessageFilter->ConnectTo(m_CombinedModality->GetTrackingDeviceDataSource());
m_TrackingToIGTLMessageFilter->SetName("Tracker Filter");
typedef itk::SimpleMemberCommand< QmitkUltrasoundCalibration > CurCommandType;
CurCommandType::Pointer newConnectionCommand = CurCommandType::New();
newConnectionCommand->SetCallbackFunction(
this, &QmitkUltrasoundCalibration::OnPlusConnected);
this->m_NewConnectionObserverTag = this->m_TrackingServer->AddObserver(
mitk::NewClientConnectionEvent(), newConnectionCommand);
//Open connections of both servers
if (m_USServer->OpenConnection())
{
MITK_INFO << "US Server opened its connection successfully";
m_USServer->StartCommunication();
}
else
{
MITK_INFO << "US Server could not open its connection";
}
if (m_TrackingServer->OpenConnection())
{
MITK_INFO << "Tracking Server opened its connection successfully";
m_TrackingServer->StartCommunication();
}
else
{
MITK_INFO << "Tracking Server could not open its connection";
}
if (m_USMessageProvider->IsCommunicating() && m_TrackingMessageProvider->IsCommunicating())
{
m_Controls.m_StartPlusCalibrationButton->setEnabled(false);
m_Controls.m_GetCalibrationFromPLUS->setEnabled(true);
m_Controls.m_StartStreaming->setEnabled(false);
m_Controls.m_SavePlusCalibration->setEnabled(false);
m_Controls.m_SetupStatus->setStyleSheet("QLabel { color : green; }");
m_Controls.m_SetupStatus->setText("Setup successfull you can now connect PLUS");
}
else
{
m_Controls.m_SetupStatus->setStyleSheet("QLabel { color : red; }");
m_Controls.m_SetupStatus->setText("Something went wrong. Please try again");
}
}
void QmitkUltrasoundCalibration::OnStopPlusCalibration()
{
//closing all server and clients when PlusCalibration is finished
if (m_USMessageProvider.IsNotNull())
{
if (m_USMessageProvider->IsStreaming())
{
m_USMessageProvider->StopStreamingOfSource(m_USImageToIGTLMessageFilter);
}
}
if (m_TrackingMessageProvider.IsNotNull())
{
if (m_TrackingMessageProvider->IsStreaming())
{
m_TrackingMessageProvider->StopStreamingOfSource(m_TrackingToIGTLMessageFilter);
}
}
if (m_USServer.IsNotNull())
{
m_USServer->CloseConnection();
}
if (m_TrackingServer.IsNotNull())
{
m_TrackingServer->CloseConnection();
}
if (m_TransformClient.IsNotNull())
{
m_TransformClient->CloseConnection();
}
m_Controls.m_GotCalibrationLabel->setText("");
m_Controls.m_ConnectionStatus->setText("");
m_Controls.m_SetupStatus->setText("");
m_Controls.m_StartPlusCalibrationButton->setEnabled(true);
m_StreamingTimer->stop();
delete m_StreamingTimer;
}
void QmitkUltrasoundCalibration::OnPlusConnected()
{
emit NewConnectionSignal();
}
void QmitkUltrasoundCalibration::OnNewConnection()
{
m_Controls.m_StartStreaming->setEnabled(true);
m_Controls.m_ConnectionStatus->setStyleSheet("QLabel { color : green; }");
m_Controls.m_ConnectionStatus->setText("Connection successfull you can now start streaming");
}
void QmitkUltrasoundCalibration::OnStreamingTimerTimeout()
{
m_USMessageProvider->Update();
m_TrackingMessageProvider->Update();
}
void QmitkUltrasoundCalibration::OnStartStreaming()
{
m_USMessageProvider->StartStreamingOfSource(m_USImageToIGTLMessageFilter, 5);
m_TrackingMessageProvider->StartStreamingOfSource(m_TrackingToIGTLMessageFilter, 5);
m_Controls.m_StartStreaming->setEnabled(false);
m_Controls.m_ConnectionStatus->setText("");
m_StreamingTimer->start((1.0 / 5.0 * 1000.0));
}
void QmitkUltrasoundCalibration::OnGetPlusCalibration()
{
m_TransformClient = mitk::IGTLClient::New(true);
m_TransformClient->SetHostname("127.0.0.1");
m_TransformClient->SetPortNumber(18946);
m_TransformDeviceSource = mitk::IGTLDeviceSource::New();
m_TransformDeviceSource->SetIGTLDevice(m_TransformClient);
m_TransformDeviceSource->Connect();
if (m_TransformDeviceSource->IsConnected())
{
MITK_INFO << "successfully connected";
m_TransformDeviceSource->StartCommunication();
if (m_TransformDeviceSource->IsCommunicating())
{
MITK_INFO << "communication started";
mitk::IGTLMessage::Pointer receivedMessage;
bool condition = false;
igtl::Matrix4x4 transformPLUS;
while (!(receivedMessage.IsNotNull() && receivedMessage->IsDataValid()))
{
std::this_thread::sleep_for(std::chrono::milliseconds(50));
m_TransformDeviceSource->Update();
receivedMessage = m_TransformDeviceSource->GetOutput();
igtl::TransformMessage::Pointer msg = dynamic_cast(m_TransformDeviceSource->GetOutput()->GetMessage().GetPointer());
if (msg == nullptr || msg.IsNull())
{
MITK_INFO << "Received message could not be casted to TransformMessage. Skipping..";
continue;
}
else
{
if (std::strcmp(msg->GetDeviceName(), "ImageToTracker") != 0)
{
MITK_INFO << "Was not Image to Tracker Transform. Skipping...";
continue;
}
else
{
msg->GetMatrix(transformPLUS);
condition = true;
break;
}
}
}
if (condition)
{
this->ProcessPlusCalibration(transformPLUS);
}
else
{
m_Controls.m_GotCalibrationLabel->setStyleSheet("QLabel { color : red; }");
m_Controls.m_GotCalibrationLabel->setText("Something went wrong. Please try again");
}
}
else
{
MITK_INFO << " no connection";
m_Controls.m_GotCalibrationLabel->setStyleSheet("QLabel { color : red; }");
m_Controls.m_GotCalibrationLabel->setText("Something went wrong. Please try again");
}
}
else
{
m_Controls.m_GotCalibrationLabel->setStyleSheet("QLabel { color : red; }");
m_Controls.m_GotCalibrationLabel->setText("Something went wrong. Please try again");
}
}
void QmitkUltrasoundCalibration::ProcessPlusCalibration(igtl::Matrix4x4& imageToTracker)
{
mitk::AffineTransform3D::Pointer imageToTrackerTransform = mitk::AffineTransform3D::New();
itk::Matrix rotationFloat = itk::Matrix();
itk::Vector translationFloat = itk::Vector();
rotationFloat[0][0] = imageToTracker[0][0];
rotationFloat[0][1] = imageToTracker[0][1];
rotationFloat[0][2] = imageToTracker[0][2];
rotationFloat[1][0] = imageToTracker[1][0];
rotationFloat[1][1] = imageToTracker[1][1];
rotationFloat[1][2] = imageToTracker[1][2];
rotationFloat[2][0] = imageToTracker[2][0];
rotationFloat[2][1] = imageToTracker[2][1];
rotationFloat[2][2] = imageToTracker[2][2];
translationFloat[0] = imageToTracker[0][3];
translationFloat[1] = imageToTracker[1][3];
translationFloat[2] = imageToTracker[2][3];
imageToTrackerTransform->SetTranslation(translationFloat);
imageToTrackerTransform->SetMatrix(rotationFloat);
m_CombinedModality->SetCalibration(imageToTrackerTransform);
m_Controls.m_ToolBox->setItemEnabled(2, true);
m_Controls.m_SavePlusCalibration->setEnabled(true);
m_Controls.m_GotCalibrationLabel->setStyleSheet("QLabel { color : green; }");
m_Controls.m_GotCalibrationLabel->setText("Recieved Calibration from PLUS you can now save it");
}
void QmitkUltrasoundCalibration::OnStopCalibrationProcess()
{
this->ClearTemporaryMembers();
m_Timer->stop();
this->GetDataStorage()->Remove(m_Node);
m_Node = 0;
this->GetDataStorage()->Remove(m_CalibNode);
m_CalibNode = 0;
this->GetDataStorage()->Remove(m_WorldNode);
m_WorldNode = 0;
m_Controls.m_ToolBox->setCurrentIndex(0);
}
void QmitkUltrasoundCalibration::OnDeviceServiceEvent(const ctkServiceEvent event)
{
if (m_CombinedModality.IsNull() || event.getType() != ctkServiceEvent::MODIFIED) { return; }
ctkServiceReference service = event.getServiceReference();
QString curDepth = service.getProperty(QString::fromStdString(mitk::USDevice::GetPropertyKeys().US_PROPKEY_BMODE_DEPTH)).toString();
if (m_CurrentDepth != curDepth)
{
m_CurrentDepth = curDepth;
this->OnReset();
}
}
void QmitkUltrasoundCalibration::OnAddCalibPoint()
{
mitk::Point3D world = this->GetRenderWindowPart()->GetSelectedPosition();
this->m_CalibPointsImage->InsertPoint(m_CalibPointsCount, world);
this->m_CalibPointsTool->InsertPoint(m_CalibPointsCount, this->m_FreezePoint);
QString text = text.number(m_CalibPointsCount + 1);
text = "Point " + text;
this->m_Controls.m_CalibPointList->addItem(text);
m_CalibPointsCount++;
SwitchFreeze();
}
void QmitkUltrasoundCalibration::OnCalibration()
{
// Compute transformation
vtkSmartPointer transform = vtkSmartPointer::New();
transform->SetSourceLandmarks(this->ConvertPointSetToVtkPolyData(m_CalibPointsImage)->GetPoints());
transform->SetTargetLandmarks(this->ConvertPointSetToVtkPolyData(m_CalibPointsTool)->GetPoints());
if( !m_CombinedModality->GetIsTrackedUltrasoundActive() )
{
if (m_Controls.m_ScaleTransform->isChecked())
{
transform->SetModeToSimilarity();
} //use affine transform
else
{
transform->SetModeToRigidBody();
} //use similarity transform: scaling is not touched
MITK_INFO << "TEST";
}
else
{
transform->SetModeToRigidBody();//use similarity transform: scaling is not touched
}
transform->Modified();
transform->Update();
// Convert from vtk to itk data types
itk::Matrix rotationFloat = itk::Matrix();
itk::Vector translationFloat = itk::Vector();
vtkSmartPointer m = transform->GetMatrix();
rotationFloat[0][0] = m->GetElement(0, 0);
rotationFloat[0][1] = m->GetElement(0, 1);
rotationFloat[0][2] = m->GetElement(0, 2);
rotationFloat[1][0] = m->GetElement(1, 0);
rotationFloat[1][1] = m->GetElement(1, 1);
rotationFloat[1][2] = m->GetElement(1, 2);
rotationFloat[2][0] = m->GetElement(2, 0);
rotationFloat[2][1] = m->GetElement(2, 1);
rotationFloat[2][2] = m->GetElement(2, 2);
translationFloat[0] = m->GetElement(0, 3);
translationFloat[1] = m->GetElement(1, 3);
translationFloat[2] = m->GetElement(2, 3);
mitk::DataNode::Pointer CalibPointsImage = mitk::DataNode::New();
CalibPointsImage->SetName("Calibration Points Image");
CalibPointsImage->SetData(m_CalibPointsImage);
this->GetDataStorage()->Add(CalibPointsImage);
mitk::DataNode::Pointer CalibPointsTracking = mitk::DataNode::New();
CalibPointsTracking->SetName("Calibration Points Tracking");
CalibPointsTracking->SetData(m_CalibPointsTool);
this->GetDataStorage()->Add(CalibPointsTracking);
mitk::PointSet::Pointer ImagePointsTransformed = m_CalibPointsImage->Clone();
this->ApplyTransformToPointSet(ImagePointsTransformed, transform);
mitk::DataNode::Pointer CalibPointsImageTransformed = mitk::DataNode::New();
CalibPointsImageTransformed->SetName("Calibration Points Image (Transformed)");
CalibPointsImageTransformed->SetData(ImagePointsTransformed);
this->GetDataStorage()->Add(CalibPointsImageTransformed);
// Set output variable
mitk::AffineTransform3D::Pointer oldUSImageTransform = m_CombinedModality->GetUltrasoundDevice()->GetOutput()->GetGeometry()->GetIndexToWorldTransform(); //including spacing!
MITK_INFO << "Old US Image transform: " << oldUSImageTransform;
mitk::AffineTransform3D::Pointer calibTransform = mitk::AffineTransform3D::New();
calibTransform->SetTranslation(translationFloat);
calibTransform->SetMatrix(rotationFloat);
MITK_INFO << "Calibration transform: " << calibTransform;
m_Transformation = mitk::AffineTransform3D::New();
if( !m_CombinedModality->GetIsTrackedUltrasoundActive() )
{
if( !m_Controls.m_ScaleTransform->isChecked() ) { m_Transformation->Compose(oldUSImageTransform); }
MITK_INFO << "Used old USImageTransform";
}
m_Transformation->Compose(calibTransform);
MITK_INFO << "New combined transform: " << m_Transformation;
mitk::SlicedGeometry3D::Pointer sliced3d = dynamic_cast (m_Node->GetData()->GetGeometry());
mitk::PlaneGeometry::Pointer plane = const_cast (sliced3d->GetPlaneGeometry(0));
plane->SetIndexToWorldTransform(m_Transformation);
// Save to US-Device
m_CombinedModality->SetCalibration(m_Transformation);
m_Controls.m_ToolBox->setItemEnabled(2, true);
// Save to NeedleProjectionFilter
m_NeedleProjectionFilter->SetTargetPlane(m_Transformation);
// Update Calibration FRE
m_CalibrationStatistics = mitk::PointSetDifferenceStatisticsCalculator::New();
mitk::PointSet::Pointer p1 = this->m_CalibPointsTool->Clone(); // We use clones to calculate statistics to avoid concurrency Problems
// Create point set with transformed image calibration points for
// calculating the difference of image calibration and tool
// calibration points in one geometry space
mitk::PointSet::Pointer p2 = mitk::PointSet::New();
int n = 0;
for (mitk::PointSet::PointsConstIterator it = m_CalibPointsImage->Begin();
it != m_CalibPointsImage->End(); ++it, ++n)
{
p2->InsertPoint(n, m_Transformation->TransformPoint(it->Value()));
}
m_CalibrationStatistics->SetPointSets(p1, p2);
//QString text = text.number(m_CalibrationStatistics->GetRMS());
QString text = QString::number(ComputeFRE(m_CalibPointsImage, m_CalibPointsTool, transform));
MITK_INFO << "Calibration FRE: " << text.toStdString().c_str();
m_Controls.m_EvalLblCalibrationFRE->setText(text);
m_Node->SetStringProperty("Calibration FRE", text.toStdString().c_str());
// Enable Button to save Calibration
m_Controls.m_CalibBtnSaveCalibration->setEnabled(true);
}
void QmitkUltrasoundCalibration::OnAddEvalTargetPoint()
{
mitk::Point3D world = this->GetRenderWindowPart()->GetSelectedPosition();
this->m_EvalPointsImage->InsertPoint(m_EvalPointsImage->GetSize(), world);
this->m_EvalPointsTool->InsertPoint(m_EvalPointsTool->GetSize(), this->m_FreezePoint);
QString text = text.number(this->m_EvalPointsTool->GetSize());
this->m_Controls.m_EvalLblNumTargetPoints->setText(text);
// Update FREs
// Update Evaluation FRE, but only if it contains more than one point (will crash otherwise)
if ((m_EvalPointsProjected->GetSize() > 1) && (m_EvalPointsTool->GetSize() > 1)) {
m_EvaluationStatistics = mitk::PointSetDifferenceStatisticsCalculator::New();
m_ProjectionStatistics = mitk::PointSetDifferenceStatisticsCalculator::New();
mitk::PointSet::Pointer p1 = this->m_EvalPointsTool->Clone(); // We use clones to calculate statistics to avoid concurrency Problems
mitk::PointSet::Pointer p2 = this->m_EvalPointsImage->Clone();
mitk::PointSet::Pointer p3 = this->m_EvalPointsProjected->Clone();
m_EvaluationStatistics->SetPointSets(p1, p2);
m_ProjectionStatistics->SetPointSets(p1, p3);
QString evalText = evalText.number(m_EvaluationStatistics->GetRMS());
QString projText = projText.number(m_ProjectionStatistics->GetRMS());
m_Controls.m_EvalLblEvaluationFRE->setText(evalText);
m_Controls.m_EvalLblProjectionFRE->setText(projText);
}
SwitchFreeze();
}
void QmitkUltrasoundCalibration::OnAddEvalProjectedPoint()
{
MITK_WARN << "Projection Evaluation may currently be inaccurate.";
// TODO: Verify correct Evaluation. Is the Point that is added really current?
mitk::Point3D projection = this->m_NeedleProjectionFilter->GetProjection()->GetPoint(1);
m_EvalPointsProjected->InsertPoint(m_EvalPointsProjected->GetSize(), projection);
QString text = text.number(this->m_EvalPointsProjected->GetSize());
this->m_Controls.m_EvalLblNumProjectionPoints->setText(text);
}
void QmitkUltrasoundCalibration::OnSaveEvaluation()
{
//Filename without suffix
QString filename = m_Controls.m_EvalFilePath->text() + "//" + m_Controls.m_EvalFilePrefix->text();
MITK_WARN << "CANNOT SAVE, ABORTING!";
/* not working any more TODO!
mitk::PointSetWriter::Pointer psWriter = mitk::PointSetWriter::New();
psWriter->SetInput(0, m_CalibPointsImage);
psWriter->SetInput(1, m_CalibPointsTool);
psWriter->SetInput(2, m_EvalPointsImage);
psWriter->SetInput(3, m_EvalPointsTool);
psWriter->SetInput(4, m_EvalPointsProjected);
psWriter->SetFileName(filename.toStdString() + ".xml");
psWriter->Write();
*/
// TODO: New writer for transformations must be implemented.
/*
mitk::TransformationFileWriter::Pointer tWriter = mitk::TransformationFileWriter::New();
tWriter->SetInput(0, m_CalibPointsImage);
tWriter->SetInput(1, m_CalibPointsTool);
tWriter->SetInput(2, m_EvalPointsImage);
tWriter->SetInput(3, m_EvalPointsTool);
tWriter->SetInput(4, m_EvalPointsProjected);
tWriter->SetOutputFilename(filename.toStdString() + ".txt");
tWriter->DoWrite(this->m_Transformation);
*/
}
void QmitkUltrasoundCalibration::OnSaveCalibration()
{
m_Controls.m_GotCalibrationLabel->setText("");
QString filename = QFileDialog::getSaveFileName(QApplication::activeWindow(),
"Save Calibration",
"",
"Calibration files *.cal");
QFile file(filename);
if (!file.open(QIODevice::WriteOnly | QIODevice::Text | QIODevice::Truncate))
{
MITK_WARN << "Cannot open file '" << filename.toStdString() << "' for writing.";
return;
}
std::string calibrationSerialization = m_CombinedModality->SerializeCalibration();
QTextStream outStream(&file);
outStream << QString::fromStdString(calibrationSerialization);
//save additional information
if (m_Controls.m_saveAdditionalCalibrationLog->isChecked())
{
mitk::SceneIO::Pointer mySceneIO = mitk::SceneIO::New();
QString filenameScene = filename + "_mitkScene.mitk";
mitk::NodePredicateNot::Pointer isNotHelperObject =
mitk::NodePredicateNot::New(mitk::NodePredicateProperty::New("helper object", mitk::BoolProperty::New(true)));
mitk::DataStorage::SetOfObjects::ConstPointer nodesToBeSaved = this->GetDataStorage()->GetSubset(isNotHelperObject);
mySceneIO->SaveScene(nodesToBeSaved, this->GetDataStorage(), filenameScene.toStdString().c_str());
}
}
void QmitkUltrasoundCalibration::OnReset()
{
this->ClearTemporaryMembers();
if (m_Transformation.IsNull()) { m_Transformation = mitk::AffineTransform3D::New(); }
m_Transformation->SetIdentity();
if (m_Node.IsNotNull() && (m_Node->GetData() != nullptr) && (m_Node->GetData()->GetGeometry() != nullptr))
{
mitk::SlicedGeometry3D::Pointer sliced3d = dynamic_cast (m_Node->GetData()->GetGeometry());
mitk::PlaneGeometry::Pointer plane = const_cast (sliced3d->GetPlaneGeometry(0));
plane->SetIndexToWorldTransform(m_Transformation);
}
QString text1 = text1.number(this->m_EvalPointsTool->GetSize());
this->m_Controls.m_EvalLblNumTargetPoints->setText(text1);
QString text2 = text2.number(this->m_EvalPointsProjected->GetSize());
this->m_Controls.m_EvalLblNumProjectionPoints->setText(text2);
}
void QmitkUltrasoundCalibration::Update()
{
//QList nodes = this->GetDataManagerSelection();
// if (nodes.empty()) return;
// Update Tracking Data
std::vector* datas = new std::vector();
datas->push_back(m_Tracker->GetOutput());
m_Controls.m_CalibTrackingStatus->SetNavigationDatas(datas);
m_Controls.m_CalibTrackingStatus->Refresh();
m_Controls.m_EvalTrackingStatus->SetNavigationDatas(datas);
m_Controls.m_EvalTrackingStatus->Refresh();
/*
if (m_Image.IsNotNull() && m_Image->IsInitialized())
{
m_Node->SetData(m_Image);
}
else
{
m_Image = m_CombinedModality->GetOutput();
m_Node->SetData(m_Image);
}*/
m_CombinedModality->Modified();
m_CombinedModality->Update();
// Update US Image
mitk::Image::Pointer image = m_CombinedModality->GetOutput();
// make sure that always the current image is set to the data node
if (image.IsNotNull() && m_Node->GetData() != image.GetPointer() && image->IsInitialized())
{
m_Node->SetData(image);
}
// Update Needle Projection
m_NeedleProjectionFilter->Update();
//only update 2d window because it is faster
this->RequestRenderWindowUpdate(mitk::RenderingManager::REQUEST_UPDATE_2DWINDOWS);
}
void QmitkUltrasoundCalibration::SwitchFreeze()
{
m_Controls.m_CalibBtnAddPoint->setEnabled(false); // generally deactivate
// We use the activity state of the timer to determine whether we are currently viewing images
if (!m_Timer->isActive()) // Activate Imaging
{
// if (m_Node) m_Node->ReleaseData();
if (m_CombinedModality.IsNull()){
m_Timer->stop();
return;
}
m_CombinedModality->Update();
m_Image = m_CombinedModality->GetOutput();
if (m_Image.IsNotNull() && m_Image->IsInitialized())
{
m_Node->SetData(m_Image);
}
std::vector datas;
datas.push_back(m_Tracker->GetOutput());
m_Controls.m_CalibTrackingStatus->SetNavigationDatas(&datas);
m_Controls.m_CalibTrackingStatus->ShowStatusLabels();
m_Controls.m_CalibTrackingStatus->Refresh();
m_Controls.m_EvalTrackingStatus->SetNavigationDatas(&datas);
m_Controls.m_EvalTrackingStatus->ShowStatusLabels();
m_Controls.m_EvalTrackingStatus->Refresh();
int interval = 40;
m_Timer->setInterval(interval);
m_Timer->start();
m_CombinedModality->SetIsFreezed(false);
}
else if (this->m_Tracker->GetOutput(0)->IsDataValid())
{
//deactivate Imaging
m_Timer->stop();
// Remember last tracking coordinates
m_FreezePoint = this->m_Tracker->GetOutput(0)->GetPosition();
m_Controls.m_CalibBtnAddPoint->setEnabled(true); // activate only, if valid point is set
m_CombinedModality->SetIsFreezed(true);
}
}
void QmitkUltrasoundCalibration::ShowNeedlePath()
{
// Init Filter
this->m_NeedleProjectionFilter->SelectInput(0);
// Create Node for Pointset
mitk::DataNode::Pointer node = this->GetDataStorage()->GetNamedNode("Needle Path");
if (node.IsNull())
{
node = mitk::DataNode::New();
node->SetName("Needle Path");
node->SetData(m_NeedleProjectionFilter->GetProjection());
node->SetBoolProperty("show contour", true);
this->GetDataStorage()->Add(node);
}
}
void QmitkUltrasoundCalibration::ClearTemporaryMembers()
{
m_CalibPointsTool->Clear();
m_CalibPointsImage->Clear();
m_CalibPointsCount = 0;
m_EvalPointsImage->Clear();
m_EvalPointsTool->Clear();
m_EvalPointsProjected->Clear();
this->m_Controls.m_CalibPointList->clear();
m_SpacingPoints->Clear();
m_Controls.m_SpacingPointsList->clear();
m_SpacingPointsCount = 0;
}
vtkSmartPointer QmitkUltrasoundCalibration::ConvertPointSetToVtkPolyData(mitk::PointSet::Pointer PointSet)
{
vtkSmartPointer returnValue = vtkSmartPointer::New();
vtkSmartPointer points = vtkSmartPointer::New();
for (int i = 0; i < PointSet->GetSize(); i++)
{
double point[3] = { PointSet->GetPoint(i)[0], PointSet->GetPoint(i)[1], PointSet->GetPoint(i)[2] };
points->InsertNextPoint(point);
}
vtkSmartPointer temp = vtkSmartPointer::New();
temp->SetPoints(points);
vtkSmartPointer vertexFilter = vtkSmartPointer::New();
vertexFilter->SetInputData(temp);
vertexFilter->Update();
returnValue->ShallowCopy(vertexFilter->GetOutput());
return returnValue;
}
double QmitkUltrasoundCalibration::ComputeFRE(mitk::PointSet::Pointer imageFiducials, mitk::PointSet::Pointer realWorldFiducials, vtkSmartPointer transform)
{
if (imageFiducials->GetSize() != realWorldFiducials->GetSize()) return -1;
double FRE = 0;
for (int i = 0; i < imageFiducials->GetSize(); ++i)
{
itk::Point current_image_fiducial_point = imageFiducials->GetPoint(i);
if (transform != nullptr)
{
current_image_fiducial_point = transform->TransformPoint(imageFiducials->GetPoint(i)[0], imageFiducials->GetPoint(i)[1], imageFiducials->GetPoint(i)[2]);
}
double cur_error_squared = current_image_fiducial_point.SquaredEuclideanDistanceTo(realWorldFiducials->GetPoint(i));
FRE += cur_error_squared;
}
FRE = sqrt(FRE / (double)imageFiducials->GetSize());
return FRE;
}
void QmitkUltrasoundCalibration::ApplyTransformToPointSet(mitk::PointSet::Pointer pointSet, vtkSmartPointer transform)
{
for (int i = 0; i < pointSet->GetSize(); ++i)
{
itk::Point current_point_transformed = itk::Point();
current_point_transformed = transform->TransformPoint(pointSet->GetPoint(i)[0], pointSet->GetPoint(i)[1], pointSet->GetPoint(i)[2]);
pointSet->SetPoint(i, current_point_transformed);
}
}
void QmitkUltrasoundCalibration::OnFreezeClicked()
{
if (m_CombinedModality->GetIsFreezed())
{
//device was already frozen so we need to delete all Spacing points because they need to be collected all at once
// no need to check if all four points are already collected, because if thats the case you can no longer click the Freeze Button
m_SpacingPoints->Clear();
m_Controls.m_SpacingPointsList->clear();
m_SpacingPointsCount = 0;
m_Controls.m_SpacingAddPoint->setEnabled(false);
}
else
{
m_Controls.m_SpacingAddPoint->setEnabled(true);
}
SwitchFreeze();
}
void QmitkUltrasoundCalibration::OnAddSpacingPoint()
{
mitk::Point3D point = this->GetRenderWindowPart()->GetSelectedPosition();
this->m_SpacingPoints->InsertPoint(m_SpacingPointsCount, point);
QString text = text.number(m_SpacingPointsCount + 1);
text = "Point " + text;
this->m_Controls.m_SpacingPointsList->addItem(text);
m_SpacingPointsCount++;
if (m_SpacingPointsCount == 4) //now we have all 4 points needed
{
m_Controls.m_SpacingAddPoint->setEnabled(false);
m_Controls.m_CalculateSpacing->setEnabled(true);
m_Controls.m_SpacingBtnFreeze->setEnabled(false);
}
}
void QmitkUltrasoundCalibration::OnCalculateSpacing()
{
mitk::Point3D horizontalOne = m_SpacingPoints->GetPoint(0);
mitk::Point3D horizontalTwo = m_SpacingPoints->GetPoint(1);
mitk::Point3D verticalOne = m_SpacingPoints->GetPoint(2);
mitk::Point3D verticalTwo = m_SpacingPoints->GetPoint(3);
//Get the distances between the points in the image
double xDistance = horizontalOne.EuclideanDistanceTo(horizontalTwo);
double yDistance = verticalOne.EuclideanDistanceTo(verticalTwo);
//Calculate the spacing of the image and fill a vector with it
double xSpacing = 30 / xDistance;
double ySpacing = 20 / yDistance;
m_CombinedModality->GetUltrasoundDevice()->SetSpacing(xSpacing, ySpacing);
//Now that the spacing is set clear all stuff and return to Calibration
m_SpacingPoints->Clear();
m_Controls.m_SpacingPointsList->clear();
m_SpacingPointsCount = 0;
m_CombinedModality->SetIsFreezed(false);
}
void QmitkUltrasoundCalibration::OnUSDepthChanged(const std::string& key, const std::string&)
{
//whenever depth of USImage is changed the spacing should no longer be overwritten
if (key == mitk::USDevice::GetPropertyKeys().US_PROPKEY_BMODE_DEPTH)
{
}
}