--- a/bindings/python/CMakeLists.txt Fri May 01 15:09:38 2009 -0700
+++ b/bindings/python/CMakeLists.txt Sat Apr 11 10:29:53 2009 -0700
@@ -12,6 +12,16 @@
static-persistence.cpp
simplex.cpp
zigzag-persistence.cpp
+ rips.cpp
)
+target_link_libraries (_dionysus ${libraries})
+
-target_link_libraries (_dionysus ${libraries})
+add_custom_target (dionysus ALL
+ ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/dionysus ${CMAKE_CURRENT_BINARY_DIR}/dionysus
+ DEPENDS dionysus/__init__.py)
+
+get_target_property (_dionysus_location _dionysus LOCATION)
+add_custom_target (dionysus-link ALL
+ ${CMAKE_COMMAND} -E create_symlink ${_dionysus_location} ${CMAKE_CURRENT_BINARY_DIR}/dionysus/_dionysus.so
+ DEPENDS _dionysus)
--- a/bindings/python/chain.cpp Fri May 01 15:09:38 2009 -0700
+++ b/bindings/python/chain.cpp Sat Apr 11 10:29:53 2009 -0700
@@ -4,7 +4,7 @@
using namespace boost::python;
-typedef SPersistence::OrderDescriptor::Chains::Chain VChain;
+typedef SPersistence::Chain VChain;
VChain::const_iterator begin(const VChain& c) { return c.begin(); }
--- a/bindings/python/dionysus.cpp Fri May 01 15:09:38 2009 -0700
+++ b/bindings/python/dionysus.cpp Sat Apr 11 10:29:53 2009 -0700
@@ -9,6 +9,8 @@
void export_chain();
void export_zigzag_persistence();
+void export_rips();
+
#ifdef LOGGING
void enable_log(std::string s)
{
@@ -25,6 +27,8 @@
export_zigzag_persistence();
+ export_rips();
+
#ifdef LOGGING
bp::def("enable_log", &enable_log);
#endif
--- a/bindings/python/dionysus/__init__.py Fri May 01 15:09:38 2009 -0700
+++ b/bindings/python/dionysus/__init__.py Sat Apr 11 10:29:53 2009 -0700
@@ -1,4 +1,5 @@
-from _dionysus import *
+from _dionysus import *
+from distances import *
#def init_with_data(s,v, d = None):
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/bindings/python/dionysus/distances.py Sat Apr 11 10:29:53 2009 -0700
@@ -0,0 +1,41 @@
+from math import sqrt
+
+def l2(p):
+ return sqrt( # take the square root
+ reduce(lambda x, y: x + y, # add them all up
+ map(lambda x: x**2, p))) # square each coordinate
+
+# Pairwise distances between the elements of `points` with respect to some `norm`
+class PairwiseDistances:
+ def __init__(self, points, norm = l2):
+ self.points = points
+ self.norm = norm
+
+ def __len__(self):
+ return len(self.points)
+
+ def __call__(self, p1, p2):
+ return self.norm([x - y for (x,y) in zip(self.points[p1], self.points[p2])])
+
+# Caches all distances specified by `distances`
+class ExplicitDistances:
+ def __init__(self, distances):
+ self.len = len(distances)
+ self.distances = []
+ for i in xrange(self.len):
+ self.distances.append([])
+ for j in xrange(self.len):
+ self.distances[-1].append(distances(i,j))
+
+ def __len__(self):
+ return self.len
+
+ def __call__(self, p1, p2):
+ return self.distances[p1][p2]
+
+# Generator of all points in a file `filename` with one point per line
+def points_file(filename):
+ fd = open(filename)
+ for line in fd.xreadlines():
+ yield map(float, line.strip().split())
+ fd.close()
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/bindings/python/python-rips.h Sat Apr 11 10:29:53 2009 -0700
@@ -0,0 +1,108 @@
+#ifndef __PYTHON_RIPS_H__
+#define __PYTHON_RIPS_H__
+
+#include <topology/rips.h>
+#include <utilities/indirect.h>
+
+#include "python-simplex.h"
+
+#include <boost/python.hpp>
+#include <boost/python/stl_iterator.hpp>
+
+
+namespace bp = boost::python;
+
+
+
+// This strange wrapper is necessary because Rips<...> stores only a const reference to the distances.
+// Something on the C++ side of things must store the actual DistancesWrapper object, and that's the
+// purpose of this class.
+class RipsWithDistances
+{
+ public:
+ class DistancesWrapper
+ {
+ public:
+ typedef unsigned IndexType;
+ typedef double DistanceType;
+
+ DistancesWrapper(bp::object distances):
+ distances_(distances) {}
+
+ DistanceType operator()(IndexType a, IndexType b) const { return bp::extract<DistanceType>(distances_(a, b)); }
+
+ IndexType size() const { return bp::len(distances_); }
+ IndexType begin() const { return 0; }
+ IndexType end() const { return size(); }
+
+ private:
+ bp::object distances_;
+ };
+
+ typedef DistancesWrapper::IndexType IndexType;
+ typedef DistancesWrapper::DistanceType DistanceType;
+
+ typedef Rips<DistancesWrapper, SimplexVD> RipsDS;
+ typedef RipsDS::Comparison Comparison;
+ typedef RipsDS::Evaluator Evaluator;
+
+ class FunctorWrapper
+ {
+ public:
+ FunctorWrapper(bp::object functor):
+ functor_(functor) {}
+
+ void operator()(const RipsDS::Simplex& s) const { functor_(s); }
+
+ private:
+ bp::object functor_;
+ };
+
+
+ RipsWithDistances(bp::object distances):
+ distances_(distances), rips_(distances_),
+ cmp_(Comparison(distances_)), eval_(distances_) {}
+
+ void generate(Dimension k, DistanceType max, bp::object functor) const
+ { rips_.generate(k, max, FunctorWrapper(functor)); }
+
+ void vertex_cofaces(IndexType v, Dimension k, DistanceType max, bp::object functor) const
+ { rips_.vertex_cofaces(v, k, max, FunctorWrapper(functor)); }
+
+ void edge_cofaces(IndexType u, IndexType v, Dimension k, DistanceType max, bp::object functor) const
+ { rips_.edge_cofaces(u, v, k, max, FunctorWrapper(functor)); }
+
+ void generate_candidates(Dimension k, DistanceType max, bp::object functor, bp::object seq) const
+ {
+ rips_.generate(k, max, FunctorWrapper(functor),
+ bp::stl_input_iterator<IndexType>(seq), bp::stl_input_iterator<IndexType>());
+ }
+
+ void vertex_cofaces_candidate(IndexType v, Dimension k, DistanceType max,
+ bp::object functor, bp::object seq) const
+ {
+ rips_.vertex_cofaces(v, k, max, FunctorWrapper(functor),
+ bp::stl_input_iterator<IndexType>(seq), bp::stl_input_iterator<IndexType>());
+ }
+
+ void edge_cofaces_candidates(IndexType u, IndexType v, Dimension k, DistanceType max,
+ bp::object functor, bp::object seq) const
+ {
+ rips_.edge_cofaces(u, v, k, max, FunctorWrapper(functor),
+ bp::stl_input_iterator<IndexType>(seq), bp::stl_input_iterator<IndexType>());
+ }
+
+ int cmp(const SimplexVD& s1, const SimplexVD& s2) const
+ { return cmp_.compare(s1, s2); }
+
+ DistanceType eval(const SimplexVD& s) const
+ { return eval_(s); }
+
+ private:
+ DistancesWrapper distances_;
+ RipsDS rips_;
+ ThreeOutcomeCompare<Comparison> cmp_; // in Python, cmp is a three outcome comparison
+ Evaluator eval_;
+};
+
+#endif
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/bindings/python/rips.cpp Sat Apr 11 10:29:53 2009 -0700
@@ -0,0 +1,31 @@
+#include <topology/rips.h>
+#include <boost/python.hpp>
+#include "python-rips.h" // defines RipsWithDistances
+
+#include <iostream>
+
+namespace bp = boost::python;
+
+/* Various wrappers for exposing Rips to Python */
+// Constructor from distances
+boost::shared_ptr<RipsWithDistances> init_from_distances(bp::object distances)
+{
+ boost::shared_ptr<RipsWithDistances> p(new RipsWithDistances(distances));
+ return p;
+}
+
+void export_rips()
+{
+ bp::class_<RipsWithDistances>("Rips", bp::no_init)
+ .def("__init__", make_constructor(&init_from_distances))
+ .def("generate", &RipsWithDistances::generate)
+ .def("generate", &RipsWithDistances::generate_candidates)
+ .def("vertex_cofaces", &RipsWithDistances::vertex_cofaces)
+ .def("vertex_cofaces", &RipsWithDistances::vertex_cofaces_candidate)
+ .def("edge_cofaces", &RipsWithDistances::edge_cofaces)
+ .def("edge_cofaces", &RipsWithDistances::edge_cofaces_candidates)
+
+ .def("cmp", &RipsWithDistances::cmp)
+ .def("eval", &RipsWithDistances::eval)
+ ;
+}
--- a/examples/rips/CMakeLists.txt Fri May 01 15:09:38 2009 -0700
+++ b/examples/rips/CMakeLists.txt Sat Apr 11 10:29:53 2009 -0700
@@ -8,3 +8,8 @@
add_executable (${t} ${t}.cpp)
target_link_libraries (${t} ${libraries} ${Boost_PROGRAM_OPTIONS_LIBRARY} ${Boost_SERIALIZATION_LIBRARY})
endforeach (t ${targets})
+
+add_custom_target (rips.py ALL
+ ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/rips.py ${CMAKE_CURRENT_BINARY_DIR}/rips.py)
+add_custom_target (rips-pairwise.py ALL
+ ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/rips-pairwise.py ${CMAKE_CURRENT_BINARY_DIR}/rips-pairwise.py)
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/examples/rips/rips-pairwise.py Sat Apr 11 10:29:53 2009 -0700
@@ -0,0 +1,40 @@
+#/usr/bin/env python
+
+
+from dionysus import Rips, PairwiseDistances, StaticPersistence, Filtration, points_file, \
+ ExplicitDistances, data_dim_cmp
+from sys import argv, exit
+import time
+
+
+if len(argv) < 4:
+ print "Usage: %s POINTS SKELETON MAX" % argv[0]
+ exit()
+
+filename = argv[1]
+skeleton = int(argv[2])
+max = float(argv[3])
+
+points = [p for p in points_file(filename)]
+distances = PairwiseDistances(points)
+#distances = ExplicitDistances(distances) # speeds up generation of the Rips complex at the expense of memory usage
+rips = Rips(distances)
+print time.asctime(), "Rips initialized"
+
+simplices = []
+rips.generate(skeleton, max, simplices.append)
+print time.asctime(), "Generated complex: %d simplices" % len(simplices)
+
+# While this step is unnecessary (Filtration below can be passed rips.cmp),
+# it greatly speeds up the running times
+for s in simplices: s.data = rips.eval(s)
+print time.asctime(), simplices[0], '...', simplices[-1]
+
+f = Filtration(simplices, data_dim_cmp) # could be rips.cmp if s.data for s in simplices is not set
+print time.asctime(), "Set up filtration"
+
+p = StaticPersistence(f)
+print time.asctime(), "Initialized StaticPersistence"
+
+p.pair_simplices()
+print time.asctime(), "Simplices paired"
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/examples/rips/rips.py Sat Apr 11 10:29:53 2009 -0700
@@ -0,0 +1,40 @@
+from math import fabs
+from dionysus import Rips, Filtration, StaticPersistence #, enable_log
+
+# Simple minded pairwise distance functor distance
+class Distances:
+ def __len__(self):
+ return 5
+
+ def __call__(self, x, y):
+ return fabs(y-x)
+
+dist = Distances()
+r = Rips(dist)
+lst = []
+lst2 = []
+
+#enable_log('rips')
+
+r.generate(1, 3, lst.append)
+r.generate(1, 3, lst2.append, [0,2,4])
+
+print "Rips complex on all vertices:", lst
+print "Rips complex on vertices [0,2,4]):", lst2
+
+print "Values:", [map(r.eval, lst)]
+print "Sorted:", sorted(lst, r.cmp)
+
+cofaces = []
+r.vertex_cofaces(2, 1, 3, cofaces.append)
+print "Cofaces of vertex 2:", cofaces
+
+cofaces = []
+r.vertex_cofaces(2, 1, 3, cofaces.append, [0,2,4])
+print "Cofaces of vertex 2 on vertices [0,2,4]:", cofaces
+
+f = Filtration(lst, r.cmp)
+p = StaticPersistence(f)
+p.pair_simplices()
+for s in p:
+ print s, s.sign()