diff --git a/Core/Code/DataManagement/mitkAbstractTransformGeometry.cpp b/Core/Code/DataManagement/mitkAbstractTransformGeometry.cpp
index db80da8ae6..5a79c7edc9 100644
--- a/Core/Code/DataManagement/mitkAbstractTransformGeometry.cpp
+++ b/Core/Code/DataManagement/mitkAbstractTransformGeometry.cpp
@@ -1,318 +1,318 @@
/*===================================================================
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 "mitkAbstractTransformGeometry.h"
#include <vtkAbstractTransform.h>
-mitk::AbstractTransformGeometry::AbstractTransformGeometry() : m_Plane(NULL), m_FrameGeometry(NULL)
+mitk::AbstractTransformGeometry::AbstractTransformGeometry() : Superclass(), m_Plane(NULL), m_FrameGeometry(NULL)
{
Initialize();
}
mitk::AbstractTransformGeometry::AbstractTransformGeometry(const AbstractTransformGeometry& other) : Superclass(other), m_ParametricBoundingBox(other.m_ParametricBoundingBox)
{
if(other.m_ParametricBoundingBox.IsNotNull())
{
m_ParametricBoundingBox = other.m_ParametricBoundingBox->DeepCopy();
this->SetParametricBounds(m_ParametricBoundingBox->GetBounds());
}
this->SetPlane(other.m_Plane);
this->SetFrameGeometry(other.m_FrameGeometry);
}
mitk::AbstractTransformGeometry::~AbstractTransformGeometry()
{
}
void mitk::AbstractTransformGeometry::PostInitialize()
{
m_ItkVtkAbstractTransform = itk::VtkAbstractTransform<ScalarType>::New();
}
vtkAbstractTransform* mitk::AbstractTransformGeometry::GetVtkAbstractTransform() const
{
return m_ItkVtkAbstractTransform->GetVtkAbstractTransform();
}
mitk::ScalarType mitk::AbstractTransformGeometry::GetParametricExtentInMM(int direction) const
{
if(m_Plane.IsNull())
{
itkExceptionMacro(<<"m_Plane is NULL.");
}
return m_Plane->GetExtentInMM(direction);
}
const itk::Transform<mitk::ScalarType, 3, 3>* mitk::AbstractTransformGeometry::GetParametricTransform() const
{
return m_ItkVtkAbstractTransform;
}
bool mitk::AbstractTransformGeometry::Project(const mitk::Point3D &pt3d_mm, mitk::Point3D &projectedPt3d_mm) const
{
assert(this->IsBoundingBoxNull()==false);
mitk::Point2D pt2d_mm;
bool isInside;
isInside = Map(pt3d_mm, pt2d_mm);
Map(pt2d_mm, projectedPt3d_mm);
return isInside;
//Point3D pt3d_units;
//pt3d_units = m_ItkVtkAbstractTransform->BackTransform(pt3d_mm);
//pt3d_units[2] = 0;
//projectedPt3d_mm = m_ItkVtkAbstractTransform->TransformPoint(pt3d_units);
//return const_cast<BoundingBox*>(m_BoundingBox.GetPointer())->IsInside(pt3d_units);
}
bool mitk::AbstractTransformGeometry::Map(const mitk::Point3D &pt3d_mm, mitk::Point2D &pt2d_mm) const
{
assert((m_ItkVtkAbstractTransform.IsNotNull()) && (m_Plane.IsNotNull()));
Point3D pt3d_units;
pt3d_units = m_ItkVtkAbstractTransform->BackTransform(pt3d_mm);
return m_Plane->Map(pt3d_units, pt2d_mm);
}
void mitk::AbstractTransformGeometry::Map(const mitk::Point2D &pt2d_mm, mitk::Point3D &pt3d_mm) const
{
assert((m_ItkVtkAbstractTransform.IsNotNull()) && (m_Plane.IsNotNull()));
m_Plane->Map(pt2d_mm, pt3d_mm);
pt3d_mm = m_ItkVtkAbstractTransform->TransformPoint(pt3d_mm);
}
bool mitk::AbstractTransformGeometry::Project(const mitk::Point3D & atPt3d_mm, const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm) const
{
itkExceptionMacro("not implemented yet - replace GetIndexToWorldTransform by m_ItkVtkAbstractTransform->GetInverseVtkAbstractTransform()");
assert(this->IsBoundingBoxNull()==false);
Vector3D vec3d_units;
vec3d_units = GetIndexToWorldTransform()->GetInverseMatrix() * vec3d_mm;
vec3d_units[2] = 0;
projectedVec3d_mm = GetIndexToWorldTransform()->TransformVector(vec3d_units);
Point3D pt3d_units;
mitk::ScalarType temp[3];
unsigned int i, j;
for (j = 0; j < 3; ++j)
temp[j] = atPt3d_mm[j] - GetIndexToWorldTransform()->GetOffset()[j];
for (i = 0; i < 3; ++i)
{
pt3d_units[i] = 0.0;
for (j = 0; j < 3; ++j)
pt3d_units[i] += GetIndexToWorldTransform()->GetInverseMatrix()[i][j] * temp[j];
}
return const_cast<BoundingBox*>(this->GetBoundingBox())->IsInside(pt3d_units);
}
bool mitk::AbstractTransformGeometry::Project(const mitk::Vector3D &/*vec3d_mm*/, mitk::Vector3D &/*projectedVec3d_mm*/) const
{
MITK_WARN << "Need additional point! No standard value defined. Please use Project(const mitk::Point3D & atPt3d_mm, const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm). Unfortunatley this one is not implemented at the moment. Sorry :(";
itkExceptionMacro("not implemented yet - replace GetIndexToWorldTransform by m_ItkVtkAbstractTransform->GetInverseVtkAbstractTransform()");
return false;
}
bool mitk::AbstractTransformGeometry::Map(const mitk::Point3D & atPt3d_mm, const mitk::Vector3D &vec3d_mm, mitk::Vector2D &vec2d_mm) const
{
assert((m_ItkVtkAbstractTransform.IsNotNull()) && (m_Plane.IsNotNull()));
ScalarType vtkpt[3], vtkvec[3];
itk2vtk(atPt3d_mm, vtkpt);
itk2vtk(vec3d_mm, vtkvec);
m_ItkVtkAbstractTransform->GetInverseVtkAbstractTransform()->TransformVectorAtPoint(vtkpt, vtkvec, vtkvec);
mitk::Vector3D vec3d_units;
vtk2itk(vtkvec, vec3d_units);
return m_Plane->Map(atPt3d_mm, vec3d_units, vec2d_mm);
}
void mitk::AbstractTransformGeometry::Map(const mitk::Point2D & atPt2d_mm, const mitk::Vector2D &vec2d_mm, mitk::Vector3D &vec3d_mm) const
{
m_Plane->Map(atPt2d_mm, vec2d_mm, vec3d_mm);
Point3D atPt3d_mm;
Map(atPt2d_mm, atPt3d_mm);
float vtkpt[3], vtkvec[3];
itk2vtk(atPt3d_mm, vtkpt);
itk2vtk(vec3d_mm, vtkvec);
m_ItkVtkAbstractTransform->GetVtkAbstractTransform()->TransformVectorAtPoint(vtkpt, vtkvec, vtkvec);
vtk2itk(vtkvec, vec3d_mm);
}
void mitk::AbstractTransformGeometry::IndexToWorld(const mitk::Point2D &pt_units, mitk::Point2D &pt_mm) const
{
m_Plane->IndexToWorld(pt_units, pt_mm);
}
void mitk::AbstractTransformGeometry::WorldToIndex(const mitk::Point2D &pt_mm, mitk::Point2D &pt_units) const
{
m_Plane->WorldToIndex(pt_mm, pt_units);
}
void mitk::AbstractTransformGeometry::IndexToWorld(const mitk::Point2D & /*atPt2d_units*/, const mitk::Vector2D &vec_units, mitk::Vector2D &vec_mm) const
{
MITK_WARN<<"Warning! Call of the deprecated function AbstractTransformGeometry::IndexToWorld(point, vec, vec). Use AbstractTransformGeometry::IndexToWorld(vec, vec) instead!";
this->IndexToWorld(vec_units, vec_mm);
}
void mitk::AbstractTransformGeometry::IndexToWorld(const mitk::Vector2D &vec_units, mitk::Vector2D &vec_mm) const
{
m_Plane->IndexToWorld(vec_units, vec_mm);
}
void mitk::AbstractTransformGeometry::WorldToIndex(const mitk::Point2D & /*atPt2d_mm*/, const mitk::Vector2D &vec_mm, mitk::Vector2D &vec_units) const
{
MITK_WARN<<"Warning! Call of the deprecated function AbstractTransformGeometry::WorldToIndex(point, vec, vec). Use AbstractTransformGeometry::WorldToIndex(vec, vec) instead!";
this->WorldToIndex(vec_mm, vec_units);
}
void mitk::AbstractTransformGeometry::WorldToIndex(const mitk::Vector2D &vec_mm, mitk::Vector2D &vec_units) const
{
m_Plane->WorldToIndex(vec_mm, vec_units);
}
bool mitk::AbstractTransformGeometry::IsAbove(const mitk::Point3D& pt3d_mm, bool considerBoundingBox) const
{
assert((m_ItkVtkAbstractTransform.IsNotNull()) && (m_Plane.IsNotNull()));
Point3D pt3d_ParametricWorld;
pt3d_ParametricWorld = m_ItkVtkAbstractTransform->BackTransform(pt3d_mm);
Point3D pt3d_ParametricUnits;
((BaseGeometry*)m_Plane)->WorldToIndex(pt3d_ParametricWorld, pt3d_ParametricUnits);
return (pt3d_ParametricUnits[2] > m_ParametricBoundingBox->GetBounds()[4]);
}
void mitk::AbstractTransformGeometry::SetVtkAbstractTransform(vtkAbstractTransform* aVtkAbstractTransform)
{
m_ItkVtkAbstractTransform->SetVtkAbstractTransform(aVtkAbstractTransform);
}
void mitk::AbstractTransformGeometry::SetPlane(const mitk::PlaneGeometry* aPlane)
{
if(aPlane!=NULL)
{
m_Plane = static_cast<mitk::PlaneGeometry*>(aPlane->Clone().GetPointer());
BoundingBox::BoundsArrayType b=m_Plane->GetBoundingBox()->GetBounds();
SetParametricBounds(b);
CalculateFrameGeometry();
}
else
{
if(m_Plane.IsNull())
return;
m_Plane=NULL;
}
Modified();
}
void mitk::AbstractTransformGeometry::CalculateFrameGeometry()
{
if((m_Plane.IsNull()) || (m_FrameGeometry.IsNotNull()))
return;
//@warning affine-transforms and bounding-box should be set by specific sub-classes!
SetBounds(m_Plane->GetBoundingBox()->GetBounds());
}
void mitk::AbstractTransformGeometry::SetFrameGeometry(const mitk::BaseGeometry* frameGeometry)
{
if((frameGeometry != NULL) && (frameGeometry->IsValid()))
{
m_FrameGeometry = static_cast<mitk::BaseGeometry*>(frameGeometry->Clone().GetPointer());
SetIndexToWorldTransform(m_FrameGeometry->GetIndexToWorldTransform());
SetBounds(m_FrameGeometry->GetBounds());
}
else
{
m_FrameGeometry = NULL;
}
}
unsigned long mitk::AbstractTransformGeometry::GetMTime() const
{
if(Superclass::GetMTime()<m_ItkVtkAbstractTransform->GetMTime())
return m_ItkVtkAbstractTransform->GetMTime();
return Superclass::GetMTime();
}
void mitk::AbstractTransformGeometry::SetOversampling(mitk::ScalarType oversampling)
{
if(m_Plane.IsNull())
{
itkExceptionMacro(<< "m_Plane is not set.");
}
mitk::BoundingBox::BoundsArrayType bounds = m_Plane->GetBounds();
bounds[1]*=oversampling; bounds[3]*=oversampling; bounds[5]*=oversampling;
SetParametricBounds(bounds);
}
itk::LightObject::Pointer mitk::AbstractTransformGeometry::InternalClone() const
{
Self::Pointer newGeometry = new AbstractTransformGeometry(*this);
newGeometry->UnRegister();
return newGeometry.GetPointer();
}
void mitk::AbstractTransformGeometry::SetParametricBounds(const BoundingBox::BoundsArrayType& bounds)
{
m_ParametricBoundingBox = BoundingBoxType::New();
BoundingBoxType::PointsContainer::Pointer pointscontainer =
BoundingBoxType::PointsContainer::New();
BoundingBoxType::PointType p;
BoundingBoxType::PointIdentifier pointid;
for(pointid=0; pointid<2;++pointid)
{
unsigned int i;
for(i=0; i<GetNDimensions(); ++i)
{
p[i] = bounds[2*i+pointid];
}
pointscontainer->InsertElement(pointid, p);
}
m_ParametricBoundingBox->SetPoints(pointscontainer);
m_ParametricBoundingBox->ComputeBoundingBox();
this->Modified();
}
const mitk::BoundingBox::BoundsArrayType& mitk::AbstractTransformGeometry::GetParametricBounds() const
{
assert(m_ParametricBoundingBox.IsNotNull());
return m_ParametricBoundingBox->GetBounds();
}
mitk::ScalarType mitk::AbstractTransformGeometry::GetParametricExtent(int direction) const
{
if (direction < 0 || direction>=3)
mitkThrow() << "Invalid direction. Must be between either 0, 1 or 2. ";
assert(m_ParametricBoundingBox.IsNotNull());
BoundingBoxType::BoundsArrayType bounds = m_ParametricBoundingBox->GetBounds();
return bounds[direction*2+1]-bounds[direction*2];
}
diff --git a/Core/Code/DataManagement/mitkBaseGeometry.cpp b/Core/Code/DataManagement/mitkBaseGeometry.cpp
index 4f10cb7607..fb51ca75e3 100644
--- a/Core/Code/DataManagement/mitkBaseGeometry.cpp
+++ b/Core/Code/DataManagement/mitkBaseGeometry.cpp
@@ -1,1095 +1,1086 @@
/*===================================================================
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 <sstream>
#include <iomanip>
#include "mitkBaseGeometry.h"
#include "mitkVector.h"
#include "mitkMatrixConvert.h"
#include <vtkMatrixToLinearTransform.h>
#include <vtkMatrix4x4.h>
#include "mitkRotationOperation.h"
#include "mitkRestorePlanePositionOperation.h"
#include "mitkApplyTransformMatrixOperation.h"
#include "mitkPointOperation.h"
#include "mitkInteractionConst.h"
#include "mitkModifiedLock.h"
mitk::BaseGeometry::BaseGeometry(): Superclass(), mitk::OperationActor(),
m_FrameOfReferenceID(0), m_IndexToWorldTransformLastModified(0), m_ImageGeometry(false), m_ModifiedLockFlag(false), m_ModifiedCalledFlag(false)
{
m_VtkMatrix = vtkMatrix4x4::New();
m_VtkIndexToWorldTransform = vtkMatrixToLinearTransform::New();
m_VtkIndexToWorldTransform->SetInput(m_VtkMatrix);
Initialize();
}
mitk::BaseGeometry::BaseGeometry(const BaseGeometry& other): Superclass(), mitk::OperationActor(), //m_TimeBounds(other.m_TimeBounds),
m_FrameOfReferenceID(other.m_FrameOfReferenceID), m_IndexToWorldTransformLastModified(other.m_IndexToWorldTransformLastModified), m_Origin(other.m_Origin),
m_ImageGeometry(other.m_ImageGeometry), m_ModifiedLockFlag(false), m_ModifiedCalledFlag(false)
{
- // DEPRECATED(m_RotationQuaternion = other.m_RotationQuaternion);
- // AffineGeometryFrame
- SetBounds(other.GetBounds());
m_VtkMatrix = vtkMatrix4x4::New();
- m_VtkMatrix->DeepCopy(other.m_VtkMatrix);
-
m_VtkIndexToWorldTransform = vtkMatrixToLinearTransform::New();
- m_VtkIndexToWorldTransform->DeepCopy(other.m_VtkIndexToWorldTransform);
m_VtkIndexToWorldTransform->SetInput(m_VtkMatrix);
other.InitializeGeometry(this);
}
mitk::BaseGeometry::~BaseGeometry()
{
m_VtkMatrix->Delete();
m_VtkIndexToWorldTransform->Delete();
}
const mitk::Point3D& mitk::BaseGeometry::GetOrigin() const
{
return m_Origin;
}
void mitk::BaseGeometry::SetOrigin(const Point3D & origin)
{
mitk::ModifiedLock lock(this);
if(origin!=GetOrigin())
{
m_Origin = origin;
m_IndexToWorldTransform->SetOffset(m_Origin.GetVectorFromOrigin());
Modified();
TransferItkToVtkTransform();
}
}
void mitk::BaseGeometry::TransferItkToVtkTransform()
{
mitk::ModifiedLock lock(this);
// copy m_IndexToWorldTransform into m_VtkIndexToWorldTransform
TransferItkTransformToVtkMatrix(m_IndexToWorldTransform.GetPointer(), m_VtkMatrix);
m_VtkIndexToWorldTransform->Modified();
}
static void CopySpacingFromTransform(mitk::AffineTransform3D* transform, mitk::Vector3D& spacing)
{
mitk::AffineTransform3D::MatrixType::InternalMatrixType vnlmatrix;
vnlmatrix = transform->GetMatrix().GetVnlMatrix();
spacing[0]=vnlmatrix.get_column(0).magnitude();
spacing[1]=vnlmatrix.get_column(1).magnitude();
spacing[2]=vnlmatrix.get_column(2).magnitude();
}
void mitk::BaseGeometry::Initialize()
{
float b[6] = {0,1,0,1,0,1};
SetFloatBounds(b);
if(m_IndexToWorldTransform.IsNull())
m_IndexToWorldTransform = TransformType::New();
else
m_IndexToWorldTransform->SetIdentity();
CopySpacingFromTransform(m_IndexToWorldTransform, m_Spacing);
vtk2itk(m_IndexToWorldTransform->GetOffset(), m_Origin);
m_VtkMatrix->Identity();
//m_TimeBounds[0]=ScalarTypeNumericTraits::NonpositiveMin(); m_TimeBounds[1]=ScalarTypeNumericTraits::max();
m_FrameOfReferenceID = 0;
m_ImageGeometry = false;
this->PostInitialize();
}
void mitk::BaseGeometry::PostInitializeGeometry(BaseGeometry * newGeometry) const
{
newGeometry->m_ImageGeometry = m_ImageGeometry;
}
void mitk::BaseGeometry::SetFloatBounds(const float bounds[6])
{
mitk::BoundingBox::BoundsArrayType b;
const float *input = bounds;
int i=0;
for(mitk::BoundingBox::BoundsArrayType::Iterator it = b.Begin(); i < 6 ;++i) *it++ = (mitk::ScalarType)*input++;
SetBounds(b);
}
void mitk::BaseGeometry::SetFloatBounds(const double bounds[6])
{
mitk::BoundingBox::BoundsArrayType b;
const double *input = bounds;
int i=0;
for(mitk::BoundingBox::BoundsArrayType::Iterator it = b.Begin(); i < 6 ;++i) *it++ = (mitk::ScalarType)*input++;
SetBounds(b);
}
/** Initialize the geometry */
void
mitk::BaseGeometry::InitializeGeometry(BaseGeometry* newGeometry) const
{
newGeometry->SetBounds(m_BoundingBox->GetBounds());
// we have to create a new transform!!
//newGeometry->SetTimeBounds(m_TimeBounds);
newGeometry->SetFrameOfReferenceID(GetFrameOfReferenceID());
if(m_IndexToWorldTransform)
{
- TransformType::Pointer indexToWorldTransform = TransformType::New();
- indexToWorldTransform->SetCenter( m_IndexToWorldTransform->GetCenter() );
- indexToWorldTransform->SetMatrix( m_IndexToWorldTransform->GetMatrix() );
- indexToWorldTransform->SetOffset( m_IndexToWorldTransform->GetOffset() );
+ TransformType::Pointer indexToWorldTransform = m_IndexToWorldTransform->Clone();
newGeometry->SetIndexToWorldTransform(indexToWorldTransform);
}
this->PostInitializeGeometry(newGeometry);
}
void mitk::BaseGeometry::PostInitialize()
{
}
/** Set the bounds */
void mitk::BaseGeometry::SetBounds(const BoundsArrayType& bounds)
{
mitk::ModifiedLock lock(this);
PreSetBounds(bounds);
m_BoundingBox = BoundingBoxType::New();
BoundingBoxType::PointsContainer::Pointer pointscontainer =
BoundingBoxType::PointsContainer::New();
BoundingBoxType::PointType p;
BoundingBoxType::PointIdentifier pointid;
for(pointid=0; pointid<2;++pointid)
{
unsigned int i;
for(i=0; i<m_NDimensions; ++i)
{
p[i] = bounds[2*i+pointid];
}
pointscontainer->InsertElement(pointid, p);
}
m_BoundingBox->SetPoints(pointscontainer);
m_BoundingBox->ComputeBoundingBox();
this->Modified();
}
void mitk::BaseGeometry::PreSetBounds(const BoundsArrayType& /*bounds*/){};
void mitk::BaseGeometry::SetIndexToWorldTransform(mitk::AffineTransform3D* transform)
{
mitk::ModifiedLock lock(this);
PreSetIndexToWorldTransform(transform);
if(m_IndexToWorldTransform.GetPointer() != transform)
{
m_IndexToWorldTransform = transform;
CopySpacingFromTransform(m_IndexToWorldTransform, m_Spacing);
vtk2itk(m_IndexToWorldTransform->GetOffset(), m_Origin);
TransferItkToVtkTransform();
Modified();
}
PostSetIndexToWorldTransform(transform);
}
void mitk::BaseGeometry::PreSetIndexToWorldTransform(mitk::AffineTransform3D* /*transform*/)
{}
void mitk::BaseGeometry::PostSetIndexToWorldTransform(mitk::AffineTransform3D* /*transform*/)
{}
const mitk::BaseGeometry::BoundsArrayType mitk::BaseGeometry::GetBounds() const
{
assert(m_BoundingBox.IsNotNull());
return m_BoundingBox->GetBounds();
}
bool mitk::BaseGeometry::IsValid() const
{
return true;
}
void mitk::BaseGeometry::SetSpacing(const mitk::Vector3D& aSpacing, bool enforceSetSpacing )
{
PreSetSpacing(aSpacing);
_SetSpacing(aSpacing, enforceSetSpacing);
}
void mitk::BaseGeometry::PreSetSpacing(const mitk::Vector3D& /*aSpacing*/)
{}
void mitk::BaseGeometry::_SetSpacing(const mitk::Vector3D& aSpacing, bool enforceSetSpacing){
if(mitk::Equal(m_Spacing, aSpacing) == false || enforceSetSpacing)
{
assert(aSpacing[0]>0 && aSpacing[1]>0 && aSpacing[2]>0);
m_Spacing = aSpacing;
AffineTransform3D::MatrixType::InternalMatrixType vnlmatrix;
vnlmatrix = m_IndexToWorldTransform->GetMatrix().GetVnlMatrix();
mitk::VnlVector col;
col = vnlmatrix.get_column(0); col.normalize(); col*=aSpacing[0]; vnlmatrix.set_column(0, col);
col = vnlmatrix.get_column(1); col.normalize(); col*=aSpacing[1]; vnlmatrix.set_column(1, col);
col = vnlmatrix.get_column(2); col.normalize(); col*=aSpacing[2]; vnlmatrix.set_column(2, col);
Matrix3D matrix;
matrix = vnlmatrix;
AffineTransform3D::Pointer transform = AffineTransform3D::New();
transform->SetMatrix(matrix);
transform->SetOffset(m_IndexToWorldTransform->GetOffset());
SetIndexToWorldTransform(transform.GetPointer());
}
}
mitk::Vector3D mitk::BaseGeometry::GetAxisVector(unsigned int direction) const
{
Vector3D frontToBack;
frontToBack.SetVnlVector(m_IndexToWorldTransform->GetMatrix().GetVnlMatrix().get_column(direction));
frontToBack *= GetExtent(direction);
return frontToBack;
}
mitk::ScalarType mitk::BaseGeometry::GetExtent(unsigned int direction) const
{
assert(m_BoundingBox.IsNotNull());
if (direction>=m_NDimensions)
mitkThrow() << "Direction is too big. This geometry is for 3D Data";
BoundsArrayType bounds = m_BoundingBox->GetBounds();
return bounds[direction*2+1]-bounds[direction*2];
}
bool mitk::BaseGeometry::Is2DConvertable()
{
bool isConvertableWithoutLoss = true;
do
{
if (this->GetSpacing()[2] != 1)
{
isConvertableWithoutLoss = false;
break;
}
if (this->GetOrigin()[2] != 0)
{
isConvertableWithoutLoss = false;
break;
}
mitk::Vector3D col0, col1, col2;
col0.SetVnlVector(this->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix().get_column(0));
col1.SetVnlVector(this->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix().get_column(1));
col2.SetVnlVector(this->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix().get_column(2));
if ((col0[2] != 0) || (col1[2] != 0) || (col2[0] != 0) || (col2[1] != 0) || (col2[2] != 1))
{
isConvertableWithoutLoss = false;
break;
}
} while (0);
return isConvertableWithoutLoss;
}
mitk::Point3D mitk::BaseGeometry::GetCenter() const
{
assert(m_BoundingBox.IsNotNull());
return m_IndexToWorldTransform->TransformPoint(m_BoundingBox->GetCenter());
}
double mitk::BaseGeometry::GetDiagonalLength2() const
{
Vector3D diagonalvector = GetCornerPoint()-GetCornerPoint(false, false, false);
return diagonalvector.GetSquaredNorm();
}
//##Documentation
//## @brief Get the length of the diagonal of the bounding-box in mm
//##
double mitk::BaseGeometry::GetDiagonalLength() const
{
return sqrt(GetDiagonalLength2());
}
mitk::Point3D mitk::BaseGeometry::GetCornerPoint(int id) const
{
assert(id >= 0);
assert(this->IsBoundingBoxNull()==false);
BoundingBox::BoundsArrayType bounds = this->GetBoundingBox()->GetBounds();
Point3D cornerpoint;
switch(id)
{
case 0: FillVector3D(cornerpoint, bounds[0],bounds[2],bounds[4]); break;
case 1: FillVector3D(cornerpoint, bounds[0],bounds[2],bounds[5]); break;
case 2: FillVector3D(cornerpoint, bounds[0],bounds[3],bounds[4]); break;
case 3: FillVector3D(cornerpoint, bounds[0],bounds[3],bounds[5]); break;
case 4: FillVector3D(cornerpoint, bounds[1],bounds[2],bounds[4]); break;
case 5: FillVector3D(cornerpoint, bounds[1],bounds[2],bounds[5]); break;
case 6: FillVector3D(cornerpoint, bounds[1],bounds[3],bounds[4]); break;
case 7: FillVector3D(cornerpoint, bounds[1],bounds[3],bounds[5]); break;
default:
{
itkExceptionMacro(<<"A cube only has 8 corners. These are labeled 0-7.");
}
}
if(m_ImageGeometry)
{
// Here i have to adjust the 0.5 offset manually, because the cornerpoint is the corner of the
// bounding box. The bounding box itself is no image, so it is corner-based
FillVector3D(cornerpoint, cornerpoint[0]-0.5, cornerpoint[1]-0.5, cornerpoint[2]-0.5);
}
return this->GetIndexToWorldTransform()->TransformPoint(cornerpoint);
}
mitk::Point3D mitk::BaseGeometry::GetCornerPoint(bool xFront, bool yFront, bool zFront) const
{
assert(this->IsBoundingBoxNull()==false);
BoundingBox::BoundsArrayType bounds = this->GetBoundingBox()->GetBounds();
Point3D cornerpoint;
cornerpoint[0] = (xFront ? bounds[0] : bounds[1]);
cornerpoint[1] = (yFront ? bounds[2] : bounds[3]);
cornerpoint[2] = (zFront ? bounds[4] : bounds[5]);
if(m_ImageGeometry)
{
// Here i have to adjust the 0.5 offset manually, because the cornerpoint is the corner of the
// bounding box. The bounding box itself is no image, so it is corner-based
FillVector3D(cornerpoint, cornerpoint[0]-0.5, cornerpoint[1]-0.5, cornerpoint[2]-0.5);
}
return this->GetIndexToWorldTransform()->TransformPoint(cornerpoint);
}
mitk::ScalarType mitk::BaseGeometry::GetExtentInMM(int direction) const
{
return m_IndexToWorldTransform->GetMatrix().GetVnlMatrix().get_column(direction).magnitude()*GetExtent(direction);
}
void mitk::BaseGeometry::SetExtentInMM(int direction, ScalarType extentInMM)
{
mitk::ModifiedLock lock(this);
ScalarType len = GetExtentInMM(direction);
if(fabs(len - extentInMM)>=mitk::eps)
{
AffineTransform3D::MatrixType::InternalMatrixType vnlmatrix;
vnlmatrix = m_IndexToWorldTransform->GetMatrix().GetVnlMatrix();
if(len>extentInMM)
vnlmatrix.set_column(direction, vnlmatrix.get_column(direction)/len*extentInMM);
else
vnlmatrix.set_column(direction, vnlmatrix.get_column(direction)*extentInMM/len);
Matrix3D matrix;
matrix = vnlmatrix;
m_IndexToWorldTransform->SetMatrix(matrix);
Modified();
}
PostSetExtentInMM(direction,extentInMM);
}
void mitk::BaseGeometry::PostSetExtentInMM(int /*direction*/, ScalarType /*extentInMM*/){};
bool mitk::BaseGeometry::IsInside(const mitk::Point3D& p) const
{
mitk::Point3D index;
WorldToIndex(p, index);
return IsIndexInside(index);
}
bool mitk::BaseGeometry::IsIndexInside(const mitk::Point3D& index) const
{
bool inside = false;
//if it is an image geometry, we need to convert the index to discrete values
//this is done by applying the rounding function also used in WorldToIndex (see line 323)
if (m_ImageGeometry)
{
mitk::Point3D discretIndex;
discretIndex[0]=itk::Math::RoundHalfIntegerUp<mitk::ScalarType>( index[0] );
discretIndex[1]=itk::Math::RoundHalfIntegerUp<mitk::ScalarType>( index[1] );
discretIndex[2]=itk::Math::RoundHalfIntegerUp<mitk::ScalarType>( index[2] );
inside = this->GetBoundingBox()->IsInside(discretIndex);
//we have to check if the index is at the upper border of each dimension,
// because the boundingbox is not centerbased
if (inside)
{
const BoundingBox::BoundsArrayType& bounds = this->GetBoundingBox()->GetBounds();
if((discretIndex[0] == bounds[1]) ||
(discretIndex[1] == bounds[3]) ||
(discretIndex[2] == bounds[5]))
inside = false;
}
}
else
inside = this->GetBoundingBox()->IsInside(index);
return inside;
}
void mitk::BaseGeometry::WorldToIndex(const mitk::Point3D &pt_mm, mitk::Point3D &pt_units) const
{
BackTransform(pt_mm, pt_units);
}
void mitk::BaseGeometry::WorldToIndex( const mitk::Vector3D &vec_mm, mitk::Vector3D &vec_units) const
{
BackTransform( vec_mm, vec_units);
}
void mitk::BaseGeometry::BackTransform(const mitk::Vector3D& in, mitk::Vector3D& out) const
{
// Get WorldToIndex transform
if (m_IndexToWorldTransformLastModified != m_IndexToWorldTransform->GetMTime())
{
m_InvertedTransform = TransformType::New();
if (!m_IndexToWorldTransform->GetInverse( m_InvertedTransform.GetPointer() ))
{
itkExceptionMacro( "Internal ITK matrix inversion error, cannot proceed." );
}
m_IndexToWorldTransformLastModified = m_IndexToWorldTransform->GetMTime();
}
// Check for valid matrix inversion
const TransformType::MatrixType& inverse = m_InvertedTransform->GetMatrix();
if(inverse.GetVnlMatrix().has_nans())
{
itkExceptionMacro( "Internal ITK matrix inversion error, cannot proceed. Matrix was: " << std::endl
<< m_IndexToWorldTransform->GetMatrix() << "Suggested inverted matrix is:" << std::endl
<< inverse );
}
// Transform vector
for (unsigned int i = 0; i < 3; i++)
{
out[i] = 0.0;
for (unsigned int j = 0; j < 3; j++)
{
out[i] += inverse[i][j]*in[j];
}
}
}
void mitk::BaseGeometry::BackTransform(const mitk::Point3D &in, mitk::Point3D& out) const
{
ScalarType temp[3];
unsigned int i, j;
const TransformType::OffsetType& offset = m_IndexToWorldTransform->GetOffset();
// Remove offset
for (j = 0; j < 3; j++)
{
temp[j] = in[j] - offset[j];
}
// Get WorldToIndex transform
if (m_IndexToWorldTransformLastModified != m_IndexToWorldTransform->GetMTime())
{
m_InvertedTransform = TransformType::New();
if (!m_IndexToWorldTransform->GetInverse( m_InvertedTransform.GetPointer() ))
{
itkExceptionMacro( "Internal ITK matrix inversion error, cannot proceed." );
}
m_IndexToWorldTransformLastModified = m_IndexToWorldTransform->GetMTime();
}
// Check for valid matrix inversion
const TransformType::MatrixType& inverse = m_InvertedTransform->GetMatrix();
if(inverse.GetVnlMatrix().has_nans())
{
itkExceptionMacro( "Internal ITK matrix inversion error, cannot proceed. Matrix was: " << std::endl
<< m_IndexToWorldTransform->GetMatrix() << "Suggested inverted matrix is:" << std::endl
<< inverse );
}
// Transform point
for (i = 0; i < 3; i++)
{
out[i] = 0.0;
for (j = 0; j < 3; j++)
{
out[i] += inverse[i][j]*temp[j];
}
}
}
mitk::VnlVector mitk::BaseGeometry::GetOriginVnl() const
{
return const_cast<Self*>(this)->m_Origin.GetVnlVector();
}
vtkLinearTransform* mitk::BaseGeometry::GetVtkTransform() const
{
return (vtkLinearTransform*)m_VtkIndexToWorldTransform;
}
void mitk::BaseGeometry::SetIdentity()
{
mitk::ModifiedLock lock(this);
m_IndexToWorldTransform->SetIdentity();
m_Origin.Fill(0);
CopySpacingFromTransform(m_IndexToWorldTransform, m_Spacing);
Modified();
TransferItkToVtkTransform();
}
void mitk::BaseGeometry::TransferVtkToItkTransform()
{
TransferVtkMatrixToItkTransform(m_VtkMatrix, m_IndexToWorldTransform.GetPointer());
CopySpacingFromTransform(m_IndexToWorldTransform, m_Spacing);
vtk2itk(m_IndexToWorldTransform->GetOffset(), m_Origin);
}
void mitk::BaseGeometry::Compose( const mitk::BaseGeometry::TransformType * other, bool pre )
{
mitk::ModifiedLock lock(this);
m_IndexToWorldTransform->Compose(other, pre);
CopySpacingFromTransform(m_IndexToWorldTransform, m_Spacing);
vtk2itk(m_IndexToWorldTransform->GetOffset(), m_Origin);
Modified();
TransferItkToVtkTransform();
}
void mitk::BaseGeometry::Compose( const vtkMatrix4x4 * vtkmatrix, bool pre )
{
mitk::BaseGeometry::TransformType::Pointer itkTransform = mitk::BaseGeometry::TransformType::New();
TransferVtkMatrixToItkTransform(vtkmatrix, itkTransform.GetPointer());
Compose(itkTransform, pre);
}
void mitk::BaseGeometry::Translate(const Vector3D & vector)
{
if((vector[0] != 0) || (vector[1] != 0) || (vector[2] != 0))
{
this->SetOrigin(m_Origin + vector);
}
}
void mitk::BaseGeometry::IndexToWorld(const mitk::Point3D &pt_units, mitk::Point3D &pt_mm) const
{
pt_mm = m_IndexToWorldTransform->TransformPoint(pt_units);
}
void mitk::BaseGeometry::IndexToWorld(const mitk::Vector3D &vec_units, mitk::Vector3D &vec_mm) const
{
vec_mm = m_IndexToWorldTransform->TransformVector(vec_units);
}
#include <vtkTransform.h>
void mitk::BaseGeometry::ExecuteOperation(Operation* operation)
{
mitk::ModifiedLock lock(this);
vtkTransform *vtktransform = vtkTransform::New();
vtktransform->SetMatrix(m_VtkMatrix);
switch (operation->GetOperationType())
{
case OpNOTHING:
break;
case OpMOVE:
{
mitk::PointOperation *pointOp = dynamic_cast<mitk::PointOperation *>(operation);
if (pointOp == NULL)
{
//mitk::StatusBar::GetInstance()->DisplayText("received wrong type of operation!See mitkAffineInteractor.cpp", 10000);
return;
}
mitk::Point3D newPos = pointOp->GetPoint();
ScalarType data[3];
vtktransform->GetPosition(data);
vtktransform->PostMultiply();
vtktransform->Translate(newPos[0], newPos[1], newPos[2]);
vtktransform->PreMultiply();
break;
}
case OpSCALE:
{
mitk::PointOperation *pointOp = dynamic_cast<mitk::PointOperation *>(operation);
if (pointOp == NULL)
{
//mitk::StatusBar::GetInstance()->DisplayText("received wrong type of operation!See mitkAffineInteractor.cpp", 10000);
return;
}
mitk::Point3D newScale = pointOp->GetPoint();
ScalarType data[3];
/* calculate new scale: newscale = oldscale * (oldscale + scaletoadd)/oldscale */
data[0] = 1 + (newScale[0] / GetMatrixColumn(0).magnitude());
data[1] = 1 + (newScale[1] / GetMatrixColumn(1).magnitude());
data[2] = 1 + (newScale[2] / GetMatrixColumn(2).magnitude());
mitk::Point3D center = const_cast<mitk::BoundingBox*>(m_BoundingBox.GetPointer())->GetCenter();
ScalarType pos[3];
vtktransform->GetPosition(pos);
vtktransform->PostMultiply();
vtktransform->Translate(-pos[0], -pos[1], -pos[2]);
vtktransform->Translate(-center[0], -center[1], -center[2]);
vtktransform->PreMultiply();
vtktransform->Scale(data[0], data[1], data[2]);
vtktransform->PostMultiply();
vtktransform->Translate(+center[0], +center[1], +center[2]);
vtktransform->Translate(pos[0], pos[1], pos[2]);
vtktransform->PreMultiply();
break;
}
case OpROTATE:
{
mitk::RotationOperation *rotateOp = dynamic_cast<mitk::RotationOperation *>(operation);
if (rotateOp == NULL)
{
//mitk::StatusBar::GetInstance()->DisplayText("received wrong type of operation!See mitkAffineInteractor.cpp", 10000);
return;
}
Vector3D rotationVector = rotateOp->GetVectorOfRotation();
Point3D center = rotateOp->GetCenterOfRotation();
ScalarType angle = rotateOp->GetAngleOfRotation();
vtktransform->PostMultiply();
vtktransform->Translate(-center[0], -center[1], -center[2]);
vtktransform->RotateWXYZ(angle, rotationVector[0], rotationVector[1], rotationVector[2]);
vtktransform->Translate(center[0], center[1], center[2]);
vtktransform->PreMultiply();
break;
}
case OpRESTOREPLANEPOSITION:
{
//Copy necessary to avoid vtk warning
vtkMatrix4x4* matrix = vtkMatrix4x4::New();
TransferItkTransformToVtkMatrix(dynamic_cast<mitk::RestorePlanePositionOperation*>(operation)->GetTransform().GetPointer(), matrix);
vtktransform->SetMatrix(matrix);
break;
}
case OpAPPLYTRANSFORMMATRIX:
{
ApplyTransformMatrixOperation *applyMatrixOp = dynamic_cast< ApplyTransformMatrixOperation* >( operation );
vtktransform->SetMatrix(applyMatrixOp->GetMatrix());
break;
}
default:
vtktransform->Delete();
return;
}
m_VtkMatrix->DeepCopy(vtktransform->GetMatrix());
TransferVtkToItkTransform();
Modified();
vtktransform->Delete();
}
mitk::VnlVector mitk::BaseGeometry::GetMatrixColumn(unsigned int direction) const
{
return m_IndexToWorldTransform->GetMatrix().GetVnlMatrix().get_column(direction);
}
mitk::BoundingBox::Pointer mitk::BaseGeometry::CalculateBoundingBoxRelativeToTransform(const mitk::AffineTransform3D* transform) const
{
mitk::BoundingBox::PointsContainer::Pointer pointscontainer=mitk::BoundingBox::PointsContainer::New();
mitk::BoundingBox::PointIdentifier pointid=0;
unsigned char i;
if(transform!=NULL)
{
mitk::AffineTransform3D::Pointer inverse = mitk::AffineTransform3D::New();
transform->GetInverse(inverse);
for(i=0; i<8; ++i)
pointscontainer->InsertElement( pointid++, inverse->TransformPoint( GetCornerPoint(i) ));
}
else
{
for(i=0; i<8; ++i)
pointscontainer->InsertElement( pointid++, GetCornerPoint(i) );
}
mitk::BoundingBox::Pointer result = mitk::BoundingBox::New();
result->SetPoints(pointscontainer);
result->ComputeBoundingBox();
return result;
}
//void mitk::BaseGeometry::SetTimeBounds(const TimeBounds& timebounds)
//{
// mitk::ModifiedLock lock(this);
//
// if(m_TimeBounds != timebounds)
// {
// m_TimeBounds = timebounds;
// Modified();
// }
// PostSetTimeBounds(timebounds);
//}
//
//void mitk::BaseGeometry::PostSetTimeBounds(const TimeBounds& timebounds)
//{}
const std::string mitk::BaseGeometry::GetTransformAsString( TransformType* transformType )
{
std::ostringstream out;
out << '[';
for( int i=0; i<3; ++i )
{
out << '[';
for( int j=0; j<3; ++j )
out << transformType->GetMatrix().GetVnlMatrix().get(i, j) << ' ';
out << ']';
}
out << "][";
for( int i=0; i<3; ++i )
out << transformType->GetOffset()[i] << ' ';
out << "]\0";
return out.str();
}
void mitk::BaseGeometry::SetIndexToWorldTransformByVtkMatrix(vtkMatrix4x4* vtkmatrix)
{
m_VtkMatrix->DeepCopy(vtkmatrix);
TransferVtkToItkTransform();
}
void mitk::BaseGeometry::WorldToIndex(const mitk::Point3D & /*atPt3d_mm*/, const mitk::Vector3D &vec_mm, mitk::Vector3D &vec_units) const
{
MITK_WARN<<"Warning! Call of the deprecated function BaseGeometry::WorldToIndex(point, vec, vec). Use BaseGeometry::WorldToIndex(vec, vec) instead!";
//BackTransform(atPt3d_mm, vec_mm, vec_units);
this->WorldToIndex(vec_mm, vec_units);
}
void mitk::BaseGeometry::IndexToWorld(const mitk::Point3D &/*atPt3d_units*/, const mitk::Vector3D &vec_units, mitk::Vector3D &vec_mm) const
{
MITK_WARN<<"Warning! Call of the deprecated function BaseGeometry::IndexToWorld(point, vec, vec). Use BaseGeometry::IndexToWorld(vec, vec) instead!";
//vec_mm = m_IndexToWorldTransform->TransformVector(vec_units);
this->IndexToWorld(vec_units, vec_mm);
}
void mitk::BaseGeometry::BackTransform(const mitk::Point3D &/*at*/, const mitk::Vector3D &in, mitk::Vector3D& out) const
{
MITK_INFO<<"Warning! Call of the deprecated function BaseGeometry::BackTransform(point, vec, vec). Use BaseGeometry::BackTransform(vec, vec) instead!";
//// Get WorldToIndex transform
//if (m_IndexToWorldTransformLastModified != m_IndexToWorldTransform->GetMTime())
//{
// m_InvertedTransform = TransformType::New();
// if (!m_IndexToWorldTransform->GetInverse( m_InvertedTransform.GetPointer() ))
// {
// itkExceptionMacro( "Internal ITK matrix inversion error, cannot proceed." );
// }
// m_IndexToWorldTransformLastModified = m_IndexToWorldTransform->GetMTime();
//}
//// Check for valid matrix inversion
//const TransformType::MatrixType& inverse = m_InvertedTransform->GetMatrix();
//if(inverse.GetVnlMatrix().has_nans())
//{
// itkExceptionMacro( "Internal ITK matrix inversion error, cannot proceed. Matrix was: " << std::endl
// << m_IndexToWorldTransform->GetMatrix() << "Suggested inverted matrix is:" << std::endl
// << inverse );
//}
//// Transform vector
//for (unsigned int i = 0; i < 3; i++)
//{
// out[i] = 0.0;
// for (unsigned int j = 0; j < 3; j++)
// {
// out[i] += inverse[i][j]*in[j];
// }
//}
this->BackTransform(in, out);
}
vtkMatrix4x4* mitk::BaseGeometry::GetVtkMatrix(){
return m_VtkMatrix;
}
bool mitk::BaseGeometry::IsBoundingBoxNull() const{
return m_BoundingBox.IsNull();
}
bool mitk::BaseGeometry::IsIndexToWorldTransformNull() const{
return m_IndexToWorldTransform.IsNull();
}
void
mitk::BaseGeometry::ChangeImageGeometryConsideringOriginOffset( const bool isAnImageGeometry )
{
// If Geometry is switched to ImageGeometry, you have to put an offset to the origin, because
// imageGeometries origins are pixel-center-based
// ... and remove the offset, if you switch an imageGeometry back to a normal geometry
// For more information please see the Geometry documentation page
if(m_ImageGeometry == isAnImageGeometry)
return;
const BoundingBox::BoundsArrayType& boundsarray =
this->GetBoundingBox()->GetBounds();
Point3D originIndex;
FillVector3D(originIndex, boundsarray[0], boundsarray[2], boundsarray[4]);
if(isAnImageGeometry == true)
FillVector3D( originIndex,
originIndex[0] + 0.5,
originIndex[1] + 0.5,
originIndex[2] + 0.5 );
else
FillVector3D( originIndex,
originIndex[0] - 0.5,
originIndex[1] - 0.5,
originIndex[2] - 0.5 );
Point3D originWorld;
originWorld = GetIndexToWorldTransform()
->TransformPoint( originIndex );
// instead could as well call IndexToWorld(originIndex,originWorld);
SetOrigin(originWorld);
this->SetImageGeometry(isAnImageGeometry);
}
//itk::LightObject::Pointer mitk::BaseGeometry::InternalClone() const
//{
// Self::Pointer newGeometry = new Self(*this);
// newGeometry->UnRegister();
// return newGeometry.GetPointer();
//}
void mitk::BaseGeometry::PrintSelf(std::ostream& os, itk::Indent indent) const
{
os << indent << " IndexToWorldTransform: ";
if(this->IsIndexToWorldTransformNull())
os << "NULL" << std::endl;
else
{
// from itk::MatrixOffsetTransformBase
unsigned int i, j;
os << std::endl;
os << indent << "Matrix: " << std::endl;
for (i = 0; i < 3; i++)
{
os << indent.GetNextIndent();
for (j = 0; j < 3; j++)
{
os << this->GetIndexToWorldTransform()->GetMatrix()[i][j] << " ";
}
os << std::endl;
}
os << indent << "Offset: " << this->GetIndexToWorldTransform()->GetOffset() << std::endl;
os << indent << "Center: " << this->GetIndexToWorldTransform()->GetCenter() << std::endl;
os << indent << "Translation: " << this->GetIndexToWorldTransform()->GetTranslation() << std::endl;
os << indent << "Inverse: " << std::endl;
for (i = 0; i < 3; i++)
{
os << indent.GetNextIndent();
for (j = 0; j < 3; j++)
{
os << this->GetIndexToWorldTransform()->GetInverseMatrix()[i][j] << " ";
}
os << std::endl;
}
// from itk::ScalableAffineTransform
os << indent << "Scale : ";
for (i = 0; i < 3; i++)
{
os << this->GetIndexToWorldTransform()->GetScale()[i] << " ";
}
os << std::endl;
}
os << indent << " BoundingBox: ";
if(this->IsBoundingBoxNull())
os << "NULL" << std::endl;
else
{
os << indent << "( ";
for (unsigned int i=0; i<3; i++)
{
os << this->GetBoundingBox()->GetBounds()[2*i] << "," << this->GetBoundingBox()->GetBounds()[2*i+1] << " ";
}
os << " )" << std::endl;
}
os << indent << " Origin: " << this->GetOrigin() << std::endl;
os << indent << " ImageGeometry: " << this->GetImageGeometry() << std::endl;
os << indent << " Spacing: " << this->GetSpacing() << std::endl;
//os << indent << " TimeBounds: " << this->GetTimeBounds() << std::endl;
}
void mitk::BaseGeometry::Modified() const{
if(!m_ModifiedLockFlag)
Superclass::Modified();
else
m_ModifiedCalledFlag = true;
}
bool mitk::Equal( const mitk::BaseGeometry::BoundingBoxType *leftHandSide, const mitk::BaseGeometry::BoundingBoxType *rightHandSide, ScalarType eps, bool verbose )
{
if(( leftHandSide == NULL) || ( rightHandSide == NULL ))
{
MITK_ERROR << "mitk::Equal( const mitk::Geometry3D::BoundingBoxType *leftHandSide, const mitk::Geometry3D::BoundingBoxType *rightHandSide, ScalarType eps, bool verbose ) does not with NULL pointer input.";
return false;
}
return Equal( *leftHandSide, *rightHandSide, eps, verbose);
}
bool mitk::Equal( const mitk::BaseGeometry::BoundingBoxType& leftHandSide, const mitk::BaseGeometry::BoundingBoxType& rightHandSide, ScalarType eps, bool verbose )
{
bool result = true;
BaseGeometry::BoundsArrayType rightBounds = rightHandSide.GetBounds();
BaseGeometry::BoundsArrayType leftBounds = leftHandSide.GetBounds();
BaseGeometry::BoundsArrayType::Iterator itLeft = leftBounds.Begin();
for( BaseGeometry::BoundsArrayType::Iterator itRight = rightBounds.Begin(); itRight != rightBounds.End(); ++itRight)
{
if(( !mitk::Equal( *itLeft, *itRight, eps )) )
{
if(verbose)
{
MITK_INFO << "[( Geometry3D::BoundingBoxType )] bounds are not equal.";
MITK_INFO << "rightHandSide is " << setprecision(12) << *itRight << " : leftHandSide is " << *itLeft << " and tolerance is " << eps;
}
result = false;
}
itLeft++;
}
return result;
}
bool mitk::Equal(const mitk::BaseGeometry *leftHandSide, const mitk::BaseGeometry *rightHandSide, ScalarType eps, bool verbose)
{
if(( leftHandSide == NULL) || ( rightHandSide == NULL ))
{
MITK_ERROR << "mitk::Equal(const mitk::Geometry3D *leftHandSide, const mitk::Geometry3D *rightHandSide, ScalarType eps, bool verbose) does not with NULL pointer input.";
return false;
}
return Equal( *leftHandSide, *rightHandSide, eps, verbose);
}
bool mitk::Equal(const mitk::BaseGeometry& leftHandSide, const mitk::BaseGeometry& rightHandSide, ScalarType eps, bool verbose)
{
bool result = true;
//Compare spacings
if( !mitk::Equal( leftHandSide.GetSpacing(), rightHandSide.GetSpacing(), eps ) )
{
if(verbose)
{
MITK_INFO << "[( Geometry3D )] Spacing differs.";
MITK_INFO << "rightHandSide is " << setprecision(12) << rightHandSide.GetSpacing() << " : leftHandSide is " << leftHandSide.GetSpacing() << " and tolerance is " << eps;
}
result = false;
}
//Compare Origins
if( !mitk::Equal( leftHandSide.GetOrigin(), rightHandSide.GetOrigin(), eps ) )
{
if(verbose)
{
MITK_INFO << "[( Geometry3D )] Origin differs.";
MITK_INFO << "rightHandSide is " << setprecision(12) << rightHandSide.GetOrigin() << " : leftHandSide is " << leftHandSide.GetOrigin() << " and tolerance is " << eps;
}
result = false;
}
//Compare Axis and Extents
for( unsigned int i=0; i<3; ++i)
{
if( !mitk::Equal( leftHandSide.GetAxisVector(i), rightHandSide.GetAxisVector(i), eps))
{
if(verbose)
{
MITK_INFO << "[( Geometry3D )] AxisVector #" << i << " differ";
MITK_INFO << "rightHandSide is " << setprecision(12) << rightHandSide.GetAxisVector(i) << " : leftHandSide is " << leftHandSide.GetAxisVector(i) << " and tolerance is " << eps;
}
result = false;
}
if( !mitk::Equal( leftHandSide.GetExtent(i), rightHandSide.GetExtent(i), eps) )
{
if(verbose)
{
MITK_INFO << "[( Geometry3D )] Extent #" << i << " differ";
MITK_INFO << "rightHandSide is " << setprecision(12) << rightHandSide.GetExtent(i) << " : leftHandSide is " << leftHandSide.GetExtent(i) << " and tolerance is " << eps;
}
result = false;
}
}
//Compare ImageGeometry Flag
if( rightHandSide.GetImageGeometry() != leftHandSide.GetImageGeometry() )
{
if(verbose)
{
MITK_INFO << "[( Geometry3D )] GetImageGeometry is different.";
MITK_INFO << "rightHandSide is " << rightHandSide.GetImageGeometry() << " : leftHandSide is " << leftHandSide.GetImageGeometry();
}
result = false;
}
//Compare BoundingBoxes
if( !mitk::Equal( *leftHandSide.GetBoundingBox(), *rightHandSide.GetBoundingBox(), eps, verbose) )
{
result = false;
}
//Compare IndexToWorldTransform Matrix
if( !mitk::Equal( *leftHandSide.GetIndexToWorldTransform(), *rightHandSide.GetIndexToWorldTransform(), eps, verbose) )
{
result = false;
}
return result;
}
bool mitk::Equal(const BaseGeometry::TransformType *leftHandSide, const BaseGeometry::TransformType *rightHandSide, ScalarType eps, bool verbose )
{
if(( leftHandSide == NULL) || ( rightHandSide == NULL ))
{
MITK_ERROR << "mitk::Equal(const Geometry3D::TransformType *leftHandSide, const Geometry3D::TransformType *rightHandSide, ScalarType eps, bool verbose ) does not with NULL pointer input.";
return false;
}
return Equal( *leftHandSide, *rightHandSide, eps, verbose);
}
bool mitk::Equal(const BaseGeometry::TransformType& leftHandSide, const BaseGeometry::TransformType& rightHandSide, ScalarType eps, bool verbose )
{
//Compare IndexToWorldTransform Matrix
if( !mitk::MatrixEqualElementWise( leftHandSide.GetMatrix(),
rightHandSide.GetMatrix() ) )
{
if(verbose)
{
MITK_INFO << "[( Geometry3D::TransformType )] Index to World Transformation matrix differs.";
MITK_INFO << "rightHandSide is " << setprecision(12) << rightHandSide.GetMatrix() << " : leftHandSide is " << leftHandSide.GetMatrix() << " and tolerance is " << eps;
}
return false;
}
return true;
}
diff --git a/Core/Code/DataManagement/mitkDisplayGeometry.cpp b/Core/Code/DataManagement/mitkDisplayGeometry.cpp
index ab93699403..7ad44c83e7 100644
--- a/Core/Code/DataManagement/mitkDisplayGeometry.cpp
+++ b/Core/Code/DataManagement/mitkDisplayGeometry.cpp
@@ -1,613 +1,614 @@
/*===================================================================
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 "mitkDisplayGeometry.h"
itk::LightObject::Pointer mitk::DisplayGeometry::InternalClone() const
{
// itkExceptionMacro(<<"calling mitk::DisplayGeometry::Clone does not make much sense.");
DisplayGeometry* returnValue = const_cast<DisplayGeometry *>(this);
return returnValue;
}
bool mitk::DisplayGeometry::IsValid() const
{
return m_WorldGeometry.IsNotNull() && m_WorldGeometry->IsValid();
}
unsigned long mitk::DisplayGeometry::GetMTime() const
{
if((m_WorldGeometry.IsNotNull()) && (PlaneGeometry::GetMTime() < m_WorldGeometry->GetMTime()))
{
Modified();
}
return PlaneGeometry::GetMTime();
}
//const mitk::TimeBounds& mitk::DisplayGeometry::GetTimeBounds() const
//{
// if(m_WorldGeometry.IsNull())
// {
// return m_TimeBounds;
// }
//
// return m_WorldGeometry->GetTimeBounds();
//}
// size definition methods
void mitk::DisplayGeometry::SetWorldGeometry(const PlaneGeometry* aWorldGeometry)
{
m_WorldGeometry = aWorldGeometry;
Modified();
}
bool mitk::DisplayGeometry::SetOriginInMM(const Vector2D& origin_mm)
{
m_OriginInMM = origin_mm;
WorldToDisplay(m_OriginInMM, m_OriginInDisplayUnits);
Modified();
return !this->RefitVisibleRect();
}
mitk::Vector2D mitk::DisplayGeometry::GetOriginInMM() const
{
return m_OriginInMM;
}
mitk::Vector2D mitk::DisplayGeometry::GetOriginInDisplayUnits() const
{
return m_OriginInDisplayUnits;
}
void mitk::DisplayGeometry::SetSizeInDisplayUnits(unsigned int width, unsigned int height, bool keepDisplayedRegion)
{
Vector2D oldSizeInMM( m_SizeInMM );
Point2D oldCenterInMM;
if(keepDisplayedRegion)
{
Point2D centerInDisplayUnits;
centerInDisplayUnits[0] = m_SizeInDisplayUnits[0]*0.5;
centerInDisplayUnits[1] = m_SizeInDisplayUnits[1]*0.5;
DisplayToWorld(centerInDisplayUnits, oldCenterInMM);
}
m_SizeInDisplayUnits[0]=width;
m_SizeInDisplayUnits[1]=height;
if(m_SizeInDisplayUnits[0] <= 0)
m_SizeInDisplayUnits[0] = 1;
if(m_SizeInDisplayUnits[1] <= 0)
m_SizeInDisplayUnits[1] = 1;
DisplayToWorld(m_SizeInDisplayUnits, m_SizeInMM);
if(keepDisplayedRegion)
{
Point2D positionOfOldCenterInCurrentDisplayUnits;
WorldToDisplay(oldCenterInMM, positionOfOldCenterInCurrentDisplayUnits);
Point2D currentNewCenterInDisplayUnits;
currentNewCenterInDisplayUnits[0] = m_SizeInDisplayUnits[0]*0.5;
currentNewCenterInDisplayUnits[1] = m_SizeInDisplayUnits[1]*0.5;
Vector2D shift;
shift=positionOfOldCenterInCurrentDisplayUnits-currentNewCenterInDisplayUnits;
MoveBy(shift);
Zoom(m_SizeInMM.GetNorm()/oldSizeInMM.GetNorm(), currentNewCenterInDisplayUnits);
}
Modified();
}
mitk::Vector2D mitk::DisplayGeometry::GetSizeInDisplayUnits() const
{
return m_SizeInDisplayUnits;
}
mitk::Vector2D mitk::DisplayGeometry::GetSizeInMM() const
{
return m_SizeInMM;
}
unsigned int mitk::DisplayGeometry::GetDisplayWidth() const
{
assert(m_SizeInDisplayUnits[0] >= 0);
return (unsigned int)m_SizeInDisplayUnits[0];
}
unsigned int mitk::DisplayGeometry::GetDisplayHeight() const
{
assert(m_SizeInDisplayUnits[1] >= 0);
return (unsigned int)m_SizeInDisplayUnits[1];
}
// zooming, panning, restriction of both
void mitk::DisplayGeometry::SetConstrainZoomingAndPanning(bool constrain)
{
m_ConstrainZoomingAndPanning = constrain;
if (m_ConstrainZoomingAndPanning)
{
this->RefitVisibleRect();
}
}
bool mitk::DisplayGeometry::GetConstrainZommingAndPanning() const
{
return m_ConstrainZoomingAndPanning;
}
bool mitk::DisplayGeometry::SetScaleFactor(ScalarType mmPerDisplayUnit)
{
if(mmPerDisplayUnit<0.0001)
{
mmPerDisplayUnit=0.0001;
}
m_ScaleFactorMMPerDisplayUnit = mmPerDisplayUnit;
assert(m_ScaleFactorMMPerDisplayUnit < itk::NumericTraits<mitk::ScalarType>::infinity());
DisplayToWorld(m_SizeInDisplayUnits, m_SizeInMM);
return !this->RefitVisibleRect();
}
mitk::ScalarType mitk::DisplayGeometry::GetScaleFactorMMPerDisplayUnit() const
{
return m_ScaleFactorMMPerDisplayUnit;
}
// Zooms with a factor (1.0=identity) around the specified center in display units
bool mitk::DisplayGeometry::Zoom(ScalarType factor, const Point2D& centerInDisplayUnits)
{
assert(factor > 0);
if ( SetScaleFactor(m_ScaleFactorMMPerDisplayUnit/factor) )
{
return SetOriginInMM(m_OriginInMM-centerInDisplayUnits.GetVectorFromOrigin()*(1-factor)*m_ScaleFactorMMPerDisplayUnit);
}
else
{
return false;
}
}
// Zooms with a factor (1.0=identity) around the specified center, but tries (if its within view contraints) to match the center in display units with the center in world coordinates.
bool mitk::DisplayGeometry::ZoomWithFixedWorldCoordinates(ScalarType factor, const Point2D& focusDisplayUnits, const Point2D& focusUnitsInMM )
{
assert(factor > 0);
SetScaleFactor(m_ScaleFactorMMPerDisplayUnit/factor);
SetOriginInMM(focusUnitsInMM.GetVectorFromOrigin()-focusDisplayUnits.GetVectorFromOrigin()*m_ScaleFactorMMPerDisplayUnit);
return true;
}
bool mitk::DisplayGeometry::MoveBy(const Vector2D& shiftInDisplayUnits)
{
SetOriginInMM(m_OriginInMM+shiftInDisplayUnits*m_ScaleFactorMMPerDisplayUnit);
Modified();
return !this->RefitVisibleRect();
}
void mitk::DisplayGeometry::Fit()
{
if((m_WorldGeometry.IsNull()) || (m_WorldGeometry->IsValid() == false)) return;
/// \FIXME: try to remove all the casts
int width=(int)m_SizeInDisplayUnits[0];
int height=(int)m_SizeInDisplayUnits[1];
ScalarType w = width;
ScalarType h = height;
const ScalarType& widthInMM = m_WorldGeometry->GetExtentInMM(0);
const ScalarType& heightInMM = m_WorldGeometry->GetExtentInMM(1);
ScalarType aspRatio=((ScalarType)widthInMM)/heightInMM;
ScalarType x = (ScalarType)w/widthInMM;
ScalarType y = (ScalarType)h/heightInMM;
if (x > y)
{
w = (int) (aspRatio*h);
}
else
{
h = (int) (w/aspRatio);
}
if(w>0)
{
SetScaleFactor(widthInMM/w);
}
Vector2D origin_display;
origin_display[0]=-(width-w)/2.0;
origin_display[1]=-(height-h)/2.0;
SetOriginInMM(origin_display*m_ScaleFactorMMPerDisplayUnit);
this->RefitVisibleRect();
Modified();
}
// conversion methods
void mitk::DisplayGeometry::DisplayToWorld(const Point2D &pt_display, Point2D &pt_mm) const
{
pt_mm[0]=m_ScaleFactorMMPerDisplayUnit*pt_display[0]+m_OriginInMM[0];
pt_mm[1]=m_ScaleFactorMMPerDisplayUnit*pt_display[1]+m_OriginInMM[1];
}
void mitk::DisplayGeometry::WorldToDisplay(const Point2D &pt_mm, Point2D &pt_display) const
{
pt_display[0]=(pt_mm[0]-m_OriginInMM[0])*(1.0/m_ScaleFactorMMPerDisplayUnit);
pt_display[1]=(pt_mm[1]-m_OriginInMM[1])*(1.0/m_ScaleFactorMMPerDisplayUnit);
}
void mitk::DisplayGeometry::DisplayToWorld(const Vector2D &vec_display, Vector2D &vec_mm) const
{
vec_mm=vec_display*m_ScaleFactorMMPerDisplayUnit;
}
void mitk::DisplayGeometry::WorldToDisplay(const Vector2D &vec_mm, Vector2D &vec_display) const
{
vec_display=vec_mm*(1.0/m_ScaleFactorMMPerDisplayUnit);
}
void mitk::DisplayGeometry::ULDisplayToMM(const Point2D &pt_ULdisplay, Point2D &pt_mm) const
{
ULDisplayToDisplay(pt_ULdisplay, pt_mm);
DisplayToWorld(pt_mm, pt_mm);
}
void mitk::DisplayGeometry::MMToULDisplay(const Point2D &pt_mm, Point2D &pt_ULdisplay) const
{
WorldToDisplay(pt_mm, pt_ULdisplay);
DisplayToULDisplay(pt_ULdisplay, pt_ULdisplay);
}
void mitk::DisplayGeometry::ULDisplayToMM(const Vector2D &vec_ULdisplay, Vector2D &vec_mm) const
{
ULDisplayToDisplay(vec_ULdisplay, vec_mm);
DisplayToWorld(vec_mm, vec_mm);
}
void mitk::DisplayGeometry::MMToULDisplay(const Vector2D &vec_mm, Vector2D &vec_ULdisplay) const
{
WorldToDisplay(vec_mm, vec_ULdisplay);
DisplayToULDisplay(vec_ULdisplay, vec_ULdisplay);
}
void mitk::DisplayGeometry::ULDisplayToDisplay(const Point2D &pt_ULdisplay, Point2D &pt_display) const
{
pt_display[0]=pt_ULdisplay[0];
pt_display[1]=GetDisplayHeight()-pt_ULdisplay[1];
}
void mitk::DisplayGeometry::DisplayToULDisplay(const Point2D &pt_display, Point2D &pt_ULdisplay) const
{
ULDisplayToDisplay(pt_display, pt_ULdisplay);
}
void mitk::DisplayGeometry::ULDisplayToDisplay(const Vector2D &vec_ULdisplay, Vector2D &vec_display) const
{
vec_display[0]= vec_ULdisplay[0];
vec_display[1]=-vec_ULdisplay[1];
}
void mitk::DisplayGeometry::DisplayToULDisplay(const Vector2D &vec_display, Vector2D &vec_ULdisplay) const
{
ULDisplayToDisplay(vec_display, vec_ULdisplay);
}
bool mitk::DisplayGeometry::Project(const Point3D &pt3d_mm, Point3D &projectedPt3d_mm) const
{
if(m_WorldGeometry.IsNotNull())
{
return m_WorldGeometry->Project(pt3d_mm, projectedPt3d_mm);
}
else
{
return false;
}
}
bool mitk::DisplayGeometry::Project(const Point3D & atPt3d_mm, const Vector3D &vec3d_mm, Vector3D &projectedVec3d_mm) const
{
if(m_WorldGeometry.IsNotNull())
{
return m_WorldGeometry->Project(atPt3d_mm, vec3d_mm, projectedVec3d_mm);
}
else
{
return false;
}
}
bool mitk::DisplayGeometry::Project(const Vector3D &vec3d_mm, Vector3D &projectedVec3d_mm) const
{
if(m_WorldGeometry.IsNotNull())
{
return m_WorldGeometry->Project(vec3d_mm, projectedVec3d_mm);
}
else
{
return false;
}
}
bool mitk::DisplayGeometry::Map(const Point3D &pt3d_mm, Point2D &pt2d_mm) const
{
if(m_WorldGeometry.IsNotNull())
{
return m_WorldGeometry->Map(pt3d_mm, pt2d_mm);
}
else
{
return false;
}
}
void mitk::DisplayGeometry::Map(const Point2D &pt2d_mm, Point3D &pt3d_mm) const
{
if(m_WorldGeometry.IsNull()) return;
m_WorldGeometry->Map(pt2d_mm, pt3d_mm);
}
bool mitk::DisplayGeometry::Map(const Point3D & atPt3d_mm, const Vector3D &vec3d_mm, Vector2D &vec2d_mm) const
{
if(m_WorldGeometry.IsNotNull())
{
return m_WorldGeometry->Map(atPt3d_mm, vec3d_mm, vec2d_mm);
}
else
{
return false;
}
}
void mitk::DisplayGeometry::Map(const Point2D & atPt2d_mm, const Vector2D &vec2d_mm, Vector3D &vec3d_mm) const
{
if(m_WorldGeometry.IsNull()) return;
m_WorldGeometry->Map(atPt2d_mm, vec2d_mm, vec3d_mm);
}
// protected methods
mitk::DisplayGeometry::DisplayGeometry()
- :m_ScaleFactorMMPerDisplayUnit(1.0)
+ : PlaneGeometry()
+ ,m_ScaleFactorMMPerDisplayUnit(1.0)
,m_WorldGeometry(NULL)
,m_ConstrainZoomingAndPanning(true)
,m_MaxWorldViewPercentage(1.0)
,m_MinWorldViewPercentage(0.1)
{
m_OriginInMM.Fill(0.0);
m_OriginInDisplayUnits.Fill(0.0);
m_SizeInMM.Fill(1.0);
m_SizeInDisplayUnits.Fill(10.0);
}
mitk::DisplayGeometry::~DisplayGeometry()
{
}
bool mitk::DisplayGeometry::RefitVisibleRect()
{
// do nothing if not asked to
if (!m_ConstrainZoomingAndPanning) return false;
// don't allow recursion (need to be fixed, singleton)
static bool inRecalculate = false;
if (inRecalculate) return false;
inRecalculate = true;
// rename some basic measures of the current viewport and world geometry (MM = milimeters Px = Pixels = display units)
float displayXMM = m_OriginInMM[0];
float displayYMM = m_OriginInMM[1];
float displayWidthPx = m_SizeInDisplayUnits[0];
float displayHeightPx = m_SizeInDisplayUnits[1];
float displayWidthMM = m_SizeInDisplayUnits[0] * m_ScaleFactorMMPerDisplayUnit;
float displayHeightMM = m_SizeInDisplayUnits[1] * m_ScaleFactorMMPerDisplayUnit;
float worldWidthMM = m_WorldGeometry->GetExtentInMM(0);
float worldHeightMM = m_WorldGeometry->GetExtentInMM(1);
// reserve variables for the correction logic to save a corrected origin and zoom factor
Vector2D newOrigin = m_OriginInMM;
bool correctPanning = false;
float newScaleFactor = m_ScaleFactorMMPerDisplayUnit;
bool correctZooming = false;
// start of the correction logic
// zoom to big means:
// at a given percentage of the world's width/height should be visible. Otherwise
// the whole screen could show only one pixel
//
// zoom to small means:
// zooming out should be limited at the point where the smaller of the world's sides is completely visible
bool zoomXtooSmall = displayWidthPx * m_ScaleFactorMMPerDisplayUnit > m_MaxWorldViewPercentage * worldWidthMM;
bool zoomXtooBig = displayWidthPx * m_ScaleFactorMMPerDisplayUnit < m_MinWorldViewPercentage * worldWidthMM;
bool zoomYtooSmall = displayHeightPx * m_ScaleFactorMMPerDisplayUnit > m_MaxWorldViewPercentage * worldHeightMM;
bool zoomYtooBig = displayHeightPx * m_ScaleFactorMMPerDisplayUnit < m_MinWorldViewPercentage * worldHeightMM;
// constrain zooming in both direction
if ( zoomXtooBig && zoomYtooBig)
{
double fx = worldWidthMM * m_MinWorldViewPercentage / displayWidthPx;
double fy = worldHeightMM * m_MinWorldViewPercentage / displayHeightPx;
newScaleFactor = fx < fy ? fx : fy;
correctZooming = true;
}
// constrain zooming in x direction
else if ( zoomXtooBig )
{
newScaleFactor = worldWidthMM * m_MinWorldViewPercentage / displayWidthPx;
correctZooming = true;
}
// constrain zooming in y direction
else if ( zoomYtooBig )
{
newScaleFactor = worldHeightMM * m_MinWorldViewPercentage / displayHeightPx;
correctZooming = true;
}
// constrain zooming out
// we stop zooming out at these situations:
//
// *** display
// --- image
//
// **********************
// * * x side maxed out
// * *
// *--------------------*
// *| |*
// *| |*
// *--------------------*
// * *
// * *
// * *
// **********************
//
// **********************
// * |------| * y side maxed out
// * | | *
// * | | *
// * | | *
// * | | *
// * | | *
// * | | *
// * | | *
// * |------| *
// **********************
//
// In both situations we center the not-maxed out direction
//
if ( zoomXtooSmall && zoomYtooSmall )
{
// determine and set the bigger scale factor
float fx = worldWidthMM * m_MaxWorldViewPercentage / displayWidthPx;
float fy = worldHeightMM * m_MaxWorldViewPercentage / displayHeightPx;
newScaleFactor = fx > fy ? fx : fy;
correctZooming = true;
}
// actually execute correction
if (correctZooming)
{
SetScaleFactor(newScaleFactor);
}
displayWidthMM = m_SizeInDisplayUnits[0] * m_ScaleFactorMMPerDisplayUnit;
displayHeightMM = m_SizeInDisplayUnits[1] * m_ScaleFactorMMPerDisplayUnit;
// constrain panning
if(worldWidthMM<displayWidthMM)
{
// zoomed out too much in x (but tolerated because y is still ok)
// --> center x
newOrigin[0] = (worldWidthMM - displayWidthMM) / 2.0;
correctPanning = true;
}
else
{
// make sure left display border inside our world
if (displayXMM < 0)
{
newOrigin[0] = 0;
correctPanning = true;
}
// make sure right display border inside our world
else if (displayXMM + displayWidthMM > worldWidthMM)
{
newOrigin[0] = worldWidthMM - displayWidthMM;
correctPanning = true;
}
}
if (worldHeightMM<displayHeightMM)
{
// zoomed out too much in y (but tolerated because x is still ok)
// --> center y
newOrigin[1] = (worldHeightMM - displayHeightMM) / 2.0;
correctPanning = true;
}
else
{
// make sure top display border inside our world
if (displayYMM + displayHeightMM > worldHeightMM)
{
newOrigin[1] = worldHeightMM - displayHeightMM;
correctPanning = true;
}
// make sure bottom display border inside our world
else
if (displayYMM < 0)
{
newOrigin[1] = 0;
correctPanning = true;
}
}
if (correctPanning)
{
SetOriginInMM( newOrigin );
}
inRecalculate = false;
if ( correctPanning || correctZooming )
{
Modified();
}
// return true if any correction has been made
return correctPanning || correctZooming;
}
void mitk::DisplayGeometry::PrintSelf(std::ostream& os, itk::Indent indent) const
{
if(m_WorldGeometry.IsNull())
{
os << indent << " WorldGeometry: " << "NULL" << std::endl;
}
else
{
m_WorldGeometry->Print(os, indent);
os << indent << " OriginInMM: " << m_OriginInMM << std::endl;
os << indent << " OriginInDisplayUnits: " << m_OriginInDisplayUnits << std::endl;
os << indent << " SizeInMM: " << m_SizeInMM << std::endl;
os << indent << " SizeInDisplayUnits: " << m_SizeInDisplayUnits << std::endl;
os << indent << " ScaleFactorMMPerDisplayUni: " << m_ScaleFactorMMPerDisplayUnit << std::endl;
}
Superclass::PrintSelf(os,indent);
}
diff --git a/Core/Code/DataManagement/mitkLandmarkProjectorBasedCurvedGeometry.cpp b/Core/Code/DataManagement/mitkLandmarkProjectorBasedCurvedGeometry.cpp
index 9e8ac72940..7aa5dc7ef1 100644
--- a/Core/Code/DataManagement/mitkLandmarkProjectorBasedCurvedGeometry.cpp
+++ b/Core/Code/DataManagement/mitkLandmarkProjectorBasedCurvedGeometry.cpp
@@ -1,83 +1,83 @@
/*===================================================================
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 "mitkLandmarkProjectorBasedCurvedGeometry.h"
#include <vtkAbstractTransform.h>
mitk::LandmarkProjectorBasedCurvedGeometry::LandmarkProjectorBasedCurvedGeometry()
- : m_LandmarkProjector(NULL), m_InterpolatingAbstractTransform(NULL), m_TargetLandmarks(NULL)
+ : AbstractTransformGeometry(), m_LandmarkProjector(NULL), m_InterpolatingAbstractTransform(NULL), m_TargetLandmarks(NULL)
{
}
mitk::LandmarkProjectorBasedCurvedGeometry::LandmarkProjectorBasedCurvedGeometry(const mitk::LandmarkProjectorBasedCurvedGeometry& other) : Superclass(other)
{
SetTargetLandmarks(other.m_TargetLandmarks);
this->SetLandmarkProjector(other.m_LandmarkProjector);
this->ComputeGeometry();
}
mitk::LandmarkProjectorBasedCurvedGeometry::~LandmarkProjectorBasedCurvedGeometry()
{
if(m_InterpolatingAbstractTransform!=NULL)
m_InterpolatingAbstractTransform->Delete();
}
void mitk::LandmarkProjectorBasedCurvedGeometry::SetLandmarkProjector(mitk::LandmarkProjector* aLandmarkProjector)
{
itkDebugMacro("setting LandmarkProjector to " << aLandmarkProjector );
if(m_LandmarkProjector != aLandmarkProjector)
{
m_LandmarkProjector = aLandmarkProjector;
if(m_LandmarkProjector.IsNotNull())
{
if(m_FrameGeometry.IsNotNull())
m_LandmarkProjector->SetFrameGeometry(m_FrameGeometry);
if(m_InterpolatingAbstractTransform == NULL)
{
itkWarningMacro(<<"m_InterpolatingAbstractTransform not set.");
}
m_LandmarkProjector->SetInterpolatingAbstractTransform(GetInterpolatingAbstractTransform());
SetVtkAbstractTransform(m_LandmarkProjector->GetCompleteAbstractTransform());
}
Modified();
}
}
void mitk::LandmarkProjectorBasedCurvedGeometry::SetFrameGeometry(const mitk::BaseGeometry* frameGeometry)
{
Superclass::SetFrameGeometry(frameGeometry);
if(m_LandmarkProjector.IsNotNull())
m_LandmarkProjector->SetFrameGeometry(frameGeometry);
}
void mitk::LandmarkProjectorBasedCurvedGeometry::ComputeGeometry()
{
if(m_LandmarkProjector.IsNull())
{
itkExceptionMacro(<< "m_LandmarkProjector is not set.");
}
m_LandmarkProjector->ProjectLandmarks(m_TargetLandmarks);
SetPlane(m_LandmarkProjector->GetParameterPlane());
}
itk::LightObject::Pointer mitk::LandmarkProjectorBasedCurvedGeometry::InternalClone() const
{
mitk::BaseGeometry::Pointer newGeometry = new LandmarkProjectorBasedCurvedGeometry(*this);
newGeometry->UnRegister();
return newGeometry.GetPointer();
}
diff --git a/Core/Code/DataManagement/mitkPlaneGeometry.cpp b/Core/Code/DataManagement/mitkPlaneGeometry.cpp
index 8ad0309caf..eca60654bd 100644
--- a/Core/Code/DataManagement/mitkPlaneGeometry.cpp
+++ b/Core/Code/DataManagement/mitkPlaneGeometry.cpp
@@ -1,916 +1,916 @@
/*===================================================================
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 "mitkPlaneGeometry.h"
#include "mitkPlaneOperation.h"
#include "mitkInteractionConst.h"
#include "mitkLine.h"
#include <vtkTransform.h>
#include <vnl/vnl_cross.h>
namespace mitk
{
PlaneGeometry::PlaneGeometry()
- : m_ScaleFactorMMPerUnitX( 1.0 ),
+ : Superclass(), m_ScaleFactorMMPerUnitX( 1.0 ),
m_ScaleFactorMMPerUnitY( 1.0 ),
m_ReferenceGeometry( NULL )
{
Initialize();
}
PlaneGeometry::~PlaneGeometry()
{
}
PlaneGeometry::PlaneGeometry(const PlaneGeometry& other)
: Superclass(other), m_ScaleFactorMMPerUnitX( other.m_ScaleFactorMMPerUnitX),
m_ScaleFactorMMPerUnitY( other.m_ScaleFactorMMPerUnitY),
m_ReferenceGeometry( other.m_ReferenceGeometry )
{
}
void
PlaneGeometry::EnsurePerpendicularNormal(mitk::AffineTransform3D *transform)
{
//ensure row(2) of transform to be perpendicular to plane, keep length.
VnlVector normal = vnl_cross_3d(
transform->GetMatrix().GetVnlMatrix().get_column(0),
transform->GetMatrix().GetVnlMatrix().get_column(1) );
normal.normalize();
ScalarType len = transform->GetMatrix()
.GetVnlMatrix().get_column(2).two_norm();
if (len==0) len = 1;
normal*=len;
Matrix3D matrix = transform->GetMatrix();
matrix.GetVnlMatrix().set_column(2, normal);
transform->SetMatrix(matrix);
}
void
PlaneGeometry::PreSetIndexToWorldTransform(mitk::AffineTransform3D *transform)
{
EnsurePerpendicularNormal(transform);
}
void
PlaneGeometry::PreSetBounds(const BoundingBox::BoundsArrayType &bounds)
{
//currently the unit rectangle must be starting at the origin [0,0]
assert(bounds[0]==0);
assert(bounds[2]==0);
//the unit rectangle must be two-dimensional
assert(bounds[1]>0);
assert(bounds[3]>0);
}
void
PlaneGeometry::IndexToWorld( const Point2D &pt_units, Point2D &pt_mm ) const
{
pt_mm[0]=m_ScaleFactorMMPerUnitX*pt_units[0];
pt_mm[1]=m_ScaleFactorMMPerUnitY*pt_units[1];
}
void PlaneGeometry::WorldToIndex( const Point2D &pt_mm, Point2D &pt_units ) const
{
pt_units[0]=pt_mm[0]*(1.0/m_ScaleFactorMMPerUnitX);
pt_units[1]=pt_mm[1]*(1.0/m_ScaleFactorMMPerUnitY);
}
void PlaneGeometry::IndexToWorld( const Point2D & /*atPt2d_units*/,
const Vector2D &vec_units, Vector2D &vec_mm) const
{
MITK_WARN<<"Warning! Call of the deprecated function PlaneGeometry::IndexToWorld(point, vec, vec). Use PlaneGeometry::IndexToWorld(vec, vec) instead!";
this->IndexToWorld(vec_units, vec_mm);
}
void PlaneGeometry::IndexToWorld(const Vector2D &vec_units, Vector2D &vec_mm) const
{
vec_mm[0] = m_ScaleFactorMMPerUnitX * vec_units[0];
vec_mm[1] = m_ScaleFactorMMPerUnitY * vec_units[1];
}
void
PlaneGeometry::WorldToIndex( const Point2D & /*atPt2d_mm*/,
const Vector2D &vec_mm, Vector2D &vec_units) const
{
MITK_WARN<<"Warning! Call of the deprecated function PlaneGeometry::WorldToIndex(point, vec, vec). Use PlaneGeometry::WorldToIndex(vec, vec) instead!";
this->WorldToIndex(vec_mm, vec_units);
}
void
PlaneGeometry::WorldToIndex( const Vector2D &vec_mm, Vector2D &vec_units) const
{
vec_units[0] = vec_mm[0] * ( 1.0 / m_ScaleFactorMMPerUnitX );
vec_units[1] = vec_mm[1] * ( 1.0 / m_ScaleFactorMMPerUnitY );
}
void
PlaneGeometry::InitializeStandardPlane( mitk::ScalarType width,
ScalarType height, const Vector3D & spacing,
PlaneGeometry::PlaneOrientation planeorientation,
ScalarType zPosition, bool frontside, bool rotated )
{
AffineTransform3D::Pointer transform;
transform = AffineTransform3D::New();
AffineTransform3D::MatrixType matrix;
AffineTransform3D::MatrixType::InternalMatrixType &vnlmatrix =
matrix.GetVnlMatrix();
vnlmatrix.set_identity();
vnlmatrix(0,0) = spacing[0];
vnlmatrix(1,1) = spacing[1];
vnlmatrix(2,2) = spacing[2];
transform->SetIdentity();
transform->SetMatrix(matrix);
InitializeStandardPlane(width, height, transform.GetPointer(),
planeorientation, zPosition, frontside, rotated);
}
void
PlaneGeometry::InitializeStandardPlane( mitk::ScalarType width,
ScalarType height, const AffineTransform3D* transform,
PlaneGeometry::PlaneOrientation planeorientation, ScalarType zPosition,
bool frontside, bool rotated )
{
Superclass::Initialize();
//construct standard view
Point3D origin;
VnlVector rightDV(3), bottomDV(3);
origin.Fill(0);
int normalDirection;
switch(planeorientation)
{
case Axial:
if(frontside)
{
if(rotated==false)
{
FillVector3D(origin, 0, 0, zPosition);
FillVector3D(rightDV, 1, 0, 0);
FillVector3D(bottomDV, 0, 1, 0);
}
else
{
FillVector3D(origin, width, height, zPosition);
FillVector3D(rightDV, -1, 0, 0);
FillVector3D(bottomDV, 0, -1, 0);
}
}
else
{
if(rotated==false)
{
FillVector3D(origin, width, 0, zPosition);
FillVector3D(rightDV, -1, 0, 0);
FillVector3D(bottomDV, 0, 1, 0);
}
else
{
FillVector3D(origin, 0, height, zPosition);
FillVector3D(rightDV, 1, 0, 0);
FillVector3D(bottomDV, 0, -1, 0);
}
}
normalDirection = 2;
break;
case Frontal:
if(frontside)
{
if(rotated==false)
{
FillVector3D(origin, 0, zPosition, 0);
FillVector3D(rightDV, 1, 0, 0);
FillVector3D(bottomDV, 0, 0, 1);
}
else
{
FillVector3D(origin, width, zPosition, height);
FillVector3D(rightDV, -1, 0, 0);
FillVector3D(bottomDV, 0, 0, -1);
}
}
else
{
if(rotated==false)
{
FillVector3D(origin, width, zPosition, 0);
FillVector3D(rightDV, -1, 0, 0);
FillVector3D(bottomDV, 0, 0, 1);
}
else
{
FillVector3D(origin, 0, zPosition, height);
FillVector3D(rightDV, 1, 0, 0);
FillVector3D(bottomDV, 0, 0, -1);
}
}
normalDirection = 1;
break;
case Sagittal:
if(frontside)
{
if(rotated==false)
{
FillVector3D(origin, zPosition, 0, 0);
FillVector3D(rightDV, 0, 1, 0);
FillVector3D(bottomDV, 0, 0, 1);
}
else
{
FillVector3D(origin, zPosition, width, height);
FillVector3D(rightDV, 0, -1, 0);
FillVector3D(bottomDV, 0, 0, -1);
}
}
else
{
if(rotated==false)
{
FillVector3D(origin, zPosition, width, 0);
FillVector3D(rightDV, 0, -1, 0);
FillVector3D(bottomDV, 0, 0, 1);
}
else
{
FillVector3D(origin, zPosition, 0, height);
FillVector3D(rightDV, 0, 1, 0);
FillVector3D(bottomDV, 0, 0, -1);
}
}
normalDirection = 0;
break;
default:
itkExceptionMacro("unknown PlaneOrientation");
}
if ( transform != NULL )
{
origin = transform->TransformPoint( origin );
rightDV = transform->TransformVector( rightDV );
bottomDV = transform->TransformVector( bottomDV );
}
ScalarType bounds[6]= { 0, width, 0, height, 0, 1 };
this->SetBounds( bounds );
if ( transform == NULL )
{
this->SetMatrixByVectors( rightDV, bottomDV );
}
else
{
this->SetMatrixByVectors(
rightDV, bottomDV,
transform->GetMatrix().GetVnlMatrix()
.get_column(normalDirection).magnitude()
);
}
this->SetOrigin(origin);
}
void
PlaneGeometry::InitializeStandardPlane( const BaseGeometry *geometry3D,
PlaneOrientation planeorientation, ScalarType zPosition,
bool frontside, bool rotated )
{
this->SetReferenceGeometry( const_cast< BaseGeometry * >( geometry3D ) );
ScalarType width, height;
const BoundingBox::BoundsArrayType& boundsarray =
geometry3D->GetBoundingBox()->GetBounds();
Vector3D originVector;
FillVector3D(originVector, boundsarray[0], boundsarray[2], boundsarray[4]);
if(geometry3D->GetImageGeometry())
{
FillVector3D( originVector,
originVector[0] - 0.5,
originVector[1] - 0.5,
originVector[2] - 0.5 );
}
switch(planeorientation)
{
case Axial:
width = geometry3D->GetExtent(0);
height = geometry3D->GetExtent(1);
break;
case Frontal:
width = geometry3D->GetExtent(0);
height = geometry3D->GetExtent(2);
break;
case Sagittal:
width = geometry3D->GetExtent(1);
height = geometry3D->GetExtent(2);
break;
default:
itkExceptionMacro("unknown PlaneOrientation");
}
InitializeStandardPlane( width, height,
geometry3D->GetIndexToWorldTransform(),
planeorientation, zPosition, frontside, rotated );
ScalarType bounds[6]= { 0, width, 0, height, 0, 1 };
this->SetBounds( bounds );
Point3D origin;
originVector = geometry3D->GetIndexToWorldTransform()
->TransformVector( originVector );
origin = GetOrigin() + originVector;
SetOrigin(origin);
}
void
PlaneGeometry::InitializeStandardPlane( const BaseGeometry *geometry3D,
bool top, PlaneOrientation planeorientation, bool frontside, bool rotated )
{
ScalarType zPosition;
switch(planeorientation)
{
case Axial:
zPosition = (top ? 0.5 : geometry3D->GetExtent(2)-1+0.5);
break;
case Frontal:
zPosition = (top ? 0.5 : geometry3D->GetExtent(1)-1+0.5);
break;
case Sagittal:
zPosition = (top ? 0.5 : geometry3D->GetExtent(0)-1+0.5);
break;
default:
itkExceptionMacro("unknown PlaneOrientation");
}
InitializeStandardPlane( geometry3D, planeorientation,
zPosition, frontside, rotated );
}
void
PlaneGeometry::InitializeStandardPlane( const Vector3D &rightVector,
const Vector3D &downVector, const Vector3D *spacing )
{
InitializeStandardPlane( rightVector.GetVnlVector(),
downVector.GetVnlVector(), spacing );
}
void
PlaneGeometry::InitializeStandardPlane( const VnlVector& rightVector,
const VnlVector &downVector, const Vector3D *spacing )
{
ScalarType width = rightVector.magnitude();
ScalarType height = downVector.magnitude();
InitializeStandardPlane( width, height, rightVector, downVector, spacing );
}
void
PlaneGeometry::InitializeStandardPlane( mitk::ScalarType width,
ScalarType height, const Vector3D &rightVector, const Vector3D &downVector,
const Vector3D *spacing )
{
InitializeStandardPlane(
width, height,
rightVector.GetVnlVector(), downVector.GetVnlVector(),
spacing );
}
void
PlaneGeometry::InitializeStandardPlane(
mitk::ScalarType width, ScalarType height,
const VnlVector &rightVector, const VnlVector &downVector,
const Vector3D *spacing )
{
assert(width > 0);
assert(height > 0);
VnlVector rightDV = rightVector; rightDV.normalize();
VnlVector downDV = downVector; downDV.normalize();
VnlVector normal = vnl_cross_3d(rightVector, downVector);
normal.normalize();
if(spacing!=NULL)
{
rightDV *= (*spacing)[0];
downDV *= (*spacing)[1];
normal *= (*spacing)[2];
}
AffineTransform3D::Pointer transform = AffineTransform3D::New();
Matrix3D matrix;
matrix.GetVnlMatrix().set_column(0, rightDV);
matrix.GetVnlMatrix().set_column(1, downDV);
matrix.GetVnlMatrix().set_column(2, normal);
transform->SetMatrix(matrix);
transform->SetOffset(this->GetIndexToWorldTransform()->GetOffset());
ScalarType bounds[6] = { 0, width, 0, height, 0, 1 };
this->SetBounds( bounds );
this->SetIndexToWorldTransform( transform );
}
void
PlaneGeometry::InitializePlane( const Point3D &origin,
const Vector3D &normal )
{
VnlVector rightVectorVnl(3), downVectorVnl;
if( Equal( normal[1], 0.0f ) == false )
{
FillVector3D( rightVectorVnl, 1.0f, -normal[0]/normal[1], 0.0f );
rightVectorVnl.normalize();
}
else
{
FillVector3D( rightVectorVnl, 0.0f, 1.0f, 0.0f );
}
downVectorVnl = vnl_cross_3d( normal.GetVnlVector(), rightVectorVnl );
downVectorVnl.normalize();
InitializeStandardPlane( rightVectorVnl, downVectorVnl );
SetOrigin(origin);
}
void
PlaneGeometry::SetMatrixByVectors( const VnlVector &rightVector,
const VnlVector &downVector, ScalarType thickness )
{
VnlVector normal = vnl_cross_3d(rightVector, downVector);
normal.normalize();
normal *= thickness;
AffineTransform3D::Pointer transform = AffineTransform3D::New();
Matrix3D matrix;
matrix.GetVnlMatrix().set_column(0, rightVector);
matrix.GetVnlMatrix().set_column(1, downVector);
matrix.GetVnlMatrix().set_column(2, normal);
transform->SetMatrix(matrix);
transform->SetOffset(this->GetIndexToWorldTransform()->GetOffset());
SetIndexToWorldTransform(transform);
}
Vector3D
PlaneGeometry::GetNormal() const
{
Vector3D frontToBack;
frontToBack.SetVnlVector( this->GetIndexToWorldTransform()
->GetMatrix().GetVnlMatrix().get_column(2) );
return frontToBack;
}
VnlVector
PlaneGeometry::GetNormalVnl() const
{
return this->GetIndexToWorldTransform()
->GetMatrix().GetVnlMatrix().get_column(2);
}
ScalarType
PlaneGeometry::DistanceFromPlane( const Point3D &pt3d_mm ) const
{
return fabs(SignedDistance( pt3d_mm ));
}
ScalarType
PlaneGeometry::SignedDistance( const Point3D &pt3d_mm ) const
{
return SignedDistanceFromPlane(pt3d_mm);
}
//Function from Geometry2D
// mitk::ScalarType
// PlaneGeometry::SignedDistance(const mitk::Point3D& pt3d_mm) const
//{
// Point3D projectedPoint;
// Project(pt3d_mm, projectedPoint);
// Vector3D direction = pt3d_mm-projectedPoint;
// ScalarType distance = direction.GetNorm();
// if(IsAbove(pt3d_mm) == false)
// distance*=-1.0;
// return distance;
//}
bool
PlaneGeometry::IsAbove( const Point3D &pt3d_mm , bool considerBoundingBox) const
{
if(considerBoundingBox)
{
Point3D pt3d_units;
BaseGeometry::WorldToIndex(pt3d_mm, pt3d_units);
return (pt3d_units[2] > this->GetBoundingBox()->GetBounds()[4]);
}
else
return SignedDistanceFromPlane(pt3d_mm) > 0;
}
bool
PlaneGeometry::IntersectionLine(
const PlaneGeometry* plane, Line3D& crossline ) const
{
Vector3D normal = this->GetNormal();
normal.Normalize();
Vector3D planeNormal = plane->GetNormal();
planeNormal.Normalize();
Vector3D direction = itk::CrossProduct( normal, planeNormal );
if ( direction.GetSquaredNorm() < eps )
return false;
crossline.SetDirection( direction );
double N1dN2 = normal * planeNormal;
double determinant = 1.0 - N1dN2 * N1dN2;
Vector3D origin = this->GetOrigin().GetVectorFromOrigin();
Vector3D planeOrigin = plane->GetOrigin().GetVectorFromOrigin();
double d1 = normal * origin;
double d2 = planeNormal * planeOrigin;
double c1 = ( d1 - d2 * N1dN2 ) / determinant;
double c2 = ( d2 - d1 * N1dN2 ) / determinant;
Vector3D p = normal * c1 + planeNormal * c2;
crossline.GetPoint().GetVnlVector() = p.GetVnlVector();
return true;
}
unsigned int
PlaneGeometry::IntersectWithPlane2D(
const PlaneGeometry* plane, Point2D& lineFrom, Point2D &lineTo ) const
{
Line3D crossline;
if ( this->IntersectionLine( plane, crossline ) == false )
return 0;
Point2D point2;
Vector2D direction2;
this->Map( crossline.GetPoint(), point2 );
this->Map( crossline.GetPoint(), crossline.GetDirection(), direction2 );
return
Line3D::RectangleLineIntersection(
0, 0, GetExtentInMM(0), GetExtentInMM(1),
point2, direction2, lineFrom, lineTo );
}
double PlaneGeometry::Angle( const PlaneGeometry *plane ) const
{
return angle(plane->GetMatrixColumn(2), GetMatrixColumn(2));
}
double PlaneGeometry::Angle( const Line3D &line ) const
{
return vnl_math::pi_over_2
- angle( line.GetDirection().GetVnlVector(), GetMatrixColumn(2) );
}
bool PlaneGeometry::IntersectionPoint(
const Line3D &line, Point3D &intersectionPoint ) const
{
Vector3D planeNormal = this->GetNormal();
planeNormal.Normalize();
Vector3D lineDirection = line.GetDirection();
lineDirection.Normalize();
double t = planeNormal * lineDirection;
if ( fabs( t ) < eps )
{
return false;
}
Vector3D diff;
diff = this->GetOrigin() - line.GetPoint();
t = ( planeNormal * diff ) / t;
intersectionPoint = line.GetPoint() + lineDirection * t;
return true;
}
bool
PlaneGeometry::IntersectionPointParam( const Line3D &line, double &t ) const
{
Vector3D planeNormal = this->GetNormal();
Vector3D lineDirection = line.GetDirection();
t = planeNormal * lineDirection;
if ( fabs( t ) < eps )
{
return false;
}
Vector3D diff;
diff = this->GetOrigin() - line.GetPoint();
t = ( planeNormal * diff ) / t;
return true;
}
bool
PlaneGeometry::IsParallel( const PlaneGeometry *plane ) const
{
return ( (Angle(plane) < 10.0 * mitk::sqrteps ) || ( Angle(plane) > ( vnl_math::pi - 10.0 * sqrteps ) ) ) ;
}
bool
PlaneGeometry::IsOnPlane( const Point3D &point ) const
{
return Distance(point) < eps;
}
bool
PlaneGeometry::IsOnPlane( const Line3D &line ) const
{
return ( (Distance( line.GetPoint() ) < eps)
&& (Distance( line.GetPoint2() ) < eps) );
}
bool
PlaneGeometry::IsOnPlane( const PlaneGeometry *plane ) const
{
return ( IsParallel( plane ) && (Distance( plane->GetOrigin() ) < eps) );
}
Point3D
PlaneGeometry::ProjectPointOntoPlane( const Point3D& pt ) const
{
ScalarType len = this->GetNormalVnl().two_norm();
return pt - this->GetNormal() * this->SignedDistanceFromPlane( pt ) / len;
}
itk::LightObject::Pointer
PlaneGeometry::InternalClone() const
{
Self::Pointer newGeometry = new PlaneGeometry(*this);
newGeometry->UnRegister();
return newGeometry.GetPointer();
}
void
PlaneGeometry::ExecuteOperation( Operation *operation )
{
vtkTransform *transform = vtkTransform::New();
transform->SetMatrix( this->GetVtkMatrix());
switch ( operation->GetOperationType() )
{
case OpORIENT:
{
mitk::PlaneOperation *planeOp = dynamic_cast< mitk::PlaneOperation * >( operation );
if ( planeOp == NULL )
{
return;
}
Point3D center = planeOp->GetPoint();
Vector3D orientationVector = planeOp->GetNormal();
Vector3D defaultVector;
FillVector3D( defaultVector, 0.0, 0.0, 1.0 );
Vector3D rotationAxis = itk::CrossProduct( orientationVector, defaultVector );
//double rotationAngle = acos( orientationVector[2] / orientationVector.GetNorm() );
double rotationAngle = atan2( (double) rotationAxis.GetNorm(), (double) (orientationVector * defaultVector) );
rotationAngle *= 180.0 / vnl_math::pi;
transform->PostMultiply();
transform->Identity();
transform->Translate( center[0], center[1], center[2] );
transform->RotateWXYZ( rotationAngle, rotationAxis[0], rotationAxis[1], rotationAxis[2] );
transform->Translate( -center[0], -center[1], -center[2] );
break;
}
case OpRESTOREPLANEPOSITION:
{
RestorePlanePositionOperation *op = dynamic_cast< mitk::RestorePlanePositionOperation* >(operation);
if(op == NULL)
{
return;
}
AffineTransform3D::Pointer transform2 = AffineTransform3D::New();
Matrix3D matrix;
matrix.GetVnlMatrix().set_column(0, op->GetTransform()->GetMatrix().GetVnlMatrix().get_column(0));
matrix.GetVnlMatrix().set_column(1, op->GetTransform()->GetMatrix().GetVnlMatrix().get_column(1));
matrix.GetVnlMatrix().set_column(2, op->GetTransform()->GetMatrix().GetVnlMatrix().get_column(2));
transform2->SetMatrix(matrix);
Vector3D offset = op->GetTransform()->GetOffset();
transform2->SetOffset(offset);
this->SetIndexToWorldTransform(transform2);
ScalarType bounds[6] = {0, op->GetWidth(), 0, op->GetHeight(), 0 ,1 };
this->SetBounds(bounds);
TransferItkToVtkTransform();
this->Modified();
transform->Delete();
return;
}
default:
Superclass::ExecuteOperation( operation );
transform->Delete();
return;
}
this->GetVtkMatrix()->DeepCopy(transform->GetMatrix());
this->TransferVtkToItkTransform();
this->Modified();
transform->Delete();
}
void PlaneGeometry::PrintSelf( std::ostream& os, itk::Indent indent ) const
{
Superclass::PrintSelf(os,indent);
os << indent << " ScaleFactorMMPerUnitX: "
<< m_ScaleFactorMMPerUnitX << std::endl;
os << indent << " ScaleFactorMMPerUnitY: "
<< m_ScaleFactorMMPerUnitY << std::endl;
os << indent << " Normal: " << GetNormal() << std::endl;
}
void PlaneGeometry::PostSetIndexToWorldTransform(
mitk::AffineTransform3D* transform)
{
m_ScaleFactorMMPerUnitX=GetExtentInMM(0)/GetExtent(0);
m_ScaleFactorMMPerUnitY=GetExtentInMM(1)/GetExtent(1);
assert(m_ScaleFactorMMPerUnitX<itk::NumericTraits<ScalarType>::infinity());
assert(m_ScaleFactorMMPerUnitY<itk::NumericTraits<ScalarType>::infinity());
}
void
PlaneGeometry::PostSetExtentInMM(int direction, ScalarType extentInMM)
{
m_ScaleFactorMMPerUnitX=GetExtentInMM(0)/GetExtent(0);
m_ScaleFactorMMPerUnitY=GetExtentInMM(1)/GetExtent(1);
assert(m_ScaleFactorMMPerUnitX<itk::NumericTraits<ScalarType>::infinity());
assert(m_ScaleFactorMMPerUnitY<itk::NumericTraits<ScalarType>::infinity());
}
bool
PlaneGeometry::Map(
const mitk::Point3D &pt3d_mm, mitk::Point2D &pt2d_mm) const
{
assert(this->IsBoundingBoxNull()==false);
Point3D pt3d_units;
BackTransform(pt3d_mm, pt3d_units);
pt2d_mm[0]=pt3d_units[0]*m_ScaleFactorMMPerUnitX;
pt2d_mm[1]=pt3d_units[1]*m_ScaleFactorMMPerUnitY;
pt3d_units[2]=0;
return const_cast<BoundingBox*>(this->GetBoundingBox())->IsInside(pt3d_units);
}
void
PlaneGeometry::Map(const mitk::Point2D &pt2d_mm, mitk::Point3D &pt3d_mm) const
{
Point3D pt3d_units;
pt3d_units[0]=pt2d_mm[0]/m_ScaleFactorMMPerUnitX;
pt3d_units[1]=pt2d_mm[1]/m_ScaleFactorMMPerUnitY;
pt3d_units[2]=0;
pt3d_mm = GetIndexToWorldTransform()->TransformPoint(pt3d_units);
}
void
PlaneGeometry::SetSizeInUnits(mitk::ScalarType width, mitk::ScalarType height)
{
ScalarType bounds[6]={0, width, 0, height, 0, 1};
ScalarType extent, newextentInMM;
if(GetExtent(0)>0)
{
extent = GetExtent(0);
if(width>extent)
newextentInMM = GetExtentInMM(0)/width*extent;
else
newextentInMM = GetExtentInMM(0)*extent/width;
SetExtentInMM(0, newextentInMM);
}
if(GetExtent(1)>0)
{
extent = GetExtent(1);
if(width>extent)
newextentInMM = GetExtentInMM(1)/height*extent;
else
newextentInMM = GetExtentInMM(1)*extent/height;
SetExtentInMM(1, newextentInMM);
}
SetBounds(bounds);
}
bool
PlaneGeometry::Project(
const mitk::Point3D &pt3d_mm, mitk::Point3D &projectedPt3d_mm) const
{
assert(this->IsBoundingBoxNull()==false);
Point3D pt3d_units;
BackTransform(pt3d_mm, pt3d_units);
pt3d_units[2] = 0;
projectedPt3d_mm = GetIndexToWorldTransform()->TransformPoint(pt3d_units);
return const_cast<BoundingBox*>(this->GetBoundingBox())->IsInside(pt3d_units);
}
bool
PlaneGeometry::Project(const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm) const
{
assert(this->IsBoundingBoxNull()==false);
Vector3D vec3d_units;
BackTransform(vec3d_mm, vec3d_units);
vec3d_units[2] = 0;
projectedVec3d_mm = GetIndexToWorldTransform()->TransformVector(vec3d_units);
return true;
}
bool
PlaneGeometry::Project(const mitk::Point3D & atPt3d_mm,
const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm) const
{
MITK_WARN << "Deprecated function! Call Project(vec3D,vec3D) instead.";
assert(this->IsBoundingBoxNull()==false);
Vector3D vec3d_units;
BackTransform(atPt3d_mm, vec3d_mm, vec3d_units);
vec3d_units[2] = 0;
projectedVec3d_mm = GetIndexToWorldTransform()->TransformVector(vec3d_units);
Point3D pt3d_units;
BackTransform(atPt3d_mm, pt3d_units);
return const_cast<BoundingBox*>(this->GetBoundingBox())->IsInside(pt3d_units);
}
bool
PlaneGeometry::Map(const mitk::Point3D & atPt3d_mm,
const mitk::Vector3D &vec3d_mm, mitk::Vector2D &vec2d_mm) const
{
Point2D pt2d_mm_start, pt2d_mm_end;
Point3D pt3d_mm_end;
bool inside=Map(atPt3d_mm, pt2d_mm_start);
pt3d_mm_end = atPt3d_mm+vec3d_mm;
inside&=Map(pt3d_mm_end, pt2d_mm_end);
vec2d_mm=pt2d_mm_end-pt2d_mm_start;
return inside;
}
void
PlaneGeometry::Map(const mitk::Point2D &/*atPt2d_mm*/,
const mitk::Vector2D &/*vec2d_mm*/, mitk::Vector3D &/*vec3d_mm*/) const
{
//@todo implement parallel to the other Map method!
assert(false);
}
void
PlaneGeometry::SetReferenceGeometry( mitk::BaseGeometry *geometry )
{
m_ReferenceGeometry = geometry;
}
mitk::BaseGeometry *
PlaneGeometry::GetReferenceGeometry() const
{
return m_ReferenceGeometry;
}
bool
PlaneGeometry::HasReferenceGeometry() const
{
return ( m_ReferenceGeometry != NULL );
}
} // namespace
diff --git a/Core/Code/DataManagement/mitkThinPlateSplineCurvedGeometry.cpp b/Core/Code/DataManagement/mitkThinPlateSplineCurvedGeometry.cpp
index 562f2d6745..33ddbac8d7 100644
--- a/Core/Code/DataManagement/mitkThinPlateSplineCurvedGeometry.cpp
+++ b/Core/Code/DataManagement/mitkThinPlateSplineCurvedGeometry.cpp
@@ -1,100 +1,101 @@
/*===================================================================
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 "mitkThinPlateSplineCurvedGeometry.h"
#include <vtkThinPlateSplineTransform.h>
#include <vtkPoints.h>
mitk::ThinPlateSplineCurvedGeometry::ThinPlateSplineCurvedGeometry()
+ : Superclass()
{
m_InterpolatingAbstractTransform = m_ThinPlateSplineTransform = vtkThinPlateSplineTransform::New();
m_VtkTargetLandmarks = vtkPoints::New();
m_VtkProjectedLandmarks = vtkPoints::New();
m_ThinPlateSplineTransform->SetInverseIterations(5000);
}
mitk::ThinPlateSplineCurvedGeometry::ThinPlateSplineCurvedGeometry(const ThinPlateSplineCurvedGeometry& other ) : Superclass(other)
{
this->SetSigma(other.GetSigma());
}
mitk::ThinPlateSplineCurvedGeometry::~ThinPlateSplineCurvedGeometry()
{
// don't need to delete m_ThinPlateSplineTransform, because it is
// the same as m_InterpolatingAbstractTransform, which will be deleted
// by the superclass.
if(m_VtkTargetLandmarks!=NULL)
m_VtkTargetLandmarks->Delete();
if(m_VtkProjectedLandmarks!=NULL)
m_VtkProjectedLandmarks->Delete();
}
bool mitk::ThinPlateSplineCurvedGeometry::IsValid() const
{
return m_TargetLandmarks.IsNotNull() && (m_TargetLandmarks->Size() >= 3) && m_LandmarkProjector.IsNotNull();
}
void mitk::ThinPlateSplineCurvedGeometry::SetSigma(double sigma)
{
m_ThinPlateSplineTransform->SetSigma(sigma);
}
double mitk::ThinPlateSplineCurvedGeometry::GetSigma() const
{
return m_ThinPlateSplineTransform->GetSigma();
}
void mitk::ThinPlateSplineCurvedGeometry::ComputeGeometry()
{
Superclass::ComputeGeometry();
const mitk::PointSet::DataType::PointsContainer *finalTargetLandmarks, *projectedTargetLandmarks;
finalTargetLandmarks = m_LandmarkProjector->GetFinalTargetLandmarks();
projectedTargetLandmarks = m_LandmarkProjector->GetProjectedLandmarks();
mitk::PointSet::DataType::PointsContainer::ConstIterator targetIt, projectedIt;
targetIt = finalTargetLandmarks->Begin();
projectedIt = projectedTargetLandmarks->Begin();
//initialize Thin-Plate-Spline
m_VtkTargetLandmarks->Reset();
m_VtkProjectedLandmarks->Reset();
vtkIdType id;
int size=finalTargetLandmarks->Size();
for(id=0; id < size; ++id, ++targetIt, ++projectedIt)
{
const mitk::PointSet::PointType& target = targetIt->Value();
m_VtkTargetLandmarks->InsertPoint(id, target[0], target[1], target[2]);
const mitk::PointSet::PointType& projected = projectedIt->Value();
m_VtkProjectedLandmarks->InsertPoint(id, projected[0], projected[1], projected[2]);
}
m_VtkTargetLandmarks->Modified();
m_VtkProjectedLandmarks->Modified();
m_ThinPlateSplineTransform->SetSourceLandmarks(m_VtkProjectedLandmarks);
m_ThinPlateSplineTransform->SetTargetLandmarks(m_VtkTargetLandmarks);
}
itk::LightObject::Pointer mitk::ThinPlateSplineCurvedGeometry::InternalClone() const
{
mitk::BaseGeometry::Pointer newGeometry = new Self(*this);
newGeometry->UnRegister();
return newGeometry.GetPointer();
}
diff --git a/Core/Code/Testing/mitkImageTest.cpp b/Core/Code/Testing/mitkImageTest.cpp
index 978fee8647..545ceb6db8 100644
--- a/Core/Code/Testing/mitkImageTest.cpp
+++ b/Core/Code/Testing/mitkImageTest.cpp
@@ -1,526 +1,517 @@
/*===================================================================
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.
===================================================================*/
// mitk includes
#include <mitkImage.h>
#include <mitkImageDataItem.h>
#include <mitkImageCast.h>
#include "mitkItkImageFileReader.h"
#include <mitkTestingMacros.h>
#include <mitkImageStatisticsHolder.h>
#include "mitkImageGenerator.h"
#include "mitkImageReadAccessor.h"
#include "mitkException.h"
#include "mitkPixelTypeMultiplex.h"
#include "mitkImagePixelReadAccessor.h"
#include "mitkImageSliceSelector.h"
// itk includes
#include <itkImage.h>
#include <itkMersenneTwisterRandomVariateGenerator.h>
// stl includes
#include <fstream>
// vtk includes
#include <vtkImageData.h>
// Checks if reference count is correct after using GetVtkImageData()
bool ImageVtkDataReferenceCheck(const char* fname) {
const std::string filename = std::string(fname);
mitk::ItkImageFileReader::Pointer imageReader = mitk::ItkImageFileReader::New();
try
{
imageReader->SetFileName(filename);
imageReader->Update();
}
catch(...) {
MITK_TEST_FAILED_MSG(<< "Could not read file for testing: " << filename);
return false;
}
{
mitk::Image::Pointer image = imageReader->GetOutput();
vtkImageData* vtk = image->GetVtkImageData();
if(vtk == NULL)
return false;
}
return true;
}
template <class T>
void TestRandomPixelAccess( const mitk::PixelType ptype, mitk::Image::Pointer image, mitk::Point3D & point, mitk::ScalarType & value )
{
// generate a random point in world coordinates
mitk::Point3D xMax, yMax, zMax, xMaxIndex, yMaxIndex, zMaxIndex;
xMaxIndex.Fill(0.0f);
yMaxIndex.Fill(0.0f);
zMaxIndex.Fill(0.0f);
xMaxIndex[0] = image->GetLargestPossibleRegion().GetSize()[0];
yMaxIndex[1] = image->GetLargestPossibleRegion().GetSize()[1];
zMaxIndex[2] = image->GetLargestPossibleRegion().GetSize()[2];
image->GetGeometry()->IndexToWorld(xMaxIndex, xMax);
image->GetGeometry()->IndexToWorld(yMaxIndex, yMax);
image->GetGeometry()->IndexToWorld(zMaxIndex, zMax);
MITK_INFO << "Origin " << image->GetGeometry()->GetOrigin()[0] << " "<< image->GetGeometry()->GetOrigin()[1] << " "<< image->GetGeometry()->GetOrigin()[2] << "";
MITK_INFO << "MaxExtend " << xMax[0] << " "<< yMax[1] << " "<< zMax[2] << "";
itk::Statistics::MersenneTwisterRandomVariateGenerator::Pointer randomGenerator = itk::Statistics::MersenneTwisterRandomVariateGenerator::New();
randomGenerator->Initialize( std::rand() ); // initialize with random value, to get sensible random points for the image
point[0] = randomGenerator->GetUniformVariate( image->GetGeometry()->GetOrigin()[0], xMax[0]);
point[1] = randomGenerator->GetUniformVariate( image->GetGeometry()->GetOrigin()[1], yMax[1]);
point[2] = randomGenerator->GetUniformVariate( image->GetGeometry()->GetOrigin()[2], zMax[2]);
MITK_INFO << "RandomPoint " << point[0] << " "<< point[1] << " "<< point[2] << "";
// test values and max/min
mitk::ScalarType imageMin = image->GetStatistics()->GetScalarValueMin();
mitk::ScalarType imageMax = image->GetStatistics()->GetScalarValueMax();
// test accessing PixelValue with coordinate leading to a negative index
const mitk::Point3D geom_origin = image->GetGeometry()->GetOrigin();
const mitk::Point3D geom_center = image->GetGeometry()->GetCenter();
// shift position from origin outside of the image ( in the opposite direction to [center-origin] vector which points in the inside)
mitk::Point3D position = geom_origin + (geom_origin - geom_center);
MITK_INFO << "Testing access outside of the image";
unsigned int dim = image->GetDimension();
if(dim == 3 || dim == 4){
mitk::ImagePixelReadAccessor<T,3> imAccess3(image,image->GetVolumeData(0));
// Comparison ?>=0 not needed since all position[i] and timestep are unsigned int
// (position[0]>=0 && position[1] >=0 && position[2]>=0 && timestep>=0)
// bug-11978 : we still need to catch index with negative values
if ( point[0] < 0 ||
point[1] < 0 ||
point[2] < 0 )
{
MITK_WARN << "Given position ("<< point << ") is out of image range, returning 0." ;
}
else {
value = static_cast<mitk::ScalarType>(imAccess3.GetPixelByWorldCoordinates(point));
MITK_TEST_CONDITION( (value >= imageMin && value <= imageMax), "Value returned is between max/min");
}
itk::Index<3> itkIndex;
image->GetGeometry()->WorldToIndex(position, itkIndex);
MITK_TEST_FOR_EXCEPTION_BEGIN(mitk::Exception);
imAccess3.GetPixelByIndexSafe(itkIndex);
MITK_TEST_FOR_EXCEPTION_END(mitk::Exception);
}
MITK_INFO << imageMin << " "<< imageMax << " "<< value << "";
}
class mitkImageTestClass
{
public:
void SetClonedGeometry_None_ClonedEqualInput()
{
mitk::Image::Pointer image = mitk::ImageGenerator::GenerateRandomImage<float>(100, 100, 100, 1, 0.2, 0.3, 0.4);
//-----------------
// geometry information for image
mitk::Point3D origin;
mitk::Vector3D right, bottom;
mitk::Vector3D spacing;
mitk::FillVector3D(origin, 17.0, 19.92, 7.83);
mitk::FillVector3D(right, 1.0, 2.0, 3.0);
mitk::FillVector3D(bottom, 0.0, -3.0, 2.0);
mitk::FillVector3D(spacing, 0.78, 0.91, 2.23);
//InitializeStandardPlane(rightVector, downVector, spacing)
mitk::PlaneGeometry::Pointer planegeometry = mitk::PlaneGeometry::New();
planegeometry->InitializeStandardPlane(100, 100, right, bottom, &spacing);
planegeometry->SetOrigin(origin);
planegeometry->ChangeImageGeometryConsideringOriginOffset(true);
image->SetClonedGeometry(planegeometry);
mitk::BaseGeometry::Pointer imageGeometry = image->GetGeometry();
itk::ScalableAffineTransform<mitk::ScalarType,3>* frameNew = imageGeometry->GetIndexToWorldTransform();
itk::ScalableAffineTransform<mitk::ScalarType,3>* frameOld = planegeometry->GetIndexToWorldTransform();
- bool matrixEqual = true;
- for (int i = 0; i < 16; ++i)
- {
- double valueNew = *(frameNew->GetMatrix()[i]);
- double valueOld = *(frameOld->GetMatrix()[i]);
- //MITK_INFO << "Index: " << i << " Old: " << valueOld << " New: " << valueNew << " Difference:" << valueOld-valueNew<< std::endl;
- matrixEqual = matrixEqual && mitk::Equal(valueNew, valueOld, mitk::eps);
- }
+ bool matrixEqual = mitk::Equal(imageGeometry, planegeometry, mitk::eps, false);
- // Disabled because this test fails on the dashboard. Does not fail on my machine.
- // See Bug 6505
- // MITK_TEST_CONDITION(matrixEqual, "Matrix elements of cloned matrix equal original matrix");
+ MITK_TEST_CONDITION(matrixEqual, "Matrix elements of cloned matrix equal original matrix");
}
};
int mitkImageTest(int argc, char* argv[])
{
MITK_TEST_BEGIN(mitkImageTest);
mitkImageTestClass tester;
tester.SetClonedGeometry_None_ClonedEqualInput();
//Create Image out of nowhere
mitk::Image::Pointer imgMem = mitk::Image::New();
mitk::PixelType pt = mitk::MakeScalarPixelType<int>();
unsigned int dim[]={100,100,20};
MITK_TEST_CONDITION_REQUIRED( imgMem.IsNotNull(), "An image was created. ");
// Initialize image
imgMem->Initialize( pt, 3, dim);
MITK_TEST_CONDITION_REQUIRED( imgMem->IsInitialized(), "Image::IsInitialized() ?");
MITK_TEST_CONDITION_REQUIRED( imgMem->GetPixelType() == pt, "PixelType was set correctly.");
int *p = NULL;
int *p2 = NULL;
try
{
mitk::ImageReadAccessor imgMemAcc(imgMem);
p = (int*)imgMemAcc.GetData();
}
catch (mitk::Exception& e)
{
MITK_ERROR << e.what();
}
MITK_TEST_CONDITION( p != NULL, "GetData() returned not-NULL pointer.");
// filling image
const unsigned int size = dim[0]*dim[1]*dim[2];
for(unsigned int i=0; i<size; ++i, ++p)
*p = (signed int)i;
// Getting it again and compare with filled values:
try
{
mitk::ImageReadAccessor imgMemAcc(imgMem);
p2 = (int*)imgMemAcc.GetData();
}
catch (mitk::Exception &e)
{
MITK_ERROR << e.what();
}
MITK_TEST_CONDITION( p2 != NULL, "GetData() returned not-NULL pointer.");
bool isEqual = true;
for(unsigned int i=0; i<size; ++i, ++p2)
{
if(*p2 != (signed int) i )
{
isEqual = false;
}
}
MITK_TEST_CONDITION( isEqual, "The values previously set as data are correct [pixelwise comparison].");
// Testing GetSliceData() and compare with filled values:
try
{
mitk::ImageReadAccessor imgMemAcc(imgMem, imgMem->GetSliceData(dim[2]/2));
p2 = (int*)imgMemAcc.GetData();
}
catch (mitk::Exception& e)
{
MITK_ERROR << e.what();
}
MITK_TEST_CONDITION_REQUIRED( p2 != NULL, "Valid slice data returned");
unsigned int xy_size = dim[0]*dim[1];
unsigned int start_mid_slice = (dim[2]/2)*xy_size;
isEqual = true;
for(unsigned int i=0; i<xy_size; ++i, ++p2)
{
if(*p2!=(signed int)(i+start_mid_slice))
{
isEqual = false;
}
}
MITK_TEST_CONDITION( isEqual, "The SliceData are correct [pixelwise comparison]. ");
imgMem = mitk::Image::New();
// testing re-initialization of test image
mitk::PixelType pType = mitk::MakePixelType<int, int, 1>();
imgMem->Initialize( pType , 3, dim);
MITK_TEST_CONDITION_REQUIRED(imgMem->GetDimension()== 3, "Testing initialization parameter dimension!");
MITK_TEST_CONDITION_REQUIRED(imgMem->GetPixelType() == pType, "Testing initialization parameter pixeltype!");
MITK_TEST_CONDITION_REQUIRED(imgMem->GetDimension(0) == dim[0] &&
imgMem->GetDimension(1)== dim[1] && imgMem->GetDimension(2)== dim[2], "Testing initialization of dimensions!");
MITK_TEST_CONDITION( imgMem->IsInitialized(), "Image is initialized.");
// Setting volume again:
try
{
mitk::ImageReadAccessor imgMemAcc(imgMem);
imgMem->SetVolume(imgMemAcc.GetData());
}
catch (mitk::Exception& e)
{
MITK_ERROR << e.what();
}
//-----------------
// geometry information for image
mitk::Point3D origin;
mitk::Vector3D right, bottom;
mitk::Vector3D spacing;
mitk::FillVector3D(origin, 17.0, 19.92, 7.83);
mitk::FillVector3D(right, 1.0, 2.0, 3.0);
mitk::FillVector3D(bottom, 0.0, -3.0, 2.0);
mitk::FillVector3D(spacing, 0.78, 0.91, 2.23);
//InitializeStandardPlane(rightVector, downVector, spacing)
mitk::PlaneGeometry::Pointer planegeometry = mitk::PlaneGeometry::New();
planegeometry->InitializeStandardPlane(100, 100, right, bottom, &spacing);
planegeometry->SetOrigin(origin);
// Testing Initialize(const mitk::PixelType& type, const mitk::Geometry3D& geometry, unsigned int slices) with PlaneGeometry and GetData(): ";
imgMem->Initialize( mitk::MakePixelType<int, int, 1>(), *planegeometry);
MITK_TEST_CONDITION_REQUIRED( imgMem->GetGeometry()->GetOrigin() == static_cast<mitk::BaseGeometry*>(planegeometry)->GetOrigin(), "Testing correct setting of geometry via initialize!");
try
{
mitk::ImageReadAccessor imgMemAcc(imgMem);
p = (int*)imgMemAcc.GetData();
}
catch (mitk::Exception& e)
{
MITK_ERROR << e.what();
}
MITK_TEST_CONDITION_REQUIRED( p!=NULL, "GetData() returned valid pointer.");
// Testing Initialize(const mitk::PixelType& type, int sDim, const mitk::PlaneGeometry& geometry) and GetData(): ";
imgMem->Initialize( mitk::MakePixelType<int, int, 1>() , 40, *planegeometry);
try
{
mitk::ImageReadAccessor imgMemAcc(imgMem);
p = (int*)imgMemAcc.GetData();
}
catch (mitk::Exception& e)
{
MITK_ERROR << e.what();
}
MITK_TEST_CONDITION_REQUIRED( p != NULL, "GetData() returned valid pointer.");
//-----------------
// testing origin information and methods
MITK_TEST_CONDITION_REQUIRED( mitk::Equal(imgMem->GetGeometry()->GetOrigin(), origin), "Testing correctness of origin via GetGeometry()->GetOrigin(): ");
// Setting origin via SetOrigin(origin): ";
mitk::FillVector3D(origin, 37.0, 17.92, 27.83); imgMem->SetOrigin(origin);
// Test origin
MITK_TEST_CONDITION_REQUIRED( mitk::Equal(imgMem->GetGeometry()->GetOrigin(), origin), "Testing correctness of changed origin via GetGeometry()->GetOrigin(): ");
MITK_TEST_CONDITION_REQUIRED( mitk::Equal(imgMem->GetSlicedGeometry()->GetPlaneGeometry(0)->GetOrigin(), origin), "Testing correctness of changed origin via GetSlicedGeometry()->GetPlaneGeometry(0)->GetOrigin(): ");
//-----------------
// testing spacing information and methodsunsigned int dim[]={100,100,20};
MITK_TEST_CONDITION_REQUIRED(mitk::Equal(imgMem->GetGeometry()->GetSpacing(), spacing), "Testing correct spacing from Geometry3D!");
mitk::FillVector3D(spacing, 7.0, 0.92, 1.83);
imgMem->SetSpacing(spacing);
MITK_TEST_CONDITION_REQUIRED( mitk::Equal(imgMem->GetGeometry()->GetSpacing(), spacing), "Testing correctness of changed spacing via GetGeometry()->GetSpacing(): ");
MITK_TEST_CONDITION_REQUIRED( mitk::Equal(imgMem->GetSlicedGeometry()->GetPlaneGeometry(0)->GetSpacing(), spacing), "Testing correctness of changed spacing via GetSlicedGeometry()->GetPlaneGeometry(0)->GetSpacing(): ");
mitk::Image::Pointer vecImg = mitk::Image::New();
try
{
mitk::ImageReadAccessor imgMemAcc(imgMem);
vecImg->Initialize( imgMem->GetPixelType(), *imgMem->GetGeometry(), 2 /* #channels */, 0 /*tDim*/ );
vecImg->SetImportChannel(const_cast<void*>(imgMemAcc.GetData()), 0, mitk::Image::CopyMemory );
vecImg->SetImportChannel(const_cast<void*>(imgMemAcc.GetData()), 1, mitk::Image::CopyMemory );
mitk::ImageReadAccessor vecImgAcc(vecImg);
mitk::ImageReadAccessor vecImgAcc0(vecImg, vecImg->GetChannelData(0));
mitk::ImageReadAccessor vecImgAcc1(vecImg, vecImg->GetChannelData(1));
MITK_TEST_CONDITION_REQUIRED(vecImgAcc0.GetData() != NULL && vecImgAcc1.GetData() != NULL, "Testing set and return of channel data!");
MITK_TEST_CONDITION_REQUIRED( vecImg->IsValidSlice(0,0,1) , "");
MITK_TEST_OUTPUT(<< " Testing whether CopyMemory worked");
MITK_TEST_CONDITION_REQUIRED(imgMemAcc.GetData() != vecImgAcc.GetData(), "");
MITK_TEST_OUTPUT(<< " Testing destruction after SetImportChannel");
vecImg = NULL;
MITK_TEST_CONDITION_REQUIRED(vecImg.IsNull() , "testing destruction!");
}
catch (mitk::Exception& e)
{
MITK_ERROR << e.what();
}
//-----------------
MITK_TEST_OUTPUT(<< "Testing initialization via vtkImageData");
MITK_TEST_OUTPUT(<< " Setting up vtkImageData");
vtkImageData* vtkimage = vtkImageData::New();
vtkimage->Initialize();
vtkimage->SetDimensions( 2, 3, 4);
double vtkorigin[] = {-350,-358.203, -1363.5};
vtkimage->SetOrigin(vtkorigin);
mitk::Point3D vtkoriginAsMitkPoint;
mitk::vtk2itk(vtkorigin, vtkoriginAsMitkPoint);
double vtkspacing[] = {1.367, 1.367, 2};
vtkimage->SetSpacing(vtkspacing);
vtkimage->AllocateScalars(VTK_SHORT,1);
std::cout<<"[PASSED]"<<std::endl;
MITK_TEST_OUTPUT(<< " Testing mitk::Image::Initialize(vtkImageData*, ...)");
mitk::Image::Pointer mitkByVtkImage = mitk::Image::New();
mitkByVtkImage ->Initialize(vtkimage);
MITK_TEST_CONDITION_REQUIRED(mitkByVtkImage->IsInitialized(), "");
vtkimage->Delete();
MITK_TEST_OUTPUT(<< " Testing whether spacing has been correctly initialized from vtkImageData");
mitk::Vector3D spacing2 = mitkByVtkImage->GetGeometry()->GetSpacing();
mitk::Vector3D vtkspacingAsMitkVector;
mitk::vtk2itk(vtkspacing, vtkspacingAsMitkVector);
MITK_TEST_CONDITION_REQUIRED(mitk::Equal(spacing2,vtkspacingAsMitkVector), "");
MITK_TEST_OUTPUT(<< " Testing whether GetSlicedGeometry(0)->GetOrigin() has been correctly initialized from vtkImageData");
mitk::Point3D origin2 = mitkByVtkImage->GetSlicedGeometry(0)->GetOrigin();
MITK_TEST_CONDITION_REQUIRED(mitk::Equal(origin2,vtkoriginAsMitkPoint), "");
MITK_TEST_OUTPUT(<< " Testing whether GetGeometry()->GetOrigin() has been correctly initialized from vtkImageData");
origin2 = mitkByVtkImage->GetGeometry()->GetOrigin();
MITK_TEST_CONDITION_REQUIRED(mitk::Equal(origin2,vtkoriginAsMitkPoint), "");
// TODO test the following initializers on channel-incorporation
// void mitk::Image::Initialize(const mitk::PixelType& type, unsigned int dimension, unsigned int *dimensions, unsigned int channels)
// void mitk::Image::Initialize(const mitk::PixelType& type, int sDim, const mitk::PlaneGeometry& geometry2d, bool flipped, unsigned int channels, int tDim )
// void mitk::Image::Initialize(const mitk::Image* image)
// void mitk::Image::Initialize(const mitkIpPicDescriptor* pic, int channels, int tDim, int sDim)
//mitk::Image::Pointer vecImg = mitk::Image::New();
//vecImg->Initialize(PixelType(typeid(float), 6, itk::ImageIOBase::SYMMETRICSECONDRANKTENSOR), *imgMem->GetGeometry(), 2 /* #channels */, 0 /*tDim*/, false /*shiftBoundingBoxMinimumToZero*/ );
//vecImg->Initialize(PixelType(typeid(itk::Vector<float,6>)), *imgMem->GetGeometry(), 2 /* #channels */, 0 /*tDim*/, false /*shiftBoundingBoxMinimumToZero*/ );
// testing access by index coordinates and by world coordinates
MITK_TEST_CONDITION_REQUIRED(argc == 2, "Check if test image is accessible!");
const std::string filename = std::string(argv[1]);
mitk::ItkImageFileReader::Pointer imageReader = mitk::ItkImageFileReader::New();
try
{
imageReader->SetFileName(filename);
imageReader->Update();
}
catch(...) {
MITK_TEST_FAILED_MSG(<< "Could not read file for testing: " << filename);
return 0;
}
mitk::Image::Pointer image = imageReader->GetOutput();
mitk::Point3D point;
mitk::ScalarType value = -1.;
mitkPixelTypeMultiplex3(TestRandomPixelAccess,image->GetImageDescriptor()->GetChannelTypeById(0),image,point,value)
{
// testing the clone method of mitk::Image
mitk::Image::Pointer cloneImage = image->Clone();
MITK_TEST_CONDITION_REQUIRED(cloneImage->GetDimension() == image->GetDimension(), "Clone (testing dimension)");
MITK_TEST_CONDITION_REQUIRED(cloneImage->GetPixelType() == image->GetPixelType(), "Clone (testing pixel type)");
// After cloning an image the geometry of both images should be equal too
MITK_TEST_CONDITION_REQUIRED(cloneImage->GetGeometry()->GetOrigin() == image->GetGeometry()->GetOrigin(), "Clone (testing origin)");
MITK_TEST_CONDITION_REQUIRED(cloneImage->GetGeometry()->GetSpacing() == image->GetGeometry()->GetSpacing(), "Clone (testing spacing)");
MITK_TEST_CONDITION_REQUIRED(mitk::MatrixEqualElementWise(cloneImage->GetGeometry()->GetIndexToWorldTransform()->GetMatrix(), image->GetGeometry()->GetIndexToWorldTransform()->GetMatrix()),
"Clone (testing transformation matrix)");
MITK_TEST_CONDITION_REQUIRED(mitk::MatrixEqualElementWise(cloneImage->GetTimeGeometry()->GetGeometryForTimeStep(cloneImage->GetDimension(3)-1)->GetIndexToWorldTransform()->GetMatrix(),
cloneImage->GetTimeGeometry()->GetGeometryForTimeStep(image->GetDimension(3)-1)->GetIndexToWorldTransform()->GetMatrix()), "Clone(testing time sliced geometry)");
for (unsigned int i = 0u; i < cloneImage->GetDimension(); ++i)
{
MITK_TEST_CONDITION_REQUIRED(cloneImage->GetDimension(i) == image->GetDimension(i), "Clone (testing dimension " << i << ")");
}
}
//access via itk
if(image->GetDimension()> 3) // CastToItk only works with 3d images so we need to check for 4d images
{
mitk::ImageTimeSelector::Pointer selector = mitk::ImageTimeSelector::New();
selector->SetTimeNr(0);
selector->SetInput(image);
selector->Update();
image = selector->GetOutput();
}
if(image->GetDimension()==3)
{
typedef itk::Image<double,3> ItkFloatImage3D;
ItkFloatImage3D::Pointer itkimage;
try
{
mitk::CastToItkImage(image, itkimage);
MITK_TEST_CONDITION_REQUIRED(itkimage.IsNotNull(), "Test conversion to itk::Image!");
}
catch (std::exception& e)
{
MITK_INFO << e.what();
}
mitk::Point3D itkPhysicalPoint;
image->GetGeometry()->WorldToItkPhysicalPoint(point, itkPhysicalPoint);
MITK_INFO << "ITKPoint " << itkPhysicalPoint[0] << " "<< itkPhysicalPoint[1] << " "<< itkPhysicalPoint[2] << "";
mitk::Point3D backTransformedPoint;
image->GetGeometry()->ItkPhysicalPointToWorld(itkPhysicalPoint, backTransformedPoint);
MITK_TEST_CONDITION_REQUIRED( mitk::Equal(point,backTransformedPoint), "Testing world->itk-physical->world consistency");
itk::Index<3> idx;
bool status = itkimage->TransformPhysicalPointToIndex(itkPhysicalPoint, idx);
MITK_INFO << "ITK Index " << idx[0] << " "<< idx[1] << " "<< idx[2] << "";
if(status && value != -1.)
{
float valByItk = itkimage->GetPixel(idx);
MITK_TEST_CONDITION_REQUIRED( mitk::Equal(valByItk, value), "Compare value of pixel returned by mitk in comparison to itk");
}
else
{
MITK_WARN<< "Index is out buffered region!";
}
}
else
{
MITK_INFO << "Image does not contain three dimensions, some test cases are skipped!";
}
// clone generated 3D image with one slice in z direction (cf. bug 11058)
unsigned int* threeDdim = new unsigned int[3];
threeDdim[0] = 100;
threeDdim[1] = 200;
threeDdim[2] = 1;
mitk::Image::Pointer threeDImage = mitk::Image::New();
threeDImage->Initialize(mitk::MakeScalarPixelType<float>(), 3, threeDdim);
mitk::Image::Pointer cloneThreeDImage = threeDImage->Clone();
// check that the clone image has the same dimensionality as the source image
MITK_TEST_CONDITION_REQUIRED( cloneThreeDImage->GetDimension() == 3, "Testing if the clone image initializes with 3D!");
MITK_TEST_CONDITION_REQUIRED( ImageVtkDataReferenceCheck(argv[1]), "Checking reference count of Image after using GetVtkImageData()");
MITK_TEST_END();
}