Context Navigation

source: trunk/python/scantable.py@ 1101

Visit:

Last change on this file since 1101 was 1099, checked in by mar637, 18 years ago
moved average_channel into average_time as an addition weight
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 57.1 KB

Rev	Line
[876]	1	from asap._asap import Scantable
[226]	2	from asap import rcParams
[976]	3	from asap import print_log, asaplog
[946]	4	from asap import selector
[189]	5	from numarray import ones,zeros
[102]	6	import sys
	7
[876]	8	class scantable(Scantable):
[102]	9	"""
	10	The ASAP container for scans
	11	"""
[710]	12
[484]	13	def __init__(self, filename, average=None, unit=None):
[102]	14	"""
	15	Create a scantable from a saved one or make a reference
	16	Parameters:
[181]	17	filename: the name of an asap table on disk
	18	or
	19	the name of a rpfits/sdfits/ms file
	20	(integrations within scans are auto averaged
	21	and the whole file is read)
	22	or
	23	[advanced] a reference to an existing
[102]	24	scantable
[484]	25	average: average all integrations withinb a scan on read.
	26	The default (True) is taken from .asaprc.
	27	unit: brightness unit; must be consistent with K or Jy.
[340]	28	Over-rides the default selected by the reader
	29	(input rpfits/sdfits/ms) or replaces the value
	30	in existing scantables
[710]	31	"""
[976]	32	if average is None:
[710]	33	average = rcParams['scantable.autoaverage']
[484]	34	varlist = vars()
[876]	35	from asap._asap import stmath
	36	self._math = stmath()
	37	if isinstance(filename, Scantable):
	38	Scantable.__init__(self, filename)
[181]	39	else:
[976]	40	if isinstance(filename,str):
	41	import os.path
	42	filename = os.path.expandvars(filename)
	43	filename = os.path.expanduser(filename)
	44	if not os.path.exists(filename):
	45	s = "File '%s' not found." % (filename)
[718]	46	if rcParams['verbose']:
[976]	47	asaplog.push(s)
	48	print asaplog.pop().strip()
[718]	49	return
[976]	50	raise IOError(s)
	51	if os.path.isdir(filename):
	52	# crude check if asap table
	53	if os.path.exists(filename+'/table.info'):
[1073]	54	Scantable.__init__(self, filename, rcParams['scantable.storage']=='disk')
[976]	55	if unit is not None:
	56	self.set_fluxunit(unit)
	57	self.set_freqframe(rcParams['scantable.freqframe'])
[718]	58	else:
[976]	59	msg = "The given file '%s'is not a valid asap table." % (filename)
	60	if rcParams['verbose']:
	61	print msg
	62	return
	63	else:
	64	raise IOError(msg)
[226]	65	else:
[976]	66	self._fill([filename],unit, average)
	67	elif (isinstance(filename,list) or isinstance(filename,tuple)) \
	68	and isinstance(filename[-1], str):
	69	self._fill(filename, unit, average)
[714]	70	print_log()
[102]	71
[876]	72	def save(self, name=None, format=None, overwrite=False):
[116]	73	"""
[710]	74	Store the scantable on disk. This can be an asap (aips++) Table, SDFITS,
[280]	75	Image FITS or MS2 format.
[116]	76	Parameters:
[1093]	77	name: the name of the outputfile. For format "ASCII"
	78	this is the root file name (data in 'name'.txt
[497]	79	and header in 'name'_header.txt)
[116]	80	format: an optional file format. Default is ASAP.
[280]	81	Allowed are - 'ASAP' (save as ASAP [aips++] Table),
[194]	82	'SDFITS' (save as SDFITS file)
[200]	83	'ASCII' (saves as ascii text file)
[226]	84	'MS2' (saves as an aips++
	85	MeasurementSet V2)
[411]	86	overwrite: If the file should be overwritten if it exists.
[256]	87	The default False is to return with warning
[411]	88	without writing the output. USE WITH CARE.
[116]	89	Example:
	90	scan.save('myscan.asap')
	91	scan.save('myscan.sdfits','SDFITS')
	92	"""
[411]	93	from os import path
[226]	94	if format is None: format = rcParams['scantable.save']
[256]	95	suffix = '.'+format.lower()
	96	if name is None or name =="":
	97	name = 'scantable'+suffix
[718]	98	from asap import asaplog
	99	msg = "No filename given. Using default name %s..." % name
	100	asaplog.push(msg)
[411]	101	name = path.expandvars(name)
[256]	102	if path.isfile(name) or path.isdir(name):
	103	if not overwrite:
[718]	104	msg = "File %s exists." % name
	105	if rcParams['verbose']:
	106	print msg
	107	return
	108	else:
	109	raise IOError(msg)
[451]	110	format2 = format.upper()
	111	if format2 == 'ASAP':
[116]	112	self._save(name)
	113	else:
[989]	114	from asap._asap import stwriter as stw
	115	w = stw(format2)
	116	w.write(self, name)
[718]	117	print_log()
[116]	118	return
	119
[102]	120	def copy(self):
	121	"""
	122	Return a copy of this scantable.
	123	Parameters:
[113]	124	none
[102]	125	Example:
	126	copiedscan = scan.copy()
	127	"""
[876]	128	sd = scantable(Scantable._copy(self))
[113]	129	return sd
	130
[1093]	131	def drop_scan(self, scanid=None):
	132	"""
	133	Return a new scantable where the specified scan number(s) has(have)
	134	been dropped.
	135	Parameters:
	136	scanid: a (list of) scan number(s)
	137	"""
	138	from asap import _is_sequence_or_number as _is_valid
	139	from asap import _to_list
	140	from asap import unique
	141	if not _is_valid(scanid):
	142	if rcParams['verbose']:
	143	print "Please specify a scanno to drop from the scantable"
	144	return
	145	else:
	146	raise RuntimeError("No scan given")
	147	try:
	148	scanid = _to_list(scanid)
	149	allscans = unique([ self.getscan(i) for i in range(self.nrow())])
	150	for sid in scanid: allscans.remove(sid)
	151	if len(allscans) == 0: raise ValueError("Can't remove all scans")
	152	except ValueError:
	153	if rcParams['verbose']:
	154	print "Couldn't find any match."
	155	return
	156	else: raise
	157	try:
	158	bsel = self.get_selection()
	159	sel = selector()
	160	sel.set_scans(allscans)
	161	self.set_selection(bsel+sel)
	162	scopy = self._copy()
	163	self.set_selection(bsel)
	164	return scantable(scopy)
	165	except RuntimeError:
	166	if rcParams['verbose']: print "Couldn't find any match."
	167	else: raise
	168
	169
[102]	170	def get_scan(self, scanid=None):
	171	"""
	172	Return a specific scan (by scanno) or collection of scans (by
	173	source name) in a new scantable.
	174	Parameters:
[513]	175	scanid: a (list of) scanno or a source name, unix-style
	176	patterns are accepted for source name matching, e.g.
	177	'*_R' gets all 'ref scans
[102]	178	Example:
[513]	179	# get all scans containing the source '323p459'
	180	newscan = scan.get_scan('323p459')
	181	# get all 'off' scans
	182	refscans = scan.get_scan('*_R')
	183	# get a susbset of scans by scanno (as listed in scan.summary())
	184	newscan = scan.get_scan([0,2,7,10])
[102]	185	"""
	186	if scanid is None:
[718]	187	if rcParams['verbose']:
	188	print "Please specify a scan no or name to retrieve from the scantable"
	189	return
	190	else:
	191	raise RuntimeError("No scan given")
	192
[102]	193	try:
[946]	194	bsel = self.get_selection()
	195	sel = selector()
[102]	196	if type(scanid) is str:
[946]	197	sel.set_name(scanid)
	198	self.set_selection(bsel+sel)
[876]	199	scopy = self._copy()
[946]	200	self.set_selection(bsel)
[876]	201	return scantable(scopy)
[102]	202	elif type(scanid) is int:
[946]	203	sel.set_scans([scanid])
	204	self.set_selection(bsel+sel)
[876]	205	scopy = self._copy()
[946]	206	self.set_selection(bsel)
[876]	207	return scantable(scopy)
[381]	208	elif type(scanid) is list:
[946]	209	sel.set_scans(scanid)
	210	self.set_selection(sel)
[876]	211	scopy = self._copy()
[946]	212	self.set_selection(bsel)
[876]	213	return scantable(scopy)
[381]	214	else:
[718]	215	msg = "Illegal scanid type, use 'int' or 'list' if ints."
	216	if rcParams['verbose']:
	217	print msg
	218	else:
	219	raise TypeError(msg)
[102]	220	except RuntimeError:
[718]	221	if rcParams['verbose']: print "Couldn't find any match."
	222	else: raise
[102]	223
	224	def __str__(self):
[876]	225	return Scantable._summary(self,True)
[102]	226
[976]	227	def summary(self, filename=None):
[102]	228	"""
	229	Print a summary of the contents of this scantable.
	230	Parameters:
	231	filename: the name of a file to write the putput to
	232	Default - no file output
[381]	233	verbose: print extra info such as the frequency table
	234	The default (False) is taken from .asaprc
[102]	235	"""
[976]	236	info = Scantable._summary(self, True)
	237	#if verbose is None: verbose = rcParams['scantable.verbosesummary']
[102]	238	if filename is not None:
[256]	239	if filename is "":
	240	filename = 'scantable_summary.txt'
[415]	241	from os.path import expandvars, isdir
[411]	242	filename = expandvars(filename)
[415]	243	if not isdir(filename):
[413]	244	data = open(filename, 'w')
	245	data.write(info)
	246	data.close()
	247	else:
[718]	248	msg = "Illegal file name '%s'." % (filename)
	249	if rcParams['verbose']:
	250	print msg
	251	else:
	252	raise IOError(msg)
	253	if rcParams['verbose']:
[794]	254	try:
	255	from IPython.genutils import page as pager
	256	except ImportError:
	257	from pydoc import pager
	258	pager(info)
[718]	259	else:
	260	return info
[710]	261
[946]	262
	263	def get_selection(self):
	264	"""
[1005]	265	Get the selection object currently set on this scantable.
	266	Parameters:
	267	none
	268	Example:
	269	sel = scan.get_selection()
	270	sel.set_ifs(0) # select IF 0
	271	scan.set_selection(sel) # apply modified selection
[946]	272	"""
	273	return selector(self._getselection())
	274
[1005]	275	def set_selection(self, selection=selector()):
[946]	276	"""
[1005]	277	Select a subset of the data. All following operations on this scantable
	278	are only applied to thi selection.
	279	Parameters:
	280	selection: a selector object (default unset the selection)
	281	Examples:
	282	sel = selector() # create a selection object
	283	self.set_scans([0,3]) # select SCANNO 0 and 3
	284	scan.set_selection(sel) # set the selection
	285	scan.summary() # will only print summary of scanno 0 an 3
	286	scan.set_selection() # unset the selection
[946]	287	"""
	288	self._setselection(selection)
	289
[553]	290	def set_cursor(self, beam=0, IF=0, pol=0):
[1006]	291	print "DEPRECATED - use set_selection"
[102]	292
[256]	293	def get_cursor(self):
[1006]	294	print "DEPRECATED - use get_selection"
[102]	295
[876]	296	def stats(self, stat='stddev', mask=None):
[102]	297	"""
[135]	298	Determine the specified statistic of the current beam/if/pol
[102]	299	Takes a 'mask' as an optional parameter to specify which
	300	channels should be excluded.
	301	Parameters:
[135]	302	stat: 'min', 'max', 'sumsq', 'sum', 'mean'
	303	'var', 'stddev', 'avdev', 'rms', 'median'
	304	mask: an optional mask specifying where the statistic
[102]	305	should be determined.
	306	Example:
[113]	307	scan.set_unit('channel')
[102]	308	msk = scan.create_mask([100,200],[500,600])
[135]	309	scan.stats(stat='mean', mask=m)
[102]	310	"""
[256]	311	from numarray import array,zeros,Float
[102]	312	if mask == None:
[876]	313	mask = []
[256]	314	axes = ['Beam','IF','Pol','Time']
[876]	315	if not self._check_ifs():
	316	raise ValueError("Cannot apply mask as the IFs have different number of channels"
	317	"Please use setselection() to select individual IFs")
[256]	318
[876]	319	statvals = self._math._stats(self, mask, stat)
	320	out = ''
	321	axes = []
	322	for i in range(self.nrow()):
	323	axis = []
	324	axis.append(self.getscan(i))
	325	axis.append(self.getbeam(i))
	326	axis.append(self.getif(i))
	327	axis.append(self.getpol(i))
	328	axis.append(self.getcycle(i))
	329	axes.append(axis)
	330	tm = self._gettime(i)
	331	src = self._getsourcename(i)
	332	out += 'Scan[%d] (%s) ' % (axis[0], src)
	333	out += 'Time[%s]:\n' % (tm)
	334	if self.nbeam(-1) > 1: out += ' Beam[%d] ' % (axis[1])
	335	if self.nif(-1) > 1: out += ' IF[%d] ' % (axis[2])
	336	if self.npol(-1) > 1: out += ' Pol[%d] ' % (axis[3])
	337	out += '= %3.3f\n' % (statvals[i])
	338	out += "--------------------------------------------------\n"
[256]	339
[876]	340	if rcParams['verbose']:
	341	print "--------------------------------------------------"
	342	print " ",stat
	343	print "--------------------------------------------------"
	344	print out
	345	retval = { 'axesnames': ['scanno','beamno','ifno','polno','cycleno'],
	346	'axes' : axes,
	347	'data': statvals}
	348	return retval
[102]	349
[876]	350	def stddev(self,mask=None):
[135]	351	"""
	352	Determine the standard deviation of the current beam/if/pol
	353	Takes a 'mask' as an optional parameter to specify which
	354	channels should be excluded.
	355	Parameters:
	356	mask: an optional mask specifying where the standard
	357	deviation should be determined.
	358
	359	Example:
	360	scan.set_unit('channel')
	361	msk = scan.create_mask([100,200],[500,600])
	362	scan.stddev(mask=m)
	363	"""
[876]	364	return self.stats(stat='stddev',mask=mask);
[135]	365
[1003]	366
	367	def column_names(self):
	368	"""
	369	Return a list of column names, which can be used for selection.
	370	"""
	371	return list(Scantable.column_names(self))
	372
[876]	373	def get_tsys(self):
[113]	374	"""
	375	Return the System temperatures.
	376	Parameters:
[876]	377
[113]	378	Returns:
[876]	379	a list of Tsys values for the current selection
[113]	380	"""
[256]	381
[876]	382	return self._row_callback(self._gettsys, "Tsys")
[256]	383
[876]	384	def _row_callback(self, callback, label):
	385	axes = []
	386	axesnames = ['scanno','beamno','ifno','polno','cycleno']
	387	out = ""
	388	outvec =[]
	389	for i in range(self.nrow()):
	390	axis = []
	391	axis.append(self.getscan(i))
	392	axis.append(self.getbeam(i))
	393	axis.append(self.getif(i))
	394	axis.append(self.getpol(i))
	395	axis.append(self.getcycle(i))
	396	axes.append(axis)
	397	tm = self._gettime(i)
	398	src = self._getsourcename(i)
	399	out += 'Scan[%d] (%s) ' % (axis[0], src)
	400	out += 'Time[%s]:\n' % (tm)
	401	if self.nbeam(-1) > 1: out += ' Beam[%d] ' % (axis[1])
	402	if self.nif(-1) > 1: out += ' IF[%d] ' % (axis[2])
	403	if self.npol(-1) > 1: out += ' Pol[%d] ' % (axis[3])
	404	outvec.append(callback(i))
	405	out += '= %3.3f\n' % (outvec[i])
	406	out += "--------------------------------------------------\n"
	407	if rcParams['verbose']:
	408	print "--------------------------------------------------"
	409	print " %s" % (label)
	410	print "--------------------------------------------------"
	411	print out
	412	retval = {'axesnames': axesnames, 'axes': axes, 'data': outvec}
	413	return retval
[256]	414
[1070]	415	def _get_column(self, callback, row=-1):
	416	"""
	417	"""
	418	if row == -1:
	419	return [callback(i) for i in range(self.nrow())]
	420	else:
	421	if 0 <= row < self.nrow():
	422	return callback(row)
[256]	423
[1070]	424
[407]	425	def get_time(self, row=-1):
[113]	426	"""
	427	Get a list of time stamps for the observations.
[181]	428	Return a string for each integration in the scantable.
[113]	429	Parameters:
[407]	430	row: row no of integration. Default -1 return all rows
[113]	431	Example:
	432	none
	433	"""
[1070]	434	return self._get_column(self._gettime, row)
[102]	435
[714]	436	def get_sourcename(self, row=-1):
	437	"""
[794]	438	Get a list source names for the observations.
[714]	439	Return a string for each integration in the scantable.
	440	Parameters:
	441	row: row no of integration. Default -1 return all rows
	442	Example:
	443	none
	444	"""
[1070]	445	return self._get_column(self._getsourcename, row)
[714]	446
[794]	447	def get_elevation(self, row=-1):
	448	"""
	449	Get a list of elevations for the observations.
	450	Return a float for each integration in the scantable.
	451	Parameters:
	452	row: row no of integration. Default -1 return all rows
	453	Example:
	454	none
	455	"""
[1070]	456	return self._get_column(self._getelevation, row)
[794]	457
	458	def get_azimuth(self, row=-1):
	459	"""
	460	Get a list of azimuths for the observations.
	461	Return a float for each integration in the scantable.
	462	Parameters:
	463	row: row no of integration. Default -1 return all rows
	464	Example:
	465	none
	466	"""
[1070]	467	return self._get_column(self._getazimuth, row)
[794]	468
	469	def get_parangle(self, row=-1):
	470	"""
	471	Get a list of parallactic angles for the observations.
	472	Return a float for each integration in the scantable.
	473	Parameters:
	474	row: row no of integration. Default -1 return all rows
	475	Example:
	476	none
	477	"""
[1070]	478	return self._get_column(self._getparangle, row)
[794]	479
[1070]	480	def get_direction(self, row=-1):
	481	"""
	482	Get a list of Positions on the sky (direction) for the observations.
	483	Return a float for each integration in the scantable.
	484	Parameters:
	485	row: row no of integration. Default -1 return all rows
	486	Example:
	487	none
	488	"""
	489	return self._get_column(self._getdirection, row)
	490
[102]	491	def set_unit(self, unit='channel'):
	492	"""
	493	Set the unit for all following operations on this scantable
	494	Parameters:
	495	unit: optional unit, default is 'channel'
[113]	496	one of '*Hz','km/s','channel', ''
[102]	497	"""
[484]	498	varlist = vars()
[113]	499	if unit in ['','pixel', 'channel']:
	500	unit = ''
	501	inf = list(self._getcoordinfo())
	502	inf[0] = unit
	503	self._setcoordinfo(inf)
[484]	504	self._add_history("set_unit",varlist)
[113]	505
[484]	506	def set_instrument(self, instr):
[358]	507	"""
	508	Set the instrument for subsequent processing
	509	Parameters:
[710]	510	instr: Select from 'ATPKSMB', 'ATPKSHOH', 'ATMOPRA',
[407]	511	'DSS-43' (Tid), 'CEDUNA', and 'HOBART'
[358]	512	"""
	513	self._setInstrument(instr)
[484]	514	self._add_history("set_instument",vars())
[718]	515	print_log()
[358]	516
[276]	517	def set_doppler(self, doppler='RADIO'):
	518	"""
	519	Set the doppler for all following operations on this scantable.
	520	Parameters:
	521	doppler: One of 'RADIO', 'OPTICAL', 'Z', 'BETA', 'GAMMA'
	522	"""
[484]	523	varlist = vars()
[276]	524	inf = list(self._getcoordinfo())
	525	inf[2] = doppler
	526	self._setcoordinfo(inf)
[484]	527	self._add_history("set_doppler",vars())
[718]	528	print_log()
[710]	529
[226]	530	def set_freqframe(self, frame=None):
[113]	531	"""
	532	Set the frame type of the Spectral Axis.
	533	Parameters:
[591]	534	frame: an optional frame type, default 'LSRK'. Valid frames are:
	535	'REST','TOPO','LSRD','LSRK','BARY',
	536	'GEO','GALACTO','LGROUP','CMB'
[113]	537	Examples:
	538	scan.set_freqframe('BARY')
	539	"""
[484]	540	if frame is None: frame = rcParams['scantable.freqframe']
	541	varlist = vars()
[113]	542	valid = ['REST','TOPO','LSRD','LSRK','BARY', \
	543	'GEO','GALACTO','LGROUP','CMB']
[591]	544
[989]	545	if frame in valid:
[113]	546	inf = list(self._getcoordinfo())
	547	inf[1] = frame
	548	self._setcoordinfo(inf)
[484]	549	self._add_history("set_freqframe",varlist)
[102]	550	else:
[718]	551	msg = "Please specify a valid freq type. Valid types are:\n",valid
	552	if rcParams['verbose']:
	553	print msg
	554	else:
	555	raise TypeError(msg)
	556	print_log()
[710]	557
[989]	558	def set_dirframe(self, frame=""):
	559	"""
	560	Set the frame type of the Direction on the sky.
	561	Parameters:
	562	frame: an optional frame type, default ''. Valid frames are:
	563	'J2000', 'B1950', 'GALACTIC'
	564	Examples:
	565	scan.set_dirframe('GALACTIC')
	566	"""
	567	varlist = vars()
	568	try:
	569	Scantable.set_dirframe(self, frame)
	570	except RuntimeError,msg:
	571	if rcParams['verbose']:
	572	print msg
	573	else:
	574	raise
	575	self._add_history("set_dirframe",varlist)
	576
[113]	577	def get_unit(self):
	578	"""
	579	Get the default unit set in this scantable
	580	Parameters:
	581	Returns:
	582	A unit string
	583	"""
	584	inf = self._getcoordinfo()
	585	unit = inf[0]
	586	if unit == '': unit = 'channel'
	587	return unit
[102]	588
[158]	589	def get_abcissa(self, rowno=0):
[102]	590	"""
[158]	591	Get the abcissa in the current coordinate setup for the currently
[113]	592	selected Beam/IF/Pol
	593	Parameters:
[226]	594	rowno: an optional row number in the scantable. Default is the
	595	first row, i.e. rowno=0
[113]	596	Returns:
[407]	597	The abcissa values and it's format string (as a dictionary)
[113]	598	"""
[256]	599	abc = self._getabcissa(rowno)
[710]	600	lbl = self._getabcissalabel(rowno)
[718]	601	print_log()
[158]	602	return abc, lbl
[113]	603
[1001]	604	def flag(self, mask=[]):
	605	"""
	606	Flag the selected data using an optional channel mask.
	607	Parameters:
	608	mask: an optional channel mask, created with create_mask. Default
	609	(no mask) is all channels.
	610	"""
	611	varlist = vars()
	612	try:
	613	self._flag(mask)
	614	except RuntimeError,msg:
	615	if rcParams['verbose']:
	616	print msg
	617	return
	618	else: raise
	619	self._add_history("flag", varlist)
	620
	621
[113]	622	def create_mask(self, args, *kwargs):
	623	"""
[102]	624	Compute and return a mask based on [min,max] windows.
[189]	625	The specified windows are to be INCLUDED, when the mask is
[113]	626	applied.
[102]	627	Parameters:
	628	[min,max],[min2,max2],...
[1024]	629	Pairs of start/end points (inclusive)specifying the regions
[102]	630	to be masked
[189]	631	invert: optional argument. If specified as True,
	632	return an inverted mask, i.e. the regions
	633	specified are EXCLUDED
[513]	634	row: create the mask using the specified row for
	635	unit conversions, default is row=0
	636	only necessary if frequency varies over rows.
[102]	637	Example:
[113]	638	scan.set_unit('channel')
	639	a)
[1024]	640	msk = scan.create_mask([400,500],[800,900])
[189]	641	# masks everything outside 400 and 500
[113]	642	# and 800 and 900 in the unit 'channel'
	643
	644	b)
[1024]	645	msk = scan.create_mask([400,500],[800,900], invert=True)
[189]	646	# masks the regions between 400 and 500
[113]	647	# and 800 and 900 in the unit 'channel'
[1024]	648	c)
	649	mask only channel 400
	650	msk = scan.create_mask([400,400])
[102]	651	"""
[513]	652	row = 0
	653	if kwargs.has_key("row"):
	654	row = kwargs.get("row")
	655	data = self._getabcissa(row)
[113]	656	u = self._getcoordinfo()[0]
[718]	657	if rcParams['verbose']:
[113]	658	if u == "": u = "channel"
[718]	659	from asap import asaplog
	660	msg = "The current mask window unit is %s" % u
[914]	661	if not self._check_ifs():
[876]	662	msg += "\nThis mask is only valid for IF=%d" % (self.getif(i))
[718]	663	asaplog.push(msg)
[102]	664	n = self.nchan()
[189]	665	msk = zeros(n)
[710]	666	# test if args is a 'list' or a 'normal *args - UGLY!!!
	667
	668	ws = (isinstance(args[-1][-1],int) or isinstance(args[-1][-1],float)) and args or args[0]
	669	for window in ws:
[102]	670	if (len(window) != 2 or window[0] > window[1] ):
[718]	671	raise TypeError("A window needs to be defined as [min,max]")
[102]	672	for i in range(n):
[1024]	673	if data[i] >= window[0] and data[i] <= window[1]:
[189]	674	msk[i] = 1
[113]	675	if kwargs.has_key('invert'):
	676	if kwargs.get('invert'):
	677	from numarray import logical_not
[710]	678	msk = logical_not(msk)
[718]	679	print_log()
[102]	680	return msk
[710]	681
[256]	682	def get_restfreqs(self):
	683	"""
	684	Get the restfrequency(s) stored in this scantable.
	685	The return value(s) are always of unit 'Hz'
	686	Parameters:
	687	none
	688	Returns:
	689	a list of doubles
	690	"""
	691	return list(self._getrestfreqs())
[102]	692
[402]	693
[931]	694	def set_restfreqs(self, freqs=None, unit='Hz'):
	695	"""
	696	Set or replace the restfrequency specified and
	697	If the 'freqs' argument holds a scalar,
	698	then that rest frequency will be applied to all the selected
	699	data. If the 'freqs' argument holds
	700	a vector, then it MUST be of equal or smaller length than
	701	the number of IFs (and the available restfrequencies will be
	702	replaced by this vector). In this case, all data have
	703	the restfrequency set per IF according
	704	to the corresponding value you give in the 'freqs' vector.
	705	E.g. 'freqs=[1e9,2e9]' would mean IF 0 gets restfreq 1e9 and
	706	IF 1 gets restfreq 2e9.
	707	You can also specify the frequencies via known line names
	708	from the built-in Lovas table.
	709	Parameters:
	710	freqs: list of rest frequency values or string idenitfiers
	711	unit: unit for rest frequency (default 'Hz')
[402]	712
[931]	713	Example:
	714	# set the given restfrequency for the whole table
	715	scan.set_restfreqs(freqs=1.4e9)
	716	# If thee number of IFs in the data is >= 2 the IF0 gets the first
	717	# value IF1 the second...
	718	scan.set_restfreqs(freqs=[1.4e9,1.67e9])
	719	#set the given restfrequency for the whole table (by name)
	720	scan.set_restfreqs(freqs="OH1667")
[391]	721
[931]	722	Note:
	723	To do more sophisticate Restfrequency setting, e.g. on a
	724	source and IF basis, use scantable.set_selection() before using
	725	this function.
	726	# provide your scantable is call scan
	727	selection = selector()
	728	selection.set_name("ORION*")
	729	selection.set_ifs([1])
	730	scan.set_selection(selection)
	731	scan.set_restfreqs(freqs=86.6e9)
	732
	733	"""
	734	varlist = vars()
	735
	736	t = type(freqs)
	737	if isinstance(freqs, int) or isinstance(freqs,float):
	738	self._setrestfreqs(freqs, unit)
	739	elif isinstance(freqs, list) or isinstance(freqs,tuple):
	740	if isinstance(freqs[-1], int) or isinstance(freqs[-1],float):
	741	sel = selector()
	742	savesel = self._getselection()
	743	for i in xrange(len(freqs)):
[946]	744	sel.set_ifs([i])
[931]	745	self._setselection(sel)
	746	self._setrestfreqs(freqs[i], unit)
	747	self._setselection(savesel)
	748	elif isinstance(freqs[-1], str):
	749	# not yet implemented
	750	pass
	751	else:
	752	return
	753	self._add_history("set_restfreqs", varlist)
	754
	755
	756
[484]	757	def history(self):
	758	hist = list(self._gethistory())
[794]	759	out = "-"*80
[484]	760	for h in hist:
[489]	761	if h.startswith("---"):
[794]	762	out += "\n"+h
[489]	763	else:
	764	items = h.split("##")
	765	date = items[0]
	766	func = items[1]
	767	items = items[2:]
[794]	768	out += "\n"+date+"\n"
	769	out += "Function: %s\n Parameters:" % (func)
[489]	770	for i in items:
	771	s = i.split("=")
[794]	772	out += "\n %s = %s" % (s[0],s[1])
	773	out += "\n"+"-"*80
	774	try:
	775	from IPython.genutils import page as pager
	776	except ImportError:
	777	from pydoc import pager
	778	pager(out)
[484]	779	return
	780
[513]	781	#
	782	# Maths business
	783	#
	784
[931]	785	def average_time(self, mask=None, scanav=False, weight='tint', align=False):
[513]	786	"""
[1070]	787	Return the (time) weighted average of a scan.
[513]	788	Note:
[1070]	789	in channels only - align if necessary
[513]	790	Parameters:
	791	one scan or comma separated scans
	792	mask: an optional mask (only used for 'var' and 'tsys'
	793	weighting)
[558]	794	scanav: True averages each scan separately
	795	False (default) averages all scans together,
[1099]	796	weight: Weighting scheme.
	797	'none' (mean no weight)
	798	'var' (1/var(spec) weighted)
	799	'tsys' (1/Tsys**2 weighted)
	800	'tint' (integration time weighted)
	801	'tintsys' (Tint/Tsys**2)
	802	'median' ( median averaging)
[535]	803	The default is 'tint'
[931]	804	align: align the spectra in velocity before averaging. It takes
	805	the time of the first spectrum as reference time.
[513]	806	Example:
	807	# time average the scantable without using a mask
[710]	808	newscan = scan.average_time()
[513]	809	"""
	810	varlist = vars()
[976]	811	if weight is None: weight = 'TINT'
[513]	812	if mask is None: mask = ()
[1099]	813	if scanav: scanav = "SCAN"
	814	else: scanav = "NONE"
[981]	815	scan = (self,)
[989]	816	try:
	817	if align:
	818	scan = (self.freq_align(insitu=False),)
[1099]	819	s = None
	820	if weight.upper() == 'MEDIAN':
	821	s = scantable(self._math._averagechannel(scan[0], 'MEDIAN', scanav))
	822	else:
	823	s = scantable(self._math._average(scan, mask, weight.upper(),
[989]	824	scanav))
	825	except RuntimeError,msg:
	826	if rcParams['verbose']:
	827	print msg
	828	return
	829	else: raise
[1099]	830	s._add_history("average_time", varlist)
[718]	831	print_log()
[513]	832	return s
[710]	833
[876]	834	def convert_flux(self, jyperk=None, eta=None, d=None, insitu=None):
[513]	835	"""
	836	Return a scan where all spectra are converted to either
	837	Jansky or Kelvin depending upon the flux units of the scan table.
	838	By default the function tries to look the values up internally.
	839	If it can't find them (or if you want to over-ride), you must
	840	specify EITHER jyperk OR eta (and D which it will try to look up
	841	also if you don't set it). jyperk takes precedence if you set both.
	842	Parameters:
	843	jyperk: the Jy / K conversion factor
	844	eta: the aperture efficiency
	845	d: the geomtric diameter (metres)
	846	insitu: if False a new scantable is returned.
	847	Otherwise, the scaling is done in-situ
	848	The default is taken from .asaprc (False)
	849	allaxes: if True apply to all spectra. Otherwise
	850	apply only to the selected (beam/pol/if)spectra only
	851	The default is taken from .asaprc (True if none)
	852	"""
	853	if insitu is None: insitu = rcParams['insitu']
[876]	854	self._math._setinsitu(insitu)
[513]	855	varlist = vars()
	856	if jyperk is None: jyperk = -1.0
	857	if d is None: d = -1.0
	858	if eta is None: eta = -1.0
[876]	859	s = scantable(self._math._convertflux(self, d, eta, jyperk))
	860	s._add_history("convert_flux", varlist)
	861	print_log()
	862	if insitu: self._assign(s)
	863	else: return s
[513]	864
[876]	865	def gain_el(self, poly=None, filename="", method="linear", insitu=None):
[513]	866	"""
	867	Return a scan after applying a gain-elevation correction.
	868	The correction can be made via either a polynomial or a
	869	table-based interpolation (and extrapolation if necessary).
	870	You specify polynomial coefficients, an ascii table or neither.
	871	If you specify neither, then a polynomial correction will be made
	872	with built in coefficients known for certain telescopes (an error
	873	will occur if the instrument is not known).
	874	The data and Tsys are divided by the scaling factors.
	875	Parameters:
	876	poly: Polynomial coefficients (default None) to compute a
	877	gain-elevation correction as a function of
	878	elevation (in degrees).
	879	filename: The name of an ascii file holding correction factors.
	880	The first row of the ascii file must give the column
	881	names and these MUST include columns
	882	"ELEVATION" (degrees) and "FACTOR" (multiply data
	883	by this) somewhere.
	884	The second row must give the data type of the
	885	column. Use 'R' for Real and 'I' for Integer.
	886	An example file would be
	887	(actual factors are arbitrary) :
	888
	889	TIME ELEVATION FACTOR
	890	R R R
	891	0.1 0 0.8
	892	0.2 20 0.85
	893	0.3 40 0.9
	894	0.4 60 0.85
	895	0.5 80 0.8
	896	0.6 90 0.75
	897	method: Interpolation method when correcting from a table.
	898	Values are "nearest", "linear" (default), "cubic"
	899	and "spline"
	900	insitu: if False a new scantable is returned.
	901	Otherwise, the scaling is done in-situ
	902	The default is taken from .asaprc (False)
	903	"""
	904
	905	if insitu is None: insitu = rcParams['insitu']
[876]	906	self._math._setinsitu(insitu)
[513]	907	varlist = vars()
	908	if poly is None:
	909	poly = ()
	910	from os.path import expandvars
	911	filename = expandvars(filename)
[876]	912	s = scantable(self._math._gainel(self, poly, filename, method))
	913	s._add_history("gain_el", varlist)
	914	print_log()
	915	if insitu: self._assign(s)
	916	else: return s
[710]	917
[931]	918	def freq_align(self, reftime=None, method='cubic', insitu=None):
[513]	919	"""
	920	Return a scan where all rows have been aligned in frequency/velocity.
	921	The alignment frequency frame (e.g. LSRK) is that set by function
	922	set_freqframe.
	923	Parameters:
	924	reftime: reference time to align at. By default, the time of
	925	the first row of data is used.
	926	method: Interpolation method for regridding the spectra.
	927	Choose from "nearest", "linear", "cubic" (default)
	928	and "spline"
	929	insitu: if False a new scantable is returned.
	930	Otherwise, the scaling is done in-situ
	931	The default is taken from .asaprc (False)
	932	"""
[931]	933	if insitu is None: insitu = rcParams["insitu"]
[876]	934	self._math._setinsitu(insitu)
[513]	935	varlist = vars()
[931]	936	if reftime is None: reftime = ""
	937	s = scantable(self._math._freq_align(self, reftime, method))
[876]	938	s._add_history("freq_align", varlist)
	939	print_log()
	940	if insitu: self._assign(s)
	941	else: return s
[513]	942
[876]	943	def opacity(self, tau, insitu=None):
[513]	944	"""
	945	Apply an opacity correction. The data
	946	and Tsys are multiplied by the correction factor.
	947	Parameters:
	948	tau: Opacity from which the correction factor is
	949	exp(tau*ZD)
	950	where ZD is the zenith-distance
	951	insitu: if False a new scantable is returned.
	952	Otherwise, the scaling is done in-situ
	953	The default is taken from .asaprc (False)
	954	"""
	955	if insitu is None: insitu = rcParams['insitu']
[876]	956	self._math._setinsitu(insitu)
[513]	957	varlist = vars()
[876]	958	s = scantable(self._math._opacity(self, tau))
	959	s._add_history("opacity", varlist)
	960	print_log()
	961	if insitu: self._assign(s)
	962	else: return s
[513]	963
	964	def bin(self, width=5, insitu=None):
	965	"""
	966	Return a scan where all spectra have been binned up.
	967	width: The bin width (default=5) in pixels
	968	insitu: if False a new scantable is returned.
	969	Otherwise, the scaling is done in-situ
	970	The default is taken from .asaprc (False)
	971	"""
	972	if insitu is None: insitu = rcParams['insitu']
[876]	973	self._math._setinsitu(insitu)
[513]	974	varlist = vars()
[876]	975	s = scantable(self._math._bin(self, width))
	976	s._add_history("bin",varlist)
	977	print_log()
	978	if insitu: self._assign(s)
	979	else: return s
[513]	980
[710]	981
[513]	982	def resample(self, width=5, method='cubic', insitu=None):
	983	"""
	984	Return a scan where all spectra have been binned up
	985	width: The bin width (default=5) in pixels
	986	method: Interpolation method when correcting from a table.
	987	Values are "nearest", "linear", "cubic" (default)
	988	and "spline"
	989	insitu: if False a new scantable is returned.
	990	Otherwise, the scaling is done in-situ
	991	The default is taken from .asaprc (False)
	992	"""
	993	if insitu is None: insitu = rcParams['insitu']
[876]	994	self._math._setinsitu(insitu)
[513]	995	varlist = vars()
[876]	996	s = scantable(self._math._resample(self, method, width))
	997	s._add_history("resample",varlist)
	998	print_log()
	999	if insitu: self._assign(s)
	1000	else: return s
[513]	1001
	1002
[946]	1003	def average_pol(self, mask=None, weight='none'):
	1004	"""
	1005	Average the Polarisations together.
	1006	Parameters:
	1007	mask: An optional mask defining the region, where the
	1008	averaging will be applied. The output will have all
	1009	specified points masked.
	1010	weight: Weighting scheme. 'none' (default), 'var' (1/var(spec)
	1011	weighted), or 'tsys' (1/Tsys**2 weighted)
	1012	"""
	1013	varlist = vars()
	1014	if mask is None:
	1015	mask = ()
[1010]	1016	s = scantable(self._math._averagepol(self, mask, weight.upper()))
[946]	1017	s._add_history("average_pol",varlist)
	1018	print_log()
[992]	1019	return s
[513]	1020
[992]	1021	def convert_pol(self, poltype=None):
	1022	"""
	1023	Convert the data to a different polarisation type.
	1024	Parameters:
	1025	poltype: The new polarisation type. Valid types are:
	1026	"linear", "stokes" and "circular"
	1027	"""
	1028	varlist = vars()
	1029	try:
	1030	s = scantable(self._math._convertpol(self, poltype))
	1031	except RuntimeError,msg:
	1032	if rcParams['verbose']:
	1033	print msg
	1034	return
	1035	else:
	1036	raise
	1037	s._add_history("convert_pol",varlist)
	1038	print_log()
	1039	return s
	1040
[876]	1041	def smooth(self, kernel="hanning", width=5.0, insitu=None):
[513]	1042	"""
	1043	Smooth the spectrum by the specified kernel (conserving flux).
	1044	Parameters:
	1045	scan: The input scan
	1046	kernel: The type of smoothing kernel. Select from
	1047	'hanning' (default), 'gaussian' and 'boxcar'.
	1048	The first three characters are sufficient.
	1049	width: The width of the kernel in pixels. For hanning this is
	1050	ignored otherwise it defauls to 5 pixels.
	1051	For 'gaussian' it is the Full Width Half
	1052	Maximum. For 'boxcar' it is the full width.
	1053	insitu: if False a new scantable is returned.
	1054	Otherwise, the scaling is done in-situ
	1055	The default is taken from .asaprc (False)
	1056	Example:
	1057	none
	1058	"""
	1059	if insitu is None: insitu = rcParams['insitu']
[876]	1060	self._math._setinsitu(insitu)
[513]	1061	varlist = vars()
[1033]	1062	s = scantable(self._math._smooth(self,kernel.lower(),width))
[876]	1063	s._add_history("smooth", varlist)
	1064	print_log()
	1065	if insitu: self._assign(s)
	1066	else: return s
[513]	1067
[876]	1068
[1061]	1069	def poly_baseline(self, mask=None, order=0, plot=False, insitu=None):
[513]	1070	"""
	1071	Return a scan which has been baselined (all rows) by a polynomial.
	1072	Parameters:
[794]	1073	scan: a scantable
	1074	mask: an optional mask
	1075	order: the order of the polynomial (default is 0)
[1061]	1076	plot: plot the fit and the residual. In this each
	1077	indivual fit has to be approved, by typing 'y'
	1078	or 'n'
[794]	1079	insitu: if False a new scantable is returned.
	1080	Otherwise, the scaling is done in-situ
	1081	The default is taken from .asaprc (False)
[513]	1082	Example:
	1083	# return a scan baselined by a third order polynomial,
	1084	# not using a mask
	1085	bscan = scan.poly_baseline(order=3)
[579]	1086	"""
[513]	1087	if insitu is None: insitu = rcParams['insitu']
	1088	varlist = vars()
	1089	if mask is None:
	1090	from numarray import ones
[876]	1091	mask = list(ones(self.nchan(-1)))
[513]	1092	from asap.asapfitter import fitter
	1093	f = fitter()
	1094	f.set_scan(self, mask)
	1095	f.set_function(poly=order)
[1061]	1096	s = f.auto_fit(insitu, plot=plot)
[876]	1097	s._add_history("poly_baseline", varlist)
[718]	1098	print_log()
[876]	1099	if insitu: self._assign(s)
	1100	else: return s
[513]	1101
[931]	1102	def auto_poly_baseline(self, mask=[], edge=(0,0), order=0,
[1061]	1103	threshold=3, plot=False, insitu=None):
[880]	1104	"""
	1105	Return a scan which has been baselined (all rows) by a polynomial.
	1106	Spectral lines are detected first using linefinder and masked out
	1107	to avoid them affecting the baseline solution.
	1108
	1109	Parameters:
	1110	mask: an optional mask retreived from scantable
	1111	edge: an optional number of channel to drop at
	1112	the edge of spectrum. If only one value is
	1113	specified, the same number will be dropped from
	1114	both sides of the spectrum. Default is to keep
[907]	1115	all channels. Nested tuples represent individual
[976]	1116	edge selection for different IFs (a number of spectral
	1117	channels can be different)
[880]	1118	order: the order of the polynomial (default is 0)
	1119	threshold: the threshold used by line finder. It is better to
	1120	keep it large as only strong lines affect the
	1121	baseline solution.
[1061]	1122	plot: plot the fit and the residual. In this each
	1123	indivual fit has to be approved, by typing 'y'
	1124	or 'n'
[880]	1125	insitu: if False a new scantable is returned.
	1126	Otherwise, the scaling is done in-situ
	1127	The default is taken from .asaprc (False)
	1128
	1129	Example:
	1130	scan2=scan.auto_poly_baseline(order=7)
	1131	"""
	1132	if insitu is None: insitu = rcParams['insitu']
	1133	varlist = vars()
	1134	from asap.asapfitter import fitter
	1135	from asap.asaplinefind import linefinder
	1136	from asap import _is_sequence_or_number as _is_valid
	1137
[976]	1138	# check whether edge is set up for each IF individually
	1139	individualEdge = False;
	1140	if len(edge)>1:
	1141	if isinstance(edge[0],list) or isinstance(edge[0],tuple):
	1142	individualEdge = True;
[907]	1143
	1144	if not _is_valid(edge, int) and not individualEdge:
[909]	1145	raise ValueError, "Parameter 'edge' has to be an integer or a \
[907]	1146	pair of integers specified as a tuple. Nested tuples are allowed \
	1147	to make individual selection for different IFs."
[919]	1148
[976]	1149	curedge = (0,0)
	1150	if individualEdge:
	1151	for edge_par in edge:
	1152	if not _is_valid(edge,int):
	1153	raise ValueError, "Each element of the 'edge' tuple has \
[907]	1154	to be a pair of integers or an integer."
	1155	else:
[976]	1156	curedge = edge;
[880]	1157
	1158	# setup fitter
	1159	f = fitter()
	1160	f.set_function(poly=order)
	1161
	1162	# setup line finder
	1163	fl=linefinder()
	1164	fl.set_options(threshold=threshold)
	1165
	1166	if not insitu:
	1167	workscan=self.copy()
	1168	else:
	1169	workscan=self
	1170
[907]	1171	fl.set_scan(workscan)
	1172
[880]	1173	rows=range(workscan.nrow())
	1174	from asap import asaplog
	1175	asaplog.push("Processing:")
	1176	for r in rows:
	1177	msg = " Scan[%d] Beam[%d] IF[%d] Pol[%d] Cycle[%d]" % (workscan.getscan(r),workscan.getbeam(r),workscan.getif(r),workscan.getpol(r), workscan.getcycle(r))
	1178	asaplog.push(msg, False)
[907]	1179
[976]	1180	# figure out edge parameter
	1181	if individualEdge:
	1182	if len(edge)>=workscan.getif(r):
	1183	raise RuntimeError, "Number of edge elements appear to be less than the number of IFs"
	1184	curedge = edge[workscan.getif(r)]
[919]	1185
[976]	1186	# setup line finder
[907]	1187	fl.find_lines(r,mask,curedge)
[880]	1188	f.set_scan(workscan, fl.get_mask())
	1189	f.x = workscan._getabcissa(r)
	1190	f.y = workscan._getspectrum(r)
	1191	f.data = None
	1192	f.fit()
	1193	x = f.get_parameters()
[1061]	1194	if plot:
	1195	f.plot(residual=True)
	1196	x = raw_input("Accept fit ( [y]/n ): ")
	1197	if x.upper() == 'N':
	1198	continue
[880]	1199	workscan._setspectrum(f.fitter.getresidual(), r)
[1061]	1200	if plot:
	1201	f._p.unmap()
	1202	f._p = None
	1203	workscan._add_history("auto_poly_baseline", varlist)
[880]	1204	if insitu:
	1205	self._assign(workscan)
	1206	else:
	1207	return workscan
	1208
[914]	1209	def rotate_linpolphase(self, angle):
	1210	"""
	1211	Rotate the phase of the complex polarization O=Q+iU correlation.
	1212	This is always done in situ in the raw data. So if you call this
	1213	function more than once then each call rotates the phase further.
	1214	Parameters:
	1215	angle: The angle (degrees) to rotate (add) by.
	1216	Examples:
	1217	scan.rotate_linpolphase(2.3)
	1218	"""
	1219	varlist = vars()
[936]	1220	self._math._rotate_linpolphase(self, angle)
[914]	1221	self._add_history("rotate_linpolphase", varlist)
	1222	print_log()
	1223	return
[710]	1224
[513]	1225
[914]	1226	def rotate_xyphase(self, angle):
	1227	"""
	1228	Rotate the phase of the XY correlation. This is always done in situ
	1229	in the data. So if you call this function more than once
	1230	then each call rotates the phase further.
	1231	Parameters:
	1232	angle: The angle (degrees) to rotate (add) by.
	1233	Examples:
	1234	scan.rotate_xyphase(2.3)
	1235	"""
	1236	varlist = vars()
[936]	1237	self._math._rotate_xyphase(self, angle)
[914]	1238	self._add_history("rotate_xyphase", varlist)
	1239	print_log()
	1240	return
	1241
	1242	def swap_linears(self):
	1243	"""
	1244	Swap the linear polarisations XX and YY
	1245	"""
	1246	varlist = vars()
[936]	1247	self._math._swap_linears(self)
[914]	1248	self._add_history("swap_linears", varlist)
	1249	print_log()
	1250	return
	1251
	1252	def invert_phase(self):
	1253	"""
	1254	Invert the phase of the complex polarisation
	1255	"""
	1256	varlist = vars()
[936]	1257	self._math._invert_phase(self)
[914]	1258	self._add_history("invert_phase", varlist)
	1259	print_log()
	1260	return
	1261
[876]	1262	def add(self, offset, insitu=None):
[513]	1263	"""
	1264	Return a scan where all spectra have the offset added
	1265	Parameters:
	1266	offset: the offset
	1267	insitu: if False a new scantable is returned.
	1268	Otherwise, the scaling is done in-situ
	1269	The default is taken from .asaprc (False)
	1270	"""
	1271	if insitu is None: insitu = rcParams['insitu']
[876]	1272	self._math._setinsitu(insitu)
[513]	1273	varlist = vars()
[876]	1274	s = scantable(self._math._unaryop(self, offset, "ADD", False))
	1275	s._add_history("add",varlist)
	1276	print_log()
	1277	if insitu:
	1278	self._assign(s)
	1279	else:
[513]	1280	return s
	1281
[876]	1282	def scale(self, factor, tsys=True, insitu=None,):
[513]	1283	"""
	1284	Return a scan where all spectra are scaled by the give 'factor'
	1285	Parameters:
	1286	factor: the scaling factor
	1287	insitu: if False a new scantable is returned.
	1288	Otherwise, the scaling is done in-situ
	1289	The default is taken from .asaprc (False)
	1290	tsys: if True (default) then apply the operation to Tsys
	1291	as well as the data
	1292	"""
	1293	if insitu is None: insitu = rcParams['insitu']
[876]	1294	self._math._setinsitu(insitu)
[513]	1295	varlist = vars()
[876]	1296	s = scantable(self._math._unaryop(self, factor, "MUL", tsys))
	1297	s._add_history("scale",varlist)
	1298	print_log()
	1299	if insitu:
	1300	self._assign(s)
	1301	else:
[513]	1302	return s
	1303
[876]	1304	def auto_quotient(self, mode='time', preserve=True):
[670]	1305	"""
	1306	This function allows to build quotients automatically.
	1307	It assumes the observation to have the same numer of
	1308	"ons" and "offs"
	1309	It will support "closest off in time" in the future
	1310	Parameters:
	1311	mode: the on/off detection mode; 'suffix' (default)
	1312	'suffix' identifies 'off' scans by the
	1313	trailing '_R' (Mopra/Parkes) or
	1314	'_e'/'_w' (Tid)
[710]	1315	preserve: you can preserve (default) the continuum or
	1316	remove it. The equations used are
[670]	1317	preserve: Output = Toff * (on/off) - Toff
[1070]	1318	remove: Output = Toff * (on/off) - Ton
[670]	1319	"""
[876]	1320	modes = ["time"]
[670]	1321	if not mode in modes:
[876]	1322	msg = "please provide valid mode. Valid modes are %s" % (modes)
	1323	raise ValueError(msg)
	1324	varlist = vars()
[1070]	1325	s = scantable(self._math._auto_quotient(self, mode, preserve))
[876]	1326	s._add_history("auto_quotient",varlist)
	1327	print_log()
	1328	return s
[710]	1329
[513]	1330
[876]	1331
	1332
[718]	1333	def freq_switch(self, insitu=None):
	1334	"""
	1335	Apply frequency switching to the data.
	1336	Parameters:
	1337	insitu: if False a new scantable is returned.
	1338	Otherwise, the swictching is done in-situ
	1339	The default is taken from .asaprc (False)
	1340	Example:
	1341	none
	1342	"""
	1343	if insitu is None: insitu = rcParams['insitu']
[876]	1344	self._math._setinsitu(insitu)
[718]	1345	varlist = vars()
[876]	1346	s = scantable(self._math._freqswitch(self))
	1347	s._add_history("freq_switch",varlist)
	1348	print_log()
	1349	if insitu: self._assign(s)
	1350	else: return s
[718]	1351
[780]	1352	def recalc_azel(self):
	1353	"""
	1354	Recalculate the azimuth and elevation for each position.
	1355	Parameters:
	1356	none
	1357	Example:
	1358	"""
	1359	varlist = vars()
[876]	1360	self._recalcazel()
[780]	1361	self._add_history("recalc_azel", varlist)
	1362	print_log()
	1363	return
	1364
[513]	1365	def __add__(self, other):
	1366	varlist = vars()
	1367	s = None
	1368	if isinstance(other, scantable):
[876]	1369	print "scantable + scantable NYI"
	1370	return
[513]	1371	elif isinstance(other, float):
[876]	1372	s = scantable(self._math._unaryop(self, other, "ADD", False))
[513]	1373	else:
[718]	1374	raise TypeError("Other input is not a scantable or float value")
[513]	1375	s._add_history("operator +", varlist)
[718]	1376	print_log()
[513]	1377	return s
	1378
	1379	def __sub__(self, other):
	1380	"""
	1381	implicit on all axes and on Tsys
	1382	"""
	1383	varlist = vars()
	1384	s = None
	1385	if isinstance(other, scantable):
[876]	1386	print "scantable - scantable NYI"
	1387	return
[513]	1388	elif isinstance(other, float):
[876]	1389	s = scantable(self._math._unaryop(self, other, "SUB", False))
[513]	1390	else:
[718]	1391	raise TypeError("Other input is not a scantable or float value")
[513]	1392	s._add_history("operator -", varlist)
[718]	1393	print_log()
[513]	1394	return s
[710]	1395
[513]	1396	def __mul__(self, other):
	1397	"""
	1398	implicit on all axes and on Tsys
	1399	"""
	1400	varlist = vars()
	1401	s = None
	1402	if isinstance(other, scantable):
[876]	1403	print "scantable * scantable NYI"
	1404	return
[513]	1405	elif isinstance(other, float):
[876]	1406	s = scantable(self._math._unaryop(self, other, "MUL", False))
[513]	1407	else:
[718]	1408	raise TypeError("Other input is not a scantable or float value")
[513]	1409	s._add_history("operator *", varlist)
[718]	1410	print_log()
[513]	1411	return s
	1412
[710]	1413
[513]	1414	def __div__(self, other):
	1415	"""
	1416	implicit on all axes and on Tsys
	1417	"""
	1418	varlist = vars()
	1419	s = None
	1420	if isinstance(other, scantable):
[876]	1421	print "scantable / scantable NYI"
	1422	return
[513]	1423	elif isinstance(other, float):
	1424	if other == 0.0:
[718]	1425	raise ZeroDivisionError("Dividing by zero is not recommended")
[876]	1426	s = scantable(self._math._unaryop(self, other, "DIV", False))
[513]	1427	else:
[718]	1428	raise TypeError("Other input is not a scantable or float value")
[513]	1429	s._add_history("operator /", varlist)
[718]	1430	print_log()
[513]	1431	return s
	1432
[530]	1433	def get_fit(self, row=0):
	1434	"""
	1435	Print or return the stored fits for a row in the scantable
	1436	Parameters:
	1437	row: the row which the fit has been applied to.
	1438	"""
	1439	if row > self.nrow():
	1440	return
[976]	1441	from asap.asapfit import asapfit
[530]	1442	fit = asapfit(self._getfit(row))
[718]	1443	if rcParams['verbose']:
[530]	1444	print fit
	1445	return
	1446	else:
	1447	return fit.as_dict()
	1448
[484]	1449	def _add_history(self, funcname, parameters):
	1450	# create date
	1451	sep = "##"
	1452	from datetime import datetime
	1453	dstr = datetime.now().strftime('%Y/%m/%d %H:%M:%S')
	1454	hist = dstr+sep
	1455	hist += funcname+sep#cdate+sep
	1456	if parameters.has_key('self'): del parameters['self']
	1457	for k,v in parameters.iteritems():
	1458	if type(v) is dict:
	1459	for k2,v2 in v.iteritems():
	1460	hist += k2
	1461	hist += "="
	1462	if isinstance(v2,scantable):
	1463	hist += 'scantable'
	1464	elif k2 == 'mask':
[513]	1465	if isinstance(v2,list) or isinstance(v2,tuple):
	1466	hist += str(self._zip_mask(v2))
	1467	else:
	1468	hist += str(v2)
[484]	1469	else:
[513]	1470	hist += str(v2)
[484]	1471	else:
	1472	hist += k
	1473	hist += "="
	1474	if isinstance(v,scantable):
	1475	hist += 'scantable'
	1476	elif k == 'mask':
[513]	1477	if isinstance(v,list) or isinstance(v,tuple):
	1478	hist += str(self._zip_mask(v))
	1479	else:
	1480	hist += str(v)
[484]	1481	else:
	1482	hist += str(v)
	1483	hist += sep
	1484	hist = hist[:-2] # remove trailing '##'
	1485	self._addhistory(hist)
	1486
[710]	1487
[484]	1488	def _zip_mask(self, mask):
	1489	mask = list(mask)
	1490	i = 0
	1491	segments = []
	1492	while mask[i:].count(1):
	1493	i += mask[i:].index(1)
	1494	if mask[i:].count(0):
	1495	j = i + mask[i:].index(0)
	1496	else:
[710]	1497	j = len(mask)
[484]	1498	segments.append([i,j])
[710]	1499	i = j
[484]	1500	return segments
[714]	1501
[626]	1502	def _get_ordinate_label(self):
	1503	fu = "("+self.get_fluxunit()+")"
	1504	import re
	1505	lbl = "Intensity"
	1506	if re.match(".K.",fu):
	1507	lbl = "Brightness Temperature "+ fu
	1508	elif re.match(".Jy.",fu):
	1509	lbl = "Flux density "+ fu
	1510	return lbl
[710]	1511
[876]	1512	def _check_ifs(self):
	1513	nchans = [self.nchan(i) for i in range(self.nif(-1))]
[889]	1514	nchans = filter(lambda t: t > 0, nchans)
[876]	1515	return (sum(nchans)/len(nchans) == nchans[0])
[976]	1516
	1517	def _fill(self, names, unit, average):
	1518	import os
	1519	varlist = vars()
	1520	from asap._asap import stfiller
	1521	first = True
	1522	fullnames = []
	1523	for name in names:
	1524	name = os.path.expandvars(name)
	1525	name = os.path.expanduser(name)
	1526	if not os.path.exists(name):
	1527	msg = "File '%s' does not exists" % (name)
	1528	if rcParams['verbose']:
	1529	asaplog.push(msg)
	1530	print asaplog.pop().strip()
	1531	return
	1532	raise IOError(msg)
	1533	fullnames.append(name)
	1534	if average:
	1535	asaplog.push('Auto averaging integrations')
[1079]	1536	stype = int(rcParams['scantable.storage'].lower() == 'disk')
[976]	1537	for name in fullnames:
[1073]	1538	tbl = Scantable(stype)
	1539	r = stfiller(tbl)
[976]	1540	msg = "Importing %s..." % (name)
	1541	asaplog.push(msg,False)
	1542	print_log()
	1543	r._open(name,-1,-1)
	1544	r._read()
[1073]	1545	#tbl = r._getdata()
[976]	1546	if average:
[981]	1547	tbl = self._math._average((tbl,),(),'NONE','SCAN')
[976]	1548	#tbl = tbl2
	1549	if not first:
	1550	tbl = self._math._merge([self, tbl])
	1551	#tbl = tbl2
	1552	Scantable.__init__(self, tbl)
	1553	r._close()
	1554	del r,tbl
	1555	first = False
	1556	if unit is not None:
	1557	self.set_fluxunit(unit)
	1558	self.set_freqframe(rcParams['scantable.freqframe'])
	1559	#self._add_history("scantable", varlist)
	1560

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

Download in other formats: