source: trunk/external-alma/atnf/PKSIO/NROReader.cc@ 2771

Last change on this file since 2771 was 2766, checked in by Takeshi Nakazato, 12 years ago

New Development: No

JIRA Issue: No

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...

Use CountedPtr for NRODataRecord object instead of raw C pointer.
File size: 22.4 KB
Line 
1//#---------------------------------------------------------------------------
2//# NROReader.cc: Base class for NRO headerdata.
3//#---------------------------------------------------------------------------
4//# Copyright (C) 2000-2006
5//# Associated Universities, Inc. Washington DC, USA.
6//#
7//# This library is free software; you can redistribute it and/or modify it
8//# under the terms of the GNU Library General Public License as published by
9//# the Free Software Foundation; either version 2 of the License, or (at your
10//# option) any later version.
11//#
12//# This library is distributed in the hope that it will be useful, but WITHOUT
13//# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14//# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
15//# License for more details.
16//#
17//# You should have received a copy of the GNU Library General Public License
18//# along with this library; if not, write to the Free Software Foundation,
19//# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
20//#
21//# Correspondence concerning AIPS++ should be addressed as follows:
22//# Internet email: aips2-request@nrao.edu.
23//# Postal address: AIPS++ Project Office
24//# National Radio Astronomy Observatory
25//# 520 Edgemont Road
26//# Charlottesville, VA 22903-2475 USA
27//#
28//# $Id$
29//#---------------------------------------------------------------------------
30//# Original: 2008/10/30, Takeshi Nakazato, NAOJ
31//#---------------------------------------------------------------------------
32
33#include <atnf/PKSIO/NROReader.h>
34#include <atnf/PKSIO/NRO45Reader.h>
35#include <atnf/PKSIO/ASTEReader.h>
36#include <atnf/PKSIO/ASTEFXReader.h>
37#include <atnf/PKSIO/NRO45FITSReader.h>
38#include <atnf/PKSIO/NROOTFDataset.h>
39#include <atnf/PKSIO/ASTEDataset.h>
40
41#include <measures/Measures/MEpoch.h>
42#include <measures/Measures/MPosition.h>
43#include <measures/Measures/MDirection.h>
44#include <measures/Measures/MCDirection.h>
45#include <measures/Measures/MeasFrame.h>
46#include <measures/Measures/MeasConvert.h>
47
48#include <casa/Logging/LogIO.h>
49#include <casa/IO/RegularFileIO.h>
50#include <casa/OS/File.h>
51#include <casa/OS/Time.h>
52
53#include <stdio.h>
54#include <string>
55#include <iomanip>
56
57using namespace std ;
58
59//
60// getNROReader
61//
62// Return an appropriate NROReader for a NRO 45m and ASTE dataset.
63//
64NROReader *getNROReader( const String filename,
65 String &datatype )
66{
67 LogIO os( LogOrigin( "", "getNROReader()", WHERE ) ) ;
68
69 // Check accessibility of the input.
70 File inFile( filename ) ;
71 if ( !inFile.exists() ) {
72 datatype = filename + " not found." ;
73 return 0 ;
74 }
75
76 if ( !inFile.isReadable() ) {
77 datatype = filename + " is not readable." ;
78 return 0 ;
79 }
80
81 // Determine the type of input.
82 NROReader *reader = 0;
83 if ( inFile.isRegular() ) {
84 FILE *file ;
85 file = fopen( filename.c_str(), "r" ) ;
86 // read LOFIL0
87 char buf[9];
88 fread( buf, 4, 1, file ) ;
89 buf[4] = '\0' ;
90 // DEBUG
91 //os << LogIO::NORMAL << "getNROReader:: buf = " << String(buf) << LogIO::POST ;
92 //
93 if ( string( buf ) == "XTEN" ) {
94 // FITS data
95 datatype = "NRO 45m FITS" ;
96 reader = new NRO45FITSReader( filename ) ;
97 }
98 else if ( string( buf ) == "RW-F") {
99 // ASTE-FX data
100 datatype = "ASTE-FX";
101 reader = new ASTEFXReader( filename );
102 } else {
103 // otherwise, read SITE0
104 NRODataset *d = new NROOTFDataset( filename ) ;
105 int size = d->getDataSize() - 188 ;
106 delete d ;
107 fseek( file, size, SEEK_SET ) ;
108 fread( buf, 8, 1, file ) ;
109 buf[8] = '\0' ;
110 // DEBUG
111 //cout << "getNROReader:: buf = " << buf << endl ;
112 //
113 if ( string( buf ) == "NRO" ) {
114 // NRO 45m data
115 datatype = "NRO 45m OTF" ;
116 reader = new NRO45Reader( filename );
117 }
118 else {
119 d = new ASTEDataset( filename ) ;
120 size = d->getDataSize() - 188 ;
121 delete d ;
122 fseek( file, size, SEEK_SET ) ;
123 fread( buf, 8, 1, file ) ;
124 buf[8] = '\0' ;
125 // DEBUG
126 //cout << "getNROReader:: buf = " << buf << endl ;
127 //
128 if ( string( buf ) == "ASTE" ) {
129 // ASTE data
130 datatype = "ASTE" ;
131 reader = new ASTEReader( filename ) ;
132 }
133 else {
134 datatype = "UNRECOGNIZED INPUT FORMAT";
135 }
136 }
137 }
138 fclose( file ) ;
139 } else {
140 datatype = "UNRECOGNIZED INPUT FORMAT";
141 }
142
143 // DEBUG
144 os << LogIO::NORMAL << "Data format of " << filename << ": " << datatype << LogIO::POST ;
145 //
146
147 // return reader if exists
148 if ( reader ) {
149 reader->read() ;
150 return reader ;
151 }
152
153 return 0 ;
154}
155
156
157//
158// getNROReader
159//
160// Search a list of directories for a NRO 45m and ASTE dataset and return an
161// appropriate NROReader for it.
162//
163NROReader* getNROReader( const String name,
164 const Vector<String> directories,
165 int &iDir,
166 String &datatype )
167{
168 int nDir = directories.size();
169 for ( iDir = 0; iDir < nDir; iDir++ ) {
170 string inName = directories[iDir] + "/" + name;
171 NROReader *reader = getNROReader( inName, datatype ) ;
172
173 if (reader != 0) {
174 return reader;
175 }
176 }
177
178 iDir = -1;
179 return 0;
180}
181
182
183//----------------------------------------------------------------------
184// constructor
185NROReader::NROReader( string name )
186 : dataset_( NULL ),
187 srcdir_( 0 ),
188 msrcdir_( 0 ),
189 freqRefFromVREF_( false ) // import frequency as REST
190{
191 // initialization
192 filename_ = name ;
193 coord_ = -1 ;
194}
195
196// destructor
197NROReader::~NROReader()
198{
199}
200
201// set frequency reference frame
202void NROReader::setFreqRefFromVREF( bool fromVREF )
203{
204 os_.origin( LogOrigin( "NROReader", "setFreqRefFromVREF", WHERE ) ) ;
205 os_ << ((fromVREF) ? "Take frequency reference frame from VREF" :
206 "Use frequency reference frame REST") << LogIO::POST ;
207
208 freqRefFromVREF_ = fromVREF ;
209}
210
211// get spectrum
212vector< vector<double> > NROReader::getSpectrum()
213{
214 return dataset_->getSpectrum() ;
215}
216
217Int NROReader::getPolarizationNum()
218{
219 return dataset_->getPolarizationNum() ;
220}
221
222double NROReader::getStartTime()
223{
224 //char *startTime = dataset_->getLOSTM() ;
225 string startTime = dataset_->getLOSTM() ;
226 //cout << "getStartTime() startTime = " << startTime << endl ;
227 return getMJD( startTime ) ;
228}
229
230double NROReader::getEndTime()
231{
232 //char *endTime = dataset_->getLOETM() ;
233 string endTime = dataset_->getLOETM() ;
234 return getMJD( endTime ) ;
235}
236
237vector<double> NROReader::getStartIntTime()
238{
239 return dataset_->getStartIntTime() ;
240}
241
242double NROReader::getMJD( char *time )
243{
244 // TODO: should be checked which time zone the time depends on
245 // 2008/11/14 Takeshi Nakazato
246 string strStartTime = string( time ) ;
247 string strYear = strStartTime.substr( 0, 4 ) ;
248 string strMonth = strStartTime.substr( 4, 2 ) ;
249 string strDay = strStartTime.substr( 6, 2 ) ;
250 string strHour = strStartTime.substr( 8, 2 ) ;
251 string strMinute = strStartTime.substr( 10, 2 ) ;
252 string strSecond = strStartTime.substr( 12, strStartTime.size() - 12 ) ;
253 uInt year = atoi( strYear.c_str() ) ;
254 uInt month = atoi( strMonth.c_str() ) ;
255 uInt day = atoi( strDay.c_str() ) ;
256 uInt hour = atoi( strHour.c_str() ) ;
257 uInt minute = atoi( strMinute.c_str() ) ;
258 double second = atof( strSecond.c_str() ) ;
259 Time t( year, month, day, hour, minute, second ) ;
260
261 return t.modifiedJulianDay() ;
262}
263
264double NROReader::getMJD( string strStartTime )
265{
266 // TODO: should be checked which time zone the time depends on
267 // 2008/11/14 Takeshi Nakazato
268 string strYear = strStartTime.substr( 0, 4 ) ;
269 string strMonth = strStartTime.substr( 4, 2 ) ;
270 string strDay = strStartTime.substr( 6, 2 ) ;
271 string strHour = strStartTime.substr( 8, 2 ) ;
272 string strMinute = strStartTime.substr( 10, 2 ) ;
273 string strSecond = strStartTime.substr( 12, strStartTime.size() - 12 ) ;
274 uInt year = atoi( strYear.c_str() ) ;
275 uInt month = atoi( strMonth.c_str() ) ;
276 uInt day = atoi( strDay.c_str() ) ;
277 uInt hour = atoi( strHour.c_str() ) ;
278 uInt minute = atoi( strMinute.c_str() ) ;
279 double second = atof( strSecond.c_str() ) ;
280 Time t( year, month, day, hour, minute, second ) ;
281
282 return t.modifiedJulianDay() ;
283}
284
285vector<Bool> NROReader::getIFs()
286{
287 return dataset_->getIFs() ;
288}
289
290vector<Bool> NROReader::getBeams()
291{
292 vector<Bool> v ;
293 vector<int> arry = dataset_->getARRY() ;
294 for ( uInt i = 0 ; i < arry.size() ; i++ ) {
295 if ( arry[i] != 0 ) {
296 v.push_back( True ) ;
297 }
298 }
299
300 // DEBUG
301 //cout << "NROReader::getBeams() number of beam is " << v.size() << endl ;
302 //
303
304 return v ;
305}
306
307// Get SRCDIRECTION in RADEC(J2000)
308Vector<Double> NROReader::getSourceDirection()
309{
310 if ( msrcdir_.nelements() == 2 )
311 return msrcdir_ ;
312
313 srcdir_.resize( 2 ) ;
314 srcdir_[0] = Double( dataset_->getRA0() ) ;
315 srcdir_[1] = Double( dataset_->getDEC0() ) ;
316 char epoch[5] ;
317 strncpy( epoch, (dataset_->getEPOCH()).c_str(), 5 ) ;
318 if ( strncmp( epoch, "B1950", 5 ) == 0 ) {
319 // convert to J2000 value
320 MDirection result =
321 MDirection::Convert( MDirection( Quantity( srcdir_[0], "rad" ),
322 Quantity( srcdir_[1], "rad" ),
323 MDirection::Ref( MDirection::B1950 ) ),
324 MDirection::Ref( MDirection::J2000 ) ) () ;
325 msrcdir_ = result.getAngle().getValue() ;
326 if ( msrcdir_[0] < 0.0 && srcdir_[0] >= 0.0 )
327 msrcdir_[0] = 2.0 * M_PI + msrcdir_[0] ;
328 //cout << "NROReader::getSourceDirection() SRCDIRECTION convert from ("
329 //<< srcra << "," << srcdec << ") B1950 to ("
330 //<< v( 0 ) << ","<< v( 1 ) << ") J2000" << endl ;
331 //LogIO os( LogOrigin( "NROReader", "getSourceDirection()", WHERE ) ) ;
332 //os << LogIO::NORMAL << "SRCDIRECTION convert from ("
333 // << srcra << "," << srcdec << ") B1950 to ("
334 // << v( 0 ) << ","<< v( 1 ) << ") J2000" << LogIO::POST ;
335 }
336 else if ( strncmp( epoch, "J2000", 5 ) == 0 ) {
337 msrcdir_.reference( srcdir_ ) ;
338 }
339
340 return msrcdir_ ;
341}
342
343// Get DIRECTION in RADEC(J2000)
344Vector<Double> NROReader::getDirection( int i )
345{
346 Vector<Double> v( 2 ) ;
347 const NRODataRecord *record = dataset_->getRecord( i ) ;
348 char epoch[5] ;
349 strncpy( epoch, (dataset_->getEPOCH()).c_str(), 5 ) ;
350 int icoord = dataset_->getSCNCD() ;
351 double scantime = dataset_->getScanTime( i ) ;
352 initConvert( icoord, scantime, epoch ) ;
353 v[0]= Double( record->SCX ) ;
354 v[1] = Double( record->SCY ) ;
355 if ( converter_.null() )
356 // no conversion
357 return v ;
358 else {
359 Vector<Double> v2 = (*converter_)( v ).getAngle().getValue() ;
360 if ( v2[0] < 0.0 && v[0] >= 0.0 )
361 v2[0] = 2.0 * M_PI + v2[0] ;
362 return v2 ;
363 }
364}
365
366void NROReader::initConvert( int icoord, double t, char *epoch )
367{
368 if ( icoord == 0 && strncmp( epoch, "J2000", 5 ) == 0 )
369 // no conversion
370 return ;
371
372 if ( converter_.null() || icoord != coord_ ) {
373 LogIO os( LogOrigin( "NROReader", "initConvert()", WHERE ) ) ;
374 coord_ = icoord ;
375 if ( coord_ == 0 ) {
376 // RADEC (B1950) -> RADEC (J2000)
377 os << "Creating converter from RADEC (B1950) to RADEC (J2000)" << LogIO::POST ;
378 converter_ = new MDirection::Convert( MDirection::B1950,
379 MDirection::J2000 ) ;
380 }
381 else if ( coord_ == 1 ) {
382 // LB -> RADEC (J2000)
383 os << "Creating converter from GALACTIC to RADEC (J2000)" << LogIO::POST ;
384 converter_ = new MDirection::Convert( MDirection::GALACTIC,
385 MDirection::J2000 ) ;
386 }
387 else {
388 // coord_ == 2
389 // AZEL -> RADEC (J2000)
390 os << "Creating converter from AZEL to RADEC (J2000)" << LogIO::POST ;
391 if ( mf_.null() ) {
392 mf_ = new MeasFrame() ;
393 vector<double> antpos = getAntennaPosition() ;
394 Vector<Quantity> qantpos( 3 ) ;
395 for ( int ip = 0 ; ip < 3 ; ip++ )
396 qantpos[ip] = Quantity( antpos[ip], "m" ) ;
397 mp_ = MPosition( MVPosition( qantpos ), MPosition::ITRF ) ;
398 mf_->set( mp_ ) ;
399 }
400 converter_ = new MDirection::Convert( MDirection::AZEL,
401 MDirection::Ref(MDirection::J2000,
402 *mf_ ) ) ;
403 }
404 }
405
406 if ( coord_ == 2 ) {
407 me_ = MEpoch( Quantity( t, "d" ), MEpoch::UTC ) ;
408 mf_->set( me_ ) ;
409 }
410}
411
412int NROReader::getHeaderInfo( Int &nchan,
413 Int &npol,
414 Int &nif,
415 Int &nbeam,
416 String &observer,
417 String &project,
418 String &obstype,
419 String &antname,
420 Vector<Double> &antpos,
421 Float &equinox,
422 String &freqref,
423 Double &reffreq,
424 Double &bw,
425 Double &utc,
426 String &fluxunit,
427 String &epoch,
428 String &poltype )
429{
430 nchan = dataset_->getNUMCH() ;
431 //cout << "nchan = " << nchan << endl ;
432 npol = getPolarizationNum() ;
433 //cout << "npol = " << npol << endl ;
434 observer = dataset_->getOBSVR() ;
435 //cout << "observer = " << observer << endl ;
436 project = dataset_->getPROJ() ;
437 //cout << "project = " << project << endl ;
438 obstype = dataset_->getSWMOD() ;
439 //cout << "obstype = " << obstype << endl ;
440 antname = dataset_->getSITE() ;
441 //cout << "antname = " << antname << endl ;
442 // TODO: should be investigated antenna position since there are
443 // no corresponding information in the header
444 // 2008/11/13 Takeshi Nakazato
445 //
446 // INFO: tentative antenna posiiton is obtained for NRO 45m from ITRF website
447 // 2008/11/26 Takeshi Nakazato
448 vector<double> pos = getAntennaPosition() ;
449 antpos = pos ;
450 //cout << "antpos = " << antpos << endl ;
451 string eq = dataset_->getEPOCH() ;
452// if ( eq.compare( 0, 5, "B1950" ) == 0 )
453// equinox = 1950.0 ;
454// else if ( eq.compare( 0, 5, "J2000" ) == 0 )
455// equinox = 2000.0 ;
456 // equinox is always 2000.0
457 equinox = 2000.0 ;
458 //cout << "equinox = " << equinox << endl ;
459 string vref = dataset_->getVREF() ;
460 if ( vref.compare( 0, 3, "LSR" ) == 0 ) {
461 if ( vref.size() == 3 ) {
462 vref.append( "K" ) ;
463 }
464 else {
465 vref[3] = 'K' ;
466 }
467 }
468 else if ( vref.compare( 0, 3, "GAL" ) == 0 ) {
469 vref = "GALACTO" ;
470 }
471 else if (vref.compare( 0, 3, "HEL" ) == 0 ) {
472 os_.origin( LogOrigin( "NROReader", "getHeaderInfo", WHERE ) ) ;
473 os_ << LogIO::WARN << "Heliocentric frame is not supported. Use Barycentric frame instead." << LogIO::POST ;
474 vref = "BARY" ;
475 }
476 //freqref = vref ;
477 //freqref = "LSRK" ;
478 //freqref = "REST" ;
479 freqref = freqRefFromVREF_ ? vref : "REST" ;
480 //cout << "freqref = " << freqref << endl ;
481 const NRODataRecord *record = dataset_->getRecord( 0 ) ;
482 reffreq = record->FREQ0 ; // rest frequency
483
484 //cout << "reffreq = " << reffreq << endl ;
485 bw = dataset_->getBEBW()[0] ;
486 //cout << "bw = " << bw << endl ;
487 utc = getStartTime() ;
488 //cout << "utc = " << utc << endl ;
489 fluxunit = "K" ;
490 //cout << "fluxunit = " << fluxunit << endl ;
491 epoch = "UTC" ;
492 //cout << "epoch = " << epoch << endl ;
493 string poltp = dataset_->getPOLTP()[0] ;
494 //cout << "poltp = '" << poltp << "'" << endl ;
495 if ( poltp == "" || poltp[0] == ' ' )
496 //poltp = "None" ;
497 poltp = "linear" ; // if no polarization type specified, set to "linear"
498 //else if ( strcmp( poltp, "LINR" ) == 0 )
499 else if ( poltp.compare( 0, 1, "LINR", 0, 1 ) == 0 )
500 poltp = "linear" ;
501 //else if ( strcmp( poltp, "CIRL" ) == 0 )
502 else if ( poltp.compare( 0, 1, "CIRL", 0, 1 ) == 0 )
503 poltp = "circular" ;
504 poltype = poltp ;
505 //cout << "poltype = " << poltype << endl ;
506
507 //vector<Bool> ifs = getIFs() ;
508 //nif = ifs.size() ;
509 nif = getNumIF() ;
510 //cout << "nif = " << nif << endl ;
511
512 //vector<Bool> beams = getBeams() ;
513 //nbeam = beams.size() ;
514 nbeam = getNumBeam() ;
515 //cout << "nbeam = " << nbeam << endl ;
516
517 return 0 ;
518}
519
520string NROReader::getScanType( int i )
521{
522 const NRODataRecord *record = dataset_->getRecord( i ) ;
523 string s = record->SCANTP ;
524
525 return s ;
526}
527
528int NROReader::getScanInfo( int irow,
529 uInt &scanno,
530 uInt &cycleno,
531 uInt &ifno,
532 uInt &beamno,
533 uInt &polno,
534 vector<double> &freqs,
535 Vector<Double> &restfreq,
536 uInt &refbeamno,
537 Double &scantime,
538 Double &interval,
539 String &srcname,
540 String &fieldname,
541 Vector<Float> &spectra,
542 Vector<uChar> &flagtra,
543 Vector<Float> &tsys,
544 Vector<Double> &direction,
545 Float &azimuth,
546 Float &elevation,
547 Float &parangle,
548 Float &opacity,
549 uInt &tcalid,
550 Int &fitid,
551 uInt &focusid,
552 Float &temperature,
553 Float &pressure,
554 Float &humidity,
555 Float &windvel,
556 Float &winddir,
557 Double &srcvel,
558 Vector<Double> &/*propermotion*/,
559 Vector<Double> &srcdir,
560 Vector<Double> &/*scanrate*/ )
561{
562 static const IPosition oneByOne( 1, 1 );
563
564 // DEBUG
565 //cout << "NROReader::getScanInfo() irow = " << irow << endl ;
566 //
567 NRODataRecord *record = dataset_->getRecord( irow ) ;
568
569 // scanno
570 scanno = (uInt)(record->ISCAN) ;
571 //cout << "scanno = " << scanno << endl ;
572
573 // cycleno
574 cycleno = 0 ;
575 //cout << "cycleno = " << cycleno << endl ;
576
577 // beamno and ifno
578 string rxname = dataset_->getRX()[0] ;
579 if ( rxname.find("MULT2") != string::npos ) {
580 string arryt = string( record->ARRYT ) ;
581 beamno = dataset_->getArrayId( arryt ) ;
582 ifno = 0 ;
583 }
584 else {
585 beamno = 0 ;
586 string arryt = string( record->ARRYT ) ;
587 ifno = dataset_->getArrayId( arryt ) ;
588 }
589 //cout << "beamno = " << beamno << endl ;
590
591 // polno
592 polno = dataset_->getPolNo( irow ) ;
593 //cout << "polno = " << polno << endl ;
594
595 // freqs (for IFNO and FREQ_ID)
596 //freqs = getFrequencies( irow ) ;
597 freqs = dataset_->getFrequencies( irow ) ;
598 //cout << "freqs = [" << freqs[0] << ", " << freqs[1] << ", " << freqs[2] << "]" << endl ;
599
600 if ( freqRefFromVREF_ ) {
601 freqs = shiftFrequency( freqs,
602 dataset_->getURVEL(),
603 dataset_->getVDEF() ) ;
604 }
605
606 // restfreq (for MOLECULE_ID)
607 restfreq.resize( oneByOne ) ;
608 restfreq[0] = record->FREQ0 ;
609 //cout << "restfreq = " << rf << endl ;
610
611 // refbeamno
612 refbeamno = 0 ;
613 //cout << "refbeamno = " << refbeamno << endl ;
614
615 // scantime
616 //scantime = Double( dataset_->getStartIntTime( irow ) ) ;
617 scantime = Double( dataset_->getScanTime( irow ) ) ;
618 //cout << "scantime = " << scantime << endl ;
619
620 // interval
621 interval = Double( dataset_->getIPTIM() ) ;
622 //cout << "interval = " << interval << endl ;
623
624 // srcname
625 srcname = String( dataset_->getOBJ() ) ;
626 //cout << "srcname = " << srcname << endl ;
627
628 // fieldname
629 fieldname = String( dataset_->getOBJ() ) ;
630 //cout << "fieldname = " << fieldname << endl ;
631
632 // spectra
633 vector<double> spec = dataset_->getSpectrum( irow ) ;
634 spectra.resize( spec.size() ) ;
635 //int index = 0 ;
636 Bool b ;
637 Float *fp = spectra.getStorage( b ) ;
638 Float *wp = fp ;
639 for ( vector<double>::iterator i = spec.begin() ;
640 i != spec.end() ; i++ ) {
641 *wp = *i ;
642 wp++ ;
643 }
644 spectra.putStorage( fp, b ) ;
645 //cout << "spec.size() = " << spec.size() << endl ;
646
647 // flagtra
648 bool setValue = !( flagtra.nelements() == spectra.nelements() ) ;
649 if ( setValue ) {
650 //cout << "flagtra resized. reset values..." << endl ;
651 flagtra.resize( spectra.nelements() ) ;
652 flagtra.set( 0 ) ;
653 }
654 //cout << "flag.size() = " << flag.size() << endl ;
655
656 // tsys
657 tsys.resize( oneByOne ) ;
658 tsys[0] = record->TSYS ;
659 //cout << "tsys[0] = " << tsys[0] << endl ;
660
661 // direction
662 direction = getDirection( irow ) ;
663 //cout << "direction = [" << direction[0] << ", " << direction[1] << "]" << endl ;
664
665 // azimuth
666 azimuth = record->RAZ ;
667 //cout << "azimuth = " << azimuth << endl ;
668
669 // elevation
670 elevation = record->REL ;
671 //cout << "elevation = " << elevation << endl ;
672
673 // parangle
674 parangle = 0.0 ;
675 //cout << "parangle = " << parangle << endl ;
676
677 // opacity
678 opacity = 0.0 ;
679 //cout << "opacity = " << opacity << endl ;
680
681 // tcalid
682 tcalid = 0 ;
683 //cout << "tcalid = " << tcalid << endl ;
684
685 // fitid
686 fitid = -1 ;
687 //cout << "fitid = " << fitid << endl ;
688
689 // focusid
690 focusid = 0 ;
691 //cout << "focusid = " << focusid << endl ;
692
693 // temperature (for WEATHER_ID)
694 temperature = Float( record->TEMP ) ;
695 //cout << "temperature = " << temperature << endl ;
696
697 // pressure (for WEATHER_ID)
698 pressure = Float( record->PATM ) ;
699 //cout << "pressure = " << pressure << endl ;
700
701 // humidity (for WEATHER_ID)
702 humidity = Float( record->PH2O ) ;
703 //cout << "humidity = " << humidity << endl ;
704
705 // windvel (for WEATHER_ID)
706 windvel = Float( record->VWIND ) ;
707 //cout << "windvel = " << windvel << endl ;
708
709 // winddir (for WEATHER_ID)
710 winddir = Float( record->DWIND ) ;
711 //cout << "winddir = " << winddir << endl ;
712
713 // srcvel
714 srcvel = dataset_->getURVEL() ;
715 //cout << "srcvel = " << srcvel << endl ;
716
717 // propermotion
718 // do nothing
719
720 // srcdir
721 srcdir = getSourceDirection() ;
722 //cout << "srcdir = [" << srcdir[0] << ", " << srcdir[1] << endl ;
723
724 // scanrate
725 // do nothing
726
727 return 0 ;
728}
729
730Int NROReader::getRowNum()
731{
732 return dataset_->getRowNum() ;
733}
734
735vector<double> NROReader::shiftFrequency( const vector<double> &f,
736 const double &v,
737 const string &vdef )
738{
739 vector<double> r( f ) ;
740 double factor = v / 2.99792458e8 ;
741 if ( vdef.compare( 0, 3, "RAD" ) == 0 ) {
742 factor = 1.0 / ( 1.0 + factor ) ;
743 r[1] *= factor ;
744 r[2] *= factor ;
745 }
746 else if ( vdef.compare( 0, 3, "OPT" ) == 0 ) {
747 factor += 1.0 ;
748 r[1] *= factor ;
749 r[2] *= factor ;
750 }
751 else {
752 cout << "vdef=" << vdef << " is not supported." << endl;
753 }
754 return r ;
755}
Note: See TracBrowser for help on using the repository browser.