00001 //# Copyright (C) 1996,1997,1998,1999,2000,2001,2003 00002 //# Associated Universities, Inc. Washington DC, USA. 00003 //# 00004 //# This library is free software; you can redistribute it and/or modify it 00005 //# under the terms of the GNU Library General Public License as published by 00006 //# the Free Software Foundation; either version 2 of the License, or (at your 00007 //# option) any later version. 00008 //# 00009 //# This library is distributed in the hope that it will be useful, but WITHOUT 00010 //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 00011 //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public 00012 //# License for more details. 00013 //# 00014 //# You should have received a copy of the GNU Library General Public License 00015 //# along with this library; if not, write to the Free Software Foundation, 00016 //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA. 00017 //# 00018 //# Correspondence concerning AIPS++ should be addressed as follows: 00019 //# Internet email: aips2-request@nrao.edu. 00020 //# Postal address: AIPS++ Project Office 00021 //# National Radio Astronomy Observatory 00022 //# 520 Edgemont Road 00023 //# Charlottesville, VA 22903-2475 USA 00024 //# 00025 //# $Id$ 00026 00027 #ifndef IMAGES_IMAGEBEAMSET_H 00028 #define IMAGES_IMAGEBEAMSET_H 00029 00030 #include <casacore/casa/aips.h> 00031 #include <casacore/casa/Arrays/Matrix.h> 00032 #include <casacore/scimath/Mathematics/GaussianBeam.h> 00033 //#include <casacore/measures/Measures/Stokes.h> 00034 //#include <map> 00035 00036 namespace casacore { 00037 00038 class SpectralCoordinate; 00039 00040 class CoordinateSystem; 00041 00042 // <summary> 00043 // Represents a set of restoring beams associated with an image. 00044 // </summary> 00045 00046 // <use visibility=export> 00047 00048 // <reviewed reviewer="" date="yyyy/mm/dd" tests="" demos=""> 00049 // </reviewed> 00050 00051 // <prerequisite> 00052 // </prerequisite> 00053 00054 // <etymology> 00055 // A Set of Beams associated with an Image. 00056 // </etymology> 00057 00058 // <synopsis> 00059 // This class represents a set of restoring beams associated with 00060 // a deconvolved image. Internally, the beams are stored in a Matrix in 00061 // which the first dimension represents the spectral axis and the second 00062 // dimension represents the polarization axis. Methods which take the number 00063 // of channels and stokes as the input parameters will accept 0 for either of 00064 // these, in the cases where the corresponding axis is absent, but internally, 00065 // the associated axis of the storage Matrix will be set to one since a Matrix 00066 // with one dimension of length 0 must be empty. If one (or both) of the axes is 00067 // of length 1, the beam associated with that position is valid for all channels or 00068 // stokes at the position. For example, if one has an image with 10 spectral channels 00069 // and 4 stokes, and one has a beam set of dimensions (1, 4) associated with that image, 00070 // all channels for a given stokes will have the same beam. Similarly, if the beam set 00071 // is of shape (10, 1) in this case, all stokes will have the same beam for a given channel. 00072 // If the axis lengths of the beam set are greater than one, they must be exactly 00073 // the same length of the corresponding axes in the associated image. 00074 // </synopsis> 00075 // 00076 // <example> 00077 00078 // </example> 00079 00080 00081 // <motivation> 00082 // Restoring beams are used many places in image analysis tasks. 00083 // </motivation> 00084 00085 // <todo> 00086 // </todo> 00087 00088 class ImageBeamSet { 00089 public: 00090 00091 typedef Array<GaussianBeam>::const_iterator BeamIter; 00092 00093 // Construct an empty beam set. 00094 ImageBeamSet(); 00095 00096 // Construct a beam set from an 2-D array of beams representing 00097 // the frequency and stokes axis. 00098 // Axis length 1 means it is valid for all channels cq. stokes. 00099 // If the image has 0 spectral channels or stokes, the corresponding 00100 // length of the axis in the provided matrix should be 1. 00101 ImageBeamSet( 00102 const Matrix<GaussianBeam>& beams 00103 ); 00104 00105 // construct an ImageBeamSet representing a single beam which is valid for 00106 // all channels and stokes 00107 ImageBeamSet(const GaussianBeam& beam); 00108 00109 // Create an ImageBeamSet of the specified shape with all 00110 // GaussianBeams initialized to <src>beam</src>. 00111 ImageBeamSet(uInt nchan, uInt nstokes, const GaussianBeam& beam=GaussianBeam::NULL_BEAM); 00112 00113 // The copy constructor (reference semantics). 00114 ImageBeamSet(const ImageBeamSet& other); 00115 00116 ~ImageBeamSet(); 00117 00118 // Assignment can change the shape (copy semantics). 00119 ImageBeamSet& operator=(const ImageBeamSet& other); 00120 00121 // Beam sets are equal if the shapes and all corresponding beams are equal. 00122 Bool operator== (const ImageBeamSet& other) const; 00123 Bool operator!= (const ImageBeamSet& other) const; 00124 00125 // Beam sets are equivalent if both have no beams or if the 00126 // expanded sets are equal. Expanded means that an axis can have 00127 // length 0 or 1 and is (virtually) expanded to the length of the matching 00128 // axis in the other beam set. 00129 Bool equivalent (const ImageBeamSet& that) const; 00130 00131 // Get the number of elements in the beam array. 00132 // <group> 00133 uInt nelements() const 00134 { return _beams.size(); } 00135 uInt size() const 00136 { return _beams.size(); } 00137 // </group> 00138 00139 Bool hasSingleBeam() const 00140 { return _beams.size() == 1; } 00141 00142 // Does this beam set contain multiple beams? 00143 Bool hasMultiBeam() const { 00144 return _beams.size() > 1; 00145 } 00146 00147 // Is the beam set empty? 00148 Bool empty() const 00149 { return _beams.empty(); } 00150 00151 // Get the shape of the beam array. The minimum value for 00152 // a component of the returned IPosition is always 1. 00153 const IPosition& shape() const 00154 { return _beams.shape(); } 00155 00156 // Get the number of channels in the beam array. NOte that this will 00157 // always return a minimum of 1, even if nchan was specified as 0 on construction. 00158 uInt nchan() const 00159 { return _beams.shape()[0]; } 00160 00161 // Get the number of stokes in the beam array. Note that this will always 00162 // return a minimum of 1, even if nstokes was specified as 0 on construction. 00163 uInt nstokes() const 00164 { return _beams.shape()[1]; } 00165 00166 // Get the single global beam. If there are multiple beams, 00167 // an exception is thrown. 00168 const GaussianBeam& getBeam() const; 00169 00170 // Get the beam at the specified location. 00171 // Note that a single channel or stokes in the beam set is valid for 00172 // all channels cq. stokes. 00173 // <group> 00174 const GaussianBeam& getBeam(Int chan, Int stokes) const; 00175 const GaussianBeam &operator()(Int chan, Int stokes) const 00176 { return getBeam (chan, stokes); } 00177 // </group> 00178 00179 // Get a beam at the given 2-dim IPosition. It should match exactly, 00180 // thus a single channel or stokes in the beam set is not valid for all. 00181 // const GaussianBeam& getBeam(const IPosition& pos) const 00182 // { return _beams(pos); } 00183 00184 // Set the beam at the given location. 00185 // The location must be within the beam set shape. 00186 // If <src>chan</src> or <src>stokes</src> is negative, then the beam applies 00187 // to all channels or stokes, respectively. If both are negative, the specified 00188 // beam becomes the global beam and the beam set is resized to (1, 1). 00189 void setBeam(Int chan, Int stokes, const GaussianBeam& beam); 00190 00191 // Resize the beam array. <src>nchan</src>=0 or <src>nstokes</src>=0 00192 // is silently changed to 1. 00193 void resize(uInt nchan, uInt nstokes); 00194 00195 // Return a subset of the beam array. 00196 // The slicer is usually the slicer used for a subimage. 00197 // The slicer can contain multiple stokes or channels, even if the 00198 // beam set has only one. 00199 ImageBeamSet subset (const Slicer& imageSlicer, 00200 const CoordinateSystem& csys) const; 00201 00202 // Get the beam array. 00203 const Matrix<GaussianBeam>& getBeams() const 00204 { return _beams; } 00205 00206 // Set the beams in this beam set. 00207 // The shape of the given array must match the beam set. 00208 // It also matches if an axis in array or beam set has length 1, which 00209 // means that it expands to the other length. 00210 void setBeams(const Matrix<GaussianBeam>& beams); 00211 00212 // Set all beams to the same value. 00213 void set(const GaussianBeam& beam); 00214 00215 // Get the beam in the set which has the smallest area. 00216 GaussianBeam getMinAreaBeam() const 00217 { return _minBeam; } 00218 00219 // Get the beam in the set which has the largest area. 00220 // Get the beam in the set which has the largest area. 00221 GaussianBeam getMaxAreaBeam() const 00222 { return _maxBeam; } 00223 00224 // Get the beam in the set which has the median area. 00225 GaussianBeam getMedianAreaBeam() const; 00226 00227 // Get the position of the beam with the minimum area. 00228 IPosition getMinAreaBeamPosition() const 00229 { return _minBeamPos; } 00230 00231 // Get the position of the beam with the maximum area. 00232 IPosition getMaxAreaBeamPosition() const 00233 { return _maxBeamPos; } 00234 00235 // Get the minimal, maximal, and median area beams and positions in the beam set matrix for 00236 // the given stokes. If the stokes axis has length 1 in the beam matrix, 00237 // it is valid for all stokes and no checking is done that <src>stokes</src> is valid; 00238 // the requested beam for the entire beam set is simply returned in this case. If the 00239 // number of stokes in the beam matrix is >1, checking is done that the specified value 00240 // of <src>stokes</src> is valid and if not, an exception is thrown. 00241 // <group> 00242 const GaussianBeam& getMinAreaBeamForPol(IPosition& pos, 00243 uInt stokes) const; 00244 00245 const GaussianBeam& getMaxAreaBeamForPol(IPosition& pos, 00246 uInt stokes) const; 00247 00248 const GaussianBeam& getMedianAreaBeamForPol(IPosition& pos, 00249 uInt stokes) const; 00250 // </group> 00251 00252 static const String& className(); 00253 00254 // Get the beam that has the smallest minor axis. If multiple beams have the smallest minor axis, 00255 // the beam in this subset with the smallest area will be returned. 00256 const GaussianBeam getSmallestMinorAxisBeam() const; 00257 00258 // convert ImageBeamSet to and from record 00259 // <group> 00260 static ImageBeamSet fromRecord(const Record& rec); 00261 Record toRecord() const; 00262 //</group> 00263 00264 // If verbose, log all beams, if not just summarize beam stats. 00265 void summarize(LogIO& log, Bool verbose, const CoordinateSystem& csys) const; 00266 00267 // Modify the beam set by rotating all beams counterclockwise through the specified angle. 00268 void rotate(const Quantity& angle); 00269 00270 private: 00271 00272 static const String _DEFAULT_AREA_UNIT; 00273 00274 Matrix<GaussianBeam> _beams; 00275 Matrix<Double> _areas; 00276 String _areaUnit; 00277 GaussianBeam _minBeam, _maxBeam; 00278 IPosition _minBeamPos, _maxBeamPos; 00279 00280 void _calculateAreas(); 00281 00282 static void _chanInfoToStream( 00283 ostream& os, const SpectralCoordinate *spCoord, 00284 const uInt chan, const uInt chanWidth, const uInt freqPrec, 00285 const uInt velWidth, const uInt velPrec 00286 ); 00287 00288 static void _beamToStream( 00289 ostream& os, const GaussianBeam& beam, 00290 const Unit& unit 00291 ); 00292 }; 00293 00294 ostream &operator<<(ostream &os, const ImageBeamSet& beamSet); 00295 00296 } 00297 00298 #endif 00299