Context Navigation

ASDMFiller.cc@ 2716

Visit:

Last change on this file since 2716 was 2710, 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...

Reference epoch for frequency frame conversion is taken from header
which is set as start time of the observation (OBSERVATION/TIME_RANGE[0]
for MS, ExecBlock/startTime for ASDM).

File size: 23.8 KB

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

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: trunk/external-alma/asdm2ASAP/ASDMFiller.cc@ 2716

Download in other formats: