source: tags/casa3.3.0asap/external-alma/asdm2ASAP/ASDMFiller.cc@ 2803

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

New Development: No

JIRA Issue: Yes CAS-1913

Ready for Test: Yes

Interface Changes: 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...

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