diff --git a/Modules/Core/include/mitkPlaneGeometry.h b/Modules/Core/include/mitkPlaneGeometry.h index 7e5155952e..29db4db6e1 100644 --- a/Modules/Core/include/mitkPlaneGeometry.h +++ b/Modules/Core/include/mitkPlaneGeometry.h @@ -1,605 +1,608 @@ /*=================================================================== 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. ===================================================================*/ /** * \brief Describes the geometry of a plane object * * Describes a two-dimensional manifold, i.e., to put it simply, * an object that can be described using a 2D coordinate-system. * * PlaneGeometry can map points between 3D world coordinates * (in mm) and the described 2D coordinate-system (in mm) by first projecting * the 3D point onto the 2D manifold and then calculating the 2D-coordinates * (in mm). These 2D-mm-coordinates can be further converted into * 2D-unit-coordinates (e.g., pixels), giving a parameter representation of * the object with parameter values inside a rectangle * (e.g., [0,0]..[width, height]), which is the bounding box (bounding range * in z-direction always [0]..[1]). * * A PlaneGeometry describes the 2D representation within a 3D object (derived from BaseGeometry). For example, * a single CT-image (slice) is 2D in the sense that you can access the * pixels using 2D-coordinates, but is also 3D, as the pixels are really * voxels, thus have an extension (thickness) in the 3rd dimension. * * * Optionally, a reference BaseGeometry can be specified, which usually would * be the geometry associated with the underlying dataset. This is currently * used for calculating the intersection of inclined / rotated planes * (represented as PlaneGeometry) with the bounding box of the associated * BaseGeometry. * * \warning The PlaneGeometry are not necessarily up-to-date and not even * initialized. As described in the previous paragraph, one of the * Generate-/Copy-/UpdateOutputInformation methods have to initialize it. * mitk::BaseData::GetPlaneGeometry() makes sure, that the PlaneGeometry is * up-to-date before returning it (by setting the update extent appropriately * and calling UpdateOutputInformation). * * Rule: everything is in mm (or ms for temporal information) if not * stated otherwise. * \ingroup Geometry */ #ifndef PLANEGEOMETRY_H_HEADER_INCLUDED_C1C68A2C #define PLANEGEOMETRY_H_HEADER_INCLUDED_C1C68A2C #include "mitkBaseGeometry.h" #include "mitkRestorePlanePositionOperation.h" #include #include namespace mitk { template class Line; typedef Line Line3D; class PlaneGeometry; /** \deprecatedSince{2014_10} This class is deprecated. Please use PlaneGeometry instead. */ DEPRECATED(typedef PlaneGeometry Geometry2D); /** * \brief Describes a two-dimensional, rectangular plane * * \ingroup Geometry */ class MITKCORE_EXPORT PlaneGeometry : public BaseGeometry { public: mitkClassMacro(PlaneGeometry, BaseGeometry); /** Method for creation through the object factory. */ itkFactorylessNewMacro(Self) itkCloneMacro(Self) enum PlaneOrientation { Axial, Sagittal, Frontal, // also known as "Coronal" in mitk. None // This defines the PlaneGeometry for the 3D renderWindow which // curiously also needs a PlaneGeometry. This should be reconsidered some time. }; virtual void IndexToWorld(const Point2D &pt_units, Point2D &pt_mm) const; virtual void WorldToIndex(const Point2D &pt_mm, Point2D &pt_units) const; //##Documentation //## @brief Convert (continuous or discrete) index coordinates of a \em vector //## \a vec_units to world coordinates (in mm) //## @deprecated First parameter (Point2D) is not used. If possible, please use void IndexToWorld(const // mitk::Vector2D& vec_units, mitk::Vector2D& vec_mm) const. //## For further information about coordinates types, please see the Geometry documentation virtual void IndexToWorld(const mitk::Point2D &atPt2d_untis, const mitk::Vector2D &vec_units, mitk::Vector2D &vec_mm) const; //##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 virtual void IndexToWorld(const mitk::Vector2D &vec_units, mitk::Vector2D &vec_mm) const; //##Documentation //## @brief Convert world coordinates (in mm) of a \em vector //## \a vec_mm to (continuous!) index coordinates. //## @deprecated First parameter (Point2D) is not used. If possible, please use void WorldToIndex(const // mitk::Vector2D& vec_mm, mitk::Vector2D& vec_units) const. //## For further information about coordinates types, please see the Geometry documentation virtual void WorldToIndex(const mitk::Point2D &atPt2d_mm, const mitk::Vector2D &vec_mm, mitk::Vector2D &vec_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 virtual void WorldToIndex(const mitk::Vector2D &vec_mm, mitk::Vector2D &vec_units) const; /** * \brief Initialize a plane with orientation \a planeorientation * (default: axial) with respect to \a BaseGeometry (default: identity). * Spacing also taken from \a BaseGeometry. * * \warning A former version of this method created a geometry with unit * spacing. For unit spacing use * * \code * // for in-plane unit spacing: * thisgeometry->SetSizeInUnits(thisgeometry->GetExtentInMM(0), * thisgeometry->GetExtentInMM(1)); * // additionally, for unit spacing in normal direction (former version * // did not do this): * thisgeometry->SetExtentInMM(2, 1.0); * \endcode */ virtual void InitializeStandardPlane(const BaseGeometry *geometry3D, PlaneOrientation planeorientation = Axial, ScalarType zPosition = 0, bool frontside = true, - bool rotated = false); + bool rotated = false, + bool top = true); /** * \brief Initialize a plane with orientation \a planeorientation * (default: axial) with respect to \a BaseGeometry (default: identity). * Spacing also taken from \a BaseGeometry. * * \param top if \a true, create plane at top, otherwise at bottom * (for PlaneOrientation Axial, for other plane locations respectively) */ virtual void InitializeStandardPlane(const BaseGeometry *geometry3D, bool top, PlaneOrientation planeorientation = Axial, bool frontside = true, bool rotated = false); /** * \brief Initialize a plane with orientation \a planeorientation * (default: axial) with respect to \a transform (default: identity) * given width and height in units. * * \a Rotated means rotated by 180 degrees (1/2 rotation) within the plane. * Rotation by 90 degrees (1/4 rotation) is not implemented as of now. * * \a Frontside/Backside: * Viewed from below = frontside in the axial case; * (radiologist's view versus neuro-surgeon's view, see: * http://www.itk.org/Wiki/images/e/ed/DICOM-OrientationDiagram-Radiologist-vs-NeuroSurgeon.png ) * Viewed from front = frontside in the coronal case; * Viewed from left = frontside in the sagittal case. * * \a Cave/Caution: Currently only RPI, LAI, LPS and RAS in the three standard planes are covered, * i.e. 12 cases of 144: 3 standard planes * 48 coordinate orientations = 144 cases. */ virtual void InitializeStandardPlane(ScalarType width, ScalarType height, const AffineTransform3D *transform = nullptr, PlaneOrientation planeorientation = Axial, ScalarType zPosition = 0, bool frontside = true, - bool rotated = false); + bool rotated = false, + bool top = true); /** * \brief Initialize plane with orientation \a planeorientation * (default: axial) given width, height and spacing. * */ virtual void InitializeStandardPlane(ScalarType width, ScalarType height, const Vector3D &spacing, PlaneOrientation planeorientation = Axial, ScalarType zPosition = 0, bool frontside = true, - bool rotated = false); + bool rotated = false, + bool top = true); /** * \brief Initialize plane by width and height in pixels, right-/down-vector * (itk) to describe orientation in world-space (vectors will be normalized) * and spacing (default: 1.0 mm in all directions). * * The vectors are normalized and multiplied by the respective spacing before * they are set in the matrix. * * This overloaded version of InitializeStandardPlane() creates only righthanded * coordinate orientations, unless spacing contains 1 or 3 negative entries. * */ virtual void InitializeStandardPlane(ScalarType width, ScalarType height, const Vector3D &rightVector, const Vector3D &downVector, const Vector3D *spacing = nullptr); /** * \brief Initialize plane by width and height in pixels, * right-/down-vector (vnl) to describe orientation in world-space (vectors * will be normalized) and spacing (default: 1.0 mm in all directions). * * The vectors are normalized and multiplied by the respective spacing * before they are set in the matrix. * * This overloaded version of InitializeStandardPlane() creates only righthanded * coordinate orientations, unless spacing contains 1 or 3 negative entries. * */ virtual void InitializeStandardPlane(ScalarType width, ScalarType height, const VnlVector &rightVector, const VnlVector &downVector, const Vector3D *spacing = nullptr); /** * \brief Initialize plane by right-/down-vector (itk) and spacing * (default: 1.0 mm in all directions). * * The length of the right-/-down-vector is used as width/height in units, * respectively. Then, the vectors are normalized and multiplied by the * respective spacing before they are set in the matrix. */ virtual void InitializeStandardPlane(const Vector3D &rightVector, const Vector3D &downVector, const Vector3D *spacing = nullptr); /** * \brief Initialize plane by right-/down-vector (vnl) and spacing * (default: 1.0 mm in all directions). * * The length of the right-/-down-vector is used as width/height in units, * respectively. Then, the vectors are normalized and multiplied by the * respective spacing before they are set in the matrix. */ virtual void InitializeStandardPlane(const VnlVector &rightVector, const VnlVector &downVector, const Vector3D *spacing = nullptr); /** * \brief Initialize plane by origin and normal (size is 1.0 mm in * all directions, direction of right-/down-vector valid but * undefined). * \warning This function can only produce righthanded coordinate orientation, not lefthanded. */ virtual void InitializePlane(const Point3D &origin, const Vector3D &normal); /** * \brief Initialize plane by right-/down-vector. * * \warning The vectors are set into the matrix as they are, * \em without normalization! * This function creates a righthanded IndexToWorldTransform, * only a negative thickness could still make it lefthanded. */ void SetMatrixByVectors(const VnlVector &rightVector, const VnlVector &downVector, ScalarType thickness = 1.0); /** * \brief Change \a transform so that the third column of the * transform-martix is perpendicular to the first two columns * */ static void EnsurePerpendicularNormal(AffineTransform3D *transform); /** * \brief Normal of the plane * */ Vector3D GetNormal() const; /** * \brief Normal of the plane as VnlVector * */ VnlVector GetNormalVnl() const; virtual ScalarType SignedDistance(const Point3D &pt3d_mm) const; /** * \brief Calculates, whether a point is below or above the plane. There are two different *calculation methods, with or without consideration of the bounding box. */ virtual bool IsAbove(const Point3D &pt3d_mm, bool considerBoundingBox = false) const; /** * \brief Distance of the point from the plane * (bounding-box \em not considered) * */ ScalarType DistanceFromPlane(const Point3D &pt3d_mm) const; /** * \brief Signed distance of the point from the plane * (bounding-box \em not considered) * * > 0 : point is in the direction of the direction vector. */ inline ScalarType SignedDistanceFromPlane(const Point3D &pt3d_mm) const { ScalarType len = GetNormalVnl().two_norm(); if (len == 0) return 0; return (pt3d_mm - GetOrigin()) * GetNormal() / len; } /** * \brief Distance of the plane from another plane * (bounding-box \em not considered) * * Result is 0 if planes are not parallel. */ ScalarType DistanceFromPlane(const PlaneGeometry *plane) const { return fabs(SignedDistanceFromPlane(plane)); } /** * \brief Signed distance of the plane from another plane * (bounding-box \em not considered) * * Result is 0 if planes are not parallel. */ inline ScalarType SignedDistanceFromPlane(const PlaneGeometry *plane) const { if (IsParallel(plane)) { return SignedDistance(plane->GetOrigin()); } return 0; } /** * \brief Calculate the intersecting line of two planes * * \return \a true planes are intersecting * \return \a false planes do not intersect */ bool IntersectionLine(const PlaneGeometry *plane, Line3D &crossline) const; /** * \brief Calculate two points where another plane intersects the border of this plane * * \return number of intersection points (0..2). First interection point (if existing) * is returned in \a lineFrom, second in \a lineTo. */ unsigned int IntersectWithPlane2D(const PlaneGeometry *plane, Point2D &lineFrom, Point2D &lineTo) const; /** * \brief Calculate the angle between two planes * * \return angle in radiants */ double Angle(const PlaneGeometry *plane) const; /** * \brief Calculate the angle between the plane and a line * * \return angle in radiants */ double Angle(const Line3D &line) const; /** * \brief Calculate intersection point between the plane and a line * * \param intersectionPoint intersection point * \return \a true if \em unique intersection exists, i.e., if line * is \em not on or parallel to the plane */ bool IntersectionPoint(const Line3D &line, Point3D &intersectionPoint) const; /** * \brief Calculate line parameter of intersection point between the * plane and a line * * \param t parameter of line: intersection point is * line.GetPoint()+t*line.GetDirection() * \return \a true if \em unique intersection exists, i.e., if line * is \em not on or parallel to the plane */ bool IntersectionPointParam(const Line3D &line, double &t) const; /** * \brief Returns whether the plane is parallel to another plane * * @return true iff the normal vectors both point to the same or exactly oposit direction */ bool IsParallel(const PlaneGeometry *plane) const; /** * \brief Returns whether the point is on the plane * (bounding-box \em not considered) */ bool IsOnPlane(const Point3D &point) const; /** * \brief Returns whether the line is on the plane * (bounding-box \em not considered) */ bool IsOnPlane(const Line3D &line) const; /** * \brief Returns whether the plane is on the plane * (bounding-box \em not considered) * * @return true iff the normal vector of the planes point to the same or the exactly oposit direction and * the distance of the planes is < eps * */ bool IsOnPlane(const PlaneGeometry *plane) const; /** * \brief Returns the lot from the point to the plane */ Point3D ProjectPointOntoPlane(const Point3D &pt) const; virtual itk::LightObject::Pointer InternalClone() const override; /** Implements operation to re-orient the plane */ virtual void ExecuteOperation(Operation *operation) override; /** * \brief Project a 3D point given in mm (\a pt3d_mm) onto the 2D * geometry. The result is a 2D point in mm (\a pt2d_mm). * * The result is a 2D point in mm (\a pt2d_mm) relative to the upper-left * corner of the geometry. To convert this point into units (e.g., pixels * in case of an image), use WorldToIndex. * \return true projection was possible * \sa Project(const mitk::Point3D &pt3d_mm, mitk::Point3D * &projectedPt3d_mm) */ virtual bool Map(const mitk::Point3D &pt3d_mm, mitk::Point2D &pt2d_mm) const; /** * \brief Converts a 2D point given in mm (\a pt2d_mm) relative to the * upper-left corner of the geometry into the corresponding * world-coordinate (a 3D point in mm, \a pt3d_mm). * * To convert a 2D point given in units (e.g., pixels in case of an * image) into a 2D point given in mm (as required by this method), use * IndexToWorld. */ virtual void Map(const mitk::Point2D &pt2d_mm, mitk::Point3D &pt3d_mm) const; /** * \brief Set the width and height of this 2D-geometry in units by calling * SetBounds. This does \a not change the extent in mm! * * For an image, this is the number of pixels in x-/y-direction. * \note In contrast to calling SetBounds directly, this does \a not change * the extent in mm! */ virtual void SetSizeInUnits(mitk::ScalarType width, mitk::ScalarType height); /** * \brief Project a 3D point given in mm (\a pt3d_mm) onto the 2D * geometry. The result is a 3D point in mm (\a projectedPt3d_mm). * * \return true projection was possible */ virtual bool Project(const mitk::Point3D &pt3d_mm, mitk::Point3D &projectedPt3d_mm) const; /** * \brief Project a 3D vector given in mm (\a vec3d_mm) onto the 2D * geometry. The result is a 2D vector in mm (\a vec2d_mm). * * The result is a 2D vector in mm (\a vec2d_mm) relative to the * upper-left * corner of the geometry. To convert this point into units (e.g., pixels * in case of an image), use WorldToIndex. * \return true projection was possible * \sa Project(const mitk::Vector3D &vec3d_mm, mitk::Vector3D * &projectedVec3d_mm) */ virtual bool Map(const mitk::Point3D &atPt3d_mm, const mitk::Vector3D &vec3d_mm, mitk::Vector2D &vec2d_mm) const; /** * \brief Converts a 2D vector given in mm (\a vec2d_mm) relative to the * upper-left corner of the geometry into the corresponding * world-coordinate (a 3D vector in mm, \a vec3d_mm). * * To convert a 2D vector given in units (e.g., pixels in case of an * image) into a 2D vector given in mm (as required by this method), use * IndexToWorld. */ virtual void Map(const mitk::Point2D &atPt2d_mm, const mitk::Vector2D &vec2d_mm, mitk::Vector3D &vec3d_mm) const; /** * \brief Project a 3D vector given in mm (\a vec3d_mm) onto the 2D * geometry. The result is a 3D vector in mm (\a projectedVec3d_mm). * * DEPRECATED. Use Project(vector,vector) instead * * \return true projection was possible */ virtual bool Project(const mitk::Point3D &atPt3d_mm, const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm) const; /** * \brief Project a 3D vector given in mm (\a vec3d_mm) onto the 2D * geometry. The result is a 3D vector in mm (\a projectedVec3d_mm). * * \return true projection was possible */ virtual bool Project(const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm) const; /** * \brief Distance of the point from the geometry * (bounding-box \em not considered) * */ inline ScalarType Distance(const Point3D &pt3d_mm) const { return fabs(SignedDistance(pt3d_mm)); } /** * \brief Set the geometrical frame of reference in which this PlaneGeometry * is placed. * * This would usually be the BaseGeometry of the underlying dataset, but * setting it is optional. */ void SetReferenceGeometry(const mitk::BaseGeometry *geometry); /** * \brief Get the geometrical frame of reference for this PlaneGeometry. */ const BaseGeometry *GetReferenceGeometry() const; bool HasReferenceGeometry() const; protected: PlaneGeometry(); PlaneGeometry(const PlaneGeometry &other); virtual ~PlaneGeometry(); virtual void PrintSelf(std::ostream &os, itk::Indent indent) const override; const mitk::BaseGeometry *m_ReferenceGeometry; //##Documentation //## @brief PreSetSpacing //## //## These virtual function allows a different beahiour in subclasses. //## Do implement them in every subclass of BaseGeometry. If not needed, use //## {Superclass::PreSetSpacing();}; virtual void PreSetSpacing(const mitk::Vector3D &aSpacing) override { Superclass::PreSetSpacing(aSpacing); }; //##Documentation //## @brief CheckBounds //## //## This function is called in SetBounds. Assertions can be implemented in this function (see PlaneGeometry.cpp). //## If you implement this function in a subclass, make sure, that all classes were your class inherits from //## have an implementation of CheckBounds //## (e.g. inheritance BaseGeometry <- A <- B. Implementation of CheckBounds in class B needs implementation in A as // well!) virtual void CheckBounds(const BoundsArrayType &bounds) override; //##Documentation //## @brief CheckIndexToWorldTransform //## //## This function is called in SetIndexToWorldTransform. Assertions can be implemented in this function (see // PlaneGeometry.cpp). //## In Subclasses of BaseGeometry, implement own conditions or call Superclass::CheckBounds(bounds);. virtual void CheckIndexToWorldTransform(mitk::AffineTransform3D *transform) override; private: /** * \brief Compares plane with another plane: \a true if IsOnPlane * (bounding-box \em not considered) */ virtual bool operator==(const PlaneGeometry *) const { return false; }; /** * \brief Compares plane with another plane: \a false if IsOnPlane * (bounding-box \em not considered) */ virtual bool operator!=(const PlaneGeometry *) const { return false; }; }; } // namespace mitk #endif /* PLANEGEOMETRY_H_HEADER_INCLUDED_C1C68A2C */ diff --git a/Modules/Core/include/mitkSlicedGeometry3D.h b/Modules/Core/include/mitkSlicedGeometry3D.h index 54fe60c0f7..04acf4d680 100644 --- a/Modules/Core/include/mitkSlicedGeometry3D.h +++ b/Modules/Core/include/mitkSlicedGeometry3D.h @@ -1,327 +1,326 @@ /*=================================================================== 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 MITKSLICEDGEOMETRY3D_H_HEADER_INCLUDED_C1EBD0AD #define MITKSLICEDGEOMETRY3D_H_HEADER_INCLUDED_C1EBD0AD #include "mitkBaseGeometry.h" #include "mitkPlaneGeometry.h" namespace mitk { class SliceNavigationController; class NavigationController; /** \brief Describes the geometry of a data object consisting of slices. * * A PlaneGeometry can be requested for each slice. In the case of * \em evenly-spaced, \em plane geometries (m_EvenlySpaced==true), * only the 2D-geometry of the first slice has to be set (to an instance of * PlaneGeometry). The 2D geometries of the other slices are calculated * by shifting the first slice in the direction m_DirectionVector by * m_Spacing.z * sliceNumber. The m_Spacing member (which is only * relevant in the case m_EvenlySpaced==true) descibes the size of a voxel * (in mm), i.e., m_Spacing.x is the voxel width in the x-direction of the * plane. It is derived from the reference geometry of this SlicedGeometry3D, * which usually would be the global geometry describing how datasets are to * be resliced. * * By default, slices are oriented in the direction of one of the main axes * (x, y, z). However, by means of rotation, it is possible to realign the * slices in any possible direction. In case of an inclined plane, the spacing * is derived as a product of the (regular) geometry spacing and the direction * vector of the plane. * * SlicedGeometry3D and the associated PlaneGeometries have to be initialized in * the method GenerateOutputInformation() of BaseProcess (or CopyInformation / * UpdateOutputInformation of BaseData, if possible, e.g., by analyzing pic * tags in Image) subclasses. See also * * \sa itk::ProcessObject::GenerateOutputInformation(), * \sa itk::DataObject::CopyInformation() and * \a itk::DataObject::UpdateOutputInformation(). * * Rule: everything is in mm (or ms for temporal information) if not * stated otherwise. * * \warning The hull (i.e., transform, bounding-box and * time-bounds) is only guaranteed to be up-to-date after calling * UpdateInformation(). * * \ingroup Geometry */ class MITKCORE_EXPORT SlicedGeometry3D : public mitk::BaseGeometry { public: mitkClassMacro(SlicedGeometry3D, BaseGeometry) /** Method for creation through the object factory. */ itkFactorylessNewMacro(Self) itkCloneMacro(Self) /** * \brief Returns the PlaneGeometry of the slice (\a s). * * If (a) m_EvenlySpaced==true, (b) we don't have a PlaneGeometry stored * for the requested slice, and (c) the first slice (s=0) * is a PlaneGeometry instance, then we calculate the geometry of the * requested as the plane of the first slice shifted by m_Spacing[3]*s * in the direction of m_DirectionVector. * * \warning The PlaneGeometries are not necessarily up-to-date and not even * initialized. * * The PlaneGeometries have to be initialized in the method * GenerateOutputInformation() of BaseProcess (or CopyInformation / * UpdateOutputInformation of BaseData, if possible, e.g., by analyzing * pic tags in Image) subclasses. See also * * \sa itk::ProcessObject::GenerateOutputInformation(), * \sa itk::DataObject::CopyInformation() and * \sa itk::DataObject::UpdateOutputInformation(). */ virtual mitk::PlaneGeometry *GetPlaneGeometry(int s) const; /** * \deprecatedSince{2014_10} Please use GetPlaneGeometry */ DEPRECATED(const PlaneGeometry *GetGeometry2D(int s)) { return GetPlaneGeometry(s); } /** * \deprecatedSince{2014_10} Please use SetPlaneGeometry */ DEPRECATED(void SetGeometry2D(PlaneGeometry *geo, int s)) { SetPlaneGeometry(geo, s); } //##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) override; // virtual void SetTimeBounds( const mitk::TimeBounds& timebounds ); virtual const mitk::BoundingBox *GetBoundingBox() const override; /** * \brief Get the number of slices */ itkGetConstMacro(Slices, unsigned int) /** * \brief Set PlaneGeometry of slice \a s. */ virtual bool SetPlaneGeometry(mitk::PlaneGeometry *geometry2D, int s); /** * \brief Check whether a slice exists */ virtual bool IsValidSlice(int s = 0) const; virtual const BaseGeometry* GetReferenceGeometry() const; virtual void SetReferenceGeometry(const BaseGeometry *referenceGeometry); bool HasReferenceGeometry() const; /** * \brief Set the SliceNavigationController corresponding to this sliced * geometry. * * The SNC needs to be informed when the number of slices in the geometry * changes, which can occur whenthe slices are re-oriented by rotation. */ virtual void SetSliceNavigationController(mitk::SliceNavigationController *snc); mitk::SliceNavigationController *GetSliceNavigationController(); /** * \brief Set/Get whether the SlicedGeometry3D is evenly-spaced * (m_EvenlySpaced) * * If (a) m_EvenlySpaced==true, (b) we don't have a PlaneGeometry stored for * the requested slice, and (c) the first slice (s=0) is a PlaneGeometry * instance, then we calculate the geometry of the requested as the plane * of the first slice shifted by m_Spacing.z * s in the direction of * m_DirectionVector. * * \sa GetPlaneGeometry */ itkGetConstMacro(EvenlySpaced, bool) virtual void SetEvenlySpaced(bool on = true); /** * \brief Set/Get the vector between slices for the evenly-spaced case * (m_EvenlySpaced==true). * * If the direction-vector is (0,0,0) (the default) and the first * 2D geometry is a PlaneGeometry, then the direction-vector will be * calculated from the plane normal. * * \sa m_DirectionVector */ virtual void SetDirectionVector(const mitk::Vector3D &directionVector); itkGetConstMacro(DirectionVector, const mitk::Vector3D &) virtual itk::LightObject::Pointer InternalClone() const override; static const std::string SLICES; const static std::string DIRECTION_VECTOR; const static std::string EVENLY_SPACED; /** * \brief Tell this instance how many PlaneGeometries it shall manage. Bounding * box and the PlaneGeometries must be set additionally by calling the respective * methods! * * \warning Bounding box and the 2D-geometries must be set additionally: use * SetBounds(), SetGeometry(). */ virtual void InitializeSlicedGeometry(unsigned int slices); /** * \brief Completely initialize this instance as evenly-spaced with slices * parallel to the provided PlaneGeometry that is used as the first slice and * for spacing calculation. * * Initializes the bounding box according to the width/height of the * PlaneGeometry and \a slices. The spacing is calculated from the PlaneGeometry. */ - virtual void InitializeEvenlySpaced(mitk::PlaneGeometry *geometry2D, unsigned int slices, bool flipped = false); + virtual void InitializeEvenlySpaced(mitk::PlaneGeometry *geometry2D, unsigned int slices); /** * \brief Completely initialize this instance as evenly-spaced with slices * parallel to the provided PlaneGeometry that is used as the first slice and * for spacing calculation (except z-spacing). * * Initializes the bounding box according to the width/height of the * PlaneGeometry and \a slices. The x-/y-spacing is calculated from the * PlaneGeometry. */ virtual void InitializeEvenlySpaced(mitk::PlaneGeometry *geometry2D, mitk::ScalarType zSpacing, - unsigned int slices, - bool flipped = false); + unsigned int slices); /** * \brief Completely initialize this instance as evenly-spaced plane slices * parallel to a side of the provided BaseGeometry and using its spacing * information. * * Initializes the bounding box according to the width/height of the * BaseGeometry and the number of slices according to * BaseGeometry::GetExtent(2). * * \param planeorientation side parallel to which the slices will be oriented * \param top if \a true, create plane at top, otherwise at bottom * (for PlaneOrientation Axial, for other plane locations respectively) * \param frontside defines the side of the plane (the definition of * front/back is somewhat arbitrary) * * \param rotate rotates the plane by 180 degree around its normal (the * definition of rotated vs not rotated is somewhat arbitrary) */ virtual void InitializePlanes(const mitk::BaseGeometry *geometry3D, mitk::PlaneGeometry::PlaneOrientation planeorientation, bool top = true, bool frontside = true, bool rotated = false); virtual void SetImageGeometry(const bool isAnImageGeometry) override; virtual void ExecuteOperation(Operation *operation) override; static double CalculateSpacing(const mitk::Vector3D &spacing, const mitk::Vector3D &d); protected: SlicedGeometry3D(); SlicedGeometry3D(const SlicedGeometry3D &other); virtual ~SlicedGeometry3D(); /** * Reinitialize plane stack after rotation. More precisely, the first plane * of the stack needs to spatially aligned, in two respects: * * 1. Re-alignment with respect to the dataset center; this is necessary * since the distance from the first plane to the center could otherwise * continuously decrease or increase. * 2. Re-alignment with respect to a given reference point; the reference * point is a location which the user wants to be exactly touched by one * plane of the plane stack. The first plane is minimally shifted to * ensure this touching. Usually, the reference point would be the * point around which the geometry is rotated. */ virtual void ReinitializePlanes(const Point3D ¢er, const Point3D &referencePoint); ScalarType GetLargestExtent(const BaseGeometry *geometry); void PrintSelf(std::ostream &os, itk::Indent indent) const override; /** Calculate "directed spacing", i.e. the spacing in directions * non-orthogonal to the coordinate axes. This is done via the * ellipsoid equation. */ double CalculateSpacing(const mitk::Vector3D &direction) const; /** The extent of the slice stack, i.e. the number of slices, depends on the * plane normal. For rotated geometries, the geometry's transform needs to * be accounted in this calculation. */ mitk::Vector3D AdjustNormal(const mitk::Vector3D &normal) const; /** * Container for the 2D-geometries contained within this SliceGeometry3D. */ mutable std::vector m_PlaneGeometries; /** * If (a) m_EvenlySpaced==true, (b) we don't have a PlaneGeometry stored * for the requested slice, and (c) the first slice (s=0) * is a PlaneGeometry instance, then we calculate the geometry of the * requested as the plane of the first slice shifted by m_Spacing.z*s * in the direction of m_DirectionVector. * * \sa GetPlaneGeometry */ bool m_EvenlySpaced; /** * Vector between slices for the evenly-spaced case (m_EvenlySpaced==true). * If the direction-vector is (0,0,0) (the default) and the first * 2D geometry is a PlaneGeometry, then the direction-vector will be * calculated from the plane normal. */ mutable mitk::Vector3D m_DirectionVector; /** Number of slices this SliceGeometry3D is descibing. */ unsigned int m_Slices; /** Underlying BaseGeometry for this SlicedGeometry */ const mitk::BaseGeometry *m_ReferenceGeometry; /** SNC correcsponding to this geometry; used to reflect changes in the * number of slices due to rotation. */ // mitk::NavigationController *m_NavigationController; mitk::SliceNavigationController *m_SliceNavigationController; //##Documentation //## @brief PreSetSpacing //## //## These virtual function allows a different beahiour in subclasses. //## Do implement them in every subclass of BaseGeometry. If not needed, use //## {Superclass::PreSetSpacing();}; virtual void PreSetSpacing(const mitk::Vector3D &aSpacing) override; }; } // namespace mitk #endif /* MITKSLICEDGEOMETRY3D_H_HEADER_INCLUDED_C1EBD0AD */ diff --git a/Modules/Core/src/Controllers/mitkSliceNavigationController.cpp b/Modules/Core/src/Controllers/mitkSliceNavigationController.cpp index a548baab5e..5d480f8eaa 100644 --- a/Modules/Core/src/Controllers/mitkSliceNavigationController.cpp +++ b/Modules/Core/src/Controllers/mitkSliceNavigationController.cpp @@ -1,633 +1,641 @@ /*=================================================================== 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 "mitkSliceNavigationController.h" #include "mitkAction.h" #include "mitkBaseRenderer.h" #include "mitkCrosshairPositionEvent.h" #include "mitkInteractionConst.h" #include "mitkOperation.h" #include "mitkOperationActor.h" #include "mitkPlaneGeometry.h" #include "mitkProportionalTimeGeometry.h" #include "mitkRenderingManager.h" #include "mitkSlicedGeometry3D.h" #include "mitkVtkPropRenderer.h" #include "mitkImage.h" #include "mitkImagePixelReadAccessor.h" #include "mitkInteractionConst.h" #include "mitkNodePredicateDataType.h" #include "mitkOperationEvent.h" #include "mitkPixelTypeMultiplex.h" #include "mitkPlaneOperation.h" #include "mitkPointOperation.h" #include "mitkStatusBar.h" #include "mitkUndoController.h" #include "mitkApplyTransformMatrixOperation.h" #include "mitkMemoryUtilities.h" #include namespace mitk { SliceNavigationController::SliceNavigationController() : BaseController(), m_InputWorldGeometry3D(NULL), m_InputWorldTimeGeometry(NULL), m_CreatedWorldGeometry(NULL), m_ViewDirection(Axial), m_DefaultViewDirection(Axial), m_RenderingManager(NULL), m_Renderer(NULL), m_Top(false), m_FrontSide(false), m_Rotated(false), m_BlockUpdate(false), m_SliceLocked(false), m_SliceRotationLocked(false), m_OldPos(0) { typedef itk::SimpleMemberCommand SNCCommandType; SNCCommandType::Pointer sliceStepperChangedCommand, timeStepperChangedCommand; sliceStepperChangedCommand = SNCCommandType::New(); timeStepperChangedCommand = SNCCommandType::New(); sliceStepperChangedCommand->SetCallbackFunction(this, &SliceNavigationController::SendSlice); timeStepperChangedCommand->SetCallbackFunction(this, &SliceNavigationController::SendTime); m_Slice->AddObserver(itk::ModifiedEvent(), sliceStepperChangedCommand); m_Time->AddObserver(itk::ModifiedEvent(), timeStepperChangedCommand); m_Slice->SetUnitName("mm"); m_Time->SetUnitName("ms"); m_Top = false; m_FrontSide = false; m_Rotated = false; } SliceNavigationController::~SliceNavigationController() {} void SliceNavigationController::SetInputWorldGeometry3D(const BaseGeometry *geometry) { if (geometry != NULL) { if (const_cast(geometry->GetBoundingBox())->GetDiagonalLength2() < eps) { itkWarningMacro("setting an empty bounding-box"); geometry = NULL; } } if (m_InputWorldGeometry3D != geometry) { m_InputWorldGeometry3D = geometry; m_InputWorldTimeGeometry = NULL; this->Modified(); } } void SliceNavigationController::SetInputWorldTimeGeometry(const TimeGeometry *geometry) { if (geometry != NULL) { if (const_cast(geometry->GetBoundingBoxInWorld())->GetDiagonalLength2() < eps) { itkWarningMacro("setting an empty bounding-box"); geometry = NULL; } } if (m_InputWorldTimeGeometry != geometry) { m_InputWorldTimeGeometry = geometry; m_InputWorldGeometry3D = NULL; this->Modified(); } } RenderingManager *SliceNavigationController::GetRenderingManager() const { mitk::RenderingManager *renderingManager = m_RenderingManager.GetPointer(); if (renderingManager != NULL) return renderingManager; if (m_Renderer != NULL) { renderingManager = m_Renderer->GetRenderingManager(); if (renderingManager != NULL) return renderingManager; } return mitk::RenderingManager::GetInstance(); } void SliceNavigationController::SetViewDirectionToDefault() { m_ViewDirection = m_DefaultViewDirection; } const char *SliceNavigationController::GetViewDirectionAsString() const { const char *viewDirectionString; switch (m_ViewDirection) { case SliceNavigationController::Axial: viewDirectionString = "Axial"; break; case SliceNavigationController::Sagittal: viewDirectionString = "Sagittal"; break; case SliceNavigationController::Frontal: viewDirectionString = "Coronal"; break; case SliceNavigationController::Original: viewDirectionString = "Original"; break; default: viewDirectionString = "No View Direction Available"; break; } return viewDirectionString; } void SliceNavigationController::Update() { if (!m_BlockUpdate) { - if (m_ViewDirection == Axial) + if (m_ViewDirection == Sagittal) + { + this->Update(Sagittal, true, true, false); + } + else if (m_ViewDirection == Frontal) + { + this->Update(Frontal, false, true, false); + } + else if (m_ViewDirection == Axial) { this->Update(Axial, false, false, true); } else { this->Update(m_ViewDirection); } } } void SliceNavigationController::Update(SliceNavigationController::ViewDirection viewDirection, bool top, bool frontside, bool rotated) { TimeGeometry::ConstPointer worldTimeGeometry = m_InputWorldTimeGeometry; if (m_BlockUpdate || (m_InputWorldTimeGeometry.IsNull() && m_InputWorldGeometry3D.IsNull()) || ((worldTimeGeometry.IsNotNull()) && (worldTimeGeometry->CountTimeSteps() == 0))) { return; } m_BlockUpdate = true; if (m_InputWorldTimeGeometry.IsNotNull() && m_LastUpdateTime < m_InputWorldTimeGeometry->GetMTime()) { Modified(); } if (m_InputWorldGeometry3D.IsNotNull() && m_LastUpdateTime < m_InputWorldGeometry3D->GetMTime()) { Modified(); } this->SetViewDirection(viewDirection); this->SetTop(top); this->SetFrontSide(frontside); this->SetRotated(rotated); if (m_LastUpdateTime < GetMTime()) { m_LastUpdateTime = GetMTime(); // initialize the viewplane SlicedGeometry3D::Pointer slicedWorldGeometry = NULL; BaseGeometry::ConstPointer currentGeometry = NULL; if (m_InputWorldTimeGeometry.IsNotNull()) if (m_InputWorldTimeGeometry->IsValidTimeStep(GetTime()->GetPos())) currentGeometry = m_InputWorldTimeGeometry->GetGeometryForTimeStep(GetTime()->GetPos()); else currentGeometry = m_InputWorldTimeGeometry->GetGeometryForTimeStep(0); else currentGeometry = m_InputWorldGeometry3D; m_CreatedWorldGeometry = NULL; switch (viewDirection) { case Original: if (worldTimeGeometry.IsNotNull()) { m_CreatedWorldGeometry = worldTimeGeometry->Clone(); worldTimeGeometry = m_CreatedWorldGeometry.GetPointer(); slicedWorldGeometry = dynamic_cast( m_CreatedWorldGeometry->GetGeometryForTimeStep(this->GetTime()->GetPos()).GetPointer()); if (slicedWorldGeometry.IsNotNull()) { break; } } else { const SlicedGeometry3D *worldSlicedGeometry = dynamic_cast(currentGeometry.GetPointer()); if (worldSlicedGeometry != NULL) { slicedWorldGeometry = static_cast(currentGeometry->Clone().GetPointer()); break; } } slicedWorldGeometry = SlicedGeometry3D::New(); slicedWorldGeometry->InitializePlanes(currentGeometry, PlaneGeometry::None, top, frontside, rotated); slicedWorldGeometry->SetSliceNavigationController(this); break; case Axial: slicedWorldGeometry = SlicedGeometry3D::New(); slicedWorldGeometry->InitializePlanes(currentGeometry, PlaneGeometry::Axial, top, frontside, rotated); slicedWorldGeometry->SetSliceNavigationController(this); break; case Frontal: slicedWorldGeometry = SlicedGeometry3D::New(); slicedWorldGeometry->InitializePlanes(currentGeometry, PlaneGeometry::Frontal, top, frontside, rotated); slicedWorldGeometry->SetSliceNavigationController(this); break; case Sagittal: slicedWorldGeometry = SlicedGeometry3D::New(); slicedWorldGeometry->InitializePlanes(currentGeometry, PlaneGeometry::Sagittal, top, frontside, rotated); slicedWorldGeometry->SetSliceNavigationController(this); break; default: itkExceptionMacro("unknown ViewDirection"); } m_Slice->SetPos(0); m_Slice->SetSteps((int)slicedWorldGeometry->GetSlices()); if (m_CreatedWorldGeometry.IsNull()) { // initialize TimeGeometry m_CreatedWorldGeometry = ProportionalTimeGeometry::New(); } if (worldTimeGeometry.IsNull()) { m_CreatedWorldGeometry = ProportionalTimeGeometry::New(); dynamic_cast(m_CreatedWorldGeometry.GetPointer()) ->Initialize(slicedWorldGeometry, 1); m_Time->SetSteps(0); m_Time->SetPos(0); m_Time->InvalidateRange(); } else { m_BlockUpdate = true; m_Time->SetSteps(worldTimeGeometry->CountTimeSteps()); m_Time->SetPos(0); const TimeBounds &timeBounds = worldTimeGeometry->GetTimeBounds(); m_Time->SetRange(timeBounds[0], timeBounds[1]); m_BlockUpdate = false; assert(worldTimeGeometry->GetGeometryForTimeStep(this->GetTime()->GetPos()).IsNotNull()); TimePointType minimumTimePoint = worldTimeGeometry->TimeStepToTimePoint(this->GetTime()->GetPos()); TimePointType stepDuration = worldTimeGeometry->TimeStepToTimePoint(this->GetTime()->GetPos() + 1) - worldTimeGeometry->TimeStepToTimePoint(this->GetTime()->GetPos()); //@todo implement for non-evenly-timed geometry! m_CreatedWorldGeometry = ProportionalTimeGeometry::New(); dynamic_cast(m_CreatedWorldGeometry.GetPointer()) ->Initialize(slicedWorldGeometry, worldTimeGeometry->CountTimeSteps()); dynamic_cast(m_CreatedWorldGeometry.GetPointer()) ->GetMinimumTimePoint(minimumTimePoint); dynamic_cast(m_CreatedWorldGeometry.GetPointer())->SetStepDuration(stepDuration); } } // unblock update; we may do this now, because if m_BlockUpdate was already // true before this method was entered, then we will never come here. m_BlockUpdate = false; // Send the geometry. Do this even if nothing was changed, because maybe // Update() was only called to re-send the old geometry and time/slice data. this->SendCreatedWorldGeometry(); this->SendSlice(); this->SendTime(); // Adjust the stepper range of slice stepper according to geometry this->AdjustSliceStepperRange(); } void SliceNavigationController::SendCreatedWorldGeometry() { // Send the geometry. Do this even if nothing was changed, because maybe // Update() was only called to re-send the old geometry. if (!m_BlockUpdate) { this->InvokeEvent(GeometrySendEvent(m_CreatedWorldGeometry, 0)); } } void SliceNavigationController::SendCreatedWorldGeometryUpdate() { if (!m_BlockUpdate) { this->InvokeEvent(GeometryUpdateEvent(m_CreatedWorldGeometry, m_Slice->GetPos())); } } void SliceNavigationController::SendSlice() { if (!m_BlockUpdate) { if (m_CreatedWorldGeometry.IsNotNull()) { this->InvokeEvent(GeometrySliceEvent(m_CreatedWorldGeometry, m_Slice->GetPos())); // send crosshair event crosshairPositionEvent.Send(); // Request rendering update for all views this->GetRenderingManager()->RequestUpdateAll(); } } } void SliceNavigationController::SendTime() { if (!m_BlockUpdate) { if (m_CreatedWorldGeometry.IsNotNull()) { this->InvokeEvent(GeometryTimeEvent(m_CreatedWorldGeometry, m_Time->GetPos())); // Request rendering update for all views this->GetRenderingManager()->RequestUpdateAll(); } } } void SliceNavigationController::SetGeometry(const itk::EventObject &) {} void SliceNavigationController::SetGeometryTime(const itk::EventObject &geometryTimeEvent) { if (m_CreatedWorldGeometry.IsNull()) { return; } const SliceNavigationController::GeometryTimeEvent *timeEvent = dynamic_cast(&geometryTimeEvent); assert(timeEvent != NULL); TimeGeometry *timeGeometry = timeEvent->GetTimeGeometry(); assert(timeGeometry != NULL); int timeStep = (int)timeEvent->GetPos(); ScalarType timeInMS; timeInMS = timeGeometry->TimeStepToTimePoint(timeStep); timeStep = m_CreatedWorldGeometry->TimePointToTimeStep(timeInMS); this->GetTime()->SetPos(timeStep); } void SliceNavigationController::SetGeometrySlice(const itk::EventObject &geometrySliceEvent) { const SliceNavigationController::GeometrySliceEvent *sliceEvent = dynamic_cast(&geometrySliceEvent); assert(sliceEvent != NULL); this->GetSlice()->SetPos(sliceEvent->GetPos()); } void SliceNavigationController::SelectSliceByPoint(const Point3D &point) { if (m_CreatedWorldGeometry.IsNull()) { return; } //@todo add time to PositionEvent and use here!! SlicedGeometry3D *slicedWorldGeometry = dynamic_cast( m_CreatedWorldGeometry->GetGeometryForTimeStep(this->GetTime()->GetPos()).GetPointer()); if (slicedWorldGeometry) { int bestSlice = -1; double bestDistance = itk::NumericTraits::max(); int s, slices; slices = slicedWorldGeometry->GetSlices(); if (slicedWorldGeometry->GetEvenlySpaced()) { mitk::PlaneGeometry *plane = slicedWorldGeometry->GetPlaneGeometry(0); const Vector3D &direction = slicedWorldGeometry->GetDirectionVector(); Point3D projectedPoint; plane->Project(point, projectedPoint); // Check whether the point is somewhere within the slice stack volume; // otherwise, the defualt slice (0) will be selected if (direction[0] * (point[0] - projectedPoint[0]) + direction[1] * (point[1] - projectedPoint[1]) + direction[2] * (point[2] - projectedPoint[2]) >= 0) { bestSlice = (int)(plane->Distance(point) / slicedWorldGeometry->GetSpacing()[2] + 0.5); } } else { Point3D projectedPoint; for (s = 0; s < slices; ++s) { slicedWorldGeometry->GetPlaneGeometry(s)->Project(point, projectedPoint); Vector3D distance = projectedPoint - point; ScalarType currentDistance = distance.GetSquaredNorm(); if (currentDistance < bestDistance) { bestDistance = currentDistance; bestSlice = s; } } } if (bestSlice >= 0) { this->GetSlice()->SetPos(bestSlice); } else { this->GetSlice()->SetPos(0); } this->SendCreatedWorldGeometryUpdate(); } } void SliceNavigationController::ReorientSlices(const Point3D &point, const Vector3D &normal) { if (m_CreatedWorldGeometry.IsNull()) { return; } PlaneOperation op(OpORIENT, point, normal); m_CreatedWorldGeometry->ExecuteOperation(&op); this->SendCreatedWorldGeometryUpdate(); } void SliceNavigationController::ReorientSlices(const mitk::Point3D &point, const mitk::Vector3D &axisVec0, const mitk::Vector3D &axisVec1) { if (m_CreatedWorldGeometry) { PlaneOperation op(OpORIENT, point, axisVec0, axisVec1); m_CreatedWorldGeometry->ExecuteOperation(&op); this->SendCreatedWorldGeometryUpdate(); } } mitk::TimeGeometry *SliceNavigationController::GetCreatedWorldGeometry() { return m_CreatedWorldGeometry; } const mitk::BaseGeometry *SliceNavigationController::GetCurrentGeometry3D() { if (m_CreatedWorldGeometry.IsNotNull()) { return m_CreatedWorldGeometry->GetGeometryForTimeStep(this->GetTime()->GetPos()); } else { return NULL; } } const mitk::PlaneGeometry *SliceNavigationController::GetCurrentPlaneGeometry() { const mitk::SlicedGeometry3D *slicedGeometry = dynamic_cast(this->GetCurrentGeometry3D()); if (slicedGeometry) { const mitk::PlaneGeometry *planeGeometry = (slicedGeometry->GetPlaneGeometry(this->GetSlice()->GetPos())); return planeGeometry; } else { return NULL; } } void SliceNavigationController::SetRenderer(BaseRenderer *renderer) { m_Renderer = renderer; } BaseRenderer *SliceNavigationController::GetRenderer() const { return m_Renderer; } void SliceNavigationController::AdjustSliceStepperRange() { const mitk::SlicedGeometry3D *slicedGeometry = dynamic_cast(this->GetCurrentGeometry3D()); const Vector3D &direction = slicedGeometry->GetDirectionVector(); int c = 0; int i, k = 0; for (i = 0; i < 3; ++i) { if (fabs(direction[i]) < 0.000000001) { ++c; } else { k = i; } } if (c == 2) { ScalarType min = slicedGeometry->GetOrigin()[k]; ScalarType max = min + slicedGeometry->GetExtentInMM(k); m_Slice->SetRange(min, max); } else { m_Slice->InvalidateRange(); } } void SliceNavigationController::ExecuteOperation(Operation *operation) { // switch on type // - select best slice for a given point // - rotate created world geometry according to Operation->SomeInfo() if (!operation || m_CreatedWorldGeometry.IsNull()) { return; } switch (operation->GetOperationType()) { case OpMOVE: // should be a point operation { if (!m_SliceLocked) // do not move the cross position { // select a slice PointOperation *po = dynamic_cast(operation); if (po && po->GetIndex() == -1) { this->SelectSliceByPoint(po->GetPoint()); } else if (po && po->GetIndex() != -1) // undo case because index != -1, index holds the old position of this slice { this->GetSlice()->SetPos(po->GetIndex()); } } break; } case OpRESTOREPLANEPOSITION: { m_CreatedWorldGeometry->ExecuteOperation(operation); this->SendCreatedWorldGeometryUpdate(); break; } case OpAPPLYTRANSFORMMATRIX: { m_CreatedWorldGeometry->ExecuteOperation(operation); this->SendCreatedWorldGeometryUpdate(); break; } default: { // do nothing break; } } } } // namespace diff --git a/Modules/Core/src/DataManagement/mitkPlaneGeometry.cpp b/Modules/Core/src/DataManagement/mitkPlaneGeometry.cpp index 931460dcfe..c9d64bffec 100644 --- a/Modules/Core/src/DataManagement/mitkPlaneGeometry.cpp +++ b/Modules/Core/src/DataManagement/mitkPlaneGeometry.cpp @@ -1,916 +1,984 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #include "mitkPlaneGeometry.h" #include "mitkInteractionConst.h" #include "mitkLine.h" #include "mitkPlaneOperation.h" +#include + #include #include namespace mitk { PlaneGeometry::PlaneGeometry() : Superclass(), m_ReferenceGeometry(nullptr) { Initialize(); } PlaneGeometry::~PlaneGeometry() {} PlaneGeometry::PlaneGeometry(const PlaneGeometry &other) : Superclass(other), m_ReferenceGeometry(other.m_ReferenceGeometry) { } void PlaneGeometry::EnsurePerpendicularNormal(mitk::AffineTransform3D *transform) { /** \brief ensure column(2) of indexToWorldTransform-matrix to be perpendicular to plane, keep length and * handedness. */ VnlVector normal = vnl_cross_3d(transform->GetMatrix().GetVnlMatrix().get_column(0), transform->GetMatrix().GetVnlMatrix().get_column(1)); normal.normalize(); // Now normal is a righthand normal unit vector, perpendicular to the plane. ScalarType len = transform->GetMatrix().GetVnlMatrix().get_column(2).two_norm(); if (len == 0) { len = 1; } normal *= len; // Get the existing normal vector zed: vnl_vector_fixed zed = transform->GetMatrix().GetVnlMatrix().get_column(2); /** If det(matrix)<0, multiply normal vector by (-1) to keep geometry lefthanded. */ if (vnl_determinant(transform->GetMatrix().GetVnlMatrix()) < 0) { MITK_DEBUG << "EnsurePerpendicularNormal(): Lefthanded geometry preserved, rh-normal: [ " << normal << " ],"; normal *= (-1.0); MITK_DEBUG << "lh-normal: [ " << normal << " ], original vector zed is: [ " << zed << " ]"; } // Now lets compare and only replace if necessary and only then warn the user: // float epsilon is precise enough here, but we need to respect numerical condition: // Higham, N., 2002, Accuracy and Stability of Numerical Algorithms, // SIAM, page 37, 2nd edition: double feps = std::numeric_limits::epsilon(); double zedsMagnitude = zed.two_norm(); feps = feps * zedsMagnitude * 2; /** Check if normal (3. column) was perpendicular: If not, replace with calculated normal vector: */ if (normal != zed) { vnl_vector_fixed parallel; for (unsigned int i = 0; i < 3; ++i) { parallel[i] = normal[i] / zed[i]; // Remember linear algebra: checking for parallelity. } // Checking if really not paralell i.e. non-colinear vectors by comparing these floating point numbers: if ((parallel[0] + feps < parallel[1] || feps + parallel[1] < parallel[0]) && (parallel[0] + feps < parallel[2] || feps + parallel[2] < parallel[0])) { MITK_WARN << "EnsurePerpendicularNormal(): Plane geometry was _/askew/_, so here it gets rectified by substituting" << " the 3rd column of the indexToWorldMatrix with an appropriate normal vector: [ " << normal << " ], original vector zed was: [ " << zed << " ]."; Matrix3D matrix = transform->GetMatrix(); matrix.GetVnlMatrix().set_column(2, normal); transform->SetMatrix(matrix); } } else { // Nothing to do, 3rd column of indexToWorldTransformMatrix already was perfectly perpendicular. } } void PlaneGeometry::CheckIndexToWorldTransform(mitk::AffineTransform3D *transform) { EnsurePerpendicularNormal(transform); } void PlaneGeometry::CheckBounds(const BoundingBox::BoundsArrayType &bounds) { // error: unused parameter 'bounds' // this happens in release mode, where the assert macro is defined empty // hence we "use" the parameter: (void)bounds; // currently the unit rectangle must be starting at the origin [0,0] assert(bounds[0] == 0); assert(bounds[2] == 0); // the unit rectangle must be two-dimensional assert(bounds[1] > 0); assert(bounds[3] > 0); } void PlaneGeometry::IndexToWorld(const Point2D &pt_units, Point2D &pt_mm) const { pt_mm[0] = GetExtentInMM(0) / GetExtent(0) * pt_units[0]; pt_mm[1] = GetExtentInMM(1) / GetExtent(1) * pt_units[1]; } void PlaneGeometry::WorldToIndex(const Point2D &pt_mm, Point2D &pt_units) const { pt_units[0] = pt_mm[0] * (1.0 / (GetExtentInMM(0) / GetExtent(0))); pt_units[1] = pt_mm[1] * (1.0 / (GetExtentInMM(1) / GetExtent(1))); } void PlaneGeometry::IndexToWorld(const Point2D & /*atPt2d_units*/, const Vector2D &vec_units, Vector2D &vec_mm) const { MITK_WARN << "Warning! Call of the deprecated function PlaneGeometry::IndexToWorld(point, vec, vec). Use " "PlaneGeometry::IndexToWorld(vec, vec) instead!"; this->IndexToWorld(vec_units, vec_mm); } void PlaneGeometry::IndexToWorld(const Vector2D &vec_units, Vector2D &vec_mm) const { vec_mm[0] = (GetExtentInMM(0) / GetExtent(0)) * vec_units[0]; vec_mm[1] = (GetExtentInMM(1) / GetExtent(1)) * vec_units[1]; } void PlaneGeometry::WorldToIndex(const Point2D & /*atPt2d_mm*/, const Vector2D &vec_mm, Vector2D &vec_units) const { MITK_WARN << "Warning! Call of the deprecated function PlaneGeometry::WorldToIndex(point, vec, vec). Use " "PlaneGeometry::WorldToIndex(vec, vec) instead!"; this->WorldToIndex(vec_mm, vec_units); } void PlaneGeometry::WorldToIndex(const Vector2D &vec_mm, Vector2D &vec_units) const { vec_units[0] = vec_mm[0] * (1.0 / (GetExtentInMM(0) / GetExtent(0))); vec_units[1] = vec_mm[1] * (1.0 / (GetExtentInMM(1) / GetExtent(1))); } void PlaneGeometry::InitializeStandardPlane(mitk::ScalarType width, ScalarType height, const Vector3D &spacing, PlaneGeometry::PlaneOrientation planeorientation, ScalarType zPosition, bool frontside, - bool rotated) + bool rotated, + bool top) { AffineTransform3D::Pointer transform; transform = AffineTransform3D::New(); AffineTransform3D::MatrixType matrix; AffineTransform3D::MatrixType::InternalMatrixType &vnlmatrix = matrix.GetVnlMatrix(); vnlmatrix.set_identity(); vnlmatrix(0, 0) = spacing[0]; vnlmatrix(1, 1) = spacing[1]; vnlmatrix(2, 2) = spacing[2]; transform->SetIdentity(); transform->SetMatrix(matrix); - InitializeStandardPlane(width, height, transform.GetPointer(), planeorientation, zPosition, frontside, rotated); + InitializeStandardPlane(width, height, transform.GetPointer(), planeorientation, zPosition, frontside, rotated, top); } void PlaneGeometry::InitializeStandardPlane(mitk::ScalarType width, mitk::ScalarType height, const AffineTransform3D *transform /* = nullptr */, PlaneGeometry::PlaneOrientation planeorientation /* = Axial */, mitk::ScalarType zPosition /* = 0 */, bool frontside /* = true */, - bool rotated /* = false */) + bool rotated /* = false */, + bool top /* = true */) { Superclass::Initialize(); /// construct standard view. // We define at the moment "frontside" as: axial from above, // coronal from front (nose), saggital from right. // TODO: Double check with medicals doctors or radiologists [ ]. // We define the orientation in patient's view, e.g. LAI is in a axial cut // (parallel to the triangle ear-ear-nose): // first axis: To the left ear of the patient // seecond axis: To the nose of the patient // third axis: To the legs of the patient. // Options are: L/R left/right; A/P anterior/posterior; I/S inferior/superior // (AKA caudal/cranial). // We note on all cases in the following switch block r.h. for right handed // or l.h. for left handed to describe the different cases. // However, which system is chosen is defined at the end of the switch block. // CAVE / be careful: the vectors right and bottom are relative to the plane // and do NOT describe e.g. the right side of the patient. Point3D origin; /** Bottom means downwards, DV means Direction Vector. Both relative to the image! */ VnlVector rightDV(3), bottomDV(3); /** Origin of this plane is by default a zero vector and implicitly in the top-left corner: */ origin.Fill(0); /** This is different to all definitions in MITK, except the QT mouse clicks. * But it is like this here and we don't want to change a running system. * Just be aware, that IN THIS FUNCTION we define the origin at the top left (e.g. your screen). */ /** NormalDirection defines which axis (i.e. column index in the transform matrix) * is perpendicular to the plane: */ int normalDirection; switch (planeorientation) // Switch through our limited choice of standard planes: { case None: /** Orientation 'None' shall be done like the axial plane orientation, * for whatever reasons. */ case Axial: if (frontside) // Radiologist's view from below. A cut along the triangle ear-ear-nose. { if (rotated == false) /** Origin in the top-left corner, x=[1; 0; 0], y=[0; 1; 0], z=[0; 0; 1], * origin=[0,0,zpos]: LAI (r.h.) * * 0---rightDV----> | * | | * | Picture of a finite, rectangular plane | * | ( insert LOLCAT-scan here ^_^ ) | * | | * v _________________________________________| * */ { FillVector3D(origin, 0, 0, zPosition); FillVector3D(rightDV, 1, 0, 0); FillVector3D(bottomDV, 0, 1, 0); } else // Origin rotated to the bottom-right corner, x=[-1; 0; 0], y=[0; -1; 0], z=[0; 0; 1], // origin=[w,h,zpos]: RPI (r.h.) { // Caveat emptor: Still using top-left as origin of index coordinate system! FillVector3D(origin, width, height, zPosition); FillVector3D(rightDV, -1, 0, 0); FillVector3D(bottomDV, 0, -1, 0); } } else // 'Backside, not frontside.' Neuro-Surgeons's view from above patient. { if (rotated == false) // x=[-1; 0; 0], y=[0; 1; 0], z=[0; 0; 1], origin=[w,0,zpos]: RAS (r.h.) { FillVector3D(origin, width, 0, zPosition); FillVector3D(rightDV, -1, 0, 0); FillVector3D(bottomDV, 0, 1, 0); } else // Origin in the bottom-left corner, x=[1; 0; 0], y=[0; -1; 0], z=[0; 0; 1], // origin=[0,h,zpos]: LPS (r.h.) { FillVector3D(origin, 0, height, zPosition); FillVector3D(rightDV, 1, 0, 0); FillVector3D(bottomDV, 0, -1, 0); } } normalDirection = 2; // That is S=Superior=z=third_axis=middlefinger in righthanded LPS-system. break; // Frontal is known as Coronal in mitk. Plane cuts through patient's ear-ear-heel-heel: case Frontal: if (frontside) { if (rotated == false) // x=[1; 0; 0], y=[0; 0; 1], z=[0; 1; 0], origin=[0,zpos,0]: LAI (r.h.) { FillVector3D(origin, 0, zPosition, 0); FillVector3D(rightDV, 1, 0, 0); FillVector3D(bottomDV, 0, 0, 1); } else // x=[-1;0;0], y=[0;0;-1], z=[0;1;0], origin=[w,zpos,h]: RAS (r.h.) { FillVector3D(origin, width, zPosition, height); FillVector3D(rightDV, -1, 0, 0); FillVector3D(bottomDV, 0, 0, -1); } } else { if (rotated == false) // x=[-1;0;0], y=[0;0;1], z=[0;1;0], origin=[w,zpos,0]: RPI (r.h.) { FillVector3D(origin, width, zPosition, 0); FillVector3D(rightDV, -1, 0, 0); FillVector3D(bottomDV, 0, 0, 1); } else // x=[1;0;0], y=[0;1;0], z=[0;0;-1], origin=[0,zpos,h]: LPS (r.h.) { FillVector3D(origin, 0, zPosition, height); FillVector3D(rightDV, 1, 0, 0); FillVector3D(bottomDV, 0, 0, -1); } } normalDirection = 1; // Normal vector = posterior direction. break; case Sagittal: // Sagittal=Medial plane, the symmetry-plane mirroring your face. if (frontside) { if (rotated == false) // x=[0;1;0], y=[0;0;1], z=[1;0;0], origin=[zpos,0,0]: LAI (r.h.) { FillVector3D(origin, zPosition, 0, 0); FillVector3D(rightDV, 0, 1, 0); FillVector3D(bottomDV, 0, 0, 1); } else // x=[0;-1;0], y=[0;0;-1], z=[1;0;0], origin=[zpos,w,h]: LPS (r.h.) { FillVector3D(origin, zPosition, width, height); FillVector3D(rightDV, 0, -1, 0); FillVector3D(bottomDV, 0, 0, -1); } } else { if (rotated == false) // x=[0;-1;0], y=[0;0;1], z=[1;0;0], origin=[zpos,w,0]: RPI (r.h.) { FillVector3D(origin, zPosition, width, 0); FillVector3D(rightDV, 0, -1, 0); FillVector3D(bottomDV, 0, 0, 1); } else // x=[0;1;0], y=[0;0;-1], z=[1;0;0], origin=[zpos,0,h]: RAS (r.h.) { FillVector3D(origin, zPosition, 0, height); FillVector3D(rightDV, 0, 1, 0); FillVector3D(bottomDV, 0, 0, -1); } } normalDirection = 0; // Normal vector = Lateral direction: Left in a LPS-system. break; default: itkExceptionMacro("unknown PlaneOrientation"); } - /// Checking if lefthanded or righthanded: - mitk::ScalarType lhOrRhSign = (+1); - if (transform != nullptr) - { - lhOrRhSign = ((0 < vnl_determinant(transform->GetMatrix().GetVnlMatrix())) ? (+1) : (-1)); - - MITK_DEBUG - << "mitk::PlaneGeometry::InitializeStandardPlane(): lhOrRhSign, normalDirection, NameOfClass, ObjectName = " - << lhOrRhSign << ", " << normalDirection << ", " << transform->GetNameOfClass() << ", " - << transform->GetObjectName() << "."; - - origin = transform->TransformPoint(origin); - rightDV = transform->TransformVector(rightDV); - bottomDV = transform->TransformVector(bottomDV); + VnlVector normal(3); + FillVector3D(normal, 0, 0, 0); + normal[normalDirection] = top ? 1 : -1; - /// Signing normal vector according to l.h./r.h. coordinate system: - this->SetMatrixByVectors( - rightDV, bottomDV, transform->GetMatrix().GetVnlMatrix().get_column(normalDirection).two_norm() * lhOrRhSign); - } - else + if ( transform != nullptr ) { - this->SetMatrixByVectors(rightDV, bottomDV); + origin = transform->TransformPoint( origin ); + rightDV = transform->TransformVector( rightDV ); + bottomDV = transform->TransformVector( bottomDV ); + normal = transform->TransformVector( normal ); } ScalarType bounds[6] = {0, width, 0, height, 0, 1}; - this->SetBounds(bounds); + AffineTransform3D::Pointer planeTransform = AffineTransform3D::New(); + Matrix3D matrix; + matrix.GetVnlMatrix().set_column(0, rightDV); + matrix.GetVnlMatrix().set_column(1, bottomDV); + matrix.GetVnlMatrix().set_column(2, normal); + planeTransform->SetMatrix(matrix); + planeTransform->SetOffset(this->GetIndexToWorldTransform()->GetOffset()); + this->SetIndexToWorldTransform(planeTransform); + this->SetOrigin(origin); } void PlaneGeometry::InitializeStandardPlane(const BaseGeometry *geometry3D, PlaneOrientation planeorientation, ScalarType zPosition, bool frontside, - bool rotated) + bool rotated, + bool top) { this->SetReferenceGeometry(geometry3D); ScalarType width, height; - const BoundingBox::BoundsArrayType &boundsarray = geometry3D->GetBoundingBox()->GetBounds(); + // Inspired by: + // http://www.na-mic.org/Wiki/index.php/Coordinate_System_Conversion_Between_ITK_and_Slicer3 + + mitk::AffineTransform3D::MatrixType matrix = geometry3D->GetIndexToWorldTransform()->GetMatrix(); - Vector3D originVector; - FillVector3D(originVector, boundsarray[0], boundsarray[2], boundsarray[4]); + matrix.GetVnlMatrix().normalize_columns(); + mitk::AffineTransform3D::MatrixType::InternalMatrixType inverseMatrix = matrix.GetInverse(); - if (geometry3D->GetImageGeometry()) + /// The index of the sagittal, coronal and axial axes in the reference geometry. + int axes[3]; + /// The direction of the sagittal, coronal and axial axes in the reference geometry. + /// +1 means that the direction is straight, i.e. greater index translates to greater + /// world coordinate. -1 means that the direction is inverted. + int directions[3]; + ScalarType extents[3]; + ScalarType spacings[3]; + for (int i = 0; i < 3; ++i) { - FillVector3D(originVector, originVector[0] - 0.5, originVector[1] - 0.5, originVector[2] - 0.5); + int dominantAxis = itk::Function::Max3( + inverseMatrix[0][i], + inverseMatrix[1][i], + inverseMatrix[2][i] + ); + axes[i] = dominantAxis; + directions[i] = itk::Function::Sign(inverseMatrix[dominantAxis][i]); + extents[i] = geometry3D->GetExtent(dominantAxis); + spacings[i] = geometry3D->GetSpacing()[dominantAxis]; } - switch (planeorientation) + + // matrix(column) = inverseTransformMatrix(row) * flippedAxes * spacing + matrix[0][0] = inverseMatrix[axes[0]][0] * directions[0] * spacings[0]; + matrix[1][0] = inverseMatrix[axes[0]][1] * directions[0] * spacings[0]; + matrix[2][0] = inverseMatrix[axes[0]][2] * directions[0] * spacings[0]; + matrix[0][1] = inverseMatrix[axes[1]][0] * directions[1] * spacings[1]; + matrix[1][1] = inverseMatrix[axes[1]][1] * directions[1] * spacings[1]; + matrix[2][1] = inverseMatrix[axes[1]][2] * directions[1] * spacings[1]; + matrix[0][2] = inverseMatrix[axes[2]][0] * directions[2] * spacings[2]; + matrix[1][2] = inverseMatrix[axes[2]][1] * directions[2] * spacings[2]; + matrix[2][2] = inverseMatrix[axes[2]][2] * directions[2] * spacings[2]; + + /// The "world origin" is the corner with the lowest physical coordinates. + /// We use it as a reference point so that we get the correct anatomical + /// orientations. + Point3D worldOrigin = geometry3D->GetOrigin(); + for (int i = 0; i < 3; ++i) { - case None: - case Axial: - width = geometry3D->GetExtent(0); - height = geometry3D->GetExtent(1); - break; - case Frontal: - width = geometry3D->GetExtent(0); - height = geometry3D->GetExtent(2); - break; - case Sagittal: - width = geometry3D->GetExtent(1); - height = geometry3D->GetExtent(2); - break; - default: - itkExceptionMacro("unknown PlaneOrientation"); + /// The distance of the plane origin from the world origin in voxels. + double offset = directions[i] > 0 ? 0.0 : extents[i]; + + if (geometry3D->GetImageGeometry()) + { + offset += directions[i] * 0.5; + } + + for (int j = 0; j < 3; ++j) + { + worldOrigin[j] -= offset * matrix[j][i]; + } } - InitializeStandardPlane( - width, height, geometry3D->GetIndexToWorldTransform(), planeorientation, zPosition, frontside, rotated); + switch(planeorientation) + { + case None: + /** Orientation 'None' shall be done like the axial plane orientation, + * for whatever reasons. */ + case Axial: + width = extents[0]; + height = extents[1]; + break; + case Frontal: + width = extents[0]; + height = extents[2]; + break; + case Sagittal: + width = extents[1]; + height = extents[2]; + break; + default: + itkExceptionMacro("unknown PlaneOrientation"); + } - ScalarType bounds[6] = {0, width, 0, height, 0, 1}; - this->SetBounds(bounds); + ScalarType bounds[6]= { 0, width, 0, height, 0, 1 }; + this->SetBounds( bounds ); - Point3D origin; - originVector = geometry3D->GetIndexToWorldTransform()->TransformVector(originVector); + AffineTransform3D::Pointer transform = AffineTransform3D::New(); + transform->SetMatrix(matrix); + transform->SetOffset(worldOrigin.GetVectorFromOrigin()); - origin = GetOrigin() + originVector; - SetOrigin(origin); + InitializeStandardPlane( + width, height, transform, planeorientation, zPosition, frontside, rotated, top); } void PlaneGeometry::InitializeStandardPlane( const BaseGeometry *geometry3D, bool top, PlaneOrientation planeorientation, bool frontside, bool rotated) { - ScalarType zPosition; - - switch (planeorientation) + /// The index of the sagittal, coronal and axial axes in world coordinate system. + int worldAxis; + switch(planeorientation) { - case Axial: - zPosition = (top ? 0.5 : geometry3D->GetExtent(2) - 0.5); - break; - case Frontal: - zPosition = (top ? 0.5 : geometry3D->GetExtent(1) - 0.5); - break; - case Sagittal: - zPosition = (top ? 0.5 : geometry3D->GetExtent(0) - 0.5); - break; - case None: - zPosition = (top ? 0 : geometry3D->GetExtent(2) - 1.0); - break; - default: - itkExceptionMacro("unknown PlaneOrientation"); + case None: + /** Orientation 'None' shall be done like the axial plane orientation, + * for whatever reasons. */ + case Axial: + worldAxis = 2; + break; + case Frontal: + worldAxis = 1; + break; + case Sagittal: + worldAxis = 0; + break; + default: + itkExceptionMacro("unknown PlaneOrientation"); } - InitializeStandardPlane(geometry3D, planeorientation, zPosition, frontside, rotated); + // Inspired by: + // http://www.na-mic.org/Wiki/index.php/Coordinate_System_Conversion_Between_ITK_and_Slicer3 + + mitk::AffineTransform3D::ConstPointer affineTransform = geometry3D->GetIndexToWorldTransform(); + mitk::AffineTransform3D::MatrixType matrix = affineTransform->GetMatrix(); + matrix.GetVnlMatrix().normalize_columns(); + mitk::AffineTransform3D::MatrixType::InternalMatrixType inverseMatrix = matrix.GetInverse(); + + /// The index of the sagittal, coronal and axial axes in the reference geometry. + int dominantAxis = itk::Function::Max3( + inverseMatrix[0][worldAxis], + inverseMatrix[1][worldAxis], + inverseMatrix[2][worldAxis]); + + ScalarType zPosition = top ? 0.5 : geometry3D->GetExtent(dominantAxis) - 0.5; + + InitializeStandardPlane(geometry3D, planeorientation, zPosition, frontside, rotated, top); } void PlaneGeometry::InitializeStandardPlane(const Vector3D &rightVector, const Vector3D &downVector, const Vector3D *spacing) { InitializeStandardPlane(rightVector.GetVnlVector(), downVector.GetVnlVector(), spacing); } void PlaneGeometry::InitializeStandardPlane(const VnlVector &rightVector, const VnlVector &downVector, const Vector3D *spacing) { ScalarType width = rightVector.two_norm(); ScalarType height = downVector.two_norm(); InitializeStandardPlane(width, height, rightVector, downVector, spacing); } void PlaneGeometry::InitializeStandardPlane(mitk::ScalarType width, ScalarType height, const Vector3D &rightVector, const Vector3D &downVector, const Vector3D *spacing) { InitializeStandardPlane(width, height, rightVector.GetVnlVector(), downVector.GetVnlVector(), spacing); } void PlaneGeometry::InitializeStandardPlane(mitk::ScalarType width, ScalarType height, const VnlVector &rightVector, const VnlVector &downVector, const Vector3D *spacing) { assert(width > 0); assert(height > 0); VnlVector rightDV = rightVector; rightDV.normalize(); VnlVector downDV = downVector; downDV.normalize(); VnlVector normal = vnl_cross_3d(rightVector, downVector); normal.normalize(); // Crossproduct vnl_cross_3d is always righthanded, but that is okay here // because in this method we create a new IndexToWorldTransform and // spacing with 1 or 3 negative components could still make it lefthanded. if (spacing != nullptr) { rightDV *= (*spacing)[0]; downDV *= (*spacing)[1]; normal *= (*spacing)[2]; } AffineTransform3D::Pointer transform = AffineTransform3D::New(); Matrix3D matrix; matrix.GetVnlMatrix().set_column(0, rightDV); matrix.GetVnlMatrix().set_column(1, downDV); matrix.GetVnlMatrix().set_column(2, normal); transform->SetMatrix(matrix); transform->SetOffset(this->GetIndexToWorldTransform()->GetOffset()); ScalarType bounds[6] = {0, width, 0, height, 0, 1}; this->SetBounds(bounds); this->SetIndexToWorldTransform(transform); } void PlaneGeometry::InitializePlane(const Point3D &origin, const Vector3D &normal) { VnlVector rightVectorVnl(3), downVectorVnl; if (Equal(normal[1], 0.0f) == false) { FillVector3D(rightVectorVnl, 1.0f, -normal[0] / normal[1], 0.0f); rightVectorVnl.normalize(); } else { FillVector3D(rightVectorVnl, 0.0f, 1.0f, 0.0f); } downVectorVnl = vnl_cross_3d(normal.GetVnlVector(), rightVectorVnl); downVectorVnl.normalize(); // Crossproduct vnl_cross_3d is always righthanded. InitializeStandardPlane(rightVectorVnl, downVectorVnl); SetOrigin(origin); } void PlaneGeometry::SetMatrixByVectors(const VnlVector &rightVector, const VnlVector &downVector, ScalarType thickness /* = 1.0 */) { VnlVector normal = vnl_cross_3d(rightVector, downVector); normal.normalize(); normal *= thickness; // Crossproduct vnl_cross_3d is always righthanded, but that is okay here // because in this method we create a new IndexToWorldTransform and // a negative thickness could still make it lefthanded. AffineTransform3D::Pointer transform = AffineTransform3D::New(); Matrix3D matrix; matrix.GetVnlMatrix().set_column(0, rightVector); matrix.GetVnlMatrix().set_column(1, downVector); matrix.GetVnlMatrix().set_column(2, normal); transform->SetMatrix(matrix); transform->SetOffset(this->GetIndexToWorldTransform()->GetOffset()); SetIndexToWorldTransform(transform); } Vector3D PlaneGeometry::GetNormal() const { Vector3D frontToBack; frontToBack.SetVnlVector(this->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix().get_column(2)); return frontToBack; } VnlVector PlaneGeometry::GetNormalVnl() const { return this->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix().get_column(2); } ScalarType PlaneGeometry::DistanceFromPlane(const Point3D &pt3d_mm) const { return fabs(SignedDistance(pt3d_mm)); } ScalarType PlaneGeometry::SignedDistance(const Point3D &pt3d_mm) const { return SignedDistanceFromPlane(pt3d_mm); } bool PlaneGeometry::IsAbove(const Point3D &pt3d_mm, bool considerBoundingBox) const { if (considerBoundingBox) { Point3D pt3d_units; BaseGeometry::WorldToIndex(pt3d_mm, pt3d_units); return (pt3d_units[2] > this->GetBoundingBox()->GetBounds()[4]); } else return SignedDistanceFromPlane(pt3d_mm) > 0; } bool PlaneGeometry::IntersectionLine(const PlaneGeometry *plane, Line3D &crossline) const { Vector3D normal = this->GetNormal(); normal.Normalize(); Vector3D planeNormal = plane->GetNormal(); planeNormal.Normalize(); Vector3D direction = itk::CrossProduct(normal, planeNormal); if (direction.GetSquaredNorm() < eps) return false; crossline.SetDirection(direction); double N1dN2 = normal * planeNormal; double determinant = 1.0 - N1dN2 * N1dN2; Vector3D origin = this->GetOrigin().GetVectorFromOrigin(); Vector3D planeOrigin = plane->GetOrigin().GetVectorFromOrigin(); double d1 = normal * origin; double d2 = planeNormal * planeOrigin; double c1 = (d1 - d2 * N1dN2) / determinant; double c2 = (d2 - d1 * N1dN2) / determinant; Vector3D p = normal * c1 + planeNormal * c2; crossline.GetPoint().GetVnlVector() = p.GetVnlVector(); return true; } unsigned int PlaneGeometry::IntersectWithPlane2D(const PlaneGeometry *plane, Point2D &lineFrom, Point2D &lineTo) const { Line3D crossline; if (this->IntersectionLine(plane, crossline) == false) return 0; Point2D point2; Vector2D direction2; this->Map(crossline.GetPoint(), point2); this->Map(crossline.GetPoint(), crossline.GetDirection(), direction2); return Line3D::RectangleLineIntersection( 0, 0, GetExtentInMM(0), GetExtentInMM(1), point2, direction2, lineFrom, lineTo); } double PlaneGeometry::Angle(const PlaneGeometry *plane) const { return angle(plane->GetMatrixColumn(2), GetMatrixColumn(2)); } double PlaneGeometry::Angle(const Line3D &line) const { return vnl_math::pi_over_2 - angle(line.GetDirection().GetVnlVector(), GetMatrixColumn(2)); } bool PlaneGeometry::IntersectionPoint(const Line3D &line, Point3D &intersectionPoint) const { Vector3D planeNormal = this->GetNormal(); planeNormal.Normalize(); Vector3D lineDirection = line.GetDirection(); lineDirection.Normalize(); double t = planeNormal * lineDirection; if (fabs(t) < eps) { return false; } Vector3D diff; diff = this->GetOrigin() - line.GetPoint(); t = (planeNormal * diff) / t; intersectionPoint = line.GetPoint() + lineDirection * t; return true; } bool PlaneGeometry::IntersectionPointParam(const Line3D &line, double &t) const { Vector3D planeNormal = this->GetNormal(); Vector3D lineDirection = line.GetDirection(); t = planeNormal * lineDirection; if (fabs(t) < eps) { return false; } Vector3D diff; diff = this->GetOrigin() - line.GetPoint(); t = (planeNormal * diff) / t; return true; } bool PlaneGeometry::IsParallel(const PlaneGeometry *plane) const { return ((Angle(plane) < 10.0 * mitk::sqrteps) || (Angle(plane) > (vnl_math::pi - 10.0 * sqrteps))); } bool PlaneGeometry::IsOnPlane(const Point3D &point) const { return Distance(point) < eps; } bool PlaneGeometry::IsOnPlane(const Line3D &line) const { return ((Distance(line.GetPoint()) < eps) && (Distance(line.GetPoint2()) < eps)); } bool PlaneGeometry::IsOnPlane(const PlaneGeometry *plane) const { return (IsParallel(plane) && (Distance(plane->GetOrigin()) < eps)); } Point3D PlaneGeometry::ProjectPointOntoPlane(const Point3D &pt) const { ScalarType len = this->GetNormalVnl().two_norm(); return pt - this->GetNormal() * this->SignedDistanceFromPlane(pt) / len; } itk::LightObject::Pointer PlaneGeometry::InternalClone() const { Self::Pointer newGeometry = new PlaneGeometry(*this); newGeometry->UnRegister(); return newGeometry.GetPointer(); } void PlaneGeometry::ExecuteOperation(Operation *operation) { vtkTransform *transform = vtkTransform::New(); transform->SetMatrix(this->GetVtkMatrix()); switch (operation->GetOperationType()) { case OpORIENT: { mitk::PlaneOperation *planeOp = dynamic_cast(operation); if (planeOp == nullptr) { return; } Point3D center = planeOp->GetPoint(); Vector3D orientationVector = planeOp->GetNormal(); Vector3D defaultVector; FillVector3D(defaultVector, 0.0, 0.0, 1.0); Vector3D rotationAxis = itk::CrossProduct(orientationVector, defaultVector); // double rotationAngle = acos( orientationVector[2] / orientationVector.GetNorm() ); double rotationAngle = atan2((double)rotationAxis.GetNorm(), (double)(orientationVector * defaultVector)); rotationAngle *= 180.0 / vnl_math::pi; transform->PostMultiply(); transform->Identity(); transform->Translate(center[0], center[1], center[2]); transform->RotateWXYZ(rotationAngle, rotationAxis[0], rotationAxis[1], rotationAxis[2]); transform->Translate(-center[0], -center[1], -center[2]); break; } case OpRESTOREPLANEPOSITION: { RestorePlanePositionOperation *op = dynamic_cast(operation); if (op == nullptr) { return; } AffineTransform3D::Pointer transform2 = AffineTransform3D::New(); Matrix3D matrix; matrix.GetVnlMatrix().set_column(0, op->GetTransform()->GetMatrix().GetVnlMatrix().get_column(0)); matrix.GetVnlMatrix().set_column(1, op->GetTransform()->GetMatrix().GetVnlMatrix().get_column(1)); matrix.GetVnlMatrix().set_column(2, op->GetTransform()->GetMatrix().GetVnlMatrix().get_column(2)); transform2->SetMatrix(matrix); Vector3D offset = op->GetTransform()->GetOffset(); transform2->SetOffset(offset); this->SetIndexToWorldTransform(transform2); ScalarType bounds[6] = {0, op->GetWidth(), 0, op->GetHeight(), 0, 1}; this->SetBounds(bounds); this->Modified(); transform->Delete(); return; } default: Superclass::ExecuteOperation(operation); transform->Delete(); return; } this->SetVtkMatrixDeepCopy(transform); this->Modified(); transform->Delete(); } void PlaneGeometry::PrintSelf(std::ostream &os, itk::Indent indent) const { Superclass::PrintSelf(os, indent); os << indent << " ScaleFactorMMPerUnitX: " << GetExtentInMM(0) / GetExtent(0) << std::endl; os << indent << " ScaleFactorMMPerUnitY: " << GetExtentInMM(1) / GetExtent(1) << std::endl; os << indent << " Normal: " << GetNormal() << std::endl; } bool PlaneGeometry::Map(const mitk::Point3D &pt3d_mm, mitk::Point2D &pt2d_mm) const { assert(this->IsBoundingBoxNull() == false); Point3D pt3d_units; Superclass::WorldToIndex(pt3d_mm, pt3d_units); pt2d_mm[0] = pt3d_units[0] * GetExtentInMM(0) / GetExtent(0); pt2d_mm[1] = pt3d_units[1] * GetExtentInMM(1) / GetExtent(1); pt3d_units[2] = 0; return const_cast(this->GetBoundingBox())->IsInside(pt3d_units); } void PlaneGeometry::Map(const mitk::Point2D &pt2d_mm, mitk::Point3D &pt3d_mm) const { // pt2d_mm is measured from the origin of the world geometry (at leats it called form BaseRendere::Mouse...Event) Point3D pt3d_units; pt3d_units[0] = pt2d_mm[0] / (GetExtentInMM(0) / GetExtent(0)); pt3d_units[1] = pt2d_mm[1] / (GetExtentInMM(1) / GetExtent(1)); pt3d_units[2] = 0; // pt3d_units is a continuos index. We divided it with the Scale Factor (= spacing in x and y) to convert it from mm // to index units. // pt3d_mm = GetIndexToWorldTransform()->TransformPoint(pt3d_units); // now we convert the 3d index to a 3D world point in mm. We could have used IndexToWorld as well as // GetITW->Transform... } void PlaneGeometry::SetSizeInUnits(mitk::ScalarType width, mitk::ScalarType height) { ScalarType bounds[6] = {0, width, 0, height, 0, 1}; ScalarType extent, newextentInMM; if (GetExtent(0) > 0) { extent = GetExtent(0); if (width > extent) newextentInMM = GetExtentInMM(0) / width * extent; else newextentInMM = GetExtentInMM(0) * extent / width; SetExtentInMM(0, newextentInMM); } if (GetExtent(1) > 0) { extent = GetExtent(1); if (width > extent) newextentInMM = GetExtentInMM(1) / height * extent; else newextentInMM = GetExtentInMM(1) * extent / height; SetExtentInMM(1, newextentInMM); } SetBounds(bounds); } bool PlaneGeometry::Project(const mitk::Point3D &pt3d_mm, mitk::Point3D &projectedPt3d_mm) const { assert(this->IsBoundingBoxNull() == false); Point3D pt3d_units; Superclass::WorldToIndex(pt3d_mm, pt3d_units); pt3d_units[2] = 0; projectedPt3d_mm = GetIndexToWorldTransform()->TransformPoint(pt3d_units); return const_cast(this->GetBoundingBox())->IsInside(pt3d_units); } bool PlaneGeometry::Project(const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm) const { assert(this->IsBoundingBoxNull() == false); Vector3D vec3d_units; Superclass::WorldToIndex(vec3d_mm, vec3d_units); vec3d_units[2] = 0; projectedVec3d_mm = GetIndexToWorldTransform()->TransformVector(vec3d_units); return true; } bool PlaneGeometry::Project(const mitk::Point3D &atPt3d_mm, const mitk::Vector3D &vec3d_mm, mitk::Vector3D &projectedVec3d_mm) const { MITK_WARN << "Deprecated function! Call Project(vec3D,vec3D) instead."; assert(this->IsBoundingBoxNull() == false); Vector3D vec3d_units; Superclass::WorldToIndex(atPt3d_mm, vec3d_mm, vec3d_units); vec3d_units[2] = 0; projectedVec3d_mm = GetIndexToWorldTransform()->TransformVector(vec3d_units); Point3D pt3d_units; Superclass::WorldToIndex(atPt3d_mm, pt3d_units); return const_cast(this->GetBoundingBox())->IsInside(pt3d_units); } bool PlaneGeometry::Map(const mitk::Point3D &atPt3d_mm, const mitk::Vector3D &vec3d_mm, mitk::Vector2D &vec2d_mm) const { Point2D pt2d_mm_start, pt2d_mm_end; Point3D pt3d_mm_end; bool inside = Map(atPt3d_mm, pt2d_mm_start); pt3d_mm_end = atPt3d_mm + vec3d_mm; inside &= Map(pt3d_mm_end, pt2d_mm_end); vec2d_mm = pt2d_mm_end - pt2d_mm_start; return inside; } void PlaneGeometry::Map(const mitk::Point2D & /*atPt2d_mm*/, const mitk::Vector2D & /*vec2d_mm*/, mitk::Vector3D & /*vec3d_mm*/) const { //@todo implement parallel to the other Map method! assert(false); } void PlaneGeometry::SetReferenceGeometry(const mitk::BaseGeometry *geometry) { m_ReferenceGeometry = geometry; } const mitk::BaseGeometry *PlaneGeometry::GetReferenceGeometry() const { return m_ReferenceGeometry; } bool PlaneGeometry::HasReferenceGeometry() const { return (m_ReferenceGeometry != nullptr); } } // namespace diff --git a/Modules/Core/src/DataManagement/mitkSlicedGeometry3D.cpp b/Modules/Core/src/DataManagement/mitkSlicedGeometry3D.cpp index 87eb07ffb0..60a3be37dc 100644 --- a/Modules/Core/src/DataManagement/mitkSlicedGeometry3D.cpp +++ b/Modules/Core/src/DataManagement/mitkSlicedGeometry3D.cpp @@ -1,975 +1,961 @@ /*=================================================================== 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 "mitkSlicedGeometry3D.h" #include "mitkAbstractTransformGeometry.h" #include "mitkApplyTransformMatrixOperation.h" #include "mitkInteractionConst.h" #include "mitkPlaneGeometry.h" #include "mitkPlaneOperation.h" #include "mitkRestorePlanePositionOperation.h" #include "mitkRotationOperation.h" #include "mitkSliceNavigationController.h" const mitk::ScalarType PI = 3.14159265359; mitk::SlicedGeometry3D::SlicedGeometry3D() : m_EvenlySpaced(true), m_Slices(0), m_ReferenceGeometry(nullptr), m_SliceNavigationController(nullptr) { m_DirectionVector.Fill(0); this->InitializeSlicedGeometry(m_Slices); } mitk::SlicedGeometry3D::SlicedGeometry3D(const SlicedGeometry3D &other) : Superclass(other), m_EvenlySpaced(other.m_EvenlySpaced), m_Slices(other.m_Slices), m_ReferenceGeometry(other.m_ReferenceGeometry), m_SliceNavigationController(other.m_SliceNavigationController) { m_DirectionVector.Fill(0); SetSpacing(other.GetSpacing()); SetDirectionVector(other.GetDirectionVector()); if (m_EvenlySpaced) { PlaneGeometry::Pointer geometry = other.m_PlaneGeometries[0]->Clone(); assert(geometry.IsNotNull()); SetPlaneGeometry(geometry, 0); } else { unsigned int s; for (s = 0; s < other.m_Slices; ++s) { if (other.m_PlaneGeometries[s].IsNull()) { assert(other.m_EvenlySpaced); m_PlaneGeometries[s] = nullptr; } else { PlaneGeometry *geometry2D = other.m_PlaneGeometries[s]->Clone(); assert(geometry2D != nullptr); SetPlaneGeometry(geometry2D, s); } } } } mitk::SlicedGeometry3D::~SlicedGeometry3D() { } mitk::PlaneGeometry *mitk::SlicedGeometry3D::GetPlaneGeometry(int s) const { mitk::PlaneGeometry::Pointer geometry2D = nullptr; if (this->IsValidSlice(s)) { geometry2D = m_PlaneGeometries[s]; // If (a) m_EvenlySpaced==true, (b) we don't have a PlaneGeometry stored // for the requested slice, and (c) the first slice (s=0) // is a PlaneGeometry instance, then we calculate the geometry of the // requested as the plane of the first slice shifted by m_Spacing[2]*s // in the direction of m_DirectionVector. if ((m_EvenlySpaced) && (geometry2D.IsNull())) { PlaneGeometry *firstSlice = m_PlaneGeometries[0]; if (firstSlice != nullptr && dynamic_cast(m_PlaneGeometries[0].GetPointer()) == nullptr) { if ((m_DirectionVector[0] == 0.0) && (m_DirectionVector[1] == 0.0) && (m_DirectionVector[2] == 0.0)) { m_DirectionVector = firstSlice->GetNormal(); m_DirectionVector.Normalize(); } Vector3D direction; direction = m_DirectionVector * this->GetSpacing()[2]; mitk::PlaneGeometry::Pointer requestedslice; requestedslice = static_cast(firstSlice->Clone().GetPointer()); requestedslice->SetOrigin(requestedslice->GetOrigin() + direction * s); geometry2D = requestedslice; m_PlaneGeometries[s] = geometry2D; } } return geometry2D; } else { return nullptr; } } const mitk::BoundingBox *mitk::SlicedGeometry3D::GetBoundingBox() const { assert(this->IsBoundingBoxNull() == false); return Superclass::GetBoundingBox(); } bool mitk::SlicedGeometry3D::SetPlaneGeometry(mitk::PlaneGeometry *geometry2D, int s) { if (this->IsValidSlice(s)) { m_PlaneGeometries[s] = geometry2D; m_PlaneGeometries[s]->SetReferenceGeometry(m_ReferenceGeometry); return true; } return false; } void mitk::SlicedGeometry3D::InitializeSlicedGeometry(unsigned int slices) { Superclass::Initialize(); m_Slices = slices; PlaneGeometry::Pointer gnull = nullptr; m_PlaneGeometries.assign(m_Slices, gnull); Vector3D spacing; spacing.Fill(1.0); this->SetSpacing(spacing); m_DirectionVector.Fill(0); } -void mitk::SlicedGeometry3D::InitializeEvenlySpaced(mitk::PlaneGeometry *geometry2D, unsigned int slices, bool flipped) +void mitk::SlicedGeometry3D::InitializeEvenlySpaced(mitk::PlaneGeometry *geometry2D, unsigned int slices) { assert(geometry2D != nullptr); - this->InitializeEvenlySpaced(geometry2D, geometry2D->GetExtentInMM(2) / geometry2D->GetExtent(2), slices, flipped); + this->InitializeEvenlySpaced(geometry2D, geometry2D->GetExtentInMM(2) / geometry2D->GetExtent(2), slices); } void mitk::SlicedGeometry3D::InitializeEvenlySpaced(mitk::PlaneGeometry *geometry2D, mitk::ScalarType zSpacing, - unsigned int slices, - bool flipped) + unsigned int slices) { assert(geometry2D != nullptr); assert(geometry2D->GetExtent(0) > 0); assert(geometry2D->GetExtent(1) > 0); geometry2D->Register(); Superclass::Initialize(); m_Slices = slices; BoundingBox::BoundsArrayType bounds = geometry2D->GetBounds(); bounds[4] = 0; bounds[5] = slices; // clear and reserve PlaneGeometry::Pointer gnull = nullptr; m_PlaneGeometries.assign(m_Slices, gnull); Vector3D directionVector = geometry2D->GetAxisVector(2); directionVector.Normalize(); directionVector *= zSpacing; - if (flipped == false) - { - // Normally we should use the following four lines to create a copy of - // the transform contrained in geometry2D, because it may not be changed - // by us. But we know that SetSpacing creates a new transform without - // changing the old (coming from geometry2D), so we can use the fifth - // line instead. We check this at (**). - // - // AffineTransform3D::Pointer transform = AffineTransform3D::New(); - // transform->SetMatrix(geometry2D->GetIndexToWorldTransform()->GetMatrix()); - // transform->SetOffset(geometry2D->GetIndexToWorldTransform()->GetOffset()); - // SetIndexToWorldTransform(transform); - - this->SetIndexToWorldTransform(const_cast(geometry2D->GetIndexToWorldTransform())); - } - else - { - directionVector *= -1.0; - this->SetIndexToWorldTransform(AffineTransform3D::New()); - this->GetIndexToWorldTransform()->SetMatrix(geometry2D->GetIndexToWorldTransform()->GetMatrix()); - - AffineTransform3D::OutputVectorType scaleVector; - FillVector3D(scaleVector, 1.0, 1.0, -1.0); - this->GetIndexToWorldTransform()->Scale(scaleVector, true); - this->GetIndexToWorldTransform()->SetOffset(geometry2D->GetIndexToWorldTransform()->GetOffset()); - } + // Normally we should use the following four lines to create a copy of + // the transform contrained in geometry2D, because it may not be changed + // by us. But we know that SetSpacing creates a new transform without + // changing the old (coming from geometry2D), so we can use the fifth + // line instead. We check this at (**). + // + // AffineTransform3D::Pointer transform = AffineTransform3D::New(); + // transform->SetMatrix(geometry2D->GetIndexToWorldTransform()->GetMatrix()); + // transform->SetOffset(geometry2D->GetIndexToWorldTransform()->GetOffset()); + // SetIndexToWorldTransform(transform); + + this->SetIndexToWorldTransform(const_cast(geometry2D->GetIndexToWorldTransform())); mitk::Vector3D spacing; FillVector3D(spacing, geometry2D->GetExtentInMM(0) / bounds[1], geometry2D->GetExtentInMM(1) / bounds[3], zSpacing); this->SetDirectionVector(directionVector); this->SetBounds(bounds); this->SetPlaneGeometry(geometry2D, 0); this->SetSpacing(spacing, true); this->SetEvenlySpaced(); // this->SetTimeBounds( geometry2D->GetTimeBounds() ); assert(this->GetIndexToWorldTransform() != geometry2D->GetIndexToWorldTransform()); // (**) see above. this->SetFrameOfReferenceID(geometry2D->GetFrameOfReferenceID()); this->SetImageGeometry(geometry2D->GetImageGeometry()); geometry2D->UnRegister(); } void mitk::SlicedGeometry3D::InitializePlanes(const mitk::BaseGeometry *geometry3D, mitk::PlaneGeometry::PlaneOrientation planeorientation, bool top, bool frontside, bool rotated) { m_ReferenceGeometry = geometry3D; PlaneGeometry::Pointer planeGeometry = mitk::PlaneGeometry::New(); planeGeometry->InitializeStandardPlane(geometry3D, top, planeorientation, frontside, rotated); - ScalarType viewSpacing = 1; - unsigned int slices = 1; - - switch (planeorientation) - { - case PlaneGeometry::None: - break; - - case PlaneGeometry::Axial: - viewSpacing = geometry3D->GetSpacing()[2]; - slices = (unsigned int)geometry3D->GetExtent(2); - break; - - case PlaneGeometry::Frontal: - viewSpacing = geometry3D->GetSpacing()[1]; - slices = (unsigned int)geometry3D->GetExtent(1); - break; - - case PlaneGeometry::Sagittal: - viewSpacing = geometry3D->GetSpacing()[0]; - slices = (unsigned int)geometry3D->GetExtent(0); - break; - - default: - itkExceptionMacro("unknown PlaneOrientation"); - } - - mitk::Vector3D normal = this->AdjustNormal(planeGeometry->GetNormal()); - - ScalarType directedExtent = std::abs(m_ReferenceGeometry->GetExtentInMM(0) * normal[0]) + - std::abs(m_ReferenceGeometry->GetExtentInMM(1) * normal[1]) + - std::abs(m_ReferenceGeometry->GetExtentInMM(2) * normal[2]); - - if (directedExtent >= viewSpacing) - { - slices = static_cast(directedExtent / viewSpacing + 0.5); - } - else - { - slices = 1; - } - - bool flipped = (top == false); - - if (frontside == false) - { - flipped = !flipped; - } - if (planeorientation == PlaneGeometry::Frontal) - { - flipped = !flipped; - } - - this->InitializeEvenlySpaced(planeGeometry, viewSpacing, slices, flipped); + int worldAxis = + planeorientation == PlaneGeometry::Sagittal ? 0 : + planeorientation == PlaneGeometry::Frontal ? 1 : 2; + + // Inspired by: + // http://www.na-mic.org/Wiki/index.php/Coordinate_System_Conversion_Between_ITK_and_Slicer3 + + mitk::AffineTransform3D::MatrixType matrix = geometry3D->GetIndexToWorldTransform()->GetMatrix(); + matrix.GetVnlMatrix().normalize_columns(); + mitk::AffineTransform3D::MatrixType::InternalMatrixType inverseMatrix = matrix.GetInverse(); + + int dominantAxis = itk::Function::Max3( + inverseMatrix[0][worldAxis], + inverseMatrix[1][worldAxis], + inverseMatrix[2][worldAxis]); + + ScalarType viewSpacing = geometry3D->GetSpacing()[dominantAxis]; + unsigned int slices = static_cast(geometry3D->GetExtent(dominantAxis)); + +#ifndef NDEBUG + int upDirection = itk::Function::Sign(inverseMatrix[dominantAxis][worldAxis]); + + /// The normal vector of an imaginary plane that points from the world origin (bottom left back + /// corner or the world, with the lowest physical coordinates) towards the inside of the volume, + /// along the renderer axis. Length is the slice thickness. + Vector3D worldPlaneNormal = inverseMatrix.get_row(dominantAxis) * (upDirection * viewSpacing); + + /// The normal of the standard plane geometry just created. + Vector3D standardPlaneNormal = planeGeometry->GetNormal(); + + /// The standard plane must be parallel to the 'world plane'. The normal of the standard plane + /// must point against the world plane if and only if 'top' is 'false'. The length of the + /// standard plane normal must be equal to the slice thickness. + assert((standardPlaneNormal - (top ? 1.0 : -1.0) * worldPlaneNormal).GetSquaredNorm() < 0.000001); +#endif + + this->InitializeEvenlySpaced(planeGeometry, viewSpacing, slices); + +#ifndef NDEBUG + /// The standard plane normal and the z axis vector of the sliced geometry must point in + /// the same direction. + Vector3D zAxisVector = this->GetAxisVector(2); + Vector3D upscaledStandardPlaneNormal = standardPlaneNormal; + upscaledStandardPlaneNormal *= slices; + assert((zAxisVector - upscaledStandardPlaneNormal).GetSquaredNorm() < 0.000001); + + /// You can use this test is to check the handedness of the coordinate system of the current + /// geometry. In principle, you can use either left- or right-handed coordinate systems, but + /// you normally want it to be consistent, that is the handedness should be the same across + /// the renderers of the same viewer. +// ScalarType det = vnl_det(this->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix()); +// MITK_DEBUG << "world axis: " << worldAxis << (det > 0 ? " ; right-handed" : " ; left-handed"); +#endif } void mitk::SlicedGeometry3D::ReinitializePlanes(const Point3D ¢er, const Point3D &referencePoint) { // Need a reference frame to align the rotated planes if (!m_ReferenceGeometry) { return; } // Get first plane of plane stack PlaneGeometry *firstPlane = m_PlaneGeometries[0]; // If plane stack is empty, exit if (!firstPlane || dynamic_cast(firstPlane)) { return; } // Calculate the "directed" spacing when taking the plane (defined by its axes // vectors and normal) as the reference coordinate frame. // // This is done by calculating the radius of the ellipsoid defined by the // original volume spacing axes, in the direction of the respective axis of the // reference frame. mitk::Vector3D axis0 = firstPlane->GetAxisVector(0); mitk::Vector3D axis1 = firstPlane->GetAxisVector(1); mitk::Vector3D normal = firstPlane->GetNormal(); normal.Normalize(); Vector3D spacing; spacing[0] = this->CalculateSpacing(axis0); spacing[1] = this->CalculateSpacing(axis1); spacing[2] = this->CalculateSpacing(normal); Superclass::SetSpacing(spacing); // Now we need to calculate the number of slices in the plane's normal // direction, so that the entire volume is covered. This is done by first // calculating the dot product between the volume diagonal (the maximum // distance inside the volume) and the normal, and dividing this value by // the directed spacing calculated above. ScalarType directedExtent = std::abs(m_ReferenceGeometry->GetExtentInMM(0) * normal[0]) + std::abs(m_ReferenceGeometry->GetExtentInMM(1) * normal[1]) + std::abs(m_ReferenceGeometry->GetExtentInMM(2) * normal[2]); if (directedExtent >= spacing[2]) { m_Slices = static_cast(directedExtent / spacing[2] + 0.5); } else { m_Slices = 1; } // The origin of our "first plane" needs to be adapted to this new extent. // To achieve this, we first calculate the current distance to the volume's // center, and then shift the origin in the direction of the normal by the // difference between this distance and half of the new extent. double centerOfRotationDistance = firstPlane->SignedDistanceFromPlane(center); if (centerOfRotationDistance > 0) { firstPlane->SetOrigin(firstPlane->GetOrigin() + normal * (centerOfRotationDistance - directedExtent / 2.0)); m_DirectionVector = normal; } else { firstPlane->SetOrigin(firstPlane->GetOrigin() + normal * (directedExtent / 2.0 + centerOfRotationDistance)); m_DirectionVector = -normal; } // Now we adjust this distance according with respect to the given reference // point: we need to make sure that the point is touched by one slice of the // new slice stack. double referencePointDistance = firstPlane->SignedDistanceFromPlane(referencePoint); int referencePointSlice = static_cast(referencePointDistance / spacing[2]); double alignmentValue = referencePointDistance / spacing[2] - referencePointSlice; firstPlane->SetOrigin(firstPlane->GetOrigin() + normal * alignmentValue * spacing[2]); // Finally, we can clear the previous geometry stack and initialize it with // our re-initialized "first plane". m_PlaneGeometries.assign(m_Slices, PlaneGeometry::Pointer(nullptr)); if (m_Slices > 0) { m_PlaneGeometries[0] = firstPlane; } // Reinitialize SNC with new number of slices m_SliceNavigationController->GetSlice()->SetSteps(m_Slices); this->Modified(); } double mitk::SlicedGeometry3D::CalculateSpacing(const mitk::Vector3D &d) const { // Need the spacing of the underlying dataset / geometry if (!m_ReferenceGeometry) { return 1.0; } const mitk::Vector3D &spacing = m_ReferenceGeometry->GetSpacing(); return SlicedGeometry3D::CalculateSpacing(spacing, d); } double mitk::SlicedGeometry3D::CalculateSpacing(const mitk::Vector3D &spacing, const mitk::Vector3D &d) { // The following can be derived from the ellipsoid equation // // 1 = x^2/a^2 + y^2/b^2 + z^2/c^2 // // where (a,b,c) = spacing of original volume (ellipsoid radii) // and (x,y,z) = scaled coordinates of vector d (according to ellipsoid) // double scaling = d[0] * d[0] / (spacing[0] * spacing[0]) + d[1] * d[1] / (spacing[1] * spacing[1]) + d[2] * d[2] / (spacing[2] * spacing[2]); scaling = sqrt(scaling); return (sqrt(d[0] * d[0] + d[1] * d[1] + d[2] * d[2]) / scaling); } mitk::Vector3D mitk::SlicedGeometry3D::AdjustNormal(const mitk::Vector3D &normal) const { TransformType::Pointer inverse = TransformType::New(); m_ReferenceGeometry->GetIndexToWorldTransform()->GetInverse(inverse); Vector3D transformedNormal = inverse->TransformVector(normal); transformedNormal.Normalize(); return transformedNormal; } void mitk::SlicedGeometry3D::SetImageGeometry(const bool isAnImageGeometry) { Superclass::SetImageGeometry(isAnImageGeometry); unsigned int s; for (s = 0; s < m_Slices; ++s) { mitk::BaseGeometry *geometry = m_PlaneGeometries[s]; if (geometry != nullptr) { geometry->SetImageGeometry(isAnImageGeometry); } } } void mitk::SlicedGeometry3D::ChangeImageGeometryConsideringOriginOffset(const bool isAnImageGeometry) { unsigned int s; for (s = 0; s < m_Slices; ++s) { mitk::BaseGeometry *geometry = m_PlaneGeometries[s]; if (geometry != nullptr) { geometry->ChangeImageGeometryConsideringOriginOffset(isAnImageGeometry); } } Superclass::ChangeImageGeometryConsideringOriginOffset(isAnImageGeometry); } bool mitk::SlicedGeometry3D::IsValidSlice(int s) const { return ((s >= 0) && (s < (int)m_Slices)); } const mitk::BaseGeometry *mitk::SlicedGeometry3D::GetReferenceGeometry() const { return m_ReferenceGeometry; } void mitk::SlicedGeometry3D::SetReferenceGeometry(const BaseGeometry *referenceGeometry) { m_ReferenceGeometry = referenceGeometry; std::vector::iterator it; for (it = m_PlaneGeometries.begin(); it != m_PlaneGeometries.end(); ++it) { (*it)->SetReferenceGeometry(referenceGeometry); } } bool mitk::SlicedGeometry3D::HasReferenceGeometry() const { return ( m_ReferenceGeometry != nullptr ); } void mitk::SlicedGeometry3D::PreSetSpacing(const mitk::Vector3D &aSpacing) { bool hasEvenlySpacedPlaneGeometry = false; mitk::Point3D origin; mitk::Vector3D rightDV, bottomDV; BoundingBox::BoundsArrayType bounds; // Check for valid spacing if (!(aSpacing[0] > 0 && aSpacing[1] > 0 && aSpacing[2] > 0)) { mitkThrow() << "You try to set a spacing with at least one element equal or " "smaller to \"0\". This might lead to a crash during rendering. Please double" " check your data!"; } // In case of evenly-spaced data: re-initialize instances of PlaneGeometry, // since the spacing influences them if ((m_EvenlySpaced) && (m_PlaneGeometries.size() > 0)) { const PlaneGeometry *planeGeometry = m_PlaneGeometries[0]; if (planeGeometry && !dynamic_cast(planeGeometry)) { this->WorldToIndex(planeGeometry->GetOrigin(), origin); this->WorldToIndex(planeGeometry->GetAxisVector(0), rightDV); this->WorldToIndex(planeGeometry->GetAxisVector(1), bottomDV); bounds = planeGeometry->GetBounds(); hasEvenlySpacedPlaneGeometry = true; } } BaseGeometry::_SetSpacing(aSpacing); mitk::PlaneGeometry::Pointer firstGeometry; // In case of evenly-spaced data: re-initialize instances of PlaneGeometry, // since the spacing influences them if (hasEvenlySpacedPlaneGeometry) { // create planeGeometry according to new spacing this->IndexToWorld(origin, origin); this->IndexToWorld(rightDV, rightDV); this->IndexToWorld(bottomDV, bottomDV); mitk::PlaneGeometry::Pointer planeGeometry = mitk::PlaneGeometry::New(); planeGeometry->SetImageGeometry(this->GetImageGeometry()); planeGeometry->SetReferenceGeometry(m_ReferenceGeometry); // Store spacing, as Initialize... needs a pointer mitk::Vector3D lokalSpacing = this->GetSpacing(); planeGeometry->InitializeStandardPlane(rightDV.GetVnlVector(), bottomDV.GetVnlVector(), &lokalSpacing); planeGeometry->SetOrigin(origin); planeGeometry->SetBounds(bounds); firstGeometry = planeGeometry; } else if ((m_EvenlySpaced) && (m_PlaneGeometries.size() > 0)) { firstGeometry = m_PlaneGeometries[0].GetPointer(); } // clear and reserve PlaneGeometry::Pointer gnull = nullptr; m_PlaneGeometries.assign(m_Slices, gnull); if (m_Slices > 0) { m_PlaneGeometries[0] = firstGeometry; } this->Modified(); } void mitk::SlicedGeometry3D::SetSliceNavigationController(SliceNavigationController *snc) { m_SliceNavigationController = snc; } mitk::SliceNavigationController *mitk::SlicedGeometry3D::GetSliceNavigationController() { return m_SliceNavigationController; } void mitk::SlicedGeometry3D::SetEvenlySpaced(bool on) { if (m_EvenlySpaced != on) { m_EvenlySpaced = on; this->Modified(); } } void mitk::SlicedGeometry3D::SetDirectionVector(const mitk::Vector3D &directionVector) { Vector3D newDir = directionVector; newDir.Normalize(); if (newDir != m_DirectionVector) { m_DirectionVector = newDir; this->Modified(); } } // void // mitk::SlicedGeometry3D::SetTimeBounds( const mitk::TimeBounds& timebounds ) //{ // Superclass::SetTimeBounds( timebounds ); // // unsigned int s; // for ( s = 0; s < m_Slices; ++s ) // { // if(m_Geometry2Ds[s].IsNotNull()) // { // m_Geometry2Ds[s]->SetTimeBounds( timebounds ); // } // } // m_TimeBounds = timebounds; //} itk::LightObject::Pointer mitk::SlicedGeometry3D::InternalClone() const { Self::Pointer newGeometry = new SlicedGeometry3D(*this); newGeometry->UnRegister(); return newGeometry.GetPointer(); } void mitk::SlicedGeometry3D::PrintSelf(std::ostream &os, itk::Indent indent) const { Superclass::PrintSelf(os, indent); os << indent << " EvenlySpaced: " << m_EvenlySpaced << std::endl; if (m_EvenlySpaced) { os << indent << " DirectionVector: " << m_DirectionVector << std::endl; } os << indent << " Slices: " << m_Slices << std::endl; os << std::endl; os << indent << " GetPlaneGeometry(0): "; if (this->GetPlaneGeometry(0) == nullptr) { os << "NULL" << std::endl; } else { this->GetPlaneGeometry(0)->Print(os, indent); } } void mitk::SlicedGeometry3D::ExecuteOperation(Operation *operation) { PlaneGeometry::Pointer geometry2D; ApplyTransformMatrixOperation *applyMatrixOp; Point3D center; switch (operation->GetOperationType()) { case OpNOTHING: break; case OpROTATE: if (m_EvenlySpaced) { // Need a reference frame to align the rotation if (m_ReferenceGeometry) { // Clear all generated geometries and then rotate only the first slice. // The other slices will be re-generated on demand // Save first slice PlaneGeometry::Pointer geometry2D = m_PlaneGeometries[0]; RotationOperation *rotOp = dynamic_cast(operation); // Generate a RotationOperation using the dataset center instead of // the supplied rotation center. This is necessary so that the rotated // zero-plane does not shift away. The supplied center is instead used // to adjust the slice stack afterwards. Point3D center = m_ReferenceGeometry->GetCenter(); RotationOperation centeredRotation( rotOp->GetOperationType(), center, rotOp->GetVectorOfRotation(), rotOp->GetAngleOfRotation()); // Rotate first slice geometry2D->ExecuteOperation(¢eredRotation); // Clear the slice stack and adjust it according to the center of // the dataset and the supplied rotation center (see documentation of // ReinitializePlanes) this->ReinitializePlanes(center, rotOp->GetCenterOfRotation()); geometry2D->SetSpacing(this->GetSpacing()); if (m_SliceNavigationController) { m_SliceNavigationController->SelectSliceByPoint(rotOp->GetCenterOfRotation()); m_SliceNavigationController->AdjustSliceStepperRange(); } BaseGeometry::ExecuteOperation(¢eredRotation); } else { // we also have to consider the case, that there is no reference geometry available. if (m_PlaneGeometries.size() > 0) { // Reach through to all slices in my container for (auto iter = m_PlaneGeometries.begin(); iter != m_PlaneGeometries.end(); ++iter) { // Test for empty slices, which can happen if evenly spaced geometry if ((*iter).IsNotNull()) { (*iter)->ExecuteOperation(operation); } } // rotate overall geometry RotationOperation *rotOp = dynamic_cast(operation); BaseGeometry::ExecuteOperation(rotOp); } } } else { // Reach through to all slices for (auto iter = m_PlaneGeometries.begin(); iter != m_PlaneGeometries.end(); ++iter) { (*iter)->ExecuteOperation(operation); } } break; case OpORIENT: if (m_EvenlySpaced) { // get operation data PlaneOperation *planeOp = dynamic_cast(operation); // Get first slice PlaneGeometry::Pointer planeGeometry = m_PlaneGeometries[0]; // Need a PlaneGeometry, a PlaneOperation and a reference frame to // carry out the re-orientation. If not all avaialble, stop here if (!m_ReferenceGeometry || (!planeGeometry || dynamic_cast(planeGeometry.GetPointer())) || !planeOp) { break; } // General Behavior: // Clear all generated geometries and then rotate only the first slice. // The other slices will be re-generated on demand // // 1st Step: Reorient Normal Vector of first plane // Point3D center = planeOp->GetPoint(); // m_ReferenceGeometry->GetCenter(); mitk::Vector3D currentNormal = planeGeometry->GetNormal(); mitk::Vector3D newNormal; if (planeOp->AreAxisDefined()) { // If planeOp was defined by one centerpoint and two axis vectors newNormal = CrossProduct(planeOp->GetAxisVec0(), planeOp->GetAxisVec1()); } else { // If planeOp was defined by one centerpoint and one normal vector newNormal = planeOp->GetNormal(); } // Get Rotation axis und angle currentNormal.Normalize(); newNormal.Normalize(); ScalarType rotationAngle = angle(currentNormal.GetVnlVector(), newNormal.GetVnlVector()); rotationAngle *= 180.0 / vnl_math::pi; // from rad to deg Vector3D rotationAxis = itk::CrossProduct(currentNormal, newNormal); if (std::abs(rotationAngle - 180) < mitk::eps) { // current Normal and desired normal are not linear independent!!(e.g 1,0,0 and -1,0,0). // Rotation Axis should be ANY vector that is 90� to current Normal mitk::Vector3D helpNormal; helpNormal = currentNormal; helpNormal[0] += 1; helpNormal[1] -= 1; helpNormal[2] += 1; helpNormal.Normalize(); rotationAxis = itk::CrossProduct(helpNormal, currentNormal); } RotationOperation centeredRotation(mitk::OpROTATE, center, rotationAxis, rotationAngle); // Rotate first slice planeGeometry->ExecuteOperation(¢eredRotation); // Reinitialize planes and select slice, if my rotations are all done. if (!planeOp->AreAxisDefined()) { // Clear the slice stack and adjust it according to the center of // rotation and plane position (see documentation of ReinitializePlanes) this->ReinitializePlanes(center, planeOp->GetPoint()); planeGeometry->SetSpacing(this->GetSpacing()); if (m_SliceNavigationController) { m_SliceNavigationController->SelectSliceByPoint(planeOp->GetPoint()); m_SliceNavigationController->AdjustSliceStepperRange(); } } // Also apply rotation on the slicedGeometry - Geometry3D (Bounding geometry) BaseGeometry::ExecuteOperation(¢eredRotation); // // 2nd step. If axis vectors were defined, rotate the plane around its normal to fit these // if (planeOp->AreAxisDefined()) { mitk::Vector3D vecAxixNew = planeOp->GetAxisVec0(); vecAxixNew.Normalize(); mitk::Vector3D VecAxisCurr = planeGeometry->GetAxisVector(0); VecAxisCurr.Normalize(); ScalarType rotationAngle = angle(VecAxisCurr.GetVnlVector(), vecAxixNew.GetVnlVector()); rotationAngle = rotationAngle * 180 / PI; // Rad to Deg // we rotate around the normal of the plane, but we do not know, if we need to rotate clockwise // or anti-clockwise. So we rotate around the crossproduct of old and new Axisvector. // Since both axis vectors lie in the plane, the crossproduct is the planes normal or the negative planes // normal rotationAxis = itk::CrossProduct(VecAxisCurr, vecAxixNew); if (std::abs(rotationAngle - 180) < mitk::eps) { // current axisVec and desired axisVec are not linear independent!!(e.g 1,0,0 and -1,0,0). // Rotation Axis can be just plane Normal. (have to rotate by 180�) rotationAxis = newNormal; } // Perfom Rotation mitk::RotationOperation op(mitk::OpROTATE, center, rotationAxis, rotationAngle); planeGeometry->ExecuteOperation(&op); // Apply changes on first slice to whole slice stack this->ReinitializePlanes(center, planeOp->GetPoint()); planeGeometry->SetSpacing(this->GetSpacing()); if (m_SliceNavigationController) { m_SliceNavigationController->SelectSliceByPoint(planeOp->GetPoint()); m_SliceNavigationController->AdjustSliceStepperRange(); } // Also apply rotation on the slicedGeometry - Geometry3D (Bounding geometry) BaseGeometry::ExecuteOperation(&op); } } else { // Reach through to all slices for (auto iter = m_PlaneGeometries.begin(); iter != m_PlaneGeometries.end(); ++iter) { (*iter)->ExecuteOperation(operation); } } break; case OpRESTOREPLANEPOSITION: if (m_EvenlySpaced) { // Save first slice PlaneGeometry::Pointer planeGeometry = m_PlaneGeometries[0]; RestorePlanePositionOperation *restorePlaneOp = dynamic_cast(operation); // Need a PlaneGeometry, a PlaneOperation and a reference frame to // carry out the re-orientation if (m_ReferenceGeometry && (planeGeometry && dynamic_cast(planeGeometry.GetPointer()) == nullptr) && restorePlaneOp) { // Clear all generated geometries and then rotate only the first slice. // The other slices will be re-generated on demand // Rotate first slice planeGeometry->ExecuteOperation(restorePlaneOp); m_DirectionVector = restorePlaneOp->GetDirectionVector(); double centerOfRotationDistance = planeGeometry->SignedDistanceFromPlane(m_ReferenceGeometry->GetCenter()); if (centerOfRotationDistance <= 0) { m_DirectionVector = -m_DirectionVector; } Vector3D spacing = restorePlaneOp->GetSpacing(); Superclass::SetSpacing(spacing); // /*Now we need to calculate the number of slices in the plane's normal // direction, so that the entire volume is covered. This is done by first // calculating the dot product between the volume diagonal (the maximum // distance inside the volume) and the normal, and dividing this value by // the directed spacing calculated above.*/ ScalarType directedExtent = std::abs(m_ReferenceGeometry->GetExtentInMM(0) * m_DirectionVector[0]) + std::abs(m_ReferenceGeometry->GetExtentInMM(1) * m_DirectionVector[1]) + std::abs(m_ReferenceGeometry->GetExtentInMM(2) * m_DirectionVector[2]); if (directedExtent >= spacing[2]) { m_Slices = static_cast(directedExtent / spacing[2] + 0.5); } else { m_Slices = 1; } m_PlaneGeometries.assign(m_Slices, PlaneGeometry::Pointer(nullptr)); if (m_Slices > 0) { m_PlaneGeometries[0] = planeGeometry; } m_SliceNavigationController->GetSlice()->SetSteps(m_Slices); this->Modified(); // End Reinitialization if (m_SliceNavigationController) { m_SliceNavigationController->GetSlice()->SetPos(restorePlaneOp->GetPos()); m_SliceNavigationController->AdjustSliceStepperRange(); } BaseGeometry::ExecuteOperation(restorePlaneOp); } } else { // Reach through to all slices for (auto iter = m_PlaneGeometries.begin(); iter != m_PlaneGeometries.end(); ++iter) { (*iter)->ExecuteOperation(operation); } } break; case OpAPPLYTRANSFORMMATRIX: // Clear all generated geometries and then transform only the first slice. // The other slices will be re-generated on demand // Save first slice geometry2D = m_PlaneGeometries[0]; applyMatrixOp = dynamic_cast(operation); // Apply transformation to first plane geometry2D->ExecuteOperation(applyMatrixOp); // Generate a ApplyTransformMatrixOperation using the dataset center instead of // the supplied rotation center. The supplied center is instead used to adjust the // slice stack afterwards (see OpROTATE). center = m_ReferenceGeometry->GetCenter(); // Clear the slice stack and adjust it according to the center of // the dataset and the supplied rotation center (see documentation of // ReinitializePlanes) this->ReinitializePlanes(center, applyMatrixOp->GetReferencePoint()); BaseGeometry::ExecuteOperation(applyMatrixOp); break; default: // let handle by base class if we don't do anything BaseGeometry::ExecuteOperation(operation); } this->Modified(); } diff --git a/Modules/Core/test/mitkPlaneGeometryTest.cpp b/Modules/Core/test/mitkPlaneGeometryTest.cpp index 9e3b222181..f1ee6e0b5e 100644 --- a/Modules/Core/test/mitkPlaneGeometryTest.cpp +++ b/Modules/Core/test/mitkPlaneGeometryTest.cpp @@ -1,1078 +1,1078 @@ /*=================================================================== 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 "mitkAffineTransform3D.h" #include "mitkBaseGeometry.h" #include "mitkGeometry3D.h" #include "mitkInteractionConst.h" #include "mitkLine.h" #include "mitkPlaneGeometry.h" #include "mitkRotationOperation.h" #include "mitkSlicedGeometry3D.h" #include "mitkThinPlateSplineCurvedGeometry.h" #include #include #include #include #include #include static const mitk::ScalarType testEps = 1E-9; // the epsilon used in this test == at least float precision. class mitkPlaneGeometryTestSuite : public mitk::TestFixture { CPPUNIT_TEST_SUITE(mitkPlaneGeometryTestSuite); MITK_TEST(TestInitializeStandardPlane); MITK_TEST(TestProjectPointOntoPlane); MITK_TEST(TestPlaneGeometryCloning); MITK_TEST(TestInheritance); MITK_TEST(TestSetExtendInMM); MITK_TEST(TestRotate); MITK_TEST(TestClone); MITK_TEST(TestPlaneComparison); MITK_TEST(TestAxialInitialization); MITK_TEST(TestFrontalInitialization); MITK_TEST(TestSaggitalInitialization); MITK_TEST(TestLefthandedCoordinateSystem); // Currently commented out, see See bug 15990 // MITK_TEST(testPlaneGeometryInitializeOrder); MITK_TEST(TestIntersectionPoint); MITK_TEST(TestCase1210); CPPUNIT_TEST_SUITE_END(); private: // private test members that are initialized by setUp() mitk::PlaneGeometry::Pointer planegeometry; mitk::Point3D origin; - mitk::Vector3D right, bottom, normal; + mitk::Vector3D right, bottom, normal, spacing; mitk::ScalarType width, height; mitk::ScalarType widthInMM, heightInMM, thicknessInMM; public: void setUp() override { planegeometry = mitk::PlaneGeometry::New(); width = 100; widthInMM = width; height = 200; heightInMM = height; thicknessInMM = 1.0; mitk::FillVector3D(origin, 4.5, 7.3, 11.2); mitk::FillVector3D(right, widthInMM, 0, 0); mitk::FillVector3D(bottom, 0, heightInMM, 0); mitk::FillVector3D(normal, 0, 0, thicknessInMM); + mitk::FillVector3D(spacing, 1.0, 1.0, thicknessInMM); - planegeometry->InitializeStandardPlane(right.GetVnlVector(), bottom.GetVnlVector()); + planegeometry->InitializeStandardPlane(right, bottom); planegeometry->SetOrigin(origin); + planegeometry->SetSpacing(spacing); } void tearDown() override {} // This test verifies inheritance behaviour, this test will fail if the behaviour changes in the future void TestInheritance() { mitk::PlaneGeometry::Pointer plane = mitk::PlaneGeometry::New(); mitk::Geometry3D::Pointer g3d = dynamic_cast(plane.GetPointer()); CPPUNIT_ASSERT_MESSAGE("Planegeometry should not be castable to Geometry 3D", g3d.IsNull()); mitk::BaseGeometry::Pointer base = dynamic_cast(plane.GetPointer()); CPPUNIT_ASSERT_MESSAGE("Planegeometry should be castable to BaseGeometry", base.IsNotNull()); g3d = mitk::Geometry3D::New(); base = dynamic_cast(g3d.GetPointer()); CPPUNIT_ASSERT_MESSAGE("Geometry3D should be castable to BaseGeometry", base.IsNotNull()); mitk::SlicedGeometry3D::Pointer sliced = mitk::SlicedGeometry3D::New(); g3d = dynamic_cast(sliced.GetPointer()); CPPUNIT_ASSERT_MESSAGE("SlicedGeometry3D should not be castable to Geometry3D", g3d.IsNull()); mitk::ThinPlateSplineCurvedGeometry::Pointer thin = mitk::ThinPlateSplineCurvedGeometry::New(); plane = dynamic_cast(thin.GetPointer()); CPPUNIT_ASSERT_MESSAGE("AbstractTransformGeometry should be castable to PlaneGeometry", plane.IsNotNull()); plane = mitk::PlaneGeometry::New(); mitk::AbstractTransformGeometry::Pointer atg = dynamic_cast(plane.GetPointer()); CPPUNIT_ASSERT_MESSAGE("PlaneGeometry should not be castable to AbstractTransofrmGeometry", atg.IsNull()); } void TestLefthandedCoordinateSystem() { /** * @brief This method tests InitializeStandardPlane() and IndexToWorld() * with a left-handed coordinate orientation or indexToWorldMatrix. * * Of course this test does not use standard Parameters, which are right-handed. * See also discussion of bug #11477: http://bugs.mitk.org/show_bug.cgi?id=11477 */ planegeometry = mitk::PlaneGeometry::New(); width = 100; widthInMM = 5; height = 200; heightInMM = 3; thicknessInMM = 1.0; mitk::FillVector3D(right, widthInMM, 0, 0); mitk::FillVector3D(bottom, 0, heightInMM, 0); // This one negative sign results in lefthanded coordinate orientation and det(matrix) < 0. mitk::FillVector3D(normal, 0, 0, -thicknessInMM); mitk::AffineTransform3D::Pointer transform = mitk::AffineTransform3D::New(); mitk::AffineTransform3D::MatrixType matrix; mitk::AffineTransform3D::MatrixType::InternalMatrixType &vnl_matrix = matrix.GetVnlMatrix(); vnl_matrix.set_column(0, right); vnl_matrix.set_column(1, bottom); vnl_matrix.set_column(2, normal); // making sure that we didn't screw up this special test case or else fail deadly: assert(vnl_determinant(vnl_matrix) < 0.0); transform->SetIdentity(); transform->SetMatrix(matrix); planegeometry->InitializeStandardPlane(width, height, transform); // Crux of the matter. CPPUNIT_ASSERT_MESSAGE( "Testing if IndexToWorldMatrix is correct after InitializeStandardPlane( width, height, transform ) ", mitk::MatrixEqualElementWise(planegeometry->GetIndexToWorldTransform()->GetMatrix(), matrix)); mitk::Point3D p_index; p_index[0] = 10.; p_index[1] = 10.; p_index[2] = 0.; mitk::Point3D p_world; mitk::Point3D p_expectedResult; p_expectedResult[0] = 50.; p_expectedResult[1] = 30.; p_expectedResult[2] = 0.; ((mitk::BaseGeometry::Pointer)planegeometry)->IndexToWorld(p_index, p_world); // Crux of the matter. CPPUNIT_ASSERT_MESSAGE("Testing if IndexToWorld(a,b) function works correctly with lefthanded matrix ", mitk::Equal(p_world, p_expectedResult, testEps)); } // See bug 1210 // Test does not use standard Parameters void TestCase1210() { mitk::PlaneGeometry::Pointer planegeometry = mitk::PlaneGeometry::New(); mitk::Point3D origin; mitk::Vector3D right, down, spacing; mitk::FillVector3D(origin, 4.5, 7.3, 11.2); mitk::FillVector3D(right, 1.015625, 1.015625, 1.1999969482421875); mitk::FillVector3D(down, 1.4012984643248170709237295832899161312802619418765e-45, 0, 0); mitk::FillVector3D(spacing, 0, 1.4713633875410579244699160624544119378442750389703e-43, 9.2806360452222355258639080851310540729807238879469e-32); std::cout << "Testing InitializeStandardPlane(rightVector, downVector, spacing = NULL): " << std::endl; CPPUNIT_ASSERT_NO_THROW(planegeometry->InitializeStandardPlane(right, down, &spacing)); /* std::cout << "Testing width, height and thickness (in units): "; if((mitk::Equal(planegeometry->GetExtent(0),width)==false) || (mitk::Equal(planegeometry->GetExtent(1),height)==false) || (mitk::Equal(planegeometry->GetExtent(2),1)==false) ) { std::cout<<"[FAILED]"<GetExtentInMM(0),widthInMM)==false) || (mitk::Equal(planegeometry->GetExtentInMM(1),heightInMM)==false) || (mitk::Equal(planegeometry->GetExtentInMM(2),thicknessInMM)==false) ) { std::cout<<"[FAILED]"< 0. * */ // Test does not use standard Parameters void TestIntersectionPoint() { // init plane with its parameter mitk::PlaneGeometry::Pointer myPlaneGeometry = mitk::PlaneGeometry::New(); mitk::Point3D origin; origin[0] = 0.0; origin[1] = 2.0; origin[2] = 0.0; mitk::Vector3D normal; normal[0] = 0.0; normal[1] = 1.0; normal[2] = 0.0; myPlaneGeometry->InitializePlane(origin, normal); // generate points and line for intersection testing // point distance of given line > 1 mitk::Point3D pointP1; pointP1[0] = 2.0; pointP1[1] = 1.0; pointP1[2] = 0.0; mitk::Point3D pointP2; pointP2[0] = 2.0; pointP2[1] = 4.0; pointP2[2] = 0.0; mitk::Vector3D lineDirection; lineDirection[0] = pointP2[0] - pointP1[0]; lineDirection[1] = pointP2[1] - pointP1[1]; lineDirection[2] = pointP2[2] - pointP1[2]; mitk::Line3D xingline(pointP1, lineDirection); mitk::Point3D calcXingPoint; myPlaneGeometry->IntersectionPoint(xingline, calcXingPoint); // point distance of given line < 1 mitk::Point3D pointP3; pointP3[0] = 2.0; pointP3[1] = 2.2; pointP3[2] = 0.0; mitk::Point3D pointP4; pointP4[0] = 2.0; pointP4[1] = 1.7; pointP4[2] = 0.0; mitk::Vector3D lineDirection2; lineDirection2[0] = pointP4[0] - pointP3[0]; lineDirection2[1] = pointP4[1] - pointP3[1]; lineDirection2[2] = pointP4[2] - pointP3[2]; mitk::Line3D xingline2(pointP3, lineDirection2); mitk::Point3D calcXingPoint2; myPlaneGeometry->IntersectionPoint(xingline2, calcXingPoint2); // intersection points must be the same CPPUNIT_ASSERT_MESSAGE("Failed to calculate Intersection Point", calcXingPoint == calcXingPoint2); } /** * @brief This method tests method ProjectPointOntoPlane. * * See also bug #3409. */ // Test does not use standard Parameters void TestProjectPointOntoPlane() { mitk::PlaneGeometry::Pointer myPlaneGeometry = mitk::PlaneGeometry::New(); // create normal mitk::Vector3D normal; normal[0] = 0.0; normal[1] = 0.0; normal[2] = 1.0; // create origin mitk::Point3D origin; origin[0] = -27.582859; origin[1] = 50; origin[2] = 200.27742; // initialize plane geometry myPlaneGeometry->InitializePlane(origin, normal); // output to descripe the test std::cout << "Testing PlaneGeometry according to bug #3409" << std::endl; std::cout << "Our normal is: " << normal << std::endl; std::cout << "So ALL projected points should have exactly the same z-value!" << std::endl; // create a number of points mitk::Point3D myPoints[5]; myPoints[0][0] = -27.582859; myPoints[0][1] = 50.00; myPoints[0][2] = 200.27742; myPoints[1][0] = -26.58662; myPoints[1][1] = 50.00; myPoints[1][2] = 200.19026; myPoints[2][0] = -26.58662; myPoints[2][1] = 50.00; myPoints[2][2] = 200.33124; myPoints[3][0] = 104.58662; myPoints[3][1] = 452.12313; myPoints[3][2] = 866.41236; myPoints[4][0] = -207.58662; myPoints[4][1] = 312.00; myPoints[4][2] = -300.12346; // project points onto plane mitk::Point3D myProjectedPoints[5]; for (unsigned int i = 0; i < 5; ++i) { myProjectedPoints[i] = myPlaneGeometry->ProjectPointOntoPlane(myPoints[i]); } // compare z-values with z-value of plane (should be equal) bool allPointsOnPlane = true; for (auto &myProjectedPoint : myProjectedPoints) { if (fabs(myProjectedPoint[2] - origin[2]) > mitk::sqrteps) { allPointsOnPlane = false; } } CPPUNIT_ASSERT_MESSAGE("All points lie not on the same plane", allPointsOnPlane); } void TestPlaneGeometryCloning() { mitk::PlaneGeometry::Pointer geometry2D = createPlaneGeometry(); try { mitk::PlaneGeometry::Pointer clone = geometry2D->Clone(); itk::Matrix matrix = clone->GetIndexToWorldTransform()->GetMatrix(); CPPUNIT_ASSERT_MESSAGE("Test if matrix element exists...", matrix[0][0] == 31); double origin = geometry2D->GetOrigin()[0]; CPPUNIT_ASSERT_MESSAGE("First Point of origin as expected...", mitk::Equal(origin, 8)); double spacing = geometry2D->GetSpacing()[0]; CPPUNIT_ASSERT_MESSAGE("First Point of spacing as expected...", mitk::Equal(spacing, 31)); } catch (...) { CPPUNIT_FAIL("Error during access on a member of cloned geometry"); } // direction [row] [coloum] MITK_TEST_OUTPUT(<< "Casting a rotated 2D ITK Image to a MITK Image and check if Geometry is still same"); } void TestPlaneGeometryInitializeOrder() { mitk::Vector3D mySpacing; mySpacing[0] = 31; mySpacing[1] = 0.1; mySpacing[2] = 5.4; mitk::Point3D myOrigin; myOrigin[0] = 8; myOrigin[1] = 9; myOrigin[2] = 10; mitk::AffineTransform3D::Pointer myTransform = mitk::AffineTransform3D::New(); itk::Matrix transMatrix; transMatrix.Fill(0); transMatrix[0][0] = 1; transMatrix[1][1] = 2; transMatrix[2][2] = 4; myTransform->SetMatrix(transMatrix); mitk::PlaneGeometry::Pointer geometry2D1 = mitk::PlaneGeometry::New(); geometry2D1->SetIndexToWorldTransform(myTransform); geometry2D1->SetSpacing(mySpacing); geometry2D1->SetOrigin(myOrigin); mitk::PlaneGeometry::Pointer geometry2D2 = mitk::PlaneGeometry::New(); geometry2D2->SetSpacing(mySpacing); geometry2D2->SetOrigin(myOrigin); geometry2D2->SetIndexToWorldTransform(myTransform); mitk::PlaneGeometry::Pointer geometry2D3 = mitk::PlaneGeometry::New(); geometry2D3->SetIndexToWorldTransform(myTransform); geometry2D3->SetSpacing(mySpacing); geometry2D3->SetOrigin(myOrigin); geometry2D3->SetIndexToWorldTransform(myTransform); CPPUNIT_ASSERT_MESSAGE("Origin of Geometry 1 matches that of Geometry 2.", mitk::Equal(geometry2D1->GetOrigin(), geometry2D2->GetOrigin())); CPPUNIT_ASSERT_MESSAGE("Origin of Geometry 1 match those of Geometry 3.", mitk::Equal(geometry2D1->GetOrigin(), geometry2D3->GetOrigin())); CPPUNIT_ASSERT_MESSAGE("Origin of Geometry 2 match those of Geometry 3.", mitk::Equal(geometry2D2->GetOrigin(), geometry2D3->GetOrigin())); CPPUNIT_ASSERT_MESSAGE("Spacing of Geometry 1 match those of Geometry 2.", mitk::Equal(geometry2D1->GetSpacing(), geometry2D2->GetSpacing())); CPPUNIT_ASSERT_MESSAGE("Spacing of Geometry 1 match those of Geometry 3.", mitk::Equal(geometry2D1->GetSpacing(), geometry2D3->GetSpacing())); CPPUNIT_ASSERT_MESSAGE("Spacing of Geometry 2 match those of Geometry 3.", mitk::Equal(geometry2D2->GetSpacing(), geometry2D3->GetSpacing())); CPPUNIT_ASSERT_MESSAGE("Transformation of Geometry 1 match those of Geometry 2.", compareMatrix(geometry2D1->GetIndexToWorldTransform()->GetMatrix(), geometry2D2->GetIndexToWorldTransform()->GetMatrix())); CPPUNIT_ASSERT_MESSAGE("Transformation of Geometry 1 match those of Geometry 3.", compareMatrix(geometry2D1->GetIndexToWorldTransform()->GetMatrix(), geometry2D3->GetIndexToWorldTransform()->GetMatrix())); CPPUNIT_ASSERT_MESSAGE("Transformation of Geometry 2 match those of Geometry 3.", compareMatrix(geometry2D2->GetIndexToWorldTransform()->GetMatrix(), geometry2D3->GetIndexToWorldTransform()->GetMatrix())); } void TestInitializeStandardPlane() { CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: width", mitk::Equal(planegeometry->GetExtent(0), width, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: height", mitk::Equal(planegeometry->GetExtent(1), height, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: depth", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: width in mm", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: heght in mm", mitk::Equal(planegeometry->GetExtentInMM(1), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: depth in mm", mitk::Equal(planegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: AxisVectorRight", mitk::Equal(planegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: AxisVectorBottom", mitk::Equal(planegeometry->GetAxisVector(1), bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with default Spacing: AxisVectorNormal", mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); mitk::Vector3D spacing; thicknessInMM = 1.5; normal.Normalize(); normal *= thicknessInMM; mitk::FillVector3D(spacing, 1.0, 1.0, thicknessInMM); planegeometry->InitializeStandardPlane(right.GetVnlVector(), bottom.GetVnlVector(), &spacing); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: width", mitk::Equal(planegeometry->GetExtent(0), width, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: height", mitk::Equal(planegeometry->GetExtent(1), height, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: depth", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: width in mm", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: height in mm", mitk::Equal(planegeometry->GetExtentInMM(1), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: depth in mm", mitk::Equal(planegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: AxisVectorRight", mitk::Equal(planegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: AxisVectorBottom", mitk::Equal(planegeometry->GetAxisVector(1), bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing correct Standard Plane initialization with custom Spacing: AxisVectorNormal", mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); ; } void TestSetExtendInMM() { normal.Normalize(); normal *= thicknessInMM; planegeometry->SetExtentInMM(2, thicknessInMM); CPPUNIT_ASSERT_MESSAGE("Testing SetExtentInMM(2, ...), querying by GetExtentInMM(2): ", mitk::Equal(planegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing SetExtentInMM(2, ...), querying by GetAxisVector(2) and comparing to normal: ", mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); planegeometry->SetOrigin(origin); CPPUNIT_ASSERT_MESSAGE("Testing SetOrigin", mitk::Equal(planegeometry->GetOrigin(), origin, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() after SetOrigin: Right", mitk::Equal(planegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() after SetOrigin: Bottom", mitk::Equal(planegeometry->GetAxisVector(1), bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() after SetOrigin: Normal", mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); mappingTests2D(planegeometry, width, height, widthInMM, heightInMM, origin, right, bottom); } void TestRotate() { // Changing the IndexToWorldTransform to a rotated version by SetIndexToWorldTransform() (keep origin): mitk::AffineTransform3D::Pointer transform = mitk::AffineTransform3D::New(); mitk::AffineTransform3D::MatrixType::InternalMatrixType vnlmatrix; vnlmatrix = planegeometry->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix(); mitk::VnlVector axis(3); mitk::FillVector3D(axis, 1.0, 1.0, 1.0); axis.normalize(); vnl_quaternion rotation(axis, 0.223); vnlmatrix = rotation.rotation_matrix_transpose() * vnlmatrix; mitk::Matrix3D matrix; matrix = vnlmatrix; transform->SetMatrix(matrix); transform->SetOffset(planegeometry->GetIndexToWorldTransform()->GetOffset()); right.SetVnlVector(rotation.rotation_matrix_transpose() * right.GetVnlVector()); bottom.SetVnlVector(rotation.rotation_matrix_transpose() * bottom.GetVnlVector()); normal.SetVnlVector(rotation.rotation_matrix_transpose() * normal.GetVnlVector()); planegeometry->SetIndexToWorldTransform(transform); // The origin changed,because m_Origin=m_IndexToWorldTransform->GetOffset()+GetAxisVector(2)*0.5 // and the AxisVector changes due to the rotation. In other words: the rotation was done around // the corner of the box, not around the planes origin. Now change it to a rotation around // the origin, simply by re-setting the origin to the original one: planegeometry->SetOrigin(origin); CPPUNIT_ASSERT_MESSAGE("Testing whether SetIndexToWorldTransform kept origin: ", mitk::Equal(planegeometry->GetOrigin(), origin, testEps)); mitk::Point2D point; point[0] = 4; point[1] = 3; mitk::Point2D dummy; planegeometry->WorldToIndex(point, dummy); planegeometry->IndexToWorld(dummy, dummy); CPPUNIT_ASSERT_MESSAGE("Testing consistency of index and world coordinates.", dummy == point); CPPUNIT_ASSERT_MESSAGE("Testing width of rotated version: ", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height of rotated version: ", mitk::Equal(planegeometry->GetExtentInMM(1), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing thickness of rotated version: ", mitk::Equal(planegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of rotated version: right ", mitk::Equal(planegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of rotated version: bottom", mitk::Equal(planegeometry->GetAxisVector(1), bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of rotated version: normal", mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing GetAxisVector(direction).GetNorm() != planegeometry->GetExtentInMM(direction) of rotated version: ", mitk::Equal(planegeometry->GetAxisVector(0).GetNorm(), planegeometry->GetExtentInMM(0), testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing GetAxisVector(direction).GetNorm() != planegeometry->GetExtentInMM(direction) of rotated version: ", mitk::Equal(planegeometry->GetAxisVector(1).GetNorm(), planegeometry->GetExtentInMM(1), testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing GetAxisVector(direction).GetNorm() != planegeometry->GetExtentInMM(direction) of rotated version: ", mitk::Equal(planegeometry->GetAxisVector(2).GetNorm(), planegeometry->GetExtentInMM(2), testEps)); mappingTests2D(planegeometry, width, height, widthInMM, heightInMM, origin, right, bottom); width *= 2; height *= 3; planegeometry->SetSizeInUnits(width, height); CPPUNIT_ASSERT_MESSAGE("Testing SetSizeInUnits() of rotated version: ", mitk::Equal(planegeometry->GetExtent(0), width, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing SetSizeInUnits() of rotated version: ", mitk::Equal(planegeometry->GetExtent(1), height, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing SetSizeInUnits() of rotated version: ", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width (in mm) of version with changed size in units: ", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height (in mm) of version with changed size in units: ", mitk::Equal(planegeometry->GetExtentInMM(1), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing thickness (in mm) of version with changed size in units: ", mitk::Equal(planegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of version with changed size in units: right ", mitk::Equal(planegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of version with changed size in units: bottom", mitk::Equal(planegeometry->GetAxisVector(1), bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of version with changed size in units: normal", mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing GetAxisVector(direction).GetNorm() != planegeometry->GetExtentInMM(direction) of rotated version: ", mitk::Equal(planegeometry->GetAxisVector(0).GetNorm(), planegeometry->GetExtentInMM(0), testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing GetAxisVector(direction).GetNorm() != planegeometry->GetExtentInMM(direction) of rotated version: ", mitk::Equal(planegeometry->GetAxisVector(1).GetNorm(), planegeometry->GetExtentInMM(1), testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing GetAxisVector(direction).GetNorm() != planegeometry->GetExtentInMM(direction) of rotated version: ", mitk::Equal(planegeometry->GetAxisVector(2).GetNorm(), planegeometry->GetExtentInMM(2), testEps)); mappingTests2D(planegeometry, width, height, widthInMM, heightInMM, origin, right, bottom); } void TestClone() { mitk::PlaneGeometry::Pointer clonedplanegeometry = dynamic_cast(planegeometry->Clone().GetPointer()); // Cave: Statement below is negated! CPPUNIT_ASSERT_MESSAGE("Testing Clone(): ", !((clonedplanegeometry.IsNull()) || (clonedplanegeometry->GetReferenceCount() != 1))); CPPUNIT_ASSERT_MESSAGE("Testing origin of cloned version: ", mitk::Equal(clonedplanegeometry->GetOrigin(), origin, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width (in units) of cloned version: ", mitk::Equal(clonedplanegeometry->GetExtent(0), width, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height (in units) of cloned version: ", mitk::Equal(clonedplanegeometry->GetExtent(1), height, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing extent (in units) of cloned version: ", mitk::Equal(clonedplanegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width (in mm) of cloned version: ", mitk::Equal(clonedplanegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height (in mm) of cloned version: ", mitk::Equal(clonedplanegeometry->GetExtentInMM(1), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing thickness (in mm) of cloned version: ", mitk::Equal(clonedplanegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of cloned version: right", mitk::Equal(clonedplanegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of cloned version: bottom", mitk::Equal(clonedplanegeometry->GetAxisVector(1), bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of cloned version: normal", mitk::Equal(clonedplanegeometry->GetAxisVector(2), normal, testEps)); mappingTests2D(clonedplanegeometry, width, height, widthInMM, heightInMM, origin, right, bottom); } void TestSaggitalInitialization() { mitk::Point3D cornerpoint0 = planegeometry->GetCornerPoint(0); - - mitk::PlaneGeometry::Pointer clonedplanegeometry = - dynamic_cast(planegeometry->Clone().GetPointer()); + mitk::PlaneGeometry::Pointer clonedplanegeometry = planegeometry->Clone(); // Testing InitializeStandardPlane(clonedplanegeometry, planeorientation = Sagittal, zPosition = 0, frontside=true): planegeometry->InitializeStandardPlane(clonedplanegeometry, mitk::PlaneGeometry::Sagittal); mitk::Vector3D newright, newbottom, newnormal; mitk::ScalarType newthicknessInMM; newright = bottom; newthicknessInMM = widthInMM / width * 1.0; // extent in normal direction is 1; newnormal = right; newnormal.Normalize(); newnormal *= newthicknessInMM; newbottom = normal; newbottom.Normalize(); newbottom *= thicknessInMM; CPPUNIT_ASSERT_MESSAGE("Testing GetCornerPoint(0) of sagitally initialized version:", mitk::Equal(planegeometry->GetCornerPoint(0), cornerpoint0, testEps)); // ok, corner was fine, so we can dare to believe the origin is ok. origin = planegeometry->GetOrigin(); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in units) of sagitally initialized version: ", mitk::Equal(planegeometry->GetExtent(0), height, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in units) of sagitally initialized version: ", mitk::Equal(planegeometry->GetExtent(1), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in units) of sagitally initialized version: ", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in mm) of sagitally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(0), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in mm) of sagitally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(1), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in mm) of sagitally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(2), newthicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of sagitally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(0), newright, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of sagitally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(1), newbottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of sagitally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(2), newnormal, testEps)); mappingTests2D(planegeometry, height, 1, heightInMM, thicknessInMM, origin, newright, newbottom); // set origin back to the one of the axial slice: origin = clonedplanegeometry->GetOrigin(); // Testing backside initialization: InitializeStandardPlane(clonedplanegeometry, planeorientation = Axial, zPosition // = 0, frontside=false, rotated=true): planegeometry->InitializeStandardPlane(clonedplanegeometry, mitk::PlaneGeometry::Axial, 0, false, true); mitk::Point3D backsideorigin; backsideorigin = origin + clonedplanegeometry->GetAxisVector(1); //+clonedplanegeometry->GetAxisVector(2); CPPUNIT_ASSERT_MESSAGE("Testing origin of backsidedly, axially initialized version: ", mitk::Equal(planegeometry->GetOrigin(), backsideorigin, testEps)); mitk::Point3D backsidecornerpoint0; backsidecornerpoint0 = cornerpoint0 + clonedplanegeometry->GetAxisVector(1); //+clonedplanegeometry->GetAxisVector(2); CPPUNIT_ASSERT_MESSAGE("Testing GetCornerPoint(0) of sagitally initialized version: ", mitk::Equal(planegeometry->GetCornerPoint(0), backsidecornerpoint0, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in units) of backsidedly, axially initialized version " "(should be same as in mm due to unit spacing, except for thickness, which is always 1): ", mitk::Equal(planegeometry->GetExtent(0), width, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in units) of backsidedly, axially initialized version " "(should be same as in mm due to unit spacing, except for thickness, which is always 1): ", mitk::Equal(planegeometry->GetExtent(1), height, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in units) of backsidedly, axially initialized version " "(should be same as in mm due to unit spacing, except for thickness, which is always 1): ", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in mm) of backsidedly, axially initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in mm) of backsidedly, axially initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(1), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width, height and thickness (in mm) of backsidedly, axially initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of backsidedly, axially initialized version: ", mitk::Equal(planegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of backsidedly, axially initialized version: ", mitk::Equal(planegeometry->GetAxisVector(1), -bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of backsidedly, axially initialized version: ", - mitk::Equal(planegeometry->GetAxisVector(2), -normal, testEps)); + mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); // T22254: Flipped sign mappingTests2D(planegeometry, width, height, widthInMM, heightInMM, backsideorigin, right, -bottom); } void TestFrontalInitialization() { mitk::Point3D cornerpoint0 = planegeometry->GetCornerPoint(0); mitk::PlaneGeometry::Pointer clonedplanegeometry = dynamic_cast(planegeometry->Clone().GetPointer()); //-------- mitk::Vector3D newright, newbottom, newnormal; mitk::ScalarType newthicknessInMM; // Testing InitializeStandardPlane(clonedplanegeometry, planeorientation = Frontal, zPosition = 0, frontside=true) planegeometry->InitializeStandardPlane(clonedplanegeometry, mitk::PlaneGeometry::Frontal); newright = right; newbottom = normal; newbottom.Normalize(); newbottom *= thicknessInMM; newthicknessInMM = heightInMM / height * 1.0 /*extent in normal direction is 1*/; newnormal = -bottom; newnormal.Normalize(); newnormal *= newthicknessInMM; CPPUNIT_ASSERT_MESSAGE("Testing GetCornerPoint(0) of frontally initialized version: ", mitk::Equal(planegeometry->GetCornerPoint(0), cornerpoint0, testEps)); // ok, corner was fine, so we can dare to believe the origin is ok. origin = planegeometry->GetOrigin(); CPPUNIT_ASSERT_MESSAGE("Testing width (in units) of frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(0), width, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height (in units) of frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(1), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing thickness (in units) of frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width (in mm) of frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height (in mm) of frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(1), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing thickness (in mm) of frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(2), newthicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of frontally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(0), newright, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of frontally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(1), newbottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of frontally initialized version: ", - mitk::Equal(planegeometry->GetAxisVector(2), newnormal, testEps)); + mitk::Equal(planegeometry->GetAxisVector(2), -newnormal, testEps)); // T22254: Flipped sign mappingTests2D(planegeometry, width, 1, widthInMM, thicknessInMM, origin, newright, newbottom); // Changing plane to in-plane unit spacing using SetSizeInUnits: planegeometry->SetSizeInUnits(planegeometry->GetExtentInMM(0), planegeometry->GetExtentInMM(1)); CPPUNIT_ASSERT_MESSAGE("Testing origin of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetOrigin(), origin, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in units) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in units) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(1), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in units) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in mm) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in mm) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(1), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in mm) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(2), newthicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(0), newright, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(1), newbottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of unit spaced, frontally initialized version: ", - mitk::Equal(planegeometry->GetAxisVector(2), newnormal, testEps)); + mitk::Equal(planegeometry->GetAxisVector(2), -newnormal, testEps)); // T22254: Flipped sign mappingTests2D(planegeometry, widthInMM, thicknessInMM, widthInMM, thicknessInMM, origin, newright, newbottom); // Changing plane to unit spacing also in normal direction using SetExtentInMM(2, 1.0): planegeometry->SetExtentInMM(2, 1.0); newnormal.Normalize(); CPPUNIT_ASSERT_MESSAGE("Testing origin of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetOrigin(), origin, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in units) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in units) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(1), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in units) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in mm) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in mm) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(1), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE( "Testing width, height and thickness (in mm) of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(2), 1.0, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(0), newright, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of unit spaced, frontally initialized version: ", mitk::Equal(planegeometry->GetAxisVector(1), newbottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of unit spaced, frontally initialized version: ", - mitk::Equal(planegeometry->GetAxisVector(2), newnormal, testEps)); + mitk::Equal(planegeometry->GetAxisVector(2), -newnormal, testEps)); // T22254: Flipped sign mappingTests2D(planegeometry, widthInMM, thicknessInMM, widthInMM, thicknessInMM, origin, newright, newbottom); } void TestAxialInitialization() { mitk::Point3D cornerpoint0 = planegeometry->GetCornerPoint(0); // Clone, move, rotate and test for 'IsParallel' and 'IsOnPlane' mitk::PlaneGeometry::Pointer clonedplanegeometry = dynamic_cast(planegeometry->Clone().GetPointer()); CPPUNIT_ASSERT_MESSAGE("Testing Clone(): ", !((clonedplanegeometry.IsNull()) || (clonedplanegeometry->GetReferenceCount() != 1))); std::cout << "Testing InitializeStandardPlane(clonedplanegeometry, planeorientation = Axial, zPosition = 0, " "frontside=true): " << std::endl; planegeometry->InitializeStandardPlane(clonedplanegeometry); CPPUNIT_ASSERT_MESSAGE("Testing origin of axially initialized version: ", mitk::Equal(planegeometry->GetOrigin(), origin)); CPPUNIT_ASSERT_MESSAGE("Testing GetCornerPoint(0) of axially initialized version: ", mitk::Equal(planegeometry->GetCornerPoint(0), cornerpoint0)); CPPUNIT_ASSERT_MESSAGE("Testing width (in units) of axially initialized version (should be same as in mm due to " "unit spacing, except for thickness, which is always 1): ", mitk::Equal(planegeometry->GetExtent(0), width, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height (in units) of axially initialized version (should be same as in mm due to " "unit spacing, except for thickness, which is always 1): ", mitk::Equal(planegeometry->GetExtent(1), height, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing thickness (in units) of axially initialized version (should be same as in mm due " "to unit spacing, except for thickness, which is always 1): ", mitk::Equal(planegeometry->GetExtent(2), 1, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing width (in mm) of axially initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(0), widthInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing height (in mm) of axially initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(1), heightInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing thickness (in mm) of axially initialized version: ", mitk::Equal(planegeometry->GetExtentInMM(2), thicknessInMM, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of axially initialized version: ", mitk::Equal(planegeometry->GetAxisVector(0), right, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of axially initialized version: ", mitk::Equal(planegeometry->GetAxisVector(1), bottom, testEps)); CPPUNIT_ASSERT_MESSAGE("Testing GetAxisVector() of axially initialized version: ", mitk::Equal(planegeometry->GetAxisVector(2), normal, testEps)); mappingTests2D(planegeometry, width, height, widthInMM, heightInMM, origin, right, bottom); } void TestPlaneComparison() { // Clone, move, rotate and test for 'IsParallel' and 'IsOnPlane' mitk::PlaneGeometry::Pointer clonedplanegeometry2 = dynamic_cast(planegeometry->Clone().GetPointer()); CPPUNIT_ASSERT_MESSAGE("Testing Clone(): ", !((clonedplanegeometry2.IsNull()) || (clonedplanegeometry2->GetReferenceCount() != 1))); CPPUNIT_ASSERT_MESSAGE("Testing wheter original and clone are at the same position", clonedplanegeometry2->IsOnPlane(planegeometry.GetPointer())); CPPUNIT_ASSERT_MESSAGE(" Asserting that origin is on the plane cloned plane:", clonedplanegeometry2->IsOnPlane(origin)); mitk::VnlVector newaxis(3); mitk::FillVector3D(newaxis, 1.0, 1.0, 1.0); newaxis.normalize(); vnl_quaternion rotation2(newaxis, 0.0); mitk::Vector3D clonednormal = clonedplanegeometry2->GetNormal(); mitk::Point3D clonedorigin = clonedplanegeometry2->GetOrigin(); auto planerot = new mitk::RotationOperation(mitk::OpROTATE, origin, clonedplanegeometry2->GetAxisVector(0), 180.0); clonedplanegeometry2->ExecuteOperation(planerot); CPPUNIT_ASSERT_MESSAGE(" Asserting that a flipped plane is still on the original plane: ", clonedplanegeometry2->IsOnPlane(planegeometry.GetPointer())); clonedorigin += clonednormal; clonedplanegeometry2->SetOrigin(clonedorigin); CPPUNIT_ASSERT_MESSAGE("Testing if the translated (cloned, flipped) plane is parallel to its origin plane: ", clonedplanegeometry2->IsParallel(planegeometry)); delete planerot; planerot = new mitk::RotationOperation(mitk::OpROTATE, origin, clonedplanegeometry2->GetAxisVector(0), 0.5); clonedplanegeometry2->ExecuteOperation(planerot); CPPUNIT_ASSERT_MESSAGE("Testing if a non-paralell plane gets recognized as not paralell [rotation +0.5 degree] : ", !clonedplanegeometry2->IsParallel(planegeometry)); delete planerot; planerot = new mitk::RotationOperation(mitk::OpROTATE, origin, clonedplanegeometry2->GetAxisVector(0), -1.0); clonedplanegeometry2->ExecuteOperation(planerot); CPPUNIT_ASSERT_MESSAGE("Testing if a non-paralell plane gets recognized as not paralell [rotation -0.5 degree] : ", !clonedplanegeometry2->IsParallel(planegeometry)); delete planerot; planerot = new mitk::RotationOperation(mitk::OpROTATE, origin, clonedplanegeometry2->GetAxisVector(0), 360.5); clonedplanegeometry2->ExecuteOperation(planerot); CPPUNIT_ASSERT_MESSAGE("Testing if a non-paralell plane gets recognized as paralell [rotation 360 degree] : ", clonedplanegeometry2->IsParallel(planegeometry)); } private: // helper Methods for the Tests mitk::PlaneGeometry::Pointer createPlaneGeometry() { mitk::Vector3D mySpacing; mySpacing[0] = 31; mySpacing[1] = 0.1; mySpacing[2] = 5.4; mitk::Point3D myOrigin; myOrigin[0] = 8; myOrigin[1] = 9; myOrigin[2] = 10; mitk::AffineTransform3D::Pointer myTransform = mitk::AffineTransform3D::New(); itk::Matrix transMatrix; transMatrix.Fill(0); transMatrix[0][0] = 1; transMatrix[1][1] = 2; transMatrix[2][2] = 4; myTransform->SetMatrix(transMatrix); mitk::PlaneGeometry::Pointer geometry2D = mitk::PlaneGeometry::New(); geometry2D->SetIndexToWorldTransform(myTransform); geometry2D->SetSpacing(mySpacing); geometry2D->SetOrigin(myOrigin); return geometry2D; } bool compareMatrix(itk::Matrix left, itk::Matrix right) { bool equal = true; for (int i = 0; i < 3; ++i) for (int j = 0; j < 3; ++j) equal &= mitk::Equal(left[i][j], right[i][j]); return equal; } /** * This function tests for correct mapping and is called several times from other tests **/ void mappingTests2D(const mitk::PlaneGeometry *planegeometry, const mitk::ScalarType &width, const mitk::ScalarType &height, const mitk::ScalarType &widthInMM, const mitk::ScalarType &heightInMM, const mitk::Point3D &origin, const mitk::Vector3D &right, const mitk::Vector3D &bottom) { std::cout << "Testing mapping Map(pt2d_mm(x=widthInMM/2.3,y=heightInMM/2.5), pt3d_mm) and compare with expected: "; mitk::Point2D pt2d_mm; mitk::Point3D pt3d_mm, expected_pt3d_mm; pt2d_mm[0] = widthInMM / 2.3; pt2d_mm[1] = heightInMM / 2.5; expected_pt3d_mm = origin + right * (pt2d_mm[0] / right.GetNorm()) + bottom * (pt2d_mm[1] / bottom.GetNorm()); planegeometry->Map(pt2d_mm, pt3d_mm); CPPUNIT_ASSERT_MESSAGE( "Testing mapping Map(pt2d_mm(x=widthInMM/2.3,y=heightInMM/2.5), pt3d_mm) and compare with expected", mitk::Equal(pt3d_mm, expected_pt3d_mm, testEps)); std::cout << "Testing mapping Map(pt3d_mm, pt2d_mm) and compare with expected: "; mitk::Point2D testpt2d_mm; planegeometry->Map(pt3d_mm, testpt2d_mm); std::cout << std::setprecision(12) << "Expected pt2d_mm " << pt2d_mm << std::endl; std::cout << std::setprecision(12) << "Result testpt2d_mm " << testpt2d_mm << std::endl; std::cout << std::setprecision(12) << "10*mitk::eps " << 10 * mitk::eps << std::endl; // This eps is temporarily set to 10*mitk::eps. See bug #15037 for details. CPPUNIT_ASSERT_MESSAGE("Testing mapping Map(pt3d_mm, pt2d_mm) and compare with expected", mitk::Equal(pt2d_mm, testpt2d_mm, 10 * mitk::eps)); std::cout << "Testing IndexToWorld(pt2d_units, pt2d_mm) and compare with expected: "; mitk::Point2D pt2d_units; pt2d_units[0] = width / 2.0; pt2d_units[1] = height / 2.0; pt2d_mm[0] = widthInMM / 2.0; pt2d_mm[1] = heightInMM / 2.0; planegeometry->IndexToWorld(pt2d_units, testpt2d_mm); std::cout << std::setprecision(12) << "Expected pt2d_mm " << pt2d_mm << std::endl; std::cout << std::setprecision(12) << "Result testpt2d_mm " << testpt2d_mm << std::endl; std::cout << std::setprecision(12) << "10*mitk::eps " << 10 * mitk::eps << std::endl; // This eps is temporarily set to 10*mitk::eps. See bug #15037 for details. CPPUNIT_ASSERT_MESSAGE("Testing IndexToWorld(pt2d_units, pt2d_mm) and compare with expected: ", mitk::Equal(pt2d_mm, testpt2d_mm, 10 * mitk::eps)); std::cout << "Testing WorldToIndex(pt2d_mm, pt2d_units) and compare with expected: "; mitk::Point2D testpt2d_units; planegeometry->WorldToIndex(pt2d_mm, testpt2d_units); std::cout << std::setprecision(12) << "Expected pt2d_units " << pt2d_units << std::endl; std::cout << std::setprecision(12) << "Result testpt2d_units " << testpt2d_units << std::endl; std::cout << std::setprecision(12) << "10*mitk::eps " << 10 * mitk::eps << std::endl; // This eps is temporarily set to 10*mitk::eps. See bug #15037 for details. CPPUNIT_ASSERT_MESSAGE("Testing WorldToIndex(pt2d_mm, pt2d_units) and compare with expected:", mitk::Equal(pt2d_units, testpt2d_units, 10 * mitk::eps)); } }; MITK_TEST_SUITE_REGISTRATION(mitkPlaneGeometry) diff --git a/Modules/Core/test/mitkSliceNavigationControllerTest.cpp b/Modules/Core/test/mitkSliceNavigationControllerTest.cpp index f0a8d0981d..78799c862a 100644 --- a/Modules/Core/test/mitkSliceNavigationControllerTest.cpp +++ b/Modules/Core/test/mitkSliceNavigationControllerTest.cpp @@ -1,649 +1,649 @@ /*=================================================================== 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 "mitkInteractionConst.h" #include "mitkPlaneGeometry.h" #include "mitkPlanePositionManager.h" #include "mitkRotationOperation.h" #include "mitkSliceNavigationController.h" #include "mitkSlicedGeometry3D.h" #include "mitkTestingMacros.h" #include "usGetModuleContext.h" #include "usModuleContext.h" #include "usServiceReference.h" #include #include #include bool operator==(const mitk::Geometry3D &left, const mitk::Geometry3D &right) { mitk::BoundingBox::BoundsArrayType leftbounds, rightbounds; leftbounds = left.GetBounds(); rightbounds = right.GetBounds(); unsigned int i; for (i = 0; i < 6; ++i) if (mitk::Equal(leftbounds[i], rightbounds[i]) == false) return false; const mitk::Geometry3D::TransformType::MatrixType &leftmatrix = left.GetIndexToWorldTransform()->GetMatrix(); const mitk::Geometry3D::TransformType::MatrixType &rightmatrix = right.GetIndexToWorldTransform()->GetMatrix(); unsigned int j; for (i = 0; i < 3; ++i) { const mitk::Geometry3D::TransformType::MatrixType::ValueType *leftvector = leftmatrix[i]; const mitk::Geometry3D::TransformType::MatrixType::ValueType *rightvector = rightmatrix[i]; for (j = 0; j < 3; ++j) if (mitk::Equal(leftvector[i], rightvector[i]) == false) return false; } const mitk::Geometry3D::TransformType::OffsetType &leftoffset = left.GetIndexToWorldTransform()->GetOffset(); const mitk::Geometry3D::TransformType::OffsetType &rightoffset = right.GetIndexToWorldTransform()->GetOffset(); for (i = 0; i < 3; ++i) if (mitk::Equal(leftoffset[i], rightoffset[i]) == false) return false; return true; } int compareGeometry(const mitk::TimeGeometry &timeGeometry, const mitk::ScalarType &width, const mitk::ScalarType &height, const mitk::ScalarType &numSlices, const mitk::ScalarType &widthInMM, const mitk::ScalarType &heightInMM, const mitk::ScalarType &thicknessInMM, const mitk::Point3D &cornerpoint0, const mitk::Vector3D &right, const mitk::Vector3D &bottom, const mitk::Vector3D &normal) { // Probleme durch umstellung von Time-SlicedGeometry auf TimeGeometry? // Eventuell gibt es keine Entsprechung mehr. const mitk::BaseGeometry::Pointer geometry = timeGeometry.GetGeometryForTimeStep(0); std::cout << "Testing width, height and thickness (in units): "; if ((mitk::Equal(geometry->GetExtent(0), width) == false) || (mitk::Equal(geometry->GetExtent(1), height) == false) || (mitk::Equal(geometry->GetExtent(2), numSlices) == false)) { std::cout << "[FAILED]" << std::endl; return EXIT_FAILURE; } std::cout << "[PASSED]" << std::endl; std::cout << "Testing width, height and thickness (in mm): "; if ((mitk::Equal(geometry->GetExtentInMM(0), widthInMM) == false) || (mitk::Equal(geometry->GetExtentInMM(1), heightInMM) == false) || (mitk::Equal(geometry->GetExtentInMM(2), thicknessInMM) == false)) { std::cout << "[FAILED]" << std::endl; return EXIT_FAILURE; } std::cout << "[PASSED]" << std::endl; std::cout << "Testing GetAxisVector(): "; std::cout << "dir=0 "; mitk::Vector3D dv; dv = right; dv.Normalize(); dv *= widthInMM; if ((mitk::Equal(geometry->GetAxisVector(0), dv) == false)) { std::cout << "[FAILED]" << std::endl; return EXIT_FAILURE; } std::cout << "[PASSED]"; std::cout << ", dir=1 "; dv = bottom; dv.Normalize(); dv *= heightInMM; if ((mitk::Equal(geometry->GetAxisVector(1), dv) == false)) { std::cout << "[FAILED]" << std::endl; return EXIT_FAILURE; } std::cout << "[PASSED]"; std::cout << ", dir=2 "; dv = normal; dv.Normalize(); dv *= thicknessInMM; if ((mitk::Equal(geometry->GetAxisVector(2), dv) == false)) { std::cout << "[FAILED]" << std::endl; return EXIT_FAILURE; } std::cout << "[PASSED]" << std::endl; std::cout << "Testing offset: "; if ((mitk::Equal(geometry->GetCornerPoint(0), cornerpoint0, (double)vnl_math::float_eps) == false)) { std::cout << "[FAILED]" << std::endl; return EXIT_FAILURE; } std::cout << "[PASSED]" << std::endl; return EXIT_SUCCESS; } int testGeometry(const mitk::Geometry3D *geometry, const mitk::ScalarType &width, const mitk::ScalarType &height, const mitk::ScalarType &numSlices, const mitk::ScalarType &widthInMM, const mitk::ScalarType &heightInMM, const mitk::ScalarType &thicknessInMM, const mitk::Point3D &cornerpoint0, const mitk::Vector3D &right, const mitk::Vector3D &bottom, const mitk::Vector3D &normal) { int result = EXIT_FAILURE; std::cout << "Comparing GetCornerPoint(0) of Geometry3D with provided cornerpoint0: "; if (mitk::Equal(geometry->GetCornerPoint(0), cornerpoint0) == false) { std::cout << "[FAILED]" << std::endl; return EXIT_FAILURE; } std::cout << "[PASSED]" << std::endl; std::cout << "Creating and initializing a SliceNavigationController with the Geometry3D: "; mitk::SliceNavigationController::Pointer sliceCtrl = mitk::SliceNavigationController::New(); sliceCtrl->SetInputWorldGeometry3D(geometry); std::cout << "[PASSED]" << std::endl; std::cout << "Testing SetViewDirection(mitk::SliceNavigationController::Axial): "; sliceCtrl->SetViewDirection(mitk::SliceNavigationController::Axial); std::cout << "[PASSED]" << std::endl; std::cout << "Testing Update(): "; sliceCtrl->Update(); std::cout << "[PASSED]" << std::endl; std::cout << "Testing result of CreatedWorldGeometry(): "; mitk::Point3D axialcornerpoint0; axialcornerpoint0 = cornerpoint0 + bottom + normal * (numSlices - 1 + 0.5); // really -1? result = compareGeometry(*sliceCtrl->GetCreatedWorldGeometry(), width, height, numSlices, widthInMM, heightInMM, thicknessInMM * numSlices, axialcornerpoint0, right, bottom * (-1.0), normal * (-1.0)); if (result != EXIT_SUCCESS) { std::cout << "[FAILED]" << std::endl; return result; } std::cout << "[PASSED]" << std::endl; std::cout << "Testing SetViewDirection(mitk::SliceNavigationController::Frontal): "; sliceCtrl->SetViewDirection(mitk::SliceNavigationController::Frontal); std::cout << "[PASSED]" << std::endl; std::cout << "Testing Update(): "; sliceCtrl->Update(); std::cout << "[PASSED]" << std::endl; std::cout << "Testing result of CreatedWorldGeometry(): "; mitk::Point3D frontalcornerpoint0; frontalcornerpoint0 = cornerpoint0 + geometry->GetAxisVector(1) * (+0.5 / geometry->GetExtent(1)); result = compareGeometry(*sliceCtrl->GetCreatedWorldGeometry(), width, numSlices, height, widthInMM, thicknessInMM * numSlices, heightInMM, frontalcornerpoint0, right, normal, bottom); if (result != EXIT_SUCCESS) { std::cout << "[FAILED]" << std::endl; return result; } std::cout << "[PASSED]" << std::endl; std::cout << "Testing SetViewDirection(mitk::SliceNavigationController::Sagittal): "; sliceCtrl->SetViewDirection(mitk::SliceNavigationController::Sagittal); std::cout << "[PASSED]" << std::endl; std::cout << "Testing Update(): " << std::endl; sliceCtrl->Update(); std::cout << "[PASSED]" << std::endl; std::cout << "Testing result of CreatedWorldGeometry(): "; mitk::Point3D sagittalcornerpoint0; sagittalcornerpoint0 = cornerpoint0 + geometry->GetAxisVector(0) * (+0.5 / geometry->GetExtent(0)); result = compareGeometry(*sliceCtrl->GetCreatedWorldGeometry(), height, numSlices, width, heightInMM, thicknessInMM * numSlices, widthInMM, sagittalcornerpoint0, bottom, normal, right); if (result != EXIT_SUCCESS) { std::cout << "[FAILED]" << std::endl; return result; } std::cout << "[PASSED]" << std::endl; return EXIT_SUCCESS; } int testReorientPlanes() { // Create PlaneGeometry mitk::PlaneGeometry::Pointer planegeometry = mitk::PlaneGeometry::New(); mitk::Point3D origin; mitk::Vector3D right, bottom, normal; mitk::ScalarType width, height; mitk::ScalarType widthInMM, heightInMM, thicknessInMM; width = 100; widthInMM = width; height = 200; heightInMM = height; thicknessInMM = 1.5; mitk::FillVector3D(origin, 4.5, 7.3, 11.2); mitk::FillVector3D(right, widthInMM, 0, 0); mitk::FillVector3D(bottom, 0, heightInMM, 0); mitk::FillVector3D(normal, 0, 0, thicknessInMM); mitk::Vector3D spacing; normal.Normalize(); normal *= thicknessInMM; mitk::FillVector3D(spacing, 1.0, 1.0, thicknessInMM); planegeometry->InitializeStandardPlane(right.GetVnlVector(), bottom.GetVnlVector(), &spacing); planegeometry->SetOrigin(origin); // Create SlicedGeometry3D out of planeGeometry mitk::SlicedGeometry3D::Pointer slicedgeometry1 = mitk::SlicedGeometry3D::New(); unsigned int numSlices = 20; - slicedgeometry1->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices, false); + slicedgeometry1->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices); // Create another slicedgeo which will be rotated mitk::SlicedGeometry3D::Pointer slicedgeometry2 = mitk::SlicedGeometry3D::New(); - slicedgeometry2->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices, false); + slicedgeometry2->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices); // Create geo3D as reference mitk::Geometry3D::Pointer geometry = mitk::Geometry3D::New(); geometry->SetBounds(slicedgeometry1->GetBounds()); geometry->SetIndexToWorldTransform(slicedgeometry1->GetIndexToWorldTransform()); // Initialize planes for (int i = 0; i < (int)numSlices; i++) { mitk::PlaneGeometry::Pointer geo2d = mitk::PlaneGeometry::New(); geo2d->Initialize(); geo2d->SetReferenceGeometry(geometry); slicedgeometry1->SetPlaneGeometry(geo2d, i); } for (int i = 0; i < (int)numSlices; i++) { mitk::PlaneGeometry::Pointer geo2d = mitk::PlaneGeometry::New(); geo2d->Initialize(); geo2d->SetReferenceGeometry(geometry); slicedgeometry2->SetPlaneGeometry(geo2d, i); } slicedgeometry1->SetReferenceGeometry(geometry); slicedgeometry2->SetReferenceGeometry(geometry); // Create SNC mitk::SliceNavigationController::Pointer sliceCtrl1 = mitk::SliceNavigationController::New(); sliceCtrl1->SetInputWorldGeometry3D(slicedgeometry1); sliceCtrl1->Update(); mitk::SliceNavigationController::Pointer sliceCtrl2 = mitk::SliceNavigationController::New(); sliceCtrl2->SetInputWorldGeometry3D(slicedgeometry2); sliceCtrl2->Update(); slicedgeometry1->SetSliceNavigationController(sliceCtrl1); slicedgeometry2->SetSliceNavigationController(sliceCtrl2); // Whats current geometry? MITK_INFO << "center: " << sliceCtrl1->GetCurrentPlaneGeometry()->GetCenter(); MITK_INFO << "normal: " << sliceCtrl1->GetCurrentPlaneGeometry()->GetNormal(); MITK_INFO << "origin: " << sliceCtrl1->GetCurrentPlaneGeometry()->GetOrigin(); MITK_INFO << "axis0 : " << sliceCtrl1->GetCurrentPlaneGeometry()->GetAxisVector(0); MITK_INFO << "aixs1 : " << sliceCtrl1->GetCurrentPlaneGeometry()->GetAxisVector(1); // // Now reorient slices (ONE POINT, ONE NORMAL) mitk::Point3D oldCenter, oldOrigin; mitk::Vector3D oldAxis0, oldAxis1; oldCenter = sliceCtrl1->GetCurrentPlaneGeometry()->GetCenter(); oldOrigin = sliceCtrl1->GetCurrentPlaneGeometry()->GetOrigin(); oldAxis0 = sliceCtrl1->GetCurrentPlaneGeometry()->GetAxisVector(0); oldAxis1 = sliceCtrl1->GetCurrentPlaneGeometry()->GetAxisVector(1); mitk::Point3D orientCenter; mitk::Vector3D orientNormal; orientCenter = oldCenter; mitk::FillVector3D(orientNormal, 0.3, 0.1, 0.8); orientNormal.Normalize(); sliceCtrl1->ReorientSlices(orientCenter, orientNormal); mitk::Point3D newCenter, newOrigin; mitk::Vector3D newNormal; newCenter = sliceCtrl1->GetCurrentPlaneGeometry()->GetCenter(); newOrigin = sliceCtrl1->GetCurrentPlaneGeometry()->GetOrigin(); newNormal = sliceCtrl1->GetCurrentPlaneGeometry()->GetNormal(); newNormal.Normalize(); itk::Index<3> orientCenterIdx; itk::Index<3> newCenterIdx; sliceCtrl1->GetCurrentGeometry3D()->WorldToIndex(orientCenter, orientCenterIdx); sliceCtrl1->GetCurrentGeometry3D()->WorldToIndex(newCenter, newCenterIdx); if ((newCenterIdx != orientCenterIdx) || (!mitk::Equal(orientNormal, newNormal))) { MITK_INFO << "Reorient Planes (1 point, 1 vector) not working as it should"; MITK_INFO << "orientCenterIdx: " << orientCenterIdx; MITK_INFO << "newCenterIdx: " << newCenterIdx; MITK_INFO << "orientNormal: " << orientNormal; MITK_INFO << "newNormal: " << newNormal; return EXIT_FAILURE; } // // Now reorient slices (center, vec0, vec1 ) mitk::Vector3D orientAxis0, orientAxis1, newAxis0, newAxis1; mitk::FillVector3D(orientAxis0, 1.0, 0.0, 0.0); mitk::FillVector3D(orientAxis1, 0.0, 1.0, 0.0); orientAxis0.Normalize(); orientAxis1.Normalize(); sliceCtrl1->ReorientSlices(orientCenter, orientAxis0, orientAxis1); newAxis0 = sliceCtrl1->GetCurrentPlaneGeometry()->GetAxisVector(0); newAxis1 = sliceCtrl1->GetCurrentPlaneGeometry()->GetAxisVector(1); newCenter = sliceCtrl1->GetCurrentPlaneGeometry()->GetCenter(); newAxis0.Normalize(); newAxis1.Normalize(); sliceCtrl1->GetCurrentGeometry3D()->WorldToIndex(orientCenter, orientCenterIdx); sliceCtrl1->GetCurrentGeometry3D()->WorldToIndex(newCenter, newCenterIdx); if ((newCenterIdx != orientCenterIdx) || (!mitk::Equal(orientAxis0, newAxis0, 1E-12)) || (!mitk::Equal(orientAxis1, newAxis1, 1E-12))) { MITK_INFO << "Reorient Planes (point, vec, vec) not working as it should"; MITK_INFO << "orientCenterIdx: " << orientCenterIdx; MITK_INFO << "newCenterIdx: " << newCenterIdx; MITK_INFO << "orientAxis0: " << orientAxis0; MITK_INFO << "newAxis0: " << newAxis0; MITK_INFO << "orientAxis1: " << orientAxis1; MITK_INFO << "newAxis1: " << newAxis1; return EXIT_FAILURE; } return EXIT_SUCCESS; } int testRestorePlanePostionOperation() { // Create PlaneGeometry mitk::PlaneGeometry::Pointer planegeometry = mitk::PlaneGeometry::New(); mitk::Point3D origin; mitk::Vector3D right, bottom, normal; mitk::ScalarType width, height; mitk::ScalarType widthInMM, heightInMM, thicknessInMM; width = 100; widthInMM = width; height = 200; heightInMM = height; thicknessInMM = 1.5; mitk::FillVector3D(origin, 4.5, 7.3, 11.2); mitk::FillVector3D(right, widthInMM, 0, 0); mitk::FillVector3D(bottom, 0, heightInMM, 0); mitk::FillVector3D(normal, 0, 0, thicknessInMM); mitk::Vector3D spacing; normal.Normalize(); normal *= thicknessInMM; mitk::FillVector3D(spacing, 1.0, 1.0, thicknessInMM); planegeometry->InitializeStandardPlane(right.GetVnlVector(), bottom.GetVnlVector(), &spacing); planegeometry->SetOrigin(origin); // Create SlicedGeometry3D out of planeGeometry mitk::SlicedGeometry3D::Pointer slicedgeometry1 = mitk::SlicedGeometry3D::New(); unsigned int numSlices = 300; - slicedgeometry1->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices, false); + slicedgeometry1->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices); // Create another slicedgeo which will be rotated mitk::SlicedGeometry3D::Pointer slicedgeometry2 = mitk::SlicedGeometry3D::New(); - slicedgeometry2->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices, false); + slicedgeometry2->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices); // Create geo3D as reference mitk::Geometry3D::Pointer geometry = mitk::Geometry3D::New(); geometry->SetBounds(slicedgeometry1->GetBounds()); geometry->SetIndexToWorldTransform(slicedgeometry1->GetIndexToWorldTransform()); // Initialize planes for (int i = 0; i < (int)numSlices; i++) { mitk::PlaneGeometry::Pointer geo2d = mitk::PlaneGeometry::New(); geo2d->Initialize(); geo2d->SetReferenceGeometry(geometry); slicedgeometry1->SetPlaneGeometry(geo2d, i); } for (int i = 0; i < (int)numSlices; i++) { mitk::PlaneGeometry::Pointer geo2d = mitk::PlaneGeometry::New(); geo2d->Initialize(); geo2d->SetReferenceGeometry(geometry); slicedgeometry2->SetPlaneGeometry(geo2d, i); } slicedgeometry1->SetReferenceGeometry(geometry); slicedgeometry2->SetReferenceGeometry(geometry); // Create SNC mitk::SliceNavigationController::Pointer sliceCtrl1 = mitk::SliceNavigationController::New(); sliceCtrl1->SetInputWorldGeometry3D(slicedgeometry1); sliceCtrl1->Update(); mitk::SliceNavigationController::Pointer sliceCtrl2 = mitk::SliceNavigationController::New(); sliceCtrl2->SetInputWorldGeometry3D(slicedgeometry2); sliceCtrl2->Update(); slicedgeometry1->SetSliceNavigationController(sliceCtrl1); slicedgeometry2->SetSliceNavigationController(sliceCtrl2); // Rotate slicedgeo2 double angle = 63.84; mitk::Vector3D rotationVector; mitk::FillVector3D(rotationVector, 0.5, 0.95, 0.23); mitk::Point3D center = slicedgeometry2->GetCenter(); mitk::RotationOperation *op = new mitk::RotationOperation(mitk::OpROTATE, center, rotationVector, angle); slicedgeometry2->ExecuteOperation(op); sliceCtrl2->Update(); us::ServiceReference serviceRef = us::GetModuleContext()->GetServiceReference(); mitk::PlanePositionManagerService *service = us::GetModuleContext()->GetService(serviceRef); service->AddNewPlanePosition(slicedgeometry2->GetPlaneGeometry(0), 178); sliceCtrl1->ExecuteOperation(service->GetPlanePosition(0)); sliceCtrl1->Update(); mitk::PlaneGeometry *planeRotated = slicedgeometry2->GetPlaneGeometry(178); mitk::PlaneGeometry *planeRestored = dynamic_cast(sliceCtrl1->GetCurrentGeometry3D())->GetPlaneGeometry(178); try { MITK_TEST_CONDITION_REQUIRED(mitk::MatrixEqualElementWise(planeRotated->GetIndexToWorldTransform()->GetMatrix(), planeRestored->GetIndexToWorldTransform()->GetMatrix()), "Testing for IndexToWorld"); MITK_TEST_CONDITION_REQUIRED(mitk::Equal(planeRotated->GetOrigin(), planeRestored->GetOrigin(), 2 * mitk::eps), "Testing for origin"); MITK_TEST_CONDITION_REQUIRED(mitk::Equal(planeRotated->GetSpacing(), planeRestored->GetSpacing()), "Testing for spacing"); MITK_TEST_CONDITION_REQUIRED( mitk::Equal( slicedgeometry2->GetDirectionVector(), dynamic_cast(sliceCtrl1->GetCurrentGeometry3D())->GetDirectionVector()), "Testing for directionvector"); MITK_TEST_CONDITION_REQUIRED( mitk::Equal(slicedgeometry2->GetSlices(), dynamic_cast(sliceCtrl1->GetCurrentGeometry3D())->GetSlices()), "Testing for numslices"); MITK_TEST_CONDITION_REQUIRED( mitk::MatrixEqualElementWise(slicedgeometry2->GetIndexToWorldTransform()->GetMatrix(), dynamic_cast(sliceCtrl1->GetCurrentGeometry3D()) ->GetIndexToWorldTransform() ->GetMatrix()), "Testing for IndexToWorld"); } catch (...) { return EXIT_FAILURE; } return EXIT_SUCCESS; } int mitkSliceNavigationControllerTest(int /*argc*/, char * /*argv*/ []) { int result = EXIT_FAILURE; std::cout << "Creating and initializing a PlaneGeometry: "; mitk::PlaneGeometry::Pointer planegeometry = mitk::PlaneGeometry::New(); mitk::Point3D origin; mitk::Vector3D right, bottom, normal; mitk::ScalarType width, height; mitk::ScalarType widthInMM, heightInMM, thicknessInMM; width = 100; widthInMM = width; height = 200; heightInMM = height; thicknessInMM = 1.5; // mitk::FillVector3D(origin, 0, 0, thicknessInMM*0.5); mitk::FillVector3D(origin, 4.5, 7.3, 11.2); mitk::FillVector3D(right, widthInMM, 0, 0); mitk::FillVector3D(bottom, 0, heightInMM, 0); mitk::FillVector3D(normal, 0, 0, thicknessInMM); mitk::Vector3D spacing; normal.Normalize(); normal *= thicknessInMM; mitk::FillVector3D(spacing, 1.0, 1.0, thicknessInMM); planegeometry->InitializeStandardPlane(right.GetVnlVector(), bottom.GetVnlVector(), &spacing); planegeometry->SetOrigin(origin); std::cout << "[PASSED]" << std::endl; std::cout << "Creating and initializing a SlicedGeometry3D with the PlaneGeometry: "; mitk::SlicedGeometry3D::Pointer slicedgeometry = mitk::SlicedGeometry3D::New(); unsigned int numSlices = 5; - slicedgeometry->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices, false); + slicedgeometry->InitializeEvenlySpaced(planegeometry, thicknessInMM, numSlices); std::cout << "[PASSED]" << std::endl; std::cout << "Creating a Geometry3D with the same extent as the SlicedGeometry3D: "; mitk::Geometry3D::Pointer geometry = mitk::Geometry3D::New(); geometry->SetBounds(slicedgeometry->GetBounds()); geometry->SetIndexToWorldTransform(slicedgeometry->GetIndexToWorldTransform()); std::cout << "[PASSED]" << std::endl; mitk::Point3D cornerpoint0; cornerpoint0 = geometry->GetCornerPoint(0); result = testGeometry( geometry, width, height, numSlices, widthInMM, heightInMM, thicknessInMM, cornerpoint0, right, bottom, normal); if (result != EXIT_SUCCESS) return result; mitk::AffineTransform3D::Pointer transform = mitk::AffineTransform3D::New(); transform->SetMatrix(geometry->GetIndexToWorldTransform()->GetMatrix()); mitk::BoundingBox::Pointer boundingbox = geometry->CalculateBoundingBoxRelativeToTransform(transform); geometry->SetBounds(boundingbox->GetBounds()); cornerpoint0 = geometry->GetCornerPoint(0); result = testGeometry( geometry, width, height, numSlices, widthInMM, heightInMM, thicknessInMM, cornerpoint0, right, bottom, normal); if (result != EXIT_SUCCESS) return result; std::cout << "Changing the IndexToWorldTransform of the geometry to a rotated version by SetIndexToWorldTransform() " "(keep cornerpoint0): "; transform = mitk::AffineTransform3D::New(); mitk::AffineTransform3D::MatrixType::InternalMatrixType vnlmatrix; vnlmatrix = planegeometry->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix(); mitk::VnlVector axis(3); mitk::FillVector3D(axis, 1.0, 1.0, 1.0); axis.normalize(); vnl_quaternion rotation(axis, 0.223); vnlmatrix = rotation.rotation_matrix_transpose() * vnlmatrix; mitk::Matrix3D matrix; matrix = vnlmatrix; transform->SetMatrix(matrix); transform->SetOffset(cornerpoint0.GetVectorFromOrigin()); right.SetVnlVector(rotation.rotation_matrix_transpose() * right.GetVnlVector()); bottom.SetVnlVector(rotation.rotation_matrix_transpose() * bottom.GetVnlVector()); normal.SetVnlVector(rotation.rotation_matrix_transpose() * normal.GetVnlVector()); geometry->SetIndexToWorldTransform(transform); std::cout << "[PASSED]" << std::endl; cornerpoint0 = geometry->GetCornerPoint(0); result = testGeometry( geometry, width, height, numSlices, widthInMM, heightInMM, thicknessInMM, cornerpoint0, right, bottom, normal); if (result != EXIT_SUCCESS) return result; // Testing Execute RestorePlanePositionOperation result = testRestorePlanePostionOperation(); if (result != EXIT_SUCCESS) return result; // Re-Orient planes not working as it should on windows // However, this might be adjusted during a geometry redesign. /* //Testing ReorientPlanes result = testReorientPlanes(); if(result!=EXIT_SUCCESS) return result; */ std::cout << "[TEST DONE]" << std::endl; return EXIT_SUCCESS; } \ No newline at end of file diff --git a/Modules/Core/test/mitkSlicedGeometry3DTest.cpp b/Modules/Core/test/mitkSlicedGeometry3DTest.cpp index b0fd527a4b..0a06d1db10 100644 --- a/Modules/Core/test/mitkSlicedGeometry3DTest.cpp +++ b/Modules/Core/test/mitkSlicedGeometry3DTest.cpp @@ -1,288 +1,184 @@ /*=================================================================== 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 "mitkImage.h" #include "mitkPlaneGeometry.h" #include "mitkSlicedGeometry3D.h" #include "mitkTestingMacros.h" #include #include #include +#include -static const mitk::ScalarType slicedGeometryEps = 1E-9; // epsilon for this testfile. Set to float precision. +static const mitk::ScalarType slicedGeometryEps = 1e-9; // Set epsilon to float precision for this test -void mitkSlicedGeometry3D_ChangeImageGeometryConsideringOriginOffset_Test() +static mitk::PlaneGeometry::Pointer createPlaneGeometry() +{ + auto planeGeometry = mitk::PlaneGeometry::New(); + planeGeometry->Initialize(); + return planeGeometry; +} + +static mitk::SlicedGeometry3D::Pointer createSlicedGeometry(const mitk::Point3D &origin, const mitk::Vector3D &spacing, int numberOfSlices) { - // Tests for Offset - MITK_TEST_OUTPUT(<< "====== NOW RUNNING: Tests for pixel-center-based offset concerns ========"); - - // create a SlicedGeometry3D - mitk::SlicedGeometry3D::Pointer slicedGeo3D = mitk::SlicedGeometry3D::New(); - int num_slices = 5; - slicedGeo3D->InitializeSlicedGeometry(num_slices); // 5 slices - mitk::Point3D newOrigin; - newOrigin[0] = 91.3; - newOrigin[1] = -13.3; - newOrigin[2] = 0; - slicedGeo3D->SetOrigin(newOrigin); - mitk::Vector3D newSpacing; - newSpacing[0] = 1.0f; - newSpacing[1] = 0.9f; - newSpacing[2] = 0.3f; - slicedGeo3D->SetSpacing(newSpacing); - // create subslices as well - for (int i = 0; i < num_slices; i++) + auto slicedGeometry = mitk::SlicedGeometry3D::New(); + slicedGeometry->InitializeSlicedGeometry(static_cast(numberOfSlices)); + slicedGeometry->SetOrigin(origin); + slicedGeometry->SetSpacing(spacing); + + for (unsigned int i = 0; i < numberOfSlices; ++i) { - mitk::PlaneGeometry::Pointer geo2d = mitk::PlaneGeometry::New(); - geo2d->Initialize(); - slicedGeo3D->SetPlaneGeometry(geo2d, i); + auto planeGeometry = createPlaneGeometry(); + slicedGeometry->SetPlaneGeometry(planeGeometry, i); } - // now run tests - - MITK_TEST_OUTPUT(<< "Testing whether slicedGeo3D->GetImageGeometry() is false by default"); - MITK_TEST_CONDITION_REQUIRED(slicedGeo3D->GetImageGeometry() == false, ""); - MITK_TEST_OUTPUT( - << "Testing whether first and last geometry in the SlicedGeometry3D have GetImageGeometry()==false by default"); - mitk::BaseGeometry *subSliceGeo2D_first = slicedGeo3D->GetPlaneGeometry(0); - mitk::BaseGeometry *subSliceGeo2D_last = slicedGeo3D->GetPlaneGeometry(num_slices - 1); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_first->GetImageGeometry() == false, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_last->GetImageGeometry() == false, ""); - - // Save some Origins and cornerpoints - mitk::Point3D OriginSlicedGeo(slicedGeo3D->GetOrigin()); - mitk::Point3D OriginFirstGeo(subSliceGeo2D_first->GetOrigin()); - mitk::Point3D OriginLastGeo(subSliceGeo2D_last->GetOrigin()); - mitk::Point3D CornerPoint0SlicedGeo(slicedGeo3D->GetCornerPoint(0)); - mitk::Point3D CornerPoint1FirstGeo(subSliceGeo2D_first->GetCornerPoint(1)); - mitk::Point3D CornerPoint2LastGeo(subSliceGeo2D_last->GetCornerPoint(2)); - - MITK_TEST_OUTPUT(<< "Calling slicedGeo3D->ChangeImageGeometryConsideringOriginOffset(true)"); - // std::cout << "vorher Origin: " << subSliceGeo2D_first->GetOrigin() << std::endl; - // std::cout << "vorher Corner: " << subSliceGeo2D_first->GetCornerPoint(0) << std::endl; - slicedGeo3D->ChangeImageGeometryConsideringOriginOffset(true); - // std::cout << "nachher Origin: " << subSliceGeo2D_first->GetOrigin() << std::endl; - // std::cout << "nachher Corner: " << subSliceGeo2D_first->GetCornerPoint(0) << std::endl; - MITK_TEST_OUTPUT(<< "Testing whether slicedGeo3D->GetImageGeometry() is now true"); - MITK_TEST_CONDITION_REQUIRED(slicedGeo3D->GetImageGeometry() == true, ""); - MITK_TEST_OUTPUT( - << "Testing whether first and last geometry in the SlicedGeometry3D have GetImageGeometry()==true now"); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_first->GetImageGeometry() == true, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_last->GetImageGeometry() == true, ""); - - MITK_TEST_OUTPUT(<< "Testing wether offset has been added to origins"); - // Manually adding Offset. - OriginSlicedGeo[0] += (slicedGeo3D->GetSpacing()[0]) / 2; - OriginSlicedGeo[1] += (slicedGeo3D->GetSpacing()[1]) / 2; - OriginSlicedGeo[2] += (slicedGeo3D->GetSpacing()[2]) / 2; - OriginFirstGeo[0] += (subSliceGeo2D_first->GetSpacing()[0]) / 2; - OriginFirstGeo[1] += (subSliceGeo2D_first->GetSpacing()[1]) / 2; - OriginFirstGeo[2] += (subSliceGeo2D_first->GetSpacing()[2]) / 2; - OriginLastGeo[0] += (subSliceGeo2D_last->GetSpacing()[0]) / 2; - OriginLastGeo[1] += (subSliceGeo2D_last->GetSpacing()[1]) / 2; - OriginLastGeo[2] += (subSliceGeo2D_last->GetSpacing()[2]) / 2; - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_first->GetCornerPoint(1) == CornerPoint1FirstGeo, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_last->GetCornerPoint(2) == CornerPoint2LastGeo, ""); - MITK_TEST_CONDITION_REQUIRED(slicedGeo3D->GetCornerPoint(0) == CornerPoint0SlicedGeo, ""); - MITK_TEST_CONDITION_REQUIRED(slicedGeo3D->GetOrigin() == OriginSlicedGeo, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_first->GetOrigin() == OriginFirstGeo, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_last->GetOrigin() == OriginLastGeo, ""); - - MITK_TEST_OUTPUT(<< "Calling slicedGeo3D->ChangeImageGeometryConsideringOriginOffset(false)"); - slicedGeo3D->ChangeImageGeometryConsideringOriginOffset(false); - MITK_TEST_OUTPUT(<< "Testing whether slicedGeo3D->GetImageGeometry() is now false"); - MITK_TEST_CONDITION_REQUIRED(slicedGeo3D->GetImageGeometry() == false, ""); - MITK_TEST_OUTPUT( - << "Testing whether first and last geometry in the SlicedGeometry3D have GetImageGeometry()==false now"); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_first->GetImageGeometry() == false, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_last->GetImageGeometry() == false, ""); - - MITK_TEST_OUTPUT(<< "Testing wether offset has been added to origins of geometry"); - - // Manually substracting Offset. - OriginSlicedGeo[0] -= (slicedGeo3D->GetSpacing()[0]) / 2; - OriginSlicedGeo[1] -= (slicedGeo3D->GetSpacing()[1]) / 2; - OriginSlicedGeo[2] -= (slicedGeo3D->GetSpacing()[2]) / 2; - OriginFirstGeo[0] -= (subSliceGeo2D_first->GetSpacing()[0]) / 2; - OriginFirstGeo[1] -= (subSliceGeo2D_first->GetSpacing()[1]) / 2; - OriginFirstGeo[2] -= (subSliceGeo2D_first->GetSpacing()[2]) / 2; - OriginLastGeo[0] -= (subSliceGeo2D_last->GetSpacing()[0]) / 2; - OriginLastGeo[1] -= (subSliceGeo2D_last->GetSpacing()[1]) / 2; - OriginLastGeo[2] -= (subSliceGeo2D_last->GetSpacing()[2]) / 2; - - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_first->GetCornerPoint(1) == CornerPoint1FirstGeo, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_last->GetCornerPoint(2) == CornerPoint2LastGeo, ""); - MITK_TEST_CONDITION_REQUIRED(slicedGeo3D->GetCornerPoint(0) == CornerPoint0SlicedGeo, ""); - MITK_TEST_CONDITION_REQUIRED(slicedGeo3D->GetOrigin() == OriginSlicedGeo, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_first->GetOrigin() == OriginFirstGeo, ""); - MITK_TEST_CONDITION_REQUIRED(subSliceGeo2D_last->GetOrigin() == OriginLastGeo, ""); - - MITK_TEST_OUTPUT(<< "ALL SUCCESSFULLY!"); + return slicedGeometry; } -int mitkSlicedGeometry3DTest(int /*argc*/, char * /*argv*/ []) +static mitk::SlicedGeometry3D::Pointer createEvenlySpacedSlicedGeometry(mitk::PlaneGeometry::Pointer planeGeometry, mitk::ScalarType spacing, int numberOfSlices) { - mitk::PlaneGeometry::Pointer planegeometry1 = mitk::PlaneGeometry::New(); - - mitk::Point3D origin; - mitk::Vector3D right, bottom, normal; - mitk::ScalarType width, height; - mitk::ScalarType widthInMM, heightInMM, thicknessInMM; - - width = 100; - widthInMM = width; - height = 200; - heightInMM = height; - thicknessInMM = 3.5; - mitk::FillVector3D(origin, 4.5, 7.3, 11.2); - mitk::FillVector3D(right, widthInMM, 0, 0); - mitk::FillVector3D(bottom, 0, heightInMM, 0); - mitk::FillVector3D(normal, 0, 0, thicknessInMM); - - std::cout << "Initializing planegeometry1 by InitializeStandardPlane(rightVector, downVector, spacing = NULL): " - << std::endl; - planegeometry1->InitializeStandardPlane(right.GetVnlVector(), bottom.GetVnlVector()); - - std::cout << "Setting planegeometry2 to a cloned version of planegeometry1: " << std::endl; - mitk::PlaneGeometry::Pointer planegeometry2; - planegeometry2 = dynamic_cast(planegeometry1->Clone().GetPointer()); - ; - - std::cout << "Changing the IndexToWorldTransform of planegeometry2 to a rotated version by " - "SetIndexToWorldTransform() (keep origin): " - << std::endl; - mitk::AffineTransform3D::Pointer transform = mitk::AffineTransform3D::New(); - mitk::AffineTransform3D::MatrixType::InternalMatrixType vnlmatrix; - vnlmatrix = planegeometry2->GetIndexToWorldTransform()->GetMatrix().GetVnlMatrix(); - mitk::VnlVector axis(3); - mitk::FillVector3D(axis, 1.0, 1.0, 1.0); - axis.normalize(); - vnl_quaternion rotation(axis, 0.123); - vnlmatrix = rotation.rotation_matrix_transpose() * vnlmatrix; - mitk::Matrix3D matrix; - matrix = vnlmatrix; - transform->SetMatrix(matrix); - transform->SetOffset(planegeometry2->GetIndexToWorldTransform()->GetOffset()); - - right.SetVnlVector(rotation.rotation_matrix_transpose() * right.GetVnlVector()); - bottom.SetVnlVector(rotation.rotation_matrix_transpose() * bottom.GetVnlVector()); - normal.SetVnlVector(rotation.rotation_matrix_transpose() * normal.GetVnlVector()); - planegeometry2->SetIndexToWorldTransform(transform); - - std::cout << "Setting planegeometry3 to the backside of planegeometry2: " << std::endl; - mitk::PlaneGeometry::Pointer planegeometry3 = mitk::PlaneGeometry::New(); - planegeometry3->InitializeStandardPlane(planegeometry2, mitk::PlaneGeometry::Axial, 0, false); - - std::cout - << "Testing SlicedGeometry3D::InitializeEvenlySpaced(planegeometry3, zSpacing = 1, slices = 5, flipped = false): " - << std::endl; - mitk::SlicedGeometry3D::Pointer slicedWorldGeometry = mitk::SlicedGeometry3D::New(); - unsigned int numSlices = 5; - slicedWorldGeometry->InitializeEvenlySpaced(planegeometry3, 1, numSlices, false); - - std::cout << "Testing availability and type (PlaneGeometry) of first geometry in the SlicedGeometry3D: "; - mitk::PlaneGeometry *accessedplanegeometry3 = - dynamic_cast(slicedWorldGeometry->GetPlaneGeometry(0)); - if (accessedplanegeometry3 == nullptr) - { - std::cout << "[FAILED]" << std::endl; - return EXIT_FAILURE; - } - std::cout << "[PASSED]" << std::endl; - - std::cout << "Testing whether the first geometry in the SlicedGeometry3D is identical to planegeometry3 by axis " - "comparison and origin: " - << std::endl; - if ((mitk::Equal(accessedplanegeometry3->GetAxisVector(0), planegeometry3->GetAxisVector(0), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3->GetAxisVector(1), planegeometry3->GetAxisVector(1), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3->GetAxisVector(2), planegeometry3->GetAxisVector(2), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3->GetOrigin(), planegeometry3->GetOrigin(), slicedGeometryEps) == false)) - { - std::cout << "[FAILED]" << std::endl; - return EXIT_FAILURE; - } - std::cout << "[PASSED]" << std::endl; - - std::cout << "Testing availability and type (PlaneGeometry) of the last geometry in the SlicedGeometry3D: "; - mitk::PlaneGeometry *accessedplanegeometry3last = - dynamic_cast(slicedWorldGeometry->GetPlaneGeometry(numSlices - 1)); - mitk::Point3D origin3last; - origin3last = planegeometry3->GetOrigin() + slicedWorldGeometry->GetDirectionVector() * (numSlices - 1); - if (accessedplanegeometry3last == nullptr) - { - std::cout << "[FAILED]" << std::endl; - return EXIT_FAILURE; - } - std::cout << "[PASSED]" << std::endl; - - std::cout - << "Testing whether the last geometry in the SlicedGeometry3D is identical to planegeometry3 by axis comparison: " - << std::endl; - if ((mitk::Equal(accessedplanegeometry3last->GetAxisVector(0), planegeometry3->GetAxisVector(0), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3last->GetAxisVector(1), planegeometry3->GetAxisVector(1), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3last->GetAxisVector(2), planegeometry3->GetAxisVector(2), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3last->GetOrigin(), origin3last, slicedGeometryEps) == false) || - (mitk::Equal(accessedplanegeometry3last->GetIndexToWorldTransform()->GetOffset(), - origin3last.GetVectorFromOrigin(), - slicedGeometryEps) == false)) - { - std::cout << "[FAILED]" << std::endl; - return EXIT_FAILURE; - } - std::cout << "[PASSED]" << std::endl; + auto slicedGeometry = mitk::SlicedGeometry3D::New(); + slicedGeometry->InitializeEvenlySpaced(planeGeometry, spacing, numberOfSlices); + return slicedGeometry; +} - std::cout << "Again for first slice - Testing availability and type (PlaneGeometry) of first geometry in the " - "SlicedGeometry3D: "; - accessedplanegeometry3 = dynamic_cast(slicedWorldGeometry->GetPlaneGeometry(0)); - if (accessedplanegeometry3 == nullptr) - { - std::cout << "[FAILED]" << std::endl; - return EXIT_FAILURE; - } - std::cout << "[PASSED]" << std::endl; - - std::cout << "Again for first slice - Testing whether the first geometry in the SlicedGeometry3D is identical to " - "planegeometry3 by axis comparison and origin: " - << std::endl; - if ((mitk::Equal(accessedplanegeometry3->GetAxisVector(0), planegeometry3->GetAxisVector(0), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3->GetAxisVector(1), planegeometry3->GetAxisVector(1), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3->GetAxisVector(2), planegeometry3->GetAxisVector(2), slicedGeometryEps) == - false) || - (mitk::Equal(accessedplanegeometry3->GetOrigin(), planegeometry3->GetOrigin(), slicedGeometryEps) == false) || - (mitk::Equal(accessedplanegeometry3->GetIndexToWorldTransform()->GetOffset(), - planegeometry3->GetOrigin().GetVectorFromOrigin(), - slicedGeometryEps) == false)) - { - std::cout << "[FAILED]" << std::endl; - return EXIT_FAILURE; - } - std::cout << "[PASSED]" << std::endl; +template static T createArray(mitk::ScalarType x, mitk::ScalarType y, mitk::ScalarType z) +{ + T array; + mitk::FillVector3D(array, x, y, z); + return array; +} + +static mitk::Point3D createPoint(mitk::ScalarType x, mitk::ScalarType y, mitk::ScalarType z) +{ + return createArray(x, y, z); +} + +static mitk::Vector3D createVector(mitk::ScalarType x, mitk::ScalarType y, mitk::ScalarType z) +{ + return createArray(x, y, z); +} + +static mitk::VnlVector createVnlVector(mitk::ScalarType x, mitk::ScalarType y, mitk::ScalarType z) +{ + mitk::VnlVector vector(3); + mitk::FillVector3D(vector, x, y, z); + return vector; +} + +void mitkSlicedGeometry3D_ChangeImageGeometryConsideringOriginOffset_Test() +{ + MITK_TEST_OUTPUT(<< "====== mitkSlicedGeometry3D_ChangeImageGeometryConsideringOriginOffset_Test() ======"); + + auto origin = createPoint(91.3, -13.3, 0); + auto spacing = createVector(1.0, 0.9, 0.3); + auto numberOfSlices = 5; + auto slicedGeometry = createSlicedGeometry(origin, spacing, numberOfSlices); + + MITK_TEST_OUTPUT(<< "SlicedGeometry3D isn't an image geometry by default"); + MITK_TEST_CONDITION_REQUIRED(slicedGeometry->GetImageGeometry() == false, ""); + + MITK_TEST_OUTPUT(<< "First and last PlaneGeometry in SlicedGeometry3D are not image geometries"); + auto firstPlaneGeometry = slicedGeometry->GetPlaneGeometry(0); + auto lastPlaneGeometry = slicedGeometry->GetPlaneGeometry(numberOfSlices - 1); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry->GetImageGeometry() == false, ""); + MITK_TEST_CONDITION_REQUIRED(lastPlaneGeometry->GetImageGeometry() == false, ""); + + auto originOfSlicedGeometry = slicedGeometry->GetOrigin(); + auto originOfFirstPlaneGeometry = firstPlaneGeometry->GetOrigin(); + auto originOfLastPlaneGeometry = lastPlaneGeometry->GetOrigin(); + auto firstCornerPointOfSlicedGeometry = slicedGeometry->GetCornerPoint(0); + auto secondCornerPointOfFirstPlaneGeometry = firstPlaneGeometry->GetCornerPoint(1); + auto thirdCornerPointOfLastPlaneGeometry = lastPlaneGeometry->GetCornerPoint(2); + + MITK_TEST_OUTPUT(<< "Calling SlicedGeometry3D::ChangeImageGeometryConsideringOriginOffset(true)"); + slicedGeometry->ChangeImageGeometryConsideringOriginOffset(true); + MITK_TEST_OUTPUT(<< "SlicedGeometry3D is an image geometry"); + MITK_TEST_CONDITION_REQUIRED(slicedGeometry->GetImageGeometry() == true, ""); + MITK_TEST_OUTPUT(<< "First and last PlaneGeometry in SlicedGeometry3D are image geometries"); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry->GetImageGeometry() == true, ""); + MITK_TEST_CONDITION_REQUIRED(lastPlaneGeometry->GetImageGeometry() == true, ""); + + MITK_TEST_OUTPUT(<< "Corner points of geometries didn't change"); + MITK_TEST_CONDITION_REQUIRED(slicedGeometry->GetCornerPoint(0) == firstCornerPointOfSlicedGeometry, ""); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry->GetCornerPoint(1) == secondCornerPointOfFirstPlaneGeometry, ""); + MITK_TEST_CONDITION_REQUIRED(lastPlaneGeometry->GetCornerPoint(2) == thirdCornerPointOfLastPlaneGeometry, ""); + + MITK_TEST_OUTPUT(<< "Offsets were added to geometry origins"); + MITK_TEST_CONDITION_REQUIRED(slicedGeometry->GetOrigin() == originOfSlicedGeometry + slicedGeometry->GetSpacing() * 0.5, ""); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry->GetOrigin() == originOfFirstPlaneGeometry + firstPlaneGeometry->GetSpacing() * 0.5, ""); + MITK_TEST_CONDITION_REQUIRED(lastPlaneGeometry->GetOrigin() == originOfLastPlaneGeometry + lastPlaneGeometry->GetSpacing() * 0.5, ""); + + MITK_TEST_OUTPUT(<< "Calling SlicedGeometry3D::ChangeImageGeometryConsideringOriginOffset(false)"); + slicedGeometry->ChangeImageGeometryConsideringOriginOffset(false); + MITK_TEST_OUTPUT(<< "SlicedGeometry3D isn't an image geometry anymore"); + MITK_TEST_CONDITION_REQUIRED(slicedGeometry->GetImageGeometry() == false, ""); + MITK_TEST_OUTPUT(<< "First and last PlaneGeometry in SlicedGeometry3D are no image geometries anymore"); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry->GetImageGeometry() == false, ""); + MITK_TEST_CONDITION_REQUIRED(lastPlaneGeometry->GetImageGeometry() == false, ""); + + MITK_TEST_OUTPUT(<< "Corner points of geometries didn't change"); + MITK_TEST_CONDITION_REQUIRED(slicedGeometry->GetCornerPoint(0) == firstCornerPointOfSlicedGeometry, ""); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry->GetCornerPoint(1) == secondCornerPointOfFirstPlaneGeometry, ""); + MITK_TEST_CONDITION_REQUIRED(lastPlaneGeometry->GetCornerPoint(2) == thirdCornerPointOfLastPlaneGeometry, ""); + + MITK_TEST_OUTPUT(<< "Offsets were subtracted from geometry origins"); + MITK_TEST_CONDITION_REQUIRED(slicedGeometry->GetOrigin() == originOfSlicedGeometry, ""); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry->GetOrigin() == originOfFirstPlaneGeometry, ""); + MITK_TEST_CONDITION_REQUIRED(lastPlaneGeometry->GetOrigin() == originOfLastPlaneGeometry, ""); +} + +int mitkSlicedGeometry3DTest(int, char *[]) +{ + mitk::ScalarType width = 100.0; + mitk::ScalarType widthInMM = width; + + mitk::ScalarType height = 200.0; + mitk::ScalarType heightInMM = height; + + mitk::ScalarType thicknessInMM = 3.0; + + auto right = createVector(widthInMM, 0.0, 0.0); + auto bottom = createVector(0.0, heightInMM, 0.0); + auto normal = createVector(0.0, 0.0, thicknessInMM); + auto spacing = createVector(1.0, 1.0, thicknessInMM); + + auto planeGeometry = mitk::PlaneGeometry::New(); + planeGeometry->InitializeStandardPlane(right, bottom, &spacing); + + auto numberOfSlices = 5; + auto slicedGeometry = createEvenlySpacedSlicedGeometry(planeGeometry, thicknessInMM, numberOfSlices); + auto firstPlaneGeometry = slicedGeometry->GetPlaneGeometry(0); + + MITK_TEST_OUTPUT(<< "Check if first PlaneGeometry of evenly spaced SlicedGeometry is valid"); + MITK_TEST_CONDITION_REQUIRED(firstPlaneGeometry != nullptr, ""); + MITK_TEST_CONDITION_REQUIRED(mitk::Equal(firstPlaneGeometry->GetAxisVector(0), planeGeometry->GetAxisVector(0), slicedGeometryEps), ""); + MITK_TEST_CONDITION_REQUIRED(mitk::Equal(firstPlaneGeometry->GetAxisVector(1), planeGeometry->GetAxisVector(1), slicedGeometryEps), ""); + MITK_TEST_CONDITION_REQUIRED(mitk::Equal(firstPlaneGeometry->GetAxisVector(2), planeGeometry->GetAxisVector(2), slicedGeometryEps), ""); + + auto lastPlaneGeometry = slicedGeometry->GetPlaneGeometry(numberOfSlices - 1); + auto expectedOriginOfLastSlice = createPoint(0.0, 0.0, thicknessInMM * (numberOfSlices - 1)); + + MITK_TEST_OUTPUT(<< "Check if origin of last PlaneGeometry is at expected location"); + MITK_TEST_CONDITION_REQUIRED(mitk::Equal(lastPlaneGeometry->GetOrigin(), expectedOriginOfLastSlice, slicedGeometryEps), ""); mitkSlicedGeometry3D_ChangeImageGeometryConsideringOriginOffset_Test(); std::cout << "[TEST DONE]" << std::endl; return EXIT_SUCCESS; } diff --git a/Modules/QtWidgetsExt/include/QmitkSliderNavigatorWidget.h b/Modules/QtWidgetsExt/include/QmitkSliderNavigatorWidget.h index a01edca994..99d5e3fd59 100644 --- a/Modules/QtWidgetsExt/include/QmitkSliderNavigatorWidget.h +++ b/Modules/QtWidgetsExt/include/QmitkSliderNavigatorWidget.h @@ -1,117 +1,126 @@ /*=================================================================== 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 QMITKSLIDERNAVIGATORWIDGET_H_ #define QMITKSLIDERNAVIGATORWIDGET_H_ #include "MitkQtWidgetsExtExports.h" #include "ui_QmitkSliderNavigator.h" #include #include #include class MITKQTWIDGETSEXT_EXPORT QmitkSliderNavigatorWidget : public QWidget, public Ui::QmitkSliderNavigator { Q_OBJECT public: QmitkSliderNavigatorWidget(QWidget *parent = nullptr, Qt::WindowFlags f = nullptr); QString GetLabelUnit(); /** * \brief Converts the passed value to a QString representation. * * If the value exceeds a certain maximum, "INF" (for "infinity") is displayed * instead. */ QString ClippedValueToString(float value); /** * \brief Returns range-minimum (displayed as label left of slider if enabled) */ QString GetMinValueLabel(); QString GetMaxValueLabel(); int GetPos(); + bool GetInverseDirection() const; + + bool GetInvertedControls() const; + public slots: /** * \brief Updates the slider with the recent changes applied to the navigator. * * Intended to be called via event mechanism, e.g. if the connected * mitk::Stepper is modified. */ void Refetch(); void SetStepper(mitk::Stepper *stepper); void ShowLabels(bool show); /** * \brief En-/disables displaying of the unit label (range will be displayed * without unit if enabled). */ void ShowLabelUnit(bool show); void SetPos(int val); void SetInverseDirection(bool inverseDirection); + void SetInvertedControls(bool invertedControls); + protected slots: - void slider_valueChanged(int); + void slider_valueChanged(double); /** * \brief Set range minimum and maximum (displayed as labels left and right * of slider if enabled) */ void SetLabelValues(float min, float max); void SetLabelValuesValid(bool minValid, bool maxValid); /** * \brief Set range unit (e.g. mm or ms) which will be displayed below range * labels if enabled. */ void SetLabelUnit(const char *unit); /** * \brief Configure slider with labels according to range and unit settings */ void SetLabels(); - void spinBox_valueChanged(int); + void spinBox_valueChanged(double); + protected: bool m_HasLabelUnit; bool m_MaxValueValid; bool m_MinValueValid; QString m_LabelUnit; mitk::Stepper::Pointer m_Stepper; bool m_InRefetch; bool m_HasLabels; float m_MinValue; float m_MaxValue; bool m_InverseDirection; + bool m_InvertedControls; + }; #endif diff --git a/Modules/QtWidgetsExt/src/QmitkSliderNavigator.ui b/Modules/QtWidgetsExt/src/QmitkSliderNavigator.ui index 458316b523..fcb08fc3da 100644 --- a/Modules/QtWidgetsExt/src/QmitkSliderNavigator.ui +++ b/Modules/QtWidgetsExt/src/QmitkSliderNavigator.ui @@ -1,191 +1,151 @@ QmitkSliderNavigator 0 0 235 137 0 0 0 0 32767 32767 QmitkSliderNavigator 2 2 2 2 0 - + - + 0 0 - - - 50 - 0 - - - - 0 - - - - - - 0 - 0 - - - - - 0 - 0 - - - - 0 - - - Qt::Horizontal + + + + true - - - - 0 0 45 0 45 32767 <p align="center">XXX</p> false - - - - true + + + + + 0 + 0 + + + + + 0 0 45 0 45 32767 <p align="center">XXX</p> false + + + ctkDoubleSlider + QWidget +
ctkDoubleSlider.h
+ 1 + + + ctkDoubleSpinBox + QWidget +
ctkDoubleSpinBox.h
+ 1 + + mitkStepper.h - - - m_Slider - valueChanged(int) - QmitkSliderNavigator - slider_valueChanged(int) - - - 20 - 20 - - - 20 - 20 - - - - - m_SpinBox - valueChanged(int) - QmitkSliderNavigator - spinBox_valueChanged(int) - - - 20 - 20 - - - 20 - 20 - - - - + diff --git a/Modules/QtWidgetsExt/src/QmitkSliderNavigatorWidget.cpp b/Modules/QtWidgetsExt/src/QmitkSliderNavigatorWidget.cpp index 9c77986ca1..c0b1918479 100644 --- a/Modules/QtWidgetsExt/src/QmitkSliderNavigatorWidget.cpp +++ b/Modules/QtWidgetsExt/src/QmitkSliderNavigatorWidget.cpp @@ -1,259 +1,311 @@ /*=================================================================== 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 "QmitkSliderNavigatorWidget.h" QmitkSliderNavigatorWidget::QmitkSliderNavigatorWidget(QWidget *parent, Qt::WindowFlags f) : QWidget(parent, f) { this->setupUi(this); + m_Slider->setOrientation(Qt::Horizontal); + m_Slider->setMinimum(0); + m_Slider->setMaximum(0); + m_Slider->setValue(0); + m_Slider->setSingleStep(1); + m_Slider->setPageStep(10); + + m_SpinBox->setMinimum( 0 ); + m_SpinBox->setMaximum(0); + m_SpinBox->setValue(0); + m_SpinBox->setDecimals(0); + m_SpinBox->setSingleStep(1); + + this->connect(m_Slider, SIGNAL(valueChanged(double)), SLOT(slider_valueChanged(double))); + this->connect(m_SpinBox, SIGNAL(valueChanged(double)), SLOT(spinBox_valueChanged(double))); + // this avoids trying to use m_Stepper until it is set to something != NULL // (additionally to the avoiding recursions during refetching) m_InRefetch = true; // Hide slider labeling -- until it is explicitly activated this->ShowLabels(false); // Set label values as invalid (N/A) this->SetLabelValuesValid(false, false); m_HasLabels = false; m_HasLabelUnit = true; m_InverseDirection = false; + m_InvertedControls = false; } void QmitkSliderNavigatorWidget::Refetch() { if (!m_InRefetch) { m_InRefetch = true; - m_Slider->setMinimum(0); - m_Slider->setMaximum(m_Stepper->GetSteps() - 1); - if (m_InverseDirection) + if (m_Stepper->GetSteps() == 0) { - m_Slider->setValue(m_Stepper->GetSteps() - 1 - m_Stepper->GetPos()); + m_Slider->setMaximum(0); + m_Slider->setValue(0); + m_SpinBox->setMaximum(0); + m_SpinBox->setValue(0); } else { - m_Slider->setValue(m_Stepper->GetPos()); - } + m_Slider->setMaximum(m_Stepper->GetSteps() - 1); + if (m_InverseDirection) + { + m_Slider->setValue(m_Stepper->GetSteps() - 1 - m_Stepper->GetPos()); + } + else + { + m_Slider->setValue(m_Stepper->GetPos()); + } - m_SpinBox->setMinimum(0); - m_SpinBox->setMaximum(m_Stepper->GetSteps() - 1); - if (m_InverseDirection) - { - m_SpinBox->setValue(m_Stepper->GetSteps() - 1 - m_Stepper->GetPos()); - } - else - { - m_SpinBox->setValue(m_Stepper->GetPos()); + m_SpinBox->setMaximum(m_Stepper->GetSteps() - 1); + if (m_InverseDirection) + { + m_SpinBox->setValue(m_Stepper->GetSteps() - 1 - m_Stepper->GetPos()); + } + else + { + m_SpinBox->setValue(m_Stepper->GetPos()); + } } if (m_Stepper->HasRange() && m_HasLabels) { // Show slider with labels according to below settings m_SliderLabelLeft->setHidden(false); m_SliderLabelRight->setHidden(false); if (m_Stepper->HasValidRange()) { this->SetLabelValuesValid(true, true); this->SetLabelValues(m_Stepper->GetRangeMin(), m_Stepper->GetRangeMax()); } else { this->SetLabelValuesValid(false, false); } if (m_Stepper->HasUnitName()) { this->SetLabelUnit(m_Stepper->GetUnitName()); } } else { // Show slider without any labels m_SliderLabelLeft->setHidden(true); m_SliderLabelRight->setHidden(true); } // Update GUI according to above settings this->SetLabels(); m_InRefetch = false; } } void QmitkSliderNavigatorWidget::SetStepper(mitk::Stepper *stepper) { m_Stepper = stepper; // this avoids trying to use m_Stepper until it is set to something != NULL // (additionally to the avoiding recursions during refetching) m_InRefetch = (stepper == nullptr); } -void QmitkSliderNavigatorWidget::slider_valueChanged(int) + +void QmitkSliderNavigatorWidget::slider_valueChanged(double) { if (!m_InRefetch) { if (m_InverseDirection) { m_Stepper->SetPos(m_Stepper->GetSteps() - 1 - m_Slider->value()); } else { m_Stepper->SetPos(m_Slider->value()); } this->Refetch(); } } void QmitkSliderNavigatorWidget::ShowLabels(bool show) { m_HasLabels = show; } void QmitkSliderNavigatorWidget::ShowLabelUnit(bool show) { m_HasLabelUnit = show; } void QmitkSliderNavigatorWidget::SetLabelValues(float min, float max) { m_MinValue = min; m_MaxValue = max; } void QmitkSliderNavigatorWidget::SetLabelValuesValid(bool minValid, bool maxValid) { m_MinValueValid = minValid; m_MaxValueValid = maxValid; } void QmitkSliderNavigatorWidget::SetLabelUnit(const char *unit) { m_LabelUnit = unit; } QString QmitkSliderNavigatorWidget::GetLabelUnit() { return m_LabelUnit; } QString QmitkSliderNavigatorWidget::ClippedValueToString(float value) { if (value < -10000000.0) { return "-INF"; } else if (value > 10000000.0) { return "+INF"; } else { return QString::number(value, 'f', 2); } } QString QmitkSliderNavigatorWidget::GetMinValueLabel() { if (m_MinValueValid) { return this->ClippedValueToString(m_MinValue); } else { return "N/A"; } } QString QmitkSliderNavigatorWidget::GetMaxValueLabel() { if (m_MaxValueValid) { return this->ClippedValueToString(m_MaxValue); } else { return "N/A"; } } void QmitkSliderNavigatorWidget::SetLabels() { QString minText = "

" + this->GetMinValueLabel(); QString maxText = "

" + this->GetMaxValueLabel(); if (m_HasLabelUnit) { minText += " " + this->GetLabelUnit(); maxText += " " + this->GetLabelUnit(); } if (m_MinValueValid) { minText += "
(pos 0)"; } if (m_MaxValueValid) { maxText += "
(pos " + QString::number(m_Stepper->GetSteps() - 1) + ")"; } minText += "

"; maxText += "

"; m_SliderLabelLeft->setText(minText); m_SliderLabelRight->setText(maxText); } -void QmitkSliderNavigatorWidget::spinBox_valueChanged(int) +void QmitkSliderNavigatorWidget::spinBox_valueChanged(double) { if (!m_InRefetch) { if (m_InverseDirection) { m_Stepper->SetPos(m_Stepper->GetSteps() - 1 - m_SpinBox->value()); } else { m_Stepper->SetPos(m_SpinBox->value()); } this->Refetch(); } } int QmitkSliderNavigatorWidget::GetPos() { return m_Stepper->GetPos(); } void QmitkSliderNavigatorWidget::SetPos(int val) { if (!m_InRefetch) { m_Stepper->SetPos(val); } } +bool QmitkSliderNavigatorWidget::GetInverseDirection() const +{ + return m_InverseDirection; +} + void QmitkSliderNavigatorWidget::SetInverseDirection(bool inverseDirection) { - m_InverseDirection = inverseDirection; + if (inverseDirection != m_InverseDirection) + { + m_InverseDirection = inverseDirection; + this->Refetch(); + } } + +bool QmitkSliderNavigatorWidget::GetInvertedControls() const +{ + return m_InvertedControls; +} + +void QmitkSliderNavigatorWidget::SetInvertedControls(bool invertedControls) +{ + if (invertedControls != m_InvertedControls) + { + m_InvertedControls = invertedControls; + m_Slider->setInvertedAppearance(invertedControls); + m_Slider->setInvertedControls(invertedControls); + m_SpinBox->setInvertedControls(invertedControls); + } +} + diff --git a/Plugins/org.mitk.gui.qt.imagenavigator/src/internal/QmitkImageNavigatorView.cpp b/Plugins/org.mitk.gui.qt.imagenavigator/src/internal/QmitkImageNavigatorView.cpp index 34abd87c69..cc645a0ce5 100644 --- a/Plugins/org.mitk.gui.qt.imagenavigator/src/internal/QmitkImageNavigatorView.cpp +++ b/Plugins/org.mitk.gui.qt.imagenavigator/src/internal/QmitkImageNavigatorView.cpp @@ -1,555 +1,633 @@ /*=================================================================== 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 "QmitkImageNavigatorView.h" +#include + #include #include #include #include #include #include #include #include #include #include #include const std::string QmitkImageNavigatorView::VIEW_ID = "org.mitk.views.imagenavigator"; static mitk::DataNode::Pointer GetTopLayerNode(const mitk::DataStorage::SetOfObjects::ConstPointer nodes, const mitk::Point3D &position, mitk::BaseRenderer* renderer) { mitk::DataNode::Pointer node; int maxlayer = -32768; if(nodes.IsNotNull()) { // find node with largest layer, that is the node shown on top in the render window for (unsigned int x = 0; x < nodes->size(); x++) { if ( (nodes->at(x)->GetData()->GetGeometry() != NULL) && nodes->at(x)->GetData()->GetGeometry()->IsInside(position) ) { int layer = 0; if(!(nodes->at(x)->GetIntProperty("layer", layer))) continue; if(layer > maxlayer) { if( static_cast(nodes->at(x))->IsVisible(renderer) ) { node = nodes->at(x); maxlayer = layer; } } } } } return node; } QmitkImageNavigatorView::QmitkImageNavigatorView() : m_AxialStepper(0) , m_SagittalStepper(0) , m_FrontalStepper(0) , m_TimeStepper(0) , m_Parent(0) , m_IRenderWindowPart(0) { } QmitkImageNavigatorView::~QmitkImageNavigatorView() { } void QmitkImageNavigatorView::CreateQtPartControl(QWidget *parent) { // create GUI widgets m_Parent = parent; m_Controls.setupUi(parent); - m_Controls.m_SliceNavigatorAxial->SetInverseDirection(true); connect(m_Controls.m_XWorldCoordinateSpinBox, SIGNAL(valueChanged(double)), this, SLOT(OnMillimetreCoordinateValueChanged())); connect(m_Controls.m_YWorldCoordinateSpinBox, SIGNAL(valueChanged(double)), this, SLOT(OnMillimetreCoordinateValueChanged())); connect(m_Controls.m_ZWorldCoordinateSpinBox, SIGNAL(valueChanged(double)), this, SLOT(OnMillimetreCoordinateValueChanged())); m_Parent->setEnabled(false); mitk::IRenderWindowPart* renderPart = this->GetRenderWindowPart(); this->RenderWindowPartActivated(renderPart); } void QmitkImageNavigatorView::SetFocus () { m_Controls.m_XWorldCoordinateSpinBox->setFocus(); } void QmitkImageNavigatorView::RenderWindowPartActivated(mitk::IRenderWindowPart* renderWindowPart) { if (this->m_IRenderWindowPart != renderWindowPart) { this->m_IRenderWindowPart = renderWindowPart; this->m_Parent->setEnabled(true); QmitkRenderWindow* renderWindow = renderWindowPart->GetQmitkRenderWindow("axial"); if (renderWindow) { if (m_AxialStepper) m_AxialStepper->deleteLater(); m_AxialStepper = new QmitkStepperAdapter(m_Controls.m_SliceNavigatorAxial, renderWindow->GetSliceNavigationController()->GetSlice(), "sliceNavigatorAxialFromSimpleExample"); m_Controls.m_SliceNavigatorAxial->setEnabled(true); m_Controls.m_AxialLabel->setEnabled(true); m_Controls.m_ZWorldCoordinateSpinBox->setEnabled(true); connect(m_AxialStepper, SIGNAL(Refetch()), this, SLOT(OnRefetch())); connect(m_AxialStepper, SIGNAL(Refetch()), this, SLOT(UpdateStatusBar())); } else { m_Controls.m_SliceNavigatorAxial->setEnabled(false); m_Controls.m_AxialLabel->setEnabled(false); m_Controls.m_ZWorldCoordinateSpinBox->setEnabled(false); } renderWindow = renderWindowPart->GetQmitkRenderWindow("sagittal"); if (renderWindow) { if (m_SagittalStepper) m_SagittalStepper->deleteLater(); m_SagittalStepper = new QmitkStepperAdapter(m_Controls.m_SliceNavigatorSagittal, renderWindow->GetSliceNavigationController()->GetSlice(), "sliceNavigatorSagittalFromSimpleExample"); m_Controls.m_SliceNavigatorSagittal->setEnabled(true); m_Controls.m_SagittalLabel->setEnabled(true); m_Controls.m_YWorldCoordinateSpinBox->setEnabled(true); connect(m_SagittalStepper, SIGNAL(Refetch()), this, SLOT(OnRefetch())); connect(m_SagittalStepper, SIGNAL(Refetch()), this, SLOT(UpdateStatusBar())); } else { m_Controls.m_SliceNavigatorSagittal->setEnabled(false); m_Controls.m_SagittalLabel->setEnabled(false); m_Controls.m_YWorldCoordinateSpinBox->setEnabled(false); } renderWindow = renderWindowPart->GetQmitkRenderWindow("coronal"); if (renderWindow) { if (m_FrontalStepper) m_FrontalStepper->deleteLater(); m_FrontalStepper = new QmitkStepperAdapter(m_Controls.m_SliceNavigatorFrontal, renderWindow->GetSliceNavigationController()->GetSlice(), "sliceNavigatorFrontalFromSimpleExample"); m_Controls.m_SliceNavigatorFrontal->setEnabled(true); m_Controls.m_CoronalLabel->setEnabled(true); m_Controls.m_XWorldCoordinateSpinBox->setEnabled(true); connect(m_FrontalStepper, SIGNAL(Refetch()), this, SLOT(OnRefetch())); connect(m_FrontalStepper, SIGNAL(Refetch()), this, SLOT(UpdateStatusBar())); } else { m_Controls.m_SliceNavigatorFrontal->setEnabled(false); m_Controls.m_CoronalLabel->setEnabled(false); m_Controls.m_XWorldCoordinateSpinBox->setEnabled(false); } mitk::SliceNavigationController* timeController = renderWindowPart->GetTimeNavigationController(); if (timeController) { if (m_TimeStepper) m_TimeStepper->deleteLater(); m_TimeStepper = new QmitkStepperAdapter(m_Controls.m_SliceNavigatorTime, timeController->GetTime(), "sliceNavigatorTimeFromSimpleExample"); m_Controls.m_SliceNavigatorTime->setEnabled(true); m_Controls.m_TimeLabel->setEnabled(true); connect(m_TimeStepper, SIGNAL(Refetch()), this, SLOT(UpdateStatusBar())); } else { m_Controls.m_SliceNavigatorTime->setEnabled(false); m_Controls.m_TimeLabel->setEnabled(false); } + + this->OnRefetch(); + this->UpdateStatusBar(); } } void QmitkImageNavigatorView::UpdateStatusBar() { if (m_IRenderWindowPart != nullptr) { mitk::Point3D position = m_IRenderWindowPart->GetSelectedPosition(); mitk::BaseRenderer::Pointer renderer = mitk::BaseRenderer::GetInstance(m_IRenderWindowPart->GetActiveQmitkRenderWindow()->GetVtkRenderWindow()); mitk::TNodePredicateDataType::Pointer isImageData = mitk::TNodePredicateDataType::New(); mitk::DataStorage::SetOfObjects::ConstPointer nodes = this->GetDataStorage()->GetSubset(isImageData).GetPointer(); if (nodes.IsNotNull()) { mitk::Image::Pointer image3D; mitk::DataNode::Pointer node; mitk::DataNode::Pointer topSourceNode; int component = 0; node = GetTopLayerNode(nodes, position, renderer); if (node.IsNotNull()) { bool isBinary(false); node->GetBoolProperty("binary", isBinary); if (isBinary) { mitk::DataStorage::SetOfObjects::ConstPointer sourcenodes = this->GetDataStorage()->GetSources(node, NULL, true); if (!sourcenodes->empty()) { topSourceNode = GetTopLayerNode(sourcenodes, position, renderer); } if (topSourceNode.IsNotNull()) { image3D = dynamic_cast(topSourceNode->GetData()); topSourceNode->GetIntProperty("Image.Displayed Component", component); } else { image3D = dynamic_cast(node->GetData()); node->GetIntProperty("Image.Displayed Component", component); } } else { image3D = dynamic_cast(node->GetData()); node->GetIntProperty("Image.Displayed Component", component); } } // get the position and pixel value from the image and build up status bar text auto statusBar = mitk::StatusBar::GetInstance(); if (image3D.IsNotNull() && statusBar != nullptr) { itk::Index<3> p; image3D->GetGeometry()->WorldToIndex(position, p); auto pixelType = image3D->GetChannelDescriptor().GetPixelType().GetPixelType(); if (pixelType == itk::ImageIOBase::RGB || pixelType == itk::ImageIOBase::RGBA) { std::string pixelValue = "Pixel RGB(A) value: "; pixelValue.append(ConvertCompositePixelValueToString(image3D, p)); statusBar->DisplayImageInfo(position, p, renderer->GetTime(), pixelValue.c_str()); } else { itk::Index<3> p; image3D->GetGeometry()->WorldToIndex(position, p); mitk::ScalarType pixelValue; mitkPixelTypeMultiplex5( mitk::FastSinglePixelAccess, image3D->GetChannelDescriptor().GetPixelType(), image3D, image3D->GetVolumeData(renderer->GetTimeStep()), p, pixelValue, component); statusBar->DisplayImageInfo(position, p, renderer->GetTime(), pixelValue); } } else { statusBar->DisplayImageInfoInvalid(); } } } } void QmitkImageNavigatorView::RenderWindowPartDeactivated(mitk::IRenderWindowPart* /*renderWindowPart*/) { m_IRenderWindowPart = 0; m_Parent->setEnabled(false); } int QmitkImageNavigatorView::GetSizeFlags(bool width) { if(!width) { return berry::Constants::MIN | berry::Constants::MAX | berry::Constants::FILL; } else { return 0; } } int QmitkImageNavigatorView::ComputePreferredSize(bool width, int /*availableParallel*/, int /*availablePerpendicular*/, int preferredResult) { if(width==false) { return 200; } else { return preferredResult; } } int QmitkImageNavigatorView::GetClosestAxisIndex(mitk::Vector3D normal) { // cos(theta) = normal . axis // cos(theta) = (a, b, c) . (d, e, f) // cos(theta) = (a, b, c) . (1, 0, 0) = a // cos(theta) = (a, b, c) . (0, 1, 0) = b // cos(theta) = (a, b, c) . (0, 0, 1) = c double absCosThetaWithAxis[3]; for (int i = 0; i < 3; i++) { absCosThetaWithAxis[i] = fabs(normal[i]); } int largestIndex = 0; double largestValue = absCosThetaWithAxis[0]; for (int i = 1; i < 3; i++) { if (absCosThetaWithAxis[i] > largestValue) { largestValue = absCosThetaWithAxis[i]; largestIndex = i; } } return largestIndex; } void QmitkImageNavigatorView::SetBorderColors() { if (m_IRenderWindowPart) { QString decoColor; QmitkRenderWindow* renderWindow = m_IRenderWindowPart->GetQmitkRenderWindow("axial"); if (renderWindow) { decoColor = GetDecorationColorOfGeometry(renderWindow); mitk::PlaneGeometry::ConstPointer geometry = renderWindow->GetSliceNavigationController()->GetCurrentPlaneGeometry(); if (geometry.IsNotNull()) { mitk::Vector3D normal = geometry->GetNormal(); int axis = this->GetClosestAxisIndex(normal); this->SetBorderColor(axis, decoColor); } } renderWindow = m_IRenderWindowPart->GetQmitkRenderWindow("sagittal"); if (renderWindow) { decoColor = GetDecorationColorOfGeometry(renderWindow); mitk::PlaneGeometry::ConstPointer geometry = renderWindow->GetSliceNavigationController()->GetCurrentPlaneGeometry(); if (geometry.IsNotNull()) { mitk::Vector3D normal = geometry->GetNormal(); int axis = this->GetClosestAxisIndex(normal); this->SetBorderColor(axis, decoColor); } } renderWindow = m_IRenderWindowPart->GetQmitkRenderWindow("coronal"); if (renderWindow) { decoColor = GetDecorationColorOfGeometry(renderWindow); mitk::PlaneGeometry::ConstPointer geometry = renderWindow->GetSliceNavigationController()->GetCurrentPlaneGeometry(); if (geometry.IsNotNull()) { mitk::Vector3D normal = geometry->GetNormal(); int axis = this->GetClosestAxisIndex(normal); this->SetBorderColor(axis, decoColor); } } } } QString QmitkImageNavigatorView::GetDecorationColorOfGeometry(QmitkRenderWindow* renderWindow) { QColor color; float rgb[3] = {1.0f, 1.0f, 1.0f}; float rgbMax = 255.0f; mitk::BaseRenderer::GetInstance(renderWindow->GetVtkRenderWindow())->GetCurrentWorldPlaneGeometryNode()->GetColor(rgb); color.setRed(static_cast(rgb[0]*rgbMax + 0.5)); color.setGreen(static_cast(rgb[1]*rgbMax + 0.5)); color.setBlue(static_cast(rgb[2]*rgbMax + 0.5)); QString colorAsString = QString(color.name()); return colorAsString; } void QmitkImageNavigatorView::SetBorderColor(int axis, QString colorAsStyleSheetString) { if (axis == 0) { this->SetBorderColor(m_Controls.m_XWorldCoordinateSpinBox, colorAsStyleSheetString); } else if (axis == 1) { this->SetBorderColor(m_Controls.m_YWorldCoordinateSpinBox, colorAsStyleSheetString); } else if (axis == 2) { this->SetBorderColor(m_Controls.m_ZWorldCoordinateSpinBox, colorAsStyleSheetString); } } void QmitkImageNavigatorView::SetBorderColor(QDoubleSpinBox *spinBox, QString colorAsStyleSheetString) { assert(spinBox); spinBox->setStyleSheet(QString("border: 2px solid ") + colorAsStyleSheetString + ";"); } void QmitkImageNavigatorView::SetStepSizes() { this->SetStepSize(0); this->SetStepSize(1); this->SetStepSize(2); } void QmitkImageNavigatorView::SetStepSize(int axis) { if (m_IRenderWindowPart) { mitk::BaseGeometry::ConstPointer geometry = m_IRenderWindowPart->GetActiveQmitkRenderWindow()->GetSliceNavigationController()->GetInputWorldGeometry3D(); if (geometry.IsNotNull()) { mitk::Point3D crossPositionInIndexCoordinates; mitk::Point3D crossPositionInIndexCoordinatesPlus1; mitk::Point3D crossPositionInMillimetresPlus1; mitk::Vector3D transformedAxisDirection; mitk::Point3D crossPositionInMillimetres = m_IRenderWindowPart->GetSelectedPosition(); geometry->WorldToIndex(crossPositionInMillimetres, crossPositionInIndexCoordinates); crossPositionInIndexCoordinatesPlus1 = crossPositionInIndexCoordinates; crossPositionInIndexCoordinatesPlus1[axis] += 1; geometry->IndexToWorld(crossPositionInIndexCoordinatesPlus1, crossPositionInMillimetresPlus1); transformedAxisDirection = crossPositionInMillimetresPlus1 - crossPositionInMillimetres; int closestAxisInMillimetreSpace = this->GetClosestAxisIndex(transformedAxisDirection); double stepSize = transformedAxisDirection.GetNorm(); this->SetStepSize(closestAxisInMillimetreSpace, stepSize); } } } void QmitkImageNavigatorView::SetStepSize(int axis, double stepSize) { if (axis == 0) { m_Controls.m_XWorldCoordinateSpinBox->setSingleStep(stepSize); } else if (axis == 1) { m_Controls.m_YWorldCoordinateSpinBox->setSingleStep(stepSize); } else if (axis == 2) { m_Controls.m_ZWorldCoordinateSpinBox->setSingleStep(stepSize); } } void QmitkImageNavigatorView::OnMillimetreCoordinateValueChanged() { if (m_IRenderWindowPart) { mitk::TimeGeometry::ConstPointer geometry = m_IRenderWindowPart->GetActiveQmitkRenderWindow()->GetSliceNavigationController()->GetInputWorldTimeGeometry(); if (geometry.IsNotNull()) { mitk::Point3D positionInWorldCoordinates; positionInWorldCoordinates[0] = m_Controls.m_XWorldCoordinateSpinBox->value(); positionInWorldCoordinates[1] = m_Controls.m_YWorldCoordinateSpinBox->value(); positionInWorldCoordinates[2] = m_Controls.m_ZWorldCoordinateSpinBox->value(); m_IRenderWindowPart->SetSelectedPosition(positionInWorldCoordinates); } } } void QmitkImageNavigatorView::OnRefetch() { if (m_IRenderWindowPart) { mitk::BaseGeometry::ConstPointer geometry = m_IRenderWindowPart->GetActiveQmitkRenderWindow()->GetSliceNavigationController()->GetInputWorldGeometry3D(); mitk::TimeGeometry::ConstPointer timeGeometry = m_IRenderWindowPart->GetActiveQmitkRenderWindow()->GetSliceNavigationController()->GetInputWorldTimeGeometry(); if (geometry.IsNull() && timeGeometry.IsNotNull()) { mitk::TimeStepType timeStep = m_IRenderWindowPart->GetActiveQmitkRenderWindow()->GetSliceNavigationController()->GetTime()->GetPos(); geometry = timeGeometry->GetGeometryForTimeStep(timeStep); } if (geometry.IsNotNull()) { mitk::BoundingBox::BoundsArrayType bounds = geometry->GetBounds(); mitk::Point3D cornerPoint1InIndexCoordinates; cornerPoint1InIndexCoordinates[0] = bounds[0]; cornerPoint1InIndexCoordinates[1] = bounds[2]; cornerPoint1InIndexCoordinates[2] = bounds[4]; mitk::Point3D cornerPoint2InIndexCoordinates; cornerPoint2InIndexCoordinates[0] = bounds[1]; cornerPoint2InIndexCoordinates[1] = bounds[3]; cornerPoint2InIndexCoordinates[2] = bounds[5]; if (!geometry->GetImageGeometry()) { cornerPoint1InIndexCoordinates[0] += 0.5; cornerPoint1InIndexCoordinates[1] += 0.5; cornerPoint1InIndexCoordinates[2] += 0.5; cornerPoint2InIndexCoordinates[0] -= 0.5; cornerPoint2InIndexCoordinates[1] -= 0.5; cornerPoint2InIndexCoordinates[2] -= 0.5; } mitk::Point3D crossPositionInWorldCoordinates = m_IRenderWindowPart->GetSelectedPosition(); mitk::Point3D cornerPoint1InWorldCoordinates; mitk::Point3D cornerPoint2InWorldCoordinates; geometry->IndexToWorld(cornerPoint1InIndexCoordinates, cornerPoint1InWorldCoordinates); geometry->IndexToWorld(cornerPoint2InIndexCoordinates, cornerPoint2InWorldCoordinates); m_Controls.m_XWorldCoordinateSpinBox->blockSignals(true); m_Controls.m_YWorldCoordinateSpinBox->blockSignals(true); m_Controls.m_ZWorldCoordinateSpinBox->blockSignals(true); m_Controls.m_XWorldCoordinateSpinBox->setMinimum(std::min(cornerPoint1InWorldCoordinates[0], cornerPoint2InWorldCoordinates[0])); m_Controls.m_YWorldCoordinateSpinBox->setMinimum(std::min(cornerPoint1InWorldCoordinates[1], cornerPoint2InWorldCoordinates[1])); m_Controls.m_ZWorldCoordinateSpinBox->setMinimum(std::min(cornerPoint1InWorldCoordinates[2], cornerPoint2InWorldCoordinates[2])); m_Controls.m_XWorldCoordinateSpinBox->setMaximum(std::max(cornerPoint1InWorldCoordinates[0], cornerPoint2InWorldCoordinates[0])); m_Controls.m_YWorldCoordinateSpinBox->setMaximum(std::max(cornerPoint1InWorldCoordinates[1], cornerPoint2InWorldCoordinates[1])); m_Controls.m_ZWorldCoordinateSpinBox->setMaximum(std::max(cornerPoint1InWorldCoordinates[2], cornerPoint2InWorldCoordinates[2])); m_Controls.m_XWorldCoordinateSpinBox->setValue(crossPositionInWorldCoordinates[0]); m_Controls.m_YWorldCoordinateSpinBox->setValue(crossPositionInWorldCoordinates[1]); m_Controls.m_ZWorldCoordinateSpinBox->setValue(crossPositionInWorldCoordinates[2]); m_Controls.m_XWorldCoordinateSpinBox->blockSignals(false); m_Controls.m_YWorldCoordinateSpinBox->blockSignals(false); m_Controls.m_ZWorldCoordinateSpinBox->blockSignals(false); + + /// Calculating 'inverse direction' property. + + mitk::AffineTransform3D::MatrixType matrix = geometry->GetIndexToWorldTransform()->GetMatrix(); + matrix.GetVnlMatrix().normalize_columns(); + mitk::AffineTransform3D::MatrixType::InternalMatrixType inverseMatrix = matrix.GetInverse(); + + for (int worldAxis = 0; worldAxis < 3; ++worldAxis) + { + QmitkRenderWindow* renderWindow = + worldAxis == 0 ? m_IRenderWindowPart->GetQmitkRenderWindow("sagittal") : + worldAxis == 1 ? m_IRenderWindowPart->GetQmitkRenderWindow("coronal") : + m_IRenderWindowPart->GetQmitkRenderWindow("axial"); + + if (renderWindow) + { + const mitk::BaseGeometry* rendererGeometry = renderWindow->GetRenderer()->GetCurrentWorldGeometry(); + + /// Because of some problems with the current way of event signalling, + /// 'Modified' events are sent out from the stepper while the renderer + /// does not have a geometry yet. Therefore, we do a nullptr check here. + /// See bug T22122. This check can be resolved after T22122 got fixed. + if (rendererGeometry) + { + int dominantAxis = itk::Function::Max3( + inverseMatrix[0][worldAxis], + inverseMatrix[1][worldAxis], + inverseMatrix[2][worldAxis]); + + bool referenceGeometryAxisInverted = inverseMatrix[dominantAxis][worldAxis] < 0; + bool rendererZAxisInverted = rendererGeometry->GetAxisVector(2)[worldAxis] < 0; + + /// `referenceGeometryAxisInverted` tells if the direction of the corresponding axis + /// of the reference geometry is flipped compared to the 'world direction' or not. + /// + /// `rendererZAxisInverted` tells if direction of the renderer geometry z axis is + /// flipped compared to the 'world direction' or not. This is the same as the indexing + /// direction in the slice navigation controller and matches the 'top' property when + /// initialising the renderer planes. (If 'top' was true then the direction is + /// inverted.) + /// + /// The world direction can be +1 ('up') that means right, anterior or superior, or + /// it can be -1 ('down') that means left, posterior or inferior, respectively. + /// + /// If these two do not match, we have to invert the index between the slice navigation + /// controller and the slider navigator widget, so that the user can see and control + /// the index according to the reference geometry, rather than the slice navigation + /// controller. The index in the slice navigation controller depends on in which way + /// the reference geometry has been resliced for the renderer, and it does not necessarily + /// match neither the world direction, nor the direction of the corresponding axis of + /// the reference geometry. Hence, it is a merely internal information that should not + /// be exposed to the GUI. + /// + /// So that one can navigate in the same world direction by dragging the slider + /// right, regardless of the direction of the corresponding axis of the reference + /// geometry, we invert the direction of the controls if the reference geometry axis + /// is inverted but the direction is not ('inversDirection' is false) or the other + /// way around. + + bool inverseDirection = referenceGeometryAxisInverted != rendererZAxisInverted; + + QmitkSliderNavigatorWidget* navigatorWidget = + worldAxis == 0 ? m_Controls.m_SliceNavigatorSagittal : + worldAxis == 1 ? m_Controls.m_SliceNavigatorFrontal : + m_Controls.m_SliceNavigatorAxial; + + navigatorWidget->SetInverseDirection(inverseDirection); + + // This should be a preference (see T22254) + // bool invertedControls = referenceGeometryAxisInverted != inverseDirection; + // navigatorWidget->SetInvertedControls(invertedControls); + } + } + } } this->SetBorderColors(); } }