Added code docs

gnuradio/gr-dji_droneid/include/gnuradio/dji_droneid/utils.h

@@ -81,6 +81,32 @@ public:
      */
     static std::pair<uint16_t, uint16_t> get_cyclic_prefix_lengths(float sample_rate);
 
+    /**
+     * Calculate the mean of a complex vector
+     * @param samples Pointer to complex vector
+     * @param sample_count Number of samples in the complex vector
+     * @return Mean of the complex vector
+     */
+    static std::complex<float> mean_fast(const std::complex<float> * samples, uint32_t sample_count);
+
+    /**
+     * Calculate the variance of a complex vector where the mean is already zero
+     * @param samples Pointer to complex vector
+     * @param sample_count Number of samples in the complex vector
+     * @return Variance of the complex vector
+     */
+    static float variance_no_mean(const std::complex<float> * samples, uint32_t sample_count);
+
+    /**
+     * Calculate the variance of a complex vector where the mean is not known to be zero
+     *
+     * This function is mostly the same as utils::variance_no_mean but does have to calculate the
+     * mean value of the input vector and subtract it from each sample, so this will take longer
+     * @param samples Pointer to complex vector
+     * @param sample_count Number of samples in the complex vector
+     * @return Variance of the complex vector
+     */
+    static float variance(const std::complex<float> * samples, uint32_t sample_count);
 private:
 };
 

gnuradio/gr-dji_droneid/lib/qa_utils.cc

@@ -113,6 +113,47 @@ struct TestFixture {
             static_cast<uint16_t>(short_cp_len[0])
         };
     }
+
+    static std::complex<float> mean(const std::vector<std::complex<float>> & samples) {
+        matlab_engine->setVariable("samples", factory.createArray({1, samples.size()}, samples.begin(), samples.end()));
+        matlab_engine->eval(u"m = mean(samples);");
+        const matlab::data::TypedArray<std::complex<float>> mean_val = matlab_engine->getVariable("m");
+
+        BOOST_REQUIRE_EQUAL(mean_val.getNumberOfElements(), 1);
+
+        matlab_engine->eval(u"clear samples m");
+
+        return mean_val[0];
+    }
+
+    static float var(const std::vector<std::complex<float>> & samples) {
+        matlab_engine->setVariable("samples", factory.createArray({1, samples.size()}, samples.begin(), samples.end()));
+        matlab_engine->eval(u"output = single(var(samples));");
+        const matlab::data::TypedArray<float> output = matlab_engine->getVariable("output");
+        matlab_engine->eval(u"clear output samples");
+        return output[0];
+    }
+
+    static std::vector<std::complex<float>> remove_mean(const std::vector<std::complex<float>> & samples) {
+        matlab_engine->setVariable("samples", factory.createArray({1, samples.size()}, samples.begin(), samples.end()));
+        matlab_engine->eval(u"samples = single(samples - mean(samples));");
+        const matlab::data::TypedArray<std::complex<float>> output_samples = matlab_engine->getVariable("samples");
+        matlab_engine->eval(u"clear samples");
+
+        BOOST_REQUIRE_EQUAL(output_samples.getNumberOfElements(), samples.size());
+
+        return {output_samples.begin(), output_samples.end()};
+    }
+
+    static std::vector<std::complex<float>> create_test_vector(const uint32_t sample_count) {
+        matlab_engine->setVariable("sample_count", factory.createScalar(static_cast<double>(sample_count)));
+        matlab_engine->eval(u"samples = single(complex(randn(1, sample_count), randn(1, sample_count)));");
+        const matlab::data::TypedArray<std::complex<float>> samples_array = matlab_engine->getVariable("samples");
+        std::vector<std::complex<float>> samples(samples_array.begin(), samples_array.end());
+        matlab_engine->eval(u"clear samples");
+
+        return samples;
+    }
 };
 
 BOOST_FIXTURE_TEST_SUITE(Utils_Test_Suite, TestFixture);
@@ -201,6 +242,57 @@ BOOST_AUTO_TEST_CASE(test_utils__get_cyclic_prefix) {
     }
 }
 
+BOOST_AUTO_TEST_CASE(test_utils__mean_fast) {
+    const uint32_t iters = 100;
+    const uint32_t length = 10000;
+
+    matlab_engine->setVariable("sample_count", factory.createScalar(static_cast<double>(length)));
+    for (uint32_t iter = 0; iter < iters; iter++) {
+        matlab_engine->eval(u"samples = single(complex(randn(1, sample_count), randn(1, sample_count)));");
+        const matlab::data::TypedArray<std::complex<float>> samples_array = matlab_engine->getVariable("samples");
+        const std::vector<std::complex<float>> samples(samples_array.begin(), samples_array.end());
+
+        const auto expected = mean(samples);
+        const auto calculated = utils::mean_fast(&samples[0], samples.size());
+
+        // MATLAB does all of its arithmetic in double precision, and the mean_fast function does everything in single
+        // This means there is going to be a large delta between the values simply due to rounding errors.  The value
+        // below is a percentage (eg 0.2 == a 0.2% delta) that the difference can be before failing the test
+        const auto accuracy = 0.2;
+
+        BOOST_REQUIRE_CLOSE(expected.real(), calculated.real(), accuracy);
+        BOOST_REQUIRE_CLOSE(expected.imag(), calculated.imag(), accuracy);
+    }
+}
+
+BOOST_AUTO_TEST_CASE(test_utils__variance) {
+    const uint32_t iters = 100;
+    const uint32_t length = 10000;
+
+    for (auto iter = decltype(iters){0}; iter < iters; iter++) {
+        const auto samples = create_test_vector(length);
+
+        const auto expected = var(samples);
+        const auto calculated = utils::variance(&samples[0], samples.size());
+
+        BOOST_REQUIRE_CLOSE(expected, calculated, 0.002);
+    }
+}
+
+BOOST_AUTO_TEST_CASE(test_utils__variance_no_mean) {
+    const uint32_t iters = 100;
+    const uint32_t length = 10000;
+
+    for (auto iter = decltype(iters){0}; iter < iters; iter++) {
+        const auto samples = remove_mean(create_test_vector(length));
+
+        const auto expected = var(samples);
+        const auto calculated = utils::variance_no_mean(&samples[0], samples.size());
+
+        BOOST_REQUIRE_CLOSE(expected, calculated, 0.002);
+    }
+}
+
 BOOST_AUTO_TEST_SUITE_END()
 
 } /* namespace dji_droneid */

gnuradio/gr-dji_droneid/lib/utils.cc

@@ -11,6 +11,7 @@
 #include <gnuradio/fft/fft.h>
 #include <gnuradio/fft/fft_shift.h>
 #include <fftw3.h>
+#include <volk/volk.h>
 
 namespace gr {
 namespace dji_droneid {
@@ -130,5 +131,46 @@ std::pair<uint16_t, uint16_t> utils::get_cyclic_prefix_lengths(const float sampl
     );
 }
 
+std::complex<float> utils::mean_fast(const std::complex<float> * const samples, const uint32_t sample_count)
+{
+    // Treat the vector of complex floats as a single vector of float values
+    auto sample_floats = reinterpret_cast<const float *>(samples);
+    float real = 0, imag = 0;
+
+    for (uint32_t idx = 0; idx < sample_count; idx++) {
+        real += *sample_floats++;
+        imag += *sample_floats++;
+    }
+
+    real = real / static_cast<float>(sample_count);
+    imag = imag / static_cast<float>(sample_count);
+
+    return {real, imag};
+}
+
+float utils::variance_no_mean(const std::complex<float> * const samples, const uint32_t sample_count)
+{
+    float total = 0;
+
+    for (auto idx = decltype(sample_count){0}; idx < sample_count; idx++) {
+        const auto & sample = samples[idx];
+        total += (sample.real() * sample.real()) + (sample.imag() * sample.imag());
+    }
+
+    return total / static_cast<float>(sample_count - 1);
+}
+float utils::variance(const std::complex<float> * const samples, const uint32_t sample_count)
+{
+    float total = 0;
+    const auto mean = mean_fast(samples, sample_count);
+
+    for (auto idx = decltype(sample_count){0}; idx < sample_count; idx++) {
+        auto sample = samples[idx] - mean;
+        total += (sample.real() * sample.real()) + (sample.imag() * sample.imag());
+    }
+
+    return total / static_cast<float>(sample_count - 1);
+}
+
 } /* namespace dji_droneid */
 } /* namespace gr */