TVI Programming Guide
Transform Visibility Iterator (TVI) development guide
Concept
Implement a class derived from TransformingVi2 overloading the filling methods of each of the columns affected by the transformation. It is not necessary to implement accessors for the columns not affected by the transformation as the TransformingVi2 base class will delegate the filling of these columns to the input ViImplementation2 object.
The input data can be accessed by the TransformingVi2::getVii() method, which returns a pointer to the input ViImplementation2 object. Also it is possible to access the 'own' transformed data via the associated VisBuffer2 object with the method TransformingVi2::getVisBufferConst(), that returns a pointer to the method VisBuffer2 object associated to each TVI. This is necessary for derived transformed data like FLAG,WEIGHT and SIGMA:
TVI life-cycle
Constructor
The TransformingVI2 class constructor merely takes as input parameter the pointer to the input (inner) Visibility Iterator Implementation object (ViImplementation2).
Classes derived from TransformingVI2 can of course implement its own constructor which takes a set of parameters defining the transformation in addition to the input ViImplementation2 object pointer. In this case the input ViImplementation2 object pointer has to be pass on to the base class TransformingVI2 constructor.
It is necessary to explicitly initialize the associated Visbuffer2 object, as this initialization does not occur automatically in the TransformingVI2 base class. This can nevertheless be done with the setVisBuffer and createAttachedVisBuffer methods which are implemented in the TransformingVI2 base class. Example (from mstransform/TVI/FreqAxisTVI.cc):
FreqAxisTVI::FreqAxisTVI( ViImplementation2 * inputVii,
const Record &configuration):
TransformingVi2 (inputVii)
{
// Parse configuration with transformation parameters
parseConfiguration(configuration);
// Initialize associated VisBuffer2 object
setVisBuffer(createAttachedVisBuffer (VbPlain,VbRekeyable));
return;
}
The reasons for the explicit VisBuffer2 initialization in the derived classes are being question CAS-8220, but for the time being it is necessary to do it in the derived classes.
The constructor of the TVI will usually receive the configuration for that TVI. There is actually no enforcement on how this configuration should be done. Some TVIs use the Record class, which allows for complex set of configurations options. In the case of using Record for the configuration, it is a good practice not to use the same fields as other TVIs are already using. This would help tools like mstransform to use the same Record to configure all the TVIs and letting the TVIs parse only the relevant configuration.
Destructor (Deprecated)
This is a cumbersome area as the ownership of the input ViImplementation2 object is not clearly defined. As a matter of fact the base class TransformingVI2 tries to delete the input ViImplementation2 object pointer at destruction time, but this can lead to segmentation faults it the input ViImplementation2 object is already deleted at the application layer.
This problematic behavior has been reported in CAS-8220, but for now my recommendation is to simply set the copy of the in ViImplementation2 object pointer to NULL at destruction time so that when the base class TransformingVI2 destructor tries to delete it there is no effect. Example (from mstransform/TVI/FreqAxisTVI.cc):
FreqAxisTVI::~FreqAxisTVI()
{
inputVii_p = NULL;
return;
}
Associated factory (Deprecated)
Each TVI should have an associated factory class derived from ViFactory. This is necessary because VisibilityIterator2 class takes a ViFactory object at construction time when it is built on top of a given input ViImplementation2 object. Otherwise VisibilityIterator2 creates a 'plain' disk-access ViImplementation2 object inside.
The factory class constructor has to take as parameters at least the same parameters that the associated TVI takes at construction time (so ViImplementation2 and transformation configuration parameters). Then the method ViFactory::createVi() has to be overloaded to return a pointer to anew object of the associated TVI class. Example (from mstransform/TVI/ChannelAverageTVI.cc)
ChannelAverageTVIFactory::ChannelAverageTVIFactory (Record &configuration,
ViImplementation2 *inputVii
{
inputVii_p = inputVii;
configuration_p = configuration;
}
vi::ViImplementation2 * ChannelAverageTVIFactory::createVi() const
{
return new ChannelAverageTVI(inputVii_p,configuration_p);
}
Methods typically overloaded for the row-invariant data transforming TVIs
Direct transformation methods
These are the data/flag/weight cubes, i.e. the actual data
void floatData (Cube<Float> & vis) const;
void visibilityObserved (Cube<Complex> & vis) const;
void visibilityCorrected (Cube<Complex> & vis) const;
void visibilityModel (Cube<Complex> & vis) const;
void flag(Cube<Bool>& flagCube) const;
void weightSpectrum(Cube<Float> &weightSp) const;
void sigmaSpectrum (Cube<Float> &sigmaSp) const;
Here it is important to take into account that:
- The input cubes have to be reshaped to the output (transformed) shape.
- Ideally the data should be filled directly in the input parameter cubes to avoid unnecessary copies.
- Each method should only transform and fill its corresponding cube, and not make assumptions about what the application layer is going to be access (e.g.: some plotms commands only require the transformed flags)
Check for examples in mstransform/TVI/ChannelAverageTVI.cc.
Derived transformation methods
These methods provide aggregated information about data-quality cubes like flag and weight/sigma, and can be used for another operations downstream that work at the row level (for instance time average). They have to be derived from the corresponding channelized transformed data:
void flagRow (Vector<Bool> & flagRow) const;
void weight (Matrix<Float> & weight) const;
void sigma (Matrix<Float> & sigma) const;
The algorithms to produce this data are common for all transformations:
- flagRow: It is calculated as the logical AND of all channelized flags
- weight: It is calculated as the average of all channelized weights per polarization, but flagged samples don't contribute to the average. Exception: When all input samples are flagged the row-level weight is the average of all flagged input samples.
- sigma: It is calculated as the average of all channelized sigmas per polarization, but flagged samples don't contribute to the average. Exception: When all input samples are flagged the row-level sigma is the average of all flagged input samples.
Actually since the above described algorithms to obtain row-level flag/weight/sigma apply to all transformation regardless of their nature, it is possible to generalize them, and as a matter of fact these algorithms are available in a separated library (mstransform/TVI/UtilsTVI.h) that can be used by all TVIs.
void accumulateWeightCube (const Cube<Float> &weightCube,
const Cube<Bool> &flags,
Matrix<Float> &result)
Using the transformed weight spectrum and channelized flags calculates the average per polarization taking into account that when all samples are flagged the result is the average of all flagged input samples.
void accumulateFlagCube (const Cube<Bool> &flagCube,
Vector<Bool> &flagRow
Performs a logical and of all channelized flags (i.e. flagRow is only set to True when all channelized flags are set to True).
In order to provide the necessary transformed data for this methods the own/associated VisBuffer2 object can be accessed with the TransformingVi2::getVisBufferConst() method. It is important to use this approach in order to benefit from the buffering mechanism (e.g.: do not recalculate transformed channelized flags after providing flagRow). Example (from mstransform/TVI/FreqAxisTVI.cc):
void FreqAxisTVI::weight (Matrix<Float> & weight) const
{
// Get flags and weightSpectrum from own VisBuffer
const Cube<Bool> &flags = getVisBufferConst()->flagCube();
const Cube<Float> &weightSpectrum = getVisBufferConst()->weightSpectrum();
// Calculate output weight
accumulateWeightCube(weightSpectrum,flags,weight);
return;
}
Shape-defining methods
In order to define the output (transformed) shape, and initialize the associated VisBuffer2 object for each sub-chunk of data TransformingVI2 resorts to the following methods to determine the number of transformed rows, channels and correlations:
Vector<Int> getChannels (Double time,
Int frameOfReference,
Int spectralWindowId,
Int msId) const
Vector<Int> getCorrelations () const
Int nRows () const
Notice the following considerations:
- Take into account that if the transformation does not change the number of elements across one axis, then the corresponding method does not have to be reimplemented, for example nRows() for row-invariant transformations.
- For the specific case of getChannels, the input parameters can be ignored as they refer to frequency relabeling and are not relevant in the context of defining the output transformed shape. Instead the recommendation is to simply return a vector of increasing integers from 0 to nOutputChan-1. Example (from mstransform/TVI/FreqAxisTVI.cc):
Vector<Int> FreqAxisTVI::getChannels (Double,Int,Int spectralWindowId,Int) const
{
Vector<Int> ret(spwOutChanNumMap_p[spectralWindowId]);
for (uInt chanIdx = 0; chanIdx<spwOutChanNumMap_p[spectralWindowId];chanIdx++)
{
ret(chanIdx) = chanIdx;
}
return ret;
}
In general the shape-defining process is a bit cumbersome and as a matter of fact its efficiency is now under question in CAS-8220, but for the time being the mechanism works and has to be implemented as described above.
Navigation methods
For row-invariant (therefore chunk/buffer invariant) transformations it is only necessary to overload the sub-chunk (buffer) navigation methods in order to guarantee that the output buffer is properly initialized to the transformed shape.
void origin()
void next()
Notice the following considerations:
In the absence of a buffering mechanism at the TVI level (as it should be) the sub-chunk/buffer methods have to merely call the method from the base class TransformingVi2::configureNewSubchunk() which in turns resorts to the shape-defining methods described in the previous section to properly defined the associated VisBuffer2 object. Example (from mstransform/TVI/FreqAxisTVI.cc):
void FreqAxisTVI::origin()
{
// Drive underlying ViImplementation2
getVii()->origin();
// Synchronize own VisBuffer
configureNewSubchunk();
return;
}
As probably notice this implementation could be directly done in the TransformingVI class, and the situation is also being analyzed under CAS-8220, but for the time being the above methods have to be implemented in the derived classes.
TVI unit testing
(Unit) testing has been prioritized for CASA 4.6 and CASA 4.7, so C++ unit tests for each separated TVI are expected (and desirable). For this purpose there are of course several strategies, but typically all of them share the need to compare CASA containers, using a given tolerance and in case of discrepancy locate the offending positions. For this purpose the convenience library mstransform/TVI/test/TestUtilsTVI.h contains the following methods:
template <class T> Bool compareVector(const Char* column,
const Vector<T> &inp,
const Vector<T> &ref,
Float tolerance = FLT_EPSILON);
template <class T> Bool compareMatrix(const Char* column,
const Matrix<T> &inp,
const Matrix<T> &ref,
Float tolerance = FLT_EPSILON);
template <class T> Bool compareCube(const Char* column,
const Cube<T> &inp,
const Cube<T> &ref,
Float tolerance = FLT_EPSILON);
When these methods compare two containers, and find a discrepancy they report the offending positions and corresponding reference / test values, e.g.:
VisibilityCubeCorrected does not match in position (row,chan,corr)=(0,0,0)
test=(-0.11878,-0.186512) reference=(-0.0421126,-0.432511)
Also another method allows comparing a custom set of accessors between two iterators
Bool compareVisibilityIterators(VisibilityIterator2 &testTVI,
VisibilityIterator2 &refTVI,
VisBufferComponents2 &columns,
Float tolerance = FLT_EPSILON);
This is useful to compare a TVI iterator vs a plain iterator pointing to an already transformed file, or maybe vs an old implementation (for instance comparing vs MSTransformIterator in mstransform/TVI/test/tChannelAverageTVI.h):
// Create MSTransformIterator pointing to reference file
MSTransformIteratorFactory refFactory(configuration);
VisibilityIterator2 refTVI(refFactory);
// Generate TVI to test
ChannelAverageTVIFactory testFactory(configuration,inputVI);
VisibilityIterator2 testTVI(testFactory);
// Determine columns to check
VisBufferComponents2 columns;
columns += FlagCube;
columns += VisibilityCubeCorrected;
columns += WeightSpectrum;
columns += SigmaSpectrum;
// Compare
res = compareVisibilityIterators(testTVI,refTVI,columns,tolerance);
Notice that in order to unit-test a TVI it is necessary to provide an input 'plain' disk-access ViImplementation2. This can be done by simple creating a VisibilityIterator2 object, and the getting the inner ViImplementation object with the getImpl() method:
ViImplementation2 * getImpl() const;
Notes regarding usage and transformation of WEIGHT/SIGMA_SPECTRUM
According to the convention established from CASA 4.4, when transforming DATA the source of weights is SIGMA_SPECTRUM, whereas when transforming CORRECTED_DATA the source of weights is WEIGHT_SPECTRUM.
Please notice that despite of using 'WEIGHT' as the name of the CORRECTED_DATA weight column, it does not mean that WEIGHT also have to be used as the source of weights for the DATA column.
Take into account the correspondence of WEIGHT/SIGMA with variance (Var) in order to define exactly how weights have to be used in each algorithm:
$Sigma = \sqrt{Var}$
$Weight = \frac{1}{Var}$
For instance, time/channel average use the inverse variance as weight to calculate a weighted average, therefore when transforming CORRECTED_DATA the weights used are directly the values contained in WEIGHT_SPECTRUM, whereas when transforming DATA the weights used are the inverse squared of the values contained in SIGMA_SPECTRUM (the so famous weight = 1 / sigma*sigma equivalence)
However it is up to each transformation algorithm to actually use weights or not. For instance the channel/time averages use weights by default and therefore are weighted averages, but other transformations may not use it if they attempt to address some specific effect (like Gibbs ringing in Hanning smooth).
Nevertheless, regardless of the particular approach of each transformation algorithm w.r.t. using weights, each TVI should produce properly the transformed weights. For this error propagation rules apply:
$Sigma = \sqrt{Var}$
$Sigma_{out} =\sqrt{\sum_{inp=1}^{inp=n}{(\frac{\partial Vis_{out}}{\partial Vis_{inp}})}^2\times Sigma_{inp}^2}$
$Var_{out} =\sum_{inp=1}^{inp=n}{(\frac{\partial Vis_{out}}{\partial Vis_{inp}})}^2\times Var_{inp}$
$Weight = \frac{1}{Var}$
$Weight_{out} =\frac{1}{\sum_{inp=1}^{inp=n}{(\frac{\partial Vis_{out}}{\partial Vis_{inp}})}^2\times \frac{1}{Weight_{inp}}}$