computation.hpp

Go to the documentation of this file.

 
#ifndef RTCTK_EXAMPLEDATATASK_COMPUTATION_HPP
#define RTCTK_EXAMPLEDATATASK_COMPUTATION_HPP
 
#include "rtctk/componentFramework/exceptions.hpp"
#include "rtctk/componentFramework/logger.hpp"
#include "rtctk/componentFramework/matrixBuffer.hpp"
 
#include <cblas.h>
#include <chrono>
#include <lapacke.h>
#include <numeric>
#include <string>
#include <vector>
 
// for python computation
#include "rtctk/componentFramework/matrixBufferPy.hpp"
 
namespace rtctk::exampleDataTask {
 
using namespace rtctk::componentFramework;
 
class Computation {
public:
    enum Algorithm { SimpleInversion, PythonInversion };
 
    void SetStaticConfig(unsigned n_slopes, unsigned n_acts) {
        m_n_slopes = n_slopes;
        m_n_acts = n_acts;
 
        // set vectors and matrix to correct size.
        m_IM.resize(m_n_acts, m_n_slopes);
        m_CM.resize(m_n_slopes, m_n_acts);
        m_ipiv.resize(m_n_acts);
 
        // some simple debug output
        LOG4CPLUS_DEBUG(GetLogger(), "m_n_slopes: " << m_n_slopes);
        LOG4CPLUS_DEBUG(GetLogger(), "m_n_acts: " << m_n_acts);
    }
 
    void SetDynamicConfig(MatrixBuffer<float>&& data) {
        m_IM = std::move(data);
 
        if (m_n_acts != (int)m_IM.GetNrows() || m_n_slopes != (int)m_IM.GetNcols()) {
            std::stringstream err_text;
            err_text << "IM wrong shape, "
                     << "expected:" << m_n_acts << " x " << m_n_slopes
                     << "received: " << m_IM.GetNrows() << " x " << m_IM.GetNcols();
            using rtctk::componentFramework::RtctkException;
            CII_THROW(RtctkException, err_text.str());
        }
    }
 
    struct Result {
        MatrixBuffer<float> const& cm;
 
        struct {
            std::chrono::duration<double> elapsed;
        } stats;
    };
 
    Result Compute(Algorithm algorithm) {
        auto time_start = std::chrono::system_clock::now();
 
        if (algorithm == SimpleInversion) {
            // invert the matrix in a very dumv way.
            memcpy((void*)m_CM.data(), (void*)m_IM.data(), m_n_slopes * m_n_acts * sizeof(float));
            LAPACKE_sgetrf(
                LAPACK_ROW_MAJOR, m_n_acts, m_n_slopes, m_CM.data(), m_n_slopes, m_ipiv.data());
            LAPACKE_sgetri(LAPACK_ROW_MAJOR, m_n_slopes, m_CM.data(), m_n_acts, m_ipiv.data());
        } else {
            // set up python and load lib (do this earlier to speed up things)
            py::scoped_interpreter py_guard;
            auto py_matrix_buffer_import = LoadMatrixBufferPyBinds();
            auto py_compute_module = py::module::import("rtctkExampleDataTaskPyLib");
            auto py_inversion = py_compute_module.attr("inversion");
            // do the computation
            py_inversion(&m_IM, &m_CM);
        }
 
        auto elapsed = std::chrono::system_clock::now() - time_start;
 
        return {m_CM, elapsed};
    }
 
private:
    // static config
    int m_n_slopes;
    int m_n_acts;
 
    // input data
    MatrixBuffer<float> m_IM;
 
    // output data
    MatrixBuffer<float> m_CM;
 
    std::vector<int> m_ipiv;
};
 
} // namespace rtctk::exampleDataTask
 
#endif  // RTCTK_EXAMPLEDATATASK_COMPUTATION_HPP