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

Last change on this file since 1653 was 1653, checked in by Kana Sugimoto, 15 years ago

New Development: No

JIRA Issue: Yes (CAS-1431)

Ready to Release: Yes

Interface Changes: Yes

What Interface Changed: A new parameter 'plot' (boolean) is added to

the method 'smooth'.

Test Programs:

Put in Release Notes: No

Module(s): sdsmooth and sdcal

Description:

A new parameter 'plot' (boolean) is added to the method scantable.smooth.
If plot=True, smoothed and unsmoothed spectra is plotted, for each row
in the input scantable, and you are asked if you accept the displayed
smoothing. If you don't, continues without smoothing the row.


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