source: trunk/src/SDAsciiWriter.cc@ 474

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

track change to SDmemTable::rowAsMaskedArray interface

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.4 KB
Line 
1//#---------------------------------------------------------------------------
2//# SDAsciiWriter.cc: ASAP class to write out single dish spectra as FITS images
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: SDAsciiWriter.cc 448 2005-02-15 07:39:54Z kil064 $
30//#---------------------------------------------------------------------------
31
32#include <casa/aips.h>
33#include <casa/Arrays/Array.h>
34#include <casa/Arrays/Vector.h>
35#include <casa/Arrays/VectorIter.h>
36#include <casa/Utilities/CountedPtr.h>
37#include <casa/Quanta/Quantum.h>
38#include <casa/Quanta/MVAngle.h>
39
40
41#include <casa/iostream.h>
42#include <casa/fstream.h>
43
44#include <coordinates/Coordinates/CoordinateUtil.h>
45#include <coordinates/Coordinates/SpectralCoordinate.h>
46#include <coordinates/Coordinates/DirectionCoordinate.h>
47#include <coordinates/Coordinates/StokesCoordinate.h>
48
49#include <measures/Measures/MEpoch.h>
50
51#include <tables/Tables/Table.h>
52#include <tables/Tables/ScalarColumn.h>
53#include <tables/Tables/ArrayColumn.h>
54
55#include "SDDefs.h"
56#include "SDContainer.h"
57#include "SDMemTable.h"
58#include "SDAsciiWriter.h"
59
60using namespace casa;
61using namespace asap;
62
63
64SDAsciiWriter::SDAsciiWriter()
65{;}
66
67SDAsciiWriter::~SDAsciiWriter()
68{;}
69
70
71Bool SDAsciiWriter::write(const SDMemTable& sdTable, const String& fileName, Bool toStokes)
72{
73
74// Get global Header from Table
75
76 SDHeader header = sdTable.getSDHeader();
77 MEpoch::Ref timeRef(MEpoch::UTC); // Should be in header
78 MDirection::Types dirRef(MDirection::J2000); // Should be in header
79
80// Column keywords
81
82 Table tab = sdTable.table();
83 ROArrayColumn<Double> dir(tab, String("DIRECTION"));
84 ROScalarColumn<Double> time(tab, "TIME");
85 ROArrayColumn<uInt> freqid(tab, "FREQID");
86 ROScalarColumn<String> src(tab, "SRCNAME");
87
88// Temps
89
90 Vector<Int> whichStokes(1,1);
91 Array<Double> whichDir;
92 Vector<Double> lonLat(2);
93 IPosition posDir(2,0);
94 const Unit RAD(String("rad"));
95
96// Open file
97
98 String fName(fileName);
99 if (fileName.length()==0) fName = String("ascii.txt");
100 ofstream of(fName.chars(), ios::trunc);
101
102// Write header
103
104 of << "row beam IF pol source longitude latitude time nchan spectrum mask"
105 << endl;
106
107// Loop over rows
108
109 const uInt nRows = sdTable.nRow();
110 for (uInt iRow=0; iRow<nRows; iRow++) {
111
112// Get data
113
114 const MaskedArray<Float>& dataIn(sdTable.rowAsMaskedArray(iRow,toStokes));
115 const Array<Float>& values = dataIn.getArray();
116 const Array<Bool>& mask = dataIn.getMask();
117
118// Epoch
119
120 Double dTmp;
121 time.get(iRow, dTmp);
122 MVEpoch tmp2(Quantum<Double>(dTmp, Unit(String("d"))));
123 MEpoch epoch(tmp2, timeRef);
124
125// Iterate through data in this row by spectra
126
127 ReadOnlyVectorIterator<Float> itData(values, asap::ChanAxis);
128 ReadOnlyVectorIterator<Bool> itMask(mask, asap::ChanAxis);
129 while (!itData.pastEnd()) {
130 const IPosition& pos = itData.pos();
131
132// FreqID
133
134 Vector<uInt> iTmp;
135 freqid.get(iRow, iTmp);
136
137// Direction
138
139 dir.get(iRow, whichDir);
140 posDir(0) = pos(asap::BeamAxis);
141 posDir(1) = 0;
142 lonLat[0] = whichDir(posDir);
143//
144 posDir(0) = pos(asap::BeamAxis);
145 posDir(1) = 1;
146 lonLat[1] = whichDir(posDir);
147
148// Write preamble
149
150 of << iRow << " " << pos(asap::BeamAxis) << " " << pos(asap::IFAxis) << " " <<
151 pos(asap::PolAxis) << " " <<
152 src(iRow) << " " << formatDirection(lonLat) << " " <<
153 dTmp << " " << itData.vector().nelements() << " ";
154// Write data
155
156 const Vector<Float>& data = itData.vector();
157 const Vector<Bool>& mask = itMask.vector();
158 const uInt n = data.nelements();
159 for (uInt i=0; i<n; i++) {
160 of << data[i] << " ";
161 }
162// Write mask
163
164 for (uInt i=0; i<n; i++) {
165 of << mask[i] << " ";
166 }
167 of << endl;
168
169// Next spectrum
170
171 itData.next();
172 itMask.next();
173 }
174 }
175//
176 of.close();
177 cout << "Wrote " << nRows << " rows into file " << fileName << endl;
178//
179 return True;
180}
181
182
183
184String SDAsciiWriter::formatDirection(const Vector<Double>& lonLat)
185{
186 MVAngle x1(lonLat(0));
187 String s1 = x1.string(MVAngle::TIME, 12);
188//
189 MVAngle x2(lonLat(1));
190 String s2 = x2.string(MVAngle::ANGLE, 12);
191//
192 String ss = s1 + String(" ") + s2;
193 return ss;
194}
195
Note: See TracBrowser for help on using the repository browser.