Line data    Source code

       1             : //# SolveDataBuffer.h: A container for data and residuals for solving
       2             : //# Copyright (C) 2008
       3             : //# Associated Universities, Inc. Washington DC, USA.
       4             : //#
       5             : //# This library is free software; you can redistribute it and/or modify it
       6             : //# under the terms of the GNU Library General Public License as published by
       7             : //# the Free Software Foundation; either version 2 of the License, or (at your
       8             : //# option) any later version.
       9             : //#
      10             : //# This library is distributed in the hope that it will be useful, but WITHOUT
      11             : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
      12             : //# FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
      13             : //# License for more details.
      14             : //#
      15             : //# You should have received a copy of the GNU Library General Public License
      16             : //# along with this library; if not, write to the Free Software Foundation,
      17             : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
      18             : //#
      19             : //# Correspondence concerning AIPS++ should be addressed as follows:
      20             : //#        Internet email: casa-feedback@nrao.edu.
      21             : //#        Postal address: AIPS++ Project Office
      22             : //#                        National Radio Astronomy Observatory
      23             : //#                        520 Edgemont Road
      24             : //#                        Charlottesville, VA 22903-2475 USA
      25             : //#
      26             : //# $Id: VisBuffer.h,v 19.14 2006/02/28 04:48:58 mvoronko Exp $
      27             : 
      28             : #ifndef SYNTHESIS_SOLVEDATABUFFER_H
      29             : #define SYNTHESIS_SOLVEDATABUFFER_H
      30             : 
      31             : #include <casacore/casa/aips.h>
      32             : #include <msvis/MSVis/VisBuffer2.h>
      33             : #include <casacore/casa/Arrays/Array.h>
      34             : #include <casacore/casa/Arrays/Cube.h>
      35             : #include <casacore/casa/Arrays/Matrix.h>
      36             : #include <casacore/casa/Containers/Block.h>
      37             : namespace casa { //# NAMESPACE CASA - BEGIN
      38             : 
      39             : //#forward
      40             : 
      41             : //<summary>SolveDataBuffer is a container for VisBuffer2 data and related residual and differentiation results related to generic calibration solving </summary>
      42             : //
      43             : // <use visibility=export>
      44             : //
      45             : // <reviewed reviewer="" date="" tests="" demos="">
      46             : 
      47             : // <prerequisite>
      48             : //   <li> <linkto class="VisBuffer">VisBuffer</linkto>
      49             : // </prerequisite>
      50             : //
      51             : // <etymology>
      52             : // SolveDataBuffer is a Buffer for Solving Data
      53             : // </etymology>
      54             : //
      55             : //<synopsis>
      56             : // 
      57             : // 
      58             : //</synopsis>
      59             : 
      60             : //<todo>
      61             : // <li> write todo list
      62             : //</todo>
      63             : 
      64             : class SolveDataBuffer 
      65             : {
      66             : public:
      67             :   // Generic ctor
      68             :   SolveDataBuffer();
      69             : 
      70             :   // Create from a VisBuffer2
      71             :   SolveDataBuffer(const vi::VisBuffer2& vb);
      72             : 
      73             :   // Copy construct
      74             :   // current iteration (or reattach).
      75             :   SolveDataBuffer(const SolveDataBuffer& sdb);
      76             : 
      77             :   // Destructor 
      78             :   ~SolveDataBuffer(); 
      79             : 
      80             :   // Assignment
      81             :   SolveDataBuffer& operator=(const SolveDataBuffer& sdb);
      82             : 
      83             :   // Does SDB contain finite weight?
      84             :   bool Ok();
      85             : 
      86             :   // Divide corrected by model
      87             :   //  NB: disable for now, may not be needed...
      88             :   //void divideCorrByModel();
      89             : 
      90             :   // Apply amp-only or phase-only to data
      91             :   void enforceAPonData(const casacore::String& apmode);
      92             : 
      93             :   // Zero flagged weights, and optionally the cross-hand weights
      94             :   void enforceSolveWeights(const casacore::Bool pHandOnly);
      95             : 
      96             :   // Set the focus channel
      97             :   //  (forms references to focus-channel flag/data/model)
      98             :   void setFocusChan(const casacore::Int focusChan=-1);
      99             : 
     100             :   // Size/init/finalize the residuals workspaces
     101             :   void sizeResiduals(const casacore::Int& nPar,const casacore::Int& nDiff);
     102             :   void initResidWithModel();
     103             :   void finalizeResiduals();
     104             : 
     105             :   // Manage working weights
     106             :   void updateWorkingFlags();
     107             :   void updateWorkingWeights(casacore::Bool, casacore::Float clamp=0.0);
     108             : 
     109             :   // Delete the workspaces
     110             :   void cleanUp();
     111             : 
     112             :   // Return pol basis
     113             :   casacore::String polBasis() const;
     114             : 
     115             :   // VB2-like data access methods (mostly const)
     116    24772171 :   casacore::Int nRows() const { return vb_->nRows(); };
     117         184 :   casacore::Int nAntennas() const { return nAnt_; };  // stored on ctor
     118       19605 :   const casacore::Vector<int>& observationId() const { return vb_->observationId(); };
     119             :   const casacore::Vector<casacore::Int>& arrayId() const { return vb_->arrayId(); };
     120    26994957 :   const casacore::Vector<casacore::Int>& antenna1() const { return vb_->antenna1(); };
     121    26991641 :   const casacore::Vector<casacore::Int>& antenna2() const { return vb_->antenna2(); };
     122             :   const casacore::Vector<casacore::Int>& dataDescriptionIds() const { return vb_->dataDescriptionIds(); };
     123      182729 :   const casacore::Vector<casacore::Int>& spectralWindow() const { return vb_->spectralWindows(); };
     124       19605 :   const casacore::Vector<casacore::Int>& scan() const { return vb_->scan(); };
     125    75832154 :   const casacore::Vector<casacore::Double>& time() const { return vb_->time(); };
     126      177957 :   const casacore::Vector<casacore::Double>& timeCentroid() const { return vb_->timeCentroid(); };
     127       21263 :   const casacore::Vector<casacore::Int>& fieldId() const { return vb_->fieldId(); };
     128       66547 :   casacore::Int nChannels() const { return vb_->nChannels(); };
     129      737200 :   const casacore::Vector<casacore::Double>& freqs() const { return freqs_; };  // stored on ctor
     130       44268 :   casacore::Double centroidFreq() const { return centroidFreq_; };        // calc/stored on ctor
     131      504312 :   casacore::Int nCorrelations() const { return vb_->nCorrelations(); };
     132     1210272 :   const casacore::Cube<casacore::Complex>& visCubeModel() const { return vb_->visCubeModel(); };
     133     2092434 :   const casacore::Cube<casacore::Complex>& visCubeCorrected() const { return vb_->visCubeCorrected(); };
     134             : 
     135             :   // The feedPa
     136        1080 :   const casacore::Vector<casacore::Float>& feedPa() const { return feedPa_; };
     137             : 
     138             :   // These are not const, because we will generally amend them:
     139    46592879 :   casacore::Vector<casacore::Bool>& flagRow() {fR_.reference(vb_->flagRow());return fR_;};
     140     6376241 :   casacore::Cube<casacore::Bool>& flagCube() {fC_.reference(vb_->flagCube());return fC_;};
     141     2969547 :   casacore::Cube<casacore::Float>& weightSpectrum() {wS_.reference(vb_->weightSpectrum());return wS_;};
     142             : 
     143             :   // <group>
     144             : 
     145             :   // Access functions
     146             :   //
     147             :   // Access to focus-channel slices of the flags, data, and model
     148     1512200 :   casacore::Cube<casacore::Bool>& infocusFlagCube() { return infocusFlagCube_p; };
     149           0 :   const casacore::Cube<casacore::Bool>& infocusFlagCube() const {return infocusFlagCube_p;};
     150             : 
     151             :   casacore::Cube<casacore::Complex>& infocusVisCube() { return infocusVisCube_p; };
     152             :   const casacore::Cube<casacore::Complex>& infocusVisCube() const {return infocusVisCube_p;};
     153             : 
     154     1009345 :   casacore::Cube<casacore::Float>& infocusWtSpec() { return infocusWtSpec_p; };
     155      419869 :   const casacore::Cube<casacore::Float>& infocusWtSpec() const {return infocusWtSpec_p;};
     156             : 
     157             :   casacore::Cube<casacore::Complex>& infocusModelVisCube() { return infocusModelVisCube_p; };
     158             :   const casacore::Cube<casacore::Complex>& infocusModelVisCube() const {return infocusModelVisCube_p;};
     159             : 
     160             :   // Working weights and flags
     161             :   //  const versions return infocus versions if working versions are empty!
     162             :   casacore::Cube<casacore::Float>& workingWtSpec() { return workingWtSpec_p;};
     163      419869 :   const casacore::Cube<casacore::Float>& const_workingWtSpec() const { return (workingWtSpec_p.nelements()==0 ? infocusWtSpec() : workingWtSpec_p); };
     164           0 :   casacore::Cube<casacore::Bool>& workingFlagCube() { return workingFlagCube_p; };
     165     1259607 :   const casacore::Cube<casacore::Bool>& const_workingFlagCube() const { return (workingFlagCube_p.nelements()==0 ? 
     166      419869 :                                                                                 (residFlagCube_p.nelements()==0 ? 
     167      839738 :                                                                                  infocusFlagCube() : residFlagCube()) : workingFlagCube_p); };
     168             :                                                                                
     169             : 
     170             : 
     171             :   // Workspace for the residual visibilities
     172      749636 :   casacore::Cube<casacore::Complex>& residuals() { return residuals_p; };
     173             :   const casacore::Cube<casacore::Complex>& residuals() const {return residuals_p;};
     174             : 
     175             :   // Workspace for flags of the residuals
     176           0 :   casacore::Cube<casacore::Bool>& residFlagCube() { return residFlagCube_p; };
     177      419869 :   const casacore::Cube<casacore::Bool>& residFlagCube() const {return residFlagCube_p;};
     178             : 
     179             :   // Workspace for the differentiated residuals
     180      419869 :   casacore::Array<casacore::Complex>& diffResiduals() { return diffResiduals_p; };
     181             :   const casacore::Array<casacore::Complex>& diffResiduals() const {return diffResiduals_p;};
     182             : 
     183             :   //</group>
     184             : 
     185             :   // Print out data, weights, flags for debugging purposes
     186             :   void reportData();
     187             : 
     188             : protected:
     189             : 
     190             :   // Handle copy from the input VB2
     191             :   void initFromVB(const vi::VisBuffer2& vb);
     192             : 
     193             : 
     194             : private:
     195             : 
     196             :   // The underlying VisBuffer2
     197             :   vi::VisBuffer2* vb_;
     198             : 
     199             :   // The number of antennas 
     200             :   casacore::Int nAnt_;
     201             : 
     202             :   // The frequencies
     203             :   //  Currently, assumed uniform over rows
     204             :   casacore::Vector<double> freqs_;
     205             :   double centroidFreq_;
     206             : 
     207             :   // The correlation types
     208             :   //  Currently, assumed uniform over rows
     209             :   //  (private; used only by polBasis, currently)
     210             :   casacore::Vector<int> corrs_;
     211             : 
     212             :   // The feedPa
     213             :   //  Currently, assumed uniform (per antenna) over rows
     214             :   casacore::Vector<float> feedPa_;
     215             : 
     216             :     // Array reference objects for things in the vb we need to mess with
     217             :   casacore::Vector<casacore::Bool> fR_;
     218             :   casacore::Cube<casacore::Bool> fC_;
     219             :   casacore::Cube<casacore::Float> wS_;
     220             :   
     221             : 
     222             :   // The current in-focus channel
     223             :   casacore::Int focusChan_p;
     224             : 
     225             :   // actual storage for the data
     226             :   casacore::Cube<casacore::Bool> infocusFlagCube_p;
     227             :   casacore::Cube<casacore::Float> infocusWtSpec_p;
     228             :   casacore::Cube<casacore::Complex> infocusVisCube_p;
     229             :   casacore::Cube<casacore::Complex> infocusModelVisCube_p;
     230             : 
     231             :   casacore::Cube<casacore::Bool> workingFlagCube_p;
     232             :   casacore::Cube<casacore::Float> workingWtSpec_p;
     233             : 
     234             :   casacore::Cube<casacore::Complex> residuals_p;
     235             :   casacore::Cube<casacore::Bool> residFlagCube_p;
     236             :   casacore::Array<casacore::Complex> diffResiduals_p;
     237             : 
     238             : 
     239             : };
     240             : 
     241             : 
     242             : class SDBList 
     243             : {
     244             : public:
     245             : 
     246             :   // Construct empty list
     247             :   SDBList();
     248             : 
     249             :   // Destructor
     250             :   ~SDBList();
     251             : 
     252             :   // How many SDBs?
     253     2471939 :   casacore::Int nSDB() const { return nSDB_; };
     254             : 
     255             :   // Generate a new SDB from an input VB2
     256             :   void add(const vi::VisBuffer2& vb);
     257             : 
     258             :   // Access an SDB by index
     259             :   SolveDataBuffer& operator()(casacore::Int i);
     260             : 
     261             :   // Aggregate meta info
     262             :   int aggregateObsId() const;
     263             :   int aggregateScan() const;
     264             :   int aggregateSpw() const;
     265             :   int aggregateFld() const;
     266             :   double aggregateTime() const;
     267             :   double aggregateTimeCentroid() const;
     268             : 
     269             :   // Return pol basis
     270             :   casacore::String polBasis() const;
     271             : 
     272             :   // How many antennas 
     273             :   //   Currently, this insists on uniformity over all SDBs
     274             :   int nAntennas() const;
     275             : 
     276             :   // How man correlations
     277             :   //   Currently, this insists on uniformity over all SDBs
     278             :   int nCorrelations() const;
     279             : 
     280             :   // How many data chans?
     281             :   //   Currently, this insists on uniformity over all SDBs
     282             :   int nChannels() const;
     283             : 
     284             :   // The frequencies
     285             :   //   Currently, this insists on uniformity over all SDBs
     286             :   const casacore::Vector<double>& freqs() const;
     287             : 
     288             :   // ~Centroid frequency over all SDBs
     289             :   casacore::Double centroidFreq() const;
     290             : 
     291             :   // Simple centroid of per-SDB centroidFreqs 
     292             :   //  NB: this differs from centroidFreq in that it is _not_ a simple average of all SDB channel freqs
     293             :   //      this matters when different SDBs have different spws with different bandwidths/channelizations
     294             :   //      Eventually, this aggregation should be weighted by aggregate bandwidth, but this is not
     295             :   //      yet available from the VB2.
     296             :   casacore::Double aggregateCentroidFreq() const;
     297             : 
     298             :   // Does the SDBList contain usable data?
     299             :   //  (at least one SDB, with non-zero net weight)
     300             :   bool Ok();
     301             : 
     302             :   // Aggregated methods
     303             :   void enforceAPonData(const casacore::String& apmode);
     304             :   void enforceSolveWeights(const casacore::Bool pHandOnly);
     305             :   void sizeResiduals(const casacore::Int& nPar, const casacore::Int& nDiff);
     306             :   void initResidWithModel();
     307             :   void finalizeResiduals();
     308             :   //  NB: disable for now, may not be needed...
     309             :   //  void divideCorrByModel();
     310             : 
     311             : 
     312             :   // Manage working flags and weights
     313             :   void updateWorkingFlags();
     314             :   void updateWorkingWeights(casacore::Bool doL1=false,casacore::Float clamp=0.0);
     315             : 
     316             :   // Print out data, weights, flags for debugging purposes
     317             :   void reportData();
     318             :   
     319             :   // Extend baseline-dependent flags to all SDBs in the list
     320             :   //  This uniformizes the baseline-dependent flags
     321             :   //  NB: Cross-hands only, for now!
     322             :   casacore::Int extendCrossHandBaselineFlags(casacore::String& message);
     323             : 
     324             : private:
     325             : 
     326             :   // How many SDBs contained herein
     327             :   casacore::Int nSDB_;
     328             : 
     329             :   // Keep SDBs as a list of pointers
     330             :   casacore::PtrBlock<SolveDataBuffer*> SDB_;
     331             : 
     332             :   // Aggregate frequency storage (so we can calculate once-ish and return a reference)
     333             :   mutable casacore::Vector<double> freqs_;
     334             :   mutable double aggCentroidFreq_;
     335             :   mutable bool aggCentroidFreqOK_;
     336             : 
     337             : };
     338             : 
     339             : 
     340             : 
     341             : } //# NAMESPACE CASA - END
     342             : 
     343             : #endif
     344             :