Context Navigation

scantable.py @ 1685

Visit:

Last change on this file since 1685 was 1685, checked in by Takeshi Nakazato, 14 years ago

New Development: No

JIRA Issue: Yes CAS-1810

Ready to Release: Yes

Interface Changes: No

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: No

Module(s): Module Names change impacts.

Description: Describe your changes here...

Help updated.

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 86.2 KB

Line
1	from asap._asap import Scantable
2	from asap import rcParams
3	from asap import print_log
4	from asap import asaplog
5	from asap import selector
6	from asap import linecatalog
7	from asap import _n_bools, mask_not, mask_and, mask_or
8
9	class scantable(Scantable):
10	"""
11	The ASAP container for scans
12	"""
13
14	def __init__(self, filename, average=None, unit=None, getpt=None, antenna=None):
15	"""
16	Create a scantable from a saved one or make a reference
17	Parameters:
18	filename: the name of an asap table on disk
19	or
20	the name of a rpfits/sdfits/ms file
21	(integrations within scans are auto averaged
22	and the whole file is read)
23	or
24	[advanced] a reference to an existing
25	scantable
26	average: average all integrations withinb a scan on read.
27	The default (True) is taken from .asaprc.
28	unit: brightness unit; must be consistent with K or Jy.
29	Over-rides the default selected by the reader
30	(input rpfits/sdfits/ms) or replaces the value
31	in existing scantables
32	getpt: for MeasurementSet input data only:
33	If True, all pointing data are filled.
34	The deafult is False, which makes time to load
35	the MS data faster in some cases.
36	antenna: Antenna selection. integer (id) or string (name
37	or id).
38	"""
39	if average is None:
40	average = rcParams['scantable.autoaverage']
41	if getpt is None:
42	getpt = True
43	if antenna is None:
44	antenna = ''
45	elif type(antenna) == int:
46	antenna = '%s'%antenna
47	elif type(antenna) == list:
48	tmpstr = ''
49	for i in range( len(antenna) ):
50	if type(antenna[i]) == int:
51	tmpstr = tmpstr + ('%s,'%(antenna[i]))
52	elif type(antenna[i]) == str:
53	tmpstr=tmpstr+antenna[i]+','
54	else:
55	asaplog.push('Bad antenna selection.')
56	print_log('ERROR')
57	return
58	antenna = tmpstr.rstrip(',')
59	varlist = vars()
60	from asap._asap import stmath
61	self._math = stmath( rcParams['insitu'] )
62	if isinstance(filename, Scantable):
63	Scantable.__init__(self, filename)
64	else:
65	if isinstance(filename, str):# or \
66	# (isinstance(filename, list) or isinstance(filename, tuple)) \
67	# and isinstance(filename[-1], str):
68	import os.path
69	filename = os.path.expandvars(filename)
70	filename = os.path.expanduser(filename)
71	if not os.path.exists(filename):
72	s = "File '%s' not found." % (filename)
73	if rcParams['verbose']:
74	asaplog.push(s)
75	#print asaplog.pop().strip()
76	print_log('ERROR')
77	return
78	raise IOError(s)
79	if os.path.isdir(filename) \
80	and not os.path.exists(filename+'/table.f1'):
81	# crude check if asap table
82	if os.path.exists(filename+'/table.info'):
83	ondisk = rcParams['scantable.storage'] == 'disk'
84	Scantable.__init__(self, filename, ondisk)
85	if unit is not None:
86	self.set_fluxunit(unit)
87	# do not reset to the default freqframe
88	#self.set_freqframe(rcParams['scantable.freqframe'])
89	else:
90	msg = "The given file '%s'is not a valid " \
91	"asap table." % (filename)
92	if rcParams['verbose']:
93	#print msg
94	asaplog.push( msg )
95	print_log( 'ERROR' )
96	return
97	else:
98	raise IOError(msg)
99	else:
100	self._fill([filename], unit, average, getpt, antenna)
101	elif (isinstance(filename, list) or isinstance(filename, tuple)) \
102	and isinstance(filename[-1], str):
103	self._fill(filename, unit, average, getpt, antenna)
104	self._add_history("scantable", varlist)
105	print_log()
106
107	def save(self, name=None, format=None, overwrite=False):
108	"""
109	Store the scantable on disk. This can be an asap (aips++) Table,
110	SDFITS or MS2 format.
111	Parameters:
112	name: the name of the outputfile. For format "ASCII"
113	this is the root file name (data in 'name'.txt
114	and header in 'name'_header.txt)
115	format: an optional file format. Default is ASAP.
116	Allowed are - 'ASAP' (save as ASAP [aips++] Table),
117	'SDFITS' (save as SDFITS file)
118	'ASCII' (saves as ascii text file)
119	'MS2' (saves as an aips++
120	MeasurementSet V2)
121	'FITS' (save as image FITS - not
122	readable by class)
123	'CLASS' (save as FITS readable by CLASS)
124	overwrite: If the file should be overwritten if it exists.
125	The default False is to return with warning
126	without writing the output. USE WITH CARE.
127	Example:
128	scan.save('myscan.asap')
129	scan.save('myscan.sdfits', 'SDFITS')
130	"""
131	from os import path
132	if format is None: format = rcParams['scantable.save']
133	suffix = '.'+format.lower()
134	if name is None or name == "":
135	name = 'scantable'+suffix
136	msg = "No filename given. Using default name %s..." % name
137	asaplog.push(msg)
138	name = path.expandvars(name)
139	if path.isfile(name) or path.isdir(name):
140	if not overwrite:
141	msg = "File %s exists." % name
142	if rcParams['verbose']:
143	#print msg
144	asaplog.push( msg )
145	print_log( 'ERROR' )
146	return
147	else:
148	raise IOError(msg)
149	format2 = format.upper()
150	if format2 == 'ASAP':
151	self._save(name)
152	else:
153	from asap._asap import stwriter as stw
154	writer = stw(format2)
155	writer.write(self, name)
156	print_log()
157	return
158
159	def copy(self):
160	"""
161	Return a copy of this scantable.
162	Note:
163	This makes a full (deep) copy. scan2 = scan1 makes a reference.
164	Parameters:
165	none
166	Example:
167	copiedscan = scan.copy()
168	"""
169	sd = scantable(Scantable._copy(self))
170	return sd
171
172	def drop_scan(self, scanid=None):
173	"""
174	Return a new scantable where the specified scan number(s) has(have)
175	been dropped.
176	Parameters:
177	scanid: a (list of) scan number(s)
178	"""
179	from asap import _is_sequence_or_number as _is_valid
180	from asap import _to_list
181	from asap import unique
182	if not _is_valid(scanid):
183	if rcParams['verbose']:
184	#print "Please specify a scanno to drop from the scantable"
185	asaplog.push( 'Please specify a scanno to drop from the scantable' )
186	print_log( 'ERROR' )
187	return
188	else:
189	raise RuntimeError("No scan given")
190	try:
191	scanid = _to_list(scanid)
192	allscans = unique([ self.getscan(i) for i in range(self.nrow())])
193	for sid in scanid: allscans.remove(sid)
194	if len(allscans) == 0:
195	raise ValueError("Can't remove all scans")
196	except ValueError:
197	if rcParams['verbose']:
198	#print "Couldn't find any match."
199	print_log()
200	asaplog.push( "Couldn't find any match." )
201	print_log( 'ERROR' )
202	return
203	else: raise
204	try:
205	bsel = self.get_selection()
206	sel = selector()
207	sel.set_scans(allscans)
208	self.set_selection(bsel+sel)
209	scopy = self._copy()
210	self.set_selection(bsel)
211	return scantable(scopy)
212	except RuntimeError:
213	if rcParams['verbose']:
214	#print "Couldn't find any match."
215	print_log()
216	asaplog.push( "Couldn't find any match." )
217	print_log( 'ERROR' )
218	else:
219	raise
220
221
222	def get_scan(self, scanid=None):
223	"""
224	Return a specific scan (by scanno) or collection of scans (by
225	source name) in a new scantable.
226	Note:
227	See scantable.drop_scan() for the inverse operation.
228	Parameters:
229	scanid: a (list of) scanno or a source name, unix-style
230	patterns are accepted for source name matching, e.g.
231	'*_R' gets all 'ref scans
232	Example:
233	# get all scans containing the source '323p459'
234	newscan = scan.get_scan('323p459')
235	# get all 'off' scans
236	refscans = scan.get_scan('*_R')
237	# get a susbset of scans by scanno (as listed in scan.summary())
238	newscan = scan.get_scan([0, 2, 7, 10])
239	"""
240	if scanid is None:
241	if rcParams['verbose']:
242	#print "Please specify a scan no or name to " \
243	# "retrieve from the scantable"
244	asaplog.push( 'Please specify a scan no or name to retrieve from the scantable' )
245	print_log( 'ERROR' )
246	return
247	else:
248	raise RuntimeError("No scan given")
249
250	try:
251	bsel = self.get_selection()
252	sel = selector()
253	if type(scanid) is str:
254	sel.set_name(scanid)
255	self.set_selection(bsel+sel)
256	scopy = self._copy()
257	self.set_selection(bsel)
258	return scantable(scopy)
259	elif type(scanid) is int:
260	sel.set_scans([scanid])
261	self.set_selection(bsel+sel)
262	scopy = self._copy()
263	self.set_selection(bsel)
264	return scantable(scopy)
265	elif type(scanid) is list:
266	sel.set_scans(scanid)
267	self.set_selection(sel)
268	scopy = self._copy()
269	self.set_selection(bsel)
270	return scantable(scopy)
271	else:
272	msg = "Illegal scanid type, use 'int' or 'list' if ints."
273	if rcParams['verbose']:
274	#print msg
275	asaplog.push( msg )
276	print_log( 'ERROR' )
277	else:
278	raise TypeError(msg)
279	except RuntimeError:
280	if rcParams['verbose']:
281	#print "Couldn't find any match."
282	print_log()
283	asaplog.push( "Couldn't find any match." )
284	print_log( 'ERROR' )
285	else: raise
286
287	def __str__(self):
288	return Scantable._summary(self, True)
289
290	def summary(self, filename=None):
291	"""
292	Print a summary of the contents of this scantable.
293	Parameters:
294	filename: the name of a file to write the putput to
295	Default - no file output
296	verbose: print extra info such as the frequency table
297	The default (False) is taken from .asaprc
298	"""
299	info = Scantable._summary(self, True)
300	#if verbose is None: verbose = rcParams['scantable.verbosesummary']
301	if filename is not None:
302	if filename is "":
303	filename = 'scantable_summary.txt'
304	from os.path import expandvars, isdir
305	filename = expandvars(filename)
306	if not isdir(filename):
307	data = open(filename, 'w')
308	data.write(info)
309	data.close()
310	else:
311	msg = "Illegal file name '%s'." % (filename)
312	if rcParams['verbose']:
313	#print msg
314	asaplog.push( msg )
315	print_log( 'ERROR' )
316	else:
317	raise IOError(msg)
318	if rcParams['verbose']:
319	try:
320	from IPython.genutils import page as pager
321	except ImportError:
322	from pydoc import pager
323	pager(info)
324	else:
325	return info
326
327	def get_spectrum(self, rowno):
328	"""Return the spectrum for the current row in the scantable as a list.
329	Parameters:
330	rowno: the row number to retrieve the spectrum from
331	"""
332	return self._getspectrum(rowno)
333
334	def get_mask(self, rowno):
335	"""Return the mask for the current row in the scantable as a list.
336	Parameters:
337	rowno: the row number to retrieve the mask from
338	"""
339	return self._getmask(rowno)
340
341	def set_spectrum(self, spec, rowno):
342	"""Return the spectrum for the current row in the scantable as a list.
343	Parameters:
344	spec: the spectrum
345	rowno: the row number to set the spectrum for
346	"""
347	assert(len(spec) == self.nchan())
348	return self._setspectrum(spec, rowno)
349
350	def get_selection(self):
351	"""
352	Get the selection object currently set on this scantable.
353	Parameters:
354	none
355	Example:
356	sel = scan.get_selection()
357	sel.set_ifs(0) # select IF 0
358	scan.set_selection(sel) # apply modified selection
359	"""
360	return selector(self._getselection())
361
362	def set_selection(self, selection=selector()):
363	"""
364	Select a subset of the data. All following operations on this scantable
365	are only applied to thi selection.
366	Parameters:
367	selection: a selector object (default unset the selection)
368	Examples:
369	sel = selector() # create a selection object
370	self.set_scans([0, 3]) # select SCANNO 0 and 3
371	scan.set_selection(sel) # set the selection
372	scan.summary() # will only print summary of scanno 0 an 3
373	scan.set_selection() # unset the selection
374	"""
375	self._setselection(selection)
376
377	def get_row(self, row=0, insitu=None):
378	"""
379	Select a row in the scantable.
380	Return a scantable with single row.
381	Parameters:
382	row: row no of integration, default is 0.
383	insitu: if False a new scantable is returned.
384	Otherwise, the scaling is done in-situ
385	The default is taken from .asaprc (False)
386	"""
387	if insitu is None: insitu = rcParams['insitu']
388	if not insitu:
389	workscan = self.copy()
390	else:
391	workscan = self
392	# Select a row
393	sel=selector()
394	sel.set_scans([workscan.getscan(row)])
395	sel.set_cycles([workscan.getcycle(row)])
396	sel.set_beams([workscan.getbeam(row)])
397	sel.set_ifs([workscan.getif(row)])
398	sel.set_polarisations([workscan.getpol(row)])
399	sel.set_name(workscan._getsourcename(row))
400	workscan.set_selection(sel)
401	if not workscan.nrow() == 1:
402	msg = "Cloud not identify single row. %d rows selected."%(workscan.nrow())
403	raise RuntimeError(msg)
404	del sel
405	if insitu:
406	self._assign(workscan)
407	else:
408	return workscan
409
410	#def stats(self, stat='stddev', mask=None):
411	def stats(self, stat='stddev', mask=None, form='3.3f'):
412	"""
413	Determine the specified statistic of the current beam/if/pol
414	Takes a 'mask' as an optional parameter to specify which
415	channels should be excluded.
416	Parameters:
417	stat: 'min', 'max', 'min_abc', 'max_abc', 'sumsq', 'sum',
418	'mean', 'var', 'stddev', 'avdev', 'rms', 'median'
419	mask: an optional mask specifying where the statistic
420	should be determined.
421	form: format string to print statistic values
422	Example:
423	scan.set_unit('channel')
424	msk = scan.create_mask([100, 200], [500, 600])
425	scan.stats(stat='mean', mask=m)
426	"""
427	if mask == None:
428	mask = []
429	axes = ['Beam', 'IF', 'Pol', 'Time']
430	if not self._check_ifs():
431	raise ValueError("Cannot apply mask as the IFs have different "
432	"number of channels. Please use setselection() "
433	"to select individual IFs")
434	rtnabc = False
435	if stat.lower().endswith('_abc'): rtnabc = True
436	getchan = False
437	if stat.lower().startswith('min') or stat.lower().startswith('max'):
438	chan = self._math._minmaxchan(self, mask, stat)
439	getchan = True
440	statvals = []
441	if not rtnabc: statvals = self._math._stats(self, mask, stat)
442
443	out = ''
444	axes = []
445	for i in range(self.nrow()):
446	axis = []
447	axis.append(self.getscan(i))
448	axis.append(self.getbeam(i))
449	axis.append(self.getif(i))
450	axis.append(self.getpol(i))
451	axis.append(self.getcycle(i))
452	axes.append(axis)
453	tm = self._gettime(i)
454	src = self._getsourcename(i)
455	refstr = ''
456	statunit= ''
457	if getchan:
458	qx, qy = self.chan2data(rowno=i, chan=chan[i])
459	if rtnabc:
460	statvals.append(qx['value'])
461	#refstr = '(value: %3.3f' % (qy['value'])+' ['+qy['unit']+'])'
462	refstr = ('(value: %'+form) % (qy['value'])+' ['+qy['unit']+'])'
463	statunit= '['+qx['unit']+']'
464	else:
465	#refstr = '(@ %3.3f' % (qx['value'])+' ['+qx['unit']+'])'
466	refstr = ('(@ %'+form) % (qx['value'])+' ['+qx['unit']+'])'
467	#statunit= ' ['+qy['unit']+']'
468	out += 'Scan[%d] (%s) ' % (axis[0], src)
469	out += 'Time[%s]:\n' % (tm)
470	if self.nbeam(-1) > 1: out += ' Beam[%d] ' % (axis[1])
471	if self.nif(-1) > 1: out += ' IF[%d] ' % (axis[2])
472	if self.npol(-1) > 1: out += ' Pol[%d] ' % (axis[3])
473	#out += '= %3.3f ' % (statvals[i]) +refstr+'\n'
474	out += ('= %'+form) % (statvals[i]) +' '+refstr+'\n'
475	out += "--------------------------------------------------\n"
476
477	if rcParams['verbose']:
478	import os
479	if os.environ.has_key( 'USER' ):
480	usr=os.environ['USER']
481	else:
482	import commands
483	usr=commands.getoutput( 'whoami' )
484	tmpfile='/tmp/tmp_'+usr+'_casapy_asap_scantable_stats'
485	f=open(tmpfile,'w')
486	print >> f, "--------------------------------------------------"
487	print >> f, " ", stat, statunit
488	print >> f, "--------------------------------------------------"
489	print >> f, out
490	f.close()
491	f=open(tmpfile,'r')
492	x=f.readlines()
493	f.close()
494	for xx in x:
495	asaplog.push( xx, False )
496	print_log()
497	#else:
498	#retval = { 'axesnames': ['scanno', 'beamno', 'ifno', 'polno', 'cycleno'],
499	# 'axes' : axes,
500	# 'data': statvals}
501	return statvals
502
503	def chan2data(self, rowno=0, chan=0):
504	"""
505	Returns channel/frequency/velocity and spectral value
506	at an arbitrary row and channel in the scantable.
507	Parameters:
508	rowno: a row number in the scantable. Default is the
509	first row, i.e. rowno=0
510	chan: a channel in the scantable. Default is the first
511	channel, i.e. pos=0
512	"""
513	if isinstance(rowno, int) and isinstance(chan, int):
514	qx = {'unit': self.get_unit(),
515	'value': self._getabcissa(rowno)[chan]}
516	qy = {'unit': self.get_fluxunit(),
517	'value': self._getspectrum(rowno)[chan]}
518	return qx, qy
519
520	def stddev(self, mask=None):
521	"""
522	Determine the standard deviation of the current beam/if/pol
523	Takes a 'mask' as an optional parameter to specify which
524	channels should be excluded.
525	Parameters:
526	mask: an optional mask specifying where the standard
527	deviation should be determined.
528
529	Example:
530	scan.set_unit('channel')
531	msk = scan.create_mask([100, 200], [500, 600])
532	scan.stddev(mask=m)
533	"""
534	return self.stats(stat='stddev', mask=mask);
535
536
537	def get_column_names(self):
538	"""
539	Return a list of column names, which can be used for selection.
540	"""
541	return list(Scantable.get_column_names(self))
542
543	def get_tsys(self):
544	"""
545	Return the System temperatures.
546	Returns:
547	a list of Tsys values for the current selection
548	"""
549
550	return self._row_callback(self._gettsys, "Tsys")
551
552	def _row_callback(self, callback, label):
553	axes = []
554	axesnames = ['scanno', 'beamno', 'ifno', 'polno', 'cycleno']
555	out = ""
556	outvec = []
557	for i in range(self.nrow()):
558	axis = []
559	axis.append(self.getscan(i))
560	axis.append(self.getbeam(i))
561	axis.append(self.getif(i))
562	axis.append(self.getpol(i))
563	axis.append(self.getcycle(i))
564	axes.append(axis)
565	tm = self._gettime(i)
566	src = self._getsourcename(i)
567	out += 'Scan[%d] (%s) ' % (axis[0], src)
568	out += 'Time[%s]:\n' % (tm)
569	if self.nbeam(-1) > 1: out += ' Beam[%d] ' % (axis[1])
570	if self.nif(-1) > 1: out += ' IF[%d] ' % (axis[2])
571	if self.npol(-1) > 1: out += ' Pol[%d] ' % (axis[3])
572	outvec.append(callback(i))
573	out += '= %3.3f\n' % (outvec[i])
574	out += "--------------------------------------------------\n"
575	if rcParams['verbose']:
576	asaplog.push("--------------------------------------------------")
577	asaplog.push(" %s" % (label))
578	asaplog.push("--------------------------------------------------")
579	asaplog.push(out)
580	print_log()
581	# disabled because the vector seems more useful
582	#retval = {'axesnames': axesnames, 'axes': axes, 'data': outvec}
583	return outvec
584
585	def _get_column(self, callback, row=-1):
586	"""
587	"""
588	if row == -1:
589	return [callback(i) for i in range(self.nrow())]
590	else:
591	if 0 <= row < self.nrow():
592	return callback(row)
593
594
595	def get_time(self, row=-1, asdatetime=False):
596	"""
597	Get a list of time stamps for the observations.
598	Return a datetime object for each integration time stamp in the scantable.
599	Parameters:
600	row: row no of integration. Default -1 return all rows
601	asdatetime: return values as datetime objects rather than strings
602	Example:
603	none
604	"""
605	from time import strptime
606	from datetime import datetime
607	times = self._get_column(self._gettime, row)
608	if not asdatetime:
609	return times
610	format = "%Y/%m/%d/%H:%M:%S"
611	if isinstance(times, list):
612	return [datetime(*strptime(i, format)[:6]) for i in times]
613	else:
614	return datetime(*strptime(times, format)[:6])
615
616
617	def get_inttime(self, row=-1):
618	"""
619	Get a list of integration times for the observations.
620	Return a time in seconds for each integration in the scantable.
621	Parameters:
622	row: row no of integration. Default -1 return all rows.
623	Example:
624	none
625	"""
626	return self._get_column(self._getinttime, row)
627
628
629	def get_sourcename(self, row=-1):
630	"""
631	Get a list source names for the observations.
632	Return a string for each integration in the scantable.
633	Parameters:
634	row: row no of integration. Default -1 return all rows.
635	Example:
636	none
637	"""
638	return self._get_column(self._getsourcename, row)
639
640	def get_elevation(self, row=-1):
641	"""
642	Get a list of elevations for the observations.
643	Return a float for each integration in the scantable.
644	Parameters:
645	row: row no of integration. Default -1 return all rows.
646	Example:
647	none
648	"""
649	return self._get_column(self._getelevation, row)
650
651	def get_azimuth(self, row=-1):
652	"""
653	Get a list of azimuths for the observations.
654	Return a float for each integration in the scantable.
655	Parameters:
656	row: row no of integration. Default -1 return all rows.
657	Example:
658	none
659	"""
660	return self._get_column(self._getazimuth, row)
661
662	def get_parangle(self, row=-1):
663	"""
664	Get a list of parallactic angles for the observations.
665	Return a float for each integration in the scantable.
666	Parameters:
667	row: row no of integration. Default -1 return all rows.
668	Example:
669	none
670	"""
671	return self._get_column(self._getparangle, row)
672
673	def get_direction(self, row=-1):
674	"""
675	Get a list of Positions on the sky (direction) for the observations.
676	Return a float for each integration in the scantable.
677	Parameters:
678	row: row no of integration. Default -1 return all rows
679	Example:
680	none
681	"""
682	return self._get_column(self._getdirection, row)
683
684	def get_directionval(self, row=-1):
685	"""
686	Get a list of Positions on the sky (direction) for the observations.
687	Return a float for each integration in the scantable.
688	Parameters:
689	row: row no of integration. Default -1 return all rows
690	Example:
691	none
692	"""
693	return self._get_column(self._getdirectionvec, row)
694
695	def set_unit(self, unit='channel'):
696	"""
697	Set the unit for all following operations on this scantable
698	Parameters:
699	unit: optional unit, default is 'channel'
700	one of '*Hz', 'km/s', 'channel', ''
701	"""
702	varlist = vars()
703	if unit in ['', 'pixel', 'channel']:
704	unit = ''
705	inf = list(self._getcoordinfo())
706	inf[0] = unit
707	self._setcoordinfo(inf)
708	self._add_history("set_unit", varlist)
709
710	def set_instrument(self, instr):
711	"""
712	Set the instrument for subsequent processing.
713	Parameters:
714	instr: Select from 'ATPKSMB', 'ATPKSHOH', 'ATMOPRA',
715	'DSS-43' (Tid), 'CEDUNA', and 'HOBART'
716	"""
717	self._setInstrument(instr)
718	self._add_history("set_instument", vars())
719	print_log()
720
721	def set_feedtype(self, feedtype):
722	"""
723	Overwrite the feed type, which might not be set correctly.
724	Parameters:
725	feedtype: 'linear' or 'circular'
726	"""
727	self._setfeedtype(feedtype)
728	self._add_history("set_feedtype", vars())
729	print_log()
730
731	def set_doppler(self, doppler='RADIO'):
732	"""
733	Set the doppler for all following operations on this scantable.
734	Parameters:
735	doppler: One of 'RADIO', 'OPTICAL', 'Z', 'BETA', 'GAMMA'
736	"""
737	varlist = vars()
738	inf = list(self._getcoordinfo())
739	inf[2] = doppler
740	self._setcoordinfo(inf)
741	self._add_history("set_doppler", vars())
742	print_log()
743
744	def set_freqframe(self, frame=None):
745	"""
746	Set the frame type of the Spectral Axis.
747	Parameters:
748	frame: an optional frame type, default 'LSRK'. Valid frames are:
749	'REST', 'TOPO', 'LSRD', 'LSRK', 'BARY',
750	'GEO', 'GALACTO', 'LGROUP', 'CMB'
751	Examples:
752	scan.set_freqframe('BARY')
753	"""
754	if frame is None: frame = rcParams['scantable.freqframe']
755	varlist = vars()
756	valid = ['REST', 'TOPO', 'LSRD', 'LSRK', 'BARY', \
757	'GEO', 'GALACTO', 'LGROUP', 'CMB']
758
759	if frame in valid:
760	inf = list(self._getcoordinfo())
761	inf[1] = frame
762	self._setcoordinfo(inf)
763	self._add_history("set_freqframe", varlist)
764	else:
765	msg = "Please specify a valid freq type. Valid types are:\n", valid
766	if rcParams['verbose']:
767	#print msg
768	asaplog.push( msg )
769	print_log( 'ERROR' )
770	else:
771	raise TypeError(msg)
772	print_log()
773
774	def set_dirframe(self, frame=""):
775	"""
776	Set the frame type of the Direction on the sky.
777	Parameters:
778	frame: an optional frame type, default ''. Valid frames are:
779	'J2000', 'B1950', 'GALACTIC'
780	Examples:
781	scan.set_dirframe('GALACTIC')
782	"""
783	varlist = vars()
784	try:
785	Scantable.set_dirframe(self, frame)
786	except RuntimeError, msg:
787	if rcParams['verbose']:
788	#print msg
789	print_log()
790	asaplog.push( str(msg) )
791	print_log( 'ERROR' )
792	else:
793	raise
794	self._add_history("set_dirframe", varlist)
795
796	def get_unit(self):
797	"""
798	Get the default unit set in this scantable
799	Returns:
800	A unit string
801	"""
802	inf = self._getcoordinfo()
803	unit = inf[0]
804	if unit == '': unit = 'channel'
805	return unit
806
807	def get_abcissa(self, rowno=0):
808	"""
809	Get the abcissa in the current coordinate setup for the currently
810	selected Beam/IF/Pol
811	Parameters:
812	rowno: an optional row number in the scantable. Default is the
813	first row, i.e. rowno=0
814	Returns:
815	The abcissa values and the format string (as a dictionary)
816	"""
817	abc = self._getabcissa(rowno)
818	lbl = self._getabcissalabel(rowno)
819	print_log()
820	return abc, lbl
821
822	def flag(self, mask=None, unflag=False):
823	"""
824	Flag the selected data using an optional channel mask.
825	Parameters:
826	mask: an optional channel mask, created with create_mask. Default
827	(no mask) is all channels.
828	unflag: if True, unflag the data
829	"""
830	varlist = vars()
831	if mask is None:
832	mask = []
833	try:
834	self._flag(mask, unflag)
835	except RuntimeError, msg:
836	if rcParams['verbose']:
837	#print msg
838	print_log()
839	asaplog.push( str(msg) )
840	print_log( 'ERROR' )
841	return
842	else: raise
843	self._add_history("flag", varlist)
844
845	def flag_row(self, rows=[], unflag=False):
846	"""
847	Flag the selected data in row-based manner.
848	Parameters:
849	rows: list of row numbers to be flagged. Default is no row (must be explicitly specified to execute row-based flagging).
850	unflag: if True, unflag the data.
851	"""
852	varlist = vars()
853	try:
854	self._flag_row(rows, unflag)
855	except RuntimeError, msg:
856	if rcParams['verbose']:
857	print_log()
858	asaplog.push( str(msg) )
859	print_log('ERROR')
860	return
861	else: raise
862	self._add_history("flag_row", varlist)
863
864
865	def lag_flag(self, frequency, width=0.0, unit="GHz", insitu=None):
866	"""
867	Flag the data in 'lag' space by providing a frequency to remove.
868	Flagged data in the scantable gets set to 0.0 before the fft.
869	No taper is applied.
870	Parameters:
871	frequency: the frequency (really a period within the bandwidth)
872	to remove
873	width: the width of the frequency to remove, to remove a
874	range of frequencies around the centre.
875	unit: the frequency unit (default "GHz")
876	Notes:
877	It is recommended to flag edges of the band or strong
878	signals beforehand.
879	"""
880	if insitu is None: insitu = rcParams['insitu']
881	self._math._setinsitu(insitu)
882	varlist = vars()
883	base = { "GHz": 1000000000., "MHz": 1000000., "kHz": 1000., "Hz": 1. }
884	if not base.has_key(unit):
885	raise ValueError("%s is not a valid unit." % unit)
886	try:
887	s = scantable(self._math._lag_flag(self, frequency*base[unit],
888	width*base[unit]))
889	except RuntimeError, msg:
890	if rcParams['verbose']:
891	#print msg
892	print_log()
893	asaplog.push( str(msg) )
894	print_log( 'ERROR' )
895	return
896	else: raise
897	s._add_history("lag_flag", varlist)
898	print_log()
899	if insitu:
900	self._assign(s)
901	else:
902	return s
903
904
905	def create_mask(self, args, *kwargs):
906	"""
907	Compute and return a mask based on [min, max] windows.
908	The specified windows are to be INCLUDED, when the mask is
909	applied.
910	Parameters:
911	[min, max], [min2, max2], ...
912	Pairs of start/end points (inclusive)specifying the regions
913	to be masked
914	invert: optional argument. If specified as True,
915	return an inverted mask, i.e. the regions
916	specified are EXCLUDED
917	row: create the mask using the specified row for
918	unit conversions, default is row=0
919	only necessary if frequency varies over rows.
920	Example:
921	scan.set_unit('channel')
922	a)
923	msk = scan.create_mask([400, 500], [800, 900])
924	# masks everything outside 400 and 500
925	# and 800 and 900 in the unit 'channel'
926
927	b)
928	msk = scan.create_mask([400, 500], [800, 900], invert=True)
929	# masks the regions between 400 and 500
930	# and 800 and 900 in the unit 'channel'
931	c)
932	mask only channel 400
933	msk = scan.create_mask([400, 400])
934	"""
935	row = 0
936	if kwargs.has_key("row"):
937	row = kwargs.get("row")
938	data = self._getabcissa(row)
939	u = self._getcoordinfo()[0]
940	if rcParams['verbose']:
941	if u == "": u = "channel"
942	msg = "The current mask window unit is %s" % u
943	i = self._check_ifs()
944	if not i:
945	msg += "\nThis mask is only valid for IF=%d" % (self.getif(i))
946	asaplog.push(msg)
947	n = self.nchan()
948	msk = _n_bools(n, False)
949	# test if args is a 'list' or a 'normal *args - UGLY!!!
950
951	ws = (isinstance(args[-1][-1], int) or isinstance(args[-1][-1], float)) \
952	and args or args[0]
953	for window in ws:
954	if (len(window) != 2 or window[0] > window[1] ):
955	raise TypeError("A window needs to be defined as [min, max]")
956	for i in range(n):
957	if data[i] >= window[0] and data[i] <= window[1]:
958	msk[i] = True
959	if kwargs.has_key('invert'):
960	if kwargs.get('invert'):
961	msk = mask_not(msk)
962	print_log()
963	return msk
964
965	def get_masklist(self, mask=None, row=0):
966	"""
967	Compute and return a list of mask windows, [min, max].
968	Parameters:
969	mask: channel mask, created with create_mask.
970	row: calcutate the masklist using the specified row
971	for unit conversions, default is row=0
972	only necessary if frequency varies over rows.
973	Returns:
974	[min, max], [min2, max2], ...
975	Pairs of start/end points (inclusive)specifying
976	the masked regions
977	"""
978	if not (isinstance(mask,list) or isinstance(mask, tuple)):
979	raise TypeError("The mask should be list or tuple.")
980	if len(mask) < 2:
981	raise TypeError("The mask elements should be > 1")
982	if self.nchan() != len(mask):
983	msg = "Number of channels in scantable != number of mask elements"
984	raise TypeError(msg)
985	data = self._getabcissa(row)
986	u = self._getcoordinfo()[0]
987	if rcParams['verbose']:
988	if u == "": u = "channel"
989	msg = "The current mask window unit is %s" % u
990	i = self._check_ifs()
991	if not i:
992	msg += "\nThis mask is only valid for IF=%d" % (self.getif(i))
993	asaplog.push(msg)
994	masklist=[]
995	ist, ien = None, None
996	ist, ien=self.get_mask_indices(mask)
997	if ist is not None and ien is not None:
998	for i in xrange(len(ist)):
999	range=[data[ist[i]],data[ien[i]]]
1000	range.sort()
1001	masklist.append([range[0],range[1]])
1002	return masklist
1003
1004	def get_mask_indices(self, mask=None):
1005	"""
1006	Compute and Return lists of mask start indices and mask end indices.
1007	Parameters:
1008	mask: channel mask, created with create_mask.
1009	Returns:
1010	List of mask start indices and that of mask end indices,
1011	i.e., [istart1,istart2,....], [iend1,iend2,....].
1012	"""
1013	if not (isinstance(mask,list) or isinstance(mask, tuple)):
1014	raise TypeError("The mask should be list or tuple.")
1015	if len(mask) < 2:
1016	raise TypeError("The mask elements should be > 1")
1017	istart=[]
1018	iend=[]
1019	if mask[0]: istart.append(0)
1020	for i in range(len(mask)-1):
1021	if not mask[i] and mask[i+1]:
1022	istart.append(i+1)
1023	elif mask[i] and not mask[i+1]:
1024	iend.append(i)
1025	if mask[len(mask)-1]: iend.append(len(mask)-1)
1026	if len(istart) != len(iend):
1027	raise RuntimeError("Numbers of mask start != mask end.")
1028	for i in range(len(istart)):
1029	if istart[i] > iend[i]:
1030	raise RuntimeError("Mask start index > mask end index")
1031	break
1032	return istart,iend
1033
1034	# def get_restfreqs(self):
1035	# """
1036	# Get the restfrequency(s) stored in this scantable.
1037	# The return value(s) are always of unit 'Hz'
1038	# Parameters:
1039	# none
1040	# Returns:
1041	# a list of doubles
1042	# """
1043	# return list(self._getrestfreqs())
1044
1045	def get_restfreqs(self, ids=None):
1046	"""
1047	Get the restfrequency(s) stored in this scantable.
1048	The return value(s) are always of unit 'Hz'
1049	Parameters:
1050	ids: (optional) a list of MOLECULE_ID for that restfrequency(s) to
1051	be retrieved
1052	Returns:
1053	dictionary containing ids and a list of doubles for each id
1054	"""
1055	if ids is None:
1056	rfreqs={}
1057	idlist = self.getmolnos()
1058	for i in idlist:
1059	rfreqs[i]=list(self._getrestfreqs(i))
1060	return rfreqs
1061	else:
1062	if type(ids)==list or type(ids)==tuple:
1063	rfreqs={}
1064	for i in ids:
1065	rfreqs[i]=list(self._getrestfreqs(i))
1066	return rfreqs
1067	else:
1068	return list(self._getrestfreqs(ids))
1069	#return list(self._getrestfreqs(ids))
1070
1071	def set_restfreqs(self, freqs=None, unit='Hz'):
1072	"""
1073	******NEED TO BE UPDATED begin**********
1074	Set or replace the restfrequency specified and
1075	If the 'freqs' argument holds a scalar,
1076	then that rest frequency will be applied to all the selected
1077	data. If the 'freqs' argument holds
1078	a vector, then it MUST be of equal or smaller length than
1079	the number of IFs (and the available restfrequencies will be
1080	replaced by this vector). In this case, all data have
1081	the restfrequency set per IF according
1082	to the corresponding value you give in the 'freqs' vector.
1083	E.g. 'freqs=[1e9, 2e9]' would mean IF 0 gets restfreq 1e9 and
1084	IF 1 gets restfreq 2e9.
1085	******NEED TO BE UPDATED end**********
1086	You can also specify the frequencies via a linecatalog.
1087
1088	Parameters:
1089	freqs: list of rest frequency values or string idenitfiers
1090	unit: unit for rest frequency (default 'Hz')
1091
1092	Example:
1093	# set the given restfrequency for the all currently selected IFs
1094	scan.set_restfreqs(freqs=1.4e9)
1095	# set multiple restfrequencies to all the selected data
1096	scan.set_restfreqs(freqs=[1.4e9, 1.41e9, 1.42e9])
1097	# If the number of IFs in the data is >= 2 the IF0 gets the first
1098	# value IF1 the second... NOTE that freqs needs to be
1099	# specified in list of list (e.g. [[],[],...] ).
1100	scan.set_restfreqs(freqs=[[1.4e9],[1.67e9]])
1101	#set the given restfrequency for the whole table (by name)
1102	scan.set_restfreqs(freqs="OH1667")
1103
1104	Note:
1105	To do more sophisticate Restfrequency setting, e.g. on a
1106	source and IF basis, use scantable.set_selection() before using
1107	this function.
1108	# provide your scantable is call scan
1109	selection = selector()
1110	selection.set_name("ORION*")
1111	selection.set_ifs([1])
1112	scan.set_selection(selection)
1113	scan.set_restfreqs(freqs=86.6e9)
1114
1115	"""
1116	varlist = vars()
1117	from asap import linecatalog
1118	# simple value
1119	if isinstance(freqs, int) or isinstance(freqs, float):
1120	# TT mod
1121	#self._setrestfreqs(freqs, "",unit)
1122	self._setrestfreqs([freqs], [""],unit)
1123	# list of values
1124	elif isinstance(freqs, list) or isinstance(freqs, tuple):
1125	# list values are scalars
1126	if isinstance(freqs[-1], int) or isinstance(freqs[-1], float):
1127	self._setrestfreqs(freqs, [""],unit)
1128	# list values are tuples, (value, name)
1129	elif isinstance(freqs[-1], dict):
1130	#sel = selector()
1131	#savesel = self._getselection()
1132	#iflist = self.getifnos()
1133	#for i in xrange(len(freqs)):
1134	# sel.set_ifs(iflist[i])
1135	# self._setselection(sel)
1136	# self._setrestfreqs(freqs[i], "",unit)
1137	#self._setselection(savesel)
1138	self._setrestfreqs(freqs["value"],
1139	freqs["name"], "MHz")
1140	elif isinstance(freqs[-1], list) or isinstance(freqs[-1], tuple):
1141	sel = selector()
1142	savesel = self._getselection()
1143	iflist = self.getifnos()
1144	if len(freqs)>len(iflist):
1145	raise ValueError("number of elements in list of list exeeds the current IF selections")
1146	for i in xrange(len(freqs)):
1147	sel.set_ifs(iflist[i])
1148	self._setselection(sel)
1149	self._setrestfreqs(freqs[i]["value"],
1150	freqs[i]["name"], "MHz")
1151	self._setselection(savesel)
1152	# freqs are to be taken from a linecatalog
1153	elif isinstance(freqs, linecatalog):
1154	sel = selector()
1155	savesel = self._getselection()
1156	for i in xrange(freqs.nrow()):
1157	sel.set_ifs(iflist[i])
1158	self._setselection(sel)
1159	self._setrestfreqs(freqs.get_frequency(i),
1160	freqs.get_name(i), "MHz")
1161	# ensure that we are not iterating past nIF
1162	if i == self.nif()-1: break
1163	self._setselection(savesel)
1164	else:
1165	return
1166	self._add_history("set_restfreqs", varlist)
1167
1168	def shift_refpix(self, delta):
1169	"""
1170	Shift the reference pixel of the Spectra Coordinate by an
1171	integer amount.
1172	Parameters:
1173	delta: the amount to shift by
1174	Note:
1175	Be careful using this with broadband data.
1176	"""
1177	Scantable.shift(self, delta)
1178
1179	def history(self, filename=None):
1180	"""
1181	Print the history. Optionally to a file.
1182	Parameters:
1183	filename: The name of the file to save the history to.
1184	"""
1185	hist = list(self._gethistory())
1186	out = "-"*80
1187	for h in hist:
1188	if h.startswith("---"):
1189	out += "\n"+h
1190	else:
1191	items = h.split("##")
1192	date = items[0]
1193	func = items[1]
1194	items = items[2:]
1195	out += "\n"+date+"\n"
1196	out += "Function: %s\n Parameters:" % (func)
1197	for i in items:
1198	s = i.split("=")
1199	out += "\n %s = %s" % (s[0], s[1])
1200	out += "\n"+"-"*80
1201	if filename is not None:
1202	if filename is "":
1203	filename = 'scantable_history.txt'
1204	import os
1205	filename = os.path.expandvars(os.path.expanduser(filename))
1206	if not os.path.isdir(filename):
1207	data = open(filename, 'w')
1208	data.write(out)
1209	data.close()
1210	else:
1211	msg = "Illegal file name '%s'." % (filename)
1212	if rcParams['verbose']:
1213	#print msg
1214	asaplog.push( msg )
1215	print_log( 'ERROR' )
1216	else:
1217	raise IOError(msg)
1218	if rcParams['verbose']:
1219	try:
1220	from IPython.genutils import page as pager
1221	except ImportError:
1222	from pydoc import pager
1223	pager(out)
1224	else:
1225	return out
1226	return
1227	#
1228	# Maths business
1229	#
1230
1231	def average_time(self, mask=None, scanav=False, weight='tint', align=False):
1232	"""
1233	Return the (time) weighted average of a scan.
1234	Note:
1235	in channels only - align if necessary
1236	Parameters:
1237	mask: an optional mask (only used for 'var' and 'tsys'
1238	weighting)
1239	scanav: True averages each scan separately
1240	False (default) averages all scans together,
1241	weight: Weighting scheme.
1242	'none' (mean no weight)
1243	'var' (1/var(spec) weighted)
1244	'tsys' (1/Tsys**2 weighted)
1245	'tint' (integration time weighted)
1246	'tintsys' (Tint/Tsys**2)
1247	'median' ( median averaging)
1248	The default is 'tint'
1249	align: align the spectra in velocity before averaging. It takes
1250	the time of the first spectrum as reference time.
1251	Example:
1252	# time average the scantable without using a mask
1253	newscan = scan.average_time()
1254	"""
1255	varlist = vars()
1256	if weight is None: weight = 'TINT'
1257	if mask is None: mask = ()
1258	if scanav: scanav = "SCAN"
1259	else: scanav = "NONE"
1260	scan = (self, )
1261	try:
1262	if align:
1263	scan = (self.freq_align(insitu=False), )
1264	s = None
1265	if weight.upper() == 'MEDIAN':
1266	s = scantable(self._math._averagechannel(scan[0], 'MEDIAN',
1267	scanav))
1268	else:
1269	s = scantable(self._math._average(scan, mask, weight.upper(),
1270	scanav))
1271	except RuntimeError, msg:
1272	if rcParams['verbose']:
1273	#print msg
1274	print_log()
1275	asaplog.push( str(msg) )
1276	print_log( 'ERROR' )
1277	return
1278	else: raise
1279	s._add_history("average_time", varlist)
1280	print_log()
1281	return s
1282
1283	def convert_flux(self, jyperk=None, eta=None, d=None, insitu=None):
1284	"""
1285	Return a scan where all spectra are converted to either
1286	Jansky or Kelvin depending upon the flux units of the scan table.
1287	By default the function tries to look the values up internally.
1288	If it can't find them (or if you want to over-ride), you must
1289	specify EITHER jyperk OR eta (and D which it will try to look up
1290	also if you don't set it). jyperk takes precedence if you set both.
1291	Parameters:
1292	jyperk: the Jy / K conversion factor
1293	eta: the aperture efficiency
1294	d: the geomtric diameter (metres)
1295	insitu: if False a new scantable is returned.
1296	Otherwise, the scaling is done in-situ
1297	The default is taken from .asaprc (False)
1298	"""
1299	if insitu is None: insitu = rcParams['insitu']
1300	self._math._setinsitu(insitu)
1301	varlist = vars()
1302	if jyperk is None: jyperk = -1.0
1303	if d is None: d = -1.0
1304	if eta is None: eta = -1.0
1305	s = scantable(self._math._convertflux(self, d, eta, jyperk))
1306	s._add_history("convert_flux", varlist)
1307	print_log()
1308	if insitu: self._assign(s)
1309	else: return s
1310
1311	def gain_el(self, poly=None, filename="", method="linear", insitu=None):
1312	"""
1313	Return a scan after applying a gain-elevation correction.
1314	The correction can be made via either a polynomial or a
1315	table-based interpolation (and extrapolation if necessary).
1316	You specify polynomial coefficients, an ascii table or neither.
1317	If you specify neither, then a polynomial correction will be made
1318	with built in coefficients known for certain telescopes (an error
1319	will occur if the instrument is not known).
1320	The data and Tsys are divided by the scaling factors.
1321	Parameters:
1322	poly: Polynomial coefficients (default None) to compute a
1323	gain-elevation correction as a function of
1324	elevation (in degrees).
1325	filename: The name of an ascii file holding correction factors.
1326	The first row of the ascii file must give the column
1327	names and these MUST include columns
1328	"ELEVATION" (degrees) and "FACTOR" (multiply data
1329	by this) somewhere.
1330	The second row must give the data type of the
1331	column. Use 'R' for Real and 'I' for Integer.
1332	An example file would be
1333	(actual factors are arbitrary) :
1334
1335	TIME ELEVATION FACTOR
1336	R R R
1337	0.1 0 0.8
1338	0.2 20 0.85
1339	0.3 40 0.9
1340	0.4 60 0.85
1341	0.5 80 0.8
1342	0.6 90 0.75
1343	method: Interpolation method when correcting from a table.
1344	Values are "nearest", "linear" (default), "cubic"
1345	and "spline"
1346	insitu: if False a new scantable is returned.
1347	Otherwise, the scaling is done in-situ
1348	The default is taken from .asaprc (False)
1349	"""
1350
1351	if insitu is None: insitu = rcParams['insitu']
1352	self._math._setinsitu(insitu)
1353	varlist = vars()
1354	if poly is None:
1355	poly = ()
1356	from os.path import expandvars
1357	filename = expandvars(filename)
1358	s = scantable(self._math._gainel(self, poly, filename, method))
1359	s._add_history("gain_el", varlist)
1360	print_log()
1361	if insitu: self._assign(s)
1362	else: return s
1363
1364	def freq_align(self, reftime=None, method='cubic', insitu=None):
1365	"""
1366	Return a scan where all rows have been aligned in frequency/velocity.
1367	The alignment frequency frame (e.g. LSRK) is that set by function
1368	set_freqframe.
1369	Parameters:
1370	reftime: reference time to align at. By default, the time of
1371	the first row of data is used.
1372	method: Interpolation method for regridding the spectra.
1373	Choose from "nearest", "linear", "cubic" (default)
1374	and "spline"
1375	insitu: if False a new scantable is returned.
1376	Otherwise, the scaling is done in-situ
1377	The default is taken from .asaprc (False)
1378	"""
1379	if insitu is None: insitu = rcParams["insitu"]
1380	self._math._setinsitu(insitu)
1381	varlist = vars()
1382	if reftime is None: reftime = ""
1383	s = scantable(self._math._freq_align(self, reftime, method))
1384	s._add_history("freq_align", varlist)
1385	print_log()
1386	if insitu: self._assign(s)
1387	else: return s
1388
1389	def opacity(self, tau, insitu=None):
1390	"""
1391	Apply an opacity correction. The data
1392	and Tsys are multiplied by the correction factor.
1393	Parameters:
1394	tau: Opacity from which the correction factor is
1395	exp(tau*ZD)
1396	where ZD is the zenith-distance
1397	insitu: if False a new scantable is returned.
1398	Otherwise, the scaling is done in-situ
1399	The default is taken from .asaprc (False)
1400	"""
1401	if insitu is None: insitu = rcParams['insitu']
1402	self._math._setinsitu(insitu)
1403	varlist = vars()
1404	s = scantable(self._math._opacity(self, tau))
1405	s._add_history("opacity", varlist)
1406	print_log()
1407	if insitu: self._assign(s)
1408	else: return s
1409
1410	def bin(self, width=5, insitu=None):
1411	"""
1412	Return a scan where all spectra have been binned up.
1413	Parameters:
1414	width: The bin width (default=5) in pixels
1415	insitu: if False a new scantable is returned.
1416	Otherwise, the scaling is done in-situ
1417	The default is taken from .asaprc (False)
1418	"""
1419	if insitu is None: insitu = rcParams['insitu']
1420	self._math._setinsitu(insitu)
1421	varlist = vars()
1422	s = scantable(self._math._bin(self, width))
1423	s._add_history("bin", varlist)
1424	print_log()
1425	if insitu: self._assign(s)
1426	else: return s
1427
1428
1429	def resample(self, width=5, method='cubic', insitu=None):
1430	"""
1431	Return a scan where all spectra have been binned up.
1432
1433	Parameters:
1434	width: The bin width (default=5) in pixels
1435	method: Interpolation method when correcting from a table.
1436	Values are "nearest", "linear", "cubic" (default)
1437	and "spline"
1438	insitu: if False a new scantable is returned.
1439	Otherwise, the scaling is done in-situ
1440	The default is taken from .asaprc (False)
1441	"""
1442	if insitu is None: insitu = rcParams['insitu']
1443	self._math._setinsitu(insitu)
1444	varlist = vars()
1445	s = scantable(self._math._resample(self, method, width))
1446	s._add_history("resample", varlist)
1447	print_log()
1448	if insitu: self._assign(s)
1449	else: return s
1450
1451
1452	def average_pol(self, mask=None, weight='none'):
1453	"""
1454	Average the Polarisations together.
1455	Parameters:
1456	mask: An optional mask defining the region, where the
1457	averaging will be applied. The output will have all
1458	specified points masked.
1459	weight: Weighting scheme. 'none' (default), 'var' (1/var(spec)
1460	weighted), or 'tsys' (1/Tsys**2 weighted)
1461	"""
1462	varlist = vars()
1463	if mask is None:
1464	mask = ()
1465	s = scantable(self._math._averagepol(self, mask, weight.upper()))
1466	s._add_history("average_pol", varlist)
1467	print_log()
1468	return s
1469
1470	def average_beam(self, mask=None, weight='none'):
1471	"""
1472	Average the Beams together.
1473	Parameters:
1474	mask: An optional mask defining the region, where the
1475	averaging will be applied. The output will have all
1476	specified points masked.
1477	weight: Weighting scheme. 'none' (default), 'var' (1/var(spec)
1478	weighted), or 'tsys' (1/Tsys**2 weighted)
1479	"""
1480	varlist = vars()
1481	if mask is None:
1482	mask = ()
1483	s = scantable(self._math._averagebeams(self, mask, weight.upper()))
1484	s._add_history("average_beam", varlist)
1485	print_log()
1486	return s
1487
1488	def convert_pol(self, poltype=None):
1489	"""
1490	Convert the data to a different polarisation type.
1491	Parameters:
1492	poltype: The new polarisation type. Valid types are:
1493	"linear", "stokes" and "circular"
1494	"""
1495	varlist = vars()
1496	try:
1497	s = scantable(self._math._convertpol(self, poltype))
1498	except RuntimeError, msg:
1499	if rcParams['verbose']:
1500	#print msg
1501	print_log()
1502	asaplog.push( str(msg) )
1503	print_log( 'ERROR' )
1504	return
1505	else:
1506	raise
1507	s._add_history("convert_pol", varlist)
1508	print_log()
1509	return s
1510
1511	#def smooth(self, kernel="hanning", width=5.0, insitu=None):
1512	def smooth(self, kernel="hanning", width=5.0, plot=False, insitu=None):
1513	"""
1514	Smooth the spectrum by the specified kernel (conserving flux).
1515	Parameters:
1516	kernel: The type of smoothing kernel. Select from
1517	'hanning' (default), 'gaussian', 'boxcar' and
1518	'rmedian'
1519	width: The width of the kernel in pixels. For hanning this is
1520	ignored otherwise it defauls to 5 pixels.
1521	For 'gaussian' it is the Full Width Half
1522	Maximum. For 'boxcar' it is the full width.
1523	For 'rmedian' it is the half width.
1524	plot: plot the original and the smoothed spectra.
1525	In this each indivual fit has to be approved, by
1526	typing 'y' or 'n'
1527	insitu: if False a new scantable is returned.
1528	Otherwise, the scaling is done in-situ
1529	The default is taken from .asaprc (False)
1530	Example:
1531	none
1532	"""
1533	if insitu is None: insitu = rcParams['insitu']
1534	self._math._setinsitu(insitu)
1535	varlist = vars()
1536
1537	if plot: orgscan = self.copy()
1538
1539	s = scantable(self._math._smooth(self, kernel.lower(), width))
1540	s._add_history("smooth", varlist)
1541
1542	if plot:
1543	if rcParams['plotter.gui']:
1544	from asap.asaplotgui import asaplotgui as asaplot
1545	else:
1546	from asap.asaplot import asaplot
1547	self._p=asaplot()
1548	self._p.set_panels()
1549	ylab=s._get_ordinate_label()
1550	#self._p.palette(0,["#777777","red"])
1551	for r in xrange(s.nrow()):
1552	xsm=s._getabcissa(r)
1553	ysm=s._getspectrum(r)
1554	xorg=orgscan._getabcissa(r)
1555	yorg=orgscan._getspectrum(r)
1556	self._p.clear()
1557	self._p.hold()
1558	self._p.set_axes('ylabel',ylab)
1559	self._p.set_axes('xlabel',s._getabcissalabel(r))
1560	self._p.set_axes('title',s._getsourcename(r))
1561	self._p.set_line(label='Original',color="#777777")
1562	self._p.plot(xorg,yorg)
1563	self._p.set_line(label='Smoothed',color="red")
1564	self._p.plot(xsm,ysm)
1565	### Ugly part for legend
1566	for i in [0,1]:
1567	self._p.subplots[0]['lines'].append([self._p.subplots[0]['axes'].lines[i]])
1568	self._p.release()
1569	### Ugly part for legend
1570	self._p.subplots[0]['lines']=[]
1571	res = raw_input("Accept smoothing ([y]/n): ")
1572	if res.upper() == 'N':
1573	s._setspectrum(yorg, r)
1574	self._p.unmap()
1575	self._p = None
1576	del orgscan
1577
1578	print_log()
1579	if insitu: self._assign(s)
1580	else: return s
1581
1582
1583	def poly_baseline(self, mask=None, order=0, plot=False, uselin=False, insitu=None):
1584	"""
1585	Return a scan which has been baselined (all rows) by a polynomial.
1586	Parameters:
1587	mask: an optional mask
1588	order: the order of the polynomial (default is 0)
1589	plot: plot the fit and the residual. In this each
1590	indivual fit has to be approved, by typing 'y'
1591	or 'n'
1592	uselin: use linear polynomial fit
1593	insitu: if False a new scantable is returned.
1594	Otherwise, the scaling is done in-situ
1595	The default is taken from .asaprc (False)
1596	Example:
1597	# return a scan baselined by a third order polynomial,
1598	# not using a mask
1599	bscan = scan.poly_baseline(order=3)
1600	"""
1601	if insitu is None: insitu = rcParams['insitu']
1602	if not insitu:
1603	workscan = self.copy()
1604	else:
1605	workscan = self
1606	varlist = vars()
1607	if mask is None:
1608	mask = [True for i in xrange(self.nchan(-1))]
1609
1610	from asap.asapfitter import fitter
1611	try:
1612	f = fitter()
1613	if uselin:
1614	f.set_function(lpoly=order)
1615	else:
1616	f.set_function(poly=order)
1617
1618	rows = range(workscan.nrow())
1619	if len(rows) > 0:
1620	self.blpars = []
1621
1622	for r in rows:
1623	# take into account flagtra info (CAS-1434)
1624	flagtra = workscan._getmask(r)
1625	actualmask = mask[:]
1626	if len(actualmask) == 0:
1627	actualmask = list(flagtra[:])
1628	else:
1629	if len(actualmask) != len(flagtra):
1630	raise RuntimeError, "Mask and flagtra have different length"
1631	else:
1632	for i in range(0, len(actualmask)):
1633	actualmask[i] = actualmask[i] and flagtra[i]
1634	f.set_scan(workscan, actualmask)
1635	f.x = workscan._getabcissa(r)
1636	f.y = workscan._getspectrum(r)
1637	f.data = None
1638	f.fit()
1639	if plot:
1640	f.plot(residual=True)
1641	x = raw_input("Accept fit ( [y]/n ): ")
1642	if x.upper() == 'N':
1643	self.blpars.append(None)
1644	continue
1645	workscan._setspectrum(f.fitter.getresidual(), r)
1646	self.blpars.append(f.get_parameters())
1647
1648	if plot:
1649	f._p.unmap()
1650	f._p = None
1651	workscan._add_history("poly_baseline", varlist)
1652	print_log()
1653	if insitu: self._assign(workscan)
1654	else: return workscan
1655	except RuntimeError:
1656	msg = "The fit failed, possibly because it didn't converge."
1657	if rcParams['verbose']:
1658	#print msg
1659	print_log()
1660	asaplog.push( str(msg) )
1661	print_log( 'ERROR' )
1662	return
1663	else:
1664	raise RuntimeError(msg)
1665
1666
1667	def auto_poly_baseline(self, mask=[], edge=(0, 0), order=0,
1668	threshold=3, chan_avg_limit=1, plot=False,
1669	insitu=None):
1670	"""
1671	Return a scan which has been baselined (all rows) by a polynomial.
1672	Spectral lines are detected first using linefinder and masked out
1673	to avoid them affecting the baseline solution.
1674
1675	Parameters:
1676	mask: an optional mask retreived from scantable
1677	edge: an optional number of channel to drop at
1678	the edge of spectrum. If only one value is
1679	specified, the same number will be dropped from
1680	both sides of the spectrum. Default is to keep
1681	all channels. Nested tuples represent individual
1682	edge selection for different IFs (a number of spectral
1683	channels can be different)
1684	order: the order of the polynomial (default is 0)
1685	threshold: the threshold used by line finder. It is better to
1686	keep it large as only strong lines affect the
1687	baseline solution.
1688	chan_avg_limit:
1689	a maximum number of consequtive spectral channels to
1690	average during the search of weak and broad lines.
1691	The default is no averaging (and no search for weak
1692	lines). If such lines can affect the fitted baseline
1693	(e.g. a high order polynomial is fitted), increase this
1694	parameter (usually values up to 8 are reasonable). Most
1695	users of this method should find the default value
1696	sufficient.
1697	plot: plot the fit and the residual. In this each
1698	indivual fit has to be approved, by typing 'y'
1699	or 'n'
1700	insitu: if False a new scantable is returned.
1701	Otherwise, the scaling is done in-situ
1702	The default is taken from .asaprc (False)
1703
1704	Example:
1705	scan2=scan.auto_poly_baseline(order=7)
1706	"""
1707	if insitu is None: insitu = rcParams['insitu']
1708	varlist = vars()
1709	from asap.asapfitter import fitter
1710	from asap.asaplinefind import linefinder
1711	from asap import _is_sequence_or_number as _is_valid
1712
1713	# check whether edge is set up for each IF individually
1714	individualedge = False;
1715	if len(edge) > 1:
1716	if isinstance(edge[0], list) or isinstance(edge[0], tuple):
1717	individualedge = True;
1718
1719	if not _is_valid(edge, int) and not individualedge:
1720	raise ValueError, "Parameter 'edge' has to be an integer or a \
1721	pair of integers specified as a tuple. Nested tuples are allowed \
1722	to make individual selection for different IFs."
1723
1724	curedge = (0, 0)
1725	if individualedge:
1726	for edgepar in edge:
1727	if not _is_valid(edgepar, int):
1728	raise ValueError, "Each element of the 'edge' tuple has \
1729	to be a pair of integers or an integer."
1730	else:
1731	curedge = edge;
1732
1733	# setup fitter
1734	f = fitter()
1735	f.set_function(poly=order)
1736
1737	# setup line finder
1738	fl = linefinder()
1739	fl.set_options(threshold=threshold,avg_limit=chan_avg_limit)
1740
1741	if not insitu:
1742	workscan = self.copy()
1743	else:
1744	workscan = self
1745
1746	fl.set_scan(workscan)
1747
1748	rows = range(workscan.nrow())
1749	# Save parameters of baseline fits & masklists as a class attribute.
1750	# NOTICE: It does not reflect changes in scantable!
1751	if len(rows) > 0:
1752	self.blpars=[]
1753	self.masklists=[]
1754	asaplog.push("Processing:")
1755	for r in rows:
1756	msg = " Scan[%d] Beam[%d] IF[%d] Pol[%d] Cycle[%d]" % \
1757	(workscan.getscan(r), workscan.getbeam(r), workscan.getif(r), \
1758	workscan.getpol(r), workscan.getcycle(r))
1759	asaplog.push(msg, False)
1760
1761	# figure out edge parameter
1762	if individualedge:
1763	if len(edge) >= workscan.getif(r):
1764	raise RuntimeError, "Number of edge elements appear to " \
1765	"be less than the number of IFs"
1766	curedge = edge[workscan.getif(r)]
1767
1768	# take into account flagtra info (CAS-1434)
1769	flagtra = workscan._getmask(r)
1770	actualmask = mask[:]
1771	if len(actualmask) == 0:
1772	actualmask = list(flagtra[:])
1773	else:
1774	if len(actualmask) != len(flagtra):
1775	raise RuntimeError, "Mask and flagtra have different length"
1776	else:
1777	for i in range(0, len(actualmask)):
1778	actualmask[i] = actualmask[i] and flagtra[i]
1779
1780	# setup line finder
1781	fl.find_lines(r, actualmask, curedge)
1782	outmask=fl.get_mask()
1783	f.set_scan(workscan, fl.get_mask())
1784	f.x = workscan._getabcissa(r)
1785	f.y = workscan._getspectrum(r)
1786	f.data = None
1787	f.fit()
1788
1789	# Show mask list
1790	masklist=workscan.get_masklist(fl.get_mask(),row=r)
1791	msg = "mask range: "+str(masklist)
1792	asaplog.push(msg, False)
1793
1794	if plot:
1795	f.plot(residual=True)
1796	x = raw_input("Accept fit ( [y]/n ): ")
1797	if x.upper() == 'N':
1798	self.blpars.append(None)
1799	self.masklists.append(None)
1800	continue
1801
1802	workscan._setspectrum(f.fitter.getresidual(), r)
1803	self.blpars.append(f.get_parameters())
1804	self.masklists.append(masklist)
1805	if plot:
1806	f._p.unmap()
1807	f._p = None
1808	workscan._add_history("auto_poly_baseline", varlist)
1809	if insitu:
1810	self._assign(workscan)
1811	else:
1812	return workscan
1813
1814	def rotate_linpolphase(self, angle):
1815	"""
1816	Rotate the phase of the complex polarization O=Q+iU correlation.
1817	This is always done in situ in the raw data. So if you call this
1818	function more than once then each call rotates the phase further.
1819	Parameters:
1820	angle: The angle (degrees) to rotate (add) by.
1821	Examples:
1822	scan.rotate_linpolphase(2.3)
1823	"""
1824	varlist = vars()
1825	self._math._rotate_linpolphase(self, angle)
1826	self._add_history("rotate_linpolphase", varlist)
1827	print_log()
1828	return
1829
1830
1831	def rotate_xyphase(self, angle):
1832	"""
1833	Rotate the phase of the XY correlation. This is always done in situ
1834	in the data. So if you call this function more than once
1835	then each call rotates the phase further.
1836	Parameters:
1837	angle: The angle (degrees) to rotate (add) by.
1838	Examples:
1839	scan.rotate_xyphase(2.3)
1840	"""
1841	varlist = vars()
1842	self._math._rotate_xyphase(self, angle)
1843	self._add_history("rotate_xyphase", varlist)
1844	print_log()
1845	return
1846
1847	def swap_linears(self):
1848	"""
1849	Swap the linear polarisations XX and YY, or better the first two
1850	polarisations as this also works for ciculars.
1851	"""
1852	varlist = vars()
1853	self._math._swap_linears(self)
1854	self._add_history("swap_linears", varlist)
1855	print_log()
1856	return
1857
1858	def invert_phase(self):
1859	"""
1860	Invert the phase of the complex polarisation
1861	"""
1862	varlist = vars()
1863	self._math._invert_phase(self)
1864	self._add_history("invert_phase", varlist)
1865	print_log()
1866	return
1867
1868	def add(self, offset, insitu=None):
1869	"""
1870	Return a scan where all spectra have the offset added
1871	Parameters:
1872	offset: the offset
1873	insitu: if False a new scantable is returned.
1874	Otherwise, the scaling is done in-situ
1875	The default is taken from .asaprc (False)
1876	"""
1877	if insitu is None: insitu = rcParams['insitu']
1878	self._math._setinsitu(insitu)
1879	varlist = vars()
1880	s = scantable(self._math._unaryop(self, offset, "ADD", False))
1881	s._add_history("add", varlist)
1882	print_log()
1883	if insitu:
1884	self._assign(s)
1885	else:
1886	return s
1887
1888	def scale(self, factor, tsys=True, insitu=None):
1889	"""
1890	Return a scan where all spectra are scaled by the give 'factor'
1891	Parameters:
1892	factor: the scaling factor (float or 1D float list)
1893	insitu: if False a new scantable is returned.
1894	Otherwise, the scaling is done in-situ
1895	The default is taken from .asaprc (False)
1896	tsys: if True (default) then apply the operation to Tsys
1897	as well as the data
1898	"""
1899	if insitu is None: insitu = rcParams['insitu']
1900	self._math._setinsitu(insitu)
1901	varlist = vars()
1902	s = None
1903	import numpy
1904	if isinstance(factor, list) or isinstance(factor, numpy.ndarray):
1905	if isinstance(factor[0], list) or isinstance(factor[0], numpy.ndarray):
1906	from asapmath import _array2dOp
1907	s = _array2dOp( self.copy(), factor, "MUL", tsys )
1908	else:
1909	s = scantable( self._math._arrayop( self.copy(), factor, "MUL", tsys ) )
1910	else:
1911	s = scantable(self._math._unaryop(self.copy(), factor, "MUL", tsys))
1912	s._add_history("scale", varlist)
1913	print_log()
1914	if insitu:
1915	self._assign(s)
1916	else:
1917	return s
1918
1919	def set_sourcetype(self, match, matchtype="pattern",
1920	sourcetype="reference"):
1921	"""
1922	Set the type of the source to be an source or reference scan
1923	using the provided pattern:
1924	Parameters:
1925	match: a Unix style pattern, regular expression or selector
1926	matchtype: 'pattern' (default) UNIX style pattern or
1927	'regex' regular expression
1928	sourcetype: the type of the source to use (source/reference)
1929	"""
1930	varlist = vars()
1931	basesel = self.get_selection()
1932	stype = -1
1933	if sourcetype.lower().startswith("r"):
1934	stype = 1
1935	elif sourcetype.lower().startswith("s"):
1936	stype = 0
1937	else:
1938	raise ValueError("Illegal sourcetype use s(ource) or r(eference)")
1939	if matchtype.lower().startswith("p"):
1940	matchtype = "pattern"
1941	elif matchtype.lower().startswith("r"):
1942	matchtype = "regex"
1943	else:
1944	raise ValueError("Illegal matchtype, use p(attern) or r(egex)")
1945	sel = selector()
1946	if isinstance(match, selector):
1947	sel = match
1948	else:
1949	sel.set_query("SRCNAME == %s('%s')" % (matchtype, match))
1950	self.set_selection(basesel+sel)
1951	self._setsourcetype(stype)
1952	self.set_selection(basesel)
1953	s._add_history("set_sourcetype", varlist)
1954
1955	def auto_quotient(self, preserve=True, mode='paired', verify=False):
1956	"""
1957	This function allows to build quotients automatically.
1958	It assumes the observation to have the same number of
1959	"ons" and "offs"
1960	Parameters:
1961	preserve: you can preserve (default) the continuum or
1962	remove it. The equations used are
1963	preserve: Output = Toff * (on/off) - Toff
1964	remove: Output = Toff * (on/off) - Ton
1965	mode: the on/off detection mode
1966	'paired' (default)
1967	identifies 'off' scans by the
1968	trailing '_R' (Mopra/Parkes) or
1969	'_e'/'_w' (Tid) and matches
1970	on/off pairs from the observing pattern
1971	'time'
1972	finds the closest off in time
1973
1974	"""
1975	modes = ["time", "paired"]
1976	if not mode in modes:
1977	msg = "please provide valid mode. Valid modes are %s" % (modes)
1978	raise ValueError(msg)
1979	varlist = vars()
1980	s = None
1981	if mode.lower() == "paired":
1982	basesel = self.get_selection()
1983	sel = selector()+basesel
1984	sel.set_query("SRCTYPE==1")
1985	self.set_selection(sel)
1986	offs = self.copy()
1987	sel.set_query("SRCTYPE==0")
1988	self.set_selection(sel)
1989	ons = self.copy()
1990	s = scantable(self._math._quotient(ons, offs, preserve))
1991	self.set_selection(basesel)
1992	elif mode.lower() == "time":
1993	s = scantable(self._math._auto_quotient(self, mode, preserve))
1994	s._add_history("auto_quotient", varlist)
1995	print_log()
1996	return s
1997
1998	def mx_quotient(self, mask = None, weight='median', preserve=True):
1999	"""
2000	Form a quotient using "off" beams when observing in "MX" mode.
2001	Parameters:
2002	mask: an optional mask to be used when weight == 'stddev'
2003	weight: How to average the off beams. Default is 'median'.
2004	preserve: you can preserve (default) the continuum or
2005	remove it. The equations used are
2006	preserve: Output = Toff * (on/off) - Toff
2007	remove: Output = Toff * (on/off) - Ton
2008	"""
2009	if mask is None: mask = ()
2010	varlist = vars()
2011	on = scantable(self._math._mx_extract(self, 'on'))
2012	preoff = scantable(self._math._mx_extract(self, 'off'))
2013	off = preoff.average_time(mask=mask, weight=weight, scanav=False)
2014	from asapmath import quotient
2015	q = quotient(on, off, preserve)
2016	q._add_history("mx_quotient", varlist)
2017	print_log()
2018	return q
2019
2020	def freq_switch(self, insitu=None):
2021	"""
2022	Apply frequency switching to the data.
2023	Parameters:
2024	insitu: if False a new scantable is returned.
2025	Otherwise, the swictching is done in-situ
2026	The default is taken from .asaprc (False)
2027	Example:
2028	none
2029	"""
2030	if insitu is None: insitu = rcParams['insitu']
2031	self._math._setinsitu(insitu)
2032	varlist = vars()
2033	s = scantable(self._math._freqswitch(self))
2034	s._add_history("freq_switch", varlist)
2035	print_log()
2036	if insitu: self._assign(s)
2037	else: return s
2038
2039	def recalc_azel(self):
2040	"""
2041	Recalculate the azimuth and elevation for each position.
2042	Parameters:
2043	none
2044	Example:
2045	"""
2046	varlist = vars()
2047	self._recalcazel()
2048	self._add_history("recalc_azel", varlist)
2049	print_log()
2050	return
2051
2052	def __add__(self, other):
2053	"""
2054	implicit on all axes and on Tsys
2055	"""
2056	return self._operation( other, "ADD" )
2057
2058	def __sub__(self, other):
2059	"""
2060	implicit on all axes and on Tsys
2061	"""
2062	return self._operation( other, 'SUB' )
2063
2064	def __mul__(self, other):
2065	"""
2066	implicit on all axes and on Tsys
2067	"""
2068	return self._operation( other, 'MUL' )
2069
2070	def __div__(self, other):
2071	"""
2072	implicit on all axes and on Tsys
2073	"""
2074	return self._operation( other, 'DIV' )
2075
2076	def get_fit(self, row=0):
2077	"""
2078	Print or return the stored fits for a row in the scantable
2079	Parameters:
2080	row: the row which the fit has been applied to.
2081	"""
2082	if row > self.nrow():
2083	return
2084	from asap.asapfit import asapfit
2085	fit = asapfit(self._getfit(row))
2086	if rcParams['verbose']:
2087	#print fit
2088	asaplog.push( '%s' %(fit) )
2089	print_log()
2090	return
2091	else:
2092	return fit.as_dict()
2093
2094	def flag_nans(self):
2095	"""
2096	Utility function to flag NaN values in the scantable.
2097	"""
2098	import numpy
2099	basesel = self.get_selection()
2100	for i in range(self.nrow()):
2101	sel = selector()+basesel
2102	sel.set_scans(self.getscan(i))
2103	sel.set_beams(self.getbeam(i))
2104	sel.set_ifs(self.getif(i))
2105	sel.set_polarisations(self.getpol(i))
2106	self.set_selection(sel)
2107	nans = numpy.isnan(self._getspectrum(0))
2108	if numpy.any(nans):
2109	bnans = [ bool(v) for v in nans]
2110	self.flag(bnans)
2111	self.set_selection(basesel)
2112
2113
2114	def _add_history(self, funcname, parameters):
2115	if not rcParams['scantable.history']:
2116	return
2117	# create date
2118	sep = "##"
2119	from datetime import datetime
2120	dstr = datetime.now().strftime('%Y/%m/%d %H:%M:%S')
2121	hist = dstr+sep
2122	hist += funcname+sep#cdate+sep
2123	if parameters.has_key('self'): del parameters['self']
2124	for k, v in parameters.iteritems():
2125	if type(v) is dict:
2126	for k2, v2 in v.iteritems():
2127	hist += k2
2128	hist += "="
2129	if isinstance(v2, scantable):
2130	hist += 'scantable'
2131	elif k2 == 'mask':
2132	if isinstance(v2, list) or isinstance(v2, tuple):
2133	hist += str(self._zip_mask(v2))
2134	else:
2135	hist += str(v2)
2136	else:
2137	hist += str(v2)
2138	else:
2139	hist += k
2140	hist += "="
2141	if isinstance(v, scantable):
2142	hist += 'scantable'
2143	elif k == 'mask':
2144	if isinstance(v, list) or isinstance(v, tuple):
2145	hist += str(self._zip_mask(v))
2146	else:
2147	hist += str(v)
2148	else:
2149	hist += str(v)
2150	hist += sep
2151	hist = hist[:-2] # remove trailing '##'
2152	self._addhistory(hist)
2153
2154
2155	def _zip_mask(self, mask):
2156	mask = list(mask)
2157	i = 0
2158	segments = []
2159	while mask[i:].count(1):
2160	i += mask[i:].index(1)
2161	if mask[i:].count(0):
2162	j = i + mask[i:].index(0)
2163	else:
2164	j = len(mask)
2165	segments.append([i, j])
2166	i = j
2167	return segments
2168
2169	def _get_ordinate_label(self):
2170	fu = "("+self.get_fluxunit()+")"
2171	import re
2172	lbl = "Intensity"
2173	if re.match(".K.", fu):
2174	lbl = "Brightness Temperature "+ fu
2175	elif re.match(".Jy.", fu):
2176	lbl = "Flux density "+ fu
2177	return lbl
2178
2179	def _check_ifs(self):
2180	nchans = [self.nchan(i) for i in range(self.nif(-1))]
2181	nchans = filter(lambda t: t > 0, nchans)
2182	return (sum(nchans)/len(nchans) == nchans[0])
2183
2184	def _fill(self, names, unit, average, getpt, antenna):
2185	import os
2186	from asap._asap import stfiller
2187	first = True
2188	fullnames = []
2189	for name in names:
2190	name = os.path.expandvars(name)
2191	name = os.path.expanduser(name)
2192	if not os.path.exists(name):
2193	msg = "File '%s' does not exists" % (name)
2194	if rcParams['verbose']:
2195	asaplog.push(msg)
2196	#print asaplog.pop().strip()
2197	print_log( 'ERROR' )
2198	return
2199	raise IOError(msg)
2200	fullnames.append(name)
2201	if average:
2202	asaplog.push('Auto averaging integrations')
2203	stype = int(rcParams['scantable.storage'].lower() == 'disk')
2204	for name in fullnames:
2205	tbl = Scantable(stype)
2206	r = stfiller(tbl)
2207	rx = rcParams['scantable.reference']
2208	r._setreferenceexpr(rx)
2209	msg = "Importing %s..." % (name)
2210	asaplog.push(msg, False)
2211	print_log()
2212	r._open(name, antenna, -1, -1, getpt)
2213	r._read()
2214	if average:
2215	tbl = self._math._average((tbl, ), (), 'NONE', 'SCAN')
2216	if not first:
2217	tbl = self._math._merge([self, tbl])
2218	Scantable.__init__(self, tbl)
2219	r._close()
2220	del r, tbl
2221	first = False
2222	if unit is not None:
2223	self.set_fluxunit(unit)
2224	#self.set_freqframe(rcParams['scantable.freqframe'])
2225
2226	def __getitem__(self, key):
2227	if key < 0:
2228	key += self.nrow()
2229	if key >= self.nrow():
2230	raise IndexError("Row index out of range.")
2231	return self._getspectrum(key)
2232
2233	def __setitem__(self, key, value):
2234	if key < 0:
2235	key += self.nrow()
2236	if key >= self.nrow():
2237	raise IndexError("Row index out of range.")
2238	if not hasattr(value, "__len__") or \
2239	len(value) > self.nchan(self.getif(key)):
2240	raise ValueError("Spectrum length doesn't match.")
2241	return self._setspectrum(value, key)
2242
2243	def __len__(self):
2244	return self.nrow()
2245
2246	def __iter__(self):
2247	for i in range(len(self)):
2248	yield self[i]
2249
2250	def _operation(self, other, opmode):
2251	varlist = vars()
2252	s = None
2253	import numpy
2254	if isinstance(other, scantable):
2255	s = scantable(self._math._binaryop(self.copy(), other, opmode))
2256	elif isinstance(other, float) or isinstance(other, int):
2257	if opmode == 'DIV' and float(other) == 0.0:
2258	raise ZeroDivisionError("Dividing by zero is not recommended")
2259	s = scantable(self._math._unaryop(self.copy(), other, opmode, False))
2260	elif isinstance(other, list) or isinstance(other, numpy.ndarray):
2261	if isinstance(other[0], list) or isinstance(other[0], numpy.ndarray):
2262	from asapmath import _array2dOp
2263	s = _array2dOp( self.copy(), other, opmode, False )
2264	else:
2265	s = scantable(self._math._arrayop(self.copy(), other, opmode, False))
2266	else:
2267	raise TypeError("Other input is not a scantable or float value or float list")
2268	opdic = {}
2269	opdic['ADD'] = '+'
2270	opdic['SUB'] = '-'
2271	opdic['MUL'] = '*'
2272	opdic['DIV'] = '/'
2273	s._add_history("operator %s" % opdic[opmode], varlist)
2274	print_log()
2275	return s
2276
2277

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

Context Navigation

source: branches/alma/python/scantable.py @ 1685

Download in other formats: