From 4764a52de827cbbe37e4d5adbb5e1c6634d5ea2f Mon Sep 17 00:00:00 2001 From: "J.A. de Jong - Redu-Sone B.V., ASCEE V.O.F" Date: Tue, 11 Oct 2022 14:50:44 +0200 Subject: [PATCH] Improved code for cubes. RtAps is not working. Don't know why. TimeBuffer code has better readability. Bugfix in output from Column to Numpy. GIL release for AvPowerSpectra::compute via Pybind. --- src/lasp/device/lasp_daq.cpp | 4 +- src/lasp/device/lasp_daqdata.cpp | 9 +++ src/lasp/device/lasp_uldaq.cpp | 6 +- src/lasp/dsp/lasp_avpowerspectra.cpp | 34 +++++----- src/lasp/dsp/lasp_avpowerspectra.h | 14 ++-- src/lasp/dsp/lasp_mathtypes.h | 3 +- src/lasp/dsp/lasp_rtaps.cpp | 76 +++++++++++++--------- src/lasp/dsp/lasp_rtaps.h | 22 +------ src/lasp/dsp/lasp_timebuffer.cpp | 46 ++++++------- src/lasp/pybind11/arma_npy.h | 3 +- src/lasp/pybind11/lasp_dsp_pybind.cpp | 38 ++++++++--- src/lasp/pybind11/lasp_pyindatahandler.cpp | 4 +- test/test_aps.py | 2 +- test/test_biquadbank.py | 2 +- test/test_fft.py | 6 +- test/test_ps.py | 2 +- 16 files changed, 157 insertions(+), 114 deletions(-) diff --git a/src/lasp/device/lasp_daq.cpp b/src/lasp/device/lasp_daq.cpp index a5d1bd5..7893929 100644 --- a/src/lasp/device/lasp_daq.cpp +++ b/src/lasp/device/lasp_daq.cpp @@ -53,11 +53,11 @@ Daq::Daq(const DeviceInfo &devinfo, const DaqConfiguration &config) if (!config.match(devinfo)) { throw rte("DaqConfiguration does not match device info"); } - if (neninchannels(false) > ninchannels) { + if (neninchannels(false) > devinfo.ninchannels) { throw rte( "Number of enabled input channels is higher than device capability"); } - if (nenoutchannels() > noutchannels) { + if (nenoutchannels() > devinfo.noutchannels) { throw rte( "Number of enabled output channels is higher than device capability"); } diff --git a/src/lasp/device/lasp_daqdata.cpp b/src/lasp/device/lasp_daqdata.cpp index dd52d0d..abae921 100644 --- a/src/lasp/device/lasp_daqdata.cpp +++ b/src/lasp/device/lasp_daqdata.cpp @@ -25,6 +25,7 @@ DaqData::DaqData(const us nframes, const us nchannels, DEBUGTRACE_ENTER; DEBUGTRACE_PRINT(sw); + assert(sw > 0 && sw <= 8); _data = new (std::align_val_t{8}) byte_t[sw * nchannels * nframes]; if (!_data) { throw rte("Could not allocate memory for DaqData!"); @@ -180,15 +181,18 @@ dmat DaqData::toFloat() const { return toFloat(); break; case (DataType::dtype_int16): + DEBUGTRACE_PRINT("Dtype = int16"); return toFloat(); break; case (DataType::dtype_int32): return toFloat(); break; case (DataType::dtype_fl32): + DEBUGTRACE_PRINT("Dtype = float32"); return toFloat(); break; case (DataType::dtype_fl64): + DEBUGTRACE_PRINT("Dtype = float64"); return toFloat(); break; default: @@ -206,5 +210,10 @@ void DaqData::print() const { cout << "First sample of first channel (as float)" << toFloat(0,0) << endl; cout << "Last sample of first channel (as float)" << toFloat(nframes-1,0) << endl; cout << "Last sample of last channel (as float)" << toFloat(nframes-1,nchannels-1) << endl; + dmat data = toFloat(); + vrd max = arma::max(data, 0); + vrd min = arma::min(data, 0); + cout << "Maximum value in buf: " << max << endl; + cout << "Minumum value in buf: " << min << endl; } diff --git a/src/lasp/device/lasp_uldaq.cpp b/src/lasp/device/lasp_uldaq.cpp index c32c6cb..7f3c567 100644 --- a/src/lasp/device/lasp_uldaq.cpp +++ b/src/lasp/device/lasp_uldaq.cpp @@ -260,9 +260,9 @@ public: us monitorOffset = monitorOutput ? 1 : 0; /* /// Put the output monitor in front */ if (monitorOutput) { - for (us sample = 0; sample < nFramesPerBlock; sample++) { - data.value(sample, 0) = - buf[totalOffset + sample * nchannels + (nchannels - 1)]; + for (us frame = 0; frame < nFramesPerBlock; frame++) { + data.value(frame, 0) = + buf[totalOffset + (frame * nchannels) + (nchannels - 1)]; } } diff --git a/src/lasp/dsp/lasp_avpowerspectra.cpp b/src/lasp/dsp/lasp_avpowerspectra.cpp index e8bc71d..7832f66 100644 --- a/src/lasp/dsp/lasp_avpowerspectra.cpp +++ b/src/lasp/dsp/lasp_avpowerspectra.cpp @@ -24,7 +24,7 @@ PowerSpectra::PowerSpectra(const vd &window) DEBUGTRACE_PRINT(win_pow); } -arma::Cube PowerSpectra::compute(const dmat &input) { +ccube PowerSpectra::compute(const dmat &input) { /// Run very often. Silence this one. /* DEBUGTRACE_ENTER; */ @@ -36,7 +36,7 @@ arma::Cube PowerSpectra::compute(const dmat &input) { cmat rfft = _fft.fft(input_tmp); - arma::cx_cube output(rfft.n_rows, input.n_cols, input.n_cols); + ccube output(rfft.n_rows, input.n_cols, input.n_cols); for (us i = 0; i < input.n_cols; i++) { /// This one can be run in parallel without any problem. Note that it is /// the inner loop that is run in parallel. @@ -63,6 +63,7 @@ AvPowerSpectra::AvPowerSpectra(const us nfft, const Window::WindowType w, } _overlap_keep = (nfft * overlap_percentage) / 100; + DEBUGTRACE_PRINT(_overlap_keep); if (_overlap_keep >= nfft) { throw rte("Overlap is too high. Results in no jump. Please " "choose a smaller overlap percentage or a higher nfft"); @@ -77,39 +78,40 @@ AvPowerSpectra::AvPowerSpectra(const us nfft, const Window::WindowType w, DEBUGTRACE_PRINT(_alpha); } } +void AvPowerSpectra::reset() { + _timeBuf.reset(); + _est.reset(); -std::optional AvPowerSpectra::compute(const dmat &timedata) { +} - /* DEBUGTRACE_ENTER; */ - /* DEBUGTRACE_PRINT(timedata.n_rows); */ - /* DEBUGTRACE_PRINT(_ps.nfft); */ +std::optional AvPowerSpectra::compute(const dmat &timedata) { + + DEBUGTRACE_ENTER; + DEBUGTRACE_PRINT(timedata.n_rows); + DEBUGTRACE_PRINT(_ps.nfft); _timeBuf.push(timedata); - if (_est.n_cols == 0) { - _est = arma::cx_cube(_ps.nfft / 2 + 1, timedata.n_cols, timedata.n_cols, - arma::fill::zeros); - } - - std::optional res; us i = 0; while (auto samples = _timeBuf.pop(_ps.nfft, _overlap_keep)) { - /* DEBUGTRACE_PRINT((int)_mode); */ + DEBUGTRACE_PRINT((int)_mode); switch (_mode) { case (Mode::Spectrogram): { + DEBUGTRACE_PRINT("Spectrogram mode"); _est = _ps.compute(samples.value()); } break; case (Mode::Averaging): { + DEBUGTRACE_PRINT("Averaging mode"); n_averages++; if (n_averages == 1) { _est = _ps.compute(samples.value()); } else { - _est = _est * (static_cast(n_averages - 1) / n_averages) + + _est = (static_cast(n_averages - 1) / n_averages) * _est + _ps.compute(samples.value()) / n_averages; } } break; case (Mode::Leaking): { - /* DEBUGTRACE_PRINT("Leaking mode"); */ + DEBUGTRACE_PRINT("Leaking mode"); if (arma::size(_est) == arma::size(0, 0, 0)) { _est = _ps.compute(samples.value()); } else { @@ -124,7 +126,7 @@ std::optional AvPowerSpectra::compute(const dmat &timedata) { } return std::nullopt; } -std::optional AvPowerSpectra::get_est() { +std::optional AvPowerSpectra::get_est() const { if (_est.n_cols > 0) return _est; return std::nullopt; diff --git a/src/lasp/dsp/lasp_avpowerspectra.h b/src/lasp/dsp/lasp_avpowerspectra.h index ed12540..95092e4 100644 --- a/src/lasp/dsp/lasp_avpowerspectra.h +++ b/src/lasp/dsp/lasp_avpowerspectra.h @@ -65,7 +65,7 @@ public: * can be accessed by: C(f, i, j). Storage is such that frequency components * are contiguous. */ - arma::Cube compute(const dmat &input); + ccube compute(const dmat &input); }; /** @@ -89,7 +89,7 @@ class AvPowerSpectra { /** * @brief Current estimate of cross-spectral density */ - arma::cx_cube _est; + ccube _est; /** * @brief Buffer of storage of time data. @@ -128,7 +128,11 @@ public: AvPowerSpectra(const AvPowerSpectra &) = delete; AvPowerSpectra &operator=(const AvPowerSpectra &) = delete; - void reset() { _est.reset(); } + /** + * @brief Reset to empty state. Clears the time buffer and sets estimator to + * empty. + */ + void reset(); /** * @brief Compute an update of the power spectra based on given time data. @@ -143,7 +147,7 @@ public: * new estimate. It can be checked by operator bool(). * */ - std::optional compute(const dmat &timedata); + std::optional compute(const dmat &timedata); /** * @brief Returns the latest estimate of cps (cross-power spectra. @@ -151,7 +155,7 @@ public: * @return Pointer (reference, not owning!) to spectral estimate date, or * nullptr, in case nothing could yet be computed. */ - std::optional get_est(); + std::optional get_est() const; /** * @brief The overlap is rounded to a certain amount of time samples. This diff --git a/src/lasp/dsp/lasp_mathtypes.h b/src/lasp/dsp/lasp_mathtypes.h index 8d3102f..4b75dd4 100644 --- a/src/lasp/dsp/lasp_mathtypes.h +++ b/src/lasp/dsp/lasp_mathtypes.h @@ -47,6 +47,7 @@ using vc = arma::Col; using vrc = arma::Row; using dmat = arma::Mat; using cmat = arma::Mat; -using cube = arma::Cube; +using ccube = arma::Cube; +using dcube = arma::Cube; const d number_pi = arma::datum::pi; diff --git a/src/lasp/dsp/lasp_rtaps.cpp b/src/lasp/dsp/lasp_rtaps.cpp index 7cfea0c..1675afe 100644 --- a/src/lasp/dsp/lasp_rtaps.cpp +++ b/src/lasp/dsp/lasp_rtaps.cpp @@ -6,42 +6,55 @@ using std::cerr; using std::endl; +RtAps::RtAps(StreamMgr &mgr, const Filter *freqWeightingFilter, + const us nfft, + const Window::WindowType w, + const d overlap_percentage, const d time_constant) + : ThreadedInDataHandler(mgr), + _ps(nfft, w, overlap_percentage, time_constant) { + + /* if (freqWeightingFilter != nullptr) { */ + /* _filterPrototype = freqWeightingFilter->clone(); */ + /* } */ + + start(); +} +RtAps::~RtAps() { + std::scoped_lock lck(_mtx); + stop(); +} bool RtAps::inCallback_threaded(const DaqData &data) { - /* DEBUGTRACE_ENTER; */ + DEBUGTRACE_ENTER; std::scoped_lock lck(_mtx); dmat fltdata = data.toFloat(); - const us nchannels = fltdata.n_cols; + data.print(); + /* const us nchannels = fltdata.n_cols; */ - if (_filterPrototype) { + /* if (_filterPrototype) { */ - // Adjust number of filters, if necessary - if (nchannels > _freqWeightingFilter.size()) { - while (nchannels > _freqWeightingFilter.size()) { - _freqWeightingFilter.emplace_back(_filterPrototype->clone()); - } + /* // Adjust number of filters, if necessary */ + /* if (nchannels > _freqWeightingFilter.size()) { */ + /* while (nchannels > _freqWeightingFilter.size()) { */ + /* _freqWeightingFilter.emplace_back(_filterPrototype->clone()); */ + /* } */ - for (auto &filter : _freqWeightingFilter) { - filter->reset(); - } - } + /* for (auto &filter : _freqWeightingFilter) { */ + /* filter->reset(); */ + /* } */ + /* } */ - // Apply filtering -#pragma omp parallel for - for (us i = 0; i < nchannels; i++) { - vd col = fltdata.col(i); - _freqWeightingFilter.at(i)->filter(col); - fltdata.col(i) = col; - } - } // End of if(_filterPrototype) + /* // Apply filtering */ + /* #pragma omp parallel for */ + /* for (us i = 0; i < nchannels; i++) { */ + /* vd col = fltdata.col(i); */ + /* _freqWeightingFilter.at(i)->filter(col); */ + /* fltdata.col(i) = col; */ + /* } */ + /* } // End of if(_filterPrototype) */ - std::optional res = _ps.compute(fltdata); - - if (res.has_value()) { - /* DEBUGTRACE_PRINT("Data ready!"); */ - _latest_est = std::make_unique(std::move(res.value())); - } + _ps.compute(fltdata); return true; } @@ -53,13 +66,18 @@ void RtAps::reset(const Daq *daq) { // Explicitly say DEBUGTRACE_ENTER; std::scoped_lock lck(_mtx); _ps.reset(); - _latest_est.reset(); } -std::unique_ptr RtAps::getCurrentValue() { +std::unique_ptr RtAps::getCurrentValue() { /* DEBUGTRACE_ENTER; */ std::scoped_lock lck(_mtx); + auto est = _ps.get_est(); + return std::make_unique(est.value_or(ccube())); - return std::move(_latest_est); + /* return std::move(_latest_est); */ + /* if (_latest_est) { */ + /* return std::make_unique(cube(*_latest_est)); */ + /* } else */ + /* return nullptr; */ } diff --git a/src/lasp/dsp/lasp_rtaps.h b/src/lasp/dsp/lasp_rtaps.h index b48f1b5..26a0510 100644 --- a/src/lasp/dsp/lasp_rtaps.h +++ b/src/lasp/dsp/lasp_rtaps.h @@ -24,8 +24,6 @@ class RtAps : public ThreadedInDataHandler { std::mutex _mtx; std::unique_ptr _filterPrototype; std::vector> _freqWeightingFilter; - bool _data_ready = false; - std::unique_ptr _latest_est; AvPowerSpectra _ps; @@ -40,29 +38,15 @@ public: */ RtAps(StreamMgr &mgr, const Filter *freqWeightingFilter, const us nfft = 2048, const Window::WindowType w = Window::WindowType::Hann, - const d overlap_percentage = 50., const d time_constant = -1) - : ThreadedInDataHandler(mgr), - _ps(nfft, w, overlap_percentage, time_constant) { - - { - std::scoped_lock lck(_mtx); - if (freqWeightingFilter != nullptr) { - _filterPrototype = freqWeightingFilter->clone(); - } - } - start(); - } - ~RtAps() { - std::scoped_lock lck(_mtx); - stop(); - } + const d overlap_percentage = 50., const d time_constant = -1); + ~RtAps(); /** * @brief Get the latest estimate of the power spectra * * @return Optionally, if available, the latest values */ - std::unique_ptr getCurrentValue(); + std::unique_ptr getCurrentValue(); bool inCallback_threaded(const DaqData &) override final; void reset(const Daq *) override final; diff --git a/src/lasp/dsp/lasp_timebuffer.cpp b/src/lasp/dsp/lasp_timebuffer.cpp index ffcf29d..f2d4fe1 100644 --- a/src/lasp/dsp/lasp_timebuffer.cpp +++ b/src/lasp/dsp/lasp_timebuffer.cpp @@ -1,6 +1,6 @@ /* #define DEBUGTRACE_ENABLED */ -#include "debugtrace.hpp" #include "lasp_timebuffer.h" +#include "debugtrace.hpp" #include #include #include @@ -23,9 +23,9 @@ public: } void push(const dmat &mat) { DEBUGTRACE_ENTER; -#if LASP_DEBUG==1 - if(!_storage.empty()) { - if(mat.n_cols != _storage.front().n_cols) { +#if LASP_DEBUG == 1 + if (!_storage.empty()) { + if (mat.n_cols != _storage.front().n_cols) { throw rte("Invalid number of channels in mat"); } } @@ -35,32 +35,36 @@ public: } } - std::optional pop(const us nsamples, us keep) { + std::optional pop(const us nframes, const us keep) { DEBUGTRACE_ENTER; - if (keep > nsamples) - throw rte("keep should be <= nsamples"); + if (keep >= nframes) { + throw rte("keep should be < nframes"); + } - if (nsamples <= n_frames()) { + if (nframes <= n_frames()) { assert(!_storage.empty()); - dmat res(nsamples, _storage.front().n_cols); + dmat res(nframes, _storage.front().n_cols); - for (us i = 0; i < nsamples; i++) { + us j = 0; + for (us i = 0; i < nframes; i++) { + + if (i + keep < nframes) { + // Just pop elements and copy over + res.row(i) = _storage.front(); + _storage.pop_front(); + + } else { - if (i + keep >= nsamples) { // Suppose keep == 0, then we never arrive here // Suppose keep == 1, then storage[0] is copyied over. // Suppose keep == 2, then storage[0[ and storage[1] is copyied over. // Etc. - res.row(i) = _storage[i + keep - nsamples]; - } else { - - // Just pop elements and copy over - res.row(i) = _storage.front(); - _storage.pop_front(); + res.row(i) = _storage[j]; + j++; } } return res; @@ -79,11 +83,9 @@ public: }; TimeBuffer::TimeBuffer() : _imp(std::make_unique()) {} -std::optional TimeBuffer::pop(const us n_rows,const us keep) { +std::optional TimeBuffer::pop(const us n_rows, const us keep) { return _imp->pop(n_rows, keep); } -TimeBuffer::~TimeBuffer(){} -void TimeBuffer::reset() { - _imp->reset(); -} +TimeBuffer::~TimeBuffer() {} +void TimeBuffer::reset() { _imp->reset(); } void TimeBuffer::push(const dmat &dat) { _imp->push(dat); } diff --git a/src/lasp/pybind11/arma_npy.h b/src/lasp/pybind11/arma_npy.h index 4fe8920..0fffd6b 100644 --- a/src/lasp/pybind11/arma_npy.h +++ b/src/lasp/pybind11/arma_npy.h @@ -9,6 +9,7 @@ namespace py = pybind11; template using pyarray = py::array_t; using dpyarray = pyarray; +using cpyarray = pyarray; /** @@ -127,6 +128,6 @@ arma::Mat NpyToCol(pyarray data) { if (data.ndim() != 1) { throw std::runtime_error("Expects a 1D array"); } - return arma::Col(data.mutable_data(0), data.size(), copy); + return arma::Col(data.mutable_data(0), data.size(), copy); } diff --git a/src/lasp/pybind11/lasp_dsp_pybind.cpp b/src/lasp/pybind11/lasp_dsp_pybind.cpp index 1e56ba8..e8cdfea 100644 --- a/src/lasp/pybind11/lasp_dsp_pybind.cpp +++ b/src/lasp/pybind11/lasp_dsp_pybind.cpp @@ -12,6 +12,7 @@ using std::cerr; using std::endl; namespace py = pybind11; +using rte = std::runtime_error; /** * \ingroup pybind @@ -29,11 +30,25 @@ void init_dsp(py::module &m) { py::class_ fft(m, "Fft"); fft.def(py::init()); - fft.def("fft", py::overload_cast(&Fft::fft)); - fft.def("fft", py::overload_cast(&Fft::fft)); + fft.def("fft", [](Fft &f, dpyarray dat) { + if (dat.ndim() == 1) { + return ColToNpy(f.fft(NpyToCol(dat))); + } else if (dat.ndim() == 2) { + return MatToNpy(f.fft(NpyToMat(dat))); + } else { + throw rte("Invalid dimensions of array"); + } + }); + fft.def("ifft", [](Fft &f, cpyarray dat) { + if (dat.ndim() == 1) { + return ColToNpy(f.ifft(NpyToCol(dat))); + } else if (dat.ndim() == 2) { + return MatToNpy(f.ifft(NpyToMat(dat))); + } else { + throw rte("Invalid dimensions of array"); + } + }); - fft.def("ifft", py::overload_cast(&Fft::ifft)); - fft.def("ifft", py::overload_cast(&Fft::ifft)); fft.def_static("load_fft_wisdom", &Fft::load_fft_wisdom); fft.def_static("store_fft_wisdom", &Fft::store_fft_wisdom); @@ -55,7 +70,7 @@ void init_dsp(py::module &m) { py::class_> sbq(m, "SeriesBiquad", filter); sbq.def(py::init([](dpyarray filter) { - return std::make_shared(NpyToCol(filter)); + return std::make_shared(NpyToCol(filter)); })); sbq.def("filter", [](SeriesBiquad &s, dpyarray input) { vd res = NpyToCol(input); @@ -89,10 +104,17 @@ void init_dsp(py::module &m) { py::arg("overlap_percentage") = 50.0, py::arg("time_constant") = -1); aps.def("compute", [](AvPowerSpectra &aps, dpyarray timedata) { - std::optional res = - aps.compute(NpyToMat(timedata)); + std::optional res; + { + py::gil_scoped_release release; + res = aps.compute(NpyToMat(timedata)); + } - return CubeToNpy(res.value_or(cube(0, 0, 0))); + return CubeToNpy(res.value_or(ccube(0, 0, 0))); + }); + aps.def("get_est", [](const AvPowerSpectra &ps) { + auto est = ps.get_est(); + return CubeToNpy(est.value_or(ccube(0, 0, 0))); }); py::class_ slm(m, "cppSLM"); diff --git a/src/lasp/pybind11/lasp_pyindatahandler.cpp b/src/lasp/pybind11/lasp_pyindatahandler.cpp index e607143..8f909c9 100644 --- a/src/lasp/pybind11/lasp_pyindatahandler.cpp +++ b/src/lasp/pybind11/lasp_pyindatahandler.cpp @@ -215,10 +215,10 @@ void init_datahandler(py::module &m) { ); rtaps.def("getCurrentValue", [](RtAps &rt) { - std::unique_ptr val = rt.getCurrentValue(); + std::unique_ptr val = rt.getCurrentValue(); if (val) { return CubeToNpy (*val); } - return CubeToNpy(cube(1,0,0)); + return CubeToNpy(ccube(1,0,0)); }); } diff --git a/test/test_aps.py b/test/test_aps.py index ed86a17..4ec248f 100644 --- a/test/test_aps.py +++ b/test/test_aps.py @@ -23,7 +23,7 @@ def test_aps1(): sig = np.random.randn(int(fs*tend)) - cps = aps.compute(sig) + cps = aps.compute(sig[:,None]) pow1 = np.sum(sig**2)/sig.size pow2 = np.sum((cps[:,0,0]).real) diff --git a/test/test_biquadbank.py b/test/test_biquadbank.py index a920610..be0a1da 100644 --- a/test/test_biquadbank.py +++ b/test/test_biquadbank.py @@ -27,7 +27,7 @@ out = bq.filter(in_) from scipy.signal import welch nfft = 48000 -freq, H1 = welch(out[:,0], +freq, H1 = welch(out, # scaling fs=fs,nfft=nfft) diff --git a/test/test_fft.py b/test/test_fft.py index 2f3b087..09e0b4e 100644 --- a/test/test_fft.py +++ b/test/test_fft.py @@ -20,7 +20,7 @@ def test_forward_fft(): sig = np.random.randn(nfft) fft_lasp = Fft(nfft) - res_lasp = fft_lasp.fft(sig)[:,0] + res_lasp = fft_lasp.fft(sig) res_npy = np.fft.rfft(sig) assert(np.isclose(np.linalg.norm(res_lasp- res_npy), 0)) @@ -38,11 +38,11 @@ def test_backward_fft(): Sig = Sigr + 1j*Sigi fft_lasp = Fft(nfft) - sig_lasp = fft_lasp.ifft(Sig)[:,0] + sig_lasp = fft_lasp.ifft(Sig) sig_py = np.fft.irfft(Sig) assert(np.isclose(np.linalg.norm(sig_py- sig_lasp), 0)) if __name__ == '__main__': test_forward_fft() - # test_backward_fft()() \ No newline at end of file + test_backward_fft() diff --git a/test/test_ps.py b/test/test_ps.py index d7aced2..b320da9 100644 --- a/test/test_ps.py +++ b/test/test_ps.py @@ -23,7 +23,7 @@ def test_ps(): # sig = 8*np.cos(omg*t) - cps = ps.compute(sig) + cps = ps.compute(sig[:,None]) pow1 = np.sum(sig**2)/sig.size pow2 = np.sum((cps[:,0,0]).real)