source: trunk/external-alma/oldasdm2ASAP/OldASDMReader.cc@ 2290

Last change on this file since 2290 was 2273, 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: List test programs

Put in Release Notes: Yes/No

Module(s): Module Names change impacts.

Description: Describe your changes here...

Support latest ASDM definition in asdm2ASAP.
Now asdm2ASAP is linked to libalma_v3.so to be able to read latest ASDM.

The importasdm task invokes asdm2ASAP when singledish is True.

Older version of asdm2ASAP is renamec as oldasdm2ASAP.
The oldasdm2ASAP and related classes are built by linking to libalma.so.
It is installed but not used by any tasks.


File size: 68.6 KB
Line 
1#include <iostream>
2#include <sstream>
3
4#include <measures/Measures/MEpoch.h>
5#include <measures/Measures/MPosition.h>
6#include <measures/Measures/MDirection.h>
7#include <measures/Measures/MCDirection.h>
8#include <measures/Measures/MeasFrame.h>
9#include <measures/Measures/MeasConvert.h>
10#include <casa/Logging/LogMessage.h>
11#include <casa/BasicSL/Constants.h>
12
13#include "OldASDMReader.h"
14#include <atnf/PKSIO/SrcType.h>
15
16using namespace std ;
17//using namespace casa ;
18using namespace asdm ;
19using namespace sdmbin ;
20
21// sec to day
22double s2d = 1.0 / 86400.0 ;
23
24OldASDMReader::OldASDMReader()
25 : asdm_(0),
26 sdmBin_(0),
27 vmsData_(0),
28 antennaId_( -1 ),
29 antennaName_( "" ),
30 row_(-1),
31 apc_(AP_CORRECTED),
32 className_("OldASDMReader")
33{
34 configDescIdList_.resize(0) ;
35 feedIdList_.resize(0) ;
36 fieldIdList_.resize(0) ;
37 mainRow_.resize(0) ;
38 ifno_.clear() ;
39 corrMode_.reset() ;
40 timeSampling_.reset() ;
41}
42
43OldASDMReader::~OldASDMReader()
44{
45 close() ;
46 logsink_ = 0 ;
47}
48
49bool OldASDMReader::open( const string &filename, const casa::Record &rec )
50{
51 casa::String funcName = "open" ;
52
53 // return value
54 bool status = true ;
55
56 // set default
57 timeSampling_.reset() ;
58 corrMode_.reset() ;
59 resolutionType_.reset() ;
60 apc_ = AP_CORRECTED ;
61
62 // parsing ASDM options
63 if ( rec.isDefined( "asdm" ) ) {
64 casa::Record asdmrec = rec.asRecord( "asdm" ) ;
65
66 // antenna
67 if ( asdmrec.isDefined( "antenna" ) ) {
68 if ( asdmrec.type( asdmrec.fieldNumber( "antenna" ) ) == casa::TpInt ) {
69 antennaId_ = asdmrec.asInt( "antenna" ) ;
70 }
71 else {
72 antennaName_ = asdmrec.asString( "antenna" ) ;
73 }
74 }
75 else {
76 antennaId_ = 0 ;
77 }
78
79 // ATM phase correction
80 if ( asdmrec.isDefined( "apc" ) ) {
81 if ( asdmrec.asBool( "apc" ) )
82 apc_ = AP_CORRECTED ;
83 else
84 apc_ = AP_UNCORRECTED ;
85 }
86
87 // time sampling
88 String timeSampling = "all" ; // take all time sampling by default
89 if ( asdmrec.isDefined( "sampling" ) ) {
90 timeSampling = asdmrec.asString( "sampling" ) ;
91 }
92 if ( timeSampling == "all" ) {
93 timeSampling_.set( INTEGRATION ) ;
94 timeSampling_.set( SUBINTEGRATION ) ;
95 }
96 else if ( timeSampling == "integration" ) {
97 timeSampling_.set( INTEGRATION ) ;
98 }
99 else if ( timeSampling == "subintegration" ) {
100 timeSampling_.set( SUBINTEGRATION ) ;
101 }
102 else {
103 //throw AipsError( "Unrecognized option for sampling" ) ;
104 logsink_->postLocally( LogMessage( "Unrecognized option for time sampling: "+String::toString(timeSampling), LogOrigin(className_,funcName,WHERE), LogMessage::WARN ) ) ;
105 status = false ;
106 }
107
108 // spectral resolution type
109 string resolutionType = "all" ;
110 if ( asdmrec.isDefined( "srt" ) ) {
111 resolutionType = string( asdmrec.asString( "srt" ) ) ;
112 }
113 string resolutionTypes[3] ;
114 Int numType = split( resolutionType, resolutionTypes, 3, "," ) ;
115 for ( Int it = 0 ; it < numType ; it++ ) {
116 if ( resolutionTypes[it] == "all" ) {
117 resolutionType_.reset() ;
118 resolutionType_.set( FULL_RESOLUTION ) ;
119 resolutionType_.set( BASEBAND_WIDE ) ;
120 resolutionType_.set( CHANNEL_AVERAGE ) ;
121 break ;
122 }
123 else if ( resolutionTypes[it] == "fr" ) {
124 resolutionType_.set( FULL_RESOLUTION ) ;
125 }
126 else if ( resolutionTypes[it] == "bw" ) {
127 resolutionType_.set( BASEBAND_WIDE ) ;
128 }
129 else if ( resolutionTypes[it] == "ca" ) {
130 resolutionType_.set( CHANNEL_AVERAGE ) ;
131 }
132 }
133 if ( resolutionType_.size() == 0 ) {
134 logsink_->postLocally( LogMessage( "Unrecognized option for spectral resolution type: "+String::toString(resolutionType), LogOrigin(className_,funcName,WHERE), LogMessage::WARN ) ) ;
135 status = false ;
136 }
137
138 // input correlation mode
139 string corrMode = "ao,ca" ;
140 if ( asdmrec.isDefined( "corr" ) ) {
141 corrMode = string( asdmrec.asString( "corr" ) ) ;
142 //logsink_->postLocally( LogMessage("corrMode = "+String(corrMode),LogOrigin(className_,funcName,WHERE)) ) ;
143 }
144 string corrModes[3] ;
145 Int numCorr = split( corrMode, corrModes, 3, "," ) ;
146 for ( Int ic = 0 ; ic < numCorr ; ic++ ) {
147 if ( corrModes[ic] == "ao" ) {
148 corrMode_.set( AUTO_ONLY ) ;
149 }
150 else if ( corrModes[ic] == "ca" ) {
151 corrMode_.set( CROSS_AND_AUTO ) ;
152 }
153 }
154 if ( corrMode_.size() == 0 ) {
155 logsink_->postLocally( LogMessage( "Invalid option for correlation mode: "+String::toString(corrMode), LogOrigin(className_,funcName,WHERE), LogMessage::WARN ) ) ;
156 status = false ;
157 }
158
159// logsink_->postLocally( LogMessage( "### asdmrec summary ###", LogOrigin(className_,funcName,WHERE) ) ) ;
160// ostringstream oss ;
161// asdmrec.print( oss ) ;
162// logsink_->postLocally( LogMessage( oss.str(), LogOrigin(className_,funcName,WHERE) ) ) ;
163// logsink_->postLocally( LogMessage( "#######################", LogOrigin(className_,funcName,WHERE) ) ) ;
164 }
165
166 // create ASDM object
167 asdm_ = new ASDM() ;
168 asdm_->setFromFile( filename ) ;
169
170 if ( antennaId_ == -1 ) {
171 AntennaTable &atab = asdm_->getAntenna() ;
172 vector<AntennaRow *> rows = atab.get() ;
173 int idx = -1 ;
174 for ( casa::uInt irow = 0 ; irow < rows.size() ; irow++ ) {
175 if ( casa::String(rows[irow]->getName()) == antennaName_ ) {
176 idx = rows[irow]->getAntennaId().getTagValue() ;
177 break ;
178 }
179 }
180 if ( idx == -1 ) {
181 close() ;
182 throw (casa::AipsError( antennaName_ + " not found." )) ;
183 }
184 else {
185 antennaId_ = idx ;
186 }
187 }
188
189 // set antenna name
190 if ( antennaName_.size() == 0 ) {
191 AntennaTable &atab = asdm_->getAntenna() ;
192 Tag tag( antennaId_, TagType::Antenna ) ;
193 AntennaRow *arow = atab.getRowByKey( tag ) ;
194 if ( arow == 0 ) {
195 close() ;
196 throw (casa::AipsError( tag.toString() + " not found." )) ;
197 }
198 else {
199 antennaName_ = casa::String( arow->getName() ) ;
200 }
201 }
202
203 // create SDMBinData object
204 sdmBin_ = new SDMBinData( asdm_, filename ) ;
205
206 // get Main rows
207 //mainRow_ = casa::Vector<MainRow *>(asdm_->getMain().get()) ;
208
209 // set up IFNO
210 setupIFNO() ;
211
212 // process Station table
213 processStation() ;
214
215 logsink_->postLocally( LogMessage( "antennaId_ = "+String::toString(antennaId_), LogOrigin(className_,funcName,WHERE) ) ) ;
216 logsink_->postLocally( LogMessage( "antennaName_ = "+antennaName_, LogOrigin(className_,funcName,WHERE) ) ) ;
217
218 return true ;
219}
220
221// void OldASDMReader::fill()
222// {
223// }
224
225void OldASDMReader::close()
226{
227 clearMainRow() ;
228
229 if ( sdmBin_ )
230 delete sdmBin_ ;
231 sdmBin_ = 0 ;
232
233 if ( asdm_ )
234 delete asdm_ ;
235 asdm_ = 0 ;
236
237 return ;
238}
239
240void OldASDMReader::fillHeader( casa::Int &nchan,
241 casa::Int &npol,
242 casa::Int &nif,
243 casa::Int &nbeam,
244 casa::String &observer,
245 casa::String &project,
246 casa::String &obstype,
247 casa::String &antennaname,
248 casa::Vector<casa::Double> &antennaposition,
249 casa::Float &equinox,
250 casa::String &freqref,
251 casa::Double &reffreq,
252 casa::Double &bandwidth,
253 casa::Double &utc,
254 casa::String &fluxunit,
255 casa::String &epoch,
256 casa::String &poltype )
257{
258 casa::String funcName = "fillHeader" ;
259
260 ExecBlockTable &ebtab = asdm_->getExecBlock() ;
261 // at the moment take first row of ExecBlock table
262 ExecBlockRow *ebrow = ebtab.get()[0] ;
263 casa::String telescopeName( ebrow->getTelescopeName() ) ;
264 AntennaTable &atab = asdm_->getAntenna() ;
265 AntennaRow *arow = atab.getRowByKey( Tag( antennaId_, TagType::Antenna ) ) ;
266 //StationTable &stab = asdm_->getStation() ;
267 //StationRow *srow = stab.getRowByKey( arow->getStationId() ) ;
268 StationRow *srow = arow->getStationUsingStationId() ;
269 casa::String stationName( srow->getName() ) ;
270
271 // antennaname
272 // <telescopeName>//<antennaName>@stationName
273 antennaname = telescopeName + "//" + antennaName_ + "@" + stationName ;
274 //logsink_->postLocally( LogMessage("antennaName = "+antennaname,LogOrigin(className_,funcName,WHERE)) ) ;
275
276 // antennaposition
277 //vector<Length> antpos = arow->getPosition() ;
278 vector<Length> antpos = srow->getPosition() ;
279 antennaposition.resize( 3 ) ;
280 for ( casa::uInt i = 0 ; i < 3 ; i++ )
281 antennaposition[i] = casa::Double( antpos[i].get() ) ;
282
283 // observer
284 observer = ebrow->getObserverName() ;
285
286 // project
287 // T.B.D. (project UID?)
288 project = "T.B.D. (" + ebrow->getProjectId().toString() + ")" ;
289
290 // utc
291 // start time of the project
292 utc = casa::Double( ebrow->getStartTime().getMJD() ) ;
293
294
295 SpectralWindowTable &spwtab = asdm_->getSpectralWindow() ;
296 vector<SpectralWindowRow *> spwrows = spwtab.get() ;
297 int nspwrow = spwrows.size() ;
298
299 // nif
300 //nif = spwrows.size() ;
301 nif = getNumIFs() ;
302
303 // nchan
304 int refidx = -1 ;
305 vector<int> nchans ;
306 for ( int irow = 0 ; irow < nspwrow ; irow++ ) {
307 nchans.push_back( spwrows[irow]->getNumChan() ) ;
308 if ( refidx == -1 && nchans[irow] != 1 && nchans[irow] != 4 )
309 refidx = irow ;
310 }
311 nchan = casa::Int( *max_element( nchans.begin(), nchans.end() ) ) ;
312
313 //logsink_->postLocally( LogMessage("refidx = "+String::toString(refidx),LogOrigin(className_,funcName,WHERE)) ) ;
314
315 // bandwidth
316 vector<double> bws ;
317 for ( int irow = 0 ; irow < nspwrow ; irow++ ) {
318 if ( nchans[irow] != 4 ) { // exclude WVR data
319 bws.push_back( spwrows[irow]->getTotBandwidth().get() ) ;
320 }
321 }
322 bandwidth = casa::Double( *max_element( bws.begin(), bws.end() ) ) ;
323
324 // reffreq
325 reffreq = casa::Double( spwrows[refidx]->getRefFreq().get() ) ;
326
327 // freqref
328 if ( spwrows[refidx]->isMeasFreqRefExists() ) {
329 string mfr = CFrequencyReferenceCode::name( spwrows[refidx]->getMeasFreqRef() ) ;
330// if (mfr == "TOPO") {
331// freqref = "TOPOCENT";
332// } else if (mfr == "GEO") {
333// freqref = "GEOCENTR";
334// } else if (mfr == "BARY") {
335// freqref = "BARYCENT";
336// } else if (mfr == "GALACTO") {
337// freqref = "GALACTOC";
338// } else if (mfr == "LGROUP") {
339// freqref = "LOCALGRP";
340// } else if (mfr == "CMB") {
341// freqref = "CMBDIPOL";
342// } else if (mfr == "REST") {
343// freqref = "SOURCE";
344// }
345 freqref = String( mfr ) ;
346 }
347 else {
348 // frequency reference is TOPOCENT by default
349 //freqref = "TOPOCENT" ;
350 freqref = "TOPO" ;
351 }
352
353
354 PolarizationTable &ptab = asdm_->getPolarization() ;
355 vector<PolarizationRow *> prows = ptab.get() ;
356 vector<int> npols ;
357 refidx = -1 ;
358 for ( unsigned int irow = 0 ; irow < prows.size() ; irow++ ) {
359 npols.push_back( prows[irow]->getNumCorr() ) ;
360 if ( refidx == -1 && npols[irow] != 1 )
361 refidx = irow ;
362 }
363 if ( refidx == -1 )
364 refidx = 0 ;
365
366 // npol
367 npol = casa::Int( *max_element( npols.begin(), npols.end() ) ) ;
368
369 // poltype
370 vector<StokesParameter> corrType = prows[refidx]->getCorrType() ;
371 if ( corrType[0] == I ||
372 corrType[0] == Q ||
373 corrType[0] == U ||
374 corrType[0] == V )
375 poltype = "stokes" ;
376 else if ( corrType[0] == RR ||
377 corrType[0] == RL ||
378 corrType[0] == LR ||
379 corrType[0] == LL )
380 poltype = "circular" ;
381 else if ( corrType[0] == XX ||
382 corrType[0] == XY ||
383 corrType[0] == YX ||
384 corrType[0] == YY )
385 poltype = "linear" ;
386 else if ( corrType[0] == PLINEAR ||
387 corrType[0] == PANGLE ) {
388 poltype = "linpol" ;
389 }
390 else {
391 poltype = "notype" ;
392 }
393
394
395 FeedTable &ftab = asdm_->getFeed() ;
396 vector<FeedRow *> frows = ftab.get() ;
397 vector<int> nbeams ;
398 for ( unsigned int irow = 0 ; irow < frows.size() ; irow++ ) {
399 if ( frows[irow]->isFeedNumExists() )
400 nbeams.push_back( frows[irow]->getFeedNum() ) ;
401 else
402 nbeams.push_back( 1 ) ;
403 }
404
405 // nbeam
406 nbeam = casa::Int( *max_element( nbeams.begin(), nbeams.end() ) ) ;
407
408 // fluxunit
409 // tentatively set 'K'? or empty?
410 fluxunit = "K" ;
411
412 // equinox
413 // tentatively set 2000.0
414 equinox = 2000.0 ;
415
416 // epoch
417 // tentatively set "UTC"
418 epoch = "UTC" ;
419
420 // obstype
421 // at the moment take observingMode attribute in SBSummary table
422 SBSummaryRow *sbrow = ebrow->getSBSummaryUsingSBSummaryId() ;
423 vector<string> obsmode = sbrow->getObservingMode() ;
424 obstype = "" ;
425 for ( unsigned int imode = 0 ; imode < obsmode.size() ; imode++ ) {
426 obstype += casa::String(obsmode[imode]) ;
427 if ( imode != obsmode.size()-1 )
428 obstype += "#" ;
429 }
430}
431
432void OldASDMReader::selectConfigDescription()
433{
434 casa::String funcName = "selectConfigDescription" ;
435
436 vector<ConfigDescriptionRow *> cdrows = asdm_->getConfigDescription().get() ;
437 vector<Tag> cdidTags ;
438 for ( unsigned int irow = 0 ; irow < cdrows.size() ; irow++ ) {
439 //logsink_->postLocally( LogMessage("correlationMode["+String::toString(irow)+"] = "+String::toString(cdrows[irow]->getCorrelationMode()),LogOrigin(className_,funcName,WHERE)) ) ;
440 if ( cdrows[irow]->getCorrelationMode() != CROSS_ONLY ) {
441 cdidTags.push_back( cdrows[irow]->getConfigDescriptionId() ) ;
442 }
443 }
444
445 configDescIdList_.resize( cdidTags.size() ) ;
446 for ( unsigned int i = 0 ; i < cdidTags.size() ; i++ ) {
447 configDescIdList_[i] = casa::uInt( cdidTags[i].getTagValue() ) ;
448 }
449}
450
451void OldASDMReader::selectFeed()
452{
453 feedIdList_.resize(0) ;
454 vector<FeedRow *> frows = asdm_->getFeed().get() ;
455 Tag atag( antennaId_, TagType::Antenna ) ;
456 for ( unsigned int irow = 0 ; irow < frows.size() ; irow++ ) {
457 casa::uInt feedId = (casa::uInt)(frows[irow]->getFeedId() ) ;
458 if ( casa::anyEQ( feedIdList_, feedId ) )
459 continue ;
460 if ( frows[irow]->getAntennaId() == atag ) {
461 unsigned int oldsize = feedIdList_.size() ;
462 feedIdList_.resize( oldsize+1, true ) ;
463 feedIdList_[oldsize] = feedId ;
464 }
465 }
466}
467
468casa::Vector<casa::uInt> OldASDMReader::getFieldIdList()
469{
470 casa::String funcName = "getFieldIdList" ;
471
472 vector<FieldRow *> frows = asdm_->getField().get() ;
473 fieldIdList_.resize( frows.size() ) ;
474 for ( unsigned int irow = 0 ; irow < frows.size() ; irow++ ) {
475 //logsink_->postLocally( LogMessage("fieldId["+String::toString(irow)+"]="+String(frows[irow]->getFieldId().toString()),LogOrigin(className_,funcName,WHERE)) ) ;
476 fieldIdList_[irow] = frows[irow]->getFieldId().getTagValue() ;
477 }
478
479 return fieldIdList_ ;
480}
481
482casa::uInt OldASDMReader::getNumMainRow()
483{
484 casa::uInt nrow = casa::uInt( mainRow_.size() ) ;
485
486 return nrow ;
487}
488
489void OldASDMReader::select()
490{
491 // selection by input CorrelationMode
492 EnumSet<CorrelationMode> esCorrs ;
493 sdmBin_->select( corrMode_ ) ;
494
495 // selection by TimeSampling
496 sdmBin_->select( timeSampling_ ) ;
497
498 // selection by SpectralResolutionType
499 sdmBin_->select( resolutionType_ ) ;
500
501 // selection by AtmPhaseCorrection and output CorrelationMode
502 EnumSet<AtmPhaseCorrection> esApcs ;
503 esApcs.set( apc_ ) ;
504 // always take only autocorrelation data
505 Enum<CorrelationMode> esCorr = AUTO_ONLY ;
506 sdmBin_->selectDataSubset( esCorr, esApcs ) ;
507
508 // select valid configDescriptionId
509 selectConfigDescription() ;
510
511 // select valid feedId
512 selectFeed() ;
513}
514
515casa::Bool OldASDMReader::setMainRow( casa::uInt irow )
516{
517 casa::Bool status = true ;
518 row_ = irow ;
519
520 unsigned int cdid = mainRow_[row_]->getConfigDescriptionId().getTagValue() ;
521 if ( (int)count(configDescIdList_.begin(),configDescIdList_.end(),cdid) == 0 )
522 status = false ;
523 else {
524 status = (casa::Bool)(sdmBin_->acceptMainRow( mainRow_[row_] )) ;
525 }
526 return status ;
527}
528
529casa::Bool OldASDMReader::setMainRow( casa::uInt configDescId, casa::uInt fieldId )
530{
531 clearMainRow() ;
532
533 Tag configDescTag( (unsigned int)configDescId, TagType::ConfigDescription ) ;
534 Tag fieldTag( (unsigned int)fieldId, TagType::Field ) ;
535 mainRow_ = casa::Vector<MainRow *>( *(asdm_->getMain().getByContext( configDescTag, fieldTag ) ) ) ;
536
537 return true ;
538}
539
540void OldASDMReader::clearMainRow()
541{
542 mainRow_.resize(0) ;
543}
544
545void OldASDMReader::setupIFNO()
546{
547 casa::String funcName = "setupIFNO" ;
548
549 vector<SpectralWindowRow *> spwrows = asdm_->getSpectralWindow().get() ;
550 unsigned int nrow = spwrows.size() ;
551 ifno_.clear() ;
552 casa::uInt idx = 0 ;
553 casa::uInt wvridx = 0 ;
554 for ( unsigned int irow = 0 ; irow < nrow ; irow++ ) {
555 casa::uInt index ;
556 if ( isWVR( spwrows[irow] ) ) {
557 //logsink_->postLocally( LogMessage(spwrows[irow]->getSpectralWindowId().toString()+" is WVR",LogOrigin(className_,funcName,WHERE)) ) ;
558 index = wvridx ;
559 }
560 else {
561 index = ++idx ;
562 }
563 ifno_.insert( pair<Tag,casa::uInt>(spwrows[irow]->getSpectralWindowId(),index) ) ;
564 //logsink_->postLocally( LogMessage(spwrows[irow]->getSpectralWindowId().toString()+": IFNO="+String::toString(index),LogOrigin(className_,funcName,WHERE)) ) ;
565 }
566}
567
568bool OldASDMReader::isWVR( SpectralWindowRow *row )
569{
570 BasebandName bbname = row->getBasebandName() ;
571 int nchan = row->getNumChan() ;
572 if ( bbname == NOBB && nchan == 4 )
573 return true ;
574 else
575 return false ;
576}
577
578// casa::Vector<casa::uInt> OldASDMReader::getDataDescIdList( casa::uInt cdid )
579// {
580// Tag cdTag( (unsigned int)cdid, TagType::ConfigDescription ) ;
581// ConfigDescriptionRow *cdrow = asdm_->getConfigDescription().getRowByKey( cdTag ) ;
582// vector<Tag> ddTags = cdrow->getDataDescriptionId() ;
583// casa::Vector<casa::uInt> ddidList( ddTags.size() ) ;
584// for ( unsigned int idd = 0 ; idd < ddTags.size() ; idd++ ) {
585// ddidList[idd] = ddTags[idd].getTagValue() ;
586// }
587// return ddidList ;
588// }
589
590// casa::Vector<casa::uInt> OldASDMReader::getSwitchCycleIdList( casa::uInt cdid )
591// {
592// Tag cdTag( (unsigned int)cdid, TagType::ConfigDescription ) ;
593// ConfigDescriptionRow *cdrow = asdm_->getConfigDescription().getRowByKey( cdTag ) ;
594// vector<Tag> scTags = cdrow->getSwitchCycleId() ;
595// casa::Vector<casa::uInt> scidList( scTags.size() ) ;
596// for ( unsigned int idd = 0 ; idd < scTags.size() ; idd++ ) {
597// scidList[idd] = scTags[idd].getTagValue() ;
598// }
599// return scidList ;
600// }
601
602// casa::Vector<casa::uInt> OldASDMReader::getFeedIdList( casa::uInt cdid )
603// {
604// casa::String funcName = "getFeedIdList" ;
605//
606// Tag cdTag( (unsigned int)cdid, TagType::ConfigDescription ) ;
607// ConfigDescriptionRow *cdrow = asdm_->getConfigDescription().getRowByKey( cdTag ) ;
608// casa::Vector<casa::uInt> feedIdList ;
609// vector<int> feedIds = cdrow->getFeedId() ;
610// for ( unsigned int ife = 0 ; ife < feedIds.size() ; ife++ ) {
611// logsink_->postLocally( LogMessage("feedIds["+String::toString(ife)+"]="+String::toString(feedIds[ife]),LogOrigin(className_,funcName,WHERE)) ) ;
612// if ( casa::anyEQ( feedIdList, casa::uInt( feedIds[ife] ) ) )
613// continue ;
614// if ( casa::anyEQ( feedIdList_, casa::uInt( feedIds[ife] ) ) ) {
615// casa::uInt oldsize = feedIdList.size() ;
616// feedIdList.resize( oldsize+1, true ) ;
617// feedIdList[oldsize] = casa::uInt( feedIds[ife] ) ;
618// }
619// }
620// logsink_->postLocally( LogMessage("feedIdList.size() = "+String::toString(feedIdList.size()),LogOrigin(className_,funcName,WHERE)) ) ;
621// return feedIdList ;
622// }
623
624casa::Bool OldASDMReader::setData()
625{
626 casa::String funcName = "setData" ;
627
628 //logsink_->postLocally( LogMessage("try to retrieve binary data",LogOrigin(className_,funcName,WHERE)) ) ;
629
630// EnumSet<AtmPhaseCorrection> esApcs ;
631// esApcs.set( apc_ ) ;
632// // always take only autocorrelation data
633// Enum<CorrelationMode> esCorr = AUTO_ONLY ;
634// vmsData_ = sdmBin_->getDataCols( esCorr, esApcs ) ;
635
636 // 2011/07/06 TN
637 // Workaround to avoid unwanted message from SDMBinData::getDataCols()
638 ostringstream oss ;
639 streambuf *buforg = cout.rdbuf(oss.rdbuf()) ;
640 vmsData_ = sdmBin_->getDataCols() ;
641 cout.rdbuf(buforg) ;
642
643// logsink_->postLocally( LogMessage("oss.str() = "+oss.str(),LogOrigin(className_,funcName,WHERE)) ) ;
644// cout << "This is test: oss.str()=" << oss.str() << endl ;
645
646
647 //logsink_->postLocally( LogMessage("processorId = "+String::toString(vmsData_->processorId),LogOrigin(className_,funcName,WHERE)) ) ;
648 //logsink_->postLocally( LogMessage("v_time.size() = "+String::toString(vmsData_->v_time.size()),LogOrigin(className_,funcName,WHERE)) ) ;
649 //logsink_->postLocally( LogMessage(" v_time[0] = "+String::toString(vmsData_->v_time[0]),LogOrigin(className_,funcName,WHERE)) ) ;
650 //logsink_->postLocally( LogMessage("v_interval.size() = "+String::toString(vmsData_->v_interval.size()),LogOrigin(className_,funcName,WHERE)) ) ;
651 //logsink_->postLocally( LogMessage(" v_interval[0] = "+String::toString(vmsData_->v_interval[0]),LogOrigin(className_,funcName,WHERE)) ) ;
652 //logsink_->postLocally( LogMessage("v_atmPhaseCorrection.size() = "+String::toString(vmsData_->v_atmPhaseCorrection.size()),LogOrigin(className_,funcName,WHERE)) ) ;
653 //logsink_->postLocally( LogMessage("binNum = "+String::toString(vmsData_->binNum),LogOrigin(className_,funcName,WHERE)) ) ;
654 //logsink_->postLocally( LogMessage("v_projectPath.size() = "+String::toString(vmsData_->v_projectPath.size()),LogOrigin(className_,funcName,WHERE)) ) ;
655 //logsink_->postLocally( LogMessage("v_antennaId1.size() = "+String::toString(vmsData_->v_antennaId1.size()),LogOrigin(className_,funcName,WHERE)) ) ;
656 //logsink_->postLocally( LogMessage("v_antennaId2.size() = "+String::toString(vmsData_->v_antennaId2.size()),LogOrigin(className_,funcName,WHERE)) ) ;
657 //logsink_->postLocally( LogMessage("v_feedId1.size() = "+String::toString(vmsData_->v_feedId1.size()),LogOrigin(className_,funcName,WHERE)) ) ;
658 //logsink_->postLocally( LogMessage("v_feedId2.size() = "+String::toString(vmsData_->v_feedId2.size()),LogOrigin(className_,funcName,WHERE)) ) ;
659 //logsink_->postLocally( LogMessage("v_dataDescId.size() = "+String::toString(vmsData_->v_dataDescId.size()),LogOrigin(className_,funcName,WHERE)) ) ;
660 //logsink_->postLocally( LogMessage("v_timeCentroid.size() = "+String::toString(vmsData_->v_timeCentroid.size()),LogOrigin(className_,funcName,WHERE)) ) ;
661 //logsink_->postLocally( LogMessage("v_exposure.size() = "+String::toString(vmsData_->v_exposure.size()),LogOrigin(className_,funcName,WHERE)) ) ;
662 //logsink_->postLocally( LogMessage("v_numData.size() = "+String::toString(vmsData_->v_numData.size()),LogOrigin(className_,funcName,WHERE)) ) ;
663 //logsink_->postLocally( LogMessage("vv_dataShape.size() = "+String::toString(vmsData_->vv_dataShape.size()),LogOrigin(className_,funcName,WHERE)) ) ;
664 //logsink_->postLocally( LogMessage("v_m_data.size() = "+String::toString(vmsData_->v_m_data.size()),LogOrigin(className_,funcName,WHERE)) ) ;
665 //logsink_->postLocally( LogMessage("v_phaseDir.size() = "+String::toString(vmsData_->v_phaseDir.size()),LogOrigin(className_,funcName,WHERE)) ) ;
666 //logsink_->postLocally( LogMessage("v_stateId.size() = "+String::toString(vmsData_->v_stateId.size()),LogOrigin(className_,funcName,WHERE)) ) ;
667 //logsink_->postLocally( LogMessage("v_msState.size() = "+String::toString(vmsData_->v_msState.size()),LogOrigin(className_,funcName,WHERE)) ) ;
668 //logsink_->postLocally( LogMessage("v_flag.size() = "+String::toString(vmsData_->v_flag.size()),LogOrigin(className_,funcName,WHERE)) ) ;
669
670 dataIdList_.clear() ;
671 unsigned int numTotalData = vmsData_->v_m_data.size() ;
672 for ( unsigned int idata = 0 ; idata < numTotalData ; idata++ ) {
673 if ( vmsData_->v_antennaId1[idata] == (int)antennaId_
674 && vmsData_->v_antennaId2[idata] == (int)antennaId_ )
675 dataIdList_.push_back( idata ) ;
676 }
677 numData_ = dataIdList_.size() ;
678 //logsink_->postLocally( LogMessage("numData_ = "+String::toString(numData_),LogOrigin(className_,funcName,WHERE)) ) ;
679 //logsink_->postLocally( LogMessage("dataSize = "+String::toString(mainRow_[row_]->getDataSize()),LogOrigin(className_,funcName,WHERE)) ) ;
680
681 return true ;
682}
683
684casa::uInt OldASDMReader::getIFNo( unsigned int idx )
685{
686 Tag ddTag( vmsData_->v_dataDescId[dataIdList_[idx]], TagType::DataDescription ) ;
687 DataDescriptionRow *ddrow = asdm_->getDataDescription().getRowByKey( ddTag ) ;
688 Tag spwid = ddrow->getSpectralWindowId() ;
689 map<Tag,casa::uInt>::iterator iter = ifno_.find( spwid ) ;
690 if ( iter != ifno_.end() )
691 return iter->second ;
692 else {
693 return 0 ;
694 }
695}
696
697int OldASDMReader::getNumPol( unsigned int idx )
698{
699 Tag ddTag( vmsData_->v_dataDescId[dataIdList_[idx]], TagType::DataDescription ) ;
700 DataDescriptionRow *ddrow = asdm_->getDataDescription().getRowByKey( ddTag ) ;
701 PolarizationRow *polrow = ddrow->getPolarizationUsingPolOrHoloId() ;
702 return polrow->getNumCorr() ;
703}
704
705void OldASDMReader::getFrequency( unsigned int idx,
706 double &refpix,
707 double &refval,
708 double &incr,
709 string &freqref )
710{
711 casa::String funcName = "getFrequency" ;
712
713 Tag ddTag( vmsData_->v_dataDescId[dataIdList_[idx]], TagType::DataDescription ) ;
714 DataDescriptionRow *ddrow = asdm_->getDataDescription().getRowByKey( ddTag ) ;
715 //Tag spwid = ddrow->getSpectralWindowId() ;
716 SpectralWindowRow *spwrow = ddrow->getSpectralWindowUsingSpectralWindowId() ;
717 int nchan = spwrow->getNumChan() ;
718 freqref = "TOPO" ;
719 if ( spwrow->isMeasFreqRefExists() )
720 freqref = CFrequencyReferenceCode::toString( spwrow->getMeasFreqRef() ) ;
721 if ( nchan == 1 ) {
722 //logsink_->postLocally( LogMessage("channel averaged data",LogOrigin(className_,funcName,WHERE)) ) ;
723 refpix = 0.0 ;
724 incr = spwrow->getTotBandwidth().get() ;
725 if ( spwrow->isChanFreqStartExists() ) {
726 refval = spwrow->getChanFreqStart().get() ;
727 }
728 else if ( spwrow->isChanFreqArrayExists() ) {
729 refval = spwrow->getChanFreqArray()[0].get() ;
730 }
731 else {
732 throw (casa::AipsError( "Either chanFreqArray or chanFreqStart must exist." )) ;
733 }
734 }
735 else if ( nchan % 2 ) {
736 // odd
737 //logsink_->postLocally( LogMessage("odd case",LogOrigin(className_,funcName,WHERE)) ) ;
738 refpix = 0.5 * ( (double)nchan - 1.0 ) ;
739 int ic = ( nchan - 1 ) / 2 ;
740 if ( spwrow->isChanWidthExists() ) {
741 incr = spwrow->getChanWidth().get() ;
742 }
743 else if ( spwrow->isChanWidthArrayExists() ) {
744 incr = spwrow->getChanWidthArray()[0].get() ;
745 }
746 else {
747 throw (casa::AipsError( "Either chanWidthArray or chanWidth must exist." )) ;
748 }
749 if ( spwrow->isChanFreqStepExists() ) {
750 if ( spwrow->getChanFreqStep().get() < 0.0 )
751 incr *= -1.0 ;
752 }
753 else if ( spwrow->isChanFreqArrayExists() ) {
754 vector<Frequency> chanFreqArr = spwrow->getChanFreqArray() ;
755 if ( chanFreqArr[0].get() > chanFreqArr[1].get() )
756 incr *= -1.0 ;
757 }
758 else {
759 throw (casa::AipsError( "Either chanFreqArray or chanFreqStep must exist." )) ;
760 }
761 if ( spwrow->isChanFreqStartExists() ) {
762 refval = spwrow->getChanFreqStart().get() + refpix * incr ;
763 }
764 else if ( spwrow->isChanFreqArrayExists() ) {
765 refval = spwrow->getChanFreqArray()[ic].get() ;
766 }
767 else {
768 throw (casa::AipsError( "Either chanFreqArray or chanFreqStart must exist." )) ;
769 }
770 }
771 else {
772 // even
773 //logsink_->postLocally( LogMessage("even case",LogOrigin(className_,funcName,WHERE)) ) ;
774 refpix = 0.5 * ( (double)nchan - 1.0 ) ;
775 int ic = nchan / 2 ;
776 if ( spwrow->isChanWidthExists() ) {
777 incr = spwrow->getChanWidth().get() ;
778 }
779 else if ( spwrow->isChanWidthArrayExists() ) {
780 incr = spwrow->getChanWidthArray()[0].get() ;
781 }
782 else {
783 throw (casa::AipsError( "Either chanWidthArray or chanWidth must exist." )) ;
784 }
785 if ( spwrow->isChanFreqStepExists() ) {
786 if ( spwrow->getChanFreqStep().get() < 0.0 )
787 incr *= -1.0 ;
788 }
789 else if ( spwrow->isChanFreqArrayExists() ) {
790 vector<Frequency> chanFreqArr = spwrow->getChanFreqArray() ;
791 if ( chanFreqArr[0].get() > chanFreqArr[1].get() )
792 incr *= -1.0 ;
793 }
794 else {
795 throw (casa::AipsError( "Either chanFreqArray or chanFreqStep must exist." )) ;
796 }
797 if ( spwrow->isChanFreqStartExists() ) {
798 refval = spwrow->getChanFreqStart().get() + refpix * incr ;
799 }
800 else if ( spwrow->isChanFreqArrayExists() ) {
801 vector<Frequency> freqs = spwrow->getChanFreqArray() ;
802 refval = 0.5 * ( freqs[ic-1].get() + freqs[ic].get() ) ;
803 }
804 else {
805 throw (casa::AipsError( "Either chanFreqArray or chanFreqStart must exist." )) ;
806 }
807 }
808}
809
810vector<double> OldASDMReader::getRestFrequency( unsigned int idx )
811{
812 vector<double> rf( 0 ) ;
813 unsigned int index = dataIdList_[idx] ;
814 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
815 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
816 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
817 Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
818 Tag spwtag = asdm_->getDataDescription().getRowByKey(ddtag)->getSpectralWindowId() ;
819 Tag ftag( vmsData_->v_fieldId[index], TagType::Field ) ;
820 FieldRow *frow = asdm_->getField().getRowByKey( ftag ) ;
821 if ( frow->isSourceIdExists() ) {
822 //logsink_->postLocally( LogMessage("sourceId exists",LogOrigin(className_,funcName,WHERE)) ) ;
823 int sid = frow->getSourceId() ;
824 SourceRow *srow = asdm_->getSource().getRowByKey( sid, tint, spwtag ) ;
825 if ( srow->isRestFrequencyExists() ) {
826 //logsink_->postLocally( LogMessage("restFrequency exists",LogOrigin(className_,funcName,WHERE)) ) ;
827 vector<Frequency> restfreq = srow->getRestFrequency() ;
828 rf.resize( restfreq.size() ) ;
829 for ( unsigned int i = 0 ; i < restfreq.size() ; i++ )
830 rf[i] = restfreq[i].get() ;
831 }
832 }
833 return rf ;
834}
835
836double OldASDMReader::getTime( unsigned int idx )
837{
838 double tsec = vmsData_->v_time[dataIdList_[idx]] ;
839 return tsec * s2d ;
840}
841
842double OldASDMReader::getInterval( unsigned int idx )
843{
844 return vmsData_->v_interval[dataIdList_[idx]] ;
845}
846
847string OldASDMReader::getSourceName( unsigned int idx )
848{
849 unsigned int index = dataIdList_[idx] ;
850 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
851 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
852 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
853 Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
854 Tag spwtag = asdm_->getDataDescription().getRowByKey(ddtag)->getSpectralWindowId() ;
855 Tag ftag( vmsData_->v_fieldId[index], TagType::Field ) ;
856 FieldRow *frow = asdm_->getField().getRowByKey( ftag ) ;
857 string srcname ;
858 if ( frow->isSourceIdExists() ) {
859 //logsink_->postLocally( LogMessage("sourceId exists",LogOrigin(className_,funcName,WHERE)) ) ;
860 int sid = frow->getSourceId() ;
861 SourceRow *srow = asdm_->getSource().getRowByKey( sid, tint, spwtag ) ;
862 srcname = srow->getSourceName() ;
863 }
864 else {
865 srcname = frow->getFieldName() ;
866 }
867 return srcname ;
868}
869
870string OldASDMReader::getFieldName( unsigned int idx )
871{
872 int fid = vmsData_->v_fieldId[dataIdList_[idx]] ;
873 Tag ftag( fid, TagType::Field ) ;
874 FieldRow *frow = asdm_->getField().getRowByKey( ftag ) ;
875 ostringstream oss ;
876 oss << frow->getFieldName() << "__" << fid ;
877 return oss.str() ;
878}
879
880int OldASDMReader::getSrcType( unsigned int scan,
881 unsigned int subscan )
882{
883 int srctype = SrcType::NOTYPE ;
884 Tag ebtag = mainRow_[row_]->getExecBlockId() ;
885 ScanRow *scanrow = asdm_->getScan().getRowByKey( ebtag, (int)scan ) ;
886 ScanIntent scanIntent = scanrow->getScanIntent()[0] ;
887 SubscanRow *subrow = asdm_->getSubscan().getRowByKey( ebtag, (int)scan, (int)subscan ) ;
888 SubscanIntent subIntent = subrow->getSubscanIntent() ;
889 SwitchingMode swmode = subrow->getSubscanMode() ;
890 if ( scanIntent == OBSERVE_TARGET ) {
891 // on sky scan
892 if ( swmode == NO_SWITCHING || swmode == POSITION_SWITCHING ) {
893 // position switching
894 // tentatively set NO_SWITCHING = POSITION_SWITCHING
895 if ( subIntent == ON_SOURCE ) {
896 srctype = SrcType::PSON ;
897 }
898 else if ( subIntent == OFF_SOURCE ) {
899 srctype = SrcType::PSOFF ;
900 }
901 }
902 else if ( swmode == FREQUENCY_SWITCHING ) {
903 // frequency switching
904 if ( subIntent == ON_SOURCE ) {
905 srctype = SrcType::FSON ;
906 }
907 else if ( subIntent == OFF_SOURCE ) {
908 srctype = SrcType::FSOFF ;
909 }
910 }
911 else if ( swmode == NUTATOR_SWITCHING ) {
912 // nutator switching
913 if ( subIntent == ON_SOURCE ) {
914 srctype = SrcType::PSON ;
915 }
916 else if ( subIntent == OFF_SOURCE ) {
917 srctype = SrcType::PSOFF ;
918 }
919 }
920 else {
921 // other switching mode
922 if ( subIntent == ON_SOURCE ) {
923 srctype = SrcType::SIG ;
924 }
925 else if ( subIntent == OFF_SOURCE ) {
926 srctype = SrcType::REF ;
927 }
928 }
929 }
930 else if ( scanIntent == CALIBRATE_ATMOSPHERE ) {
931 if ( swmode == NO_SWITCHING || swmode == POSITION_SWITCHING ) {
932 // position switching
933 // tentatively set NO_SWITCHING = POSITION_SWITCHING
934 if ( subIntent == ON_SOURCE ) {
935 srctype = SrcType::PONCAL ;
936 }
937 else if ( subIntent == OFF_SOURCE ) {
938 srctype = SrcType::POFFCAL ;
939 }
940 }
941 else if ( swmode == FREQUENCY_SWITCHING ) {
942 // frequency switching
943 if ( subIntent == ON_SOURCE ) {
944 srctype = SrcType::FONCAL ;
945 }
946 else if ( subIntent == OFF_SOURCE ) {
947 srctype = SrcType::FOFFCAL ;
948 }
949 }
950 else if ( swmode == NUTATOR_SWITCHING ) {
951 // nutator switching
952 if ( subIntent == ON_SOURCE ) {
953 srctype = SrcType::PONCAL ;
954 }
955 else if ( subIntent == OFF_SOURCE ) {
956 srctype = SrcType::POFFCAL ;
957 }
958 }
959 else {
960 // other switching mode
961 if ( subIntent == ON_SOURCE ) {
962 srctype = SrcType::CAL ;
963 }
964 else if ( subIntent == OFF_SOURCE ) {
965 srctype = SrcType::CAL ;
966 }
967 }
968 }
969 else {
970 // off sky (e.g. calibrator device) scan
971 if ( subIntent == ON_SOURCE ) {
972 srctype = SrcType::SIG ;
973 }
974 else if ( subIntent == OFF_SOURCE ) {
975 srctype = SrcType::REF ;
976 }
977 else if ( subIntent == HOT ) {
978 srctype = SrcType::HOT ;
979 }
980 else if ( subIntent == AMBIENT ) {
981 srctype = SrcType::SKY ;
982 }
983 else {
984 srctype = SrcType::CAL ;
985 }
986 }
987
988 return srctype ;
989}
990
991unsigned int OldASDMReader::getSubscanNo( unsigned int idx )
992{
993 //logsink_->postLocally( LogMessage("subscan"+String::toString(vmsData_->v_msState[dataIdList_[idx]].subscanNum)+": obsmode="+String::toString(vmsData_->v_msState[dataIdList_[idx]].obsMode),LogOrigin(className_,funcName,WHERE)) ) ;
994 return vmsData_->v_msState[dataIdList_[idx]].subscanNum ;
995}
996
997vector<double> OldASDMReader::getSourceDirection( unsigned int idx )
998{
999 vector<double> dir( 2, 0.0 ) ;
1000 unsigned int index = dataIdList_[idx] ;
1001 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1002 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1003 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1004 Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
1005 Tag spwtag = asdm_->getDataDescription().getRowByKey(ddtag)->getSpectralWindowId() ;
1006 Tag ftag( vmsData_->v_fieldId[index], TagType::Field ) ;
1007 FieldRow *frow = asdm_->getField().getRowByKey( ftag ) ;
1008 string srcname ;
1009 if ( frow->isSourceIdExists() ) {
1010 //logsink_->postLocally( LogMessage("sourceId exists",LogOrigin(className_,funcName,WHERE)) ) ;
1011 int sid = frow->getSourceId() ;
1012 SourceRow *srow = asdm_->getSource().getRowByKey( sid, tint, spwtag ) ;
1013 vector<Angle> srcdir = srow->getDirection() ;
1014 dir[0] = limitedAngle( srcdir[0].get() ) ;
1015 dir[1] = limitedAngle( srcdir[1].get() ) ;
1016 if ( srow->isDirectionCodeExists() ) {
1017 DirectionReferenceCode dircode = srow->getDirectionCode() ;
1018 //logsink_->postLocally( LogMessage("dircode="+CDirectionReferenceCode::toString(dircode),LogOrigin(className_,funcName,WHERE)) ) ;
1019 if ( dircode != J2000 ) {
1020 // if not J2000, need direction conversion
1021 string ref = CDirectionReferenceCode::toString( dircode ) ;
1022 double mjd = vmsData_->v_time[index] * s2d ;
1023 Tag atag( antennaId_, TagType::Antenna ) ;
1024 AntennaRow *arow = asdm_->getAntenna().getRowByKey( atag ) ;
1025 StationRow *srow = arow->getStationUsingStationId() ;
1026 vector<Length> antposL = srow->getPosition() ;
1027 casa::Vector<casa::Double> antpos( 3 ) ;
1028 for ( int i = 0 ; i < 3 ; i++ )
1029 antpos[i] = antposL[i].get() ;
1030 dir = toJ2000( dir, ref, mjd, antpos ) ;
1031 }
1032 }
1033 }
1034 return dir ;
1035}
1036
1037void OldASDMReader::getSourceDirection( unsigned int idx,
1038 vector<double> &dir,
1039 string &ref )
1040{
1041 dir.resize( 2 ) ;
1042 dir[0] = 0.0 ;
1043 dir[1] = 0.0 ;
1044 ref = "J2000" ;
1045 unsigned int index = dataIdList_[idx] ;
1046 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1047 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1048 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1049 Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
1050 Tag spwtag = asdm_->getDataDescription().getRowByKey(ddtag)->getSpectralWindowId() ;
1051 Tag ftag( vmsData_->v_fieldId[index], TagType::Field ) ;
1052 FieldRow *frow = asdm_->getField().getRowByKey( ftag ) ;
1053 string srcname ;
1054 if ( frow->isSourceIdExists() ) {
1055 //logsink_->postLocally( LogMessage("sourceId exists",LogOrigin(className_,funcName,WHERE)) ) ;
1056 int sid = frow->getSourceId() ;
1057 SourceRow *srow = asdm_->getSource().getRowByKey( sid, tint, spwtag ) ;
1058 vector<Angle> srcdir = srow->getDirection() ;
1059 if ( srow->isDirectionCodeExists() ) {
1060 ref = CDirectionReferenceCode::toString( srow->getDirectionCode() ) ;
1061 }
1062 dir[0] = limitedAngle( srcdir[0].get() ) ;
1063 dir[1] = limitedAngle( srcdir[1].get() ) ;
1064 }
1065}
1066
1067void OldASDMReader::getSourceDirection( vector<double> &dir, string &ref )
1068{
1069 dir.resize( 2 ) ;
1070 ref = "J2000" ; // default is J2000
1071 SourceTable &tab = asdm_->getSource() ;
1072 SourceRow *row = tab.get()[0] ;
1073 vector<Angle> dirA = row->getDirection() ;
1074 dir[0] = dirA[0].get() ;
1075 dir[1] = dirA[1].get() ;
1076 if ( row->isDirectionCodeExists() ) {
1077 ref = CDirectionReferenceCode::toString( row->getDirectionCode() ) ;
1078 }
1079}
1080
1081vector<double> OldASDMReader::getSourceProperMotion( unsigned int idx )
1082{
1083 vector<double> pm( 2, 0.0 ) ;
1084 unsigned int index = dataIdList_[idx] ;
1085 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1086 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1087 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1088 Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
1089 Tag spwtag = asdm_->getDataDescription().getRowByKey(ddtag)->getSpectralWindowId() ;
1090 Tag ftag( vmsData_->v_fieldId[index], TagType::Field ) ;
1091 FieldRow *frow = asdm_->getField().getRowByKey( ftag ) ;
1092 string srcname ;
1093 if ( frow->isSourceIdExists() ) {
1094 //logsink_->postLocally( LogMessage("sourceId exists",LogOrigin(className_,funcName,WHERE)) ) ;
1095 int sid = frow->getSourceId() ;
1096 SourceRow *srow = asdm_->getSource().getRowByKey( sid, tint, spwtag ) ;
1097 vector<AngularRate> srcpm = srow->getProperMotion() ;
1098 pm[0] = srcpm[0].get() ;
1099 pm[1] = srcpm[1].get() ;
1100 }
1101 return pm ;
1102}
1103
1104double OldASDMReader::getSysVel( unsigned int idx )
1105{
1106 double sysvel = 0.0 ;
1107 unsigned int index = dataIdList_[idx] ;
1108 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1109 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1110 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1111 Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
1112 Tag spwtag = asdm_->getDataDescription().getRowByKey(ddtag)->getSpectralWindowId() ;
1113 Tag ftag( vmsData_->v_fieldId[index], TagType::Field ) ;
1114 FieldRow *frow = asdm_->getField().getRowByKey( ftag ) ;
1115 string srcname ;
1116 if ( frow->isSourceIdExists() ) {
1117 //logsink_->postLocally( LogMessage("sourceId exists",LogOrigin(className_,funcName,WHERE)) ) ;
1118 int sid = frow->getSourceId() ;
1119 SourceRow *srow = asdm_->getSource().getRowByKey( sid, tint, spwtag ) ;
1120 if ( srow->isSysVelExists() ) {
1121 vector<Speed> sysvelV = srow->getSysVel() ;
1122 if ( sysvelV.size() > 0 )
1123 sysvel = sysvelV[0].get() ;
1124 }
1125 }
1126 return sysvel ;
1127}
1128
1129unsigned int OldASDMReader::getFlagRow( unsigned int idx )
1130{
1131 return vmsData_->v_flag[dataIdList_[idx]] ;
1132}
1133
1134vector<unsigned int> OldASDMReader::getDataShape( unsigned int idx )
1135{
1136 return vmsData_->vv_dataShape[dataIdList_[idx]] ;
1137}
1138
1139float * OldASDMReader::getSpectrum( unsigned int idx )
1140{
1141 map<AtmPhaseCorrection, float*> data = vmsData_->v_m_data[dataIdList_[idx]] ;
1142 //map<AtmPhaseCorrection, float*>::iterator iter = data.find(AP_UNCORRECTED) ;
1143 map<AtmPhaseCorrection, float*>::iterator iter = data.find(apc_) ;
1144 float *autoCorr = iter->second ;
1145 return autoCorr ;
1146}
1147
1148// bool * OldASDMReader::getFlagChannel( unsigned int idx )
1149// {
1150// return 0 ;
1151// }
1152
1153vector< vector<float> > OldASDMReader::getTsys( unsigned int idx )
1154{
1155 vector< vector<float> > defaultTsys( 1, vector<float>( 1, 1.0 ) ) ;
1156 SysCalRow *scrow = getSysCalRow( idx ) ;
1157 if ( scrow != 0 && scrow->isTsysSpectrumExists() ) {
1158 vector< vector<Temperature> > tsysSpec = scrow->getTsysSpectrum() ;
1159 unsigned int numReceptor = tsysSpec.size() ;
1160 unsigned int numChan = tsysSpec[0].size() ;
1161 vector< vector<float> > tsys( numReceptor, vector<float>( numChan, 1.0 ) ) ;
1162 for ( unsigned int ir = 0 ; ir < numReceptor ; ir++ ) {
1163 for ( unsigned int ic = 0 ; ic < numChan ; ic++ ) {
1164 tsys[ir][ic] = tsysSpec[ir][ic].get() ;
1165 }
1166 }
1167 return tsys ;
1168 }
1169// else if ( scrow->isTsysExists() ) {
1170// vector<Temperature> tsysScalar = scrow->getTsys() ;
1171// unsigned int numReceptor = tsysScalar.size() ;
1172// vector< vector<float> > tsys( numReceptor ) ;
1173// for ( unsigned int ir = 0 ; ir < numReceptor ; ir++ )
1174// tsys[ir] = vector<float>( 1, tsysScalar[ir].get() ) ;
1175// return tsys ;
1176// }
1177 else {
1178 return defaultTsys ;
1179 }
1180}
1181
1182vector< vector<float> > OldASDMReader::getTcal( unsigned int idx )
1183{
1184 vector< vector<float> > defaultTcal( 1, vector<float>( 1, 1.0 ) ) ;
1185 SysCalRow *scrow = getSysCalRow( idx ) ;
1186 if ( scrow != 0 && scrow->isTcalSpectrumExists() ) {
1187 vector< vector<Temperature> > tcalSpec = scrow->getTcalSpectrum() ;
1188 unsigned int numReceptor = tcalSpec.size() ;
1189 unsigned int numChan = tcalSpec[0].size() ;
1190 vector< vector<float> > tcal( numReceptor, vector<float>( numChan, 1.0 ) ) ;
1191 for ( unsigned int ir = 0 ; ir < numReceptor ; ir++ ) {
1192 for ( unsigned int ic = 0 ; ic < numChan ; ic++ ) {
1193 tcal[ir][ic] = tcalSpec[ir][ic].get() ;
1194 }
1195 }
1196 return tcal ;
1197 }
1198// else if ( scrow->isTcalExists() ) {
1199// vector<Temperature> tcalScalar = scrow->getTcal() ;
1200// unsigned int numReceptor = tcalScalar.size() ;
1201// vector< vector<float> > tcal( numReceptor, vector<float>( 1, 1.0 ) ) ;
1202// for ( unsigned int ir = 0 ; ir < numReceptor ; ir++ )
1203// tcal[ir][0] = tcalScalar[ir][0].get() ;
1204// return tcal ;
1205// }
1206 else {
1207 return defaultTcal ;
1208 }
1209}
1210
1211void OldASDMReader::getTcalAndTsys( unsigned int idx,
1212 vector< vector<float> > &tcal,
1213 vector< vector<float> > &tsys )
1214{
1215 String funcName = "getTcalAndTsys" ;
1216
1217 vector< vector<float> > defaultT( 1, vector<float>( 1, 1.0 ) ) ;
1218 SysCalRow *scrow = getSysCalRow( idx ) ;
1219 if ( scrow == 0 ) {
1220 tcal = defaultT ;
1221 tsys = defaultT ;
1222 }
1223 else {
1224 if ( scrow->isTsysSpectrumExists() ) {
1225 vector< vector<Temperature> > tsysSpec = scrow->getTsysSpectrum() ;
1226 unsigned int numReceptor = tsysSpec.size() ;
1227 unsigned int numChan = tsysSpec[0].size() ;
1228 tsys.resize( numReceptor ) ;
1229 for ( unsigned int ir = 0 ; ir < numReceptor ; ir++ ) {
1230 tsys[ir].resize( numChan ) ;
1231 for ( unsigned int ic = 0 ; ic < numChan ; ic++ ) {
1232 tsys[ir][ic] = tsysSpec[ir][ic].get() ;
1233 }
1234 }
1235 }
1236 else {
1237 tsys = defaultT ;
1238 }
1239 if ( scrow->isTcalSpectrumExists() ) {
1240 vector< vector<Temperature> > tcalSpec = scrow->getTcalSpectrum() ;
1241 unsigned int numReceptor = tcalSpec.size() ;
1242 unsigned int numChan = tcalSpec[0].size() ;
1243 tcal.resize( numReceptor ) ;
1244 for ( unsigned int ir = 0 ; ir < numReceptor ; ir++ ) {
1245 tcal[ir].resize( numReceptor ) ;
1246 for ( unsigned int ic = 0 ; ic < numChan ; ic++ ) {
1247 tcal[ir][ic] = tcalSpec[ir][ic].get() ;
1248 }
1249 }
1250 }
1251 else {
1252 tcal = defaultT ;
1253 }
1254 }
1255}
1256
1257vector<float> OldASDMReader::getOpacity( unsigned int idx )
1258{
1259 vector<float> tau(0) ;
1260 CalAtmosphereTable &atmtab = asdm_->getCalAtmosphere() ;
1261 unsigned int nrow = atmtab.size() ;
1262 if ( nrow > 0 ) {
1263 unsigned int index = dataIdList_[idx] ;
1264 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1265 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1266 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1267 //int feedid = vmsData_->v_feedId1[index] ;
1268 //Tag atag( antennaId_, TagType::Antenna ) ;
1269 //Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
1270 //DataDescriptionRow *ddrow = asdm_->getDataDescription().getRowByKey(ddtag) ;
1271 //Tag spwtag = ddrow->getSpectralWindowId() ;
1272 //SpectralWindowRow *spwrow = ddrow->getSpectralWindowUsingSpectralWindowId() ;
1273 //BasebandName bbname = spwrow->getBasebandName() ;
1274 //FeedRow *frow = asdm_->getFeed().getRowByKey( atag, spwtag, tint, feedid ) ;
1275 //int nfeed = frow->getNumReceptor() ;
1276 //vector<ReceiverRow *> rrows = frow->getReceivers() ;
1277 vector<CalAtmosphereRow *> atmrows = atmtab.get() ;
1278 //ReceiverBand rb = rrows[0]->getFrequencyBand() ;
1279 int row0 = -1 ;
1280 double eps = tint.getStart().getMJD() ;
1281 for ( unsigned int irow = 0 ; irow < nrow ; irow++ ) {
1282 CalAtmosphereRow *atmrow = atmrows[irow] ;
1283 if ( casa::String(atmrow->getAntennaName()) != antennaName_
1284 //|| atmrow->getReceiverBand() != rb
1285 //|| atmrow->getBasebandName() != bbname
1286 || atmrow->getCalDataUsingCalDataId()->getCalType() != CAL_ATMOSPHERE )
1287 continue ;
1288 else {
1289 double dt = tint.getStart().getMJD() - atmrow->getEndValidTime().getMJD() ;
1290 if ( dt >= 0 && dt < eps ) {
1291 eps = dt ;
1292 row0 = (int)irow ;
1293 }
1294 }
1295 }
1296 if ( row0 != -1 ) {
1297 CalAtmosphereRow *atmrow = atmrows[row0] ;
1298 tau = atmrow->getTau() ;
1299 }
1300 }
1301 else {
1302 tau.resize( 1 ) ;
1303 tau[0] = 0.0 ;
1304 }
1305 return tau ;
1306}
1307
1308void OldASDMReader::getWeatherInfo( unsigned int idx,
1309 float &temperature,
1310 float &pressure,
1311 float &humidity,
1312 float &windspeed,
1313 float &windaz )
1314{
1315 casa::String funcName = "getWeatherInfo" ;
1316
1317 temperature = 0.0 ;
1318 pressure = 0.0 ;
1319 humidity = 0.0 ;
1320 windspeed = 0.0 ;
1321 windaz = 0.0 ;
1322
1323 //logsink_->postLocally( LogMessage("weatherStationId_ = "+String::toString(weatherStationId_),LogOrigin(className_,funcName,WHERE)) ) ;
1324
1325 WeatherTable &wtab = asdm_->getWeather() ;
1326 if ( wtab.size() == 0 || weatherStationId_ == -1 )
1327 return ;
1328
1329 unsigned int index = dataIdList_[idx] ;
1330 //Tag anttag( antennaId_, TagType::Antenna ) ;
1331 //Tag sttag = (asdm_->getAntenna().getRowByKey( anttag ))->getStationId() ;
1332 Tag sttag( (unsigned int)weatherStationId_, TagType::Station ) ;
1333 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1334 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1335 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1336 //WeatherRow *wrow = wtab.getRowByKey( sttag, tint ) ;
1337 vector<WeatherRow *> *wrows = wtab.getByContext( sttag ) ;
1338 WeatherRow *wrow = (*wrows)[0] ;
1339 unsigned int nrow = wrows->size() ;
1340 //logsink_->postLocally( LogMessage("There are "+String::toString(nrow)+" rows for given context: stationId "+String::toString(weatherStationId_),LogOrigin(className_,funcName,WHERE)) ) ;
1341 ArrayTime startTime = getMidTime( tint ) ;
1342 if ( startTime < (*wrows)[0]->getTimeInterval().getStart() ) {
1343 temperature = (*wrows)[0]->getTemperature().get() ;
1344 pressure = (*wrows)[0]->getPressure().get() ;
1345 humidity = (*wrows)[0]->getRelHumidity().get() ;
1346 windspeed = (*wrows)[0]->getWindSpeed().get() ;
1347 windaz = (*wrows)[0]->getWindDirection().get() ;
1348 }
1349 else if ( startTime > getEndTime( (*wrows)[nrow-1]->getTimeInterval() ) ) {
1350 temperature = (*wrows)[nrow-1]->getTemperature().get() ;
1351 pressure = (*wrows)[nrow-1]->getPressure().get() ;
1352 humidity = (*wrows)[nrow-1]->getRelHumidity().get() ;
1353 windspeed = (*wrows)[nrow-1]->getWindSpeed().get() ;
1354 windaz = (*wrows)[nrow-1]->getWindDirection().get() ;
1355 }
1356 else {
1357 for ( unsigned int irow = 1 ; irow < wrows->size() ; irow++ ) {
1358 wrow = (*wrows)[irow-1] ;
1359 ArrayTime tStart = wrow->getTimeInterval().getStart() ;
1360 ArrayTime tEnd = (*wrows)[irow]->getTimeInterval().getStart() ;
1361 if ( startTime >= tStart && startTime <= tEnd ) {
1362 temperature = wrow->getTemperature().get() ;
1363 pressure = wrow->getPressure().get() ;
1364 humidity = wrow->getRelHumidity().get() ;
1365 windspeed = wrow->getWindSpeed().get() ;
1366 windaz = wrow->getWindDirection().get() ;
1367 break ;
1368 }
1369 }
1370 }
1371
1372 // Pa -> hPa
1373 pressure /= 100.0 ;
1374
1375 return ;
1376}
1377
1378void OldASDMReader::processStation()
1379{
1380 antennaPad_.resize(0) ;
1381 weatherStation_.resize(0) ;
1382 StationTable &stab = asdm_->getStation() ;
1383 vector<StationRow *> srows = stab.get() ;
1384 for ( unsigned int irow = 0 ; irow < srows.size() ; irow++ ) {
1385 StationType stype = srows[irow]->getType() ;
1386 Tag stag = srows[irow]->getStationId() ;
1387 if ( stype == ANTENNA_PAD )
1388 antennaPad_.push_back( stag ) ;
1389 else if ( stype == WEATHER_STATION )
1390 weatherStation_.push_back( stag ) ;
1391 }
1392
1393 weatherStationId_ = getClosestWeatherStation() ;
1394}
1395
1396int OldASDMReader::getClosestWeatherStation()
1397{
1398 if ( weatherStation_.size() == 0 )
1399 return -1 ;
1400
1401 Tag atag( antennaId_, TagType::Antenna ) ;
1402 Tag stag = (asdm_->getAntenna().getRowByKey( atag ))->getStationId() ;
1403 vector<double> apos( 3 ) ;
1404 StationTable &stab = asdm_->getStation() ;
1405 StationRow *srow = stab.getRowByKey( stag ) ;
1406 vector<Length> pos = srow->getPosition() ;
1407 apos[0] = pos[0].get() ;
1408 apos[1] = pos[1].get() ;
1409 apos[2] = pos[2].get() ;
1410
1411 double eps = 1.0e20 ;
1412 int retval = -1 ;
1413 for ( unsigned int ir = 0 ; ir < weatherStation_.size() ; ir++ ) {
1414 srow = stab.getRowByKey( weatherStation_[ir] ) ;
1415 vector<Length> wpos = srow->getPosition() ;
1416 double dist = (apos[0]-wpos[0].get())*(apos[0]-wpos[0].get())
1417 + (apos[1]-wpos[1].get())*(apos[1]-wpos[1].get())
1418 + (apos[2]-wpos[2].get())*(apos[2]-wpos[2].get()) ;
1419 if ( dist < eps ) {
1420 retval = (int)(weatherStation_[ir].getTagValue()) ;
1421 }
1422 }
1423
1424 return retval ;
1425}
1426
1427void OldASDMReader::getPointingInfo( unsigned int idx,
1428 vector<double> &dir,
1429 double &az,
1430 double &el,
1431 vector<double> &srate )
1432{
1433 String funcName = "getPointingInfo" ;
1434
1435 dir.resize(0) ;
1436 az = -1.0 ;
1437 el = -1.0 ;
1438 srate.resize(0) ;
1439
1440 Tag atag( antennaId_, TagType::Antenna ) ;
1441 unsigned int index = dataIdList_[idx] ;
1442 vector<PointingRow *> *prows = asdm_->getPointing().getByContext( atag ) ;
1443
1444 if ( prows == 0 )
1445 return ;
1446
1447 PointingRow *prow ;
1448 PointingRow *qrow ;
1449 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1450 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1451 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1452
1453 unsigned int nrow = prows->size() ;
1454 //logsink_->postLocally( LogMessage("There are " << nrow << " rows for given context: antennaId "+String::toString(antennaId_),LogOrigin(className_,funcName,WHERE)) ) ;
1455
1456// for ( unsigned int irow = 0 ; irow < nrow ; irow++ ) {
1457// prow = (*prows)[irow] ;
1458// ArrayTimeInterval ati = prow->getTimeInterval() ;
1459// ArrayTime pst = ati.getStart() ;
1460// ArrayTime pet( ati.getStartInMJD()+ati.getDurationInDays() ) ;
1461// logsink_->postLocally( LogMessage("start: "+pst.toFITS(),LogOrigin(className_,funcName,WHERE)) ) ;
1462// logsink_->postLocally( LogMessage("end: "+pet.toFITS(),LogOrigin(className_,funcName,WHERE)) ) ;
1463// }
1464
1465 srate.resize(2) ;
1466 srate[0] = 0.0 ;
1467 srate[1] = 0.0 ;
1468 az = 0.0 ;
1469 el = 0.0 ;
1470 //double tcen = 0.0 ;
1471 double tcen = getMidTime( tint ).getMJD() ;
1472
1473 //
1474 // shape of pointingDirection is (numSample,2) if usePolynomial = False, while it is
1475 // (numTerm,2) if usePolynomial = True.
1476 //
1477 // In the former case, typical sampled time interval is 48msec, which is very small
1478 // compared with typical integration time (~a few sec).
1479 // Scan rate for this case is always [0,0] (or get slope?).
1480 //
1481 // In the later case, scan rate is (pointingDirection[1][0],pointingDirection[1][1])
1482 //
1483 // PointingDirection seems to store AZEL
1484 //
1485 ArrayTimeInterval pTime0 = (*prows)[0]->getTimeInterval() ;
1486 ArrayTimeInterval pTime1 = (*prows)[nrow-1]->getTimeInterval() ;
1487 //if ( tint.getStartInMJD()+tint.getDurationInDays() < pTime0.getStartInMJD() ) {
1488 if ( getEndTime( tint ) < getStartTime( pTime0 ) ) {
1489 logsink_->postLocally( LogMessage( "ArrayTimeInterval out of bounds: no data for given position (tint < ptime)", LogOrigin(className_,funcName,WHERE), LogMessage::WARN ) ) ;
1490 prow = (*prows)[0] ;
1491 vector< vector<double> > dirA = pointingDir( prow ) ;
1492 az = dirA[0][0] ;
1493 el = dirA[0][1] ;
1494 if ( prow->getUsePolynomials() && prow->getNumTerm() > 1 ) {
1495 srate[0] = dirA[1][0] ;
1496 srate[1] = dirA[1][1] ;
1497 }
1498 }
1499 //else if ( tint.getStartInMJD() > pTime1.getStartInMJD()+pTime1.getDurationInDays() ) {
1500 else if ( getStartTime( tint ) > getEndTime( pTime1 ) ) {
1501 logsink_->postLocally( LogMessage( "ArrayTimeInterval out of bounds: no data for given position (tint > ptime)", LogOrigin(className_,funcName,WHERE), LogMessage::WARN ) ) ;
1502 prow = (*prows)[nrow-1] ;
1503 int numSample = prow->getNumSample() ;
1504 vector< vector<double> > dirA = pointingDir( prow ) ;
1505 if ( prow->getUsePolynomials() ) {
1506 az = dirA[0][0] ;
1507 el = dirA[0][1] ;
1508 if ( prow->getNumTerm() > 1 ) {
1509 srate[0] = dirA[1][0] ;
1510 srate[1] = dirA[1][1] ;
1511 }
1512 }
1513 else {
1514 az = dirA[numSample-1][0] ;
1515 el = dirA[numSample-1][1] ;
1516 }
1517 }
1518 else {
1519 ArrayTime startTime = tint.getStart() ;
1520 ArrayTime endTime = getEndTime( tint ) ;
1521 int row0 = -1 ;
1522 int row1 = -1 ;
1523 int row2 = -1 ;
1524 double dt0 = getMidTime( tint ).getMJD() ;
1525 for ( unsigned int irow = 0 ; irow < nrow ; irow++ ) {
1526 prow = (*prows)[irow] ;
1527 double dt = getMidTime( tint ).getMJD() - getMidTime( prow->getTimeInterval() ).getMJD() ;
1528 if ( dt > 0 && dt < dt0 ) {
1529 dt0 = dt ;
1530 row2 = irow ;
1531 }
1532 if ( prow->getTimeInterval().contains( startTime ) )
1533 row0 = irow ;
1534 else if ( prow->getTimeInterval().contains( endTime ) )
1535 row1 = irow ;
1536 if ( row0 != -1 && row1 != -1 )
1537 break ;
1538 }
1539 //logsink_->postLocally( LogMessage("row0 = "+String::toString(row0)+", row1 = "+String::toString(row1)+", row2 = "+String::toString(row2),LogOrigin(className_,funcName,WHERE)) ) ;
1540 if ( row0 == -1 && row1 == -1 ) {
1541 prow = (*prows)[row2] ;
1542 qrow = (*prows)[row2+1] ;
1543 double t0 = getEndTime( prow->getTimeInterval() ).getMJD() ;
1544 double t1 = qrow->getTimeInterval().getStartInMJD() ;
1545 vector< vector<double> > dirA = pointingDir( prow ) ;
1546 vector< vector<double> > dirB = pointingDir( qrow ) ;
1547 double da0 = dirA[0][0] ;
1548 double db0 = dirB[0][0] ;
1549 double da1 = dirA[0][1] ;
1550 double db1 = dirB[0][1] ;
1551 if ( prow->getUsePolynomials() && qrow->getUsePolynomials() ) {
1552 double dt = ( tcen - t0 ) / ( t1 - t0 ) ;
1553 az = da0 + ( db0 - da0 ) * dt ;
1554 el = da1 + ( db1 - da1 ) * dt ;
1555 if ( prow->getNumTerm() > 0 && qrow->getNumTerm() > 1 ) {
1556 double ra0 = dirA[1][0] ;
1557 double rb0 = dirB[1][0] ;
1558 double ra1 = dirA[1][1] ;
1559 double rb1 = dirB[1][1] ;
1560 srate[0] = ra0 + ( rb0 - ra0 ) * dt ;
1561 srate[1] = ra1 + ( rb1 - ra1 ) * dt ;
1562 }
1563 }
1564 else if ( !(qrow->getUsePolynomials()) ) {
1565 double dt = ( tcen - t0 ) / ( t1 - t0 ) ;
1566 az = da0 + ( db0 - da0 ) * dt ;
1567 el = da1 + ( db1 - da1 ) * dt ;
1568 }
1569 //else {
1570 // nothing to do
1571 //}
1572 }
1573 else {
1574 int ndir = 0 ;
1575 if ( row0 == -1 ) {
1576 row0 = row1 ;
1577 }
1578 else if ( row1 == -1 ) {
1579 row1 = row0 ;
1580 }
1581 prow = (*prows)[row0] ;
1582 qrow = (*prows)[row1] ;
1583 if ( prow->getUsePolynomials() && qrow->getUsePolynomials() ) {
1584 //logsink_->postLocally( LogMessage("usePolynomial = True",LogOrigin(className_,funcName,WHERE)) ) ;
1585 if ( row0 == row1 ) {
1586 prow = (*prows)[row0] ;
1587 vector< vector<double> > dirA = pointingDir( prow ) ;
1588 az = dirA[0][0] ;
1589 el = dirA[0][1] ;
1590 if ( prow->getNumTerm() > 1 ) {
1591 srate[0] = dirA[1][0] ;
1592 srate[1] = dirA[1][1] ;
1593 }
1594 }
1595 else {
1596 prow = (*prows)[row0] ;
1597 qrow = (*prows)[row1] ;
1598 vector< vector<double> > dirA = pointingDir( prow ) ;
1599 vector< vector<double> > dirB = pointingDir( qrow ) ;
1600 double t0 = getMidTime( prow->getTimeInterval() ).getMJD() ;
1601 double t1 = getMidTime( qrow->getTimeInterval() ).getMJD() ;
1602 double dt = ( tcen - t0 ) / ( t1 - t0 ) ;
1603 double da0 = dirA[0][0] ;
1604 double db0 = dirB[0][0] ;
1605 double da1 = dirA[0][1] ;
1606 double db1 = dirB[0][1] ;
1607 az = da0 + ( db0 - da0 ) * dt ;
1608 el = da1 + ( db1 - db0 ) * dt ;
1609 if ( prow->getNumTerm() > 0 && qrow->getNumTerm() > 1 ) {
1610 double ra0 = dirA[1][0] ;
1611 double rb0 = dirB[1][0] ;
1612 double ra1 = dirA[1][1] ;
1613 double rb1 = dirB[1][1] ;
1614 srate[0] = ra0 + ( rb0 - ra0 ) * dt ;
1615 srate[1] = ra1 + ( rb1 - ra1 ) * dt ;
1616 }
1617 }
1618 }
1619 else if ( prow->getUsePolynomials() == qrow->getUsePolynomials() ) {
1620 //logsink_->postLocally( LogMessage("numSample == numTerm",LogOrigin(className_,funcName,WHERE)) ) ;
1621 tcen = 0.0 ;
1622 for ( int irow = row0 ; irow <= row1 ; irow++ ) {
1623 prow = (*prows)[irow] ;
1624 int numSample = prow->getNumSample() ;
1625 //logsink_->postLocally( LogMessage("numSample = "+String::toString(numSample),LogOrigin(className_,funcName,WHERE)) ) ;
1626 vector< vector<double> > dirA = pointingDir( prow ) ;
1627 if ( prow->isSampledTimeIntervalExists() ) {
1628 //logsink_->postLocally( LogMessage("sampledTimeIntervalExists",LogOrigin(className_,funcName,WHERE)) ) ;
1629 vector<ArrayTimeInterval> stime = prow->getSampledTimeInterval() ;
1630 for ( int isam = 0 ; isam < numSample ; isam++ ) {
1631 //if ( tint.overlaps( stime[isam] ) ) {
1632 if ( tint.contains( stime[isam] ) ) {
1633 az += dirA[isam][0] ;
1634 el += dirA[isam][1] ;
1635 tcen += getMidTime( stime[isam] ).getMJD() ;
1636 ndir++ ;
1637 }
1638 }
1639 }
1640 else {
1641 double sampleStart = prow->getTimeInterval().getStartInMJD() ;
1642 double sampleInterval = prow->getTimeInterval().getDurationInDays() / (double)numSample ;
1643 //logsink_->postLocally( LogMessage("sampleStart = "+String::toString(sampleStart),LogOrigin(className_,funcName,WHERE)) )
1644 //logsink_->postLocally( LogMessage("sampleInterval = "+String::toString(sampleInterval),LogOrigin(className_,funcName,WHERE)) ) ;
1645 //logsink_->postLocally( LogMessage("tint = "+tint.toString(),LogOrigin(className_,funcName,WHERE)) ) ;
1646 for ( int isam = 0 ; isam < numSample ; isam++ ) {
1647 ArrayTimeInterval stime( sampleStart+isam*sampleInterval, sampleInterval ) ;
1648 //if ( tint.overlaps( stime ) ) {
1649 if ( tint.contains( stime ) ) {
1650 az += dirA[isam][0] ;
1651 el += dirA[isam][1] ;
1652 tcen += getMidTime( stime ).getMJD() ;
1653 ndir++ ;
1654 }
1655 }
1656 }
1657 }
1658 if ( ndir > 1 ) {
1659 az /= (double)ndir ;
1660 el /= (double)ndir ;
1661 tcen /= (double)ndir ;
1662 }
1663 }
1664 //else {
1665 // nothing to do
1666 //}
1667 }
1668
1669 AntennaRow *arow = asdm_->getAntenna().getRowByKey( Tag( antennaId_, TagType::Antenna ) ) ;
1670 StationRow *srow = arow->getStationUsingStationId() ;
1671 vector<Length> antposL = srow->getPosition() ;
1672 casa::Vector<casa::Double> antpos( 3 ) ;
1673 for ( int i = 0 ; i < 3 ; i++ )
1674 antpos[i] = antposL[i].get() ;
1675 //logsink_->postLocally( LogMessage("tcen = "+String::toString(tcen),LogOrigin(className_,funcName,WHERE)) ) ;
1676 //logsink_->postLocally( LogMessage("antpos = "+String::toString(antpos),LogOrigin(className_,funcName,WHERE)) ) ;
1677 toJ2000( dir, az, el, tcen, antpos ) ;
1678
1679 }
1680
1681 return ;
1682}
1683
1684ArrayTime OldASDMReader::getMidTime( const ArrayTimeInterval &t )
1685{
1686 return ArrayTime( t.getStartInMJD() + 0.5 * t.getDurationInDays() ) ;
1687}
1688
1689ArrayTime OldASDMReader::getEndTime( const ArrayTimeInterval &t )
1690{
1691 return ArrayTime( t.getStartInMJD() + t.getDurationInDays() ) ;
1692}
1693
1694ArrayTime OldASDMReader::getStartTime( const ArrayTimeInterval &t )
1695{
1696 return ArrayTime( t.getStartInMJD() ) ;
1697}
1698
1699void OldASDMReader::toJ2000( vector<double> &dir,
1700 double az,
1701 double el,
1702 double mjd,
1703 casa::Vector<casa::Double> antpos )
1704{
1705 String funcName = "toJ2000" ;
1706
1707 vector<double> azel( 2 ) ;
1708 azel[0] = az ;
1709 azel[1] = el ;
1710 dir = toJ2000( azel, "AZELGEO", mjd, antpos ) ;
1711}
1712
1713vector<double> OldASDMReader::toJ2000( vector<double> dir,
1714 casa::String dirref,
1715 double mjd,
1716 casa::Vector<casa::Double> antpos )
1717{
1718 casa::String funcName = "toJ2000" ;
1719
1720 vector<double> newd( dir ) ;
1721 if ( dirref != "J2000" ) {
1722 casa::MEpoch me( casa::Quantity( (casa::Double)mjd, "d" ), casa::MEpoch::UTC ) ;
1723 casa::Vector<casa::Quantity> qantpos( 3 ) ;
1724 qantpos[0] = casa::Quantity( antpos[0], "m" ) ;
1725 qantpos[1] = casa::Quantity( antpos[1], "m" ) ;
1726 qantpos[2] = casa::Quantity( antpos[2], "m" ) ;
1727 casa::MPosition mp( casa::MVPosition( qantpos ),
1728 casa::MPosition::ITRF ) ;
1729 //ostringstream oss ;
1730 //mp.print( oss ) ;
1731 //logsink_->postLocally( LogMessage(oss.str(),LogOrigin(className_,funcName,WHERE)) ) ;
1732
1733 casa::MeasFrame mf( me, mp ) ;
1734 casa::MDirection::Types dirtype ;
1735 casa::Bool b = casa::MDirection::getType( dirtype, dirref ) ;
1736 if ( b ) {
1737 casa::MDirection::Convert toj2000( dirtype,
1738 casa::MDirection::Ref( casa::MDirection::J2000, mf ) ) ;
1739 casa::Vector<casa::Double> cdir = toj2000( dir ).getAngle( "rad" ).getValue() ;
1740 //logsink_->postLocally( LogMessage("cdir = "+String::toString(cdir),LogOrigin(className_,funcName,WHERE)) ) ;
1741 newd[0] = (double)(cdir[0]) ;
1742 newd[1] = (double)(cdir[1]) ;
1743 }
1744 }
1745 return newd ;
1746}
1747
1748void OldASDMReader::setLogger( CountedPtr<LogSinkInterface> &logsink )
1749{
1750 logsink_ = logsink ;
1751}
1752
1753string OldASDMReader::getFrame()
1754{
1755 casa::String funcName = "getFrame" ;
1756
1757 // default is TOPO
1758 string frame = "TOPO" ;
1759
1760 SpectralWindowTable &spwtab = asdm_->getSpectralWindow() ;
1761 vector<SpectralWindowRow *> rows = spwtab.get() ;
1762 vector<FrequencyReferenceCode> measFreqRef( rows.size() ) ;
1763 int nref = 0 ;
1764 for ( unsigned int irow = 0 ; irow < rows.size() ; irow++ ) {
1765 int nchan = rows[irow]->getNumChan() ;
1766 if ( nchan != 4 ) {
1767 if ( rows[irow]->isMeasFreqRefExists() ) {
1768 measFreqRef[nref] = rows[irow]->getMeasFreqRef() ;
1769 nref++ ;
1770 }
1771 }
1772 }
1773 if ( nref != 0 ) {
1774 frame = CFrequencyReferenceCode::toString( measFreqRef[0] ) ;
1775 }
1776
1777 //logsink_->postLocally( LogMessage("frame = "+String::toString(frame),LogOrigin(className_,funcName,WHERE)) ) ;
1778
1779 return frame ;
1780}
1781
1782int OldASDMReader::getNumIFs()
1783{
1784 casa::String funcName = "getNumIFs" ;
1785
1786 int nif = 0 ;
1787 vector<SpectralWindowRow *> rows = asdm_->getSpectralWindow().get() ;
1788 unsigned int nrow = rows.size() ;
1789 // check if all rows have freqGroup attribute
1790 bool freqGroupExists = true ;
1791 bool countedWvr = false ;
1792 for ( unsigned int irow = 0 ; irow < nrow ; irow++ ) {
1793 freqGroupExists &= rows[irow]->isFreqGroupExists() ;
1794 if ( rows[irow]->getNumChan() == 4 ) {
1795 if ( !countedWvr ) {
1796 countedWvr = true ;
1797 nif++ ;
1798 }
1799 }
1800 else {
1801 nif++ ;
1802 }
1803 }
1804
1805 if ( freqGroupExists ) {
1806 vector<int> freqGroup(0) ;
1807 for ( unsigned int irow = 0 ; irow < nrow ; irow++ ) {
1808 int fg = rows[irow]->getFreqGroup() ;
1809 if ( (int)count( freqGroup.begin(), freqGroup.end(), fg ) == 0 ) {
1810 freqGroup.push_back( fg ) ;
1811 }
1812 }
1813 nif = freqGroup.size() ;
1814 }
1815
1816 //logsink_->postLocally( LogMessage("nif = "+String::toString(nif),LogOrigin(className_,funcName,WHERE)) ) ;
1817
1818 return nif ;
1819}
1820
1821SysCalRow *OldASDMReader::getSysCalRow( unsigned int idx )
1822{
1823 String funcName = "getSysCalRow" ;
1824
1825 SysCalRow *row = 0 ;
1826 unsigned int index = dataIdList_[idx] ;
1827 Tag anttag( antennaId_, TagType::Antenna ) ;
1828 int feedid = vmsData_->v_feedId1[index] ;
1829 //ArrayTimeInterval tint( vmsData_->v_time[index]*s2d, vmsData_->v_interval[index]*s2d ) ;
1830 double startSec = vmsData_->v_time[index] - 0.5 * vmsData_->v_interval[index] ;
1831 ArrayTimeInterval tint( startSec*s2d, vmsData_->v_interval[index]*s2d ) ;
1832 Tag ddtag( vmsData_->v_dataDescId[index], TagType::DataDescription ) ;
1833 Tag spwtag = asdm_->getDataDescription().getRowByKey(ddtag)->getSpectralWindowId() ;
1834 //int nchan = asdm_->getSpectralWindow().getRowByKey(spwtag)->getNumChan() ;
1835 vector< vector<float> > defaultTsys( 1, vector<float>( 1, 1.0 ) ) ;
1836 SysCalTable &sctab = asdm_->getSysCal() ;
1837 //vector<SysCalRow *> rows = sctab.get() ;
1838 vector<SysCalRow *> *rows = sctab.getByContext( anttag, spwtag, feedid ) ;
1839 //if ( nrow == 0 ) {
1840 if ( rows == 0 ) {
1841 //logsink_->postLocally( LogMessage("no rows in SysCal table",LogOrigin(className_,funcName,WHERE)) ) ;
1842 row = 0 ;
1843 }
1844 else {
1845 unsigned int nrow = rows->size() ;
1846 //logsink_->postLocally( LogMessage("nrow = "+String::toString(nrow),LogOrigin(className_,funcName,WHERE)) ) ;
1847 int scindex = -1 ;
1848 if ( nrow == 1 ) {
1849 scindex = 0 ;
1850 }
1851 else if ( getEndTime( tint ) <= getStartTime( (*rows)[0]->getTimeInterval() ) )
1852 scindex = 0 ;
1853 else {
1854 for ( unsigned int irow = 0 ; irow < nrow-1 ; irow++ ) {
1855 ArrayTime t = getMidTime( tint ) ;
1856 if ( t > getStartTime( (*rows)[irow]->getTimeInterval() )
1857 && t <= getStartTime( (*rows)[irow+1]->getTimeInterval() ) ) {
1858 scindex = irow ;
1859 break ;
1860 }
1861 }
1862 if ( scindex == -1 )
1863 scindex = nrow-1 ;
1864 }
1865 //logsink_->postLocally( LogMessage("irow = "+String::toString(scindex),LogOrigin(className_,funcName,WHERE)) ) ;
1866 row = (*rows)[scindex] ;
1867 }
1868 return row ;
1869}
1870
1871double OldASDMReader::limitedAngle( double angle )
1872{
1873 if ( angle > C::pi )
1874 while ( angle > C::pi ) angle -= C::_2pi ;
1875 else if ( angle < -C::pi )
1876 while ( angle < -C::pi ) angle += C::_2pi ;
1877 return angle ;
1878}
1879
1880vector< vector<double> > OldASDMReader::pointingDir( PointingRow *row )
1881{
1882 vector< vector<Angle> > aTar = row->getTarget() ;
1883 vector< vector<Angle> > aOff = row->getOffset() ;
1884 vector< vector<Angle> > aDir = row->getPointingDirection() ;
1885 vector< vector<Angle> > aEnc = row->getEncoder() ;
1886 unsigned int n = aTar.size() ;
1887 vector< vector<double> > dir( n ) ;
1888 double factor = 1.0 / cos( aTar[0][1].get() ) ;
1889 for ( unsigned int i = 0 ; i < n ; i++ ) {
1890 dir[i].resize( 2 ) ;
1891 /**
1892 * This is approximate way to add offset taking tracking error
1893 * into account
1894 *
1895 * az = dir[0][0] = target[0][0] + offset[0][0] / cos(el)
1896 * + encorder[0][0] - direction[0][0]
1897 * el = dir[0][1] = target[0][1] + offset[0][1]
1898 * + encorder[0][1] - direction[0][1]
1899 **/
1900 dir[i][0] = aTar[i][0].get() + factor * aOff[i][0].get()
1901 + aEnc[i][0].get() - aDir[i][0].get() ;
1902 dir[i][1] = aTar[i][1].get() + aOff[i][1].get()
1903 + aEnc[i][1].get() - aDir[i][1].get() ;
1904 }
1905 return dir ;
1906}
Note: See TracBrowser for help on using the repository browser.