diff --git a/steps/StationAdder.cc b/steps/StationAdder.cc
index 1793a45e6046833e40b959e2589ae96b3a050d17..9a1cce166791f8dd9402c07d1e01518734703226 100644
--- a/steps/StationAdder.cc
+++ b/steps/StationAdder.cc
@@ -28,7 +28,6 @@
#include "../base/FlagCounter.h"
using casacore::ArrayColumn;
-using casacore::IPosition;
using casacore::MPosition;
using casacore::MVPosition;
using casacore::Regex;
@@ -142,7 +141,7 @@ void StationAdder::updateInfo(const DPInfo& infoIn) {
}
// The superstation position is the average of its parts.
newNames.push_back(iter->first);
- newPosition *= 1. / parts.size();
+ newPosition *= 1.0 / parts.size();
newPoss.push_back(MPosition(newPosition, MPosition::ITRF));
itsParts.push_back(casacore::Vector<int>(parts));
// Set the diameter of the new station by determining the
@@ -154,7 +153,7 @@ void StationAdder::updateInfo(const DPInfo& infoIn) {
newPosition -
MPosition::Convert(antennaPos[inx], MPosition::ITRF)().getValue();
const casacore::Vector<double>& diff = mvdiff.getValue();
- double dist = std::sqrt(std::accumulate(diff.cbegin(), diff.cend(), 0.,
+ double dist = std::sqrt(std::accumulate(diff.cbegin(), diff.cend(), 0.0,
casacore::SumSqr<double>()));
// Add the radius of the station used.
maxdist = std::max(maxdist, dist + 0.5 * antennaDiam[inx]);
@@ -244,13 +243,6 @@ void StationAdder::updateInfo(const DPInfo& infoIn) {
// Setup the UVW calculator (for new baselines).
itsUVWCalc = std::make_unique<base::UVWCalculator>(
infoIn.phaseCenter(), infoIn.arrayPos(), antennaPos);
- // Size the buffer to cater for the new baselines.
- const std::array<std::size_t, 3> shape{getInfo().nbaselines(),
- getInfo().nchan(), getInfo().ncorr()};
- itsBuf.ResizeData(shape);
- itsBuf.ResizeFlags(shape);
- itsBuf.ResizeWeights(shape);
- itsBuf.ResizeUvw(getInfo().nbaselines());
}
void StationAdder::show(std::ostream& os) const {
@@ -279,38 +271,51 @@ void StationAdder::showTimings(std::ostream& os, double duration) const {
os << " StationAdder " << itsName << '\n';
}
-bool StationAdder::process(const DPBuffer& buf) {
+bool StationAdder::process(std::unique_ptr<base::DPBuffer> buffer) {
itsTimer.start();
- // Get the various data arrays.
- const casacore::Array<casacore::Complex>& data = buf.GetCasacoreData();
- const casacore::Array<bool>& flags = buf.GetCasacoreFlags();
- const casacore::Array<float>& weights = buf.GetCasacoreWeights();
- const casacore::Array<double>& uvws = buf.GetCasacoreUvw();
- // Copy the data; only the first baselines will be filled.
- std::copy(data.data(), data.data() + data.size(), itsBuf.GetData().data());
- std::copy(flags.data(), flags.data() + flags.size(),
- itsBuf.GetFlags().data());
- std::copy(weights.data(), weights.data() + weights.size(),
- itsBuf.GetWeights().data());
- std::copy(uvws.data(), uvws.data() + uvws.size(), itsBuf.GetUvw().data());
+
+ // Resize of data buffers can be destructive, therefore we must:
+ // 1. Store the data locally.
+ // TODO(AST-1330): Avoid copying the data, e.g., using 'swap'.
+ const xt::xtensor<std::complex<float>, 3> input_data = buffer->GetData();
+ const xt::xtensor<bool, 3> input_flags = buffer->GetFlags();
+ const xt::xtensor<float, 3> input_weights = buffer->GetWeights();
+ const xt::xtensor<double, 2> input_uvws = buffer->GetUvw();
+ // 2. Resize the buffer to the new sizes that were set in our info object
+ const std::array<std::size_t, 3> new_shape{
+ getInfo().nbaselines(), getInfo().nchan(), getInfo().ncorr()};
+ buffer->ResizeData(new_shape);
+ buffer->ResizeFlags(new_shape);
+ buffer->ResizeWeights(new_shape);
+ buffer->ResizeUvw(getInfo().nbaselines());
+ // 3. Copy the data back into the resized buffer; only the existing baselines
+ // at the start will be filled the additional new baselines at the end will be
+ // empty.
+ std::copy_n(input_data.data(), input_data.size(), buffer->GetData().data());
+ std::copy_n(input_flags.data(), input_flags.size(),
+ buffer->GetFlags().data());
+ std::copy_n(input_weights.data(), input_weights.size(),
+ buffer->GetWeights().data());
+ std::copy_n(input_uvws.data(), input_uvws.size(), buffer->GetUvw().data());
+
// Now calculate the data pointers of the new baselines.
- unsigned int nrOldBL = data.shape()[2];
- unsigned int nrcc = data.shape()[0] * data.shape()[1];
- casacore::Complex* dataPtr = itsBuf.GetData().data() + data.size();
- bool* flagPtr = itsBuf.GetFlags().data() + data.size();
- float* wghtPtr = itsBuf.GetWeights().data() + data.size();
- double* uvwPtr = itsBuf.GetUvw().data() + uvws.size();
+ unsigned int nrOldBL = input_data.shape()[0];
+ unsigned int nrcc = input_data.shape()[2] * input_data.shape()[1];
+ std::complex<float>* dataPtr = buffer->GetData().data() + input_data.size();
+ bool* flagPtr = buffer->GetFlags().data() + input_data.size();
+ float* wghtPtr = buffer->GetWeights().data() + input_data.size();
+ double* uvwPtr = buffer->GetUvw().data() + input_uvws.size();
std::vector<unsigned int> npoints(nrcc);
- std::vector<casacore::Complex> dataFlg(nrcc);
+ std::vector<std::complex<float>> dataFlg(nrcc);
std::vector<float> wghtFlg(nrcc);
// Loop over all new baselines.
for (unsigned int i = 0; i < itsBufRows.size(); ++i) {
// Clear the data for the new baseline.
for (unsigned int k = 0; k < nrcc; ++k) {
- dataPtr[k] = casacore::Complex();
+ dataPtr[k] = std::complex<float>();
wghtPtr[k] = 0.;
npoints[k] = 0;
- dataFlg[k] = casacore::Complex();
+ dataFlg[k] = std::complex<float>();
wghtFlg[k] = 0.;
}
@@ -331,11 +336,11 @@ bool StationAdder::process(const DPBuffer& buf) {
}
blnr--; // decrement because blnr+1 is stored in itsBufRows
// Get pointers to the input baseline data.
- const casacore::Complex* inDataPtr =
- (itsBuf.GetData().data() + blnr * nrcc);
- const bool* inFlagPtr = (itsBuf.GetFlags().data() + blnr * nrcc);
- const float* inWghtPtr = (itsBuf.GetWeights().data() + blnr * nrcc);
- const double* inUvwPtr = (itsBuf.GetUvw().data() + blnr * 3);
+ const std::complex<float>* inDataPtr = &buffer->GetData()(blnr, 0, 0);
+ const bool* inFlagPtr = &buffer->GetFlags()(blnr, 0, 0);
+ const float* inWghtPtr = &buffer->GetWeights()(blnr, 0, 0);
+ const double* inUvwPtr = &buffer->GetUvw()(blnr, 0);
+
// Add the data, uvw, and weights if not flagged.
// Write 4 loops to avoid having to test inside the loop.
// Count the flagged points separately, so it can be used
@@ -428,8 +433,9 @@ bool StationAdder::process(const DPBuffer& buf) {
}
} else {
unsigned int blnr = nrOldBL + i;
- const std::array<double, 3> uvws = itsUVWCalc->getUVW(
- getInfo().getAnt1()[blnr], getInfo().getAnt2()[blnr], buf.getTime());
+ const std::array<double, 3> uvws =
+ itsUVWCalc->getUVW(getInfo().getAnt1()[blnr],
+ getInfo().getAnt2()[blnr], buffer->getTime());
uvwPtr[0] = uvws[0];
uvwPtr[1] = uvws[1];
uvwPtr[2] = uvws[2];
@@ -440,10 +446,8 @@ bool StationAdder::process(const DPBuffer& buf) {
wghtPtr += nrcc;
uvwPtr += 3;
}
- itsBuf.setTime(buf.getTime());
- itsBuf.setExposure(buf.getExposure());
itsTimer.stop();
- getNextStep()->process(itsBuf);
+ getNextStep()->process(std::move(buffer));
return true;
}
@@ -486,7 +490,7 @@ void StationAdder::addToMS(const string& msName) {
stat = statCol(0);
}
}
- casacore::Vector<double> offset(3, 0.);
+ casacore::Vector<double> offset(3, 0.0);
// Put the data for each new antenna.
for (unsigned int i = origNant; i < getInfo().antennaNames().size(); ++i) {
antTab.addRow();
diff --git a/steps/StationAdder.h b/steps/StationAdder.h
index c12f1bde2d27bc6afd0bf8bc9ab0c6eb03407326..4d78f8968c99c8337c95c5b03b553075093a8408 100644
--- a/steps/StationAdder.h
+++ b/steps/StationAdder.h
@@ -61,7 +61,7 @@ class StationAdder : 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> buffer) override;
/// Finish the processing of this step and subsequent steps.
void finish() override;
@@ -94,7 +94,6 @@ class StationAdder : public Step {
casacore::Table& antTab);
std::string itsName;
- base::DPBuffer itsBuf;
common::ParameterRecord itsStatRec; ///< stations definitions
std::vector<casacore::Vector<int>>
itsParts; ///< the stations in each superstation
diff --git a/steps/test/unit/tStationAdder.cc b/steps/test/unit/tStationAdder.cc
index bfb530bd859fc7b9a0cb28dd3d63a5f7a3b2790e..50478405a04d9d06af371acd263310fe7b39a1bf 100644
--- a/steps/test/unit/tStationAdder.cc
+++ b/steps/test/unit/tStationAdder.cc
@@ -8,6 +8,10 @@
#include <casacore/casa/Arrays/ArrayLogical.h>
#include <casacore/casa/BasicMath/Math.h>
+#include <xtensor/xcomplex.hpp>
+#include <xtensor/xio.hpp>
+#include <xtensor/xview.hpp>
+
#include <boost/test/unit_test.hpp>
#include <boost/test/data/test_case.hpp>
@@ -33,7 +37,8 @@ BOOST_AUTO_TEST_SUITE(stationadder)
// It can be used with different nr of times, channels, etc.
class TestInput : public dp3::steps::MockInput {
public:
- TestInput(int ntime, int nbl, int nchan, int ncorr)
+ TestInput(std::size_t ntime, std::size_t nbl, std::size_t nchan,
+ std::size_t ncorr)
: itsCount(0),
itsNTime(ntime),
itsNBl(nbl),
@@ -47,7 +52,7 @@ class TestInput : public dp3::steps::MockInput {
vector<int> ant2(nbl);
int st1 = 0;
int st2 = 0;
- for (int i = 0; i < nbl; ++i) {
+ for (std::size_t i = 0; i < nbl; ++i) {
ant1[i] = st1;
ant2[i] = st2;
if (++st2 == 4) {
@@ -85,44 +90,47 @@ class TestInput : public dp3::steps::MockInput {
antPos[3] = casacore::MPosition(
casacore::Quantum<casacore::Vector<double>>(vals, "m"),
casacore::MPosition::ITRF);
- vector<double> antDiam(4, 70.);
+ vector<double> antDiam(4, 70.0);
info().setAntennas(antNames, antDiam, antPos, ant1, ant2);
// Define the frequencies.
- std::vector<double> chanWidth(nchan, 1000000.);
+ std::vector<double> chanWidth(nchan, 1000000.0);
std::vector<double> chanFreqs;
- for (int i = 0; i < nchan; i++) {
- chanFreqs.push_back(10500000. + i * 1000000.);
+ for (std::size_t i = 0; i < nchan; i++) {
+ chanFreqs.push_back(10500000.0 + i * 1000000.0);
}
info().setChannels(std::move(chanFreqs), std::move(chanWidth));
}
private:
- bool process(const DPBuffer&) override {
+ bool process(std::unique_ptr<DPBuffer>) override {
// Stop when all times are done.
if (itsCount == itsNTime) {
return false;
}
- casacore::Cube<casacore::Complex> data(itsNCorr, itsNChan, itsNBl);
+
+ const std::array<std::size_t, 3> data_shape{itsNBl, itsNChan, itsNCorr};
+
+ xt::xtensor<std::complex<float>, 3> data(data_shape);
for (int i = 0; i < int(data.size()); ++i) {
- data.data()[i] =
- casacore::Complex(i + itsCount * 10, i - 10 + itsCount * 6);
+ data[i] = std::complex<float>(i + itsCount * 10.0f,
+ i - 10.0f + itsCount * 6.0f);
}
casacore::Cube<float> weights(itsNCorr, itsNChan, itsNBl);
indgen(weights, 0.5f, 0.01f);
casacore::Matrix<double> uvw(3, itsNBl);
- for (int i = 0; i < itsNBl; ++i) {
+ for (std::size_t i = 0; i < itsNBl; ++i) {
uvw(0, i) = 1 + itsCount + i;
uvw(1, i) = 2 + itsCount + i;
uvw(2, i) = 3 + itsCount + i;
}
DPBuffer buf;
buf.setTime(itsCount * 30 + 4472025740.0);
- buf.setData(data);
+ buf.ResizeData(data_shape);
+ buf.GetData() = data;
buf.setWeights(weights);
buf.setUVW(uvw);
- casacore::Cube<bool> flags(data.shape());
- flags = false;
- buf.setFlags(flags);
+ buf.ResizeFlags(data_shape);
+ buf.GetFlags().fill(false);
getNextStep()->process(buf);
++itsCount;
return true;
@@ -133,13 +141,14 @@ 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;
};
// Class to check result of TestInput run by test1.
class TestOutput : public dp3::steps::test::ThrowStep {
public:
- TestOutput(int ntime, int nbl, int nchan, int ncorr, bool sumauto)
+ TestOutput(int ntime, std::size_t nbl, std::size_t nchan, std::size_t ncorr,
+ bool sumauto)
: itsCount(0),
itsNTime(ntime),
itsNBl(nbl),
@@ -148,27 +157,27 @@ class TestOutput : public dp3::steps::test::ThrowStep {
itsSumAuto(sumauto) {}
private:
- void addData(casacore::Cube<casacore::Complex>& to,
- const casacore::Cube<casacore::Complex>& from, int bl) {
- to += from(casacore::IPosition(3, 0, 0, bl),
- casacore::IPosition(3, to.nrow() - 1, to.ncolumn() - 1, bl));
+ void addData(xt::xtensor<std::complex<float>, 2>& to,
+ const xt::xtensor<std::complex<float>, 3>& from, int bl) {
+ to += xt::view(from, bl, xt::all(), xt::all());
}
- void addConjData(casacore::Cube<casacore::Complex>& to,
- const casacore::Cube<casacore::Complex>& from, int bl) {
- to +=
- conj(from(casacore::IPosition(3, 0, 0, bl),
- casacore::IPosition(3, to.nrow() - 1, to.ncolumn() - 1, bl)));
+ void addConjData(xt::xtensor<std::complex<float>, 2>& to,
+ const xt::xtensor<std::complex<float>, 3>& from, int bl) {
+ to += xt::conj(xt::view(from, bl, xt::all(), xt::all()));
}
- bool process(const DPBuffer& buf) override {
- casacore::Cube<casacore::Complex> data(itsNCorr, itsNChan, itsNBl);
+ bool process(std::unique_ptr<DPBuffer> buffer) override {
+ const std::array<std::size_t, 3> shape{itsNBl, itsNChan, itsNCorr};
+ xt::xtensor<std::complex<float>, 3> data(shape);
+
for (int i = 0; i < int(data.size()); ++i) {
- data.data()[i] =
- casacore::Complex(i + itsCount * 10, i - 10 + itsCount * 6);
+ data[i] = std::complex<float>(i + itsCount * 10, i - 10 + itsCount * 6);
}
casacore::Cube<float> weights(itsNCorr, itsNChan, itsNBl);
indgen(weights, 0.5f, 0.01f);
- casacore::Cube<casacore::Complex> databl0(itsNCorr, itsNChan, 1);
- casacore::Cube<casacore::Complex> databl1(itsNCorr, itsNChan, 1);
+ xt::xtensor<std::complex<float>, 2> databl0(
+ std::array<std::size_t, 2>{itsNChan, itsNCorr});
+ xt::xtensor<std::complex<float>, 2> databl1(
+ std::array<std::size_t, 2>{itsNChan, itsNCorr});
// "{ns:[rs01.s01, rs02.s01, cs01.s02]}" was given resulting in 2 new
// baselines (ns-ns and cs01.s01-ns).
// Thus adding the baselines below.
@@ -195,39 +204,43 @@ class TestOutput : public dp3::steps::test::ThrowStep {
addConjData(databl1, data, 2);
addConjData(databl1, data, 6);
addConjData(databl1, data, 14);
- casacore::Matrix<double> uvw(3, itsNBl);
- for (int i = 0; i < itsNBl; ++i) {
- uvw(0, i) = 1 + itsCount + i;
- uvw(1, i) = 2 + itsCount + i;
- uvw(2, i) = 3 + itsCount + i;
+ xt::xtensor<double, 2> uvw(std::array<std::size_t, 2>{itsNBl, 3});
+ for (std::size_t i = 0; i < itsNBl; ++i) {
+ uvw(i, 0) = 1 + itsCount + i;
+ uvw(i, 1) = 2 + itsCount + i;
+ uvw(i, 2) = 3 + itsCount + i;
}
- casacore::IPosition end(3, itsNCorr - 1, itsNChan - 1, itsNBl - 1);
- BOOST_CHECK(
- allEQ(buf.GetCasacoreData()(casacore::IPosition(3, 0), end), data));
- BOOST_CHECK_EQUAL(buf.GetCasacoreFlags().shape(),
- casacore::IPosition(3, itsNCorr, itsNChan, itsNBl + 2));
- BOOST_CHECK(allEQ(buf.GetCasacoreFlags(), false));
- BOOST_CHECK(allEQ(buf.GetCasacoreWeights()(casacore::IPosition(3, 0), end),
- weights));
- BOOST_CHECK(
- allEQ(buf.GetCasacoreUvw()(casacore::IPosition(2, 0),
- casacore::IPosition(2, 2, itsNBl - 1)),
- uvw));
+
+ auto xt_weights = aocommon::xt::CreateSpan(
+ weights.data(),
+ std::array{static_cast<std::size_t>(weights.shape()[2]),
+ static_cast<std::size_t>(weights.shape()[1]),
+ static_cast<std::size_t>(weights.shape()[0])});
+
+ BOOST_TEST(data == xt::view(buffer->GetData(), xt::range(0, itsNBl),
+ xt::all(), xt::all()));
+ BOOST_TEST(buffer->GetFlags().shape() ==
+ (std::array<std::size_t, 3>{itsNBl + 2, itsNChan, itsNCorr}),
+ boost::test_tools::per_element());
+ BOOST_CHECK(!xt::any(buffer->GetFlags()));
+ BOOST_TEST(xt_weights == xt::view(buffer->GetWeights(),
+ xt::range(0, itsNBl), xt::all(),
+ xt::all()));
+ BOOST_TEST(uvw ==
+ xt::view(buffer->GetUvw(), xt::range(0, itsNBl), xt::all()));
+
// Now check data of new baselines.
- end[2] = itsNBl;
- BOOST_CHECK(allNear(
- buf.GetCasacoreData()(casacore::IPosition(3, 0, 0, itsNBl), end),
- databl0 / weight, 1e-5));
- BOOST_CHECK(allNear(
- buf.GetCasacoreWeights()(casacore::IPosition(3, 0, 0, itsNBl), end),
- weight, 1e-5));
- end[2] = itsNBl + 1;
- BOOST_CHECK(allNear(
- buf.GetCasacoreData()(casacore::IPosition(3, 0, 0, itsNBl + 1), end),
- databl1 / 6.f, 1e-5));
- BOOST_CHECK(allNear(
- buf.GetCasacoreWeights()(casacore::IPosition(3, 0, 0, itsNBl + 1), end),
- 6.f, 1e-5));
+ BOOST_CHECK(
+ xt::allclose(xt::view(buffer->GetData(), itsNBl, xt::all(), xt::all()),
+ (databl0 / weight)));
+ BOOST_CHECK(xt::allclose(
+ xt::view(buffer->GetWeights(), itsNBl, xt::all(), xt::all()), weight));
+ BOOST_CHECK(xt::allclose(
+ xt::view(buffer->GetData(), itsNBl + 1, xt::all(), xt::all()),
+ (databl1 / 6.0f)));
+ BOOST_CHECK(xt::allclose(
+ xt::view(buffer->GetWeights(), itsNBl + 1, xt::all(), xt::all()),
+ 6.0f));
itsCount++;
return true;
}
@@ -247,9 +260,11 @@ class TestOutput : public dp3::steps::test::ThrowStep {
BOOST_CHECK_EQUAL(int(infoIn.antennaPos().size()), 5);
BOOST_CHECK_EQUAL(infoIn.antennaNames()[4], "ns");
casacore::Vector<double> pos1(infoIn.antennaPos()[4].getValue().getValue());
- BOOST_CHECK(casacore::near(pos1[0], (3828763. + 3828746. + 3828713.) / 3));
- BOOST_CHECK(casacore::near(pos1[1], (442449. + 442592. + 442878.) / 3));
- BOOST_CHECK(casacore::near(pos1[2], (5064923. + 5064924. + 5064926.) / 3));
+ BOOST_CHECK(
+ casacore::near(pos1[0], (3828763.0 + 3828746.0 + 3828713.0) / 3));
+ BOOST_CHECK(casacore::near(pos1[1], (442449.0 + 442592.0 + 442878.0) / 3));
+ BOOST_CHECK(
+ casacore::near(pos1[2], (5064923.0 + 5064924.0 + 5064926.0) / 3));
// Check diam.
double d1 = sqrt((pos1[0] - 3828763) * (pos1[0] - 3828763) +
(pos1[1] - 442449) * (pos1[1] - 442449) +
@@ -265,7 +280,7 @@ class TestOutput : public dp3::steps::test::ThrowStep {
}
int itsCount;
- int itsNTime, itsNBl, itsNChan, itsNCorr;
+ std::size_t itsNTime, itsNBl, itsNChan, itsNCorr;
bool itsSumAuto;
};