CORDIC_Rotate_APFX/sources/CCordicRotateRom/CCordicRotateRom.hpp.in

/*
 *
 * Copyright 2022 Camille "DrasLorus" Monière.
 *
 * This file is part of CORDIC_Rotate_APFX.
 *
 * This program is free software: you can redistribute it and/or modify it under the terms of the GNU
 * Lesser General Public License as published by the Free Software Foundation, either version 3 of
 * the License, or (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without
 * even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
 * Lesser General Public License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public License along with this program.
 * If not, see <https://www.gnu.org/licenses/>.
 *
 */

#ifndef C_CORDIC_ROTATE_ROM_@CORDIC_W@_@CORDIC_STAGES@_@CORDIC_Q@_@CORDIC_DIVIDER@
#define C_CORDIC_ROTATE_ROM_@CORDIC_W@_@CORDIC_STAGES@_@CORDIC_Q@_@CORDIC_DIVIDER@

#include <climits>
#include <cmath>
#include <cstdint>
#include <cstdlib>

#include <complex>

#include <ap_fixed.h>
#include <ap_int.h>

#include "CCordicRotateRomTemplate.hpp"
#include "CordicRoms/cordic_rom_@ROM_TYPE@_@CORDIC_W@_@CORDIC_STAGES@_@CORDIC_Q@_@CORDIC_DIVIDER@.hpp"
#include "definitions.hpp"

#ifndef KN_STATIC_TABLE_DEFINED
#define KN_STATIC_TABLE_DEFINED 1
// ``` GNU Octave
// kn_values(X) = prod(1 ./ abs(1 + 1j * 2.^ (-(0:X))))
// ```
static constexpr double kn_values[7] = {
    0.70710678118655, 0.632455532033680, 0.613571991077900,
    0.608833912517750, 0.607648256256170, 0.607351770141300, 0.607277644093530};
#endif // KN_STATIC_TABLE_DEFINED

template <unsigned TIn_I>
class CCordicRotateRom<TIn_I, @ROM_TYPE@, @CORDIC_W@, @CORDIC_STAGES@, @CORDIC_Q@, @CORDIC_DIVIDER@> {
    static_assert(@CORDIC_W@ > 0, "Inputs can't be on zero bits.");
    static_assert(@CORDIC_STAGES@ < 8, "7 stages of CORDIC is the maximum supported.");
    static_assert(@CORDIC_STAGES@ > 1, "2 stages of CORDIC is the minimum.");
    static_assert(is_pow_2(@CORDIC_DIVIDER@), "divider must be a power of 2.");

public:
    static constexpr unsigned In_W      = @CORDIC_W@;
    static constexpr unsigned In_I      = TIn_I;
    static constexpr unsigned Out_W     = In_W + 2;
    static constexpr unsigned Out_I     = In_I + 2;
    static constexpr unsigned nb_stages = @CORDIC_STAGES@;
    static constexpr unsigned q         = @CORDIC_Q@;

    static constexpr uint64_t kn_i             = uint64_t(kn_values[nb_stages - 1] * double(1U << 4)); // 4 bits are enough
    static constexpr uint64_t in_scale_factor  = uint64_t(1U << (In_W - In_I));
    static constexpr uint64_t out_scale_factor = uint64_t(1U << (Out_W - Out_I));

    static constexpr double   rotation    = pi / @CORDIC_DIVIDER@;
    static constexpr unsigned addr_length = cordic_roms::@ROM_TYPE@_@CORDIC_W@_@CORDIC_STAGES@_@CORDIC_Q@_@CORDIC_DIVIDER@_size;

    static constexpr int64_t scale_cordic(int64_t in) {
        return in * kn_i / 16U;
    }

#if !defined(__SYNTHESIS__) && defined(SOFTWARE)
    static constexpr std::complex<int64_t> cordic(std::complex<int64_t> x_in,
                                                  uint8_t               counter) {

        int64_t A = x_in.real();
        int64_t B = x_in.imag();

        const uint8_t R    = cordic_roms::@ROM_TYPE@_@CORDIC_W@_@CORDIC_STAGES@_@CORDIC_Q@_@CORDIC_DIVIDER@[counter];
        uint8_t       mask = 0x01;
        if ((R & mask) == mask) {
            A = -A;
            B = -B;
        }

        for (uint8_t u = 1; u < nb_stages + 1; u++) {
            mask = mask << 1;

            const int64_t Ri = (R & mask) == mask ? 1 : -1;

            const int64_t I = A + Ri * (B / int64_t(1U << (u - 1)));
            B               = B - Ri * (A / int64_t(1U << (u - 1)));
            A               = I;
        }

        return {(A), (B)};
    }

    static constexpr double scale_cordic(double in) {
        return in * kn_values[nb_stages - 1];
    }

    static constexpr std::complex<double> cordic(std::complex<double> x_in,
                                                 uint8_t              counter) {
        const std::complex<int64_t> fx_x_in(int64_t(x_in.real() * double(in_scale_factor)),
                                            int64_t(x_in.imag() * double(in_scale_factor)));

        const std::complex<int64_t> fx_out = cordic(fx_x_in, counter);
        return {scale_cordic(double(fx_out.real())) / double(out_scale_factor), scale_cordic(double(fx_out.imag())) / double(out_scale_factor)};
    }

#endif

    static ap_int<Out_W> scale_cordic(const ap_int<Out_W> & in) {
        const ap_int<Out_W + 4> tmp = in * ap_uint<4>(kn_i);
        return ap_int<Out_W>(tmp >> 4);
    }

    static void cordic(const ap_int<In_W> & re_in, const ap_int<In_W> & im_in,
                       const ap_uint<8> & counter,
                       ap_int<Out_W> & re_out, ap_int<Out_W> & im_out) {

        const ap_uint<nb_stages + 1> R = cordic_roms::@ROM_TYPE@_@CORDIC_W@_@CORDIC_STAGES@_@CORDIC_Q@_@CORDIC_DIVIDER@[counter];

        ap_int<Out_W> A = bool(R[0]) ? ap_int<In_W>(-re_in) : re_in;
        ap_int<Out_W> B = bool(R[0]) ? ap_int<In_W>(-im_in) : im_in;

        for (uint8_t u = 1; u < nb_stages + 1; u++) { // nb_stages stages

            const bool Ri = bool(R[u]);

            // Results in (X / 2^(u - 1)), meaning only the
            // Out_W - u LSBs are meaninfull in shifted_X
            // Can't use range access since 11111111 (-1) would become 00001111 (15).
            // Would be possible if the loop is manually unrolled, to predict bitsize,
            // thus directly put 1111 into 4 bits (so still -1).
            const ap_int<Out_W> shifted_A = A >> (u - 1); // A(Out_W - 1, u - 1);
            const ap_int<Out_W> shifted_B = B >> (u - 1); // B(Out_W - 1, u - 1);

            const ap_int<Out_W> arc_step_A
                = Ri
                    ? ap_int<Out_W>(-shifted_A)
                    : shifted_A;
            const ap_int<Out_W> arc_step_B
                = Ri
                    ? shifted_B
                    : ap_int<Out_W>(-shifted_B);

            const ap_int<Out_W + 1> I = A + arc_step_B;
            B                         = B + arc_step_A;
            A                         = I;
        }

        re_out = A;
        im_out = B;
    }

    constexpr CCordicRotateRom() = default;
};

#endif // C_CORDIC_ROTATE_ROM_W_STAGES_Q_DIVIDER