Context Navigation

source: trunk/python/scantable.py @ 714

Visit:

Last change on this file since 714 was 714, checked in by mar637, 19 years ago
updated logging; added public get_sourcename
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 58.1 KB

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

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

Download in other formats: