" << std::endl;
#ifdef USE_COEFFICIENTS
std::cerr << " --modulus compute homology with coefficients in the prime field Z/
Z"
<< std::endl;
#endif
exit(exit_code);
}
int main(int argc, char** argv) {
if (argc < 2) print_usage_and_exit(-1);
const char* filename = nullptr;
index_t dim_max = 1;
value_t threshold = std::numeric_limits::max();
#ifdef USE_COEFFICIENTS
coefficient_t modulus = 2;
#else
const coefficient_t modulus = 2;
#endif
for (index_t i = 1; i < argc; ++i) {
const std::string arg(argv[i]);
if (arg == "--help") {
print_usage_and_exit(0);
} else if (arg == "--dim") {
std::string parameter = std::string(argv[++i]);
size_t next_pos;
dim_max = std::stol(parameter, &next_pos);
if (next_pos != parameter.size()) print_usage_and_exit(-1);
} else if (arg == "--threshold") {
std::string parameter = std::string(argv[++i]);
size_t next_pos;
threshold = std::stof(parameter, &next_pos);
if (next_pos != parameter.size()) print_usage_and_exit(-1);
#ifdef USE_COEFFICIENTS
} else if (arg == "--modulus") {
std::string parameter = std::string(argv[++i]);
size_t next_pos;
modulus = std::stol(parameter, &next_pos);
if (next_pos != parameter.size() || !is_prime(modulus)) print_usage_and_exit(-1);
#endif
} else {
if (filename) { print_usage_and_exit(-1); }
filename = argv[i];
}
}
std::ifstream input_stream(filename);
if (input_stream.fail()) {
std::cerr << "couldn't open file " << filename << std::endl;
exit(-1);
}
#ifdef FILE_FORMAT_POINT_CLOUD
std::vector> points;
std::string line;
value_t value;
while (std::getline(input_stream, line)) {
std::vector point;
std::istringstream s(line);
while (s >> value) point.push_back(value);
if (!point.empty()) points.push_back(point);
}
euclidean_distance_matrix eucl_dist(std::move(points));
index_t n = eucl_dist.size();
std::cout << "point cloud with " << n << " points" << std::endl;
std::vector distances;
for (int i = 0; i < n; ++i)
for (int j = 0; j < i; ++j)
if (i > j) distances.push_back(eucl_dist(i,j));
compressed_lower_distance_matrix dist(std::move(distances));
std::cout << "distance matrix with " << n << " points" << std::endl;
#endif
#ifdef FILE_FORMAT_DISTANCE_MATRIX
std::vector distances;
std::string line;
value_t value;
for (int i = 0; std::getline(input_stream, line); ++i) {
std::istringstream s(line);
for (int j = 0; j < i && s >> value; ++j) distances.push_back(value);
}
compressed_lower_distance_matrix dist(std::move(distances));
index_t n = dist.size();
std::cout << "distance matrix with " << n << " points" << std::endl;
#endif
#ifdef FILE_FORMAT_LOWER_DISTANCE_MATRIX
std::vector distances;
value_t value;
while (input_stream >> value) {
distances.push_back(value);
input_stream.ignore();
}
compressed_lower_distance_matrix dist(std::move(distances));
index_t n = dist.size();
std::cout << "distance matrix with " << n << " points" << std::endl;
#endif
#ifdef FILE_FORMAT_UPPER_DISTANCE_MATRIX
std::vector distances;
value_t value;
while (input_stream >> value) {
distances.push_back(value);
input_stream.ignore();
}
compressed_upper_distance_matrix dist(std::move(distances));
index_t n = dist.size();
std::cout << "distance matrix with " << n << " points" << std::endl;
#endif
#ifdef FILE_FORMAT_DIPHA
if (read(input_stream) != 8067171840) {
std::cerr << filename << " is not a Dipha file (magic number: 8067171840)" << std::endl;
exit(-1);
}
if (read(input_stream) != 7) {
std::cerr << filename << " is not a Dipha distance matrix (file type: 7)" << std::endl;
exit(-1);
}
index_t n = read(input_stream);
std::vector distances;
for (int i = 0; i < n; ++i)
for (int j = 0; j < n; ++j)
if (i > j)
distances.push_back(read(input_stream));
else
read(input_stream);
compressed_lower_distance_matrix dist(std::move(distances));
std::cout << "distance matrix with " << n << " points" << std::endl;
#endif
auto result = std::minmax_element(dist.distances.begin(), dist.distances.end());
std::cout << "value range: [" << *result.first << "," << *result.second << "]" << std::endl;
dim_max = std::min(dim_max, n - 2);
binomial_coeff_table binomial_coeff(n, dim_max + 2);
std::vector multiplicative_inverse(multiplicative_inverse_vector(modulus));
std::vector columns_to_reduce;
for (index_t index = n; index-- > 0;) columns_to_reduce.push_back(diameter_index_t(0, index));
for (index_t dim = 0; dim <= dim_max; ++dim) {
rips_filtration_comparator comp(dist, dim + 1, binomial_coeff);
rips_filtration_comparator comp_prev(dist, dim, binomial_coeff);
hash_map pivot_column_index;
pivot_column_index.reserve(columns_to_reduce.size());
compute_pairs(columns_to_reduce, pivot_column_index, dist, comp, comp_prev, dim, n,
threshold, modulus, multiplicative_inverse, binomial_coeff);
if (dim < dim_max)
assemble_columns_to_reduce(columns_to_reduce, pivot_column_index, comp, dim, n,
threshold, binomial_coeff);
}
}