[353] | 1 | //#--------------------------------------------------------------------------- |
---|
[878] | 2 | //# STAttr.cc: A collection of attributes for different telescopes |
---|
[353] | 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 <casa/aips.h> |
---|
| 33 | #include <casa/Arrays/Vector.h> |
---|
| 34 | #include <casa/Arrays/ArrayMath.h> |
---|
| 35 | #include <casa/Exceptions.h> |
---|
| 36 | #include <casa/Quanta/QC.h> |
---|
| 37 | #include <casa/Quanta/Quantum.h> |
---|
[394] | 38 | #include <casa/Quanta/MVTime.h> |
---|
[353] | 39 | |
---|
| 40 | #include <measures/Measures/MEpoch.h> |
---|
| 41 | |
---|
| 42 | #include <scimath/Mathematics/InterpolateArray1D.h> |
---|
| 43 | |
---|
[878] | 44 | #include "STAttr.h" |
---|
[717] | 45 | |
---|
[353] | 46 | using namespace casa; |
---|
| 47 | using namespace asap; |
---|
| 48 | |
---|
[878] | 49 | STAttr::STAttr() |
---|
[353] | 50 | { |
---|
[361] | 51 | initData(); |
---|
[353] | 52 | } |
---|
| 53 | |
---|
[2163] | 54 | STAttr::STAttr(const STAttr& other): |
---|
| 55 | Logger() |
---|
[353] | 56 | { |
---|
[2163] | 57 | (void) other; //suppress unused warning |
---|
[717] | 58 | initData(); // state just private 'static' data |
---|
[353] | 59 | } |
---|
| 60 | |
---|
[1188] | 61 | STAttr& STAttr::operator=(const STAttr& other) |
---|
[353] | 62 | { |
---|
| 63 | if (this != &other) { |
---|
[717] | 64 | ; // state just private 'static' data |
---|
[353] | 65 | } |
---|
| 66 | return *this; |
---|
| 67 | } |
---|
| 68 | |
---|
[878] | 69 | STAttr::~STAttr() |
---|
[353] | 70 | {;} |
---|
| 71 | |
---|
| 72 | |
---|
[878] | 73 | Float STAttr::diameter(Instrument inst) const |
---|
[353] | 74 | { |
---|
| 75 | Float D = 1.0; |
---|
[362] | 76 | switch (inst) { |
---|
[1391] | 77 | case ALMA: |
---|
| 78 | { |
---|
| 79 | D = 12.0; |
---|
| 80 | } |
---|
| 81 | break; |
---|
| 82 | |
---|
[362] | 83 | case ATMOPRA: |
---|
| 84 | { |
---|
| 85 | D = 22.0; |
---|
| 86 | } |
---|
| 87 | break; |
---|
| 88 | case ATPKSMB: |
---|
| 89 | case ATPKSHOH: |
---|
| 90 | { |
---|
| 91 | D = 64.0; |
---|
| 92 | } |
---|
| 93 | break; |
---|
| 94 | case TIDBINBILLA: |
---|
| 95 | { |
---|
| 96 | D = 70.0; |
---|
| 97 | } |
---|
| 98 | break; |
---|
| 99 | case CEDUNA: |
---|
| 100 | { |
---|
| 101 | D = 30.0; |
---|
| 102 | } |
---|
| 103 | break; |
---|
[1391] | 104 | case GBT: |
---|
| 105 | { |
---|
| 106 | D = 104.9; |
---|
| 107 | } |
---|
| 108 | break; |
---|
[362] | 109 | case HOBART: |
---|
| 110 | { |
---|
| 111 | D = 26.0; |
---|
| 112 | } |
---|
| 113 | break; |
---|
| 114 | default: |
---|
| 115 | { |
---|
| 116 | throw(AipsError("Unknown instrument")); |
---|
| 117 | } |
---|
[353] | 118 | } |
---|
| 119 | return D; |
---|
| 120 | } |
---|
| 121 | |
---|
[1188] | 122 | Vector<Float> STAttr::beamEfficiency(Instrument inst, const MEpoch& dateObs, |
---|
[717] | 123 | const Vector<Float>& freqs) const |
---|
[353] | 124 | { |
---|
| 125 | |
---|
[717] | 126 | // Look at date where appropriate |
---|
[394] | 127 | MVTime t(dateObs.getValue()); |
---|
| 128 | uInt year = t.year(); |
---|
[717] | 129 | |
---|
[353] | 130 | Vector<Float> facs(freqs.nelements(),1.0); |
---|
[362] | 131 | switch (inst) { |
---|
[717] | 132 | case ATMOPRA: |
---|
[362] | 133 | { |
---|
[717] | 134 | if (year<2003) { |
---|
| 135 | pushLog("There is no beam efficiency data from before 2003" |
---|
| 136 | " - using 2003 data"); |
---|
[1188] | 137 | facs = interp(freqs/1.0e9f, MopEtaBeamX_, MopEtaBeam2003Y_); |
---|
[717] | 138 | } else if (year==2003) { |
---|
| 139 | pushLog("Using beam efficiency data from 2003"); |
---|
[1188] | 140 | facs = interp(freqs/1.0e9f, MopEtaBeamX_, MopEtaBeam2003Y_); |
---|
[717] | 141 | } else { |
---|
| 142 | pushLog("Using beam efficiency data from 2004"); |
---|
[1188] | 143 | facs = interp(freqs/1.0e9f, MopEtaBeamX_, MopEtaBeam2004Y_); |
---|
[717] | 144 | } |
---|
[362] | 145 | } |
---|
| 146 | break; |
---|
[717] | 147 | default: |
---|
| 148 | { |
---|
| 149 | pushLog("No beam efficiency data for this instrument - assuming unity"); |
---|
| 150 | } |
---|
[353] | 151 | } |
---|
| 152 | return facs; |
---|
| 153 | } |
---|
| 154 | |
---|
[1188] | 155 | Vector<Float> STAttr::apertureEfficiency(Instrument inst, |
---|
| 156 | const MEpoch& dateObs, |
---|
[717] | 157 | const Vector<Float>& freqs) const |
---|
[353] | 158 | { |
---|
[1188] | 159 | |
---|
[717] | 160 | // Look at date where appropriate |
---|
| 161 | MVTime t(dateObs.getValue()); |
---|
| 162 | uInt year = t.year(); |
---|
[1188] | 163 | |
---|
[717] | 164 | Vector<Float> facs(freqs.nelements(),1.0); |
---|
| 165 | switch (inst) { |
---|
| 166 | case ATMOPRA: |
---|
| 167 | { |
---|
| 168 | if (year<2004) { |
---|
| 169 | pushLog("There is no aperture efficiency data from before 2004" |
---|
| 170 | " - using 2004 data"); |
---|
| 171 | facs = interp(freqs/1.0e9f, MopEtaApX_, MopEtaAp2004Y_); |
---|
| 172 | } else { |
---|
| 173 | pushLog("Using aperture efficiency data from 2004"); |
---|
| 174 | facs = interp(freqs/1.0e9f, MopEtaApX_, MopEtaAp2004Y_); |
---|
| 175 | } |
---|
| 176 | } |
---|
| 177 | break; |
---|
| 178 | case TIDBINBILLA: |
---|
| 179 | { |
---|
| 180 | facs = interp(freqs/1.0e9f, TidEtaApX_, TidEtaApY_); |
---|
| 181 | } |
---|
| 182 | break; |
---|
| 183 | default: |
---|
| 184 | { |
---|
| 185 | pushLog("No aperture efficiency data for this instrument" |
---|
| 186 | " - assuming unity"); |
---|
| 187 | } |
---|
| 188 | } |
---|
| 189 | return facs; |
---|
[353] | 190 | } |
---|
| 191 | |
---|
[878] | 192 | Vector<Float> STAttr::JyPerK(Instrument inst, const MEpoch& dateObs, |
---|
[717] | 193 | const Vector<Float>& freqs) const |
---|
[353] | 194 | { |
---|
[1188] | 195 | |
---|
[717] | 196 | // Find what we need |
---|
| 197 | Vector<Float> etaAp = apertureEfficiency(inst, dateObs, freqs); |
---|
| 198 | Float D = diameter(inst); |
---|
| 199 | // Compute it |
---|
[353] | 200 | Vector<Float> facs(freqs.nelements(),1.0); |
---|
| 201 | for (uInt i=0; i<freqs.nelements(); i++) { |
---|
[878] | 202 | facs(i) = STAttr::findJyPerK(etaAp(i), D); |
---|
[353] | 203 | } |
---|
| 204 | return facs; |
---|
| 205 | } |
---|
| 206 | |
---|
| 207 | |
---|
[878] | 208 | Float STAttr::findJyPerK(Float etaAp, Float D) |
---|
[717] | 209 | // |
---|
| 210 | // Converts K -> Jy |
---|
| 211 | // D in m |
---|
| 212 | // |
---|
[353] | 213 | { |
---|
[717] | 214 | Double kb = QC::k.getValue(Unit(String("erg/K"))); |
---|
| 215 | Float gA = C::pi * D * D / 4.0; |
---|
| 216 | return (2.0 * 1.0e19 * kb / etaAp / gA); |
---|
[353] | 217 | } |
---|
| 218 | |
---|
| 219 | |
---|
[878] | 220 | Vector<Float> STAttr::gainElevationPoly(Instrument inst) const |
---|
[362] | 221 | { |
---|
[353] | 222 | |
---|
[717] | 223 | // Look at date where appropriate |
---|
| 224 | switch (inst) { |
---|
| 225 | case TIDBINBILLA: |
---|
| 226 | { |
---|
| 227 | return TidGainElPoly_.copy(); |
---|
| 228 | } |
---|
| 229 | break; |
---|
[1345] | 230 | // assume HOH for K-band |
---|
| 231 | case ATPKSHOH: |
---|
| 232 | { |
---|
| 233 | return ParkesGainElPoly_.copy(); |
---|
| 234 | } |
---|
| 235 | break; |
---|
[717] | 236 | default: |
---|
| 237 | { |
---|
| 238 | Vector<Float> t; |
---|
| 239 | return t.copy(); |
---|
| 240 | } |
---|
| 241 | } |
---|
[362] | 242 | } |
---|
| 243 | |
---|
[1188] | 244 | std::string STAttr::feedPolType(Instrument inst) const |
---|
[507] | 245 | { |
---|
[1188] | 246 | std::string type; |
---|
[717] | 247 | switch (inst) { |
---|
| 248 | case ATMOPRA: |
---|
| 249 | case ATPKSMB: |
---|
| 250 | case ATPKSHOH: |
---|
| 251 | { |
---|
[1188] | 252 | type = "linear"; |
---|
[717] | 253 | } |
---|
| 254 | break; |
---|
| 255 | case TIDBINBILLA: |
---|
| 256 | { |
---|
[1188] | 257 | type = "circular"; |
---|
[717] | 258 | } |
---|
| 259 | break; |
---|
| 260 | default: |
---|
| 261 | { |
---|
[1188] | 262 | type = "linear"; |
---|
[717] | 263 | } |
---|
| 264 | } |
---|
| 265 | return type; |
---|
[507] | 266 | } |
---|
[362] | 267 | |
---|
| 268 | |
---|
[353] | 269 | // Private |
---|
[878] | 270 | Vector<Float> STAttr::interp(const Vector<Float>& xOut, |
---|
[1188] | 271 | const Vector<Float>& xIn, |
---|
[717] | 272 | const Vector<Float>& yIn) const |
---|
| 273 | { |
---|
| 274 | Int method = 1; // Linear |
---|
| 275 | Vector<Float> yOut; |
---|
| 276 | Vector<Bool> mOut; |
---|
[353] | 277 | |
---|
[717] | 278 | Vector<Bool> mIn(xIn.nelements(),True); |
---|
[1188] | 279 | |
---|
[717] | 280 | InterpolateArray1D<Float,Float>::interpolate(yOut, mOut, xOut, |
---|
| 281 | xIn, yIn, mIn, |
---|
| 282 | method, True, True); |
---|
[1188] | 283 | |
---|
[717] | 284 | return yOut; |
---|
[353] | 285 | } |
---|
| 286 | |
---|
[1188] | 287 | void STAttr::initData() |
---|
[717] | 288 | // |
---|
| 289 | // Mopra data from Mopra web page |
---|
| 290 | // Tid data from Tid web page |
---|
| 291 | // X in GHz |
---|
| 292 | // |
---|
[362] | 293 | { |
---|
| 294 | |
---|
[717] | 295 | // Beam efficiency |
---|
[353] | 296 | MopEtaBeamX_.resize(3); |
---|
| 297 | MopEtaBeamX_(0) = 86.0; |
---|
| 298 | MopEtaBeamX_(1) = 100.0; |
---|
| 299 | MopEtaBeamX_(2) = 115.0; |
---|
[717] | 300 | |
---|
[353] | 301 | MopEtaBeam2003Y_.resize(3); |
---|
| 302 | MopEtaBeam2003Y_(0) = 0.39; |
---|
| 303 | MopEtaBeam2003Y_(1) = 0.37; |
---|
[717] | 304 | MopEtaBeam2003Y_(2) = 0.37; // replicated from (1) |
---|
[1188] | 305 | |
---|
[353] | 306 | MopEtaBeam2004Y_.resize(3); |
---|
| 307 | MopEtaBeam2004Y_(0) = 0.49; |
---|
| 308 | MopEtaBeam2004Y_(1) = 0.44; |
---|
[361] | 309 | MopEtaBeam2004Y_(2) = 0.42; |
---|
[362] | 310 | |
---|
[717] | 311 | // Aperture efficiency |
---|
[353] | 312 | MopEtaApX_.resize(2); |
---|
| 313 | MopEtaApX_(0) = 86.0; |
---|
| 314 | MopEtaApX_(1) = 115.0; |
---|
[717] | 315 | |
---|
[353] | 316 | MopEtaAp2004Y_.resize(2); |
---|
| 317 | MopEtaAp2004Y_(0) = 0.33; |
---|
| 318 | MopEtaAp2004Y_(1) = 0.24; |
---|
[717] | 319 | |
---|
[467] | 320 | TidEtaApX_.resize(2); |
---|
| 321 | TidEtaApY_.resize(2); |
---|
| 322 | TidEtaApX_(0) = 18.0; |
---|
| 323 | TidEtaApX_(1) = 26.5; |
---|
| 324 | TidEtaApY_(0) = 0.4848; |
---|
| 325 | TidEtaApY_(1) = 0.4848; |
---|
[362] | 326 | |
---|
[717] | 327 | // Gain elevation correction polynomial coefficients (for elevation |
---|
| 328 | // in degrees) |
---|
[362] | 329 | |
---|
| 330 | TidGainElPoly_.resize(3); |
---|
| 331 | TidGainElPoly_(0) = 3.58788e-1; |
---|
| 332 | TidGainElPoly_(1) = 2.87243e-2; |
---|
| 333 | TidGainElPoly_(2) = -3.219093e-4; |
---|
[1345] | 334 | |
---|
[1618] | 335 | // 2009-09-15 - 13mm (22.2GHz) receiver |
---|
[1345] | 336 | ParkesGainElPoly_.resize(3); |
---|
[1618] | 337 | ParkesGainElPoly_(0) = -0.194031; |
---|
| 338 | ParkesGainElPoly_(1) = 0.457724e-1; |
---|
| 339 | ParkesGainElPoly_(2) = -0.438659e-3; |
---|
[353] | 340 | } |
---|
[475] | 341 | |
---|
| 342 | |
---|
[878] | 343 | Instrument STAttr::convertInstrument(const String& instrument, |
---|
[475] | 344 | Bool throwIt) |
---|
| 345 | { |
---|
[717] | 346 | String t(instrument); |
---|
| 347 | t.upcase(); |
---|
[1188] | 348 | |
---|
[878] | 349 | // The strings are what STReader returns, after cunning |
---|
[717] | 350 | // interrogation of station names... :-( |
---|
| 351 | Instrument inst = asap::UNKNOWNINST; |
---|
| 352 | if (t==String("DSS-43")) { |
---|
| 353 | inst = TIDBINBILLA; |
---|
[1391] | 354 | } else if (t==String("ALMA")) { |
---|
| 355 | inst = ALMA; |
---|
[717] | 356 | } else if (t==String("ATPKSMB")) { |
---|
| 357 | inst = ATPKSMB; |
---|
| 358 | } else if (t==String("ATPKSHOH")) { |
---|
| 359 | inst = ATPKSHOH; |
---|
| 360 | } else if (t==String("ATMOPRA")) { |
---|
| 361 | inst = ATMOPRA; |
---|
| 362 | } else if (t==String("CEDUNA")) { |
---|
| 363 | inst = CEDUNA; |
---|
[1391] | 364 | } else if (t==String("GBT")) { |
---|
| 365 | inst = GBT; |
---|
[717] | 366 | } else if (t==String("HOBART")) { |
---|
| 367 | inst = HOBART; |
---|
| 368 | } else { |
---|
| 369 | if (throwIt) { |
---|
| 370 | throw AipsError("Unrecognized instrument" |
---|
| 371 | " - use function scan.set_instrument to set"); |
---|
| 372 | } |
---|
| 373 | } |
---|
| 374 | return inst; |
---|
[475] | 375 | } |
---|