Line data    Source code

       1             : //# AspMatrixCleaner.h: Minor Cycle for Asp deconvolution
       2             : //# Copyright (C) 1996,1997,1998,1999,2000,2001,2002,2003
       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: aips2-request@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             : //#
      27             : //# $Id: AspMatrixCleaner.h Genie H. 2020-04-06 <mhsieh@nrao.edu $
      28             : 
      29             : #ifndef SYNTHESIS_ASPMATRIXCLEANER_H
      30             : #define SYNTHESIS_ASPMATRIXCLEANER_H
      31             : 
      32             : //# Includes
      33             : #include <casacore/scimath/Mathematics/FFTServer.h>
      34             : #include <synthesis/MeasurementEquations/MatrixCleaner.h>
      35             : #include <deque>
      36             : #include <casacore/coordinates/Coordinates/CoordinateSystem.h>
      37             : #include <casacore/coordinates/Coordinates/SpectralCoordinate.h>
      38             : 
      39             : namespace casa { //# NAMESPACE CASA - BEGIN
      40             : 
      41             : class AspMatrixCleaner : public MatrixCleaner
      42             : {
      43             : public:
      44             :   // Create a cleaner : default constructor
      45             :   AspMatrixCleaner();
      46             : 
      47             :   // The destructor does nothing special.
      48             :   ~AspMatrixCleaner();
      49             : 
      50             :   casacore::Bool setaspcontrol(const casacore::Int niter,
      51             :       const casacore::Float gain, const casacore::Quantity& aThreshold,
      52             :       const casacore::Quantity& fThreshold);
      53             : 
      54             : 
      55             :   // Clean an image.
      56             :   //return value gives you a hint of what's happening
      57             :   //  1 = converged
      58             :   //  0 = not converged but behaving normally
      59             :   // -1 = not converged and stopped on cleaning consecutive smallest scale
      60             :   // -2 = not converged and either large scale hit negative or diverging
      61             :   // -3 = clean is diverging rather than converging
      62             :   casacore::Int aspclean(casacore::Matrix<casacore::Float> & model, casacore::Bool doPlotProgress=false);
      63             : 
      64             :   // helper functions for ASP
      65             :   float getPsfGaussianWidth(casacore::ImageInterface<casacore::Float>& psf);
      66             :   void getLargestScaleSize(casacore::ImageInterface<casacore::Float>& psf);
      67             : 
      68             :   // Make an image of the specified scale by Gaussian
      69             :   void makeInitScaleImage(casacore::Matrix<casacore::Float>& iscale, const casacore::Float& scaleSize);
      70             :   void makeScaleImage(casacore::Matrix<casacore::Float>& iscale, const casacore::Float& scaleSize, const casacore::Float& amp, const casacore::IPosition& center);
      71             : 
      72             :   void setInitScales();
      73             :   void setInitScaleXfrs(const casacore::Float width);
      74             : 
      75             :   // calculate the convolutions of the psf with the initial scales
      76             :   void setInitScalePsfs();
      77             : 
      78             :   casacore::Bool setInitScaleMasks(const casacore::Array<casacore::Float> arrmask, const casacore::Float& maskThreshold = 0.99);
      79             : 
      80             :   void maxDirtyConvInitScales(float& strengthOptimum, int& optimumScale, casacore::IPosition& positionOptimum);
      81             : 
      82             :   // returns the active-set aspen for cleaning
      83             :   std::vector<casacore::Float> getActiveSetAspen();
      84             : 
      85             :   // Juat define the active-set aspen scales
      86             :   void defineAspScales(std::vector<casacore::Float>& scaleSizes);
      87             : 
      88             :   void switchedToHogbom(bool runlong= false);
      89             :   void setOrigDirty(const casacore::Matrix<casacore::Float>& dirty);
      90          40 :   void setFusedThreshold(const casacore::Float fusedThreshold = 0.0) { itsFusedThreshold = fusedThreshold; }
      91          40 :   void setUserLargestScale(const casacore::Int largestScale = -1) { itsUserLargestScale = float(largestScale); }
      92             : 
      93             :   // setter/getter
      94             :   float getterPsfWidth() { return itsPsfWidth; }
      95             :   bool getterSwitchedHogbom() { return itsSwitchedToHogbom; }
      96          40 :   casacore::Matrix<casacore::Float>  getterResidual() { return (*itsDirty); }
      97          40 :   float getterPeakResidual() { return itsPeakResidual; }
      98             : 
      99             :   void setBinSizeForSumFlux(const casacore::Int binSize = 4) { itsBinSizeForSumFlux = binSize ; } ;
     100             :   void getFluxByBins(const std::vector<casacore::Float>& scaleSizes,const std::vector<casacore::Float>& optimum, casacore::Int binSize, std::vector<casacore::Float>&  sumFluxByBins, std::vector<casacore::Float>&  rangeFluxByBins);
     101             : 
     102             : 
     103             : //protected:
     104             : private:
     105             : 
     106             :   using MatrixCleaner::findMaxAbs;
     107             :   using MatrixCleaner::findMaxAbsMask;
     108             :   using MatrixCleaner::itsGain;
     109             :   using MatrixCleaner::itsMaxNiter;
     110             :   using MatrixCleaner::itsThreshold;
     111             :   using MatrixCleaner::itsMask;
     112             :   using MatrixCleaner::itsScalesValid;
     113             :   using MatrixCleaner::itsNscales;
     114             :   using MatrixCleaner::itsMaskThreshold;
     115             :   using MatrixCleaner::itsDirty;
     116             :   using MatrixCleaner::itsXfr;
     117             :   using MatrixCleaner::itsScaleSizes;
     118             :   casacore::Block<casacore::Matrix<casacore::Float> > itsInitScales;
     119             :   casacore::Block<casacore::Matrix<casacore::Complex> > itsInitScaleXfrs;
     120             :   casacore::Block<casacore::Matrix<casacore::Float> > itsDirtyConvInitScales;
     121             :   casacore::Block<casacore::Matrix<casacore::Float> > itsInitScaleMasks;
     122             :   casacore::Block<casacore::Matrix<casacore::Float> > itsPsfConvInitScales;
     123             : 
     124             :   using MatrixCleaner::itsIteration;
     125             :   using MatrixCleaner::itsStartingIter;
     126             :   using MatrixCleaner::itsFracThreshold;
     127             :   using MatrixCleaner::itsMaximumResidual;
     128             :   using MatrixCleaner::itsStrengthOptimum;
     129             :   using MatrixCleaner::itsTotalFlux;
     130             :   using MatrixCleaner::index;
     131             : 
     132             :   using MatrixCleaner::destroyScales;
     133             :   casacore::Bool destroyAspScales();
     134             :   casacore::Bool destroyInitScales();
     135             :   using MatrixCleaner::destroyMasks;
     136             :   casacore::Bool destroyInitMasks();
     137             :   casacore::Float computeThreshold() const;
     138             : 
     139             :   using MatrixCleaner::itsIgnoreCenterBox;
     140             :   using MatrixCleaner::itsStopAtLargeScaleNegative;
     141             :   using MatrixCleaner::itsStopPointMode;
     142             :   using MatrixCleaner::itsDidStopPointMode;
     143             :   using MatrixCleaner::psfShape_p;
     144             :   using MatrixCleaner::noClean_p;
     145             : 
     146             :   // FFTServer
     147             :   casacore::FFTServer<casacore::Float,casacore::Complex> fft;
     148             : 
     149             :   // set to 0, 1*, 2*, 4*, and 8* width for initial scales in Asp
     150             :   std::vector<casacore::Float> itsInitScaleSizes;
     151             :   //std::vector<casacore::Float> itsAspScaleSizes; // No longer needed. permanent list for making model image
     152             :   //std::vector<casacore::Float> itsAspAmplitude;
     153             :   //std::vector<casacore::IPosition> itsAspCenter;
     154             :   casacore::Int itsNInitScales;
     155             :   std::deque<int> itsNumIterNoGoodAspen;
     156             :   float itsPsfWidth;
     157             :   bool itsUseZhang;
     158             :   bool itsSwitchedToHogbom;
     159             :   unsigned int itsNumHogbomIter;
     160             :   unsigned int itsNthHogbom;
     161             :   std::vector<casacore::Float> itsGoodAspActiveSet; // avtice-set of aspens (updated)
     162             :   std::vector<casacore::Float> itsGoodAspAmplitude; // amplitude of avtice-set of aspens (updated)
     163             :   std::vector<casacore::IPosition> itsGoodAspCenter; // center of aspens in active-set
     164             :   //std::vector<casacore::Float> itsPrevAspActiveSet; // No longer needed. avtice-set of aspens (before bfgs)
     165             :   //std::vector<casacore::Float> itsPrevAspAmplitude; // No longer needed. amplitude of avtice-set of aspens (before bfgs)
     166             :   casacore::Int itsOptimumScale;
     167             :   casacore::IPosition itsPositionOptimum;
     168             :   casacore::Float itsOptimumScaleSize;
     169             :   double itsUsedMemoryMB;
     170             :   float itsPeakResidual;
     171             :   float itsPrevPeakResidual;
     172             :   casacore::CountedPtr<casacore::Matrix<casacore::Float> > itsOrigDirty;
     173             : 
     174             :   const casacore::Int itsDefaultNorm = 1;
     175             :   casacore::Int itsNormMethod;
     176             :   casacore::Float itsFusedThreshold;
     177             :   casacore::Int itsBinSizeForSumFlux ;   // number of bins for histogram of the sum of Flux
     178             :   unsigned int itsNumNoChange; // number of times peakres rarely changes
     179             :   float itsLargestInitScale; // estimated largest initial scale
     180             :   float itsUserLargestScale; // user-specified largest initial scale
     181             : };
     182             : 
     183             : } //# NAMESPACE CASA - END
     184             : 
     185             : #endif