source: trunk/src/SDAsciiWriter.cc@ 323

Last change on this file since 323 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: 5.7 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 215 2005-01-19 04:03:10Z mar637 $
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)
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));
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 cerr << "Wrote " << nRows << " rows into file " << fileName << endl;
178//
179 return True;
180}
181
182
183Int SDAsciiWriter::convertStokes(Int val)
184{
185 Stokes::StokesTypes stokes = Stokes::RR;
186 if (val==0) {
187 stokes = Stokes::RR;
188 } else if (val==1) {
189 stokes = Stokes::LL;
190 } else if (val==2) {
191 stokes = Stokes::RL;
192 } else if (val==3) {
193 stokes = Stokes::LR;
194 } else {
195 stokes = Stokes::Undefined;
196 }
197//
198 return Int(stokes);
199}
200
201
202String SDAsciiWriter::formatDirection(const Vector<Double>& lonLat)
203{
204 MVAngle x1(lonLat(0));
205 String s1 = x1.string(MVAngle::TIME, 12);
206//
207 MVAngle x2(lonLat(1));
208 String s2 = x2.string(MVAngle::ANGLE, 12);
209//
210 String ss = s1 + String(" ") + s2;
211 return ss;
212}
213
Note: See TracBrowser for help on using the repository browser.