diff --git a/Surface_mesh_approximation/include/CGAL/vsa_approximation.h b/Surface_mesh_approximation/include/CGAL/vsa_approximation.h
index 79b44d779fc..89e185de3ab 100644
--- a/Surface_mesh_approximation/include/CGAL/vsa_approximation.h
+++ b/Surface_mesh_approximation/include/CGAL/vsa_approximation.h
@@ -850,6 +850,22 @@ public:
       facet_proxy_map[f] = fproxy_map[f];
   }
 
+  /*!
+   * @brief Get the facet region of the specified proxy.
+   * @tparam OutputIterator output iterator with `boost::graph_traits<TriangleMesh>::%face_descriptor` as value type
+   * @param px_idx proxy index
+   * @param out_itr output iterator
+   */
+  template <typename OutputIterator>
+  void get_proxy_region(const std::size_t px_idx, OutputIterator out_itr) const {
+    if (px_idx >= proxies.size())
+      return;
+
+    BOOST_FOREACH(face_descriptor f, faces(*m_pmesh))
+      if (fproxy_map[f] == px_idx)
+        *out_itr++ = f;
+  }
+
   /*!
    * @brief Get the proxies.
    * @tparam OutputIterator output iterator with Proxy as value type
diff --git a/Surface_mesh_approximation/test/Surface_mesh_approximation/vsa_class_interface_test.cpp b/Surface_mesh_approximation/test/Surface_mesh_approximation/vsa_class_interface_test.cpp
index d071de7e178..48a856584d0 100644
--- a/Surface_mesh_approximation/test/Surface_mesh_approximation/vsa_class_interface_test.cpp
+++ b/Surface_mesh_approximation/test/Surface_mesh_approximation/vsa_class_interface_test.cpp
@@ -99,6 +99,11 @@ int main()
   std::cout << "get outputs" << std::endl;
   approx.get_proxy_map(proxy_pmap);
 
+  for (std::size_t i = 0; i < approx.get_proxies_size(); ++i) {
+    std::list<Facet_handle> patch;
+    approx.get_proxy_region(i, std::back_inserter(patch));
+  }
+
   std::vector<Plane_proxy> proxies;
   approx.get_proxies(std::back_inserter(proxies));