diff --git a/Core/Code/DataManagement/mitkBaseData.cpp b/Core/Code/DataManagement/mitkBaseData.cpp index 113522b2bd..e2a7530b25 100644 --- a/Core/Code/DataManagement/mitkBaseData.cpp +++ b/Core/Code/DataManagement/mitkBaseData.cpp @@ -1,296 +1,298 @@ /*=================================================================== 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 "mitkBaseData.h" #include #include #include +#include mitk::BaseData::BaseData() : m_RequestedRegionInitialized(false), m_SourceOutputIndexDuplicate(0), m_Initialized(true) { m_TimeGeometry = mitk::ProportionalTimeGeometry::New(); m_PropertyList = PropertyList::New(); } mitk::BaseData::BaseData( const BaseData &other ): itk::DataObject(), mitk::OperationActor(), m_RequestedRegionInitialized(other.m_RequestedRegionInitialized), m_SourceOutputIndexDuplicate(other.m_SourceOutputIndexDuplicate), m_Initialized(other.m_Initialized) { m_TimeGeometry = dynamic_cast(other.m_TimeGeometry->Clone().GetPointer()); m_PropertyList = other.m_PropertyList->Clone(); } mitk::BaseData::~BaseData() { } void mitk::BaseData::InitializeTimeGeometry(unsigned int timeSteps) { - mitk::BaseGeometry::Pointer g3d = mitk::BaseGeometry::New(); - g3d->Initialize(); + mitk::Geometry3D::Pointer geo3D = mitk::Geometry3D::New(); + mitk::BaseGeometry::Pointer baseGeo = dynamic_cast(geo3D.GetPointer()); + baseGeo->Initialize(); if ( timeSteps > 1 ) { mitk::ScalarType timeBounds[] = {0.0, 1.0}; - g3d->SetTimeBounds( timeBounds ); + baseGeo->SetTimeBounds( timeBounds ); } // The geometry is propagated automatically to the other items, // if EvenlyTimed is true... //Old timeGeometry->InitializeEvenlyTimed( g3d.GetPointer(), timeSteps ); TimeGeometry::Pointer timeGeometry = this->GetTimeGeometry(); timeGeometry->Initialize(); timeGeometry->Expand(timeSteps); for (TimeStepType step = 0; step < timeSteps; ++step) { - timeGeometry->SetTimeStepGeometry(g3d.GetPointer(),step); + timeGeometry->SetTimeStepGeometry(baseGeo.GetPointer(),step); } } void mitk::BaseData::UpdateOutputInformation() { if ( this->GetSource() ) { this->GetSource()->UpdateOutputInformation(); } if (m_TimeGeometry.IsNotNull()) { m_TimeGeometry->UpdateBoundingBox(); } } const mitk::TimeGeometry* mitk::BaseData::GetUpdatedTimeGeometry() { SetRequestedRegionToLargestPossibleRegion(); UpdateOutputInformation(); return GetTimeGeometry(); } void mitk::BaseData::Expand( unsigned int timeSteps ) { if (m_TimeGeometry.IsNotNull() ) { ProportionalTimeGeometry * propTimeGeometry = dynamic_cast (m_TimeGeometry.GetPointer()); if (propTimeGeometry) { propTimeGeometry->Expand(timeSteps); return; } mitkThrow() << "TimeGeometry is of an unkown Type. Could not expand it. "; } else { this->InitializeTimeGeometry(timeSteps); } } const mitk::BaseGeometry* mitk::BaseData::GetUpdatedGeometry(int t) { SetRequestedRegionToLargestPossibleRegion(); UpdateOutputInformation(); return GetGeometry(t); } void mitk::BaseData::SetGeometry(BaseGeometry* geometry) { ProportionalTimeGeometry::Pointer timeGeometry = ProportionalTimeGeometry::New(); if(geometry!=NULL) { timeGeometry->Initialize(geometry, 1); } SetTimeGeometry(timeGeometry); return; } void mitk::BaseData::SetTimeGeometry(TimeGeometry* geometry) { m_TimeGeometry = geometry; this->Modified(); } void mitk::BaseData::SetClonedGeometry(const BaseGeometry* aGeometry3D) { SetGeometry(static_cast(aGeometry3D->Clone().GetPointer())); } void mitk::BaseData::SetClonedTimeGeometry(const TimeGeometry* geometry) { TimeGeometry::Pointer clonedGeometry = geometry->Clone(); SetTimeGeometry(clonedGeometry.GetPointer()); } void mitk::BaseData::SetClonedGeometry(const BaseGeometry* aGeometry3D, unsigned int time) { if (m_TimeGeometry) { m_TimeGeometry->SetTimeStepGeometry(static_cast(aGeometry3D->Clone().GetPointer()),time); } } bool mitk::BaseData::IsEmptyTimeStep(unsigned int) const { return IsInitialized() == false; } bool mitk::BaseData::IsEmpty() const { if(IsInitialized() == false) return true; const TimeGeometry* timeGeometry = const_cast(this)->GetUpdatedTimeGeometry(); if(timeGeometry == NULL) return true; unsigned int timeSteps = timeGeometry->CountTimeSteps(); for ( unsigned int t = 0 ; t < timeSteps ; ++t ) { if(IsEmptyTimeStep(t) == false) return false; } return true; } itk::SmartPointer mitk::BaseData::GetSource() const { return static_cast(Superclass::GetSource().GetPointer()); } mitk::PropertyList::Pointer mitk::BaseData::GetPropertyList() const { return m_PropertyList; } mitk::BaseProperty::Pointer mitk::BaseData::GetProperty(const char *propertyKey) const { return m_PropertyList->GetProperty(propertyKey); } void mitk::BaseData::SetProperty(const char *propertyKey, BaseProperty* propertyValue) { m_PropertyList->SetProperty(propertyKey, propertyValue); } void mitk::BaseData::SetPropertyList(PropertyList *pList) { m_PropertyList = pList; } void mitk::BaseData::SetOrigin(const mitk::Point3D& origin) { TimeGeometry* timeGeom = GetTimeGeometry(); assert (timeGeom != NULL); BaseGeometry* geometry; TimeStepType steps = timeGeom->CountTimeSteps(); for (TimeStepType timestep = 0; timestep < steps; ++timestep) { geometry = GetGeometry(timestep); if (geometry != NULL) { geometry->SetOrigin(origin); } } } unsigned long mitk::BaseData::GetMTime() const { unsigned long time = Superclass::GetMTime(); if(m_TimeGeometry.IsNotNull()) { if((time < m_TimeGeometry->GetMTime())) { Modified(); return Superclass::GetMTime(); } } return time; } void mitk::BaseData::Graft(const itk::DataObject*) { itkExceptionMacro(<< "Graft not implemented for mitk::BaseData subclass " << this->GetNameOfClass()) } void mitk::BaseData::CopyInformation( const itk::DataObject* data ) { const Self* bd = dynamic_cast(data); if (bd != NULL) { m_PropertyList = bd->GetPropertyList()->Clone(); if (bd->GetTimeGeometry()!=NULL) { m_TimeGeometry = bd->GetTimeGeometry()->Clone(); } } else { // pointer could not be cast back down; this can be the case if your filters input // and output objects differ in type; then you have to write your own GenerateOutputInformation method itkExceptionMacro(<< "mitk::BaseData::CopyInformation() cannot cast " << typeid(data).name() << " to " << typeid(Self*).name() ); } } bool mitk::BaseData::IsInitialized() const { return m_Initialized; } void mitk::BaseData::Clear() { this->ClearData(); this->InitializeEmpty(); } void mitk::BaseData::ClearData() { if(m_Initialized) { ReleaseData(); m_Initialized = false; } } void mitk::BaseData::ExecuteOperation(mitk::Operation* /*operation*/) { //empty by default. override if needed! } void mitk::BaseData::PrintSelf(std::ostream& os, itk::Indent indent) const { os << std::endl; os << indent << " TimeGeometry: "; if(GetTimeGeometry() == NULL) os << "NULL" << std::endl; else GetTimeGeometry()->Print(os, indent); } diff --git a/Core/Code/DataManagement/mitkBaseGeometry.cpp b/Core/Code/DataManagement/mitkBaseGeometry.cpp index 1d6f1e7662..8b1c9166b4 100644 --- a/Core/Code/DataManagement/mitkBaseGeometry.cpp +++ b/Core/Code/DataManagement/mitkBaseGeometry.cpp @@ -1,1067 +1,1068 @@ /*=================================================================== 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 #include #include "mitkBaseGeometry.h" #include "mitkvector.h" #include "mitkMatrixConvert.h" #include #include #include "mitkRotationOperation.h" #include "mitkRestorePlanePositionOperation.h" #include "mitkApplyTransformMatrixOperation.h" #include "mitkPointOperation.h" #include "mitkInteractionConst.h" mitk::BaseGeometry::BaseGeometry(): Superclass(), mitk::OperationActor(), m_FrameOfReferenceID(0), m_IndexToWorldTransformLastModified(0), m_ImageGeometry(false) { m_VtkMatrix = vtkMatrix4x4::New(); m_VtkIndexToWorldTransform = vtkMatrixToLinearTransform::New(); m_VtkIndexToWorldTransform->SetInput(m_VtkMatrix); Initialize(); } mitk::BaseGeometry::BaseGeometry(const BaseGeometry& other): Superclass(), 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) { // 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) { if(origin!=GetOrigin()) { m_Origin = origin; m_IndexToWorldTransform->SetOffset(m_Origin.GetVectorFromOrigin()); Modified(); TransferItkToVtkTransform(); } } void mitk::BaseGeometry::TransferItkToVtkTransform() { // 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() ); newGeometry->SetIndexToWorldTransform(indexToWorldTransform); } this->PostInitializeGeometry(newGeometry); } void mitk::BaseGeometry::PostInitialize() { - m_ImageGeometry = false; } /** Set the bounds */ void mitk::BaseGeometry::SetBounds(const BoundsArrayType& bounds) { 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; iInsertElement(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) { 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; } 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 << "[( BaseGeometry::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; } bool mitk::Equal( const mitk::BaseGeometry::BoundingBoxType *leftHandSide, const mitk::BaseGeometry::BoundingBoxType *rightHandSide, ScalarType eps, bool verbose ) { bool result = true; if( rightHandSide == NULL ) { if(verbose) MITK_INFO << "[( BaseGeometry::BoundingBoxType )] rightHandSide NULL."; return false; } if( leftHandSide == NULL ) { if(verbose) MITK_INFO << "[( BaseGeometry::BoundingBoxType )] leftHandSide NULL."; return false; } 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 << "[( BaseGeometry::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) { bool result = true; if( rightHandSide == NULL ) { if(verbose) MITK_INFO << "[( Geometry3D )] rightHandSide NULL."; return false; } if( leftHandSide == NULL) { if(verbose) MITK_INFO << "[( Geometry3D )] leftHandSide NULL."; return false; } //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; } 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) { 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( index[0] ); discretIndex[1]=itk::Math::RoundHalfIntegerUp( index[1] ); discretIndex[2]=itk::Math::RoundHalfIntegerUp( 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(this)->m_Origin.GetVnlVector(); } vtkLinearTransform* mitk::BaseGeometry::GetVtkTransform() const { return (vtkLinearTransform*)m_VtkIndexToWorldTransform; } void mitk::BaseGeometry::SetIdentity() { 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 ) { 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 void mitk::BaseGeometry::ExecuteOperation(Operation* operation) { vtkTransform *vtktransform = vtkTransform::New(); vtktransform->SetMatrix(m_VtkMatrix); switch (operation->GetOperationType()) { case OpNOTHING: break; case OpMOVE: { mitk::PointOperation *pointOp = dynamic_cast(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(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(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(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(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) { 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(); -} +//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; } diff --git a/Core/Code/DataManagement/mitkBaseGeometry.h b/Core/Code/DataManagement/mitkBaseGeometry.h index d36b2f1c50..bd726e9e54 100644 --- a/Core/Code/DataManagement/mitkBaseGeometry.h +++ b/Core/Code/DataManagement/mitkBaseGeometry.h @@ -1,631 +1,630 @@ /*=================================================================== 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. ===================================================================*/ #ifndef BaseGeometry_H_HEADER_INCLUDED #define BaseGeometry_H_HEADER_INCLUDED #include #include #include "mitkoperationactor.h" #include #include "mitkvector.h" #include #include #include "itkScalableAffineTransform.h" #include class vtkMatrix4x4; class vtkMatrixToLinearTransform; class vtkLinearTransform; namespace mitk { //##Documentation //## @brief Standard 3D-BoundingBox typedef //## //## Standard 3D-BoundingBox typedef to get rid of template arguments (3D, type). typedef itk::BoundingBox BoundingBox; //##Documentation //## @brief Standard typedef for time-bounds typedef itk::FixedArray TimeBounds; typedef itk::FixedArray FixedArrayType; //##Documentation //## @brief BaseGeometry xxxxxxxxxxxxxx //## //## xxxxxxxxxxx //## //## Rule: everything is in mm (ms) if not stated otherwise. //## @ingroup Geometry class MITK_CORE_EXPORT BaseGeometry : public itk::Object, public OperationActor { public: mitkClassMacro(BaseGeometry, itk::Object); - /** Method for creation through the object factory. */ - itkNewMacro(Self); - mitkNewMacro1Param(Self,Self); - //##Documentation + + //##Documentation //## @brief clones the geometry //## //## Overwrite in all sub-classes. //## Normally looks like: //## \code //## Self::Pointer newGeometry = new Self(*this); //## newGeometry->UnRegister(); //## return newGeometry.GetPointer(); //## \endcode - virtual itk::LightObject::Pointer InternalClone() const; + virtual itk::LightObject::Pointer InternalClone() const =0; // ********************************** TypeDef ********************************** typedef itk::QuaternionRigidTransform< ScalarType > QuaternionTransformType; typedef QuaternionTransformType::VnlQuaternionType VnlQuaternionType; typedef itk::ScalableAffineTransform TransformType; typedef itk::BoundingBox BoundingBoxType; typedef BoundingBoxType::BoundsArrayType BoundsArrayType; typedef BoundingBoxType::Pointer BoundingBoxPointer; // ********************************** Origin, Spacing ********************************** //##Documentation //## @brief Get the origin, e.g. the upper-left corner of the plane const Point3D& GetOrigin() const; //##Documentation //## @brief Set the origin, i.e. the upper-left corner of the plane //## void SetOrigin(const Point3D& origin); //##Documentation //## @brief Get the spacing (size of a pixel). //## itkGetConstReferenceMacro(Spacing, mitk::Vector3D); //##Documentation //## @brief Set the spacing (m_Spacing) void SetSpacing(const mitk::Vector3D& aSpacing, bool enforceSetSpacing = false); //##Documentation //## @brief Get the origin as VnlVector //## //## \sa GetOrigin VnlVector GetOriginVnl() const; // ********************************** other functions ********************************** //##Documentation //## @brief Get the DICOM FrameOfReferenceID referring to the //## used world coordinate system itkGetConstMacro(FrameOfReferenceID, unsigned int); //##Documentation //## @brief Set the DICOM FrameOfReferenceID referring to the //## used world coordinate system itkSetMacro(FrameOfReferenceID, unsigned int); itkGetConstMacro(IndexToWorldTransformLastModified, unsigned long); //##Documentation //## @brief Is this Geometry3D in a state that is valid? virtual bool IsValid() const; // ********************************** Initialize ********************************** //##Documentation //## @brief Initialize the Geometry3D void Initialize(); void InitializeGeometry(Self * newGeometry) const; // ********************************** Transformations Set/Get ********************************** // a bit of a misuse, but we want only doxygen to see the following: #ifdef DOXYGEN_SKIP //##Documentation //## @brief Get the transformation used to convert from index //## to world coordinates itkGetObjectMacro(IndexToWorldTransform, AffineTransform3D); #endif //## @brief Set the transformation used to convert from index //## to world coordinates void SetIndexToWorldTransform(mitk::AffineTransform3D* transform); //##Documentation //## @brief Convenience method for setting the ITK transform //## (m_IndexToWorldTransform) via an vtkMatrix4x4 //## \sa SetIndexToWorldTransform virtual void SetIndexToWorldTransformByVtkMatrix(vtkMatrix4x4* vtkmatrix); /** Set/Get the IndexToWorldTransform */ itkGetConstObjectMacro(IndexToWorldTransform, AffineTransform3D); itkGetObjectMacro(IndexToWorldTransform, AffineTransform3D); vtkMatrix4x4* GetVtkMatrix(); //##Documentation //## @brief Get the m_IndexToWorldTransform as a vtkLinearTransform vtkLinearTransform* GetVtkTransform() const; //##Documentation //## @brief Set the transform to identity and origin to 0 //## virtual void SetIdentity(); // ********************************** Transformations ********************************** //##Documentation //## @brief Compose new IndexToWorldTransform with a given transform. //## //## This method composes m_IndexToWorldTransform with another transform, //## modifying self to be the composition of self and other. //## If the argument pre is true, then other is precomposed with self; //## that is, the resulting transformation consists of first applying //## other to the source, followed by self. If pre is false or omitted, //## then other is post-composed with self; that is the resulting //## transformation consists of first applying self to the source, //## followed by other. void Compose( const BaseGeometry::TransformType * other, bool pre = 0 ); //##Documentation //## @brief Compose new IndexToWorldTransform with a given vtkMatrix4x4. //## //## Converts the vtkMatrix4x4 into a itk-transform and calls the previous method. void Compose( const vtkMatrix4x4 * vtkmatrix, bool pre = 0 ); //##Documentation //## @brief Translate the origin by a vector //## void Translate(const Vector3D& vector); //##Documentation //##@brief executes affine operations (translate, rotate, scale) virtual void ExecuteOperation(Operation* operation); //##Documentation //## @brief Convert world coordinates (in mm) of a \em point to (continuous!) index coordinates //## \warning If you need (discrete) integer index coordinates (e.g., for iterating easily over an image), //## use WorldToIndex(const mitk::Point3D& pt_mm, itk::Index &index). //## For further information about coordinates types, please see the Geometry documentation void WorldToIndex(const mitk::Point3D& pt_mm, mitk::Point3D& pt_units) const; //##Documentation //## @brief Convert world coordinates (in mm) of a \em vector //## \a vec_mm to (continuous!) index coordinates. //## For further information about coordinates types, please see the Geometry documentation void WorldToIndex(const mitk::Vector3D& vec_mm, mitk::Vector3D& vec_units) const; //##Documentation //## @brief Convert world coordinates (in mm) of a \em point to (discrete!) index coordinates. //## This method rounds to integer indices! //## For further information about coordinates types, please see the Geometry documentation template void WorldToIndex(const mitk::Point3D& pt_mm, itk::Index &index) const { typedef itk::Index IndexType; mitk::Point3D pt_units; this->WorldToIndex(pt_mm, pt_units); int i, dim=index.GetIndexDimension(); if(dim>3) { index.Fill(0); dim=3; } for(i=0;i( pt_units[i] ); } } //##Documentation //## @brief Convert (continuous or discrete) index coordinates of a \em vector //## \a vec_units to world coordinates (in mm) //## For further information about coordinates types, please see the Geometry documentation void IndexToWorld(const mitk::Vector3D& vec_units, mitk::Vector3D& vec_mm) const; //##Documentation //## @brief Convert (continuous or discrete) index coordinates of a \em point to world coordinates (in mm) //## For further information about coordinates types, please see the Geometry documentation void IndexToWorld(const mitk::Point3D& pt_units, mitk::Point3D& pt_mm) const; //##Documentation //## @brief Convert (discrete) index coordinates of a \em point to world coordinates (in mm) //## For further information about coordinates types, please see the Geometry documentation template void IndexToWorld(const itk::Index &index, mitk::Point3D& pt_mm ) const { mitk::Point3D pt_units; pt_units.Fill(0); int i, dim=index.GetIndexDimension(); if(dim>3) { dim=3; } for(i=0;i void ItkPhysicalPointToWorld(const itk::Point& itkPhysicalPoint, mitk::Point3D& pt_mm) const { mitk::vtk2itk(itkPhysicalPoint, pt_mm); } //##Documentation //## @brief Deprecated for use with ITK version 3.10 or newer. //## Convert world coordinates (in mm) of a \em point to //## ITK physical coordinates (in mm, but without a possible rotation) //## //## This method is useful if you have want to access an mitk::Image //## via an itk::Image. ITK v3.8 and older did not support rotated (tilted) //## images, i.e., ITK images are always parallel to the coordinate axes. //## When accessing a (possibly rotated) mitk::Image via an itk::Image //## the rotational part of the transformation in the Geometry3D is //## simply discarded; in other word: only the origin and spacing is //## used by ITK, not the complete matrix available in MITK. //## With WorldToItkPhysicalPoint you can convert an MITK world //## coordinate (including the rotation) into a coordinate that //## can be used with the ITK image as a ITK physical coordinate //## (excluding the rotation). template void WorldToItkPhysicalPoint(const mitk::Point3D& pt_mm, itk::Point& itkPhysicalPoint) const { mitk::vtk2itk(pt_mm, itkPhysicalPoint); } // ********************************** BoundingBox ********************************** /** Get the bounding box */ itkGetConstObjectMacro(BoundingBox, BoundingBoxType); //##Documentation //## @brief Get the time bounds (in ms) itkGetConstReferenceMacro(TimeBounds, TimeBounds); // a bit of a misuse, but we want only doxygen to see the following: #ifdef DOXYGEN_SKIP //##Documentation //## @brief Get bounding box (in index/unit coordinates) itkGetConstObjectMacro(BoundingBox, BoundingBoxType); //##Documentation //## @brief Get bounding box (in index/unit coordinates) as a BoundsArrayType const BoundsArrayType GetBounds() const; #endif const BoundsArrayType GetBounds() const; //##Documentation //## \brief Set the bounding box (in index/unit coordinates) //## //## Only possible via the BoundsArray to make clear that a //## copy of the bounding-box is stored, not a reference to it. void SetBounds(const BoundsArrayType& bounds); //##Documentation //## @brief Set the bounding box (in index/unit coordinates) via a float array void SetFloatBounds(const float bounds[6]); //##Documentation //## @brief Set the bounding box (in index/unit coordinates) via a double array void SetFloatBounds(const double bounds[6]); //##Documentation //## @brief Get a VnlVector along bounding-box in the specified //## @a direction, length is spacing //## //## \sa GetAxisVector VnlVector GetMatrixColumn(unsigned int direction) const; //##Documentation //## @brief Calculates a bounding-box around the geometry relative //## to a coordinate system defined by a transform //## mitk::BoundingBox::Pointer CalculateBoundingBoxRelativeToTransform(const mitk::AffineTransform3D* transform) const; //##Documentation //## @brief Set the time bounds (in ms) void SetTimeBounds(const TimeBounds& timebounds); // ********************************** Geometry ********************************** #ifdef DOXYGEN_SKIP //##Documentation //## @brief Get the extent of the bounding box (in index/unit coordinates) //## //## To access the extent in mm use GetExtentInMM ScalarType GetExtent(unsigned int direction) const; #endif /** Get the extent of the bounding box */ ScalarType GetExtent(unsigned int direction) const; //##Documentation //## @brief Get the extent of the bounding-box in the specified @a direction in mm //## //## Equals length of GetAxisVector(direction). ScalarType GetExtentInMM(int direction) const; //##Documentation //## @brief Get vector along bounding-box in the specified @a direction in mm //## //## The length of the vector is the size of the bounding-box in the //## specified @a direction in mm //## \sa GetMatrixColumn Vector3D GetAxisVector(unsigned int direction) const; //##Documentation //## @brief Checks, if the given geometry can be converted to 2D without information loss //## e.g. when a 2D image is saved, the matrix is usually cropped to 2x2, and when you load it back to MITK //## it will be filled with standard values. This function checks, if information would be lost during this //## procedure virtual bool Is2DConvertable(); //##Documentation //## @brief Get the center of the bounding-box in mm //## Point3D GetCenter() const; //##Documentation //## @brief Get the squared length of the diagonal of the bounding-box in mm //## double GetDiagonalLength2() const; //##Documentation //## @brief Get the length of the diagonal of the bounding-box in mm //## double GetDiagonalLength() const; //##Documentation //## @brief Get the position of the corner number \a id (in world coordinates) //## //## See SetImageGeometry for how a corner is defined on images. Point3D GetCornerPoint(int id) const; //##Documentation //## @brief Get the position of a corner (in world coordinates) //## //## See SetImageGeometry for how a corner is defined on images. Point3D GetCornerPoint(bool xFront=true, bool yFront=true, bool zFront=true) const; //##Documentation //## @brief Set the extent of the bounding-box in the specified @a direction in mm //## //## @note This changes the matrix in the transform, @a not the bounds, which are given in units! void SetExtentInMM(int direction, ScalarType extentInMM); //##Documentation //## @brief Test whether the point \a p (world coordinates in mm) is //## inside the bounding box bool IsInside(const mitk::Point3D& p) const; //##Documentation //## @brief Test whether the point \a p ((continous!)index coordinates in units) is //## inside the bounding box bool IsIndexInside(const mitk::Point3D& index) const; //##Documentation //## @brief Convenience method for working with ITK indices template bool IsIndexInside(const itk::Index &index) const { int i, dim=index.GetIndexDimension(); Point3D pt_index; pt_index.Fill(0); for ( i = 0; i < dim; ++i ) { pt_index[i] = index[i]; } return IsIndexInside(pt_index); } // ********************************* Image Geometry ******************************** //##Documentation //## @brief When switching from an Image Geometry to a normal Geometry (and the other way around), you have to change the origin as well (See Geometry Documentation)! This function will change the "isImageGeometry" bool flag and changes the origin respectively. virtual void ChangeImageGeometryConsideringOriginOffset( const bool isAnImageGeometry ); //##Documentation //## @brief Is this an ImageGeometry? //## //## For more information, see SetImageGeometry itkGetConstMacro(ImageGeometry, bool); //##Documentation //## @brief Define that this Geometry3D is refering to an Image //## //## A geometry referring to an Image needs a slightly different //## definition of the position of the corners (see GetCornerPoint). //## The position of a voxel is defined by the position of its center. //## If we would use the origin (position of the (center of) the first //## voxel) as a corner and display this point, it would seem to be //## \em not at the corner but a bit within the image. Even worse for //## the opposite corner of the image: here the corner would appear //## outside the image (by half of the voxel diameter). Thus, we have //## to correct for this and to be able to do that, we need to know //## that the Geometry3D is referring to an Image. itkSetMacro(ImageGeometry, bool); itkBooleanMacro(ImageGeometry); protected: // ********************************** Constructor ********************************** BaseGeometry(); BaseGeometry(const BaseGeometry& other); virtual ~BaseGeometry(); virtual void PrintSelf(std::ostream& os, itk::Indent indent) const; void BackTransform(const mitk::Point3D& in, mitk::Point3D& out) const; //Without redundant parameter Point3D void BackTransform(const mitk::Vector3D& in, mitk::Vector3D& out) const; //##Documentation //## @brief Deprecated void BackTransform(const mitk::Point3D& at, const mitk::Vector3D& in, mitk::Vector3D& out) const; //##Documentation //## @brief Copy the ITK transform //## (m_IndexToWorldTransform) to the VTK transform //## \sa SetIndexToWorldTransform void TransferItkToVtkTransform(); //##Documentation //## @brief Copy the VTK transform //## to the ITK transform (m_IndexToWorldTransform) //## \sa SetIndexToWorldTransform void TransferVtkToItkTransform(); static const std::string GetTransformAsString( TransformType* transformType ); + virtual void PreSetBounds(const BoundsArrayType& bounds); + virtual void PostInitialize(); virtual void PostInitializeGeometry(Self * newGeometry) const; - virtual void PreSetBounds(const BoundsArrayType& bounds); virtual void PostSetExtentInMM(int direction, ScalarType extentInMM); virtual void PostSetTimeBounds(const TimeBounds& timebounds); virtual void PreSetIndexToWorldTransform(mitk::AffineTransform3D* transform); virtual void PostSetIndexToWorldTransform(mitk::AffineTransform3D* transform); virtual void PreSetSpacing(const mitk::Vector3D& aSpacing); void _SetSpacing(const mitk::Vector3D& aSpacing, bool enforceSetSpacing = false); itkGetConstMacro(NDimensions, unsigned int); bool IsBoundingBoxNull() const; bool IsIndexToWorldTransformNull() const; private: // ********************************** Variables ********************************** mitk::Vector3D m_Spacing; AffineTransform3D::Pointer m_IndexToWorldTransform; BoundingBoxPointer m_BoundingBox; vtkMatrixToLinearTransform* m_VtkIndexToWorldTransform; vtkMatrix4x4* m_VtkMatrix; unsigned int m_FrameOfReferenceID; mitk::TimeBounds m_TimeBounds; //##Documentation //## @brief Origin, i.e. upper-left corner of the plane //## Point3D m_Origin; static const unsigned int m_NDimensions = 3; mutable TransformType::Pointer m_InvertedTransform; mutable unsigned long m_IndexToWorldTransformLastModified; bool m_ImageGeometry; // DEPRECATED(VnlQuaternionType m_RotationQuaternion); }; // ********************************** Equal Functions ********************************** // // Static compare functions mainly for testing // /** * @brief Equal A function comparing two bounding boxes (BoundingBoxType) for beeing identical. * * @ingroup MITKTestingAPI * * The function compares the bounds (elementwise). * * The parameter eps is a tolarence value for all methods which are internally used for comparion. * @param rightHandSide Compare this against leftHandSide. * @param leftHandSide Compare this against rightHandSide. * @param eps Tolarence for comparison. You can use mitk::eps in most cases. * @param verbose Flag indicating if the user wants detailed console output or not. * @return True, if all comparison are true. False in any other case. */ MITK_CORE_EXPORT bool Equal( const mitk::BaseGeometry::BoundingBoxType *leftHandSide, const mitk::BaseGeometry::BoundingBoxType *rightHandSide, mitk::ScalarType eps, bool verbose); //ToDo // // Static compare functions mainly for testing // /** * @brief Equal A function comparing two geometries for beeing identical. * * @ingroup MITKTestingAPI * * The function compares the spacing, origin, axisvectors, extents, the matrix of the * IndexToWorldTransform (elementwise), the bounding (elementwise) and the ImageGeometry flag. * * The parameter eps is a tolarence value for all methods which are internally used for comparion. * If you want to use different tolarance values for different parts of the geometry, feel free to use * the other comparison methods and write your own implementation of Equal. * @param rightHandSide Compare this against leftHandSide. * @param leftHandSide Compare this against rightHandSide. * @param eps Tolarence for comparison. You can use mitk::eps in most cases. * @param verbose Flag indicating if the user wants detailed console output or not. * @return True, if all comparison are true. False in any other case. */ MITK_CORE_EXPORT bool Equal(const mitk::BaseGeometry* leftHandSide, const mitk::BaseGeometry* rightHandSide, mitk::ScalarType eps, bool verbose); //ToDo /** * @brief Equal A function comparing two transforms (TransformType) for beeing identical. * * @ingroup MITKTestingAPI * * The function compares the IndexToWorldTransform (elementwise). * * The parameter eps is a tolarence value for all methods which are internally used for comparion. * @param rightHandSide Compare this against leftHandSide. * @param leftHandSide Compare this against rightHandSide. * @param eps Tolarence for comparison. You can use mitk::eps in most cases. * @param verbose Flag indicating if the user wants detailed console output or not. * @return True, if all comparison are true. False in any other case. */ MITK_CORE_EXPORT bool Equal(const mitk::BaseGeometry::TransformType *leftHandSide, const mitk::BaseGeometry::TransformType *rightHandSide, mitk::ScalarType eps, bool verbose); //ToDo } // namespace mitk #endif /* BaseGeometry_H_HEADER_INCLUDED */ diff --git a/Core/Code/DataManagement/mitkImage.cpp b/Core/Code/DataManagement/mitkImage.cpp index de2fede13f..3bcfe03708 100644 --- a/Core/Code/DataManagement/mitkImage.cpp +++ b/Core/Code/DataManagement/mitkImage.cpp @@ -1,1384 +1,1384 @@ /*=================================================================== 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 #include "mitkImage.h" #include "mitkImageStatisticsHolder.h" #include "mitkPixelTypeMultiplex.h" #include #include "mitkCompareImageDataFilter.h" //VTK #include //Other #include #define FILL_C_ARRAY( _arr, _size, _value) for(unsigned int i=0u; i<_size; i++) \ { _arr[i] = _value; } mitk::Image::Image() : m_Dimension(0), m_Dimensions(NULL), m_ImageDescriptor(NULL), m_OffsetTable(NULL), m_CompleteData(NULL), m_ImageStatistics(NULL) { m_Dimensions = new unsigned int[MAX_IMAGE_DIMENSIONS]; FILL_C_ARRAY( m_Dimensions, MAX_IMAGE_DIMENSIONS, 0u); m_Initialized = false; } mitk::Image::Image(const Image &other) : SlicedData(other), m_Dimension(0), m_Dimensions(NULL), m_ImageDescriptor(NULL), m_OffsetTable(NULL), m_CompleteData(NULL), m_ImageStatistics(NULL) { m_Dimensions = new unsigned int[MAX_IMAGE_DIMENSIONS]; FILL_C_ARRAY( m_Dimensions, MAX_IMAGE_DIMENSIONS, 0u); this->Initialize( other.GetPixelType(), other.GetDimension(), other.GetDimensions()); //Since the above called "Initialize" method doesn't take the geometry into account we need to set it //here manually TimeGeometry::Pointer cloned = other.GetTimeGeometry()->Clone(); this->SetTimeGeometry(cloned.GetPointer()); if (this->GetDimension() > 3) { const unsigned int time_steps = this->GetDimension(3); for (unsigned int i = 0u; i < time_steps; ++i) { ImageDataItemPointer volume = const_cast(other).GetVolumeData(i); this->SetVolume(volume->GetData(), i); } } else { ImageDataItemPointer volume = const_cast(other).GetVolumeData(0); this->SetVolume(volume->GetData(), 0); } } mitk::Image::~Image() { Clear(); m_ReferenceCountLock.Lock(); m_ReferenceCount = 3; m_ReferenceCountLock.Unlock(); m_ReferenceCountLock.Lock(); m_ReferenceCount = 0; m_ReferenceCountLock.Unlock(); if(m_OffsetTable != NULL) delete [] m_OffsetTable; if(m_ImageStatistics != NULL) delete m_ImageStatistics; } const mitk::PixelType mitk::Image::GetPixelType(int n) const { return this->m_ImageDescriptor->GetChannelTypeById(n); } unsigned int mitk::Image::GetDimension() const { return m_Dimension; } unsigned int mitk::Image::GetDimension(int i) const { if((i>=0) && (i<(int)m_Dimension)) return m_Dimensions[i]; return 1; } void* mitk::Image::GetData() { if(m_Initialized==false) { if(GetSource().IsNull()) return NULL; if(GetSource()->Updating()==false) GetSource()->UpdateOutputInformation(); } m_CompleteData=GetChannelData(); // update channel's data // if data was not available at creation point, the m_Data of channel descriptor is NULL // if data present, it won't be overwritten m_ImageDescriptor->GetChannelDescriptor(0).SetData(m_CompleteData->GetData()); return m_CompleteData->GetData(); } template void AccessPixel( const mitk::PixelType ptype, void* data, const unsigned int offset, double& value ) { value = 0.0; if( data == NULL ) return; if(ptype.GetBpe() != 24) { value = (double) (((T*) data)[ offset ]); } else { const unsigned int rgboffset = 3 * offset; double returnvalue = (((T*) data)[rgboffset ]); returnvalue += (((T*) data)[rgboffset + 1]); returnvalue += (((T*) data)[rgboffset + 2]); value = returnvalue; } } double mitk::Image::GetPixelValueByIndex(const mitk::Index3D &position, unsigned int timestep) { double value = 0; if (this->GetTimeSteps() < timestep) { timestep = this->GetTimeSteps(); } value = 0.0; const unsigned int* imageDims = this->m_ImageDescriptor->GetDimensions(); const mitk::PixelType ptype = this->m_ImageDescriptor->GetChannelTypeById(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 ( position[0] < 0 || position[1] < 0 || position[2] < 0 ) { MITK_WARN << "Given position ("<< position << ") is out of image range, returning 0." ; } // check if the given position is inside the index range of the image, the 3rd dimension needs to be compared only if the dimension is not 0 else if ( (unsigned int)position[0] >= imageDims[0] || (unsigned int)position[1] >= imageDims[1] || ( imageDims[2] && (unsigned int)position[2] >= imageDims[2] )) { MITK_WARN << "Given position ("<< position << ") is out of image range, returning 0." ; } else { const unsigned int offset = position[0] + position[1]*imageDims[0] + position[2]*imageDims[0]*imageDims[1] + timestep*imageDims[0]*imageDims[1]*imageDims[2]; mitkPixelTypeMultiplex3( AccessPixel, ptype, this->GetData(), offset, value ); } return value; } double mitk::Image::GetPixelValueByWorldCoordinate(const mitk::Point3D& position, unsigned int timestep) { double value = 0.0; if (this->GetTimeSteps() < timestep) { timestep = this->GetTimeSteps(); } Index3D itkIndex; this->GetGeometry()->WorldToIndex(position, itkIndex); value = this->GetPixelValueByIndex( itkIndex, timestep); return value; } mitk::ImageVtkAccessor* mitk::Image::GetVtkImageData(int t, int n) { if(m_Initialized==false) { if(GetSource().IsNull()) return NULL; if(GetSource()->Updating()==false) GetSource()->UpdateOutputInformation(); } ImageDataItemPointer volume=GetVolumeData(t, n); if(volume.GetPointer()==NULL || volume->GetVtkImageData(this) == NULL) return NULL; SlicedGeometry3D* geom3d = GetSlicedGeometry(t); const mitk::Vector3D vspacing = (geom3d->GetSpacing()); double dspacing[3] = {vspacing[0],vspacing[1],vspacing[2]}; volume->GetVtkImageData(this)->SetSpacing( dspacing ); return volume->GetVtkImageData(this); } mitk::Image::ImageDataItemPointer mitk::Image::GetSliceData(int s, int t, int n, void *data, ImportMemoryManagementType importMemoryManagement) { if(IsValidSlice(s,t,n)==false) return NULL; const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); // slice directly available? int pos=GetSliceIndex(s,t,n); if(m_Slices[pos].GetPointer()!=NULL) return m_Slices[pos]; // is slice available as part of a volume that is available? ImageDataItemPointer sl, ch, vol; vol=m_Volumes[GetVolumeIndex(t,n)]; if((vol.GetPointer()!=NULL) && (vol->IsComplete())) { sl=new ImageDataItem(*vol, m_ImageDescriptor, 2, data, importMemoryManagement == ManageMemory, ((size_t) s)*m_OffsetTable[2]*(ptypeSize)); sl->SetComplete(true); return m_Slices[pos]=sl; } // is slice available as part of a channel that is available? ch=m_Channels[n]; if((ch.GetPointer()!=NULL) && (ch->IsComplete())) { sl=new ImageDataItem(*ch, m_ImageDescriptor, 2, data, importMemoryManagement == ManageMemory, (((size_t) s)*m_OffsetTable[2]+((size_t) t)*m_OffsetTable[3])*(ptypeSize)); sl->SetComplete(true); return m_Slices[pos]=sl; } // slice is unavailable. Can we calculate it? if((GetSource().IsNotNull()) && (GetSource()->Updating()==false)) { // ... wir mussen rechnen!!! .... m_RequestedRegion.SetIndex(0, 0); m_RequestedRegion.SetIndex(1, 0); m_RequestedRegion.SetIndex(2, s); m_RequestedRegion.SetIndex(3, t); m_RequestedRegion.SetIndex(4, n); m_RequestedRegion.SetSize(0, m_Dimensions[0]); m_RequestedRegion.SetSize(1, m_Dimensions[1]); m_RequestedRegion.SetSize(2, 1); m_RequestedRegion.SetSize(3, 1); m_RequestedRegion.SetSize(4, 1); m_RequestedRegionInitialized=true; GetSource()->Update(); if(IsSliceSet(s,t,n)) //yes: now we can call ourselves without the risk of a endless loop (see "if" above) return GetSliceData(s,t,n,data,importMemoryManagement); else return NULL; } else { ImageDataItemPointer item = AllocateSliceData(s,t,n,data,importMemoryManagement); item->SetComplete(true); return item; } } mitk::Image::ImageDataItemPointer mitk::Image::GetVolumeData(int t, int n, void *data, ImportMemoryManagementType importMemoryManagement) { if(IsValidVolume(t,n)==false) return NULL; ImageDataItemPointer ch, vol; // volume directly available? int pos=GetVolumeIndex(t,n); vol=m_Volumes[pos]; if((vol.GetPointer()!=NULL) && (vol->IsComplete())) return vol; const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); // is volume available as part of a channel that is available? ch=m_Channels[n]; if((ch.GetPointer()!=NULL) && (ch->IsComplete())) { vol=new ImageDataItem(*ch, m_ImageDescriptor, 3, data, importMemoryManagement == ManageMemory, (((size_t) t)*m_OffsetTable[3])*(ptypeSize)); vol->SetComplete(true); return m_Volumes[pos]=vol; } // let's see if all slices of the volume are set, so that we can (could) combine them to a volume bool complete=true; unsigned int s; for(s=0;sSetComplete(true); } else { mitk::PixelType chPixelType = this->m_ImageDescriptor->GetChannelTypeById(n); vol=m_Volumes[pos]; // ok, let's combine the slices! if(vol.GetPointer()==NULL) vol=new ImageDataItem( chPixelType, 3, m_Dimensions, NULL, true); vol->SetComplete(true); size_t size=m_OffsetTable[2]*(ptypeSize); for(s=0;sGetParent()!=vol) { // copy data of slices in volume size_t offset = ((size_t) s)*size; std::memcpy(static_cast(vol->GetData())+offset, sl->GetData(), size); // FIXME mitkIpPicDescriptor * pic = sl->GetPicDescriptor(); // replace old slice with reference to volume sl=new ImageDataItem(*vol, m_ImageDescriptor, 2, data, importMemoryManagement == ManageMemory, ((size_t) s)*size); sl->SetComplete(true); //mitkIpFuncCopyTags(sl->GetPicDescriptor(), pic); m_Slices[posSl]=sl; } } //if(vol->GetPicDescriptor()->info->tags_head==NULL) // mitkIpFuncCopyTags(vol->GetPicDescriptor(), m_Slices[GetSliceIndex(0,t,n)]->GetPicDescriptor()); } return m_Volumes[pos]=vol; } // volume is unavailable. Can we calculate it? if((GetSource().IsNotNull()) && (GetSource()->Updating()==false)) { // ... wir muessen rechnen!!! .... m_RequestedRegion.SetIndex(0, 0); m_RequestedRegion.SetIndex(1, 0); m_RequestedRegion.SetIndex(2, 0); m_RequestedRegion.SetIndex(3, t); m_RequestedRegion.SetIndex(4, n); m_RequestedRegion.SetSize(0, m_Dimensions[0]); m_RequestedRegion.SetSize(1, m_Dimensions[1]); m_RequestedRegion.SetSize(2, m_Dimensions[2]); m_RequestedRegion.SetSize(3, 1); m_RequestedRegion.SetSize(4, 1); m_RequestedRegionInitialized=true; GetSource()->Update(); if(IsVolumeSet(t,n)) //yes: now we can call ourselves without the risk of a endless loop (see "if" above) return GetVolumeData(t,n,data,importMemoryManagement); else return NULL; } else { ImageDataItemPointer item = AllocateVolumeData(t,n,data,importMemoryManagement); item->SetComplete(true); return item; } } mitk::Image::ImageDataItemPointer mitk::Image::GetChannelData(int n, void *data, ImportMemoryManagementType importMemoryManagement) { if(IsValidChannel(n)==false) return NULL; ImageDataItemPointer ch, vol; ch=m_Channels[n]; if((ch.GetPointer()!=NULL) && (ch->IsComplete())) return ch; // let's see if all volumes are set, so that we can (could) combine them to a channel if(IsChannelSet(n)) { // if there is only one time frame we do not need to combine anything if(m_Dimensions[3]<=1) { vol=GetVolumeData(0,n,data,importMemoryManagement); ch=new ImageDataItem(*vol, m_ImageDescriptor, m_ImageDescriptor->GetNumberOfDimensions(), data, importMemoryManagement == ManageMemory); ch->SetComplete(true); } else { const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); ch=m_Channels[n]; // ok, let's combine the volumes! if(ch.GetPointer()==NULL) ch=new ImageDataItem(this->m_ImageDescriptor, NULL, true); ch->SetComplete(true); size_t size=m_OffsetTable[m_Dimension-1]*(ptypeSize); unsigned int t; ImageDataItemPointerArray::iterator slicesIt = m_Slices.begin()+n*m_Dimensions[2]*m_Dimensions[3]; for(t=0;tGetParent()!=ch) { // copy data of volume in channel size_t offset = ((size_t) t)*m_OffsetTable[3]*(ptypeSize); std::memcpy(static_cast(ch->GetData())+offset, vol->GetData(), size); // REVEIW FIX mitkIpPicDescriptor * pic = vol->GetPicDescriptor(); // replace old volume with reference to channel vol=new ImageDataItem(*ch, m_ImageDescriptor, 3, data, importMemoryManagement == ManageMemory, offset); vol->SetComplete(true); //mitkIpFuncCopyTags(vol->GetPicDescriptor(), pic); m_Volumes[posVol]=vol; // get rid of slices - they may point to old volume ImageDataItemPointer dnull=NULL; for(unsigned int i = 0; i < m_Dimensions[2]; ++i, ++slicesIt) { assert(slicesIt != m_Slices.end()); *slicesIt = dnull; } } } // REVIEW FIX // if(ch->GetPicDescriptor()->info->tags_head==NULL) // mitkIpFuncCopyTags(ch->GetPicDescriptor(), m_Volumes[GetVolumeIndex(0,n)]->GetPicDescriptor()); } return m_Channels[n]=ch; } // channel is unavailable. Can we calculate it? if((GetSource().IsNotNull()) && (GetSource()->Updating()==false)) { // ... wir muessen rechnen!!! .... m_RequestedRegion.SetIndex(0, 0); m_RequestedRegion.SetIndex(1, 0); m_RequestedRegion.SetIndex(2, 0); m_RequestedRegion.SetIndex(3, 0); m_RequestedRegion.SetIndex(4, n); m_RequestedRegion.SetSize(0, m_Dimensions[0]); m_RequestedRegion.SetSize(1, m_Dimensions[1]); m_RequestedRegion.SetSize(2, m_Dimensions[2]); m_RequestedRegion.SetSize(3, m_Dimensions[3]); m_RequestedRegion.SetSize(4, 1); m_RequestedRegionInitialized=true; GetSource()->Update(); // did it work? if(IsChannelSet(n)) //yes: now we can call ourselves without the risk of a endless loop (see "if" above) return GetChannelData(n,data,importMemoryManagement); else return NULL; } else { ImageDataItemPointer item = AllocateChannelData(n,data,importMemoryManagement); item->SetComplete(true); return item; } } bool mitk::Image::IsSliceSet(int s, int t, int n) const { if(IsValidSlice(s,t,n)==false) return false; if(m_Slices[GetSliceIndex(s,t,n)].GetPointer()!=NULL) return true; ImageDataItemPointer ch, vol; vol=m_Volumes[GetVolumeIndex(t,n)]; if((vol.GetPointer()!=NULL) && (vol->IsComplete())) return true; ch=m_Channels[n]; if((ch.GetPointer()!=NULL) && (ch->IsComplete())) return true; return false; } bool mitk::Image::IsVolumeSet(int t, int n) const { if(IsValidVolume(t,n)==false) return false; ImageDataItemPointer ch, vol; // volume directly available? vol=m_Volumes[GetVolumeIndex(t,n)]; if((vol.GetPointer()!=NULL) && (vol->IsComplete())) return true; // is volume available as part of a channel that is available? ch=m_Channels[n]; if((ch.GetPointer()!=NULL) && (ch->IsComplete())) return true; // let's see if all slices of the volume are set, so that we can (could) combine them to a volume unsigned int s; for(s=0;sIsComplete())) return true; // let's see if all volumes are set, so that we can (could) combine them to a channel unsigned int t; for(t=0;t(data), s, t, n, CopyMemory); } bool mitk::Image::SetVolume(const void *data, int t, int n) { // const_cast is no risk for ImportMemoryManagementType == CopyMemory return SetImportVolume(const_cast(data), t, n, CopyMemory); } bool mitk::Image::SetChannel(const void *data, int n) { // const_cast is no risk for ImportMemoryManagementType == CopyMemory return SetImportChannel(const_cast(data), n, CopyMemory); } bool mitk::Image::SetImportSlice(void *data, int s, int t, int n, ImportMemoryManagementType importMemoryManagement) { if(IsValidSlice(s,t,n)==false) return false; ImageDataItemPointer sl; const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); if(IsSliceSet(s,t,n)) { sl=GetSliceData(s,t,n,data,importMemoryManagement); if(sl->GetManageMemory()==false) { sl=AllocateSliceData(s,t,n,data,importMemoryManagement); if(sl.GetPointer()==NULL) return false; } if ( sl->GetData() != data ) std::memcpy(sl->GetData(), data, m_OffsetTable[2]*(ptypeSize)); sl->Modified(); //we have changed the data: call Modified()! Modified(); } else { sl=AllocateSliceData(s,t,n,data,importMemoryManagement); if(sl.GetPointer()==NULL) return false; if ( sl->GetData() != data ) std::memcpy(sl->GetData(), data, m_OffsetTable[2]*(ptypeSize)); //we just added a missing slice, which is not regarded as modification. //Therefore, we do not call Modified()! } return true; } bool mitk::Image::SetImportVolume(void *data, int t, int n, ImportMemoryManagementType importMemoryManagement) { if(IsValidVolume(t,n)==false) return false; const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); ImageDataItemPointer vol; if(IsVolumeSet(t,n)) { vol=GetVolumeData(t,n,data,importMemoryManagement); if(vol->GetManageMemory()==false) { vol=AllocateVolumeData(t,n,data,importMemoryManagement); if(vol.GetPointer()==NULL) return false; } if ( vol->GetData() != data ) std::memcpy(vol->GetData(), data, m_OffsetTable[3]*(ptypeSize)); vol->Modified(); vol->SetComplete(true); //we have changed the data: call Modified()! Modified(); } else { vol=AllocateVolumeData(t,n,data,importMemoryManagement); if(vol.GetPointer()==NULL) return false; if ( vol->GetData() != data ) { std::memcpy(vol->GetData(), data, m_OffsetTable[3]*(ptypeSize)); } vol->SetComplete(true); this->m_ImageDescriptor->GetChannelDescriptor(n).SetData( vol->GetData() ); //we just added a missing Volume, which is not regarded as modification. //Therefore, we do not call Modified()! } return true; } bool mitk::Image::SetImportChannel(void *data, int n, ImportMemoryManagementType importMemoryManagement) { if(IsValidChannel(n)==false) return false; // channel descriptor const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); ImageDataItemPointer ch; if(IsChannelSet(n)) { ch=GetChannelData(n,data,importMemoryManagement); if(ch->GetManageMemory()==false) { ch=AllocateChannelData(n,data,importMemoryManagement); if(ch.GetPointer()==NULL) return false; } if ( ch->GetData() != data ) std::memcpy(ch->GetData(), data, m_OffsetTable[4]*(ptypeSize)); ch->Modified(); ch->SetComplete(true); //we have changed the data: call Modified()! Modified(); } else { ch=AllocateChannelData(n,data,importMemoryManagement); if(ch.GetPointer()==NULL) return false; if ( ch->GetData() != data ) std::memcpy(ch->GetData(), data, m_OffsetTable[4]*(ptypeSize)); ch->SetComplete(true); this->m_ImageDescriptor->GetChannelDescriptor(n).SetData( ch->GetData() ); //we just added a missing Channel, which is not regarded as modification. //Therefore, we do not call Modified()! } return true; } void mitk::Image::Initialize() { ImageDataItemPointerArray::iterator it, end; for( it=m_Slices.begin(), end=m_Slices.end(); it!=end; ++it ) { (*it)=NULL; } for( it=m_Volumes.begin(), end=m_Volumes.end(); it!=end; ++it ) { (*it)=NULL; } for( it=m_Channels.begin(), end=m_Channels.end(); it!=end; ++it ) { (*it)=NULL; } m_CompleteData = NULL; if( m_ImageStatistics == NULL) { m_ImageStatistics = new mitk::ImageStatisticsHolder( this ); } SetRequestedRegionToLargestPossibleRegion(); } void mitk::Image::Initialize(const mitk::ImageDescriptor::Pointer inDesc) { // store the descriptor this->m_ImageDescriptor = inDesc; // initialize image this->Initialize( inDesc->GetChannelDescriptor(0).GetPixelType(), inDesc->GetNumberOfDimensions(), inDesc->GetDimensions(), 1 ); } void mitk::Image::Initialize(const mitk::PixelType& type, unsigned int dimension, const unsigned int *dimensions, unsigned int channels) { Clear(); m_Dimension=dimension; if(!dimensions) itkExceptionMacro(<< "invalid zero dimension image"); unsigned int i; for(i=0;im_ImageDescriptor = mitk::ImageDescriptor::New(); this->m_ImageDescriptor->Initialize( this->m_Dimensions, this->m_Dimension ); for(i=0;i<4;++i) { m_LargestPossibleRegion.SetIndex(i, 0); m_LargestPossibleRegion.SetSize (i, m_Dimensions[i]); } m_LargestPossibleRegion.SetIndex(i, 0); m_LargestPossibleRegion.SetSize(i, channels); if(m_LargestPossibleRegion.GetNumberOfPixels()==0) { delete [] m_Dimensions; m_Dimensions = NULL; return; } for( unsigned int i=0u; im_ImageDescriptor->AddNewChannel( type ); } PlaneGeometry::Pointer planegeometry = PlaneGeometry::New(); planegeometry->InitializeStandardPlane(m_Dimensions[0], m_Dimensions[1]); SlicedGeometry3D::Pointer slicedGeometry = SlicedGeometry3D::New(); slicedGeometry->InitializeEvenlySpaced(planegeometry, m_Dimensions[2]); if(dimension>=4) { TimeBounds timebounds; timebounds[0] = 0.0; timebounds[1] = 1.0; slicedGeometry->SetTimeBounds(timebounds); } ProportionalTimeGeometry::Pointer timeGeometry = ProportionalTimeGeometry::New(); timeGeometry->Initialize(slicedGeometry, m_Dimensions[3]); for (TimeStepType step = 0; step < timeGeometry->CountTimeSteps(); ++step) { timeGeometry->GetGeometryForTimeStep(step)->ImageGeometryOn(); } SetTimeGeometry(timeGeometry); ImageDataItemPointer dnull=NULL; m_Channels.assign(GetNumberOfChannels(), dnull); m_Volumes.assign(GetNumberOfChannels()*m_Dimensions[3], dnull); m_Slices.assign(GetNumberOfChannels()*m_Dimensions[3]*m_Dimensions[2], dnull); ComputeOffsetTable(); Initialize(); m_Initialized = true; } void mitk::Image::Initialize(const mitk::PixelType& type, const mitk::BaseGeometry& geometry, unsigned int channels, int tDim ) { mitk::ProportionalTimeGeometry::Pointer timeGeometry = ProportionalTimeGeometry::New(); - BaseGeometry::Pointer geometry3D = geometry.Clone(); - timeGeometry->Initialize(geometry3D.GetPointer(), tDim); + itk::LightObject::Pointer lopointer = geometry.Clone(); + timeGeometry->Initialize(dynamic_cast(lopointer.GetPointer()), tDim); this->Initialize(type, *timeGeometry, channels, tDim); } void mitk::Image::Initialize(const mitk::PixelType& type, const mitk::TimeGeometry& geometry, unsigned int channels, int tDim ) { unsigned int dimensions[5]; dimensions[0] = (unsigned int)(geometry.GetGeometryForTimeStep(0)->GetExtent(0)+0.5); dimensions[1] = (unsigned int)(geometry.GetGeometryForTimeStep(0)->GetExtent(1)+0.5); dimensions[2] = (unsigned int)(geometry.GetGeometryForTimeStep(0)->GetExtent(2)+0.5); dimensions[3] = (tDim > 0) ? tDim : geometry.CountTimeSteps(); dimensions[4] = 0; unsigned int dimension = 2; if ( dimensions[2] > 1 ) dimension = 3; if ( dimensions[3] > 1 ) dimension = 4; Initialize( type, dimension, dimensions, channels ); if (geometry.CountTimeSteps() > 1) { TimeGeometry::Pointer cloned = geometry.Clone(); SetTimeGeometry(cloned.GetPointer()); } else Superclass::SetGeometry(geometry.GetGeometryForTimeStep(0)); /* //Old //TODO_GOETZ Really necessary? mitk::BoundingBox::BoundsArrayType bounds = geometry.GetBoundingBoxInWorld()->GetBounds(); if( (bounds[0] != 0.0) || (bounds[2] != 0.0) || (bounds[4] != 0.0) ) { SlicedGeometry3D* slicedGeometry = GetSlicedGeometry(0); mitk::Point3D origin; origin.Fill(0.0); slicedGeometry->IndexToWorld(origin, origin); bounds[1]-=bounds[0]; bounds[3]-=bounds[2]; bounds[5]-=bounds[4]; bounds[0] = 0.0; bounds[2] = 0.0; bounds[4] = 0.0; this->m_ImageDescriptor->Initialize( this->m_Dimensions, this->m_Dimension ); slicedGeometry->SetBounds(bounds); slicedGeometry->GetIndexToWorldTransform()->SetOffset(origin.GetVnlVector().data_block()); ProportionalTimeGeometry::Pointer timeGeometry = ProportionalTimeGeometry::New(); timeGeometry->Initialize(slicedGeometry, m_Dimensions[3]); SetTimeGeometry(timeGeometry); }*/ } void mitk::Image::Initialize(const mitk::PixelType& type, int sDim, const mitk::PlaneGeometry& geometry2d, bool flipped, unsigned int channels, int tDim ) { SlicedGeometry3D::Pointer slicedGeometry = SlicedGeometry3D::New(); slicedGeometry->InitializeEvenlySpaced(static_cast(geometry2d.Clone().GetPointer()), sDim, flipped); Initialize(type, *slicedGeometry, channels, tDim); } void mitk::Image::Initialize(const mitk::Image* image) { Initialize(image->GetPixelType(), *image->GetTimeGeometry()); } void mitk::Image::Initialize(vtkImageData* vtkimagedata, int channels, int tDim, int sDim, int pDim) { if(vtkimagedata==NULL) return; m_Dimension=vtkimagedata->GetDataDimension(); unsigned int i, *tmpDimensions=new unsigned int[m_Dimension>4?m_Dimension:4]; for(i=0;iGetDimensions()[i]; if(m_Dimension<4) { unsigned int *p; for(i=0,p=tmpDimensions+m_Dimension;i<4-m_Dimension;++i, ++p) *p=1; } if(pDim>=0) { tmpDimensions[1]=pDim; if(m_Dimension < 2) m_Dimension = 2; } if(sDim>=0) { tmpDimensions[2]=sDim; if(m_Dimension < 3) m_Dimension = 3; } if(tDim>=0) { tmpDimensions[3]=tDim; if(m_Dimension < 4) m_Dimension = 4; } switch ( vtkimagedata->GetScalarType() ) { case VTK_BIT: case VTK_CHAR: //pixelType.Initialize(typeid(char), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_UNSIGNED_CHAR: //pixelType.Initialize(typeid(unsigned char), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_SHORT: //pixelType.Initialize(typeid(short), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_UNSIGNED_SHORT: //pixelType.Initialize(typeid(unsigned short), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_INT: //pixelType.Initialize(typeid(int), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_UNSIGNED_INT: //pixelType.Initialize(typeid(unsigned int), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_LONG: //pixelType.Initialize(typeid(long), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_UNSIGNED_LONG: //pixelType.Initialize(typeid(unsigned long), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_FLOAT: //pixelType.Initialize(typeid(float), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; case VTK_DOUBLE: //pixelType.Initialize(typeid(double), vtkimagedata->GetNumberOfScalarComponents()); Initialize(mitk::MakeScalarPixelType(), m_Dimension, tmpDimensions, channels); break; default: break; } /* Initialize(pixelType, m_Dimension, tmpDimensions, channels); */ const double *spacinglist = vtkimagedata->GetSpacing(); Vector3D spacing; FillVector3D(spacing, spacinglist[0], 1.0, 1.0); if(m_Dimension>=2) spacing[1]=spacinglist[1]; if(m_Dimension>=3) spacing[2]=spacinglist[2]; // access origin of vtkImage Point3D origin; double vtkorigin[3]; vtkimagedata->GetOrigin(vtkorigin); FillVector3D(origin, vtkorigin[0], 0.0, 0.0); if(m_Dimension>=2) origin[1]=vtkorigin[1]; if(m_Dimension>=3) origin[2]=vtkorigin[2]; SlicedGeometry3D* slicedGeometry = GetSlicedGeometry(0); // re-initialize PlaneGeometry with origin and direction PlaneGeometry* planeGeometry = static_cast(slicedGeometry->GetGeometry2D(0)); planeGeometry->SetOrigin(origin); // re-initialize SlicedGeometry3D slicedGeometry->SetOrigin(origin); slicedGeometry->SetSpacing(spacing); ProportionalTimeGeometry::Pointer timeGeometry = ProportionalTimeGeometry::New(); timeGeometry->Initialize(slicedGeometry, m_Dimensions[3]); SetTimeGeometry(timeGeometry); delete [] tmpDimensions; } bool mitk::Image::IsValidSlice(int s, int t, int n) const { if(m_Initialized) return ((s>=0) && (s<(int)m_Dimensions[2]) && (t>=0) && (t< (int) m_Dimensions[3]) && (n>=0) && (n< (int)GetNumberOfChannels())); else return false; } bool mitk::Image::IsValidVolume(int t, int n) const { if(m_Initialized) return IsValidSlice(0, t, n); else return false; } bool mitk::Image::IsValidChannel(int n) const { if(m_Initialized) return IsValidSlice(0, 0, n); else return false; } void mitk::Image::ComputeOffsetTable() { if(m_OffsetTable!=NULL) delete [] m_OffsetTable; m_OffsetTable=new size_t[m_Dimension>4 ? m_Dimension+1 : 4+1]; unsigned int i; size_t num=1; m_OffsetTable[0] = 1; for (i=0; i < m_Dimension; ++i) { num *= m_Dimensions[i]; m_OffsetTable[i+1] = num; } for (;i < 4; ++i) m_OffsetTable[i+1] = num; } bool mitk::Image::IsValidTimeStep(int t) const { return ( ( m_Dimension >= 4 && t <= (int)m_Dimensions[3] && t > 0 ) || (t == 0) ); } void mitk::Image::Expand(unsigned int timeSteps) { if(timeSteps < 1) itkExceptionMacro(<< "Invalid timestep in Image!"); Superclass::Expand(timeSteps); } int mitk::Image::GetSliceIndex(int s, int t, int n) const { if(IsValidSlice(s,t,n)==false) return false; return ((size_t)s)+((size_t) t)*m_Dimensions[2]+((size_t) n)*m_Dimensions[3]*m_Dimensions[2]; //?? } int mitk::Image::GetVolumeIndex(int t, int n) const { if(IsValidVolume(t,n)==false) return false; return ((size_t)t)+((size_t) n)*m_Dimensions[3]; //?? } mitk::Image::ImageDataItemPointer mitk::Image::AllocateSliceData(int s, int t, int n, void *data, ImportMemoryManagementType importMemoryManagement) { int pos; pos=GetSliceIndex(s,t,n); const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); // is slice available as part of a volume that is available? ImageDataItemPointer sl, ch, vol; vol=m_Volumes[GetVolumeIndex(t,n)]; if(vol.GetPointer()!=NULL) { sl=new ImageDataItem(*vol, m_ImageDescriptor, 2, data, importMemoryManagement == ManageMemory, ((size_t) s)*m_OffsetTable[2]*(ptypeSize)); sl->SetComplete(true); return m_Slices[pos]=sl; } // is slice available as part of a channel that is available? ch=m_Channels[n]; if(ch.GetPointer()!=NULL) { sl=new ImageDataItem(*ch, m_ImageDescriptor, 2, data, importMemoryManagement == ManageMemory, (((size_t) s)*m_OffsetTable[2]+((size_t) t)*m_OffsetTable[3])*(ptypeSize)); sl->SetComplete(true); return m_Slices[pos]=sl; } // allocate new volume (instead of a single slice to keep data together!) m_Volumes[GetVolumeIndex(t,n)]=vol=AllocateVolumeData(t,n,NULL,importMemoryManagement); sl=new ImageDataItem(*vol, m_ImageDescriptor, 2, data, importMemoryManagement == ManageMemory, ((size_t) s)*m_OffsetTable[2]*(ptypeSize)); sl->SetComplete(true); return m_Slices[pos]=sl; ////ALTERNATIVE: //// allocate new slice //sl=new ImageDataItem(*m_PixelType, 2, m_Dimensions); //m_Slices[pos]=sl; //return vol; } mitk::Image::ImageDataItemPointer mitk::Image::AllocateVolumeData(int t, int n, void *data, ImportMemoryManagementType importMemoryManagement) { int pos; pos=GetVolumeIndex(t,n); const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); // is volume available as part of a channel that is available? ImageDataItemPointer ch, vol; ch=m_Channels[n]; if(ch.GetPointer()!=NULL) { vol=new ImageDataItem(*ch, m_ImageDescriptor, 3, data,importMemoryManagement == ManageMemory, (((size_t) t)*m_OffsetTable[3])*(ptypeSize)); return m_Volumes[pos]=vol; } mitk::PixelType chPixelType = this->m_ImageDescriptor->GetChannelTypeById(n); // allocate new volume if(importMemoryManagement == CopyMemory) { vol=new ImageDataItem( chPixelType, 3, m_Dimensions, NULL, true); if(data != NULL) std::memcpy(vol->GetData(), data, m_OffsetTable[3]*(ptypeSize)); } else { vol=new ImageDataItem( chPixelType, 3, m_Dimensions, data, importMemoryManagement == ManageMemory); } m_Volumes[pos]=vol; return vol; } mitk::Image::ImageDataItemPointer mitk::Image::AllocateChannelData(int n, void *data, ImportMemoryManagementType importMemoryManagement) { ImageDataItemPointer ch; // allocate new channel if(importMemoryManagement == CopyMemory) { const size_t ptypeSize = this->m_ImageDescriptor->GetChannelTypeById(n).GetSize(); ch=new ImageDataItem(this->m_ImageDescriptor, NULL, true); if(data != NULL) std::memcpy(ch->GetData(), data, m_OffsetTable[4]*(ptypeSize)); } else { ch=new ImageDataItem(this->m_ImageDescriptor, data, importMemoryManagement == ManageMemory); } m_Channels[n]=ch; return ch; } unsigned int* mitk::Image::GetDimensions() const { return m_Dimensions; } void mitk::Image::Clear() { Superclass::Clear(); delete [] m_Dimensions; m_Dimensions = NULL; } void mitk::Image::SetGeometry(BaseGeometry* aGeometry3D) { // Please be aware of the 0.5 offset/pixel-center issue! See Geometry documentation for further information if(aGeometry3D->GetImageGeometry()==false) { MITK_INFO << "WARNING: Applied a non-image geometry onto an image. Please be SURE that this geometry is pixel-center-based! If it is not, you need to call Geometry3D->ChangeImageGeometryConsideringOriginOffset(true) before calling image->setGeometry(..)\n"; } Superclass::SetGeometry(aGeometry3D); for (TimeStepType step = 0; step < GetTimeGeometry()->CountTimeSteps(); ++step) GetTimeGeometry()->GetGeometryForTimeStep(step)->ImageGeometryOn(); } void mitk::Image::PrintSelf(std::ostream& os, itk::Indent indent) const { unsigned char i; if(m_Initialized) { os << indent << " Dimension: " << m_Dimension << std::endl; os << indent << " Dimensions: "; for(i=0; i < m_Dimension; ++i) os << GetDimension(i) << " "; os << std::endl; for(unsigned int ch=0; ch < this->m_ImageDescriptor->GetNumberOfChannels(); ch++) { mitk::PixelType chPixelType = this->m_ImageDescriptor->GetChannelTypeById(ch); os << indent << " Channel: " << this->m_ImageDescriptor->GetChannelName(ch) << std::endl; os << indent << " PixelType: " << chPixelType.GetPixelTypeAsString() << std::endl; os << indent << " BytesPerElement: " << chPixelType.GetSize() << std::endl; os << indent << " ComponentType: " << chPixelType.GetComponentTypeAsString() << std::endl; os << indent << " NumberOfComponents: " << chPixelType.GetNumberOfComponents() << std::endl; os << indent << " BitsPerComponent: " << chPixelType.GetBitsPerComponent() << std::endl; } } else { os << indent << " Image not initialized: m_Initialized: false" << std::endl; } Superclass::PrintSelf(os,indent); } bool mitk::Image::IsRotated() const { const mitk::BaseGeometry* geo = this->GetGeometry(); bool ret = false; if(geo) { const vnl_matrix_fixed & mx = geo->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix(); mitk::ScalarType ref = 0; for(short k = 0; k < 3; ++k) ref += mx[k][k]; ref/=1000; // Arbitrary value; if a non-diagonal (nd) element is bigger then this, matrix is considered nd. for(short i = 0; i < 3; ++i) { for(short j = 0; j < 3; ++j) { if(i != j) { if(std::abs(mx[i][j]) > ref) // matrix is nd ret = true; } } } } return ret; } mitk::ScalarType mitk::Image::GetScalarValueMin(int t) const { return m_ImageStatistics->GetScalarValueMin(t); } //## \brief Get the maximum for scalar images mitk::ScalarType mitk::Image::GetScalarValueMax(int t) const { return m_ImageStatistics->GetScalarValueMax(t); } //## \brief Get the second smallest value for scalar images mitk::ScalarType mitk::Image::GetScalarValue2ndMin(int t) const { return m_ImageStatistics->GetScalarValue2ndMin(t); } mitk::ScalarType mitk::Image::GetScalarValueMinNoRecompute( unsigned int t ) const { return m_ImageStatistics->GetScalarValueMinNoRecompute(t); } mitk::ScalarType mitk::Image::GetScalarValue2ndMinNoRecompute( unsigned int t ) const { return m_ImageStatistics->GetScalarValue2ndMinNoRecompute(t); } mitk::ScalarType mitk::Image::GetScalarValue2ndMax(int t) const { return m_ImageStatistics->GetScalarValue2ndMax(t); } mitk::ScalarType mitk::Image::GetScalarValueMaxNoRecompute( unsigned int t) const { return m_ImageStatistics->GetScalarValueMaxNoRecompute(t); } mitk::ScalarType mitk::Image::GetScalarValue2ndMaxNoRecompute( unsigned int t ) const { return m_ImageStatistics->GetScalarValue2ndMaxNoRecompute(t); } mitk::ScalarType mitk::Image::GetCountOfMinValuedVoxels(int t ) const { return m_ImageStatistics->GetCountOfMinValuedVoxels(t); } mitk::ScalarType mitk::Image::GetCountOfMaxValuedVoxels(int t) const { return m_ImageStatistics->GetCountOfMaxValuedVoxels(t); } unsigned int mitk::Image::GetCountOfMaxValuedVoxelsNoRecompute( unsigned int t ) const { return m_ImageStatistics->GetCountOfMaxValuedVoxelsNoRecompute(t); } unsigned int mitk::Image::GetCountOfMinValuedVoxelsNoRecompute( unsigned int t ) const { return m_ImageStatistics->GetCountOfMinValuedVoxelsNoRecompute(t); } bool mitk::Equal(const mitk::Image* rightHandSide, const mitk::Image* leftHandSide, ScalarType eps, bool verbose) { bool returnValue = true; if( rightHandSide == NULL ) { if(verbose) MITK_INFO << "[( Image )] rightHandSide is NULL."; return false; } if( leftHandSide == NULL ) { if(verbose) MITK_INFO << "[( Image )] leftHandSide is NULL."; return false; } // Dimensionality if( rightHandSide->GetDimension() != leftHandSide->GetDimension() ) { if(verbose) { MITK_INFO << "[( Image )] Dimensionality differs."; MITK_INFO << "rightHandSide is " << rightHandSide->GetDimension() << "leftHandSide is " << leftHandSide->GetDimension(); } returnValue = false; } // Pair-wise dimension (size) comparison unsigned int minDimensionality = std::min(rightHandSide->GetDimension(),leftHandSide->GetDimension()); for( unsigned int i=0; i< minDimensionality; ++i) { if( rightHandSide->GetDimension(i) != leftHandSide->GetDimension(i) ) { returnValue = false; if(verbose) { MITK_INFO << "[( Image )] dimension differs."; MITK_INFO << "rightHandSide->GetDimension("<GetDimension(i) << "leftHandSide->GetDimension("<GetDimension(i); } } } // Pixeltype mitk::PixelType pixelTypeRightHandSide = rightHandSide->GetPixelType(); mitk::PixelType pixelTypeLeftHandSide = leftHandSide->GetPixelType(); if( !( pixelTypeRightHandSide == pixelTypeLeftHandSide ) ) { if(verbose) { MITK_INFO << "[( Image )] PixelType differs."; MITK_INFO << "rightHandSide is " << pixelTypeRightHandSide.GetTypeAsString() << "leftHandSide is " << pixelTypeLeftHandSide.GetTypeAsString(); } returnValue = false; } // Geometries if( !mitk::Equal( leftHandSide->GetGeometry(), rightHandSide->GetGeometry(), eps, verbose) ) { if(verbose) { MITK_INFO << "[( Image )] Geometries differ."; } returnValue = false; } // Pixel values - default mode [ 0 threshold in difference ] // compare only if all previous checks were successfull, otherwise the ITK filter will throw an exception if( returnValue ) { mitk::CompareImageDataFilter::Pointer compareFilter = mitk::CompareImageDataFilter::New(); compareFilter->SetInput(0, rightHandSide); compareFilter->SetInput(1, leftHandSide); compareFilter->SetTolerance(eps); compareFilter->Update(); if(( !compareFilter->GetResult() ) ) { returnValue = false; if(verbose) { MITK_INFO << "[(Image)] Pixel values differ: "; compareFilter->GetCompareResults().PrintSelf(); } } } return returnValue; } diff --git a/Core/Code/DataManagement/mitkProportionalTimeGeometry.cpp b/Core/Code/DataManagement/mitkProportionalTimeGeometry.cpp index f466c2dfcc..92b5f95c5c 100644 --- a/Core/Code/DataManagement/mitkProportionalTimeGeometry.cpp +++ b/Core/Code/DataManagement/mitkProportionalTimeGeometry.cpp @@ -1,237 +1,241 @@ /*=================================================================== 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 #include mitk::ProportionalTimeGeometry::ProportionalTimeGeometry() : m_FirstTimePoint(0.0), m_StepDuration(1.0) { } mitk::ProportionalTimeGeometry::~ProportionalTimeGeometry() { } void mitk::ProportionalTimeGeometry::Initialize() { m_FirstTimePoint = 0.0; m_StepDuration = 1.0; m_GeometryVector.resize(1); } mitk::TimeStepType mitk::ProportionalTimeGeometry::CountTimeSteps () const { return static_cast(m_GeometryVector.size() ); } mitk::TimePointType mitk::ProportionalTimeGeometry::GetMinimumTimePoint () const { return m_FirstTimePoint; } mitk::TimePointType mitk::ProportionalTimeGeometry::GetMaximumTimePoint () const { TimePointType timePoint = m_FirstTimePoint + m_StepDuration * CountTimeSteps(); if (timePoint >std::numeric_limits().max()) timePoint = std::numeric_limits().max(); return timePoint; } mitk::TimeBounds mitk::ProportionalTimeGeometry::GetTimeBounds () const { TimeBounds bounds; bounds[0] = this->GetMinimumTimePoint(); bounds[1] = this->GetMaximumTimePoint(); return bounds; } bool mitk::ProportionalTimeGeometry::IsValidTimePoint (TimePointType timePoint) const { return this->GetMinimumTimePoint() <= timePoint && timePoint < this->GetMaximumTimePoint(); } bool mitk::ProportionalTimeGeometry::IsValidTimeStep (TimeStepType timeStep) const { return timeStep < this->CountTimeSteps(); } mitk::TimePointType mitk::ProportionalTimeGeometry::TimeStepToTimePoint( TimeStepType timeStep) const { if (m_FirstTimePoint <= std::numeric_limits::min() || m_FirstTimePoint >= std::numeric_limits::max() || m_StepDuration <= std::numeric_limits::min() || m_StepDuration >= std::numeric_limits::max()) { return static_cast(timeStep); } return m_FirstTimePoint + timeStep * m_StepDuration; } mitk::TimeStepType mitk::ProportionalTimeGeometry::TimePointToTimeStep( TimePointType timePoint) const { if (m_FirstTimePoint <= timePoint) return static_cast((timePoint -m_FirstTimePoint) / m_StepDuration); else return 0; } mitk::BaseGeometry::Pointer mitk::ProportionalTimeGeometry::GetGeometryForTimeStep( TimeStepType timeStep) const { if (IsValidTimeStep(timeStep)) { return dynamic_cast(m_GeometryVector[timeStep].GetPointer()); } else { return 0; } } mitk::BaseGeometry::Pointer mitk::ProportionalTimeGeometry::GetGeometryForTimePoint(TimePointType timePoint) const { if (this->IsValidTimePoint(timePoint)) { TimeStepType timeStep = this->TimePointToTimeStep(timePoint); return this->GetGeometryForTimeStep(timeStep); } else { return 0; } } mitk::BaseGeometry::Pointer mitk::ProportionalTimeGeometry::GetGeometryCloneForTimeStep( TimeStepType timeStep) const { if (timeStep > m_GeometryVector.size()) return 0; - return m_GeometryVector[timeStep]->Clone(); + itk::LightObject::Pointer lopointer = m_GeometryVector[timeStep]->Clone(); + return dynamic_cast(lopointer.GetPointer()); } bool mitk::ProportionalTimeGeometry::IsValid() const { bool isValid = true; isValid &= m_GeometryVector.size() > 0; isValid &= m_StepDuration > 0; return isValid; } void mitk::ProportionalTimeGeometry::ClearAllGeometries() { m_GeometryVector.clear(); } void mitk::ProportionalTimeGeometry::ReserveSpaceForGeometries(TimeStepType numberOfGeometries) { m_GeometryVector.reserve(numberOfGeometries); } void mitk::ProportionalTimeGeometry::Expand(mitk::TimeStepType size) { m_GeometryVector.reserve(size); while (m_GeometryVector.size() < size) { - m_GeometryVector.push_back(BaseGeometry::New()); + Geometry3D::Pointer geo3D = Geometry3D::New(); + m_GeometryVector.push_back(dynamic_cast(geo3D.GetPointer())); } } void mitk::ProportionalTimeGeometry::SetTimeStepGeometry(BaseGeometry* geometry, TimeStepType timeStep) { assert(timeStep<=m_GeometryVector.size()); if (timeStep == m_GeometryVector.size()) m_GeometryVector.push_back(geometry); m_GeometryVector[timeStep] = geometry; } itk::LightObject::Pointer mitk::ProportionalTimeGeometry::InternalClone() const { itk::LightObject::Pointer parent = Superclass::InternalClone(); ProportionalTimeGeometry::Pointer newTimeGeometry = dynamic_cast (parent.GetPointer()); newTimeGeometry->m_FirstTimePoint = this->m_FirstTimePoint; newTimeGeometry->m_StepDuration = this->m_StepDuration; newTimeGeometry->m_GeometryVector.clear(); newTimeGeometry->Expand(this->CountTimeSteps()); for (TimeStepType i =0; i < CountTimeSteps(); ++i) { - BaseGeometry::Pointer tempGeometry = GetGeometryForTimeStep(i)->Clone(); + itk::LightObject::Pointer lopointer=GetGeometryForTimeStep(i)->Clone(); + BaseGeometry::Pointer tempGeometry = dynamic_cast(lopointer.GetPointer()); newTimeGeometry->SetTimeStepGeometry(tempGeometry.GetPointer(),i); } return parent; } void mitk::ProportionalTimeGeometry::Initialize (BaseGeometry* geometry, TimeStepType timeSteps) { timeSteps = (timeSteps > 0) ? timeSteps : 1; m_FirstTimePoint = geometry->GetTimeBounds()[0]; m_StepDuration = geometry->GetTimeBounds()[1] - geometry->GetTimeBounds()[0]; this->ReserveSpaceForGeometries(timeSteps); try{ for (TimeStepType currentStep = 0; currentStep < timeSteps; ++currentStep) { mitk::TimeBounds timeBounds; if (timeSteps > 1) { timeBounds[0] = m_FirstTimePoint + currentStep * m_StepDuration; timeBounds[1] = m_FirstTimePoint + (currentStep+1) * m_StepDuration; } else { timeBounds = geometry->GetTimeBounds(); } - - BaseGeometry::Pointer clonedGeometry = geometry->Clone(); + itk::LightObject::Pointer lopointer=geometry->Clone(); + BaseGeometry::Pointer clonedGeometry = dynamic_cast(lopointer.GetPointer()); this->SetTimeStepGeometry(clonedGeometry.GetPointer(), currentStep); GetGeometryForTimeStep(currentStep)->SetTimeBounds(timeBounds); } } catch (...) { MITK_INFO << "Cloning of geometry produced an error!"; } Update(); } void mitk::ProportionalTimeGeometry::Initialize (TimeStepType timeSteps) { - mitk::BaseGeometry::Pointer geometry = mitk::BaseGeometry::New(); + mitk::Geometry3D::Pointer geo3D = Geometry3D::New(); + mitk::BaseGeometry::Pointer geometry = dynamic_cast(geo3D.GetPointer()); geometry->Initialize(); if ( timeSteps > 1 ) { mitk::ScalarType timeBounds[] = {0.0, 1.0}; geometry->SetTimeBounds( timeBounds ); } this->Initialize(geometry.GetPointer(), timeSteps); } void mitk::ProportionalTimeGeometry::PrintSelf(std::ostream& os, itk::Indent indent) const { os << indent << " TimeSteps: " << this->CountTimeSteps() << std::endl; os << indent << " FirstTimePoint: " << this->GetFirstTimePoint() << std::endl; os << indent << " StepDuration: " << this->GetStepDuration() << " ms" << std::endl; os << std::endl; os << indent << " GetGeometryForTimeStep(0): "; if(GetGeometryForTimeStep(0).IsNull()) os << "NULL" << std::endl; else GetGeometryForTimeStep(0)->Print(os, indent); } diff --git a/Core/Code/DataManagement/mitkProportionalTimeGeometry.h b/Core/Code/DataManagement/mitkProportionalTimeGeometry.h index 60c2d0b449..a76e300754 100644 --- a/Core/Code/DataManagement/mitkProportionalTimeGeometry.h +++ b/Core/Code/DataManagement/mitkProportionalTimeGeometry.h @@ -1,206 +1,207 @@ /*=================================================================== 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. ===================================================================*/ #ifndef ProportionalTimeGeometry_h #define ProportionalTimeGeometry_h //MITK #include +#include #include #include namespace mitk { /** * \brief Organizes geometries over proportional time steps * * For this TimeGeometry implementation it is assumed that * the durations of the time steps are equidistant, e.g. * the durations of the time steps in one ProportionalTimeGeometry * are the same. The geometries of the time steps are independent, * and not linked to each other. Since the timeBounds of the * geometries are different for each time step it is not possible * to set the same geometry to different time steps. Instead * copies should be used. * * \addtogroup geometry */ class MITK_CORE_EXPORT ProportionalTimeGeometry : public TimeGeometry { public: mitkClassMacro(ProportionalTimeGeometry, TimeGeometry); ProportionalTimeGeometry(); typedef ProportionalTimeGeometry self; itkNewMacro(self); /** * \brief Returns the number of time steps. * * Returns the number of time steps for which * geometries are saved. The number of time steps * is also the upper bound of the time steps. The * minimum time steps is always 0. */ virtual TimeStepType CountTimeSteps() const; /** * \brief Returns the first time point for which the object is valid. * * Returns the first valid time point for this geometry. If only one * time steps available it usually goes from -max to +max. The time point * is given in ms. */ virtual TimePointType GetMinimumTimePoint () const; /** * \brief Returns the last time point for which the object is valid * * Gives the last time point for which a valid geometrie is saved in * this time geometry. The time point is given in ms. */ virtual TimePointType GetMaximumTimePoint () const; /** * \brief Get the time bounds (in ms) */ virtual TimeBounds GetTimeBounds( ) const; /** * \brief Tests if a given time point is covered by this object * * Returns true if a geometry can be returned for the given time * point and falls if not. The time point must be given in ms. */ virtual bool IsValidTimePoint (TimePointType timePoint) const; /** * \brief Test for the given time step if a geometry is availible * * Returns true if a geometry is defined for the given time step. * Otherwise false is returned. * The time step is defined as positiv number. */ virtual bool IsValidTimeStep (TimeStepType timeStep) const; /** * \brief Converts a time step to a time point * * Converts a time step to a time point in a way that * the new time point indicates the same geometry as the time step. * If the original time steps does not point to a valid geometry, * a time point is calculated that also does not point to a valid * geometry, but no exception is raised. */ virtual TimePointType TimeStepToTimePoint (TimeStepType timeStep) const; /** * \brief Converts a time point to the corresponding time step * * Converts a time point to a time step in a way that * the new time step indicates the same geometry as the time point. * If a negativ invalid time point is given always time step 0 is * returned. If an positiv invalid time step is given an invalid * time step will be returned. */ virtual TimeStepType TimePointToTimeStep (TimePointType timePoint) const; /** * \brief Returns the geometry which corresponds to the given time step * * Returns a clone of the geometry which defines the given time step. If * the given time step is invalid an null-pointer is returned. */ virtual BaseGeometry::Pointer GetGeometryCloneForTimeStep( TimeStepType timeStep) const; /** * \brief Returns the geometry which corresponds to the given time point * * Returns the geometry which defines the given time point. If * the given time point is invalid an null-pointer is returned. * * If the returned geometry is changed this will affect the saved * geometry. */ virtual BaseGeometry::Pointer GetGeometryForTimePoint ( TimePointType timePoint) const; /** * \brief Returns the geometry which corresponds to the given time step * * Returns the geometry which defines the given time step. If * the given time step is invalid an null-pointer is returned. * * If the returned geometry is changed this will affect the saved * geometry. */ virtual BaseGeometry::Pointer GetGeometryForTimeStep ( TimeStepType timeStep) const; /** * \brief Tests if all necessary informations are set and the object is valid */ virtual bool IsValid () const; /** * \brief Initilizes a new object with one time steps which contains an empty geometry. */ virtual void Initialize(); /** * \brief Expands the time geometry to the given number of time steps. * * Initializes the new time steps with empty geometries. * Shrinking is not supported. */ virtual void Expand(TimeStepType size); /** * \brief Sets the geometry for the given time step * * This method does not afflict other time steps, since the geometry for * each time step is saved individually. */ virtual void SetTimeStepGeometry(BaseGeometry* geometry, TimeStepType timeStep); /** * \brief Makes a deep copy of the current object */ virtual itk::LightObject::Pointer InternalClone () const; itkGetConstMacro(FirstTimePoint, TimePointType); itkSetMacro(FirstTimePoint, TimePointType); itkGetConstMacro(StepDuration, TimePointType); itkSetMacro(StepDuration, TimePointType); // void SetGeometryForTimeStep(TimeStepType timeStep, BaseGeometry& geometry); void ClearAllGeometries (); // void AddGeometry(BaseGeometry geometry); void ReserveSpaceForGeometries (TimeStepType numberOfGeometries); /** * \brief Initializes the TimeGeometry with equally time Step geometries * * Saves a copy for each time step. */ void Initialize (BaseGeometry* geometry, TimeStepType timeSteps); /** * \brief Initialize the TimeGeometry with empty BaseGeometry */ void Initialize (TimeStepType timeSteps); virtual void PrintSelf(std::ostream& os, itk::Indent indent) const; protected: virtual ~ProportionalTimeGeometry(); std::vector m_GeometryVector; TimePointType m_FirstTimePoint; TimePointType m_StepDuration; }; // end class ProportialTimeGeometry } // end namespace MITK #endif // ProportionalTimeGeometry_h diff --git a/Core/Code/Testing/mitkBaseDataTest.cpp b/Core/Code/Testing/mitkBaseDataTest.cpp index d6c16af0c1..c0a8a3f3a9 100644 --- a/Core/Code/Testing/mitkBaseDataTest.cpp +++ b/Core/Code/Testing/mitkBaseDataTest.cpp @@ -1,120 +1,121 @@ /*=================================================================== 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 "mitkBaseDataTestImplementation.h" #include "mitkStringProperty.h" #include "mitkTestingMacros.h" #include #include #include "itkImage.h" int mitkBaseDataTest(int /*argc*/, char* /*argv*/[]) { MITK_TEST_BEGIN("BaseData") //Create a BaseData implementation MITK_INFO << "Creating a base data instance..."; mitk::BaseDataTestImplementation::Pointer baseDataImpl = mitk::BaseDataTestImplementation::New(); MITK_TEST_CONDITION_REQUIRED(baseDataImpl.IsNotNull(),"Testing instantiation"); MITK_TEST_CONDITION(baseDataImpl->IsInitialized(), "BaseDataTestImplementation is initialized"); MITK_TEST_CONDITION(baseDataImpl->IsEmpty(), "BaseDataTestImplementation is initialized and empty"); mitk::BaseDataTestImplementation::Pointer cloneBaseData = baseDataImpl->Clone(); MITK_TEST_CONDITION_REQUIRED(cloneBaseData.IsNotNull(),"Testing instantiation of base data clone"); MITK_TEST_CONDITION(cloneBaseData->IsInitialized(), "Clone of BaseDataTestImplementation is initialized"); MITK_TEST_CONDITION(cloneBaseData->IsEmpty(), "Clone of BaseDataTestImplementation is initialized and empty"); MITK_INFO << "Testing setter and getter for geometries..."; //test method GetTimeGeometry() MITK_TEST_CONDITION(baseDataImpl->GetTimeGeometry(), "Testing creation of TimeGeometry"); mitk::TimeGeometry* geo = NULL; baseDataImpl->SetTimeGeometry(geo); MITK_TEST_CONDITION(baseDataImpl->GetTimeGeometry() == NULL, "Reset Geometry"); mitk::ProportionalTimeGeometry::Pointer geo2 = mitk::ProportionalTimeGeometry::New(); baseDataImpl->SetTimeGeometry(geo2); geo2->Initialize(2); MITK_TEST_CONDITION(baseDataImpl->GetTimeGeometry() == geo2.GetPointer(), "Correct Reinit of TimeGeometry"); //test method GetGeometry(int timeStep) MITK_TEST_CONDITION(baseDataImpl->GetGeometry(1) != NULL, "... and single Geometries"); //test method Expand(unsigned int timeSteps) baseDataImpl->Expand(5); MITK_TEST_CONDITION(baseDataImpl->GetTimeSteps() == 5, "Expand the geometry to further time slices!"); //test method GetUpdatedGeometry(int timeStep); - mitk::BaseGeometry::Pointer geo3 = mitk::BaseGeometry::New(); + mitk::Geometry3D::Pointer geometry3D = mitk::Geometry3D::New(); + mitk::BaseGeometry::Pointer geo3 = dynamic_cast(geometry3D.GetPointer()); mitk::ProportionalTimeGeometry::Pointer timeGeometry = dynamic_cast(baseDataImpl->GetTimeGeometry()); if (timeGeometry.IsNotNull() ) { timeGeometry->SetTimeStepGeometry(geo3,1); } MITK_TEST_CONDITION(baseDataImpl->GetUpdatedGeometry(1) == geo3, "Set Geometry for time step 1"); MITK_TEST_CONDITION(baseDataImpl->GetMTime()!= 0, "Check if modified time is set"); baseDataImpl->SetClonedGeometry(geo3, 1); mitk::ScalarType x[3]; x[0] = 2; x[1] = 4; x[2] = 6; mitk::Point3D p3d(x); baseDataImpl->SetOrigin(p3d); geo3->SetOrigin(p3d); MITK_TEST_CONDITION(baseDataImpl->GetGeometry(1)->GetOrigin() == geo3->GetOrigin(), "Testing Origin set"); cloneBaseData = baseDataImpl->Clone(); MITK_TEST_CONDITION(cloneBaseData->GetGeometry(1)->GetOrigin() == geo3->GetOrigin(), "Testing origin set in clone!"); MITK_TEST_CONDITION(!baseDataImpl->IsEmptyTimeStep(1), "Is not empty before clear()!"); baseDataImpl->Clear(); MITK_TEST_CONDITION(baseDataImpl->IsEmptyTimeStep(1), "...but afterwards!"); //test method Set-/GetProperty() baseDataImpl->SetProperty("property38", mitk::StringProperty::New("testproperty")); //baseDataImpl->SetProperty("visibility", mitk::BoolProperty::New()); MITK_TEST_CONDITION(baseDataImpl->GetProperty("property38")->GetValueAsString() == "testproperty","Check if base property is set correctly!"); cloneBaseData = baseDataImpl->Clone(); MITK_TEST_CONDITION(cloneBaseData->GetProperty("property38")->GetValueAsString() == "testproperty", "Testing origin set in clone!"); //test method Set-/GetPropertyList mitk::PropertyList::Pointer propertyList = mitk::PropertyList::New(); propertyList->SetFloatProperty("floatProperty1", 123.45); propertyList->SetBoolProperty("visibility",true); propertyList->SetStringProperty("nameXY","propertyName"); baseDataImpl->SetPropertyList(propertyList); bool value = false; MITK_TEST_CONDITION(baseDataImpl->GetPropertyList() == propertyList, "Check if base property list is set correctly!"); MITK_TEST_CONDITION(baseDataImpl->GetPropertyList()->GetBoolProperty("visibility", value) == true, "Check if base property is set correctly in the property list!"); //test method UpdateOutputInformation() baseDataImpl->UpdateOutputInformation(); MITK_TEST_CONDITION(baseDataImpl->GetUpdatedTimeGeometry() == geo2, "TimeGeometry update!"); //Test method CopyInformation() mitk::BaseDataTestImplementation::Pointer newBaseData = mitk::BaseDataTestImplementation::New(); newBaseData->CopyInformation(baseDataImpl); MITK_TEST_CONDITION_REQUIRED( newBaseData->GetTimeGeometry()->CountTimeSteps() == 5, "Check copying of of Basedata Data Object!"); MITK_TEST_END() } diff --git a/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidget.cpp b/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidget.cpp index 03f81028b5..eb10925ff4 100644 --- a/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidget.cpp +++ b/Modules/IGTUI/Qmitk/QmitkInteractiveTransformationWidget.cpp @@ -1,288 +1,289 @@ /*=================================================================== 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 "QmitkInteractiveTransformationWidget.h" // mitk includes #include "mitkRenderingManager.h" #include "mitkNavigationData.h" // vtk includes #include "vtkMatrix4x4.h" #include "vtkLinearTransform.h" const std::string QmitkInteractiveTransformationWidget::VIEW_ID = "org.mitk.views.interactivetransformationwidget"; QmitkInteractiveTransformationWidget::QmitkInteractiveTransformationWidget(QWidget* parent, Qt::WindowFlags f) : QWidget(parent, f), m_Geometry(NULL), m_ResetGeometry(NULL), m_Controls(NULL) { CreateQtPartControl(this); CreateConnections(); m_TranslationVector.Fill(0.0f); m_RotateSliderPos.Fill(0.0f); } QmitkInteractiveTransformationWidget::~QmitkInteractiveTransformationWidget() { } void QmitkInteractiveTransformationWidget::CreateQtPartControl(QWidget *parent) { if (!m_Controls) { // create GUI widgets m_Controls = new Ui::QmitkInteractiveTransformationWidgetControls; m_Controls->setupUi(parent); } } void QmitkInteractiveTransformationWidget::CreateConnections() { if ( m_Controls ) { // translations connect( (QObject*)(m_Controls->m_XTransSlider), SIGNAL(valueChanged(int)), this, SLOT(OnXTranslationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_XTransSpinBox), SIGNAL(valueChanged(int)), this, SLOT(OnXTranslationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_YTransSlider), SIGNAL(valueChanged(int)), this, SLOT(OnYTranslationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_YTransSpinBox), SIGNAL(valueChanged(int)), this, SLOT(OnYTranslationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_ZTransSlider), SIGNAL(valueChanged(int)), this, SLOT(OnZTranslationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_ZTransSpinBox), SIGNAL(valueChanged(int)), this, SLOT(OnZTranslationValueChanged(int)) ); // rotations connect( (QObject*)(m_Controls->m_XRotSlider), SIGNAL(valueChanged(int)), this, SLOT(OnXRotationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_XRotSpinBox), SIGNAL(valueChanged(int)), this, SLOT(OnXRotationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_YRotSlider), SIGNAL(valueChanged(int)), this, SLOT(OnYRotationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_YRotSpinBox), SIGNAL(valueChanged(int)), this, SLOT(OnYRotationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_ZRotSlider), SIGNAL(valueChanged(int)), this, SLOT(OnZRotationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_ZRotSpinBox), SIGNAL(valueChanged(int)), this, SLOT(OnZRotationValueChanged(int)) ); connect( (QObject*)(m_Controls->m_ResetPB), SIGNAL(clicked()), this, SLOT(OnResetGeometry()) ); connect( (QObject*)(m_Controls->m_UseManipulatedToolTipPB), SIGNAL(clicked()), this, SLOT(OnApplyManipulatedToolTip()) ); } } void QmitkInteractiveTransformationWidget::SetGeometry( mitk::BaseGeometry::Pointer geometry, mitk::BaseGeometry::Pointer defaultValues ) { m_Geometry = geometry; - m_ResetGeometry = geometry->Clone(); + itk::LightObject::Pointer lopointer = geometry->Clone(); + m_ResetGeometry = dynamic_cast(lopointer.GetPointer()); //set default values if (defaultValues.IsNotNull()) { //first: some conversion mitk::NavigationData::Pointer transformConversionHelper = mitk::NavigationData::New(defaultValues->GetIndexToWorldTransform()); double eulerAlphaDegrees = transformConversionHelper->GetOrientation().rotation_euler_angles()[0] / vnl_math::pi * 180; double eulerBetaDegrees = transformConversionHelper->GetOrientation().rotation_euler_angles()[1] / vnl_math::pi * 180; double eulerGammaDegrees = transformConversionHelper->GetOrientation().rotation_euler_angles()[2] / vnl_math::pi * 180; //set translation OnXTranslationValueChanged(defaultValues->GetIndexToWorldTransform()->GetOffset()[0]); OnYTranslationValueChanged(defaultValues->GetIndexToWorldTransform()->GetOffset()[1]); OnZTranslationValueChanged(defaultValues->GetIndexToWorldTransform()->GetOffset()[2]); //set rotation OnXRotationValueChanged(eulerAlphaDegrees); OnYRotationValueChanged(eulerBetaDegrees); OnZRotationValueChanged(eulerGammaDegrees); } else { //reset everything OnXTranslationValueChanged(0); OnYTranslationValueChanged(0); OnZTranslationValueChanged(0); OnXRotationValueChanged(0); OnYRotationValueChanged(0); OnZRotationValueChanged(0); } } mitk::BaseGeometry::Pointer QmitkInteractiveTransformationWidget::GetGeometry() { return m_Geometry; } ///////////////////////////////////////////////////////////////////////////////////////////// // Section to allow interactive positioning of the moving surface ///////////////////////////////////////////////////////////////////////////////////////////// void QmitkInteractiveTransformationWidget::OnXTranslationValueChanged( int v ) { mitk::Vector3D translationParams; translationParams[0] = v; translationParams[1] = m_Controls->m_YTransSlider->value(); translationParams[2] = m_Controls->m_ZTransSlider->value(); SetSliderX(v); this->Translate(translationParams); } void QmitkInteractiveTransformationWidget::SetSliderX(int v) { m_Controls->m_XTransSlider->setValue(v); m_Controls->m_XTransSpinBox->setValue(v); } void QmitkInteractiveTransformationWidget::OnYTranslationValueChanged( int v ) { mitk::Vector3D translationParams; translationParams[0] = m_Controls->m_XTransSlider->value(); translationParams[1] = v; translationParams[2] = m_Controls->m_ZTransSlider->value(); SetSliderY(v); this->Translate(translationParams); } void QmitkInteractiveTransformationWidget::SetSliderY(int v) { m_Controls->m_YTransSlider->setValue(v); m_Controls->m_YTransSpinBox->setValue(v); } void QmitkInteractiveTransformationWidget::OnZTranslationValueChanged( int v ) { mitk::Vector3D translationParams; translationParams[0] = m_Controls->m_XTransSlider->value(); translationParams[1] = m_Controls->m_YTransSlider->value(); translationParams[2] = v; SetSliderZ(v); this->Translate(translationParams); } void QmitkInteractiveTransformationWidget::SetSliderZ(int v) { m_Controls->m_ZTransSlider->setValue(v); m_Controls->m_ZTransSpinBox->setValue(v); } void QmitkInteractiveTransformationWidget::Translate( mitk::Vector3D translateVector) { mitk::Vector3D translateVec; // transform the translation vector translateVec[0] = translateVector[0] - m_TranslationVector[0]; translateVec[1] = translateVector[1] - m_TranslationVector[1]; translateVec[2] = translateVector[2] - m_TranslationVector[2]; // set the new translation vector to member variable m_TranslationVector[0] = translateVector[0]; m_TranslationVector[1] = translateVector[1]; m_TranslationVector[2] = translateVector[2]; m_Geometry->Translate( translateVec ); qApp->processEvents(); mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkInteractiveTransformationWidget::OnXRotationValueChanged( int v ) { mitk::Vector3D rotationParams; rotationParams[0] = v; rotationParams[1] = m_Controls->m_YRotSlider->value(); rotationParams[2] = m_Controls->m_ZRotSlider->value(); m_Controls->m_XRotSlider->setValue(v); m_Controls->m_XRotSpinBox->setValue(v); this->Rotate(rotationParams); } void QmitkInteractiveTransformationWidget::OnYRotationValueChanged( int v ) { mitk::Vector3D rotationParams; rotationParams[0] = m_Controls->m_XRotSlider->value(); rotationParams[1] = v; rotationParams[2] = m_Controls->m_ZRotSlider->value(); m_Controls->m_YRotSlider->setValue(v); m_Controls->m_YRotSpinBox->setValue(v); this->Rotate(rotationParams); } void QmitkInteractiveTransformationWidget::OnZRotationValueChanged( int v ) { mitk::Vector3D rotationParams; rotationParams[0]=m_Controls->m_XRotSlider->value(); rotationParams[1]=m_Controls->m_YRotSlider->value(); rotationParams[2]=v; m_Controls->m_ZRotSlider->setValue(v); m_Controls->m_ZRotSpinBox->setValue(v); this->Rotate(rotationParams); } void QmitkInteractiveTransformationWidget::Rotate(mitk::Vector3D rotateVector) { //0: from degrees to radians double radianX = rotateVector[0] * vnl_math::pi / 180; double radianY = rotateVector[1] * vnl_math::pi / 180; double radianZ = rotateVector[2] * vnl_math::pi / 180; //1: from euler angles to quaternion mitk::Quaternion rotation(radianX,radianY,radianZ); //2: Conversion to navigation data / transform mitk::NavigationData::Pointer rotationTransform = mitk::NavigationData::New(); rotationTransform->SetOrientation(rotation); //3: Apply transform //also remember old transform, but without rotation, because rotation is completely stored in the sliders mitk::NavigationData::Pointer oldTransform = mitk::NavigationData::New(m_Geometry->GetIndexToWorldTransform()); mitk::Quaternion identity(0,0,0,1); oldTransform->SetOrientation(identity); //compose old transform with the new one rotationTransform->Compose(oldTransform); //and apply it... m_Geometry->SetIndexToWorldTransform(rotationTransform->GetAffineTransform3D()); qApp->processEvents(); mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkInteractiveTransformationWidget::OnResetGeometry() { m_Controls->m_XRotSlider->setValue(0); m_Controls->m_YRotSlider->setValue(0); m_Controls->m_ZRotSlider->setValue(0); m_Controls->m_XRotSpinBox->setValue(0); m_Controls->m_YRotSpinBox->setValue(0); m_Controls->m_ZRotSpinBox->setValue(0); m_Controls->m_XTransSlider->setValue(0); m_Controls->m_YTransSlider->setValue(0); m_Controls->m_ZTransSlider->setValue(0); m_Controls->m_XTransSpinBox->setValue(0); m_Controls->m_YTransSpinBox->setValue(0); m_Controls->m_ZTransSpinBox->setValue(0); qApp->processEvents(); // reset the input to its initial state. //???xxx??? m_ResetGeometry->TransferItkToVtkTransform(); m_Geometry->SetIdentity(); m_Geometry->Compose(m_ResetGeometry->GetVtkTransform()->GetMatrix()); mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkInteractiveTransformationWidget::OnApplyManipulatedToolTip() { emit ApplyManipulatedToolTip(); } \ No newline at end of file diff --git a/Modules/MitkExt/Algorithms/mitkProbeFilter.cpp b/Modules/MitkExt/Algorithms/mitkProbeFilter.cpp index 8b0ba5a2bf..27f6851664 100644 --- a/Modules/MitkExt/Algorithms/mitkProbeFilter.cpp +++ b/Modules/MitkExt/Algorithms/mitkProbeFilter.cpp @@ -1,213 +1,214 @@ /*=================================================================== 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 "mitkProbeFilter.h" #include "mitkSurface.h" #include "mitkImage.h" #include #include #include #include #include mitk::ProbeFilter::ProbeFilter() { } mitk::ProbeFilter::~ProbeFilter() { } const mitk::Surface *mitk::ProbeFilter::GetInput(void) { if (this->GetNumberOfInputs() < 1) { return 0; } return static_cast< const mitk::Surface * >(this->ProcessObject::GetInput(0) ); } const mitk::Image *mitk::ProbeFilter::GetSource(void) { return static_cast< const mitk::Image * >(this->ProcessObject::GetInput(1)); } void mitk::ProbeFilter::SetInput(const mitk::Surface *input) { this->ProcessObject::SetNthInput( 0, const_cast< mitk::Surface * >( input ) ); } void mitk::ProbeFilter::SetSource(const mitk::Image *source) { this->ProcessObject::SetNthInput( 1, const_cast< mitk::Image * >( source ) ); } void mitk::ProbeFilter::GenerateOutputInformation() { mitk::Surface::ConstPointer input = this->GetInput(); mitk::Image::ConstPointer source = this->GetSource(); mitk::Surface::Pointer output = this->GetOutput(); if(input.IsNull()) return; if(source.IsNull()) return; if(input->GetGeometry()==NULL) return; if(source->GetGeometry()==NULL) return; if( (input->GetTimeGeometry()->CountTimeSteps()==1) && (source->GetTimeGeometry()->CountTimeSteps()>1) ) { - BaseGeometry::Pointer geometry3D = BaseGeometry::New(); + Geometry3D::Pointer geo3D = Geometry3D::New(); + BaseGeometry::Pointer geometry3D = dynamic_cast(geo3D.GetPointer()); geometry3D->Initialize(); geometry3D->SetBounds(source->GetTimeGeometry()->GetBoundsInWorld()); geometry3D->SetTimeBounds(source->GetTimeGeometry()->GetGeometryForTimeStep(0)->GetTimeBounds()); ProportionalTimeGeometry::Pointer outputTimeGeometry = ProportionalTimeGeometry::New(); outputTimeGeometry->Initialize(geometry3D, source->GetTimeGeometry()->CountTimeSteps()); output->Expand(outputTimeGeometry->CountTimeSteps()); output->SetTimeGeometry( outputTimeGeometry ); } else output->SetGeometry( static_cast(input->GetGeometry()->Clone().GetPointer()) ); itkDebugMacro(<<"GenerateOutputInformation()"); } void mitk::ProbeFilter::GenerateData() { mitk::Surface *input = const_cast< mitk::Surface * >(this->GetInput()); mitk::Image *source = const_cast< mitk::Image * >(this->GetSource()); mitk::Surface::Pointer output = this->GetOutput(); itkDebugMacro(<<"Generating Data"); if(output.IsNull()) { itkDebugMacro(<<"Output is NULL."); return; } mitk::Surface::RegionType outputRegion = output->GetRequestedRegion(); const TimeGeometry *outputTimeGeometry = output->GetTimeGeometry(); const TimeGeometry *inputTimeGeometry = input->GetTimeGeometry(); const TimeGeometry *sourceTimeGeometry = source->GetTimeGeometry(); TimePointType timeInMS; int timestep=0; int tstart, tmax; tstart=outputRegion.GetIndex(3); tmax=tstart+outputRegion.GetSize(3); int t; for(t=tstart;tTimeStepToTimePoint( t ); vtkProbeFilter* probe = vtkProbeFilter::New(); timestep = inputTimeGeometry->TimePointToTimeStep( timeInMS ); probe->SetInputData( input->GetVtkPolyData(timestep) ); timestep = sourceTimeGeometry->TimePointToTimeStep( timeInMS ); probe->SetSourceData( source->GetVtkImageData(timestep) ); output->SetVtkPolyData( probe->GetPolyDataOutput(), t ); probe->Update(); probe->Delete(); } } void mitk::ProbeFilter::GenerateInputRequestedRegion() { Superclass::GenerateInputRequestedRegion(); mitk::Surface *input = const_cast< mitk::Surface * >( this->GetInput() ); mitk::Image *source = const_cast< mitk::Image * >( this->GetSource() ); if(input==NULL) return; if(source==NULL) return; mitk::Surface::Pointer output = this->GetOutput(); mitk::Surface::RegionType outputRegion = output->GetRequestedRegion(); const TimeGeometry *outputTimeGeometry = output->GetTimeGeometry(); mitk::Surface::RegionType inputSurfaceRegion = outputRegion; Image::RegionType sourceImageRegion = source->GetLargestPossibleRegion(); if(outputRegion.GetSize(3)<1) { mitk::Surface::RegionType::SizeType surfacesize; surfacesize.Fill(0); inputSurfaceRegion.SetSize(surfacesize); input->SetRequestedRegion( &inputSurfaceRegion ); mitk::Image::RegionType::SizeType imagesize; imagesize.Fill(0); sourceImageRegion.SetSize(imagesize); source->SetRequestedRegion( &sourceImageRegion ); return; } //create and set input requested region for the input surface const TimeGeometry *inputTimeGeometry = input->GetTimeGeometry(); ScalarType timeInMS; int timestep=0; // convert the start-index-time of output in start-index-time of input via millisecond-time timeInMS = outputTimeGeometry->TimeStepToTimePoint(outputRegion.GetIndex(3)); timestep = inputTimeGeometry->TimePointToTimeStep( timeInMS ); if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( inputTimeGeometry->IsValidTimeStep( timestep ) ) ) inputSurfaceRegion.SetIndex( 3, timestep ); else inputSurfaceRegion.SetIndex( 3, 0 ); // convert the end-index-time of output in end-index-time of input via millisecond-time timeInMS = outputTimeGeometry->TimeStepToTimePoint(outputRegion.GetIndex(3)+outputRegion.GetSize(3)-1); timestep = inputTimeGeometry->TimePointToTimeStep( timeInMS ); if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( outputTimeGeometry->IsValidTimeStep( timestep ) ) ) inputSurfaceRegion.SetSize( 3, timestep - inputSurfaceRegion.GetIndex(3) + 1 ); else inputSurfaceRegion.SetSize( 3, 1 ); input->SetRequestedRegion( &inputSurfaceRegion ); //create and set input requested region for the source image const TimeGeometry *sourceTimeGeometry = source->GetTimeGeometry(); // convert the start-index-time of output in start-index-time of source via millisecond-time timeInMS = outputTimeGeometry->TimeStepToTimePoint(outputRegion.GetIndex(3)); timestep = sourceTimeGeometry->TimePointToTimeStep( timeInMS ); if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( sourceTimeGeometry->IsValidTimeStep( timestep ) ) ) sourceImageRegion.SetIndex( 3, timestep ); else sourceImageRegion.SetIndex( 3, 0 ); // convert the end-index-time of output in end-index-time of source via millisecond-time timeInMS = outputTimeGeometry->TimeStepToTimePoint(outputRegion.GetIndex(3)+outputRegion.GetSize(3)-1); timestep = sourceTimeGeometry->TimePointToTimeStep( timeInMS ); if( ( timeInMS > ScalarTypeNumericTraits::NonpositiveMin() ) && ( outputTimeGeometry->IsValidTimeStep( timestep ) ) ) sourceImageRegion.SetSize( 3, timestep - sourceImageRegion.GetIndex(3) + 1 ); else sourceImageRegion.SetSize( 3, 1 ); sourceImageRegion.SetIndex( 4, 0 ); sourceImageRegion.SetSize( 4, 1 ); source->SetRequestedRegion( &sourceImageRegion ); } diff --git a/Modules/MitkExt/Interactions/mitkAffineInteractor3D.cpp b/Modules/MitkExt/Interactions/mitkAffineInteractor3D.cpp index 5e8478cb1c..70aa724730 100644 --- a/Modules/MitkExt/Interactions/mitkAffineInteractor3D.cpp +++ b/Modules/MitkExt/Interactions/mitkAffineInteractor3D.cpp @@ -1,488 +1,489 @@ /*=================================================================== 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 "mitkAffineInteractor3D.h" #include "mitkPointOperation.h" #include "mitkPositionEvent.h" #include "mitkStatusBar.h" #include "mitkDataNode.h" #include "mitkInteractionConst.h" #include "mitkAction.h" #include "mitkStateEvent.h" #include "mitkOperationEvent.h" #include "mitkUndoController.h" #include "mitkStateMachineFactory.h" #include "mitkStateTransitionOperation.h" #include "mitkBaseRenderer.h" #include "mitkRenderingManager.h" #include "mitkRotationOperation.h" #include #include #include #include #include #include #include #include namespace mitk { //how precise must the user pick the point //default value AffineInteractor3D ::AffineInteractor3D(const char * type, DataNode* dataNode, int /* n */ ) : Interactor( type, dataNode ), m_Precision( 6.5 ), m_InteractionMode( INTERACTION_MODE_TRANSLATION ) { - m_OriginalGeometry = BaseGeometry::New(); + Geometry3D::Pointer geo3D = Geometry3D::New(); + m_OriginalGeometry = dynamic_cast(geo3D.GetPointer()); // Initialize vector arithmetic m_ObjectNormal[0] = 0.0; m_ObjectNormal[1] = 0.0; m_ObjectNormal[2] = 1.0; } AffineInteractor3D::~AffineInteractor3D() { } void AffineInteractor3D::SetInteractionMode( unsigned int interactionMode ) { m_InteractionMode = interactionMode; } void AffineInteractor3D::SetInteractionModeToTranslation() { m_InteractionMode = INTERACTION_MODE_TRANSLATION; } void AffineInteractor3D::SetInteractionModeToRotation() { m_InteractionMode = INTERACTION_MODE_ROTATION; } unsigned int AffineInteractor3D::GetInteractionMode() const { return m_InteractionMode; } void AffineInteractor3D::SetPrecision( ScalarType precision ) { m_Precision = precision; } // Overwritten since this class can handle it better! float AffineInteractor3D ::CanHandleEvent(StateEvent const* stateEvent) const { float returnValue = 0.5; // If it is a key event that can be handled in the current state, // then return 0.5 DisplayPositionEvent const *disPosEvent = dynamic_cast (stateEvent->GetEvent()); // Key event handling: if (disPosEvent == NULL) { // Check if the current state has a transition waiting for that key event. if (this->GetCurrentState()->GetTransition(stateEvent->GetId())!=NULL) { return 0.5; } else { return 0.0; } } //on MouseMove do nothing! //if (stateEvent->GetEvent()->GetType() == Type_MouseMove) //{ // return 0.0; //} //if the event can be understood and if there is a transition waiting for that event if (this->GetCurrentState()->GetTransition(stateEvent->GetId())!=NULL) { returnValue = 0.5;//it can be understood } //int timeStep = disPosEvent->GetSender()->GetTimeStep(); //CurveModel *curveModel = dynamic_cast( // m_DataNode->GetData() ); //if ( curveModel != NULL ) //{ // // Get the PlaneGeometry of the window the user interacts with (for 2D point // // projection) // BaseRenderer *renderer = stateEvent->GetEvent()->GetSender(); // const PlaneGeometry *projectionPlane = renderer->GetCurrentWorldGeometry2D(); // // For reading on the points, Ids etc // //CurveModel::PointSetType *pointSet = curveModel->GetPointSet( timeStep ); // //if ( pointSet == NULL ) // //{ // // return 0.0; // //} //} return returnValue; } bool AffineInteractor3D ::ExecuteAction( Action *action, StateEvent const *stateEvent ) { bool ok = false; // Get data object BaseData *data = m_DataNode->GetData(); if ( data == NULL ) { MITK_ERROR << "No data object present!"; return ok; } // Get Event and extract renderer const Event *event = stateEvent->GetEvent(); BaseRenderer *renderer = NULL; vtkRenderWindow *renderWindow = NULL; vtkRenderWindowInteractor *renderWindowInteractor = NULL; vtkRenderer *currentVtkRenderer = NULL; vtkCamera *camera = NULL; if ( event != NULL ) { renderer = event->GetSender(); if ( renderer != NULL ) { renderWindow = renderer->GetRenderWindow(); if ( renderWindow != NULL ) { renderWindowInteractor = renderWindow->GetInteractor(); if ( renderWindowInteractor != NULL ) { currentVtkRenderer = renderWindowInteractor ->GetInteractorStyle()->GetCurrentRenderer(); if ( currentVtkRenderer != NULL ) { camera = currentVtkRenderer->GetActiveCamera(); } } } } } // Check if we have a DisplayPositionEvent const DisplayPositionEvent *dpe = dynamic_cast< const DisplayPositionEvent * >( stateEvent->GetEvent() ); if ( dpe != NULL ) { m_CurrentPickedPoint = dpe->GetWorldPosition(); m_CurrentPickedDisplayPoint = dpe->GetDisplayPosition(); } // Get the timestep to also support 3D+t int timeStep = 0; ScalarType timeInMS = 0.0; if ( renderer != NULL ) { timeStep = renderer->GetTimeStep( data ); timeInMS = renderer->GetTime(); } // If data is an mitk::Surface, extract it Surface *surface = dynamic_cast< Surface * >( data ); vtkPolyData *polyData = NULL; if ( surface != NULL ) { polyData = surface->GetVtkPolyData( timeStep ); // Extract surface normal from surface (if existent, otherwise use default) vtkPointData *pointData = polyData->GetPointData(); if ( pointData != NULL ) { vtkDataArray *normal = polyData->GetPointData()->GetVectors( "planeNormal" ); if ( normal != NULL ) { m_ObjectNormal[0] = normal->GetComponent( 0, 0 ); m_ObjectNormal[1] = normal->GetComponent( 0, 1 ); m_ObjectNormal[2] = normal->GetComponent( 0, 2 ); } } } // Get geometry object m_Geometry = data->GetGeometry( timeStep ); // Make sure that the data (if time-resolved) has enough entries; // if not, create the required extra ones (empty) data->Expand( timeStep+1 ); switch (action->GetActionId()) { case AcDONOTHING: ok = true; break; case AcCHECKOBJECT: { // Re-enable VTK interactor (may have been disabled previously) if ( renderWindowInteractor != NULL ) { renderWindowInteractor->Enable(); } // Check if we have a DisplayPositionEvent const DisplayPositionEvent *dpe = dynamic_cast< const DisplayPositionEvent * >( stateEvent->GetEvent() ); if ( dpe == NULL ) { ok = true; break; } // Check if an object is present at the current mouse position DataNode *pickedNode = dpe->GetPickedObjectNode(); StateEvent *newStateEvent; if ( pickedNode == m_DataNode ) { // Yes: object will be selected newStateEvent = new StateEvent( EIDYES ); } else { // No: back to start state newStateEvent = new StateEvent( EIDNO ); } this->HandleEvent( newStateEvent ); ok = true; break; } case AcDESELECTOBJECT: { // Color object white m_DataNode->SetColor( 1.0, 1.0, 1.0 ); RenderingManager::GetInstance()->RequestUpdateAll(); // Colorize surface / wireframe as inactive this->ColorizeSurface( polyData, m_CurrentPickedPoint, -1.0 ); ok = true; break; } case AcSELECTPICKEDOBJECT: { // Color object red m_DataNode->SetColor( 1.0, 0.0, 0.0 ); RenderingManager::GetInstance()->RequestUpdateAll(); // Colorize surface / wireframe dependend on distance from picked point this->ColorizeSurface( polyData, m_CurrentPickedPoint, 0.0 ); ok = true; break; } case AcINITMOVE: { // Disable VTK interactor until MITK interaction has been completed if ( renderWindowInteractor != NULL ) { renderWindowInteractor->Disable(); } // Check if we have a DisplayPositionEvent const DisplayPositionEvent *dpe = dynamic_cast< const DisplayPositionEvent * >( stateEvent->GetEvent() ); if ( dpe == NULL ) { ok = true; break; } //DataNode *pickedNode = dpe->GetPickedObjectNode(); m_InitialPickedPoint = m_CurrentPickedPoint; m_InitialPickedDisplayPoint = m_CurrentPickedDisplayPoint; if ( currentVtkRenderer != NULL ) { vtkInteractorObserver::ComputeDisplayToWorld( currentVtkRenderer, m_InitialPickedDisplayPoint[0], m_InitialPickedDisplayPoint[1], 0.0, //m_InitialInteractionPickedPoint[2], m_InitialPickedPointWorld ); } // Make deep copy of current Geometry3D of the plane data->UpdateOutputInformation(); // make sure that the Geometry is up-to-date m_OriginalGeometry = static_cast< BaseGeometry * >( data->GetGeometry( timeStep )->Clone().GetPointer() ); ok = true; break; } case AcMOVE: { // Check if we have a DisplayPositionEvent const DisplayPositionEvent *dpe = dynamic_cast< const DisplayPositionEvent * >( stateEvent->GetEvent() ); if ( dpe == NULL ) { ok = true; break; } if ( currentVtkRenderer != NULL ) { vtkInteractorObserver::ComputeDisplayToWorld( currentVtkRenderer, m_CurrentPickedDisplayPoint[0], m_CurrentPickedDisplayPoint[1], 0.0, //m_InitialInteractionPickedPoint[2], m_CurrentPickedPointWorld ); } Vector3D interactionMove; interactionMove[0] = m_CurrentPickedPointWorld[0] - m_InitialPickedPointWorld[0]; interactionMove[1] = m_CurrentPickedPointWorld[1] - m_InitialPickedPointWorld[1]; interactionMove[2] = m_CurrentPickedPointWorld[2] - m_InitialPickedPointWorld[2]; if ( m_InteractionMode == INTERACTION_MODE_TRANSLATION ) { Point3D origin = m_OriginalGeometry->GetOrigin(); Vector3D transformedObjectNormal; data->GetGeometry( timeStep )->IndexToWorld( m_ObjectNormal, transformedObjectNormal ); data->GetGeometry( timeStep )->SetOrigin( origin + transformedObjectNormal * (interactionMove * transformedObjectNormal) ); } else if ( m_InteractionMode == INTERACTION_MODE_ROTATION ) { if ( camera ) { double vpn[3]; camera->GetViewPlaneNormal( vpn ); Vector3D viewPlaneNormal; viewPlaneNormal[0] = vpn[0]; viewPlaneNormal[1] = vpn[1]; viewPlaneNormal[2] = vpn[2]; Vector3D rotationAxis = itk::CrossProduct( viewPlaneNormal, interactionMove ); rotationAxis.Normalize(); int *size = currentVtkRenderer->GetSize(); double l2 = (m_CurrentPickedDisplayPoint[0] - m_InitialPickedDisplayPoint[0]) * (m_CurrentPickedDisplayPoint[0] - m_InitialPickedDisplayPoint[0]) + (m_CurrentPickedDisplayPoint[1] - m_InitialPickedDisplayPoint[1]) * (m_CurrentPickedDisplayPoint[1] - m_InitialPickedDisplayPoint[1]); double rotationAngle = 360.0 * sqrt(l2/(size[0]*size[0]+size[1]*size[1])); // Use center of data bounding box as center of rotation Point3D rotationCenter = m_OriginalGeometry->GetCenter();; // Reset current Geometry3D to original state (pre-interaction) and // apply rotation RotationOperation op( OpROTATE, rotationCenter, rotationAxis, rotationAngle ); BaseGeometry::Pointer newGeometry = static_cast< BaseGeometry * >( m_OriginalGeometry->Clone().GetPointer() ); newGeometry->ExecuteOperation( &op ); data->SetClonedGeometry(newGeometry, timeStep); } } RenderingManager::GetInstance()->RequestUpdateAll(); ok = true; break; } default: return Superclass::ExecuteAction( action, stateEvent ); } return ok; } bool AffineInteractor3D::ColorizeSurface( vtkPolyData *polyData, const Point3D & /*pickedPoint*/, double scalar ) { if ( polyData == NULL ) { return false; } //vtkPoints *points = polyData->GetPoints(); vtkPointData *pointData = polyData->GetPointData(); if ( pointData == NULL ) { return false; } vtkDataArray *scalars = pointData->GetScalars(); if ( scalars == NULL ) { return false; } for ( unsigned int i = 0; i < pointData->GetNumberOfTuples(); ++i ) { scalars->SetComponent( i, 0, scalar ); } polyData->Modified(); pointData->Update(); return true; } } // namespace diff --git a/Modules/MitkExt/Testing/mitkPlaneFitTest.cpp b/Modules/MitkExt/Testing/mitkPlaneFitTest.cpp index be5dc4a2a6..2b038c967c 100644 --- a/Modules/MitkExt/Testing/mitkPlaneFitTest.cpp +++ b/Modules/MitkExt/Testing/mitkPlaneFitTest.cpp @@ -1,97 +1,98 @@ /*=================================================================== 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 #include #include #include #include #include int mitkPlaneFitTest(int, char*[] ) { //float bounds[]={0.0f,10.0f,0.0f,10.0f,0.0f,5.0f}; mitk::PlaneFit::Pointer PlaneFit = mitk::PlaneFit::New(); mitk::PointSet::Pointer PointSet = mitk::PointSet::New(); - mitk::BaseGeometry::Pointer BaseGeometry = mitk::BaseGeometry::New(); + mitk::PlaneGeometry::Pointer planeGeo= mitk::PlaneGeometry::New(); + mitk::BaseGeometry::Pointer BaseGeometry = dynamic_cast(planeGeo.GetPointer()); mitk::Point3D Point; //first without any point, then incrementally add points within thre points there will be a plane geometry std::cout <<"Start PlaneFitTest "<GetPointSet()->GetPoints()->InsertElement(position, Point); } //Set Input PlaneFit->SetInput(PointSet); const mitk::PointSet* testPointSet = PlaneFit->GetInput(); std::cout<<" Size test of Input Method: "; if( testPointSet->GetSize() == PointSet->GetSize() ) { std::cout<<"[PASSED]"<Update(); const mitk::Point3D ¢roid = PlaneFit->GetCentroid(); mitk::Point3D expectedCentroid; expectedCentroid[0]=2.5; expectedCentroid[1]=3.75; expectedCentroid[2]=2.5; if ( centroid == expectedCentroid ) { std::cout<<"[PASSED]"<( PlaneFit->GetOutput()->GetGeometry()); if( PlaneGeometry ) { std::cout<<"[PASSED]"< mitk::DiffSliceOperation::DiffSliceOperation():Operation(1) { m_TimeStep = 0; m_Slice = NULL; m_Image = NULL; m_WorldGeometry = NULL; m_SliceGeometry = NULL; m_ImageIsValid = false; } mitk::DiffSliceOperation::DiffSliceOperation(mitk::Image* imageVolume, vtkImageData* slice, - BaseGeometry* sliceGeometry, + SlicedGeometry3D* sliceGeometry, unsigned int timestep, BaseGeometry* currentWorldGeometry):Operation(1) { - m_WorldGeometry = currentWorldGeometry->Clone(); + itk::LightObject::Pointer lopointer = currentWorldGeometry->Clone(); + m_WorldGeometry = dynamic_cast(lopointer.GetPointer()); /* Quick fix for bug 12338. Guard object - fix this when clone method of PlaneGeometry is cloning the reference geometry (see bug 13392)*/ - m_GuardReferenceGeometry = mitk::BaseGeometry::New(); + //xxxx m_GuardReferenceGeometry = mitk::BaseGeometry::New(); m_GuardReferenceGeometry = dynamic_cast(m_WorldGeometry.GetPointer())->GetReferenceGeometry(); /*---------------------------------------------------------------------------------------------------*/ m_SliceGeometry = sliceGeometry->Clone(); m_TimeStep = timestep; /*m_zlibSliceContainer = CompressedImageContainer::New(); m_zlibSliceContainer->SetImage( slice );*/ m_Slice = vtkSmartPointer::New(); m_Slice->DeepCopy(slice); m_Image = imageVolume; if ( m_Image) { /*add an observer to listen to the delete event of the image, this is necessary because the operation is then invalid*/ itk::SimpleMemberCommand< DiffSliceOperation >::Pointer command = itk::SimpleMemberCommand< DiffSliceOperation >::New(); command->SetCallbackFunction( this, &DiffSliceOperation::OnImageDeleted ); //get the id of the observer, used to remove it later on m_DeleteObserverTag = imageVolume->AddObserver( itk::DeleteEvent(), command ); m_ImageIsValid = true; } else m_ImageIsValid = false; } mitk::DiffSliceOperation::~DiffSliceOperation() { m_Slice = NULL; m_WorldGeometry = NULL; //m_zlibSliceContainer = NULL; if (m_ImageIsValid) { //if the image is still there, we have to remove the observer from it m_Image->RemoveObserver( m_DeleteObserverTag ); } m_Image = NULL; } vtkImageData* mitk::DiffSliceOperation::GetSlice() { //Image::ConstPointer image = m_zlibSliceContainer->GetImage().GetPointer(); return m_Slice; } bool mitk::DiffSliceOperation::IsValid() { return m_ImageIsValid && (m_Slice.GetPointer() != NULL) && (m_WorldGeometry.IsNotNull());//TODO improve } void mitk::DiffSliceOperation::OnImageDeleted() { //if our imageVolume is removed e.g. from the datastorage the operation is no lnger valid m_ImageIsValid = false; } \ No newline at end of file diff --git a/Modules/Segmentation/Algorithms/mitkDiffSliceOperation.h b/Modules/Segmentation/Algorithms/mitkDiffSliceOperation.h index c3f01ef5c7..1a8bcfd28a 100644 --- a/Modules/Segmentation/Algorithms/mitkDiffSliceOperation.h +++ b/Modules/Segmentation/Algorithms/mitkDiffSliceOperation.h @@ -1,139 +1,139 @@ /*=================================================================== 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. ===================================================================*/ #ifndef mitkDiffSliceOperation_h_Included #define mitkDiffSliceOperation_h_Included #include "SegmentationExports.h" #include "mitkCommon.h" #include //#include "mitkCompressedImageContainer.h" #include #include #include //DEPRECATED #include namespace mitk { /** \brief An Operation for applying an edited slice to the volume. \sa DiffSliceOperationApplier The information for the operation is specified by properties: imageVolume the volume where the slice was extracted from. slice the slice to be applied. timestep the timestep in an 4D image. currentWorldGeometry specifies the axis where the slice has to be applied in the volume. This Operation can be used to realize undo-redo functionality for e.g. segmentation purposes. */ class Segmentation_EXPORT DiffSliceOperation : public Operation { public: mitkClassMacro(DiffSliceOperation, OperationActor); //itkNewMacro(DiffSliceOperation); //mitkNewMacro4Param(DiffSliceOperation,mitk::Image,mitk::Image,unsigned int, mitk::PlaneGeometry); /** \brief Creates an empty instance. Note that it is not valid yet. The properties of the object have to be set. */ DiffSliceOperation(); /** \brief */ - DiffSliceOperation( mitk::Image* imageVolume, vtkImageData* slice, BaseGeometry* sliceGeometry, unsigned int timestep, BaseGeometry* currentWorldGeometry); + DiffSliceOperation( mitk::Image* imageVolume, vtkImageData* slice, SlicedGeometry3D* sliceGeometry, unsigned int timestep, BaseGeometry* currentWorldGeometry); /** \brief * * \deprecatedSince{2013_09} Please use TimeGeometry instead of TimeSlicedGeometry. For more information see http://www.mitk.org/Development/Refactoring%20of%20the%20Geometry%20Classes%20-%20Part%201 */ DEPRECATED(DiffSliceOperation( mitk::Image* imageVolume, vtkImageData* slice, TimeSlicedGeometry* sliceGeometry, unsigned int timestep, BaseGeometry* currentWorldGeometry)); /** \brief Check if it is a valid operation.*/ bool IsValid(); /** \brief Set the image volume.*/ void SetImage(mitk::Image* image){ this->m_Image = image;} /** \brief Get th image volume.*/ mitk::Image* GetImage(){return this->m_Image;} /** \brief Set thee slice to be applied.*/ void SetImage(vtkImageData* slice){ this->m_Slice = slice;} /** \brief Get the slice that is applied in the operation.*/ vtkImageData* GetSlice(); /** \brief Get timeStep.*/ void SetTimeStep(unsigned int timestep){this->m_TimeStep = timestep;} /** \brief Set timeStep*/ unsigned int GetTimeStep(){return this->m_TimeStep;} /** \brief Set the axis where the slice has to be applied in the volume.*/ - void SetSliceGeometry(BaseGeometry* sliceGeometry){this->m_SliceGeometry = sliceGeometry;} + void SetSliceGeometry(SlicedGeometry3D* sliceGeometry){this->m_SliceGeometry = sliceGeometry;} /** \brief Get the axis where the slice has to be applied in the volume.*/ - BaseGeometry* GetSliceGeometry(){return this->m_SliceGeometry;} + SlicedGeometry3D* GetSliceGeometry(){return this->m_SliceGeometry;} /** \brief Set the axis where the slice has to be applied in the volume. * \deprecatedSince{2013_09} Please use TimeGeometry instead of TimeSlicedGeometry. For more information see http://www.mitk.org/Development/Refactoring%20of%20the%20Geometry%20Classes%20-%20Part%201 */ void SetSliceGeometry(TimeSlicedGeometry* sliceGeometry); /** \brief Set the axis where the slice has to be applied in the volume.*/ void SetCurrentWorldGeometry(BaseGeometry* worldGeometry){this->m_WorldGeometry = worldGeometry;} /** \brief Get the axis where the slice has to be applied in the volume.*/ BaseGeometry* GetWorldGeometry(){return this->m_WorldGeometry;} /** \brief Set the axis where the slice has to be applied in the volume. * \deprecatedSince{2013_09} Please use TimeGeometry instead of TimeSlicedGeometry. For more information see http://www.mitk.org/Development/Refactoring%20of%20the%20Geometry%20Classes%20-%20Part%201 */ void SetCurrentWorldGeometry(TimeSlicedGeometry* worldGeometry); protected: virtual ~DiffSliceOperation(); /** \brief Callback for image observer.*/ void OnImageDeleted(); //CompressedImageContainer::Pointer m_zlibSliceContainer; mitk::Image* m_Image; vtkSmartPointer m_Slice; - BaseGeometry::Pointer m_SliceGeometry; + SlicedGeometry3D::Pointer m_SliceGeometry; unsigned int m_TimeStep; BaseGeometry::Pointer m_WorldGeometry; bool m_ImageIsValid; unsigned long m_DeleteObserverTag; mitk::BaseGeometry::Pointer m_GuardReferenceGeometry; }; } #endif diff --git a/Modules/Segmentation/Interactions/mitkSegTool2D.cpp b/Modules/Segmentation/Interactions/mitkSegTool2D.cpp index 4e8dbed00e..8c0dc1af2d 100644 --- a/Modules/Segmentation/Interactions/mitkSegTool2D.cpp +++ b/Modules/Segmentation/Interactions/mitkSegTool2D.cpp @@ -1,402 +1,402 @@ /*=================================================================== 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 "mitkSegTool2D.h" #include "mitkToolManager.h" #include "mitkDataStorage.h" #include "mitkBaseRenderer.h" #include "mitkPlaneGeometry.h" #include "mitkExtractImageFilter.h" #include "mitkExtractDirectedPlaneImageFilter.h" //Include of the new ImageExtractor #include "mitkExtractDirectedPlaneImageFilterNew.h" #include "mitkPlanarCircle.h" #include "mitkOverwriteSliceImageFilter.h" #include "mitkOverwriteDirectedPlaneImageFilter.h" #include "usGetModuleContext.h" //Includes for 3DSurfaceInterpolation #include "mitkImageToContourFilter.h" #include "mitkSurfaceInterpolationController.h" //includes for resling and overwriting #include #include #include #include #include #include "mitkOperationEvent.h" #include "mitkUndoController.h" #define ROUND(a) ((a)>0 ? (int)((a)+0.5) : -(int)(0.5-(a))) mitk::SegTool2D::SegTool2D(const char* type) :Tool(type), m_LastEventSender(NULL), m_LastEventSlice(0), m_Contourmarkername ("Position"), m_ShowMarkerNodes (false), m_3DInterpolationEnabled(true) { } mitk::SegTool2D::~SegTool2D() { } float mitk::SegTool2D::CanHandleEvent( StateEvent const *stateEvent) const { const PositionEvent* positionEvent = dynamic_cast(stateEvent->GetEvent()); if (!positionEvent) return 0.0; if ( positionEvent->GetSender()->GetMapperID() != BaseRenderer::Standard2D ) return 0.0; // we don't want anything but 2D //This are the mouse event that are used by the statemachine patterns for zooming and panning. This must be possible although a tool is activ if (stateEvent->GetId() == EIDRIGHTMOUSEBTN || stateEvent->GetId() == EIDMIDDLEMOUSEBTN || stateEvent->GetId() == EIDRIGHTMOUSEBTNANDCTRL || stateEvent->GetId() == EIDMIDDLEMOUSERELEASE || stateEvent->GetId() == EIDRIGHTMOUSERELEASE || stateEvent->GetId() == EIDRIGHTMOUSEBTNANDMOUSEMOVE || stateEvent->GetId() == EIDMIDDLEMOUSEBTNANDMOUSEMOVE || stateEvent->GetId() == EIDCTRLANDRIGHTMOUSEBTNANDMOUSEMOVE || stateEvent->GetId() == EIDCTRLANDRIGHTMOUSEBTNRELEASE ) { //Since the usual segmentation tools currently do not need right click interaction but the mitkDisplayVectorInteractor return 0.0; } else { return 1.0; } } bool mitk::SegTool2D::DetermineAffectedImageSlice( const Image* image, const PlaneGeometry* plane, int& affectedDimension, int& affectedSlice ) { assert(image); assert(plane); // compare normal of plane to the three axis vectors of the image Vector3D normal = plane->GetNormal(); Vector3D imageNormal0 = image->GetSlicedGeometry()->GetAxisVector(0); Vector3D imageNormal1 = image->GetSlicedGeometry()->GetAxisVector(1); Vector3D imageNormal2 = image->GetSlicedGeometry()->GetAxisVector(2); normal.Normalize(); imageNormal0.Normalize(); imageNormal1.Normalize(); imageNormal2.Normalize(); imageNormal0.SetVnlVector( vnl_cross_3d(normal.GetVnlVector(),imageNormal0.GetVnlVector()) ); imageNormal1.SetVnlVector( vnl_cross_3d(normal.GetVnlVector(),imageNormal1.GetVnlVector()) ); imageNormal2.SetVnlVector( vnl_cross_3d(normal.GetVnlVector(),imageNormal2.GetVnlVector()) ); double eps( 0.00001 ); // axial if ( imageNormal2.GetNorm() <= eps ) { affectedDimension = 2; } // sagittal else if ( imageNormal1.GetNorm() <= eps ) { affectedDimension = 1; } // frontal else if ( imageNormal0.GetNorm() <= eps ) { affectedDimension = 0; } else { affectedDimension = -1; // no idea return false; } // determine slice number in image BaseGeometry* imageGeometry = image->GetGeometry(0); Point3D testPoint = imageGeometry->GetCenter(); Point3D projectedPoint; plane->Project( testPoint, projectedPoint ); Point3D indexPoint; imageGeometry->WorldToIndex( projectedPoint, indexPoint ); affectedSlice = ROUND( indexPoint[affectedDimension] ); MITK_DEBUG << "indexPoint " << indexPoint << " affectedDimension " << affectedDimension << " affectedSlice " << affectedSlice; // check if this index is still within the image if ( affectedSlice < 0 || affectedSlice >= static_cast(image->GetDimension(affectedDimension)) ) return false; return true; } mitk::Image::Pointer mitk::SegTool2D::GetAffectedImageSliceAs2DImage(const PositionEvent* positionEvent, const Image* image) { if (!positionEvent) return NULL; assert( positionEvent->GetSender() ); // sure, right? unsigned int timeStep = positionEvent->GetSender()->GetTimeStep( image ); // get the timestep of the visible part (time-wise) of the image // first, we determine, which slice is affected const PlaneGeometry* planeGeometry( dynamic_cast (positionEvent->GetSender()->GetCurrentWorldGeometry2D() ) ); return this->GetAffectedImageSliceAs2DImage(planeGeometry, image, timeStep); } mitk::Image::Pointer mitk::SegTool2D::GetAffectedImageSliceAs2DImage(const PlaneGeometry* planeGeometry, const Image* image, unsigned int timeStep) { if ( !image || !planeGeometry ) return NULL; //Make sure that for reslicing and overwriting the same alogrithm is used. We can specify the mode of the vtk reslicer vtkSmartPointer reslice = vtkSmartPointer::New(); //set to false to extract a slice reslice->SetOverwriteMode(false); reslice->Modified(); //use ExtractSliceFilter with our specific vtkImageReslice for overwriting and extracting mitk::ExtractSliceFilter::Pointer extractor = mitk::ExtractSliceFilter::New(reslice); extractor->SetInput( image ); extractor->SetTimeStep( timeStep ); extractor->SetWorldGeometry( planeGeometry ); extractor->SetVtkOutputRequest(false); extractor->SetResliceTransformByGeometry( image->GetTimeGeometry()->GetGeometryForTimeStep( timeStep ) ); extractor->Modified(); extractor->Update(); Image::Pointer slice = extractor->GetOutput(); /*============= BEGIN undo feature block ========================*/ //specify the undo operation with the non edited slice - m_undoOperation = new DiffSliceOperation(const_cast(image), extractor->GetVtkOutput(), slice->GetGeometry(), timeStep, const_cast(planeGeometry)); + m_undoOperation = new DiffSliceOperation(const_cast(image), extractor->GetVtkOutput(), dynamic_cast(slice->GetGeometry()), timeStep, const_cast(planeGeometry)); /*============= END undo feature block ========================*/ return slice; } mitk::Image::Pointer mitk::SegTool2D::GetAffectedWorkingSlice(const PositionEvent* positionEvent) { DataNode* workingNode( m_ToolManager->GetWorkingData(0) ); if ( !workingNode ) return NULL; Image* workingImage = dynamic_cast(workingNode->GetData()); if ( !workingImage ) return NULL; return GetAffectedImageSliceAs2DImage( positionEvent, workingImage ); } mitk::Image::Pointer mitk::SegTool2D::GetAffectedReferenceSlice(const PositionEvent* positionEvent) { DataNode* referenceNode( m_ToolManager->GetReferenceData(0) ); if ( !referenceNode ) return NULL; Image* referenceImage = dynamic_cast(referenceNode->GetData()); if ( !referenceImage ) return NULL; return GetAffectedImageSliceAs2DImage( positionEvent, referenceImage ); } void mitk::SegTool2D::WriteBackSegmentationResult (const PositionEvent* positionEvent, Image* slice) { if(!positionEvent) return; const PlaneGeometry* planeGeometry( dynamic_cast (positionEvent->GetSender()->GetCurrentWorldGeometry2D() ) ); if( planeGeometry && slice) { DataNode* workingNode( m_ToolManager->GetWorkingData(0) ); Image* image = dynamic_cast(workingNode->GetData()); unsigned int timeStep = positionEvent->GetSender()->GetTimeStep( image ); this->WriteBackSegmentationResult(planeGeometry, slice, timeStep); slice->DisconnectPipeline(); ImageToContourFilter::Pointer contourExtractor = ImageToContourFilter::New(); contourExtractor->SetInput(slice); contourExtractor->Update(); mitk::Surface::Pointer contour = contourExtractor->GetOutput(); if (m_3DInterpolationEnabled && contour->GetVtkPolyData()->GetNumberOfPoints() > 0 && image->GetDimension() == 3) { unsigned int pos = this->AddContourmarker(positionEvent); us::ServiceReference serviceRef = us::GetModuleContext()->GetServiceReference(); PlanePositionManagerService* service = us::GetModuleContext()->GetService(serviceRef); mitk::SurfaceInterpolationController::GetInstance()->AddNewContour( contour, service->GetPlanePosition(pos)); contour->DisconnectPipeline(); } } } void mitk::SegTool2D::WriteBackSegmentationResult (const PlaneGeometry* planeGeometry, Image* slice, unsigned int timeStep) { if(!planeGeometry || !slice) return; DataNode* workingNode( m_ToolManager->GetWorkingData(0) ); Image* image = dynamic_cast(workingNode->GetData()); //Make sure that for reslicing and overwriting the same alogrithm is used. We can specify the mode of the vtk reslicer vtkSmartPointer reslice = vtkSmartPointer::New(); //Set the slice as 'input' reslice->SetInputSlice(slice->GetVtkImageData()); //set overwrite mode to true to write back to the image volume reslice->SetOverwriteMode(true); reslice->Modified(); mitk::ExtractSliceFilter::Pointer extractor = mitk::ExtractSliceFilter::New(reslice); extractor->SetInput( image ); extractor->SetTimeStep( timeStep ); extractor->SetWorldGeometry( planeGeometry ); extractor->SetVtkOutputRequest(true); extractor->SetResliceTransformByGeometry( image->GetGeometry( timeStep ) ); extractor->Modified(); extractor->Update(); //the image was modified within the pipeline, but not marked so image->Modified(); image->GetVtkImageData()->Modified(); /*============= BEGIN undo feature block ========================*/ //specify the undo operation with the edited slice - m_doOperation = new DiffSliceOperation(image, extractor->GetVtkOutput(),slice->GetGeometry(), timeStep, const_cast(planeGeometry)); + m_doOperation = new DiffSliceOperation(image, extractor->GetVtkOutput(),dynamic_cast(slice->GetGeometry()), timeStep, const_cast(planeGeometry)); //create an operation event for the undo stack OperationEvent* undoStackItem = new OperationEvent( DiffSliceOperationApplier::GetInstance(), m_doOperation, m_undoOperation, "Segmentation" ); //add it to the undo controller UndoController::GetCurrentUndoModel()->SetOperationEvent( undoStackItem ); //clear the pointers as the operation are stored in the undocontroller and also deleted from there m_undoOperation = NULL; m_doOperation = NULL; /*============= END undo feature block ========================*/ mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void mitk::SegTool2D::SetShowMarkerNodes(bool status) { m_ShowMarkerNodes = status; } void mitk::SegTool2D::SetEnable3DInterpolation(bool enabled) { m_3DInterpolationEnabled = enabled; } unsigned int mitk::SegTool2D::AddContourmarker ( const PositionEvent* positionEvent ) { const mitk::PlaneGeometry* plane = dynamic_cast (dynamic_cast< const mitk::SlicedGeometry3D*>( positionEvent->GetSender()->GetSliceNavigationController()->GetCurrentGeometry3D())->GetGeometry2D(0)); us::ServiceReference serviceRef = us::GetModuleContext()->GetServiceReference(); PlanePositionManagerService* service = us::GetModuleContext()->GetService(serviceRef); unsigned int size = service->GetNumberOfPlanePositions(); unsigned int id = service->AddNewPlanePosition(plane, positionEvent->GetSender()->GetSliceNavigationController()->GetSlice()->GetPos()); mitk::PlanarCircle::Pointer contourMarker = mitk::PlanarCircle::New(); mitk::Point2D p1; plane->Map(plane->GetCenter(), p1); mitk::Point2D p2 = p1; p2[0] -= plane->GetSpacing()[0]; p2[1] -= plane->GetSpacing()[1]; contourMarker->PlaceFigure( p1 ); contourMarker->SetCurrentControlPoint( p1 ); contourMarker->SetGeometry2D( const_cast(plane)); std::stringstream markerStream; mitk::DataNode* workingNode (m_ToolManager->GetWorkingData(0)); markerStream << m_Contourmarkername ; markerStream << " "; markerStream << id+1; DataNode::Pointer rotatedContourNode = DataNode::New(); rotatedContourNode->SetData(contourMarker); rotatedContourNode->SetProperty( "name", StringProperty::New(markerStream.str()) ); rotatedContourNode->SetProperty( "isContourMarker", BoolProperty::New(true)); rotatedContourNode->SetBoolProperty( "PlanarFigureInitializedWindow", true, positionEvent->GetSender() ); rotatedContourNode->SetProperty( "includeInBoundingBox", BoolProperty::New(false)); rotatedContourNode->SetProperty( "helper object", mitk::BoolProperty::New(!m_ShowMarkerNodes)); rotatedContourNode->SetProperty( "planarfigure.drawcontrolpoints", BoolProperty::New(false)); rotatedContourNode->SetProperty( "planarfigure.drawname", BoolProperty::New(false)); rotatedContourNode->SetProperty( "planarfigure.drawoutline", BoolProperty::New(false)); rotatedContourNode->SetProperty( "planarfigure.drawshadow", BoolProperty::New(false)); if (plane) { if ( id == size ) { m_ToolManager->GetDataStorage()->Add(rotatedContourNode, workingNode); } else { mitk::NodePredicateProperty::Pointer isMarker = mitk::NodePredicateProperty::New("isContourMarker", mitk::BoolProperty::New(true)); mitk::DataStorage::SetOfObjects::ConstPointer markers = m_ToolManager->GetDataStorage()->GetDerivations(workingNode,isMarker); for ( mitk::DataStorage::SetOfObjects::const_iterator iter = markers->begin(); iter != markers->end(); ++iter) { std::string nodeName = (*iter)->GetName(); unsigned int t = nodeName.find_last_of(" "); unsigned int markerId = atof(nodeName.substr(t+1).c_str())-1; if(id == markerId) { return id; } } m_ToolManager->GetDataStorage()->Add(rotatedContourNode, workingNode); } } return id; } void mitk::SegTool2D::InteractiveSegmentationBugMessage( const std::string& message ) { MITK_ERROR << "********************************************************************************" << std::endl << " " << message << std::endl << "********************************************************************************" << std::endl << " " << std::endl << " If your image is rotated or the 2D views don't really contain the patient image, try to press the button next to the image selection. " << std::endl << " " << std::endl << " Please file a BUG REPORT: " << std::endl << " http://bugs.mitk.org" << std::endl << " Contain the following information:" << std::endl << " - What image were you working on?" << std::endl << " - Which region of the image?" << std::endl << " - Which tool did you use?" << std::endl << " - What did you do?" << std::endl << " - What happened (not)? What did you expect?" << std::endl; }