source: branches/hpc34/external-alma/asdm2ASAP/ASDMFiller.cc@ 3098

Last change on this file since 3098 was 2407, checked in by Takeshi Nakazato, 13 years ago

New Development: No

JIRA Issue: No

Ready for Test: Yes

Interface Changes: Yes/No

What Interface Changed: Please list interface changes

Test Programs: test_importasdm_sd

Put in Release Notes: Yes/No

Module(s): Module Names change impacts.

Description: Describe your changes here...

some clean up

File size: 23.7 KB
Line 
1#include<iostream>
2
3#include <STHeader.h>
4
5#include <measures/Measures/MEpoch.h>
6#include <measures/Measures/MPosition.h>
7#include <measures/Measures/MDirection.h>
8#include <measures/Measures/MCDirection.h>
9#include <measures/Measures/MFrequency.h>
10#include <measures/Measures/MCFrequency.h>
11#include <measures/Measures/MeasFrame.h>
12#include <measures/Measures/MeasConvert.h>
13#include <casa/Arrays/Vector.h>
14#include <casa/Arrays/Matrix.h>
15#include <casa/Quanta/MVTime.h>
16#include <casa/Logging/LogMessage.h>
17
18#include "ASDMFiller.h"
19
20using namespace std ;
21using namespace casa ;
22using namespace asap ;
23
24ASDMFiller::ASDMFiller( CountedPtr<Scantable> stable )
25 : FillerBase( stable ),
26 antennaId_( -1 ),
27 antennaName_( "" ),
28 className_("ASDMFiller")
29{
30 reader_ = new ASDMReader() ;
31}
32
33ASDMFiller::~ASDMFiller()
34{
35 // nothing to do?
36 logsink_ = 0 ;
37}
38
39void ASDMFiller::setLogger( CountedPtr<LogSinkInterface> &logsink )
40{
41 logsink_ = logsink ;
42 if ( !(reader_.null()) ) {
43 reader_->setLogger( logsink ) ;
44 }
45}
46
47bool ASDMFiller::open( const string &filename, const Record &rec )
48{
49 String funcName = "open" ;
50 bool status = reader_->open( filename, rec ) ;
51
52 antennaId_ = reader_->getAntennaId() ;
53 antennaName_ = reader_->getAntennaName() ;
54
55 //logsink_->postLocally( LogMessage("antennaId_ = "+String::toString(antennaId_),LogOrigin(className_,funcName,WHERE)) ) ;
56 //logsink_->postLocally( LogMessage("antennaName_ = "+antennaName_,LogOrigin(className_,funcName,WHERE)) ) ;
57
58 return status ;
59}
60
61void ASDMFiller::fill()
62{
63 String funcName = "fill" ;
64
65 // header
66 fillHeader() ;
67
68 // set Frame for FREQUENCIES table
69 string sFreqFrame = reader_->getFrame() ;
70 //MFrequency::Types freqFrame = toFrameType( sFreqFrame ) ;
71 MFrequency::Types freqFrame = MFrequency::LSRK ;
72 table_->frequencies().setFrame( freqFrame, false ) ;
73 table_->frequencies().setFrame( freqFrame, true ) ;
74 //logsink_->postLocally( LogMessage("sFreqFrame = "+sFreqFrame,LogOrigin(className_,funcName,WHERE)) ) ;
75
76 Vector<casa::Double> antpos = table_->getHeader().antennaposition ;
77
78 // data selection
79 reader_->select() ;
80
81 // pick up valid configDescriptionId
82 Vector<uInt> configDescIdList = reader_->getConfigDescriptionIdList() ;
83 uInt numConfigDescId = configDescIdList.size() ;
84
85 //logsink_->postLocally( LogMessage("configDescIdList = "+String::toString(configDescIdList),LogOrigin(className_,funcName,WHERE)) ) ;
86
87 // get field list
88 Vector<uInt> fieldIdList = reader_->getFieldIdList() ;
89 uInt numFieldId = fieldIdList.size() ;
90
91 //logsink_->postLocally( LogMessage("fieldIdList = "+String::toString(fieldIdList),LogOrigin(className_,funcName,WHERE)) ) ;
92
93 // BEAMNO is always 0 since ALMA antenna is single beam
94 uInt beamno = 0 ;
95
96 // REFBEAMNO is -1
97 setReferenceBeam() ;
98
99 // fill FOCUS_ID and add FOCUS row if necessary
100 setFocus() ;
101
102 // CYCLENO
103 map< unsigned int, unsigned int > cycleno ;
104 map< unsigned int, unsigned int >::iterator citer ;
105
106 for ( unsigned int ifield = 0 ; ifield < numFieldId ; ifield++ ) {
107 for ( uInt icon = 0 ; icon < numConfigDescId ; icon++ ) {
108 //logsink_->postLocally( LogMessage("start configDescId "+String::toString(configDescIdList[icon])+" fieldId "+String::toString(fieldIdList[ifield]),LogOrigin(className_,funcName,WHERE)) ) ;
109
110 if ( !(reader_->setMainRow( configDescIdList[icon], fieldIdList[ifield] )) ) {
111 //logsink_->postLocally( LogMessage("skip configDescId "+String::toString(configDescIdList[icon])+" fieldId "+String::toString(fieldIdList[ifield]),LogOrigin(className_,funcName,WHERE)) ) ;
112 continue ;
113 }
114
115 // number of rows
116 uInt nrow = reader_->getNumMainRow() ;
117
118 //logsink_->postLocally( LogMessage("There are "+String::toString(nrow)+" rows in Main table corresponding to configDescId "+String::toString(configDescIdList[icon])+" fieldId "+String::toString(fieldIdList[ifield]),LogOrigin(className_,funcName,WHERE)) ) ;
119
120 for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
121
122 // set main row
123 if ( !(reader_->setMainRow( irow )) ) {
124 // skip row since the row doesn't have valid configDescId
125 //logsink_->postLocally( LogMessage("skip "+String::toString(irow)+" since ASDMReader::setMainrow() returns False",LogOrigin(className_,funcName,WHERE)) ) ;
126 continue ;
127 }
128
129 // scan and subscan
130 unsigned int scanno = reader_->getScanNoOfCurrentRow() ;
131 //logsink_->postLocally( LogMessage("scanno = "+String::toString(scanno),LogOrigin(className_,funcName,WHERE)) ) ;
132 citer = cycleno.find( scanno ) ;
133 if ( citer == cycleno.end() )
134 cycleno[scanno] = 0 ;
135
136 // set data
137 if ( !(reader_->setData()) ) {
138 // skip row since reader failed to retrieve data
139 //logsink_->postLocally( LogMessage("skip "+String::toString(irow)+" since ASDMReader::setData() returns False",LogOrigin(className_,funcName,WHERE)) ) ;
140 continue ;
141 }
142
143 unsigned int numData = reader_->getNumData() ;
144 double refpix = 0.0 ;
145 double refval = 0.0 ;
146 double incr = 0.0 ;
147 string freqref = "" ;
148
149 //logsink_->postLocally( LogMessage("numData = "+String::toString(numData),LogOrigin(className_,funcName,WHERE)) ) ;
150 for ( unsigned int idata = 0 ; idata < numData ; idata++ ) {
151 // prepare to extract binary data
152 //logsink_->postLocally( LogMessage("prepare data...",LogOrigin(className_,funcName,WHERE)) ) ;
153 reader_->prepareData( idata ) ;
154
155 // subscan number
156 unsigned int subscanno = reader_->getSubscanNo() ;
157 //logsink_->postLocally( LogMessage("subscanno = "+String::toString(subscanno),LogOrigin(className_,funcName,WHERE)) ) ;
158
159 // IFNO
160 uInt ifno = reader_->getIFNo() ;
161 //logsink_->postLocally( LogMessage("ifno = "+String::toString(ifno),LogOrigin(className_,funcName,WHERE)) ) ;
162 // source spec
163 int srctype = reader_->getSrcType( scanno, subscanno ) ;
164 //logsink_->postLocally( LogMessage("srctype = "+String::toString(srctype),LogOrigin(className_,funcName,WHERE)) ) ;
165 string srcname ;
166 string fieldname ;
167 vector<double> srcDirection ;
168 vector<double> srcProperMotion ;
169 double sysVel ;
170 vector<double> rf ;
171 reader_->getSourceProperty( srcname,
172 fieldname,
173 srcDirection,
174 srcProperMotion,
175 sysVel,
176 rf ) ;
177 //logsink_->postLocally( LogMessage("srcname = "+String::toString(srcname),LogOrigin(className_,funcName,WHERE)) ) ;
178
179 // fill MOLECULE_ID and add MOLECULES row if necessary
180 Vector<casa::Double> restFreqs( rf.size() ) ;
181 for ( uInt i = 0 ; i < rf.size() ; i++ )
182 restFreqs[i] = (casa::Double)(rf[i]) ;
183 setMolecule( restFreqs ) ;
184
185 // time and interval
186 casa::Double mjd = (casa::Double)(reader_->getTime()) ;
187 casa::Double interval = (casa::Double)(reader_->getInterval()) ;
188 //logsink_->postLocally( LogMessage("mjd = "+String::toString(mjd),LogOrigin(className_,funcName,WHERE)) ) ;
189
190 // fill TIME and INTERVAL
191 setTime( mjd, interval ) ;
192
193 // fill SRCNAME, SRCTYPE, FIELDNAME, SRCDIRECTION, SRCPROPERMOTION, and SRCVELOCITY
194 Vector<casa::Double> srcDir( 2 ) ;
195 srcDir[0] = (casa::Double)(srcDirection[0]) ;
196 srcDir[1] = (casa::Double)(srcDirection[1]) ;
197 Vector<casa::Double> srcPM( 2 ) ;
198 srcPM[0] = (casa::Double)(srcProperMotion[0]) ;
199 srcPM[1] = (casa::Double)(srcProperMotion[1]) ;
200 setSource( srcname, srctype, fieldname, srcDir, srcPM, (casa::Double)sysVel ) ;
201
202 // fill FLAGROW
203 unsigned int flagrow = reader_->getFlagRow() ;
204 setFlagrow( (uInt)flagrow ) ;
205 //logsink_->postLocally( LogMessage("flagrow = "+String::toString(flagrow),LogOrigin(className_,funcName,WHERE)) ) ;
206 // fill WEATHER_ID and add WEATHER row if necessary
207 float temperature ;
208 float pressure ;
209 float humidity ;
210 float windspeed ;
211 float windaz ;
212 reader_->getWeatherInfo( temperature,
213 pressure,
214 humidity,
215 windspeed,
216 windaz ) ;
217 setWeather2( (casa::Float)temperature,
218 (casa::Float)pressure,
219 (casa::Float)humidity,
220 (casa::Float)windspeed,
221 (casa::Float)windaz ) ;
222 //logsink_->postLocally( LogMessage("temperature = "+String::toString(temperature),LogOrigin(className_,funcName,WHERE)) ) ;
223 // fill AZIMUTH, ELEVATION, DIRECTION and SCANRATE
224 vector<double> dir ;
225 double az ;
226 double el ;
227 vector<double> srate ;
228 reader_->getPointingInfo( dir,
229 az,
230 el,
231 srate ) ;
232 Vector<casa::Double> scanRate( 2, 0.0 ) ;
233 Vector<casa::Double> direction( 2, 0.0 ) ;
234 if ( srate.size() > 0 ) {
235 scanRate[0] = (casa::Double)(srate[0]) ;
236 scanRate[1] = (casa::Double)(srate[1]) ;
237 }
238 setScanRate( scanRate ) ;
239 if ( dir.size() > 0 ) {
240 direction[0] = (casa::Double)(dir[0]) ;
241 direction[1] = (casa::Double)(dir[1]) ;
242 }
243 else {
244 toJ2000( direction, az, el, mjd, antpos ) ;
245 }
246 //logsink_->postLocally( LogMessage("direction = "+String::toString(direction),LogOrigin(className_,funcName,WHERE)) ) ;
247 setDirection( direction, (casa::Float)az, (casa::Float)el ) ;
248
249 // REFPIX, REFVAL, INCREMENT
250 String ifkey = getIFKey( ifno ) ;
251 if ( ifrec_.isDefined( ifkey ) ) {
252 getFrequencyRec( ifkey, refpix, refval, incr ) ;
253 }
254 else {
255 reader_->getFrequency( refpix, refval, incr, freqref ) ;
256 refval = (double)toLSRK( casa::Double(refval),
257 String(freqref),
258 mjd,
259 antpos,
260 //direction,
261 srcDir,
262 "J2000" ) ;
263 setFrequencyRec( ifkey, refpix, refval, incr ) ;
264 }
265
266 // fill FREQ_ID and add FREQUENCIES row if necessary
267 setFrequency( (casa::Double)refpix, (casa::Double)refval, (casa::Double)incr ) ;
268
269 // loop on polarization
270 vector<unsigned int> dataShape = reader_->getDataShape() ;
271// ostringstream oss ;
272// for ( unsigned int i = 0 ; i < dataShape.size() ; i++ ) {
273// if ( i == 0 )
274// oss << "dataShape=[" << dataShape[i] << ", " ;
275// else if ( i == dataShape.size()-1 )
276// oss << dataShape[i] << "]" ;
277// else
278// oss << dataShape[i] << ", " ;
279// }
280// logsink_->postLocally( LogMessage(oss.str(),LogOrigin(className_,funcName,WHERE)) ) ;
281
282 unsigned int numPol = dataShape[0] ;
283 unsigned int numChan = dataShape[1] ;
284
285 //logsink_->postLocally( LogMessage("numPol = "+String::toString(numPol),LogOrigin(className_,funcName,WHERE)) ) ;
286
287 // OPACITY
288 vector<float> tau = reader_->getOpacity() ;
289 Vector<casa::Float> opacity = toVector( tau, numPol ) ;
290
291 // SPECTRA, FLAGTRA, TSYS, TCAL
292 float *sp = reader_->getSpectrum() ;
293 //logsink_->postLocally( LogMessage("sp[0] = "+String::toString(sp[0]),LogOrigin(className_,funcName,WHERE)) ) ;
294 vector< vector<float> > ts ;
295 vector< vector<float> > tc ;
296 reader_->getTcalAndTsys( tc, ts ) ;
297 Matrix<casa::Float> spectra = toMatrix( sp, numPol, numChan ) ;
298 //logsink_->postLocally( LogMessage("spectra(0,0) = "+String::toString(spectra(0,0)),LogOrigin(className_,funcName,WHERE)) ) ;
299 Vector<uChar> flagtra( numChan, 0 ) ;
300 Matrix<casa::Float> tsys = toMatrix( ts, numPol, numChan ) ;
301 Matrix<casa::Float> tcal = toMatrix( tc, numPol, numChan ) ;
302 //logsink_->postLocally( LogMessage("tsys(0,0) = "+String::toString(tsys(0,0)),LogOrigin(className_,funcName,WHERE)) ) ;
303// String caltime = "" ;
304// if ( anyNE( tcal, (casa::Float)1.0 ) )
305// caltime = toTcalTime( mjd ) ;
306 String caltime = toTcalTime( mjd ) ;
307
308 for ( unsigned int ipol = 0 ; ipol < numPol ; ipol++ ) {
309
310 // fill SCANNO, CYCLENO, IFNO, POLNO, and BEAMNO
311 setIndex( (uInt)scanno-1, (uInt)cycleno[scanno], ifno, ipol, beamno ) ;
312
313 // fill SPECTRA, FLAGTRA, TSYS
314 setSpectrum( spectra.row(ipol), flagtra, tsys.row(ipol) ) ;
315
316 // fill TCAL_ID and add TCAL row if necessary
317 setTcal2( caltime, tcal.row(ipol) ) ;
318
319 // fill OPACITY
320 setOpacity( opacity[ipol] ) ;
321
322 // commit row
323 commitRow() ;
324 }
325
326 // increment CYCLENO
327 cycleno[scanno]++ ;
328 }
329 }
330 }
331 }
332
333 return ;
334}
335
336void ASDMFiller::close()
337{
338 reader_->close() ;
339 reader_ = 0 ;
340
341 return ;
342}
343
344void ASDMFiller::fillHeader()
345{
346 STHeader hdr ;
347
348 reader_->fillHeader( hdr.nchan,
349 hdr.npol,
350 hdr.nif,
351 hdr.nbeam,
352 hdr.observer,
353 hdr.project,
354 hdr.obstype,
355 hdr.antennaname,
356 hdr.antennaposition,
357 hdr.equinox,
358 hdr.freqref,
359 hdr.reffreq,
360 hdr.bandwidth,
361 hdr.utc,
362 hdr.fluxunit,
363 hdr.epoch,
364 hdr.poltype ) ;
365
366 if ( hdr.freqref != "LSRK" ) {
367// String freqref ;
368// if (hdr.freqref == "TOPOCENT") {
369// freqref = "TOPO";
370// } else if (hdr.freqref == "GEOCENTR") {
371// freqref = "GEO";
372// } else if (hdr.freqref == "BARYCENT") {
373// freqref = "BARY";
374// } else if (hdr.freqref == "GALACTOC") {
375// freqref = "GALACTO";
376// } else if (hdr.freqref == "LOCALGRP") {
377// freqref = "LGROUP";
378// } else if (hdr.freqref == "CMBDIPOL") {
379// freqref = "CMB";
380// } else if (hdr.freqref == "SOURCE") {
381// freqref = "REST";
382// }
383 vector<double> sdir ;
384 string ref ;
385 reader_->getSourceDirection( sdir, ref ) ;
386 Vector<casa::Double> sourceDir( sdir ) ;
387 hdr.reffreq = toLSRK( hdr.reffreq, hdr.freqref, hdr.utc, hdr.antennaposition, sdir, String(ref) ) ;
388 hdr.freqref = "LSRK" ;
389 }
390
391 setHeader( hdr ) ;
392}
393
394String ASDMFiller::getIFKey( uInt ifno )
395{
396 return "IFNO"+String::toString( ifno ) ;
397}
398
399void ASDMFiller::getFrequencyRec( String key,
400 double &refpix,
401 double &refval,
402 double &incr )
403{
404 Record frec = ifrec_.asRecord( key ) ;
405 refpix = frec.asdouble( "REFPIX" ) ;
406 refval = frec.asdouble( "REFVAL" ) ;
407 incr = frec.asdouble( "INCREMENT" ) ;
408}
409
410void ASDMFiller::setFrequencyRec( String key,
411 double refpix,
412 double refval,
413 double incr )
414{
415 Record frec ;
416 frec.define( "REFPIX", refpix ) ;
417 frec.define( "REFVAL", refval ) ;
418 frec.define( "INCREMENT", incr ) ;
419 ifrec_.defineRecord( key, frec ) ;
420}
421
422Matrix<casa::Float> ASDMFiller::toMatrix( float *sp,
423 unsigned int npol,
424 unsigned int nchan )
425{
426 Matrix<casa::Float> mSp( npol, nchan ) ;
427 if ( npol <= 2 ) {
428 // 1 or 2 polarization case
429 for ( unsigned int ich = 0 ; ich < nchan ; ich++ ) {
430 for ( unsigned int ipol = 0 ; ipol < npol ; ipol++ ) {
431 mSp(ipol,ich) = (casa::Float)(sp[npol*ich+ipol]) ;
432 }
433 }
434 }
435 else {
436 // 4 polarization case
437 for ( unsigned int ich = 0 ; ich < nchan ; ich++ ) {
438 mSp(0,ich) = (casa::Float)(sp[4*ich]) ; // Re(XX)
439 mSp(1,ich) = (casa::Float)(sp[4*ich+4]) ; // Re(YY)
440 mSp(2,ich) = (casa::Float)(sp[4*ich+2]) ; // Re(XY)
441 mSp(3,ich) = (casa::Float)(sp[4*ich+3]) ; // Im(XY)
442 }
443 }
444 return mSp ;
445}
446
447Matrix<casa::Float> ASDMFiller::toMatrix( vector< vector<float> > &tsys,
448 unsigned int npol,
449 unsigned int nchan )
450{
451 unsigned int numRec = tsys.size() ;
452 unsigned int numChan = tsys[0].size() ;
453 Matrix<casa::Float> ret ;
454 if ( npol == numRec && nchan == numChan ) {
455 ret.resize( npol, nchan ) ;
456 for ( unsigned int ip = 0 ; ip < npol ; ip++ )
457 for ( unsigned int ic = 0 ; ic < nchan ; ic++ )
458 ret( ip, ic ) = (casa::Float)(tsys[ip][ic]) ;
459 }
460 else if ( npol == numRec && numChan == 1 ) {
461 ret.resize( npol, 1 ) ;
462 for ( unsigned int ip = 0 ; ip < npol ; ip++ )
463 ret( ip, 0 ) = (casa::Float)(tsys[0][0]) ;
464 }
465 else if ( numRec == 1 && nchan == numChan ) {
466 ret.resize( npol, nchan ) ;
467 for ( unsigned int ip = 0 ; ip < npol ; ip++ )
468 for ( unsigned int ic = 0 ; ic < nchan ; ic++ )
469 ret( ip, ic ) = (casa::Float)(tsys[0][ic]) ;
470 }
471 else if ( numRec == 1 && numChan == 1 ) {
472 ret.resize( npol, 1 ) ;
473 for ( unsigned int ip = 0 ; ip < npol ; ip++ )
474 ret( ip, 0 ) = (casa::Float)(tsys[0][0]) ;
475 }
476 else if ( numRec == 2 && npol == 4 && numChan == nchan ) {
477 // TODO: How to determine Tsys for XY?
478 // at the moment Tsys[XY] = 0.5*(Tsys[X]+Tsys[Y])
479 ret.resize( npol, nchan ) ;
480 for ( unsigned int ic = 0 ; ic < nchan ; ic++ ) {
481 casa::Float tsysxy = (casa::Float)(0.5*(tsys[0][ic]+tsys[1][ic])) ;
482 ret( 0, ic ) = (casa::Float)(tsys[0][ic]) ;
483 ret( 1, ic ) = (casa::Float)(tsys[1][ic]) ;
484 ret( 2, ic ) = tsysxy ;
485 ret( 3, ic ) = tsysxy ;
486 }
487 }
488 else if ( numRec == 2 && npol == 4 && numChan == 1 ) {
489 // TODO: How to determine Tsys for XY?
490 // at the moment Tsys[XY] = 0.5*(Tsys[X]+Tsys[Y])
491 ret.resize( npol, 1 ) ;
492 casa::Float tsysxy = (casa::Float)(0.5*(tsys[0][0]+tsys[1][0])) ;
493 ret( 0, 0 ) = (casa::Float)(tsys[0][0]) ;
494 ret( 1, 0 ) = (casa::Float)(tsys[1][0]) ;
495 ret( 2, 0 ) = tsysxy ;
496 ret( 3, 0 ) = tsysxy ;
497 }
498 else {
499 // I don't know how to handle ...
500 for ( unsigned int ip = 0 ; ip < npol ; ip++ )
501 for ( unsigned int ic = 0 ; ic < nchan ; ic++ )
502 ret( ip, ic ) = (casa::Float)(tsys[0][ic]) ;
503 }
504 return ret ;
505}
506
507Vector<casa::Float> ASDMFiller::toVector( vector<float> &tau,
508 unsigned int npol )
509{
510 String funcName = "toVector" ;
511
512 Vector<casa::Float> ret( npol ) ;
513 //logsink_->postLocally( LogMessage("tau0="+String::toString(tau[0]),LogOrigin(className_,funcName,WHERE)) ) ;
514 if ( npol == 4 ) {
515 ret[0] = (casa::Float)tau[0] ;
516 ret[1] = (casa::Float)tau[1] ;
517 ret[2] = 0.5 * ( ret[0] + ret[1] ) ;
518 ret[3] = ret[2] ;
519 }
520 else if ( npol == tau.size() ) {
521 for ( unsigned int ipol = 0 ; ipol < npol ; ipol++ )
522 ret[ipol] = (casa::Float)tau[ipol] ;
523 }
524 else {
525 // I don't know how to handle...
526 for ( unsigned int ipol = 0 ; ipol < npol ; ipol++ )
527 ret[ipol] = (casa::Float)tau[0] ;
528 }
529 //logsink_->postLocally( LogMessage("tau="+String::toString(ret),LogOrigin(className_,funcName,WHERE)) ) ;
530 return ret ;
531}
532
533String ASDMFiller::toTcalTime( casa::Double mjd )
534{
535 return MVTime( mjd ).string( MVTime::YMD ) ;
536}
537
538void ASDMFiller::toJ2000( Vector<casa::Double> &dir,
539 double az,
540 double el,
541 casa::Double mjd,
542 Vector<casa::Double> antpos )
543{
544 String funcName = "toJ2000" ;
545
546 Vector<casa::Double> azel( 2 ) ;
547 azel[0] = az ;
548 azel[1] = el ;
549// MEpoch me( Quantity( mjd, "d" ), MEpoch::UTC ) ;
550// Vector<Quantity> qantpos( 3 ) ;
551// qantpos[0] = Quantity( antpos[0], "m" ) ;
552// qantpos[1] = Quantity( antpos[1], "m" ) ;
553// qantpos[2] = Quantity( antpos[2], "m" ) ;
554// MPosition mp( MVPosition( qantpos ),
555// MPosition::ITRF ) ;
556// // mp.print( os_.output() ) ;
557// MeasFrame mf( me, mp ) ;
558// MDirection::Convert toj2000( MDirection::AZELGEO,
559// MDirection::Ref( MDirection::J2000, mf ) ) ;
560// dir = toj2000( azel ).getAngle( "rad" ).getValue() ;
561 dir = toJ2000( azel, "AZELGEO", mjd, antpos ) ;
562 //logsink_->postLocally( LogMessage("dir = "+String::toString(dir),LogOrigin(className_,funcName,WHERE)) ) ;
563}
564
565Vector<casa::Double> ASDMFiller::toJ2000( Vector<casa::Double> dir,
566 String dirref,
567 casa::Double mjd,
568 Vector<casa::Double> antpos )
569{
570 Vector<casa::Double> newd( dir ) ;
571 if ( dirref != "J2000" ) {
572 MEpoch me( Quantity( mjd, "d" ), MEpoch::UTC ) ;
573 Vector<Quantity> qantpos( 3 ) ;
574 qantpos[0] = Quantity( antpos[0], "m" ) ;
575 qantpos[1] = Quantity( antpos[1], "m" ) ;
576 qantpos[2] = Quantity( antpos[2], "m" ) ;
577 MPosition mp( MVPosition( qantpos ),
578 MPosition::ITRF ) ;
579 // mp.print( os_.output() ) ;
580 MeasFrame mf( me, mp ) ;
581 MDirection::Types dirtype ;
582 Bool b = MDirection::getType( dirtype, dirref ) ;
583 if ( b ) {
584 MDirection::Convert toj2000( dirtype,
585 MDirection::Ref( MDirection::J2000, mf ) ) ;
586 newd = toj2000( dir ).getAngle( "rad" ).getValue() ;
587 }
588 }
589 return newd ;
590}
591
592MFrequency::Types ASDMFiller::toFrameType( string &s )
593{
594 MFrequency::Types ftype = MFrequency::DEFAULT ;
595 if ( s == "LABREST" )
596 ftype = MFrequency::REST ;
597 else {
598 Bool b = MFrequency::getType( ftype, String(s) ) ;
599 if (!b)
600 ftype = MFrequency::DEFAULT ;
601 }
602 return ftype ;
603}
604
605casa::Double ASDMFiller::toLSRK( casa::Double freq,
606 String freqref,
607 casa::Double utc,
608 Vector<casa::Double> antpos,
609 Vector<casa::Double> dir,
610 String dirref )
611{
612 String funcName = "toLSRK" ;
613
614 //logsink_->postLocally( LogMessage("freqref = "+freqref,LogOrigin(className_,funcName,WHERE)) ) ;
615 casa::Double newf = freq ;
616 if ( freqref != "LSRK" ) {
617 MEpoch me( Quantum<casa::Double>( utc, Unit("d") ), MEpoch::UTC ) ;
618 Vector< Quantum<casa::Double> > antposQ( 3 ) ;
619 for ( int i = 0 ; i < 3 ; i++ )
620 antposQ[i] = Quantum<casa::Double>( antpos[i], Unit("m") ) ;
621 MPosition mp( antposQ, MPosition::ITRF ) ;
622 MDirection::Types dirtype ;
623 Bool b = MDirection::getType( dirtype, dirref ) ;
624 if ( !b )
625 dirtype = MDirection::J2000 ;
626 MDirection md( Quantum<casa::Double>( dir[0], Unit("rad") ),
627 Quantum<casa::Double>( dir[1], Unit("rad") ),
628 dirtype ) ;
629 MeasFrame mf( me, mp, md ) ;
630 MFrequency::Types freqtype ;
631 b = MFrequency::getType( freqtype, freqref ) ;
632 if ( !b )
633 freqtype = MFrequency::TOPO ;
634 MFrequency::Convert tolsr( freqtype,
635 MFrequency::Ref( MFrequency::LSRK, mf ) ) ;
636 newf = tolsr( Quantum<casa::Double>( freq, Unit("Hz") ) ).get( "Hz" ).getValue() ;
637 //logsink_->postLocally( LogMessage("freq = "+String::toString(freq)+", newf = "+String::toString(newf),LogOrigin(className_,funcName,WHERE)) ) ;
638 }
639 return newf ;
640}
Note: See TracBrowser for help on using the repository browser.