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

Last change on this file since 2940 was 2940, checked in by Malte Marquarding, 10 years ago

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