mnncorrect.hpp

Go to the documentation of this file.

#ifndef MNNCORRECT_HPP
#define MNNCORRECT_HPP
 
#include <algorithm>
#include <vector>
#include <numeric>
#include <stdexcept>
#include <cstddef>
 
#include "knncolle/knncolle.hpp"
#include "sanisizer/sanisizer.hpp"
 
#include "AutomaticOrder.hpp"
#include "restore_input_order.hpp"
#include "utils.hpp"
 
namespace mnncorrect {
 
template<typename Index_, typename Float_, class Matrix_ = knncolle::Matrix<Index_, Float_> >
struct Options {
    int num_neighbors = 15;
 
    int num_steps = 1;
 
    std::shared_ptr<knncolle::Builder<Index_, Float_, Float_, Matrix_> > builder;
 
    MergePolicy merge_policy = MergePolicy::RSS;
 
    int num_threads = 1;
};
 
namespace internal {
 
template<typename Index_, typename Float_, class Matrix_>
void compute(const std::size_t num_dim, const std::vector<Index_>& num_obs, const std::vector<const Float_*>& batches, Float_* const output, const Options<Index_, Float_, Matrix_>& options) {
    auto builder = options.builder;
    if (!builder) {
        typedef knncolle::EuclideanDistance<Float_, Float_> Euclidean;
        builder.reset(new knncolle::VptreeBuilder<Index_, Float_, Float_, Matrix_, Euclidean>(std::make_shared<Euclidean>()));
    }
 
    AutomaticOrder<Index_, Float_, Matrix_> runner(
        num_dim,
        num_obs,
        batches,
        output,
        *builder,
        options.num_neighbors,
        options.num_steps,
        options.merge_policy,
        options.num_threads
    );
 
    runner.merge();
}
 
}
template<typename Index_, typename Float_, class Matrix_>
void compute(const std::size_t num_dim, const std::vector<Index_>& num_obs, const std::vector<const Float_*>& batches, Float_* const output, const Options<Index_, Float_, Matrix_>& options) {
    internal::compute(num_dim, num_obs, batches, output, options);
}
 
template<typename Index_, typename Float_, class Matrix_>
void compute(const std::size_t num_dim, const std::vector<Index_>& num_obs, const Float_* const input, Float_* const output, const Options<Index_, Float_, Matrix_>& options) {
    std::vector<const Float_*> batches;
    batches.reserve(num_obs.size());
 
    Index_ accumulated = 0;
    for (const auto n : num_obs) {
        batches.push_back(input + sanisizer::product_unsafe<std::size_t>(accumulated, num_dim));
 
        // After this check, all internal functions may assume that the total number of observations fits in an Index_.
        accumulated = sanisizer::sum<decltype(I(accumulated))>(accumulated, n);
    }
 
    compute(num_dim, num_obs, batches, output, options);
}
 
template<typename Index_, typename Float_, typename Batch_, class Matrix_>
void compute(const std::size_t num_dim, const Index_ num_obs, const Float_* const input, const Batch_* const batch, Float_* const output, const Options<Index_, Float_, Matrix_>& options) {
    const BatchIndex nbatches = (num_obs ? sanisizer::sum<BatchIndex>(*std::max_element(batch, batch + num_obs), 1) : static_cast<BatchIndex>(0));
    auto sizes = sanisizer::create<std::vector<Index_> >(nbatches);
    for (Index_ o = 0; o < num_obs; ++o) {
        ++sizes[batch[o]];
    }
 
    // Avoiding the need to allocate a temporary buffer
    // if we're already dealing with contiguous batches.
    bool already_sorted = true;
    for (Index_ o = 1; o < num_obs; ++o) {
       if (batch[o] < batch[o-1]) {
           already_sorted = false;
           break;
       }
    }
    if (already_sorted) {
        compute(num_dim, sizes, input, output, options);
        return;
    }
 
    Index_ accumulated = 0;
    auto offsets = sanisizer::create<std::vector<Index_> >(nbatches);
    std::vector<Float_> tmp(sanisizer::product<typename std::vector<Float_>::size_type>(num_dim, num_obs));
    auto ptrs = sanisizer::create<std::vector<const Float_*> >(nbatches);
    for (BatchIndex b = 0; b < nbatches; ++b) {
        ptrs[b] = tmp.data() + sanisizer::product_unsafe<std::size_t>(accumulated, num_dim);
        offsets[b] = accumulated;
        accumulated += sizes[b]; // this won't overflow as know that num_obs fits in an Index_.
    }
 
    for (Index_ o = 0; o < num_obs; ++o) {
        auto& offset = offsets[batch[o]];
        std::copy_n(
            input + sanisizer::product_unsafe<std::size_t>(o, num_dim),
            num_dim,
            tmp.data() + sanisizer::product_unsafe<std::size_t>(offset, num_dim)
        );
        ++offset;
    }
 
    internal::compute(num_dim, sizes, ptrs, output, options);
    internal::restore_input_order(num_dim, sizes, batch, output);
}
 
}
 
#endif