source: tags/Release2.1.0b/python/scantable.py @ 1250

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

disabled returning the value form stat when in verbose

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 62.6 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._setrestfreqs(freqs[i]["value"],
815                                       freqs[i]["name"], "MHz")
816                    self._setselection(sel)
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        f = fitter()
1207        f.set_scan(self, mask)
1208        f.set_function(poly=order)
1209        s = f.auto_fit(insitu, plot=plot)
1210        s._add_history("poly_baseline", varlist)
1211        print_log()
1212        if insitu: self._assign(s)
1213        else: return s
1214
1215    def auto_poly_baseline(self, mask=[], edge=(0, 0), order=0,
1216                           threshold=3, plot=False, insitu=None):
1217        """
1218        Return a scan which has been baselined (all rows) by a polynomial.
1219        Spectral lines are detected first using linefinder and masked out
1220        to avoid them affecting the baseline solution.
1221
1222        Parameters:
1223            mask:       an optional mask retreived from scantable
1224            edge:       an optional number of channel to drop at
1225                        the edge of spectrum. If only one value is
1226                        specified, the same number will be dropped from
1227                        both sides of the spectrum. Default is to keep
1228                        all channels. Nested tuples represent individual
1229                        edge selection for different IFs (a number of spectral
1230                        channels can be different)
1231            order:      the order of the polynomial (default is 0)
1232            threshold:  the threshold used by line finder. It is better to
1233                        keep it large as only strong lines affect the
1234                        baseline solution.
1235            plot:       plot the fit and the residual. In this each
1236                        indivual fit has to be approved, by typing 'y'
1237                        or 'n'
1238            insitu:     if False a new scantable is returned.
1239                        Otherwise, the scaling is done in-situ
1240                        The default is taken from .asaprc (False)
1241
1242        Example:
1243            scan2=scan.auto_poly_baseline(order=7)
1244        """
1245        if insitu is None: insitu = rcParams['insitu']
1246        varlist = vars()
1247        from asap.asapfitter import fitter
1248        from asap.asaplinefind import linefinder
1249        from asap import _is_sequence_or_number as _is_valid
1250
1251        # check whether edge is set up for each IF individually
1252        individualedge = False;
1253        if len(edge) > 1:
1254            if isinstance(edge[0], list) or isinstance(edge[0], tuple):
1255                individualedge = True;
1256
1257        if not _is_valid(edge, int) and not individualedge:
1258            raise ValueError, "Parameter 'edge' has to be an integer or a \
1259            pair of integers specified as a tuple. Nested tuples are allowed \
1260            to make individual selection for different IFs."
1261
1262        curedge = (0, 0)
1263        if individualedge:
1264            for edgepar in edge:
1265                if not _is_valid(edgepar, int):
1266                    raise ValueError, "Each element of the 'edge' tuple has \
1267                                       to be a pair of integers or an integer."
1268        else:
1269            curedge = edge;
1270
1271        # setup fitter
1272        f = fitter()
1273        f.set_function(poly=order)
1274
1275        # setup line finder
1276        fl = linefinder()
1277        fl.set_options(threshold=threshold)
1278
1279        if not insitu:
1280            workscan = self.copy()
1281        else:
1282            workscan = self
1283
1284        fl.set_scan(workscan)
1285
1286        rows = range(workscan.nrow())
1287        asaplog.push("Processing:")
1288        for r in rows:
1289            msg = " Scan[%d] Beam[%d] IF[%d] Pol[%d] Cycle[%d]" % \
1290                (workscan.getscan(r), workscan.getbeam(r), workscan.getif(r), \
1291                 workscan.getpol(r), workscan.getcycle(r))
1292            asaplog.push(msg, False)
1293
1294            # figure out edge parameter
1295            if individualedge:
1296                if len(edge) >= workscan.getif(r):
1297                    raise RuntimeError, "Number of edge elements appear to " \
1298                                        "be less than the number of IFs"
1299                    curedge = edge[workscan.getif(r)]
1300
1301            # setup line finder
1302            fl.find_lines(r, mask, curedge)
1303            f.set_scan(workscan, fl.get_mask())
1304            f.x = workscan._getabcissa(r)
1305            f.y = workscan._getspectrum(r)
1306            f.data = None
1307            f.fit()
1308            x = f.get_parameters()
1309            if plot:
1310                f.plot(residual=True)
1311                x = raw_input("Accept fit ( [y]/n ): ")
1312                if x.upper() == 'N':
1313                    continue
1314            workscan._setspectrum(f.fitter.getresidual(), r)
1315        if plot:
1316            f._p.unmap()
1317            f._p = None
1318        workscan._add_history("auto_poly_baseline", varlist)
1319        if insitu:
1320            self._assign(workscan)
1321        else:
1322            return workscan
1323
1324    def rotate_linpolphase(self, angle):
1325        """
1326        Rotate the phase of the complex polarization O=Q+iU correlation.
1327        This is always done in situ in the raw data.  So if you call this
1328        function more than once then each call rotates the phase further.
1329        Parameters:
1330            angle:   The angle (degrees) to rotate (add) by.
1331        Examples:
1332            scan.rotate_linpolphase(2.3)
1333        """
1334        varlist = vars()
1335        self._math._rotate_linpolphase(self, angle)
1336        self._add_history("rotate_linpolphase", varlist)
1337        print_log()
1338        return
1339
1340
1341    def rotate_xyphase(self, angle):
1342        """
1343        Rotate the phase of the XY correlation.  This is always done in situ
1344        in the data.  So if you call this function more than once
1345        then each call rotates the phase further.
1346        Parameters:
1347            angle:   The angle (degrees) to rotate (add) by.
1348        Examples:
1349            scan.rotate_xyphase(2.3)
1350        """
1351        varlist = vars()
1352        self._math._rotate_xyphase(self, angle)
1353        self._add_history("rotate_xyphase", varlist)
1354        print_log()
1355        return
1356
1357    def swap_linears(self):
1358        """
1359        Swap the linear polarisations XX and YY
1360        """
1361        varlist = vars()
1362        self._math._swap_linears(self)
1363        self._add_history("swap_linears", varlist)
1364        print_log()
1365        return
1366
1367    def invert_phase(self):
1368        """
1369        Invert the phase of the complex polarisation
1370        """
1371        varlist = vars()
1372        self._math._invert_phase(self)
1373        self._add_history("invert_phase", varlist)
1374        print_log()
1375        return
1376
1377    def add(self, offset, insitu=None):
1378        """
1379        Return a scan where all spectra have the offset added
1380        Parameters:
1381            offset:      the offset
1382            insitu:      if False a new scantable is returned.
1383                         Otherwise, the scaling is done in-situ
1384                         The default is taken from .asaprc (False)
1385        """
1386        if insitu is None: insitu = rcParams['insitu']
1387        self._math._setinsitu(insitu)
1388        varlist = vars()
1389        s = scantable(self._math._unaryop(self, offset, "ADD", False))
1390        s._add_history("add", varlist)
1391        print_log()
1392        if insitu:
1393            self._assign(s)
1394        else:
1395            return s
1396
1397    def scale(self, factor, tsys=True, insitu=None, ):
1398        """
1399        Return a scan where all spectra are scaled by the give 'factor'
1400        Parameters:
1401            factor:      the scaling factor
1402            insitu:      if False a new scantable is returned.
1403                         Otherwise, the scaling is done in-situ
1404                         The default is taken from .asaprc (False)
1405            tsys:        if True (default) then apply the operation to Tsys
1406                         as well as the data
1407        """
1408        if insitu is None: insitu = rcParams['insitu']
1409        self._math._setinsitu(insitu)
1410        varlist = vars()
1411        s = scantable(self._math._unaryop(self, factor, "MUL", tsys))
1412        s._add_history("scale", varlist)
1413        print_log()
1414        if insitu:
1415            self._assign(s)
1416        else:
1417            return s
1418
1419    def auto_quotient(self, mode='time', preserve=True):
1420        """
1421        This function allows to build quotients automatically.
1422        It assumes the observation to have the same numer of
1423        "ons" and "offs"
1424        It will support "closest off in time" in the future
1425        Parameters:
1426            mode:           the on/off detection mode; 'suffix' (default)
1427                            'suffix' identifies 'off' scans by the
1428                            trailing '_R' (Mopra/Parkes) or
1429                            '_e'/'_w' (Tid)
1430            preserve:       you can preserve (default) the continuum or
1431                            remove it.  The equations used are
1432                            preserve: Output = Toff * (on/off) - Toff
1433                            remove:   Output = Toff * (on/off) - Ton
1434        """
1435        modes = ["time"]
1436        if not mode in modes:
1437            msg = "please provide valid mode. Valid modes are %s" % (modes)
1438            raise ValueError(msg)
1439        varlist = vars()
1440        s = scantable(self._math._auto_quotient(self, mode, preserve))
1441        s._add_history("auto_quotient", varlist)
1442        print_log()
1443        return s
1444
1445    def mx_quotient(self, mask = None, weight='median', preserve=True):
1446        """
1447        Form a quotient using "off" beams when observing in "MX" mode.
1448        Parameters:
1449            mask:           an optional mask to be used when weight == 'stddev'
1450            weight:         How to average the off beams.  Default is 'median'.
1451            preserve:       you can preserve (default) the continuum or
1452                            remove it.  The equations used are
1453                            preserve: Output = Toff * (on/off) - Toff
1454                            remove:   Output = Toff * (on/off) - Ton
1455        """
1456        if mask is None: mask = ()
1457        varlist = vars()
1458        on = scantable(self._math._mx_extract(self, 'on'))
1459        preoff = scantable(self._math._mx_extract(self, 'off'))
1460        off = preoff.average_time(mask=mask, weight=weight, scanav=False)
1461        from asapmath  import quotient
1462        q = quotient(on, off, preserve)
1463        q._add_history("mx_quotient", varlist)
1464        print_log()
1465        return q
1466
1467    def freq_switch(self, insitu=None):
1468        """
1469        Apply frequency switching to the data.
1470        Parameters:
1471            insitu:      if False a new scantable is returned.
1472                         Otherwise, the swictching is done in-situ
1473                         The default is taken from .asaprc (False)
1474        Example:
1475            none
1476        """
1477        if insitu is None: insitu = rcParams['insitu']
1478        self._math._setinsitu(insitu)
1479        varlist = vars()
1480        s = scantable(self._math._freqswitch(self))
1481        s._add_history("freq_switch", varlist)
1482        print_log()
1483        if insitu: self._assign(s)
1484        else: return s
1485
1486    def recalc_azel(self):
1487        """
1488        Recalculate the azimuth and elevation for each position.
1489        Parameters:
1490            none
1491        Example:
1492        """
1493        varlist = vars()
1494        self._recalcazel()
1495        self._add_history("recalc_azel", varlist)
1496        print_log()
1497        return
1498
1499    def __add__(self, other):
1500        varlist = vars()
1501        s = None
1502        if isinstance(other, scantable):
1503            print "scantable + scantable NYI"
1504            return
1505        elif isinstance(other, float):
1506            s = scantable(self._math._unaryop(self, other, "ADD", False))
1507        else:
1508            raise TypeError("Other input is not a scantable or float value")
1509        s._add_history("operator +", varlist)
1510        print_log()
1511        return s
1512
1513    def __sub__(self, other):
1514        """
1515        implicit on all axes and on Tsys
1516        """
1517        varlist = vars()
1518        s = None
1519        if isinstance(other, scantable):
1520            print "scantable - scantable NYI"
1521            return
1522        elif isinstance(other, float):
1523            s = scantable(self._math._unaryop(self, other, "SUB", False))
1524        else:
1525            raise TypeError("Other input is not a scantable or float value")
1526        s._add_history("operator -", varlist)
1527        print_log()
1528        return s
1529
1530    def __mul__(self, other):
1531        """
1532        implicit on all axes and on Tsys
1533        """
1534        varlist = vars()
1535        s = None
1536        if isinstance(other, scantable):
1537            print "scantable * scantable NYI"
1538            return
1539        elif isinstance(other, float):
1540            s = scantable(self._math._unaryop(self, other, "MUL", False))
1541        else:
1542            raise TypeError("Other input is not a scantable or float value")
1543        s._add_history("operator *", varlist)
1544        print_log()
1545        return s
1546
1547
1548    def __div__(self, other):
1549        """
1550        implicit on all axes and on Tsys
1551        """
1552        varlist = vars()
1553        s = None
1554        if isinstance(other, scantable):
1555            print "scantable / scantable NYI"
1556            return
1557        elif isinstance(other, float):
1558            if other == 0.0:
1559                raise ZeroDivisionError("Dividing by zero is not recommended")
1560            s = scantable(self._math._unaryop(self, other, "DIV", False))
1561        else:
1562            raise TypeError("Other input is not a scantable or float value")
1563        s._add_history("operator /", varlist)
1564        print_log()
1565        return s
1566
1567    def get_fit(self, row=0):
1568        """
1569        Print or return the stored fits for a row in the scantable
1570        Parameters:
1571            row:    the row which the fit has been applied to.
1572        """
1573        if row > self.nrow():
1574            return
1575        from asap.asapfit import asapfit
1576        fit = asapfit(self._getfit(row))
1577        if rcParams['verbose']:
1578            print fit
1579            return
1580        else:
1581            return fit.as_dict()
1582
1583    def _add_history(self, funcname, parameters):
1584        # create date
1585        sep = "##"
1586        from datetime import datetime
1587        dstr = datetime.now().strftime('%Y/%m/%d %H:%M:%S')
1588        hist = dstr+sep
1589        hist += funcname+sep#cdate+sep
1590        if parameters.has_key('self'): del parameters['self']
1591        for k, v in parameters.iteritems():
1592            if type(v) is dict:
1593                for k2, v2 in v.iteritems():
1594                    hist += k2
1595                    hist += "="
1596                    if isinstance(v2, scantable):
1597                        hist += 'scantable'
1598                    elif k2 == 'mask':
1599                        if isinstance(v2, list) or isinstance(v2, tuple):
1600                            hist += str(self._zip_mask(v2))
1601                        else:
1602                            hist += str(v2)
1603                    else:
1604                        hist += str(v2)
1605            else:
1606                hist += k
1607                hist += "="
1608                if isinstance(v, scantable):
1609                    hist += 'scantable'
1610                elif k == 'mask':
1611                    if isinstance(v, list) or isinstance(v, tuple):
1612                        hist += str(self._zip_mask(v))
1613                    else:
1614                        hist += str(v)
1615                else:
1616                    hist += str(v)
1617            hist += sep
1618        hist = hist[:-2] # remove trailing '##'
1619        self._addhistory(hist)
1620
1621
1622    def _zip_mask(self, mask):
1623        mask = list(mask)
1624        i = 0
1625        segments = []
1626        while mask[i:].count(1):
1627            i += mask[i:].index(1)
1628            if mask[i:].count(0):
1629                j = i + mask[i:].index(0)
1630            else:
1631                j = len(mask)
1632            segments.append([i, j])
1633            i = j
1634        return segments
1635
1636    def _get_ordinate_label(self):
1637        fu = "("+self.get_fluxunit()+")"
1638        import re
1639        lbl = "Intensity"
1640        if re.match(".K.", fu):
1641            lbl = "Brightness Temperature "+ fu
1642        elif re.match(".Jy.", fu):
1643            lbl = "Flux density "+ fu
1644        return lbl
1645
1646    def _check_ifs(self):
1647        nchans = [self.nchan(i) for i in range(self.nif(-1))]
1648        nchans = filter(lambda t: t > 0, nchans)
1649        return (sum(nchans)/len(nchans) == nchans[0])
1650
1651    def _fill(self, names, unit, average):
1652        import os
1653        from asap._asap import stfiller
1654        first = True
1655        fullnames = []
1656        for name in names:
1657            name = os.path.expandvars(name)
1658            name = os.path.expanduser(name)
1659            if not os.path.exists(name):
1660                msg = "File '%s' does not exists" % (name)
1661                if rcParams['verbose']:
1662                    asaplog.push(msg)
1663                    print asaplog.pop().strip()
1664                    return
1665                raise IOError(msg)
1666            fullnames.append(name)
1667        if average:
1668            asaplog.push('Auto averaging integrations')
1669        stype = int(rcParams['scantable.storage'].lower() == 'disk')
1670        for name in fullnames:
1671            tbl = Scantable(stype)
1672            r = stfiller(tbl)
1673            msg = "Importing %s..." % (name)
1674            asaplog.push(msg, False)
1675            print_log()
1676            r._open(name, -1, -1)
1677            r._read()
1678            #tbl = r._getdata()
1679            if average:
1680                tbl = self._math._average((tbl, ), (), 'NONE', 'SCAN')
1681                #tbl = tbl2
1682            if not first:
1683                tbl = self._math._merge([self, tbl])
1684                #tbl = tbl2
1685            Scantable.__init__(self, tbl)
1686            r._close()
1687            del r, tbl
1688            first = False
1689        if unit is not None:
1690            self.set_fluxunit(unit)
1691        self.set_freqframe(rcParams['scantable.freqframe'])
1692
Note: See TracBrowser for help on using the repository browser.