diff --git a/CEP/DP3/DPPP/include/DPPP/BaselineSelection.h b/CEP/DP3/DPPP/include/DPPP/BaselineSelection.h
index de02f2b113341986d7fb8fa434015d534573bf6d..7f83ace2753098399f3f912ccfaf8c97aa3d1b0a 100644
--- a/CEP/DP3/DPPP/include/DPPP/BaselineSelection.h
+++ b/CEP/DP3/DPPP/include/DPPP/BaselineSelection.h
@@ -54,7 +54,8 @@ namespace LOFAR {
// <li> maxbl: maximum baseline length (in m); only if minmax=true
// </ul>
BaselineSelection (const ParSet&, const string& prefix,
- bool minmax=false);
+ bool minmax=false,
+ const string& defaultCorrType=string());
// Is there any selection?
bool hasSelection() const;
diff --git a/CEP/DP3/DPPP/include/DPPP/DPInfo.h b/CEP/DP3/DPPP/include/DPPP/DPInfo.h
index cc3efd17f87b93cda0f5c351907bdc55b097f998..f14b5ac0c51b86ecd35eda861b178636b4aaef6d 100644
--- a/CEP/DP3/DPPP/include/DPPP/DPInfo.h
+++ b/CEP/DP3/DPPP/include/DPPP/DPInfo.h
@@ -82,10 +82,6 @@ namespace LOFAR {
const casa::Vector<casa::Int>& ant1,
const casa::Vector<casa::Int>& ant2);
- // Set the info for the given baselines.
- void set (const casa::Vector<casa::Int>& ant1,
- const casa::Vector<casa::Int>& ant2);
-
// Update the info for the given average factors.
// If chanAvg is higher than the actual nr of channels, it is reset.
// The same is true for timeAvg.
@@ -153,6 +149,17 @@ namespace LOFAR {
double refFreq() const
{ return itsRefFreq; }
+ // Get the antenna numbers actually used in the (selected) baselines.
+ // E.g. [0,2,5,6]
+ const vector<int>& antennaUsed() const
+ { return itsAntUsed; }
+
+ // Get the indices of all antennae in the used antenna vector above.
+ // -1 means that the antenna is not used.
+ // E.g. [0,-1,1,-1,-1,2,3] for the example above.
+ const vector<int>& antennaMap() const
+ { return itsAntMap; }
+
// Are the visibility data needed?
bool needVisData() const
{ return itsNeedVisData; }
@@ -175,6 +182,10 @@ namespace LOFAR {
const vector<double>& getBaselineLengths() const;
private:
+ // Set which antennae are actually used.
+ void setAntUsed();
+
+ //# Data members.
bool itsNeedVisData; //# Are the visibility data needed?
bool itsNeedWrite; //# Does the last step need to write?
string itsMSName;
@@ -200,6 +211,8 @@ namespace LOFAR {
double itsRefFreq;
casa::Vector<casa::String> itsAntNames;
vector<casa::MPosition> itsAntPos;
+ vector<int> itsAntUsed;
+ vector<int> itsAntMap;
casa::Vector<casa::Int> itsAnt1; //# ant1 of all baselines
casa::Vector<casa::Int> itsAnt2; //# ant2 of all baselines
mutable vector<double> itsBLength; //# baseline lengths
diff --git a/CEP/DP3/DPPP/include/DPPP/DPStep.h b/CEP/DP3/DPPP/include/DPPP/DPStep.h
index 6c6b1d900bc1bda8f2f598cabd6c6d58201fbd7e..fc56fadc0fdec33917b9b981462b959f27906e8d 100644
--- a/CEP/DP3/DPPP/include/DPPP/DPStep.h
+++ b/CEP/DP3/DPPP/include/DPPP/DPStep.h
@@ -118,7 +118,7 @@ namespace LOFAR {
// The default implementation copies the info.
virtual void updateInfo (const DPInfo&);
- //# Data mamabers.
+ //# Data members.
DPStep::ShPtr itsNextStep;
DPInfo itsInfo;
};
diff --git a/CEP/DP3/DPPP/include/DPPP/Demixer.h b/CEP/DP3/DPPP/include/DPPP/Demixer.h
index 857a089bcf6ec00f100952487c55140fba1a77fc..6fc90536be860b4692c9fc5267affa75a7f3bbbf 100644
--- a/CEP/DP3/DPPP/include/DPPP/Demixer.h
+++ b/CEP/DP3/DPPP/include/DPPP/Demixer.h
@@ -32,6 +32,7 @@
#include <DPPP/DPBuffer.h>
#include <DPPP/Patch.h>
#include <DPPP/PhaseShift.h>
+#include <DPPP/Filter.h>
#include <casa/Arrays/Cube.h>
#include <casa/Quanta/Quantum.h>
@@ -102,6 +103,9 @@ namespace LOFAR {
uint nChanOut,
uint nChanAvg);
+ // Do the demixing.
+ void handleDemix();
+
// Deproject the sources without a model.
void deproject (casa::Array<casa::DComplex>& factors,
vector<MultiResultStep*> avgResults,
@@ -113,18 +117,27 @@ namespace LOFAR {
// Export the solutions to a ParmDB.
void dumpSolutions();
+ // Merge the data of the selected baselines from the subtract buffer
+ // into the full buffer.
+ void mergeSubtractResult();
+
//# Data members.
DPInput* itsInput;
string itsName;
string itsSkyName;
string itsInstrumentName;
+ BaselineSelection itsSelBL;
+ Filter itsFilter;
vector<PhaseShift*> itsPhaseShifts;
- //# Phase shift and average steps.
+ //# Phase shift and average steps for demix.
vector<DPStep::ShPtr> itsFirstSteps;
//# Result of phase shifting and averaging the directions of interest
//# at the demix resolution.
vector<MultiResultStep*> itsAvgResults;
+ DPStep::ShPtr itsAvgStepSubtr;
+ Filter* itsFilterSubtr;
//# Result of averaging the target at the subtract resolution.
+ MultiResultStep* itsAvgResultFull;
MultiResultStep* itsAvgResultSubtr;
//# Ignore target in demixing?
bool itsIgnoreTarget;
diff --git a/CEP/DP3/DPPP/include/DPPP/Filter.h b/CEP/DP3/DPPP/include/DPPP/Filter.h
index 9b0acb82bc5cba0fbb390e0862b5bf32d2c087bc..ee63f38274933208b913bcfcabb5ad47766a4634 100644
--- a/CEP/DP3/DPPP/include/DPPP/Filter.h
+++ b/CEP/DP3/DPPP/include/DPPP/Filter.h
@@ -130,6 +130,9 @@ namespace LOFAR {
// Parameters are obtained from the parset using the given prefix.
Filter (DPInput* input, const ParSet&, const string& prefix);
+ // Construct the object for the given MS and baseline selection.
+ Filter (DPInput* input, const BaselineSelection&);
+
virtual ~Filter();
// Process the next data chunk.
@@ -148,6 +151,18 @@ namespace LOFAR {
// Show the timings.
virtual void showTimings (std::ostream&, double duration) const;
+ // Does the filter step has an actual selection?
+ bool hasSelection() const
+ { return itsDoSelect; }
+
+ // Get the indices of the selected baselines.
+ const vector<uint>& getIndicesBL() const
+ { return itsSelBL; }
+
+ // Get the buffer.
+ const DPBuffer& getBuffer() const
+ { return itsBuf; }
+
private:
//# Data members.
DPInput* itsInput;
@@ -158,6 +173,7 @@ namespace LOFAR {
BaselineSelection itsBaselines;
uint itsStartChan;
vector<uint> itsSelBL; //# Index of baselines to select
+ bool itsDoSelect; //# Any selection?
NSTimer itsTimer;
};
diff --git a/CEP/DP3/DPPP/include/DPPP/MSReader.h b/CEP/DP3/DPPP/include/DPPP/MSReader.h
index b50f16d485885d2fe9754dbed183f8b4dfd3b84b..62ba12d5494b726e9871a30e6f068ed07aa6bc81 100644
--- a/CEP/DP3/DPPP/include/DPPP/MSReader.h
+++ b/CEP/DP3/DPPP/include/DPPP/MSReader.h
@@ -31,7 +31,6 @@
#include <DPPP/DPBuffer.h>
#include <DPPP/UVWCalculator.h>
#include <DPPP/FlagCounter.h>
-#include <DPPP/BaselineSelection.h>
#include <tables/Tables/TableIter.h>
#include <tables/Tables/RefRows.h>
#include <casa/Arrays/Slicer.h>
diff --git a/CEP/DP3/DPPP/src/BaselineSelection.cc b/CEP/DP3/DPPP/src/BaselineSelection.cc
index efa40f0baeaf4668585e1bf0721ab0e8f5201a03..d5d3f45f5c15a820b04775f493d6f066c8cff951 100644
--- a/CEP/DP3/DPPP/src/BaselineSelection.cc
+++ b/CEP/DP3/DPPP/src/BaselineSelection.cc
@@ -39,9 +39,10 @@ namespace LOFAR {
BaselineSelection::BaselineSelection (const ParSet& parset,
const string& prefix,
- bool minmax)
+ bool minmax,
+ const string& defaultCorrType)
: itsStrBL (parset.getString (prefix + "baseline", "")),
- itsCorrType (parset.getString (prefix + "corrtype", "")),
+ itsCorrType (parset.getString (prefix + "corrtype", defaultCorrType)),
itsRangeBL (parset.getDoubleVector (prefix + "blrange",
vector<double>()))
{
@@ -69,9 +70,10 @@ namespace LOFAR {
void BaselineSelection::show (ostream& os) const
{
- os << " baseline: " << itsStrBL << std::endl;
- os << " corrtype: " << itsCorrType << std::endl;
- os << " blrange: " << itsRangeBL << std::endl;
+ os << " Baseline selection:" << std::endl;
+ os << " baseline: " << itsStrBL << std::endl;
+ os << " corrtype: " << itsCorrType << std::endl;
+ os << " blrange: " << itsRangeBL << std::endl;
}
Matrix<bool> BaselineSelection::apply (const DPInfo& info) const
diff --git a/CEP/DP3/DPPP/src/DPInfo.cc b/CEP/DP3/DPPP/src/DPInfo.cc
index 43e7efa192d586a315c1ebd123ef07ed4cd965ad..3901be35654abb9bd7eae13c23d688d07156e9c5 100644
--- a/CEP/DP3/DPPP/src/DPInfo.cc
+++ b/CEP/DP3/DPPP/src/DPInfo.cc
@@ -111,13 +111,28 @@ namespace LOFAR {
itsAntPos = antPos;
itsAnt1.reference (ant1);
itsAnt2.reference (ant2);
+ // Set which antennae are used.
+ setAntUsed();
}
- void DPInfo::set (const Vector<Int>& ant1,
- const Vector<Int>& ant2)
+ void DPInfo::setAntUsed()
{
- itsAnt1.reference (ant1);
- itsAnt2.reference (ant2);
+ itsAntUsed.clear();
+ itsAntMap.resize (itsAntNames.size());
+ std::fill (itsAntMap.begin(), itsAntMap.end(), -1);
+ for (uint i=0; i<itsAnt1.size(); ++i) {
+ ASSERT (itsAnt1[i] >= 0 && itsAnt1[i] < int(itsAntMap.size()) &&
+ itsAnt2[i] >= 0 && itsAnt2[i] < int(itsAntMap.size()));
+ itsAntMap[itsAnt1[i]] = 0;
+ itsAntMap[itsAnt2[i]] = 0;
+ }
+ itsAntUsed.reserve (itsAntNames.size());
+ for (uint i=0; i<itsAntMap.size(); ++i) {
+ if (itsAntMap[i] == 0) {
+ itsAntMap[i] = itsAntUsed.size();
+ itsAntUsed.push_back (i);
+ }
+ }
}
uint DPInfo::update (uint chanAvg, uint timeAvg)
@@ -171,8 +186,8 @@ namespace LOFAR {
Vector<Int> ant1 (baselines.size());
Vector<Int> ant2 (baselines.size());
for (uint i=0; i<baselines.size(); ++i) {
- ant1 = itsAnt1[i];
- ant2 = itsAnt2[i];
+ ant1[i] = itsAnt1[baselines[i]];
+ ant2[i] = itsAnt2[baselines[i]];
}
itsAnt1.reference (ant1);
itsAnt2.reference (ant2);
@@ -180,6 +195,7 @@ namespace LOFAR {
itsBLength.resize (0);
itsAutoCorrIndex.resize (0);
}
+ setAntUsed();
}
const vector<double>& DPInfo::getBaselineLengths() const
diff --git a/CEP/DP3/DPPP/src/DPStep.cc b/CEP/DP3/DPPP/src/DPStep.cc
index 19317fe0374f6e71db140e946749c8971afb65eb..f04b5d47e0034a6a9360eed7642e2a39699e1767 100644
--- a/CEP/DP3/DPPP/src/DPStep.cc
+++ b/CEP/DP3/DPPP/src/DPStep.cc
@@ -78,6 +78,7 @@ namespace LOFAR {
bool ResultStep::process (const DPBuffer& buf)
{
itsBuffer = buf;
+ getNextStep()->process (buf);
return true;
}
@@ -100,6 +101,7 @@ namespace LOFAR {
bool MultiResultStep::process (const DPBuffer& buf)
{
itsBuffers.push_back (buf);
+ getNextStep()->process (buf);
return true;
}
diff --git a/CEP/DP3/DPPP/src/Demixer.cc b/CEP/DP3/DPPP/src/Demixer.cc
index cc04916e0ae1764713f2070a16da643515c749f4..390664d89184392a0ceae60b6ca7380768e1ba74 100644
--- a/CEP/DP3/DPPP/src/Demixer.cc
+++ b/CEP/DP3/DPPP/src/Demixer.cc
@@ -75,7 +75,9 @@ namespace LOFAR {
itsName (prefix),
itsSkyName (parset.getString(prefix+"skymodel", "sky")),
itsInstrumentName (parset.getString(prefix+"instrumentmodel",
- "instrument")),
+ "instrument")),
+ itsSelBL (parset, prefix, false, "cross"),
+ itsFilter (input, itsSelBL),
itsAvgResultSubtr (0),
itsIgnoreTarget (parset.getBool (prefix+"ignoretarget", false)),
itsTargetSource (parset.getString(prefix+"targetsource", string())),
@@ -137,6 +139,9 @@ namespace LOFAR {
"An empty name is given for the sky and/or instrument model");
ASSERTSTR (!itsIgnoreTarget || itsTargetSource.empty(),
"Target source name cannot be given if ignoretarget=true");
+ // Add a null step as last step in the filter.
+ DPStep::ShPtr nullStep(new NullStep());
+ itsFilter.setNextStep (nullStep);
// Default nr of time chunks is maximum number of threads.
if (itsNTimeChunk == 0) {
itsNTimeChunk = OpenMP::maxThreads();
@@ -178,18 +183,26 @@ namespace LOFAR {
itsFactors.resize (itsNTimeChunk);
itsFactorsSubtr.resize (itsNTimeChunkSubtr);
itsPhaseShifts.reserve (itsNDir-1); // not needed for target direction
- itsFirstSteps.reserve (itsNDir+1); // one extra for itsAvgSubtr
+ itsFirstSteps.reserve (itsNDir);
itsAvgResults.reserve (itsNDir);
- // Create the steps for the sources to be removed.
- // Demixing consists of the following steps:
- // - phaseshift data to each demix source
- // - average data in each direction, also for original phasecenter.
- // - determine demix factors for all directions
- // - use BBS to predict and solve in each direction. It is possible to
- // predict more directions than to solve (for strong sources in field).
- // - use BBS to subtract the solved sources using the demix factors.
- // The averaging used here can be smaller than used when solving.
+ // Create the solve and subtract steps for the sources to be removed.
+ // Solving consists of the following steps:
+ // - select the requested baselines (longer baselines may need no demix)
+ // - phaseshift selected data to each demix source
+ // - average selected data in each direction, also original phasecenter.
+ // - determine and average demix factors for all directions
+ // - predict and solve in each direction. It is possible to predict
+ // more directions than to solve (for strong sources in field).
+ // Subtract consists of the following steps:
+ // - average all data (possibly different averaging than used in solve)
+ // - determine and average demix factors (using select output in solve)
+ // - select the requested baselines
+ // - subtract sources for selected data
+ // - merge subtract result into averaged data. This is not needed if
+ // no selection is done.
+ // Note that multiple time chunks are handled jointly, so a
+ // MultiResultStep is used to catch the results of all time chunks.
for (uint i=0; i<itsNDir-1; ++i) {
// First make the phaseshift and average steps for each demix source.
// The resultstep gets the result.
@@ -224,12 +237,18 @@ namespace LOFAR {
itsAvgResults.push_back (targetAvgRes);
// Create the data average step for the subtract.
- DPStep::ShPtr targetAvgSubtr(new Averager(input, prefix,
- itsNChanAvgSubtr,
- itsNTimeAvgSubtr));
+ // The entire average result is needed for the next NDPPP step.
+ // Only the selected baselines need to be subtracted, so add a
+ // filter step as the last one.
+ itsAvgStepSubtr = DPStep::ShPtr(new Averager(input, prefix,
+ itsNChanAvgSubtr,
+ itsNTimeAvgSubtr));
+ itsAvgResultFull = new MultiResultStep(itsNTimeChunkSubtr);
+ itsFilterSubtr = new Filter(input, itsSelBL);
itsAvgResultSubtr = new MultiResultStep(itsNTimeChunkSubtr);
- targetAvgSubtr->setNextStep (DPStep::ShPtr(itsAvgResultSubtr));
- itsFirstSteps.push_back (targetAvgSubtr);
+ itsAvgStepSubtr->setNextStep (DPStep::ShPtr(itsAvgResultFull));
+ itsAvgResultFull->setNextStep (DPStep::ShPtr(itsFilterSubtr));
+ itsFilterSubtr->setNextStep (DPStep::ShPtr(itsAvgResultSubtr));
// while(itsCutOffs.size() < itsNModel) {
// itsCutOffs.push_back(0.0);
@@ -242,37 +261,43 @@ namespace LOFAR {
info() = infoIn;
// Get size info.
- itsNStation = infoIn.antennaNames().size();
- itsNChanIn = infoIn.nchan();
- itsNBl = infoIn.nbaselines();
- itsNCorr = infoIn.ncorr();
+ itsNChanIn = infoIn.nchan();
+ itsNCorr = infoIn.ncorr();
ASSERTSTR (itsNCorr==4, "Demixing requires data with 4 polarizations");
+
+ // Handle possible data selection.
+ itsFilter.setInfo (infoIn);
+ const DPInfo& infoSel = itsFilter.getInfo();
+ itsNBl = infoSel.nbaselines();
itsFactorBuf.resize (IPosition(4, itsNCorr, itsNChanIn, itsNBl,
itsNDir*(itsNDir-1)/2));
itsFactorBufSubtr.resize (IPosition(4, itsNCorr, itsNChanIn, itsNBl,
itsNDir*(itsNDir-1)/2));
-
+ for (uint i=0; i<itsNBl; ++i) {
+ itsBaselines.push_back (Baseline(infoSel.getAnt1()[i],
+ infoSel.getAnt2()[i]));
+ }
+ // Get nr of stations actually used.
+ ///itsNStation = infoSel.antennaUsed().size();
+ itsNStation = infoSel.antennaNames().size();
// Adapt averaging to available nr of channels and times.
// Use a copy of the DPInfo, otherwise it is updated multiple times.
- DPInfo infoDemix(infoIn);
- itsNTimeAvg = std::min (itsNTimeAvg, infoIn.ntime());
+ DPInfo infoDemix(infoSel);
+ itsNTimeAvg = std::min (itsNTimeAvg, infoSel.ntime());
itsNChanAvg = infoDemix.update (itsNChanAvg, itsNTimeAvg);
itsNChanOut = infoDemix.nchan();
itsTimeIntervalAvg = infoDemix.timeInterval();
itsNTimeDemix = infoDemix.ntime();
- for (size_t i=0; i<infoIn.getAnt1().size(); ++i) {
- itsBaselines.push_back (Baseline(infoIn.getAnt1()[i],
- infoIn.getAnt2()[i]));
- }
// Let the internal steps update their data.
for (uint i=0; i<itsFirstSteps.size(); ++i) {
- itsFirstSteps[i]->setInfo (infoIn);
+ itsFirstSteps[i]->setInfo (infoSel);
}
+ itsAvgStepSubtr->setInfo (infoIn);
// Update the info of this object.
info().setNeedVisData();
info().setNeedWrite();
- itsNTimeAvgSubtr = std::min (itsNTimeAvgSubtr, infoIn.ntime());
+ itsNTimeAvgSubtr = std::min (itsNTimeAvgSubtr, infoSel.ntime());
itsNChanAvgSubtr = info().update (itsNChanAvgSubtr, itsNTimeAvgSubtr);
itsNChanOutSubtr = info().nchan();
ASSERTSTR (itsNChanAvg % itsNChanAvgSubtr == 0,
@@ -317,6 +342,14 @@ namespace LOFAR {
os << "Demixer " << itsName << std::endl;
os << " skymodel: " << itsSkyName << std::endl;
os << " instrumentmodel: " << itsInstrumentName << std::endl;
+ itsSelBL.show (os);
+ if (itsSelBL.hasSelection()) {
+ os << " demixing " << itsFilter.getInfo().nbaselines()
+ << " out of " << getInfo().nbaselines() << " baselines ("
+ << itsFilter.getInfo().antennaUsed().size()
+ << " out of " << getInfo().antennaUsed().size()
+ << " stations)" << std::endl;
+ }
os << " targetsource: " << itsTargetSource << std::endl;
os << " subtractsources: " << itsSubtrSources << std::endl;
os << " modelsources: " << itsModelSources << std::endl;
@@ -388,16 +421,22 @@ namespace LOFAR {
itsTimer));
}
- // Do the initial steps (phaseshift and average).
+ // Do the filter step first.
+ itsFilter.process (newBuf);
+ const DPBuffer& selBuf = itsFilter.getBuffer();
+ // Do the next steps (phaseshift and average) on the filter output.
itsTimerPhaseShift.start();
for (int i=0; i<int(itsFirstSteps.size()); ++i) {
- itsFirstSteps[i]->process(newBuf);
+ itsFirstSteps[i]->process(selBuf);
}
+ // Do the average and filter step for the output for all data.
+ itsAvgStepSubtr->process (newBuf);
itsTimerPhaseShift.stop();
- // For each itsNTimeAvg times, calculate the phase rotation per direction.
+ // For each itsNTimeAvg times, calculate the phase rotation per direction
+ // for the selected data.
itsTimerDemix.start();
- addFactors (newBuf, itsFactorBuf);
+ addFactors (selBuf, itsFactorBuf);
if (itsNTimeIn % itsNTimeAvg == 0) {
makeFactors (itsFactorBuf, itsFactors[itsNTimeOut],
itsAvgResults[0]->get()[itsNTimeOut].getWeights(),
@@ -409,8 +448,8 @@ namespace LOFAR {
itsNTimeOut++;
}
// Subtract is done with different averaging parameters, so calculate the
- // factors for it.
- addFactors (newBuf, itsFactorBufSubtr);
+ // factors for it (again for selected data only).
+ addFactors (selBuf, itsFactorBufSubtr);
if (itsNTimeIn % itsNTimeAvgSubtr == 0) {
makeFactors (itsFactorBufSubtr, itsFactorsSubtr[itsNTimeOutSubtr],
itsAvgResultSubtr->get()[itsNTimeOutSubtr].getWeights(),
@@ -424,36 +463,49 @@ namespace LOFAR {
// Estimate gains and subtract source contributions when sufficient time
// slots have been collected.
if (itsNTimeOut == itsNTimeChunk) {
- if(itsNModel > 0) {
- itsTimerSolve.start();
- demix();
- itsTimerSolve.stop();
- }
-
- // Clear the input buffers.
- for (size_t i=0; i<itsAvgResults.size(); ++i) {
- itsAvgResults[i]->clear();
- }
-
- // Let the next step process the data.
- for (uint i=0; i<itsNTimeOutSubtr; ++i) {
- itsTimer.stop();
- getNextStep()->process (itsAvgResultSubtr->get()[i]);
- itsTimer.start();
- }
+ handleDemix();
+ }
+ itsTimer.stop();
+ return true;
+ }
- // Clear the output buffer.
- itsAvgResultSubtr->clear();
+ void Demixer::handleDemix()
+ {
+ if(itsNModel > 0) {
+ itsTimerSolve.start();
+ demix();
+ itsTimerSolve.stop();
+ // If selection was done, merge the subtract results back into the
+ // buffer.
+ if (itsSelBL.hasSelection()) {
+ mergeSubtractResult();
+ }
+ }
- // Reset counters.
- itsNTimeIn = 0;
- itsNTimeOut = 0;
- itsNTimeOutSubtr = 0;
- itsTimeIndex += itsNTimeChunk;
+ // Clear the input buffers.
+ for (size_t i=0; i<itsAvgResults.size(); ++i) {
+ itsAvgResults[i]->clear();
+ }
+ // Let the next step process the data.
+ for (uint i=0; i<itsNTimeOutSubtr; ++i) {
+ itsTimer.stop();
+ if (itsSelBL.hasSelection()) {
+ getNextStep()->process (itsAvgResultFull->get()[i]);
+ } else {
+ getNextStep()->process (itsAvgResultSubtr->get()[i]);
+ }
+ itsTimer.start();
}
- itsTimer.stop();
- return true;
+ // Clear the output buffer.
+ itsAvgResultFull->clear();
+ itsAvgResultSubtr->clear();
+
+ // Reset counters.
+ itsNTimeIn = 0;
+ itsNTimeOut = 0;
+ itsNTimeOutSubtr = 0;
+ itsTimeIndex += itsNTimeChunk;
}
void Demixer::finish()
@@ -492,26 +544,7 @@ namespace LOFAR {
itsFactorsSubtr.resize(itsNTimeOutSubtr);
// Demix the source directions.
- if(itsNModel > 0) {
- itsTimerSolve.start();
- demix();
- itsTimerSolve.stop();
- }
-
- // Clear the input buffers.
- for (size_t i=0; i<itsAvgResults.size(); ++i) {
- itsAvgResults[i]->clear();
- }
-
- // Let the next step process the data.
- for (uint i=0; i<itsNTimeOutSubtr; ++i) {
- itsTimer.stop();
- getNextStep()->process (itsAvgResultSubtr->get()[i]);
- itsTimer.start();
- }
-
- // Clear the output buffer.
- itsAvgResultSubtr->clear();
+ handleDemix();
}
// Write solutions to disk in ParmDB format.
@@ -525,6 +558,22 @@ namespace LOFAR {
getNextStep()->finish();
}
+ void Demixer::mergeSubtractResult()
+ {
+ // Merge the selected baselines from the subtract buffer into the
+ // full buffer. Do it for all timestamps.
+ for (uint i=0; i<itsNTimeOutSubtr; ++i) {
+ const Array<Complex>& arr = itsAvgResultSubtr->get()[i].getData();
+ size_t nr = arr.shape()[0] * arr.shape()[1];
+ const Complex* in = arr.data();
+ Complex* out = itsAvgResultFull->get()[i].getData().data();
+ for (size_t j=0; j<itsFilter.getIndicesBL().size(); ++j) {
+ size_t inx = itsFilter.getIndicesBL()[j];
+ memcpy (out+inx*nr, in+j*nr, nr*sizeof(Complex));
+ }
+ }
+ }
+
void Demixer::addFactors (const DPBuffer& newBuf,
Array<DComplex>& factorBuf)
{
@@ -775,15 +824,24 @@ namespace LOFAR {
const size_t nCr = 4;
const size_t nSamples = nBl * nCh * nCr;
+ ///for (uint j=0; j<nTimeSubtr; ++j) {
+ ///cout << itsAvgResultSubtr->get()[j].getData()<<endl;
+ /// cout << itsAvgResultSubtr->get()[j].getFlags()<<endl;
+ /// cout << itsAvgResultSubtr->get()[j].getWeights()<<endl;
+ /// cout << itsAvgResultSubtr->get()[j].getUVW()<<endl;
+ ///}
+ ///for (uint i=0; i<itsFactors.size(); ++i) {
+ /// cout << itsFactors[i] << endl;
+ ///}
+ ///for (uint i=0; i<itsFactorsSubtr.size(); ++i) {
+ /// cout << itsFactorsSubtr[i] << endl;
+ ///}
+
vector<ThreadPrivateStorage> threadStorage(nThread);
for(vector<ThreadPrivateStorage>::iterator it = threadStorage.begin(),
end = threadStorage.end(); it != end; ++it)
{
initThreadPrivateStorage(*it, nDr, nSt, nBl, nCh, nChSubtr);
-
- // Copy solutions from global solution array to thread private solution
- // array (solution propagation between chunks).
- copy(itsLastKnowns.begin(), itsLastKnowns.end(), it->unknowns.begin());
}
const_cursor<double> cr_freq = casa_const_cursor(itsFreqDemix);
@@ -795,6 +853,11 @@ namespace LOFAR {
{
const size_t thread = OpenMP::threadNum();
ThreadPrivateStorage &storage = threadStorage[thread];
+ // Copy solutions from global solution array to thread private solution
+ // array (solution propagation between chunks).
+ if (!itsPropagateSolutions || ts==0) {
+ copy(itsLastKnowns.begin(), itsLastKnowns.end(), storage.unknowns.begin());
+ }
// Simulate.
//
@@ -978,42 +1041,30 @@ namespace LOFAR {
resolution[1] = itsTimeIntervalAvg;
parmDB.setDefaultSteps(resolution);
- // Convert station names from casa::String to std::string.
- ASSERT(getInfo().antennaNames().size() == itsNStation);
- vector<string> stations(itsNStation);
- copy(getInfo().antennaNames().begin(), getInfo().antennaNames().end(),
- stations.begin());
-
vector<BBS::Parm> parms;
for(size_t dr = 0; dr < itsNModel; ++dr) {
- for(size_t st = 0; st < itsNStation; ++st) {
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:0:0:Real:" + stations[st] + ":"
- + itsAllSources[dr])));
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:0:0:Imag:" + stations[st] + ":"
- + itsAllSources[dr])));
-
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:0:1:Real:" + stations[st] + ":"
- + itsAllSources[dr])));
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:0:1:Imag:" + stations[st] + ":"
- + itsAllSources[dr])));
-
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:1:0:Real:" + stations[st] + ":"
- + itsAllSources[dr])));
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:1:0:Imag:" + stations[st] + ":"
- + itsAllSources[dr])));
-
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:1:1:Real:" + stations[st] + ":"
- + itsAllSources[dr])));
- parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
- "DirectionalGain:1:1:Imag:" + stations[st] + ":"
- + itsAllSources[dr])));
+ for(size_t st = 0; st < getInfo().antennaNames().size(); ++st) {
+ string suffix (string(getInfo().antennaNames()[st]) + ':' +
+ itsAllSources[dr]);
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:0:0:Real:" + suffix)));
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:0:0:Imag:" + suffix)));
+
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:0:1:Real:" + suffix)));
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:0:1:Imag:" + suffix)));
+
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:1:0:Real:" + suffix)));
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:1:0:Imag:" + suffix)));
+
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:1:1:Real:" + suffix)));
+ parms.push_back(BBS::Parm(parmCache, parmSet.addParm(parmDB,
+ "DirectionalGain:1:1:Imag:" + suffix)));
}
}
diff --git a/CEP/DP3/DPPP/src/Filter.cc b/CEP/DP3/DPPP/src/Filter.cc
index 1021a26fd94e4efc086c936da8e56653cd1ee9f5..9ab9a6274772e99ddc9750aea4493f02ffdc8812 100644
--- a/CEP/DP3/DPPP/src/Filter.cc
+++ b/CEP/DP3/DPPP/src/Filter.cc
@@ -43,7 +43,16 @@ namespace LOFAR {
itsName (prefix),
itsStartChanStr (parset.getString(prefix+"startchan", "0")),
itsNrChanStr (parset.getString(prefix+"nchan", "0")),
- itsBaselines (parset, prefix)
+ itsBaselines (parset, prefix),
+ itsDoSelect (false)
+ {}
+
+ Filter::Filter (DPInput* input, const BaselineSelection& baselines)
+ : itsInput (input),
+ itsStartChanStr ("0"),
+ itsNrChanStr ("0"),
+ itsBaselines (baselines),
+ itsDoSelect (false)
{}
Filter::~Filter()
@@ -76,6 +85,7 @@ namespace LOFAR {
} else {
nrChan = std::min (nrChan, maxNrChan);
}
+ itsDoSelect = itsStartChan>0 || nrChan<maxNrChan;
// Handle possible baseline selection.
if (itsBaselines.hasSelection()) {
Matrix<bool> selbl(itsBaselines.apply (infoIn));
@@ -87,14 +97,19 @@ namespace LOFAR {
itsSelBL.push_back (i);
}
}
+ if (itsSelBL.size() < ant1.size()) {
+ itsDoSelect = true;
+ }
+ }
+ if (itsDoSelect) {
+ // Update the DPInfo object.
+ info().update (itsStartChan, nrChan, itsSelBL);
+ // Shape the arrays in the buffer.
+ IPosition shape (3, infoIn.ncorr(), nrChan, getInfo().nbaselines());
+ itsBuf.getData().resize (shape);
+ itsBuf.getFlags().resize (shape);
+ itsBuf.getWeights().resize (shape);
}
- // Update the DPInfo object.
- info().update (itsStartChan, nrChan, itsSelBL);
- // Shape the arrays in the buffer.
- IPosition shape (3, infoIn.ncorr(), nrChan, getInfo().nbaselines());
- itsBuf.getData().resize (shape);
- itsBuf.getFlags().resize (shape);
- itsBuf.getWeights().resize (shape);
}
void Filter::show (std::ostream& os) const
@@ -117,6 +132,17 @@ namespace LOFAR {
bool Filter::process (const DPBuffer& buf)
{
itsTimer.start();
+ if (!itsDoSelect) {
+ itsBuf = buf; // uses reference semantics
+ itsTimer.stop();
+ getNextStep()->process (itsBuf);
+ return true;
+ }
+ // Make sure no other object references the DATA and UVW arrays.
+ itsBuf.getData().unique();
+ itsBuf.getFlags().unique();
+ itsBuf.getWeights().unique();
+ itsBuf.getFullResFlags().unique();
// Get the various data arrays.
RefRows rowNrs(buf.getRowNrs());
const Array<Complex>& data = buf.getData();
@@ -154,6 +180,7 @@ namespace LOFAR {
} else {
// Copy the data of the selected baselines and channels.
itsBuf.getUVW().resize (IPosition(2, 3, getInfo().nbaselines()));
+ itsBuf.getUVW().unique();
Complex* toData = itsBuf.getData().data();
Bool* toFlag = itsBuf.getFlags().data();
Float* toWeight = itsBuf.getWeights().data();