Context Navigation

source: branches/Release12/python/scantable.py @ 798

Visit:

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

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

Download in other formats: