source: trunk/src/SDAttr.cc@ 407

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

Detect Year now in selecting beam/aperture efficiency data
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 6.9 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 default:
164 {
165 cerr << "No aperture efficiency data for this instrument - assuming unity" << endl;
166 }
167 }
168 return facs;
169}
170
171Vector<Float> SDAttr::JyPerK (Instrument inst, const MEpoch& dateObs, const Vector<Float>& freqs) const
172{
173
174// FInd what we need
175
176 Vector<Float> etaAp = apertureEfficiency (inst, dateObs, freqs);
177 Float D = diameter(inst);
178
179// Compute it
180
181 Vector<Float> facs(freqs.nelements(),1.0);
182 for (uInt i=0; i<freqs.nelements(); i++) {
183 facs(i) = SDAttr::findJyPerKFac (etaAp(i), D);
184 }
185//
186 return facs;
187}
188
189
190Float SDAttr::findJyPerKFac (Float etaAp, Float D)
191//
192// Converts K -> Jy
193// D in m
194//
195{
196 Double kb = QC::k.getValue(Unit(String("erg/K")));
197 Float gA = C::pi * D * D / 4.0;
198 return (2.0 * 1.0e19 * kb / etaAp / gA);
199}
200
201
202Vector<Float> SDAttr::gainElevationPoly (Instrument inst) const
203{
204
205// Look at date where appropriate
206
207 switch (inst) {
208 case TIDBINBILLA:
209 {
210 return TidGainElPoly_.copy();
211 }
212 break;
213 default:
214 {
215 Vector<Float> t;
216 return t.copy();
217 }
218 }
219}
220
221
222
223
224// Private
225
226Vector<Float> SDAttr::interp (const Vector<Float>& xOut, const Vector<Float>& xIn,
227 const Vector<Float>& yIn) const
228{
229 Int method = 1; // Linear
230 Vector<Float> yOut;
231 Vector<Bool> mOut;
232//
233 Vector<Bool> mIn(xIn.nelements(),True);
234//
235 InterpolateArray1D<Float,Float>::interpolate(yOut, mOut, xOut,
236 xIn, yIn, mIn,
237 method, True, True);
238//
239 return yOut;
240}
241
242void SDAttr::initData ()
243//
244// Mopra data from online Mopra guide.
245//
246{
247
248// Beam efficiency
249
250 MopEtaBeamX_.resize(3);
251 MopEtaBeamX_(0) = 86.0;
252 MopEtaBeamX_(1) = 100.0;
253 MopEtaBeamX_(2) = 115.0;
254//
255 MopEtaBeam2003Y_.resize(3);
256 MopEtaBeam2003Y_(0) = 0.39;
257 MopEtaBeam2003Y_(1) = 0.37;
258 MopEtaBeam2003Y_(2) = 0.37; // replicated from (1)
259//
260 MopEtaBeam2004Y_.resize(3);
261 MopEtaBeam2004Y_(0) = 0.49;
262 MopEtaBeam2004Y_(1) = 0.44;
263 MopEtaBeam2004Y_(2) = 0.42;
264
265// Aperture efficiency
266
267 MopEtaApX_.resize(2);
268 MopEtaApX_(0) = 86.0;
269 MopEtaApX_(1) = 115.0;
270//
271 MopEtaAp2004Y_.resize(2);
272 MopEtaAp2004Y_(0) = 0.33;
273 MopEtaAp2004Y_(1) = 0.24;
274
275// Gain elevation correction polynomial coefficients (for elevation in degrees)
276
277 TidGainElPoly_.resize(3);
278 TidGainElPoly_(0) = 3.58788e-1;
279 TidGainElPoly_(1) = 2.87243e-2;
280 TidGainElPoly_(2) = -3.219093e-4;
281}
Note: See TracBrowser for help on using the repository browser.