source: trunk/src/SDMemTable.cc@ 385

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

have another go at summary function. this time
rework with list of scan boundaries
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 38.1 KB
Line 
1//#---------------------------------------------------------------------------
2//# SDMemTable.cc: A MemoryTable container 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
32#include <map>
33
34#include <casa/aips.h>
35#include <casa/iostream.h>
36#include <casa/iomanip.h>
37#include <casa/Arrays/Array.h>
38#include <casa/Arrays/ArrayMath.h>
39#include <casa/Arrays/MaskArrMath.h>
40#include <casa/Arrays/ArrayLogical.h>
41#include <casa/Arrays/ArrayAccessor.h>
42#include <casa/Arrays/Vector.h>
43#include <casa/Quanta/MVAngle.h>
44
45#include <tables/Tables/TableParse.h>
46#include <tables/Tables/TableDesc.h>
47#include <tables/Tables/SetupNewTab.h>
48#include <tables/Tables/ScaColDesc.h>
49#include <tables/Tables/ArrColDesc.h>
50
51#include <tables/Tables/ExprNode.h>
52#include <tables/Tables/TableRecord.h>
53#include <measures/Measures/MFrequency.h>
54#include <measures/Measures/MeasTable.h>
55#include <coordinates/Coordinates/CoordinateUtil.h>
56#include <casa/Quanta/MVTime.h>
57#include <casa/Quanta/MVAngle.h>
58
59#include "SDDefs.h"
60#include "SDMemTable.h"
61#include "SDContainer.h"
62#include "MathUtils.h"
63
64
65using namespace casa;
66using namespace asap;
67
68SDMemTable::SDMemTable() :
69 IFSel_(0),
70 beamSel_(0),
71 polSel_(0)
72{
73 setup();
74 attach();
75}
76
77SDMemTable::SDMemTable(const std::string& name) :
78 IFSel_(0),
79 beamSel_(0),
80 polSel_(0)
81{
82 Table tab(name);
83 table_ = tab.copyToMemoryTable("dummy");
84 //cerr << "hello from C SDMemTable @ " << this << endl;
85 attach();
86}
87
88SDMemTable::SDMemTable(const SDMemTable& other, Bool clear)
89{
90 IFSel_= other.IFSel_;
91 beamSel_= other.beamSel_;
92 polSel_= other.polSel_;
93 chanMask_ = other.chanMask_;
94 table_ = other.table_.copyToMemoryTable(String("dummy"));
95 // clear all rows()
96 if (clear) {
97 table_.removeRow(this->table_.rowNumbers());
98 } else {
99 IFSel_ = other.IFSel_;
100 beamSel_ = other.beamSel_;
101 polSel_ = other.polSel_;
102 }
103//
104 attach();
105 //cerr << "hello from CC SDMemTable @ " << this << endl;
106}
107
108SDMemTable::SDMemTable(const Table& tab, const std::string& exprs) :
109 IFSel_(0),
110 beamSel_(0),
111 polSel_(0)
112{
113 cout << exprs << endl;
114 Table t = tableCommand(exprs,tab);
115 if (t.nrow() == 0)
116 throw(AipsError("Query unsuccessful."));
117 table_ = t.copyToMemoryTable("dummy");
118 attach();
119 renumber();
120}
121
122SDMemTable::~SDMemTable()
123{
124 //cerr << "goodbye from SDMemTable @ " << this << endl;
125}
126
127SDMemTable SDMemTable::getScan(Int scanID) const
128{
129 String cond("SELECT * from $1 WHERE SCANID == ");
130 cond += String::toString(scanID);
131 return SDMemTable(table_, cond);
132}
133
134SDMemTable &SDMemTable::operator=(const SDMemTable& other)
135{
136 if (this != &other) {
137 IFSel_= other.IFSel_;
138 beamSel_= other.beamSel_;
139 polSel_= other.polSel_;
140 chanMask_.resize(0);
141 chanMask_ = other.chanMask_;
142 table_ = other.table_.copyToMemoryTable(String("dummy"));
143 attach();
144 }
145 //cerr << "hello from ASS SDMemTable @ " << this << endl;
146 return *this;
147}
148
149SDMemTable SDMemTable::getSource(const std::string& source) const
150{
151 String cond("SELECT * from $1 WHERE SRCNAME == ");
152 cond += source;
153 return SDMemTable(table_, cond);
154}
155
156void SDMemTable::setup()
157{
158 TableDesc td("", "1", TableDesc::Scratch);
159 td.comment() = "A SDMemTable";
160//
161 td.addColumn(ScalarColumnDesc<Double>("TIME"));
162 td.addColumn(ScalarColumnDesc<String>("SRCNAME"));
163 td.addColumn(ArrayColumnDesc<Float>("SPECTRA"));
164 td.addColumn(ArrayColumnDesc<uChar>("FLAGTRA"));
165 td.addColumn(ArrayColumnDesc<Float>("TSYS"));
166 td.addColumn(ScalarColumnDesc<Int>("SCANID"));
167 td.addColumn(ScalarColumnDesc<Double>("INTERVAL"));
168 td.addColumn(ArrayColumnDesc<uInt>("FREQID"));
169 td.addColumn(ArrayColumnDesc<Double>("DIRECTION"));
170 td.addColumn(ScalarColumnDesc<String>("FIELDNAME"));
171 td.addColumn(ScalarColumnDesc<String>("TCALTIME"));
172 td.addColumn(ArrayColumnDesc<Float>("TCAL"));
173 td.addColumn(ScalarColumnDesc<Float>("AZIMUTH"));
174 td.addColumn(ScalarColumnDesc<Float>("ELEVATION"));
175 td.addColumn(ScalarColumnDesc<Float>("PARANGLE"));
176 td.addColumn(ScalarColumnDesc<Int>("REFBEAM"));
177 td.addColumn(ArrayColumnDesc<String>("HISTORY"));
178
179 // Now create a new table from the description.
180
181 SetupNewTable aNewTab("dummy", td, Table::New);
182 table_ = Table(aNewTab, Table::Memory, 0);
183}
184
185void SDMemTable::attach()
186{
187 timeCol_.attach(table_, "TIME");
188 srcnCol_.attach(table_, "SRCNAME");
189 specCol_.attach(table_, "SPECTRA");
190 flagsCol_.attach(table_, "FLAGTRA");
191 tsCol_.attach(table_, "TSYS");
192 scanCol_.attach(table_, "SCANID");
193 integrCol_.attach(table_, "INTERVAL");
194 freqidCol_.attach(table_, "FREQID");
195 dirCol_.attach(table_, "DIRECTION");
196 fldnCol_.attach(table_, "FIELDNAME");
197 tcaltCol_.attach(table_, "TCALTIME");
198 tcalCol_.attach(table_, "TCAL");
199 azCol_.attach(table_, "AZIMUTH");
200 elCol_.attach(table_, "ELEVATION");
201 paraCol_.attach(table_, "PARANGLE");
202 rbeamCol_.attach(table_, "REFBEAM");
203 histCol_.attach(table_, "HISTORY");
204}
205
206
207std::string SDMemTable::getSourceName(Int whichRow) const
208{
209 String name;
210 srcnCol_.get(whichRow, name);
211 return name;
212}
213
214std::string SDMemTable::getTime(Int whichRow, Bool showDate) const
215{
216 Double tm;
217 if (whichRow > -1) {
218 timeCol_.get(whichRow, tm);
219 } else {
220 table_.keywordSet().get("UTC",tm);
221 }
222 MVTime mvt(tm);
223 if (showDate)
224 mvt.setFormat(MVTime::YMD);
225 else
226 mvt.setFormat(MVTime::TIME);
227 ostringstream oss;
228 oss << mvt;
229 return String(oss);
230}
231
232double SDMemTable::getInterval(Int whichRow) const
233{
234 Double intval;
235 integrCol_.get(whichRow, intval);
236 return intval;
237}
238
239bool SDMemTable::setIF(Int whichIF)
240{
241 if ( whichIF >= 0 && whichIF < nIF()) {
242 IFSel_ = whichIF;
243 return true;
244 }
245 return false;
246}
247
248bool SDMemTable::setBeam(Int whichBeam)
249{
250 if ( whichBeam >= 0 && whichBeam < nBeam()) {
251 beamSel_ = whichBeam;
252 return true;
253 }
254 return false;
255}
256
257bool SDMemTable::setPol(Int whichPol)
258{
259 if ( whichPol >= 0 && whichPol < nPol()) {
260 polSel_ = whichPol;
261 return true;
262 }
263 return false;
264}
265
266void SDMemTable::resetCursor()
267{
268 polSel_ = 0;
269 IFSel_ = 0;
270 beamSel_ = 0;
271}
272
273bool SDMemTable::setMask(std::vector<int> whichChans)
274{
275 std::vector<int>::iterator it;
276 uInt n = flagsCol_.shape(0)(3);
277 if (whichChans.empty()) {
278 chanMask_ = std::vector<bool>(n,true);
279 return true;
280 }
281 chanMask_.resize(n,true);
282 for (it = whichChans.begin(); it != whichChans.end(); ++it) {
283 if (*it < n) {
284 chanMask_[*it] = false;
285 }
286 }
287 return true;
288}
289
290std::vector<bool> SDMemTable::getMask(Int whichRow) const {
291 std::vector<bool> mask;
292 Array<uChar> arr;
293 flagsCol_.get(whichRow, arr);
294 ArrayAccessor<uChar, Axis<asap::BeamAxis> > aa0(arr);
295 aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
296 ArrayAccessor<uChar, Axis<asap::IFAxis> > aa1(aa0);
297 aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
298 ArrayAccessor<uChar, Axis<asap::PolAxis> > aa2(aa1);
299 aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
300
301 Bool useUserMask = ( chanMask_.size() == arr.shape()(3) );
302
303 std::vector<bool> tmp;
304 tmp = chanMask_; // WHY the fxxx do I have to make a copy here
305 std::vector<bool>::iterator miter;
306 miter = tmp.begin();
307
308 for (ArrayAccessor<uChar, Axis<asap::ChanAxis> > i(aa2); i != i.end(); ++i) {
309 bool out =!static_cast<bool>(*i);
310 if (useUserMask) {
311 out = out && (*miter);
312 miter++;
313 }
314 mask.push_back(out);
315 }
316 return mask;
317}
318std::vector<float> SDMemTable::getSpectrum(Int whichRow) const
319{
320 std::vector<float> spectrum;
321 Array<Float> arr;
322 specCol_.get(whichRow, arr);
323 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(arr);
324 aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
325 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
326 aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
327 ArrayAccessor<Float, Axis<asap::PolAxis> > aa2(aa1);
328 aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
329 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > i(aa2); i != i.end(); ++i) {
330 spectrum.push_back(*i);
331 }
332 return spectrum;
333}
334std::vector<string> SDMemTable::getCoordInfo() const
335{
336 String un;
337 Table t = table_.keywordSet().asTable("FREQUENCIES");
338 String sunit;
339 t.keywordSet().get("UNIT",sunit);
340 String dpl;
341 t.keywordSet().get("DOPPLER",dpl);
342 if (dpl == "") dpl = "RADIO";
343 String rfrm;
344 t.keywordSet().get("REFFRAME",rfrm);
345 std::vector<string> inf;
346 inf.push_back(sunit);
347 inf.push_back(rfrm);
348 inf.push_back(dpl);
349 t.keywordSet().get("BASEREFFRAME",rfrm);
350 inf.push_back(rfrm);
351 return inf;
352}
353
354void SDMemTable::setCoordInfo(std::vector<string> theinfo)
355{
356 std::vector<string>::iterator it;
357 String un,rfrm, brfrm,dpl;
358 un = theinfo[0]; // Abcissa unit
359 rfrm = theinfo[1]; // Active (or conversion) frame
360 dpl = theinfo[2]; // Doppler
361 brfrm = theinfo[3]; // Base frame
362//
363 Table t = table_.rwKeywordSet().asTable("FREQUENCIES");
364//
365 Vector<Double> rstf;
366 t.keywordSet().get("RESTFREQS",rstf);
367//
368 Bool canDo = True;
369 Unit u1("km/s");Unit u2("Hz");
370 if (Unit(un) == u1) {
371 Vector<Double> rstf;
372 t.keywordSet().get("RESTFREQS",rstf);
373 if (rstf.nelements() == 0) {
374 throw(AipsError("Can't set unit to km/s if no restfrequencies are specified"));
375 }
376 } else if (Unit(un) != u2 && un != "") {
377 throw(AipsError("Unit not conformant with Spectral Coordinates"));
378 }
379 t.rwKeywordSet().define("UNIT", un);
380//
381 MFrequency::Types mdr;
382 if (!MFrequency::getType(mdr, rfrm)) {
383
384 Int a,b;const uInt* c;
385 const String* valid = MFrequency::allMyTypes(a, b, c);
386 String pfix = "Please specify a legal frame type. Types are\n";
387 throw(AipsError(pfix+(*valid)));
388 } else {
389 t.rwKeywordSet().define("REFFRAME",rfrm);
390 }
391//
392 MDoppler::Types dtype;
393 dpl.upcase();
394 if (!MDoppler::getType(dtype, dpl)) {
395 throw(AipsError("Doppler type unknown"));
396 } else {
397 t.rwKeywordSet().define("DOPPLER",dpl);
398 }
399//
400 if (!MFrequency::getType(mdr, brfrm)) {
401 Int a,b;const uInt* c;
402 const String* valid = MFrequency::allMyTypes(a, b, c);
403 String pfix = "Please specify a legal frame type. Types are\n";
404 throw(AipsError(pfix+(*valid)));
405 } else {
406 t.rwKeywordSet().define("BASEREFFRAME",brfrm);
407 }
408}
409
410
411std::vector<double> SDMemTable::getAbcissa(Int whichRow) const
412{
413 std::vector<double> abc(nChan());
414
415// Get header units keyword
416
417 Table t = table_.keywordSet().asTable("FREQUENCIES");
418 String sunit;
419 t.keywordSet().get("UNIT",sunit);
420 if (sunit == "") sunit = "pixel";
421 Unit u(sunit);
422
423// Easy if just wanting pixels
424
425 if (sunit==String("pixel")) {
426 // assume channels/pixels
427 std::vector<double>::iterator it;
428 uInt i=0;
429 for (it = abc.begin(); it != abc.end(); ++it) {
430 (*it) = Double(i++);
431 }
432//
433 return abc;
434 }
435
436// Continue with km/s or Hz. Get FreqID
437
438 Vector<uInt> freqIDs;
439 freqidCol_.get(whichRow, freqIDs);
440 uInt freqID = freqIDs(IFSel_);
441
442// Get SpectralCoordinate, set reference frame conversion,
443// velocity conversion, and rest freq state
444
445 SpectralCoordinate spc = getSpectralCoordinate(freqID, whichRow);
446 Vector<Double> pixel(nChan());
447 indgen(pixel);
448//
449 if (u == Unit("km/s")) {
450 Vector<Double> world;
451 spc.pixelToVelocity(world,pixel);
452 std::vector<double>::iterator it;
453 uInt i = 0;
454 for (it = abc.begin(); it != abc.end(); ++it) {
455 (*it) = world[i];
456 i++;
457 }
458 } else if (u == Unit("Hz")) {
459
460// Set world axis units
461
462 Vector<String> wau(1); wau = u.getName();
463 spc.setWorldAxisUnits(wau);
464//
465 std::vector<double>::iterator it;
466 Double tmp;
467 uInt i = 0;
468 for (it = abc.begin(); it != abc.end(); ++it) {
469 spc.toWorld(tmp,pixel[i]);
470 (*it) = tmp;
471 i++;
472 }
473 }
474 return abc;
475}
476
477std::string SDMemTable::getAbcissaString(Int whichRow) const
478{
479 Table t = table_.keywordSet().asTable("FREQUENCIES");
480//
481 String sunit;
482 t.keywordSet().get("UNIT",sunit);
483 if (sunit == "") sunit = "pixel";
484 Unit u(sunit);
485//
486 Vector<uInt> freqIDs;
487 freqidCol_.get(whichRow, freqIDs);
488 uInt freqID = freqIDs(IFSel_);
489
490// Get SpectralCoordinate, with frame, velocity, rest freq state set
491
492 SpectralCoordinate spc = getSpectralCoordinate(freqID, whichRow);
493//
494 String s = "Channel";
495 if (u == Unit("km/s")) {
496 s = CoordinateUtil::axisLabel(spc,0,True,True,True);
497 } else if (u == Unit("Hz")) {
498 Vector<String> wau(1);wau = u.getName();
499 spc.setWorldAxisUnits(wau);
500//
501 s = CoordinateUtil::axisLabel(spc,0,True,True,False);
502 }
503 return s;
504}
505
506void SDMemTable::setSpectrum(std::vector<float> spectrum, int whichRow)
507{
508 Array<Float> arr;
509 specCol_.get(whichRow, arr);
510 if (spectrum.size() != arr.shape()(3)) {
511 throw(AipsError("Attempting to set spectrum with incorrect length."));
512 }
513
514 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(arr);
515 aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
516 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
517 aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
518 ArrayAccessor<Float, Axis<asap::PolAxis> > aa2(aa1);
519 aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
520
521 std::vector<float>::iterator it = spectrum.begin();
522 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > i(aa2); i != i.end(); ++i) {
523 (*i) = Float(*it);
524 it++;
525 }
526 specCol_.put(whichRow, arr);
527}
528
529void SDMemTable::getSpectrum(Vector<Float>& spectrum, Int whichRow) const
530{
531 Array<Float> arr;
532 specCol_.get(whichRow, arr);
533 spectrum.resize(arr.shape()(3));
534 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(arr);
535 aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
536 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
537 aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
538 ArrayAccessor<Float, Axis<asap::PolAxis> > aa2(aa1);
539 aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
540
541 ArrayAccessor<Float, Axis<asap::BeamAxis> > va(spectrum);
542 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > i(aa2); i != i.end(); ++i) {
543 (*va) = (*i);
544 va++;
545 }
546}
547/*
548void SDMemTable::getMask(Vector<Bool>& mask, Int whichRow) const {
549 Array<uChar> arr;
550 flagsCol_.get(whichRow, arr);
551 mask.resize(arr.shape()(3));
552
553 ArrayAccessor<uChar, Axis<asap::BeamAxis> > aa0(arr);
554 aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
555 ArrayAccessor<uChar, Axis<asap::IFAxis> > aa1(aa0);
556 aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
557 ArrayAccessor<uChar, Axis<asap::PolAxis> > aa2(aa1);
558 aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
559
560 Bool useUserMask = ( chanMask_.size() == arr.shape()(3) );
561
562 ArrayAccessor<Bool, Axis<asap::BeamAxis> > va(mask);
563 std::vector<bool> tmp;
564 tmp = chanMask_; // WHY the fxxx do I have to make a copy here. The
565 // iterator should work on chanMask_??
566 std::vector<bool>::iterator miter;
567 miter = tmp.begin();
568
569 for (ArrayAccessor<uChar, Axis<asap::ChanAxis> > i(aa2); i != i.end(); ++i) {
570 bool out =!static_cast<bool>(*i);
571 if (useUserMask) {
572 out = out && (*miter);
573 miter++;
574 }
575 (*va) = out;
576 va++;
577 }
578}
579*/
580MaskedArray<Float> SDMemTable::rowAsMaskedArray(uInt whichRow,
581 Bool useSelection) const
582{
583 Array<Float> arr;
584 Array<uChar> farr;
585 specCol_.get(whichRow, arr);
586 flagsCol_.get(whichRow, farr);
587 Array<Bool> barr(farr.shape());convertArray(barr, farr);
588 MaskedArray<Float> marr;
589 if (useSelection) {
590 ArrayAccessor<Float, Axis<asap::BeamAxis> > aa0(arr);
591 aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
592 ArrayAccessor<Float, Axis<asap::IFAxis> > aa1(aa0);
593 aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
594 ArrayAccessor<Float, Axis<asap::PolAxis> > aa2(aa1);
595 aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
596
597 ArrayAccessor<Bool, Axis<asap::BeamAxis> > baa0(barr);
598 baa0.reset(baa0.begin(uInt(beamSel_)));//go to beam
599 ArrayAccessor<Bool, Axis<asap::IFAxis> > baa1(baa0);
600 baa1.reset(baa1.begin(uInt(IFSel_)));// go to IF
601 ArrayAccessor<Bool, Axis<asap::PolAxis> > baa2(baa1);
602 baa2.reset(baa2.begin(uInt(polSel_)));// go to pol
603
604 Vector<Float> a(arr.shape()(3));
605 Vector<Bool> b(barr.shape()(3));
606 ArrayAccessor<Float, Axis<asap::BeamAxis> > a0(a);
607 ArrayAccessor<Bool, Axis<asap::BeamAxis> > b0(b);
608
609 ArrayAccessor<Bool, Axis<asap::ChanAxis> > j(baa2);
610 for (ArrayAccessor<Float, Axis<asap::ChanAxis> > i(aa2);
611 i != i.end(); ++i) {
612 (*a0) = (*i);
613 (*b0) = !(*j);
614 j++;
615 a0++;
616 b0++;
617 }
618 marr.setData(a,b);
619 } else {
620 marr.setData(arr,!barr);
621 }
622 return marr;
623}
624
625Float SDMemTable::getTsys(Int whichRow) const
626{
627 Array<Float> arr;
628 tsCol_.get(whichRow, arr);
629 Float out;
630 IPosition ip(arr.shape());
631 ip(0) = beamSel_;ip(1) = IFSel_;ip(2) = polSel_;ip(3)=0;
632 out = arr(ip);
633 return out;
634}
635
636MDirection SDMemTable::getDirection(Int whichRow, Bool refBeam) const
637{
638 MDirection::Types mdr = getDirectionReference();
639 Array<Double> posit;
640 dirCol_.get(whichRow,posit);
641 Vector<Double> wpos(2);
642 Int rb;
643 rbeamCol_.get(whichRow,rb);
644 wpos[0] = posit(IPosition(2,beamSel_,0));
645 wpos[1] = posit(IPosition(2,beamSel_,1));
646 if (refBeam && rb > -1) { // use refBeam instead if it exists
647 wpos[0] = posit(IPosition(2,rb,0));
648 wpos[1] = posit(IPosition(2,rb,1));
649 }
650
651 Quantum<Double> lon(wpos[0],Unit(String("rad")));
652 Quantum<Double> lat(wpos[1],Unit(String("rad")));
653 return MDirection(lon, lat, mdr);
654}
655
656MEpoch SDMemTable::getEpoch(Int whichRow) const
657{
658 MEpoch::Types met = getTimeReference();
659//
660 Double obstime;
661 timeCol_.get(whichRow,obstime);
662 MVEpoch tm2(Quantum<Double>(obstime, Unit(String("d"))));
663 return MEpoch(tm2, met);
664}
665
666MPosition SDMemTable::getAntennaPosition () const
667{
668 Vector<Double> antpos;
669 table_.keywordSet().get("AntennaPosition", antpos);
670 MVPosition mvpos(antpos(0),antpos(1),antpos(2));
671 return MPosition(mvpos);
672}
673
674
675SpectralCoordinate SDMemTable::getSpectralCoordinate(uInt freqID) const
676{
677
678 Table t = table_.keywordSet().asTable("FREQUENCIES");
679 if (freqID> t.nrow() ) {
680 throw(AipsError("SDMemTable::getSpectralCoordinate - freqID out of range"));
681 }
682
683 Double rp,rv,inc;
684 String rf;
685 ROScalarColumn<Double> rpc(t, "REFPIX");
686 ROScalarColumn<Double> rvc(t, "REFVAL");
687 ROScalarColumn<Double> incc(t, "INCREMENT");
688 t.keywordSet().get("BASEREFFRAME",rf);
689
690// Create SpectralCoordinate (units Hz)
691
692 MFrequency::Types mft;
693 if (!MFrequency::getType(mft, rf)) {
694 cerr << "Frequency type unknown assuming TOPO" << endl;
695 mft = MFrequency::TOPO;
696 }
697 rpc.get(freqID, rp);
698 rvc.get(freqID, rv);
699 incc.get(freqID, inc);
700//
701 SpectralCoordinate spec(mft,rv,inc,rp);
702//
703 return spec;
704}
705
706
707SpectralCoordinate SDMemTable::getSpectralCoordinate(uInt freqID,
708 uInt whichRow) const
709{
710
711// Create basic SC
712
713 SpectralCoordinate spec = getSpectralCoordinate (freqID);
714//
715 Table t = table_.keywordSet().asTable("FREQUENCIES");
716
717// Get rest frequencies from table, one per IF ???
718
719 Vector<Double> vec;
720 t.keywordSet().get("RESTFREQS",vec);
721 if (vec.nelements() > 0) {
722 spec.setRestFrequencies(vec);
723
724// Select rest freq
725
726 if (IFSel_ < vec.nelements()) {
727 spec.selectRestFrequency(uInt(IFSel_));
728 } else {
729 cerr << "There is no rest frequency for this IF; selecting rest freq for IF=0" << endl;
730 spec.selectRestFrequency(0u);
731 }
732 }
733
734// Set up frame conversion layer
735
736 String frm;
737 t.keywordSet().get("REFFRAME",frm);
738 if (frm == "") frm = "TOPO";
739 MFrequency::Types mtype;
740 if (!MFrequency::getType(mtype, frm)) {
741 cerr << "Frequency type unknown assuming TOPO" << endl; // SHould never happen
742 mtype = MFrequency::TOPO;
743 }
744
745// Set reference frame conversion (requires row)
746
747 MDirection direct = getDirection(whichRow);
748 MEpoch epoch = getEpoch(whichRow);
749 MPosition pos = getAntennaPosition();
750 if (!spec.setReferenceConversion(mtype,epoch,pos,direct)) {
751 throw(AipsError("Couldn't convert frequency frame."));
752 }
753
754// Now velocity conversion if appropriate
755
756 String unitStr;
757 t.keywordSet().get("UNIT",unitStr);
758//
759 String dpl;
760 t.keywordSet().get("DOPPLER",dpl);
761 MDoppler::Types dtype;
762 MDoppler::getType(dtype, dpl);
763
764// Only set velocity unit if non-blank and non-Hz
765
766 if (!unitStr.empty()) {
767 Unit unitU(unitStr);
768 if (unitU==Unit("Hz")) {
769 } else {
770 spec.setVelocity(unitStr, dtype);
771 }
772 }
773//
774 return spec;
775}
776
777
778Bool SDMemTable::setCoordinate(const SpectralCoordinate& speccord,
779 uInt freqID) {
780 Table t = table_.rwKeywordSet().asTable("FREQUENCIES");
781 if (freqID > t.nrow() ) {
782 throw(AipsError("SDMemTable::setCoordinate - coord no out of range"));
783 }
784 ScalarColumn<Double> rpc(t, "REFPIX");
785 ScalarColumn<Double> rvc(t, "REFVAL");
786 ScalarColumn<Double> incc(t, "INCREMENT");
787
788 rpc.put(freqID, speccord.referencePixel()[0]);
789 rvc.put(freqID, speccord.referenceValue()[0]);
790 incc.put(freqID, speccord.increment()[0]);
791
792 return True;
793}
794
795Int SDMemTable::nCoordinates() const
796{
797 return table_.keywordSet().asTable("FREQUENCIES").nrow();
798}
799
800void SDMemTable::setRestFreqs(std::vector<double> freqs,
801 const std::string& theunit)
802{
803 Vector<Double> tvec(freqs);
804 Quantum<Vector<Double> > q(tvec, String(theunit));
805 tvec.resize();
806 tvec = q.getValue("Hz");
807 Table t = table_.keywordSet().asTable("FREQUENCIES");
808 t.rwKeywordSet().define("RESTFREQS",tvec);
809}
810
811std::vector<double> SDMemTable::getRestFreqs() const
812{
813 Table t = table_.keywordSet().asTable("FREQUENCIES");
814 Vector<Double> tvec;
815 t.keywordSet().get("RESTFREQS",tvec);
816 std::vector<double> stlout;
817 tvec.tovector(stlout);
818 return stlout;
819}
820
821bool SDMemTable::putSDFreqTable(const SDFrequencyTable& sdft)
822{
823 TableDesc td("", "1", TableDesc::Scratch);
824 td.addColumn(ScalarColumnDesc<Double>("REFPIX"));
825 td.addColumn(ScalarColumnDesc<Double>("REFVAL"));
826 td.addColumn(ScalarColumnDesc<Double>("INCREMENT"));
827 SetupNewTable aNewTab("freqs", td, Table::New);
828 Table aTable (aNewTab, Table::Memory, sdft.length());
829 ScalarColumn<Double> sc0(aTable, "REFPIX");
830 ScalarColumn<Double> sc1(aTable, "REFVAL");
831 ScalarColumn<Double> sc2(aTable, "INCREMENT");
832 for (uInt i=0; i < sdft.length(); ++i) {
833 sc0.put(i,sdft.referencePixel(i));
834 sc1.put(i,sdft.referenceValue(i));
835 sc2.put(i,sdft.increment(i));
836 }
837 String rf = sdft.refFrame();
838 if (rf.contains("TOPO")) rf = "TOPO";
839
840 aTable.rwKeywordSet().define("BASEREFFRAME", rf);
841 aTable.rwKeywordSet().define("REFFRAME", rf);
842 aTable.rwKeywordSet().define("EQUINOX", sdft.equinox());
843 aTable.rwKeywordSet().define("UNIT", String(""));
844 aTable.rwKeywordSet().define("DOPPLER", String("RADIO"));
845 Vector<Double> rfvec;
846 String rfunit;
847 sdft.restFrequencies(rfvec,rfunit);
848 Quantum<Vector<Double> > q(rfvec, rfunit);
849 rfvec.resize();
850 rfvec = q.getValue("Hz");
851 aTable.rwKeywordSet().define("RESTFREQS", rfvec);
852 table_.rwKeywordSet().defineTable ("FREQUENCIES", aTable);
853 return True;
854}
855
856SDFrequencyTable SDMemTable::getSDFreqTable() const
857{
858 const Table& t = table_.keywordSet().asTable("FREQUENCIES");
859 SDFrequencyTable sdft;
860
861// Add refpix/refval/incr. What are the units ? Hz I suppose
862// but it's nowhere described...
863
864 Vector<Double> refPix, refVal, incr;
865 ScalarColumn<Double> refPixCol(t, "REFPIX");
866 ScalarColumn<Double> refValCol(t, "REFVAL");
867 ScalarColumn<Double> incrCol(t, "INCREMENT");
868 refPix = refPixCol.getColumn();
869 refVal = refValCol.getColumn();
870 incr = incrCol.getColumn();
871//
872 uInt n = refPix.nelements();
873 for (uInt i=0; i<n; i++) {
874 sdft.addFrequency(refPix[i], refVal[i], incr[i]);
875 }
876
877// Frequency reference frame. I don't know if this
878// is the correct frame. It might be 'REFFRAME'
879// rather than 'BASEREFFRAME' ?
880
881 String baseFrame;
882 t.keywordSet().get("BASEREFFRAME",baseFrame);
883 sdft.setRefFrame(baseFrame);
884
885// Equinox
886
887 Float equinox;
888 t.keywordSet().get("EQUINOX", equinox);
889 sdft.setEquinox(equinox);
890
891// Rest Frequency
892
893 Vector<Double> restFreqs;
894 t.keywordSet().get("RESTFREQS", restFreqs);
895 for (uInt i=0; i<restFreqs.nelements(); i++) {
896 sdft.addRestFrequency(restFreqs[i]);
897 }
898 sdft.setRestFrequencyUnit(String("Hz"));
899//
900 return sdft;
901}
902
903bool SDMemTable::putSDContainer(const SDContainer& sdc)
904{
905 uInt rno = table_.nrow();
906 table_.addRow();
907
908// mjd.put(rno, sdc.timestamp);
909 timeCol_.put(rno, sdc.timestamp);
910 srcnCol_.put(rno, sdc.sourcename);
911 fldnCol_.put(rno, sdc.fieldname);
912 specCol_.put(rno, sdc.getSpectrum());
913 flagsCol_.put(rno, sdc.getFlags());
914 tsCol_.put(rno, sdc.getTsys());
915 scanCol_.put(rno, sdc.scanid);
916 integrCol_.put(rno, sdc.interval);
917 freqidCol_.put(rno, sdc.getFreqMap());
918 dirCol_.put(rno, sdc.getDirection());
919 rbeamCol_.put(rno, sdc.refbeam);
920 tcalCol_.put(rno, sdc.tcal);
921 tcaltCol_.put(rno, sdc.tcaltime);
922 azCol_.put(rno, sdc.azimuth);
923 elCol_.put(rno, sdc.elevation);
924 paraCol_.put(rno, sdc.parangle);
925 histCol_.put(rno, sdc.getHistory());
926
927 return true;
928}
929
930SDContainer SDMemTable::getSDContainer(uInt whichRow) const
931{
932 SDContainer sdc(nBeam(),nIF(),nPol(),nChan());
933 timeCol_.get(whichRow, sdc.timestamp);
934 srcnCol_.get(whichRow, sdc.sourcename);
935 integrCol_.get(whichRow, sdc.interval);
936 scanCol_.get(whichRow, sdc.scanid);
937 fldnCol_.get(whichRow, sdc.fieldname);
938 rbeamCol_.get(whichRow, sdc.refbeam);
939 azCol_.get(whichRow, sdc.azimuth);
940 elCol_.get(whichRow, sdc.elevation);
941 paraCol_.get(whichRow, sdc.parangle);
942 Vector<Float> tc;
943 tcalCol_.get(whichRow, tc);
944 sdc.tcal[0] = tc[0];sdc.tcal[1] = tc[1];
945 tcaltCol_.get(whichRow, sdc.tcaltime);
946 Array<Float> spectrum;
947 Array<Float> tsys;
948 Array<uChar> flagtrum;
949 Vector<uInt> fmap;
950 Array<Double> direction;
951 Vector<String> histo;
952 specCol_.get(whichRow, spectrum);
953 sdc.putSpectrum(spectrum);
954 flagsCol_.get(whichRow, flagtrum);
955 sdc.putFlags(flagtrum);
956 tsCol_.get(whichRow, tsys);
957 sdc.putTsys(tsys);
958 freqidCol_.get(whichRow, fmap);
959 sdc.putFreqMap(fmap);
960 dirCol_.get(whichRow, direction);
961 sdc.putDirection(direction);
962 histCol_.get(whichRow, histo);
963 sdc.putHistory(histo);
964 return sdc;
965}
966
967bool SDMemTable::putSDHeader(const SDHeader& sdh)
968{
969 table_.rwKeywordSet().define("nIF", sdh.nif);
970 table_.rwKeywordSet().define("nBeam", sdh.nbeam);
971 table_.rwKeywordSet().define("nPol", sdh.npol);
972 table_.rwKeywordSet().define("nChan", sdh.nchan);
973 table_.rwKeywordSet().define("Observer", sdh.observer);
974 table_.rwKeywordSet().define("Project", sdh.project);
975 table_.rwKeywordSet().define("Obstype", sdh.obstype);
976 table_.rwKeywordSet().define("AntennaName", sdh.antennaname);
977 table_.rwKeywordSet().define("AntennaPosition", sdh.antennaposition);
978 table_.rwKeywordSet().define("Equinox", sdh.equinox);
979 table_.rwKeywordSet().define("FreqRefFrame", sdh.freqref);
980 table_.rwKeywordSet().define("FreqRefVal", sdh.reffreq);
981 table_.rwKeywordSet().define("Bandwidth", sdh.bandwidth);
982 table_.rwKeywordSet().define("UTC", sdh.utc);
983 table_.rwKeywordSet().define("FluxUnit", sdh.fluxunit);
984 table_.rwKeywordSet().define("Epoch", sdh.epoch);
985 return true;
986}
987
988SDHeader SDMemTable::getSDHeader() const
989{
990 SDHeader sdh;
991 table_.keywordSet().get("nBeam",sdh.nbeam);
992 table_.keywordSet().get("nIF",sdh.nif);
993 table_.keywordSet().get("nPol",sdh.npol);
994 table_.keywordSet().get("nChan",sdh.nchan);
995 table_.keywordSet().get("Observer", sdh.observer);
996 table_.keywordSet().get("Project", sdh.project);
997 table_.keywordSet().get("Obstype", sdh.obstype);
998 table_.keywordSet().get("AntennaName", sdh.antennaname);
999 table_.keywordSet().get("AntennaPosition", sdh.antennaposition);
1000 table_.keywordSet().get("Equinox", sdh.equinox);
1001 table_.keywordSet().get("FreqRefFrame", sdh.freqref);
1002 table_.keywordSet().get("FreqRefVal", sdh.reffreq);
1003 table_.keywordSet().get("Bandwidth", sdh.bandwidth);
1004 table_.keywordSet().get("UTC", sdh.utc);
1005 table_.keywordSet().get("FluxUnit", sdh.fluxunit);
1006 table_.keywordSet().get("Epoch", sdh.epoch);
1007 return sdh;
1008}
1009void SDMemTable::makePersistent(const std::string& filename)
1010{
1011 table_.deepCopy(filename,Table::New);
1012}
1013
1014Int SDMemTable::nScan() const {
1015 Int n = 0;
1016 Int previous = -1;Int current=0;
1017 for (uInt i=0; i< scanCol_.nrow();i++) {
1018 scanCol_.getScalar(i,current);
1019 if (previous != current) {
1020 previous = current;
1021 n++;
1022 }
1023 }
1024 return n;
1025}
1026
1027String SDMemTable::formatSec(Double x) const
1028{
1029 Double xcop = x;
1030 MVTime mvt(xcop/24./3600.); // make days
1031
1032 if (x < 59.95)
1033 return String(" ") + mvt.string(MVTime::TIME_CLEAN_NO_HM, 7)+"s";
1034 else if (x < 3599.95)
1035 return String(" ") + mvt.string(MVTime::TIME_CLEAN_NO_H,7)+" ";
1036 else {
1037 ostringstream oss;
1038 oss << setw(2) << std::right << setprecision(1) << mvt.hour();
1039 oss << ":" << mvt.string(MVTime::TIME_CLEAN_NO_H,7) << " ";
1040 return String(oss);
1041 }
1042};
1043
1044String SDMemTable::formatDirection(const MDirection& md) const
1045{
1046 Vector<Double> t = md.getAngle(Unit(String("rad"))).getValue();
1047 Int prec = 7;
1048
1049 MVAngle mvLon(t[0]);
1050 String sLon = mvLon.string(MVAngle::TIME,prec);
1051 MVAngle mvLat(t[1]);
1052 String sLat = mvLat.string(MVAngle::ANGLE+MVAngle::DIG2,prec);
1053 return sLon + String(" ") + sLat;
1054}
1055
1056
1057std::string SDMemTable::getFluxUnit() const
1058{
1059 String tmp;
1060 table_.keywordSet().get("FluxUnit", tmp);
1061 return tmp;
1062}
1063
1064void SDMemTable::setFluxUnit(const std::string& unit)
1065{
1066 String tmp(unit);
1067 Unit tU(tmp);
1068 if (tU==Unit("K") || tU==Unit("Jy")) {
1069 table_.rwKeywordSet().define(String("FluxUnit"), tmp);
1070 } else {
1071 throw AipsError("Illegal unit - must be compatible with Jy or K");
1072 }
1073}
1074
1075
1076void SDMemTable::setInstrument(const std::string& name)
1077{
1078 Bool throwIt = True;
1079 Instrument ins = convertInstrument (name, throwIt);
1080 String nameU(name);
1081 nameU.upcase();
1082 table_.rwKeywordSet().define(String("AntennaName"), nameU);
1083}
1084
1085std::string SDMemTable::summary(bool verbose) const {
1086
1087// Format header info
1088
1089 ostringstream oss;
1090 oss << endl;
1091 oss << "--------------------------------------------------------------------------------" << endl;
1092 oss << " Scan Table Summary" << endl;
1093 oss << "--------------------------------------------------------------------------------" << endl;
1094 oss.flags(std::ios_base::left);
1095 oss << setw(15) << "Beams:" << setw(4) << nBeam() << endl
1096 << setw(15) << "IFs:" << setw(4) << nIF() << endl
1097 << setw(15) << "Polarisations:" << setw(4) << nPol() << endl
1098 << setw(15) << "Channels:" << setw(4) << nChan() << endl;
1099 oss << endl;
1100 String tmp;
1101 table_.keywordSet().get("Observer", tmp);
1102 oss << setw(15) << "Observer:" << tmp << endl;
1103 oss << setw(15) << "Obs Date:" << getTime(-1,True) << endl;
1104 table_.keywordSet().get("Project", tmp);
1105 oss << setw(15) << "Project:" << tmp << endl;
1106 table_.keywordSet().get("Obstype", tmp);
1107 oss << setw(15) << "Obs. Type:" << tmp << endl;
1108 table_.keywordSet().get("AntennaName", tmp);
1109 oss << setw(15) << "Antenna Name:" << tmp << endl;
1110 table_.keywordSet().get("FluxUnit", tmp);
1111 oss << setw(15) << "Flux Unit:" << tmp << endl;
1112 Table t = table_.keywordSet().asTable("FREQUENCIES");
1113 Vector<Double> vec;
1114 t.keywordSet().get("RESTFREQS",vec);
1115 oss << setw(15) << "Rest Freqs:";
1116 if (vec.nelements() > 0) {
1117 oss << setprecision(0) << vec << " [Hz]" << endl;
1118 } else {
1119 oss << "None set" << endl;
1120 }
1121 oss << setw(15) << "Abcissa:" << getAbcissaString() << endl;
1122 oss << setw(15) << "Cursor:" << "Beam[" << getBeam() << "] "
1123 << "IF[" << getIF() << "] " << "Pol[" << getPol() << "]" << endl;
1124 oss << endl;
1125//
1126 String dirtype ="Position ("+ MDirection::showType(getDirectionReference()) + ")";
1127 oss << setw(5) << "Scan"
1128 << setw(15) << "Source"
1129 << setw(24) << dirtype
1130 << setw(10) << "Time"
1131 << setw(18) << "Integration"
1132 << setw(7) << "FreqIDs" << endl;
1133 oss << "--------------------------------------------------------------------------------" << endl;
1134
1135// Generate list of scan start and end integrations
1136
1137 Vector<Int> scanIDs = scanCol_.getColumn();
1138 Vector<uInt> startInt, endInt;
1139 mathutil::scanBoundaries(startInt, endInt, scanIDs);
1140//
1141 const uInt nScans = startInt.nelements();
1142 String name;
1143 Vector<uInt> freqIDs, listFQ;
1144 uInt nInt;
1145//
1146 for (uInt i=0; i<nScans; i++) {
1147
1148// Get things from first integration of scan
1149
1150 String time = getTime(startInt(i),False);
1151 String tInt = formatSec(Double(getInterval(startInt(i))));
1152 String posit = formatDirection(getDirection(startInt(i),True));
1153 srcnCol_.getScalar(startInt(i),name);
1154
1155// Find all the FreqIDs in this scan
1156
1157 listFQ.resize(0);
1158 for (uInt j=startInt(i); j<endInt(i)+1; j++) {
1159 freqidCol_.get(j, freqIDs);
1160 for (uInt k=0; k<freqIDs.nelements(); k++) {
1161 mathutil::addEntry(listFQ, freqIDs(k));
1162 }
1163 }
1164//
1165 nInt = endInt(i) - startInt(i) + 1;
1166 oss << setw(3) << std::right << i << std::left << setw(2) << " "
1167 << setw(15) << name
1168 << setw(24) << posit
1169 << setw(10) << time
1170 << setw(3) << std::right << nInt << setw(3) << " x " << std::left
1171 << setw(6) << tInt
1172 << " " << listFQ << endl;
1173 }
1174 oss << endl;
1175 oss << "Table contains " << table_.nrow() << " integration(s) in " << nScans << " scan(s)." << endl;
1176
1177// Frequency Table
1178
1179 if (verbose) {
1180 std::vector<string> info = getCoordInfo();
1181 SDFrequencyTable sdft = getSDFreqTable();
1182 oss << endl << endl;
1183 oss << "FreqID Frame RefFreq(Hz) RefPix Increment(Hz)" << endl;
1184 oss << "--------------------------------------------------------------------------------" << endl;
1185 for (uInt i=0; i<sdft.length(); i++) {
1186 oss << setw(8) << i << setw(8)
1187 << info[3] << setw(16) << setprecision(8)
1188 << sdft.referenceValue(i) << setw(10)
1189 << sdft.referencePixel(i) << setw(12)
1190 << sdft.increment(i) << endl;
1191 }
1192 oss << "--------------------------------------------------------------------------------" << endl;
1193 }
1194 return String(oss);
1195}
1196
1197Int SDMemTable::nBeam() const
1198{
1199 Int n;
1200 table_.keywordSet().get("nBeam",n);
1201 return n;
1202}
1203Int SDMemTable::nIF() const {
1204 Int n;
1205 table_.keywordSet().get("nIF",n);
1206 return n;
1207}
1208Int SDMemTable::nPol() const {
1209 Int n;
1210 table_.keywordSet().get("nPol",n);
1211 return n;
1212}
1213Int SDMemTable::nChan() const {
1214 Int n;
1215 table_.keywordSet().get("nChan",n);
1216 return n;
1217}
1218bool SDMemTable::appendHistory(const std::string& hist, int whichRow)
1219{
1220 Vector<String> history;
1221 histCol_.get(whichRow, history);
1222 history.resize(history.nelements()+1,True);
1223 history[history.nelements()-1] = hist;
1224 histCol_.put(whichRow, history);
1225}
1226
1227std::vector<std::string> SDMemTable::history(int whichRow) const
1228{
1229 Vector<String> history;
1230 histCol_.get(whichRow, history);
1231 std::vector<std::string> stlout;
1232 // there is no Array<String>.tovector(std::vector<std::string>), so
1233 // do it by hand
1234 for (uInt i=0; i<history.nelements(); ++i) {
1235 stlout.push_back(history[i]);
1236 }
1237 return stlout;
1238}
1239/*
1240void SDMemTable::maskChannels(const std::vector<Int>& whichChans ) {
1241
1242 std::vector<int>::iterator it;
1243 ArrayAccessor<uChar, Axis<asap::PolAxis> > j(flags_);
1244 for (it = whichChans.begin(); it != whichChans.end(); it++) {
1245 j.reset(j.begin(uInt(*it)));
1246 for (ArrayAccessor<uChar, Axis<asap::BeamAxis> > i(j); i != i.end(); ++i) {
1247 for (ArrayAccessor<uChar, Axis<asap::IFAxis> > ii(i); ii != ii.end(); ++ii) {
1248 for (ArrayAccessor<uChar, Axis<asap::ChanAxis> > iii(ii);
1249 iii != iii.end(); ++iii) {
1250 (*iii) =
1251 }
1252 }
1253 }
1254 }
1255
1256}
1257*/
1258void SDMemTable::flag(int whichRow)
1259{
1260 Array<uChar> arr;
1261 flagsCol_.get(whichRow, arr);
1262
1263 ArrayAccessor<uChar, Axis<asap::BeamAxis> > aa0(arr);
1264 aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
1265 ArrayAccessor<uChar, Axis<asap::IFAxis> > aa1(aa0);
1266 aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
1267 ArrayAccessor<uChar, Axis<asap::PolAxis> > aa2(aa1);
1268 aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
1269
1270 for (ArrayAccessor<uChar, Axis<asap::ChanAxis> > i(aa2); i != i.end(); ++i) {
1271 (*i) = uChar(True);
1272 }
1273
1274 flagsCol_.put(whichRow, arr);
1275}
1276
1277MDirection::Types SDMemTable::getDirectionReference() const
1278{
1279 Float eq;
1280 table_.keywordSet().get("Equinox",eq);
1281 std::map<float,string> mp;
1282 mp[2000.0] = "J2000";
1283 mp[1950.0] = "B1950";
1284 MDirection::Types mdr;
1285 if (!MDirection::getType(mdr, mp[eq])) {
1286 mdr = MDirection::J2000;
1287 cerr << "Unknown equinox using J2000" << endl;
1288 }
1289//
1290 return mdr;
1291}
1292
1293MEpoch::Types SDMemTable::getTimeReference() const
1294{
1295 MEpoch::Types met;
1296 String ep;
1297 table_.keywordSet().get("Epoch",ep);
1298 if (!MEpoch::getType(met, ep)) {
1299 cerr << "Epoch type uknown - using UTC" << endl;
1300 met = MEpoch::UTC;
1301 }
1302//
1303 return met;
1304}
1305
1306
1307Instrument SDMemTable::convertInstrument(const String& instrument,
1308 Bool throwIt)
1309{
1310 String t(instrument);
1311 t.upcase();
1312
1313// The strings are what SDReader returns, after cunning interrogation
1314// of station names... :-(
1315
1316 Instrument inst = asap::UNKNOWN;
1317 if (t==String("DSS-43")) {
1318 inst = TIDBINBILLA;
1319 } else if (t==String("ATPKSMB")) {
1320 inst = ATPKSMB;
1321 } else if (t==String("ATPKSHOH")) {
1322 inst = ATPKSHOH;
1323 } else if (t==String("ATMOPRA")) {
1324 inst = ATMOPRA;
1325 } else if (t==String("CEDUNA")) {
1326 inst = CEDUNA;
1327 } else if (t==String("HOBART")) {
1328 inst = HOBART;
1329 } else {
1330 if (throwIt) {
1331 throw AipsError("Unrecognized instrument - use function scan.set_instrument to set");
1332 }
1333 }
1334 return inst;
1335}
1336
1337void SDMemTable::renumber()
1338{
1339 uInt nRow = scanCol_.nrow();
1340 Int newscanid = 0;
1341 Int cIdx;// the current scanid
1342 // get the first scanid
1343 scanCol_.getScalar(0,cIdx);
1344 Int pIdx = cIdx;// the scanid of the previous row
1345 for (uInt i=0; i<nRow;++i) {
1346 scanCol_.getScalar(i,cIdx);
1347 if (pIdx == cIdx) {
1348 // renumber
1349 scanCol_.put(i,newscanid);
1350 } else {
1351 ++newscanid;
1352 pIdx = cIdx; // store scanid
1353 --i; // don't increment next loop
1354 }
1355 }
1356}
Note: See TracBrowser for help on using the repository browser.