source: trunk/python/scantable.py@ 1261

Last change on this file since 1261 was 1259, checked in by mar637, 18 years ago

Merge from Release2.1.0b tag

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