[casacore] r21578 committed - Merged in CASA changes

[casa...@googlecode.com]

Revision: 21578
Author: gervandiepen
Date: Wed Mar 18 15:01:43 2015 UTC
Log: Merged in CASA changes

https://code.google.com/p/casacore/source/detail?r=21578

Modified:
/branches/nov14/casa/BasicSL/Complex.h
/branches/nov14/casa/Containers/ValueHolderRep.cc
/branches/nov14/casa/Logging/LogOrigin.h
/branches/nov14/casa/Quanta/QMath.h
/branches/nov14/casa/Utilities/CountedPtr.h
/branches/nov14/changescripts/updateinclude
/branches/nov14/images/Images/ImageStatistics.h
/branches/nov14/images/Images/ImageStatistics.tcc
/branches/nov14/images/Images/test/dImageStatistics.cc
/branches/nov14/images/Images/test/tImageStatistics.cc
/branches/nov14/lattices/LatticeMath/LatticeHistograms.h
/branches/nov14/lattices/LatticeMath/LatticeStatistics.h
/branches/nov14/lattices/LatticeMath/LatticeStatistics.tcc
/branches/nov14/lattices/LatticeMath/LatticeStatsBase.h
/branches/nov14/ms/MeasurementSets/MSKeys.h
/branches/nov14/ms/MeasurementSets/MSMetaData.cc
/branches/nov14/ms/MeasurementSets/MSMetaData.h
/branches/nov14/ms/MeasurementSets/MSSummary.cc
/branches/nov14/ms/MeasurementSets/test/tMSMetaData.cc

=======================================
--- /branches/nov14/casa/BasicSL/Complex.h Wed Dec 10 08:06:42 2014 UTC
+++ /branches/nov14/casa/BasicSL/Complex.h Wed Mar 18 15:01:43 2015 UTC
@@ -344,7 +344,7 @@
namespace std {
inline float conj(float x) { return x; }
inline double conj(double x) { return x; }
-#ifndef AIPS_CXX11
+#if !(defined(AIPS_CXX11) || (defined(__APPLE_CC__) && __APPLE_CC__ >
5621))
inline float real(float x) { return x; }
inline double real(double x) { return x; }
inline float imag(float ) { return 0; }
=======================================
--- /branches/nov14/casa/Containers/ValueHolderRep.cc Wed Dec 10 08:06:42
2014 UTC
+++ /branches/nov14/casa/Containers/ValueHolderRep.cc Wed Mar 18 15:01:43
2015 UTC
@@ -1047,6 +1047,9 @@
case TpArrayInt:
rec.define (id, *static_cast<Array<Int>*>(itsPtr));
break;
+ case TpArrayUInt:
+ rec.define (id, *static_cast<Array<uInt>*>(itsPtr));
+ break;
case TpArrayInt64:
rec.define (id, *static_cast<Array<Int64>*>(itsPtr));
break;
=======================================
--- /branches/nov14/casa/Logging/LogOrigin.h Tue Mar 3 12:22:11 2015 UTC
+++ /branches/nov14/casa/Logging/LogOrigin.h Wed Mar 18 15:01:43 2015 UTC
@@ -156,8 +156,8 @@
// <src>where</src> will be defined with the <src>WHERE</src> macro.
LogOrigin &sourceLocation(const SourceLocation *where);

- // Returns <src>class::function</src> for a member function, or
- // <src>::function</src> for a global function.
+ // Returns <src>class\::function</src> for a member function, or
+ // <src>\::function</src> for a global function.
String fullName() const;

// Turn the entire origin into a String.
=======================================
--- /branches/nov14/casa/Quanta/QMath.h Wed Dec 10 08:06:42 2014 UTC
+++ /branches/nov14/casa/Quanta/QMath.h Wed Mar 18 15:01:43 2015 UTC
@@ -229,8 +229,10 @@
// <group name="foreign">
Int ceil(const Int &val);
Int floor(const Int &val);
+#if !(defined(AIPS_CXX11) || (defined(__APPLE_CC__) && __APPLE_CC__ >
5621))
Float real(const Float &val);
-Double real(const Double&val);
+Double real(const Double &val);
+#endif
Array<Complex> operator *(const Array<Complex> &in, Double f);
Array<Complex> operator /(const Array<Complex> &in, Double f);
Array<DComplex> operator *(const Array<DComplex> &in, Double f);
=======================================
--- /branches/nov14/casa/Utilities/CountedPtr.h Thu Mar 12 10:32:41 2015 UTC
+++ /branches/nov14/casa/Utilities/CountedPtr.h Wed Mar 18 15:01:43 2015 UTC
@@ -30,12 +30,18 @@

#include <casacore/casa/aips.h>

-#if defined AIPS_CXX11
-# include <memory>
-# define SHARED_PTR std::shared_ptr
+#if (defined(AIPS_CXX11) || (defined(__APPLE_CC__) && __APPLE_CC__ > 5621))
+#include <memory>
+#define SHARED_PTR std::shared_ptr
+#define DYNAMIC_POINTER_CAST std::dynamic_pointer_cast
+#define CONST_POINTER_CAST std::const_pointer_cast
+#define STATIC_POINTER_CAST std::static_pointer_cast
#else
-# include <tr1/memory>
-# define SHARED_PTR std::tr1::shared_ptr
+#include <tr1/memory>
+#define SHARED_PTR std::tr1::shared_ptr
+#define DYNAMIC_POINTER_CAST std::tr1::dynamic_pointer_cast
+#define CONST_POINTER_CAST std::tr1::const_pointer_cast
+#define STATIC_POINTER_CAST std::tr1::static_pointer_cast
#endif

namespace casacore { //#Begin casa namespace
@@ -79,12 +85,6 @@
template<class t>
class CountedPtr
{
-#ifdef AIPS_CXX11
- typedef std::shared_ptr<t> PointerRep;
-#else
- typedef std::tr1::shared_ptr<t> PointerRep;
-#endif
-

protected:
// Helper class to make deletion of object optional.
@@ -99,6 +99,9 @@


public:
+
+
+
// This constructor allows for the creation of a null
// <src>CountedPtr</src>. The assignment operator can be used
// to assign a null <src>CountedPtr</src> from another
@@ -120,17 +123,17 @@
: pointerRep_p (val, Deleter<t> (delit))
{}

- // Create from a shared_ptr.
- CountedPtr (const PointerRep& rep)
- : pointerRep_p (rep)
- {}
-
// This copy constructor allows <src>CountedPtr</src>s to be
// initialized from other <src>CountedPtr</src>s for which the pointer
TP*
// is convertible to T*.
template<typename TP>
CountedPtr(const CountedPtr<TP>& that)
- : pointerRep_p (that.pointerRep_p)
+ : pointerRep_p(that.pointerRep_p)
+ {}
+
+ // Create from a shared_ptr.
+ CountedPtr (const SHARED_PTR<t>& rep)
+ : pointerRep_p (rep)
{}

// This destructor only deletes the really stored data when it was
@@ -220,25 +223,13 @@
// <group>
template<typename U>
CountedPtr<U> static_ptr_cast() const
-#ifdef AIPS_CXX11
- { return CountedPtr<U> (std::static_pointer_cast<U> (pointerRep_p));
}
-#else
- { return CountedPtr<U> (std::tr1::static_pointer_cast<U>
(pointerRep_p)); }
-#endif
+ { return CountedPtr<U> (STATIC_POINTER_CAST<U> (pointerRep_p)); }
template<typename U>
CountedPtr<U> const_ptr_cast() const
-#ifdef AIPS_CXX11
- { return CountedPtr<U> (std::const_pointer_cast<U> (pointerRep_p)); }
-#else
- { return CountedPtr<U> (std::tr1::const_pointer_cast<U>
(pointerRep_p)); }
-#endif
+ { return CountedPtr<U> (CONST_POINTER_CAST<U> (pointerRep_p)); }
template<typename U>
CountedPtr<U> dynamic_ptr_cast() const
-#ifdef AIPS_CXX11
- { return CountedPtr<U> (std::dynamic_pointer_cast<U>
(pointerRep_p)); }
-#else
- { return CountedPtr<U> (std::tr1::dynamic_pointer_cast<U>
(pointerRep_p)); }
-#endif
+ { return CountedPtr<U> (DYNAMIC_POINTER_CAST<U> (pointerRep_p)); }
// </group>

// Sometimes it is useful to know if there is more than one
@@ -260,6 +251,7 @@
// Make all types of CountedPtr a friend for the templated operator=.
template<typename TP> friend class CountedPtr;

+ typedef SHARED_PTR<t> PointerRep;

PointerRep pointerRep_p;
};
=======================================
--- /branches/nov14/changescripts/updateinclude Wed Jan 28 09:59:04 2015 UTC
+++ /branches/nov14/changescripts/updateinclude Wed Mar 18 15:01:43 2015 UTC
@@ -13,9 +13,11 @@
sed -f $pgmpath/latticemath.sed "$FILE.n2" > "$FILE.n3"
sed -f $pgmpath/lel.sed "$FILE.n3" > "$FILE.n4"
sed -f $pgmpath/regions.sed "$FILE.n4" > "$FILE.n5"
- if ! diff "$FILE.n5" "$FILE" >& /dev/null; then
- mv "$FILE.n5" "$FILE"
+ sed -f $pgmpath/mssel.sed "$FILE.n5" > "$FILE.n6"
+ sed -f $pgmpath/msoper.sed "$FILE.n6" > "$FILE.n7"
+ if ! diff "$FILE.n7" "$FILE" >& /dev/null; then
+ mv "$FILE.n7" "$FILE"
# echo " changed #includes in $FILE"
fi
- rm -f "$FILE".n[12345]
+ rm -f "$FILE".n[1234567]
done
=======================================
--- /branches/nov14/images/Images/ImageStatistics.h Wed Feb 11 14:00:31
2015 UTC
+++ /branches/nov14/images/Images/ImageStatistics.h Wed Mar 18 15:01:43
2015 UTC
@@ -162,6 +162,16 @@
AccumType sum, Double beamAreaInPixels
) const;

+ virtual Bool _computeFlux(
+ Array<AccumType>& flux, const Array<AccumType>& npts,
+ const Array<AccumType>& sum
+ );
+
+ virtual Bool _computeFlux(
+ AccumType& flux, AccumType sum, const IPosition& pos,
+ Bool posInLattice
+ );
+
private:
// Data

@@ -181,7 +191,7 @@
// Get beam area in pixels if possible. Return False if the beam area
could not be
// calculated.
virtual Bool _getBeamArea(
- Array<Double>& beamArea
+ Array<Double>& beamArea, String& msg
) const;

// List min and max with world coordinates
@@ -202,6 +212,7 @@
);


+
//# Make members of parent class known.
protected:
using LatticeStatistics<T>::locInLattice;
=======================================
--- /branches/nov14/images/Images/ImageStatistics.tcc Mon Feb 16 07:05:15
2015 UTC
+++ /branches/nov14/images/Images/ImageStatistics.tcc Wed Mar 18 15:01:43
2015 UTC
@@ -155,7 +155,7 @@
}

template <class T> Bool ImageStatistics<T>::_getBeamArea(
- Array<Double>& beamArea
+ Array<Double>& beamArea, String& msg
) const {
ImageInfo ii = pInImage_p->imageInfo();
Bool hasMultiBeams = ii.hasMultipleBeams();
@@ -164,19 +164,19 @@

// use contains() not == so moment maps are dealt with nicely
if (! hasMultiBeams && ! hasSingleBeam ) {
- _messages.push_back(String("Image has no beam, cannot compute flux"));
+ msg = "Image has no beam";
return False;
}
else if (! cSys.hasDirectionCoordinate()) {
- _messages.push_back(String("Image does not have a direction coordinate,
cannot computer flux"));
+ msg = "Image does not have a direction coordinate";
return False;
}
else {
String imageUnits = pInImage_p->units().getName();
imageUnits.upcase();
if (! imageUnits.contains("JY/BEAM")) {
- _messages.push_back(String("Image brightness units not conformant with
Jy/beam, cannot compute flux"));
- return False;
+ msg = "Image brightness units not conformant with Jy/beam";
+ return False;
}
}
DirectionCoordinate dCoord = cSys.directionCoordinate();
@@ -199,14 +199,20 @@
);
return True;
}
+
// per plane beams
- // ensure both the spectral and polarization axes are display axes
- Bool foundSpec = !cSys.hasSpectralAxis() || False;
- Bool foundPol = !cSys.hasPolarizationCoordinate() || False;
+ // ensure both the spectral and polarization axes are display axes,
+ // a degenerate axis is considered not to be a cursor axis for
+ // this purpose since no aggregation along that axis actually occurs
+ IPosition shape = pInImage_p->shape();
+
+ Bool foundSpec = ! cSys.hasSpectralAxis()
+ || shape[cSys.spectralAxisNumber(False)] == 1;
+ Bool foundPol = ! cSys.hasPolarizationCoordinate()
+ || shape[cSys.polarizationAxisNumber(False)] == 1;
Int specAxis = foundSpec ? -1 : cSys.spectralAxisNumber();
Int polAxis = foundPol ? -1 : cSys.polarizationAxisNumber();
Bool found = False;
- const ImageBeamSet& beams = ii.getBeamSet();
Int storageSpecAxis = -1;
Int storagePolAxis = -1;

@@ -227,24 +233,29 @@
if (! found) {
// if per-plane beams, either the spectral axis and/or the
// polarization axis is not a display axis
- // or else the image has no beam
+ msg = "One or both of the spectral or polarization axes is "
+ "not a display axis, not degenerate, and the image has multiple beams";
return False;
}
-
+ const ImageBeamSet& beams = ii.getBeamSet();
IPosition beamsShape = beams.shape();
+ /*
if (cSys.hasSpectralAxis()) {
AlwaysAssert(
beamsShape[0] == beamAreaShape[storageSpecAxis],
AipsError
);
}
+ */
Int beamPolAxis = -1;
if (cSys.hasPolarizationCoordinate()) {
beamPolAxis = specAxis < 0 ? 0 : 1;
+ /*
AlwaysAssert(
beamsShape[beamPolAxis] == beamAreaShape[storagePolAxis],
AipsError
);
+ */
}
IPosition curPos(beamAreaShape.nelements(), 0);
GaussianBeam curBeam;
@@ -468,16 +479,24 @@
ostringstream oss;
if (_canDoFlux()) {
Array<Double> beamArea;
- Bool hasBeam = _getBeamArea(beamArea);
- Quantum<AccumType> qFlux = _flux(
- sum, hasBeam ? *(beamArea.begin()) : 0
- );
- AccumType val = qFlux.getValue();
- String unit = qFlux.getFullUnit().getName();
- oss << " -- flux density [flux]: "
+ String msg;
+ Bool hasBeam = _getBeamArea(beamArea, msg);
+ String unit;
+ if (! hasBeam) {
+ unit = pInImage_p->units().getName();
+ unit.downcase();
+ }
+ if (hasBeam || ! unit.contains("/beam")) {
+ Quantum<AccumType> qFlux = _flux(
+ sum, hasBeam ? *(beamArea.begin()) : 0
+ );
+ AccumType val = qFlux.getValue();
+ String unit = qFlux.getFullUnit().getName();
+ oss << " -- flux density [flux]: "
<< val << " " << unit;
- messages.push_back(oss.str());
- oss.str("");
+ messages.push_back(oss.str());
+ oss.str("");
+ }
}

if (LattStatsSpecialize::hasSomePoints(nPts)) {
@@ -609,6 +628,70 @@
}
return Quantum<AccumType>(flux, fUnit);
}
+
+template <class T> Bool ImageStatistics<T>::_computeFlux(
+ Array<AccumType>& flux,const Array<AccumType>& npts, const
Array<AccumType>& sum
+) {
+ Array<Double> beamArea;
+ String msg;
+ Bool gotBeamArea = _getBeamArea(beamArea, msg);
+ if (! gotBeamArea) {
+ String unit = pInImage_p->units().getName();
+ unit.downcase();
+ if (unit.contains("/beam")) {
+ os_p << LogIO::WARN << "Unable to compute flux density: "
+ << msg << LogIO::POST;
+ return False;
+ }
+ }
+ ReadOnlyVectorIterator<AccumType> sumIt(sum);
+ ReadOnlyVectorIterator<AccumType> nPtsIt(npts);
+ VectorIterator<AccumType> fluxIt(flux);
+ PtrHolder<ReadOnlyVectorIterator<Double> > beamAreaIter(
+ gotBeamArea ? new ReadOnlyVectorIterator<Double>(beamArea) : 0
+ );
+ uInt n1 = nPtsIt.vector().nelements();
+ while (!nPtsIt.pastEnd()) {
+ for (uInt i=0; i<n1; ++i) {
+ if (nPtsIt.vector()(i) > 0.5) {
+ fluxIt.vector()(i) = _flux(
+ sumIt.vector()(i), gotBeamArea ? beamAreaIter->vector()(i) : 0
+ ).getValue();
+ }
+ }
+ nPtsIt.next();
+ sumIt.next();
+ fluxIt.next();
+ if (gotBeamArea) {
+ beamAreaIter->next();
+ }
+ }
+ return True;
+}
+
+template <class T> Bool ImageStatistics<T>::_computeFlux(
+ AccumType& flux, AccumType sum, const IPosition& pos,
+ Bool posInLattice
+) {
+ Array<Double> beamArea;
+ String msg;
+ if (_getBeamArea(beamArea, msg)) {
+ IPosition beamPos = pos;
+ if (posInLattice) {
+ this->_latticePosToStoragePos(beamPos, pos);
+ }
+ flux = _flux(sum, beamArea(beamPos)).getValue();
+ }
+ else {
+ String unit = pInImage_p->units().getName();
+ unit.downcase();
+ if (unit.contains("/beam")) {
+ return False;
+ }
+ flux = _flux(sum, 0).getValue();
+ }
+ return True;
+}

template <class T> Bool ImageStatistics<T>::_canDoFlux() const {
String unit = pInImage_p->units().getName();
=======================================
--- /branches/nov14/images/Images/test/dImageStatistics.cc Wed Jan 28
10:01:12 2015 UTC
+++ /branches/nov14/images/Images/test/dImageStatistics.cc Wed Mar 18
15:01:43 2015 UTC
@@ -316,10 +316,12 @@
}
if (validInputs(PLOTTING)) {
PGPlotter plotter(device);
+ /*
if (!stats.setPlotting(plotter, statisticTypes, nxy)) {
os << stats.errorMessage() << LogIO::POST;
return 1;
}
+ */
}

// Recover things
@@ -430,6 +432,7 @@
os << stats.errorMessage() << LogIO::POST;
return 1;
}
+ /*
if (validInputs(PLOTTING)) {
PGPlotter plotter(device);
if (!stats.setPlotting(plotter, statisticTypes, nxy)) {
@@ -437,6 +440,7 @@
return 1;
}
}
+ */

// Display statistics

=======================================
--- /branches/nov14/images/Images/test/tImageStatistics.cc Fri Nov 21
12:31:42 2014 UTC
+++ /branches/nov14/images/Images/test/tImageStatistics.cc Wed Mar 18
15:01:43 2015 UTC
@@ -146,6 +146,74 @@
}
}
}
+ {
+ writeTestString(
+ "test getConvertedStatistic(FLUX) and getStats() for FLUX"
+ );
+ CoordinateSystem csys = CoordinateUtil::defaultCoords3D();
+ IPosition shape(3, 10, 15, 20);
+ TempImage<Float> tim(TiledShape(shape), csys);
+ Array<Float> arr(shape);
+ indgen(arr);
+ tim.put(arr);
+ tim.setUnits("Jy/beam");
+ GaussianBeam beam (
+ Quantity(3, "arcmin"),
+ Quantity(2.5, "arcmin"),
+ Quantity(30, "deg")
+ );
+ ImageInfo info = tim.imageInfo();
+ info.setRestoringBeam(beam);
+ tim.setImageInfo(info);
+ ImageStatistics<Float> stats(tim);
+ Array<Float> flux;
+ AlwaysAssert(
+ stats.getConvertedStatistic (flux, LatticeStatsBase::FLUX), AipsError
+ );
+ AlwaysAssert(near(*flux.begin(), 529348.9), AipsError);
+ Vector<Int> axes(2, 0);
+ axes[1] = 1;
+ stats.setAxes(axes);
+ AlwaysAssert(
+ stats.getConvertedStatistic (flux, LatticeStatsBase::FLUX), AipsError
+ );
+ AlwaysAssert(flux.shape() == IPosition(1, 20), AipsError);
+ Vector<Double> statVals;
+ Double area = beam.getArea("arcmin2");
+ for (uInt i=0; i<20; ++i) {
+ Float expFlux = sum(arr(IPosition(3,0,0,i), IPosition(3, 9, 14,
i)))/area;
+ AlwaysAssert(near(flux(IPosition(1,i)), expFlux), AipsError);
+ AlwaysAssert(
+ stats.getStats(
+ statVals, IPosition(1, i), False
+ ), AipsError
+ );
+ AlwaysAssert(near(statVals[LatticeStatsBase::FLUX], expFlux),
AipsError);
+ }
+ tim.setUnits("K");
+ stats = ImageStatistics<Float>(tim);
+ AlwaysAssert(
+ stats.getConvertedStatistic (flux, LatticeStatsBase::FLUX), AipsError
+ );
+ AlwaysAssert(near(*flux.begin(), 1.61946e10), AipsError);
+ axes = Vector<Int>(2, 0);
+ axes[1] = 1;
+ stats.setAxes(axes);
+ AlwaysAssert(
+ stats.getConvertedStatistic (flux, LatticeStatsBase::FLUX), AipsError
+ );
+ AlwaysAssert(flux.shape() == IPosition(1, 20), AipsError);
+ for (uInt i=0; i<20; ++i) {
+ Float expFlux = sum(arr(IPosition(3,0,0,i), IPosition(3, 9, 14,
i)))*3600;
+ AlwaysAssert(near(flux(IPosition(1,i)), expFlux), AipsError);
+ AlwaysAssert(
+ stats.getStats(
+ statVals, IPosition(1, i), False
+ ), AipsError
+ );
+ AlwaysAssert(near(statVals[LatticeStatsBase::FLUX], expFlux),
AipsError);
+ }
+ }
cout << "ok" << endl;
}
catch (AipsError x) {
=======================================
--- /branches/nov14/lattices/LatticeMath/LatticeHistograms.h Wed Jan 28
10:01:12 2015 UTC
+++ /branches/nov14/lattices/LatticeMath/LatticeHistograms.h Wed Mar 18
15:01:43 2015 UTC
@@ -35,6 +35,7 @@
#include <casacore/lattices/LatticeMath/LatticeStatistics.h>
#include <casacore/lattices/LatticeMath/LatticeProgress.h>
#include <casacore/casa/Logging/LogIO.h>
+#include <casacore/casa/System/PGPlotter.h>
#include <casacore/scimath/Mathematics/NumericTraits.h>
#include <casacore/casa/iosfwd.h>

=======================================
--- /branches/nov14/lattices/LatticeMath/LatticeStatistics.h Tue Mar 3
12:22:11 2015 UTC
+++ /branches/nov14/lattices/LatticeMath/LatticeStatistics.h Wed Mar 18
15:01:43 2015 UTC
@@ -40,7 +40,6 @@
#include <casacore/lattices/LatticeMath/LatticeStatsDataProvider.h>
#include <casacore/lattices/LatticeMath/MaskedLatticeStatsDataProvider.h>
#include <casacore/scimath/Mathematics/NumericTraits.h>
-#include <casacore/casa/System/PGPlotter.h>
#include <casacore/casa/Utilities/DataType.h>
#include <casacore/casa/BasicSL/String.h>
#include <casacore/casa/Logging/LogIO.h>
@@ -284,21 +283,6 @@
// not list the output when making a plot.
Bool setList(const Bool& doList);

-// This functions enable you to specify which statistics you would like to
-// plot, sets the name of the PGPLOT plotting device and the number of
-// subplots in x and y per page. If you set <src>device</src>
-// but offer a zero length array for <src>nxy</src> then <src>nxy</src> is
-// set to [1,1]. Similarly, the default for <src>statsToPlot</src> is
-// to plot the mean and standard deviation. Use the helper function
-// <src>LatticeStatsBase::toStatisticTypes(String& stats)</src> to convert
-// a <src>String</src> to the desired <src>Vector<Int> statsToPlot</src>.
-// A return value of <src>False</src> indicates invalid plotting arguments
-// or that the internal state of the class is bad. If you don't call this
function,
-// the default state of the class is to not make plots.
- Bool setPlotting(PGPlotter& plotter,
- const Vector<Int>& statsToPlot,
- const Vector<Int>& nxy);
-
// Display the statistics by listing and/or plotting them. If you don't
call
// this function then you won't see anything ! A return value of
<src>False</src>
// indicates an invalid plotting device, or that the internal state of the
class is bad.
@@ -314,10 +298,6 @@
Int zAxis=-1, Int zLayer=-1,
Int hAxis=-1, Int hLayer=-1);

-// CLose plotter
- void closePlotting();
-
-
// Return the display axes. The returned vector will be valid only if
<src>setAxes</src>
// has been called, or if one of the active "display" or "get*" methods
has been called.
Vector<Int> displayAxes() const {return displayAxes_p;}
@@ -436,15 +416,6 @@
String error_p;
//
// Virtual Functions. See implementation to figure it all out !
-//
-// Get beam volume if possible. Your lattice needs to be
-// an ImageInterface for you to be able to do this.
-// See for example, class ImageStatistics. When you provide
-// the beam, then the Flux statistic, if requested, can be
-// computed. Returns False if beam not available, else True.
-// The implementation here returns False.
-
- virtual Bool _getBeamArea (Array<Double>& beamArea) const;

// FIXME The indirect dependence of this class on ImageInterface related
// issues (eg flux density) breaks encapsulation. All the
ImageInterface related code should be
@@ -498,6 +469,26 @@

// get the storage lattice shape
inline IPosition _storageLatticeShape() const { return
pStoreLattice_p->shape(); }
+
+ virtual Bool _computeFlux(
+ Array<AccumType>& flux, const Array<AccumType>& npts, const
Array<AccumType>& sum
+ );
+
+ virtual Bool _computeFlux(
+ AccumType& flux, AccumType sum, const IPosition& pos,
+ Bool posInLattice
+ );
+
+ // convert a position in the input lattice to the corresponding
+ // position in the stats storage lattice. The number of elements
+ // in storagePos will not be changed and only the first N elements
+ // will be modified where N = the number of elements in latticePos.
+ // <src>storagePos</src> must therefore have at least as many elements
+ // as <src>latticePos</src>. Returns False if
+ //<src>latticePos</src> is inconsistent with the input lattice.
+ void _latticePosToStoragePos(
+ IPosition& storagePos, const IPosition& latticePos
+ );

private:

@@ -505,7 +496,6 @@
CountedPtr<TempLattice<AccumType> > pStoreLattice_p;
Vector<Int> nxy_p, statsToPlot_p;
Vector<T> range_p;
- PGPlotter plotter_p;
Bool noInclude_p, noExclude_p;

Bool needStorageLattice_p, doneSomeGoodPoints_p, someGoodPointsValue_p;
@@ -544,71 +534,24 @@
LatticeStatsBase::StatisticsTypes type,
Bool dropDeg);

-// Convert a AccumType to a Float for plotting.
- static Float convertATtoF (AccumType value) {return
Float(std::real(value));};
-
-// Find the next good or bad point in an array
- Bool findNextDatum (uInt& iFound,
- const uInt& n,
- const Vector<AccumType>& mask,
- const uInt& iStart,
- const Bool& findGood) const;
-
-// Find the next label in a list of comma delimitered labels
- Bool findNextLabel (String& subLabel,
- Int& iLab,
- String& label) const;
-
// Find the median per cursorAxes chunk
void generateRobust ();

// Create a new storage lattice
Bool generateStorageLattice ();

-// Examine an array and determine how many segments of good points it
consists
-// of. A good point occurs if the array value is greater than zero.
- void lineSegments (uInt& nSeg,
- Vector<uInt>& start,
- Vector<uInt>& nPts,
- const Vector<AccumType>& mask) const;
-
// Given a location in the lattice and a statistic type, work
// out where to put it in the storage lattice
IPosition locInStorageLattice(const IPosition& latticePosition,
LatticeStatsBase::StatisticsTypes type)
const;

-// Draw each Y-axis sublabel in a string with a different colour
- void multiColourYLabel (String& label,
- PGPlotter& plotter,
- const String& LRLoc,
- const Vector<uInt>& colours,
- const Int& nLabs) const;
-
-// Plot an array which may have some blanked points.
-// Thus we plot it in segments
- void multiPlot (PGPlotter& plotter,
- const Vector<AccumType>& x,
- const Vector<AccumType>& y,
- const Vector<AccumType>& n) const;
-
// Find min and max of good data in arrays specified by pointers
void minMax (Bool& none, AccumType& dMin, AccumType& dMax,
const Vector<AccumType>& d,
const Vector<AccumType>& n) const;

-// Find the next nice PGPLOT colour index
- Int niceColour (Bool& initColours) const;
-
- // FIXME AFAIK, CASA no longer uses plotStats since we no longer use
pgplot. Can we
- // remove this method?
-// Plot the statistics
-
- Bool plotStats (const IPosition& dPos,
- const Matrix<AccumType>& ord,
- PGPlotter& plotter);
-
// Retrieve a statistic from the storage lattice and return in an array
- Bool retrieveStorageStatistic (Array<AccumType>& slice,
+ Bool retrieveStorageStatistic (Array<AccumType>& slice,
const LatticeStatsBase::StatisticsTypes
type,
const Bool dropDeg);

@@ -625,21 +568,9 @@
// See if there were some valid points found in the storage lattice
Bool someGoodPoints ();

-
// Stretch min and max by 5%
void stretchMinMax (AccumType& dMin, AccumType& dMax) const;

- // convert a position in the input lattice to the corresponding
- // position in the stats storage lattice. The number of elements
- // in storagePos will not be changed and only the first N elements
- // will be modified where N = the number of elements in latticePos.
- // <src>storagePos</src> must therefore have at least as many elements
- // as <src>latticePos</src>. Returns False if
- //<src>latticePos</src> is inconsistent with the input lattice.
- void _latticePosToStoragePos(
- IPosition& storagePos, const IPosition& latticePos
- );
-
CountedPtr<StatisticsAlgorithm<AccumType, const T*, const Bool*> >
_createStatsAlgorithm() const;

void _configureDataProviders(
@@ -650,128 +581,7 @@
void _doStatsLoop(uInt nsets, CountedPtr<LattStatsProgress>
progressMeter);
};

-// <summary> Generate statistics, tile by tile, from a masked lattice
</summary>
-//
-// <use visibility=export>
-//
-// <reviewed reviewer="" date="yyyy/mm/dd" tests="" demos="">
-// </reviewed>
-//
-// <prerequisite>
-// <li> <linkto class=LatticeApply>LatticeApply</linkto>
-// <li> <linkto class=TiledCollapser>TiledCollapser</linkto>
-// </prerequisite>
-//
-// <etymology>
-// This class is used by <src>LatticeStatistics</src> to generate
-// statistical sum from an input <src>MaskedLattice</src>.
-// The input lattice is iterated through in tile-sized chunks
-// and fed to an object of this class.
-// </etymology>
-//
-// <synopsis>
-// <src>StatsTiledCollapser</src> is derived from
<src>TiledCollapser</src> which
-// is a base class used to define methods. Objects of this base class are
-// used by <src>LatticeApply</src> functions. In this particular case,
-// we are interested in <src>LatticeApply::tiledApply</src>. This
function iterates
-// through a <src>MaskedLattice</src> and allows you to collapse one or
more
-// axes, computing some values from it, and placing those values into
-// an output <src>MaskedLattice</src>. It iterates through the input
-// lattice in optimal tile-sized chunks. <src>LatticeStatistics</src>
-// uses a <src>StatsTiledCollapser</src> object which it gives to
-// <src>LatticeApply::tiledApply</src> for digestion. After it has
-// done its work, <src>LatticeStatistics</src> then accesses the output
-// <src>Lattice</src> that it made.
-// </synopsis>
-//
-// <example>
-// <srcblock>
-//// Create collapser. Control information is passed in via the constructor
-//
-// StatsTiledCollapser<T> collapser(range_p, noInclude_p, noExclude_p,
-// fixedMinMax_p, blcParent_p);
-//
-//// This is the first output axis getting collapsed values. In
LatticeStatistics
-//// this is the last axis of the output lattice
-//
-// Int newOutAxis = outLattice.ndim()-1;
-//
-//// tiledApply does the work by passing the collapser data in chunks
-//// and by writing the results into the output lattice
-//
-// LatticeApply<T>::tiledApply(outLattice, inLattice,
-// collapser, collapseAxes,
-// newOutAxis);
-//
-// </srcblock>
-// In this example, a collapser is made and passed to LatticeApply.
-// Afterwards, the output Lattice is available for use.
-// The Lattices must all be the correct shapes on input to tiledApply
-// </example>
-//
-// <motivation>
-// The LatticeApply classes enable the ugly details of optimal
-// Lattice iteration to be hidden from the user.
-// </motivation>
-//
-// <todo asof="1998/05/10">
-// <li>
-// </todo>
-
-/*
-template <class T, class U=T>
-class StatsTiledCollapser : public TiledCollapser<T,U>
-{
-public:
-// Constructor provides pixel selection range and whether that
-// range is an inclusion or exclusion range. If
<src>fixedMinMax=True</src>
-// and an inclusion range is given, the min and max is set to
-// that inclusion range.
- StatsTiledCollapser(const Vector<T>& pixelRange, Bool noInclude,
- Bool noExclude, Bool fixedMinMax);
-
-// Initialize process, making some checks
- virtual void init (uInt nOutPixelsPerCollapse);
-
-// Initialiaze the accumulator
- virtual void initAccumulator (uInt n1, uInt n3);
-
-// Process the data in the current chunk.
- virtual void process (
- uInt accumIndex1, uInt accumIndex3,
- const T* inData, const Bool* inMask,
- uInt dataIncr, uInt maskIncr,
- uInt nrval, const IPosition& startPos,
- const IPosition& shape
- );
-
-// End the accumulation process and return the result arrays
- virtual void endAccumulator(Array<U>& result,
- Array<Bool>& resultMask,
- const IPosition& shape);
-
-// Can handle null mask
- virtual Bool canHandleNullMask() const {return True;};
-
-// Find the location of the minimum and maximum data values
-// in the input lattice.
- void minMaxPos(IPosition& minPos, IPosition& maxPos);
-
-private:
-
- Bool noInclude_p, noExclude_p, fixedMinMax_p;
- IPosition minPos_p, maxPos_p;
-
- uInt n1_p;
- uInt n3_p;
- vector<ClassicalStatistics<U, const T*, const Bool*> > _cs;
- vector<std::pair<U, U> > _ranges;
- DataType _type;
- vector<Bool> _first;
-};
-*/
-
-} //# NAMESPACE CASACORE - END
+} //# NAMESPACE CASA - END

#ifndef CASACORE_NO_AUTO_TEMPLATES
#include <casacore/lattices/LatticeMath/LatticeStatistics.tcc>
=======================================
--- /branches/nov14/lattices/LatticeMath/LatticeStatistics.tcc Tue Mar 3
12:22:11 2015 UTC
+++ /branches/nov14/lattices/LatticeMath/LatticeStatistics.tcc Wed Mar 18
15:01:43 2015 UTC
@@ -215,7 +215,6 @@
statsToPlot_p = other.statsToPlot_p;
range_p.resize(other.range_p.size());
range_p = other.range_p;
- plotter_p = other.plotter_p;
doList_p = other.doList_p;
noInclude_p = other.noInclude_p;
noExclude_p = other.noExclude_p;
@@ -252,8 +251,6 @@
LatticeStatistics<T>::~LatticeStatistics() {
delete pInLattice_p;
pInLattice_p = 0;
- //delete pStoreLattice_p;
- // pStoreLattice_p = 0;
}

template <class T>
@@ -383,97 +380,6 @@
return True;
}

-template <class T>
-Bool LatticeStatistics<T>::setPlotting(PGPlotter& plotter,
- const Vector<Int>& statsToPlot,
- const Vector<Int>& nxy)
-//
-// Assign the desired PGPLOT device name and number
-// of subplots
-//
-{
- if (!goodParameterStatus_p) {
- return False;
- }
-
-// Is new plotter attached ?
-
- if (!plotter.isAttached()) {
- if (haveLogger_p) {
- error_p = "Input plotter is not attached";
- }
- goodParameterStatus_p = False;
- return False;
- }
-
-// Don't reattach to the same plotter. The assignment will
-// close the previous device
-
- if (plotter_p.isAttached()) {
- if (plotter_p.qid() != plotter.qid()) plotter_p = plotter;
- } else {
- plotter_p = plotter;
- }
-
-
-// Make sure requested statistics are valid
-// Set need robust statistics flag here as well
-
- statsToPlot_p.resize(0);
- statsToPlot_p = statsToPlot;
- for (uInt i=0; i<statsToPlot_p.nelements(); i++) {
- if (statsToPlot_p(i) < 0 || statsToPlot_p(i) > NSTATS-1) {
- error_p = "Invalid statistic requested for display";
- goodParameterStatus_p = False;
- return False;
- }
-
-// If the user wants robust stats, signal this and if they
-// did not previously ask for them, signify we need to
-// regenerate the storage lattice as well - the robust
-// stats are just written directly into the storage lattice
-
- if (
- statsToPlot_p(i)==Int(LatticeStatsBase::MEDIAN)
- || statsToPlot_p(i)==Int(LatticeStatsBase::MEDABSDEVMED)
- || statsToPlot_p(i)==Int(LatticeStatsBase::QUARTILE)
- || statsToPlot_p(i)==Int(LatticeStatsBase::Q1)
- || statsToPlot_p(i)==Int(LatticeStatsBase::Q3)
- ) {
- if (!doRobust_p) {
- needStorageLattice_p = True;
- }
- doRobust_p = True;
- }
- }
-
-
-// Plotting device and subplots. nxy_p is set to [1,1] if zero length
-
- nxy_p.resize(0);
- nxy_p = nxy;
- ostringstream os;
- if (!LatticeStatsBase::setNxy(nxy_p, os)) {
- error_p = "Invalid number of subplots";
- goodParameterStatus_p = False;
- return False;
- }
-
-
-// Set mean and sigma if no statistics requested
-
- if (statsToPlot_p.nelements()==0) {
- error_p = "No plot statistics requested, setting mean and std dev";
- statsToPlot_p.resize(2);
- statsToPlot_p(0) = MEAN;
- statsToPlot_p(1) = SIGMA;
- }
-
- return True;
-}
-
-
-
template <class T>
Bool LatticeStatistics<T>::setNewLattice(const MaskedLattice<T>& lattice)
{
@@ -506,7 +412,6 @@
needStorageLattice_p = True;
return True;
}
-

template <class T>
Bool LatticeStatistics<T>::getConvertedStatistic (Array<T>& stats,
@@ -572,7 +477,6 @@
}
return True;
}
-

template <class T>
Bool LatticeStatistics<T>::getStats(
@@ -604,22 +508,12 @@
stats.resize(0);
return True;
}
- //stats(MEAN) = LattStatsSpecialize::getMean(stats(SUM), n);
- //stats(VARIANCE) =
LattStatsSpecialize::getVariance(stats(SUM),
stats(SUMSQ), n);
stats(SIGMA) = LattStatsSpecialize::getSigma(stats(VARIANCE));
stats(RMS) = LattStatsSpecialize::getRms(stats(SUMSQ), n);
stats(FLUX) = 0;
if (_canDoFlux()) {
- Array<Double> beamArea;
- if (_getBeamArea(beamArea)) {
- IPosition beamPos = pos;
- if (posInLattice) {
- this->_latticePosToStoragePos(beamPos, pos);
- }
- stats(FLUX) = _flux(stats(SUM), beamArea(beamPos)).getValue();
- }
- else {
- stats(FLUX) = _flux(stats(SUM), 0).getValue();
+ if (! _computeFlux(stats(FLUX), stats(SUM), pos, posInLattice)) {
+ return False;
}
}
return True;
@@ -666,11 +560,23 @@
//
return (maxFull_p > minFull_p);
}
-
-

// Private functions

+template <class T>
+Bool LatticeStatistics<T>::_computeFlux(
+ Array<AccumType>&, const Array<AccumType>&, const Array<AccumType>&
+) {
+ ThrowCc("This object does not support computing fluxes");
+}
+
+template <class T>
+Bool LatticeStatistics<T>::_computeFlux(
+ AccumType&, AccumType, const IPosition&,
+ Bool
+) {
+ ThrowCc("This object does not support computing fluxes");
+}

template <class T>
Bool LatticeStatistics<T>::calculateStatistic (Array<AccumType>& slice,
@@ -723,34 +629,15 @@
}
}
else if (type==FLUX) {
- if (! _canDoFlux()) {
+ if (_canDoFlux()) {
+ retrieveStorageStatistic (sum, SUM, dropDeg);
+ return _computeFlux(slice, nPts, sum);
+ }
+ else {
slice.resize(IPosition(0,0));
return False;
}
- Array<Double> beamArea;
- Bool gotBeamArea = _getBeamArea(beamArea);
- retrieveStorageStatistic (sum, SUM, dropDeg);
- ReadOnlyVectorIterator<AccumType> sumIt(sum);
- PtrHolder<ReadOnlyVectorIterator<Double> > beamAreaIter(
- gotBeamArea ? new ReadOnlyVectorIterator<Double>(beamArea) : 0
- );
- while (!nPtsIt.pastEnd()) {
- for (uInt i=0; i<n1; i++) {
- if (LattStatsSpecialize::hasSomePoints(nPtsIt.vector()(i))) {
- //sliceIt.vector()(i) = sumIt.vector()(i) /
beamAreaIter.vector()(i);
- sliceIt.vector()(i) = _flux(
- sumIt.vector()(i), gotBeamArea ? beamAreaIter->vector()(i) :
0
- ).getValue();
- }
- }
- nPtsIt.next();
- sumIt.next();
- sliceIt.next();
- if (gotBeamArea) {
- beamAreaIter->next();
- }
- }
- }
+ }
else if (type==SIGMA) {
Array<AccumType> variance;
retrieveStorageStatistic (variance, VARIANCE, dropDeg);
@@ -920,7 +807,6 @@
_algConf.algorithm == StatisticsData::CLASSICAL
&& timeOld < timeNew
) {
- //cout << "using old method" << endl;
// use older method for higher performance in the large loop count
// regime
//timer.mark();
@@ -941,12 +827,9 @@
collapser.minMaxPos(minPos_p, maxPos_p);
}
else {
- //cout << "using new method" << endl;
- //timer.mark();
_doStatsLoop(nsets, pProgressMeter);
}
pProgressMeter = NULL;
- //cout << "time " << timer.real() << endl;

// Do robust statistics separately as required.
generateRobust();
@@ -1020,11 +903,6 @@
uInt nSublatticeSteps = nelem > 4096*4096
? nelem/subLat.advisedMaxPixels()
: 1;
- /*
- cout << "nelem " << nelem << endl;
- cout << "nice cursor shape " << subLat.niceCursorShape() << endl;
- cout << "nsteps " << nsets*nSublatticeSteps << endl;
- */
progressMeter->init(nsets*nSublatticeSteps);
}
if(subLat.isMasked()) {
@@ -1536,23 +1414,6 @@
setStream(os, oPrec);
os << setw(oDWidth)
<< dMax;
-
- /*
- if (!fixedMinMax_p) {
- os << "Minimum value ";
- os << setw(oWidth) << String(os6);
- if (type==TpFloat) {
- os << " at " << blcParent_p + minPos_p+1;
- }
-
- os << "Maximum value ";
- os << setw(oWidth) << String(os7);
- if (type==TpFloat) {
- os << " at " << blcParent_p + maxPos_p+1 << endl;
- }
- }
- */
-
stats += os.str();
stats += '\n';
return True;
@@ -1620,10 +1481,6 @@
if (_canDoFlux()) {
ord(i,FLUX) = _flux(matrix(i,SUM), area).getValue();
}
- /*
- ord(i,VARIANCE) = LattStatsSpecialize::getVariance(
- matrix(i,SUM), matrix(i,SUMSQ), nPts);
- */
ord(i,SIGMA) = LattStatsSpecialize::getSigma(
matrix(i,VARIANCE));
ord(i,RMS) = LattStatsSpecialize::getRms(
@@ -1767,24 +1624,6 @@
Int layer = 0;
for ( pixelIterator.reset(); !pixelIterator.atEnd(); pixelIterator++ )
{
IPosition dPos = pixelIterator.position();
- /*
- if (displayAxes_p.nelements() == 2) {
- if (zAx == 1) {
- if (dPos[1] != zLayer)
- continue;
- else
- layer = hLayer;
- }
- if (hAx == 1) {
- if (dPos[1] != hLayer)
- continue;
- else
- layer = zLayer;
- }
- }
- if (displayAxes_p.nelements() == 1)
- layer = zLayer;
- */
if (displayAxes_p.nelements() == 2) {
if (zAx == 1) {
if (dPos[1] != zLayer) {
@@ -1914,27 +1753,6 @@
IPosition blc(pInLattice_p->ndim(),0);
IPosition trc(pInLattice_p->shape()-1);

- /*
- Int len0;
- if (nStatsAxes == 1) {
- len0 = 8;
- os << setw(len0) << "Profile ";
- }
- else if (nStatsAxes == 2) {
- len0 = 6;
- os << setw(len0) << "Plane ";
- }
- else if (nStatsAxes == 3) {
- len0 = 5;
- os << setw(len0) << "Cube ";
- }
- else {
- len0 = 11;
- os << setw(len0) << "Hyper-cube ";
- }
- */
-
-
os << setw(10) << "Npts";
os << setw(oDWidth) << "Sum";
if (_canDoFlux()) {
@@ -2099,18 +1917,11 @@

// Do we have anything to do

- if (!doList_p && !plotter_p.isAttached() && haveLogger_p) {
+ if (!doList_p && haveLogger_p) {
os_p << LogIO::NORMAL1 << "There is nothing to plot or list" <<
LogIO::POST;
return True;
}
-// Set up some plotting things

- if (plotter_p.isAttached()) {
- plotter_p.subp(nxy_p(0), nxy_p(1));
- plotter_p.ask(True);
- plotter_p.sch (1.2);
- plotter_p.svp(0.1,0.9,0.1,0.9);
- }
// Generate storage lattice if required

if (needStorageLattice_p) {
@@ -2159,12 +1970,6 @@
const AccumType& nPts = matrix(i,NPTS);
if (LattStatsSpecialize::hasSomePoints(nPts)) {
ord(i,MEAN) = LattStatsSpecialize::getMean(matrix(i,SUM),
nPts);
- /*
- * see 2012may23 comment above
- if (hasBeam) {
- ord(i,FLUX) = matrix(i,SUM) / beamArea;
- }
- */
ord(i,SIGMA) = LattStatsSpecialize::getSigma(ord(i,VARIANCE));
ord(i,RMS) = LattStatsSpecialize::getRms(matrix(i,SUMSQ),
nPts);
}
@@ -2178,13 +1983,6 @@
ord(j,i) = matrix(j,i);
}
}
-
-// Plot statistics
-
- if (plotter_p.isAttached()) {
- if (!plotStats (pixelIterator.position(), ord, plotter_p)) return
False;
- }
-

// List statistics

@@ -2192,805 +1990,11 @@
if (!listStats(hasBeam, pixelIterator.position(), ord)) return
False;
}
}
-
-
-
-
-// Finish up
-
- if (plotter_p.isAttached()) {
- plotter_p.updt();
- }
return True;
}
-
-
-
-template <class T>
-void LatticeStatistics<T>::lineSegments (uInt& nSeg,
- Vector<uInt>& start,
- Vector<uInt>& nPts,
- const Vector<AccumType>& mask) const
-//
-// Examine an array and determine how many segments
-// of good points it consists of. A good point
-// occurs if the array value is greater than zero.
-//
-// Inputs:
-// mask The array. Note that even if <T> is complex, only
-// the real part of this array contains the information.
-// Outputs:
-// nSeg Number of segments
-// start Indices of start of each segment
-// nPts Number of points in segment
-//
-{
- Bool finish = False;
- nSeg = 0;
- uInt iGood, iBad;
- const uInt n = mask.nelements();
- start.resize(n);
- nPts.resize(n);
-
- for (uInt i=0; !finish;) {
- if (!findNextDatum (iGood, n, mask, i, True)) {
- finish = True;
- } else {
- nSeg++;
- start(nSeg-1) = iGood;
-
- if (!findNextDatum (iBad, n, mask, iGood, False)) {
- nPts(nSeg-1) = n - start(nSeg-1);
- finish = True;
- } else {
- nPts(nSeg-1) = iBad - start(nSeg-1);
- i = iBad + 1;
- }
- }
- }
- start.resize(nSeg,True);
- nPts.resize(nSeg,True);
-}
-
-
-template <class T>
-void LatticeStatistics<T>::multiColourYLabel (String& label,
- PGPlotter& plotter,
- const String& LRLoc,
- const Vector<uInt>& colours,
- const Int& nLabs) const
-
-//
-// Draw each Y-axis sublabel in a string with a different colour
-//
-{
-// Get attributes
-
-
- Vector<Float> result= plotter.qwin();
- Float y1 = result(2);
- Float y2 = result(3);
- Int sci = plotter.qci();
-
-
-// Find y-location of start of string as fraction of window
-
- result.resize(0);
- result = plotter.qtxt(0.0, 0.0, 90.0, 0.0, label);
- Vector<Float> xb = result(Slice(0,4));
- Vector<Float> yb = result(Slice(4,4));
- Float dy = yb(2)-yb(0);
- Float yLoc = abs(0.5*(y2-y1-dy)/(y2-y1));
-
-
-// Loop over number of sub-labels and write them in colour
-
- String subLabel;
- Float just = 0.0;
- Float disp = 2.5;
- if (LRLoc == "R") disp = 3.0;
- for (Int iLab=0; iLab<nLabs; iLab++) {
-
-// Fish out next sub label
-
- if (!findNextLabel (subLabel, iLab, label)) {
- plotter.sci (sci);
- return;
- }
-
-
-// Write it
-
- if (iLab < nLabs-1) subLabel = subLabel + ",";
- if (iLab > 0) subLabel.prepend(" ");
- plotter.sci (colours(iLab));
- plotter.mtxt (LRLoc, disp, yLoc, just, subLabel);
-
-
-// Increment y location. pgqtxt won't count a leading blank so
-// replace it with a character for length counting purposes. These
-// stupid string classes make this very hard work.
-
- String s2;
- if (iLab > 0) {
- String s(subLabel(1,subLabel.length()-1));
- s2 = "x" + s;
- } else
- s2 = subLabel;
- result.resize(0);
- result = plotter.qtxt (0.0, 0.0, 90.0, 0.0, s2.chars());
- xb = result(Slice(0,4));
- yb = result(Slice(4,4));
- dy = abs((yb(2)-yb(0))/(y2-y1));
- yLoc += dy;
- }
-
-// Set colour back to what it was
-
- plotter.sci (sci);
- return;
-}
-
-
-
-
-template <class T>
-void LatticeStatistics<T>::multiPlot (PGPlotter& plotter,
- const Vector<AccumType>& x,
- const Vector<AccumType>& y,
- const Vector<AccumType>& mask) const
-//
-// Plot an array which may have some blanked points.
-// Thus we plot it in segments
-//
-// Currently, only the real part of the <AccumType> data is
-// plotted
-//
-// Inputs:
-// x,y,mask Abcissa, ordinate, and "masking" array
-// (if > 0 plot it)
-{
-
-// Find number of segments in this array
-
- uInt nSeg = 0;
- Vector<uInt> start;
- Vector<uInt> nPts;
- lineSegments (nSeg, start, nPts, mask);
-
-// Loop over segments and plot them
-
- Vector<Float> xF, yF;
- for (uInt i=0; i<nSeg; i++) {
- const uInt ip = start(i);
- if (nPts(i) == 1) {
- xF.resize(1);
- yF.resize(1);
- xF(0) = convertATtoF(x(ip));
- yF(0) = convertATtoF(y(ip));
- plotter.pt (xF, yF, 1);
- } else {
- xF.resize(nPts(i));
- yF.resize(nPts(i));
- for (uInt j=0; j<nPts(i); j++) {
- xF(j) = convertATtoF(x(start(i)+j));
- yF(j) = convertATtoF(y(start(i)+j));
- }
- plotter.line (xF, yF);
- }
- }
-}
-
-template <class T>
-Int LatticeStatistics<T>::niceColour (Bool& initColours) const
-{
- static Int colourIndex = 1;
- if (initColours) {
- colourIndex = 1;
- initColours = False;
- }
-
- colourIndex++;
- if (colourIndex == 4 || colourIndex == 14) colourIndex++;
- return colourIndex;
-}
-
-
-template <class T>
-Bool LatticeStatistics<T>::plotStats (
- const IPosition& dPos,
- const Matrix<AccumType>& stats,
- PGPlotter& plotter)
-//
-// Plot the desired statistics.
-//
-// Inputs:
-// dPos The location of the start of the cursor in the
-// storage lattice for this line
-// stats Statistics matrix
-//
-{
-// The plotting for Complex just take the real part which is
-// not very useful. Until I do someting better, stub it out
-
- T* dummy = 0;
- DataType type = whatType(dummy);
- if (type!=TpFloat) {
- os_p << LogIO::WARN << "Plotting not yet available for complex
images " << LogIO::POST;
- return True;
- }
-
-// Work out what we are plotting
-
- const uInt n = statsToPlot_p.nelements();
- Bool doMean, doSigma, doVar, doRms, doSum, doSumSq;
- Bool doMin, doMax, doNPts, doFlux, doMedian;
- Bool doMedAbsDevMed, doQuartile;
- linearSearch(doMean, statsToPlot_p, Int(MEAN), n);
- linearSearch(doMedian, statsToPlot_p, Int(MEDIAN), n);
- linearSearch(doMedAbsDevMed, statsToPlot_p, Int(MEDABSDEVMED), n);
- linearSearch(doQuartile, statsToPlot_p, Int(QUARTILE), n);
- linearSearch(doSigma, statsToPlot_p, Int(SIGMA), n);
- linearSearch(doVar, statsToPlot_p, Int(VARIANCE), n);
- linearSearch(doRms, statsToPlot_p, Int(RMS), n);
- linearSearch(doSum, statsToPlot_p, Int(SUM), n);
- linearSearch(doSumSq, statsToPlot_p, Int(SUMSQ), n);
- linearSearch(doMin, statsToPlot_p, Int(MIN), n);
- linearSearch(doMax, statsToPlot_p, Int(MAX), n);
- linearSearch(doNPts, statsToPlot_p, Int(NPTS), n);
- linearSearch(doFlux, statsToPlot_p, Int(FLUX), n);
- if (! _canDoFlux() ) {
- doFlux = False;
- }
-//
- Bool none;
- Bool first = True;
- Int nL = 0;
- Int nR = 0;
-
-// Generate abcissa. Note that T(n) where T is COmplex results in (n+0i)
-
- const Int n1 = stats.shape()(0);
- Vector<AccumType> abc(n1);
- for (Int j=0; j<n1; j++) abc(j) = AccumType(j+1);
-
-// Find extrema. Return if there were no valid points to plot
-
- AccumType yMin, yMax, xMin, xMax, yLMin, yLMax, yRMin, yRMax;
- // avoid compiler warning by initializing yMin, yMax
- yMin = 0;
- yMax = 0;
- minMax(none, xMin, xMax, abc, stats.column(NPTS));
- if (none) return True;
-
-// Left hand y axis
-
- if (doMean) {
- minMax(none, yLMin, yLMax, stats.column(MEAN), stats.column(NPTS));
- first = False;
- nL++;
- }
- if (doMedian) {
- minMax(none, yMin, yMax, stats.column(MEDIAN), stats.column(NPTS));
- if (first) {
- yLMin = yMin;
- yLMax = yMax;
- } else {
- yLMin = min(yLMin,yMin);
- yLMax = max(yLMax,yMax);
- }
- first = False;
- nL++;
- }
- if (doFlux) {
- minMax(none, yMin, yMax, stats.column(FLUX), stats.column(NPTS));
- if (first) {
- yLMin = yMin;
- yLMax = yMax;
- } else {
- yLMin = min(yLMin,yMin);
- yLMax = max(yLMax,yMax);
- }
- first = False;
- nL++;
- }
- if (doSum) {
- minMax(none, yMin, yMax, stats.column(SUM), stats.column(NPTS));
- if (first) {
- yLMin = yMin;
- yLMax = yMax;
- } else {
- yLMin = min(yLMin,yMin);
- yLMax = max(yLMax,yMax);
- }
- first = False;
- nL++;
- }
- if (doSumSq) {
- minMax(none, yMin, yMax, stats.column(SUMSQ), stats.column(NPTS));
- if (first) {
- yLMin = yMin;
- yLMax = yMax;
- } else {
- yLMin = min(yLMin,yMin);
- yLMax = max(yLMax,yMax);
- }
- first = False;
- nL++;
- }
- if (doMin) {
- minMax(none, yMin, yMax, stats.column(MIN), stats.column(NPTS));
- if (first) {
- yLMin = yMin;
- yLMax = yMax;
- } else {
- yLMin = min(yLMin,yMin);
- yLMax = max(yLMax,yMax);
- }
- first = False;
- nL++;
- }
- if (doMax) {
- minMax(none, yMin, yMax, stats.column(MAX), stats.column(NPTS));
- if (first) {
- yLMin = yMin;
- yLMax = yMax;
- } else {
- yLMin = min(yLMin,yMin);
- yLMax = max(yLMax,yMax);
- }
- first = False;
- nL++;
- }
- if (doNPts) {
- minMax(none, yMin, yMax, stats.column(NPTS), stats.column(NPTS));
- if (first) {
- yLMin = yMin;
- yLMax = yMax;
- } else {
- yLMin = min(yLMin,yMin);
- yLMax = max(yLMax,yMax);
- }
- first = False;
- nL++;
- }
-
-
-// Right hand y axis
-
- first = True;
- if (doSigma) {
- minMax(none, yRMin, yRMax, stats.column(SIGMA), stats.column(NPTS));
- first = False;
- nR++;
- }
- if (doVar) {
- minMax(none, yMin, yMax, stats.column(VARIANCE), stats.column(NPTS));
- if (first) {
- yRMin = yMin;
- yRMax = yMax;
- } else {
- yRMin = min(yRMin,yMin);
- yRMax = max(yRMax,yMax);
- }
- first = False;
- nR++;
- }
- if (doRms) {
- minMax(none, yMin, yMax, stats.column(RMS), stats.column(NPTS));
- if (first) {
- yRMin = yMin;
- yRMax = yMax;
- } else {
- yRMin = min(yRMin,yMin);
- yRMax = max(yRMax,yMax);
- }
- nR++;
- }
- if (doMedAbsDevMed) {
- minMax(none, yMin, yMax, stats.column(MEDABSDEVMED),
stats.column(NPTS));
- if (first) {
- yRMin = yMin;
- yRMax = yMax;
- } else {
- yRMin = min(yRMin,yMin);
- yRMax = max(yRMax,yMax);
- }
- nR++;
- }
- if (doQuartile) {
- minMax(none, yMin, yMax, stats.column(QUARTILE), stats.column(NPTS));
- if (first) {
- yRMin = yMin;
- yRMax = yMax;
- } else {
- yRMin = min(yRMin,yMin);
- yRMax = max(yRMax,yMax);
- }
- nR++;
- }
-//
- stretchMinMax(xMin, xMax);
- if (nL>0) stretchMinMax(yLMin, yLMax);
- if (nR>0) stretchMinMax(yRMin, yRMax);
-
-// Set labels.
-
- String hLabel, xLabel;
- getLabels(hLabel, xLabel, dPos);
-//
- String yLLabel = "";
- String yRLabel = "";
-
- Int nLLabs = 0;
- if (nL>0) {
- if (doMean) {
- yLLabel += "Mean,";
- nLLabs++;
- }
- if (doMedian) {
- yLLabel += "Median,";
- nLLabs++;
- }
- if (doFlux) {
- yLLabel += "Flux,";
- nLLabs++;
- }
- if (doSum) {
- yLLabel += "Sum,";
- nLLabs++;
- }
- if (doSumSq) {
- yLLabel += "Sum Squared,";
- nLLabs++;
- }
- if (doMin) {
- yLLabel += "Min,";
- nLLabs++;
- }
- if (doMax) {
- yLLabel += "Max,";
- nLLabs++;
- }
- if (doNPts) {
- yLLabel += "nPts,";
- nLLabs++;
- }
- yLLabel.del(Int(yLLabel.length()-1),1);
- }
-
- Int nRLabs = 0;
- if (nR>0) {
- if (doSigma) {
- yRLabel += "Std dev,";
- nRLabs++;
- }
- if (doVar) {
- yRLabel += "Variance,";
- nRLabs++;
- }
- if (doRms) {
- yRLabel += "Rms,";
- nRLabs++;
- }
- if (doMedAbsDevMed) {
- yRLabel += "MedAbsDevMed,";
- nRLabs++;
- }
- if (doQuartile) {
- yRLabel += "Quartile,";
- nRLabs++;
- }
- yRLabel.del(Int(yRLabel.length()-1),1);
- }
-
-// Do plots. Here we convert to real for now. To properly
-// make this deal with Complex I will have to be cleverer
-
- Vector<uInt> lCols(nL);
- Vector<uInt> rCols(nR);
- Int ls = 0;
- Int i = -1;
- Bool initColours = True;
- plotter.page();
-
- if (nL>0) {
- plotter.swin(real(xMin), real(xMax), real(yLMin), real(yLMax));
- if (nR>0) {
- plotter.box("BCNST", 0.0, 0, "BNST", 0.0, 0);
- } else {
- plotter.box("BCNST", 0.0, 0, "BCNST", 0.0, 0);
- }
- plotter.lab(xLabel, "", "");
-
- if (doMean) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(MEAN), stats.column(NPTS));
- }
- if (doMedian) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(MEDIAN), stats.column(NPTS));
- }
- if (doFlux) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(FLUX), stats.column(NPTS));
- }
- if (doSum) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(SUM), stats.column(NPTS));
- }
- if (doSumSq) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(SUMSQ), stats.column(NPTS));
- }
- if (doMin) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(MIN), stats.column(NPTS));
- }
- if (doMax) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(MAX), stats.column(NPTS));
- }
- if (doNPts) {
- if (++ls > 5) ls = 1;
- plotter.sls (ls);
-
- lCols(++i) = niceColour (initColours);
- plotter.sci (lCols(i));
-
- multiPlot(plotter, abc, stats.column(NPTS), stats.column(NPTS));
- }
-
-// Y label
-
- multiColourYLabel (yLLabel, plotter_p, "L", lCols, nLLabs);
-
- }
- plotter.sls (1);
- plotter.sci (1);
-
-
- i = -1;
- if (nR>0) {
- plotter.swin(real(xMin), real(xMax), real(yRMin), real(yRMax));
- plotter.sci (1);
- if (nL>0)
- plotter.box("", 0.0, 0, "CMST", 0.0, 0);
- else {
- plotter.box("BCNST", 0.0, 0, "BCMST", 0.0, 0);
- plotter.lab(xLabel, "", "");
- }
-
***The diff for this file has been truncated for email.***
=======================================
--- /branches/nov14/lattices/LatticeMath/LatticeStatsBase.h Thu Jan 29
09:15:47 2015 UTC
+++ /branches/nov14/lattices/LatticeMath/LatticeStatsBase.h Wed Mar 18
15:01:43 2015 UTC
@@ -62,19 +62,11 @@
// </synopsis>
//
// <example>
-// <srcBlock>
+// <srcblock>
// Vector<Int> statsToPlot =
LatticeStatsBase::toStatisticTypes("mean,rms,sigma");
-// </srcBlock>
+// </srcblock>
// </example>
//
-// <motivation>
-// My sensibilities were offended at having to say
-//
-// <src>Vector<Int> statsToPlot =
LatticeStatistics<Float>::toStatisticTypes("mean,rms,sigma");</src>
-//
-// when the <src><Float></src> was meaningless.
-// </motivation>
-//
// <todo asof="yyyy/mm/dd">
// </todo>

=======================================
--- /branches/nov14/ms/MeasurementSets/MSKeys.h Thu Mar 12 10:32:41 2015 UTC
+++ /branches/nov14/ms/MeasurementSets/MSKeys.h Wed Mar 18 15:01:43 2015 UTC
@@ -22,7 +22,7 @@
//# 520 Edgemont Road
//# Charlottesville, VA 22903-2475 USA
//#
-//# $Id: MSTable.h 21521 2014-12-10 08:06:42Z gervandiepen $
+//# $Id$

#ifndef MS_MSKEYS_H
#define MS_MSKEYS_H
=======================================
--- /branches/nov14/ms/MeasurementSets/MSMetaData.cc Thu Mar 12 10:32:41
2015 UTC
+++ /branches/nov14/ms/MeasurementSets/MSMetaData.cc Wed Mar 18 15:01:43
2015 UTC
@@ -840,6 +840,19 @@
}
return spws;
}
+
+std::vector<std::set<uInt> > MSMetaData::getSpwToDataDescriptionIDMap()
const {
+ // TODO perhaps cache the result, but atm doesn't seem worth doing
+ std::map<std::pair<uInt, uInt>, uInt> spwPolToDDID =
getSpwIDPolIDToDataDescIDMap();
+ std::vector<std::set<uInt> > mymap(nSpw(True));
+ std::map<std::pair<uInt, uInt>, uInt>::const_iterator iter =
spwPolToDDID.begin();
+ std::map<std::pair<uInt, uInt>, uInt>::const_iterator end =
spwPolToDDID.end();
+ while (iter != end) {
+ mymap[iter->first.first].insert(iter->second);
+ ++iter;
+ }
+ return mymap;
+}

uInt MSMetaData::nSpw(Bool includewvr) const {
if (_nSpw > 0) {
@@ -904,7 +917,7 @@
Vector<Int>::const_iterator curField = allFieldIDs->begin();
fieldToSpwMap.clear();
spwToFieldMap.resize(nSpw(True));
- std::map<Int, uInt> ddidToSpwMap = _getDataDescIDToSpwMap();
+ vector<uInt> ddidToSpwMap = _getDataDescIDToSpwMap();
for (
Vector<Int>::const_iterator curDDID=allDDIDs->begin();
curDDID!=endDDID; ++curDDID, ++curField
@@ -1151,7 +1164,7 @@
std::map<ScanKey, std::set<uInt> > scanToDDIDMap;
vector<std::set<ScanKey> > ddIDToScanMap;
_getScansAndDDIDMaps(scanToDDIDMap, ddIDToScanMap);
- std::map<Int, uInt> ddToSpw = _getDataDescIDToSpwMap();
+ vector<uInt> ddToSpw = _getDataDescIDToSpwMap();
std::map<ScanKey, std::set<uInt> >::const_iterator iter =
scanToDDIDMap.begin();
std::map<ScanKey, std::set<uInt> >::const_iterator end =
scanToDDIDMap.end();
std::set<uInt>::const_iterator dIter;
@@ -1902,7 +1915,7 @@
Vector<Double>::const_iterator iIter = intervals.begin();
scanSpwToAverageIntervalMap.clear();
std::map<ScanKey, std::map<uInt, uInt> > counts;
- std::map<Int, uInt> dataDesIDToSpwMap = _getDataDescIDToSpwMap();
+ vector<uInt> dataDesIDToSpwMap = _getDataDescIDToSpwMap();
ScanKey scanKey;
//std::pair<ScanKey, uInt> scanSpwToAverageIntervalKey;
uInt mysize = 0;
@@ -1918,7 +1931,7 @@
scanToTimeRangeMap[scanKey].first =
min(scanToTimeRangeMap[scanKey].first, *tIter-half);
scanToTimeRangeMap[scanKey].second =
max(scanToTimeRangeMap[scanKey].second, *tIter+half);
}
- uInt spw = dataDesIDToSpwMap.find(*dIter)->second;
+ uInt spw = dataDesIDToSpwMap[*dIter];
if (
(
scanSpwToAverageIntervalMap.find(scanKey)
@@ -2207,18 +2220,20 @@
Record polTable;
polTable.define("n correlations", Vector<Int>(_getNumCorrs()));
Record dataDescTable;
- std::map<Int, uInt> ddToSpw = _getDataDescIDToSpwMap();
- std::map<Int, uInt> ddToPolID = _getDataDescIDToPolIDMap();
- std::map<Int, uInt>::const_iterator siter = ddToSpw.begin();
- std::map<Int, uInt>::const_iterator send = ddToSpw.end();
- std::map<Int, uInt>::const_iterator piter = ddToPolID.begin();
+ vector<uInt> ddToSpw = _getDataDescIDToSpwMap();
+ vector<uInt> ddToPolID = getDataDescIDToPolIDMap();
+ vector<uInt>::const_iterator siter = ddToSpw.begin();
+ vector<uInt>::const_iterator send = ddToSpw.end();
+ vector<uInt>::const_iterator piter = ddToPolID.begin();
vector<Int> spws(ddToSpw.size());
vector<Int> polids(ddToPolID.size());
+ uInt ddid = 0;
while (siter != send) {
- spws[siter->first] = siter->second;
- polids[piter->first] = piter->second;
+ spws[ddid] = *siter;
+ polids[ddid] = *piter;
++siter;
++piter;
+ ++ddid;
}
dataDescTable.define("spectral windows", Vector<Int>(spws));
dataDescTable.define("polarization ids", Vector<Int>(polids));
@@ -2771,14 +2786,16 @@
Vector<Double> times = ScalarColumn<Double>(result, "TIME").getColumn();
// each row represents a unique baseline, data description ID, and time
combination
uInt nrows = result.nrow();
- std::map<Int, uInt> dataDescToSpwIdMap = _getDataDescIDToSpwMap();
+ vector<uInt> dataDescToSpwIdMap = _getDataDescIDToSpwMap();
std::set<uInt> avgSpw, tdmSpw, fdmSpw, wvrSpw, sqldSpw;
vector<SpwProperties> spwInfo = _getSpwInfo(avgSpw, tdmSpw, fdmSpw,
wvrSpw, sqldSpw);
std::map<Double, Double> timeToBWMap = _getTimeToTotalBWMap(
times, ddIDs
);
for (uInt i=0; i<nrows; ++i) {
- Quantum<Vector<Double> > channelWidths =
spwInfo[dataDescToSpwIdMap.find(ddIDs[i])->second].chanwidths;
+ uInt ddID = ddIDs[i];
+ uInt spw = dataDescToSpwIdMap[ddID];
+ Quantum<Vector<Double> > channelWidths = spwInfo[spw].chanwidths;
Matrix<Bool> flagsMatrix(ArrayColumn<Bool>(result, "FLAG").get(i));
uInt nCorrelations = flagsMatrix.nrow();
Double denom =
(timeToBWMap.find(times[i])->second)*maxNBaselines*nCorrelations;
@@ -2892,7 +2909,7 @@
++dIter;
}
std::map<Double, std::set<uInt> >::const_iterator end1 =
timeToDDIDMap.end();
- std::map<Int, uInt> dataDescIDToSpwMap = _getDataDescIDToSpwMap();
+ vector<uInt> dataDescIDToSpwMap = _getDataDescIDToSpwMap();
std::set<uInt> avgSpw, tdmSpw, fdmSpw, wvrSpw, sqldSpw;
vector<SpwProperties> spwInfo = _getSpwInfo(avgSpw, tdmSpw, fdmSpw,
wvrSpw, sqldSpw);
for (
@@ -2906,7 +2923,8 @@
std::set<uInt>::const_iterator dIter=ddIDs.begin();
dIter!=end2; ++dIter
) {
- timeToBWMap[iter->first] +=
spwInfo[dataDescIDToSpwMap.find(*dIter)->second].bandwidth;
+ uInt spw = dataDescIDToSpwMap[*dIter];
+ timeToBWMap[iter->first] += spwInfo[spw].bandwidth;
}
}
return timeToBWMap;
@@ -2995,12 +3013,13 @@
uInt i = 0;
//uInt64 count = 0;
// a flag value of True means the datum is bad (flagged), so False =>
unflagged
- std::map<Int, uInt> dataDescIDToSpwMap = _getDataDescIDToSpwMap();
+ vector<uInt> dataDescIDToSpwMap = _getDataDescIDToSpwMap();
std::set<uInt> a, b, c, d, e;
vector<SpwProperties> spwInfo = _getSpwInfo(a, b, c, d, e);
CountedPtr<ArrayColumn<Bool> > flags = _getFlags();
while (a1Iter != aEnd) {
- SpwProperties spwProp = spwInfo[dataDescIDToSpwMap.find(*dIter)->second];
+ uInt spw = dataDescIDToSpwMap[*dIter];
+ SpwProperties spwProp = spwInfo[spw];
Vector<Double> channelWidths(
Vector<Double>(spwProp.chanwidths.getValue("Hz"))
);
@@ -3114,7 +3133,7 @@
Vector<Int>::const_iterator endDDID = dataDescIDs->end();
CountedPtr<Vector<Int> > states = _getStateIDs();
Vector<Int>::const_iterator curState = states->begin();
- std::map<Int, uInt> dataDescToSpwMap = _getDataDescIDToSpwMap();
+ vector<uInt> dataDescToSpwMap = _getDataDescIDToSpwMap();
while (curDDID!=endDDID) {
uInt spw = dataDescToSpwMap[*curDDID];
std::set<String> intents = stateToIntentsMap[*curState];
@@ -3285,21 +3304,22 @@
}
}

-std::map<std::pair<uInt, uInt>, Int>
MSMetaData::getSpwIDPolIDToDataDescIDMap() {
+std::map<std::pair<uInt, uInt>, uInt>
MSMetaData::getSpwIDPolIDToDataDescIDMap() const {
if (! _spwPolIDToDataDescIDMap.empty()) {
return _spwPolIDToDataDescIDMap;
}
- std::map<Int, uInt> dataDescIDToSpwMap = _getDataDescIDToSpwMap();
- std::map<Int, uInt>::const_iterator i1 = dataDescIDToSpwMap.begin();
- std::map<Int, uInt>::const_iterator end = dataDescIDToSpwMap.end();
- std::map<std::pair<uInt, uInt>, Int> spwPolIDToDataDescIDMap;
- std::map<Int, uInt> dataDescIDToPolIDMap = _getDataDescIDToPolIDMap();
+ vector<uInt> dataDescIDToSpwMap = _getDataDescIDToSpwMap();
+ vector<uInt>::const_iterator i1 = dataDescIDToSpwMap.begin();
+ vector<uInt>::const_iterator end = dataDescIDToSpwMap.end();
+ std::map<std::pair<uInt, uInt>, uInt> spwPolIDToDataDescIDMap;
+ vector<uInt> dataDescIDToPolIDMap = getDataDescIDToPolIDMap();
+ uInt dataDesc = 0;
while (i1 != end) {
- Int dataDesc = i1->first;
- uInt spw = i1->second;
- uInt polID = dataDescIDToPolIDMap.at(dataDesc);
+ uInt spw = *i1;
+ uInt polID = dataDescIDToPolIDMap[dataDesc];
spwPolIDToDataDescIDMap[std::make_pair(spw, polID)] = dataDesc;
++i1;
+ ++dataDesc;
}
uInt mysize = 2*sizeof(Int)*spwPolIDToDataDescIDMap.size();
if (_cacheUpdated(mysize)) {
@@ -3341,7 +3361,7 @@
else {
dir2 = _getInterpolatedDirection(pCols, pidx2, time);
}
- return std::make_pair<MDirection, MDirection>(dir1, dir2);
+ return std::make_pair(dir1, dir2);
}

MDirection MSMetaData::_getInterpolatedDirection(
@@ -3388,13 +3408,14 @@
return MDirection(qlon, qlat, rf);
}

-std::map<Int, uInt> MSMetaData::_getDataDescIDToSpwMap() const {
+vector<uInt> MSMetaData::_getDataDescIDToSpwMap() const {
if (! _dataDescIDToSpwMap.empty()) {
return _dataDescIDToSpwMap;
}
String spwColName =
MSDataDescription::columnName(MSDataDescriptionEnums::SPECTRAL_WINDOW_ID);
ROScalarColumn<Int> spwCol(_ms->dataDescription(), spwColName);
- std::map<Int, uInt> dataDescToSpwMap = _toUIntMap(spwCol.getColumn());
+ Vector<Int> spws = spwCol.getColumn();
+ vector<uInt> dataDescToSpwMap(spws.begin(), spws.end());
uInt mysize = sizeof(Int) * dataDescToSpwMap.size();
if (_cacheUpdated(mysize)) {
_dataDescIDToSpwMap = dataDescToSpwMap;
@@ -3413,8 +3434,8 @@
if (! _scanSpwToPolIDMap.empty()) {
return _scanSpwToPolIDMap.find(std::pair<ScanKey, uInt>(scanKey,
spwid))->second;
}
- std::map<Int, uInt> ddToPolMap = _getDataDescIDToPolIDMap();
- std::map<Int, uInt> ddToSpwMap = _getDataDescIDToSpwMap();
+ vector<uInt> ddToPolMap = getDataDescIDToPolIDMap();
+ vector<uInt> ddToSpwMap = _getDataDescIDToSpwMap();
std::map<ScanKey, std::set<uInt> > scanToDDIDMap;
vector<std::set<ScanKey> > ddIDToScanMap;
_getScansAndDDIDMaps(scanToDDIDMap, ddIDToScanMap);
@@ -3451,13 +3472,14 @@
return size;
}

-std::map<Int, uInt> MSMetaData::_getDataDescIDToPolIDMap() const {
+vector<uInt> MSMetaData::getDataDescIDToPolIDMap() const {
if (! _dataDescIDToPolIDMap.empty()) {
return _dataDescIDToPolIDMap;
}
- String spwColName =
MSDataDescription::columnName(MSDataDescriptionEnums::POLARIZATION_ID);
- ROScalarColumn<Int> spwCol(_ms->dataDescription(), spwColName);
- std::map<Int, uInt> dataDescToPolIDMap = _toUIntMap(spwCol.getColumn());
+ String polColName =
MSDataDescription::columnName(MSDataDescriptionEnums::POLARIZATION_ID);
+ ROScalarColumn<Int> polCol(_ms->dataDescription(), polColName);
+ Vector<Int> pols = polCol.getColumn();
+ vector<uInt> dataDescToPolIDMap(pols.begin(), pols.end());
uInt mysize = sizeof(Int) * dataDescToPolIDMap.size();
if (_cacheUpdated(mysize)) {
_dataDescIDToPolIDMap = dataDescToPolIDMap;
=======================================
--- /branches/nov14/ms/MeasurementSets/MSMetaData.h Thu Mar 12 10:32:41
2015 UTC
+++ /branches/nov14/ms/MeasurementSets/MSMetaData.h Wed Mar 18 15:01:43
2015 UTC
@@ -97,6 +97,8 @@
// get the antenna stations for the specified antenna names
vector<String> getAntennaStations(const vector<String>& antennaNames);

+ vector<uInt> getDataDescIDToPolIDMap() const;
+
std::map<String, std::set<Int> > getIntentToFieldsMap();

std::map<String, std::set<ScanKey> > getIntentToScansMap();
@@ -146,6 +148,9 @@
// get the number of fields.
uInt nFields() const;

+ // get a mapping of spectral window ID to data descrption IDs
+ std::vector<std::set<uInt> > getSpwToDataDescriptionIDMap() const;
+
// get a set of spectral windows corresponding to the specified fieldID
std::set<uInt> getSpwsForField(const Int fieldID);

@@ -158,17 +163,6 @@
// get the set of field names corresponding to the specified spectral
window.
std::set<String> getFieldNamesForSpw(const uInt spw);

- // get the number of scans in the dataset
- uInt nScans();
-
- // get the number of observations (from the OBSERVATIONS table) in the
dataset
- uInt nObservations() const;
-
- // get the number of arrays (from the ARRAY table) in the dataset
- uInt nArrays();
-
- // get the number of data description IDs (from the DATA_DESCRIPTION
table)
- uInt nDataDescriptions() const;

// get the set of spectral windows for the specified scan number.
std::set<uInt> getSpwsForScan(const ScanKey& scan) const;
@@ -179,8 +173,6 @@
// get the number of antennas in the ANTENNA table
uInt nAntennas() const;

-
-
// ALMA-specific. get set of spectral windows used for TDM. These are
windows that have
// 64, 128, or 256 channels
std::set<uInt> getTDMSpw();
@@ -297,7 +289,7 @@
std::set<Int> getScansForField(const String& field, Int obsID, Int
arrayID) const;

// The first value of the pair is spw, the second is polarization ID.
- std::map<std::pair<uInt, uInt>, Int> getSpwIDPolIDToDataDescIDMap();
+ std::map<std::pair<uInt, uInt>, uInt> getSpwIDPolIDToDataDescIDMap()
const;

vector<String> getSpwNames() const;

@@ -327,6 +319,18 @@
// get the effective total exposure time. This is the effective time
spent collecting unflagged data.
Quantity getEffectiveTotalExposureTime();

+ // get the number of scans in the dataset
+ uInt nScans();
+
+ // get the number of observations (from the OBSERVATIONS table) in the
dataset
+ uInt nObservations() const;
+
+ // get the number of arrays (from the ARRAY table) in the dataset
+ uInt nArrays();
+
+ // get the number of data description IDs (from the DATA_DESCRIPTION
table)
+ uInt nDataDescriptions() const;
+
// get the number of unflagged rows
Double nUnflaggedRows() const;

@@ -417,12 +421,12 @@
mutable uInt _nStates, _nACRows, _nXCRows, _nSpw, _nFields, _nAntennas,
_nObservations, _nScans, _nArrays, _nrows, _nPol, _nDataDescIDs;
mutable std::map<ScanKey, std::set<uInt> > _scanToSpwsMap,
_scanToDDIDsMap;
- mutable std::map<Int, uInt> _dataDescIDToSpwMap, _dataDescIDToPolIDMap
/*, _scanToNRowsMap */;
+ mutable vector<uInt> _dataDescIDToSpwMap, _dataDescIDToPolIDMap /*,
_scanToNRowsMap */;
std::map<Int, std::set<uInt> > _fieldToSpwMap;
mutable std::map<ScanKey, std::set<Int> > _scanToStatesMap,
_scanToFieldsMap;
mutable std::map<Int, std::set<Int> > _fieldToStatesMap,
_stateToFieldsMap, _sourceToFieldsMap /*,
_arrayToScansMap, _scanToAntennasMap */;
- std::map<std::pair<uInt, uInt>, Int> _spwPolIDToDataDescIDMap;
+ mutable std::map<std::pair<uInt, uInt>, uInt> _spwPolIDToDataDescIDMap;
std::map<String, uInt> _antennaNameToIDMap;
//mutable std::map<ArrayKey, ArrayProperties> _arrayProperties;
//mutable std::map<ScanKey, ScanProperties> _scanProperties;
@@ -550,9 +554,7 @@

//std::map<Int, std::set<Int> > _getArrayIDToScansMap() const;

- std::map<Int, uInt> _getDataDescIDToSpwMap() const;
-
- std::map<Int, uInt> _getDataDescIDToPolIDMap() const;
+ vector<uInt> _getDataDescIDToSpwMap() const;

void _getFieldsAndIntentsMaps(
vector<std::set<String> >& fieldToIntentsMap,
=======================================
--- /branches/nov14/ms/MeasurementSets/MSSummary.cc Tue Mar 3 12:22:11
2015 UTC
+++ /branches/nov14/ms/MeasurementSets/MSSummary.cc Wed Mar 18 15:01:43
2015 UTC
@@ -1591,23 +1591,14 @@
// Line is (V1): RefFreq RestFreq Molecule Trans'n Resol BW Numch Correls
// V2 subtable: SOURCE SOURCE POLARIZ'N

- // Define the column widths
Int widthLead = 2;
Int widthFreq = 12;
Int widthFrqNum = 7;
- //Int widthMol = 10;
- //Int widthTrans = 12;
Int widthNumChan = 8;
- // Int widthCorrTypes = msPolC.corrType()(0).nelements()*4;
- // Int widthCorrType = 4;

- // Write the column headers
os.output().setf(ios::left, ios::adjustfield);
os.output().width(widthLead); os << " ";
os.output().width(widthFreq); os << "Ref.Freq";
- //os.output().width(widthFreq); os << "RestFreq";
- //os.output().width(widthMol); os << "Molecule";
- //os.output().width(widthTrans); os << "Transition";
os.output().width(widthNumChan); os << "#Chans";
os.output().width(widthFreq); os << "Resolution";
os.output().width(widthFreq); os << "TotalBW";
@@ -1774,7 +1765,6 @@

// determine the data_desc_ids present in the main table
MSRange msr(*pMS);
- //ant2 = msr.range(MSS::ANTENNA2).asArrayInt(RecordFieldId(0));
Vector<Int> ddId =
msr.range(MSS::DATA_DESC_ID).asArrayInt(RecordFieldId(0));
Vector<uInt> uddId(ddId.nelements());
for (uInt i=0; i<ddId.nelements(); i++) uddId(i)=ddId(i);
@@ -1783,8 +1773,8 @@
Vector<Int> polIds = msDDC.polarizationId().getColumnCells(uddId);
const Int option=Sort::HeapSort | Sort::NoDuplicates;
const Sort::Order order=Sort::Ascending;
- Int nSpw=GenSort<Int>::sort (spwIds, order, option);
- Int nPol=GenSort<Int>::sort (polIds, order, option);
+ uInt nSpw = GenSort<Int>::sort (spwIds, order, option);
+ Int nPol = GenSort<Int>::sort (polIds, order, option);

if (ddId.nelements()>0) {
os << "Spectral Windows: ";
@@ -1811,8 +1801,6 @@
Int widthCorrType = 4;
uInt widthBBCNo = 8;

-
-
// Write the column headers
os.output().setf(ios::left, ios::adjustfield);
os.output().width(widthLead); os << " ";
@@ -1843,56 +1831,73 @@
vector<uInt> bbcNo = hasBBCNo ? _msmd->getBBCNos() : vector<uInt>();

os.output().precision(9);
- // For each row of the DataDesc subtable, write the info
- for (uInt i=0; i<ddId.nelements(); i++) {
- Int dd=ddId(i);
- Int spw = msDDC.spectralWindowId()(dd);
- Int pol = msDDC.polarizationId()(dd);
- os.output().setf(ios::left, ios::adjustfield);
- os.output().width(widthLead); os << " ";
- // 1th column: Spectral Window Id
- os.output().width(widthSpwId); os << (spw);
- // 2nd column: SPW name
- os.output().width(widthName);
- os << names[spw] << " ";
+ // order by spwid, not ddid, CAS-7376
+ Vector<Int>::const_iterator iter = spwIds.begin();
+ Vector<Int>::const_iterator end = spwIds.end();
+ std::vector<std::set<uInt> > spwToDDID =
_msmd->getSpwToDataDescriptionIDMap();
+ std::set<uInt> uDDIDSet(uddId.begin(), uddId.end());
+ vector<uInt> ddToPolID = _msmd->getDataDescIDToPolIDMap();
+ while (iter != end) {
+ // for (uInt i=0; i<ddId.nelements(); i++) {
+ Int spw = *iter;
+ std::set<uInt> ddids = spwToDDID[spw];
+ std::set<uInt>::const_iterator diter = ddids.begin();
+ std::set<uInt>::const_iterator dend = ddids.end();
+ while (diter != dend) {
+ uInt dd = *diter;
+ if (uDDIDSet.find(dd) == uDDIDSet.end()) {
+ // data description ID not in main table, so not reported here
+ continue;
+ }
+ Int pol = ddToPolID[dd];
+ os.output().setf(ios::left, ios::adjustfield);
+ os.output().width(widthLead); os << " ";
+ // 1th column: Spectral Window Id
+ os.output().width(widthSpwId); os << (spw);
+ // 2nd column: SPW name
+ os.output().width(widthName);
+ os << names[spw] << " ";

- // 3rd column: number of channels in the spectral window
- os.output().setf(ios::right, ios::adjustfield);
- os.output().width(widthNumChan);
- os << nChans[spw] << " ";
- // 4th column: Reference Frame info
- // os.output().setf(ios::left, ios::adjustfield);
- os.output().width(widthFrame);
- os<< msSWC.refFrequencyMeas()(spw).getRefString();
- // 5th column: Chan 1 freq (may be at high freq end of band!)
- os.output().setf(ios::fixed);
- os.output().precision(3);
- os.output().width(widthFrqNum);
- os<< chanFreqs[spw].getValue("MHz")[0];
- // 6th column: channel resolution
- os.output().width(widthFrqNum+2);
- os << chanWidths[spw].getValue("kHz")[0];
- // 7th column: total bandwidth of the spectral window
- os.output().width(widthFrqNum);
- os.output().precision(1);
- os<< bandwidths[spw]/1000;
- os.output().width(widthFrqNum);
- os.output().precision(4);
- os << centerFreqs[spw].getValue("MHz") << " ";
- if (hasBBCNo) {
- os.output().width(widthBBCNo);
- os<< bbcNo[spw];
- }
- // 8th column: reference frequency
- // os.output().width(widthFrqNum);
- // os<< msSWC.refFrequency()(spw)/1.0e6;
- // 9th column: the correlation type(s)
- for (uInt j=0; j<msPolC.corrType()(pol).nelements(); j++) {
- os.output().width(widthCorrType);
- Int index = msPolC.corrType()(pol)(IPosition(1,j));
- os << Stokes::name(Stokes::type(index));
+ // 3rd column: number of channels in the spectral window
+ os.output().setf(ios::right, ios::adjustfield);
+ os.output().width(widthNumChan);
+ os << nChans[spw] << " ";
+ // 4th column: Reference Frame info
+ // os.output().setf(ios::left, ios::adjustfield);
+ os.output().width(widthFrame);
+ os<< msSWC.refFrequencyMeas()(spw).getRefString();
+ // 5th column: Chan 1 freq (may be at high freq end of band!)
+ os.output().setf(ios::fixed);
+ os.output().precision(3);
+ os.output().width(widthFrqNum);
+ os<< chanFreqs[spw].getValue("MHz")[0];
+ // 6th column: channel resolution
+ os.output().width(widthFrqNum+2);
+ os << chanWidths[spw].getValue("kHz")[0];
+ // 7th column: total bandwidth of the spectral window
+ os.output().width(widthFrqNum);
+ os.output().precision(1);
+ os<< bandwidths[spw]/1000;
+ os.output().width(widthFrqNum);
+ os.output().precision(4);
+ os << centerFreqs[spw].getValue("MHz") << " ";
+ if (hasBBCNo) {
+ os.output().width(widthBBCNo);
+ os<< bbcNo[spw];
+ }
+ // 8th column: reference frequency
+ // os.output().width(widthFrqNum);
+ // os<< msSWC.refFrequency()(spw)/1.0e6;
+ // 9th column: the correlation type(s)
+ for (uInt j=0; j<msPolC.corrType()(pol).nelements(); j++) {
+ os.output().width(widthCorrType);
+ Int index = msPolC.corrType()(pol)(IPosition(1,j));
+ os << Stokes::name(Stokes::type(index));
+ }
+ ++diter;
}
os << endl;
+ ++iter;
}
}
os << LogIO::POST;
=======================================
--- /branches/nov14/ms/MeasurementSets/test/tMSMetaData.cc Thu Mar 12
10:32:41 2015 UTC
+++ /branches/nov14/ms/MeasurementSets/test/tMSMetaData.cc Wed Mar 18
15:01:43 2015 UTC
@@ -47,7 +47,7 @@
const std::set<uInt>::const_iterator end = set.end();
for (
std::set<uInt>::const_iterator iter = set.begin();
- iter!=end; iter++
+ iter!=end; ++iter
) {
if (iter!=set.begin()) {
cout << ", ";
@@ -61,7 +61,7 @@
const std::set<String>::const_iterator end = set.end();
for (
std::set<String>::const_iterator iter = set.begin();
- iter!=end; iter++
+ iter!=end; ++iter
) {
if (iter!=set.begin()) {
cout << ", ";
@@ -80,7 +80,7 @@
cout << "*** cache size " << md.getCache() << endl;

cout << "*** test getScansForState()" << endl;
- for (uInt stateID=0; stateID<md.nStates(); stateID++) {
+ for (uInt stateID=0; stateID<md.nStates(); ++stateID) {
std::set<Int> scans = md.getScansForState(stateID, 0, 0);
std::set<Int> expec;
if (stateID < 5) {
@@ -142,7 +142,7 @@
scanKey.arrayID = 0;
for (
std::set<Int>::const_iterator scanNum = scans.begin();
- scanNum!=scans.end(); scanNum++
+ scanNum!=scans.end(); ++scanNum
) {
scanKey.scan = *scanNum;
std::set<String> intents = md.getIntentsForScan(scanKey);
@@ -222,14 +222,14 @@
cout << "*** test getSpwsForIntent()" << endl;
for (
std::set<String>::const_iterator intent=uniqueIntents.begin();
- intent!=uniqueIntents.end(); intent++
+ intent!=uniqueIntents.end(); ++intent
) {
std::set<uInt> exp;
if (
- *intent == "CALIBRATE_AMPLI#ON_SOURCE"
- || *intent == "CALIBRATE_BANDPASS#ON_SOURCE"
- || *intent == "CALIBRATE_PHASE#ON_SOURCE"
- || *intent == "OBSERVE_TARGET#ON_SOURCE"
+ *intent == "CALIBRATE_AMPLI#ON_SOURCE"
+ || *intent == "CALIBRATE_BANDPASS#ON_SOURCE"
+ || *intent == "CALIBRATE_PHASE#ON_SOURCE"
+ || *intent == "OBSERVE_TARGET#ON_SOURCE"
) {
uInt myints[] = {0, 17, 18, 19, 20, 21, 22, 23, 24};
exp.insert(myints, myints+9);
@@ -267,7 +267,7 @@
AlwaysAssert(nSpw == 40, AipsError);
AlwaysAssert(md.nSpw(False) == 24, AipsError);
cout << "*** test getIntentsForSpw()" << endl;
- for (uInt spw=0; spw<nSpw; spw++) {
+ for (uInt spw=0; spw<nSpw; ++spw) {
std::set<String> exp;
if (spw == 0) {
String mystr[] = {
@@ -324,7 +324,7 @@
"3C279", "J1337-129", "Titan",
"J1625-254", "V866 Sco", "RNO 90"
};
- for (uInt i=0; i<nFields; i++) {
+ for (uInt i=0; i<nFields; ++i) {
std::set<uInt> exp;
if (i==0 || i==3) {
uInt myints[] = {
@@ -354,7 +354,7 @@

}
cout << "*** test getFieldIDsForSpw()" << endl;
- for (uInt i=0; i<md.nSpw(True); i++) {
+ for (uInt i=0; i<md.nSpw(True); ++i) {
std::set<Int> exp;
std::set<String> expNames;
if (i==0) {
@@ -403,7 +403,7 @@
scanKey.arrayID = 0;
for (
std::set<Int>::const_iterator scan=scanNumbers.begin();
- scan!=scanNumbers.end(); scan++
+ scan!=scanNumbers.end(); ++scan
) {
std::set<uInt> exp;
if (*scan == 1 || *scan==5 || *scan==8) {
@@ -438,7 +438,7 @@
AlwaysAssert(md.getSpwsForScan(scanKey) == exp, AipsError);
for (
std::set<uInt>::const_iterator spw=exp.begin();
- spw!=exp.end(); spw++
+ spw!=exp.end(); ++spw
) {
std::set<uInt> exppols;
std::set<uInt> pols = md.getPolarizationIDs(0, 0, *scan, *spw);
@@ -452,7 +452,7 @@
}
}
cout << "*** test getScansForSpw()" << endl;
- for (uInt i=0; i<md.nSpw(True); i++) {
+ for (uInt i=0; i<md.nSpw(True); ++i) {
std::set<Int> exp;
if (i==0) {
Int myints[] = {
@@ -496,7 +496,7 @@
"DV07", "DV08", "DV10", "DV12", "DV13",
"DV14", "DV15", "DV16", "DV17", "DV18"
};
- for (uInt i=0; i<md.nAntennas(); i++) {
+ for (uInt i=0; i<md.nAntennas(); ++i) {
vector<uInt> ids(1);
ids[0] = i;
std::map<String, uInt> mymap;
@@ -507,7 +507,7 @@
}
cout << "*** test getAntennaID()" << endl;
std::map<String, uInt> mymap;
- for (uInt i=0; i<md.nAntennas(); i++) {
+ for (uInt i=0; i<md.nAntennas(); ++i) {
vector<uInt> ids(1);
ids[0] = i;
AlwaysAssert(
@@ -626,7 +626,7 @@
std::set<Int> scanNumbers = md.getScanNumbers(0, 0);
for (
std::set<Int>::const_iterator curScan=scanNumbers.begin();
- curScan!=scanNumbers.end(); curScan++
+ curScan!=scanNumbers.end(); ++curScan
) {
expec.clear();
if (*curScan == 1 || *curScan == 5 || *curScan == 8) {
@@ -688,7 +688,7 @@
std::set<String> intents = md.getIntents();
for (
std::set<String>::const_iterator intent=intents.begin();
- intent!=intents.end(); intent++
+ intent!=intents.end(); ++intent
) {
std::set<Int> expec;
if (
@@ -754,7 +754,7 @@
{
cout << "*** test getScansForFieldID()" << endl;
std::set<Int> expec;
- for (uInt i=0; i<6; i++) {
+ for (uInt i=0; i<6; ++i) {
expec.clear();
switch(i) {
case 0:
@@ -800,7 +800,7 @@
}
{
cout << "*** test getFieldIDsForField()" << endl;
- for (uInt i=0; i<6; i++) {
+ for (uInt i=0; i<6; ++i) {
std::set<Int> expec;
expec.insert(i);
String name = i == 0 ? "3C279"
@@ -817,7 +817,7 @@
}
{
cout << "*** test getScansForField()" << endl;
- for (uInt i=0; i<6; i++) {
+ for (uInt i=0; i<6; ++i) {
std::set<Int> expec;
String name;
switch(i) {
@@ -876,7 +876,7 @@
std::set<Int> curScanSet;
for (
std::set<Int>::const_iterator curScan=scans.begin();
- curScan!=scans.end(); curScan++
+ curScan!=scans.end(); ++curScan
) {
std::set<Int> expec;
curScanSet.insert(*curScan);
@@ -935,7 +935,7 @@
std::set<String> intents = md.getIntents();
for (
std::set<String>::const_iterator intent=intents.begin();
- intent!=intents.end(); intent++
+ intent!=intents.end(); ++intent
) {
std::set<Int> expec;
if (
@@ -987,7 +987,7 @@
}
{
cout << "*** test getFieldNamesForFieldIDs()" << endl;
- for (uInt i=0; i<md.nFields(); i++) {
+ for (uInt i=0; i<md.nFields(); ++i) {
String name;
switch(i) {
case 0:
@@ -1035,15 +1035,15 @@
{
cout << "*** test getTimesForField()" << endl;
uInt nfields = md.nFields();
- for (uInt i=0; i< nfields; i++) {
+ for (uInt i=0; i< nfields; ++i) {
std::set<Double> times = md.getTimesForField(i);
uInt expec = i == 0 ? 818
- : i == 1 ? 81
- : i == 2 ? 248
- : i == 3 ? 402
- : i == 4 ? 963
- : i == 5 ? 965
- : 0;
+ : i == 1 ? 81
+ : i == 2 ? 248
+ : i == 3 ? 402
+ : i == 4 ? 963
+ : i == 5 ? 965
+ : 0;
AlwaysAssert(md.getTimesForField(i).size() == expec, AipsError);
}
}
@@ -1110,7 +1110,7 @@

for (
std::map<Double, Double>::const_iterator
iter=exposuresForTimes.begin();
- iter!=exposuresForTimes.end(); iter++
+ iter!=exposuresForTimes.end(); ++iter
) {
cout << std::setprecision(11) << iter->first << ", " << iter->second
<< endl;
}
@@ -1120,7 +1120,7 @@
{
cout << "*** Test getIntentsForField()" << endl;
uInt nFields = md.nFields();
- for (uInt i=0; i<nFields; i++) {
+ for (uInt i=0; i<nFields; ++i) {
std::set<String> expec;
switch (i) {
case 0: {
@@ -1197,14 +1197,14 @@
}
{
cout << "*** test BBCNosToSpwMap()" << endl;
- for (uInt i=0; i<3; i++) {
+ for (uInt i=0; i<3; ++i) {
MSMetaData::SQLDSwitch sqldSwitch = i == 0 ? MSMetaData::SQLD_INCLUDE
: i == 1 ? MSMetaData::SQLD_EXCLUDE : MSMetaData::SQLD_ONLY;
std::map<uInt, std::set<uInt> > got = md.getBBCNosToSpwMap(sqldSwitch);
std::map<uInt, std::set<uInt> >::const_iterator end = got.end();
for (
std::map<uInt, std::set<uInt> >::const_iterator iter=got.begin();
- iter!=end; iter++
+ iter!=end; ++iter
) {
std::set<uInt> expec;
switch(iter->first) {
@@ -1275,12 +1275,12 @@
}
{
cout << "*** test getSpwIDPolIDToDataDescIDMap()" << endl;
- std::map<std::pair<uInt, uInt>, Int> dataDescToPolID =
md.getSpwIDPolIDToDataDescIDMap();
- std::map<std::pair<uInt, uInt>, Int>::const_iterator iter;
- std::map<std::pair<uInt, uInt>, Int>::const_iterator begin =
dataDescToPolID.begin();
- std::map<std::pair<uInt, uInt>, Int>::const_iterator end =
dataDescToPolID.end();
+ std::map<std::pair<uInt, uInt>, uInt> dataDescToPolID =
md.getSpwIDPolIDToDataDescIDMap();
+ std::map<std::pair<uInt, uInt>, uInt>::const_iterator iter;
+ std::map<std::pair<uInt, uInt>, uInt>::const_iterator begin =
dataDescToPolID.begin();
+ std::map<std::pair<uInt, uInt>, uInt>::const_iterator end =
dataDescToPolID.end();
for(
- iter=begin; iter!=end; iter++
+ iter=begin; iter!=end; ++iter
) {
std::pair<uInt, uInt> mypair = iter->first;
uInt spw = mypair.first;
@@ -1311,7 +1311,7 @@
{
cout << "*** test getUniqueFiedIDs()" << endl;
std::set<Int> expec;
- for (Int i=0; i<6; i++) {
+ for (Int i=0; i<6; ++i) {
expec.insert(i);
}
AlwaysAssert(md.getUniqueFiedIDs() == expec, AipsError);
@@ -1343,7 +1343,7 @@
187550000000.0, 187550000000.0
};
vector<Double> expec(mine, mine + 39);
- for (uInt i=0; i<40; i++) {
+ for (uInt i=0; i<40; ++i) {
AlwaysAssert(abs(centers[i].getValue("Hz")/mine[i] - 1) < 1e-8,
AipsError);
}
}
@@ -1357,7 +1357,7 @@
};
AlwaysAssert(res.size() == 6, AipsError);
AlwaysAssert(res2.size() == 6, AipsError);
- for (Int i=0; i<6; i++) {
+ for (Int i=0; i<6; ++i) {
AlwaysAssert(res[i].size() == 1 && *(res[i].begin()) == i, AipsError);
AlwaysAssert(
res2[i].size() == 1 && *(res2[i].begin()) == names[i], AipsError
@@ -1440,7 +1440,7 @@
exp = 1868;
}
AlwaysAssert(nTimes == exp, AipsError);
- intent++;
+ ++intent;
}
{
cout << "*** test getSummary()" << endl;