- Timestamp:
- 07/06/12 18:30:00 (12 years ago)
- Location:
- branches/hpc34/external-alma/atnf/PKSIO
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/hpc34/external-alma/atnf/PKSIO/NROReader.cc
r2489 r2588 182 182 //---------------------------------------------------------------------- 183 183 // constructor 184 NROReader::NROReader( string name ) { 184 NROReader::NROReader( string name ) 185 : dataset_( NULL ), 186 srcdir_( 0 ), 187 msrcdir_( 0 ) 188 { 185 189 // initialization 186 190 filename_ = name ; 187 dataset_ = NULL;191 coord_ = -1 ; 188 192 } 189 193 … … 298 302 LogIO os( LogOrigin( "NROReader", "getSourceDirection()", WHERE ) ) ; 299 303 300 Vector<Double> v ; 301 Double srcra = Double( dataset_->getRA0() ) ; 302 Double srcdec = Double( dataset_->getDEC0() ) ; 304 if ( msrcdir_.nelements() == 2 ) 305 return msrcdir_ ; 306 307 srcdir_.resize( 2 ) ; 308 srcdir_[0] = Double( dataset_->getRA0() ) ; 309 srcdir_[1] = Double( dataset_->getDEC0() ) ; 303 310 char epoch[5] ; 304 311 strncpy( epoch, (dataset_->getEPOCH()).c_str(), 5 ) ; … … 306 313 // convert to J2000 value 307 314 MDirection result = 308 MDirection::Convert( MDirection( Quantity( src ra, "rad" ),309 Quantity( srcd ec, "rad" ),315 MDirection::Convert( MDirection( Quantity( srcdir_[0], "rad" ), 316 Quantity( srcdir_[1], "rad" ), 310 317 MDirection::Ref( MDirection::B1950 ) ), 311 318 MDirection::Ref( MDirection::J2000 ) ) () ; 312 v = result.getAngle().getValue() ; 313 Double srcra2 = v( 0 ) ; 314 if ( srcra2 < 0.0 && srcra >= 0.0 ) 315 v( 0 ) = 2.0 * M_PI + srcra2 ; 319 msrcdir_ = result.getAngle().getValue() ; 320 if ( msrcdir_[0] < 0.0 && srcdir_[0] >= 0.0 ) 321 msrcdir_[0] = 2.0 * M_PI + msrcdir_[0] ; 316 322 //cout << "NROReader::getSourceDirection() SRCDIRECTION convert from (" 317 323 //<< srcra << "," << srcdec << ") B1950 to (" … … 322 328 } 323 329 else if ( strncmp( epoch, "J2000", 5 ) == 0 ) { 324 v.resize( 2 ) ; 325 v( 0 ) = srcra ; 326 v( 1 ) = srcdec ; 327 } 328 329 return v ; 330 msrcdir_.reference( srcdir_ ) ; 331 } 332 333 return msrcdir_ ; 330 334 } 331 335 … … 335 339 LogIO os( LogOrigin( "NROReader", "getDirection()", WHERE ) ) ; 336 340 337 Vector<Double> v ;341 Vector<Double> v( 2 ) ; 338 342 NRODataRecord *record = dataset_->getRecord( i ) ; 339 343 char epoch[5] ; 340 344 strncpy( epoch, (dataset_->getEPOCH()).c_str(), 5 ) ; 341 345 int icoord = dataset_->getSCNCD() ; 342 Double dirx = Double( record->SCX ) ; 343 Double diry = Double( record->SCY ) ; 344 if ( icoord == 1 ) { 345 // convert from LB to RADEC 346 MDirection result = 347 MDirection::Convert( MDirection( Quantity( dirx, "rad" ), 348 Quantity( diry, "rad" ), 349 MDirection::Ref( MDirection::GALACTIC ) ), 350 MDirection::Ref( MDirection::J2000 ) ) () ; 351 v = result.getAngle().getValue() ; 352 Double dirx2 = v( 0 ) ; 353 if ( dirx2 < 0.0 && dirx >= 0.0 ) 354 v( 0 ) = 2.0 * M_PI + dirx2 ; 355 //cout << "NROReader::getDirection() DIRECTION convert from (" 356 //<< dirx << "," << diry << ") LB to (" 357 //<< v( 0 ) << ","<< v( 1 ) << ") RADEC" << endl ; 358 //os << LogIO::NORMAL << "DIRECTION convert from (" 359 // << dirx << "," << diry << ") LB to (" 360 // << v( 0 ) << ","<< v( 1 ) << ") RADEC" << LogIO::POST ; 361 } 362 else if ( icoord == 2 ) { 363 // convert from AZEL to RADEC 364 vector<double> antpos = getAntennaPosition() ; 365 Vector<Quantity> qantpos( 3 ) ; 366 for ( int ip = 0 ; ip < 3 ; ip++ ) 367 qantpos[ip] = Quantity( antpos[ip], "m" ) ; 368 Double scantime = Double( dataset_->getScanTime( i ) ) ; 369 MEpoch me( Quantity( scantime, "d" ), MEpoch::UTC ) ; 370 MPosition mp( MVPosition( qantpos ), MPosition::ITRF ) ; 371 MeasFrame mf( me, mp ) ; 372 MDirection result = 373 MDirection::Convert( MDirection( Quantity( dirx, "rad" ), 374 Quantity( diry, "rad" ), 375 MDirection::Ref( MDirection::AZEL ) ), 376 MDirection::Ref( MDirection::J2000, mf ) ) () ; 377 v = result.getAngle().getValue() ; 378 //cout << "NROReader::getDirection() DIRECTION convert from (" 379 //<< dirx << "," << diry << ") AZEL to (" 380 //<< v( 0 ) << ","<< v( 1 ) << ") RADEC" << endl ; 381 //os << LogIO::NORMAL << "DIRECTION convert from (" 382 // << dirx << "," << diry << ") AZEL to (" 383 // << v( 0 ) << ","<< v( 1 ) << ") RADEC" << LogIO::POST ; 384 } 385 else if ( icoord == 0 ) { 386 if ( strncmp( epoch, "B1950", 5 ) == 0 ) { 387 // convert to J2000 value 388 MDirection result = 389 MDirection::Convert( MDirection( Quantity( dirx, "rad" ), 390 Quantity( diry, "rad" ), 391 MDirection::Ref( MDirection::B1950 ) ), 392 MDirection::Ref( MDirection::J2000 ) ) () ; 393 v = result.getAngle().getValue() ; 394 Double dirx2 = v( 0 ) ; 395 if ( dirx2 < 0.0 && dirx >= 0.0 ) 396 v( 0 ) = 2.0 * M_PI + dirx2 ; 397 //cout << "STFiller::readNRO() DIRECTION convert from (" 398 //<< dirx << "," << diry << ") B1950 to (" 399 //<< v( 0 ) << ","<< v( 1 ) << ") J2000" << endl ; 400 //os << LogIO::NORMAL << "DIRECTION convert from (" 401 // << dirx << "," << diry << ") B1950 to (" 402 // << v( 0 ) << ","<< v( 1 ) << ") J2000" << LogIO::POST ; 403 } 404 else if ( strncmp( epoch, "J2000", 5 ) == 0 ) { 405 v.resize( 2 ) ; 406 v( 0 ) = dirx ; 407 v( 1 ) = diry ; 408 } 409 } 410 411 return v ; 346 double scantime = dataset_->getScanTime( i ) ; 347 initConvert( icoord, scantime, epoch ) ; 348 v[0]= Double( record->SCX ) ; 349 v[1] = Double( record->SCY ) ; 350 if ( converter_.null() ) 351 // no conversion 352 return v ; 353 else { 354 Vector<Double> v2 = (*converter_)( v ).getAngle().getValue() ; 355 if ( v2[0] < 0.0 && v[0] >= 0.0 ) 356 v2[0] = 2.0 * M_PI + v2[0] ; 357 return v2 ; 358 } 359 } 360 361 void NROReader::initConvert( int icoord, double t, char *epoch ) 362 { 363 if ( icoord == 0 && strncmp( epoch, "J2000", 5 ) == 0 ) 364 // no conversion 365 return ; 366 367 if ( converter_.null() || icoord != coord_ ) { 368 LogIO os( LogOrigin( "NROReader", "initConvert()", WHERE ) ) ; 369 coord_ = icoord ; 370 if ( coord_ == 0 ) { 371 // RADEC (B1950) -> RADEC (J2000) 372 os << "Creating converter from RADEC (B1950) to RADEC (J2000)" << LogIO::POST ; 373 converter_ = new MDirection::Convert( MDirection::B1950, 374 MDirection::J2000 ) ; 375 } 376 else if ( coord_ == 1 ) { 377 // LB -> RADEC (J2000) 378 os << "Creating converter from GALACTIC to RADEC (J2000)" << LogIO::POST ; 379 converter_ = new MDirection::Convert( MDirection::GALACTIC, 380 MDirection::J2000 ) ; 381 } 382 else { 383 // coord_ == 2 384 // AZEL -> RADEC (J2000) 385 os << "Creating converter from AZEL to RADEC (J2000)" << LogIO::POST ; 386 if ( mf_.null() ) { 387 mf_ = new MeasFrame() ; 388 vector<double> antpos = getAntennaPosition() ; 389 Vector<Quantity> qantpos( 3 ) ; 390 for ( int ip = 0 ; ip < 3 ; ip++ ) 391 qantpos[ip] = Quantity( antpos[ip], "m" ) ; 392 MPosition mp( MVPosition( qantpos ), MPosition::ITRF ) ; 393 mf_->set( mp ) ; 394 } 395 converter_ = new MDirection::Convert( MDirection::AZEL, 396 MDirection::Ref(MDirection::J2000, 397 *mf_ ) ) ; 398 } 399 } 400 401 if ( coord_ == 2 ) { 402 MEpoch me( Quantity( t, "d" ), MEpoch::UTC ) ; 403 mf_->set( me ) ; 404 } 412 405 } 413 406 -
branches/hpc34/external-alma/atnf/PKSIO/NROReader.h
r2289 r2588 48 48 #include <measures/Measures/MeasFrame.h> 49 49 #include <casa/Logging/LogIO.h> 50 #include <casa/Utilities/CountedPtr.h> 50 51 51 52 //#include <fitsio.h> … … 206 207 // Get DIRECTION in RADEC(J2000) 207 208 virtual Vector<Double> getDirection( int i ) ; 209 virtual void initConvert( int icoord, double t, char *epoch ) ; 208 210 209 211 // filename … … 212 214 // dataset 213 215 NRODataset *dataset_ ; 216 217 // source direction 218 Vector<Double> srcdir_ ; 219 Vector<Double> msrcdir_ ; 220 221 // for direction conversion 222 CountedPtr<MDirection::Convert> converter_ ; 223 CountedPtr<MeasFrame> mf_ ; 224 int coord_ ; 214 225 215 226 // Logger
Note:
See TracChangeset
for help on using the changeset viewer.