size_t count = 0;
// As long as we did not either end, look for common elements
- while (first1!=last1 && first2!=last2)
+ while (first1 != last1 && first2 != last2)
{
if (*first1 < *first2)
++first1;
--- /dev/null
+#ifndef image_h_uiHiKahYooroomai
+#define image_h_uiHiKahYooroomai
+
+#include "simplicialcomplex.h"
+
+
+namespace libstick {
+
+/** We consider an image, i.e. a rectangular arrangement of pixels, of given
+ * height and width. Every pixel has a specific value (e.g., grey value) given
+ * by the 'values' array as concatenation of all rows. We add every pixel as
+ * vertex with its value as simplex value. We further triangulate the resulting
+ * rectangular latice as a grid graph together with south-west to north-east
+ * edges for each of the small squares of the grid graph. Finally the resulting
+ * triangles are added. For all non-vertex simplices the maximum of its
+ * vertex-values is used as its value. The first pixel is considered to be at
+ * bottom-left. */
+template<class VT, class IT>
+void add_image_to_complex(simplicial_complex<2, IT, VT> &s, const VT *values, size_t width, size_t height) {
+ // Add the vertices
+ const IT v1st = s.size();
+ for (unsigned i=0; i < height; ++i)
+ for (unsigned j=0; j < width; ++j)
+ s.add_simplex(0, NULL, values[i*width + j]);
+
+ // row i+1 o-------o
+ // | /|
+ // | / |
+ // | 2/ |
+ // 0 | / |
+ // | / |
+ // | / |
+ // |/ 1 |
+ // row i o-------o
+ // ^
+ // column j
+
+ IT faces[3];
+
+ // Adding horizontal edges
+ const IT he1st = s.size();
+ for (unsigned i=0; i < height; ++i)
+ for (unsigned j=0; j < width-1; ++j) {
+ faces[0] = v1st + i*width + j;
+ faces[1] = faces[0] + 1;
+ s.add_simplex(1, faces);
+ }
+
+ // Adding vertical edges
+ const IT ve1st = s.size();
+ for (unsigned i=0; i < height-1; ++i)
+ for (unsigned j=0; j < width; ++j) {
+ faces[0] = v1st + i*width + j;
+ faces[1] = faces[0] + width;
+ s.add_simplex(1, faces);
+ }
+
+ // Adding diagonal edges
+ const IT de1st = s.size();
+ for (unsigned i=0; i < height-1; ++i)
+ for (unsigned j=0; j < width-1; ++j) {
+ faces[0] = v1st + i*width + j;
+ faces[1] = faces[0] + width + 1;
+ s.add_simplex(1, faces);
+ }
+
+ // Add triangles
+ for (unsigned i=0; i < height-1; ++i)
+ for (unsigned j=0; j < width-1; ++j) {
+ faces[0] = he1st + i*(width-1) + j;
+ faces[1] = ve1st + i*width + j + 1;
+ faces[2] = de1st + i*(width-1) + j;
+ s.add_simplex(2, faces);
+
+ faces[0] = he1st + (i+1)*(width-1) + j;
+ faces[1] = ve1st + i*width + j;
+ faces[2] = de1st + i*(width-1) + j;
+ s.add_simplex(2, faces);
+ }
+}
+
+
+}
+
+#endif
+
}
}
#ifndef NDEBUG
- // Shakespeare principle: A simplex either gaves birth or kills, be
+ // Shakespeare principle: A simplex either gives birth or kills, be
// or not be, tertium non datur.
for (unsigned c=0; c < lowestones.size(); ++c) {
if (lowestones.get_column(c).size() != 0)
* Each 0-dimensional simplex automatically has this simplex as its face.
* Consequently, the innner class simplex_order gives the extended boundary
* matrix. */
-template<int MAXDIM, class IT=uint32_t, class VT=double>
+template<int MAXDIM, class IT, class VT>
class simplicial_complex {
public:
}
/** Get i-th simplex in the simplex order. */
- const Simplex& get_simplex(size_t i) const {
+ const Simplex& get_simplex(index_type i) const {
assert(i < size());
return c.simplices[order.at(i)];
}
}
/** Add a simplex to the complex. The dimension of the faces must be
- * dim-1, and they must already be part of the complex. */
- void add_simplex(int dim, index_type* faces, value_type value) {
- add_simplex(Simplex::create(dim, faces, value));
+ * dim-1, and they must already be part of the complex. Returns the
+ * index of the added simplex. */
+ index_type add_simplex(int dim, index_type* faces, value_type value) {
+ return add_simplex(Simplex::create(dim, faces, value));
+ }
+
+ /** Add a simplex to the complex of at least dimension 1. The dimension
+ * of the faces must be dim-1, and they must already be part of the
+ * complex. The value of the simplex is set to the maximum value of its
+ * faces. Returns the index of the added simplex. */
+ index_type add_simplex(int dim, index_type* faces) {
+ assert(dim >= 1);
+
+ // Get max value of its faces
+ VT value = simplices[faces[0]].value;
+ for (size_t i=0; i < Simplex::face_count_bydim(dim); ++i)
+ value = std::max(value, simplices[faces[i]].value);
+
+ return add_simplex(dim, faces, value);
}
/** Add a simplex to the complex. The dimension of the faces must be
- * dim-1, and they must already be part of the complex. */
- void add_simplex(Simplex s) {
+ * dim-1, and they must already be part of the complex. Returns the
+ * index of the added simplex. */
+ index_type add_simplex(Simplex s) {
// Check requirements for faces
for (unsigned i=0; i < s.face_count(); ++i) {
// Faces are already in complex.
assert(simplices[s.faces[i]].dim == s.dim-1);
}
+ // index_type must be large enough
+ assert(simplices.size() < std::numeric_limits<IT>::max());
+
+ index_type idx = simplices.size();
simplices.push_back(s);
+ return idx;
}
/** Add an array of simplices */
add_executable(tests ${tests_SRC})
target_link_libraries(tests
#stick
+ m
${libcpptest_LIBRARIES})
--- /dev/null
+#ifndef image_h_PooBeibuceushabo
+#define image_h_PooBeibuceushabo
+
+#include <cpptest.h>
+#include <cpptest-suite.h>
+
+#include <libstick-0.1/image.h>
+
+#include <cassert>
+#include <math.h>
+
+
+using namespace libstick;
+
+
+class image_TestSuite: public Test::Suite {
+
+ private:
+
+ public:
+ image_TestSuite()
+ {
+ TEST_ADD(image_TestSuite::test_incidences);
+ }
+
+ protected:
+ virtual void setup() {
+ }
+
+ virtual void tear_down() {
+ }
+
+ void test_incidences() {
+ const size_t height=128, width=64;
+
+ //Create "some" image
+ float image[height*width];
+ for (unsigned i=0; i < height; ++i)
+ for (unsigned j=0; j < width; ++j)
+ image[i*width + j] = (sin(j/10.0) + 0.5*sin(j/4.0))*cos(i/7.0) + i*j*0.5*1e-3;
+
+ // Create complex and add image
+ typedef simplicial_complex<2, uint32_t, float> scomplex;
+ scomplex s;
+ add_image_to_complex(s, image, width, height);
+ assert(s.is_monotone());
+
+ typedef scomplex::simplex_order sorder;
+ sorder so(s);
+ assert(so.is_filtration());
+
+ typedef sorder::boundary_matrix bmatrix;
+ bmatrix bm = so.get_boundary_matrix();
+
+ // Check for the right vertex incidences
+ for (unsigned i=0; i < height; ++i) {
+ for (unsigned j=0; j < width; ++j) {
+ // Get matrix row of this vertex. Its size is the number of edges
+ // incident to this vertex
+ const bmatrix::row_type &row = bm.get_row(1 + i*width + j);
+
+ // Somewhere at the boundary of the image
+ if (i == 0 || i == height-1 || j == 0 || j == width-1) {
+ // top-left and bottom-right corner
+ if ((i == 0 && j == 0) || (i == height - 1 && j == width - 1)) {
+ TEST_ASSERT(row.size() == 3);
+ } else if ((i == 0 && j == width - 1) || (i == height - 1 && j == 0)) {
+ TEST_ASSERT(row.size() == 2);
+ } else {
+ TEST_ASSERT(row.size() == 4);
+ }
+ } else {
+ TEST_ASSERT(row.size() == 6);
+ }
+ }
+ }
+ }
+};
+
+#endif
#include "booleanmatrix.h"
#include "simplicialcomplex.h"
#include "persistence.h"
+#include "image.h"
using namespace std;
ts.add(auto_ptr<Test::Suite>(new boolean_matrix_TestSuite));
ts.add(auto_ptr<Test::Suite>(new simplicial_complex_TestSuite));
ts.add(auto_ptr<Test::Suite>(new persistence_TestSuite));
+ ts.add(auto_ptr<Test::Suite>(new image_TestSuite));
Test::TextOutput output(Test::TextOutput::Verbose);
if (ts.run(output))
projects / libstick.git / commitdiff