Ignore:
Location:
trunk/python
Files:
1 added
13 edited

Legend:

Unmodified
Added
Removed
  • trunk/python/CMakeLists.txt

    r2704 r3112  
    3535     ${PYTHONDIR}/opacity.py
    3636     ${PYTHONDIR}/parameters.py
    37 #     ${PYTHONDIR}/plotter2.py
     37     ${PYTHONDIR}/plotter2.py
    3838     ${PYTHONDIR}/sbseparator.py
    3939     ${PYTHONDIR}/scantable.py
  • trunk/python/__init__.py

    r2704 r3112  
    5454from asapgrid import asapgrid, asapgrid2
    5555from edgemarker import edgemarker
    56 #from plotter2 import plotter2
     56if is_casapy():
     57    from plotter2 import plotter2
    5758from sbseparator import sbseparator
    5859from _asap import srctype
    5960
    60 __date__ = '$Date$'.split()[1]
    61 __version__  = 'trunk'
     61__date__ = get_asap_revdate()
     62__version__  = '4.0.0a'
    6263__revision__ = get_revision()
    6364
  • trunk/python/asapfitter.py

    r2704 r3112  
    157157
    158158        if self.data is not None:
     159            if self.data._getflagrow(row):
     160                raise RuntimeError,"Can not fit flagged row."
    159161            self.x = self.data._getabcissa(row)
    160162            self.y = self.data._getspectrum(row)
  • trunk/python/asapgrid.py

    r2704 r3112  
    486486                    #print irow
    487487        # show image
    488         extent=[self.trc[0]+0.5*self.cellx,
    489                 self.blc[0]-0.5*self.cellx,
     488        extent=[self.blc[0]-0.5*self.cellx,
     489                self.trc[0]+0.5*self.cellx,
    490490                self.blc[1]-0.5*self.celly,
    491491                self.trc[1]+0.5*self.celly]
  • trunk/python/asaplotbase.py

    r2704 r3112  
    99from matplotlib.figure import Figure, Text
    1010from matplotlib.font_manager import FontProperties as FP
    11 from numpy import sqrt
     11from numpy import sqrt, ceil
    1212from matplotlib import rc, rcParams
    1313from matplotlib.ticker import OldScalarFormatter
     
    100100        return False
    101101
     102    def _subplotsOk(self, rows, cols, npanel=0):
     103        """
     104        Check if the axes in subplots are actually the ones plotted on
     105        the figure. Returns a bool.
     106        This method is to detect manual layout changes using mpl methods.
     107        """
     108        # compare with user defined layout
     109        if (rows is not None) and (rows != self.rows):
     110            return False
     111        if (cols is not None) and (cols != self.cols):
     112            return False
     113        # check number of subplots
     114        figaxes = self.figure.get_axes()
     115        np = self.rows*self.cols
     116        if npanel > np:
     117            return False
     118        if len(figaxes) != np:
     119            return False
     120        if len(self.subplots) != len(figaxes):
     121            return False
     122        # compare axes instance in this class and on the plotter
     123        ok = True
     124        for ip in range(np):
     125            if self.subplots[ip]['axes'] != figaxes[ip]:
     126                ok = False
     127                break
     128        return ok
    102129
    103130    ### Delete artists ###
  • trunk/python/asapmath.py

    r2704 r3112  
     1import re
    12from asap.scantable import scantable
    23from asap.parameters import rcParams
     
    9091    else:
    9192        #s = scantable(stm._average(alignedlst, mask, weight.upper(), scanav))
    92         s = scantable(stm._new_average(alignedlst, compel, mask, weight.upper(), scanav))
     93        s = scantable(stm._new_average(alignedlst, compel, mask,
     94                                       weight.upper(), scanav))
    9395    s._add_history("average_time",varlist)
    9496
     
    391393            p.quit()
    392394            del p
    393             return scabtab
     395            return scantab
    394396        p.quit()
    395397        del p
     
    611613            p.quit()
    612614            del p
    613             return scabtab
     615            return scantab
    614616        p.quit()
    615617        del p
     
    814816            p.quit()
    815817            del p
    816             return scabtab
     818            return scantab
    817819        p.quit()
    818820        del p
     
    822824
    823825@asaplog_post_dec
    824 def merge(*args):
     826def merge(*args, **kwargs):
    825827    """
    826828    Merge a list of scanatables, or comma-sperated scantables into one
     
    828830    Parameters:
    829831        A list [scan1, scan2] or scan1, scan2.
     832        freq_tol: frequency tolerance for merging IFs. numeric values
     833                  in units of Hz (1.0e6 -> 1MHz) and string ('1MHz')
     834                  is allowed.
    830835    Example:
    831836        myscans = [scan1, scan2]
     
    833838        # or equivalent
    834839        sameallscans = merge(scan1, scan2)
     840        # with freqtol
     841        allscans = merge(scan1, scan2, freq_tol=1.0e6)
     842        # or equivalently
     843        allscans = merge(scan1, scan2, freq_tol='1MHz')
    835844    """
    836845    varlist = vars()
     
    841850    else:
    842851        lst = tuple(args)
     852    if kwargs.has_key('freq_tol'):
     853        freq_tol = str(kwargs['freq_tol'])
     854        if len(freq_tol) > 0 and re.match('.+[GMk]Hz$', freq_tol) is None:
     855            freq_tol += 'Hz'
     856    else:
     857        freq_tol = ''
    843858    varlist["args"] = "%d scantables" % len(lst)
    844859    # need special formatting her for history...
     
    849864            msg = "Please give a list of scantables"
    850865            raise TypeError(msg)
    851     s = scantable(stm._merge(lst))
     866    s = scantable(stm._merge(lst, freq_tol))
    852867    s._add_history("merge", varlist)
    853868    return s
     
    949964
    950965@asaplog_post_dec
    951 def splitant(filename, outprefix='',overwrite=False):
     966def splitant(filename, outprefix='',overwrite=False, getpt=True):
    952967    """
    953968    Split Measurement set by antenna name, save data as a scantables,
    954     and return a list of filename.
     969    and return a list of filename. Note that frequency reference frame
     970    is imported as it is in Measurement set.
    955971    Notice this method can only be available from CASA.
    956972    Prameter
     
    963979                    The default False is to return with warning
    964980                    without writing the output. USE WITH CARE.
    965 
     981       getpt        Whether to import direction from MS/POINTING
     982                    table or not. Default is True (import direction).
    966983    """
    967984    # Import the table toolkit from CASA
    968     from casac import casac
     985    from taskinit import gentools
    969986    from asap.scantable import is_ms
    970     tb = casac.table()
     987    tb = gentools(['tb'])[0]
    971988    # Check the input filename
    972989    if isinstance(filename, str):
     
    978995            raise IOError(s)
    979996        # check if input file is MS
    980         #if not os.path.isdir(filename) \
    981         #       or not os.path.exists(filename+'/ANTENNA') \
    982         #       or not os.path.exists(filename+'/table.f1'):
    983997        if not is_ms(filename):
    984998            s = "File '%s' is not a Measurement set." % (filename)
     
    9961010    tb.open(tablename=filename,nomodify=True)
    9971011    ant1=tb.getcol('ANTENNA1',0,-1,1)
    998     #anttab=tb.getkeyword('ANTENNA').split()[-1]
    9991012    anttab=tb.getkeyword('ANTENNA').lstrip('Table: ')
    10001013    tb.close()
    1001     #tb.open(tablename=filename+'/ANTENNA',nomodify=True)
    10021014    tb.open(tablename=anttab,nomodify=True)
    10031015    nant=tb.nrows()
    10041016    antnames=tb.getcol('NAME',0,nant,1)
    10051017    tb.close()
    1006     tmpname='asapmath.splitant.tmp'
    10071018    for antid in set(ant1):
    1008         tb.open(tablename=filename,nomodify=True)
    1009         tbsel=tb.query('ANTENNA1 == %s && ANTENNA2 == %s'%(antid,antid),tmpname)
    1010         scan=scantable(tmpname,average=False,getpt=True,antenna=int(antid))
     1019        scan=scantable(filename,average=False,antenna=int(antid),getpt=getpt)
    10111020        outname=prefix+antnames[antid]+'.asap'
    10121021        scan.save(outname,format='ASAP',overwrite=overwrite)
    1013         tbsel.close()
    1014         tb.close()
    1015         del tbsel
    10161022        del scan
    10171023        outfiles.append(outname)
    1018         os.system('rm -rf '+tmpname)
    1019     del tb
    10201024    return outfiles
    10211025
    10221026@asaplog_post_dec
    1023 def _array2dOp( scan, value, mode="ADD", tsys=False, insitu=None):
     1027def _array2dOp( scan, value, mode="ADD", tsys=False, insitu=None, skip_flaggedrow=False):
    10241028    """
    10251029    This function is workaround on the basic operation of scantable
     
    10321036    insitu:  if False, a new scantable is returned.
    10331037             Otherwise, the array operation is done in-sitsu.
     1038    skip_flaggedrow: skip operation for row-flagged spectra.
    10341039    """
    10351040    if insitu is None: insitu = rcParams['insitu']
     
    10401045    stm._setinsitu(insitu)
    10411046    if len( value ) == 1:
    1042         s = scantable( stm._arrayop( scan, value[0], mode, tsys ) )
     1047        s = scantable( stm._arrayop( scan, value[0], mode, tsys, skip_flaggedrow ) )
    10431048    elif len( value ) != nrow:
    10441049        raise ValueError( 'len(value) must be 1 or conform to scan.nrow()' )
     
    10581063            s.set_selection( sel )
    10591064            if len( value[irow] ) == 1:
    1060                 stm._unaryop( s, value[irow][0], mode, tsys )
     1065                stm._unaryop( s, value[irow][0], mode, tsys, skip_flaggedrow )
    10611066            else:
    10621067                #stm._arrayop( s, value[irow], mode, tsys, 'channel' )
    1063                 stm._arrayop( s, value[irow], mode, tsys )
     1068                stm._arrayop( s, value[irow], mode, tsys, skip_flaggedrow )
    10641069        s.set_selection(basesel)
    10651070    return s
  • trunk/python/asapplotter.py

    r2704 r3112  
    112112        self._plotter.legend(self._legendloc)
    113113
     114    ### TODO: it's probably better to define following two methods in
     115    ###       backend dependent class.
    114116    def _new_custombar(self):
    115117        backend=matplotlib.get_backend()
     
    130132            return True
    131133        return False
     134    ### end of TODO
    132135
    133136    def _assert_plotter(self,action="status",errmsg=None):
     
    276279
    277280
    278     ### Forwards to matplotlib axes ###
     281    ### Forwards to methods in matplotlib axes ###
    279282    def text(self, *args, **kwargs):
    280283        self._assert_plotter(action="reload")
     
    415418                del self._data
    416419                msg = "A new scantable is set to the plotter. "\
    417                       "The masks and data selections are reset."
    418                 asaplog.push( msg )
     420                      "The masks, data selections, and labels are reset."
     421                asaplog.push(msg)
    419422            self._data = scan
    420423            # reset
     
    514517        self._rows = rows
    515518        self._cols = cols
     519        # new layout is set. need to reset counters for multi page plotting
     520        self._reset_counters()
    516521        if refresh and self._data: self.plot(self._data)
    517522        return
     
    905910            msg = "Can only set mask after a first call to plot()"
    906911            raise RuntimeError(msg)
    907         if len(mask):
     912        if (mask is not None) and len(mask):
    908913            if isinstance(mask, list) or isinstance(mask, tuple):
    909914                self._usermask = array(mask)
     
    926931    ### Reset methods ###
    927932    def _reset(self):
    928         self._usermask = []
    929         self._usermaskspectra = None
     933        """Reset method called when new data is set"""
     934        # reset selections and masks
     935        self.set_selection(None, False)
     936        self.set_mask(None, None, False)
     937        # reset offset
    930938        self._offset = None
    931         self.set_selection(None, False)
     939        # reset header
    932940        self._reset_header()
     941        # reset labels
     942        self._lmap = None # related to stack
     943        self.set_title(None, None, False)
     944        self.set_ordinate(None, None, False)
     945        self.set_abcissa(None, None, False)
    933946
    934947    def _reset_header(self):
     
    938951        self._startrow = 0
    939952        self._ipanel = -1
     953        self._panelrows = []
    940954        self._reset_header()
    941         self._panelrows = []
    942955        if self.casabar_exists():
    943956            self._plotter.figmgr.casabar.set_pagecounter(1)
     
    10311044        if self._panelling == 'i':
    10321045            ganged = False
    1033         if not firstpage:
    1034             # not the first page just clear the axis
     1046        if (not firstpage) and \
     1047               self._plotter._subplotsOk(self._rows, self._cols, n):
     1048            # Not the first page and subplot number is ok.
     1049            # Just clear the axis
    10351050            nx = self._plotter.cols
    10361051            ipaxx = n - nx - 1 #the max panel id to supress x-label
     
    13021317            return start,end
    13031318
     1319    def _get_date_axis_setup(self, dates, axlim=None):
     1320        """
     1321        Returns proper axis title and formatters for a list of dates
     1322        Input
     1323            dates : a list of datetime objects returned by,
     1324                    e.g. scantable.get_time(asdatetime=True)
     1325            axlim : a tuple of min and max day range in the plot axis.
     1326                    if defined, the values are taken into account.
     1327        Output
     1328            a set of
     1329            * date axis title string
     1330            * formatter of date axis
     1331            * major axis locator
     1332            * minor axis locator
     1333        """
     1334        from matplotlib import pylab as PL
     1335        from matplotlib.dates import DateFormatter
     1336        from matplotlib.dates import HourLocator, MinuteLocator,SecondLocator, DayLocator, YearLocator, MonthLocator
     1337        t = PL.date2num(dates)
     1338        tmin = min(t)
     1339        tmax = max(t)
     1340        if axlim is not None:
     1341            tmin = min(tmin, min(axlim))
     1342            tmax = max(tmax, max(axlim))
     1343        tdel = tmax - tmin # interval in day
     1344        dstr = dates[0].strftime('%Y/%m/%d')
     1345        if tdel > 365.0: # >1year (also the case for single or very small time range)
     1346            majloc = YearLocator()
     1347            minloc = MonthLocator(range(1,12,6))
     1348            timefmt = DateFormatter('%Y/%m/%d')
     1349        elif tdel > 1.0: # >1day
     1350            dstr2 = dates[len(dates)-1].strftime('%Y/%m/%d')
     1351            dstr = dstr + " - " + dstr2
     1352            majloc = DayLocator()
     1353            minloc = HourLocator(range(0,23,12))
     1354            timefmt = DateFormatter("%b%d")
     1355        elif tdel > 24./60.: # 9.6h - 1day
     1356            timefmt = DateFormatter('%H:%M')
     1357            majloc = HourLocator()
     1358            minloc = MinuteLocator(range(0,60,30))
     1359        elif tdel > 2./24.: # 2h-9.6h
     1360            timefmt = DateFormatter('%H:%M')
     1361            majloc = HourLocator()
     1362            minloc = MinuteLocator(range(0,60,10))
     1363        elif tdel> 10./24./60.: # 10min-2h
     1364            timefmt = DateFormatter('%H:%M')
     1365            majloc = MinuteLocator(range(0,60,10))
     1366            minloc = MinuteLocator()
     1367        else: # <10min
     1368            timefmt = DateFormatter('%H:%M')
     1369            majloc = MinuteLocator()
     1370            minloc = SecondLocator(30)
     1371        return (dstr, timefmt, majloc, minloc)
     1372
    13041373    def plotazel(self, scan=None, outfile=None):
    13051374        """
     
    13091378        visible = rcParams['plotter.gui']
    13101379        from matplotlib import pylab as PL
    1311         from matplotlib.dates import DateFormatter
    13121380        from pytz import timezone
    1313         from matplotlib.dates import HourLocator, MinuteLocator,SecondLocator, DayLocator
    13141381        from matplotlib.ticker import MultipleLocator
    1315         from numpy import array, pi
     1382        from numpy import array, pi, ma
    13161383        if self._plotter and (PL.gcf() == self._plotter.figure):
    13171384            # the current figure is ASAP plotter. Use mpl plotter
     
    13251392        self._data = scan
    13261393        dates = self._data.get_time(asdatetime=True)
     1394        # for flag handling
     1395        mask = [ self._data._is_all_chan_flagged(i) for i in range(self._data.nrow())]
    13271396        t = PL.date2num(dates)
    13281397        tz = timezone('UTC')
     
    13371406                                 wspace=wsp,hspace=hsp)
    13381407
    1339         tdel = max(t) - min(t)
     1408        tdel = max(t) - min(t) # interval in day
    13401409        ax = PL.subplot(2,1,1)
    1341         el = array(self._data.get_elevation())*180./pi
     1410        el = ma.masked_array(array(self._data.get_elevation())*180./pi, mask)
    13421411        PL.ylabel('El [deg.]')
    1343         dstr = dates[0].strftime('%Y/%m/%d')
    1344         if tdel > 1.0:
    1345             dstr2 = dates[len(dates)-1].strftime('%Y/%m/%d')
    1346             dstr = dstr + " - " + dstr2
    1347             majloc = DayLocator()
    1348             minloc = HourLocator(range(0,23,12))
    1349             timefmt = DateFormatter("%b%d")
    1350         elif tdel > 24./60.:
    1351             timefmt = DateFormatter('%H:%M')
    1352             majloc = HourLocator()
    1353             minloc = MinuteLocator(30)
    1354         else:
    1355             timefmt = DateFormatter('%H:%M')
    1356             majloc = MinuteLocator(interval=5)
    1357             minloc = SecondLocator(30)
    1358 
     1412        (dstr, timefmt, majloc, minloc) = self._get_date_axis_setup(dates)
     1413       
    13591414        PL.title(dstr)
    13601415        if tdel == 0.0:
     
    13631418        else:
    13641419            PL.plot_date(t,el,'o', markersize=2, markerfacecolor='b', markeredgecolor='b',tz=tz)
    1365             #ax.grid(True)
    1366             ax.xaxis.set_major_formatter(timefmt)
    1367             ax.xaxis.set_major_locator(majloc)
    1368             ax.xaxis.set_minor_locator(minloc)
     1420            #ax.xaxis.set_major_formatter(timefmt)
     1421            #ax.xaxis.set_major_locator(majloc)
     1422            #ax.xaxis.set_minor_locator(minloc)
    13691423        ax.yaxis.grid(True)
    1370         yloc = MultipleLocator(30)
    13711424        ax.set_ylim(0,90)
    1372         ax.yaxis.set_major_locator(yloc)
     1425        #yloc = MultipleLocator(30)
     1426        #ax.yaxis.set_major_locator(yloc)
    13731427        if tdel > 1.0:
    13741428            labels = ax.get_xticklabels()
     
    13771431
    13781432        # Az plot
    1379         az = array(self._data.get_azimuth())*180./pi
     1433        az = ma.masked_array(array(self._data.get_azimuth())*180./pi, mask)
    13801434        if min(az) < 0:
    13811435            for irow in range(len(az)):
     
    13891443        else:
    13901444            PL.plot_date(t,az,'o', markersize=2,markeredgecolor='b',markerfacecolor='b',tz=tz)
    1391             ax2.xaxis.set_major_formatter(timefmt)
    1392             ax2.xaxis.set_major_locator(majloc)
    1393             ax2.xaxis.set_minor_locator(minloc)
    1394         #ax2.grid(True)
     1445            #ax2.xaxis.set_major_formatter(timefmt)
     1446            #ax2.xaxis.set_major_locator(majloc)
     1447            #ax2.xaxis.set_minor_locator(minloc)
    13951448        ax2.set_ylim(0,360)
    13961449        ax2.yaxis.grid(True)
    13971450        #hfmt = DateFormatter('%H')
    13981451        #hloc = HourLocator()
    1399         yloc = MultipleLocator(60)
    1400         ax2.yaxis.set_major_locator(yloc)
     1452        #yloc = MultipleLocator(60)
     1453        #ax2.yaxis.set_major_locator(yloc)
    14011454        if tdel > 1.0:
    14021455            labels = ax2.get_xticklabels()
    14031456            PL.setp(labels, fontsize=10)
    1404             PL.xlabel('Time (UT [day])')
    1405         else:
    1406             PL.xlabel('Time (UT [hour])')
     1457        #    PL.xlabel('Time (UT [day])')
     1458        #else:
     1459        #    PL.xlabel('Time (UT [hour])')
     1460        PL.xlabel('Time (UT)')
    14071461
    14081462        PL.ion()
     
    14171471        plot telescope pointings
    14181472        Parameters:
    1419             infile  : input filename or scantable instance
     1473            scan    : input scantable instance
    14201474            colorby : change color by either
    14211475                      'type'(source type)|'scan'|'if'|'pol'|'beam'
     
    14251479        """
    14261480        self._plotmode = "pointing"
    1427         from numpy import array, pi
     1481        from numpy import array, pi, ma
    14281482        from asap import scantable
    14291483        # check for scantable
     
    15011555                self._data.set_selection(basesel)
    15021556                continue
    1503             print "Plotting direction of %s = %s" % (colorby, str(idx))
     1557            #print "Plotting direction of %s = %s" % (colorby, str(idx))
    15041558            # getting data to plot
    15051559            dir = array(self._data.get_directionval()).transpose()
     1560            # for flag handling
     1561            mask = [ self._data._is_all_chan_flagged(i) for i in range(self._data.nrow())]
    15061562            ra = dir[0]*180./pi
    1507             dec = dir[1]*180./pi
     1563            dec = ma.masked_array(dir[1]*180./pi, mask)
    15081564            # actual plot
    15091565            self._plotter.set_line(label=(sellab+str(idx)))
     
    15311587        # reverse x-axis
    15321588        xmin, xmax = self.gca().get_xlim()
    1533         self._plotter.set_limits(xlim=[xmax,xmin])
     1589        ymin, ymax = self.gca().get_ylim()
     1590        # expand plotrange if xmin==xmax or ymin==ymax
     1591        if abs(ymax-ymin) < 1.e-3: #~4arcsec
     1592            delx = 0.5*abs(xmax - xmin)
     1593            if delx < 5.e-4:
     1594                dxy = 5.e-4 #~2arcsec
     1595                (ymin, ymax) = (ymin-dxy, ymax+dxy)
     1596                (xmin, xmax) = (xmin-dxy, xmax+dxy)
     1597            (ymin, ymax) = (ymin-delx, ymax+delx)
     1598        elif abs(xmax-xmin) < 1.e-3:
     1599            dely = 0.5*abs(ymax - ymin)
     1600            (xmin, xmax) = (xmin-dely, xmax+dely)
     1601        self._plotter.set_limits(xlim=[xmax,xmin], ylim=[ymin, ymax])
    15341602
    15351603        self._plotter.release()
     
    15871655    # plotting in time is not yet implemented..
    15881656    @asaplog_post_dec
    1589     def plottp(self, scan=None):
     1657    def plottp(self, scan=None, colorby=''):
     1658        """
     1659        Plot averaged spectra (total power) in time or in row ID (colorby='')
     1660        Parameters:
     1661            scan    : input scantable instance
     1662            colorby : change color by either
     1663                      'type'(source type)|'scan'|'if'|'pol'|'beam'|''
     1664        """
    15901665        self._plotmode = "totalpower"
    15911666        from asap import scantable
     
    16181693            left=lef,bottom=bot,right=rig,top=top,wspace=wsp,hspace=hsp)
    16191694        if self.casabar_exists(): self._plotter.figmgr.casabar.disable_button()
    1620         self._plottp(self._data)
     1695        if len(colorby) == 0:
     1696            self._plottp(self._data)
     1697        else:
     1698            self._plottp_in_time(self._data,colorby)
    16211699        if self._minmaxy is not None:
    16221700            self._plotter.set_limits(ylim=self._minmaxy)
     
    16251703        self._plotter.show(hardrefresh=False)
    16261704        return
     1705
     1706    def _plottp_in_time(self,scan,colorby):
     1707        """
     1708        private method for plotting total power data in time
     1709        Parameters:
     1710            scan    : input scantable instance
     1711            colorby : change color by either
     1712                      'type'(source type)|'scan'|'if'|'pol'|'beam'
     1713        """
     1714        from numpy import ma, array, arange, logical_not
     1715        r=0
     1716        nr = scan.nrow()
     1717        a0,b0 = -1,-1
     1718        allxlim = []
     1719        allylim = []
     1720        y=[]
     1721        self._plotter.set_panels()
     1722        self._plotter.palette(0)
     1723        # check of overlay settings
     1724        time_types = ['type','scan'] # time dependent meta-data
     1725        misc_types = ['if','pol','beam'] # time independent meta-data
     1726        validtypes=time_types + misc_types
     1727        stype = None
     1728        col_msg = "Invalid choice of 'colorby' (choices: %s)" % str(validtypes)
     1729        colorby = colorby.lower()
     1730        if (colorby in validtypes):
     1731            stype = colorby[0]
     1732        elif len(colorby) > 0:
     1733            raise ValueError(col_msg)
     1734        if not stype:
     1735            raise ValueError(col_msg)
     1736        # Selection and sort order
     1737        basesel = scan.get_selection()
     1738        if colorby in misc_types: misc_types.pop(misc_types.index(colorby))
     1739        sel_lbl = ""
     1740        for meta in misc_types:
     1741            idx = getattr(scan,'get'+meta+'nos')()
     1742            if len(idx) > 1: getattr(basesel, 'set_'+meta+'s')([idx[0]])
     1743            sel_lbl += ("%s%d, " % (meta.upper(), idx[0]))
     1744        sel_lbl = sel_lbl.rstrip(', ')
     1745        scan.set_selection(basesel)
     1746        if len(sel_lbl) > 0:
     1747            asaplog.push("Selection contains multiple IFs/Pols/Beams. Plotting the first ones: %s" % sel_lbl)
     1748            asaplog.post("WARN")
     1749        if stype == 't':
     1750            selIds = range(15)
     1751            sellab = "src type "
     1752        else:
     1753            selIds = getattr(scan,'get'+colorby+'nos')()
     1754            sellab = colorby.upper()
     1755        selFunc = "set_"+colorby+"s"
     1756        basesel.set_order(["TIME"])
     1757        # define axes labels
     1758        xlab = self._abcissa or 'Time (UTC)'
     1759        ylab = self._ordinate or scan._get_ordinate_label()
     1760        self._plotter.set_axes('xlabel',xlab)
     1761        self._plotter.set_axes('ylabel',ylab)
     1762        # define the panel title
     1763        if len(sel_lbl) > 0: lbl = sel_lbl
     1764        else: lbl = self._get_label(scan, r, 's', self._title)
     1765        if isinstance(lbl, list) or isinstance(lbl, tuple):
     1766            # get default label
     1767            lbl = self._get_label(scan, r, self._panelling, None)
     1768        self._plotter.set_axes('title',lbl)
     1769        # linestyle
     1770        lstyle = '' if colorby in time_types else ':'
     1771        alldates = []
     1772        for idx in selIds:
     1773            sel = selector() + basesel
     1774            bid = getattr(basesel,'get_'+colorby+"s")()
     1775            if (len(bid) > 0) and (not idx in bid):
     1776                # base selection doesn't contain idx
     1777                # Note summation of selector is logical sum if
     1778                continue
     1779            getattr(sel, selFunc)([idx])
     1780            if not sel.is_empty():
     1781                try:
     1782                    scan.set_selection(sel)
     1783                except RuntimeError, instance:
     1784                    if stype == 't' and str(instance).startswith("Selection contains no data."):
     1785                        continue
     1786                    else:
     1787                        scan.set_selection(basesel)
     1788                        raise RuntimeError, instance
     1789            if scan.nrow() == 0:
     1790                scan.set_selection(basesel)
     1791                continue
     1792            y=array(scan._get_column(scan._getspectrum,-1))
     1793            m = array(scan._get_column(scan._getmask,-1))
     1794            y = ma.masked_array(y,mask=logical_not(array(m,copy=False)))
     1795            # try to handle spectral data somewhat...
     1796            try:
     1797                l,m = y.shape
     1798            except ValueError, e:
     1799                raise ValueError(str(e)+" This error usually occurs when you select multiple spws with different number of channels. Try selecting single spw and retry.")
     1800            if m > 1:
     1801                y=y.mean(axis=1)
     1802            # flag handling
     1803            m = [ scan._is_all_chan_flagged(i) for i in range(scan.nrow()) ]
     1804            y = ma.masked_array(y,mask=m)
     1805            if len(y) == 0: continue
     1806            # line label
     1807            llbl=sellab+str(idx)
     1808            from matplotlib.dates import date2num
     1809            from pytz import timezone
     1810            dates = self._data.get_time(asdatetime=True)
     1811            alldates += list(dates)
     1812            x = date2num(dates)
     1813            tz = timezone('UTC')
     1814            # get color
     1815            lc = self._plotter.colormap[self._plotter.color]
     1816            self._plotter.palette( (self._plotter.color+1) % len(self._plotter.colormap) )
     1817            # actual plotting
     1818            self._plotter.axes.plot_date(x,y,tz=tz,label=llbl,linestyle=lstyle,color=lc,
     1819                                         marker='o',markersize=3,markeredgewidth=0)
     1820
     1821        # legend and axis formatting
     1822        ax = self.gca()
     1823        (dstr, timefmt, majloc, minloc) = self._get_date_axis_setup(alldates, ax.get_xlim())
     1824        ax.xaxis.set_major_formatter(timefmt)
     1825        ax.xaxis.set_major_locator(majloc)
     1826        ax.xaxis.set_minor_locator(minloc)
     1827        self._plotter.axes.legend(loc=self._legendloc)
    16271828
    16281829    def _plottp(self,scan):
     
    16611862        x = arange(len(y))
    16621863        # try to handle spectral data somewhat...
    1663         l,m = y.shape
     1864        try:
     1865            l,m = y.shape
     1866        except ValueError, e:
     1867                raise ValueError(str(e)+" This error usually occurs when you select multiple spws with different number of channels. Try selecting single spw and retry.")
    16641868        if m > 1:
    16651869            y=y.mean(axis=1)
     1870        # flag handling
     1871        m = [ scan._is_all_chan_flagged(i) for i in range(scan.nrow()) ]
     1872        y = ma.masked_array(y,mask=m)
    16661873        plotit = self._plotter.plot
    16671874        llbl = self._get_label(scan, r, self._stacking, None)
     
    17011908            selstr += '\n'
    17021909            self._headtext['selstr'] = selstr
    1703         ssel=(selstr+self._data.get_selection().__str__()+self._selection.__str__() or 'none')
    1704         headstr.append('***Selections***\n'+ssel)
     1910        #ssel=(selstr+self._data.get_selection().__str__()+self._selection.__str__() or 'none')
     1911        curr_selstr = selstr+self._data.get_selection().__str__() or "none"
     1912        ssel=(curr_selstr+"\n" +self._selection.__str__())
     1913        headstr.append('\n\n***Selections***\n'+ssel.replace('$','\$'))
    17051914
    17061915        if plot:
     
    17521961    # plot spectra by pointing
    17531962    @asaplog_post_dec
    1754     def plotgrid(self, scan=None,center=None,spacing=None,rows=None,cols=None):
     1963    def plotgrid(self, scan=None,center="",spacing=[],rows=None,cols=None):
    17551964        """
    17561965        Plot spectra based on direction.
     
    17581967        Parameters:
    17591968            scan:      a scantable to plot
    1760             center:    the grid center direction (a list) of plots in the
    1761                        unit of DIRECTION column.
     1969            center:    the grid center direction (a string)
    17621970                       (default) the center of map region
     1971                       (example) 'J2000 19h30m00s -25d00m00s'
    17631972            spacing:   a list of horizontal (R.A.) and vertical (Dec.)
    1764                        spacing in the unit of DIRECTION column.
     1973                       spacing.
    17651974                       (default) Calculated by the extent of map region and
     1975                       (example) ['1arcmin', '1arcmin']
    17661976                       the number of rows and cols to cover
    17671977            rows:      number of panels (grid points) in horizontal direction
     
    17871997
    17881998        # Rows and cols
    1789         if rows:
    1790             self._rows = int(rows)
    1791         else:
    1792             msg = "Number of rows to plot are not specified. "
    1793             if self._rows:
    1794                 msg += "Using previous value = %d" % (self._rows)
    1795                 asaplog.push(msg)
    1796             else:
    1797                 self._rows = 1
    1798                 msg += "Setting rows = %d" % (self._rows)
    1799                 asaplog.post()
    1800                 asaplog.push(msg)
    1801                 asaplog.post("WARN")
    1802         if cols:
    1803             self._cols = int(cols)
    1804         else:
    1805             msg = "Number of cols to plot are not specified. "
    1806             if self._cols:
    1807                 msg += "Using previous value = %d" % (self._cols)
    1808                 asaplog.push(msg)
    1809             else:
    1810                 self._cols = 1
    1811                 msg += "Setting cols = %d" % (self._cols)
    1812                 asaplog.post()
    1813                 asaplog.push(msg)
    1814                 asaplog.post("WARN")
    1815 
    1816         # Center and spacing
    1817         dirarr = array(self._data.get_directionval()).transpose()
    1818         print "Pointing range: (x, y) = (%f - %f, %f - %f)" %\
    1819               (dirarr[0].min(),dirarr[0].max(),dirarr[1].min(),dirarr[1].max())
    1820         dircent = [0.5*(dirarr[0].max() + dirarr[0].min()),
    1821                    0.5*(dirarr[1].max() + dirarr[1].min())]
    1822         del dirarr
    1823         if center is None:
    1824             #asaplog.post()
    1825             asaplog.push("Grid center is not specified. Automatically calculated from pointing center.")
    1826             #asaplog.post("WARN")
    1827             #center = [dirarr[0].mean(), dirarr[1].mean()]
    1828             center = dircent
    1829         elif (type(center) in (list, tuple)) and len(center) > 1:
    1830             from numpy import pi
    1831             # make sure center_x is in +-pi of pointing center
    1832             # (assumes dirs are in rad)
    1833             rotnum = round(abs(center[0] - dircent[0])/(2*pi))
    1834             if center[0] < dircent[0]: rotnum *= -1
    1835             cenx = center[0] - rotnum*2*pi
    1836             center = [cenx, center[1]]
    1837         else:
    1838             msg = "Direction of grid center should be a list of float (R.A., Dec.)"
    1839             raise ValueError, msg
    1840         asaplog.push("Grid center: (%f, %f) " % (center[0],center[1]))
    1841 
    1842         if spacing is None:
    1843             #asaplog.post()
    1844             asaplog.push("Grid spacing not specified. Automatically calculated from map coverage")
    1845             #asaplog.post("WARN")
    1846             # automatically get spacing
    1847             dirarr = array(self._data.get_directionval()).transpose()
    1848             wx = 2. * max(abs(dirarr[0].max()-center[0]),
    1849                           abs(dirarr[0].min()-center[0]))
    1850             wy = 2. * max(abs(dirarr[1].max()-center[1]),
    1851                           abs(dirarr[1].min()-center[1]))
    1852             ## slightly expand area to plot the edges
    1853             #wx *= 1.1
    1854             #wy *= 1.1
    1855             xgrid = wx/max(self._cols-1.,1.)
    1856             #xgrid = wx/float(max(self._cols,1.))
    1857             xgrid *= cos(center[1])
    1858             ygrid = wy/max(self._rows-1.,1.)
    1859             #ygrid = wy/float(max(self._rows,1.))
    1860             # single pointing (identical R.A. and/or Dec. for all spectra.)
    1861             if xgrid == 0:
    1862                 xgrid = 1.
    1863             if ygrid == 0:
    1864                 ygrid = 1.
    1865             # spacing should be negative to transpose plot
    1866             spacing = [- xgrid, - ygrid]
    1867             del dirarr, xgrid, ygrid
    1868         #elif isinstance(spacing, str):
    1869         #    # spacing is a quantity
    1870         elif (type(spacing) in (list, tuple)) and len(spacing) > 1:
    1871             for i in xrange(2):
    1872                 val = spacing[i]
    1873                 if not isinstance(val, float):
    1874                     raise TypeError("spacing should be a list of float")
    1875                 if val > 0.:
    1876                     spacing[i] = -val
    1877             spacing = spacing[0:2]
    1878         else:
    1879             msg = "Invalid spacing."
    1880             raise TypeError(msg)
    1881         asaplog.push("Spacing: (%f, %f) (projected)" % (spacing[0],spacing[1]))
    1882 
     1999        if (self._rows is None):
     2000            rows = max(1, rows)
     2001        if (self._cols is None):
     2002            cols = max(1, cols)
     2003        self.set_layout(rows,cols,False)
     2004
     2005        # Select the first IF, POL, and BEAM for plotting
    18832006        ntotpl = self._rows * self._cols
    18842007        ifs = self._data.getifnos()
     
    19072030            asaplog.push(msg)
    19082031            asaplog.post("WARN")
    1909        
     2032
     2033        # Prepare plotter
    19102034        self._assert_plotter(action="reload")
    19112035        self._plotter.hold()
     
    19232047        from asap._asap import plothelper as plhelper
    19242048        ph = plhelper(self._data)
    1925         ph.set_gridval(self._cols, self._rows, spacing[0], spacing[1],
    1926                           center[0], center[1], epoch="J2000", projname="SIN")
     2049        #ph.set_gridval(self._cols, self._rows, spacing[0], spacing[1],
     2050        #                  center[0], center[1], epoch="J2000", projname="SIN")
     2051        if type(spacing) in (list, tuple, array):
     2052            if len(spacing) == 0:
     2053                spacing = ["", ""]
     2054            elif len(spacing) == 1:
     2055                spacing = [spacing[0], spacing[0]]
     2056        else:
     2057            spacing = [spacing, spacing]
     2058        ph.set_grid(self._cols, self._rows, spacing[0], spacing[1], \
     2059                    center, projname="SIN")
     2060
    19272061        # Actual plot
    19282062        npl = 0
  • trunk/python/customgui_base.py

    r2704 r3112  
    11import os
     2import weakref
    23import matplotlib, numpy
    34from asap.logging import asaplog, asaplog_post_dec
     
    1314class CustomToolbarCommon:
    1415    def __init__(self,parent):
    15         self.plotter = parent
     16        self.plotter = weakref.ref(parent)
    1617        #self.figmgr=self.plotter._plotter.figmgr
     18
     19    def _get_plotter(self):
     20        # check for the validity of the plotter and
     21        # return the plotter class instance if its valid.
     22        if self.plotter() is None:
     23            raise RuntimeError, "Internal Error. The plotter has been destroyed."
     24        else:
     25            return self.plotter()
    1726
    1827    ### select the nearest spectrum in pick radius
     
    2938        if event.button != 1:
    3039            return
     40
     41        # check for the validity of plotter and get the plotter
     42        theplotter = self._get_plotter()
    3143
    3244        xclick = event.xdata
     
    6274        del pind, inds, xlin, ylin
    6375        # Spectra are Picked
    64         theplot = self.plotter._plotter
     76        theplot = theplotter._plotter
    6577        thetoolbar = self.figmgr.toolbar
    6678        thecanvas = self.figmgr.canvas
     
    154166            return
    155167
     168        # check for the validity of plotter and get the plotter
     169        theplotter = self._get_plotter()
     170
    156171        self._thisregion = {'axes': event.inaxes,'xs': event.x,
    157172                            'worldx': [event.xdata,event.xdata],
     
    160175        self.xdataold = event.xdata
    161176
    162         self.plotter._plotter.register('button_press',None)
    163         self.plotter._plotter.register('motion_notify', self._xspan_draw)
    164         self.plotter._plotter.register('button_press', self._xspan_end)
     177        theplotter._plotter.register('button_press',None)
     178        theplotter._plotter.register('motion_notify', self._xspan_draw)
     179        theplotter._plotter.register('button_press', self._xspan_end)
    165180
    166181    def _xspan_draw(self,event):
     
    203218            xdataend = self.xdataold
    204219
     220        # check for the validity of plotter and get the plotter
     221        theplotter = self._get_plotter()
     222
    205223        self._thisregion['worldx'][1] = xdataend
    206224        # print statistics of spectra in subplot
     
    208226       
    209227        # release event
    210         self.plotter._plotter.register('button_press',None)
    211         self.plotter._plotter.register('motion_notify',None)
     228        theplotter._plotter.register('button_press',None)
     229        theplotter._plotter.register('motion_notify',None)
    212230        # Clear up region selection
    213231        self._thisregion = None
     
    215233        self.xold = None
    216234        # finally recover region selection event
    217         self.plotter._plotter.register('button_press',self._single_mask)
     235        theplotter._plotter.register('button_press',self._single_mask)
    218236
    219237    def _subplot_stats(self,selection):
     
    321339    ### actual plotting of the new page
    322340    def _new_page(self,goback=False):
     341        # check for the validity of plotter and get the plotter
     342        theplotter = self._get_plotter()
     343
    323344        top = None
    324         header = self.plotter._headtext
     345        header = theplotter._headtext
    325346        reset = False
    326347        doheader = (isinstance(header['textobj'],list) and \
    327348                    len(header['textobj']) > 0)
    328349        if doheader:
    329             top = self.plotter._plotter.figure.subplotpars.top
     350            top = theplotter._plotter.figure.subplotpars.top
    330351            fontsize = header['textobj'][0].get_fontproperties().get_size()
    331         if self.plotter._startrow <= 0:
     352        if theplotter._startrow <= 0:
    332353            msg = "The page counter is reset due to chages of plot settings. "
    333354            msg += "Plotting from the first page."
     
    342363                if header.has_key('selstr'):
    343364                    selstr = header['selstr']
    344             self.plotter._reset_header()
    345 
    346         self.plotter._plotter.hold()
     365            theplotter._reset_header()
     366
     367        theplotter._plotter.hold()
    347368        if goback:
    348369            self._set_prevpage_counter()
    349         #self.plotter._plotter.clear()
    350         self.plotter._plot(self.plotter._data)
     370        #theplotter._plotter.clear()
     371        theplotter._plot(theplotter._data)
    351372        pagenum = self._get_pagenum()
    352373        self.set_pagecounter(pagenum)
     
    354375        #if header['textobj']:
    355376        if doheader and pagenum == 1:
    356             if top and top != self.plotter._margins[3]:
     377            if top and top != theplotter._margins[3]:
    357378                # work around for sdplot in CASA. complete checking in future?
    358                 self.plotter._plotter.figure.subplots_adjust(top=top)
     379                theplotter._plotter.figure.subplots_adjust(top=top)
    359380            if reset:
    360                 self.plotter.print_header(plot=True,fontsize=fontsize,selstr=selstr, extrastr=extrastr)
     381                theplotter.print_header(plot=True,fontsize=fontsize,selstr=selstr, extrastr=extrastr)
    361382            else:
    362                 self.plotter._header_plot(header['string'],fontsize=fontsize)
    363         self.plotter._plotter.release()
    364         self.plotter._plotter.tidy()
    365         self.plotter._plotter.show(hardrefresh=False)
     383                theplotter._header_plot(header['string'],fontsize=fontsize)
     384        theplotter._plotter.release()
     385        theplotter._plotter.tidy()
     386        theplotter._plotter.show(hardrefresh=False)
    366387        del top
    367388
    368389    ### calculate the panel ID and start row to plot the previous page
    369390    def _set_prevpage_counter(self):
     391        # check for the validity of plotter and get the plotter
     392        theplotter = self._get_plotter()
     393
    370394        # set row and panel counters to those of the 1st panel of previous page
    371395        maxpanel = 16
    372396        # the ID of the last panel in current plot
    373         lastpanel = self.plotter._ipanel
     397        lastpanel = theplotter._ipanel
    374398        # the number of current subplots
    375         currpnum = len(self.plotter._plotter.subplots)
     399        currpnum = len(theplotter._plotter.subplots)
    376400        # the nuber of previous subplots
    377401        prevpnum = None
    378         if self.plotter._rows and self.plotter._cols:
     402        if theplotter._rows and theplotter._cols:
    379403            # when user set layout
    380             prevpnum = self.plotter._rows*self.plotter._cols
     404            prevpnum = theplotter._rows*theplotter._cols
    381405        else:
    382406            # no user specification
     
    385409        start_ipanel = max(lastpanel-currpnum-prevpnum+1, 0)
    386410        # set the pannel ID of the last panel of prev-prev page
    387         self.plotter._ipanel = start_ipanel-1
    388         if self.plotter._panelling == 'r':
    389             self.plotter._startrow = start_ipanel
     411        theplotter._ipanel = start_ipanel-1
     412        if theplotter._panelling == 'r':
     413            theplotter._startrow = start_ipanel
    390414        else:
    391415            # the start row number of the next panel
    392             self.plotter._startrow = self.plotter._panelrows[start_ipanel]
     416            theplotter._startrow = theplotter._panelrows[start_ipanel]
    393417        del lastpanel,currpnum,prevpnum,start_ipanel
    394418
     
    405429
    406430    def _get_pagenum(self):
     431        # check for the validity of plotter and get the plotter
     432        theplotter = self._get_plotter()
     433       
    407434        # get the ID of last panel in the current page
    408         idlastpanel = self.plotter._ipanel
     435        idlastpanel = theplotter._ipanel
    409436        # max panels in a page
    410         ppp = self.plotter._plotter.rows*self.plotter._plotter.cols
     437        ppp = theplotter._plotter.rows*theplotter._plotter.cols
    411438        return int(idlastpanel/ppp)+1
    412439
     
    683710class CustomFlagToolbarCommon:
    684711    def __init__(self,parent):
    685         self.plotter=parent
     712        self.plotter=weakref.ref(parent)
    686713        #self.figmgr=self.plotter._plotter.figmgr
    687714        self._selregions = {}
     
    691718        self.xdataold=None
    692719
     720    def _get_plotter(self):
     721        # check for the validity of the plotter and
     722        # return the plotter class instance if its valid.
     723        if self.plotter() is None:
     724            raise RuntimeError, "Internal Error. The plotter has been destroyed."
     725        else:
     726            return self.plotter()
     727
    693728    ### select the nearest spectrum in pick radius
    694729    ###    and display spectral value on the toolbar.
     
    704739        if event.button != 1:
    705740            return
     741
     742        # check for the validity of plotter and get the plotter
     743        theplotter = self._get_plotter()
    706744
    707745        xclick = event.xdata
     
    737775        del pind, inds, xlin, ylin
    738776        # Spectra are Picked
    739         theplot = self.plotter._plotter
     777        theplot = theplotter._plotter
    740778        thetoolbar = self.figmgr.toolbar
    741779        thecanvas = self.figmgr.canvas
     
    819857        if event.button != 1 or event.inaxes == None:
    820858            return
     859        # check for the validity of plotter and get the plotter
     860        theplotter = self._get_plotter()
    821861        # this row resolution assumes row panelling
    822862        irow = int(self._getrownum(event.inaxes))
     
    829869        self._thisregion = {'axes': event.inaxes,'xs': event.x,
    830870                            'worldx': [event.xdata,event.xdata]}
    831         self.plotter._plotter.register('button_press',None)
     871        theplotter._plotter.register('button_press',None)
    832872        self.xold = event.x
    833873        self.xdataold = event.xdata
    834         self.plotter._plotter.register('motion_notify', self._xspan_draw)
    835         self.plotter._plotter.register('button_press', self._xspan_end)
     874        theplotter._plotter.register('motion_notify', self._xspan_draw)
     875        theplotter._plotter.register('button_press', self._xspan_end)
    836876
    837877    def _xspan_draw(self,event):
     
    882922        self._thisregion['axes'].set_xlim(axlimx)
    883923       
    884         self.plotter._plotter.canvas.draw()
     924        # check for the validity of plotter and get the plotter
     925        theplotter = self._get_plotter()
     926
     927        theplotter._plotter.canvas.draw()
    885928        self._polygons.append(pregion)
    886929        srow = self._getrownum(self._thisregion['axes'])
     
    895938
    896939        # release event
    897         self.plotter._plotter.register('button_press',None)
    898         self.plotter._plotter.register('motion_notify',None)
     940        theplotter._plotter.register('button_press',None)
     941        theplotter._plotter.register('motion_notify',None)
    899942        # Clear up region selection
    900943        self._thisregion = None
     
    902945        self.xold = None
    903946        # finally recover region selection event
    904         self.plotter._plotter.register('button_press',self._add_region)
     947        theplotter._plotter.register('button_press',self._add_region)
    905948
    906949    ### add panels to selections
     
    911954        if event.button != 1 or event.inaxes == None:
    912955            return
     956        # check for the validity of plotter and get the plotter
     957        theplotter = self._get_plotter()
     958
    913959        selax = event.inaxes
    914960        # this row resolution assumes row panelling
     
    919965        shadow = Rectangle((0,0),1,1,facecolor='0.7',transform=selax.transAxes,visible=True)
    920966        self._polygons.append(selax.add_patch(shadow))
    921         #self.plotter._plotter.show(False)
    922         self.plotter._plotter.canvas.draw()
     967        #theplotter._plotter.show(False)
     968        theplotter._plotter.canvas.draw()
    923969        asaplog.push("row "+str(irow)+" is selected")
    924970        ## check for region selection of the spectra and overwrite it.
     
    9561002            asaplog.push("Invalid panel specification")
    9571003            asaplog.post('ERROR')
    958         strow = self._getrownum(self.plotter._plotter.subplots[0]['axes'])
    959         enrow = self._getrownum(self.plotter._plotter.subplots[-1]['axes'])
     1004
     1005        # check for the validity of plotter and get the plotter
     1006        theplotter = self._get_plotter()
     1007
     1008        strow = self._getrownum(theplotter._plotter.subplots[0]['axes'])
     1009        enrow = self._getrownum(theplotter._plotter.subplots[-1]['axes'])
    9601010        for irow in range(int(strow),int(enrow)+1):
    9611011            if regions.has_key(str(irow)):
    962                 ax = self.plotter._plotter.subplots[irow - int(strow)]['axes']
     1012                ax = theplotter._plotter.subplots[irow - int(strow)]['axes']
    9631013                mlist = regions.pop(str(irow))
    9641014                # WORKAROUND for the issue axvspan started to reset xlim.
     
    9701020                del ax,mlist,axlimx
    9711021            if irow in panels:
    972                 ax = self.plotter._plotter.subplots[irow - int(strow)]['axes']
     1022                ax = theplotter._plotter.subplots[irow - int(strow)]['axes']
    9731023                shadow = Rectangle((0,0),1,1,facecolor='0.7',
    9741024                                   transform=ax.transAxes,visible=True)
    9751025                self._polygons.append(ax.add_patch(shadow))
    9761026                del ax,shadow
    977         self.plotter._plotter.canvas.draw()
     1027        theplotter._plotter.canvas.draw()
    9781028        del regions,panels,strow,enrow
    9791029
     
    9831033            for shadow in self._polygons:
    9841034                shadow.remove()
    985             if refresh: self.plotter._plotter.canvas.draw()
     1035            if refresh:
     1036                # check for the validity of plotter and get the plotter
     1037                theplotter = self._get_plotter()
     1038                theplotter._plotter.canvas.draw()
    9861039        self._polygons = []
    9871040
     
    10071060            asaplog.post('WARN')
    10081061            return
     1062
    10091063        self._pause_buttons(operation="start",msg="Flagging data...")
    10101064        self._flag_operation(rows=self._selpanels,
     
    10151069        asaplog.push(sout)
    10161070        del sout
    1017         self.plotter._ismodified = True
     1071        # check for the validity of plotter and get the plotter
     1072        theplotter = self._get_plotter()
     1073
     1074        theplotter._ismodified = True
    10181075        self._clearup_selections(refresh=False)
    10191076        self._plot_page(pagemode="current")
     
    10361093        asaplog.push(sout)
    10371094        del sout
    1038         self.plotter._ismodified = True
     1095
     1096        # check for the validity of plotter and get the plotter
     1097        theplotter = self._get_plotter()
     1098        theplotter._ismodified = True
    10391099        self._clearup_selections(refresh=False)
    10401100        self._plot_page(pagemode="current")
     
    10441104    @asaplog_post_dec
    10451105    def _flag_operation(self,rows=None,regions=None,unflag=False):
    1046         scan = self.plotter._data
     1106        # check for the validity of plotter and get the plotter
     1107        theplotter = self._get_plotter()
     1108
     1109        scan = theplotter._data
    10471110        if not scan:
    10481111            asaplog.post()
     
    10791142    @asaplog_post_dec
    10801143    def _selected_stats(self,rows=None,regions=None):
    1081         scan = self.plotter._data
     1144        # check for the validity of plotter and get the plotter
     1145        theplotter = self._get_plotter()
     1146
     1147        scan = theplotter._data
    10821148        if not scan:
    10831149            asaplog.post()
     
    11641230    ### actual plotting of the new page
    11651231    def _plot_page(self,pagemode="next"):
    1166         if self.plotter._startrow <= 0:
     1232        # check for the validity of plotter and get the plotter
     1233        theplotter = self._get_plotter()
     1234        if theplotter._startrow <= 0:
    11671235            msg = "The page counter is reset due to chages of plot settings. "
    11681236            msg += "Plotting from the first page."
     
    11721240            goback = False
    11731241
    1174         self.plotter._plotter.hold()
    1175         #self.plotter._plotter.legend(1)
     1242        theplotter._plotter.hold()
     1243        #theplotter._plotter.legend(1)
    11761244        self._set_plot_counter(pagemode)
    1177         self.plotter._plot(self.plotter._data)
     1245        theplotter._plot(theplotter._data)
    11781246        self.set_pagecounter(self._get_pagenum())
    1179         self.plotter._plotter.release()
    1180         self.plotter._plotter.tidy()
    1181         self.plotter._plotter.show(hardrefresh=False)
     1247        theplotter._plotter.release()
     1248        theplotter._plotter.tidy()
     1249        theplotter._plotter.show(hardrefresh=False)
    11821250
    11831251    ### calculate the panel ID and start row to plot a page
     
    11941262            # nothing necessary to plot the next page
    11951263            return
     1264
     1265        # check for the validity of plotter and get the plotter
     1266        theplotter = self._get_plotter()
     1267
    11961268        # set row and panel counters to those of the 1st panel of previous page
    11971269        maxpanel = 25
    11981270        # the ID of the last panel in current plot
    1199         lastpanel = self.plotter._ipanel
     1271        lastpanel = theplotter._ipanel
    12001272        # the number of current subplots
    1201         currpnum = len(self.plotter._plotter.subplots)
     1273        currpnum = len(theplotter._plotter.subplots)
    12021274
    12031275        # the nuber of previous subplots
     
    12081280            ## previous page
    12091281            prevpnum = None
    1210             if self.plotter._rows and self.plotter._cols:
     1282            if theplotter._rows and theplotter._cols:
    12111283                # when user set layout
    1212                 prevpnum = self.plotter._rows*self.plotter._cols
     1284                prevpnum = theplotter._rows*theplotter._cols
    12131285            else:
    12141286                # no user specification
     
    12181290
    12191291        # set the pannel ID of the last panel of the prev(-prev) page
    1220         self.plotter._ipanel = start_ipanel-1
    1221         if self.plotter._panelling == 'r':
    1222             self.plotter._startrow = start_ipanel
     1292        theplotter._ipanel = start_ipanel-1
     1293        if theplotter._panelling == 'r':
     1294            theplotter._startrow = start_ipanel
    12231295        else:
    12241296            # the start row number of the next panel
    1225             self.plotter._startrow = self.plotter._panelrows[start_ipanel]
     1297            theplotter._startrow = theplotter._panelrows[start_ipanel]
    12261298        del lastpanel,currpnum,start_ipanel
    12271299
     
    12381310
    12391311    def _get_pagenum(self):
     1312        # check for the validity of plotter and get the plotter
     1313        theplotter = self._get_plotter()
    12401314        # get the ID of last panel in the current page
    1241         idlastpanel = self.plotter._ipanel
     1315        idlastpanel = theplotter._ipanel
    12421316        # max panels in a page
    1243         ppp = self.plotter._plotter.rows*self.plotter._plotter.cols
     1317        ppp = theplotter._plotter.rows*theplotter._plotter.cols
    12441318        return int(idlastpanel/ppp)+1
    12451319
  • trunk/python/env.py

    r2704 r3112  
    22"""
    33__all__ = ["is_casapy", "is_ipython", "setup_env", "get_revision",
    4            "is_asap_cli"]
     4           "is_asap_cli", "get_asap_revdate"]
    55
    66import sys
     
    5252        os.environ["CASAPATH"] = "%s %s somwhere" % ( asapdata, plf)
    5353
    54 def get_revision():
     54def get_revinfo_file():
    5555    """Get the revision of the software. Only useful within casapy."""
    5656    if not is_casapy:
     
    5959    versioninfo = sys.version_info
    6060    pyversion = '%s.%s'%(versioninfo[0],versioninfo[1])
     61    revinfo = None
    6162    if os.path.isdir(casapath[0]+'/'+casapath[1]+'/python/%s/asap'%(pyversion)):
    6263        # for casa developer environment (linux or darwin)
     
    6869        else:
    6970            revinfo=casapath[0]+'/lib/python%s/asap/svninfo.txt'%(pyversion)
     71    return revinfo
     72
     73def get_revision():
     74    revinfo=get_revinfo_file()
    7075    if os.path.isfile(revinfo):
    7176        f = file(revinfo)
     
    7580        return revsionno.rstrip()
    7681    return ' unknown '
     82
     83
     84def get_asap_revdate():
     85    revinfo=get_revinfo_file()
     86    if os.path.isfile(revinfo):
     87        f = file(revinfo)
     88        f.readline()
     89        f.readline()
     90        revdate=f.readline()
     91        return revdate.rstrip().lstrip()
     92    return 'unknown'
     93
     94
  • trunk/python/sbseparator.py

    r2704 r3112  
    99from asap.selector import selector
    1010from asap.asapgrid import asapgrid2
    11 #from asap._asap import sidebandsep
     11from asap._asap import SBSeparator
    1212
    1313class sbseparator:
     
    2121      information in scantable in sideband sparation. Frequency
    2222      setting of SIGNAL sideband is stored in output table for now.
    23     - Flag information (stored in FLAGTRA) is ignored.
    2423
    2524    Example:
     
    3837    """
    3938    def __init__(self, infiles):
    40         self.intables = None
    41         self.signalShift = []
    42         self.imageShift = []
    43         self.dsbmode = True
    44         self.getboth = False
    45         self.rejlimit = 0.2
    46         self.baseif = -1
    47         self.freqtol = 10.
    48         self.freqframe = ""
    49         self.solveother = False
    50         self.dirtol = [1.e-5, 1.e-5] # direction tolerance in rad (2 arcsec)
    51 
    52         self.tables = []
    53         self.nshift = -1
    54         self.nchan = -1
    55 
    56         self.set_data(infiles)
    57        
    58         #self.separator = sidebandsep()
    59 
    60     @asaplog_post_dec
    61     def set_data(self, infiles):
    62         """
    63         Set data to be processed.
    64 
    65         infiles  : a list of filenames or scantables
    66         """
    67         if not (type(infiles) in (list, tuple, numpy.ndarray)):
    68             infiles = [infiles]
    69         if isinstance(infiles[0], scantable):
    70             # a list of scantable
    71             for stab in infiles:
    72                 if not isinstance(stab, scantable):
    73                     asaplog.post()
    74                     raise TypeError, "Input data is not a list of scantables."
    75 
    76             #self.separator._setdata(infiles)
    77             self._reset_data()
    78             self.intables = infiles
    79         else:
    80             # a list of filenames
    81             for name in infiles:
    82                 if not os.path.exists(name):
    83                     asaplog.post()
    84                     raise ValueError, "Could not find input file '%s'" % name
    85            
    86             #self.separator._setdataname(infiles)
    87             self._reset_data()
    88             self.intables = infiles
    89 
    90         asaplog.push("%d files are set to process" % len(self.intables))
     39        self._separator = SBSeparator(infiles)
    9140
    9241
    93     def _reset_data(self):
    94         del self.intables
    95         self.intables = None
    96         self.signalShift = []
    97         #self.imageShift = []
    98         self.tables = []
    99         self.nshift = -1
    100         self.nchan = -1
    101 
    102     @asaplog_post_dec
    10342    def set_frequency(self, baseif, freqtol, frame=""):
    10443        """
     
    10746        Parameters:
    10847          - reference IFNO to process in the first table in the list
    109           - frequency tolerance from reference IF to select data
     48          - frequency tolerance from reference IF to select data (string)
    11049          frame  : frequency frame to select IF
    11150        """
    112         self._reset_if()
    113         self.baseif = baseif
    114         if isinstance(freqtol,dict) and freqtol["unit"] == "Hz":
    115             if freqtol['value'] > 0.:
    116                 self.freqtol = freqtol
    117             else:
    118                 asaplog.post()
    119                 asaplog.push("Frequency tolerance should be positive value.")
    120                 asaplog.post("ERROR")
    121                 return
    122         else:
    123             # torelance in channel unit
    124             if freqtol > 0:
    125                 self.freqtol = float(freqtol)
    126             else:
    127                 asaplog.post()
    128                 asaplog.push("Frequency tolerance should be positive value.")
    129                 asaplog.post("ERROR")
    130                 return
    131         self.freqframe = frame
     51        if type(freqtol) in (float, int):
     52            freqtol = str(freqtol)
     53        elif isinstance(freqtol, dict):
     54            try:
     55                freqtol = str(freqtol['value']) + freqtol['unit']
     56            except:
     57                raise ValueError, "Invalid frequency tolerance."
     58        self._separator.set_freq(baseif, freqtol, frame)
    13259
    133     def _reset_if(self):
    134         self.baseif = -1
    135         self.freqtol = 0
    136         self.freqframe = ""
    137         self.signalShift = []
    138         #self.imageShift = []
    139         self.tables = []
    140         self.nshift = 0
    141         self.nchan = -1
    14260
    143     @asaplog_post_dec
    144     def set_dirtol(self, dirtol=[1.e-5,1.e-5]):
     61    def set_dirtol(self, dirtol=["2arcsec", "2arcsec"]):
    14562        """
    14663        Set tolerance of direction to select data
    14764        """
    148         # direction tolerance in rad
    149         if not (type(dirtol) in [list, tuple, numpy.ndarray]):
    150             dirtol = [dirtol, dirtol]
    151         if len(dirtol) == 1:
    152             dirtol = [dirtol[0], dirtol[0]]
    153         if len(dirtol) > 1:
    154             self.dirtol = dirtol[0:2]
    155         else:
    156             asaplog.post()
    157             asaplog.push("Invalid direction tolerance. Should be a list of float in unit radian")
    158             asaplog.post("ERROR")
    159             return
    160         asaplog.post("Set direction tolerance [%f, %f] (rad)" % \
    161                      (self.dirtol[0], self.dirtol[1]))
     65        if isinstance(dirtol, str):
     66            dirtol = [dirtol]
    16267
    163     @asaplog_post_dec
    164     def set_shift(self, mode="DSB", imageshift=None):
     68        self._separator.set_dirtol(dirtol)
     69   
     70           
     71    def set_shift(self, imageshift=[]):
    16572        """
    16673        Set shift mode and channel shift of image band.
    16774
    168         mode       : shift mode ['DSB'|'SSB']
    169                      When mode='DSB', imageshift is assumed to be equal
    170                      to the shift of signal sideband but in opposite direction.
    171         imageshift : a list of number of channel shift in image band of
    172                      each scantable. valid only mode='SSB'
     75        imageshift : a list of number of channels shifted in image
     76                     side band of each scantable.
     77                     If the shifts are not set, they are assumed to be
     78                     equal to those of signal side band, but in opposite
     79                     direction as usual by LO1 offsetting of DSB receivers.
    17380        """
    174         if mode.upper().startswith("S"):
    175             if not imageshift:
    176                 raise ValueError, "Need to set shift value of image sideband"
    177             self.dsbmode = False
    178             self.imageShift = imageshift
    179             asaplog.push("Image sideband shift is set manually: %s" % str(self.imageShift))
    180         else:
    181             # DSB mode
    182             self.dsbmode = True
    183             self.imageShift = []
     81        if not imageshift:
     82            imageshift = []
     83        self._separator.set_shift(imageshift)
     84
    18485
    18586    @asaplog_post_dec
     
    18889        Resolve both image and signal sideband when True.
    18990        """
    190         self.getboth = flag
    191         if self.getboth:
    192             asaplog.push("Both signal and image sidebands are solved and output as separate tables.")
     91        self._separator.solve_both(flag)
     92        if flag:
     93            asaplog.push("Both signal and image sidebands will be solved and stored in separate tables.")
    19394        else:
    194             asaplog.push("Only signal sideband is solved and output as an table.")
     95            asaplog.push("Only signal sideband will be solved and stored in an table.")
    19596
    19697    @asaplog_post_dec
     
    199100        Set rejection limit of solution.
    200101        """
    201         #self.separator._setlimit(abs(threshold))
    202         self.rejlimit = threshold
    203         asaplog.push("The threshold of rejection is set to %f" % self.rejlimit)
     102        self._separator.set_limit(threshold)
    204103
    205104
     
    210109        when True.
    211110        """
    212         self.solveother = flag
     111        self._separator.subtract_other(flag)
    213112        if flag:
    214113            asaplog.push("Expert mode: solution are obtained by subtraction of the other sideband.")
    215114
     115
     116    def set_lo1(self, lo1, frame="TOPO", reftime=-1, refdir=""):
     117        """
     118        Set LO1 frequency to calculate frequency of image sideband.
     119
     120        lo1     : LO1 frequency
     121        frame   : the frequency frame of LO1
     122        reftime : the reference time used in frequency frame conversion.
     123        refdir  : the reference direction used in frequency frame conversion.
     124        """
     125        self._separator.set_lo1(lo1, frame, reftime, refdir)
     126
     127
     128    def set_lo1root(self, name):
     129        """
     130        Set MS name which stores LO1 frequency of signal side band.
     131        It is used to calculate frequency of image sideband.
     132
     133        name : MS name which contains 'ASDM_SPECTRALWINDOW' and
     134               'ASDM_RECEIVER' tables.
     135        """
     136        self._separator.set_lo1root(name)
    216137
    217138    @asaplog_post_dec
     
    223144        overwrite : overwrite existing table
    224145        """
    225         # List up valid scantables and IFNOs to convolve.
    226         #self.separator._separate()
    227         self._setup_shift()
    228         #self._preprocess_tables()
    229 
    230         nshift = len(self.tables)
    231         signaltab = self._grid_outtable(self.tables[0].copy())
    232         if self.getboth:
    233             imagetab = signaltab.copy()
    234 
    235         rejrow = []
    236         for irow in xrange(signaltab.nrow()):
    237             currpol = signaltab.getpol(irow)
    238             currbeam = signaltab.getbeam(irow)
    239             currdir = signaltab.get_directionval(irow)
    240             spec_array, tabidx = self._get_specarray(polid=currpol,\
    241                                                      beamid=currbeam,\
    242                                                      dir=currdir)
    243             #if not spec_array:
    244             if len(tabidx) == 0:
    245                 asaplog.post()
    246                 asaplog.push("skipping row %d" % irow)
    247                 rejrow.append(irow)
    248                 continue
    249             signal = self._solve_signal(spec_array, tabidx)
    250             signaltab.set_spectrum(signal, irow)
    251             if self.getboth:
    252                 image = self._solve_image(spec_array, tabidx)
    253                 imagetab.set_spectrum(image, irow)
    254        
    255         # TODO: Need to remove rejrow form scantables here
    256         signaltab.flag_row(rejrow)
    257         if self.getboth:
    258             imagetab.flag_row(rejrow)
    259        
    260         if outname == "":
    261             outname = "sbsepareted.asap"
    262         signalname = outname + ".signalband"
    263         if os.path.exists(signalname):
    264             if not overwrite:
    265                 raise Exception, "Output file '%s' exists." % signalname
    266             else:
    267                 shutil.rmtree(signalname)
    268         signaltab.save(signalname)
    269         if self.getboth:
    270             # Warnings
     146        out_default = "sbseparated.asap"
     147        if len(outname) == 0:
     148            outname = out_default
    271149            asaplog.post()
    272             asaplog.push("Saving IMAGE sideband.")
    273             asaplog.push("Note, frequency information of IMAGE sideband cannot be properly filled so far. (future development)")
    274             asaplog.push("Storing frequency setting of SIGNAL sideband in FREQUENCIES table for now.")
     150            asaplog.push("The output file name is not specified.")
     151            asaplog.push("Using default name '%s'" % outname)
    275152            asaplog.post("WARN")
    276153
    277             imagename = outname + ".imageband"
    278             if os.path.exists(imagename):
    279                 if not overwrite:
    280                     raise Exception, "Output file '%s' exists." % imagename
    281                 else:
    282                     shutil.rmtree(imagename)
    283             imagetab.save(imagename)
     154        if os.path.exists(outname):
     155            if overwrite:
     156                asaplog.push("removing the old file '%s'" % outname)
     157                shutil.rmtree(outname)
     158            else:
     159                asaplog.post()
     160                asaplog.push("Output file '%s' already exists." % outname)
     161                asaplog.post("ERROR")
     162                return False
    284163
     164        self._separator.separate(outname)
    285165
    286     def _solve_signal(self, data, tabidx=None):
    287         if not tabidx:
    288             tabidx = range(len(data))
    289 
    290         tempshift = []
    291         dshift = []
    292         if self.solveother:
    293             for idx in tabidx:
    294                 tempshift.append(-self.imageShift[idx])
    295                 dshift.append(self.signalShift[idx] - self.imageShift[idx])
    296         else:
    297             for idx in tabidx:
    298                 tempshift.append(-self.signalShift[idx])
    299                 dshift.append(self.imageShift[idx] - self.signalShift[idx])
    300 
    301         shiftdata = numpy.zeros(data.shape, numpy.float)
    302         for i in range(len(data)):
    303             shiftdata[i] = self._shiftSpectrum(data[i], tempshift[i])
    304         ifftdata = self._Deconvolution(shiftdata, dshift, self.rejlimit)
    305         result_image = self._combineResult(ifftdata)
    306         if not self.solveother:
    307             return result_image
    308         result_signal = self._subtractOtherSide(shiftdata, dshift, result_image)
    309         return result_signal
    310 
    311 
    312     def _solve_image(self, data, tabidx=None):
    313         if not tabidx:
    314             tabidx = range(len(data))
    315 
    316         tempshift = []
    317         dshift = []
    318         if self.solveother:
    319             for idx in tabidx:
    320                 tempshift.append(-self.signalShift[idx])
    321                 dshift.append(self.imageShift[idx] - self.signalShift[idx])
    322         else:
    323             for idx in tabidx:
    324                 tempshift.append(-self.imageShift[idx])
    325                 dshift.append(self.signalShift[idx] - self.imageShift[idx])
    326 
    327         shiftdata = numpy.zeros(data.shape, numpy.float)
    328         for i in range(len(data)):
    329             shiftdata[i] = self._shiftSpectrum(data[i], tempshift[i])
    330         ifftdata = self._Deconvolution(shiftdata, dshift, self.rejlimit)
    331         result_image = self._combineResult(ifftdata)
    332         if not self.solveother:
    333             return result_image
    334         result_signal = self._subtractOtherSide(shiftdata, dshift, result_image)
    335         return result_signal
    336 
    337     @asaplog_post_dec
    338     def _grid_outtable(self, table):
    339         # Generate gridded table for output (Just to get rows)
    340         gridder = asapgrid2(table)
    341         gridder.setIF(self.baseif)
    342        
    343         cellx = str(self.dirtol[0])+"rad"
    344         celly = str(self.dirtol[1])+"rad"
    345         dirarr = numpy.array(table.get_directionval()).transpose()
    346         mapx = dirarr[0].max() - dirarr[0].min()
    347         mapy = dirarr[1].max() - dirarr[1].min()
    348         nx = max(1, numpy.ceil(mapx/self.dirtol[0]))
    349         ny = max(1, numpy.ceil(mapy/self.dirtol[0]))
    350 
    351         asaplog.push("Regrid output scantable with cell = [%s, %s]" % \
    352                      (cellx, celly))
    353         gridder.defineImage(nx=nx, ny=ny, cellx=cellx, celly=celly)
    354         gridder.setFunc(func='box', convsupport=1)
    355         gridder.setWeight(weightType='uniform')
    356         gridder.grid()
    357         return gridder.getResult()
    358 
    359     @asaplog_post_dec
    360     def _get_specarray(self, polid=None, beamid=None, dir=None):
    361         ntable = len(self.tables)
    362         spec_array = numpy.zeros((ntable, self.nchan), numpy.float)
    363         nspec = 0
    364         asaplog.push("Start data selection by POL=%d, BEAM=%d, direction=[%f, %f]" % (polid, beamid, dir[0], dir[1]))
    365         tabidx = []
    366         for itab in range(ntable):
    367             tab = self.tables[itab]
    368             # Select rows by POLNO and BEAMNO
    369             try:
    370                 tab.set_selection(pols=[polid], beams=[beamid])
    371                 if tab.nrow() > 0: tabidx.append(itab)
    372             except: # no selection
    373                 asaplog.post()
    374                 asaplog.push("table %d - No spectrum ....skipping the table" % (itab))
    375                 asaplog.post("WARN")
    376                 continue
    377 
    378             # Select rows by direction
    379             spec = numpy.zeros(self.nchan, numpy.float)
    380             selrow = []
    381             for irow in range(tab.nrow()):
    382                 currdir = tab.get_directionval(irow)
    383                 if (abs(currdir[0] - dir[0]) > self.dirtol[0]) or \
    384                    (abs(currdir[1] - dir[1]) > self.dirtol[1]):
    385                     continue
    386                 selrow.append(irow)
    387             if len(selrow) == 0:
    388                 asaplog.post()
    389                 asaplog.push("table %d - No spectrum ....skipping the table" % (itab))
    390                 asaplog.post("WARN")
    391                 continue
    392             else:
    393                 seltab = tab.copy()
    394                 seltab.set_selection(selector(rows=selrow))
    395            
    396             if tab.nrow() > 1:
    397                 asaplog.push("table %d - More than a spectrum selected. averaging rows..." % (itab))
    398                 tab = seltab.average_time(scanav=False, weight="tintsys")
    399             else:
    400                 tab = seltab
    401 
    402             spec_array[nspec] = tab._getspectrum()
    403             nspec += 1
    404 
    405         if nspec != ntable:
    406             asaplog.post()
    407             #asaplog.push("Some tables has no spectrum with POL=%d BEAM=%d. averaging rows..." % (polid, beamid))
    408             asaplog.push("Could not find corresponding rows in some tables.")
    409             asaplog.push("Number of spectra selected = %d (table: %d)" % (nspec, ntable))
    410             if nspec < 2:
    411                 asaplog.push("At least 2 spectra are necessary for convolution")
    412                 asaplog.post("ERROR")
    413                 return False, tabidx
    414 
    415         return spec_array[0:nspec], tabidx
    416            
    417 
    418     @asaplog_post_dec
    419     def _setup_shift(self):
    420         ### define self.tables, self.signalShift, and self.imageShift
    421         if not self.intables:
    422             asaplog.post()
    423             raise RuntimeError, "Input data is not defined."
    424         #if self.baseif < 0:
    425         #    asaplog.post()
    426         #    raise RuntimeError, "Reference IFNO is not defined."
    427        
    428         byname = False
    429         #if not self.intables:
    430         if isinstance(self.intables[0], str):
    431             # A list of file name is given
    432             if not os.path.exists(self.intables[0]):
    433                 asaplog.post()
    434                 raise RuntimeError, "Could not find '%s'" % self.intables[0]
    435            
    436             stab = scantable(self.intables[0],average=False)
    437             ntab = len(self.intables)
    438             byname = True
    439         else:
    440             stab = self.intables[0]
    441             ntab = len(self.intables)
    442 
    443         if len(stab.getbeamnos()) > 1:
    444             asaplog.post()
    445             asaplog.push("Mult-beam data is not supported by this module.")
    446             asaplog.post("ERROR")
    447             return
    448 
    449         valid_ifs = stab.getifnos()
    450         if self.baseif < 0:
    451             self.baseif = valid_ifs[0]
    452             asaplog.post()
    453             asaplog.push("IFNO is not selected. Using the first IF in the first scantable. Reference IFNO = %d" % (self.baseif))
    454        
    455         if not (self.baseif in valid_ifs):
    456             asaplog.post()
    457             errmsg = "IF%d does not exist in the first scantable" %  \
    458                      self.baseif
    459             raise RuntimeError, errmsg
    460 
    461         asaplog.push("Start selecting tables and IFNOs to solve.")
    462         asaplog.push("Cheching frequency of the reference IF")
    463         unit_org = stab.get_unit()
    464         coord = stab._getcoordinfo()
    465         frame_org = coord[1]
    466         stab.set_unit("Hz")
    467         if len(self.freqframe) > 0:
    468             stab.set_freqframe(self.freqframe)
    469         stab.set_selection(ifs=[self.baseif])
    470         spx = stab._getabcissa()
    471         stab.set_selection()
    472         basech0 = spx[0]
    473         baseinc = spx[1]-spx[0]
    474         self.nchan = len(spx)
    475         if isinstance(self.freqtol, float):
    476             vftol = abs(baseinc * self.freqtol)
    477             self.freqtol = dict(value=vftol, unit="Hz")
    478         else:
    479             vftol = abs(self.freqtol['value'])
    480         inctol = abs(baseinc/float(self.nchan))
    481         asaplog.push("Reference frequency setup (Table = 0, IFNO = %d):  nchan = %d, chan0 = %f Hz, incr = %f Hz" % (self.baseif, self.nchan, basech0, baseinc))
    482         asaplog.push("Allowed frequency tolerance = %f Hz ( %f channels)" % (vftol, vftol/baseinc))
    483         poltype0 = stab.poltype()
    484        
    485         self.tables = []
    486         self.signalShift = []
    487         if self.dsbmode:
    488             self.imageShift = []
    489 
    490         for itab in range(ntab):
    491             asaplog.push("Table %d:" % itab)
    492             tab_selected = False
    493             if itab > 0:
    494                 if byname:
    495                     stab = scantable(self.intables[itab],average=False)
    496                 else:
    497                     stab = self.intables[itab]
    498                 unit_org = stab.get_unit()
    499                 coord = stab._getcoordinfo()
    500                 frame_org = coord[1]
    501                 stab.set_unit("Hz")
    502                 if len(self.freqframe) > 0:
    503                     stab.set_freqframe(self.freqframe)
    504 
    505             # Check POLTYPE should be equal to the first table.
    506             if stab.poltype() != poltype0:
    507                 asaplog.post()
    508                 raise Exception, "POLTYPE should be equal to the first table."
    509             # Multiple beam data may not handled properly
    510             if len(stab.getbeamnos()) > 1:
    511                 asaplog.post()
    512                 asaplog.push("table contains multiple beams. It may not be handled properly.")
    513                 asaplog.push("WARN")
    514 
    515             for ifno in stab.getifnos():
    516                 stab.set_selection(ifs=[ifno])
    517                 spx = stab._getabcissa()
    518                 if (len(spx) != self.nchan) or \
    519                    (abs(spx[0]-basech0) > vftol) or \
    520                    (abs(spx[1]-spx[0]-baseinc) > inctol):
    521                     continue
    522                 tab_selected = True
    523                 seltab = stab.copy()
    524                 seltab.set_unit(unit_org)
    525                 seltab.set_freqframe(frame_org)
    526                 self.tables.append(seltab)
    527                 self.signalShift.append((spx[0]-basech0)/baseinc)
    528                 if self.dsbmode:
    529                     self.imageShift.append(-self.signalShift[-1])
    530                 asaplog.push("- IF%d selected: sideband shift = %16.12e channels" % (ifno, self.signalShift[-1]))
    531             stab.set_selection()
    532             stab.set_unit(unit_org)
    533             stab.set_freqframe(frame_org)
    534             if not tab_selected:
    535                 asaplog.post()
    536                 asaplog.push("No data selected in Table %d" % itab)
    537                 asaplog.post("WARN")
    538 
    539         asaplog.push("Total number of IFs selected = %d" % len(self.tables))
    540         if len(self.tables) < 2:
    541             asaplog.post()
    542             raise RuntimeError, "At least 2 IFs are necessary for convolution!"
    543 
    544         if not self.dsbmode and len(self.imageShift) != len(self.signalShift):
    545             asaplog.post()
    546             errmsg = "User defined channel shift of image sideband has %d elements, while selected IFNOs are %d" % (len(self.imageShift), len(self.signalShift))
    547             raise RuntimeError, errmsg
    548 
    549         self.signalShift = numpy.array(self.signalShift)
    550         self.imageShift = numpy.array(self.imageShift)
    551         self.nshift = len(self.tables)
    552 
    553     @asaplog_post_dec
    554     def _preprocess_tables(self):
    555         ### temporary method to preprocess data
    556         ### Do time averaging for now.
    557         for itab in range(len(self.tables)):
    558             self.tables[itab] = self.tables[itab].average_time(scanav=False, weight="tintsys")
    559        
    560 
    561 #     def save(self, outfile, outform="ASAP", overwrite=False):
    562 #         if not overwrite and os.path.exists(outfile):
    563 #             raise RuntimeError, "Output file '%s' already exists" % outfile
    564 #
    565 #         #self.separator._save(outfile, outform)
    566 
    567 #     def done(self):
    568 #         self.close()
    569 
    570 #     def close(self):
    571 #         pass
    572 #         #del self.separator
    573    
    574 
    575 
    576 ########################################################################
    577     def _Deconvolution(self, data_array, shift, threshold=0.00000001):
    578         FObs = []
    579         Reject = 0
    580         nshift, nchan = data_array.shape
    581         nspec = nshift*(nshift-1)/2
    582         ifftObs  = numpy.zeros((nspec, nchan), numpy.float)
    583         for i in range(nshift):
    584            F = FFT.fft(data_array[i])
    585            FObs.append(F)
    586         z = 0
    587         for i in range(nshift):
    588             for j in range(i+1, nshift):
    589                 Fobs = (FObs[i]+FObs[j])/2.0
    590                 dX = (shift[j]-shift[i])*2.0*math.pi/float(self.nchan)
    591                 #print 'dX,i,j=',dX,i,j
    592                 for k in range(1,self.nchan):
    593                     if math.fabs(math.sin(dX*k)) > threshold:
    594                         Fobs[k] += ((FObs[i][k]-FObs[j][k])/2.0/(1.0-math.cos(dX*k))*math.sin(dX*k))*1.0j
    595                     else: Reject += 1
    596                 ifftObs[z] = FFT.ifft(Fobs)
    597                 z += 1
    598         print 'Threshold=%s Reject=%d' % (threshold, Reject)
    599         return ifftObs
    600 
    601     def _combineResult(self, ifftObs):
    602         nspec = len(ifftObs)
    603         sum = ifftObs[0]
    604         for i in range(1,nspec):
    605             sum += ifftObs[i]
    606         return(sum/float(nspec))
    607 
    608     def _subtractOtherSide(self, data_array, shift, Data):
    609         sum = numpy.zeros(len(Data), numpy.float)
    610         numSP = len(data_array)
    611         for i in range(numSP):
    612             SPsub = data_array[i] - Data
    613             sum += self._shiftSpectrum(SPsub, -shift[i])
    614         return(sum/float(numSP))
    615 
    616     def _shiftSpectrum(self, data, Shift):
    617         Out = numpy.zeros(self.nchan, numpy.float)
    618         w2 = Shift % 1
    619         w1 = 1.0 - w2
    620         for i in range(self.nchan):
    621             c1 = int((Shift + i) % self.nchan)
    622             c2 = (c1 + 1) % self.nchan
    623             Out[c1] += data[i] * w1
    624             Out[c2] += data[i] * w2
    625         return Out.copy()
  • trunk/python/scantable.py

    r2704 r3112  
    22
    33import os
     4import re
    45import tempfile
    56import numpy
     
    4546        l=f.readline()
    4647        f.close()
    47         if ( l.find('Scantable') != -1 ):
    48             return True
    49         elif ( l.find('Measurement Set') == -1 and
    50                l.find('Image') == -1 ):
     48        match_pattern = '^Type = (Scantable)? *$'
     49        if re.match(match_pattern,l):
    5150            return True
    5251        else:
     
    175174    return str(showprogress).lower() + ',' + str(minnrow)
    176175
     176def pack_blinfo(blinfo, maxirow):
     177    """\
     178    convert a dictionary or a list of dictionaries of baseline info
     179    into a list of comma-separated strings.
     180    """
     181    if isinstance(blinfo, dict):
     182        res = do_pack_blinfo(blinfo, maxirow)
     183        return [res] if res != '' else []
     184    elif isinstance(blinfo, list) or isinstance(blinfo, tuple):
     185        res = []
     186        for i in xrange(len(blinfo)):
     187            resi = do_pack_blinfo(blinfo[i], maxirow)
     188            if resi != '':
     189                res.append(resi)
     190    return res
     191
     192def do_pack_blinfo(blinfo, maxirow):
     193    """\
     194    convert a dictionary of baseline info for a spectrum into
     195    a comma-separated string.
     196    """
     197    dinfo = {}
     198    for key in ['row', 'blfunc', 'masklist']:
     199        if blinfo.has_key(key):
     200            val = blinfo[key]
     201            if key == 'row':
     202                irow = val
     203            if isinstance(val, list) or isinstance(val, tuple):
     204                slval = []
     205                for i in xrange(len(val)):
     206                    if isinstance(val[i], list) or isinstance(val[i], tuple):
     207                        for j in xrange(len(val[i])):
     208                            slval.append(str(val[i][j]))
     209                    else:
     210                        slval.append(str(val[i]))
     211                sval = ",".join(slval)
     212            else:
     213                sval = str(val)
     214
     215            dinfo[key] = sval
     216        else:
     217            raise ValueError("'"+key+"' is missing in blinfo.")
     218
     219    if irow >= maxirow: return ''
     220   
     221    for key in ['order', 'npiece', 'nwave']:
     222        if blinfo.has_key(key):
     223            val = blinfo[key]
     224            if isinstance(val, list) or isinstance(val, tuple):
     225                slval = []
     226                for i in xrange(len(val)):
     227                    slval.append(str(val[i]))
     228                sval = ",".join(slval)
     229            else:
     230                sval = str(val)
     231            dinfo[key] = sval
     232
     233    fspec_keys = {'poly': 'order', 'chebyshev': 'order', 'cspline': 'npiece', 'sinusoid': 'nwave'}
     234
     235    fspec_key = fspec_keys[dinfo['blfunc']]
     236    if not blinfo.has_key(fspec_key):
     237        raise ValueError("'"+fspec_key+"' is missing in blinfo.")
     238
     239    clip_params_n = 0
     240    for key in ['clipthresh', 'clipniter']:
     241        if blinfo.has_key(key):
     242            clip_params_n += 1
     243            dinfo[key] = str(blinfo[key])
     244
     245    if clip_params_n == 0:
     246        dinfo['clipthresh'] = '0.0'
     247        dinfo['clipniter']  = '0'
     248    elif clip_params_n != 2:
     249        raise ValueError("both 'clipthresh' and 'clipniter' must be given for n-sigma clipping.")
     250
     251    lf_params_n = 0
     252    for key in ['thresh', 'edge', 'chan_avg_limit']:
     253        if blinfo.has_key(key):
     254            lf_params_n += 1
     255            val = blinfo[key]
     256            if isinstance(val, list) or isinstance(val, tuple):
     257                slval = []
     258                for i in xrange(len(val)):
     259                    slval.append(str(val[i]))
     260                sval = ",".join(slval)
     261            else:
     262                sval = str(val)
     263            dinfo[key] = sval
     264
     265    if lf_params_n == 3:
     266        dinfo['use_linefinder'] = 'true'
     267    elif lf_params_n == 0:
     268        dinfo['use_linefinder'] = 'false'
     269        dinfo['thresh']         = ''
     270        dinfo['edge']           = ''
     271        dinfo['chan_avg_limit'] = ''
     272    else:
     273        raise ValueError("all of 'thresh', 'edge' and 'chan_avg_limit' must be given to use linefinder.")
     274   
     275    slblinfo = [dinfo['row'], dinfo['blfunc'], dinfo[fspec_key], dinfo['masklist'], \
     276                dinfo['clipthresh'], dinfo['clipniter'], \
     277                dinfo['use_linefinder'], dinfo['thresh'], dinfo['edge'], dinfo['chan_avg_limit']]
     278   
     279    return ":".join(slblinfo)
     280
     281def parse_fitresult(sres):
     282    """\
     283    Parse the returned value of apply_bltable() or sub_baseline() and
     284    extract row number, the best-fit coefficients and rms, then pack
     285    them into a dictionary.
     286    The input value is generated by Scantable::packFittingResults() and
     287    formatted as 'row:coeff[0],coeff[1],..,coeff[n-1]:rms'.
     288    """
     289    res = []
     290    for i in xrange(len(sres)):
     291        (srow, scoeff, srms) = sres[i].split(":")
     292        row = int(srow)
     293        rms = float(srms)
     294        lscoeff = scoeff.split(",")
     295        coeff = []
     296        for j in xrange(len(lscoeff)):
     297            coeff.append(float(lscoeff[j]))
     298        res.append({'row': row, 'coeff': coeff, 'rms': rms})
     299
     300    return res
     301
     302def is_number(s):
     303    s = s.strip()
     304    res = True
     305    try:
     306        a = float(s)
     307        res = True
     308    except:
     309        res = False
     310    finally:
     311        return res
     312
     313def is_frequency(s):
     314    s = s.strip()
     315    return (s[-2:].lower() == "hz")
     316
     317def get_freq_by_string(s1, s2):
     318    if not (is_number(s1) and is_frequency(s2)):
     319        raise RuntimeError("Invalid input string.")
     320   
     321    prefix_list = ["a", "f", "p", "n", "u", "m", ".", "k", "M", "G", "T", "P", "E"]
     322    factor_list = [1e-18, 1e-15, 1e-12, 1e-9, 1e-6, 1e-3, 1.0, 1e+3, 1e+6, 1e+9, 1e+12, 1e+15, 1e+18]
     323
     324    s1 = s1.strip()
     325    s2 = s2.strip()
     326   
     327    prefix = s2[-3:-2]
     328    if is_number(prefix):
     329        res1 = float(s1)
     330        res2 = float(s2[:-2])
     331    else:
     332        factor = factor_list[prefix_list.index(prefix)]
     333        res1 = float(s1) * factor
     334        res2 = float(s2[:-3]) * factor
     335
     336    return (res1, res2)
     337
     338def is_velocity(s):
     339    s = s.strip()
     340    return (s[-3:].lower() == "m/s")
     341
     342def get_velocity_by_string(s1, s2):
     343    if not (is_number(s1) and is_velocity(s2)):
     344        raise RuntimeError("Invalid input string.")
     345
     346    # note that the default velocity unit is km/s
     347    prefix_list = [".", "k"]
     348    factor_list = [1e-3, 1.0]
     349
     350    s1 = s1.strip()
     351    s2 = s2.strip()
     352
     353    prefix = s2[-4:-3]
     354    if is_number(prefix): # in case velocity unit m/s
     355        res1 = float(s1) * 1e-3
     356        res2 = float(s2[:-3]) * 1e-3
     357    else:
     358        factor = factor_list[prefix_list.index(prefix)]
     359        res1 = float(s1) * factor
     360        res2 = float(s2[:-4]) * factor
     361
     362    return (res1, res2)
     363
     364def get_frequency_by_velocity(restfreq, vel, doppler):
     365    # vel is in unit of km/s
     366
     367    # speed of light
     368    vel_c = 299792.458
     369
     370    import math
     371    r = vel / vel_c
     372
     373    if doppler.lower() == 'radio':
     374        return restfreq * (1.0 - r)
     375    if doppler.lower() == 'optical':
     376        return restfreq / (1.0 + r)
     377    else:
     378        return restfreq * math.sqrt((1.0 - r) / (1.0 + r))
     379
     380def get_restfreq_in_Hz(s_restfreq):
     381    value = 0.0
     382    unit = ""
     383    s = s_restfreq.replace(" ","")
     384
     385    for i in range(len(s))[::-1]:
     386        if s[i].isalpha():
     387            unit = s[i] + unit
     388        else:
     389            value = float(s[0:i+1])
     390            break
     391
     392    if (unit == "") or (unit.lower() == "hz"):
     393        return value
     394    elif (len(unit) == 3) and (unit[1:3].lower() == "hz"):
     395        unitprefix = unit[0]
     396        factor = 1.0
     397
     398        prefix_list = ["a", "f", "p", "n", "u", "m", ".", "k", "M", "G", "T", "P", "E"]
     399        factor_list = [1e-18, 1e-15, 1e-12, 1e-9, 1e-6, 1e-3, 1.0, 1e+3, 1e+6, 1e+9, 1e+12, 1e+15, 1e+18]
     400        factor = factor_list[prefix_list.index(unitprefix)]
     401        """
     402        if (unitprefix == 'a'):
     403            factor = 1.0e-18
     404        elif (unitprefix == 'f'):
     405            factor = 1.0e-15
     406        elif (unitprefix == 'p'):
     407            factor = 1.0e-12
     408        elif (unitprefix == 'n'):
     409            factor = 1.0e-9
     410        elif (unitprefix == 'u'):
     411            factor = 1.0e-6
     412        elif (unitprefix == 'm'):
     413            factor = 1.0e-3
     414        elif (unitprefix == 'k'):
     415            factor = 1.0e+3
     416        elif (unitprefix == 'M'):
     417            factor = 1.0e+6
     418        elif (unitprefix == 'G'):
     419            factor = 1.0e+9
     420        elif (unitprefix == 'T'):
     421            factor = 1.0e+12
     422        elif (unitprefix == 'P'):
     423            factor = 1.0e+15
     424        elif (unitprefix == 'E'):
     425            factor = 1.0e+18
     426        """
     427        return value*factor
     428    else:
     429        mesg = "wrong unit of restfreq."
     430        raise Exception, mesg
     431
     432def normalise_restfreq(in_restfreq):
     433    if isinstance(in_restfreq, float):
     434        return in_restfreq
     435    elif isinstance(in_restfreq, int) or isinstance(in_restfreq, long):
     436        return float(in_restfreq)
     437    elif isinstance(in_restfreq, str):
     438        return get_restfreq_in_Hz(in_restfreq)
     439    elif isinstance(in_restfreq, list) or isinstance(in_restfreq, numpy.ndarray):
     440        if isinstance(in_restfreq, numpy.ndarray):
     441            if len(in_restfreq.shape) > 1:
     442                mesg = "given in numpy.ndarray, in_restfreq must be 1-D."
     443                raise Exception, mesg
     444       
     445        res = []
     446        for i in xrange(len(in_restfreq)):
     447            elem = in_restfreq[i]
     448            if isinstance(elem, float):
     449                res.append(elem)
     450            elif isinstance(elem, int) or isinstance(elem, long):
     451                res.append(float(elem))
     452            elif isinstance(elem, str):
     453                res.append(get_restfreq_in_Hz(elem))
     454            elif isinstance(elem, dict):
     455                if isinstance(elem["value"], float):
     456                    res.append(elem)
     457                elif isinstance(elem["value"], int):
     458                    dictelem = {}
     459                    dictelem["name"]  = elem["name"]
     460                    dictelem["value"] = float(elem["value"])
     461                    res.append(dictelem)
     462                elif isinstance(elem["value"], str):
     463                    dictelem = {}
     464                    dictelem["name"]  = elem["name"]
     465                    dictelem["value"] = get_restfreq_in_Hz(elem["value"])
     466                    res.append(dictelem)
     467            else:
     468                mesg = "restfreq elements must be float, int, or string."
     469                raise Exception, mesg
     470        return res
     471    else:
     472        mesg = "wrong type of restfreq given."
     473        raise Exception, mesg
     474
     475def set_restfreq(s, restfreq):
     476    rfset = (restfreq != '') and (restfreq != [])
     477    if rfset:
     478        s.set_restfreqs(normalise_restfreq(restfreq))
     479
    177480class scantable(Scantable):
    178481    """\
     
    210513                          Default (false) is taken from rc file.
    211514
     515            getpt:        Whether to import direction from MS/POINTING
     516                          table properly or not.
     517                          This is effective only when filename is MS.
     518                          The default (True) is to import direction
     519                          from MS/POINTING.
    212520        """
    213521        if average is None:
     
    244552                    self._fill([filename], unit, average, opts)
    245553                elif os.path.isfile(filename):
    246                     self._fill([filename], unit, average)
     554                    opts={'nro': {}}
     555                    nrokeys=['freqref']
     556                    for key in nrokeys:
     557                        if key in args.keys():
     558                            opts['nro'][key] = args[key]
     559                    self._fill([filename], unit, average, opts)
    247560                    # only apply to new data not "copy constructor"
    248561                    self.parallactify(parallactify)
     
    578891
    579892    @asaplog_post_dec
    580     def stats(self, stat='stddev', mask=None, form='3.3f', row=None):
     893    def stats(self, stat='stddev', mask=None, form='3.3f', row=None, skip_flaggedrow=False):
    581894        """\
    582895        Determine the specified statistic of the current beam/if/pol
     
    596909            row:     row number of spectrum to process.
    597910                     (default is None: for all rows)
     911
     912            skip_flaggedrow: if True, skip outputting text for flagged
     913                             spectra. default is False.
    598914
    599915        Example:
     
    608924                             "number of channels. Please use setselection() "
    609925                             "to select individual IFs")
    610         rtnabc = False
    611         if stat.lower().endswith('_abc'): rtnabc = True
    612926        getchan = False
    613927        if stat.lower().startswith('min') or stat.lower().startswith('max'):
     
    615929            getchan = True
    616930            statvals = []
    617         if not rtnabc:
     931
     932        rtnabc = False
     933        if stat.lower().endswith('_abc'):
     934            rtnabc = True
     935        else:
    618936            if row == None:
    619937                statvals = self._math._stats(self, mask, stat)
     
    641959            statunit= ''
    642960            if getchan:
    643                 qx, qy = self.chan2data(rowno=i, chan=chan[i])
    644                 if rtnabc:
    645                     statvals.append(qx['value'])
    646                     refstr = ('(value: %'+form) % (qy['value'])+' ['+qy['unit']+'])'
    647                     statunit= '['+qx['unit']+']'
     961                if self._is_all_chan_flagged(i):
     962                    if rtnabc:
     963                        statvals.append(None)
    648964                else:
    649                     refstr = ('(@ %'+form) % (qx['value'])+' ['+qx['unit']+'])'
     965                    qx, qy = self.chan2data(rowno=i, chan=chan[i])
     966                    if rtnabc:
     967                        statvals.append(qx['value'])
     968                        refstr = ('(value: %'+form) % (qy['value'])+' ['+qy['unit']+'])'
     969                        statunit= '['+qx['unit']+']'
     970                    else:
     971                        refstr = ('(@ %'+form) % (qx['value'])+' ['+qx['unit']+'])'
     972
     973            if self._is_all_chan_flagged(i):
     974                if not rtnabc:
     975                    statvals[i] = None
     976                if skip_flaggedrow:
     977                    continue
    650978
    651979            tm = self._gettime(i)
     
    659987            if len(rows) > 1:
    660988                # out += ('= %'+form) % (outvec[i]) +'   '+refstr+'\n'
    661                 out += ('= %'+form) % (statvals[i]) +'   '+refstr+'\n'
     989                if statvals[i] is None:
     990                    out += ('= None(flagged)') + '   '+refstr+'\n'
     991                else:
     992                    out += ('= %'+form) % (statvals[i]) +'   '+refstr+'\n'
    662993            else:
    663994                # out += ('= %'+form) % (outvec[0]) +'   '+refstr+'\n'
    664                 out += ('= %'+form) % (statvals[0]) +'   '+refstr+'\n'
     995                if statvals[0] is None:
     996                    out += ('= None(flagged)') + '   '+refstr+'\n'
     997                else:
     998                    out += ('= %'+form) % (statvals[0]) +'   '+refstr+'\n'
    665999            out +=  sep+"\n"
    6661000
     
    6831017        asaplog.push(''.join(x), False)
    6841018
     1019        if skip_flaggedrow:
     1020            nstatvals = len(statvals)
     1021            for i in reversed(xrange(nstatvals)):
     1022                if statvals[i] is None:
     1023                    del statvals[i]
    6851024        return statvals
    6861025
     
    7651104        return self._get_column( self._gettsysspectrum, row )
    7661105
     1106    def set_tsys(self, values, row=-1):
     1107        """\
     1108        Set the Tsys value(s) of the given 'row' or the whole scantable
     1109        (selection).
     1110       
     1111        Parameters:
     1112
     1113            values:    a scalar or list (if Tsys is a vector) of Tsys value(s)
     1114            row:       the row number to apply Tsys values to.
     1115                       (default all rows)
     1116                       
     1117        """
     1118           
     1119        if not hasattr(values, "__len__"):
     1120            values = [values]
     1121        self._settsys(values, row)
     1122
    7671123    def get_weather(self, row=-1):
    7681124        """\
    769         Return the weather informations.
     1125        Return the weather information.
    7701126
    7711127        Parameters:
     
    7781134
    7791135        """
    780 
     1136        if row >= len(self):
     1137            raise IndexError("row out of range")
    7811138        values = self._get_column(self._get_weather, row)
    7821139        if row > -1:
     
    7881145            out = []
    7891146            for r in values:
    790 
    7911147                out.append({'temperature': r[0],
    7921148                            'pressure': r[1], 'humidity' : r[2],
     
    10051361        self._add_history("set_feedtype", vars())
    10061362
     1363    @asaplog_post_dec
     1364    def get_doppler(self):
     1365        """\
     1366        Get the doppler.
     1367        """
     1368        return self._getcoordinfo()[2]
     1369   
    10071370    @asaplog_post_dec
    10081371    def set_doppler(self, doppler='RADIO'):
     
    14711834
    14721835    @asaplog_post_dec
     1836    def parse_spw_selection(self, selectstring, restfreq=None, frame=None, doppler=None):
     1837        """
     1838        Parse MS type spw/channel selection syntax.
     1839
     1840        Parameters:
     1841            selectstring : A string expression of spw and channel selection.
     1842                           Comma-separated expressions mean different spw -
     1843                           channel combinations. Spws and channel selections
     1844                           are partitioned by a colon ':'. In a single
     1845                           selection expression, you can put multiple values
     1846                           separated by semicolons ';'. Both for spw and
     1847                           channel selection, allowed cases include single
     1848                           value, blank('') or asterisk('*') to specify all
     1849                           available values, two values connected with a
     1850                           tilde ('~') to specify an inclusive range. Unit
     1851                           strings for frequency or velocity can be added to
     1852                           the tilde-connected values. For channel selection
     1853                           expression, placing a '<' or a '>' is possible to
     1854                           specify a semi-infinite interval as well.
     1855
     1856                     examples:
     1857                           '' or '*'   = all spws (all channels)
     1858                           '<2,4~6,9'  = Spws 0,1,4,5,6,9 (all channels)
     1859                           '3:3~45;60' = channels 3 to 45 and 60 in spw 3
     1860                           '0~1:2~6,8' = channels 2 to 6 in spws 0,1, and
     1861                                         all channels in spw8
     1862                           '1.3~1.5GHz' = all spws whose central frequency
     1863                                          falls in frequency range between
     1864                                          1.3GHz and 1.5GHz.
     1865                           '1.3~1.5GHz:1.3~1.5GHz' = channels which fall
     1866                                                     between the specified
     1867                                                     frequency range in spws
     1868                                                     whose central frequency
     1869                                                     falls in the specified
     1870                                                     frequency range.
     1871                           '1:-200~250km/s' = channels that fall between the
     1872                                              specified velocity range in
     1873                                              spw 1.
     1874            restfreq: the rest frequency.
     1875                     examples: '115.2712GHz', 115271201800.0
     1876            frame:   an optional frame type, default 'LSRK'. Valid frames are:
     1877                     'TOPO', 'LSRD', 'LSRK', 'BARY',
     1878                     'GEO', 'GALACTO', 'LGROUP', 'CMB'
     1879            doppler: one of 'RADIO', 'OPTICAL', 'Z', 'BETA', 'GAMMA'
     1880        Returns:
     1881        A dictionary of selected (valid) spw and masklist pairs,
     1882        e.g. {'0': [[50,250],[350,462]], '2': [[100,400],[550,974]]}
     1883        """
     1884        if not isinstance(selectstring, str):
     1885            asaplog.post()
     1886            asaplog.push("Expression of spw/channel selection must be a string.")
     1887            asaplog.post("ERROR")
     1888
     1889        orig_unit = self.get_unit()
     1890        self.set_unit('channel')
     1891       
     1892        if restfreq is not None:
     1893            orig_molids = self._getmolidcol_list()
     1894            set_restfreq(self, restfreq)
     1895
     1896        orig_coord = self._getcoordinfo()
     1897
     1898        if frame is not None:
     1899            orig_frame = orig_coord[1]
     1900            self.set_freqframe(frame)
     1901
     1902        if doppler is not None:
     1903            orig_doppler = orig_coord[2]
     1904            self.set_doppler(doppler)
     1905       
     1906        valid_ifs = self.getifnos()
     1907
     1908        comma_sep = selectstring.split(",")
     1909        res = {}
     1910
     1911        for cms_elem in comma_sep:
     1912            colon_sep = cms_elem.split(":")
     1913           
     1914            if (len(colon_sep) > 2):
     1915                raise RuntimeError("Invalid selection expression: more than two colons!")
     1916           
     1917            # parse spw expression and store result in spw_list.
     1918            # allowed cases include '', '*', 'a', '<a', '>a', 'a~b',
     1919            # 'a~b*Hz' (where * can be '', 'k', 'M', 'G' etc.),
     1920            # 'a~b*m/s' (where * can be '' or 'k') and also
     1921            # several of the above expressions connected with ';'.
     1922           
     1923            spw_list = []
     1924
     1925            semicolon_sep = colon_sep[0].split(";")
     1926           
     1927            for scs_elem in semicolon_sep:
     1928                scs_elem = scs_elem.strip()
     1929               
     1930                lt_sep = scs_elem.split("<")
     1931                gt_sep = scs_elem.split(">")
     1932                ti_sep = scs_elem.split("~")
     1933               
     1934                lt_sep_length = len(lt_sep)
     1935                gt_sep_length = len(gt_sep)
     1936                ti_sep_length = len(ti_sep)
     1937               
     1938                len_product = lt_sep_length * gt_sep_length * ti_sep_length
     1939
     1940                if (len_product > 2):
     1941                    # '<', '>' and '~' must not coexist in a single spw expression
     1942                   
     1943                    raise RuntimeError("Invalid spw selection.")
     1944               
     1945                elif (len_product == 1):
     1946                    # '', '*', or single spw number.
     1947                   
     1948                    if (scs_elem == "") or (scs_elem == "*"):
     1949                        spw_list = valid_ifs[:] # deep copy
     1950                   
     1951                    else: # single number
     1952                        expr = int(scs_elem)
     1953                        spw_list.append(expr)
     1954                        if expr not in valid_ifs:
     1955                            asaplog.push("Invalid spw given. Ignored.")
     1956                   
     1957                else: # (len_product == 2)
     1958                    # namely, one of '<', '>' or '~' appears just once.
     1959                   
     1960                    if (lt_sep_length == 2): # '<a'
     1961                        if is_number(lt_sep[1]):
     1962                            no_valid_spw = True
     1963                            for i in valid_ifs:
     1964                                if (i < float(lt_sep[1])):
     1965                                    spw_list.append(i)
     1966                                    no_valid_spw = False
     1967
     1968                            if no_valid_spw:
     1969                                raise ValueError("Invalid spw selection ('<" + str(lt_sep[1]) + "').")
     1970                       
     1971                        else:
     1972                            raise RuntimeError("Invalid spw selection.")
     1973                       
     1974                    elif (gt_sep_length == 2): # '>a'
     1975                        if is_number(gt_sep[1]):
     1976                            no_valid_spw = True
     1977                            for i in valid_ifs:
     1978                                if (i > float(gt_sep[1])):
     1979                                    spw_list.append(i)
     1980                                    no_valid_spw = False
     1981                           
     1982                            if no_valid_spw:
     1983                                raise ValueError("Invalid spw selection ('>" + str(gt_sep[1]) + "').")
     1984                       
     1985                        else:
     1986                            raise RuntimeError("Invalid spw selection.")
     1987                       
     1988                    else: # (ti_sep_length == 2) where both boundaries inclusive
     1989                        expr0 = ti_sep[0].strip()
     1990                        expr1 = ti_sep[1].strip()
     1991
     1992                        if is_number(expr0) and is_number(expr1):
     1993                            # 'a~b'
     1994                            expr_pmin = min(float(expr0), float(expr1))
     1995                            expr_pmax = max(float(expr0), float(expr1))
     1996                            has_invalid_spw = False
     1997                            no_valid_spw = True
     1998                           
     1999                            for i in valid_ifs:
     2000                                if (expr_pmin <= i) and (i <= expr_pmax):
     2001                                    spw_list.append(i)
     2002                                    no_valid_spw = False
     2003                                else:
     2004                                    has_invalid_spw = True
     2005
     2006                            if has_invalid_spw:
     2007                                msg = "Invalid spw is given. Ignored."
     2008                                asaplog.push(msg)
     2009                                asaplog.post()
     2010
     2011                            if no_valid_spw:
     2012                                raise ValueError("No valid spw in range ('" + str(expr_pmin) + "~" + str(expr_pmax) + "').")
     2013                       
     2014                        elif is_number(expr0) and is_frequency(expr1):
     2015                            # 'a~b*Hz'
     2016                            (expr_f0, expr_f1) = get_freq_by_string(expr0, expr1)
     2017                            expr_fmin = min(expr_f0, expr_f1)
     2018                            expr_fmax = max(expr_f0, expr_f1)
     2019                            no_valid_spw = True
     2020                           
     2021                            for coord in self._get_coordinate_list():
     2022                                spw = coord['if']
     2023                               
     2024                                """
     2025                                expr_p0 = coord['coord'].to_pixel(expr_f0)
     2026                                expr_p1 = coord['coord'].to_pixel(expr_f1)
     2027                                expr_pmin = min(expr_p0, expr_p1)
     2028                                expr_pmax = max(expr_p0, expr_p1)
     2029                               
     2030                                pmin = 0.0
     2031                                pmax = float(self.nchan(spw) - 1)
     2032
     2033                                if ((expr_pmax - pmin)*(expr_pmin - pmax) <= 0.0):
     2034                                    spw_list.append(spw)
     2035                                    no_valid_spw = False
     2036                                """
     2037
     2038                                crd = coord['coord']
     2039                                fhead = crd.to_frequency(0)
     2040                                ftail = crd.to_frequency(self.nchan(spw) - 1)
     2041                                fcen  = (fhead + ftail) / 2.0
     2042
     2043                                if ((expr_fmin <= fcen) and (fcen <= expr_fmax)):
     2044                                    spw_list.append(spw)
     2045                                    no_valid_spw = False
     2046                               
     2047                            if no_valid_spw:
     2048                                raise ValueError("No valid spw in range ('" + str(expr0) + "~" + str(expr1) + "').")
     2049                               
     2050                        elif is_number(expr0) and is_velocity(expr1):
     2051                            # 'a~b*m/s'
     2052                            (expr_v0, expr_v1) = get_velocity_by_string(expr0, expr1)
     2053                            expr_vmin = min(expr_v0, expr_v1)
     2054                            expr_vmax = max(expr_v0, expr_v1)
     2055                            no_valid_spw = True
     2056                           
     2057                            for coord in self._get_coordinate_list():
     2058                                spw = coord['if']
     2059
     2060                                """
     2061                                pmin = 0.0
     2062                                pmax = float(self.nchan(spw) - 1)
     2063                               
     2064                                vel0 = coord['coord'].to_velocity(pmin)
     2065                                vel1 = coord['coord'].to_velocity(pmax)
     2066                               
     2067                                vmin = min(vel0, vel1)
     2068                                vmax = max(vel0, vel1)
     2069
     2070                                if ((expr_vmax - vmin)*(expr_vmin - vmax) <= 0.0):
     2071                                    spw_list.append(spw)
     2072                                    no_valid_spw = False
     2073                                """
     2074
     2075                                crd = coord['coord']
     2076                                fhead = crd.to_frequency(0)
     2077                                ftail = crd.to_frequency(self.nchan(spw) - 1)
     2078                                fcen  = (fhead + ftail) / 2.0
     2079                                vcen  = crd.to_velocity(crd.to_pixel(fcen))
     2080
     2081                                if ((expr_vmin <= vcen) and (vcen <= expr_vmax)):
     2082                                    spw_list.append(spw)
     2083                                    no_valid_spw = False
     2084
     2085                            if no_valid_spw:
     2086                                raise ValueError("No valid spw in range ('" + str(expr0) + "~" + str(expr1) + "').")
     2087                           
     2088                        else:
     2089                            # cases such as 'aGHz~bkm/s' are not allowed now
     2090                            raise RuntimeError("Invalid spw selection.")
     2091
     2092            # check spw list and remove invalid ones.
     2093            # if no valid spw left, emit ValueError.
     2094            if len(spw_list) == 0:
     2095                raise ValueError("No valid spw in given range.")
     2096           
     2097            # parse channel expression and store the result in crange_list.
     2098            # allowed cases include '', 'a~b', 'a*Hz~b*Hz' (where * can be
     2099            # '', 'k', 'M', 'G' etc.), 'a*m/s~b*m/s' (where * can be '' or 'k')
     2100            # and also several of the above expressions connected with ';'.
     2101           
     2102            for spw in spw_list:
     2103                pmin = 0.0
     2104                pmax = float(self.nchan(spw) - 1)
     2105
     2106                molid = self._getmolidcol_list()[self.get_first_rowno_by_if(spw)]
     2107               
     2108                if (len(colon_sep) == 1):
     2109                    # no expression for channel selection,
     2110                    # which means all channels are to be selected.
     2111                    crange_list = [[pmin, pmax]]
     2112               
     2113                else: # (len(colon_sep) == 2)
     2114                    crange_list = []
     2115                   
     2116                    found = False
     2117                    for i in self._get_coordinate_list():
     2118                        if (i['if'] == spw):
     2119                            coord = i['coord']
     2120                            found = True
     2121                            break
     2122
     2123                    if found:
     2124                        semicolon_sep = colon_sep[1].split(";")
     2125                        for scs_elem in semicolon_sep:
     2126                            scs_elem = scs_elem.strip()
     2127
     2128                            ti_sep = scs_elem.split("~")
     2129                            ti_sep_length = len(ti_sep)
     2130
     2131                            if (ti_sep_length > 2):
     2132                                raise RuntimeError("Invalid channel selection.")
     2133                       
     2134                            elif (ti_sep_length == 1):
     2135                                if (scs_elem == "") or (scs_elem == "*"):
     2136                                    # '' and '*' for all channels
     2137                                    crange_list = [[pmin, pmax]]
     2138                                    break
     2139                                elif (is_number(scs_elem)):
     2140                                    # single channel given
     2141                                    crange_list.append([float(scs_elem), float(scs_elem)])
     2142                                else:
     2143                                    raise RuntimeError("Invalid channel selection.")
     2144
     2145                            else: #(ti_sep_length == 2)
     2146                                expr0 = ti_sep[0].strip()
     2147                                expr1 = ti_sep[1].strip()
     2148
     2149                                if is_number(expr0) and is_number(expr1):
     2150                                    # 'a~b'
     2151                                    expr_pmin = min(float(expr0), float(expr1))
     2152                                    expr_pmax = max(float(expr0), float(expr1))
     2153
     2154                                elif is_number(expr0) and is_frequency(expr1):
     2155                                    # 'a~b*Hz'
     2156                                    (expr_f0, expr_f1) = get_freq_by_string(expr0, expr1)
     2157                                    expr_p0 = coord.to_pixel(expr_f0)
     2158                                    expr_p1 = coord.to_pixel(expr_f1)
     2159                                    expr_pmin = min(expr_p0, expr_p1)
     2160                                    expr_pmax = max(expr_p0, expr_p1)
     2161
     2162                                elif is_number(expr0) and is_velocity(expr1):
     2163                                    # 'a~b*m/s'
     2164                                    restf = self.get_restfreqs()[molid][0]
     2165                                    (expr_v0, expr_v1) = get_velocity_by_string(expr0, expr1)
     2166                                    dppl = self.get_doppler()
     2167                                    expr_f0 = get_frequency_by_velocity(restf, expr_v0, dppl)
     2168                                    expr_f1 = get_frequency_by_velocity(restf, expr_v1, dppl)
     2169                                    expr_p0 = coord.to_pixel(expr_f0)
     2170                                    expr_p1 = coord.to_pixel(expr_f1)
     2171                                    expr_pmin = min(expr_p0, expr_p1)
     2172                                    expr_pmax = max(expr_p0, expr_p1)
     2173                           
     2174                                else:
     2175                                    # cases such as 'aGHz~bkm/s' are not allowed now
     2176                                    raise RuntimeError("Invalid channel selection.")
     2177
     2178                                cmin = max(pmin, expr_pmin)
     2179                                cmax = min(pmax, expr_pmax)
     2180                                # if the given range of channel selection has overwrap with
     2181                                # that of current spw, output the overwrap area.
     2182                                if (cmin <= cmax):
     2183                                    cmin = float(int(cmin + 0.5))
     2184                                    cmax = float(int(cmax + 0.5))
     2185                                    crange_list.append([cmin, cmax])
     2186
     2187                    if (len(crange_list) == 0):
     2188                        crange_list.append([])
     2189
     2190                if (len(crange_list[0]) > 0):
     2191                    if res.has_key(spw):
     2192                        res[spw].extend(crange_list)
     2193                    else:
     2194                        res[spw] = crange_list
     2195
     2196        for spw in res.keys():
     2197            if spw in valid_ifs:
     2198                # remove duplicated channel ranges
     2199                for i in reversed(xrange(len(res[spw]))):
     2200                    for j in xrange(i):
     2201                        if ((res[spw][i][0]-res[spw][j][1])*(res[spw][i][1]-res[spw][j][0]) <= 0) or \
     2202                            (min(abs(res[spw][i][0]-res[spw][j][1]),abs(res[spw][j][0]-res[spw][i][1])) == 1):
     2203                            asaplog.post()
     2204                            merge_warn_mesg = "Spw " + str(spw) + ": overwrapping channel ranges are merged."
     2205                            asaplog.push(merge_warn_mesg)
     2206                            asaplog.post('WARN')
     2207                            res[spw][j][0] = min(res[spw][i][0], res[spw][j][0])
     2208                            res[spw][j][1] = max(res[spw][i][1], res[spw][j][1])
     2209                            res[spw].pop(i)
     2210                            break
     2211            else:
     2212                del res[spw]
     2213
     2214        if len(res) == 0:
     2215            raise RuntimeError("No valid spw.")
     2216       
     2217        # restore original values
     2218        self.set_unit(orig_unit)
     2219        if restfreq is not None:
     2220            self._setmolidcol_list(orig_molids)
     2221        if frame is not None:
     2222            self.set_freqframe(orig_frame)
     2223        if doppler is not None:
     2224            self.set_doppler(orig_doppler)
     2225       
     2226        return res
     2227   
     2228    @asaplog_post_dec
     2229    def get_first_rowno_by_if(self, ifno):
     2230        found = False
     2231        for irow in xrange(self.nrow()):
     2232            if (self.getif(irow) == ifno):
     2233                res = irow
     2234                found = True
     2235                break
     2236
     2237        if not found: raise RuntimeError("No valid spw.")
     2238       
     2239        return res
     2240
     2241    @asaplog_post_dec
     2242    def _get_coordinate_list(self):
     2243        res = []
     2244        spws = self.getifnos()
     2245        for spw in spws:
     2246            elem = {}
     2247            elem['if']    = spw
     2248            elem['coord'] = self.get_coordinate(self.get_first_rowno_by_if(spw))
     2249            res.append(elem)
     2250
     2251        return res
     2252   
     2253    @asaplog_post_dec
    14732254    def parse_maskexpr(self, maskstring):
    14742255        """
     
    15012282            maskstring = str(valid_ifs)[1:-1]
    15022283        ## split each selection "IF range[:CHAN range]"
    1503         sellist = maskstring.split(',')
     2284        # split maskstring by "<spaces>,<spaces>"
     2285        comma_sep = re.compile('\s*,\s*')
     2286        sellist = comma_sep.split(maskstring)
     2287        # separator by "<spaces>:<spaces>"
     2288        collon_sep = re.compile('\s*:\s*')
    15042289        for currselstr in sellist:
    1505             selset = currselstr.split(':')
     2290            selset = collon_sep.split(currselstr)
    15062291            # spw and mask string (may include ~, < or >)
    15072292            spwmasklist = self._parse_selection(selset[0], typestr='integer',
     
    16372422        if smode == 'r':
    16382423            minidx = 0
    1639             maxidx = self.nrow()
     2424            maxidx = self.nrow()-1
    16402425        else:
    16412426            idx = getattr(self,"get"+mode+"nos")()
     
    16432428            maxidx = max(idx)
    16442429            del idx
    1645         sellist = selexpr.split(',')
     2430        # split selexpr by "<spaces>,<spaces>"
     2431        comma_sep = re.compile('\s*,\s*')
     2432        sellist = comma_sep.split(selexpr)
    16462433        idxlist = []
    16472434        for currselstr in sellist:
     
    16512438            for thelist in currlist:
    16522439                idxlist += range(thelist[0],thelist[1]+1)
     2440        # remove duplicated elements after first ones
     2441        for i in reversed(xrange(len(idxlist))):
     2442            if idxlist.index(idxlist[i]) < i:
     2443                idxlist.pop(i)
     2444
     2445        # remove elements outside range [minidx, maxidx] for smode='r'
     2446        if smode == 'r':
     2447            for i in reversed(xrange(len(idxlist))):
     2448                if (idxlist[i] < minidx) or (idxlist[i] > maxidx):
     2449                    idxlist.pop(i)
     2450       
    16532451        msg = "Selected %s: %s" % (mode.upper()+"NO", str(idxlist))
    16542452        asaplog.push(msg)
     
    16762474            --> returns [[0.,2.5],[5.0,7.0],[9.,9.]]
    16772475        """
    1678         selgroups = selstr.split(';')
     2476        # split selstr by '<spaces>;<spaces>'
     2477        semi_sep = re.compile('\s*;\s*')
     2478        selgroups = semi_sep.split(selstr)
    16792479        sellists = []
    16802480        if typestr.lower().startswith('int'):
     
    16852485       
    16862486        for currsel in  selgroups:
     2487            if currsel.strip() == '*' or len(currsel.strip()) == 0:
     2488                minsel = minval
     2489                maxsel = maxval
    16872490            if currsel.find('~') > 0:
    16882491                # val0 <= x <= val1
    16892492                minsel = formatfunc(currsel.split('~')[0].strip())
    1690                 maxsel = formatfunc(currsel.split('~')[1].strip()) 
     2493                maxsel = formatfunc(currsel.split('~')[1].strip())
    16912494            elif currsel.strip().find('<=') > -1:
    16922495                bound = currsel.split('<=')
     
    19072710
    19082711    @asaplog_post_dec
    1909     def history(self, filename=None):
     2712    def history(self, filename=None, nrows=-1, start=0):
    19102713        """\
    19112714        Print the history. Optionally to a file.
     
    19162719
    19172720        """
    1918         hist = list(self._gethistory())
     2721        n = self._historylength()
     2722        if nrows == -1:
     2723            nrows = n
     2724        if start+nrows > n:
     2725            nrows = nrows-start
     2726        if n > 1000 and nrows == n:
     2727            nrows = 1000
     2728            start = n-1000
     2729            asaplog.push("Warning: History has {0} entries. Displaying last "
     2730                         "1000".format(n))
     2731        hist = list(self._gethistory(nrows, start))
    19192732        out = "-"*80
    19202733        for h in hist:
    1921             if h.startswith("---"):
    1922                 out = "\n".join([out, h])
     2734            if not h.strip():
     2735                continue
     2736            if h.find("---") >-1:
     2737                continue
    19232738            else:
    19242739                items = h.split("##")
     
    19332748                    s = i.split("=")
    19342749                    out += "\n   %s = %s" % (s[0], s[1])
    1935                 out = "\n".join([out, "-"*80])
     2750                out = "\n".join([out, "*"*80])
    19362751        if filename is not None:
    19372752            if filename is "":
    19382753                filename = 'scantable_history.txt'
    1939             import os
    19402754            filename = os.path.expandvars(os.path.expanduser(filename))
    19412755            if not os.path.isdir(filename):
     
    19522766    #
    19532767    @asaplog_post_dec
    1954     def average_time(self, mask=None, scanav=False, weight='tint', align=False):
     2768    def average_time(self, mask=None, scanav=False, weight='tint', align=False,
     2769                     avmode="NONE"):
    19552770        """\
    19562771        Return the (time) weighted average of a scan. Scans will be averaged
     
    19802795            align:    align the spectra in velocity before averaging. It takes
    19812796                      the time of the first spectrum as reference time.
     2797            avmode:   'SOURCE' - also select by source name -  or
     2798                      'NONE' (default). Not applicable for scanav=True or
     2799                      weight=median
    19822800
    19832801        Example::
     
    19902808        weight = weight or 'TINT'
    19912809        mask = mask or ()
    1992         scanav = (scanav and 'SCAN') or 'NONE'
     2810        scanav = (scanav and 'SCAN') or avmode.upper()
    19932811        scan = (self, )
    19942812
    19952813        if align:
    19962814            scan = (self.freq_align(insitu=False), )
     2815            asaplog.push("Note: Alignment is don on a source-by-source basis")
     2816            asaplog.push("Note: Averaging (by default) is not")
     2817            # we need to set it to SOURCE averaging here           
    19972818        s = None
    19982819        if weight.upper() == 'MEDIAN':
     
    25393360            raise RuntimeError(msg)
    25403361
    2541 
    2542     @asaplog_post_dec
    2543     def sinusoid_baseline(self, insitu=None, mask=None, applyfft=None,
     3362    @asaplog_post_dec
     3363    def apply_bltable(self, insitu=None, retfitres=None, inbltable=None, outbltable=None, overwrite=None):
     3364        """\
     3365        Subtract baseline based on parameters written in Baseline Table.
     3366
     3367        Parameters:
     3368            insitu:        if True, baseline fitting/subtraction is done
     3369                           in-situ. If False, a new scantable with
     3370                           baseline subtracted is returned. Actually,
     3371                           format of the returned value depends on both
     3372                           insitu and retfitres (see below).
     3373                           The default is taken from .asaprc (False)
     3374            retfitres:     if True, the results of baseline fitting (i.e.,
     3375                           coefficients and rms) are returned.
     3376                           default is False.
     3377                           The format of the returned value of this
     3378                           function varies as follows:
     3379                           (1) in case insitu=True and retfitres=True:
     3380                               fitting result.
     3381                           (2) in case insitu=True and retfitres=False:
     3382                               None.
     3383                           (3) in case insitu=False and retfitres=True:
     3384                               a dictionary containing a new scantable
     3385                               (with baseline subtracted) and the fitting
     3386                               results.
     3387                           (4) in case insitu=False and retfitres=False:
     3388                               a new scantable (with baseline subtracted).
     3389            inbltable:     name of input baseline table. The row number of
     3390                           scantable and that of inbltable must be
     3391                           identical.
     3392            outbltable:    name of output baseline table where baseline
     3393                           parameters and fitting results recorded.
     3394                           default is ''(no output).
     3395            overwrite:     if True when an existing baseline table is
     3396                           specified for outbltable, overwrites it.
     3397                           Otherwise there is no harm.
     3398                           default is False.
     3399        """
     3400
     3401        try:
     3402            varlist = vars()
     3403            if retfitres      is None: retfitres      = False
     3404            if inbltable      is None: raise ValueError("bltable missing.")
     3405            if outbltable     is None: outbltable     = ''
     3406            if overwrite      is None: overwrite      = False
     3407
     3408            if insitu is None: insitu = rcParams['insitu']
     3409            if insitu:
     3410                workscan = self
     3411            else:
     3412                workscan = self.copy()
     3413
     3414            sres = workscan._apply_bltable(inbltable,
     3415                                           retfitres,
     3416                                           outbltable,
     3417                                           os.path.exists(outbltable),
     3418                                           overwrite)
     3419            if retfitres: res = parse_fitresult(sres)
     3420
     3421            workscan._add_history('apply_bltable', varlist)
     3422
     3423            if insitu:
     3424                self._assign(workscan)
     3425                if retfitres:
     3426                    return res
     3427                else:
     3428                    return None
     3429            else:
     3430                if retfitres:
     3431                    return {'scantable': workscan, 'fitresults': res}
     3432                else:
     3433                    return workscan
     3434       
     3435        except RuntimeError, e:
     3436            raise_fitting_failure_exception(e)
     3437
     3438    @asaplog_post_dec
     3439    def sub_baseline(self, insitu=None, retfitres=None, blinfo=None, bltable=None, overwrite=None):
     3440        """\
     3441        Subtract baseline based on parameters written in the input list.
     3442
     3443        Parameters:
     3444            insitu:        if True, baseline fitting/subtraction is done
     3445                           in-situ. If False, a new scantable with
     3446                           baseline subtracted is returned. Actually,
     3447                           format of the returned value depends on both
     3448                           insitu and retfitres (see below).
     3449                           The default is taken from .asaprc (False)
     3450            retfitres:     if True, the results of baseline fitting (i.e.,
     3451                           coefficients and rms) are returned.
     3452                           default is False.
     3453                           The format of the returned value of this
     3454                           function varies as follows:
     3455                           (1) in case insitu=True and retfitres=True:
     3456                               fitting result.
     3457                           (2) in case insitu=True and retfitres=False:
     3458                               None.
     3459                           (3) in case insitu=False and retfitres=True:
     3460                               a dictionary containing a new scantable
     3461                               (with baseline subtracted) and the fitting
     3462                               results.
     3463                           (4) in case insitu=False and retfitres=False:
     3464                               a new scantable (with baseline subtracted).
     3465            blinfo:        baseline parameter set stored in a dictionary
     3466                           or a list of dictionary. Each dictionary
     3467                           corresponds to each spectrum and must contain
     3468                           the following keys and values:
     3469                             'row': row number,
     3470                             'blfunc': function name. available ones include
     3471                                       'poly', 'chebyshev', 'cspline' and
     3472                                       'sinusoid',
     3473                             'order': maximum order of polynomial. needed
     3474                                      if blfunc='poly' or 'chebyshev',
     3475                             'npiece': number or piecewise polynomial.
     3476                                       needed if blfunc='cspline',
     3477                             'nwave': a list of sinusoidal wave numbers.
     3478                                      needed if blfunc='sinusoid', and
     3479                             'masklist': min-max windows for channel mask.
     3480                                         the specified ranges will be used
     3481                                         for fitting.
     3482            bltable:       name of output baseline table where baseline
     3483                           parameters and fitting results recorded.
     3484                           default is ''(no output).
     3485            overwrite:     if True when an existing baseline table is
     3486                           specified for bltable, overwrites it.
     3487                           Otherwise there is no harm.
     3488                           default is False.
     3489                           
     3490        Example:
     3491            sub_baseline(blinfo=[{'row':0, 'blfunc':'poly', 'order':5,
     3492                                  'masklist':[[10,350],[352,510]]},
     3493                                 {'row':1, 'blfunc':'cspline', 'npiece':3,
     3494                                  'masklist':[[3,16],[19,404],[407,511]]}
     3495                                  ])
     3496
     3497                the first spectrum (row=0) will be fitted with polynomial
     3498                of order=5 and the next one (row=1) will be fitted with cubic
     3499                spline consisting of 3 pieces.
     3500        """
     3501
     3502        try:
     3503            varlist = vars()
     3504            if retfitres      is None: retfitres      = False
     3505            if blinfo         is None: blinfo         = []
     3506            if bltable        is None: bltable        = ''
     3507            if overwrite      is None: overwrite      = False
     3508
     3509            if insitu is None: insitu = rcParams['insitu']
     3510            if insitu:
     3511                workscan = self
     3512            else:
     3513                workscan = self.copy()
     3514
     3515            nrow = workscan.nrow()
     3516
     3517            in_blinfo = pack_blinfo(blinfo=blinfo, maxirow=nrow)
     3518
     3519            sres = workscan._sub_baseline(in_blinfo,
     3520                                          retfitres,
     3521                                          bltable,
     3522                                          os.path.exists(bltable),
     3523                                          overwrite)
     3524            if retfitres: res = parse_fitresult(sres)
     3525           
     3526            workscan._add_history('sub_baseline', varlist)
     3527
     3528            if insitu:
     3529                self._assign(workscan)
     3530                if retfitres:
     3531                    return res
     3532                else:
     3533                    return None
     3534            else:
     3535                if retfitres:
     3536                    return {'scantable': workscan, 'fitresults': res}
     3537                else:
     3538                    return workscan
     3539
     3540        except RuntimeError, e:
     3541            raise_fitting_failure_exception(e)
     3542
     3543    @asaplog_post_dec
     3544    def calc_aic(self, value=None, blfunc=None, order=None, mask=None,
     3545                 whichrow=None, uselinefinder=None, edge=None,
     3546                 threshold=None, chan_avg_limit=None):
     3547        """\
     3548        Calculates and returns model selection criteria for a specified
     3549        baseline model and a given spectrum data.
     3550        Available values include Akaike Information Criterion (AIC), the
     3551        corrected Akaike Information Criterion (AICc) by Sugiura(1978),
     3552        Bayesian Information Criterion (BIC) and the Generalised Cross
     3553        Validation (GCV).
     3554
     3555        Parameters:
     3556            value:         name of model selection criteria to calculate.
     3557                           available ones include 'aic', 'aicc', 'bic' and
     3558                           'gcv'. default is 'aicc'.
     3559            blfunc:        baseline function name. available ones include
     3560                           'chebyshev', 'cspline' and 'sinusoid'.
     3561                           default is 'chebyshev'.
     3562            order:         parameter for basline function. actually stands for
     3563                           order of polynomial (order) for 'chebyshev',
     3564                           number of spline pieces (npiece) for 'cspline' and
     3565                           maximum wave number for 'sinusoid', respectively.
     3566                           default is 5 (which is also the default order value
     3567                           for [auto_]chebyshev_baseline()).
     3568            mask:          an optional mask. default is [].
     3569            whichrow:      row number. default is 0 (the first row)
     3570            uselinefinder: use sd.linefinder() to flag out line regions
     3571                           default is True.
     3572            edge:           an optional number of channel to drop at
     3573                            the edge of spectrum. If only one value is
     3574                            specified, the same number will be dropped
     3575                            from both sides of the spectrum. Default
     3576                            is to keep all channels. Nested tuples
     3577                            represent individual edge selection for
     3578                            different IFs (a number of spectral channels
     3579                            can be different)
     3580                            default is (0, 0).
     3581            threshold:      the threshold used by line finder. It is
     3582                            better to keep it large as only strong lines
     3583                            affect the baseline solution.
     3584                            default is 3.
     3585            chan_avg_limit: a maximum number of consequtive spectral
     3586                            channels to average during the search of
     3587                            weak and broad lines. The default is no
     3588                            averaging (and no search for weak lines).
     3589                            If such lines can affect the fitted baseline
     3590                            (e.g. a high order polynomial is fitted),
     3591                            increase this parameter (usually values up
     3592                            to 8 are reasonable). Most users of this
     3593                            method should find the default value sufficient.
     3594                            default is 1.
     3595
     3596        Example:
     3597            aic = scan.calc_aic(blfunc='chebyshev', order=5, whichrow=0)
     3598        """
     3599
     3600        try:
     3601            varlist = vars()
     3602
     3603            if value          is None: value          = 'aicc'
     3604            if blfunc         is None: blfunc         = 'chebyshev'
     3605            if order          is None: order          = 5
     3606            if mask           is None: mask           = []
     3607            if whichrow       is None: whichrow       = 0
     3608            if uselinefinder  is None: uselinefinder  = True
     3609            if edge           is None: edge           = (0, 0)
     3610            if threshold      is None: threshold      = 3
     3611            if chan_avg_limit is None: chan_avg_limit = 1
     3612
     3613            return self._calc_aic(value, blfunc, order, mask,
     3614                                  whichrow, uselinefinder, edge,
     3615                                  threshold, chan_avg_limit)
     3616           
     3617        except RuntimeError, e:
     3618            raise_fitting_failure_exception(e)
     3619
     3620    @asaplog_post_dec
     3621    def sinusoid_baseline(self, mask=None, applyfft=None,
    25443622                          fftmethod=None, fftthresh=None,
    2545                           addwn=None, rejwn=None, clipthresh=None,
    2546                           clipniter=None, plot=None,
    2547                           getresidual=None, showprogress=None,
    2548                           minnrow=None, outlog=None, blfile=None, csvformat=None):
     3623                          addwn=None, rejwn=None,
     3624                          insitu=None,
     3625                          clipthresh=None, clipniter=None,
     3626                          plot=None,
     3627                          getresidual=None,
     3628                          showprogress=None, minnrow=None,
     3629                          outlog=None,
     3630                          blfile=None, csvformat=None,
     3631                          bltable=None):
    25493632        """\
    25503633        Return a scan which has been baselined (all rows) with sinusoidal
     
    25523635
    25533636        Parameters:
    2554             insitu:        if False a new scantable is returned.
    2555                            Otherwise, the scaling is done in-situ
    2556                            The default is taken from .asaprc (False)
    25573637            mask:          an optional mask
    25583638            applyfft:      if True use some method, such as FFT, to find
     
    25863666                           wave numbers which are specified both in addwn
    25873667                           and rejwn will NOT be used. default is [].
     3668            insitu:        if False a new scantable is returned.
     3669                           Otherwise, the scaling is done in-situ
     3670                           The default is taken from .asaprc (False)
    25883671            clipthresh:    Clipping threshold. (default is 3.0, unit: sigma)
    25893672            clipniter:     maximum number of iteration of 'clipthresh'-sigma
     
    26053688                           (default is '': no file/logger output)
    26063689            csvformat:     if True blfile is csv-formatted, default is False.
     3690            bltable:       name of a baseline table where fitting results
     3691                           (coefficients, rms, etc.) are to be written.
     3692                           if given, fitting results will NOT be output to
     3693                           scantable (insitu=True) or None will be
     3694                           returned (insitu=False).
     3695                           (default is "": no table output)
    26073696
    26083697        Example:
     
    26263715                workscan = self.copy()
    26273716           
    2628             #if mask          is None: mask          = [True for i in xrange(workscan.nchan())]
    26293717            if mask          is None: mask          = []
    26303718            if applyfft      is None: applyfft      = True
     
    26423730            if blfile        is None: blfile        = ''
    26433731            if csvformat     is None: csvformat     = False
    2644 
    2645             if csvformat:
    2646                 scsvformat = "T"
    2647             else:
    2648                 scsvformat = "F"
     3732            if bltable       is None: bltable       = ''
     3733
     3734            sapplyfft = 'true' if applyfft else 'false'
     3735            fftinfo = ','.join([sapplyfft, fftmethod.lower(), str(fftthresh).lower()])
     3736
     3737            scsvformat = 'T' if csvformat else 'F'
    26493738
    26503739            #CURRENTLY, PLOT=true is UNAVAILABLE UNTIL sinusoidal fitting is implemented as a fitter method.
    2651             workscan._sinusoid_baseline(mask, applyfft, fftmethod.lower(),
    2652                                         str(fftthresh).lower(),
     3740            workscan._sinusoid_baseline(mask,
     3741                                        fftinfo,
     3742                                        #applyfft, fftmethod.lower(),
     3743                                        #str(fftthresh).lower(),
    26533744                                        workscan._parse_wn(addwn),
    26543745                                        workscan._parse_wn(rejwn),
     
    26573748                                        pack_progress_params(showprogress,
    26583749                                                             minnrow),
    2659                                         outlog, scsvformat+blfile)
     3750                                        outlog, scsvformat+blfile,
     3751                                        bltable)
    26603752            workscan._add_history('sinusoid_baseline', varlist)
    2661            
    2662             if insitu:
    2663                 self._assign(workscan)
     3753
     3754            if bltable == '':
     3755                if insitu:
     3756                    self._assign(workscan)
     3757                else:
     3758                    return workscan
    26643759            else:
    2665                 return workscan
     3760                if not insitu:
     3761                    return None
    26663762           
    26673763        except RuntimeError, e:
     
    26703766
    26713767    @asaplog_post_dec
    2672     def auto_sinusoid_baseline(self, insitu=None, mask=None, applyfft=None,
     3768    def auto_sinusoid_baseline(self, mask=None, applyfft=None,
    26733769                               fftmethod=None, fftthresh=None,
    2674                                addwn=None, rejwn=None, clipthresh=None,
    2675                                clipniter=None, edge=None, threshold=None,
    2676                                chan_avg_limit=None, plot=None,
    2677                                getresidual=None, showprogress=None,
    2678                                minnrow=None, outlog=None, blfile=None, csvformat=None):
     3770                               addwn=None, rejwn=None,
     3771                               insitu=None,
     3772                               clipthresh=None, clipniter=None,
     3773                               edge=None, threshold=None, chan_avg_limit=None,
     3774                               plot=None,
     3775                               getresidual=None,
     3776                               showprogress=None, minnrow=None,
     3777                               outlog=None,
     3778                               blfile=None, csvformat=None,
     3779                               bltable=None):
    26793780        """\
    26803781        Return a scan which has been baselined (all rows) with sinusoidal
     
    26843785
    26853786        Parameters:
    2686             insitu:         if False a new scantable is returned.
    2687                             Otherwise, the scaling is done in-situ
    2688                             The default is taken from .asaprc (False)
    26893787            mask:           an optional mask retreived from scantable
    26903788            applyfft:       if True use some method, such as FFT, to find
     
    27183816                            wave numbers which are specified both in addwn
    27193817                            and rejwn will NOT be used. default is [].
     3818            insitu:         if False a new scantable is returned.
     3819                            Otherwise, the scaling is done in-situ
     3820                            The default is taken from .asaprc (False)
    27203821            clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
    27213822            clipniter:      maximum number of iteration of 'clipthresh'-sigma
     
    27573858                            (default is "": no file/logger output)
    27583859            csvformat:      if True blfile is csv-formatted, default is False.
     3860            bltable:        name of a baseline table where fitting results
     3861                            (coefficients, rms, etc.) are to be written.
     3862                            if given, fitting results will NOT be output to
     3863                            scantable (insitu=True) or None will be
     3864                            returned (insitu=False).
     3865                            (default is "": no table output)
    27593866
    27603867        Example:
     
    27753882                workscan = self.copy()
    27763883           
    2777             #if mask           is None: mask           = [True for i in xrange(workscan.nchan())]
    27783884            if mask           is None: mask           = []
    27793885            if applyfft       is None: applyfft       = True
     
    27943900            if blfile         is None: blfile         = ''
    27953901            if csvformat      is None: csvformat      = False
    2796 
    2797             if csvformat:
    2798                 scsvformat = "T"
    2799             else:
    2800                 scsvformat = "F"
     3902            if bltable        is None: bltable        = ''
     3903
     3904            sapplyfft = 'true' if applyfft else 'false'
     3905            fftinfo = ','.join([sapplyfft, fftmethod.lower(), str(fftthresh).lower()])
     3906
     3907            scsvformat = 'T' if csvformat else 'F'
    28013908
    28023909            #CURRENTLY, PLOT=true is UNAVAILABLE UNTIL sinusoidal fitting is implemented as a fitter method.
    2803             workscan._auto_sinusoid_baseline(mask, applyfft,
    2804                                              fftmethod.lower(),
    2805                                              str(fftthresh).lower(),
     3910            workscan._auto_sinusoid_baseline(mask,
     3911                                             fftinfo,
    28063912                                             workscan._parse_wn(addwn),
    28073913                                             workscan._parse_wn(rejwn),
     
    28123918                                             pack_progress_params(showprogress,
    28133919                                                                  minnrow),
    2814                                              outlog, scsvformat+blfile)
     3920                                             outlog, scsvformat+blfile, bltable)
    28153921            workscan._add_history("auto_sinusoid_baseline", varlist)
    2816            
    2817             if insitu:
    2818                 self._assign(workscan)
     3922
     3923            if bltable == '':
     3924                if insitu:
     3925                    self._assign(workscan)
     3926                else:
     3927                    return workscan
    28193928            else:
    2820                 return workscan
     3929                if not insitu:
     3930                    return None
    28213931           
    28223932        except RuntimeError, e:
     
    28243934
    28253935    @asaplog_post_dec
    2826     def cspline_baseline(self, insitu=None, mask=None, npiece=None,
     3936    def cspline_baseline(self, mask=None, npiece=None, insitu=None,
    28273937                         clipthresh=None, clipniter=None, plot=None,
    28283938                         getresidual=None, showprogress=None, minnrow=None,
    2829                          outlog=None, blfile=None, csvformat=None):
     3939                         outlog=None, blfile=None, csvformat=None,
     3940                         bltable=None):
    28303941        """\
    28313942        Return a scan which has been baselined (all rows) by cubic spline
     
    28333944
    28343945        Parameters:
     3946            mask:         An optional mask
     3947            npiece:       Number of pieces. (default is 2)
    28353948            insitu:       If False a new scantable is returned.
    28363949                          Otherwise, the scaling is done in-situ
    28373950                          The default is taken from .asaprc (False)
    2838             mask:         An optional mask
    2839             npiece:       Number of pieces. (default is 2)
    28403951            clipthresh:   Clipping threshold. (default is 3.0, unit: sigma)
    28413952            clipniter:    maximum number of iteration of 'clipthresh'-sigma
     
    28573968                          (default is "": no file/logger output)
    28583969            csvformat:    if True blfile is csv-formatted, default is False.
     3970            bltable:      name of a baseline table where fitting results
     3971                          (coefficients, rms, etc.) are to be written.
     3972                          if given, fitting results will NOT be output to
     3973                          scantable (insitu=True) or None will be
     3974                          returned (insitu=False).
     3975                          (default is "": no table output)
    28593976
    28603977        Example:
     
    28783995                workscan = self.copy()
    28793996
    2880             #if mask         is None: mask         = [True for i in xrange(workscan.nchan())]
    28813997            if mask         is None: mask         = []
    28823998            if npiece       is None: npiece       = 2
     
    28894005            if outlog       is None: outlog       = False
    28904006            if blfile       is None: blfile       = ''
    2891             if csvformat     is None: csvformat     = False
    2892 
    2893             if csvformat:
    2894                 scsvformat = "T"
    2895             else:
    2896                 scsvformat = "F"
     4007            if csvformat    is None: csvformat    = False
     4008            if bltable      is None: bltable      = ''
     4009
     4010            scsvformat = 'T' if csvformat else 'F'
    28974011
    28984012            #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method.
    2899             workscan._cspline_baseline(mask, npiece, clipthresh, clipniter,
     4013            workscan._cspline_baseline(mask, npiece,
     4014                                       clipthresh, clipniter,
    29004015                                       getresidual,
    29014016                                       pack_progress_params(showprogress,
    29024017                                                            minnrow),
    2903                                        outlog, scsvformat+blfile)
     4018                                       outlog, scsvformat+blfile,
     4019                                       bltable)
    29044020            workscan._add_history("cspline_baseline", varlist)
    2905            
    2906             if insitu:
    2907                 self._assign(workscan)
     4021
     4022            if bltable == '':
     4023                if insitu:
     4024                    self._assign(workscan)
     4025                else:
     4026                    return workscan
    29084027            else:
    2909                 return workscan
     4028                if not insitu:
     4029                    return None
    29104030           
    29114031        except RuntimeError, e:
     
    29134033
    29144034    @asaplog_post_dec
    2915     def auto_cspline_baseline(self, insitu=None, mask=None, npiece=None,
     4035    def auto_cspline_baseline(self, mask=None, npiece=None, insitu=None,
    29164036                              clipthresh=None, clipniter=None,
    29174037                              edge=None, threshold=None, chan_avg_limit=None,
    29184038                              getresidual=None, plot=None,
    29194039                              showprogress=None, minnrow=None, outlog=None,
    2920                               blfile=None, csvformat=None):
     4040                              blfile=None, csvformat=None, bltable=None):
    29214041        """\
    29224042        Return a scan which has been baselined (all rows) by cubic spline
     
    29264046
    29274047        Parameters:
     4048            mask:           an optional mask retreived from scantable
     4049            npiece:         Number of pieces. (default is 2)
    29284050            insitu:         if False a new scantable is returned.
    29294051                            Otherwise, the scaling is done in-situ
    29304052                            The default is taken from .asaprc (False)
    2931             mask:           an optional mask retreived from scantable
    2932             npiece:         Number of pieces. (default is 2)
    29334053            clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
    29344054            clipniter:      maximum number of iteration of 'clipthresh'-sigma
     
    29704090                            (default is "": no file/logger output)
    29714091            csvformat:      if True blfile is csv-formatted, default is False.
     4092            bltable:        name of a baseline table where fitting results
     4093                            (coefficients, rms, etc.) are to be written.
     4094                            if given, fitting results will NOT be output to
     4095                            scantable (insitu=True) or None will be
     4096                            returned (insitu=False).
     4097                            (default is "": no table output)
    29724098
    29734099        Example:
     
    30034129            if blfile         is None: blfile         = ''
    30044130            if csvformat      is None: csvformat      = False
    3005 
    3006             if csvformat:
    3007                 scsvformat = "T"
    3008             else:
    3009                 scsvformat = "F"
     4131            if bltable        is None: bltable        = ''
     4132
     4133            scsvformat = 'T' if csvformat else 'F'
    30104134
    30114135            #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method.
    3012             workscan._auto_cspline_baseline(mask, npiece, clipthresh,
    3013                                             clipniter,
     4136            workscan._auto_cspline_baseline(mask, npiece,
     4137                                            clipthresh, clipniter,
    30144138                                            normalise_edge_param(edge),
    30154139                                            threshold,
     
    30174141                                            pack_progress_params(showprogress,
    30184142                                                                 minnrow),
    3019                                             outlog, scsvformat+blfile)
     4143                                            outlog,
     4144                                            scsvformat+blfile,
     4145                                            bltable)
    30204146            workscan._add_history("auto_cspline_baseline", varlist)
    3021            
    3022             if insitu:
    3023                 self._assign(workscan)
     4147
     4148            if bltable == '':
     4149                if insitu:
     4150                    self._assign(workscan)
     4151                else:
     4152                    return workscan
    30244153            else:
    3025                 return workscan
     4154                if not insitu:
     4155                    return None
    30264156           
    30274157        except RuntimeError, e:
     
    30294159
    30304160    @asaplog_post_dec
    3031     def chebyshev_baseline(self, insitu=None, mask=None, order=None,
     4161    def chebyshev_baseline(self, mask=None, order=None, insitu=None,
    30324162                           clipthresh=None, clipniter=None, plot=None,
    30334163                           getresidual=None, showprogress=None, minnrow=None,
    3034                            outlog=None, blfile=None, csvformat=None):
     4164                           outlog=None, blfile=None, csvformat=None,
     4165                           bltable=None):
    30354166        """\
    30364167        Return a scan which has been baselined (all rows) by Chebyshev polynomials.
    30374168
    30384169        Parameters:
     4170            mask:         An optional mask
     4171            order:        the maximum order of Chebyshev polynomial (default is 5)
    30394172            insitu:       If False a new scantable is returned.
    30404173                          Otherwise, the scaling is done in-situ
    30414174                          The default is taken from .asaprc (False)
    3042             mask:         An optional mask
    3043             order:        the maximum order of Chebyshev polynomial (default is 5)
    30444175            clipthresh:   Clipping threshold. (default is 3.0, unit: sigma)
    30454176            clipniter:    maximum number of iteration of 'clipthresh'-sigma
     
    30614192                          (default is "": no file/logger output)
    30624193            csvformat:    if True blfile is csv-formatted, default is False.
     4194            bltable:      name of a baseline table where fitting results
     4195                          (coefficients, rms, etc.) are to be written.
     4196                          if given, fitting results will NOT be output to
     4197                          scantable (insitu=True) or None will be
     4198                          returned (insitu=False).
     4199                          (default is "": no table output)
    30634200
    30644201        Example:
     
    30824219                workscan = self.copy()
    30834220
    3084             #if mask         is None: mask         = [True for i in xrange(workscan.nchan())]
    30854221            if mask         is None: mask         = []
    30864222            if order        is None: order        = 5
     
    30934229            if outlog       is None: outlog       = False
    30944230            if blfile       is None: blfile       = ''
    3095             if csvformat     is None: csvformat     = False
    3096 
    3097             if csvformat:
    3098                 scsvformat = "T"
    3099             else:
    3100                 scsvformat = "F"
     4231            if csvformat    is None: csvformat    = False
     4232            if bltable      is None: bltable      = ''
     4233
     4234            scsvformat = 'T' if csvformat else 'F'
    31014235
    31024236            #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method.
    3103             workscan._chebyshev_baseline(mask, order, clipthresh, clipniter,
     4237            workscan._chebyshev_baseline(mask, order,
     4238                                         clipthresh, clipniter,
    31044239                                         getresidual,
    31054240                                         pack_progress_params(showprogress,
    31064241                                                              minnrow),
    3107                                          outlog, scsvformat+blfile)
     4242                                         outlog, scsvformat+blfile,
     4243                                         bltable)
    31084244            workscan._add_history("chebyshev_baseline", varlist)
    3109            
    3110             if insitu:
    3111                 self._assign(workscan)
     4245
     4246            if bltable == '':
     4247                if insitu:
     4248                    self._assign(workscan)
     4249                else:
     4250                    return workscan
    31124251            else:
    3113                 return workscan
     4252                if not insitu:
     4253                    return None
    31144254           
    31154255        except RuntimeError, e:
     
    31174257
    31184258    @asaplog_post_dec
    3119     def auto_chebyshev_baseline(self, insitu=None, mask=None, order=None,
     4259    def auto_chebyshev_baseline(self, mask=None, order=None, insitu=None,
    31204260                              clipthresh=None, clipniter=None,
    31214261                              edge=None, threshold=None, chan_avg_limit=None,
    31224262                              getresidual=None, plot=None,
    31234263                              showprogress=None, minnrow=None, outlog=None,
    3124                               blfile=None, csvformat=None):
     4264                              blfile=None, csvformat=None, bltable=None):
    31254265        """\
    31264266        Return a scan which has been baselined (all rows) by Chebyshev polynomials.
     
    31294269
    31304270        Parameters:
     4271            mask:           an optional mask retreived from scantable
     4272            order:          the maximum order of Chebyshev polynomial (default is 5)
    31314273            insitu:         if False a new scantable is returned.
    31324274                            Otherwise, the scaling is done in-situ
    31334275                            The default is taken from .asaprc (False)
    3134             mask:           an optional mask retreived from scantable
    3135             order:          the maximum order of Chebyshev polynomial (default is 5)
    31364276            clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
    31374277            clipniter:      maximum number of iteration of 'clipthresh'-sigma
     
    31734313                            (default is "": no file/logger output)
    31744314            csvformat:      if True blfile is csv-formatted, default is False.
     4315            bltable:        name of a baseline table where fitting results
     4316                            (coefficients, rms, etc.) are to be written.
     4317                            if given, fitting results will NOT be output to
     4318                            scantable (insitu=True) or None will be
     4319                            returned (insitu=False).
     4320                            (default is "": no table output)
    31754321
    31764322        Example:
     
    31914337                workscan = self.copy()
    31924338           
    3193             #if mask           is None: mask           = [True for i in xrange(workscan.nchan())]
    31944339            if mask           is None: mask           = []
    31954340            if order          is None: order          = 5
     
    32064351            if blfile         is None: blfile         = ''
    32074352            if csvformat      is None: csvformat      = False
    3208 
    3209             if csvformat:
    3210                 scsvformat = "T"
    3211             else:
    3212                 scsvformat = "F"
     4353            if bltable        is None: bltable        = ''
     4354
     4355            scsvformat = 'T' if csvformat else 'F'
    32134356
    32144357            #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method.
    3215             workscan._auto_chebyshev_baseline(mask, order, clipthresh,
    3216                                               clipniter,
     4358            workscan._auto_chebyshev_baseline(mask, order,
     4359                                              clipthresh, clipniter,
    32174360                                              normalise_edge_param(edge),
    32184361                                              threshold,
     
    32204363                                              pack_progress_params(showprogress,
    32214364                                                                   minnrow),
    3222                                               outlog, scsvformat+blfile)
     4365                                              outlog, scsvformat+blfile,
     4366                                              bltable)
    32234367            workscan._add_history("auto_chebyshev_baseline", varlist)
    3224            
    3225             if insitu:
    3226                 self._assign(workscan)
     4368
     4369            if bltable == '':
     4370                if insitu:
     4371                    self._assign(workscan)
     4372                else:
     4373                    return workscan
    32274374            else:
    3228                 return workscan
     4375                if not insitu:
     4376                    return None
    32294377           
    32304378        except RuntimeError, e:
     
    32324380
    32334381    @asaplog_post_dec
    3234     def poly_baseline(self, mask=None, order=None, insitu=None, plot=None,
     4382    def poly_baseline(self, mask=None, order=None, insitu=None,
     4383                      clipthresh=None, clipniter=None, plot=None,
    32354384                      getresidual=None, showprogress=None, minnrow=None,
    3236                       outlog=None, blfile=None, csvformat=None):
     4385                      outlog=None, blfile=None, csvformat=None,
     4386                      bltable=None):
    32374387        """\
    32384388        Return a scan which has been baselined (all rows) by a polynomial.
    32394389        Parameters:
     4390            mask:         an optional mask
     4391            order:        the order of the polynomial (default is 0)
    32404392            insitu:       if False a new scantable is returned.
    32414393                          Otherwise, the scaling is done in-situ
    32424394                          The default is taken from .asaprc (False)
    3243             mask:         an optional mask
    3244             order:        the order of the polynomial (default is 0)
     4395            clipthresh:   Clipping threshold. (default is 3.0, unit: sigma)
     4396            clipniter:    maximum number of iteration of 'clipthresh'-sigma
     4397                          clipping (default is 0)
    32454398            plot:         plot the fit and the residual. In this each
    32464399                          indivual fit has to be approved, by typing 'y'
     
    32584411                          (default is "": no file/logger output)
    32594412            csvformat:    if True blfile is csv-formatted, default is False.
     4413            bltable:      name of a baseline table where fitting results
     4414                          (coefficients, rms, etc.) are to be written.
     4415                          if given, fitting results will NOT be output to
     4416                          scantable (insitu=True) or None will be
     4417                          returned (insitu=False).
     4418                          (default is "": no table output)
    32604419
    32614420        Example:
     
    32754434                workscan = self.copy()
    32764435
    3277             #if mask         is None: mask         = [True for i in \
    3278             #                                           xrange(workscan.nchan())]
    32794436            if mask         is None: mask         = []
    32804437            if order        is None: order        = 0
     4438            if clipthresh   is None: clipthresh   = 3.0
     4439            if clipniter    is None: clipniter    = 0
    32814440            if plot         is None: plot         = False
    32824441            if getresidual  is None: getresidual  = True
     
    32844443            if minnrow      is None: minnrow      = 1000
    32854444            if outlog       is None: outlog       = False
    3286             if blfile       is None: blfile       = ""
     4445            if blfile       is None: blfile       = ''
    32874446            if csvformat    is None: csvformat    = False
    3288 
    3289             if csvformat:
    3290                 scsvformat = "T"
    3291             else:
    3292                 scsvformat = "F"
     4447            if bltable      is None: bltable      = ''
     4448
     4449            scsvformat = 'T' if csvformat else 'F'
    32934450
    32944451            if plot:
     
    33544511                    blf.close()
    33554512            else:
    3356                 workscan._poly_baseline(mask, order, getresidual,
     4513                workscan._poly_baseline(mask, order,
     4514                                        clipthresh, clipniter, #
     4515                                        getresidual,
    33574516                                        pack_progress_params(showprogress,
    33584517                                                             minnrow),
    3359                                         outlog, scsvformat+blfile)
     4518                                        outlog, scsvformat+blfile,
     4519                                        bltable)  #
    33604520           
    33614521            workscan._add_history("poly_baseline", varlist)
     
    33704530
    33714531    @asaplog_post_dec
    3372     def auto_poly_baseline(self, mask=None, order=None, edge=None,
    3373                            threshold=None, chan_avg_limit=None,
    3374                            plot=None, insitu=None,
    3375                            getresidual=None, showprogress=None,
    3376                            minnrow=None, outlog=None, blfile=None, csvformat=None):
     4532    def auto_poly_baseline(self, mask=None, order=None, insitu=None,
     4533                           clipthresh=None, clipniter=None,
     4534                           edge=None, threshold=None, chan_avg_limit=None,
     4535                           getresidual=None, plot=None,
     4536                           showprogress=None, minnrow=None, outlog=None,
     4537                           blfile=None, csvformat=None, bltable=None):
    33774538        """\
    33784539        Return a scan which has been baselined (all rows) by a polynomial.
     
    33814542
    33824543        Parameters:
     4544            mask:           an optional mask retreived from scantable
     4545            order:          the order of the polynomial (default is 0)
    33834546            insitu:         if False a new scantable is returned.
    33844547                            Otherwise, the scaling is done in-situ
    33854548                            The default is taken from .asaprc (False)
    3386             mask:           an optional mask retreived from scantable
    3387             order:          the order of the polynomial (default is 0)
     4549            clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
     4550            clipniter:      maximum number of iteration of 'clipthresh'-sigma
     4551                            clipping (default is 0)
    33884552            edge:           an optional number of channel to drop at
    33894553                            the edge of spectrum. If only one value is
     
    34214585                            (default is "": no file/logger output)
    34224586            csvformat:      if True blfile is csv-formatted, default is False.
     4587            bltable:        name of a baseline table where fitting results
     4588                            (coefficients, rms, etc.) are to be written.
     4589                            if given, fitting results will NOT be output to
     4590                            scantable (insitu=True) or None will be
     4591                            returned (insitu=False).
     4592                            (default is "": no table output)
    34234593
    34244594        Example:
     
    34364606                workscan = self.copy()
    34374607
    3438             #if mask           is None: mask           = [True for i in xrange(workscan.nchan())]
    34394608            if mask           is None: mask           = []
    34404609            if order          is None: order          = 0
     4610            if clipthresh     is None: clipthresh     = 3.0
     4611            if clipniter      is None: clipniter      = 0
    34414612            if edge           is None: edge           = (0, 0)
    34424613            if threshold      is None: threshold      = 3
     
    34494620            if blfile         is None: blfile         = ''
    34504621            if csvformat      is None: csvformat      = False
    3451 
    3452             if csvformat:
    3453                 scsvformat = "T"
    3454             else:
    3455                 scsvformat = "F"
     4622            if bltable        is None: bltable        = ''
     4623
     4624            scsvformat = 'T' if csvformat else 'F'
    34564625
    34574626            edge = normalise_edge_param(edge)
     
    35234692                if outblfile: blf.close()
    35244693            else:
    3525                 workscan._auto_poly_baseline(mask, order, edge, threshold,
     4694                workscan._auto_poly_baseline(mask, order,
     4695                                             clipthresh, clipniter,
     4696                                             edge, threshold,
    35264697                                             chan_avg_limit, getresidual,
    35274698                                             pack_progress_params(showprogress,
    35284699                                                                  minnrow),
    3529                                              outlog, scsvformat+blfile)
    3530 
     4700                                             outlog, scsvformat+blfile,
     4701                                             bltable)
    35314702            workscan._add_history("auto_poly_baseline", varlist)
    3532            
    3533             if insitu:
    3534                 self._assign(workscan)
     4703
     4704            if bltable == '':
     4705                if insitu:
     4706                    self._assign(workscan)
     4707                else:
     4708                    return workscan
    35354709            else:
    3536                 return workscan
     4710                if not insitu:
     4711                    return None
    35374712           
    35384713        except RuntimeError, e:
     
    36444819        self._math._setinsitu(insitu)
    36454820        varlist = vars()
    3646         s = scantable(self._math._unaryop(self, offset, "ADD", False))
     4821        s = scantable(self._math._unaryop(self, offset, "ADD", False, False))
    36474822        s._add_history("add", varlist)
    36484823        if insitu:
     
    36674842            tsys:        if True (default) then apply the operation to Tsys
    36684843                         as well as the data
    3669 
    36704844        """
    36714845        if insitu is None: insitu = rcParams['insitu']
     
    36784852                                                         numpy.ndarray):
    36794853                from asapmath import _array2dOp
    3680                 s = _array2dOp( self, factor, "MUL", tsys, insitu )
     4854                s = _array2dOp( self, factor, "MUL", tsys, insitu, True )
    36814855            else:
    36824856                s = scantable( self._math._arrayop( self, factor,
    3683                                                     "MUL", tsys ) )
     4857                                                    "MUL", tsys, True ) )
    36844858        else:
    3685             s = scantable(self._math._unaryop(self, factor, "MUL", tsys))
     4859            s = scantable(self._math._unaryop(self, factor, "MUL", tsys, True ))
    36864860        s._add_history("scale", varlist)
    36874861        if insitu:
     
    37304904        self._add_history("set_sourcetype", varlist)
    37314905
     4906
     4907    def set_sourcename(self, name):
     4908        varlist = vars()
     4909        self._setsourcename(name)
     4910        self._add_history("set_sourcename", varlist)
     4911
    37324912    @asaplog_post_dec
    37334913    @preserve_selection
     
    37664946        s = None
    37674947        if mode.lower() == "paired":
     4948            from asap._asap import srctype
    37684949            sel = self.get_selection()
    3769             sel.set_query("SRCTYPE==psoff")
     4950            #sel.set_query("SRCTYPE==psoff")
     4951            sel.set_types(srctype.psoff)
    37704952            self.set_selection(sel)
    37714953            offs = self.copy()
    3772             sel.set_query("SRCTYPE==pson")
     4954            #sel.set_query("SRCTYPE==pson")
     4955            sel.set_types(srctype.pson)
    37734956            self.set_selection(sel)
    37744957            ons = self.copy()
     
    38875070            if other == 0.0:
    38885071                raise ZeroDivisionError("Dividing by zero is not recommended")
    3889             s = scantable(self._math._unaryop(self, other, mode, False))
     5072            s = scantable(self._math._unaryop(self, other, mode, False, True))
    38905073        elif isinstance(other, list) or isinstance(other, numpy.ndarray):
    38915074            if isinstance(other[0], list) \
    38925075                    or isinstance(other[0], numpy.ndarray):
    38935076                from asapmath import _array2dOp
    3894                 s = _array2dOp( self, other, mode, False )
     5077                s = _array2dOp(self, other, mode, False)
    38955078            else:
    3896                 s = scantable( self._math._arrayop( self, other,
    3897                                                     mode, False ) )
     5079                s = scantable(self._math._arrayop(self, other, mode, False, True))
    38985080        else:
    38995081            raise TypeError("Other input is not a scantable or float value")
     
    39285110            self.set_selection(basesel+sel)
    39295111            nans = numpy.isnan(self._getspectrum(0))
    3930         if numpy.any(nans):
    3931             bnans = [ bool(v) for v in nans]
    3932             self.flag(bnans)
     5112            if numpy.any(nans):
     5113                bnans = [ bool(v) for v in nans]
     5114                self.flag(bnans)
     5115       
     5116        self.set_selection(basesel)
    39335117
    39345118    def get_row_selector(self, rowno):
  • trunk/python/selector.py

    r2704 r3112  
    11import re
     2import math
     3import string
    24from asap._asap import selector as _selector, srctype
    35from asap.utils import unique, _to_list
     
    200202            raise TypeError('Unknown row number type. Use lists of integers.')
    201203
     204    def set_msselection_field(self, selection):
     205        """
     206        Set a field selection in msselection syntax. The msselection
     207        suppports the following syntax:
     208
     209        pattern match:
     210            - UNIX style pattern match for source name using '*'
     211              (compatible with set_name)
     212
     213        field id selection:
     214            - simple number in string ('0', '1', etc.)
     215            - range specification using '~' ('0~1', etc.)
     216            - range specification using '>' or '<' in combination
     217              with '=' ('>=1', '<3', etc.)
     218
     219        comma separated multiple selection:
     220            - selections can be combined by using ',' ('0,>1',
     221              'mysource*,2~4', etc.)
     222        """
     223        selection_list =  map(string.strip, selection.split(','))
     224        query_list = list(self.generate_query(selection_list))
     225        if len(query_list) > 0:
     226            original_query = self.get_query()
     227            if len(original_query) == 0 or re.match('.*(SRC|FIELD)NAME.*',original_query):
     228                query = 'SELECT FROM $1 WHERE ' + ' || '.join(query_list)
     229            else:
     230                query = 'SELECT FROM $1 WHERE (' + original_query + ') && (' + ' || '.join(query_list) + ')'
     231            self._settaql(query)
     232
     233    def generate_query(self, selection_list):
     234        for s in selection_list:
     235            if s.isdigit() or re.match('^[<>]=?[0-9]*$', s) or \
     236                    re.match('^[0-9]+~[0-9]+$', s):
     237                #print '"%s" is ID selection using < or <='%(s)
     238                a = FieldIdRegexGenerator(s)
     239                yield '(%s)'%(a.get_regex())
     240            elif len(s) > 0:
     241                #print '"%s" is UNIX style pattern match'%(s)
     242                yield '(SRCNAME == pattern(\'%s\'))'%(s)
     243       
    202244    def get_scans(self):
    203245        return list(self._getscans())
     
    275317                union.set_query(qs)
    276318        return union
     319
     320class FieldIdRegexGenerator(object):
     321    def __init__(self, pattern):
     322        if pattern.isdigit():
     323            self.regex = 'FIELDNAME == regex(\'.+__%s$\')'%(pattern)
     324        else:
     325            self.regex = None
     326            ineq = None
     327            if pattern.find('<') >= 0:
     328                ineq = '<'
     329                s = pattern.strip().lstrip(ineq).lstrip('=')
     330                if not s.isdigit():
     331                    raise RuntimeError('Invalid syntax: %s'%(pattern))
     332                self.id = int(s) + (-1 if pattern.find('=') < 0 else 0)
     333                self.template = string.Template('FIELDNAME == regex(\'.+__${reg}$\')')
     334            elif pattern.find('>') >= 0:
     335                ineq = '>'
     336                s = pattern.strip().lstrip(ineq).lstrip('=')
     337                if not s.isdigit():
     338                    raise RuntimeError('Invalid syntax: %s'%(pattern))
     339                self.id = int(s) + (-1 if pattern.find('=') >= 0 else 0)
     340                self.template = string.Template('FIELDNAME == regex(\'.+__[0-9]+$\') && FIELDNAME != regex(\'.+__${reg}$\')')
     341            elif pattern.find('~') >= 0:
     342                s = map(string.strip, pattern.split('~'))
     343                if len(s) == 2 and s[0].isdigit() and s[1].isdigit():
     344                    id0 = int(s[0])
     345                    id1 = int(s[1])
     346                    if id0 == 0:
     347                        self.id = id1
     348                        self.template = string.Template('FIELDNAME == regex(\'.+__${reg}$\')')
     349                    else:
     350                        self.id = [id0-1,id1]
     351                        self.template = string.Template('FIELDNAME == regex(\'.+__${reg}$\') && FIELDNAME != regex(\'.+__${optreg}$\')')
     352                else:
     353                    raise RuntimeError('Invalid syntax: %s'%(pattern))
     354            else:
     355                raise RuntimeError('Invalid syntax: %s'%(pattern))
     356            #print 'self.id=',self.id
     357
     358    def get_regex(self):
     359        if self.regex is not None:
     360            # 'X'
     361            return self.regex
     362        elif isinstance(self.id, list):
     363            # 'X~Y'
     364            return self.template.safe_substitute(reg=self.__compile(self.id[1]),
     365                                                 optreg=self.__compile(self.id[0]))
     366        else:
     367            # '<(=)X' or '>(=)X'
     368            return self.template.safe_substitute(reg=self.__compile(self.id))
     369
     370    def __compile(self, idx):
     371        pattern = ''
     372        if idx >= 0:
     373            numerics = map(int,list(str(idx)))
     374            #numerics.reverse()
     375            num_digits = len(numerics)
     376            #print 'numerics=',numerics
     377            if num_digits == 1:
     378                if numerics[0] == 0:
     379                    pattern = '0'
     380                else:
     381                    pattern = '[0-%s]'%(numerics[0])
     382            elif num_digits == 2:
     383                pattern = '(%s)'%('|'.join(
     384                        list(self.__gen_two_digit_pattern(numerics))))
     385            elif num_digits == 3:
     386                pattern = '(%s)'%('|'.join(
     387                        list(self.__gen_three_digit_pattern(numerics))))
     388            else:
     389                raise RuntimeError('ID > 999 is not supported')
     390        else:
     391            raise RuntimeError('ID must be >= 0')
     392        return pattern
     393
     394    def __gen_two_digit_pattern(self, numerics):
     395        assert len(numerics) == 2
     396        yield '[0-9]'
     397        if numerics[0] == 2:
     398            yield '1[0-9]'
     399        elif numerics[0] > 2:
     400            yield '[1-%s][0-9]'%(numerics[0]-1)
     401        if numerics[1] == 0:
     402            yield '%s%s'%(numerics[0],numerics[1])
     403        else:
     404            yield '%s[0-%s]'%(numerics[0],numerics[1])
     405
     406    def __gen_three_digit_pattern(self, numerics):
     407        assert len(numerics) == 3
     408        yield '[0-9]'
     409        yield '[1-9][0-9]'
     410        if numerics[0] == 2:
     411            yield '1[0-9][0-9]'
     412        elif numerics[0] > 2:
     413            yield '[1-%s][0-9][0-9]'%(numerics[0]-1)
     414        if numerics[1] == 0:
     415            if numerics[2] == 0:
     416                yield '%s00'%(numerics[0])
     417            else:
     418                yield '%s0[0-%s]'%(numerics[0],numerics[2])
     419        else:
     420            if numerics[1] > 1:
     421                yield '%s[0-%s][0-9]'%(numerics[0],numerics[1]-1)
     422            elif numerics[1] == 1:
     423                yield '%s0[0-9]'%(numerics[0])
     424            if numerics[0] == 0:
     425                yield '%s%s%s'%(numerics[0],numerics[1],numerics[2])
     426            else:
     427                yield '%s%s[0-%s]'%(numerics[0],numerics[1],numerics[2])
  • trunk/python/utils.py

    r2704 r3112  
    1616
    1717def _is_sequence_or_number(param, ptype=int):
    18     if isinstance(param,tuple) or isinstance(param,list):
     18    import numpy
     19    #if isinstance(param,tuple) or isinstance(param,list):
     20    if type(param) in (tuple, list, numpy.ndarray):
    1921        if len(param) == 0: return True # empty list
    2022        out = True
     
    3133        else: return [param]
    3234    if _is_sequence_or_number(param, ptype):
    33         return param
     35        return list(param)
    3436    return None
    3537
Note: See TracChangeset for help on using the changeset viewer.