source: trunk/src/SDContainer.cc@ 240

Last change on this file since 240 was 215, checked in by mar637, 20 years ago

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