Changed implementation of WeightedRips to store simplex values (max distance between simplices' vertices) as an invisible layer on top of each simplex object, so that the data() field of WeightedRips has been freed for use by the users again.
#include <algorithm>
#include <utility>
#include <boost/utility.hpp>
#include <iostream>
#include <utilities/log.h>
#include <utilities/counter.h>
#include <utilities/indirect.h>
#include <boost/iterator/counting_iterator.hpp>
#include <functional>
#ifdef LOGGING
static rlog::RLogChannel* rlRips = DEF_CHANNEL("rips/info", rlog::Log_Debug);
static rlog::RLogChannel* rlRipsDebug = DEF_CHANNEL("rips/debug", rlog::Log_Debug);
#endif // LOGGING
#ifdef COUNTERS
static Counter* cClique = GetCounter("rips/clique");
#endif // COUNTERS
template<class D, class S>
template<class Functor, class Iterator>
void
Rips<D,S>::
generate(Dimension k, DistanceType max, const Functor& f, Iterator bg, Iterator end) const
{
rLog(rlRipsDebug, "Entered generate with %d indices", distances().size());
WithinDistance neighbor(distances(), max);
// current = empty
// candidates = everything
VertexContainer current;
VertexContainer candidates(bg, end);
bron_kerbosch(current, candidates, boost::prior(candidates.begin()), k, neighbor, f);
}
template<class D, class S>
template<class Functor, class Iterator>
void
Rips<D,S>::
vertex_cofaces(IndexType v, Dimension k, DistanceType max, const Functor& f, Iterator bg, Iterator end) const
{
WithinDistance neighbor(distances(), max);
// current = [v]
// candidates = everything - [v]
VertexContainer current; current.push_back(v);
VertexContainer candidates;
for (Iterator cur = bg; cur != end; ++cur)
if (*cur != v && neighbor(v, *cur))
candidates.push_back(*cur);
bron_kerbosch(current, candidates, boost::prior(candidates.begin()), k, neighbor, f);
}
template<class D, class S>
template<class Functor, class Iterator>
void
Rips<D,S>::
edge_cofaces(IndexType u, IndexType v, Dimension k, DistanceType max, const Functor& f, Iterator bg, Iterator end) const
{
rLog(rlRipsDebug, "In edge_cofaces(%d, %d)", u, v);
WithinDistance neighbor(distances(), max);
AssertMsg(neighbor(u,v), "The new edge must be in the complex");
// current = [u,v]
// candidates = everything - [u,v]
VertexContainer current; current.push_back(u); current.push_back(v);
VertexContainer candidates;
for (Iterator cur = bg; cur != end; ++cur)
if (*cur != u && *cur != v && neighbor(v,*cur) && neighbor(u,*cur))
{
candidates.push_back(*cur);
rLog(rlRipsDebug, " added candidate: %d", *cur);
}
bron_kerbosch(current, candidates, boost::prior(candidates.begin()), k, neighbor, f);
}
template<class D, class S>
template<class Functor, class Iterator>
void
Rips<D,S>::
cofaces(const Simplex& s, Dimension k, DistanceType max, const Functor& f, Iterator bg, Iterator end) const
{
rLog(rlRipsDebug, "In cofaces(%s)", tostring(s).c_str());
WithinDistance neighbor(distances(), max);
// current = s.vertices()
VertexContainer current(s.vertices().begin(), s.vertices().end());
// candidates = everything - s.vertices() that is a neighbor() of every vertex in the simplex
VertexContainer candidates;
typedef difference_iterator<Iterator,
typename VertexContainer::const_iterator,
std::less<Vertex> > DifferenceIterator;
for (DifferenceIterator cur = DifferenceIterator(bg, end, s.vertices().begin(), s.vertices().end());
cur != DifferenceIterator(end, end, s.vertices().end(), s.vertices().end());
++cur)
{
bool nghbr = true;
for (typename VertexContainer::const_iterator v = s.vertices().begin(); v != s.vertices().end(); ++v)
if (!neighbor(*v, *cur)) { nghbr = false; break; }
if (nghbr)
{
candidates.push_back(*cur);
rLog(rlRipsDebug, " added candidate: %d", *cur);
}
}
bron_kerbosch(current, candidates, boost::prior(candidates.begin()), k, neighbor, f, false);
}
template<class D, class S>
template<class Functor, class NeighborTest>
void
Rips<D,S>::
bron_kerbosch(VertexContainer& current,
const VertexContainer& candidates,
typename VertexContainer::const_iterator excluded,
Dimension max_dim,
const NeighborTest& neighbor,
const Functor& functor,
bool check_initial) const
{
rLog(rlRipsDebug, "Entered bron_kerbosch");
if (check_initial && !current.empty())
{
Simplex s(current);
rLog(rlRipsDebug, "Reporting simplex: %s", tostring(s).c_str());
functor(s);
}
if (current.size() == max_dim + 1)
return;
rLog(rlRipsDebug, "Traversing %d vertices", candidates.end() - boost::next(excluded));
for (typename VertexContainer::const_iterator cur = boost::next(excluded); cur != candidates.end(); ++cur)
{
current.push_back(*cur);
rLog(rlRipsDebug, " current.size() = %d, current.back() = %d", current.size(), current.back());
VertexContainer new_candidates;
for (typename VertexContainer::const_iterator ccur = candidates.begin(); ccur != cur; ++ccur)
if (neighbor(*ccur, *cur))
new_candidates.push_back(*ccur);
size_t ex = new_candidates.size();
for (typename VertexContainer::const_iterator ccur = boost::next(cur); ccur != candidates.end(); ++ccur)
if (neighbor(*ccur, *cur))
new_candidates.push_back(*ccur);
excluded = new_candidates.begin() + (ex - 1);
bron_kerbosch(current, new_candidates, excluded, max_dim, neighbor, functor);
current.pop_back();
}
}
template<class Distances_, class Simplex_>
typename Rips<Distances_, Simplex_>::DistanceType
Rips<Distances_, Simplex_>::
distance(const Simplex& s1, const Simplex& s2) const
{
DistanceType mx = 0;
for (typename Simplex::VertexContainer::const_iterator a = s1.vertices().begin(); a != s1.vertices().end(); ++a)
for (typename Simplex::VertexContainer::const_iterator b = s2.vertices().begin(); b != s2.vertices().end(); ++b)
mx = std::max(mx, distances_(*a,*b));
return mx;
}
template<class Distances_, class Simplex_>
typename Rips<Distances_, Simplex_>::DistanceType
Rips<Distances_, Simplex_>::
max_distance() const
{
DistanceType mx = 0;
for (IndexType a = distances_.begin(); a != distances_.end(); ++a)
for (IndexType b = boost::next(a); b != distances_.end(); ++b)
mx = std::max(mx, distances_(a,b));
return mx;
}
template<class Distances_, class Simplex_>
typename Rips<Distances_, Simplex_>::DistanceType
Rips<Distances_, Simplex_>::Evaluator::
operator()(const Simplex& s) const
{
DistanceType mx = 0;
for (typename Simplex::VertexContainer::const_iterator a = s.vertices().begin(); a != s.vertices().end(); ++a)
for (typename Simplex::VertexContainer::const_iterator b = boost::next(a); b != s.vertices().end(); ++b)
mx = std::max(mx, distances_(*a,*b));
return mx;
}