source: trunk/src/SDContainer.cc@ 314

Last change on this file since 314 was 308, checked in by kil064, 20 years ago

use function 'near' in function SDFreqTable:addFrequency

to test for pre-exitsing entgries

fix error with refPix being Int instead of double in SDFreqTable

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 10.7 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/Exceptions.h>
33#include <tables/Tables/Table.h>
34#include <casa/Arrays/IPosition.h>
35#include <casa/Arrays/Matrix.h>
36#include <casa/Arrays/ArrayAccessor.h>
37#include <casa/BasicMath/Math.h>
38#include <casa/Quanta/MVTime.h>
39
40
41#include "SDDefs.h"
42#include "SDContainer.h"
43
44using namespace casa;
45using namespace asap;
46
47void SDHeader::print() const {
48 MVTime mvt(this->utc);
49 mvt.setFormat(MVTime::YMD);
50 cout << "Observer: " << this->observer << endl
51 << "Project: " << this->project << endl
52 << "Obstype: " << this->obstype << endl
53 << "Antenna: " << this->antennaname << endl
54 << "Ant. Position: " << this->antennaposition << endl
55 << "Equinox: " << this->equinox << endl
56 << "Freq. ref.: " << this->freqref << endl
57 << "Ref. frequency: " << this->reffreq << endl
58 << "Bandwidth: " << this->bandwidth << endl
59 << "Time (utc): "
60 << mvt
61 << endl;
62 //setprecision(10) << this->utc << endl;
63}
64
65
66SDContainer::SDContainer(uInt nBeam, uInt nIF, uInt nPol, uInt nChan)
67 : nBeam_(nBeam),
68 nIF_(nIF),
69 nPol_(nPol),
70 nChan_(nChan),
71 spectrum_(IPosition(4,nBeam,nIF,nPol,nChan)),
72 flags_(IPosition(4,nBeam,nIF,nPol,nChan)),
73 tsys_(IPosition(4,nBeam,nIF,nPol,nChan)),
74 freqidx_(nIF),
75 direction_(IPosition(2,nBeam,2)) {
76 uChar x = 0;
77 flags_ = ~x;
78 tcal.resize(2);
79}
80
81SDContainer::SDContainer(IPosition shp)
82 : nBeam_(shp(0)),
83 nIF_(shp(1)),
84 nPol_(shp(2)),
85 nChan_(shp(3)),
86 spectrum_(shp),
87 flags_(shp),
88 tsys_(shp),
89 freqidx_(shp(1)) {
90 IPosition ip(2,shp(0),2);
91 direction_.resize(ip);
92 uChar x = 0;
93 flags_ = ~x;
94 tcal.resize(2);
95}
96
97SDContainer::~SDContainer() {
98}
99
100Bool SDContainer::resize(IPosition shp) {
101 nBeam_ = shp(0);
102 nIF_ = shp(1);
103 nPol_ = shp(2);
104 nChan_ = shp(3);
105 spectrum_.resize(shp);
106 flags_.resize(shp);
107 tsys_.resize(shp);
108 freqidx_.resize(shp(1));
109 IPosition ip(2,shp(0),2);
110 direction_.resize(ip);
111}
112
113Bool SDContainer::putSpectrum(const Array<Float>& spec) {
114 spectrum_ = spec;
115}
116Bool SDContainer::putFlags(const Array<uChar>& flag) {
117 flags_ = flag;
118}
119Bool SDContainer::putTsys(const Array<Float>& tsys) {
120 tsys_ = tsys;
121}
122
123Bool SDContainer::setSpectrum(const Matrix<Float>& spec,
124 uInt whichBeam, uInt whichIF) {
125
126 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(spectrum_);
127 aa0.reset(aa0.begin(whichBeam));
128 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
129 aa1.reset(aa1.begin(whichIF));
130
131 //Vector<Float> pols(nPol);
132 ArrayAccessor<Float, Axis<asap::IFAxis> > j(spec);
133 IPosition shp0 = spectrum_.shape();
134 IPosition shp1 = spec.shape();
135 if ( (shp0(2) != shp1(1)) || (shp0(3) != shp1(0)) ) {
136 throw(AipsError("Arrays not conformant"));
137 return False;
138 }
139 // assert dimensions are the same....
140 for (ArrayAccessor<Float, Axis<asap::PolAxis> > i(aa1);i != i.end(); ++i) {
141 ArrayAccessor<Float, Axis<asap::BeamAxis> > jj(j);
142 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > ii(i);ii != ii.end(); ++ii) {
143 (*ii) = (*jj);
144 jj++;
145 }
146 j++;
147 }
148 // unset flags for this spectrum, they might be set again by the
149 // setFlags method
150
151 IPosition shp = flags_.shape();
152 IPosition start(4,whichBeam,whichIF,0,0);
153 IPosition end(4,whichBeam,whichIF,shp(2)-1,shp(3)-1);
154 Array<uChar> arr(flags_(start,end));
155 arr = uChar(0);
156}
157
158Bool SDContainer::setFlags(const Matrix<uChar>& flag,
159 uInt whichBeam, uInt whichIF) {
160
161 ArrayAccessor<uChar, Axis<asap::BeamAxis> > aa0(flags_);
162 aa0.reset(aa0.begin(whichBeam));
163 ArrayAccessor<uChar, Axis<asap::IFAxis> > aa1(aa0);
164 aa1.reset(aa1.begin(whichIF));
165
166 ArrayAccessor<uChar, Axis<asap::IFAxis> > j(flag);
167 IPosition shp0 = flags_.shape();
168 IPosition shp1 = flag.shape();
169 if ( (shp0(2) != shp1(1)) || (shp0(3) != shp1(0)) ) {
170 cerr << "Arrays not conformant" << endl;
171 return False;
172 }
173
174 // assert dimensions are the same....
175 for (ArrayAccessor<uChar, Axis<asap::PolAxis> > i(aa1);i != i.end(); ++i) {
176 ArrayAccessor<uChar, Axis<asap::BeamAxis> > jj(j);
177 for (ArrayAccessor<uChar, Axis<asap::ChanAxis> > ii(i);ii != ii.end(); ++ii) {
178 (*ii) = (*jj);
179 jj++;
180 }
181 j++;
182 }
183 return True;
184}
185
186Bool SDContainer::setTsys(const Vector<Float>& tsys,
187 uInt whichBeam, uInt whichIF) {
188 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(tsys_);
189 aa0.reset(aa0.begin(whichBeam));
190 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
191 aa1.reset(aa1.begin(whichIF));
192 // assert dimensions are the same....
193 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > i(aa1);i != i.end(); ++i) {
194 ArrayAccessor<Float, Axis<asap::BeamAxis> > j(tsys);
195 for (ArrayAccessor<Float, Axis<asap::PolAxis> > ii(i);ii != ii.end(); ++ii) {
196 (*ii) = (*j);
197 j++;
198 }
199 }
200}
201
202Array<Float> SDContainer::getSpectrum(uInt whichBeam, uInt whichIF) const
203{
204 Matrix<Float> spectra(nChan_, nPol_);
205
206 // Beam.
207 ArrayAccessor<Float, Axis<asap::BeamAxis> > i0(spectrum_);
208 i0.reset(i0.begin(whichBeam));
209
210 // IF.
211 ArrayAccessor<Float, Axis<asap::IFAxis> > i1(i0);
212 i1.reset(i1.begin(whichIF));
213
214 // Polarization.
215 ArrayAccessor<Float, Axis<asap::PolAxis> > i2(i1);
216 ArrayAccessor<Float, Axis<asap::IFAxis> > o1(spectra);
217
218 while (i2 != i2.end()) {
219 // Channel.
220 ArrayAccessor<Float, Axis<asap::ChanAxis> > i3(i2);
221 ArrayAccessor<Float, Axis<asap::BeamAxis> > o0(o1);
222
223 while (i3 != i3.end()) {
224 *o0 = *i3;
225
226 i3++;
227 o0++;
228 }
229
230 i2++;
231 o1++;
232 }
233
234 return spectra.copy();
235}
236
237
238Array<uChar> SDContainer::getFlags(uInt whichBeam, uInt whichIF) const
239{
240 Matrix<uChar> flagtra(nChan_, nPol_);
241
242 // Beam.
243 ArrayAccessor<uChar, Axis<asap::BeamAxis> > i0(flags_);
244 i0.reset(i0.begin(whichBeam));
245
246 // IF.
247 ArrayAccessor<uChar, Axis<asap::IFAxis> > i1(i0);
248 i1.reset(i1.begin(whichIF));
249
250 // Polarization.
251 ArrayAccessor<uChar, Axis<asap::PolAxis> > i2(i1);
252 ArrayAccessor<uChar, Axis<asap::IFAxis> > o1(flagtra);
253
254 while (i2 != i2.end()) {
255 // Channel.
256 ArrayAccessor<uChar, Axis<asap::ChanAxis> > i3(i2);
257 ArrayAccessor<uChar, Axis<asap::BeamAxis> > o0(o1);
258
259 while (i3 != i3.end()) {
260 *o0 = *i3;
261
262 i3++;
263 o0++;
264 }
265
266 i2++;
267 o1++;
268 }
269
270 return flagtra.copy();
271}
272
273Array<Float> SDContainer::getTsys(uInt whichBeam, uInt whichIF) const
274{
275 Vector<Float> tsys(nPol_);
276
277 // Beam.
278 ArrayAccessor<Float, Axis<asap::BeamAxis> > i0(tsys_);
279 i0.reset(i0.begin(whichBeam));
280
281 // IF.
282 ArrayAccessor<Float, Axis<asap::IFAxis> > i1(i0);
283 i1.reset(i1.begin(whichIF));
284
285 // Channel.
286 ArrayAccessor<Float, Axis<asap::ChanAxis> > i3(i1);
287
288 // Polarization.
289 ArrayAccessor<Float, Axis<asap::PolAxis> > i2(i3);
290 ArrayAccessor<Float, Axis<asap::BeamAxis> > o0(tsys);
291
292 while (i2 != i2.end()) {
293 *o0 = *i2;
294
295 i2++;
296 o0++;
297 }
298 return tsys.copy();
299}
300
301Array<Double> SDContainer::getDirection(uInt whichBeam) const {
302 Vector<Double> direct(2);
303 ArrayAccessor<Double, Axis<asap::BeamAxis> > i0(direction_);
304 i0.reset(i0.begin(whichBeam));
305 ArrayAccessor<Double, Axis<asap::BeamAxis> > o0(direct);
306 ArrayAccessor<Double, Axis<asap::IFAxis> > i1(i0);
307 while (i1 != i1.end()) {
308 *o0 = *i1;
309 i1++;
310 o0++;
311 }
312 return direct.copy();
313}
314
315
316Bool SDContainer::setFrequencyMap(uInt freqslot, uInt whichIF) {
317 freqidx_[whichIF] = freqslot;
318 return True;
319}
320
321Bool SDContainer::putFreqMap(const Vector<uInt>& freqs) {
322 freqidx_.resize();
323 freqidx_ = freqs;
324 return True;
325}
326
327Bool SDContainer::setDirection(const Vector<Double>& point, uInt whichBeam) {
328 if (point.nelements() != 2) return False;
329 ArrayAccessor<Double, Axis<asap::BeamAxis> > aa0(direction_);
330 aa0.reset(aa0.begin(whichBeam));
331 ArrayAccessor<Double, Axis<asap::BeamAxis> > jj(point);
332 for (ArrayAccessor<Double, Axis<asap::IFAxis> > i(aa0);i != i.end(); ++i) {
333
334 (*i) = (*jj);
335 jj++;
336 }
337 return True;
338}
339
340Bool SDContainer::putDirection(const Array<Double>& dir) {
341 direction_.resize();
342 direction_ = dir;
343 return True;
344}
345
346Bool SDContainer::appendHistory(const String& hist)
347{
348 history_.resize(history_.nelements()+1,True);
349 history_[history_.nelements()-1] = hist;
350 return True;
351}
352Bool SDContainer::putHistory(const Vector<String>& hist)
353{
354 history_.resize();
355 history_ = hist;
356 return True;
357}
358
359Int SDFrequencyTable::addFrequency(Double refPix, Double refVal, Double inc) {
360 Int idx = -1;
361 Bool addit = False;
362 if (length() > 0) {
363 for (uInt i=0; i< length();++i) {
364 if ( near(refVal,refVal_[i]) ) {
365 if (near(refPix,refPix_[i]) )
366 if ( near(inc,increment_[i]) )
367 idx = Int(i);
368 }
369 }
370 if (idx >= 0) {
371 return idx;
372 }
373 }
374 nFreq_ += 1;
375 refPix_.resize(nFreq_,True);
376 refVal_.resize(nFreq_,True);
377 increment_.resize(nFreq_,True);
378 refPix_[nFreq_-1] = refPix;
379 refVal_[nFreq_-1] = refVal;
380 increment_[nFreq_-1] = inc;
381 idx = nFreq_-1;
382 return idx;
383}
384
385void SDFrequencyTable::addRestFrequency(Double val)
386{
387 if (restFreqs_.nelements() == 0) {
388 restFreqs_.resize(1);
389 restFreqs_[0] = val;
390 } else {
391 Bool found = False;
392 for (uInt i=0;i<restFreqs_.nelements();++i) {
393 if (restFreqs_[i] == val) {
394 found = True;
395 return;
396 }
397 }
398 if (!found) {
399 restFreqs_.resize(restFreqs_.nelements()+1,True);
400 restFreqs_[restFreqs_.nelements()-1] = val;
401 }
402 }
403}
404void SDFrequencyTable::restFrequencies(Vector<Double>& rfs,
405 String& rfunit ) const
406{
407 rfs.resize(restFreqs_.nelements());
408 rfs = restFreqs_;
409 rfunit = restFreqUnit_;
410}
Note: See TracBrowser for help on using the repository browser.