Line data Source code
1 : //# Copyright (C) 1996,1997,1998,1999,2000,2001,2002
2 : //# Associated Universities, Inc. Washington DC, USA.
3 : //#
4 : //# This library is free software; you can redistribute it and/or modify it
5 : //# under the terms of the GNU Library General Public License as published by
6 : //# the Free Software Foundation; either version 2 of the License, or (at your
7 : //# option) any later version.
8 : //#
9 : //# This library is distributed in the hope that it will be useful, but WITHOUT
10 : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11 : //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
12 : //# License for more details.
13 : //#
14 : //# You should have received a copy of the GNU Library General Public License
15 : //# along with this library; if not, write to the Free Software Foundation,
16 : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
17 : //#
18 : //# Correspondence concerning AIPS++ should be addressed as follows:
19 : //# Internet email: aips2-request@nrao.edu.
20 : //# Postal address: AIPS++ Project Office
21 : //# National Radio Astronomy Observatory
22 : //# 520 Edgemont Road
23 : //# Charlottesville, VA 22903-2475 USA
24 : //#
25 :
26 : #ifndef IMAGEANALYSIS_IMAGEPOLTASK_H
27 : #define IMAGEANALYSIS_IMAGEPOLTASK_H
28 :
29 : #include <imageanalysis/ImageAnalysis/ImageTask.h>
30 :
31 : namespace casa {
32 :
33 : // <summary>
34 : // Begin moving polarization tasks to ImageTask framework.
35 : // </summary>
36 :
37 : // <use visibility=export>
38 :
39 : // <prerequisite>
40 : // <li> <linkto class=casacore::ImageExpr>ImageExpr</linkto>
41 : // <li> <linkto class=casacore::ImageInterface>ImageInterface</linkto>
42 : // </prerequisite>
43 :
44 : // <etymology>
45 : // Polarimetric analysis of Images
46 : // </etymology>
47 :
48 : // <synopsis>
49 : // This class provides polarimetric image analysis capability.
50 : // It takes an image with a casacore::Stokes axis (some combination
51 : // of IQUV is needed) as its input.
52 : //
53 : // Many functions return casacore::ImageExpr objects. These are
54 : // read-only images.
55 : //
56 : // Sometimes the standard deviation of the noise is needed.
57 : // This is for debiasing polarized intensity images or working out
58 : // error images. By default it is worked out for you with a
59 : // clipped mean algorithm. However, you can provide sigma if you
60 : // know it accurately. It should be the standard deviation of the noise in
61 : // the absence of signal. You won't measure that very well from
62 : // casacore::Stokes I if it is dynamic range limited. Better to get it from
63 : // V or Q or U. When this class needs the standard deviation of
64 : // the noise, it will try and get it from V or Q and U and finally I.
65 : //
66 : // However, note that the functions sigmaStokes{I,Q,U,V} DO return the standard
67 : // deviation of the noise for that specific casacore::Stokes type.
68 : //
69 : // The casacore::ImageExpr objects returned have the brightness units and ImageInfo
70 : // set. The MiscInfo (a permanent record) and logSink are not set.
71 : //
72 : // </synopsis>
73 : //
74 : // <motivation>
75 : // Basic image analysis capability
76 : // </motivation>
77 :
78 : // <todo asof="1999/11/01">
79 : // <li> plotting for function rotationMeasure
80 : // <li> some assessment of the curvature of pa-l**2
81 : // </todo>
82 :
83 : class ImagePolTask : public ImageTask<casacore::Float> {
84 : public:
85 :
86 : enum StokesTypes {I, Q, U, V};
87 :
88 : ImagePolTask() = delete;
89 :
90 : virtual ~ImagePolTask();
91 :
92 : virtual casacore::String getClass() const;
93 :
94 : /*
95 : // Get the linearly polarized intensity image and its
96 : // standard deviation. If wish to debias the image, you
97 : // can either provide <src>sigma</src> (the standard
98 : // deviation of the termal noise ) or if <src>sigma</src> is non-positive,
99 : // it will be worked out for you with a clipped mean algorithm.
100 : // <group>
101 : casacore::ImageExpr<casacore::Float> linPolInt(
102 : casacore::Bool debias, casacore::Float clip=10.0, casacore::Float sigma=-1.0
103 : );
104 : */
105 :
106 : casacore::Float sigmaLinPolInt(
107 : casacore::Float clip=10.0, casacore::Float sigma=-1.0
108 : );
109 : // </group>
110 :
111 : protected:
112 :
113 : ImagePolTask(
114 : const SPCIIF image, const casacore::String& outname,
115 : casacore::Bool overwrite
116 : );
117 :
118 : casacore::Bool _checkQUBeams(
119 : casacore::Bool requireChannelEquality, casacore::Bool throws=casacore::True
120 : ) const;
121 :
122 : // Change the casacore::Stokes casacore::Coordinate for the given
123 : // complex image to be of the specified casacore::Stokes type
124 : void _fiddleStokesCoordinate(
125 : casacore::ImageInterface<casacore::Float>& ie,
126 : casacore::Stokes::StokesTypes type
127 : ) const;
128 :
129 0 : std::vector<casacore::Coordinate::Type> _getNecessaryCoordinates() const {
130 0 : static const std::vector<casacore::Coordinate::Type> v { casacore::Coordinate::STOKES };
131 0 : return v;
132 : }
133 :
134 0 : CasacRegionManager::StokesControl _getStokesControl() const {
135 0 : return CasacRegionManager::USE_ALL_STOKES;
136 : }
137 :
138 : SPCIIF _getStokesImage(StokesTypes type) const;
139 :
140 : // Make a LEN for the give types of polarized intensity
141 : casacore::LatticeExprNode _makePolIntNode(
142 : casacore::Bool debias, casacore::Float clip, casacore::Float sigma,
143 : casacore::Bool doLin, casacore::Bool doCirc
144 : );
145 :
146 : void _maskAndZeroNaNs(SPIIF out);
147 :
148 : void _setDoLinDoCirc(
149 : casacore::Bool& doLin, casacore::Bool& doCirc,
150 : casacore::Bool requireI
151 : ) const;
152 :
153 : void _setInfo(
154 : casacore::ImageInterface<Float>& im, const StokesTypes stokes
155 : ) const;
156 :
157 : private:
158 :
159 : // These blocks are always size 4, with IQUV in slots 0,1,2,3
160 : // If your image is IV only, they still use slots 0 and 3
161 : std::vector<SPIIF> _stokesImage = std::vector<SPIIF>(4);
162 : std::vector<std::shared_ptr<casacore::LatticeStatistics<casacore::Float>>> _stokesStats
163 : = std::vector<std::shared_ptr<casacore::LatticeStatistics<casacore::Float>>>(4);
164 : casacore::Matrix<casacore::Bool> _beamsEqMat
165 : = casacore::Matrix<casacore::Bool>(4, 4, casacore::False);
166 : Float _oldClip = 0.0;
167 :
168 : static const casacore::String CLASS_NAME;
169 :
170 : casacore::Bool _checkBeams(
171 : const std::vector<StokesTypes>& stokes,
172 : casacore::Bool requireChannelEquality,
173 : casacore::Bool throws=true
174 : ) const;
175 :
176 : casacore::Bool _checkIQUBeams(
177 : casacore::Bool requireChannelEquality, casacore::Bool throws=casacore::True
178 : ) const;
179 :
180 : casacore::Bool _checkIVBeams(
181 : casacore::Bool requireChannelEquality, casacore::Bool throws=casacore::True
182 : ) const;
183 :
184 : void _createBeamsEqMat();
185 :
186 : // Find the casacore::Stokes in the construction image and assign pointers
187 : void _findStokes();
188 :
189 : void _fiddleStokesCoordinate(
190 : casacore::CoordinateSystem& cSys, casacore::Stokes::StokesTypes type
191 : ) const;
192 :
193 :
194 : // Make a casacore::SubImage from the construction image for the specified pixel
195 : // along the specified pixel axis
196 : SPIIF _makeSubImage(
197 : casacore::IPosition& blc, casacore::IPosition& trc,
198 : casacore::Int axis, casacore::Int pix
199 : ) const;
200 :
201 : // Get the best estimate of the statistical noise. This gives you
202 : // the standard deviation with outliers from the mean
203 : // clipped first. The idea is to not be confused by source or dynamic range issues.
204 : // Generally casacore::Stokes V is empty of sources (not always), then Q and U are generally
205 : // less bright than I. So this function first tries V, then Q and U
206 : // and lastly I to give you its noise estimate
207 : casacore::Float _sigma (casacore::Float clip=10.0);
208 :
209 : // Find the standard deviation for the casacore::Stokes image specified by the integer index
210 : casacore::Float _sigma(StokesTypes index, casacore::Float clip);
211 :
212 : // Return I, Q, U or V for specified integer index (0-3)
213 : casacore::String _stokesName (StokesTypes index) const;
214 :
215 : };
216 :
217 : }
218 :
219 : #endif
|