source: trunk/src/SDContainer.cc@ 335

Last change on this file since 335 was 327, checked in by kil064, 20 years ago

format better in SDDataDesc::summary
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 11.4 KB
Line 
1//#---------------------------------------------------------------------------
2//# SDContainer.cc: A container class for single dish integrations
3//#---------------------------------------------------------------------------
4//# Copyright (C) 2004
5//# ATNF
6//#
7//# This program is free software; you can redistribute it and/or modify it
8//# under the terms of the GNU General Public License as published by the Free
9//# Software Foundation; either version 2 of the License, or (at your option)
10//# any later version.
11//#
12//# This program is distributed in the hope that it will be useful, but
13//# WITHOUT ANY WARRANTY; without even the implied warranty of
14//# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
15//# Public License for more details.
16//#
17//# You should have received a copy of the GNU General Public License along
18//# with this program; if not, write to the Free Software Foundation, Inc.,
19//# 675 Massachusetts Ave, Cambridge, MA 02139, USA.
20//#
21//# Correspondence concerning this software should be addressed as follows:
22//# Internet email: Malte.Marquarding@csiro.au
23//# Postal address: Malte Marquarding,
24//# Australia Telescope National Facility,
25//# P.O. Box 76,
26//# Epping, NSW, 2121,
27//# AUSTRALIA
28//#
29//# $Id:
30//#---------------------------------------------------------------------------
31#include <casa/aips.h>
32#include <casa/iostream.h>
33#include <casa/iomanip.h>
34#include <casa/Exceptions.h>
35#include <tables/Tables/Table.h>
36#include <casa/Arrays/IPosition.h>
37#include <casa/Arrays/Matrix.h>
38#include <casa/Arrays/ArrayAccessor.h>
39#include <casa/BasicMath/Math.h>
40#include <casa/Quanta/MVTime.h>
41
42
43#include "SDDefs.h"
44#include "SDContainer.h"
45
46using namespace casa;
47using namespace asap;
48
49void SDHeader::print() const {
50 MVTime mvt(this->utc);
51 mvt.setFormat(MVTime::YMD);
52 cout << "Observer: " << this->observer << endl
53 << "Project: " << this->project << endl
54 << "Obstype: " << this->obstype << endl
55 << "Antenna: " << this->antennaname << endl
56 << "Ant. Position: " << this->antennaposition << endl
57 << "Equinox: " << this->equinox << endl
58 << "Freq. ref.: " << this->freqref << endl
59 << "Ref. frequency: " << this->reffreq << endl
60 << "Bandwidth: " << this->bandwidth << endl
61 << "Time (utc): "
62 << mvt
63 << endl;
64 //setprecision(10) << this->utc << endl;
65}
66
67
68SDContainer::SDContainer(uInt nBeam, uInt nIF, uInt nPol, uInt nChan)
69 : nBeam_(nBeam),
70 nIF_(nIF),
71 nPol_(nPol),
72 nChan_(nChan),
73 spectrum_(IPosition(4,nBeam,nIF,nPol,nChan)),
74 flags_(IPosition(4,nBeam,nIF,nPol,nChan)),
75 tsys_(IPosition(4,nBeam,nIF,nPol,nChan)),
76 freqidx_(nIF),
77 direction_(IPosition(2,nBeam,2)) {
78 uChar x = 0;
79 flags_ = ~x;
80 tcal.resize(2);
81}
82
83SDContainer::SDContainer(IPosition shp)
84 : nBeam_(shp(0)),
85 nIF_(shp(1)),
86 nPol_(shp(2)),
87 nChan_(shp(3)),
88 spectrum_(shp),
89 flags_(shp),
90 tsys_(shp),
91 freqidx_(shp(1)) {
92 IPosition ip(2,shp(0),2);
93 direction_.resize(ip);
94 uChar x = 0;
95 flags_ = ~x;
96 tcal.resize(2);
97}
98
99SDContainer::~SDContainer() {
100}
101
102Bool SDContainer::resize(IPosition shp) {
103 nBeam_ = shp(0);
104 nIF_ = shp(1);
105 nPol_ = shp(2);
106 nChan_ = shp(3);
107 spectrum_.resize(shp);
108 flags_.resize(shp);
109 tsys_.resize(shp);
110 freqidx_.resize(shp(1));
111 IPosition ip(2,shp(0),2);
112 direction_.resize(ip);
113}
114
115Bool SDContainer::putSpectrum(const Array<Float>& spec) {
116 spectrum_ = spec;
117}
118Bool SDContainer::putFlags(const Array<uChar>& flag) {
119 flags_ = flag;
120}
121Bool SDContainer::putTsys(const Array<Float>& tsys) {
122 tsys_ = tsys;
123}
124
125Bool SDContainer::setSpectrum(const Matrix<Float>& spec,
126 uInt whichBeam, uInt whichIF) {
127
128 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(spectrum_);
129 aa0.reset(aa0.begin(whichBeam));
130 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
131 aa1.reset(aa1.begin(whichIF));
132
133 //Vector<Float> pols(nPol);
134 ArrayAccessor<Float, Axis<asap::IFAxis> > j(spec);
135 IPosition shp0 = spectrum_.shape();
136 IPosition shp1 = spec.shape();
137 if ( (shp0(2) != shp1(1)) || (shp0(3) != shp1(0)) ) {
138 throw(AipsError("Arrays not conformant"));
139 return False;
140 }
141 // assert dimensions are the same....
142 for (ArrayAccessor<Float, Axis<asap::PolAxis> > i(aa1);i != i.end(); ++i) {
143 ArrayAccessor<Float, Axis<asap::BeamAxis> > jj(j);
144 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > ii(i);ii != ii.end(); ++ii) {
145 (*ii) = (*jj);
146 jj++;
147 }
148 j++;
149 }
150 // unset flags for this spectrum, they might be set again by the
151 // setFlags method
152
153 IPosition shp = flags_.shape();
154 IPosition start(4,whichBeam,whichIF,0,0);
155 IPosition end(4,whichBeam,whichIF,shp(2)-1,shp(3)-1);
156 Array<uChar> arr(flags_(start,end));
157 arr = uChar(0);
158}
159
160Bool SDContainer::setFlags(const Matrix<uChar>& flag,
161 uInt whichBeam, uInt whichIF) {
162
163 ArrayAccessor<uChar, Axis<asap::BeamAxis> > aa0(flags_);
164 aa0.reset(aa0.begin(whichBeam));
165 ArrayAccessor<uChar, Axis<asap::IFAxis> > aa1(aa0);
166 aa1.reset(aa1.begin(whichIF));
167
168 ArrayAccessor<uChar, Axis<asap::IFAxis> > j(flag);
169 IPosition shp0 = flags_.shape();
170 IPosition shp1 = flag.shape();
171 if ( (shp0(2) != shp1(1)) || (shp0(3) != shp1(0)) ) {
172 cerr << "Arrays not conformant" << endl;
173 return False;
174 }
175
176 // assert dimensions are the same....
177 for (ArrayAccessor<uChar, Axis<asap::PolAxis> > i(aa1);i != i.end(); ++i) {
178 ArrayAccessor<uChar, Axis<asap::BeamAxis> > jj(j);
179 for (ArrayAccessor<uChar, Axis<asap::ChanAxis> > ii(i);ii != ii.end(); ++ii) {
180 (*ii) = (*jj);
181 jj++;
182 }
183 j++;
184 }
185 return True;
186}
187
188Bool SDContainer::setTsys(const Vector<Float>& tsys,
189 uInt whichBeam, uInt whichIF) {
190 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(tsys_);
191 aa0.reset(aa0.begin(whichBeam));
192 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
193 aa1.reset(aa1.begin(whichIF));
194 // assert dimensions are the same....
195 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > i(aa1);i != i.end(); ++i) {
196 ArrayAccessor<Float, Axis<asap::BeamAxis> > j(tsys);
197 for (ArrayAccessor<Float, Axis<asap::PolAxis> > ii(i);ii != ii.end(); ++ii) {
198 (*ii) = (*j);
199 j++;
200 }
201 }
202}
203
204Array<Float> SDContainer::getSpectrum(uInt whichBeam, uInt whichIF) const
205{
206 Matrix<Float> spectra(nChan_, nPol_);
207
208 // Beam.
209 ArrayAccessor<Float, Axis<asap::BeamAxis> > i0(spectrum_);
210 i0.reset(i0.begin(whichBeam));
211
212 // IF.
213 ArrayAccessor<Float, Axis<asap::IFAxis> > i1(i0);
214 i1.reset(i1.begin(whichIF));
215
216 // Polarization.
217 ArrayAccessor<Float, Axis<asap::PolAxis> > i2(i1);
218 ArrayAccessor<Float, Axis<asap::IFAxis> > o1(spectra);
219
220 while (i2 != i2.end()) {
221 // Channel.
222 ArrayAccessor<Float, Axis<asap::ChanAxis> > i3(i2);
223 ArrayAccessor<Float, Axis<asap::BeamAxis> > o0(o1);
224
225 while (i3 != i3.end()) {
226 *o0 = *i3;
227
228 i3++;
229 o0++;
230 }
231
232 i2++;
233 o1++;
234 }
235
236 return spectra.copy();
237}
238
239
240Array<uChar> SDContainer::getFlags(uInt whichBeam, uInt whichIF) const
241{
242 Matrix<uChar> flagtra(nChan_, nPol_);
243
244 // Beam.
245 ArrayAccessor<uChar, Axis<asap::BeamAxis> > i0(flags_);
246 i0.reset(i0.begin(whichBeam));
247
248 // IF.
249 ArrayAccessor<uChar, Axis<asap::IFAxis> > i1(i0);
250 i1.reset(i1.begin(whichIF));
251
252 // Polarization.
253 ArrayAccessor<uChar, Axis<asap::PolAxis> > i2(i1);
254 ArrayAccessor<uChar, Axis<asap::IFAxis> > o1(flagtra);
255
256 while (i2 != i2.end()) {
257 // Channel.
258 ArrayAccessor<uChar, Axis<asap::ChanAxis> > i3(i2);
259 ArrayAccessor<uChar, Axis<asap::BeamAxis> > o0(o1);
260
261 while (i3 != i3.end()) {
262 *o0 = *i3;
263
264 i3++;
265 o0++;
266 }
267
268 i2++;
269 o1++;
270 }
271
272 return flagtra.copy();
273}
274
275Array<Float> SDContainer::getTsys(uInt whichBeam, uInt whichIF) const
276{
277 Vector<Float> tsys(nPol_);
278
279 // Beam.
280 ArrayAccessor<Float, Axis<asap::BeamAxis> > i0(tsys_);
281 i0.reset(i0.begin(whichBeam));
282
283 // IF.
284 ArrayAccessor<Float, Axis<asap::IFAxis> > i1(i0);
285 i1.reset(i1.begin(whichIF));
286
287 // Channel.
288 ArrayAccessor<Float, Axis<asap::ChanAxis> > i3(i1);
289
290 // Polarization.
291 ArrayAccessor<Float, Axis<asap::PolAxis> > i2(i3);
292 ArrayAccessor<Float, Axis<asap::BeamAxis> > o0(tsys);
293
294 while (i2 != i2.end()) {
295 *o0 = *i2;
296
297 i2++;
298 o0++;
299 }
300 return tsys.copy();
301}
302
303Array<Double> SDContainer::getDirection(uInt whichBeam) const {
304 Vector<Double> direct(2);
305 ArrayAccessor<Double, Axis<asap::BeamAxis> > i0(direction_);
306 i0.reset(i0.begin(whichBeam));
307 ArrayAccessor<Double, Axis<asap::BeamAxis> > o0(direct);
308 ArrayAccessor<Double, Axis<asap::IFAxis> > i1(i0);
309 while (i1 != i1.end()) {
310 *o0 = *i1;
311 i1++;
312 o0++;
313 }
314 return direct.copy();
315}
316
317
318Bool SDContainer::setFrequencyMap(uInt freqslot, uInt whichIF) {
319 freqidx_[whichIF] = freqslot;
320 return True;
321}
322
323Bool SDContainer::putFreqMap(const Vector<uInt>& freqs) {
324 freqidx_.resize();
325 freqidx_ = freqs;
326 return True;
327}
328
329Bool SDContainer::setDirection(const Vector<Double>& point, uInt whichBeam) {
330 if (point.nelements() != 2) return False;
331 ArrayAccessor<Double, Axis<asap::BeamAxis> > aa0(direction_);
332 aa0.reset(aa0.begin(whichBeam));
333 ArrayAccessor<Double, Axis<asap::BeamAxis> > jj(point);
334 for (ArrayAccessor<Double, Axis<asap::IFAxis> > i(aa0);i != i.end(); ++i) {
335
336 (*i) = (*jj);
337 jj++;
338 }
339 return True;
340}
341
342Bool SDContainer::putDirection(const Array<Double>& dir) {
343 direction_.resize();
344 direction_ = dir;
345 return True;
346}
347
348Bool SDContainer::appendHistory(const String& hist)
349{
350 history_.resize(history_.nelements()+1,True);
351 history_[history_.nelements()-1] = hist;
352 return True;
353}
354Bool SDContainer::putHistory(const Vector<String>& hist)
355{
356 history_.resize();
357 history_ = hist;
358 return True;
359}
360
361Int SDFrequencyTable::addFrequency(Double refPix, Double refVal, Double inc) {
362 Int idx = -1;
363 Bool addit = False;
364 if (length() > 0) {
365 for (uInt i=0; i< length();++i) {
366 if ( near(refVal,refVal_[i]) ) {
367 if (near(refPix,refPix_[i]) )
368 if ( near(inc,increment_[i]) )
369 idx = Int(i);
370 }
371 }
372 if (idx >= 0) {
373 return idx;
374 }
375 }
376 nFreq_ += 1;
377 refPix_.resize(nFreq_,True);
378 refVal_.resize(nFreq_,True);
379 increment_.resize(nFreq_,True);
380 refPix_[nFreq_-1] = refPix;
381 refVal_[nFreq_-1] = refVal;
382 increment_[nFreq_-1] = inc;
383 idx = nFreq_-1;
384 return idx;
385}
386
387void SDFrequencyTable::addRestFrequency(Double val)
388{
389 if (restFreqs_.nelements() == 0) {
390 restFreqs_.resize(1);
391 restFreqs_[0] = val;
392 } else {
393 Bool found = False;
394 for (uInt i=0;i<restFreqs_.nelements();++i) {
395 if (restFreqs_[i] == val) {
396 found = True;
397 return;
398 }
399 }
400 if (!found) {
401 restFreqs_.resize(restFreqs_.nelements()+1,True);
402 restFreqs_[restFreqs_.nelements()-1] = val;
403 }
404 }
405}
406void SDFrequencyTable::restFrequencies(Vector<Double>& rfs,
407 String& rfunit ) const
408{
409 rfs.resize(restFreqs_.nelements());
410 rfs = restFreqs_;
411 rfunit = restFreqUnit_;
412}
413
414
415uInt SDDataDesc::addEntry (const String& source, uInt freqID, const MDirection& dir)
416{
417
418// See if already exists
419
420 if (n_ > 0) {
421 for (uInt i=0; i<n_; i++) {
422 if (source==source_[i] && freqID==freqID_[i]) {
423 return i;
424 }
425 }
426 }
427
428// Not found - add it
429
430 n_ += 1;
431 source_.resize(n_,True);
432 freqID_.resize(n_,True);
433 dir_.resize(n_,True,True);
434//
435 source_[n_-1] = source;
436 freqID_[n_-1] = freqID;
437 dir_[n_-1] = dir;
438//
439 return n_-1;
440}
441
442
443void SDDataDesc::summary() const
444{
445 cerr << "Source FreqID" << endl;
446 for (uInt i=0; i<n_; i++) {
447 cerr << setw(11) << source_(i) << freqID_(i) << endl;
448 }
449}
450
Note: See TracBrowser for help on using the repository browser.