Line data Source code
1 : //# SkyCal.h: this defines <ClassName>, which ...
2 : //# Copyright (C) 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: 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 : //#
27 : //# $Id$
28 :
29 : #ifndef _SYNTHESIS_SKY_CAL_H_
30 : #define _SYNTHESIS_SKY_CAL_H_
31 :
32 : #include <casacore/casa/Arrays/Matrix.h>
33 : #include <casacore/casa/Exceptions/Error.h>
34 :
35 : namespace casa { //# NAMESPACE CASA - BEGIN
36 :
37 : // <summary>
38 : //#! A one line summary of the class. This summary (shortened a bit
39 : //#! if necessary so that it fits along with the "ClassFileName.h" in 75
40 : //#! characters) should also appear as the first line of this file.
41 : //#! Be sure to use the word "abstract" here if this class is, indeed,
42 : //#! an abstract base class.
43 : // </summary>
44 :
45 : // <use visibility=local> or <use visibility=export>
46 : //#! If a class is intended for use by application programmers, or
47 : //#! people writing other libraries, then specify that this class
48 : //#! has an "export" visibility: it defines an interface that
49 : //#! will be seen outside of its module. On the other hand, if the
50 : //#! class has a "local" visibility, then it is known and used only
51 : //#! within its module.
52 :
53 : // <reviewed reviewer="" date="yyyy/mm/dd" tests="" demos="">
54 : //#! for example:
55 : //#! <reviewed reviewer="pshannon@nrao.edu" date="1994/10/10" tests="tMyClass, t1MyClass" demos="dMyClass, d1MyClass">
56 : //#! This is a well-designed class, without any obvious problems.
57 : //#! However, the addition of several more demo programs would
58 : //#! go a *long* way towards making it more usable.
59 : //#! </reviewed>
60 : //#!
61 : //#! (In time, the documentation extractor will be able handle reviewed
62 : //#! attributes spread over more than one line.)
63 : //#!
64 : //#! see "Coding Standards and Guidelines" (AIPS++ note 167) and
65 : //#! "AIPS++ Code Review Process" (note 170) for a full explanation
66 : //#! It is up to the class author (the programmer) to fill in these fields:
67 : //#! tests, demos
68 : //#! The reviewer fills in
69 : //#! reviewer, date
70 : // </reviewed>
71 :
72 : // <prerequisite>
73 : //#! Classes or concepts you should understand before using this class.
74 : // <li> SomeClass
75 : // <li> SomeOtherClass
76 : // <li> some concept
77 : // </prerequisite>
78 : //
79 : // <etymology>
80 : //#! Except when it is obvious (e.g., "casacore::Array") explain how the class name
81 : //#! expresses the role of this class. Example: casacore::IPosition is short for
82 : //#! "Integral Position" - a specialized integer vector for specifying
83 : //#! array dimensions and indices.
84 : // </etymology>
85 : //
86 : // <synopsis>
87 : //#! What does the class do? How? For whom? This should include code
88 : //#! fragments as appropriate to support text. Code fragments shall be
89 : //#! delimited by <srcblock> </srcblock> tags. The synopsis section will
90 : //#! usually be dozens of lines long.
91 : // </synopsis>
92 : //
93 : // <example>
94 : //#! One or two concise (~10-20 lines) examples, with a modest amount of
95 : //#! text to support code fragments. Use <srcblock> and </srcblock> to
96 : //#! delimit example code.
97 : // </example>
98 : //
99 : // <motivation>
100 : //#! Insight into a class is often provided by a description of
101 : //#! the circumstances that led to its conception and design. Describe
102 : //#! them here.
103 : // </motivation>
104 : //
105 : // <templating arg=T>
106 : //#! A list of member functions that must be provided by classes that
107 : //#! appear as actual template arguments. For example: imagine that you
108 : //#! are writing a templated sort class, which does a quicksort on a
109 : //#! list of arbitrary objects. Anybody who uses your templated class
110 : //#! must make sure that the actual argument class (say, casacore::Int or
111 : //#! casacore::String or casacore::Matrix) has comparison operators defined.
112 : //#! This tag must be repeated for each template formal argument in the
113 : //#! template class definition -- that's why this tag has the "arg" attribute.
114 : //#! (Most templated classes, however, are templated on only a single
115 : //#! argument.)
116 : // <li>
117 : // <li>
118 : // </templating>
119 : //
120 : // <thrown>
121 : //#! A list of exceptions thrown if errors are discovered in the function.
122 : //#! This tag will appear in the body of the header file, preceding the
123 : //#! declaration of each function which throws an exception.
124 : // <li>
125 : // <li>
126 : // </thrown>
127 : //
128 : // <todo asof="yyyy/mm/dd">
129 : //#! A casacore::List of bugs, limitations, extensions or planned refinements.
130 : //#! The programmer should fill in a date in the "asof" field, which
131 : //#! will usually be the date at which the class is submitted for review.
132 : //#! If, during the review, new "todo" items come up, then the "asof"
133 : //#! date should be changed to the end of the review period.
134 : // <li> add this feature
135 : // <li> fix this bug
136 : // <li> start discussion of this possible extension
137 : // </todo>
138 :
139 : template<class DataType, class CalDataType>
140 : class SkyCal
141 : {
142 : // all the member functions are inline
143 : public:
144 759 : SkyCal()
145 759 : : npol_(2),
146 759 : nchan_(1),
147 759 : m0_(NULL),
148 759 : ok0_(NULL),
149 759 : m_(NULL),
150 759 : mi_(NULL),
151 759 : ok_(NULL),
152 759 : oki_(NULL),
153 759 : cOne_(1.0),
154 759 : cZero_(0.0),
155 759 : scalardata_(false)
156 759 : {}
157 :
158 1518 : virtual ~SkyCal() {}
159 :
160 : // stride
161 : // In Mueller series, typesize is defined as "number of polarizations"
162 : // while SkyCal definition is "number of elements in the DATA cell",
163 : // which is npol * nchan in practice.
164 5076 : casacore::uInt typesize() const { return npol_ * nchan_; }
165 5076 : void setNumChannel(casacore::uInt n) { nchan_ = n; }
166 5076 : void setNumPolarization(casacore::uInt n) { npol_ = n; }
167 :
168 : // Set scalardata_
169 : // TBD: Handle this better; for now, we need to set this from
170 : // an external call so we handle single-corr data properly
171 : // when setting non-corr-dep flags
172 0 : void setScalarData(casacore::Bool scalardata) const { scalardata_=scalardata; }
173 :
174 : // Synchronize with leading element in external array
175 0 : void sync(CalDataType& mat) { m0_=&mat; origin(); }
176 5490 : void sync(CalDataType& mat, casacore::Bool& ok) { m0_=&mat; ok0_=&ok; origin(); }
177 :
178 : // Reset to origin
179 5490 : void origin() {m_=m0_;ok_=ok0_;}
180 :
181 : // Increment to next vector (according to len)
182 : // In practice, this operator increments row index
183 0 : void operator++() { m_+=typesize(); if (ok_) ok_+=typesize();}
184 0 : void operator++(int) { m_+=typesize(); if (ok_) ok_+=typesize();}
185 :
186 : // Advance step matrices forward (according to len)
187 0 : void advance(const casacore::Int& step) { m_+=(step*typesize()); if (ok_) ok_+=(step*typesize());}
188 :
189 : // In-place invert
190 0 : void invert() {}
191 :
192 : // Set matrix elements according to ok flag
193 : // (so we don't have to check ok flags atomically in apply)
194 0 : void setMatByOk()
195 : {
196 0 : throw(casacore::AipsError("Illegal use of SkyCal::setMatByOk"));
197 : }
198 :
199 : // In-place multiply onto a data with flag information
200 : // apply implements position switch calibration: (ON - OFF)/OFF
201 : //
202 : // This processes the data corresponding to each DATA cell
203 : // (npol * nchan) together in contrast to Mueller series, which
204 : // processes one casacore::Stokes vector, i.e., process each channel
205 : // individually.
206 5490 : void apply(casacore::Matrix<DataType>& v, casacore::Matrix<casacore::Bool> &f)
207 : {
208 5490 : if (f.shape() == v.shape()) {
209 5490 : flag(f);
210 : }
211 :
212 5490 : constexpr size_t nPolCal = 2;
213 :
214 : casacore::Bool deleteIt;
215 5490 : DataType *data = v.getStorage(deleteIt);
216 3393382 : for (size_t i = 0; i < nchan_; ++i) {
217 3387892 : auto const iData = i * npol_;
218 3387892 : auto const iMat = i * nPolCal;
219 10158676 : for (size_t j = 0; j < npol_; ++j) {
220 : // (ON - OFF) / OFF
221 6770784 : data[iData+j] = (data[iData+j] - m_[iMat+j]) / m_[iMat+j];
222 : }
223 : }
224 5490 : v.putStorage(data, deleteIt);
225 5490 : }
226 :
227 0 : void apply(casacore::Matrix<DataType>& v, casacore::Matrix<casacore::Bool> &f, casacore::Vector<casacore::Bool>& vflag)
228 : {
229 0 : if (!ok_) throw(casacore::AipsError("Illegal use of SkyCal::applyFlag."));
230 :
231 0 : applyFlag(vflag);
232 0 : apply(v, f);
233 0 : }
234 :
235 : // Apply only flags according to cal flags
236 : //
237 : // Similar to apply, flagging also processes each DATA cell together.
238 0 : void applyFlag(casacore::Vector<casacore::Bool>& vflag)
239 : {
240 0 : if (!ok_) throw(casacore::AipsError("Illegal use of SkyCal::applyFlag(vflag)."));
241 :
242 0 : constexpr size_t nPolCal = 2;
243 :
244 0 : if (scalardata_) {
245 0 : for (size_t i = 0; i < nchan_; ++i) {
246 0 : vflag[i] |= (!ok_[i*nPolCal]);
247 : }
248 : }
249 : else {
250 0 : for (size_t i = 0; i < nchan_; ++i) {
251 0 : for (size_t j = 0; j < npol_; ++j) {
252 0 : vflag[i] |= !(ok_[i*nPolCal + j]);
253 : }
254 : }
255 : }
256 0 : }
257 :
258 5490 : void flag(casacore::Matrix<casacore::Bool>& v)
259 : {
260 5490 : constexpr size_t nPolCal = 2;
261 :
262 : casacore::Bool deleteIt;
263 5490 : casacore::Bool *data = v.getStorage(deleteIt);
264 3393382 : for (size_t i = 0; i < nchan_; ++i) {
265 3387892 : auto const iData = i * npol_;
266 3387892 : auto const iMat = i * nPolCal;
267 10158676 : for (size_t j = 0; j < npol_; ++j) {
268 6770784 : data[iData+j] |= (!ok_[iMat+j]); // data: false is valid, ok_: true is valid
269 : }
270 : }
271 5490 : v.putStorage(data, deleteIt);
272 5490 : }
273 :
274 : // Multiply onto a vis VisVector, preserving input (copy then in-place apply)
275 0 : void apply(casacore::Matrix<DataType>& out, casacore::Matrix<casacore::Bool> &outFlag,
276 : const casacore::Matrix<DataType>& in, const casacore::Matrix<casacore::Bool> &inFlag)
277 : {
278 0 : out = in;
279 0 : outFlag = inFlag;
280 0 : apply(out, outFlag);
281 0 : }
282 :
283 : // print it out
284 : // unused
285 : /*
286 : friend std::ostream& operator<<(std::ostream& os, const SkyCal<DataType, CalDataType>& mat)
287 : {
288 : return os;
289 : }
290 : */
291 :
292 : protected:
293 :
294 : private:
295 : casacore::uInt npol_;
296 : casacore::uInt nchan_;
297 :
298 : // Pointer to origin
299 : CalDataType *m0_;
300 : casacore::Bool *ok0_;
301 :
302 : // Moving pointer
303 : CalDataType *m_, *mi_;
304 : casacore::Bool *ok_, *oki_;
305 :
306 : // casacore::Complex unity, zero (for use in invert and similar methods)
307 : const CalDataType cOne_,cZero_;
308 :
309 : mutable casacore::Bool scalardata_;
310 :
311 : // Copy ctor protected
312 : SkyCal(const SkyCal<DataType, CalDataType>& mat);
313 : };
314 :
315 : } //# NAMESPACE CASA - END
316 :
317 : //#include <synthesis/MeasurementComponents/SkyCal.tcc>
318 :
319 : #endif /* _SYNTHESIS_SKY_CAL_H_ */
320 :
321 :
322 :
|
