source: trunk/src/SDAttr.cc@ 474

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

add Tid aperture efficiency data
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.2 KB
Line 
1//#---------------------------------------------------------------------------
2//# SDAttr.cc: A collection of attributes for different telescopes
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 "SDAttr.h"
33#include <casa/aips.h>
34#include <casa/Arrays/Vector.h>
35#include <casa/Arrays/ArrayMath.h>
36#include <casa/Exceptions.h>
37#include <casa/Quanta/QC.h>
38#include <casa/Quanta/Quantum.h>
39#include <casa/Quanta/MVTime.h>
40
41#include <measures/Measures/MEpoch.h>
42
43#include <scimath/Mathematics/InterpolateArray1D.h>
44
45using namespace casa;
46using namespace asap;
47
48SDAttr::SDAttr ()
49{
50 initData();
51}
52
53SDAttr::SDAttr(const SDAttr& other)
54{
55 initData(); // state just private 'static' data
56}
57
58SDAttr& SDAttr::operator=(const SDAttr& other)
59{
60 if (this != &other) {
61 ; // state just private 'static' data
62 }
63 return *this;
64}
65
66SDAttr::~SDAttr()
67{;}
68
69
70Float SDAttr::diameter (Instrument inst) const
71{
72 Float D = 1.0;
73 switch (inst) {
74 case ATMOPRA:
75 {
76 D = 22.0;
77 }
78 break;
79 case ATPKSMB:
80 case ATPKSHOH:
81 {
82 D = 64.0;
83 }
84 break;
85 case TIDBINBILLA:
86 {
87 D = 70.0;
88 }
89 break;
90 case CEDUNA:
91 {
92 D = 30.0;
93 }
94 break;
95 case HOBART:
96 {
97 D = 26.0;
98 }
99 break;
100 default:
101 {
102 throw(AipsError("Unknown instrument"));
103 }
104 }
105//
106 return D;
107}
108
109Vector<Float> SDAttr::beamEfficiency (Instrument inst, const MEpoch& dateObs, const Vector<Float>& freqs) const
110{
111
112// Look at date where appropriate
113
114 MVTime t(dateObs.getValue());
115 uInt year = t.year();
116//
117 Vector<Float> facs(freqs.nelements(),1.0);
118 switch (inst) {
119 case ATMOPRA:
120 {
121 if (year<2003) {
122 cerr << "There is no beam efficiency data from before 2003 - using 2003 data" << endl;
123 facs = interp (freqs/1.0e9f, MopEtaBeamX_, MopEtaBeam2003Y_);
124 } else if (year==2003) {
125 cerr << "Using beam efficiency data from 2003" << endl;
126 facs = interp (freqs/1.0e9f, MopEtaBeamX_, MopEtaBeam2003Y_);
127 } else {
128 cerr << "Using beam efficiency data from 2004" << endl;
129 facs = interp (freqs/1.0e9f, MopEtaBeamX_, MopEtaBeam2004Y_);
130 }
131 }
132 break;
133 default:
134 {
135 cerr << "No beam efficiency data for this instrument - assuming unity" << endl;
136 }
137 }
138//
139 return facs;
140}
141
142Vector<Float> SDAttr::apertureEfficiency (Instrument inst, const MEpoch& dateObs, const Vector<Float>& freqs) const
143{
144
145// Look at date where appropriate
146
147 MVTime t(dateObs.getValue());
148 uInt year = t.year();
149//
150 Vector<Float> facs(freqs.nelements(),1.0);
151 switch (inst) {
152 case ATMOPRA:
153 {
154 if (year<2004) {
155 cerr << "There is no aperture efficiency data from before 2004 - using 2004 data" << endl;
156 facs = interp (freqs/1.0e9f, MopEtaApX_, MopEtaAp2004Y_);
157 } else {
158 cerr << "Using aperture efficiency data from 2004" << endl;
159 facs = interp (freqs/1.0e9f, MopEtaApX_, MopEtaAp2004Y_);
160 }
161 }
162 break;
163 case TIDBINBILLA:
164 {
165 facs = interp (freqs/1.0e9f, TidEtaApX_, TidEtaApY_);
166 }
167 default:
168 {
169 cerr << "No aperture efficiency data for this instrument - assuming unity" << endl;
170 }
171 }
172 return facs;
173}
174
175Vector<Float> SDAttr::JyPerK (Instrument inst, const MEpoch& dateObs, const Vector<Float>& freqs) const
176{
177
178// FInd what we need
179
180 Vector<Float> etaAp = apertureEfficiency (inst, dateObs, freqs);
181 Float D = diameter(inst);
182
183// Compute it
184
185 Vector<Float> facs(freqs.nelements(),1.0);
186 for (uInt i=0; i<freqs.nelements(); i++) {
187 facs(i) = SDAttr::findJyPerKFac (etaAp(i), D);
188 }
189//
190 return facs;
191}
192
193
194Float SDAttr::findJyPerKFac (Float etaAp, Float D)
195//
196// Converts K -> Jy
197// D in m
198//
199{
200 Double kb = QC::k.getValue(Unit(String("erg/K")));
201 Float gA = C::pi * D * D / 4.0;
202 return (2.0 * 1.0e19 * kb / etaAp / gA);
203}
204
205
206Vector<Float> SDAttr::gainElevationPoly (Instrument inst) const
207{
208
209// Look at date where appropriate
210
211 switch (inst) {
212 case TIDBINBILLA:
213 {
214 return TidGainElPoly_.copy();
215 }
216 break;
217 default:
218 {
219 Vector<Float> t;
220 return t.copy();
221 }
222 }
223}
224
225
226
227
228// Private
229
230Vector<Float> SDAttr::interp (const Vector<Float>& xOut, const Vector<Float>& xIn,
231 const Vector<Float>& yIn) const
232{
233 Int method = 1; // Linear
234 Vector<Float> yOut;
235 Vector<Bool> mOut;
236//
237 Vector<Bool> mIn(xIn.nelements(),True);
238//
239 InterpolateArray1D<Float,Float>::interpolate(yOut, mOut, xOut,
240 xIn, yIn, mIn,
241 method, True, True);
242//
243 return yOut;
244}
245
246void SDAttr::initData ()
247//
248// Mopra data from Mopra web page
249// Tid data from Tid web page
250// X in GHz
251//
252{
253
254// Beam efficiency
255
256 MopEtaBeamX_.resize(3);
257 MopEtaBeamX_(0) = 86.0;
258 MopEtaBeamX_(1) = 100.0;
259 MopEtaBeamX_(2) = 115.0;
260//
261 MopEtaBeam2003Y_.resize(3);
262 MopEtaBeam2003Y_(0) = 0.39;
263 MopEtaBeam2003Y_(1) = 0.37;
264 MopEtaBeam2003Y_(2) = 0.37; // replicated from (1)
265//
266 MopEtaBeam2004Y_.resize(3);
267 MopEtaBeam2004Y_(0) = 0.49;
268 MopEtaBeam2004Y_(1) = 0.44;
269 MopEtaBeam2004Y_(2) = 0.42;
270
271// Aperture efficiency
272
273 MopEtaApX_.resize(2);
274 MopEtaApX_(0) = 86.0;
275 MopEtaApX_(1) = 115.0;
276//
277 MopEtaAp2004Y_.resize(2);
278 MopEtaAp2004Y_(0) = 0.33;
279 MopEtaAp2004Y_(1) = 0.24;
280//
281 TidEtaApX_.resize(2);
282 TidEtaApY_.resize(2);
283 TidEtaApX_(0) = 18.0;
284 TidEtaApX_(1) = 26.5;
285 TidEtaApY_(0) = 0.4848;
286 TidEtaApY_(1) = 0.4848;
287
288// Gain elevation correction polynomial coefficients (for elevation in degrees)
289
290 TidGainElPoly_.resize(3);
291 TidGainElPoly_(0) = 3.58788e-1;
292 TidGainElPoly_(1) = 2.87243e-2;
293 TidGainElPoly_(2) = -3.219093e-4;
294}
Note: See TracBrowser for help on using the repository browser.