source: trunk/src/SDContainer.cc @ 83

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

chnaged namesapce form "atnf_sd" to "asap"

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