Initial version of the sample calculator

libpsk/include/psk/spectrum_calculator.hpp

@@ -0,0 +1,74 @@
+/*
+ * spectrum_calculator.hpp
+ *
+ *  Created on: 15-05-2013
+ *      Author: lucck
+ */
+/*----------------------------------------------------------*/
+#ifndef LIBPSK_SPECTRUM_CALCULATOR_HPP_
+#define LIBPSK_SPECTRUM_CALCULATOR_HPP_
+/*----------------------------------------------------------*/
+#include <cstddef>
+#include <complex>
+#include <cstring>
+/*----------------------------------------------------------*/
+namespace ham {
+namespace psk {
+
+/*----------------------------------------------------------*/
+//Class spectrum calculator calculate the current spectrum
+class spectrum_calculator {
+	//Number of samples
+	static constexpr auto MLOG2 = 9;
+	static constexpr auto WIDTH = 1<<MLOG2;
+	//Make object noncopyable
+	spectrum_calculator(const spectrum_calculator&) = delete;
+	spectrum_calculator& operator=(const spectrum_calculator&) = delete;
+public:
+	//Constructor
+	spectrum_calculator()
+	{}
+	//Spectrum data type
+	typedef short pow_t;
+public:
+	//Copy sample buffer
+	void copy_samples( const pow_t* samples )
+	{
+		std::memcpy( m_real, samples, sizeof(pow_t) * WIDTH );
+		m_energy_calculated = false;
+	}
+	const pow_t& operator[]( size_t idx )
+	{
+		if( !m_energy_calculated )
+		{
+			m_energy_calculated = true;
+			calculate_samples();
+		}
+		//TODO: Range checking
+		return m_real[idx];
+	}
+private:
+	//Calculate samplees
+	void calculate_samples();
+private:
+	//Private pointer for store data
+	std::complex<pow_t>m_cplx[WIDTH];
+	/*
+	    For any complex number z, reinterpret_cast<T(&)[2]>(z)[0] is the real part of z
+	    and reinterpret_cast<T(&)[2]>(z)[1] is the imaginary part of z.
+		For any pointer to an element of an array of complex numbers p
+		and any valid array index i, reinterpret_cast<T*>(p)[2*i] is the real part
+		of the complex number p[i], and reinterpret_cast<T*>(p)[2*i + 1] is the
+		imaginary part of the complex number p[i]
+	 */
+	pow_t* const m_real { reinterpret_cast<pow_t* const>(m_cplx) };
+	bool m_energy_calculated { false };
+};
+
+/*----------------------------------------------------------*/
+
+} /* namespace psk */
+} /* namespace ham */
+/*----------------------------------------------------------*/
+#endif /* LIBPSK_SPECTRUM_CALCULATOR_HPP_ */
+/*----------------------------------------------------------*/

libpsk/src/spectrum_calculator.cpp

@@ -0,0 +1,44 @@
+/*
+ * spectrum_calculator.cpp
+ *
+ *  Created on: 15-05-2013
+ *      Author: lucck
+ */
+
+/*----------------------------------------------------------*/
+#include "psk/spectrum_calculator.hpp"
+#include "dsp/fft.h"
+#include "dsp/sqrt_int.hpp"
+
+/*----------------------------------------------------------*/
+namespace ham {
+namespace psk {
+/*----------------------------------------------------------*/
+namespace
+{
+	//Inline POW2
+	constexpr inline unsigned pow2( unsigned u )
+	{
+		return u * u;
+	}
+}
+/*----------------------------------------------------------*/
+//Calculate samplees
+void spectrum_calculator::calculate_samples()
+{
+	namespace fft = dsp::refft;
+	namespace dint = dsp::integer;
+	//Calculate FFT real transform
+	fft::fft_real( m_cplx, m_real, MLOG2 );
+	//Calculate energies
+	//TODO: FInaly calculate only to WITDTH/2
+	for( int i=0; i<WIDTH; i++ )
+	{
+		m_real[i] = dint::sqrt( pow2(m_cplx[i].real()) + pow2(m_cplx[i].imag()) );
+	}
+}
+/*----------------------------------------------------------*/
+
+} /* namespace psk */
+} /* namespace ham */
+/*----------------------------------------------------------*/

tests/fftsamples/Makefile

@@ -34,7 +34,7 @@ LISTING = n
 DEBUG ?=  y
 
 #Source C++ files
-CPPSRC += $(wildcard *.cpp)
+CPPSRC += $(wildcard *.cpp) ../../libpsk/src/spectrum_calculator.cpp
 
 include ../../unix.mk
 

tests/fftsamples/test.cpp

@@ -5,6 +5,7 @@
 #include <iostream>
 #include <fstream>
 #include "dsp/fft.h"
+#include "psk/spectrum_calculator.hpp"
 
 /* ---------------------------------------------------------------- */
 template <typename T, typename K >
@@ -59,7 +60,7 @@ float cplx_abs( std::complex<short> v )
 	for(int s=0; s< N_SAMPL; s++ )
 	  abst[s] = 20.0*log10((abst[s]+1)/MAX);
  */
-constexpr auto N_SAMPL = 1024*1;
+constexpr auto N_SAMPL = 512;
 constexpr auto FREQUENCY = N_SAMPL / 10;
 
 int main()
@@ -80,7 +81,14 @@ int main()
         t[s] = s/(float)N_SAMPL;
     }
 	plot( t, y, N_SAMPL , "SI");
-    if( 0 )
+
+#if 1
+	ham::psk::spectrum_calculator spc;
+	spc.copy_samples( y );
+	plot( t, &spc[0], N_SAMPL, "FINAL");
+
+#else
+	if( 0 )
     {
 		std::complex<short> x[N_SAMPL];
 		bzero(x, sizeof x);
@@ -97,6 +105,7 @@ int main()
     		abst[s] = cplx_abs( x[s] );
     	plot( t, abst, N_SAMPL , "FFT");
     }
+#endif
 }