changeset 1766:ae18b5cfbd3e

VolumeImageGeometry::DetectProjection() reports if normal is opposite
author Sebastien Jodogne <s.jodogne@gmail.com>
date Tue, 11 May 2021 11:31:34 +0200
parents 177f93c8371c
children e7b4479dea6f
files OrthancStone/Sources/Toolbox/SlicesSorter.cpp OrthancStone/Sources/Toolbox/SlicesSorter.h OrthancStone/Sources/Volumes/VolumeImageGeometry.cpp OrthancStone/Sources/Volumes/VolumeImageGeometry.h UnitTestsSources/GeometryToolboxTests.cpp
diffstat 5 files changed, 36 insertions(+), 54 deletions(-) [+]
line wrap: on
line diff
--- a/OrthancStone/Sources/Toolbox/SlicesSorter.cpp	Tue May 11 11:01:22 2021 +0200
+++ b/OrthancStone/Sources/Toolbox/SlicesSorter.cpp	Tue May 11 11:31:34 2021 +0200
@@ -259,38 +259,6 @@
   }
 
 
-  bool SlicesSorter::LookupClosestSlice(size_t& index,
-                                        double& distance,
-                                        const CoordinateSystem3D& slice) const
-  {
-    // TODO Turn this linear-time lookup into a log-time lookup,
-    // keeping track of whether the slices are sorted along the normal
-
-    bool found = false;
-    
-    distance = std::numeric_limits<double>::infinity();
-    
-    for (size_t i = 0; i < slices_.size(); i++)
-    {
-      assert(slices_[i] != NULL);
-
-      double tmp;
-      if (CoordinateSystem3D::ComputeDistance(tmp, slices_[i]->GetGeometry(), slice))
-      {
-        if (!found ||
-            tmp < distance)
-        {
-          index = i;
-          distance = tmp;
-          found = true;
-        }
-      }
-    }
-
-    return found;
-  }
-
-
   bool SlicesSorter::ComputeSpacingBetweenSlices(double& spacing /* out */) const
   {
     if (GetSlicesCount() <= 1)
--- a/OrthancStone/Sources/Toolbox/SlicesSorter.h	Tue May 11 11:01:22 2021 +0200
+++ b/OrthancStone/Sources/Toolbox/SlicesSorter.h	Tue May 11 11:31:34 2021 +0200
@@ -86,11 +86,6 @@
     // parallel to the reference normal that will be selected.
     bool Sort();
 
-    // TODO - Remove this
-    bool LookupClosestSlice(size_t& index,
-                            double& distance,
-                            const CoordinateSystem3D& slice) const;
-
     // WARNING - The slices must have been sorted before calling this method
     bool ComputeSpacingBetweenSlices(double& spacing /* out */) const;
 
--- a/OrthancStone/Sources/Volumes/VolumeImageGeometry.cpp	Tue May 11 11:01:22 2021 +0200
+++ b/OrthancStone/Sources/Volumes/VolumeImageGeometry.cpp	Tue May 11 11:31:34 2021 +0200
@@ -229,19 +229,20 @@
 
 
   bool VolumeImageGeometry::DetectProjection(VolumeProjection& projection,
+                                             bool& isOpposite,
                                              const Vector& planeNormal) const
   {
-    if (GeometryToolbox::IsParallel(planeNormal, axialGeometry_.GetNormal()))
+    if (GeometryToolbox::IsParallelOrOpposite(isOpposite, planeNormal, axialGeometry_.GetNormal()))
     {
       projection = VolumeProjection_Axial;
       return true;
     }
-    else if (GeometryToolbox::IsParallel(planeNormal, coronalGeometry_.GetNormal()))
+    else if (GeometryToolbox::IsParallelOrOpposite(isOpposite, planeNormal, coronalGeometry_.GetNormal()))
     {
       projection = VolumeProjection_Coronal;
       return true;
     }
-    else if (GeometryToolbox::IsParallel(planeNormal, sagittalGeometry_.GetNormal()))
+    else if (GeometryToolbox::IsParallelOrOpposite(isOpposite, planeNormal, sagittalGeometry_.GetNormal()))
     {
       projection = VolumeProjection_Sagittal;
       return true;
@@ -257,7 +258,8 @@
                                         unsigned int& slice,
                                         const CoordinateSystem3D& plane) const
   {
-    if (!DetectProjection(projection, plane.GetNormal()))
+    bool isOpposite;
+    if (!DetectProjection(projection, isOpposite, plane.GetNormal()))
     {
       return false;
     }
--- a/OrthancStone/Sources/Volumes/VolumeImageGeometry.h	Tue May 11 11:01:22 2021 +0200
+++ b/OrthancStone/Sources/Volumes/VolumeImageGeometry.h	Tue May 11 11:31:34 2021 +0200
@@ -92,8 +92,8 @@
     }
 
     void SetSizeInVoxels(unsigned int width,
-                 unsigned int height,
-                 unsigned int depth);
+                         unsigned int height,
+                         unsigned int depth);
 
     // Set the geometry of the first axial slice (i.e. the one whose
     // depth == 0)
@@ -118,6 +118,7 @@
                           float z) const;
 
     bool DetectProjection(VolumeProjection& projection,
+                          bool& isOpposite,
                           const Vector& planeNormal) const;
 
     /**
--- a/UnitTestsSources/GeometryToolboxTests.cpp	Tue May 11 11:01:22 2021 +0200
+++ b/UnitTestsSources/GeometryToolboxTests.cpp	Tue May 11 11:31:34 2021 +0200
@@ -666,12 +666,16 @@
   ASSERT_DOUBLE_EQ(-3.0 / 2.0 + 30.0 * 3.0, p[2]);
 
   VolumeProjection proj;
-  ASSERT_TRUE(g.DetectProjection(proj, g.GetAxialGeometry().GetNormal()));
+  bool isOpposite;
+  ASSERT_TRUE(g.DetectProjection(proj, isOpposite, g.GetAxialGeometry().GetNormal()));
   ASSERT_EQ(VolumeProjection_Axial, proj);
-  ASSERT_TRUE(g.DetectProjection(proj, g.GetCoronalGeometry().GetNormal()));
+  ASSERT_FALSE(isOpposite);
+  ASSERT_TRUE(g.DetectProjection(proj, isOpposite, g.GetCoronalGeometry().GetNormal()));
   ASSERT_EQ(VolumeProjection_Coronal, proj);
-  ASSERT_TRUE(g.DetectProjection(proj, g.GetSagittalGeometry().GetNormal()));
+  ASSERT_FALSE(isOpposite);
+  ASSERT_TRUE(g.DetectProjection(proj, isOpposite, g.GetSagittalGeometry().GetNormal()));
   ASSERT_EQ(VolumeProjection_Sagittal, proj);
+  ASSERT_FALSE(isOpposite);
 
   ASSERT_EQ(10u, g.GetProjectionWidth(VolumeProjection_Axial));
   ASSERT_EQ(20u, g.GetProjectionHeight(VolumeProjection_Axial));
@@ -1034,6 +1038,7 @@
     double spacingZ;
     ASSERT_TRUE(slices.ComputeSpacingBetweenSlices(spacingZ));
     ASSERT_FLOAT_EQ(1, spacingZ);
+    ASSERT_TRUE(slices.AreAllSlicesDistinct());
 
     ASSERT_FLOAT_EQ(300.0, slices.GetSliceGeometry(0).GetOrigin() [0]);
     ASSERT_FLOAT_EQ(302.5, slices.GetSliceGeometry(0).GetOrigin() [1]);
@@ -1056,20 +1061,31 @@
     geometry.SetVoxelDimensions(spacingXY, spacingXY, 1);
     geometry.SetAxialGeometry(slices.GetSliceGeometry(0));
 
-    OrthancStone::Vector p, q;
-
-    p = OrthancStone::LinearAlgebra::CreateVector(0, 0, 0, 1);
-    q = OrthancStone::LinearAlgebra::Product(geometry.GetTransform(), p);
+    OrthancStone::Vector q = geometry.GetCoordinates(0, 0, 0);
+    ASSERT_EQ(3u, q.size());
     ASSERT_FLOAT_EQ(300 + spacingXY / 2.0, q[0]);
     ASSERT_FLOAT_EQ(302.5 + spacingXY / 2.0, q[1]);
     ASSERT_FLOAT_EQ(323.11 - 0.5, q[2]);
-    ASSERT_FLOAT_EQ(1, q[3]);
 
-    p = OrthancStone::LinearAlgebra::CreateVector(1, 1, 1, 1);
-    q = OrthancStone::LinearAlgebra::Product(geometry.GetTransform(), p);
+    q = geometry.GetCoordinates(1, 1, 1);
+    ASSERT_EQ(3u, q.size());
     ASSERT_FLOAT_EQ(300 - double(width) * spacingXY + spacingXY / 2.0, q[0]);
     ASSERT_FLOAT_EQ(302.5 - double(height) * spacingXY + spacingXY / 2.0, q[1]);
     ASSERT_FLOAT_EQ(323.11 + 368.0 - 0.5, q[2]);
-    ASSERT_FLOAT_EQ(1, q[3]);
+
+    OrthancStone::VolumeProjection projection;
+    unsigned int slice;
+    ASSERT_TRUE(geometry.DetectSlice(projection, slice, OrthancStone::CoordinateSystem3D("300\\302.5\\690.11", "-1\\0\\0\\0\\-1\\0")));
+    ASSERT_EQ(OrthancStone::VolumeProjection_Axial, projection);
+    ASSERT_EQ(367u, slice);
+    ASSERT_TRUE(geometry.DetectSlice(projection, slice, OrthancStone::CoordinateSystem3D("300\\302.5\\323.11", "-1\\0\\0\\0\\-1\\0")));
+    ASSERT_EQ(OrthancStone::VolumeProjection_Axial, projection);
+    ASSERT_EQ(0u, slice);
+
+    // DOSE instance: RD1.2.752.243.1.1.20210202150624529.3790.85357_DoseTPS.dcm
+    ASSERT_TRUE(geometry.DetectSlice(projection, slice,
+                                     OrthancStone::CoordinateSystem3D("-217.0492\\-161.4141\\376.61", "1\\0\\0\\0\\1\\0")));
+    ASSERT_EQ(OrthancStone::VolumeProjection_Axial, projection);
+    ASSERT_EQ(376u - 323u, slice);
   }
 }