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

Last change on this file since 787 was 787, checked in by phi196, 18 years ago

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