Context Navigation

MSFiller.cpp@ 1987

Visit:

Last change on this file since 1987 was 1987, checked in by Takeshi Nakazato, 15 years ago

New Development: No

JIRA Issue: Yes CAS-2718

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

Improved MSFiller/MSWriter.
Changed Scantable::get_antennaName() to fit with MSFiller/MSWriter.

File size: 59.7 KB

Rev	Line
[1974]	1	//
	2	// C++ Interface: MSFiller
	3	//
	4	// Description:
	5	//
	6	// This class is specific filler for MS format
	7	//
	8	// Takeshi Nakazato <takeshi.nakazato@nao.ac.jp>, (C) 2010
	9	//
	10	// Copyright: See COPYING file that comes with this distribution
	11	//
	12	//
	13
	14	#include <iostream>
	15	#include <map>
	16
	17	#include <tables/Tables/ExprNode.h>
	18	#include <tables/Tables/TableIter.h>
[1987]	19	#include <tables/Tables/TableColumn.h>
[1974]	20	#include <tables/Tables/ScalarColumn.h>
	21	#include <tables/Tables/ArrayColumn.h>
	22	#include <tables/Tables/RefRows.h>
[1987]	23	#include <tables/Tables/TableParse.h>
	24	#include <tables/Tables/RefRows.h>
[1974]	25
	26	#include <casa/Containers/Block.h>
	27	#include <casa/Logging/LogIO.h>
	28	#include <casa/Arrays/Slicer.h>
	29	#include <casa/Quanta/MVTime.h>
[1987]	30	#include <casa/OS/Path.h>
[1974]	31
	32	#include <measures/Measures/Stokes.h>
	33	#include <measures/Measures/MEpoch.h>
[1987]	34	#include <measures/Measures/MCEpoch.h>
	35	#include <measures/Measures/MFrequency.h>
	36	#include <measures/Measures/MCFrequency.h>
	37	#include <measures/Measures/MPosition.h>
	38	#include <measures/Measures/MCPosition.h>
	39	#include <measures/Measures/MDirection.h>
	40	#include <measures/Measures/MCDirection.h>
	41	#include <measures/Measures/MeasConvert.h>
[1974]	42	#include <measures/TableMeasures/ScalarMeasColumn.h>
[1987]	43	#include <measures/TableMeasures/ArrayMeasColumn.h>
	44	#include <measures/TableMeasures/ScalarQuantColumn.h>
	45	#include <measures/TableMeasures/ArrayQuantColumn.h>
[1974]	46
	47	#include <atnf/PKSIO/SrcType.h>
	48
	49	#include "MSFiller.h"
	50	#include "STHeader.h"
	51
	52	using namespace casa ;
	53	using namespace std ;
	54
	55	namespace asap {
	56	MSFiller::MSFiller( casa::CountedPtr<Scantable> stable )
	57	: table_( stable ),
[1987]	58	tablename_( "" ),
[1974]	59	antenna_( -1 ),
	60	getPt_( False ),
	61	isFloatData_( False ),
	62	isData_( False ),
	63	isDoppler_( False ),
	64	isFlagCmd_( False ),
	65	isFreqOffset_( False ),
	66	isHistory_( False ),
	67	isProcessor_( False ),
	68	isSysCal_( False ),
	69	isWeather_( False )
	70	{
	71	os_ = LogIO() ;
	72	os_.origin( LogOrigin( "MSFiller", "MSFiller()", WHERE ) ) ;
	73	}
	74
	75	MSFiller::~MSFiller()
	76	{
	77	os_.origin( LogOrigin( "MSFiller", "~MSFiller()", WHERE ) ) ;
	78	}
	79
	80	bool MSFiller::open( const std::string &filename, const casa::Record &rec )
	81	{
	82	os_.origin( LogOrigin( "MSFiller", "open()", WHERE ) ) ;
[1987]	83	//os_ << "start MSFiller::open()" << LogIO::POST ;
[1974]	84	//os_ << " filename = " << filename << endl ;
	85
	86	// parsing MS options
	87	if ( rec.isDefined( "ms" ) ) {
	88	Record msrec = rec.asRecord( "ms" ) ;
	89	if ( msrec.isDefined( "getpt" ) ) {
	90	getPt_ = msrec.asBool( "getpt" ) ;
	91	}
	92	if ( msrec.isDefined( "antenna" ) ) {
	93	if ( msrec.type( msrec.fieldNumber( "antenna" ) ) == TpInt ) {
	94	antenna_ = msrec.asInt( "antenna" ) ;
	95	}
	96	else {
	97	antenna_ = atoi( msrec.asString( "antenna" ).c_str() ) ;
	98	}
	99	}
	100	else {
	101	antenna_ = 0 ;
	102	}
	103	}
	104
	105	os_ << "Parsing MS options" << endl ;
	106	os_ << " getPt = " << getPt_ << endl ;
	107	os_ << " antenna = " << antenna_ << LogIO::POST ;
	108
[1987]	109	MeasurementSet *tmpMS = new MeasurementSet( filename, Table::Old ) ;
	110	//mstable_ = (*tmpMS)( tmpMS->col("ANTENNA1") == antenna_
	111	// && tmpMS->col("ANTENNA1") == tmpMS->col("ANTENNA2") ) ;
	112	tablename_ = tmpMS->tableName() ;
	113	mstable_ = MeasurementSet( (*tmpMS)( tmpMS->col("ANTENNA1") == antenna_
	114	&& tmpMS->col("ANTENNA1") == tmpMS->col("ANTENNA2") ) ) ;
	115	// stringstream ss ;
	116	// ss << "SELECT FROM $1 WHERE ANTENNA1 == ANTENNA2 && ANTENNA1 == " << antenna_ ;
	117	// String taql( ss.str() ) ;
	118	// mstable_ = MeasurementSet( tableCommand( taql, *tmpMS ) ) ;
	119	delete tmpMS ;
[1974]	120
	121	// check which data column exists
[1987]	122	isFloatData_ = mstable_.tableDesc().isColumn( "FLOAT_DATA" ) ;
	123	isData_ = mstable_.tableDesc().isColumn( "DATA" ) ;
[1974]	124
[1987]	125	//os_ << "end MSFiller::open()" << LogIO::POST ;
[1974]	126	return true ;
	127	}
	128
	129	void MSFiller::fill()
	130	{
	131	os_.origin( LogOrigin( "MSFiller", "fill()", WHERE ) ) ;
[1987]	132	//os_ << "start MSFiller::fill()" << LogIO::POST ;
[1974]	133
	134	// Initialize header
	135	STHeader sdh ;
	136	sdh.nchan = 0 ;
	137	sdh.npol = 0 ;
	138	sdh.nif = 0 ;
	139	sdh.nbeam = 0 ;
	140	sdh.observer = "" ;
	141	sdh.project = "" ;
	142	sdh.obstype = "" ;
	143	sdh.antennaname = "" ;
	144	sdh.antennaposition.resize( 0 ) ;
	145	sdh.equinox = 0.0 ;
	146	sdh.freqref = "" ;
	147	sdh.reffreq = -1.0 ;
	148	sdh.bandwidth = 0.0 ;
	149	sdh.utc = 0.0 ;
	150	sdh.fluxunit = "" ;
	151	sdh.epoch = "" ;
	152	sdh.poltype = "" ;
	153
	154	// check if optional table exists
[1987]	155	//const TableRecord msrec = tablesel_.keywordSet() ;
	156	const TableRecord msrec = mstable_.keywordSet() ;
[1974]	157	isDoppler_ = msrec.isDefined( "DOPPLER" ) ;
	158	isFlagCmd_ = msrec.isDefined( "FLAG_CMD" ) ;
	159	isFreqOffset_ = msrec.isDefined( "FREQ_OFFSET" ) ;
	160	isHistory_ = msrec.isDefined( "HISTORY" ) ;
	161	isProcessor_ = msrec.isDefined( "PROCESSOR" ) ;
	162	isSysCal_ = msrec.isDefined( "SYSCAL" ) ;
	163	isWeather_ = msrec.isDefined( "WEATHER" ) ;
	164
	165	// Access to MS subtables
[1987]	166	MSField fieldtab = mstable_.field() ;
	167	MSPolarization poltab = mstable_.polarization() ;
	168	MSDataDescription ddtab = mstable_.dataDescription() ;
	169	MSObservation obstab = mstable_.observation() ;
	170	MSSource srctab = mstable_.source() ;
	171	MSSpectralWindow spwtab = mstable_.spectralWindow() ;
	172	MSSysCal caltab = mstable_.sysCal() ;
[1974]	173	if ( caltab.nrow() == 0 )
	174	isSysCal_ = False ;
[1987]	175	MSPointing pointtab = mstable_.pointing() ;
	176	if ( mstable_.weather().nrow() == 0 )
[1974]	177	isWeather_ = False ;
[1987]	178	MSState stattab = mstable_.state() ;
	179	MSAntenna anttab = mstable_.antenna() ;
[1974]	180
[1987]	181	// TEST
	182	// memory allocation by boost::object_pool
	183	boost::object_pool<ROTableColumn> *tpoolr = new boost::object_pool<ROTableColumn> ;
	184	boost::object_pool<TableColumn> *tpoolw = new boost::object_pool<TableColumn> ;
	185	//
[1974]	186
[1987]	187	// SUBTABLES: FREQUENCIES
[1974]	188	table_->frequencies().setFrame( "LSRK" ) ;
	189	table_->frequencies().setFrame( "LSRK", True ) ;
	190
	191	// SUBTABLES: WEATHER
	192	if ( isWeather_ )
	193	fillWeather() ;
	194
	195	// SUBTABLES: FOCUS
	196	fillFocus() ;
	197
	198	// SUBTABLES: TCAL
	199	if ( isSysCal_ )
[1987]	200	fillTcal( tpoolr ) ;
[1974]	201
	202	// SUBTABLES: FIT
	203	//fillFit() ;
	204
	205	// SUBTABLES: HISTORY
	206	//fillHistory() ;
	207
[1987]	208	// shared pointers
	209	ROTableColumn *tcolr ;
	210	TableColumn *tcolw ;
	211
	212	// Scantable columns
	213	Table stab = table_->table() ;
	214	TableColumn *scannoCol = tpoolw->construct( stab, "SCANNO" ) ;
	215	TableColumn *cyclenoCol = tpoolw->construct( stab, "CYCLENO" ) ;
	216	TableColumn *beamnoCol = tpoolw->construct( stab, "BEAMNO" ) ;
	217	TableColumn *ifnoCol = tpoolw->construct( stab, "IFNO" ) ;
	218	TableColumn *polnoCol = tpoolw->construct( stab, "POLNO" ) ;
	219	TableColumn *freqidCol = tpoolw->construct( stab, "FREQ_ID" ) ;
	220	TableColumn *molidCol = tpoolw->construct( stab, "MOLECULE_ID" ) ;
	221	TableColumn *flagrowCol = tpoolw->construct( stab, "FLAGROW" ) ;
	222	ScalarMeasColumn<MEpoch> *timeCol = new ScalarMeasColumn<MEpoch>( stab, "TIME" ) ;
	223	TableColumn *intervalCol = tpoolw->construct( stab, "INTERVAL" ) ;
	224	TableColumn *srcnameCol = tpoolw->construct( stab, "SRCNAME" ) ;
	225	TableColumn *srctypeCol = tpoolw->construct( stab, "SRCTYPE" ) ;
	226	TableColumn *fieldnameCol = tpoolw->construct( stab, "SRCNAME" ) ;
	227	ArrayColumn<Float> *spCol = new ArrayColumn<Float>( stab, "SPECTRA" ) ;
	228	ArrayColumn<uChar> *flCol = new ArrayColumn<uChar>( stab, "FLAGTRA" ) ;
	229	ArrayColumn<Float> *tsysCol = new ArrayColumn<Float>( stab, "TSYS" ) ;
	230	ArrayColumn<Double> *dirCol = new ArrayColumn<Double>( stab, "DIRECTION" ) ;
	231	TableColumn *azCol = tpoolw->construct( stab, "AZIMUTH" ) ;
	232	TableColumn *elCol = tpoolw->construct( stab, "ELEVATION" ) ;
	233	TableColumn *tcalidCol = tpoolw->construct( stab, "TCAL_ID" ) ;
	234	TableColumn *focusidCol = tpoolw->construct( stab, "FOCUS_ID" ) ;
	235	TableColumn *weatheridCol = tpoolw->construct( stab, "WEATHER_ID" ) ;
	236	ArrayColumn<Double> *srcpmCol = new ArrayColumn<Double>( stab, "SRCPROPERMOTION" ) ;
	237	ArrayColumn<Double> *srcdirCol = new ArrayColumn<Double>( stab, "SRCDIRECTION" ) ;
	238	TableColumn *srcvelCol = tpoolw->construct( stab, "SRCVELOCITY" ) ;
	239	ArrayColumn<Double> *scanrateCol = new ArrayColumn<Double>( stab, "SCANRATE" ) ;
	240
[1974]	241	// MAIN
	242	// Iterate over several ids
	243	Int oldnr = table_->nrow() ;
	244	map<Int, uInt> ifmap ; // (IFNO, FREQ_ID) pair
[1987]	245	ROArrayQuantColumn<Double> *sharedQDArrCol = new ROArrayQuantColumn<Double>( anttab, "POSITION" ) ;
	246	Vector< Quantum<Double> > antpos = (*sharedQDArrCol)( antenna_ ) ;
	247	delete sharedQDArrCol ;
[1974]	248	MPosition mp( MVPosition( antpos ), MPosition::ITRF ) ;
[1987]	249	//pointtab = pointtab( pointtab.col("ANTENNA_ID")==antenna_ ).sort("TIME") ;
	250	pointtab = MSPointing( pointtab( pointtab.col("ANTENNA_ID")==antenna_ ).sort("TIME") ) ;
	251	tcolr = tpoolr->construct( anttab, "STATION" ) ;
	252	String stationName = tcolr->asString( antenna_ ) ;
	253	tpoolr->destroy( tcolr ) ;
	254	tcolr = tpoolr->construct( anttab, "NAME" ) ;
	255	String antennaName = tcolr->asString( antenna_ ) ;
	256	tpoolr->destroy( tcolr ) ;
	257	sdh.antennaposition.resize( 3 ) ;
	258	for ( int i = 0 ; i < 3 ; i++ )
	259	sdh.antennaposition[i] = antpos[i].getValue( "m" ) ;
	260	String telescopeName = "" ;
	261
[1978]	262	//
	263	// ITERATION: OBSERVATION_ID
	264	//
	265	Int added0 = 0 ;
	266	Int current0 = table_->nrow() ;
[1987]	267	TableIterator iter0( mstable_, "OBSERVATION_ID" ) ;
[1974]	268	while( !iter0.pastEnd() ) {
[1987]	269	Table t0 = iter0.table() ;
	270	tcolr = tpoolr->construct( t0, "OBSERVATION_ID" ) ;
	271	Int obsId = tcolr->asInt( 0 ) ;
	272	tpoolr->destroy( tcolr ) ;
	273	if ( sdh.observer == "" ) {
	274	tcolr = tpoolr->construct( obstab, "OBSERVER" ) ;
	275	sdh.observer = tcolr->asString( obsId ) ;
	276	tpoolr->destroy( tcolr ) ;
	277	}
	278	if ( sdh.project == "" ) {
	279	tcolr = tpoolr->construct( obstab, "PROJECT" ) ;
	280	sdh.observer = tcolr->asString( obsId ) ;
	281	tpoolr->destroy( tcolr ) ;
	282	}
	283	ROArrayMeasColumn<MEpoch> *tmpMeasCol = new ROArrayMeasColumn<MEpoch>( obstab, "TIME_RANGE" ) ;
	284	MEpoch me = (*tmpMeasCol)( obsId )( IPosition(1,0) ) ;
	285	delete tmpMeasCol ;
[1974]	286	if ( sdh.utc == 0.0 ) {
[1987]	287	sdh.utc = me.get( "s" ).getValue() ;
[1974]	288	}
[1987]	289	if ( telescopeName == "" ) {
	290	tcolr = tpoolr->construct( obstab, "TELESCOPE_NAME" ) ;
	291	sdh.observer = tcolr->asString( obsId ) ;
	292	tpoolr->destroy( tcolr ) ;
	293	}
[1978]	294	Int nbeam = 0 ;
[1974]	295	//
	296	// ITERATION: FEED1
	297	//
[1978]	298	Int added1 = 0 ;
	299	Int current1 = table_->nrow() ;
[1974]	300	TableIterator iter1( t0, "FEED1" ) ;
	301	while( !iter1.pastEnd() ) {
[1987]	302	Table t1 = iter1.table() ;
	303	// assume FEED1 == FEED2
	304	tcolr = tpoolr->construct( t1, "FEED1" ) ;
	305	Int feedId = tcolr->asInt( 0 ) ;
	306	tpoolr->destroy( tcolr ) ;
[1974]	307	nbeam++ ;
	308	//
	309	// ITERATION: FIELD_ID
	310	//
[1978]	311	Int added2 = 0 ;
	312	Int current2 = table_->nrow() ;
[1974]	313	TableIterator iter2( t1, "FIELD_ID" ) ;
	314	while( !iter2.pastEnd() ) {
[1987]	315	Table t2 = iter2.table() ;
	316	tcolr = tpoolr->construct( t2, "FIELD_ID" ) ;
	317	Int fieldId = tcolr->asInt( 0 ) ;
	318	tpoolr->destroy( tcolr ) ;
	319	tcolr = tpoolr->construct( fieldtab, "SOURCE_ID" ) ;
	320	Int srcId = tcolr->asInt( fieldId ) ;
	321	tpoolr->destroy( tcolr ) ;
	322	tcolr = tpoolr->construct( fieldtab, "NAME" ) ;
	323	String fieldName = tcolr->asString( fieldId ) + "__" + String::toString(fieldId) ;
	324	tpoolr->destroy( tcolr ) ;
[1974]	325	//
	326	// ITERATION: DATA_DESC_ID
	327	//
[1978]	328	Int added3 = 0 ;
	329	Int current3 = table_->nrow() ;
[1974]	330	TableIterator iter3( t2, "DATA_DESC_ID" ) ;
	331	while( !iter3.pastEnd() ) {
[1987]	332	Table t3 = iter3.table() ;
	333	tcolr = tpoolr->construct( t3, "DATA_DESC_ID" ) ;
	334	Int ddId = tcolr->asInt( 0 ) ;
	335	tpoolr->destroy( tcolr ) ;
	336	tcolr = tpoolr->construct( ddtab, "POLARIZATION_ID" ) ;
	337	Int polId = tcolr->asInt( ddId ) ;
	338	tpoolr->destroy( tcolr ) ;
	339	tcolr = tpoolr->construct( ddtab, "SPECTRAL_WINDOW_ID" ) ;
	340	Int spwId = tcolr->asInt( ddId ) ;
	341	tpoolr->destroy( tcolr ) ;
[1974]	342	// polarization information
[1987]	343	tcolr = tpoolr->construct( poltab, "NUM_CORR" ) ;
	344	Int npol = tcolr->asInt( polId ) ;
	345	tpoolr->destroy( tcolr ) ;
	346	ROArrayColumn<Int> *roArrICol = new ROArrayColumn<Int>( poltab, "CORR_TYPE" ) ;
	347	Vector<Int> corrtype = (*roArrICol)( polId ) ;
	348	delete roArrICol ;
[1974]	349	//os_ << "npol = " << npol << LogIO::POST ;
	350	//os_ << "corrtype = " << corrtype << LogIO::POST ;
	351	if ( sdh.npol < npol ) sdh.npol = npol ;
	352	if ( sdh.poltype == "" ) sdh.poltype = getPolType( corrtype[0] ) ;
	353	// source information
	354	MSSource srctabSel( srctab( srctab.col("SOURCE_ID") == srcId && srctab.col("SPECTRAL_WINDOW_ID") == spwId ) ) ;
	355	//os_ << "srcId = " << srcId << " spwId = " << spwId << " nrow = " << srctabSel.nrow() << LogIO::POST ;
[1987]	356	tcolr = tpoolr->construct( srctabSel, "NAME" ) ;
	357	String srcName = tcolr->asString( 0 ) ;
	358	tpoolr->destroy( tcolr ) ;
[1974]	359	//os_ << "srcName = " << srcName << LogIO::POST ;
[1987]	360	ROArrayColumn<Double> *roArrDCol = new ROArrayColumn<Double>( srctabSel, "PROPER_MOTION" ) ;
	361	Array<Double> srcPM = (*roArrDCol)( 0 ) ;
	362	delete roArrDCol ;
[1974]	363	//os_ << "srcPM = " << srcPM << LogIO::POST ;
[1987]	364	roArrDCol = new ROArrayColumn<Double>( srctabSel, "DIRECTION" ) ;
	365	Array<Double> srcDir = (*roArrDCol)( 0 ) ;
	366	delete roArrDCol ;
[1974]	367	//os_ << "srcDir = " << srcDir << LogIO::POST ;
[1987]	368	Array<Double> sysVels ;
[1974]	369	Double sysVel = 0.0 ;
[1987]	370	if ( srctabSel.tableDesc().isColumn( "SYSVEL" ) ) {
	371	roArrDCol = new ROArrayColumn<Double>( srctabSel, "SYSVEL" ) ;
	372	sysVels = (*roArrDCol)( 0 ) ;
	373	delete roArrDCol ;
	374	}
[1974]	375	if ( !sysVels.empty() ) {
	376	//os_ << "sysVels.shape() = " << sysVels.shape() << LogIO::POST ;
[1975]	377	// NB: assume all SYSVEL values are the same
[1987]	378	sysVel = sysVels( IPosition(1,0) ) ;
[1974]	379	}
[1987]	380	//delete tmpArrCol ;
[1974]	381	//os_ << "sysVel = " << sysVel << LogIO::POST ;
[1987]	382	ROScalarMeasColumn<MDirection> *tmpMeasCol = new ROScalarMeasColumn<MDirection>( srctabSel, "DIRECTION" ) ;
	383	MDirection md = (*tmpMeasCol)( 0 ) ;
	384	delete tmpMeasCol ;
[1974]	385	// for MOLECULES subtable
[1987]	386	tcolr = tpoolr->construct( srctabSel, "NUM_LINES" ) ;
	387	Int numLines = tcolr->asInt( 0 ) ;
	388	tpoolr->destroy( tcolr ) ;
[1974]	389	//os_ << "numLines = " << numLines << LogIO::POST ;
	390	Vector<Double> restFreqs( numLines, 0.0 ) ;
	391	Vector<String> transitionName( numLines, "" ) ;
	392	if ( numLines != 0 ) {
[1987]	393	if ( srctabSel.tableDesc().isColumn( "REST_FREQUENCY" ) ) {
	394	sharedQDArrCol = new ROArrayQuantColumn<Double>( srctabSel, "REST_FREQUENCY" ) ;
	395	Array< Quantum<Double> > qRestFreqs = (*sharedQDArrCol)( 0 ) ;
	396	delete sharedQDArrCol ;
	397	for ( int i = 0 ; i < numLines ; i++ ) {
	398	restFreqs[i] = qRestFreqs( IPosition( 1, i ) ).getValue( "Hz" ) ;
	399	}
[1974]	400	}
	401	//os_ << "restFreqs = " << restFreqs << LogIO::POST ;
	402	if ( srctabSel.tableDesc().isColumn( "TRANSITION" ) ) {
[1987]	403	tcolr = tpoolr->construct( srctabSel, "TRANSITION" ) ;
	404	transitionName = tcolr->asString( 0 ) ;
	405	tpoolr->destroy( tcolr ) ;
[1974]	406	//os_ << "transitionNameCol.nrow() = " << transitionNameCol.nrow() << LogIO::POST ;
	407	}
	408	}
	409	uInt molId = table_->molecules().addEntry( restFreqs, transitionName, transitionName ) ;
	410	// spectral setup
	411	MeasFrame mf( me, mp, md ) ;
[1987]	412	tcolr = tpoolr->construct( spwtab, "MEAS_FREQ_REF" ) ;
	413	MFrequency::Types freqRef = MFrequency::castType((uInt)(tcolr->asInt(spwId))) ;
	414	tpoolr->destroy( tcolr ) ;
	415	tcolr = tpoolr->construct( spwtab, "NUM_CHAN" ) ;
	416	Int nchan = tcolr->asInt( spwId ) ;
	417	tpoolr->destroy( tcolr ) ;
[1974]	418	Bool even = False ;
	419	if ( (nchan/2)*2 == nchan ) even = True ;
	420	if ( sdh.nchan < nchan ) sdh.nchan = nchan ;
[1987]	421	ROScalarQuantColumn<Double> *tmpQuantCol = new ROScalarQuantColumn<Double>( spwtab, "TOTAL_BANDWIDTH" ) ;
	422	Double totbw = (*tmpQuantCol)( spwId ).getValue( "Hz" ) ;
	423	delete tmpQuantCol ;
[1974]	424	if ( sdh.bandwidth < totbw ) sdh.bandwidth = totbw ;
	425	if ( sdh.freqref == "" )
	426	//sdh.freqref = MFrequency::showType( freqRef ) ;
	427	sdh.freqref = "LSRK" ;
	428	if ( sdh.reffreq == -1.0 ) {
[1987]	429	tmpQuantCol = new ROScalarQuantColumn<Double>( spwtab, "REF_FREQUENCY" ) ;
	430	Quantum<Double> qreffreq = (*tmpQuantCol)( spwId ) ;
	431	delete tmpQuantCol ;
[1974]	432	if ( freqRef == MFrequency::LSRK ) {
	433	sdh.reffreq = qreffreq.getValue("Hz") ;
	434	}
	435	else {
	436	MFrequency::Convert tolsr( freqRef, MFrequency::Ref( MFrequency::LSRK, mf ) ) ;
	437	sdh.reffreq = tolsr( qreffreq ).get("Hz").getValue() ;
	438	}
	439	}
	440	Int refchan = nchan / 2 ;
	441	IPosition refip( 1, refchan ) ;
	442	Double refpix = 0.5*(nchan-1) ;
	443	Double refval = 0.0 ;
[1987]	444	sharedQDArrCol = new ROArrayQuantColumn<Double>( spwtab, "CHAN_WIDTH" ) ;
	445	Double increment = (*sharedQDArrCol)( spwId )( refip ).getValue( "Hz" ) ;
	446	delete sharedQDArrCol ;
[1974]	447	//os_ << "nchan = " << nchan << " refchan = " << refchan << "(even=" << even << ") refpix = " << refpix << LogIO::POST ;
[1987]	448	sharedQDArrCol = new ROArrayQuantColumn<Double>( spwtab, "CHAN_FREQ" ) ;
	449	Vector< Quantum<Double> > chanFreqs = (*sharedQDArrCol)( spwId ) ;
	450	delete sharedQDArrCol ;
[1974]	451	if ( freqRef == MFrequency::LSRK ) {
	452	if ( even ) {
	453	IPosition refip0( 1, refchan-1 ) ;
	454	Double refval0 = chanFreqs(refip0).getValue("Hz") ;
	455	Double refval1 = chanFreqs(refip).getValue("Hz") ;
	456	refval = 0.5 * ( refval0 + refval1 ) ;
	457	}
	458	else {
	459	refval = chanFreqs(refip).getValue("Hz") ;
	460	}
	461	}
	462	else {
	463	MFrequency::Convert tolsr( freqRef, MFrequency::Ref( MFrequency::LSRK, mf ) ) ;
	464	if ( even ) {
	465	IPosition refip0( 1, refchan-1 ) ;
	466	Double refval0 = chanFreqs(refip0).getValue("Hz") ;
	467	Double refval1 = chanFreqs(refip).getValue("Hz") ;
	468	refval = 0.5 * ( refval0 + refval1 ) ;
	469	refval = tolsr( refval ).get("Hz").getValue() ;
	470	}
	471	else {
	472	refval = tolsr( chanFreqs(refip) ).get("Hz").getValue() ;
	473	}
	474	}
	475	uInt freqId = table_->frequencies().addEntry( refpix, refval, increment ) ;
	476	if ( ifmap.find( spwId ) == ifmap.end() ) {
	477	ifmap.insert( pair<Int, uInt>(spwId,freqId) ) ;
	478	//os_ << "added to ifmap: (" << spwId << "," << freqId << ")" << LogIO::POST ;
	479	}
	480	// for TSYS and TCAL
	481	MSSysCal caltabsel( caltab( caltab.col("ANTENNA_ID") == antenna_ && caltab.col("FEED_ID") == feedId && caltab.col("SPECTRAL_WINDOW_ID") == spwId ).sort("TIME") ) ;
	482	//
	483	// ITERATION: SCAN_NUMBER
	484	//
[1978]	485	Int added4 = 0 ;
	486	Int current4 = table_->nrow() ;
[1974]	487	TableIterator iter4( t3, "SCAN_NUMBER" ) ;
	488	while( !iter4.pastEnd() ) {
[1987]	489	Table t4 = iter4.table() ;
	490	tcolr = tpoolr->construct( t4, "SCAN_NUMBER" ) ;
	491	Int scanNum = tcolr->asInt( 0 ) ;
	492	tpoolr->destroy( tcolr ) ;
[1974]	493	uInt cycle = 0 ;
	494	//
	495	// ITERATION: STATE_ID
	496	//
[1978]	497	Int added5 = 0 ;
	498	Int current5 = table_->nrow() ;
[1974]	499	TableIterator iter5( t4, "STATE_ID" ) ;
	500	while( !iter5.pastEnd() ) {
[1987]	501	Table t5 = iter5.table() ;
	502	tcolr = tpoolr->construct( t5, "STATE_ID" ) ;
	503	Int stateId = tcolr->asInt( 0 ) ;
	504	tpoolr->destroy( tcolr ) ;
	505	tcolr = tpoolr->construct( stattab, "OBS_MODE" ) ;
	506	String obstype = tcolr->asString( stateId ) ;
	507	tpoolr->destroy( tcolr ) ;
[1974]	508	if ( sdh.obstype == "" ) sdh.obstype = obstype ;
	509
	510	Int nrow = t5.nrow() ;
	511	Int prevnr = oldnr ;
	512	Int addednr = 0 ;
[1987]	513	Int nloop = 0 ;
[1974]	514
	515	// SPECTRA, FLAG
	516	ROArrayColumn<Bool> mFlagCol( t5, "FLAG" ) ;
	517	if ( isFloatData_ ) {
	518	//os_ << "FLOAT_DATA exists" << LogIO::POST ;
	519	ROArrayColumn<Float> mFloatDataCol( t5, "FLOAT_DATA" ) ;
[1987]	520	addednr = nrow*npol ;
	521	oldnr += addednr ;
	522	table_->table().addRow( addednr ) ;
	523	nloop = npol ;
	524	for ( Int irow = 0 ; irow < nrow ; irow++ ) {
	525	Matrix<Float> sp = mFloatDataCol( irow ) ;
	526	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
	527	spCol->put( prevnr+ipol*nrow+irow, sp.row(ipol) ) ;
	528	}
[1974]	529	}
[1987]	530	for ( Int irow = 0 ; irow < nrow ; irow++ ) {
	531	Matrix<Bool> flb = mFlagCol( irow ) ;
	532	Matrix<uChar> fl( flb.shape() ) ;
	533	convertArray( fl, flb ) ;
	534	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
	535	flCol->put( prevnr+ipol*nrow+irow, fl.row(ipol) ) ;
	536	}
	537	}
[1974]	538	if ( sdh.fluxunit == "" ) {
	539	const TableRecord dataColKeys = mFloatDataCol.keywordSet() ;
	540	if ( dataColKeys.isDefined( "UNIT" ) )
	541	sdh.fluxunit = dataColKeys.asString( "UNIT" ) ;
	542	}
	543	}
	544	else if ( isData_ ) {
	545	//os_ << "DATA exists" << LogIO::POST ;
	546	ROArrayColumn<Complex> mDataCol( t5, "DATA" ) ;
[1987]	547	addednr = nrow*npol ;
	548	oldnr += addednr ;
	549	table_->table().addRow( addednr ) ;
	550	nloop = npol ;
	551	for ( Int irow = 0 ; irow < nrow ; irow++ ) {
	552	Bool crossOK = False ;
	553	Matrix<Complex> sp = mDataCol( irow ) ;
	554	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
	555	if ( corrtype[ipol] == Stokes::XY \|\| corrtype[ipol] == Stokes::YX
	556	\|\| corrtype[ipol] == Stokes::RL \|\| corrtype[ipol] == Stokes::LR ) {
	557	if ( !crossOK ) {
	558	crossOK = True ;
	559	Int pidx = prevnr + ipol * nrow + irow ;
	560	spCol->put( pidx, real( sp.row(ipol) ) ) ;
	561	if ( corrtype[ipol] == Stokes::XY \|\| corrtype[ipol] == Stokes::RL )
	562	spCol->put( pidx+nrow, imag( sp.row(ipol) ) ) ;
	563	else
	564	spCol->put( pidx+nrow, imag( conj(sp.row(ipol)) ) ) ;
[1974]	565	}
	566	}
[1987]	567	else {
	568	spCol->put( prevnr+ipol*nrow+irow, real( sp.row(ipol) ) ) ;
	569	}
[1974]	570	}
[1987]	571	}
	572	for ( Int irow = 0 ; irow < nrow ; irow++ ) {
	573	Bool crossOK = False ;
	574	Matrix<Bool> flb = mFlagCol( irow ) ;
	575	Matrix<uChar> fl( flb.shape() ) ;
	576	convertArray( fl, flb ) ;
	577	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
	578	if ( corrtype[ipol] == Stokes::XY \|\| corrtype[ipol] == Stokes::YX
	579	\|\| corrtype[ipol] == Stokes::RL \|\| corrtype[ipol] == Stokes::LR ) {
	580	if ( !crossOK ) {
	581	crossOK = True ;
	582	Int pidx = prevnr + ipol * nrow + irow ;
	583	flCol->put( pidx, fl.row(ipol) ) ;
	584	flCol->put( pidx+nrow, fl.row(ipol+1) ) ;
	585	}
	586	}
	587	else {
	588	flCol->put( prevnr+ipol*nrow+irow, fl.row(ipol) ) ;
	589	}
[1974]	590	}
	591	}
	592	if ( sdh.fluxunit == "" ) {
	593	const TableRecord dataColKeys = mDataCol.keywordSet() ;
	594	if ( dataColKeys.isDefined( "UNIT" ) )
	595	sdh.fluxunit = dataColKeys.asString( "UNIT" ) ;
	596	}
	597	}
	598
	599	// number of rows added in this cycle
	600	//os_ << "prevnr = " << prevnr << LogIO::POST ;
	601	//os_ << "addednr = " << addednr << LogIO::POST ;
	602	RefRows rows( prevnr, prevnr+addednr-1 ) ;
	603
	604	// TIME
[1987]	605	ROScalarMeasColumn<MEpoch> *mTimeCol = new ROScalarMeasColumn<MEpoch>( t5, "TIME" ) ;
[1974]	606	Int tidx = prevnr ;
[1987]	607	for ( Int i = 0 ; i < nloop ; i++ ) {
[1974]	608	for ( Int j = 0 ; j < nrow ; j++ ) {
[1987]	609	timeCol->put( tidx++, (*mTimeCol)( j ) ) ;
[1974]	610	}
	611	}
	612
	613	// TSYS
	614	Vector<Double> sysCalTime ;
	615	if ( isSysCal_ ) {
[1987]	616	sysCalTime = getSysCalTime( caltabsel, *mTimeCol ) ;
[1974]	617	tidx = prevnr ;
	618	uInt calidx = 0 ;
	619	for ( Int i = 0 ; i < nrow ; i++ ) {
[1987]	620	Matrix<Float> tsys ;
[1974]	621	calidx = getTsys( calidx, tsys, caltabsel, sysCalTime(i) ) ;
	622	//os_ << "tsys = " << tsys << LogIO::POST ;
[1987]	623	uInt ncol = tsys.ncolumn() ;
	624	if ( ncol == 0 ) {
	625	IPosition cShape = IPosition( 2, npol, 1 ) ;
[1974]	626	tsys.resize( cShape ) ;
	627	tsys = 1.0 ;
	628	}
	629	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	630	//floatArrCol->put( prevnr+i+nrow*ipol, tsys.row( ipol ) ) ;
	631	tsysCol->put( prevnr+i+nrow*ipol, tsys.row( ipol ) ) ;
[1974]	632	}
	633	}
	634	}
	635	else {
[1987]	636	Vector<Float> tsys( 1, 1.0 ) ;
	637	for ( Int i = prevnr ; i < prevnr+addednr ; i++ )
	638	tsysCol->put( i, tsys ) ;
[1974]	639	}
	640
[1987]	641
	642	// INTERVAL
	643	tcolr = tpoolr->construct( t5, "INTERVAL" ) ;
	644	//Vector<Double> integ = mIntervalCol->getColumn() ;
	645	for ( int i = 0 ; i < nloop ; i++ ) {
	646	//Int startrow = prevnr + i ;
	647	//Int endrow = startrow + nrow - 1 ;
	648	//RefRows prows( startrow, endrow ) ;
	649	//intervalCol->putColumnCells( prows, integ ) ;
	650	for ( int j = 0 ; j < nrow ; j++ ) {
	651	intervalCol->putScalar( prevnr+i*nrow+j, (Double)(tcolr->asdouble( j )) ) ;
	652	}
	653	}
	654	tpoolr->destroy( tcolr ) ;
	655
	656	// SRCTYPE
	657	Int srcType = getSrcType( stateId, tpoolr ) ;
	658	for ( int i = 0 ; i < addednr ; i++ ) {
	659	srctypeCol->putScalar( prevnr+i, srcType ) ;
	660	}
	661	//Vector<Int> *srcType = new Vector<Int>( addednr, getSrcType( stateId ) ) ;
	662	//srcTypeCol->putColumnCells( rows, *srcType ) ;
	663
[1974]	664	// DIRECTION, AZIMUTH, ELEVATION, SCANRATE
	665	Vector<Double> defaultScanrate( 2, 0.0 ) ;
	666	uInt diridx = 0 ;
	667	MDirection::Types dirType ;
	668	if ( getPt_ ) {
	669	for ( Int i = 0 ; i < nrow ; i++ ) {
	670	Vector<Double> dir ;
	671	Vector<Double> scanrate ;
	672	String refString ;
[1987]	673	diridx = getDirection( diridx, dir, scanrate, refString, pointtab, (*mTimeCol)(i).get("s").getValue() ) ;
[1974]	674	//os_ << "diridx = " << diridx << " dmTimeCol(" << i << ") = " << mTimeCol(i).get("s").getValue()-mTimeCol(0).get("s").getValue() << LogIO::POST ;
	675	//os_ << "dir = " << dir << LogIO::POST ;
	676	//os_ << "scanrate = " << scanrate << LogIO::POST ;
	677	//os_ << "refString = " << refString << LogIO::POST ;
	678	MDirection::getType( dirType, refString ) ;
	679	//os_ << "dirType = " << dirType << LogIO::POST ;
[1987]	680	mf.resetEpoch( (*mTimeCol)(i) ) ;
[1974]	681	mf.resetDirection( MDirection( MVDirection(dir), dirType ) ) ;
	682	if ( refString == "J2000" ) {
	683	//os_ << "J2000" << LogIO::POST ;
	684	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	685	dirCol->put( prevnr+i+nrow*ipol, dir ) ;
[1974]	686	}
	687	MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ;
	688	Vector<Double> azel = toazel( dir ).getAngle("rad").getValue() ;
	689	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	690	azCol->putScalar( prevnr+i+nrow*ipol, (Float)azel(0) ) ;
	691	elCol->putScalar( prevnr+i+nrow*ipol, (Float)azel(1) ) ;
[1974]	692	}
	693	}
	694	else if ( refString(0,4) == "AZEL" ) {
	695	//os_ << "AZEL" << LogIO::POST ;
	696	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	697	azCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(0) ) ;
	698	elCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(1) ) ;
[1974]	699	}
	700	MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ;
	701	Vector<Double> newdir = toj2000( dir ).getAngle("rad").getValue() ;
	702	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	703	dirCol->put( prevnr+i+nrow*ipol, newdir ) ;
[1974]	704	}
	705	}
	706	else {
	707	//os_ << "OTHER: " << refString << LogIO::POST ;
	708	MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ;
	709	Vector<Double> azel = toazel( dir ).getAngle("rad").getValue() ;
	710	MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ;
	711	Vector<Double> newdir = toj2000( dir ).getAngle("rad").getValue() ;
	712	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	713	dirCol->put( prevnr+i+nrow*ipol, newdir ) ;
	714	azCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(0) ) ;
	715	elCol->putScalar( prevnr+i+nrow*ipol, (Float)dir(1) ) ;
[1974]	716	}
	717	}
	718	if ( scanrate.size() != 0 ) {
	719	//os_ << "scanrate.size() = " << scanrate.size() << LogIO::POST ;
	720	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	721	scanrateCol->put( prevnr+i+nrow*ipol, scanrate ) ;
[1974]	722	}
	723	}
	724	else {
	725	//os_ << "scanrate.size() = " << scanrate.size() << LogIO::POST ;
	726	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	727	scanrateCol->put( prevnr+i+nrow*ipol, defaultScanrate ) ;
[1974]	728	}
	729	}
	730	}
	731	}
	732	else {
	733	// All directions are set to source direction
[1987]	734	ROArrayMeasColumn<MDirection> dmcol( pointtab, "DIRECTION" ) ;
	735	ROArrayColumn<Double> dcol( pointtab, "DIRECTION" ) ;
[1974]	736	IPosition ip( dmcol(0).shape().nelements(), 0 ) ;
	737	IPosition outp( 1, 2 ) ;
	738	String ref = dmcol(0)(ip).getRefString() ;
	739	Slice ds( 0, 2, 1 ) ;
	740	Slice ds0( 0, 1, 1 ) ;
	741	Slicer dslice0( ds, ds0 ) ;
	742	Vector<Double> defaultDir = dcol(0)(dslice0).reform(outp) ;
	743	MDirection::getType( dirType, "J2000" ) ;
	744	mf.resetDirection( MDirection( MVDirection(srcDir), dirType ) ) ;
	745	if ( ref != "J2000" ) {
[1987]	746	ROScalarMeasColumn<MEpoch> tmCol( pointtab, "TIME" ) ;
	747	mf.resetEpoch( tmCol( 0 ) ) ;
[1974]	748	MDirection::Convert toj2000( dirType, MDirection::Ref( MDirection::J2000, mf ) ) ;
	749	defaultDir = toj2000( defaultDir ).getAngle("rad").getValue() ;
	750	}
	751	for ( Int i = 0 ; i < nrow ; i++ ) {
[1987]	752	mf.resetEpoch( (*mTimeCol)(i) ) ;
[1974]	753	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
	754	Int localidx = prevnr+i+nrow*ipol ;
	755	MDirection::Convert toazel( dirType, MDirection::Ref( MDirection::AZEL, mf ) ) ;
	756	Vector<Double> azel = toazel( defaultDir ).getAngle("rad").getValue() ;
[1987]	757	azCol->putScalar( localidx, (Float)azel(0) ) ;
	758	elCol->putScalar( localidx, (Float)azel(1) ) ;
	759	dirCol->put( localidx, defaultDir ) ;
	760	scanrateCol->put( localidx, defaultScanrate ) ;
[1974]	761	}
	762	}
	763	}
	764
[1987]	765	// CYCLENO
	766	for ( int i = 0 ; i < nloop ; i++ ) {
	767	for ( int j = 0 ; j < nrow ; j++ ) {
	768	cyclenoCol->putScalar( prevnr+nrow*i+j, cycle+j ) ;
	769	}
	770	}
	771	cycle += nrow ;
	772
	773	// POLNO
	774	Int pidx = 0 ;
	775	Bool crossOK = False ;
	776	for ( int i = 0 ; i < npol ; i++ ) {
	777	Vector<uInt> polnos = getPolNo( corrtype[i] ) ;
	778	if ( polnos.size() > 1 ) {
	779	if ( crossOK ) continue ;
	780	else crossOK = True ;
	781	}
	782	for ( uInt j = 0 ; j < polnos.size() ; j++ ) {
	783	for ( Int irow = 0 ; irow < nrow ; irow++ ) {
	784	polnoCol->putScalar( prevnr+pidx*nrow+irow, polnos[j] ) ;
	785	}
	786	pidx++ ;
	787	}
	788	}
	789
	790	// FLAGROW
	791	tcolr = tpoolr->construct( t5, "FLAG_ROW" ) ;
	792	for ( int i = 0 ; i < nloop ; i++ ) {
	793	for ( int j = 0 ; j < nrow ; j++ ) {
	794	flagrowCol->putScalar( prevnr+nrow*i+j, (uInt)(tcolr->asBool( j )) ) ;
	795	}
	796	}
	797	tpoolr->destroy( tcolr ) ;
	798
[1974]	799	// TCAL_ID
	800	if ( isSysCal_ ) {
	801	for( Int irow = 0 ; irow < nrow ; irow++ ) {
	802	Vector<uInt> tcalids = getTcalId( feedId, spwId, sysCalTime[irow] ) ;
	803	if ( tcalids.size() == 0 ) {
	804	tcalids.resize( npol ) ;
	805	tcalids = 0 ;
	806	}
	807	for ( Int ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	808	tcalidCol->putScalar( prevnr+irow+nrow*ipol, tcalids[ipol] ) ;
[1974]	809	}
	810	}
	811	}
	812	else {
[1987]	813	//Vector<uInt> tcalid( addednr, 0 ) ;
	814	//uIntCol->putColumnCells( rows, tcalid ) ;
	815	uInt tcalid = 0 ;
	816	for ( Int irow = 0 ; irow < addednr ; irow++ )
	817	tcalidCol->putScalar( prevnr+irow, tcalid ) ;
[1974]	818	}
	819
	820	// WEATHER_ID
	821	if ( isWeather_ ) {
[1987]	822	uInt wid = 0 ;
	823	for ( Int irow = 0 ; irow < nrow ; irow++ ) {
	824	wid = getWeatherId( wid, (*mTimeCol)(irow).get("s").getValue() ) ;
	825	for ( Int ipol = 0 ; ipol < nloop ; ipol++ ) {
	826	weatheridCol->putScalar( prevnr+ipol*nrow+irow, wid ) ;
	827	}
[1974]	828	}
	829	}
	830
[1987]	831	delete mTimeCol ;
[1974]	832
	833	//os_ << "field: " << fieldId << " scan: " << scanNum << " obs: " << obsId << " state: " << stateId << " ddid: " << ddId << endl ;
[1987]	834	//os_ << "addednr = " << addednr << endl ;
[1978]	835	added5 += addednr ;
[1974]	836	iter5.next() ;
	837	}
[1978]	838
	839	// SCANNO
[1987]	840	// MS: 1-base
	841	// Scantable: 0-base
	842	Int dest5 = current5 + added5 ;
	843	scanNum -= 1 ;
	844	for ( Int irow = current5 ; irow < dest5 ; irow++ )
	845	scannoCol->putScalar( irow, (uInt)scanNum ) ;
[1978]	846
[1987]	847	//os_ << "added5 = " << added5 << endl ;
[1978]	848	added4 += added5 ;
[1974]	849	iter4.next() ;
	850	}
[1978]	851
	852	// IFNO
[1987]	853	Int dest4 = current4 + added4 ;
	854	for ( Int irow = current4 ; irow < dest4 ; irow++ )
	855	ifnoCol->putScalar( irow, (uInt)spwId ) ;
[1978]	856
	857	// FREQ_ID
[1987]	858	uInt fid = ifmap[spwId] ;
	859	for ( Int irow = current4 ; irow < dest4 ; irow++ )
	860	freqidCol->putScalar( irow, fid ) ;
[1978]	861
	862	// MOLECULE_ID
[1987]	863	for ( Int irow = current4 ; irow < dest4 ; irow++ )
	864	molidCol->putScalar( irow, molId ) ;
[1978]	865
	866	// SRCNAME
[1987]	867	for ( Int irow = current4 ; irow < dest4 ; irow++ )
	868	srcnameCol->putScalar( irow, srcName ) ;
[1978]	869
	870	// SRCVELOCITY, SRCPROPERMOTION and SRCDIRECTION
	871	// no reference conversion for direction at the moment (assume J2000)
	872	// no reference conversion for velocity at the moment (assume LSRK)
[1987]	873	for ( Int irow = current4 ; irow < dest4 ; irow++ )
	874	srcpmCol->put( irow, srcPM ) ;
	875	for ( Int irow = current4 ; irow < dest4 ; irow++ )
	876	srcdirCol->put( irow, srcDir ) ;
	877	for ( Int irow = current4 ; irow < dest4 ; irow++ )
	878	srcvelCol->putScalar( irow, sysVel ) ;
[1978]	879
[1987]	880	//os_ << "added4 = " << added4 << endl ;
[1978]	881	added3 += added4 ;
[1974]	882	iter3.next() ;
	883	}
[1978]	884
	885	// FIELDNAME
[1987]	886	Int dest3 = current3 + added3 ;
	887	for ( Int irow = current3 ; irow < dest3 ; irow++ )
	888	fieldnameCol->putScalar( irow, fieldName ) ;
[1978]	889
[1987]	890	//os_ << "added3 = " << added3 << endl ;
[1978]	891	added2 += added3 ;
[1974]	892	iter2.next() ;
	893	}
[1978]	894
	895	// BEAMNO
[1987]	896	Int dest2 = current2 + added2 ;
	897	for ( Int irow = current2 ; irow < dest2 ; irow++ )
	898	beamnoCol->putScalar( irow, (uInt)feedId ) ;
[1978]	899
	900	// FOCUS_ID
	901	// tentative
[1987]	902	uInt focusId = 0 ;
	903	for ( Int irow = current2 ; irow < dest2 ; irow++ )
	904	focusidCol->putScalar( irow, focusId ) ;
	905
	906	//os_ << "added2 = " << added2 << endl ;
[1978]	907	added1 += added2 ;
[1974]	908	iter1.next() ;
	909	}
	910	if ( sdh.nbeam < nbeam ) sdh.nbeam = nbeam ;
[1978]	911
[1987]	912	//os_ << "added1 = " << added1 << endl ;
[1978]	913	added0 += added1 ;
[1974]	914	iter0.next() ;
	915	}
	916
[1987]	917	//os_ << "added0 = " << added0 << endl ;
	918
[1978]	919	// REFBEAMNO
	920	// set 0 at the moment
[1987]	921	tcolw = tpoolw->construct( table_->table(), "REFBEAMNO" ) ;
	922	for ( Int irow = current0 ; irow < added0 ; irow++ )
	923	tcolw->putScalar( irow, 0 ) ;
	924	tpoolw->destroy( tcolw ) ;
[1978]	925
	926	// FIT_ID
	927	// nothing to do
[1987]	928	tcolw = tpoolw->construct( table_->table(), "FIT_ID" ) ;
	929	for ( Int irow = current0 ; irow < added0 ; irow++ )
	930	tcolw->putScalar( irow, -1 ) ;
	931	tpoolw->destroy( tcolw ) ;
[1978]	932
[1987]	933	// OPACITY
	934	// not used?
	935	tcolw = tpoolw->construct( table_->table(), "OPACITY" ) ;
	936	for ( Int irow = current0 ; irow < added0 ; irow++ )
	937	tcolw->putScalar( irow, 0.0 ) ;
	938	tpoolw->destroy( tcolw ) ;
[1978]	939
[1987]	940	// delete Scantable columns
	941	tpoolw->destroy( scannoCol ) ;
	942	tpoolw->destroy( cyclenoCol ) ;
	943	tpoolw->destroy( beamnoCol ) ;
	944	tpoolw->destroy( ifnoCol ) ;
	945	tpoolw->destroy( polnoCol ) ;
	946	tpoolw->destroy( freqidCol ) ;
	947	tpoolw->destroy( molidCol ) ;
	948	tpoolw->destroy( flagrowCol ) ;
	949	delete timeCol ;
	950	tpoolw->destroy( intervalCol ) ;
	951	tpoolw->destroy( srcnameCol ) ;
	952	tpoolw->destroy( srctypeCol ) ;
	953	tpoolw->destroy( fieldnameCol ) ;
	954	delete spCol ;
	955	delete flCol ;
	956	delete tsysCol ;
	957	delete dirCol ;
	958	tpoolw->destroy( azCol ) ;
	959	tpoolw->destroy( elCol ) ;
	960	tpoolw->destroy( tcalidCol ) ;
	961	tpoolw->destroy( focusidCol ) ;
	962	tpoolw->destroy( weatheridCol ) ;
	963	delete srcpmCol ;
	964	delete srcdirCol ;
	965	tpoolw->destroy( srcvelCol ) ;
	966	delete scanrateCol ;
	967
	968	delete tpoolr ;
	969	delete tpoolw ;
	970
	971
[1978]	972	// Table Keywords
[1974]	973	sdh.nif = ifmap.size() ;
[1987]	974	if ( ( telescopeName == "" ) \|\| ( antennaName == telescopeName ) ) {
[1974]	975	sdh.antennaname = antennaName ;
	976	}
	977	else {
	978	sdh.antennaname = telescopeName + "//" + antennaName ;
	979	}
	980	if ( stationName != "" ) {
	981	sdh.antennaname += "@" + stationName ;
	982	}
[1987]	983	ROArrayColumn<Double> pdirCol( pointtab, "DIRECTION" ) ;
	984	String dirref = pdirCol.keywordSet().asRecord("MEASINFO").asString("Ref") ;
[1974]	985	if ( dirref == "AZELGEO" \|\| dirref == "AZEL" ) {
	986	dirref = "J2000" ;
	987	}
	988	sscanf( dirref.chars()+1, "%f", &sdh.equinox ) ;
	989	sdh.epoch = "UTC" ;
	990	if (sdh.freqref == "TOPO") {
	991	sdh.freqref = "TOPOCENT";
	992	} else if (sdh.freqref == "GEO") {
	993	sdh.freqref = "GEOCENTR";
	994	} else if (sdh.freqref == "BARY") {
	995	sdh.freqref = "BARYCENT";
	996	} else if (sdh.freqref == "GALACTO") {
	997	sdh.freqref = "GALACTOC";
	998	} else if (sdh.freqref == "LGROUP") {
	999	sdh.freqref = "LOCALGRP";
	1000	} else if (sdh.freqref == "CMB") {
	1001	sdh.freqref = "CMBDIPOL";
	1002	} else if (sdh.freqref == "REST") {
	1003	sdh.freqref = "SOURCE";
	1004	}
	1005	table_->setHeader( sdh ) ;
[1987]	1006
	1007	// save path to POINTING table
	1008	//Path datapath(mstable_.tableName()) ;
	1009	Path datapath( tablename_ ) ;
	1010	String pTabName = datapath.absoluteName() + "/POINTING" ;
	1011	table_->table().rwKeywordSet().define( "POINTING", pTabName ) ;
	1012
	1013	// for GBT
	1014	if ( antennaName == "GBT" ) {
	1015	String goTabName = datapath.absoluteName() + "/GBT_GO" ;
	1016	table_->table().rwKeywordSet().define( "GBT_GO", goTabName ) ;
	1017	}
	1018	//os_ << "end MSFiller::fill()" << LogIO::POST ;
[1974]	1019	}
	1020
	1021	void MSFiller::close()
	1022	{
[1987]	1023	//tablesel_.closeSubTables() ;
[1974]	1024	mstable_.closeSubTables() ;
[1987]	1025	//tablesel_.unlock() ;
[1974]	1026	mstable_.unlock() ;
	1027	}
	1028
[1987]	1029	Int MSFiller::getSrcType( Int stateId, boost::object_pool<ROTableColumn> *tpool )
[1974]	1030	{
[1987]	1031	//os_ << "start MSFiller::getSrcType()" << LogIO::POST ;
[1974]	1032
	1033	MSState statetab = mstable_.state() ;
[1987]	1034	ROTableColumn *sharedCol ;
	1035	sharedCol = tpool->construct( statetab, "OBS_MODE" ) ;
	1036	String obsMode = sharedCol->asString( stateId ) ;
	1037	tpool->destroy( sharedCol ) ;
	1038	sharedCol = tpool->construct( statetab, "SIG" ) ;
	1039	Bool sig = sharedCol->asBool( stateId ) ;
	1040	tpool->destroy( sharedCol ) ;
	1041	sharedCol = tpool->construct( statetab, "REF" ) ;
	1042	Bool ref = sharedCol->asBool( stateId ) ;
	1043	tpool->destroy( sharedCol ) ;
	1044	sharedCol = tpool->construct( statetab, "CAL" ) ;
	1045	Double cal = (Double)(sharedCol->asdouble( stateId )) ;
	1046	tpool->destroy( sharedCol ) ;
[1974]	1047	//os_ << "OBS_MODE = " << obsMode << LogIO::POST ;
	1048
	1049	// determine separator
	1050	String sep = "" ;
	1051	if ( obsMode.find( ":" ) != String::npos ) {
	1052	sep = ":" ;
	1053	}
	1054	else if ( obsMode.find( "." ) != String::npos ) {
	1055	sep = "." ;
	1056	}
	1057
	1058	// determine SRCTYPE
	1059	Int srcType = SrcType::NOTYPE ;
	1060	if ( sep == ":" ) {
	1061	// sep == ":"
	1062	//
	1063	// GBT case
	1064	//
	1065	// obsMode1=Nod
	1066	// NOD
	1067	// obsMode1=OffOn
	1068	// obsMode2=PSWITCHON: PSON
	1069	// obsMode2=PSWITCHOFF: PSOFF
	1070	// obsMode1=??
	1071	// obsMode2=FSWITCH:
	1072	// SIG=1: FSON
	1073	// REF=1: FSOFF
[1987]	1074	// Calibration scan if CAL != 0
[1974]	1075	Int epos = obsMode.find_first_of( sep ) ;
	1076	Int nextpos = obsMode.find_first_of( sep, epos+1 ) ;
	1077	String obsMode1 = obsMode.substr( 0, epos ) ;
	1078	String obsMode2 = obsMode.substr( epos+1, nextpos-epos-1 ) ;
	1079	if ( obsMode1 == "Nod" ) {
	1080	srcType = SrcType::NOD ;
	1081	}
	1082	else if ( obsMode1 == "OffOn" ) {
	1083	if ( obsMode2 == "PSWITCHON" ) srcType = SrcType::PSON ;
	1084	if ( obsMode2 == "PSWITCHOFF" ) srcType = SrcType::PSOFF ;
	1085	}
	1086	else {
	1087	if ( obsMode2 == "FSWITCH" ) {
	1088	if ( sig ) srcType = SrcType::FSON ;
	1089	if ( ref ) srcType = SrcType::FSOFF ;
	1090	}
	1091	}
[1987]	1092	if ( cal > 0.0 ) {
	1093	if ( srcType == SrcType::NOD )
	1094	srcType = SrcType::NODCAL ;
	1095	else if ( srcType == SrcType::PSON )
	1096	srcType = SrcType::PONCAL ;
	1097	else if ( srcType == SrcType::PSOFF )
	1098	srcType = SrcType::POFFCAL ;
	1099	else if ( srcType == SrcType::FSON )
	1100	srcType = SrcType::FONCAL ;
	1101	else if ( srcType == SrcType::FSOFF )
	1102	srcType = SrcType::FOFFCAL ;
	1103	else
	1104	srcType = SrcType::CAL ;
	1105	}
[1974]	1106	}
	1107	else if ( sep == "." ) {
	1108	// sep == "."
	1109	//
	1110	// ALMA & EVLA case (MS via ASDM)
	1111	//
	1112	// obsMode1=CALIBRATE_*
	1113	// obsMode2=ON_SOURCE: PONCAL
	1114	// obsMode2=OFF_SOURCE: POFFCAL
	1115	// obsMode1=OBSERVE_TARGET
	1116	// obsMode2=ON_SOURCE: PON
	1117	// obsMode2=OFF_SOURCE: POFF
	1118	Int epos = obsMode.find_first_of( sep ) ;
	1119	Int nextpos = obsMode.find_first_of( sep, epos+1 ) ;
	1120	String obsMode1 = obsMode.substr( 0, epos ) ;
	1121	String obsMode2 = obsMode.substr( epos+1, nextpos-epos-1 ) ;
	1122	if ( obsMode1.find( "CALIBRATE_" ) == 0 ) {
	1123	if ( obsMode2 == "ON_SOURCE" ) srcType = SrcType::PONCAL ;
	1124	if ( obsMode2 == "OFF_SOURCE" ) srcType = SrcType::POFFCAL ;
	1125	}
	1126	else if ( obsMode1 == "OBSERVE_TARGET" ) {
	1127	if ( obsMode2 == "ON_SOURCE" ) srcType = SrcType::PSON ;
	1128	if ( obsMode2 == "OFF_SOURCE" ) srcType = SrcType::PSOFF ;
	1129	}
	1130	}
	1131	else {
	1132	if ( sig ) srcType = SrcType::SIG ;
	1133	if ( ref ) srcType = SrcType::REF ;
	1134	}
	1135
	1136	//os_ << "srcType = " << srcType << LogIO::POST ;
[1987]	1137	//os_ << "end MSFiller::getSrcType()" << LogIO::POST ;
[1974]	1138	return srcType ;
	1139	}
	1140
	1141	Vector<uInt> MSFiller::getPolNo( Int corrType )
	1142	{
[1987]	1143	//os_ << "start MSFiller::getPolNo()" << LogIO::POST ;
[1974]	1144	Vector<uInt> polno( 1 ) ;
	1145
	1146	if ( corrType == Stokes::I \|\| corrType == Stokes::RR \|\| corrType == Stokes::XX ) {
	1147	polno = 0 ;
	1148	}
	1149	else if ( corrType == Stokes::Q \|\| corrType == Stokes::LL \|\| corrType == Stokes::YY ) {
	1150	polno = 1 ;
	1151	}
	1152	else if ( corrType == Stokes::U ) {
	1153	polno = 2 ;
	1154	}
	1155	else if ( corrType == Stokes::V ) {
	1156	polno = 3 ;
	1157	}
	1158	else if ( corrType == Stokes::RL \|\| corrType == Stokes::XY \|\| corrType == Stokes::LR \|\| corrType == Stokes::RL ) {
	1159	polno.resize( 2 ) ;
	1160	polno[0] = 2 ;
	1161	polno[1] = 3 ;
	1162	}
	1163	else if ( corrType == Stokes::Plinear ) {
	1164	polno[0] = 1 ;
	1165	}
	1166	else if ( corrType == Stokes::Pangle ) {
	1167	polno[0] = 2 ;
	1168	}
	1169	else {
	1170	polno = 99 ;
	1171	}
	1172	//os_ << "polno = " << polno << LogIO::POST ;
[1987]	1173	//os_ << "end MSFiller::getPolNo()" << LogIO::POST ;
[1974]	1174
	1175	return polno ;
	1176	}
	1177
	1178	String MSFiller::getPolType( Int corrType )
	1179	{
[1987]	1180	//os_ << "start MSFiller::getPolType()" << LogIO::POST ;
[1974]	1181	String poltype = "" ;
	1182
	1183	if ( corrType == Stokes::I \|\| corrType == Stokes::Q \|\| corrType == Stokes::U \|\| corrType == Stokes::V )
	1184	poltype = "stokes" ;
	1185	else if ( corrType == Stokes::XX \|\| corrType == Stokes::YY \|\| corrType == Stokes::XY \|\| corrType == Stokes::YX )
	1186	poltype = "linear" ;
	1187	else if ( corrType == Stokes::RR \|\| corrType == Stokes::LL \|\| corrType == Stokes::RL \|\| corrType == Stokes::LR )
	1188	poltype = "circular" ;
	1189	else if ( corrType == Stokes::Plinear \|\| corrType == Stokes::Pangle )
	1190	poltype = "linpol" ;
	1191
[1987]	1192	//os_ << "end MSFiller::getPolType()" << LogIO::POST ;
[1974]	1193	return poltype ;
	1194	}
	1195
	1196	void MSFiller::fillWeather()
	1197	{
[1987]	1198	//os_ << "start MSFiller::fillWeather()" << LogIO::POST ;
	1199	Table mWeather = mstable_.weather() ;
	1200	//Table mWeatherSel = mWeather( mWeather.col("ANTENNA_ID") == antenna_ ).sort("TIME") ;
	1201	Table mWeatherSel( mWeather( mWeather.col("ANTENNA_ID") == antenna_ ).sort("TIME") ) ;
[1974]	1202	//os_ << "mWeatherSel.nrow() = " << mWeatherSel.nrow() << LogIO::POST ;
	1203	if ( mWeatherSel.nrow() == 0 ) {
	1204	os_ << "No rows with ANTENNA_ID = " << antenna_ << ", Try -1..." << LogIO::POST ;
[1987]	1205	mWeatherSel = Table( MSWeather( mWeather( mWeather.col("ANTENNA_ID") == -1 ) ) ) ;
[1974]	1206	if ( mWeatherSel.nrow() == 0 ) {
	1207	os_ << "No rows in WEATHER table" << LogIO::POST ;
	1208	}
	1209	}
[1987]	1210	uInt wnrow = mWeatherSel.nrow() ;
[1974]	1211	//os_ << "wnrow = " << wnrow << LogIO::POST ;
	1212
	1213	if ( wnrow == 0 )
	1214	return ;
	1215
	1216	Table wtab = table_->weather().table() ;
	1217	wtab.addRow( wnrow ) ;
	1218
[1987]	1219	ScalarColumn<Float> *fCol ;
	1220	ROScalarColumn<Float> *sharedFloatCol ;
	1221	if ( mWeatherSel.tableDesc().isColumn( "TEMPERATURE" ) ) {
	1222	fCol = new ScalarColumn<Float>( wtab, "TEMPERATURE" ) ;
	1223	sharedFloatCol = new ROScalarColumn<Float>( mWeatherSel, "TEMPERATURE" ) ;
	1224	fCol->putColumn( *sharedFloatCol ) ;
	1225	delete sharedFloatCol ;
	1226	delete fCol ;
	1227	}
	1228	if ( mWeatherSel.tableDesc().isColumn( "PRESSURE" ) ) {
	1229	fCol = new ScalarColumn<Float>( wtab, "PRESSURE" ) ;
	1230	sharedFloatCol = new ROScalarColumn<Float>( mWeatherSel, "PRESSURE" ) ;
	1231	fCol->putColumn( *sharedFloatCol ) ;
	1232	delete sharedFloatCol ;
	1233	delete fCol ;
	1234	}
	1235	if ( mWeatherSel.tableDesc().isColumn( "REL_HUMIDITY" ) ) {
	1236	fCol = new ScalarColumn<Float>( wtab, "HUMIDITY" ) ;
	1237	sharedFloatCol = new ROScalarColumn<Float>( mWeatherSel, "REL_HUMIDITY" ) ;
	1238	fCol->putColumn( *sharedFloatCol ) ;
	1239	delete sharedFloatCol ;
	1240	delete fCol ;
	1241	}
	1242	if ( mWeatherSel.tableDesc().isColumn( "WIND_SPEED" ) ) {
	1243	fCol = new ScalarColumn<Float>( wtab, "WINDSPEED" ) ;
	1244	sharedFloatCol = new ROScalarColumn<Float>( mWeatherSel, "WIND_SPEED" ) ;
	1245	fCol->putColumn( *sharedFloatCol ) ;
	1246	delete sharedFloatCol ;
	1247	delete fCol ;
	1248	}
	1249	if ( mWeatherSel.tableDesc().isColumn( "WIND_DIRECTION" ) ) {
	1250	fCol = new ScalarColumn<Float>( wtab, "WINDAZ" ) ;
	1251	sharedFloatCol = new ROScalarColumn<Float>( mWeatherSel, "WIND_DIRECTION" ) ;
	1252	fCol->putColumn( *sharedFloatCol ) ;
	1253	delete sharedFloatCol ;
	1254	delete fCol ;
	1255	}
[1974]	1256	ScalarColumn<uInt> idCol( wtab, "ID" ) ;
[1987]	1257	for ( uInt irow = 0 ; irow < wnrow ; irow++ )
	1258	idCol.put( irow, irow ) ;
[1974]	1259
[1987]	1260	ROScalarQuantColumn<Double> tqCol( mWeatherSel, "TIME" ) ;
	1261	ROScalarColumn<Double> tCol( mWeatherSel, "TIME" ) ;
	1262	String tUnit = tqCol.getUnits() ;
	1263	mwTime_ = tCol.getColumn() ;
[1974]	1264	if ( tUnit == "d" )
	1265	mwTime_ *= 86400.0 ;
[1987]	1266	tqCol.attach( mWeatherSel, "INTERVAL" ) ;
	1267	tCol.attach( mWeatherSel, "INTERVAL" ) ;
	1268	String iUnit = tqCol.getUnits() ;
	1269	mwInterval_ = tCol.getColumn() ;
[1974]	1270	if ( iUnit == "d" )
	1271	mwInterval_ *= 86400.0 ;
	1272	//os_ << "mwTime[0] = " << mwTime_[0] << " mwInterval[0] = " << mwInterval_[0] << LogIO::POST ;
[1987]	1273	//os_ << "end MSFiller::fillWeather()" << LogIO::POST ;
[1974]	1274	}
	1275
	1276	void MSFiller::fillFocus()
	1277	{
[1987]	1278	//os_ << "start MSFiller::fillFocus()" << LogIO::POST ;
[1974]	1279	// tentative
	1280	Table tab = table_->focus().table() ;
	1281	tab.addRow( 1 ) ;
	1282	ScalarColumn<uInt> idCol( tab, "ID" ) ;
	1283	idCol.put( 0, 0 ) ;
[1987]	1284	//os_ << "end MSFiller::fillFocus()" << LogIO::POST ;
[1974]	1285	}
	1286
[1987]	1287	void MSFiller::fillTcal( boost::object_pool<ROTableColumn> *tpool )
[1974]	1288	{
[1987]	1289	//os_ << "start MSFiller::fillTcal()" << LogIO::POST ;
	1290
	1291	//MSSysCal sctab = mstable_.sysCal() ;
	1292	Table sctab = mstable_.sysCal() ;
[1974]	1293	if ( sctab.nrow() == 0 ) {
	1294	os_ << "No SysCal rows" << LogIO::POST ;
	1295	return ;
	1296	}
	1297	Bool isSp = sctab.tableDesc().isColumn( "TCAL_SPECTRUM" ) ;
[1987]	1298	//MSSysCal sctabsel = sctab( sctab.col("ANTENNA_ID") == antenna_ ) ;
	1299	//Table sctabsel = sctab( sctab.col("ANTENNA_ID") == antenna_ ) ;
	1300	Table sctabsel( sctab( sctab.col("ANTENNA_ID") == antenna_ ) ) ;
[1974]	1301	if ( sctabsel.nrow() == 0 ) {
	1302	os_ << "No SysCal rows" << LogIO::POST ;
	1303	return ;
	1304	}
[1987]	1305	ROArrayColumn<Float> *tmpTcalCol = new ROArrayColumn<Float>( sctabsel, "TCAL" ) ;
	1306	uInt npol = tmpTcalCol->shape( 0 )(0) ;
	1307	delete tmpTcalCol ;
[1974]	1308	//os_ << "fillTcal(): npol = " << npol << LogIO::POST ;
	1309	Table tab = table_->tcal().table() ;
	1310	ScalarColumn<uInt> idCol( tab, "ID" ) ;
	1311	ScalarColumn<String> timeCol( tab, "TIME" ) ;
	1312	ArrayColumn<Float> tcalCol( tab, "TCAL" ) ;
[1987]	1313	ROTableColumn *sharedCol ;
[1974]	1314	uInt oldnr = 0 ;
	1315	uInt newnr = 0 ;
	1316	TableIterator iter0( sctabsel, "FEED_ID" ) ;
	1317	while( !iter0.pastEnd() ) {
[1987]	1318	Table t0 = iter0.table() ;
	1319	sharedCol = tpool->construct( t0, "FEED_ID" ) ;
	1320	Int feedId = sharedCol->asInt( 0 ) ;
	1321	tpool->destroy( sharedCol ) ;
[1974]	1322	String ffield = "FEED" + String::toString( feedId ) ;
	1323	Record rec ;
	1324	TableIterator iter1( t0, "SPECTRAL_WINDOW_ID" ) ;
	1325	while( !iter1.pastEnd() ) {
	1326	MSSysCal t1( iter1.table().sort("TIME") ) ;
	1327	uInt nrow = t1.nrow() ;
[1987]	1328	sharedCol = tpool->construct( t1, "SPECTRAL_WINDOW_ID" ) ;
	1329	Int spwId = sharedCol->asInt( 0 ) ;
	1330	tpool->destroy( sharedCol ) ;
[1974]	1331	String spwfield = "SPW" + String::toString( spwId ) ;
[1987]	1332	ROScalarQuantColumn<Double> scTimeCol( t1, "TIME" ) ;
[1974]	1333	ROArrayColumn<Float> scTcalCol ;
	1334	IPosition newShape( 2, 1, nrow ) ;
	1335	if ( isSp ) {
[1987]	1336	scTcalCol.attach( t1, "TCAL_SPECTRUM" ) ;
[1974]	1337	newShape[0] = scTcalCol.shape(0)(1) ;
	1338	}
	1339	else {
[1987]	1340	scTcalCol.attach( t1, "TCAL" ) ;
[1974]	1341	}
[1987]	1342	tab.addRow( nrow*npol ) ;
	1343	newnr += nrow*npol ;
[1974]	1344	for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
[1987]	1345	String sTime = MVTime( scTimeCol(irow) ).string( MVTime::YMD ) ;
	1346	for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
	1347	timeCol.put( oldnr+irow+ipol*nrow, sTime ) ;
	1348	}
[1974]	1349	}
[1987]	1350	uInt idx = oldnr ;
[1974]	1351	for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
[1987]	1352	for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
	1353	idCol.put( oldnr+ipol*nrow+irow, idx++ ) ;
[1974]	1354	}
[1987]	1355	}
	1356	Vector<uInt> idminmax( 2, oldnr ) ;
	1357	for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
	1358	Matrix<Float> subtcal = scTcalCol( irow ) ;
	1359	for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) {
	1360	tcalCol.put( oldnr+ipol*nrow+irow, subtcal.row( ipol ) ) ;
[1974]	1361	}
	1362	}
	1363	idminmax[1] = newnr - 1 ;
[1987]	1364	oldnr = newnr ;
[1974]	1365	rec.define( spwfield, idminmax ) ;
	1366	iter1++ ;
	1367	}
	1368	tcalrec_.defineRecord( ffield, rec ) ;
	1369	iter0++ ;
	1370	}
	1371
[1987]	1372
[1974]	1373	//tcalrec_.print( std::cout ) ;
[1987]	1374	//os_ << "end MSFiller::fillFocus()" << LogIO::POST ;
[1974]	1375	}
	1376
	1377	uInt MSFiller::getWeatherId( uInt idx, Double wtime )
	1378	{
[1987]	1379	//os_ << "start MSFiller::getWeatherId()" << LogIO::POST ;
[1974]	1380	uInt nrow = mwTime_.size() ;
	1381	if ( nrow == 0 )
	1382	return 0 ;
	1383	uInt wid = nrow ;
	1384	for ( uInt i = idx ; i < nrow-1 ; i++ ) {
	1385	Double tStart = mwTime_[i]-0.5*mwInterval_[i] ;
	1386	// use of INTERVAL column is problematic
	1387	// since there are "blank" time of weather monitoring
	1388	//Double tEnd = tStart + mwInterval_[i] ;
	1389	Double tEnd = mwTime_[i+1]-0.5*mwInterval_[i+1] ;
	1390	//os_ << "tStart = " << tStart << " dtEnd = " << tEnd-tStart << " dwtime = " << wtime-tStart << LogIO::POST ;
	1391	if ( wtime >= tStart && wtime <= tEnd ) {
	1392	wid = i ;
	1393	break ;
	1394	}
	1395	}
	1396	if ( wid == nrow ) {
	1397	uInt i = nrow - 1 ;
	1398	Double tStart = mwTime_[i-1]+0.5*mwInterval_[i-1] ;
	1399	Double tEnd = mwTime_[i]+0.5*mwInterval_[i] ;
	1400	//os_ << "tStart = " << tStart << " dtEnd = " << tEnd-tStart << " dwtime = " << wtime-tStart << LogIO::POST ;
	1401	if ( wtime >= tStart && wtime <= tEnd )
	1402	wid = i ;
	1403	}
	1404
	1405	//if ( wid == nrow )
	1406	//os_ << LogIO::WARN << "Couldn't find correct WEATHER_ID for time " << wtime << LogIO::POST ;
	1407
[1987]	1408	//os_ << "end MSFiller::getWeatherId()" << LogIO::POST ;
[1974]	1409	return wid ;
	1410	}
	1411
	1412	Vector<Double> MSFiller::getSysCalTime( MSSysCal &tab, MEpoch::ROScalarColumn &tcol )
	1413	{
[1987]	1414	//os_ << "start MSFiller::getSysCalTime()" << LogIO::POST ;
[1974]	1415	uInt nrow = tcol.table().nrow() ;
	1416	Vector<Double> tstr( nrow, -1.0 ) ;
	1417	if ( tab.nrow() == 0 )
	1418	return tstr ;
	1419	uInt scnrow = tab.nrow() ;
[1987]	1420	ROScalarMeasColumn<MEpoch> scTimeCol( tab, "TIME" ) ;
	1421	ROScalarQuantColumn<Double> scIntervalCol( tab, "INTERVAL" ) ;
[1974]	1422	uInt idx = 0 ;
	1423	const Double half = 0.5e0 ;
	1424	for ( uInt i = 0 ; i < nrow ; i++ ) {
	1425	Double t = tcol( i ).get( "s" ).getValue() ;
	1426	for ( uInt j = idx ; j < scnrow-1 ; j++ ) {
	1427	Double tsc1 = scTimeCol( j ).get( "s" ).getValue() ;
	1428	Double dt1 = scIntervalCol( j ).getValue("s") ;
	1429	Double tsc2 = scTimeCol( j+1 ).get( "s" ).getValue() ;
	1430	Double dt2 = scIntervalCol( j+1 ).getValue("s") ;
	1431	if ( t > tsc1-halfdt1 && t <= tsc2-halfdt2 ) {
	1432	tstr[i] = tsc1 ;
	1433	idx = j ;
	1434	break ;
	1435	}
	1436	}
	1437	if ( tstr[i] == -1.0 ) {
	1438	Double tsc = scTimeCol( scnrow-1 ).get( "s" ).getValue() ;
	1439	Double dt = scIntervalCol( scnrow-1 ).getValue( "s" ) ;
	1440	if ( t <= tsc+0.5*dt )
	1441	tstr[i] = tsc ;
	1442	}
	1443	}
[1987]	1444	//os_ << "end MSFiller::getSysCalTime()" << LogIO::POST ;
[1974]	1445	return tstr ;
	1446	}
	1447
[1987]	1448	uInt MSFiller::getTsys( uInt idx, Matrix<Float> &tsys, MSSysCal &tab, Double t )
[1974]	1449	{
[1987]	1450	//os_ << "start MSFiller::getTsys()" << LogIO::POST ;
[1974]	1451	uInt nrow = tab.nrow() ;
	1452	if ( nrow == 0 ) {
	1453	os_ << "No SysCal rows" << LogIO::POST ;
	1454	tsys.resize( IPosition(0) ) ;
	1455	return 0 ;
	1456	}
	1457	Bool isSp = tab.tableDesc().isColumn( "TSYS_SPECTRUM" ) ;
[1987]	1458	ROScalarMeasColumn<MEpoch> scTimeCol( tab, "TIME" ) ;
[1974]	1459	ROArrayColumn<Float> mTsysCol ;
	1460	if ( isSp ) {
[1987]	1461	mTsysCol.attach( tab, "TSYS_SPECTRUM" ) ;
[1974]	1462	}
	1463	else {
[1987]	1464	mTsysCol.attach( tab, "TSYS" ) ;
[1974]	1465	}
	1466	for ( uInt i = idx ; i < nrow ; i++ ) {
	1467	Double tref = scTimeCol( i ).get( "s" ).getValue() ;
	1468	if ( t == tref ) {
	1469	tsys.reference( mTsysCol( i ) ) ;
	1470	idx = i ;
	1471	break ;
	1472	}
	1473	}
[1987]	1474	//os_ << "MSFiller::getTsys() idx = " << idx << " tsys = " << tsys << LogIO::POST ;
	1475	//os_ << "end MSFiller::getTsys()" << LogIO::POST ;
[1974]	1476	return idx ;
	1477	}
	1478
	1479	Vector<uInt> MSFiller::getTcalId( Int fid, Int spwid, Double t )
	1480	{
[1987]	1481	//os_ << "start MSFiller::getTcalId()" << LogIO::POST ;
[1974]	1482	String feed = "FEED" + String::toString(fid) ;
	1483	String spw = "SPW" + String::toString(spwid) ;
	1484	String sctime = MVTime( Quantum<Double>(t,"s") ).string(MVTime::YMD) ;
	1485	Table ttab = table_->tcal().table() ;
	1486	if ( ttab.nrow() == 0 ) {
	1487	os_ << "No TCAL rows" << LogIO::POST ;
	1488	Vector<uInt> tcalids( 0 ) ;
	1489	return tcalids ;
	1490	}
	1491	Vector<uInt> ids = tcalrec_.asRecord(feed).asArrayuInt(spw) ;
[1987]	1492	//Table ttabsel = ttab( ttab.col("TIME") == sctime && ttab.col("ID") >= ids[0] && ttab.col("ID") <= ids[1] ).sort("ID") ;
	1493	Table ttabsel( ttab( ttab.col("TIME") == sctime && ttab.col("ID") >= ids[0] && ttab.col("ID") <= ids[1] ).sort("ID")) ;
[1974]	1494	uInt nrow = ttabsel.nrow() ;
	1495	Vector<uInt> tcalids( nrow ) ;
	1496	if ( nrow == 0 ) {
	1497	os_ << "No TCAL rows" << LogIO::POST ;
	1498	return tcalids ;
	1499	}
	1500	ROScalarColumn<uInt> idCol( ttabsel, "ID" ) ;
	1501	tcalids[0] = idCol(0) ;
	1502	if ( nrow == 2 ) {
	1503	tcalids[1] = idCol(1) ;
	1504	}
	1505	else if ( nrow == 3 ) {
	1506	tcalids[1] = idCol(2) ;
	1507	tcalids[2] = idCol(1) ;
	1508	}
	1509	else if ( nrow == 4 ) {
	1510	tcalids[1] = idCol(3) ;
	1511	tcalids[2] = idCol(1) ;
	1512	tcalids[3] = idCol(2) ;
	1513	}
	1514
[1987]	1515	//os_ << "end MSFiller::getTcalId()" << LogIO::POST ;
[1974]	1516	return tcalids ;
	1517	}
	1518
[1987]	1519	uInt MSFiller::getDirection( uInt idx, Vector<Double> &dir, Vector<Double> &srate, String &ref, MSPointing &tab, Double t )
[1974]	1520	{
[1987]	1521	//os_ << "start MSFiller::getDirection()" << LogIO::POST ;
[1974]	1522	// assume that cols is sorted by TIME
	1523	Bool doInterp = False ;
[1987]	1524	//uInt nrow = cols.nrow() ;
	1525	uInt nrow = tab.nrow() ;
[1974]	1526	if ( nrow == 0 )
	1527	return 0 ;
[1987]	1528	ROScalarMeasColumn<MEpoch> tcol( tab, "TIME" ) ;
	1529	ROArrayMeasColumn<MDirection> dmcol( tab, "DIRECTION" ) ;
	1530	ROArrayColumn<Double> dcol( tab, "DIRECTION" ) ;
[1974]	1531	// ensure that tcol(idx) < t
	1532	//os_ << "tcol(idx) = " << tcol(idx).get("s").getValue() << " t = " << t << " diff = " << tcol(idx).get("s").getValue()-t << endl ;
	1533	while ( tcol(idx).get("s").getValue() > t && idx > 0 )
	1534	idx-- ;
	1535	//os_ << "idx = " << idx << LogIO::POST ;
	1536
	1537	// index search
	1538	for ( uInt i = idx ; i < nrow ; i++ ) {
	1539	Double tref = tcol( i ).get( "s" ).getValue() ;
	1540	if ( tref == t ) {
	1541	idx = i ;
	1542	break ;
	1543	}
	1544	else if ( tref > t ) {
	1545	if ( i == 0 ) {
	1546	idx = i ;
	1547	}
	1548	else {
	1549	idx = i-1 ;
	1550	doInterp = True ;
	1551	}
	1552	break ;
	1553	}
	1554	else {
	1555	idx = nrow - 1 ;
	1556	}
	1557	}
	1558	//os_ << "searched idx = " << idx << LogIO::POST ;
	1559
	1560	//os_ << "dmcol(idx).shape() = " << dmcol(idx).shape() << LogIO::POST ;
	1561	IPosition ip( dmcol(idx).shape().nelements(), 0 ) ;
	1562	//os_ << "ip = " << ip << LogIO::POST ;
	1563	ref = dmcol(idx)(ip).getRefString() ;
	1564	//os_ << "ref = " << ref << LogIO::POST ;
	1565	if ( doInterp ) {
	1566	//os_ << "do interpolation" << LogIO::POST ;
	1567	//os_ << "dcol(idx).shape() = " << dcol(idx).shape() << LogIO::POST ;
	1568	Double tref0 = tcol(idx).get("s").getValue() ;
	1569	Double tref1 = tcol(idx+1).get("s").getValue() ;
[1987]	1570	Matrix<Double> mdir0 = dcol( idx ) ;
	1571	Matrix<Double> mdir1 = dcol( idx+1 ) ;
	1572	Vector<Double> dir0 = mdir0.column( 0 ) ;
[1974]	1573	//os_ << "dir0 = " << dir0 << LogIO::POST ;
[1987]	1574	Vector<Double> dir1 = mdir1.column( 0 ) ;
[1974]	1575	//os_ << "dir1 = " << dir1 << LogIO::POST ;
	1576	Double dt0 = t - tref0 ;
	1577	Double dt1 = tref1 - t ;
	1578	dir.reference( (dt0dir1+dt1dir0)/(dt0+dt1) ) ;
[1987]	1579	if ( mdir0.ncolumn() > 1 ) {
	1580	if ( dt0 >= dt1 )
	1581	srate.reference( mdir0.column( 1 ) ) ;
	1582	else
	1583	srate.reference( mdir1.column( 1 ) ) ;
[1974]	1584	}
	1585	//os_ << "dir = " << dir << LogIO::POST ;
	1586	}
	1587	else {
	1588	//os_ << "no interpolation" << LogIO::POST ;
[1987]	1589	Matrix<Double> mdir0 = dcol( idx ) ;
	1590	dir.reference( mdir0.column( 0 ) ) ;
	1591	if ( mdir0.ncolumn() > 1 )
	1592	srate.reference( mdir0.column( 1 ) ) ;
[1974]	1593	}
	1594
[1987]	1595	//os_ << "end MSFiller::getDirection()" << LogIO::POST ;
[1974]	1596	return idx ;
	1597	}
	1598
	1599	} ;
	1600

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

Context Navigation

source: trunk/src/MSFiller.cpp@ 1987

Download in other formats: