Mercurial > hg > orthanc-stone

changeset 1874:08f2476e8f5e

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
added classes OrientedIntegerLine2D and RectanglesIntegerProjection
author Sebastien Jodogne <s.jodogne@gmail.com>
date Tue, 11 Jan 2022 18:58:37 +0100
parents e0966648ebd0
children b896f20d24ca
files OrthancStone/Resources/CMake/OrthancStoneConfiguration.cmake OrthancStone/Sources/Toolbox/Internals/OrientedIntegerLine2D.cpp OrthancStone/Sources/Toolbox/Internals/OrientedIntegerLine2D.h OrthancStone/Sources/Toolbox/Internals/RectanglesIntegerProjection.cpp OrthancStone/Sources/Toolbox/Internals/RectanglesIntegerProjection.h UnitTestsSources/ComputationalGeometryTests.cpp
diffstat 6 files changed, 611 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
--- a/OrthancStone/Resources/CMake/OrthancStoneConfiguration.cmake	Tue Jan 11 15:36:04 2022 +0100
+++ b/OrthancStone/Resources/CMake/OrthancStoneConfiguration.cmake	Tue Jan 11 18:58:37 2022 +0100
@@ -424,6 +424,10 @@
   ${ORTHANC_STONE_ROOT}/Toolbox/ImageGeometry.h
   ${ORTHANC_STONE_ROOT}/Toolbox/ImageToolbox.cpp
   ${ORTHANC_STONE_ROOT}/Toolbox/ImageToolbox.h
+  ${ORTHANC_STONE_ROOT}/Toolbox/Internals/OrientedIntegerLine2D.cpp
+  ${ORTHANC_STONE_ROOT}/Toolbox/Internals/OrientedIntegerLine2D.h
+  ${ORTHANC_STONE_ROOT}/Toolbox/Internals/RectanglesIntegerProjection.cpp
+  ${ORTHANC_STONE_ROOT}/Toolbox/Internals/RectanglesIntegerProjection.h
   ${ORTHANC_STONE_ROOT}/Toolbox/LinearAlgebra.cpp
   ${ORTHANC_STONE_ROOT}/Toolbox/LinearAlgebra.h
   ${ORTHANC_STONE_ROOT}/Toolbox/PixelTestPatterns.h
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/OrthancStone/Sources/Toolbox/Internals/OrientedIntegerLine2D.cpp	Tue Jan 11 18:58:37 2022 +0100
@@ -0,0 +1,99 @@
+/**
+ * Stone of Orthanc
+ * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics
+ * Department, University Hospital of Liege, Belgium
+ * Copyright (C) 2017-2022 Osimis S.A., Belgium
+ * Copyright (C) 2021-2022 Sebastien Jodogne, ICTEAM UCLouvain, Belgium
+ *
+ * This program is free software: you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation, either version 3 of
+ * the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this program. If not, see
+ * <http://www.gnu.org/licenses/>.
+ **/
+
+
+#include "OrientedIntegerLine2D.h"
+
+#include <cassert>
+#include <map>
+
+
+namespace OrthancStone
+{
+  namespace Internals
+  {
+    OrientedIntegerLine2D::OrientedIntegerLine2D(size_t x1,
+                                                 size_t y1,
+                                                 size_t x2,
+                                                 size_t y2) :
+      x1_(x1),
+      y1_(y1),
+      x2_(x2),
+      y2_(y2)
+    {
+    }
+
+    void OrientedIntegerLine2D::ExtractChains(std::list<Chain>& chains,
+                                              const std::vector<OrientedIntegerLine2D>& edges)
+    {
+      chains.clear();
+      
+      typedef std::map< std::pair<size_t, size_t>, std::list<size_t> > Index;
+      Index index;
+
+      for (size_t i = 0; i < edges.size(); i++)
+      {
+        index[std::make_pair(edges[i].GetX1(), edges[i].GetY1())].push_back(i);
+      }
+
+      std::vector<bool> visited(edges.size(), false);
+
+      for (size_t i = 0; i < edges.size(); i++)
+      {
+        if (!visited[i])
+        {
+          Chain chain;
+          
+          size_t current = i;
+          
+          chain.push_back(std::make_pair(edges[current].GetX1(), edges[current].GetY1()));
+          
+          for (;;)
+          {
+            visited[current] = true;
+
+            std::pair<size_t, size_t> start(edges[current].GetX1(), edges[current].GetY1());
+            std::pair<size_t, size_t> end(edges[current].GetX2(), edges[current].GetY2());
+            
+            assert(index.find(start) != index.end());
+            index[start].pop_back();
+            
+            chain.push_back(end);
+            Index::iterator next = index.find(end);
+
+            if (next == index.end() /* non-closed chain */ ||
+                next->second.empty())
+            {
+              break;
+            }
+            else
+            {
+              current = next->second.back();
+            }
+          }
+          
+          chains.push_back(chain);
+        }
+      }
+    }
+  }
+}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/OrthancStone/Sources/Toolbox/Internals/OrientedIntegerLine2D.h	Tue Jan 11 18:58:37 2022 +0100
@@ -0,0 +1,98 @@
+/**
+ * Stone of Orthanc
+ * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics
+ * Department, University Hospital of Liege, Belgium
+ * Copyright (C) 2017-2022 Osimis S.A., Belgium
+ * Copyright (C) 2021-2022 Sebastien Jodogne, ICTEAM UCLouvain, Belgium
+ *
+ * This program is free software: you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation, either version 3 of
+ * the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this program. If not, see
+ * <http://www.gnu.org/licenses/>.
+ **/
+
+
+#pragma once
+
+#include <vector>
+#include <list>
+
+
+namespace OrthancStone
+{
+  namespace Internals
+  {
+    /**
+     * This is an oriented 2D line with unsigned integer coordinates.
+     **/
+    class OrientedIntegerLine2D
+    {
+    private:
+      size_t  x1_;
+      size_t  y1_;
+      size_t  x2_;
+      size_t  y2_;
+
+    public:
+      OrientedIntegerLine2D(size_t x1,
+                            size_t y1,
+                            size_t x2,
+                            size_t y2);
+
+      bool IsVertical() const
+      {
+        return (x1_ == x2_);
+      }
+
+      bool IsHorizontal() const
+      {
+        return (y1_ == y2_);
+      }
+
+      bool IsEmpty() const
+      {
+        return (x1_ == x2_ &&
+                y1_ == y2_);
+      }
+
+      size_t GetX1() const
+      {
+        return x1_;
+      }
+
+      size_t GetY1() const
+      {
+        return y1_;
+      }
+
+      size_t GetX2() const
+      {
+        return x2_;
+      }
+
+      size_t GetY2() const
+      {
+        return y2_;
+      }
+
+      bool IsDownward() const
+      {
+        return (y1_ < y2_);
+      }
+
+      typedef std::list< std::pair<size_t, size_t> >  Chain;
+      
+      static void ExtractChains(std::list<Chain>& chains,
+                                const std::vector<OrientedIntegerLine2D>& edges);
+    };
+  }
+}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/OrthancStone/Sources/Toolbox/Internals/RectanglesIntegerProjection.cpp	Tue Jan 11 18:58:37 2022 +0100
@@ -0,0 +1,194 @@
+/**
+ * Stone of Orthanc
+ * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics
+ * Department, University Hospital of Liege, Belgium
+ * Copyright (C) 2017-2022 Osimis S.A., Belgium
+ * Copyright (C) 2021-2022 Sebastien Jodogne, ICTEAM UCLouvain, Belgium
+ *
+ * This program is free software: you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation, either version 3 of
+ * the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this program. If not, see
+ * <http://www.gnu.org/licenses/>.
+ **/
+
+
+#include "RectanglesIntegerProjection.h"
+
+#include <OrthancException.h>
+
+#include <algorithm>
+#include <cassert>
+
+
+namespace OrthancStone
+{
+  namespace Internals
+  {
+    class RectanglesIntegerProjection::Endpoint
+    {
+    private:
+      size_t  intervalIndex_;
+      double  value_;
+      bool    isLow_;
+
+    public:
+      Endpoint(size_t intervalIndex,
+               double value,
+               bool isLow) :
+        intervalIndex_(intervalIndex),
+        value_(value),
+        isLow_(isLow)
+      {
+      }
+
+      bool operator< (const Endpoint& other) const
+      {
+        if (value_ < other.value_)
+        {
+          return true;
+        }
+        else if (value_ == other.value_)
+        {
+          return isLow_;
+        }
+        else
+        {
+          return false;
+        }
+      }
+
+      size_t GetIntervalIndex() const
+      {
+        return intervalIndex_;
+      }
+
+      double GetValue() const
+      {
+        return value_;
+      }
+
+      bool IsLow() const
+      {
+        return isLow_;
+      }
+    };
+
+
+    RectanglesIntegerProjection::RectanglesIntegerProjection(const std::list<Extent2D>& rectangles,
+                                                             bool isHorizontal)
+    {
+      std::vector<Endpoint> endpoints;
+      endpoints.reserve(2 * rectangles.size());
+      
+      size_t count = 0;
+      for (std::list<Extent2D>::const_iterator it = rectangles.begin(); it != rectangles.end(); ++it)
+      {
+        if (!it->IsEmpty())
+        {
+          if (isHorizontal)
+          {
+            assert(it->GetX1() < it->GetX2());
+            endpoints.push_back(Endpoint(count, it->GetX1(), true));
+            endpoints.push_back(Endpoint(count, it->GetX2(), false));
+          }
+          else
+          {
+            assert(it->GetY1() < it->GetY2());
+            endpoints.push_back(Endpoint(count, it->GetY1(), true));
+            endpoints.push_back(Endpoint(count, it->GetY2(), false));
+          }
+
+          assert(endpoints[endpoints.size() - 2] < endpoints[endpoints.size() - 1]);
+          
+          count++;
+        }
+      }
+
+      std::sort(endpoints.begin(), endpoints.end());
+
+      intervalsLow_.resize(count);
+      intervalsHigh_.resize(count);
+      endpointsToDouble_.reserve(count);
+
+      for (size_t i = 0; i < endpoints.size(); i++)
+      {
+        if (endpointsToDouble_.empty() ||
+            endpointsToDouble_.back() < endpoints[i].GetValue())
+        {
+          endpointsToDouble_.push_back(endpoints[i].GetValue());
+        }
+
+        size_t intervalIndex = endpoints[i].GetIntervalIndex();
+        
+        if (endpoints[i].IsLow())
+        {
+          intervalsLow_[intervalIndex] = endpointsToDouble_.size() - 1;
+        }
+        else
+        {
+          intervalsHigh_[intervalIndex] = endpointsToDouble_.size() - 1;
+        }
+      }
+
+      for (size_t i = 0; i < intervalsLow_.size(); i++)
+      {
+        assert(intervalsLow_[i] < intervalsHigh_[i]);
+      }
+    }
+
+
+    double RectanglesIntegerProjection::GetEndpointCoordinate(size_t index) const
+    {
+      if (index >= endpointsToDouble_.size())
+      {
+        throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange);
+      }
+      else
+      {
+        return endpointsToDouble_[index];
+      }
+    }
+
+    
+    size_t RectanglesIntegerProjection::GetProjectedRectanglesCount() const
+    {
+      assert(intervalsLow_.size() == intervalsHigh_.size());
+      return intervalsLow_.size();
+    }
+
+    
+    size_t RectanglesIntegerProjection::GetProjectedRectangleLow(size_t index) const
+    {
+      if (index >= GetProjectedRectanglesCount())
+      {
+        throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange);
+      }
+      else
+      {
+        return intervalsLow_[index];
+      }
+    }
+
+    
+    size_t RectanglesIntegerProjection::GetProjectedRectangleHigh(size_t index) const
+    {
+      if (index >= GetProjectedRectanglesCount())
+      {
+        throw Orthanc::OrthancException(Orthanc::ErrorCode_ParameterOutOfRange);
+      }
+      else
+      {
+        return intervalsHigh_[index];
+      }
+    }
+  }
+}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/OrthancStone/Sources/Toolbox/Internals/RectanglesIntegerProjection.h	Tue Jan 11 18:58:37 2022 +0100
@@ -0,0 +1,63 @@
+/**
+ * Stone of Orthanc
+ * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics
+ * Department, University Hospital of Liege, Belgium
+ * Copyright (C) 2017-2022 Osimis S.A., Belgium
+ * Copyright (C) 2021-2022 Sebastien Jodogne, ICTEAM UCLouvain, Belgium
+ *
+ * This program is free software: you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation, either version 3 of
+ * the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this program. If not, see
+ * <http://www.gnu.org/licenses/>.
+ **/
+
+
+#pragma once
+
+#include "../Extent2D.h"
+
+#include <boost/noncopyable.hpp>
+#include <list>
+#include <vector>
+
+namespace OrthancStone
+{
+  namespace Internals
+  {
+    class RectanglesIntegerProjection : public boost::noncopyable
+    {
+    private:
+      class Endpoint;
+
+      std::vector<double>  endpointsToDouble_;
+      std::vector<size_t>  intervalsLow_;
+      std::vector<size_t>  intervalsHigh_;
+
+    public:
+      RectanglesIntegerProjection(const std::list<Extent2D>& rectangles,
+                                  bool isHorizontal);
+      
+      size_t GetEndpointsCount() const
+      {
+        return endpointsToDouble_.size();
+      }
+
+      double GetEndpointCoordinate(size_t index) const;
+      
+      size_t GetProjectedRectanglesCount() const;
+      
+      size_t GetProjectedRectangleLow(size_t index) const;
+      
+      size_t GetProjectedRectangleHigh(size_t index) const;
+    };
+  }
+}
--- a/UnitTestsSources/ComputationalGeometryTests.cpp	Tue Jan 11 15:36:04 2022 +0100
+++ b/UnitTestsSources/ComputationalGeometryTests.cpp	Tue Jan 11 18:58:37 2022 +0100
@@ -22,6 +22,8 @@
 
 #include <gtest/gtest.h>
 
+#include "../OrthancStone/Sources/Toolbox/Internals/OrientedIntegerLine2D.h"
+#include "../OrthancStone/Sources/Toolbox/Internals/RectanglesIntegerProjection.h"
 #include "../OrthancStone/Sources/Toolbox/SegmentTree.h"
 
 #include <Logging.h>
@@ -326,3 +328,154 @@
   root.VisitSegment(8, 9, minus);
   ASSERT_TRUE(CheckCounter(root, 0));
 }
+
+
+TEST(UnionOfRectangles, RectanglesIntegerProjection)
+{
+  std::list<OrthancStone::Extent2D> rectangles;
+  rectangles.push_back(OrthancStone::Extent2D(10, 20, 30, 40));
+
+  {
+    OrthancStone::Internals::RectanglesIntegerProjection h(rectangles, true);
+    ASSERT_EQ(2u, h.GetEndpointsCount());
+    ASSERT_EQ(10, h.GetEndpointCoordinate(0));
+    ASSERT_EQ(30, h.GetEndpointCoordinate(1));
+    ASSERT_EQ(1u, h.GetProjectedRectanglesCount());
+    ASSERT_EQ(0u, h.GetProjectedRectangleLow(0));
+    ASSERT_EQ(1u, h.GetProjectedRectangleHigh(0));
+
+    ASSERT_THROW(h.GetEndpointCoordinate(2), Orthanc::OrthancException);
+    ASSERT_THROW(h.GetProjectedRectangleLow(1), Orthanc::OrthancException);
+    ASSERT_THROW(h.GetProjectedRectangleHigh(1), Orthanc::OrthancException);
+  }
+
+  {
+    OrthancStone::Internals::RectanglesIntegerProjection h(rectangles, false);
+    ASSERT_EQ(2u, h.GetEndpointsCount());
+    ASSERT_EQ(20, h.GetEndpointCoordinate(0));
+    ASSERT_EQ(40, h.GetEndpointCoordinate(1));
+    ASSERT_EQ(1u, h.GetProjectedRectanglesCount());
+    ASSERT_EQ(0u, h.GetProjectedRectangleLow(0));
+    ASSERT_EQ(1u, h.GetProjectedRectangleHigh(0));
+  }
+
+  rectangles.push_back(OrthancStone::Extent2D(20, 30, 40, 50));
+
+  {
+    OrthancStone::Internals::RectanglesIntegerProjection h(rectangles, true);
+    ASSERT_EQ(4u, h.GetEndpointsCount());
+    ASSERT_EQ(10, h.GetEndpointCoordinate(0));
+    ASSERT_EQ(20, h.GetEndpointCoordinate(1));
+    ASSERT_EQ(30, h.GetEndpointCoordinate(2));
+    ASSERT_EQ(40, h.GetEndpointCoordinate(3));
+    ASSERT_EQ(2u, h.GetProjectedRectanglesCount());
+    ASSERT_EQ(0u, h.GetProjectedRectangleLow(0));
+    ASSERT_EQ(2u, h.GetProjectedRectangleHigh(0));
+    ASSERT_EQ(1u, h.GetProjectedRectangleLow(1));
+    ASSERT_EQ(3u, h.GetProjectedRectangleHigh(1));
+  }
+
+  {
+    OrthancStone::Internals::RectanglesIntegerProjection h(rectangles, false);
+    ASSERT_EQ(4u, h.GetEndpointsCount());
+    ASSERT_EQ(20, h.GetEndpointCoordinate(0));
+    ASSERT_EQ(30, h.GetEndpointCoordinate(1));
+    ASSERT_EQ(40, h.GetEndpointCoordinate(2));
+    ASSERT_EQ(50, h.GetEndpointCoordinate(3));
+    ASSERT_EQ(2u, h.GetProjectedRectanglesCount());
+    ASSERT_EQ(0u, h.GetProjectedRectangleLow(0));
+    ASSERT_EQ(2u, h.GetProjectedRectangleHigh(0));
+    ASSERT_EQ(1u, h.GetProjectedRectangleLow(1));
+    ASSERT_EQ(3u, h.GetProjectedRectangleHigh(1));
+  }
+}
+
+
+static void Convert(std::vector<size_t>& horizontal,
+                    std::vector<size_t>& vertical,
+                    const OrthancStone::Internals::OrientedIntegerLine2D::Chain& chain)
+{
+  horizontal.clear();
+  vertical.clear();
+
+  for (OrthancStone::Internals::OrientedIntegerLine2D::Chain::const_iterator
+         it = chain.begin(); it != chain.end(); ++it)
+  {
+    horizontal.push_back(it->first);
+    vertical.push_back(it->second);
+  }
+}
+
+
+TEST(UnionOfRectangles, ExtractChains)
+{
+  std::vector<OrthancStone::Internals::OrientedIntegerLine2D> edges;
+  edges.push_back(OrthancStone::Internals::OrientedIntegerLine2D(0, 0, 10, 0));
+  edges.push_back(OrthancStone::Internals::OrientedIntegerLine2D(10, 0, 10, 20));
+  edges.push_back(OrthancStone::Internals::OrientedIntegerLine2D(10, 20, 0, 20));
+
+  std::list<OrthancStone::Internals::OrientedIntegerLine2D::Chain> chains;
+  OrthancStone::Internals::OrientedIntegerLine2D::ExtractChains(chains, edges);
+
+  std::vector<size_t> h, v;
+
+  ASSERT_EQ(1u, chains.size());
+
+  Convert(h, v, chains.front());
+  ASSERT_EQ(4u, h.size());
+  ASSERT_EQ(0u, h[0]);
+  ASSERT_EQ(10u, h[1]);
+  ASSERT_EQ(10u, h[2]);
+  ASSERT_EQ(0u, h[3]);
+  ASSERT_EQ(4u, v.size());
+  ASSERT_EQ(0u, v[0]);
+  ASSERT_EQ(0u, v[1]);
+  ASSERT_EQ(20u, v[2]);
+  ASSERT_EQ(20u, v[3]);
+
+  edges.push_back(OrthancStone::Internals::OrientedIntegerLine2D(5, 5, 10, 5));
+  OrthancStone::Internals::OrientedIntegerLine2D::ExtractChains(chains, edges);
+
+  ASSERT_EQ(2u, chains.size());
+  
+  Convert(h, v, chains.front());
+  ASSERT_EQ(4u, h.size());
+  ASSERT_EQ(0u, h[0]);
+  ASSERT_EQ(10u, h[1]);
+  ASSERT_EQ(10u, h[2]);
+  ASSERT_EQ(0u, h[3]);
+  ASSERT_EQ(4u, v.size());
+  ASSERT_EQ(0u, v[0]);
+  ASSERT_EQ(0u, v[1]);
+  ASSERT_EQ(20u, v[2]);
+  ASSERT_EQ(20u, v[3]);
+  
+  Convert(h, v, chains.back());
+  ASSERT_EQ(2u, h.size());
+  ASSERT_EQ(5u, h[0]);
+  ASSERT_EQ(10u, h[1]);
+  ASSERT_EQ(2u, v.size());
+  ASSERT_EQ(5u, v[0]);
+  ASSERT_EQ(5u, v[1]);
+
+  edges.push_back(OrthancStone::Internals::OrientedIntegerLine2D(0, 20, 5, 5));
+  OrthancStone::Internals::OrientedIntegerLine2D::ExtractChains(chains, edges);
+
+  ASSERT_EQ(1u, chains.size());
+  
+  Convert(h, v, chains.front());
+  ASSERT_EQ(6u, h.size());
+  ASSERT_EQ(0u, h[0]);
+  ASSERT_EQ(10u, h[1]);
+  ASSERT_EQ(10u, h[2]);
+  ASSERT_EQ(0u, h[3]);
+  ASSERT_EQ(5u, h[4]);
+  ASSERT_EQ(10u, h[5]);
+  ASSERT_EQ(6u, v.size());
+  ASSERT_EQ(0u, v[0]);
+  ASSERT_EQ(0u, v[1]);
+  ASSERT_EQ(20u, v[2]);
+  ASSERT_EQ(20u, v[3]);
+  ASSERT_EQ(5u, v[4]);
+  ASSERT_EQ(5u, v[5]);
+}