diff --git a/steps/Interpolate.cc b/steps/Interpolate.cc index 3029cb0dd91046704a34e07912fb76dca423fb07..44f456ccff1187693ccf97f7ad2f5c1afe6ad6df 100644 --- a/steps/Interpolate.cc +++ b/steps/Interpolate.cc @@ -11,14 +11,10 @@ #include "../common/ParameterSet.h" #include "../common/StringTools.h" -#include <casacore/casa/Arrays/ArrayMath.h> - #include <iostream> #include <iomanip> #include <thread> -using casacore::IPosition; - using dp3::base::DPBuffer; using dp3::base::DPInfo; using dp3::base::FlagCounter; @@ -60,11 +56,10 @@ void Interpolate::showTimings(std::ostream& os, double duration) const { os << " Interpolate " << _name << '\n'; } -bool Interpolate::process(const DPBuffer& buf) { +bool Interpolate::process(std::unique_ptr<base::DPBuffer> input_buffer) { _timer.start(); // Collect the data in buffers. - _buffers.emplace_back(); - _buffers.back().copy(buf); + _buffers.emplace_back(std::move(input_buffer)); // If we have a full window of data, interpolate everything // up to the middle of the window if (_buffers.size() >= _windowSize) { @@ -84,11 +79,12 @@ bool Interpolate::process(const DPBuffer& buf) { } void Interpolate::sendFrontBufferToNextStep() { - casacore::IPosition shp = _buffers.front().GetCasacoreData().shape(); - size_t nPol = shp[0], nChan = shp[1], nBl = shp[2], n = nPol * nChan * nBl; + auto front_buffer = std::move(_buffers.front()); + const auto& shp = front_buffer->GetData().shape(); + size_t nPol = shp[2], nChan = shp[1], nBl = shp[0], n = nPol * nChan * nBl; // Set all flags to false - bool* flags = _buffers.front().GetFlags().data(); - casacore::Complex* data = _buffers.front().GetData().data(); + bool* flags = front_buffer->GetFlags().data(); + std::complex<float>* data = front_buffer->GetData().data(); std::fill(flags, flags + n, false); // Flag NaN values (values for which the entire window was flagged on input) for (size_t i = 0; i != n; ++i) { @@ -101,7 +97,7 @@ void Interpolate::sendFrontBufferToNextStep() { } _timer.stop(); - getNextStep()->process(_buffers.front()); + getNextStep()->process(std::move(front_buffer)); _timer.start(); _buffers.pop_front(); @@ -128,9 +124,9 @@ void Interpolate::finish() { #define BUFFER_SIZE 1024 void Interpolate::interpolateTimestep(size_t index) { - const IPosition shp = _buffers.front().GetCasacoreData().shape(); - const size_t nPol = shp[0], nChan = shp[1], nPerBl = nPol * nChan, - nBl = shp[2]; + const auto& shp = _buffers.front()->GetData().shape(); + const size_t nPol = shp[2], nChan = shp[1], nPerBl = nPol * nChan, + nBl = shp[0]; std::vector<std::thread> threads; size_t nthreads = std::min<size_t>(getInfo().nThreads(), 8); @@ -141,7 +137,7 @@ void Interpolate::interpolateTimestep(size_t index) { threads.emplace_back(&Interpolate::interpolationThread, this); for (size_t bl = 0; bl < nBl; ++bl) { - bool* flags = _buffers[index].GetFlags().data() + bl * nPerBl; + bool* flags = _buffers[index]->GetFlags().data() + bl * nPerBl; for (size_t ch = 0; ch != nChan; ++ch) { for (size_t p = 0; p != nPol; ++p) { if (*flags) { @@ -167,8 +163,8 @@ void Interpolate::interpolationThread() { void Interpolate::interpolateSample(size_t timestep, size_t baseline, size_t channel, size_t pol) { - const casacore::IPosition shp = _buffers.front().GetCasacoreData().shape(); - const size_t nPol = shp[0], nChan = shp[1], + const auto& shp = _buffers.front()->GetData().shape(); + const size_t nPol = shp[2], nChan = shp[1], timestepBegin = (timestep > _windowSize / 2) ? (timestep - _windowSize / 2) : 0, @@ -183,9 +179,9 @@ void Interpolate::interpolateSample(size_t timestep, size_t baseline, float windowSum = 0.0; for (size_t t = timestepBegin; t != timestepEnd; ++t) { - casacore::Complex* data = _buffers[t].GetData().data() + - (baseline * nChan + channelBegin) * nPol + pol; - const bool* flags = _buffers[t].GetFlags().data() + + std::complex<float>* data = _buffers[t]->GetData().data() + + (baseline * nChan + channelBegin) * nPol + pol; + const bool* flags = _buffers[t]->GetFlags().data() + (baseline * nChan + channelBegin) * nPol + pol; const float* row = &_kernelLookup[_windowSize * (t + int(_windowSize / 2) - timestep)]; @@ -203,15 +199,15 @@ void Interpolate::interpolateSample(size_t timestep, size_t baseline, } // This write is multithreaded, but is allowed because this value is never // read from in the loops above (because flagged values are skipped). - casacore::Complex& value = + std::complex<float>& value = _buffers[timestep] - .GetCasacoreData() + ->GetData() .data()[(baseline * nChan + channel) * nPol + pol]; if (windowSum != 0.0) value = valueSum / windowSum; else - value = casacore::Complex(std::numeric_limits<float>::quiet_NaN(), - std::numeric_limits<float>::quiet_NaN()); + value = std::complex<float>(std::numeric_limits<float>::quiet_NaN(), + std::numeric_limits<float>::quiet_NaN()); } } // namespace steps diff --git a/steps/Interpolate.h b/steps/Interpolate.h index bfba897cbf09ccf62414b6db7e95cb01118990cb..d96f6908bc85cf60d8aa36a3d308c30bbcc8eff4 100644 --- a/steps/Interpolate.h +++ b/steps/Interpolate.h @@ -15,8 +15,6 @@ #include <aocommon/lane.h> -#include <casacore/casa/Arrays/Cube.h> - namespace dp3 { namespace steps { @@ -39,7 +37,7 @@ class Interpolate : public Step { /// Process the data. /// It keeps the data. /// When processed, it invokes the process function of the next step. - bool process(const base::DPBuffer&) override; + bool process(std::unique_ptr<base::DPBuffer> input_buffer) override; /// Finish the processing of this step and subsequent steps. void finish() override; @@ -72,7 +70,7 @@ class Interpolate : public Step { std::string _name; size_t _interpolatedPos; - std::deque<base::DPBuffer> _buffers; + std::deque<std::unique_ptr<base::DPBuffer>> _buffers; size_t _windowSize; common::NSTimer _timer; aocommon::Lane<Sample> _lane; diff --git a/steps/test/unit/tInterpolate.cc b/steps/test/unit/tInterpolate.cc index e283b14163c22d4d913a31ef5fe8ee10f2adc4c7..cc014dd25898a54764dfbffb170a67a7084155e2 100644 --- a/steps/test/unit/tInterpolate.cc +++ b/steps/test/unit/tInterpolate.cc @@ -12,16 +12,11 @@ #include "../../../common/ParameterSet.h" #include "../../../common/StringTools.h" -#include <casacore/casa/Arrays/ArrayMath.h> -#include <casacore/casa/Arrays/ArrayLogical.h> -#include <casacore/casa/BasicMath/Math.h> - #include <boost/test/unit_test.hpp> #include <string> #include <vector> -using dp3::base::DPBuffer; using dp3::base::DPInfo; using dp3::steps::Step; @@ -33,7 +28,8 @@ BOOST_AUTO_TEST_SUITE(interpolate) // It can be used with different nr of times, channels, etc. class TestInput : public dp3::steps::MockInput { public: - TestInput(int ntime, int nant, int nchan, int ncorr, bool flag) + TestInput(std::size_t ntime, std::size_t nant, std::size_t nchan, + std::size_t ncorr, bool flag) : itsCount(0), itsNTime(ntime), itsNBl(nant * (nant + 1) / 2), @@ -52,7 +48,7 @@ class TestInput : public dp3::steps::MockInput { std::vector<int> ant2(itsNBl); int st1 = 0; int st2 = 0; - for (int i = 0; i < itsNBl; ++i) { + for (std::size_t i = 0; i < itsNBl; ++i) { ant1[i] = st1; ant2[i] = st2; if (++st2 == 4) { @@ -87,35 +83,37 @@ class TestInput : public dp3::steps::MockInput { // Define the frequencies. std::vector<double> chanFreqs; std::vector<double> chanWidth(nchan, 100000.); - for (int i = 0; i < nchan; i++) { + for (std::size_t i = 0; i < nchan; i++) { chanFreqs.push_back(1050000. + i * 100000.); } info().setChannels(std::move(chanFreqs), std::move(chanWidth)); } private: - bool process(const DPBuffer&) override { + bool process(std::unique_ptr<dp3::base::DPBuffer>) override { // Stop when all times are done. if (itsCount == itsNTime) { return false; } - casacore::Cube<std::complex<float>> data(itsNCorr, itsNChan, itsNBl); - for (int i = 0; i < int(data.size()); ++i) { + const std::array<std::size_t, 3> data_shape{itsNBl, itsNChan, itsNCorr}; + xt::xtensor<std::complex<float>, 3> data(data_shape); + for (std::size_t i = 0; i < data.size(); ++i) { data.data()[i] = std::complex<float>(1.6, 0.9); } if (itsCount == 5) { data += std::complex<float>(10., 10.); } - DPBuffer buf; - buf.setTime(itsCount * 5 + 2); // same interval as in updateAveragInfo - buf.setData(data); - casacore::Cube<float> weights(data.shape()); - weights = 1.; - buf.setWeights(weights); - casacore::Cube<bool> flags(data.shape()); - flags = itsFlag; - buf.setFlags(flags); - getNextStep()->process(buf); + auto buffer = std::make_unique<dp3::base::DPBuffer>(); + buffer->setTime(itsCount * 5 + 2); // same interval as in updateAveragInfo + buffer->ResizeData(data_shape); + buffer->GetData() = data; + xt::xtensor<float, 3> weights(data_shape, 1.0); + buffer->ResizeWeights(data_shape); + buffer->GetWeights() = weights; + xt::xtensor<bool, 3> flags(data_shape, itsFlag); + buffer->ResizeFlags(data_shape); + buffer->GetFlags() = flags; + getNextStep()->process(std::move(buffer)); ++itsCount; return true; } @@ -125,14 +123,15 @@ class TestInput : public dp3::steps::MockInput { // Do nothing / keep the info set in the constructor. } - int itsCount, itsNTime, itsNBl, itsNChan, itsNCorr; + std::size_t itsCount, itsNTime, itsNBl, itsNChan, itsNCorr; bool itsFlag; }; // Class to check result. class TestOutput : public dp3::steps::test::ThrowStep { public: - TestOutput(int ntime, int nant, int nchan, int ncorr) + TestOutput(std::size_t ntime, std::size_t nant, std::size_t nchan, + std::size_t ncorr) : itsCount(0), itsNTime(ntime), itsNBl(nant * (nant + 1) / 2), @@ -140,38 +139,39 @@ class TestOutput : public dp3::steps::test::ThrowStep { itsNCorr(ncorr) {} private: - bool process(const DPBuffer& buf) override { + bool process( + const std::unique_ptr<dp3::base::DPBuffer> input_buffer) override { // Fill expected result in similar way as TestInput. - casacore::Cube<casacore::Complex> result(itsNCorr, itsNChan, itsNBl); + const std::array<std::size_t, 3> data_shape{itsNBl, itsNChan, itsNCorr}; + xt::xtensor<std::complex<float>, 3> result(data_shape); for (std::size_t i = 0; i < result.size(); ++i) { - result.data()[i] = casacore::Complex(1.6, 0.9); + result.data()[i] = std::complex<float>(1.6, 0.9); } if (itsCount == 5) { - result += casacore::Complex(10., 10.); + result += std::complex<float>(10.0, 10.0); } // Check the result. - BOOST_CHECK(casacore::allNear(buf.GetCasacoreData(), result, 1e-10)); - BOOST_CHECK(casacore::near(buf.getTime(), 2 + 5. * itsCount)); + BOOST_CHECK(xt::allclose(input_buffer->GetData(), result, 1e-10)); + BOOST_CHECK(xt::allclose(input_buffer->getTime(), 2 + 5.0 * itsCount)); ++itsCount; return true; } void finish() override {} void updateInfo(const DPInfo& info) override { - BOOST_CHECK_EQUAL(int(info.origNChan()), itsNChan); - BOOST_CHECK_EQUAL(int(info.nchan()), itsNChan); - BOOST_CHECK_EQUAL(int(info.ntime()), itsNTime); + BOOST_CHECK_EQUAL(info.origNChan(), itsNChan); + BOOST_CHECK_EQUAL(info.nchan(), itsNChan); + BOOST_CHECK_EQUAL(info.ntime(), itsNTime); BOOST_CHECK_EQUAL(info.firstTime(), 100.0); BOOST_CHECK_EQUAL(info.timeInterval(), 5.0); - BOOST_CHECK_EQUAL(int(info.nchanAvg()), 1); - BOOST_CHECK_EQUAL(int(info.ntimeAvg()), 1); - BOOST_CHECK_EQUAL(int(info.chanFreqs().size()), itsNChan); - BOOST_CHECK_EQUAL(int(info.chanWidths().size()), itsNChan); + BOOST_CHECK_EQUAL(info.nchanAvg(), 1); + BOOST_CHECK_EQUAL(info.ntimeAvg(), 1); + BOOST_CHECK_EQUAL(info.chanFreqs().size(), itsNChan); + BOOST_CHECK_EQUAL(info.chanWidths().size(), itsNChan); BOOST_CHECK(info.msName().empty()); } - int itsCount; - int itsNTime, itsNBl, itsNChan, itsNCorr; + std::size_t itsCount, itsNTime, itsNBl, itsNChan, itsNCorr; }; BOOST_AUTO_TEST_CASE(test1) {