1 | //#---------------------------------------------------------------------------
|
---|
2 | //# PKSSDwriter.cc: Class to write Parkes multibeam data to an SDFITS file.
|
---|
3 | //#---------------------------------------------------------------------------
|
---|
4 | //# Copyright (C) 2000-2006
|
---|
5 | //# Associated Universities, Inc. Washington DC, USA.
|
---|
6 | //#
|
---|
7 | //# This library is free software; you can redistribute it and/or modify it
|
---|
8 | //# under the terms of the GNU Library General Public License as published by
|
---|
9 | //# the Free Software Foundation; either version 2 of the License, or (at your
|
---|
10 | //# option) any later version.
|
---|
11 | //#
|
---|
12 | //# This library is distributed in the hope that it will be useful, but
|
---|
13 | //# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
|
---|
14 | //# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
|
---|
15 | //# License for more details.
|
---|
16 | //#
|
---|
17 | //# You should have received a copy of the GNU Library General Public License
|
---|
18 | //# along with this library; if not, write to the Free Software Foundation,
|
---|
19 | //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
|
---|
20 | //#
|
---|
21 | //# Correspondence concerning AIPS++ should be addressed as follows:
|
---|
22 | //# Internet email: aips2-request@nrao.edu.
|
---|
23 | //# Postal address: AIPS++ Project Office
|
---|
24 | //# National Radio Astronomy Observatory
|
---|
25 | //# 520 Edgemont Road
|
---|
26 | //# Charlottesville, VA 22903-2475 USA
|
---|
27 | //#
|
---|
28 | //# $Id: PKSSDwriter.cc,v 19.11 2006/07/05 05:41:12 mcalabre Exp $
|
---|
29 | //#---------------------------------------------------------------------------
|
---|
30 |
|
---|
31 | #include <atnf/PKSIO/PKSMBrecord.h>
|
---|
32 | #include <atnf/PKSIO/PKSSDwriter.h>
|
---|
33 |
|
---|
34 | #include <casa/Quanta/MVTime.h>
|
---|
35 |
|
---|
36 |
|
---|
37 | //--------------------------------------------------- PKSSDwriter::PKSSDwriter
|
---|
38 |
|
---|
39 | // Default constructor.
|
---|
40 |
|
---|
41 | PKSSDwriter::PKSSDwriter()
|
---|
42 | {
|
---|
43 | }
|
---|
44 |
|
---|
45 | //-------------------------------------------------- PKSSDwriter::~PKSSDwriter
|
---|
46 |
|
---|
47 | // Destructor.
|
---|
48 |
|
---|
49 | PKSSDwriter::~PKSSDwriter()
|
---|
50 | {
|
---|
51 | close();
|
---|
52 | }
|
---|
53 |
|
---|
54 | //-------------------------------------------------------- PKSSDwriter::create
|
---|
55 |
|
---|
56 | // Create the SDFITS file and and write static data.
|
---|
57 |
|
---|
58 | Int PKSSDwriter::create(
|
---|
59 | const String sdName,
|
---|
60 | const String observer,
|
---|
61 | const String project,
|
---|
62 | const String antName,
|
---|
63 | const Vector<Double> antPosition,
|
---|
64 | const String obsMode,
|
---|
65 | const Float equinox,
|
---|
66 | const String dopplerFrame,
|
---|
67 | const Vector<uInt> nChan,
|
---|
68 | const Vector<uInt> nPol,
|
---|
69 | const Vector<Bool> haveXPol,
|
---|
70 | const Bool haveBase)
|
---|
71 | {
|
---|
72 | double antPos[3];
|
---|
73 | antPos[0] = antPosition(0);
|
---|
74 | antPos[1] = antPosition(1);
|
---|
75 | antPos[2] = antPosition(2);
|
---|
76 |
|
---|
77 | cNIF = nChan.nelements();
|
---|
78 | if (nPol.nelements() != cNIF || haveXPol.nelements() != cNIF) {
|
---|
79 | cerr << "PKSSDwriter::create: "
|
---|
80 | << "Inconsistent number of IFs for nChan, nPol, and/or haveXPol."
|
---|
81 | << endl;
|
---|
82 | return 1;
|
---|
83 | }
|
---|
84 |
|
---|
85 | cNChan.assign(nChan);
|
---|
86 | cNPol.assign(nPol);
|
---|
87 |
|
---|
88 | cHaveXPol.resize(cNIF);
|
---|
89 | for (uInt iIF = 0; iIF < cNIF; iIF++) {
|
---|
90 | // Convert Bool -> uInt.
|
---|
91 | cHaveXPol(iIF) = haveXPol(iIF) ? 1 : 0;
|
---|
92 | }
|
---|
93 |
|
---|
94 | cHaveBase = haveBase;
|
---|
95 |
|
---|
96 | // Storage in the trivial cNChan, cNPol, and cHaveXPol arrays should always
|
---|
97 | // be contiguous so the pointer returned by getStorage() shouldn't need to
|
---|
98 | // be deleted via freeStorage() (i.e. deleteIt always returned False). This
|
---|
99 | // storage will, of course, be deleted when the PKSwriter object is deleted.
|
---|
100 | Bool deleteIt;
|
---|
101 | Int status = cSDwriter.create((char *)sdName.chars(),
|
---|
102 | (char *)observer.chars(), (char *)project.chars(),
|
---|
103 | (char *)antName.chars(), antPos, (char *)obsMode.chars(), equinox,
|
---|
104 | (char *)dopplerFrame.chars(), cNIF,
|
---|
105 | (int *)cNChan.getStorage(deleteIt),
|
---|
106 | (int *)cNPol.getStorage(deleteIt),
|
---|
107 | (int *)cHaveXPol.getStorage(deleteIt), (int)cHaveBase, 1);
|
---|
108 | if (status) {
|
---|
109 | cSDwriter.reportError();
|
---|
110 | cSDwriter.deleteFile();
|
---|
111 | close();
|
---|
112 | }
|
---|
113 |
|
---|
114 | return status;
|
---|
115 | }
|
---|
116 |
|
---|
117 | //--------------------------------------------------------- PKSSDwriter::write
|
---|
118 |
|
---|
119 | // Write the next data record.
|
---|
120 |
|
---|
121 | Int PKSSDwriter::write(
|
---|
122 | const Int scanNo,
|
---|
123 | const Int cycleNo,
|
---|
124 | const Double mjd,
|
---|
125 | const Double interval,
|
---|
126 | const String fieldName,
|
---|
127 | const String srcName,
|
---|
128 | const Vector<Double> srcDir,
|
---|
129 | const Vector<Double> srcPM,
|
---|
130 | const Double srcVel,
|
---|
131 | const String obsMode,
|
---|
132 | const Int IFno,
|
---|
133 | const Double refFreq,
|
---|
134 | const Double bandwidth,
|
---|
135 | const Double freqInc,
|
---|
136 | const Double restFreq,
|
---|
137 | const Vector<Float> tcal,
|
---|
138 | const String tcalTime,
|
---|
139 | const Float azimuth,
|
---|
140 | const Float elevation,
|
---|
141 | const Float parAngle,
|
---|
142 | const Float focusAxi,
|
---|
143 | const Float focusTan,
|
---|
144 | const Float focusRot,
|
---|
145 | const Float temperature,
|
---|
146 | const Float pressure,
|
---|
147 | const Float humidity,
|
---|
148 | const Float windSpeed,
|
---|
149 | const Float windAz,
|
---|
150 | const Int refBeam,
|
---|
151 | const Int beamNo,
|
---|
152 | const Vector<Double> direction,
|
---|
153 | const Vector<Double> scanRate,
|
---|
154 | const Vector<Float> tsys,
|
---|
155 | const Vector<Float> sigma,
|
---|
156 | const Vector<Float> calFctr,
|
---|
157 | const Matrix<Float> baseLin,
|
---|
158 | const Matrix<Float> baseSub,
|
---|
159 | const Matrix<Float> &spectra,
|
---|
160 | const Matrix<uChar> &flagged,
|
---|
161 | const Complex xCalFctr,
|
---|
162 | const Vector<Complex> &xPol)
|
---|
163 | {
|
---|
164 | // Do basic checks.
|
---|
165 | uInt iIF = IFno - 1;
|
---|
166 | if (IFno < 1 || Int(cNIF) < IFno) {
|
---|
167 | cerr << "PKSDwriter::write: "
|
---|
168 | << "Invalid IF number " << IFno
|
---|
169 | << " (maximum " << cNIF << ")." << endl;
|
---|
170 | return 1;
|
---|
171 | }
|
---|
172 |
|
---|
173 | uInt nChan = spectra.nrow();
|
---|
174 | if (nChan != cNChan(iIF)) {
|
---|
175 | cerr << "PKSDwriter::write: "
|
---|
176 | << "Wrong number of channels for IF " << IFno << "," << endl
|
---|
177 | << " "
|
---|
178 | << "got " << nChan << " should be " << cNChan(iIF) << "." << endl;
|
---|
179 | return 1;
|
---|
180 | }
|
---|
181 |
|
---|
182 | uInt nPol = spectra.ncolumn();
|
---|
183 | if (nPol != cNPol(iIF)) {
|
---|
184 | cerr << "PKSDwriter::write: "
|
---|
185 | << "Wrong number of polarizations for IF " << IFno << "," << endl
|
---|
186 | << " "
|
---|
187 | << "got " << nPol << " should be " << cNPol(iIF) << "." << endl;
|
---|
188 | return 1;
|
---|
189 | }
|
---|
190 |
|
---|
191 | // Extract calendar information from mjd.
|
---|
192 | MVTime time(mjd);
|
---|
193 | Int year = time.year();
|
---|
194 | Int month = time.month();
|
---|
195 | Int day = time.monthday();
|
---|
196 |
|
---|
197 | // Transfer data to a single-IF PKSMBrecord.
|
---|
198 | PKSMBrecord mbrec(1);
|
---|
199 |
|
---|
200 | // Start with basic beam- and IF-independent bookkeeping information.
|
---|
201 | mbrec.scanNo = scanNo;
|
---|
202 | mbrec.cycleNo = cycleNo;
|
---|
203 |
|
---|
204 | sprintf(mbrec.datobs, "%4.4d-%2.2d-%2.2d", year, month, day);
|
---|
205 | mbrec.utc = fmod(mjd, 1.0) * 86400.0;
|
---|
206 | mbrec.exposure = float(interval);
|
---|
207 |
|
---|
208 | strncpy(mbrec.srcName, (char *)srcName.chars(), 17);
|
---|
209 | mbrec.srcRA = srcDir(0);
|
---|
210 | mbrec.srcDec = srcDir(1);
|
---|
211 |
|
---|
212 | mbrec.restFreq = restFreq;
|
---|
213 |
|
---|
214 | strncpy(mbrec.obsType, (char *)obsMode.chars(), 16);
|
---|
215 |
|
---|
216 | // Now beam-dependent parameters.
|
---|
217 | mbrec.beamNo = beamNo;
|
---|
218 | mbrec.ra = direction(0);
|
---|
219 | mbrec.dec = direction(1);
|
---|
220 | mbrec.raRate = scanRate(0);
|
---|
221 | mbrec.decRate = scanRate(1);
|
---|
222 |
|
---|
223 | // Now IF-dependent parameters.
|
---|
224 | mbrec.nIF = 1;
|
---|
225 | mbrec.IFno[0] = IFno;
|
---|
226 | mbrec.nChan[0] = nChan;
|
---|
227 | mbrec.nPol[0] = nPol;
|
---|
228 |
|
---|
229 | mbrec.fqRefPix[0] = (nChan/2) + 1;
|
---|
230 | mbrec.fqRefVal[0] = refFreq;
|
---|
231 | mbrec.fqDelt[0] = freqInc;
|
---|
232 |
|
---|
233 | // Now the data itself.
|
---|
234 | for (uInt i = 0; i < tsys.nelements(); i++) {
|
---|
235 | mbrec.tsys[0][i] = tsys(i);
|
---|
236 | }
|
---|
237 |
|
---|
238 | for (uInt ipol = 0; ipol < nPol; ipol++) {
|
---|
239 | mbrec.calfctr[0][ipol] = calFctr(ipol);
|
---|
240 | }
|
---|
241 |
|
---|
242 | if (cHaveXPol(iIF)) {
|
---|
243 | mbrec.xcalfctr[0][0] = xCalFctr.real();
|
---|
244 | mbrec.xcalfctr[0][1] = xCalFctr.imag();
|
---|
245 | } else {
|
---|
246 | mbrec.xcalfctr[0][0] = 0.0f;
|
---|
247 | mbrec.xcalfctr[0][1] = 0.0f;
|
---|
248 | }
|
---|
249 |
|
---|
250 | if (cHaveBase) {
|
---|
251 | mbrec.haveBase = 1;
|
---|
252 |
|
---|
253 | for (uInt ipol = 0; ipol < nPol; ipol++) {
|
---|
254 | mbrec.baseLin[0][ipol][0] = baseLin(0,ipol);
|
---|
255 | mbrec.baseLin[0][ipol][1] = baseLin(1,ipol);
|
---|
256 |
|
---|
257 | for (uInt j = 0; j < baseSub.nrow(); j++) {
|
---|
258 | mbrec.baseSub[0][ipol][j] = baseSub(j,ipol);
|
---|
259 | }
|
---|
260 | for (uInt j = baseSub.nrow(); j < 9; j++) {
|
---|
261 | mbrec.baseSub[0][ipol][j] = 0.0f;
|
---|
262 | }
|
---|
263 | }
|
---|
264 |
|
---|
265 | } else {
|
---|
266 | mbrec.haveBase = 0;
|
---|
267 | }
|
---|
268 |
|
---|
269 | Bool delSpectra = False;
|
---|
270 | const Float *specstor = spectra.getStorage(delSpectra);
|
---|
271 | mbrec.spectra[0] = (float *)specstor;
|
---|
272 |
|
---|
273 | Bool delFlagged = False;
|
---|
274 | const uChar *flagstor = flagged.getStorage(delFlagged);
|
---|
275 | mbrec.flagged[0] = (unsigned char *)flagstor;
|
---|
276 |
|
---|
277 | Bool delXPol = False;
|
---|
278 | const Complex *xpolstor;
|
---|
279 | if (cHaveXPol(iIF)) {
|
---|
280 | xpolstor = xPol.getStorage(delXPol);
|
---|
281 | } else {
|
---|
282 | xpolstor = 0;
|
---|
283 | }
|
---|
284 | mbrec.xpol[0] = (float *)xpolstor;
|
---|
285 |
|
---|
286 | // Finish off with system calibration parameters.
|
---|
287 | mbrec.extraSysCal = 1;
|
---|
288 | mbrec.refBeam = refBeam;
|
---|
289 | for (uInt i = 0; i < tcal.nelements(); i++) {
|
---|
290 | mbrec.tcal[0][i] = tcal(i);
|
---|
291 | }
|
---|
292 | strncpy(mbrec.tcalTime, (char *)tcalTime.chars(), 16);
|
---|
293 | mbrec.azimuth = azimuth;
|
---|
294 | mbrec.elevation = elevation;
|
---|
295 | mbrec.parAngle = parAngle;
|
---|
296 | mbrec.focusAxi = focusAxi;
|
---|
297 | mbrec.focusTan = focusTan;
|
---|
298 | mbrec.focusRot = focusRot;
|
---|
299 | mbrec.temp = temperature;
|
---|
300 | mbrec.pressure = pressure;
|
---|
301 | mbrec.humidity = humidity;
|
---|
302 | mbrec.windSpeed = windSpeed;
|
---|
303 | mbrec.windAz = windAz;
|
---|
304 |
|
---|
305 | Int status = cSDwriter.write(mbrec);
|
---|
306 | if (status) {
|
---|
307 | cSDwriter.reportError();
|
---|
308 | status = 1;
|
---|
309 | }
|
---|
310 |
|
---|
311 | spectra.freeStorage(specstor, delSpectra);
|
---|
312 | flagged.freeStorage(flagstor, delFlagged);
|
---|
313 | xPol.freeStorage(xpolstor, delXPol);
|
---|
314 |
|
---|
315 | return status;
|
---|
316 | }
|
---|
317 |
|
---|
318 | //--------------------------------------------------------- PKSSDwriter::close
|
---|
319 |
|
---|
320 | // Close the SDFITS file.
|
---|
321 |
|
---|
322 | void PKSSDwriter::close()
|
---|
323 | {
|
---|
324 | cSDwriter.close();
|
---|
325 | }
|
---|